Epetra_SerialDenseMatrix Class Reference

Epetra_SerialDenseMatrix: A class for constructing and using real double precision general dense matrices. More...

#include <Epetra_SerialDenseMatrix.h>

Inheritance diagram for Epetra_SerialDenseMatrix:
Inheritance graph
[legend]

List of all members.

Protected Member Functions

void CopyMat (double *Source, int Source_LDA, int NumRows, int NumCols, double *Target, int Target_LDA, bool add=false)
void CleanupData ()

Protected Attributes

int M_
int N_
bool A_Copied_
Epetra_DataAccess CV_
int LDA_
double * A_
bool UseTranspose_

Friends

class Epetra_VbrMatrix

Constructor/Destructor Methods



 Epetra_SerialDenseMatrix (bool set_object_label=true)
 Default constructor; defines a zero size object.
 Epetra_SerialDenseMatrix (int NumRows, int NumCols, bool set_object_label=true)
 Shaped constructor; defines a variable-sized object.
 Epetra_SerialDenseMatrix (Epetra_DataAccess CV, double *A_in, int LDA_in, int NumRows, int NumCols, bool set_object_label=true)
 Set object values from two-dimensional array.
 Epetra_SerialDenseMatrix (const Epetra_SerialDenseMatrix &Source)
 Epetra_SerialDenseMatrix copy constructor.
virtual ~Epetra_SerialDenseMatrix ()
 Epetra_SerialDenseMatrix destructor.

Shaping/sizing Methods



int Shape (int NumRows, int NumCols)
 Set dimensions of a Epetra_SerialDenseMatrix object; init values to zero.
int Reshape (int NumRows, int NumCols)
 Reshape a Epetra_SerialDenseMatrix object.

Mathematical methods



int Multiply (char TransA, char TransB, double ScalarAB, const Epetra_SerialDenseMatrix &A, const Epetra_SerialDenseMatrix &B, double ScalarThis)
 Matrix-Matrix multiplication, this = ScalarThis*this + ScalarAB*A*B.
int Multiply (bool transA, const Epetra_SerialDenseMatrix &x, Epetra_SerialDenseMatrix &y)
 Matrix-Vector multiplication, y = A*x, where 'this' == A.
int Multiply (char SideA, double ScalarAB, const Epetra_SerialSymDenseMatrix &A, const Epetra_SerialDenseMatrix &B, double ScalarThis)
 Matrix-Matrix multiplication with a symmetric matrix A.
int Scale (double ScalarA)
 Inplace scalar-matrix product A = a A.
virtual double NormOne () const
 Computes the 1-Norm of the this matrix.
virtual double NormInf () const
 Computes the Infinity-Norm of the this matrix.

Data Accessor methods



Epetra_SerialDenseMatrixoperator= (const Epetra_SerialDenseMatrix &Source)
 Value copy from one matrix to another.
bool operator== (const Epetra_SerialDenseMatrix &rhs) const
 Comparison operator.
bool operator!= (const Epetra_SerialDenseMatrix &rhs) const
 Inequality operator.
Epetra_SerialDenseMatrixoperator+= (const Epetra_SerialDenseMatrix &Source)
 Add one matrix to another.
double & operator() (int RowIndex, int ColIndex)
 Element access function.
const double & operator() (int RowIndex, int ColIndex) const
 Element access function.
double * operator[] (int ColIndex)
 Column access function.
const double * operator[] (int ColIndex) const
 Column access function.
int Random ()
 Set matrix values to random numbers.
int M () const
 Returns row dimension of system.
int N () const
 Returns column dimension of system.
double * A () const
 Returns pointer to the this matrix.
double * A ()
 Returns pointer to the this matrix.
int LDA () const
 Returns the leading dimension of the this matrix.
Epetra_DataAccess CV () const
 Returns the data access mode of the this matrix.

I/O methods



virtual void Print (ostream &os) const
 Print service methods; defines behavior of ostream << operator.

Deprecated methods (will be removed in later versions of this class)



virtual double OneNorm () const
 Computes the 1-Norm of the this matrix (identical to NormOne() method).
virtual double InfNorm () const
 Computes the Infinity-Norm of the this matrix (identical to NormInf() method).

