Epetra_BasicRowMatrix: A class for simplifying the development of Epetra_RowMatrix adapters. More...
#include <Epetra_BasicRowMatrix.h>
Public Member Functions | |
Constructor/Destructor | |
Epetra_BasicRowMatrix (const Epetra_Comm &Comm) | |
Epetra_BasicRowMatrix constructor. More... | |
virtual | ~Epetra_BasicRowMatrix () |
Epetra_BasicRowMatrix Destructor. | |
Setup functions | |
void | SetMaps (const Epetra_Map &RowMap, const Epetra_Map &ColMap) |
Set maps (Version 1); call this function or the next, but not both. More... | |
void | SetMaps (const Epetra_Map &RowMap, const Epetra_Map &ColMap, const Epetra_Map &DomainMap, const Epetra_Map &RangeMap) |
Set maps (Version 2); call this function or the previous, but not both. More... | |
User-required implementation methods | |
virtual int | ExtractMyRowCopy (int MyRow, int Length, int &NumEntries, double *Values, int *Indices) const =0 |
Returns a copy of the specified local row in user-provided arrays. More... | |
virtual int | ExtractMyEntryView (int CurEntry, double *&Value, int &RowIndex, int &ColIndex)=0 |
Returns a reference to the ith entry in the matrix, along with its row and column index. More... | |
virtual int | ExtractMyEntryView (int CurEntry, double const *&Value, int &RowIndex, int &ColIndex) const =0 |
Returns a const reference to the ith entry in the matrix, along with its row and column index. More... | |
virtual int | NumMyRowEntries (int MyRow, int &NumEntries) const =0 |
Return the current number of values stored for the specified local row. More... | |
Computational methods | |
virtual int | Multiply (bool TransA, const Epetra_MultiVector &X, Epetra_MultiVector &Y) const |
Returns the result of a Epetra_BasicRowMatrix multiplied by a Epetra_MultiVector X in Y. More... | |
virtual int | Solve (bool Upper, bool Trans, bool UnitDiagonal, const Epetra_MultiVector &X, Epetra_MultiVector &Y) const |
Returns the result of a Epetra_BasicRowMatrix solve with a Epetra_MultiVector X in Y (not implemented). More... | |
virtual int | ExtractDiagonalCopy (Epetra_Vector &Diagonal) const |
Returns a copy of the main diagonal in a user-provided vector. More... | |
virtual int | InvRowSums (Epetra_Vector &x) const |
Computes the sum of absolute values of the rows of the Epetra_BasicRowMatrix, results returned in x. More... | |
virtual int | LeftScale (const Epetra_Vector &x) |
Scales the Epetra_BasicRowMatrix on the left with a Epetra_Vector x. More... | |
virtual int | InvColSums (Epetra_Vector &x) const |
Computes the sum of absolute values of the columns of the Epetra_BasicRowMatrix, results returned in x. More... | |
virtual int | RightScale (const Epetra_Vector &x) |
Scales the Epetra_BasicRowMatrix on the right with a Epetra_Vector x. More... | |
Matrix Properties Query Methods | |
virtual bool | Filled () const |
If FillComplete() has been called, this query returns true, otherwise it returns false, presently always returns true. | |
bool | LowerTriangular () const |
If matrix is lower triangular, this query returns true, otherwise it returns false. | |
virtual bool | UpperTriangular () const |
If matrix is upper triangular, this query returns true, otherwise it returns false. | |
Attribute access functions | |
virtual double | NormInf () const |
Returns the infinity norm of the global matrix. More... | |
virtual double | NormOne () const |
Returns the one norm of the global matrix. More... | |
virtual int | NumGlobalNonzeros () const |
Returns the number of nonzero entries in the global matrix. More... | |
virtual long long | NumGlobalNonzeros64 () const |
virtual int | NumGlobalRows () const |
Returns the number of global matrix rows. | |
virtual long long | NumGlobalRows64 () const |
virtual int | NumGlobalCols () const |
Returns the number of global matrix columns. | |
virtual long long | NumGlobalCols64 () const |
virtual int | NumGlobalDiagonals () const |
Returns the number of global nonzero diagonal entries. | |
virtual long long | NumGlobalDiagonals64 () const |
virtual int | NumMyNonzeros () const |
Returns the number of nonzero entries in the calling processor's portion of the matrix. | |
virtual int | NumMyRows () const |
Returns the number of matrix rows owned by the calling processor. | |
virtual int | NumMyCols () const |
Returns the number of matrix columns owned by the calling processor. | |
virtual int | NumMyDiagonals () const |
Returns the number of local nonzero diagonal entries. | |
virtual int | MaxNumEntries () const |
Returns the maximum number of nonzero entries across all rows on this processor. | |
virtual const Epetra_Map & | OperatorDomainMap () const |
Returns the Epetra_Map object associated with the domain of this operator. | |
virtual const Epetra_Map & | OperatorRangeMap () const |
Returns the Epetra_Map object associated with the range of this operator (same as domain). | |
virtual const Epetra_BlockMap & | Map () const |
Implement the Epetra_SrcDistObjec::Map() function. | |
virtual const Epetra_Map & | RowMatrixRowMap () const |
Returns the Row Map object needed for implementing Epetra_RowMatrix. | |
virtual const Epetra_Map & | RowMatrixColMap () const |
Returns the Column Map object needed for implementing Epetra_RowMatrix. | |
virtual const Epetra_Import * | RowMatrixImporter () const |
Returns the Epetra_Import object that contains the import operations for distributed operations. | |
virtual const Epetra_Comm & | Comm () const |
Returns a pointer to the Epetra_Comm communicator associated with this matrix. | |
I/O Methods | |
virtual void | Print (std::ostream &os) const |
Print method. | |
Additional methods required to support the Epetra_RowMatrix interface | |
virtual int | SetUseTranspose (bool use_transpose) |
If set true, transpose of this operator will be applied. More... | |
virtual const char * | Label () const |
Returns a character string describing the operator. | |
virtual int | Apply (const Epetra_MultiVector &X, Epetra_MultiVector &Y) const |
Returns the result of a Epetra_RowMatrix applied to a Epetra_MultiVector X in Y. More... | |
virtual int | ApplyInverse (const Epetra_MultiVector &X, Epetra_MultiVector &Y) const |
Returns the result of a Epetra_RowMatrix inverse applied to an Epetra_MultiVector X in Y. More... | |
bool | HasNormInf () const |
Returns true because this class can compute an Inf-norm. | |
virtual bool | UseTranspose () const |
Returns the current UseTranspose setting. | |
Additional accessor methods | |
virtual const Epetra_Import * | Importer () const |
Returns the Epetra_Import object that contains the import operations for distributed operations, returns zero if none. More... | |
virtual const Epetra_Export * | Exporter () const |
Returns the Epetra_Export object that contains the export operations for distributed operations, returns zero if none. More... | |
Public Member Functions inherited from Epetra_CompObject | |
Epetra_CompObject & | operator= (const Epetra_CompObject &src) |
Epetra_CompObject () | |
Basic Epetra_CompObject constuctor. | |
Epetra_CompObject (const Epetra_CompObject &Source) | |
Epetra_CompObject copy constructor. | |
virtual | ~Epetra_CompObject () |
Epetra_CompObject destructor. | |
void | SetFlopCounter (const Epetra_Flops &FlopCounter_in) |
Set the internal Epetra_Flops() pointer. | |
void | SetFlopCounter (const Epetra_CompObject &CompObject) |
Set the internal Epetra_Flops() pointer to the flop counter of another Epetra_CompObject. | |
void | UnsetFlopCounter () |
Set the internal Epetra_Flops() pointer to 0 (no flops counted). | |
Epetra_Flops * | GetFlopCounter () const |
Get the pointer to the Epetra_Flops() object associated with this object, returns 0 if none. | |
void | ResetFlops () const |
Resets the number of floating point operations to zero for this multi-vector. | |
double | Flops () const |
Returns the number of floating point operations with this multi-vector. | |
void | UpdateFlops (int Flops_in) const |
Increment Flop count for this object. | |
void | UpdateFlops (long int Flops_in) const |
Increment Flop count for this object. | |
void | UpdateFlops (long long Flops_in) const |
Increment Flop count for this object. | |
void | UpdateFlops (double Flops_in) const |
Increment Flop count for this object. | |
void | UpdateFlops (float Flops_in) const |
Increment Flop count for this object. | |
Public Member Functions inherited from Epetra_Object | |
Epetra_Object (int TracebackModeIn=-1, bool set_label=true) | |
Epetra_Object Constructor. More... | |
Epetra_Object (const char *const Label, int TracebackModeIn=-1) | |
Epetra_Object Constructor. More... | |
Epetra_Object (const Epetra_Object &Object) | |
Epetra_Object Copy Constructor. More... | |
virtual | ~Epetra_Object () |
Epetra_Object Destructor. More... | |
virtual int | ReportError (const std::string Message, int ErrorCode) const |
Error reporting method. | |
virtual void | SetLabel (const char *const Label) |
Epetra_Object Label definition using char *. More... | |
Public Member Functions inherited from Epetra_RowMatrix | |
virtual | ~Epetra_RowMatrix () |
Destructor. | |
Public Member Functions inherited from Epetra_Operator | |
virtual | ~Epetra_Operator () |
Destructor. | |
Public Member Functions inherited from Epetra_SrcDistObject | |
virtual | ~Epetra_SrcDistObject () |
Epetra_SrcDistObject destructor. | |
Protected Member Functions | |
void | Setup () |
void | UpdateImportVector (int NumVectors) const |
void | UpdateExportVector (int NumVectors) const |
void | SetImportExport () |
Post-construction modifications | |
virtual void | ComputeStructureConstants () const |
Update the constants associated with the structure of the matrix: Call only if structure changes from the initial RowMatrix. | |
virtual void | ComputeNumericConstants () const |
Update the constants associated with the values of the matrix: Call only if values changes from the initial RowMatrix. | |
Protected Member Functions inherited from Epetra_Object | |
std::string | toString (const int &x) const |
std::string | toString (const long long &x) const |
std::string | toString (const double &x) const |
Protected Attributes | |
Epetra_Comm * | Comm_ |
Epetra_Map | OperatorDomainMap_ |
Epetra_Map | OperatorRangeMap_ |
Epetra_Map | RowMatrixRowMap_ |
Epetra_Map | RowMatrixColMap_ |
int | NumMyNonzeros_ |
long long | NumGlobalNonzeros_ |
int | MaxNumEntries_ |
double | NormInf_ |
double | NormOne_ |
int | NumMyRows_ |
int | NumMyCols_ |
bool | UseTranspose_ |
bool | HasNormInf_ |
bool | LowerTriangular_ |
bool | UpperTriangular_ |
bool | HaveStructureConstants_ |
bool | HaveNumericConstants_ |
bool | HaveMaps_ |
Epetra_MultiVector * | ImportVector_ |
Epetra_MultiVector * | ExportVector_ |
Epetra_Import * | Importer_ |
Epetra_Export * | Exporter_ |
Protected Attributes inherited from Epetra_CompObject | |
Epetra_Flops * | FlopCounter_ |
Additional Inherited Members | |
Static Public Member Functions inherited from Epetra_Object | |
static void | SetTracebackMode (int TracebackModeValue) |
Set the value of the Epetra_Object error traceback report mode. More... | |
static int | GetTracebackMode () |
Get the value of the Epetra_Object error report mode. | |
static std::ostream & | GetTracebackStream () |
Get the output stream for error reporting. | |
Static Public Attributes inherited from Epetra_Object | |
static int | TracebackMode |
Epetra_BasicRowMatrix: A class for simplifying the development of Epetra_RowMatrix adapters.
The Epetra_BasicRowMatrix is an adapter class for Epetra_RowMatrix that implements most of the Epetra_RowMatrix methods using reasonable default implementations. The Epetra_RowMatrix class has 39 pure virtual methods, requiring the adapter class to implement all of them. Epetra_BasicRowMatrix has only 4 pure virtual methods that must be implemented (See Epetra_JadMatrix for an example):
An alternative is possible if you do not want to provide a non-trivial implementation of the ExtraMyEntryView methods (See Epetra_VbrRowMatrix for an example):
In addition, most adapters will probably re-implement the Multiply() method and perhaps the Solve() method, although one or the other may be implemented to return -1, signaling that there is no valid implementation. By default, the Multiply() method is implemented using ExtractMyRowCopy, which can usual be improved upon. By default Solve() and ApplyInverse() are implemented to return -1 (not implemented).
All other implemented methods in Epetra_BasicRowMatrix should not exhibit a signficant performance degradation, either because they are relatively small and fast, or because they are not a significant portion of the runtime for most codes. All methods are virtual, so they can be re-implemented by the adapter.
In addition to implementing the above methods, an adapter must inherit the Epetra_BasicRowMatrix interface and call the Epetra_BasicRowMatrix constructor as part of the adapter constructor. There are two constructors. The first requires the user to pass in the RowMap and ColMap, both of which are Epetra_Map objects. On each processor the RowMap (ColMap) must contain the global IDs (GIDs) of the rows (columns) that the processor cares about. The first constructor requires only these two maps, assuming that the RowMap will also serve as the DomainMap and RangeMap. In this case, the RowMap must be 1-to-1, meaning that if a global ID appears on one processor, it appears only once on that processor and does not appear on any other processor. For many sparse matrix data structures, it is the case that a given row is completely owned by one processor and that the global matrix is square. The first constructor is for this situation.
The second constructor allows the caller to specify all four maps. In this case the DomainMap, the layout of multivectors/vectors that are in the domain of the matrix (the x vector if computing y = A*x), must be 1-to-1. Also, the RangeMap, the layout of y must be 1-to-1. The RowMap and ColMap do not need to be 1-to-1, but the GIDs must be found in the RangeMap and DomainMap, respectively.
Note that Epetra_Operator is a base class for Epetra_RowMatrix, so any adapter for Epetra_BasicRowMatrix (or Epetra_RowMatrix) is also an adapter for Epetra_Operator.
An example of how to provide an adapter for Epetra_BasicRowMatrix can be found by looking at Epetra_JadMatrix.
Epetra_BasicRowMatrix::Epetra_BasicRowMatrix | ( | const Epetra_Comm & | Comm | ) |
Epetra_BasicRowMatrix constructor.
This constructor requires a valid Epetra_Comm object as its only argument. The constructor will use Comm to build Epetra_Maps objects: RowMap, ColMap, DomainMap and RangeMap. However, these will be zero-length (trivial) maps that must be reset by calling one of the two SetMap() methods listed below.
[in] | Comm | An Epetra_Comm containing a valid Comm object. |
|
inlinevirtual |
Returns the result of a Epetra_RowMatrix applied to a Epetra_MultiVector X in Y.
X | (In) - A Epetra_MultiVector of dimension NumVectors to multiply with matrix. |
Y | (Out) - A Epetra_MultiVector of dimension NumVectors containing result. |
Implements Epetra_Operator.
References Epetra_RowMatrix::Multiply(), and UseTranspose().
|
inlinevirtual |
Returns the result of a Epetra_RowMatrix inverse applied to an Epetra_MultiVector X in Y.
X | (In) - A Epetra_MultiVector of dimension NumVectors to solve for. |
Y | (Out) - A Epetra_MultiVector of dimension NumVectors containing result. |
Implements Epetra_Operator.
|
inlinevirtual |
Returns the Epetra_Export object that contains the export operations for distributed operations, returns zero if none.
If RowMatrixRowMap!=OperatorRangeMap, then this method returns a pointer to an Epetra_Export object that exports objects from an RowMatrixRowMap layout to a OperatorRangeMap layout. This operation is needed for sparse matrix-vector
multiplication, y = Ax, to scatter-add y elements generated during local multiplication operations.
If RowMatrixRowMap==OperatorRangeMap, then the pointer will be returned as 0. For a typical Epetra_RowMatrix object, this pointer will be zero since it is often the case that RowMatrixRowMap==OperatorRangeMap.
|
virtual |
Returns a copy of the main diagonal in a user-provided vector.
Diagonal | (Out) - Extracted main diagonal. |
Implements Epetra_RowMatrix.
|
pure virtual |
Returns a reference to the ith entry in the matrix, along with its row and column index.
CurEntry | (In) - Index of local entry (from 0 to NumMyNonzeros()-1) to extract. |
Value | (Out) - Extracted reference to current values. |
RowIndex | (Out) - Row index for current entry. |
ColIndex | (Out) - Column index for current entry. |
Implemented in Epetra_VbrRowMatrix, and Epetra_JadMatrix.
|
pure virtual |
Returns a const reference to the ith entry in the matrix, along with its row and column index.
CurEntry | (In) - Index of local entry (from 0 to NumMyNonzeros()-1) to extract. |
Value | (Out) - Extracted reference to current values. |
RowIndex | (Out) - Row index for current entry. |
ColIndex | (Out) - Column index for current entry. |
Implemented in Epetra_VbrRowMatrix, and Epetra_JadMatrix.
|
pure virtual |
Returns a copy of the specified local row in user-provided arrays.
MyRow | (In) - Local row to extract. |
Length | (In) - Length of Values and Indices. |
NumEntries | (Out) - Number of nonzero entries extracted. |
Values | (Out) - Extracted values for this row. |
Indices | (Out) - Extracted global column indices for the corresponding values. |
Implements Epetra_RowMatrix.
Implemented in Epetra_VbrRowMatrix, and Epetra_JadMatrix.
|
inlinevirtual |
Returns the Epetra_Import object that contains the import operations for distributed operations, returns zero if none.
If RowMatrixColMap!=OperatorDomainMap, then this method returns a pointer to an Epetra_Import object that imports objects from an OperatorDomainMap layout to a RowMatrixColMap layout. This operation is needed for sparse matrix-vector
multiplication, y = Ax, to gather x elements for local multiplication operations.
If RowMatrixColMap==OperatorDomainMap, then the pointer will be returned as 0.
|
virtual |
Computes the sum of absolute values of the columns of the Epetra_BasicRowMatrix, results returned in x.
The vector x will return such that x[j] will contain the inverse of sum of the absolute values of the \e this matrix will be sca such that A(i,j) = x(j)*A(i,j) where i denotes the global row number of A and j denotes the global column number of A. Using the resulting vector from this function as input to
RighttScale() will make the one norm of the resulting matrix exactly 1.
x | (Out) - An Epetra_Vector containing the column sums of the this matrix. |
Implements Epetra_RowMatrix.
|
virtual |
Computes the sum of absolute values of the rows of the Epetra_BasicRowMatrix, results returned in x.
The vector x will return such that x[i] will contain the inverse of sum of the absolute values of the \e this matrix will be scaled such that A(i,j) = x(i)*A(i,j) where i denotes the global row number of A and j denotes the global column number of A. Using the resulting vector from this function as input to LeftScale()
will make the infinity norm of the resulting matrix exactly 1.
x | (Out) - An Epetra_Vector containing the row sums of the this matrix. |
Implements Epetra_RowMatrix.
|
virtual |
Scales the Epetra_BasicRowMatrix on the left with a Epetra_Vector x.
The \e this matrix will be scaled such that A(i,j) = x(i)*A(i,j) where i denotes the row number of A and j denotes the column number of A.
x | (In) - An Epetra_Vector to solve for. |
Implements Epetra_RowMatrix.
Reimplemented in Epetra_VbrRowMatrix.
|
virtual |
Returns the result of a Epetra_BasicRowMatrix multiplied by a Epetra_MultiVector X in Y.
TransA | (In) - If true, multiply by the transpose of matrix, otherwise just use matrix. |
X | (Out) - An Epetra_MultiVector of dimension NumVectors to multiply with matrix. |
Y | (Out) - An Epetra_MultiVector of dimension NumVectorscontaining result. |
Implements Epetra_RowMatrix.
Reimplemented in Epetra_VbrRowMatrix, and Epetra_JadMatrix.
|
inlinevirtual |
Returns the infinity norm of the global matrix.
Returns the quantity such that
.
Implements Epetra_RowMatrix.
|
inlinevirtual |
Returns the one norm of the global matrix.
Returns the quantity such that
.
Implements Epetra_RowMatrix.
|
inlinevirtual |
Returns the number of nonzero entries in the global matrix.
Note that if the data decomposition is defined such that some nonzeros appear on multiple processors, then those nonzeros will be counted multiple times.
Implements Epetra_RowMatrix.
References Epetra_Operator::OperatorDomainMap(), and Epetra_Operator::OperatorRangeMap().
Referenced by Epetra_VbrRowMatrix::LeftScale(), and Epetra_VbrRowMatrix::RightScale().
|
pure virtual |
Return the current number of values stored for the specified local row.
Similar to NumMyEntries() except NumEntries is returned as an argument and error checking is done on the input value MyRow.
MyRow | (In) - Local row. |
NumEntries | (Out) - Number of nonzero values. |
Implements Epetra_RowMatrix.
Implemented in Epetra_VbrRowMatrix, and Epetra_JadMatrix.
|
virtual |
Scales the Epetra_BasicRowMatrix on the right with a Epetra_Vector x.
The \e this matrix will be scaled such that A(i,j) = x(j)*A(i,j) where i denotes the global row number of A and j denotes the global column number of A.
x | (In) - The Epetra_Vector used for scaling this. |
Implements Epetra_RowMatrix.
Reimplemented in Epetra_VbrRowMatrix.
void Epetra_BasicRowMatrix::SetMaps | ( | const Epetra_Map & | RowMap, |
const Epetra_Map & | ColMap | ||
) |
Set maps (Version 1); call this function or the next, but not both.
This method takes a row and column map. On each processor these maps describe the global rows and columns, resp, that the processor will care about. Note that the ColMap does not have to be one-to-one. In other words, a column ID can appear on more than one processor. The RowMap must be 1-to-1.
[in] | RowMap | An Epetra_Map containing on each processor a list of GIDs of rows that the processor cares about. |
[in] | ColMap | An Epetra_Map containing on each processor a list of GIDs of columns that the processor cares about. |
In this method, the domain and range maps are assumed to be the same as the row map. Note that this requires that the global matrix be square. If the matrix is not square, or the domain vectors or range vectors do not have the same layout as the rows, then the second constructor should be called.
Referenced by Epetra_VbrRowMatrix::Epetra_VbrRowMatrix().
void Epetra_BasicRowMatrix::SetMaps | ( | const Epetra_Map & | RowMap, |
const Epetra_Map & | ColMap, | ||
const Epetra_Map & | DomainMap, | ||
const Epetra_Map & | RangeMap | ||
) |
Set maps (Version 2); call this function or the previous, but not both.
This constructor takes a row, column, domain and range map. On each processor these maps describe the global rows, columns, domain and range, resp, that the processor will care about. The domain and range maps must be one-to-one, but note that the row and column maps do not have to be one-to-one. In other words, a row ID can appear on more than one processor, as can a column ID.
[in] | RowMap | An Epetra_Map containing on each processor a list of GIDs of rows that the processor cares about. |
[in] | ColMap | An Epetra_Map containing on each processor a list of GIDs of columns that the processor cares about. |
[in] | DomainMap | An Epetra_Map describing the distribution of domain vectors and multivectors. |
[in] | RangeMap | An Epetra_Map describing the distribution of range vectors and multivectors. |
|
inlinevirtual |
If set true, transpose of this operator will be applied.
This flag allows the transpose of the given operator to be used implicitly. Setting this flag affects only the Apply() and ApplyInverse() methods. If the implementation of this interface
does not support transpose use, this method should return a value of -1.
use_transpose | (In) - If true, multiply by the transpose of operator, otherwise just use operator. |
Implements Epetra_Operator.
|
inlinevirtual |
Returns the result of a Epetra_BasicRowMatrix solve with a Epetra_MultiVector X in Y (not implemented).
Upper | (In) - If true, solve Ux = y, otherwise solve Lx = y. |
Trans | (In) - If true, solve transpose problem. |
UnitDiagonal | (In) - If true, assume diagonal is unit (whether it's stored or not). |
X | (In) - An Epetra_MultiVector of dimension NumVectors to solve for. |
Y | (Out) - An Epetra_MultiVector of dimension NumVectors containing result. |
Implements Epetra_RowMatrix.
Reimplemented in Epetra_VbrRowMatrix, and Epetra_JadMatrix.