Komplex
Development
|
#include <Komplex_RowMatrix.h>
Public Member Functions | |
Komplex_RowMatrix (Epetra_DataAccess CV, const Epetra_BlockMap &RowMap, int NumEntriesPerRow, Komplex_KForms KForm=K1) | |
Default Komplex_RowMatrix constuctor with fixed number of indices per row. More... | |
Komplex_RowMatrix (Epetra_DataAccess CV, const Epetra_BlockMap &RowMap, int *NumEntriesPerRow, Komplex_KForms KForm=K1) | |
Default Komplex_RowMatrix constructor with variable number of indices per row. More... | |
Komplex_RowMatrix (Epetra_DataAccess CV, const Epetra_BlockMap &RowMap, const Epetra_BlockMap &ColMap, int NumEntriesPerRow, Komplex_KForms KForm=K1) | |
Constructor with fixed number of indices per row and both row and column maps. More... | |
Komplex_RowMatrix (Epetra_DataAccess CV, const Epetra_BlockMap &RowMap, const Epetra_BlockMap &ColMap, int *NumEntriesPerRow, Komplex_KForms KForm=K1) | |
Constructor with variable number of indices per row and both row and column maps. More... | |
Komplex_RowMatrix (Epetra_DataAccess CV, const Epetra_CrsGraph &Graph, Komplex_KForms KForm=K1) | |
Constructor that uses an existing Epetra_CrsGraph object. More... | |
Komplex_RowMatrix (double cr, double ci, Epetra_RowMatrix *A, Komplex_KForms KForm=K1) | |
General Komplex_RowMatrix constructor taking one Epetra_RowMatrix object. More... | |
Komplex_RowMatrix (double c0r, double c0i, Epetra_RowMatrix *A0, double c1r, double c1i, Epetra_RowMatrix *A1, Komplex_KForms KForm=K1) | |
General Komplex_RowMatrix constructor. More... | |
Komplex_RowMatrix (const Komplex_RowMatrix &Matrix) | |
virtual | ~Komplex_RowMatrix () |
Komplex_RowMatrix Destructor. | |
Komplex_RowMatrix & | operator= (const Komplex_RowMatrix &src) |
bool | Filled () const |
If this matrix has been filled, this query returns true; otherwise it returns false. | |
bool | StorageOptimized () const |
If OptimizeStorage() has been called, this query returns true; otherwise it returns false. | |
bool | IndicesAreGlobal () const |
If matrix indices have not been transformed to local, this query returns true; otherwise it returns false. | |
bool | IndicesAreLocal () const |
If matrix indices have been transformed to local, this query returns true; otherwise it returns false. | |
bool | IndicesAreContiguous () const |
If matrix indices are packed into single array (done in OptimizeStorage()) return true, otherwise false. | |
bool | NoDiagonal () const |
If matrix has no diagonal entries in global index space, this query returns true; otherwise it returns false. | |
int | IndexBase () const |
Returns the index base for row and column indices for this graph. | |
const Epetra_CrsGraph & | Graph () const |
Returns a reference to the Epetra_CrsGraph object associated with this matrix. | |
const Epetra_Map & | RowMap () const |
Returns the Epetra_Map object associated with the rows of this matrix. | |
const Epetra_Map & | ColMap () const |
Returns the Epetra_Map object that describes the set of column-indices that appear in each processor's locally owned matrix rows. More... | |
const Epetra_Map & | DomainMap () const |
Returns the Epetra_Map object associated with the domain of this matrix operator. More... | |
const Epetra_Map & | RangeMap () const |
Returns the Epetra_Map object associated with the range of this matrix operator. More... | |
const Epetra_Import * | Importer () const |
Returns the Epetra_Import object that contains the import operations for distributed operations. | |
const Epetra_Export * | Exporter () const |
Returns the Epetra_Export object that contains the export operations for distributed operations. | |
const Epetra_Comm & | Comm () const |
Returns a pointer to the Epetra_Comm communicator associated with this matrix. | |
int | PutScalar (double ScalarConstant) |
Initialize all values in the matrix with constant value. More... | |
int | Scale (double ScalarConstant) |
Multiply all values in the matrix by a constant value (in place: A <- ScalarConstant * A). More... | |
int | ReplaceDiagonalValues (const Epetra_Vector &Diagonal) |
Replaces diagonal values of the matrix with those in the user-provided vector. More... | |
int | FillComplete () |
Signal that data entry is complete. Perform transformations to local index space. | |
int | FillComplete (const Epetra_Map &DomainMap, const Epetra_Map &RangeMap) |
Signal that data entry is complete. Perform transformations to local index space. | |
int | Multiply (bool TransA, const Epetra_MultiVector &X, Epetra_MultiVector &Y) const |
Returns the result of a Komplex_RowMatrix multiplied by a Epetra_MultiVector X in Y. More... | |
int | Solve (bool Upper, bool Trans, bool UnitDiagonal, const Epetra_Vector &x, Epetra_Vector &y) const |
Returns the result of a solve using the Komplex_RowMatrix on a Epetra_Vector x in y. More... | |
int | Solve (bool Upper, bool Trans, bool UnitDiagonal, const Epetra_MultiVector &X, Epetra_MultiVector &Y) const |
Returns result of a local-only solve using a triangular Komplex_RowMatrix with Epetra_MultiVectors X and Y. More... | |
int | InvRowSums (Epetra_Vector &x) const |
Computes the sum of absolute values of the rows of the Komplex_RowMatrix, results returned in x. More... | |
int | LeftScale (const Epetra_Vector &x) |
Scales the Komplex_RowMatrix on the left with a Epetra_Vector x. More... | |
int | InvColSums (Epetra_Vector &x) const |
Computes the sum of absolute values of the columns of the Komplex_RowMatrix, results returned in x. More... | |
int | RightScale (const Epetra_Vector &x) |
Scales the Komplex_RowMatrix on the right with a Epetra_Vector x. More... | |
bool | LowerTriangular () const |
If matrix is lower triangular in local index space, this query returns true; otherwise it returns false. | |
bool | UpperTriangular () const |
If matrix is upper triangular in local index space, this query returns true; otherwise it returns false. | |
double | NormInf () const |
Returns the infinity norm of the global matrix. | |
double | NormOne () const |
Returns the one norm of the global matrix. | |
int | NumMyNonzeros () const |
Returns the number of nonzero entries owned by the calling processor. | |
int | NumMyRows () const |
Returns the number of matrix rows owned by the calling processor. | |
int | NumMyCols () const |
Returns the number of matrix columns owned by the calling processor. | |
int | NumGlobalRows () const |
Returns the number of global matrix rows. | |
int | NumGlobalCols () const |
Returns the number of global matrix columns. | |
int | NumGlobalNonzeros () const |
Returns the number of nonzero entries in the global matrix. | |
int | NumMyDiagonals () const |
Returns the number of local nonzero diagonal entries, based on global row/column index comparisons. | |
int | NumGlobalDiagonals () const |
Returns the number of global nonzero diagonal entries, based on global row/column index comparisons. | |
int | ExtractGlobalRowCopy (int GlobalRow, int Length, int &NumEntries, double *Values, int *Indices) const |
Returns a copy of the specified global row in user-provided arrays. More... | |
int | ExtractGlobalRowCopy (int GlobalRow, int Length, int &NumEntries, double *Values) const |
Returns a copy of the specified global row values in user-provided array. More... | |
int | ExtractMyRowCopy (int MyRow, int Length, int &NumEntries, double *Values, int *Indices) const |
Returns a copy of the specified local row in user-provided arrays. More... | |
int | NumMyRowEntries (int MyRow, int &NumEntries) const |
Return the current number of values stored for the specified local row. More... | |
int | ExtractDiagonalCopy (Epetra_Vector &Diagonal) const |
Returns a copy of the main diagonal in a user-provided vector. More... | |
int | MaxNumEntries () const |
Returns the maximum of NumMyRowEntries() over all rows. | |
const Epetra_Map & | RowMatrixRowMap () const |
Returns the Epetra_Map object associated with the rows of this matrix. | |
const Epetra_Map & | RowMatrixColMap () const |
Returns the Epetra_Map object associated with columns of this matrix. | |
const Epetra_Import * | RowMatrixImporter () const |
Returns the Epetra_Import object that contains the import operations for distributed operations. | |
void | Print (ostream &os) const |
Print method. | |
const char * | Label () const |
Returns a character string describing the operator. | |
int | SetUseTranspose (bool UseTranspose) |
If set true, transpose of this operator will be applied. More... | |
int | Apply (const Epetra_MultiVector &X, Epetra_MultiVector &Y) const |
Returns the result of a Epetra_Operator applied to a Epetra_MultiVector X in Y. More... | |
int | ApplyInverse (const Epetra_MultiVector &X, Epetra_MultiVector &Y) const |
Returns the result of a Epetra_Operator inverse applied to an Epetra_MultiVector X in Y. More... | |
bool | HasNormInf () const |
Returns true because this class can compute an Inf-norm. | |
bool | UseTranspose () const |
Returns the current UseTranspose setting. | |
const Epetra_Map & | OperatorDomainMap () const |
Returns the Epetra_Map object associated with the domain of this matrix operator. | |
const Epetra_Map & | OperatorRangeMap () const |
Returns the Epetra_Map object associated with the range of this matrix operator. | |
Public Member Functions inherited from Epetra_RowMatrix | |
virtual const Epetra_BlockMap & | Map () const =0 |
Komplex_RowMatrix: A class for the construction and use of complex-valued matrices stored in the equivalent real Komplex format. This class implements the pure virtual Epetra_RowMatrix class.
Constructing Komplex_RowMatrix objects
Constructing the Komplex_RowMatrix: The user defines the complex-valued matrix C = (c0r+i*c0i)*A0 +(c1r+i*c1i)*A1. Using this general expression for the complex matrix allows easy formulation of a variety of common complex problems. A0 and A1 are stored in Epetra_VbrMatrix or Epetra_CrsMatrix format.
The different K forms (K1, K2, K3, K4, K14, and K23) can easily convert back and forth by going from one K form to the canonical form to another K form. The Komplex_Ordering that each Komplex_RowMatrix object has is what determines the conversions. Let Kanon stand for the canonical form of a complex matrix in equivalent real formulation. Then any K form is equivalent to: P_l * Kanon * P_r * D_r, where P_l, P_r are specific permutation matrices and D_r is a specific right diagonal scaling matrix.
This is helpful because certain forms are advantageous in certain conditions. To be able to convert back and forth during preconditioning and solving should allow for faster, more accurate solutions.
Komplex_RowMatrix::Komplex_RowMatrix | ( | Epetra_DataAccess | CV, |
const Epetra_BlockMap & | RowMap, | ||
int | NumEntriesPerRow, | ||
Komplex_KForms | KForm = K1 |
||
) |
Default Komplex_RowMatrix constuctor with fixed number of indices per row.
Creates a Komplex_RowMatrix object and allocates storage.
DataAccess | (In) Copy or view. |
RowMap | (In) A Epetra_LocalMap, Epetra_Map or Epetra_BlockMap. |
NumEntriesPerRow | (In) Integer giving the approximate number of entries per row. |
KForm | (In) The Komplex_KForms to use for this RowMatrix; by default, it is set to K1. |
Komplex_RowMatrix::Komplex_RowMatrix | ( | Epetra_DataAccess | CV, |
const Epetra_BlockMap & | RowMap, | ||
int * | NumEntriesPerRow, | ||
Komplex_KForms | KForm = K1 |
||
) |
Default Komplex_RowMatrix constructor with variable number of indices per row.
Creates a Komplex_RowMatrix object and allocates storage.
DataAccess | (In) Copy or view. |
RowMap | (In) A Epetra_LocalMap, Epetra_Map or Epetra_BlockMap. |
NumEntriesPerRow | (In) Integer array giving the approximate number of entries per row. |
KForm | (In) The Komplex_KForms to use for this RowMatrix; by default, it is set to K1. |
Komplex_RowMatrix::Komplex_RowMatrix | ( | Epetra_DataAccess | CV, |
const Epetra_BlockMap & | RowMap, | ||
const Epetra_BlockMap & | ColMap, | ||
int | NumEntriesPerRow, | ||
Komplex_KForms | KForm = K1 |
||
) |
Constructor with fixed number of indices per row and both row and column maps.
Creates a Komplex_RowMatrix object and allocates storage.
DataAccess | (In) Copy or view. |
RowMap | (In) A Epetra_LocalMap, Epetra_Map or Epetra_BlockMap describing the layout of the rows. |
ColMap | (In) A Epetra_LocalMap, Epetra_Map or Epetra_BlockMap describing the layout of the columns. |
NumEntriesPerRow | (In) Integer giving the approximate number of entries per row. |
KForm | (In) The Komplex_KForms to use for this RowMatrix; by default, it is set to K1. |
Komplex_RowMatrix::Komplex_RowMatrix | ( | Epetra_DataAccess | CV, |
const Epetra_BlockMap & | RowMap, | ||
const Epetra_BlockMap & | ColMap, | ||
int * | NumEntriesPerRow, | ||
Komplex_KForms | KForm = K1 |
||
) |
Constructor with variable number of indices per row and both row and column maps.
Creates a Komplex_RowMatrix object and allocates storage.
DataAccess | (In) Copy or view. |
RowMap | (In) A Epetra_LocalMap, Epetra_Map or Epetra_BlockMap describing the layout of the rows. |
ColMap | (In) A Epetra_LocalMap, Epetra_Map or Epetra_BlockMap describing the layout of the columns. |
NumEntriesPerRow | (In) Integer array giving the approximate number of entries per row. |
KForm | (In) The Komplex_KForms to use for this RowMatrix; by default, it is set to K1. |
Komplex_RowMatrix::Komplex_RowMatrix | ( | Epetra_DataAccess | CV, |
const Epetra_CrsGraph & | Graph, | ||
Komplex_KForms | KForm = K1 |
||
) |
Constructor that uses an existing Epetra_CrsGraph object.
Creates a Komplex_RowMatrix object and allocates storage.
DataAccess | (In) Copy or view. |
Graph | (In) A Epetra_CrsGraph. |
KForm | (In) The Komplex_KForms to use for this RowMatrix; by default, it is set to K1. |
Komplex_RowMatrix::Komplex_RowMatrix | ( | double | cr, |
double | ci, | ||
Epetra_RowMatrix * | A, | ||
Komplex_KForms | KForm = K1 |
||
) |
General Komplex_RowMatrix constructor taking one Epetra_RowMatrix object.
Creates a Komplex_RowMatrix object and fills it.
cr | (In) The real part of the complex coefficient multiplying A. |
ci | (In) The imag part of the complex coefficient multiplying A. |
A | (In) An Epetra_RowMatrix that is used to define the true complex matrix, real and imaginary values are interleaved. |
KForm | (In) The Komplex_KForms to use for this RowMatrix; by default, it is set to K1. |
Komplex_RowMatrix::Komplex_RowMatrix | ( | double | c0r, |
double | c0i, | ||
Epetra_RowMatrix * | A0, | ||
double | c1r, | ||
double | c1i, | ||
Epetra_RowMatrix * | A1, | ||
Komplex_KForms | KForm = K1 |
||
) |
General Komplex_RowMatrix constructor.
Creates a Komplex_RowMatrix object and fills it.
c0r | (In) The real part of the complex coefficient multiplying A0. |
c0i | (In) The imag part of the complex coefficient multiplying A0. |
A0 | (In) An Epetra_RowMatrix that is one of the matrices used to define the true complex matrix. |
c1r | (In) The real part of the complex coefficient multiplying A1. |
c1i | (In) The imag part of the complex coefficient multiplying A1. |
A1 | (In) An Epetra_RowMatrix that is the second of the matrices used to define the true complex matrix. |
KForm | (In) The Komplex_KForms to use for this operator; by default, it is set to K1. |
Komplex_RowMatrix::Komplex_RowMatrix | ( | const Komplex_RowMatrix & | Matrix | ) |
Copy constructor.
|
virtual |
Returns the result of a Epetra_Operator applied to a Epetra_MultiVector X in Y.
X | (In) An Epetra_MultiVector of dimension NumVectors to multiply with matrix. |
Y | (Out) An Epetra_MultiVector of dimension NumVectors containing result. |
Implements Epetra_Operator.
|
virtual |
Returns the result of a Epetra_Operator inverse applied to an Epetra_MultiVector X in Y.
In this implementation, we use several existing attributes to determine how virtual method ApplyInverse() should call the concrete method Solve(). We pass in the UpperTriangular(), the Epetra_CrsMatrix::UseTranspose(), and NoDiagonal() methods. The most notable warning is that if a matrix has no diagonal values we assume that there is an implicit unit diagonal that should be accounted for when doing a triangular solve.
X | (In) An Epetra_MultiVector of dimension NumVectors to solve for. |
Y | (Out) An Epetra_MultiVector of dimension NumVectors containing result. |
Implements Epetra_Operator.
const Epetra_Map & Komplex_RowMatrix::ColMap | ( | ) | const |
Returns the Epetra_Map object that describes the set of column-indices that appear in each processor's locally owned matrix rows.
Note that if the matrix was constructed with only a row-map, then until FillComplete() is called, this method returns a column-map that is a copy of the row-map. That 'initial' column-map is replaced with a computed column-map (that contains the set of column-indices appearing in each processor's local portion of the matrix) when FillComplete() is called.
const Epetra_Map & Komplex_RowMatrix::DomainMap | ( | ) | const |
Returns the Epetra_Map object associated with the domain of this matrix operator.
|
virtual |
Returns a copy of the main diagonal in a user-provided vector.
Diagonal | (Out) Extracted main diagonal. |
Implements Epetra_RowMatrix.
int Komplex_RowMatrix::ExtractGlobalRowCopy | ( | int | GlobalRow, |
int | Length, | ||
int & | NumEntries, | ||
double * | Values, | ||
int * | Indices | ||
) | const |
Returns a copy of the specified global row in user-provided arrays.
GlobalRow | (In) Global row to extract. |
ILength | (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. |
int Komplex_RowMatrix::ExtractGlobalRowCopy | ( | int | GlobalRow, |
int | Length, | ||
int & | NumEntries, | ||
double * | Values | ||
) | const |
Returns a copy of the specified global row values in user-provided array.
GlobalRow | (In) Global row to extract. |
Length | (In) Length of Values. |
NumEntries | (Out) Number of nonzero entries extracted. |
Values | (Out) Extracted values for this row. |
|
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 local column indices for the corresponding values. |
Implements Epetra_RowMatrix.
|
virtual |
Computes the sum of absolute values of the columns of the Komplex_RowMatrix, 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 this matrix and 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. Using the resulting vector from this function as input to RightScale() will make the one norm of the resulting matrix exactly 1.
x | (Out) A 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 Komplex_RowMatrix, 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 this matrix and 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) A Epetra_Vector containing the row sums of the this matrix. |
Implements Epetra_RowMatrix.
|
virtual |
Scales the Komplex_RowMatrix on the left with a Epetra_Vector x.
The 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) A Epetra_Vector to use to scale. |
Implements Epetra_RowMatrix.
|
virtual |
Returns the result of a Komplex_RowMatrix multiplied by a Epetra_MultiVector X in Y.
TransA | (In) If true, multiply by the transpose of the matrix, otherwise just use the matrix. |
X | (In) A Epetra_MultiVector of dimension NumVectors to multiply with the matrix. |
Y | (Out) A Epetra_MultiVector of dimension NumVectors containing the result. |
Implements Epetra_RowMatrix.
|
virtual |
Return the current number of values stored for the specified local row.
MyRow | (In) Local row. |
NumEntries | (Out) Number of nonzero values. |
Implements Epetra_RowMatrix.
int Komplex_RowMatrix::PutScalar | ( | double | ScalarConstant | ) |
Initialize all values in the matrix with constant value.
ScalarConstant | (In) Value to use. |
const Epetra_Map & Komplex_RowMatrix::RangeMap | ( | ) | const |
Returns the Epetra_Map object associated with the range of this matrix operator.
int Komplex_RowMatrix::ReplaceDiagonalValues | ( | const Epetra_Vector & | Diagonal | ) |
Replaces diagonal values of the matrix with those in the user-provided vector.
This routine is meant to allow replacement of { existing} diagonal values. If a diagonal value does not exist for a given row, the corresponding value in the input Epetra_Vector will be ignored and the return code will be set to 1.
The Epetra_Map associated with the input Epetra_Vector must be compatible with the RowMap of the matrix.
Diagonal | (In) New values to be placed in the main diagonal. |
|
virtual |
Scales the Komplex_RowMatrix on the right with a Epetra_Vector x.
The 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.
int Komplex_RowMatrix::Scale | ( | double | ScalarConstant | ) |
Multiply all values in the matrix by a constant value (in place: A <- ScalarConstant * A).
ScalarConstant | (In) Value to use. |
|
virtual |
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.
UseTranspose | (In) If true, multiply by the transpose of operator, otherwise just use operator. |
Implements Epetra_Operator.
int Komplex_RowMatrix::Solve | ( | bool | Upper, |
bool | Trans, | ||
bool | UnitDiagonal, | ||
const Epetra_Vector & | x, | ||
Epetra_Vector & | y | ||
) | const |
Returns the result of a solve using the Komplex_RowMatrix on a Epetra_Vector x in y.
Upper | (In) If true, solve Ux = y, otherwise solve Lx = y. |
Trans | (In) If true, solve the transpose problem. |
UnitDiagonal | (In) If true, assume diagonal is unit (whether it's stored or not). |
x | (In) A Epetra_Vector to solve for. |
y | (Out) A Epetra_Vector containing the result. |
|
virtual |
Returns result of a local-only solve using a triangular Komplex_RowMatrix with Epetra_MultiVectors X and Y.
This method will perform a triangular solve independently on each processor of the parallel machine. No communication is performed.
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) A Epetra_MultiVector of dimension NumVectors to solve for. |
Y | (Out) A Epetra_MultiVector of dimension NumVectors containing result. |
Implements Epetra_RowMatrix.