Additional methods to support Epetra_SerialDenseOperator interface



virtual int SetUseTranspose (bool UseTranspose_in)
 If set true, transpose of this operator will be applied.
virtual int Apply (const Epetra_SerialDenseMatrix &X, Epetra_SerialDenseMatrix &Y)
 Returns the result of a Epetra_SerialDenseOperator applied to a Epetra_SerialDenseMatrix X in Y.
virtual int ApplyInverse (const Epetra_SerialDenseMatrix &X, Epetra_SerialDenseMatrix &Y)
 Returns the result of a Epetra_SerialDenseOperator inverse applied to an Epetra_SerialDenseMatrix X in Y.
virtual const char * Label () const
 Returns a character string describing the operator.
virtual bool UseTranspose () const
 Returns the current UseTranspose setting.
virtual bool HasNormInf () const
 Returns true if the this object can provide an approximate Inf-norm, false otherwise.
virtual int RowDim () const
 Returns the row dimension of operator.
virtual int ColDim () const
 Returns the column dimension of operator.

Detailed Description

Epetra_SerialDenseMatrix: A class for constructing and using real double precision general dense matrices.

The Epetra_SerialDenseMatrix class enables the construction and use of real-valued, general, double-precision dense matrices. It is built on the BLAS, and derives from the Epetra_BLAS.

The Epetra_SerialDenseMatrix class is intended to provide very basic support for dense rectangular matrices.

Constructing Epetra_SerialDenseMatrix Objects

There are four Epetra_SerialDenseMatrix constructors. The first constructs a zero-sized object which should be made to appropriate length using the Shape() or Reshape() functions and then filled with the [] or () operators. The second constructs an object sized to the dimensions specified, which should be filled with the [] or () operators. The third is a constructor that accepts user data as a 2D array, and the fourth is a copy constructor. The third constructor has two data access modes (specified by the Epetra_DataAccess argument):

  1. Copy mode - Allocates memory and makes a copy of the user-provided data. In this case, the user data is not needed after construction.
  2. View mode - Creates a "view" of the user data. In this case, the user data is required to remain intact for the life of the object.
Warning:
View mode is extremely dangerous from a data hiding perspective. Therefore, we strongly encourage users to develop code using Copy mode first and only use the View mode in a secondary optimization phase.

Extracting Data from Epetra_SerialDenseMatrix Objects

Once a Epetra_SerialDenseMatrix is constructed, it is possible to view the data via access functions.

Warning:
Use of these access functions cam be extremely dangerous from a data hiding perspective.

Vector and Utility Functions

Once a Epetra_SerialDenseMatrix is constructed, several mathematical functions can be applied to the object. Specifically:

Counting floating point operations The Epetra_SerialDenseMatrix class has Epetra_CompObject as a base class. Thus, floating point operations are counted and accumulated in the Epetra_Flop object (if any) that was set using the SetFlopCounter() method in the Epetra_CompObject base class.

Definition at line 95 of file Epetra_SerialDenseMatrix.h.


Constructor & Destructor Documentation

Epetra_SerialDenseMatrix::Epetra_SerialDenseMatrix ( bool  set_object_label = true  ) 

Default constructor; defines a zero size object.

Epetra_SerialDenseMatrix objects defined by the default constructor should be sized with the Shape() or Reshape functions. Values should be defined by using the [] or () operators.

Definition at line 36 of file Epetra_SerialDenseMatrix.cpp.

Epetra_SerialDenseMatrix::Epetra_SerialDenseMatrix ( int  NumRows,
int  NumCols,
bool  set_object_label = true 
)

Shaped constructor; defines a variable-sized object.

Parameters:
In NumRows - Number of rows in object.
In NumCols - Number of columns in object.

Epetra_SerialDenseMatrix objects defined by the shaped constructor are already shaped to the dimensions given as a parameters. All values are initialized to 0. Calling this constructor is equivalent to using the default constructor, and then calling the Shape function on it. Values should be defined by using the [] or () operators.

Definition at line 53 of file Epetra_SerialDenseMatrix.cpp.

Epetra_SerialDenseMatrix::Epetra_SerialDenseMatrix ( Epetra_DataAccess  CV,
double *  A_in,
int  LDA_in,
int  NumRows,
int  NumCols,
bool  set_object_label = true 
)

