GEOS icon indicating copy to clipboard operation
GEOS copied to clipboard

[Cleanup] Redefine LAI set/insert/add interface functions...or delete them entirely

Open rrsettgast opened this issue 2 years ago • 2 comments

Describe the issue The LAI objects have set/insert/add interface functions that take values, pointers, and arrayslice. This either requires the interface function to copy the data and move it to the correct memory space, or it requires the caller to manage the memory motion rather than allowing the interface functions to manage the memory motion as part of their responsibility. For example, MatrixBase contains the following interface functions:

  virtual void add( globalIndex const rowIndex,
                    globalIndex const colIndex,
                    real64 const value ) = 0;

  virtual void set( globalIndex const rowIndex,
                    globalIndex const colIndex,
                    real64 const value ) = 0;

  virtual void insert( globalIndex const rowIndex,
                       globalIndex const colIndex,
                       real64 const value ) = 0;

These are single point interface functions. They are very inefficient, and require that the interface function copy the data, then move it to the correct memory space. These should be removed.

virtual void add( globalIndex const rowIndex,
                    globalIndex const * colIndices,
                    real64 const * values,
                    localIndex const size ) = 0;

  virtual void set( globalIndex const rowIndex,
                    globalIndex const * colIndices,
                    real64 const * values,
                    localIndex const size ) = 0;

  virtual void insert( globalIndex const rowIndex,
                       globalIndex const * colIndices,
                       real64 const * values,
                       localIndex const size ) = 0;

  virtual void add( globalIndex const rowIndex,
                    arraySlice1d< globalIndex const > const & colIndices,
                    arraySlice1d< real64 const > const & values ) = 0;

  virtual void set( globalIndex const rowIndex,
                    arraySlice1d< globalIndex const > const & colIndices,
                    arraySlice1d< real64 const > const & values ) = 0;

  virtual void insert( globalIndex const rowIndex,
                       arraySlice1d< globalIndex const > const & colIndices,
                       arraySlice1d< real64 const > const & values ) = 0;

These are data insertions for a single row. Pointers/arraySlice make this interface hard to manage memory motion.

  virtual void add( arraySlice1d< globalIndex const > const & rowIndices,
                    arraySlice1d< globalIndex const > const & colIndices,
                    arraySlice2d< real64 const > const & values ) override;

  virtual void set( arraySlice1d< globalIndex const > const & rowIndices,
                    arraySlice1d< globalIndex const > const & colIndices,
                    arraySlice2d< real64 const > const & values ) override;

  virtual void insert( arraySlice1d< globalIndex const > const & rowIndices,
                       arraySlice1d< globalIndex const > const & colIndices,
                       arraySlice2d< real64 const > const & values ) override;

  virtual void add( globalIndex const * rowIndices,
                    globalIndex const * colIndices,
                    real64 const * values,
                    localIndex const numRows,
                    localIndex const numCols ) override;

  virtual void set( globalIndex const * rowIndices,
                    globalIndex const * colIndices,
                    real64 const * values,
                    localIndex const numRows,
                    localIndex const numCols ) override;

  virtual void insert( globalIndex const * rowIndices,
                       globalIndex const * colIndices,
                       real64 const * values,
                       localIndex const numRows,
                       localIndex const numCols ) override;

Dense data insertions. These are not as bad since they should be launched inside a kernel...but they are not __device__ callable.

Proposed cleanup Remove all of these functions. We should just be using CRSMatrix objects and creating these objects from CRSMatrix. We can then focus on getting a shallow copy from CRSMatrix to the derived classes of MatrixBase.

rrsettgast avatar Nov 02 '21 16:11 rrsettgast

It would be nice to better encapsulate our current assembly pattern, which is pretty raw, e.g.

https://github.com/GEOSX/GEOSX/blob/e1f643dc63297bde35e2f0eb3e2f3e394b4629ce/src/coreComponents/physicsSolvers/solidMechanics/SolidMechanicsSmallStrainQuasiStaticKernel.hpp#L290-L305

https://github.com/GEOSX/GEOSX/blob/e1f643dc63297bde35e2f0eb3e2f3e394b4629ce/src/coreComponents/physicsSolvers/multiphysics/CompositionalMultiphaseReservoir.cpp#L279-L289

joshua-white avatar Nov 02 '21 22:11 joshua-white

Proposed cleanup Remove all of these functions. We should just be using CRSMatrix objects and creating these objects from CRSMatrix. We can then focus on getting a shallow copy from CRSMatrix to the derived classes of MatrixBase.

For the hypre interface I would push the cleanup even further. The proposed changes are removing the need for the IJ layer. I think we can redesign the interface manipulating directly the low level hypre_ParCSRMatrix object - which is what we already do in plenty of methods.

This method in mfem (link) is a good example for what has to change in the current create() method:

https://github.com/GEOSX/GEOSX/blob/2baf4927d912262c60bb03c36d73a13f2564473e/src/coreComponents/linearAlgebra/interfaces/hypre/HypreMatrix.hpp#L131-L132

This effort may require a dedicated PR.

castelletto1 avatar Nov 04 '21 03:11 castelletto1