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 
123  virtual void createWithLocalSize( localIndex const localRows,
124  localIndex const localCols,
125  localIndex const maxEntriesPerRow,
126  MPI_Comm const & comm ) override;
127 
128  virtual void createWithGlobalSize( globalIndex const globalRows,
129  globalIndex const globalCols,
130  localIndex const maxEntriesPerRow,
131  MPI_Comm const & comm ) override;
132 
136  virtual bool created() const override;
137 
138  virtual void reset() override;
139 
140  virtual void set( real64 const value ) override;
141 
142  virtual void zero() override;
143 
144  virtual void open() override;
145 
146  virtual void close() override;
147 
148  virtual void add( globalIndex const rowIndex,
149  globalIndex const colIndex,
150  real64 const value ) override;
151 
152  virtual void set( globalIndex const rowIndex,
153  globalIndex const colIndex,
154  real64 const value ) override;
155 
156  virtual void insert( globalIndex const rowIndex,
157  globalIndex const colIndex,
158  real64 const value ) override;
159 
160  virtual void add( globalIndex const rowIndex,
161  globalIndex const * colIndices,
162  real64 const * values,
163  localIndex const size ) override;
164 
165  virtual void set( globalIndex const rowIndex,
166  globalIndex const * colIndices,
167  real64 const * values,
168  localIndex const size ) override;
169 
170  virtual void insert( globalIndex const rowIndex,
171  globalIndex const * colIndices,
172  real64 const * values,
173  localIndex const size ) override;
174 
175  virtual void add( globalIndex const rowIndex,
176  arraySlice1d< globalIndex const > const & colIndices,
177  arraySlice1d< real64 const > const & values ) override;
178 
179  virtual void set( globalIndex const rowIndex,
180  arraySlice1d< globalIndex const > const & colIndices,
181  arraySlice1d< real64 const > const & values ) override;
182 
183  virtual void insert( globalIndex const rowIndex,
184  arraySlice1d< globalIndex const > const & colIndices,
185  arraySlice1d< real64 const > const & values ) override;
186 
187  virtual void add( arraySlice1d< globalIndex const > const & rowIndices,
188  arraySlice1d< globalIndex const > const & colIndices,
189  arraySlice2d< real64 const > const & values ) override;
190 
191  virtual void set( arraySlice1d< globalIndex const > const & rowIndices,
192  arraySlice1d< globalIndex const > const & colIndices,
193  arraySlice2d< real64 const > const & values ) override;
194 
195  virtual void insert( arraySlice1d< globalIndex const > const & rowIndices,
196  arraySlice1d< globalIndex const > const & colIndices,
197  arraySlice2d< real64 const > const & values ) override;
198 
199  virtual void add( globalIndex const * rowIndices,
200  globalIndex const * colIndices,
201  real64 const * values,
202  localIndex const numRows,
203  localIndex const numCols ) override;
204 
205  virtual void set( globalIndex const * rowIndices,
206  globalIndex const * colIndices,
207  real64 const * values,
208  localIndex const numRows,
209  localIndex const numCols ) override;
210 
211  virtual void insert( globalIndex const * rowIndices,
212  globalIndex const * colIndices,
213  real64 const * values,
214  localIndex const numRows,
215  localIndex const numCols ) override;
216 
217  virtual void insert( arrayView1d< globalIndex const > const & rowIndices,
218  arrayView1d< globalIndex const > const & colIndices,
219  arrayView1d< real64 const > const & values ) override;
220 
221  virtual void apply( PetscVector const & src,
222  PetscVector & dst ) const override;
223 
224  virtual void multiply( PetscMatrix const & src,
225  PetscMatrix & dst ) const override;
226 
227  virtual void applyTranspose( Vector const & src,
228  Vector & dst ) const override;
229 
230  virtual void leftMultiplyTranspose( PetscMatrix const & src,
231  PetscMatrix & dst ) const override;
232 
233  virtual void rightMultiplyTranspose( PetscMatrix const & src,
234  PetscMatrix & dst ) const override;
235 
236  virtual void multiplyRAP( PetscMatrix const & R,
237  PetscMatrix const & P,
238  PetscMatrix & dst ) const override;
239 
240  virtual void multiplyPtAP( PetscMatrix const & P,
241  PetscMatrix & dst ) const override;
242 
243  virtual void gemv( real64 const alpha,
244  PetscVector const & x,
245  real64 const beta,
246  PetscVector & y,
247  bool useTranspose = false ) const override;
248 
249  virtual void scale( real64 const scalingFactor ) override;
250 
251  virtual void leftScale( PetscVector const & vec ) override;
252 
253  virtual void rightScale( PetscVector const & vec ) override;
254 
255  virtual void leftRightScale( PetscVector const & vecLeft,
256  PetscVector const & vecRight ) override;
257 
258  virtual void rescaleRows( arrayView1d< globalIndex const > const & rowIndices,
259  RowSumType const rowSumType ) override;
260 
261  virtual void transpose( PetscMatrix & dst ) const override;
262 
264  integer const dofsPerNode ) const override;
265 
266  virtual real64 clearRow( globalIndex const row,
267  bool const keepDiag = false,
268  real64 const diagValue = 0.0 ) override;
269 
270  virtual void addEntries( PetscMatrix const & src,
271  MatrixPatternOp const op,
272  real64 const scale = 1.0 ) override;
273 
274  virtual void addDiagonal( PetscVector const & src,
275  real64 const scale ) override;
276 
277  virtual void clampEntries( real64 const lo,
278  real64 const hi,
279  bool const excludeDiag ) override;
280 
284  virtual localIndex maxRowLength() const override;
285 
286  virtual localIndex rowLength( globalIndex const globalRowIndex ) const override;
287 
288  virtual void getRowLengths( arrayView1d< localIndex > const & lengths ) const override;
289 
290  virtual void extractDiagonal( PetscVector & dst ) const override;
291 
292  virtual void getRowSums( PetscVector & dst,
293  RowSumType const rowSumType ) const override;
294 
295  virtual void getRowCopy( globalIndex globalRow,
296  arraySlice1d< globalIndex > const & colIndices,
297  arraySlice1d< real64 > const & values ) const override;
298 
302  virtual globalIndex numGlobalRows() const override;
303 
307  virtual globalIndex numGlobalCols() const override;
308 
312  virtual localIndex numLocalRows() const override;
313 
317  virtual localIndex numLocalCols() const override;
318 
322  virtual globalIndex ilower() const override;
323 
327  virtual globalIndex iupper() const override;
328 
332  virtual globalIndex jlower() const override;
333 
337  virtual globalIndex jupper() const override;
338 
342  virtual localIndex numLocalNonzeros() const override;
343 
347  virtual globalIndex numGlobalNonzeros() const override;
348 
352  virtual real64 normInf() const override;
353 
357  virtual real64 norm1() const override;
358 
362  virtual real64 normFrobenius() const override;
363 
367  virtual real64 normMax() const override;
368 
369  virtual real64 normMax( arrayView1d< globalIndex const > const & m ) const override;
370 
371  virtual localIndex getLocalRowID( globalIndex const index ) const override;
372 
373  virtual globalIndex getGlobalRowID( localIndex const index ) const override;
374 
378  virtual MPI_Comm comm() const override;
379 
380  virtual void print( std::ostream & os = std::cout ) const override;
381 
382  virtual void write( string const & filename,
383  LAIOutputFormat const format = LAIOutputFormat::MATRIX_MARKET ) const override;
384 
386 
391  const Mat & unwrapped() const;
392 
398 
399 private:
400 
402  Mat m_mat{};
403 
405  array1d< globalIndex > m_rowsToClear;
406 
408  array1d< real64 > m_diagValues;
409 
410 };
411 
412 } // namespace geos
413 
414 #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 residual(Vector const &x, Vector const &b, Vector &r) const override
Compute residual r = b - A * x.
Definition: MatrixBase.hpp:551
void setDofManager(DofManager const *const dofManager)
Associate a DofManager with this matrix.
Definition: MatrixBase.hpp:157
bool assembled() const
Query matrix assembled status.
Definition: MatrixBase.hpp:108
bool closed() const
Query matrix closed status.
Definition: MatrixBase.hpp:102
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:204
bool insertable() const
Query matrix status.
Definition: MatrixBase.hpp:130
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:186
DofManager const * dofManager() const
Definition: MatrixBase.hpp:165
bool ready() const
Query matrix ready status.
Definition: MatrixBase.hpp:115
bool modifiable() const
Query matrix status.
Definition: MatrixBase.hpp:123
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:248
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
Get the number of local columns.
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
Returns the total number of nonzeros in the matrix.
virtual localIndex numLocalRows() const override
Get the number of local rows.
PetscMatrix & operator=(PetscMatrix const &src)
Copy assignment.
virtual globalIndex numGlobalRows() const override
Get the number of global rows.
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 real64 normMax(arrayView1d< globalIndex const > const &m) const override
Returns the max norm of the matrix on a subset of rows.
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
Get the number of global columns.
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 localIndex maxRowLength() const override
Returns the number of nonzero entries in the longest row of the matrix.
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 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 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
Get the number of global rows.
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 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 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
Returns the number of nonzeros in the local portion of the matrix.
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:180
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:88
LAIOutputFormat
Definition: common.hpp:142
ArraySlice< T, 2, USD > arraySlice2d
Alias for 2D array slice.
Definition: DataTypes.hpp:200
double real64
64-bit floating point type.
Definition: DataTypes.hpp:99
GEOS_LOCALINDEX_TYPE localIndex
Local index type (for indexing objects within an MPI partition).
Definition: DataTypes.hpp:85
ArraySlice< T, 1, USD > arraySlice1d
Alias for 1D array slice.
Definition: DataTypes.hpp:184
std::int32_t integer
Signed integer type.
Definition: DataTypes.hpp:82
Array< T, 1 > array1d
Alias for 1D array.
Definition: DataTypes.hpp:176
MatrixPatternOp
Describes relationship between and treatment of nonzero patterns of arguments in matrix functions lik...
Definition: MatrixBase.hpp:48