Set object values from two-dimensional array.

Parameters:
In Epetra_DataAccess - Enumerated type set to Copy or View.
In A - Pointer to an array of double precision numbers. The first vector starts at A. The second vector starts at A+LDA, the third at A+2*LDA, and so on.
In LDA - The "Leading Dimension", or stride between vectors in memory.
In NumRows - Number of rows in object.
In NumCols - Number of columns in object.

See Detailed Description section for further discussion.

Definition at line 79 of file Epetra_SerialDenseMatrix.cpp.

Epetra_SerialDenseMatrix::Epetra_SerialDenseMatrix ( const Epetra_SerialDenseMatrix Source  ) 

Epetra_SerialDenseMatrix copy constructor.

Definition at line 118 of file Epetra_SerialDenseMatrix.cpp.

Epetra_SerialDenseMatrix::~Epetra_SerialDenseMatrix (  )  [virtual]

Epetra_SerialDenseMatrix destructor.

Definition at line 192 of file Epetra_SerialDenseMatrix.cpp.


Member Function Documentation

int Epetra_SerialDenseMatrix::Shape ( int  NumRows,
int  NumCols 
)

Set dimensions of a Epetra_SerialDenseMatrix object; init values to zero.

Parameters:
In NumRows - Number of rows in object.
In NumCols - Number of columns in object.

Allows user to define the dimensions of a Epetra_SerialDenseMatrix at any point. This function can be called at any point after construction. Any values that were previously in this object are destroyed and the resized matrix starts off with all zero values.

Returns:
Integer error code, set to 0 if successful.

Definition at line 173 of file Epetra_SerialDenseMatrix.cpp.

int Epetra_SerialDenseMatrix::Reshape ( int  NumRows,
int  NumCols 
)

Reshape a Epetra_SerialDenseMatrix object.

Parameters:
In NumRows - Number of rows in object.
In NumCols - Number of columns in object.

Allows user to define the dimensions of a Epetra_SerialDenseMatrix at any point. This function can be called at any point after construction. Any values that were previously in this object are copied into the new shape. If the new shape is smaller than the original, the upper left portion of the original matrix (the principal submatrix) is copied to the new matrix.

Returns:
Integer error code, set to 0 if successful.

Definition at line 144 of file Epetra_SerialDenseMatrix.cpp.

int Epetra_SerialDenseMatrix::Multiply ( char  TransA,
char  TransB,
double  ScalarAB,
const Epetra_SerialDenseMatrix A,
const Epetra_SerialDenseMatrix B,
double  ScalarThis 
)

Matrix-Matrix multiplication, this = ScalarThis*this + ScalarAB*A*B.

This function performs a variety of matrix-matrix multiply operations.

Parameters:
In TransA - Operate with the transpose of A if = 'T', else no transpose if = 'N'.
In TransB - Operate with the transpose of B if = 'T', else no transpose if = 'N'.
In ScalarAB - Scalar to multiply with A*B.
In A - Dense Matrix.
In B - Dense Matrix.
In ScalarThis - Scalar to multiply with this.
Returns:
Integer error code, set to 0 if successful.

Definition at line 407 of file Epetra_SerialDenseMatrix.cpp.

int Epetra_SerialDenseMatrix::Multiply ( bool  transA,
const Epetra_SerialDenseMatrix x,
Epetra_SerialDenseMatrix y 
)

Matrix-Vector multiplication, y = A*x, where 'this' == A.

Definition at line 440 of file Epetra_SerialDenseMatrix.cpp.

int Epetra_SerialDenseMatrix::Multiply ( char  SideA,
double  ScalarAB,
const Epetra_SerialSymDenseMatrix A,
const Epetra_SerialDenseMatrix B,
double  ScalarThis 
)

Matrix-Matrix multiplication with a symmetric matrix A.

If SideA = 'L', compute this = ScalarThis*this + ScalarAB*A*B. If SideA = 'R', compute this = ScalarThis*this + ScalarAB*B*A.

This function performs a variety of matrix-matrix multiply operations.

