MFEM
v3.0
|
Data type dense matrix. More...
#include <densemat.hpp>
Public Member Functions | |
DenseMatrix () | |
DenseMatrix (const DenseMatrix &) | |
Copy constructor. More... | |
DenseMatrix (int s) | |
Creates square matrix of size s. More... | |
DenseMatrix (int m, int n) | |
Creates rectangular matrix of size m x n. More... | |
DenseMatrix (const DenseMatrix &mat, char ch) | |
Creates rectangular matrix equal to the transpose of mat. More... | |
DenseMatrix (double *d, int h, int w) | |
void | UseExternalData (double *d, int h, int w) |
void | ClearExternalData () |
int | Size () const |
For backward compatibility define Size to be synonym of Width() More... | |
void | SetSize (int s) |
If the matrix is not a square matrix of size s then recreate it. More... | |
void | SetSize (int h, int w) |
If the matrix is not a matrix of size (h x w) then recreate it. More... | |
double * | Data () const |
Returns vector of the elements. More... | |
double & | operator() (int i, int j) |
Returns reference to a_{ij}. More... | |
const double & | operator() (int i, int j) const |
Returns constant reference to a_{ij}. More... | |
double | operator* (const DenseMatrix &m) const |
Matrix inner product: tr(A^t B) More... | |
double | Trace () const |
Trace of a square matrix. More... | |
virtual double & | Elem (int i, int j) |
Returns reference to a_{ij}. More... | |
virtual const double & | Elem (int i, int j) const |
Returns constant reference to a_{ij}. More... | |
void | Mult (const double *x, double *y) const |
Matrix vector multiplication. More... | |
virtual void | Mult (const Vector &x, Vector &y) const |
Matrix vector multiplication. More... | |
void | MultTranspose (const double *x, double *y) const |
Multiply a vector with the transpose matrix. More... | |
virtual void | MultTranspose (const Vector &x, Vector &y) const |
Multiply a vector with the transpose matrix. More... | |
void | AddMult (const Vector &x, Vector &y) const |
y += A.x More... | |
double | InnerProduct (const double *x, const double *y) const |
Compute y^t A x. More... | |
void | LeftScaling (const Vector &s) |
LeftScaling this = diag(s) * this. More... | |
void | InvLeftScaling (const Vector &s) |
InvLeftScaling this = diag(1./s) * this. More... | |
void | RightScaling (const Vector &s) |
RightScaling: this = this * diag(s);. More... | |
void | InvRightScaling (const Vector &s) |
InvRightScaling: this = this * diag(1./s);. More... | |
void | SymmetricScaling (const Vector &s) |
SymmetricScaling this = diag(sqrt(s)) * this * diag(sqrt(s)) More... | |
void | InvSymmetricScaling (const Vector &s) |
InvSymmetricScaling this = diag(sqrt(1./s)) * this * diag(sqrt(1./s)) More... | |
double | InnerProduct (const Vector &x, const Vector &y) const |
Compute y^t A x. More... | |
virtual MatrixInverse * | Inverse () const |
Returns a pointer to the inverse matrix. More... | |
void | Invert () |
Replaces the current matrix with its inverse. More... | |
double | Det () const |
Calculates the determinant of the matrix (for 2x2 or 3x3 matrices) More... | |
double | Weight () const |
void | Add (const double c, const DenseMatrix &A) |
Adds the matrix A multiplied by the number c to the matrix. More... | |
DenseMatrix & | operator= (double c) |
Sets the matrix elements equal to constant c. More... | |
DenseMatrix & | operator= (const double *d) |
Copy the matrix entries from the given array. More... | |
DenseMatrix & | operator= (const DenseMatrix &m) |
Sets the matrix size and elements equal to those of m. More... | |
DenseMatrix & | operator+= (DenseMatrix &m) |
DenseMatrix & | operator-= (DenseMatrix &m) |
DenseMatrix & | operator*= (double c) |
void | Neg () |
(*this) = -(*this) More... | |
void | Norm2 (double *v) const |
Take the 2-norm of the columns of A and store in v. More... | |
double | MaxMaxNorm () const |
Compute the norm ||A|| = max_{ij} |A_{ij}|. More... | |
double | FNorm () const |
Compute the Frobenius norm of the matrix. More... | |
void | Eigenvalues (Vector &ev) |
void | Eigenvalues (Vector &ev, DenseMatrix &evect) |
void | Eigensystem (Vector &ev, DenseMatrix &evect) |
void | SingularValues (Vector &sv) const |
int | Rank (double tol) const |
double | CalcSingularvalue (const int i) const |
Return the i-th singular value (decreasing order) of NxN matrix, N=1,2,3. More... | |
void | CalcEigenvalues (double *lambda, double *vec) const |
void | GetColumn (int c, Vector &col) |
void | GetColumnReference (int c, Vector &col) |
void | GetDiag (Vector &d) |
Returns the diagonal of the matrix. More... | |
void | Getl1Diag (Vector &l) |
Returns the l1 norm of the rows of the matrix v_i = sum_j |a_ij|. More... | |
void | Diag (double c, int n) |
Creates n x n diagonal matrix with diagonal elements c. More... | |
void | Diag (double *diag, int n) |
Creates n x n diagonal matrix with diagonal given by diag. More... | |
void | Transpose () |
(*this) = (*this)^t More... | |
void | Transpose (DenseMatrix &A) |
(*this) = A^t More... | |
void | Symmetrize () |
(*this) = 1/2 ((*this) + (*this)^t) More... | |
void | Lump () |
void | GradToCurl (DenseMatrix &curl) |
void | GradToDiv (Vector &div) |
void | CopyRows (DenseMatrix &A, int row1, int row2) |
Copy rows row1 through row2 from A to *this. More... | |
void | CopyCols (DenseMatrix &A, int col1, int col2) |
Copy columns col1 through col2 from A to *this. More... | |
void | CopyMN (DenseMatrix &A, int m, int n, int Aro, int Aco) |
Copy the m x n submatrix of A at (Aro)ffset, (Aco)loffset to *this. More... | |
void | CopyMN (DenseMatrix &A, int row_offset, int col_offset) |
Copy matrix A to the location in *this at row_offset, col_offset. More... | |
void | CopyMNt (DenseMatrix &A, int row_offset, int col_offset) |
Copy matrix A^t to the location in *this at row_offset, col_offset. More... | |
void | CopyMNDiag (double c, int n, int row_offset, int col_offset) |
Copy c on the diagonal of size n to *this at row_offset, col_offset. More... | |
void | CopyMNDiag (double *diag, int n, int row_offset, int col_offset) |
Copy diag on the diagonal of size n to *this at row_offset, col_offset. More... | |
void | AddMatrix (DenseMatrix &A, int ro, int co) |
Perform (ro+i,co+j)+=A(i,j) for 0<=i<A.Height, 0<=j<A.Width. More... | |
void | AddMatrix (double a, DenseMatrix &A, int ro, int co) |
Perform (ro+i,co+j)+=a*A(i,j) for 0<=i<A.Height, 0<=j<A.Width. More... | |
void | AddToVector (int offset, Vector &v) const |
Add the matrix 'data' to the Vector 'v' at the given 'offset'. More... | |
void | GetFromVector (int offset, const Vector &v) |
Get the matrix 'data' from the Vector 'v' at the given 'offset'. More... | |
void | AdjustDofDirection (Array< int > &dofs) |
void | SetRow (int row, double value) |
Set all entries of a row to the specified value. More... | |
void | SetCol (int col, double value) |
Set all entries of a column to the specified value. More... | |
int | CheckFinite () const |
virtual void | Print (std::ostream &out=std::cout, int width_=4) const |
Prints matrix to stream out. More... | |
virtual void | PrintMatlab (std::ostream &out=std::cout) const |
virtual void | PrintT (std::ostream &out=std::cout, int width_=4) const |
Prints the transpose matrix to stream out. More... | |
void | TestInversion () |
Invert and print the numerical conditioning of the inversion. More... | |
virtual | ~DenseMatrix () |
Destroys dense matrix. More... | |
Public Member Functions inherited from mfem::Matrix | |
Matrix (int s) | |
Creates a square matrix of size s. More... | |
Matrix (int h, int w) | |
Creates a matrix of the given height and width. More... | |
virtual void | Finalize (int) |
Finalizes the matrix initialization. More... | |
virtual | ~Matrix () |
Destroys matrix. More... | |
Public Member Functions inherited from mfem::Operator | |
Operator (int s=0) | |
Construct a square Operator with given size s (default 0) More... | |
Operator (int h, int w) | |
int | Height () const |
Get the height (size of output) of the Operator. Synonym with NumRows. More... | |
int | NumRows () const |
int | Width () const |
Get the width (size of input) of the Operator. Synonym with NumCols. More... | |
int | NumCols () const |
virtual Operator & | GetGradient (const Vector &x) const |
Evaluate the gradient operator at the point x. More... | |
void | PrintMatlab (std::ostream &out, int n=0, int m=0) |
Prints operator with input size n and output size m in matlab format. More... | |
virtual | ~Operator () |
Friends | |
class | DenseTensor |
class | DenseMatrixInverse |
void | Mult (const DenseMatrix &b, const DenseMatrix &c, DenseMatrix &a) |
Matrix matrix multiplication. A = B * C. More... | |
Additional Inherited Members | |
Protected Attributes inherited from mfem::Operator | |
int | height |
int | width |
Data type dense matrix.
Definition at line 22 of file densemat.hpp.
mfem::DenseMatrix::DenseMatrix | ( | ) |
Default constructor for DenseMatrix. Sets data = NULL and height = width = 0.
Definition at line 32 of file densemat.cpp.
mfem::DenseMatrix::DenseMatrix | ( | const DenseMatrix & | m | ) |
Copy constructor.
Definition at line 37 of file densemat.cpp.
|
explicit |
Creates square matrix of size s.
Definition at line 45 of file densemat.cpp.
mfem::DenseMatrix::DenseMatrix | ( | int | m, |
int | n | ||
) |
Creates rectangular matrix of size m x n.
Definition at line 54 of file densemat.cpp.
mfem::DenseMatrix::DenseMatrix | ( | const DenseMatrix & | mat, |
char | ch | ||
) |
Creates rectangular matrix equal to the transpose of mat.
Definition at line 63 of file densemat.cpp.
|
inline |
Definition at line 53 of file densemat.hpp.
|
virtual |
Destroys dense matrix.
Definition at line 2422 of file densemat.cpp.
void mfem::DenseMatrix::Add | ( | const double | c, |
const DenseMatrix & | A | ||
) |
Adds the matrix A multiplied by the number c to the matrix.
Definition at line 375 of file densemat.cpp.
void mfem::DenseMatrix::AddMatrix | ( | DenseMatrix & | A, |
int | ro, | ||
int | co | ||
) |
Perform (ro+i,co+j)+=A(i,j) for 0<=i<A.Height, 0<=j<A.Width.
Definition at line 2246 of file densemat.cpp.
void mfem::DenseMatrix::AddMatrix | ( | double | a, |
DenseMatrix & | A, | ||
int | ro, | ||
int | co | ||
) |
Perform (ro+i,co+j)+=a*A(i,j) for 0<=i<A.Height, 0<=j<A.Width.
Definition at line 2272 of file densemat.cpp.
y += A.x
Definition at line 184 of file densemat.cpp.
void mfem::DenseMatrix::AddToVector | ( | int | offset, |
Vector & | v | ||
) | const |
Add the matrix 'data' to the Vector 'v' at the given 'offset'.
Definition at line 2298 of file densemat.cpp.
void mfem::DenseMatrix::AdjustDofDirection | ( | Array< int > & | dofs | ) |
If (dofs[i] < 0 and dofs[j] >= 0) or (dofs[i] >= 0 and dofs[j] < 0) then (*this)(i,j) = -(*this)(i,j).
Definition at line 2316 of file densemat.cpp.
void mfem::DenseMatrix::CalcEigenvalues | ( | double * | lambda, |
double * | vec | ||
) | const |
Return the eigenvalues (in increasing order) and eigenvectors of a 2x2 or 3x3 symmetric matrix.
Definition at line 1739 of file densemat.cpp.
double mfem::DenseMatrix::CalcSingularvalue | ( | const int | i | ) | const |
Return the i-th singular value (decreasing order) of NxN matrix, N=1,2,3.
Definition at line 1451 of file densemat.cpp.
|
inline |
Count the number of entries in the matrix for which isfinite is false, i.e. the entry is a NaN or +/-Inf.
Definition at line 254 of file densemat.hpp.
|
inline |
Definition at line 57 of file densemat.hpp.
void mfem::DenseMatrix::CopyCols | ( | DenseMatrix & | A, |
int | col1, | ||
int | col2 | ||
) |
Copy columns col1 through col2 from A to *this.
Definition at line 2179 of file densemat.cpp.
void mfem::DenseMatrix::CopyMN | ( | DenseMatrix & | A, |
int | m, | ||
int | n, | ||
int | Aro, | ||
int | Aco | ||
) |
Copy the m x n submatrix of A at (Aro)ffset, (Aco)loffset to *this.
Definition at line 2188 of file densemat.cpp.
void mfem::DenseMatrix::CopyMN | ( | DenseMatrix & | A, |
int | row_offset, | ||
int | col_offset | ||
) |
Copy matrix A to the location in *this at row_offset, col_offset.
Definition at line 2199 of file densemat.cpp.
void mfem::DenseMatrix::CopyMNDiag | ( | double | c, |
int | n, | ||
int | row_offset, | ||
int | col_offset | ||
) |
Copy c on the diagonal of size n to *this at row_offset, col_offset.
Definition at line 2219 of file densemat.cpp.
void mfem::DenseMatrix::CopyMNDiag | ( | double * | diag, |
int | n, | ||
int | row_offset, | ||
int | col_offset | ||
) |
Copy diag on the diagonal of size n to *this at row_offset, col_offset.
Definition at line 2232 of file densemat.cpp.
void mfem::DenseMatrix::CopyMNt | ( | DenseMatrix & | A, |
int | row_offset, | ||
int | col_offset | ||
) |
Copy matrix A^t to the location in *this at row_offset, col_offset.
Definition at line 2209 of file densemat.cpp.
void mfem::DenseMatrix::CopyRows | ( | DenseMatrix & | A, |
int | row1, | ||
int | row2 | ||
) |
Copy rows row1 through row2 from A to *this.
Definition at line 2170 of file densemat.cpp.
|
inline |
Returns vector of the elements.
Definition at line 69 of file densemat.hpp.
double mfem::DenseMatrix::Det | ( | ) | const |
Calculates the determinant of the matrix (for 2x2 or 3x3 matrices)
Definition at line 321 of file densemat.cpp.
void mfem::DenseMatrix::Diag | ( | double | c, |
int | n | ||
) |
Creates n x n diagonal matrix with diagonal elements c.
Definition at line 2014 of file densemat.cpp.
void mfem::DenseMatrix::Diag | ( | double * | diag, |
int | n | ||
) |
Creates n x n diagonal matrix with diagonal given by diag.
Definition at line 2025 of file densemat.cpp.
|
inline |
Definition at line 171 of file densemat.hpp.
|
inline |
Definition at line 165 of file densemat.hpp.
|
inline |
Definition at line 168 of file densemat.hpp.
|
virtual |
|
virtual |
Returns constant reference to a_{ij}.
Implements mfem::Matrix.
Definition at line 115 of file densemat.cpp.
double mfem::DenseMatrix::FNorm | ( | ) | const |
Compute the Frobenius norm of the matrix.
Definition at line 578 of file densemat.cpp.
void mfem::DenseMatrix::GetColumn | ( | int | c, |
Vector & | col | ||
) |
Definition at line 1977 of file densemat.cpp.
|
inline |
Definition at line 186 of file densemat.hpp.
void mfem::DenseMatrix::GetDiag | ( | Vector & | d | ) |
Returns the diagonal of the matrix.
Definition at line 1991 of file densemat.cpp.
void mfem::DenseMatrix::GetFromVector | ( | int | offset, |
const Vector & | v | ||
) |
Get the matrix 'data' from the Vector 'v' at the given 'offset'.
Definition at line 2307 of file densemat.cpp.
void mfem::DenseMatrix::Getl1Diag | ( | Vector & | l | ) |
Returns the l1 norm of the rows of the matrix v_i = sum_j |a_ij|.
Definition at line 2001 of file densemat.cpp.
void mfem::DenseMatrix::GradToCurl | ( | DenseMatrix & | curl | ) |
Given a DShape matrix (from a scalar FE), stored in this, returns the CurlShape matrix. If *this is a N by D matrix, then curl is a D*N by D(D-1)/2 matrix. The size of curl must be set outside. The dimension D can be either 2 or 3.
Definition at line 2096 of file densemat.cpp.
void mfem::DenseMatrix::GradToDiv | ( | Vector & | div | ) |
Given a DShape matrix (from a scalar FE), stored in *this, returns the DivShape vector. If *this is a N by dim matrix, then div is a dim*N vector. The size of div must be set outside.
Definition at line 2153 of file densemat.cpp.
double mfem::DenseMatrix::InnerProduct | ( | const double * | x, |
const double * | y | ||
) | const |
Compute y^t A x.
Definition at line 202 of file densemat.cpp.
Compute y^t A x.
Definition at line 121 of file densemat.hpp.
|
virtual |
Returns a pointer to the inverse matrix.
Implements mfem::Matrix.
Definition at line 316 of file densemat.cpp.
void mfem::DenseMatrix::Invert | ( | ) |
Replaces the current matrix with its inverse.
Definition at line 465 of file densemat.cpp.
void mfem::DenseMatrix::InvLeftScaling | ( | const Vector & | s | ) |
InvLeftScaling this = diag(1./s) * this.
Definition at line 227 of file densemat.cpp.
void mfem::DenseMatrix::InvRightScaling | ( | const Vector & | s | ) |
InvRightScaling: this = this * diag(1./s);.
Definition at line 249 of file densemat.cpp.
void mfem::DenseMatrix::InvSymmetricScaling | ( | const Vector & | s | ) |
InvSymmetricScaling this = diag(sqrt(1./s)) * this * diag(sqrt(1./s))
Definition at line 282 of file densemat.cpp.
void mfem::DenseMatrix::LeftScaling | ( | const Vector & | s | ) |
LeftScaling this = diag(s) * this.
Definition at line 218 of file densemat.cpp.
void mfem::DenseMatrix::Lump | ( | ) |
Definition at line 2082 of file densemat.cpp.
double mfem::DenseMatrix::MaxMaxNorm | ( | ) | const |
Compute the norm ||A|| = max_{ij} |A_{ij}|.
Definition at line 562 of file densemat.cpp.
void mfem::DenseMatrix::Mult | ( | const double * | x, |
double * | y | ||
) | const |
Matrix vector multiplication.
Definition at line 120 of file densemat.cpp.
Matrix vector multiplication.
Implements mfem::Operator.
Definition at line 136 of file densemat.cpp.
void mfem::DenseMatrix::MultTranspose | ( | const double * | x, |
double * | y | ||
) | const |
Multiply a vector with the transpose matrix.
Definition at line 161 of file densemat.cpp.
Multiply a vector with the transpose matrix.
Reimplemented from mfem::Operator.
Definition at line 174 of file densemat.cpp.
void mfem::DenseMatrix::Neg | ( | ) |
(*this) = -(*this)
Definition at line 447 of file densemat.cpp.
void mfem::DenseMatrix::Norm2 | ( | double * | v | ) | const |
Take the 2-norm of the columns of A and store in v.
Definition at line 551 of file densemat.cpp.
|
inline |
Returns reference to a_{ij}.
Definition at line 475 of file densemat.hpp.
|
inline |
Returns constant reference to a_{ij}.
Definition at line 485 of file densemat.hpp.
double mfem::DenseMatrix::operator* | ( | const DenseMatrix & | m | ) | const |
Matrix inner product: tr(A^t B)
Definition at line 146 of file densemat.cpp.
DenseMatrix & mfem::DenseMatrix::operator*= | ( | double | c | ) |
Definition at line 438 of file densemat.cpp.
DenseMatrix & mfem::DenseMatrix::operator+= | ( | DenseMatrix & | m | ) |
Definition at line 415 of file densemat.cpp.
DenseMatrix & mfem::DenseMatrix::operator-= | ( | DenseMatrix & | m | ) |
Definition at line 427 of file densemat.cpp.
DenseMatrix & mfem::DenseMatrix::operator= | ( | double | c | ) |
Sets the matrix elements equal to constant c.
Definition at line 382 of file densemat.cpp.
DenseMatrix & mfem::DenseMatrix::operator= | ( | const double * | d | ) |
Copy the matrix entries from the given array.
Definition at line 391 of file densemat.cpp.
DenseMatrix & mfem::DenseMatrix::operator= | ( | const DenseMatrix & | m | ) |
Sets the matrix size and elements equal to those of m.
Definition at line 402 of file densemat.cpp.
|
virtual |
Prints matrix to stream out.
Reimplemented from mfem::Matrix.
Definition at line 2353 of file densemat.cpp.
|
virtual |
Definition at line 2371 of file densemat.cpp.
|
virtual |
Prints the transpose matrix to stream out.
Definition at line 2386 of file densemat.cpp.
int mfem::DenseMatrix::Rank | ( | double | tol | ) | const |
Definition at line 868 of file densemat.cpp.
void mfem::DenseMatrix::RightScaling | ( | const Vector & | s | ) |
RightScaling: this = this * diag(s);.
Definition at line 236 of file densemat.cpp.
void mfem::DenseMatrix::SetCol | ( | int | col, |
double | value | ||
) |
Set all entries of a column to the specified value.
Definition at line 2347 of file densemat.cpp.
void mfem::DenseMatrix::SetRow | ( | int | row, |
double | value | ||
) |
Set all entries of a row to the specified value.
Definition at line 2341 of file densemat.cpp.
void mfem::DenseMatrix::SetSize | ( | int | s | ) |
If the matrix is not a square matrix of size s then recreate it.
Definition at line 73 of file densemat.cpp.
void mfem::DenseMatrix::SetSize | ( | int | h, |
int | w | ||
) |
If the matrix is not a matrix of size (h x w) then recreate it.
Definition at line 91 of file densemat.cpp.
void mfem::DenseMatrix::SingularValues | ( | Vector & | sv | ) | const |
Definition at line 829 of file densemat.cpp.
|
inline |
For backward compatibility define Size to be synonym of Width()
Definition at line 60 of file densemat.hpp.
void mfem::DenseMatrix::SymmetricScaling | ( | const Vector & | s | ) |
SymmetricScaling this = diag(sqrt(s)) * this * diag(sqrt(s))
Definition at line 262 of file densemat.cpp.
void mfem::DenseMatrix::Symmetrize | ( | ) |
(*this) = 1/2 ((*this) + (*this)^t)
Definition at line 2067 of file densemat.cpp.
void mfem::DenseMatrix::TestInversion | ( | ) |
Invert and print the numerical conditioning of the inversion.
Definition at line 2404 of file densemat.cpp.
double mfem::DenseMatrix::Trace | ( | ) | const |
Trace of a square matrix.
Definition at line 301 of file densemat.cpp.
void mfem::DenseMatrix::Transpose | ( | ) |
(*this) = (*this)^t
Definition at line 2036 of file densemat.cpp.
void mfem::DenseMatrix::Transpose | ( | DenseMatrix & | A | ) |
(*this) = A^t
Definition at line 2058 of file densemat.cpp.
|
inline |
Definition at line 54 of file densemat.hpp.
double mfem::DenseMatrix::Weight | ( | ) | const |
Definition at line 348 of file densemat.cpp.
|
friend |
Definition at line 29 of file densemat.hpp.
|
friend |
Definition at line 24 of file densemat.hpp.
|
friend |
Matrix matrix multiplication. A = B * C.
Definition at line 2445 of file densemat.cpp.