GEOS
PetscMatrix.hpp
Go to the documentation of this file.
1 /*
2  * ------------------------------------------------------------------------------------------------------------
3  * SPDX-License-Identifier: LGPL-2.1-only
4  *
5  * Copyright (c) 2016-2024 Lawrence Livermore National Security LLC
6  * Copyright (c) 2018-2024 TotalEnergies
7  * Copyright (c) 2018-2024 The Board of Trustees of the Leland Stanford Junior University
8  * Copyright (c) 2023-2024 Chevron
9  * Copyright (c) 2019- GEOS/GEOSX Contributors
10  * All rights reserved
11  *
12  * See top level LICENSE, COPYRIGHT, CONTRIBUTORS, NOTICE, and ACKNOWLEDGEMENTS files for details.
13  * ------------------------------------------------------------------------------------------------------------
14  */
15 
20 #ifndef GEOS_LINEARALGEBRA_INTERFACES_PETSCSPARSEMATRIX_HPP_
21 #define GEOS_LINEARALGEBRA_INTERFACES_PETSCSPARSEMATRIX_HPP_
22 
23 #include "common/DataTypes.hpp"
28 
36 
38 extern "C" struct _p_Mat;
39 
41 
42 namespace geos
43 {
44 
49 class PetscMatrix final : public virtual LinearOperator< PetscVector >,
50  private MatrixBase< PetscMatrix, PetscVector >
51 {
52 public:
53 
56 
59 
61  using Mat = struct _p_Mat *;
62 
67 
72 
77  PetscMatrix( PetscMatrix const & src );
78 
83  PetscMatrix( PetscMatrix && src ) noexcept;
84 
90  PetscMatrix & operator=( PetscMatrix const & src );
91 
97  PetscMatrix & operator=( PetscMatrix && src ) noexcept;
98 
102  ~PetscMatrix() override;
103 
105 
110 
113  using MatrixBase::create;
114  using MatrixBase::closed;
115  using MatrixBase::assembled;
118  using MatrixBase::ready;
119  using MatrixBase::residual;
122  using MatrixBase::extract;
126 
127  virtual void createWithLocalSize( localIndex const localRows,
128  localIndex const localCols,
129  localIndex const maxEntriesPerRow,
130  MPI_Comm const & comm ) override;
131 
132  virtual void createWithGlobalSize( globalIndex const globalRows,
133  globalIndex const globalCols,
134  localIndex const maxEntriesPerRow,
135  MPI_Comm const & comm ) override;
136 
140  virtual bool created() const override;
141 
142  virtual void reset() override;
143 
144  virtual void set( real64 const value ) override;
145 
146  virtual void zero() override;
147 
148  virtual void open() override;
149 
150  virtual void close() override;
151 
152  virtual void add( globalIndex const rowIndex,
153  globalIndex const colIndex,
154  real64 const value ) override;
155 
156  virtual void set( globalIndex const rowIndex,
157  globalIndex const colIndex,
158  real64 const value ) override;
159 
160  virtual void insert( globalIndex const rowIndex,
161  globalIndex const colIndex,
162  real64 const value ) override;
163 
164  virtual void add( globalIndex const rowIndex,
165  globalIndex const * colIndices,
166  real64 const * values,
167  localIndex const size ) override;
168 
169  virtual void set( globalIndex const rowIndex,
170  globalIndex const * colIndices,
171  real64 const * values,
172  localIndex const size ) override;
173 
174  virtual void insert( globalIndex const rowIndex,
175  globalIndex const * colIndices,
176  real64 const * values,
177  localIndex const size ) override;
178 
179  virtual void add( globalIndex const rowIndex,
180  arraySlice1d< globalIndex const > const & colIndices,
181  arraySlice1d< real64 const > const & values ) override;
182 
183  virtual void set( globalIndex const rowIndex,
184  arraySlice1d< globalIndex const > const & colIndices,
185  arraySlice1d< real64 const > const & values ) override;
186 
187  virtual void insert( globalIndex const rowIndex,
188  arraySlice1d< globalIndex const > const & colIndices,
189  arraySlice1d< real64 const > const & values ) override;
190 
191  virtual void add( arraySlice1d< globalIndex const > const & rowIndices,
192  arraySlice1d< globalIndex const > const & colIndices,
193  arraySlice2d< real64 const > const & values ) override;
194 
195  virtual void set( arraySlice1d< globalIndex const > const & rowIndices,
196  arraySlice1d< globalIndex const > const & colIndices,
197  arraySlice2d< real64 const > const & values ) override;
198 
199  virtual void insert( arraySlice1d< globalIndex const > const & rowIndices,
200  arraySlice1d< globalIndex const > const & colIndices,
201  arraySlice2d< real64 const > const & values ) override;
202 
203  virtual void add( globalIndex const * rowIndices,
204  globalIndex const * colIndices,
205  real64 const * values,
206  localIndex const numRows,
207  localIndex const numCols ) override;
208 
209  virtual void set( globalIndex const * rowIndices,
210  globalIndex const * colIndices,
211  real64 const * values,
212  localIndex const numRows,
213  localIndex const numCols ) override;
214 
215  virtual void insert( globalIndex const * rowIndices,
216  globalIndex const * colIndices,
217  real64 const * values,
218  localIndex const numRows,
219  localIndex const numCols ) override;
220 
221  virtual void insert( arrayView1d< globalIndex const > const & rowIndices,
222  arrayView1d< globalIndex const > const & colIndices,
223  arrayView1d< real64 const > const & values ) override;
224 
225  virtual void apply( PetscVector const & src,
226  PetscVector & dst ) const override;
227 
228  virtual void multiply( PetscMatrix const & src,
229  PetscMatrix & dst ) const override;
230 
231  virtual void applyTranspose( Vector const & src,
232  Vector & dst ) const override;
233 
234  virtual void leftMultiplyTranspose( PetscMatrix const & src,
235  PetscMatrix & dst ) const override;
236 
237  virtual void rightMultiplyTranspose( PetscMatrix const & src,
238  PetscMatrix & dst ) const override;
239 
240  virtual void multiplyRAP( PetscMatrix const & R,
241  PetscMatrix const & P,
242  PetscMatrix & dst ) const override;
243 
244  virtual void multiplyPtAP( PetscMatrix const & P,
245  PetscMatrix & dst ) const override;
246 
247  virtual void gemv( real64 const alpha,
248  PetscVector const & x,
249  real64 const beta,
250  PetscVector & y,
251  bool useTranspose = false ) const override;
252 
253  virtual void scale( real64 const scalingFactor ) override;
254 
255  virtual void leftScale( PetscVector const & vec ) override;
256 
257  virtual void rightScale( PetscVector const & vec ) override;
258 
259  virtual void leftRightScale( PetscVector const & vecLeft,
260  PetscVector const & vecRight ) override;
261 
262  virtual void computeScalingVector( PetscVector & scaling ) const override;
263 
264  virtual void rescaleRows( arrayView1d< globalIndex const > const & rowIndices,
265  RowSumType const rowSumType ) override;
266 
267  virtual void transpose( PetscMatrix & dst ) const override;
268 
270  integer const dofsPerNode ) const override;
271 
272  virtual real64 clearRow( globalIndex const row,
273  bool const keepDiag = false,
274  real64 const diagValue = 0.0 ) override;
275 
276  virtual void addEntries( PetscMatrix const & src,
277  MatrixPatternOp const op,
278  real64 const scale = 1.0 ) override;
279 
280  virtual void addDiagonal( PetscVector const & src,
281  real64 const scale ) override;
282 
283  virtual void clampEntries( real64 const lo,
284  real64 const hi,
285  bool const excludeDiag ) override;
286 
290  virtual localIndex maxRowLengthLocal() const override;
291 
292  virtual localIndex rowLength( globalIndex const globalRowIndex ) const override;
293 
294  virtual void getRowLengths( arrayView1d< localIndex > const & lengths ) const override;
295 
296  virtual void getRowLocalLengths( arrayView1d< localIndex > const & lengths ) const override;
297 
298  virtual void extractDiagonal( PetscVector & dst ) const override;
299 
300  virtual void extract( CRSMatrixView< real64, globalIndex > const & localMat ) const override;
301 
302  virtual void extract( CRSMatrixView< real64, globalIndex const > const & localMat ) const override;
303 
304  virtual void extractLocal( CRSMatrixView< real64, localIndex > const & localMat ) const override;
305 
306  virtual void getRowSums( PetscVector & dst,
307  RowSumType const rowSumType ) const override;
308 
309  virtual void getRowCopy( globalIndex globalRow,
310  arraySlice1d< globalIndex > const & colIndices,
311  arraySlice1d< real64 > const & values ) const override;
312 
316  virtual globalIndex numGlobalRows() const override;
317 
321  virtual globalIndex numGlobalCols() const override;
322 
326  virtual localIndex numLocalRows() const override;
327 
331  virtual localIndex numLocalCols() const override;
332 
336  virtual globalIndex ilower() const override;
337 
341  virtual globalIndex iupper() const override;
342 
346  virtual globalIndex jlower() const override;
347 
351  virtual globalIndex jupper() const override;
352 
356  virtual localIndex numLocalNonzeros() const override;
357 
361  virtual globalIndex numGlobalNonzeros() const override;
362 
366  virtual real64 normInf() const override;
367 
371  virtual real64 norm1() const override;
372 
376  virtual real64 normFrobenius() const override;
377 
381  virtual real64 normMax() const override;
382 
383  virtual real64 normMax( arrayView1d< globalIndex const > const & m ) const override;
384 
385  virtual localIndex getLocalRowID( globalIndex const index ) const override;
386 
387  virtual globalIndex getGlobalRowID( localIndex const index ) const override;
388 
392  virtual MPI_Comm comm() const override;
393 
394  virtual void print( std::ostream & os = std::cout ) const override;
395 
396  virtual void write( string const & filename,
397  LAIOutputFormat const format = LAIOutputFormat::MATRIX_MARKET ) const override;
398 
400 
405  const Mat & unwrapped() const;
406 
412 
413 private:
414 
416  Mat m_mat{};
417 
419  array1d< globalIndex > m_rowsToClear;
420 
422  array1d< real64 > m_diagValues;
423 
424 };
425 
426 } // namespace geos
427 
428 #endif /*GEOS_LINEARALGEBRA_INTERFACES_PETSCSPARSEMATRIX_HPP_*/
Abstract base class for linear operators.
Common base template for all matrix wrapper types.
Definition: MatrixBase.hpp:76
virtual void multiplyPtAP(Matrix const &P, Matrix &dst) const
Compute the triple product dst = P^T * this * P
Definition: MatrixBase.hpp:648
virtual void residual(Vector const &x, Vector const &b, Vector &r) const override
Compute residual r = b - A * x.
Definition: MatrixBase.hpp:553
void setDofManager(DofManager const *const dofManager)
Associate a DofManager with this matrix.
Definition: MatrixBase.hpp:159
bool assembled() const
Query matrix assembled status.
Definition: MatrixBase.hpp:110
bool closed() const
Query matrix closed status.
Definition: MatrixBase.hpp:104
virtual void createWithGlobalSize(globalIndex const globalSize, localIndex const maxEntriesPerRow, MPI_Comm const &comm)
Create a square matrix from global number of rows.
Definition: MatrixBase.hpp:206
bool insertable() const
Query matrix status.
Definition: MatrixBase.hpp:132
virtual void createWithLocalSize(localIndex const localSize, localIndex const maxEntriesPerRow, MPI_Comm const &comm)
Create a square matrix from local number of rows.
Definition: MatrixBase.hpp:188
CRSMatrix< real64, globalIndex > extract() const
Extract local part of the matrix.
Definition: MatrixBase.hpp:875
DofManager const * dofManager() const
Definition: MatrixBase.hpp:167
bool ready() const
Query matrix ready status.
Definition: MatrixBase.hpp:117
bool modifiable() const
Query matrix status.
Definition: MatrixBase.hpp:125
virtual void multiplyP1tAP2(Matrix const &P1, Matrix const &P2, Matrix &dst) const
Compute the triple product dst = P1^T * this * P2
Definition: MatrixBase.hpp:628
CRSMatrix< real64, localIndex > extractLocal() const
Extract block of the matrix corresponding to locally owned columns.
Definition: MatrixBase.hpp:898
virtual void create(CRSMatrixView< real64 const, globalIndex const > const &localMatrix, localIndex const numLocalColumns, MPI_Comm const &comm)
Create parallel matrix from a local CRS matrix.
Definition: MatrixBase.hpp:250
Facilitates exporting Petsc matrix and associated vector objects (either in parallel or serial).
Definition: PetscExport.hpp:44
This class creates and provides basic support for the Mat matrix object type used in PETSc.
Definition: PetscMatrix.hpp:51
virtual localIndex numLocalCols() const override
virtual void set(globalIndex const *rowIndices, globalIndex const *colIndices, real64 const *values, localIndex const numRows, localIndex const numCols) override
Set a dense block of values.
virtual globalIndex numGlobalNonzeros() const override
virtual localIndex numLocalRows() const override
PetscMatrix & operator=(PetscMatrix const &src)
Copy assignment.
virtual globalIndex numGlobalRows() const override
struct _p_Mat * Mat
Alias for PETSc matrix struct pointer.
Definition: PetscMatrix.hpp:61
virtual void clampEntries(real64 const lo, real64 const hi, bool const excludeDiag) override
Clamp each matrix value between values of lo and hi.
PetscMatrix(PetscMatrix const &src)
Copy constructor.
virtual void getRowLocalLengths(arrayView1d< localIndex > const &lengths) const override
Get the local row lengths of every local row (number of nonzeros in local columns)
virtual real64 normMax(arrayView1d< globalIndex const > const &m) const override
Returns the max norm of the matrix on a subset of rows.
virtual void computeScalingVector(PetscVector &scaling) const override
Compute left and right scaling vectors for diagonal scaling.
virtual void transpose(PetscMatrix &dst) const override
Matrix transposition.
virtual void set(globalIndex const rowIndex, globalIndex const *colIndices, real64 const *values, localIndex const size) override
Set elements to one row using c-style arrays.
virtual void insert(globalIndex const rowIndex, arraySlice1d< globalIndex const > const &colIndices, arraySlice1d< real64 const > const &values) override
Insert elements of one row using array1d.
virtual void zero() override
Set all elements to zero.
virtual globalIndex numGlobalCols() const override
Mat & unwrapped()
Returns a non-const pointer to the underlying matrix.
virtual void scale(real64 const scalingFactor) override
Multiply all elements by scalingFactor.
virtual void close() override
Assemble and compress the matrix.
PetscMatrix()
Empty matrix constructor.
virtual void leftMultiplyTranspose(PetscMatrix const &src, PetscMatrix &dst) const override
Matrix/Matrix transpose multiplication.
virtual globalIndex jupper() const override
Returns index one past the last global col owned by that processor.
virtual void insert(globalIndex const rowIndex, globalIndex const *colIndices, real64 const *values, localIndex const size) override
Insert elements to one row using c-style arrays.
virtual void getRowSums(PetscVector &dst, RowSumType const rowSumType) const override
Populate a vector with row sums of this.
virtual void rescaleRows(arrayView1d< globalIndex const > const &rowIndices, RowSumType const rowSumType) override
Rescales selected rows of matrix using row sum reciprocal as a factor.
virtual void add(globalIndex const *rowIndices, globalIndex const *colIndices, real64 const *values, localIndex const numRows, localIndex const numCols) override
Add a dense block of values.
virtual void apply(PetscVector const &src, PetscVector &dst) const override
Apply operator to a vector, dst = this(src).
virtual void add(globalIndex const rowIndex, arraySlice1d< globalIndex const > const &colIndices, arraySlice1d< real64 const > const &values) override
Add elements to one row using array1d.
virtual void set(globalIndex const rowIndex, arraySlice1d< globalIndex const > const &colIndices, arraySlice1d< real64 const > const &values) override
Set elements of one row using array1d.
virtual void open() override
Open matrix for adding new entries.
virtual void getRowLengths(arrayView1d< localIndex > const &lengths) const override
Get the row lengths of every local row.
virtual void set(arraySlice1d< globalIndex const > const &rowIndices, arraySlice1d< globalIndex const > const &colIndices, arraySlice2d< real64 const > const &values) override
Set a dense block of values.
virtual void leftScale(PetscVector const &vec) override
Pre-multiplies (left) with diagonal matrix consisting of the values in vec.
PetscMatrix(PetscMatrix &&src) noexcept
Move constructor.
virtual void createWithLocalSize(localIndex const localRows, localIndex const localCols, localIndex const maxEntriesPerRow, MPI_Comm const &comm) override
Create a rectangular matrix from number of rows/columns.
virtual void print(std::ostream &os=std::cout) const override
Print the matrix in Trilinos format to a stream.
virtual void separateComponentFilter(PetscMatrix &dst, integer const dofsPerNode) const override
Apply a separate component approximation (filter) to this matrix.
virtual void rightScale(PetscVector const &vec) override
Post-multiplies (right) with diagonal matrix consisting of the values in vec.
virtual MPI_Comm comm() const override
Get the MPI communicator the matrix was created with.
virtual void extract(CRSMatrixView< real64, globalIndex > const &localMat) const override
Extract local part of the matrix using previously allocated storage.
virtual real64 clearRow(globalIndex const row, bool const keepDiag=false, real64 const diagValue=0.0) override
Clear a row, and optionally set diagonal element to diagValue.
virtual void multiplyPtAP(PetscMatrix const &P, PetscMatrix &dst) const override
Compute the triple product dst = P^T * this * P
virtual localIndex maxRowLengthLocal() const override
Returns the number of nonzero entries in the longest row of the matrix.
virtual real64 normFrobenius() const override
Returns the Frobenius norm of the matrix.
virtual void leftRightScale(PetscVector const &vecLeft, PetscVector const &vecRight) override
Post-multiplies (right) with diagonal matrix consisting of the values in vecRight and pre-multiplies ...
virtual void insert(globalIndex const *rowIndices, globalIndex const *colIndices, real64 const *values, localIndex const numRows, localIndex const numCols) override
Insert dense matrix.
virtual void insert(arrayView1d< globalIndex const > const &rowIndices, arrayView1d< globalIndex const > const &colIndices, arrayView1d< real64 const > const &values) override
Insert values stored in 3 linear vectors.
virtual localIndex rowLength(globalIndex const globalRowIndex) const override
Get row length via global row index.
virtual void multiplyRAP(PetscMatrix const &R, PetscMatrix const &P, PetscMatrix &dst) const override
Compute the triple product dst = R * this * P
~PetscMatrix() override
Destructor.
virtual void add(globalIndex const rowIndex, globalIndex const colIndex, real64 const value) override
Add to one element.
virtual globalIndex ilower() const override
Returns the index of the first global row owned by that processor.
virtual real64 normMax() const override
Returns the max norm of the matrix (the largest absolute element value).
virtual void getRowCopy(globalIndex globalRow, arraySlice1d< globalIndex > const &colIndices, arraySlice1d< real64 > const &values) const override
Returns a copy of the data in row globalRow.
virtual globalIndex jlower() const override
Returns the index of the first global col owned by that processor.
const Mat & unwrapped() const
Returns a const pointer to the underlying matrix.
virtual void applyTranspose(Vector const &src, Vector &dst) const override
Apply transpose of the matrix to a vector.
virtual void insert(globalIndex const rowIndex, globalIndex const colIndex, real64 const value) override
Insert one element.
virtual void add(arraySlice1d< globalIndex const > const &rowIndices, arraySlice1d< globalIndex const > const &colIndices, arraySlice2d< real64 const > const &values) override
Add a dense block of values.
virtual real64 normInf() const override
Returns the infinity norm of the matrix.
virtual void createWithGlobalSize(globalIndex const globalRows, globalIndex const globalCols, localIndex const maxEntriesPerRow, MPI_Comm const &comm) override
Create a rectangular matrix from number of rows/columns.
virtual globalIndex iupper() const override
Returns index one past the last global row owned by that processor.
virtual bool created() const override
virtual void insert(arraySlice1d< globalIndex const > const &rowIndices, arraySlice1d< globalIndex const > const &colIndices, arraySlice2d< real64 const > const &values) override
Insert a dense block of values.
virtual real64 norm1() const override
Returns the one norm of the matrix.
virtual void gemv(real64 const alpha, PetscVector const &x, real64 const beta, PetscVector &y, bool useTranspose=false) const override
Compute gemv y = alpha*A*x + beta*y.
virtual void extractLocal(CRSMatrixView< real64, localIndex > const &localMat) const override
Extract block of the matrix corresponding to locally owned columns using previously allocated storage...
virtual void rightMultiplyTranspose(PetscMatrix const &src, PetscMatrix &dst) const override
Matrix/Matrix transpose multiplication.
virtual localIndex getLocalRowID(globalIndex const index) const override
Map a global row index to local row index.
virtual void extract(CRSMatrixView< real64, globalIndex const > const &localMat) const override
Extract local part of the matrix using previously allocated storage and structure.
virtual void addEntries(PetscMatrix const &src, MatrixPatternOp const op, real64 const scale=1.0) override
Add entries of another matrix to this.
virtual localIndex numLocalNonzeros() const override
virtual void set(globalIndex const rowIndex, globalIndex const colIndex, real64 const value) override
Set one element.
virtual void add(globalIndex const rowIndex, globalIndex const *colIndices, real64 const *values, localIndex const size) override
Add elements to one row using c-style arrays.
virtual globalIndex getGlobalRowID(localIndex const index) const override
Map a local row index to global row index.
virtual void write(string const &filename, LAIOutputFormat const format=LAIOutputFormat::MATRIX_MARKET) const override
Write the matrix to filename in a matlab-compatible format.
PetscMatrix & operator=(PetscMatrix &&src) noexcept
Move assignment.
virtual void multiply(PetscMatrix const &src, PetscMatrix &dst) const override
Matrix/Matrix multiplication.
virtual void addDiagonal(PetscVector const &src, real64 const scale) override
Add (scaled) entries of a vector to the diagonal of this matrix.
virtual void reset() override
Reset the matrix to default state.
virtual void set(real64 const value) override
Set all non-zero elements to a value.
virtual void extractDiagonal(PetscVector &dst) const override
Extract diagonal values into a vector.
This class creates and provides basic support for Vec vector object type used in PETSc.
Definition: PetscVector.hpp:46
ArrayView< T, 1 > arrayView1d
Alias for 1D array view.
Definition: DataTypes.hpp:179
RowSumType
Type of row sum to compute.
Definition: MatrixBase.hpp:36
GEOS_GLOBALINDEX_TYPE globalIndex
Global index type (for indexing objects across MPI partitions).
Definition: DataTypes.hpp:87
LAIOutputFormat
Definition: common.hpp:142
ArraySlice< T, 2, USD > arraySlice2d
Alias for 2D array slice.
Definition: DataTypes.hpp:199
double real64
64-bit floating point type.
Definition: DataTypes.hpp:98
GEOS_LOCALINDEX_TYPE localIndex
Local index type (for indexing objects within an MPI partition).
Definition: DataTypes.hpp:84
ArraySlice< T, 1, USD > arraySlice1d
Alias for 1D array slice.
Definition: DataTypes.hpp:183
std::int32_t integer
Signed integer type.
Definition: DataTypes.hpp:81
LvArray::CRSMatrixView< T, COL_INDEX, INDEX_TYPE const, LvArray::ChaiBuffer > CRSMatrixView
Alias for CRS Matrix View.
Definition: DataTypes.hpp:309
Array< T, 1 > array1d
Alias for 1D array.
Definition: DataTypes.hpp:175
MatrixPatternOp
Describes relationship between and treatment of nonzero patterns of arguments in matrix functions lik...
Definition: MatrixBase.hpp:48