Parameters:
In SideA - Specifies order of A relative to B.
In ScalarAB - Scalar to multiply with A*B.
In A - Symmetric Dense Matrix, either upper or lower triangle will be used depending on value of A.Upper().
In B - Dense Matrix.
In ScalarThis - Scalar to multiply with this.
Returns:
Integer error code, set to 0 if successful.

Definition at line 475 of file Epetra_SerialDenseMatrix.cpp.

int Epetra_SerialDenseMatrix::Scale ( double  ScalarA  ) 

Inplace scalar-matrix product A = a A.

Scale a matrix, entry-by-entry using the value ScalarA.

Parameters:
ScalarA (In) Scalar to multiply with A.
Returns:
Integer error code, set to 0 if successful.

Reimplemented in Epetra_SerialSymDenseMatrix.

Definition at line 393 of file Epetra_SerialDenseMatrix.cpp.

double Epetra_SerialDenseMatrix::NormOne ( void   )  const [virtual]

Computes the 1-Norm of the this matrix.

Returns:
Integer error code, set to 0 if successful.

Reimplemented in Epetra_SerialSymDenseMatrix.

Definition at line 355 of file Epetra_SerialDenseMatrix.cpp.

double Epetra_SerialDenseMatrix::NormInf ( void   )  const [virtual]

Computes the Infinity-Norm of the this matrix.

Implements Epetra_SerialDenseOperator.

Reimplemented in Epetra_SerialDenseVector, and Epetra_SerialSymDenseMatrix.

Definition at line 371 of file Epetra_SerialDenseMatrix.cpp.

Epetra_SerialDenseMatrix & Epetra_SerialDenseMatrix::operator= ( const Epetra_SerialDenseMatrix Source  ) 

Value copy from one matrix to another.

The operator= allows one to copy the values from one existing SerialDenseMatrix to another, as long as there is enough room in the target to hold the source.

Returns:
Values of the left hand side matrix are modified by the values of the right hand side matrix.

Definition at line 208 of file Epetra_SerialDenseMatrix.cpp.

bool Epetra_SerialDenseMatrix::operator== ( const Epetra_SerialDenseMatrix rhs  )  const

Comparison operator.

operator== compares two Epetra_SerialDenseMatrix objects, returns false if sizes are different, or if any coefficients differ by an amount greater than Epetra_MinDouble.

Definition at line 269 of file Epetra_SerialDenseMatrix.cpp.

bool Epetra_SerialDenseMatrix::operator!= ( const Epetra_SerialDenseMatrix rhs  )  const [inline]

Inequality operator.

operator!= simply returns the negation of operator==.

Definition at line 292 of file Epetra_SerialDenseMatrix.h.

Epetra_SerialDenseMatrix & Epetra_SerialDenseMatrix::operator+= ( const Epetra_SerialDenseMatrix Source  ) 

Add one matrix to another.

The operator+= allows one to add the values from one existin SerialDenseMatrix to another, as long as there is enough room in the target to hold the source.

Returns:
Values of the left hand side matrix are modified by the addition of the values of the right hand side matrix.

Definition at line 290 of file Epetra_SerialDenseMatrix.cpp.

double & Epetra_SerialDenseMatrix::operator() ( int  RowIndex,
int  ColIndex 
) [inline]

Element access function.

The parentheses operator returns the element in the ith row and jth column if A(i,j) is specified, the expression A[j][i] (note that i and j are reversed) will return the same element. Thus, A(i,j) = A[j][i] for all valid i and j.

Returns:
Element from the specified row and column.
Warning:
No bounds checking is done unless Epetra is compiled with HAVE_EPETRA_ARRAY_BOUNDS_CHECK.

Definition at line 487 of file Epetra_SerialDenseMatrix.h.

const double & Epetra_SerialDenseMatrix::operator() ( int  RowIndex,
int  ColIndex 
) const [inline]

Element access function.

The parentheses operator returns the element in the ith row and jth column if A(i,j) is specified, the expression A[j][i] (note that i and j are reversed) will return the same element. Thus, A(i,j) = A[j][i] for all valid i and j.

Returns:
Element from the specified row and column.
Warning:
No bounds checking is done unless Epetra is compiled with HAVE_EPETRA_ARRAY_BOUNDS_CHECK.

