GEOS
GEOS copied to clipboard
[Cleanup] Redefine LAI set/insert/add interface functions...or delete them entirely
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
.
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
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.