Definition at line 499 of file Epetra_SerialDenseMatrix.h.

double * Epetra_SerialDenseMatrix::operator[] ( int  ColIndex  )  [inline]

Column access function.

The parentheses operator returns the element in the ith row and jth column if A(i,j) is specified, the expression A[j][i] (note that i and j are reversed) will return the same element. Thus, A(i,j) = A[j][i] for all valid i and j.

Returns:
Pointer to address of specified column.
Warning:
No bounds checking can be done for the index i in the expression A[j][i].
No bounds checking is done unless Epetra is compiled with HAVE_EPETRA_ARRAY_BOUNDS_CHECK.

Reimplemented in Epetra_SerialDenseVector.

Definition at line 511 of file Epetra_SerialDenseMatrix.h.

const double * Epetra_SerialDenseMatrix::operator[] ( int  ColIndex  )  const [inline]

Column access function.

The parentheses operator returns the element in the ith row and jth column if A(i,j) is specified, the expression A[j][i] (note that i and j are reversed) will return the same element. Thus, A(i,j) = A[j][i] for all valid i and j.

Returns:
Pointer to address of specified column.
Warning:
No bounds checking can be done for the index i in the expression A[j][i].
No bounds checking is done unless Epetra is compiled with HAVE_EPETRA_ARRAY_BOUNDS_CHECK.

Reimplemented in Epetra_SerialDenseVector.

Definition at line 520 of file Epetra_SerialDenseMatrix.h.

int Epetra_SerialDenseMatrix::Random (  ) 

Set matrix values to random numbers.

SerialDenseMatrix uses the random number generator provided by Epetra_Util. The matrix values will be set to random values on the interval (-1.0, 1.0).

Returns:
Integer error code, set to 0 if successful.

Reimplemented in Epetra_SerialDenseVector.

Definition at line 532 of file Epetra_SerialDenseMatrix.cpp.

int Epetra_SerialDenseMatrix::M (  )  const [inline]

Returns row dimension of system.

Definition at line 365 of file Epetra_SerialDenseMatrix.h.

int Epetra_SerialDenseMatrix::N (  )  const [inline]

Returns column dimension of system.

Definition at line 368 of file Epetra_SerialDenseMatrix.h.

double* Epetra_SerialDenseMatrix::A (  )  const [inline]

Returns pointer to the this matrix.

Definition at line 371 of file Epetra_SerialDenseMatrix.h.

double* Epetra_SerialDenseMatrix::A (  )  [inline]

Returns pointer to the this matrix.

Definition at line 374 of file Epetra_SerialDenseMatrix.h.

int Epetra_SerialDenseMatrix::LDA (  )  const [inline]

Returns the leading dimension of the this matrix.

Definition at line 377 of file Epetra_SerialDenseMatrix.h.

Epetra_DataAccess Epetra_SerialDenseMatrix::CV (  )  const [inline]

Returns the data access mode of the this matrix.

Reimplemented in Epetra_SerialDenseVector.

Definition at line 380 of file Epetra_SerialDenseMatrix.h.

void Epetra_SerialDenseMatrix::Print ( ostream &  os  )  const [virtual]

Print service methods; defines behavior of ostream << operator.

Reimplemented from Epetra_Object.

Reimplemented in Epetra_SerialDenseVector.

Definition at line 508 of file Epetra_SerialDenseMatrix.cpp.

virtual double Epetra_SerialDenseMatrix::OneNorm (  )  const [inline, virtual]

Computes the 1-Norm of the this matrix (identical to NormOne() method).

Returns:
Integer error code, set to 0 if successful.

Reimplemented in Epetra_SerialSymDenseMatrix.

Definition at line 396 of file Epetra_SerialDenseMatrix.h.

virtual double Epetra_SerialDenseMatrix::InfNorm (  )  const [inline, virtual]

Computes the Infinity-Norm of the this matrix (identical to NormInf() method).

Reimplemented in Epetra_SerialSymDenseMatrix.

Definition at line 399 of file Epetra_SerialDenseMatrix.h.

virtual int Epetra_SerialDenseMatrix::SetUseTranspose ( bool  UseTranspose_in  )  [inline, 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.

Parameters:
In UseTranspose -If true, multiply by the transpose of operator, otherwise just use operator.
Returns:
Integer error code, set to 0 if successful. Set to -1 if this implementation does not support transpose.

Implements Epetra_SerialDenseOperator.

Definition at line 415 of file Epetra_SerialDenseMatrix.h.

int Epetra_SerialDenseMatrix::Apply ( const Epetra_SerialDenseMatrix X,
Epetra_SerialDenseMatrix Y 
) [virtual]

Returns the result of a Epetra_SerialDenseOperator applied to a Epetra_SerialDenseMatrix X in Y.

Parameters:
In X - A Epetra_SerialDenseMatrix to multiply with operator.
Out Y -A Epetra_SerialDenseMatrix containing result.
Returns:
Integer error code, set to 0 if successful.

Implements Epetra_SerialDenseOperator.

Definition at line 547 of file Epetra_SerialDenseMatrix.cpp.

virtual int Epetra_SerialDenseMatrix::ApplyInverse ( const Epetra_SerialDenseMatrix X,
Epetra_SerialDenseMatrix Y 
) [inline, virtual]

Returns the result of a Epetra_SerialDenseOperator inverse applied to an Epetra_SerialDenseMatrix X in Y.

Parameters:
In X - A Epetra_SerialDenseMatrix to solve for.
Out Y -A Epetra_SerialDenseMatrix containing result.
Returns:
Integer error code, set to 0 if successful.

Implements Epetra_SerialDenseOperator.

Definition at line 438 of file Epetra_SerialDenseMatrix.h.

virtual const char* Epetra_SerialDenseMatrix::Label (  )  const [inline, virtual]

Returns a character string describing the operator.

Reimplemented from Epetra_Object.

Definition at line 446 of file Epetra_SerialDenseMatrix.h.

virtual bool Epetra_SerialDenseMatrix::UseTranspose (  )  const [inline, virtual]

Returns the current UseTranspose setting.

Implements Epetra_SerialDenseOperator.

Definition at line 449 of file Epetra_SerialDenseMatrix.h.

virtual bool Epetra_SerialDenseMatrix::HasNormInf (  )  const [inline, virtual]

Returns true if the this object can provide an approximate Inf-norm, false otherwise.

Implements Epetra_SerialDenseOperator.

Definition at line 452 of file Epetra_SerialDenseMatrix.h.

virtual int Epetra_SerialDenseMatrix::RowDim (  )  const [inline, virtual]

Returns the row dimension of operator.

Implements Epetra_SerialDenseOperator.

Definition at line 455 of file Epetra_SerialDenseMatrix.h.

virtual int Epetra_SerialDenseMatrix::ColDim (  )  const [inline, virtual]

Returns the column dimension of operator.

Implements Epetra_SerialDenseOperator.

Definition at line 458 of file Epetra_SerialDenseMatrix.h.

void Epetra_SerialDenseMatrix::CopyMat ( double *  Source,
int  Source_LDA,
int  NumRows,
int  NumCols,
double *  Target,
int  Target_LDA,
bool  add = false 
) [protected]

Definition at line 302 of file Epetra_SerialDenseMatrix.cpp.

void Epetra_SerialDenseMatrix::CleanupData (  )  [protected]

Definition at line 197 of file Epetra_SerialDenseMatrix.cpp.


Friends And Related Function Documentation

friend class Epetra_VbrMatrix [friend]

Definition at line 477 of file Epetra_SerialDenseMatrix.h.


Member Data Documentation

Definition at line 467 of file Epetra_SerialDenseMatrix.h.

Definition at line 468 of file Epetra_SerialDenseMatrix.h.

Definition at line 469 of file Epetra_SerialDenseMatrix.h.

Definition at line 470 of file Epetra_SerialDenseMatrix.h.

Definition at line 479 of file Epetra_SerialDenseMatrix.h.

double* Epetra_SerialDenseMatrix::A_ [protected]

Definition at line 480 of file Epetra_SerialDenseMatrix.h.

Definition at line 482 of file Epetra_SerialDenseMatrix.h.


The documentation for this class was generated from the following files:
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines
Generated on Wed Apr 13 09:58:46 2011 for Epetra Package Browser (Single Doxygen Collection) by  doxygen 1.6.3