MFEM v4.9.0
Finite element discretization library
Loading...
Searching...
No Matches
mfem::DenseMatrix Class Reference

Data type dense matrix using column-major storage. More...

#include <densemat.hpp>

Inheritance diagram for mfem::DenseMatrix:
[legend]
Collaboration diagram for mfem::DenseMatrix:
[legend]

Public Member Functions

 DenseMatrix ()
 
 DenseMatrix (int s)
 Creates square matrix of size s.
 
 DenseMatrix (int m, int n)
 Creates rectangular matrix of size m x n.
 
 DenseMatrix (const DenseMatrix &mat, char ch)
 Creates rectangular matrix equal to the transpose of mat.
 
 DenseMatrix (real_t *d, int h, int w)
 Construct a DenseMatrix using an existing data array.
 
 DenseMatrix (const DenseMatrix &)=default
 Copy constructor (deep copy).
 
 DenseMatrix (DenseMatrix &&)=default
 Move constructor.
 
DenseMatrixoperator= (const DenseMatrix &)=default
 Copy assignment (deep copy).
 
DenseMatrixoperator= (DenseMatrix &&)=default
 Move assignment.
 
template<int M, int N, typename T = real_t>
 DenseMatrix (const T(&values)[M][N])
 
void UseExternalData (real_t *d, int h, int w)
 Change the data array and the size of the DenseMatrix.
 
void Reset (real_t *d, int h, int w)
 Change the data array and the size of the DenseMatrix.
 
void ClearExternalData ()
 
void Clear ()
 Delete the matrix data array (if owned) and reset the matrix state.
 
int Size () const
 For backward compatibility define Size to be synonym of Width()
 
int TotalSize () const
 Total size = width*height.
 
void SetSize (int s)
 Change the size of the DenseMatrix to s x s.
 
void SetSize (int h, int w)
 Change the size of the DenseMatrix to h x w.
 
real_tData () const
 Returns the matrix data array. Warning: this method casts away constness.
 
real_tGetData () const
 Returns the matrix data array. Warning: this method casts away constness.
 
Memory< real_t > & GetMemory ()
 
const Memory< real_t > & GetMemory () const
 
bool OwnsData () const
 Return the DenseMatrix data (host pointer) ownership flag.
 
real_toperator() (int i, int j)
 Returns reference to a_{ij}.
 
const real_toperator() (int i, int j) const
 Returns constant reference to a_{ij}.
 
real_t operator* (const DenseMatrix &m) const
 Matrix inner product: tr(A^t B)
 
real_t Trace () const
 Trace of a square matrix.
 
real_tElem (int i, int j) override
 Returns reference to a_{ij}.
 
const real_tElem (int i, int j) const override
 Returns constant reference to a_{ij}.
 
void Mult (const real_t *x, real_t *y) const
 Matrix vector multiplication.
 
void Mult (const real_t *x, Vector &y) const
 Matrix vector multiplication.
 
void Mult (const Vector &x, real_t *y) const
 Matrix vector multiplication.
 
void Mult (const Vector &x, Vector &y) const override
 Matrix vector multiplication.
 
void AbsMult (const Vector &x, Vector &y) const override
 Absolute-value matrix vector multiplication.
 
void MultTranspose (const real_t *x, real_t *y) const
 Multiply a vector with the transpose matrix.
 
void MultTranspose (const real_t *x, Vector &y) const
 Multiply a vector with the transpose matrix.
 
void MultTranspose (const Vector &x, real_t *y) const
 Multiply a vector with the transpose matrix.
 
void MultTranspose (const Vector &x, Vector &y) const override
 Multiply a vector with the transpose matrix.
 
void AbsMultTranspose (const Vector &x, Vector &y) const override
 Multiply a vector with the absolute-value transpose matrix.
 
void AddMult (const Vector &x, Vector &y, const real_t a=1.0) const override
 y += a * A.x
 
void AddMultTranspose (const Vector &x, Vector &y, const real_t a=1.0) const override
 y += a * A^t x
 
void AddMult_a (real_t a, const Vector &x, Vector &y) const
 y += a * A.x
 
void AddMultTranspose_a (real_t a, const Vector &x, Vector &y) const
 y += a * A^t x
 
real_t InnerProduct (const real_t *x, const real_t *y) const
 Compute y^t A x.
 
void LeftScaling (const Vector &s)
 LeftScaling this = diag(s) * this.
 
void InvLeftScaling (const Vector &s)
 InvLeftScaling this = diag(1./s) * this.
 
void RightScaling (const Vector &s)
 RightScaling: this = this * diag(s);.
 
void InvRightScaling (const Vector &s)
 InvRightScaling: this = this * diag(1./s);.
 
void SymmetricScaling (const Vector &s)
 SymmetricScaling this = diag(sqrt(s)) * this * diag(sqrt(s))
 
void InvSymmetricScaling (const Vector &s)
 InvSymmetricScaling this = diag(sqrt(1./s)) * this * diag(sqrt(1./s))
 
real_t InnerProduct (const Vector &x, const Vector &y) const
 Compute y^t A x.
 
MatrixInverseInverse () const override
 Returns a pointer to the inverse matrix.
 
void Invert ()
 Replaces the current matrix with its inverse.
 
void SquareRootInverse ()
 Replaces the current matrix with its square root inverse.
 
void Exponential ()
 
real_t Det () const
 
real_t Weight () const
 
void Set (real_t alpha, const real_t *A)
 Set the matrix to alpha * A, assuming that A has the same dimensions as the matrix and uses column-major layout.
 
void Set (real_t alpha, const DenseMatrix &A)
 Set the matrix to alpha * A.
 
void Add (const real_t c, const DenseMatrix &A)
 Adds the matrix A multiplied by the number c to the matrix.
 
void Add (const real_t c, const real_t *A)
 
DenseMatrixoperator= (real_t c)
 Sets the matrix elements equal to constant c.
 
DenseMatrixoperator= (const real_t *d)
 Copy the matrix entries from the given array.
 
DenseMatrixoperator+= (const real_t *m)
 
DenseMatrixoperator+= (const DenseMatrix &m)
 
DenseMatrixoperator-= (const DenseMatrix &m)
 
DenseMatrixoperator*= (real_t c)
 
void Neg ()
 (*this) = -(*this)
 
void Norm2 (real_t *v) const
 Take the 2-norm of the columns of A and store in v.
 
void Norm2 (Vector &v) const
 Take the 2-norm of the columns of A and store in v.
 
real_t MaxMaxNorm () const
 Compute the norm ||A|| = max_{ij} |A_{ij}|.
 
real_t FNorm () const
 Compute the Frobenius norm of the matrix.
 
real_t FNorm2 () const
 Compute the square of the Frobenius norm of the matrix.
 
void Eigenvalues (Vector &ev)
 Compute eigenvalues of A x = ev x where A = *this.
 
void Eigenvalues (Vector &ev, DenseMatrix &evect)
 Compute eigenvalues and eigenvectors of A x = ev x where A = *this.
 
void Eigensystem (Vector &ev, DenseMatrix &evect)
 Compute eigenvalues and eigenvectors of A x = ev x where A = *this.
 
void Eigenvalues (DenseMatrix &b, Vector &ev)
 
void Eigenvalues (DenseMatrix &b, Vector &ev, DenseMatrix &evect)
 Compute generalized eigenvalues of A x = ev B x, where A = *this.
 
void Eigensystem (DenseMatrix &b, Vector &ev, DenseMatrix &evect)
 
void SingularValues (Vector &sv) const
 
int Rank (real_t tol) const
 
real_t CalcSingularvalue (const int i) const
 Return the i-th singular value (decreasing order) of NxN matrix, N=1,2,3.
 
void CalcEigenvalues (real_t *lambda, real_t *vec) const
 
void GetRow (int r, Vector &row) const
 
void GetColumn (int c, Vector &col) const
 
real_tGetColumn (int col)
 
const real_tGetColumn (int col) const
 
void GetColumnReference (int c, Vector &col)
 
void SetRow (int r, const real_t *row)
 
void SetRow (int r, const Vector &row)
 
void SetCol (int c, const real_t *col)
 
void SetCol (int c, const Vector &col)
 
void SetRow (int row, real_t value)
 Set all entries of a row to the specified value.
 
void SetCol (int col, real_t value)
 Set all entries of a column to the specified value.
 
void GetDiag (Vector &d) const
 Returns the diagonal of the matrix.
 
void Getl1Diag (Vector &l) const
 Returns the l1 norm of the rows of the matrix v_i = sum_j |a_ij|.
 
void GetRowSums (Vector &l) const
 Compute the row sums of the DenseMatrix.
 
void Diag (real_t c, int n)
 Creates n x n diagonal matrix with diagonal elements c.
 
void Diag (real_t *diag, int n)
 Creates n x n diagonal matrix with diagonal given by diag.
 
void Transpose ()
 (*this) = (*this)^t
 
void Transpose (const DenseMatrix &A)
 (*this) = A^t
 
void Symmetrize ()
 (*this) = 1/2 ((*this) + (*this)^t)
 
void Lump ()
 
void GradToCurl (DenseMatrix &curl)
 
void GradToVectorCurl2D (DenseMatrix &curl)
 
void GradToDiv (Vector &div)
 
void CopyRows (const DenseMatrix &A, int row1, int row2)
 Copy rows row1 through row2 from A to *this.
 
void CopyCols (const DenseMatrix &A, int col1, int col2)
 Copy columns col1 through col2 from A to *this.
 
void CopyMN (const DenseMatrix &A, int m, int n, int Aro, int Aco)
 Copy the m x n submatrix of A at row/col offsets Aro/Aco to *this.
 
void CopyMN (const DenseMatrix &A, int row_offset, int col_offset)
 Copy matrix A to the location in *this at row_offset, col_offset.
 
void CopyMNt (const DenseMatrix &A, int row_offset, int col_offset)
 Copy matrix A^t to the location in *this at row_offset, col_offset.
 
void CopyMN (const DenseMatrix &A, int m, int n, int Aro, int Aco, int row_offset, int col_offset)
 
void CopyMNDiag (real_t c, int n, int row_offset, int col_offset)
 Copy c on the diagonal of size n to *this at row_offset, col_offset.
 
void CopyMNDiag (real_t *diag, int n, int row_offset, int col_offset)
 Copy diag on the diagonal of size n to *this at row_offset, col_offset.
 
void CopyExceptMN (const DenseMatrix &A, int m, int n)
 Copy All rows and columns except m and n from A.
 
void AddMatrix (DenseMatrix &A, int ro, int co)
 Perform (ro+i,co+j)+=A(i,j) for 0<=i.
 
void AddMatrix (real_t a, const DenseMatrix &A, int ro, int co)
 Perform (ro+i,co+j)+=a*A(i,j) for 0<=i.
 
void GetSubMatrix (const Array< int > &idx, DenseMatrix &A) const
 
void GetSubMatrix (const Array< int > &idx_i, const Array< int > &idx_j, DenseMatrix &A) const
 
void GetSubMatrix (int ibeg, int iend, DenseMatrix &A)
 
void GetSubMatrix (int ibeg, int iend, int jbeg, int jend, DenseMatrix &A)
 
void SetSubMatrix (const Array< int > &idx, const DenseMatrix &A)
 Set (*this)(idx[i],idx[j]) = A(i,j)
 
void SetSubMatrix (const Array< int > &idx_i, const Array< int > &idx_j, const DenseMatrix &A)
 Set (*this)(idx_i[i],idx_j[j]) = A(i,j)
 
void SetSubMatrix (int ibeg, const DenseMatrix &A)
 
void SetSubMatrix (int ibeg, int jbeg, const DenseMatrix &A)
 
void AddSubMatrix (const Array< int > &idx, const DenseMatrix &A)
 (*this)(idx[i],idx[j]) += A(i,j)
 
void AddSubMatrix (const Array< int > &idx_i, const Array< int > &idx_j, const DenseMatrix &A)
 (*this)(idx_i[i],idx_j[j]) += A(i,j)
 
void AddSubMatrix (int ibeg, const DenseMatrix &A)
 
void AddSubMatrix (int ibeg, int jbeg, const DenseMatrix &A)
 
void AddToVector (int offset, Vector &v) const
 Add the matrix 'data' to the Vector 'v' at the given 'offset'.
 
void GetFromVector (int offset, const Vector &v)
 Get the matrix 'data' from the Vector 'v' at the given 'offset'.
 
void AdjustDofDirection (Array< int > &dofs)
 
void Threshold (real_t eps)
 Replace small entries, abs(a_ij) <= eps, with zero.
 
int CheckFinite () const
 
void Print (std::ostream &out=mfem::out, int width_=4) const override
 Prints matrix to stream out.
 
void PrintMatlab (std::ostream &out=mfem::out) const override
 Prints operator in Matlab format.
 
virtual void PrintMathematica (std::ostream &out=mfem::out) const
 
virtual void PrintT (std::ostream &out=mfem::out, int width_=4) const
 Prints the transpose matrix to stream out.
 
void TestInversion ()
 Invert and print the numerical conditioning of the inversion.
 
std::size_t MemoryUsage () const
 
const real_tRead (bool on_dev=true) const
 Shortcut for mfem::Read( GetMemory(), TotalSize(), on_dev).
 
const real_tHostRead () const
 Shortcut for mfem::Read(GetMemory(), TotalSize(), false).
 
real_tWrite (bool on_dev=true)
 Shortcut for mfem::Write(GetMemory(), TotalSize(), on_dev).
 
real_tHostWrite ()
 Shortcut for mfem::Write(GetMemory(), TotalSize(), false).
 
real_tReadWrite (bool on_dev=true)
 Shortcut for mfem::ReadWrite(GetMemory(), TotalSize(), on_dev).
 
real_tHostReadWrite ()
 Shortcut for mfem::ReadWrite(GetMemory(), TotalSize(), false).
 
void Swap (DenseMatrix &other)
 
- Public Member Functions inherited from mfem::Matrix
 Matrix (int s)
 Creates a square matrix of size s.
 
 Matrix (int h, int w)
 Creates a matrix of the given height and width.
 
bool IsSquare () const
 Returns whether the matrix is a square matrix.
 
virtual void Finalize (int)
 Finalizes the matrix initialization.
 
virtual ~Matrix ()
 Destroys matrix.
 
- Public Member Functions inherited from mfem::Operator
void InitTVectors (const Operator *Po, const Operator *Ri, const Operator *Pi, Vector &x, Vector &b, Vector &X, Vector &B) const
 Initializes memory for true vectors of linear system.
 
 Operator (int s=0)
 Construct a square Operator with given size s (default 0).
 
 Operator (int h, int w)
 Construct an Operator with the given height (output size) and width (input size).
 
int Height () const
 Get the height (size of output) of the Operator. Synonym with NumRows().
 
int NumRows () const
 Get the number of rows (size of output) of the Operator. Synonym with Height().
 
int Width () const
 Get the width (size of input) of the Operator. Synonym with NumCols().
 
int NumCols () const
 Get the number of columns (size of input) of the Operator. Synonym with Width().
 
virtual MemoryClass GetMemoryClass () const
 Return the MemoryClass preferred by the Operator.
 
virtual void ArrayMult (const Array< const Vector * > &X, Array< Vector * > &Y) const
 Operator application on a matrix: Y=A(X).
 
virtual void ArrayMultTranspose (const Array< const Vector * > &X, Array< Vector * > &Y) const
 Action of the transpose operator on a matrix: Y=A^t(X).
 
virtual void ArrayAddMult (const Array< const Vector * > &X, Array< Vector * > &Y, const real_t a=1.0) const
 Operator application on a matrix: Y+=A(X) (default) or Y+=a*A(X).
 
virtual void ArrayAddMultTranspose (const Array< const Vector * > &X, Array< Vector * > &Y, const real_t a=1.0) const
 Operator transpose application on a matrix: Y+=A^t(X) (default) or Y+=a*A^t(X).
 
virtual OperatorGetGradient (const Vector &x) const
 Evaluate the gradient operator at the point x. The default behavior in class Operator is to generate an error.
 
virtual void AssembleDiagonal (Vector &diag) const
 Computes the diagonal entries into diag. Typically, this operation only makes sense for linear Operators. In some cases, only an approximation of the diagonal is computed.
 
virtual const OperatorGetProlongation () const
 Prolongation operator from linear algebra (linear system) vectors, to input vectors for the operator. NULL means identity.
 
virtual const OperatorGetRestriction () const
 Restriction operator from input vectors for the operator to linear algebra (linear system) vectors. NULL means identity.
 
virtual const OperatorGetOutputProlongation () const
 Prolongation operator from linear algebra (linear system) vectors, to output vectors for the operator. NULL means identity.
 
virtual const OperatorGetOutputRestrictionTranspose () const
 Transpose of GetOutputRestriction, directly available in this form to facilitate matrix-free RAP-type operators.
 
virtual const OperatorGetOutputRestriction () const
 Restriction operator from output vectors for the operator to linear algebra (linear system) vectors. NULL means identity.
 
void FormLinearSystem (const Array< int > &ess_tdof_list, Vector &x, Vector &b, Operator *&A, Vector &X, Vector &B, int copy_interior=0)
 Form a constrained linear system using a matrix-free approach.
 
void FormRectangularLinearSystem (const Array< int > &trial_tdof_list, const Array< int > &test_tdof_list, Vector &x, Vector &b, Operator *&A, Vector &X, Vector &B)
 Form a column-constrained linear system using a matrix-free approach.
 
virtual void RecoverFEMSolution (const Vector &X, const Vector &b, Vector &x)
 Reconstruct a solution vector x (e.g. a GridFunction) from the solution X of a constrained linear system obtained from Operator::FormLinearSystem() or Operator::FormRectangularLinearSystem().
 
void FormSystemOperator (const Array< int > &ess_tdof_list, Operator *&A)
 Return in A a parallel (on truedofs) version of this square operator.
 
void FormRectangularSystemOperator (const Array< int > &trial_tdof_list, const Array< int > &test_tdof_list, Operator *&A)
 Return in A a parallel (on truedofs) version of this rectangular operator (including constraints).
 
void FormDiscreteOperator (Operator *&A)
 Return in A a parallel (on truedofs) version of this rectangular operator.
 
void PrintMatlab (std::ostream &out, int n, int m=0) const
 Prints operator with input size n and output size m in Matlab format.
 
virtual ~Operator ()
 Virtual destructor.
 
Type GetType () const
 Return the type ID of the Operator class.
 

Friends

class DenseTensor
 
class DenseMatrixInverse
 

Additional Inherited Members

- Public Types inherited from mfem::Operator
enum  DiagonalPolicy { DIAG_ZERO , DIAG_ONE , DIAG_KEEP }
 Defines operator diagonal policy upon elimination of rows and/or columns. More...
 
enum  Type {
  ANY_TYPE , MFEM_SPARSEMAT , Hypre_ParCSR , PETSC_MATAIJ ,
  PETSC_MATIS , PETSC_MATSHELL , PETSC_MATNEST , PETSC_MATHYPRE ,
  PETSC_MATGENERIC , Complex_Operator , MFEM_ComplexSparseMat , Complex_Hypre_ParCSR ,
  Complex_DenseMat , MFEM_Block_Matrix , MFEM_Block_Operator
}
 Enumeration defining IDs for some classes derived from Operator. More...
 
- Protected Member Functions inherited from mfem::Operator
void FormConstrainedSystemOperator (const Array< int > &ess_tdof_list, ConstrainedOperator *&Aout)
 see FormSystemOperator()
 
void FormRectangularConstrainedSystemOperator (const Array< int > &trial_tdof_list, const Array< int > &test_tdof_list, RectangularConstrainedOperator *&Aout)
 see FormRectangularSystemOperator()
 
OperatorSetupRAP (const Operator *Pi, const Operator *Po)
 Returns RAP Operator of this, using input/output Prolongation matrices Pi corresponds to "P", Po corresponds to "Rt".
 
- Protected Attributes inherited from mfem::Operator
int height
 Dimension of the output / number of rows in the matrix.
 
int width
 Dimension of the input / number of columns in the matrix.
 

Detailed Description

Data type dense matrix using column-major storage.

Definition at line 23 of file densemat.hpp.

Constructor & Destructor Documentation

◆ DenseMatrix() [1/8]

mfem::DenseMatrix::DenseMatrix ( )

Default constructor for DenseMatrix. Sets data = NULL and height = width = 0.

Definition at line 42 of file densemat.cpp.

◆ DenseMatrix() [2/8]

mfem::DenseMatrix::DenseMatrix ( int s)
explicit

Creates square matrix of size s.

Definition at line 44 of file densemat.cpp.

◆ DenseMatrix() [3/8]

mfem::DenseMatrix::DenseMatrix ( int m,
int n )

Creates rectangular matrix of size m x n.

Definition at line 54 of file densemat.cpp.

◆ DenseMatrix() [4/8]

mfem::DenseMatrix::DenseMatrix ( const DenseMatrix & mat,
char ch )

Creates rectangular matrix equal to the transpose of mat.

Definition at line 66 of file densemat.cpp.

◆ DenseMatrix() [5/8]

mfem::DenseMatrix::DenseMatrix ( real_t * d,
int h,
int w )
inline

Construct a DenseMatrix using an existing data array.

The DenseMatrix does not assume ownership of the data array, i.e. it will not delete the array.

Definition at line 55 of file densemat.hpp.

◆ DenseMatrix() [6/8]

mfem::DenseMatrix::DenseMatrix ( const DenseMatrix & )
default

Copy constructor (deep copy).

◆ DenseMatrix() [7/8]

mfem::DenseMatrix::DenseMatrix ( DenseMatrix && )
default

Move constructor.

◆ DenseMatrix() [8/8]

template<int M, int N, typename T = real_t>
mfem::DenseMatrix::DenseMatrix ( const T(&) values[M][N])
inlineexplicit

Create a dense matrix using a braced initializer list The inner lists correspond to rows of the matrix

Definition at line 73 of file densemat.hpp.

Member Function Documentation

◆ AbsMult()

void mfem::DenseMatrix::AbsMult ( const Vector & x,
Vector & y ) const
overridevirtual

Absolute-value matrix vector multiplication.

Reimplemented from mfem::Operator.

Definition at line 135 of file densemat.cpp.

◆ AbsMultTranspose()

void mfem::DenseMatrix::AbsMultTranspose ( const Vector & x,
Vector & y ) const
overridevirtual

Multiply a vector with the absolute-value transpose matrix.

Reimplemented from mfem::Operator.

Definition at line 185 of file densemat.cpp.

◆ Add() [1/2]

void mfem::DenseMatrix::Add ( const real_t c,
const DenseMatrix & A )

Adds the matrix A multiplied by the number c to the matrix.

Definition at line 589 of file densemat.cpp.

◆ Add() [2/2]

void mfem::DenseMatrix::Add ( const real_t c,
const real_t * A )

Adds the matrix A multiplied by the number c to the matrix, assuming A has the same dimensions and uses column-major layout.

Definition at line 600 of file densemat.cpp.

◆ AddMatrix() [1/2]

void mfem::DenseMatrix::AddMatrix ( DenseMatrix & A,
int ro,
int co )

Perform (ro+i,co+j)+=A(i,j) for 0<=i.

Definition at line 1720 of file densemat.cpp.

◆ AddMatrix() [2/2]

void mfem::DenseMatrix::AddMatrix ( real_t a,
const DenseMatrix & A,
int ro,
int co )

Perform (ro+i,co+j)+=a*A(i,j) for 0<=i.

Definition at line 1750 of file densemat.cpp.

◆ AddMult()

void mfem::DenseMatrix::AddMult ( const Vector & x,
Vector & y,
const real_t a = 1.0 ) const
overridevirtual

y += a * A.x

Reimplemented from mfem::Operator.

Definition at line 194 of file densemat.cpp.

◆ AddMult_a()

void mfem::DenseMatrix::AddMult_a ( real_t a,
const Vector & x,
Vector & y ) const

y += a * A.x

Definition at line 241 of file densemat.cpp.

◆ AddMultTranspose()

void mfem::DenseMatrix::AddMultTranspose ( const Vector & x,
Vector & y,
const real_t a = 1.0 ) const
overridevirtual

y += a * A^t x

Reimplemented from mfem::Operator.

Definition at line 217 of file densemat.cpp.

◆ AddMultTranspose_a()

void mfem::DenseMatrix::AddMultTranspose_a ( real_t a,
const Vector & x,
Vector & y ) const

y += a * A^t x

Definition at line 262 of file densemat.cpp.

◆ AddSubMatrix() [1/4]

void mfem::DenseMatrix::AddSubMatrix ( const Array< int > & idx,
const DenseMatrix & A )

(*this)(idx[i],idx[j]) += A(i,j)

Definition at line 1997 of file densemat.cpp.

◆ AddSubMatrix() [2/4]

void mfem::DenseMatrix::AddSubMatrix ( const Array< int > & idx_i,
const Array< int > & idx_j,
const DenseMatrix & A )

(*this)(idx_i[i],idx_j[j]) += A(i,j)

Definition at line 2025 of file densemat.cpp.

◆ AddSubMatrix() [3/4]

void mfem::DenseMatrix::AddSubMatrix ( int ibeg,
const DenseMatrix & A )

Add the submatrix A to this with row and column offset ibeg

Definition at line 2056 of file densemat.cpp.

◆ AddSubMatrix() [4/4]

void mfem::DenseMatrix::AddSubMatrix ( int ibeg,
int jbeg,
const DenseMatrix & A )

Add the submatrix A to this with row and column offsets ibeg and jbeg respectively

Definition at line 2082 of file densemat.cpp.

◆ AddToVector()

void mfem::DenseMatrix::AddToVector ( int offset,
Vector & v ) const

Add the matrix 'data' to the Vector 'v' at the given 'offset'.

Definition at line 2110 of file densemat.cpp.

◆ AdjustDofDirection()

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 2132 of file densemat.cpp.

◆ CalcEigenvalues()

void mfem::DenseMatrix::CalcEigenvalues ( real_t * lambda,
real_t * vec ) const

Return the eigenvalues (in increasing order) and eigenvectors of a 2x2 or 3x3 symmetric matrix.

Definition at line 1290 of file densemat.cpp.

◆ CalcSingularvalue()

real_t 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 1265 of file densemat.cpp.

◆ CheckFinite()

int mfem::DenseMatrix::CheckFinite ( ) const
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 471 of file densemat.hpp.

◆ Clear()

void mfem::DenseMatrix::Clear ( )
inline

Delete the matrix data array (if owned) and reset the matrix state.

Definition at line 106 of file densemat.hpp.

◆ ClearExternalData()

void mfem::DenseMatrix::ClearExternalData ( )
inline

Clear the data array and the dimensions of the DenseMatrix. This method should not be used with DenseMatrix that owns its current data array.

Definition at line 103 of file densemat.hpp.

◆ CopyCols()

void mfem::DenseMatrix::CopyCols ( const DenseMatrix & A,
int col1,
int col2 )

Copy columns col1 through col2 from A to *this.

Definition at line 1576 of file densemat.cpp.

◆ CopyExceptMN()

void mfem::DenseMatrix::CopyExceptMN ( const DenseMatrix & A,
int m,
int n )

Copy All rows and columns except m and n from A.

Definition at line 1694 of file densemat.cpp.

◆ CopyMN() [1/3]

void mfem::DenseMatrix::CopyMN ( const DenseMatrix & A,
int m,
int n,
int Aro,
int Aco )

Copy the m x n submatrix of A at row/col offsets Aro/Aco to *this.

Definition at line 1589 of file densemat.cpp.

◆ CopyMN() [2/3]

void mfem::DenseMatrix::CopyMN ( const DenseMatrix & A,
int m,
int n,
int Aro,
int Aco,
int row_offset,
int col_offset )

Copy the m x n submatrix of A at row/col offsets Aro/Aco to *this at row_offset, col_offset

Definition at line 1628 of file densemat.cpp.

◆ CopyMN() [3/3]

void mfem::DenseMatrix::CopyMN ( const DenseMatrix & A,
int row_offset,
int col_offset )

Copy matrix A to the location in *this at row_offset, col_offset.

Definition at line 1602 of file densemat.cpp.

◆ CopyMNDiag() [1/2]

void mfem::DenseMatrix::CopyMNDiag ( real_t * 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 1676 of file densemat.cpp.

◆ CopyMNDiag() [2/2]

void mfem::DenseMatrix::CopyMNDiag ( real_t 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 1659 of file densemat.cpp.

◆ CopyMNt()

void mfem::DenseMatrix::CopyMNt ( const 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 1615 of file densemat.cpp.

◆ CopyRows()

void mfem::DenseMatrix::CopyRows ( const DenseMatrix & A,
int row1,
int row2 )

Copy rows row1 through row2 from A to *this.

Definition at line 1563 of file densemat.cpp.

◆ Data()

real_t * mfem::DenseMatrix::Data ( ) const
inline

Returns the matrix data array. Warning: this method casts away constness.

Definition at line 122 of file densemat.hpp.

◆ Det()

real_t mfem::DenseMatrix::Det ( ) const

Calculates the determinant of the matrix (optimized for 2x2, 3x3, and 4x4 matrices)

Definition at line 496 of file densemat.cpp.

◆ Diag() [1/2]

void mfem::DenseMatrix::Diag ( real_t * diag,
int n )

Creates n x n diagonal matrix with diagonal given by diag.

Definition at line 1402 of file densemat.cpp.

◆ Diag() [2/2]

void mfem::DenseMatrix::Diag ( real_t c,
int n )

Creates n x n diagonal matrix with diagonal elements c.

Definition at line 1387 of file densemat.cpp.

◆ Eigensystem() [1/2]

void mfem::DenseMatrix::Eigensystem ( DenseMatrix & b,
Vector & ev,
DenseMatrix & evect )
inline

Compute generalized eigenvalues and eigenvectors of A x = ev B x, where A = *this

Definition at line 313 of file densemat.hpp.

◆ Eigensystem() [2/2]

void mfem::DenseMatrix::Eigensystem ( Vector & ev,
DenseMatrix & evect )
inline

Compute eigenvalues and eigenvectors of A x = ev x where A = *this.

Definition at line 299 of file densemat.hpp.

◆ Eigenvalues() [1/4]

void mfem::DenseMatrix::Eigenvalues ( DenseMatrix & b,
Vector & ev )
inline

Compute generalized eigenvalues and eigenvectors of A x = ev B x, where A = *this

Definition at line 304 of file densemat.hpp.

◆ Eigenvalues() [2/4]

void mfem::DenseMatrix::Eigenvalues ( DenseMatrix & b,
Vector & ev,
DenseMatrix & evect )
inline

Compute generalized eigenvalues of A x = ev B x, where A = *this.

Definition at line 308 of file densemat.hpp.

◆ Eigenvalues() [3/4]

void mfem::DenseMatrix::Eigenvalues ( Vector & ev)
inline

Compute eigenvalues of A x = ev x where A = *this.

Definition at line 291 of file densemat.hpp.

◆ Eigenvalues() [4/4]

void mfem::DenseMatrix::Eigenvalues ( Vector & ev,
DenseMatrix & evect )
inline

Compute eigenvalues and eigenvectors of A x = ev x where A = *this.

Definition at line 295 of file densemat.hpp.

◆ Elem() [1/2]

const real_t & mfem::DenseMatrix::Elem ( int i,
int j ) const
overridevirtual

Returns constant reference to a_{ij}.

Implements mfem::Matrix.

Definition at line 103 of file densemat.cpp.

◆ Elem() [2/2]

real_t & mfem::DenseMatrix::Elem ( int i,
int j )
overridevirtual

Returns reference to a_{ij}.

Implements mfem::Matrix.

Definition at line 98 of file densemat.cpp.

◆ Exponential()

void mfem::DenseMatrix::Exponential ( )

Replaces the current matrix with its exponential (currently only supports 2x2 matrices)

Formulas from Corollary 2.4 of doi:10.1109/9.233156 Note typo in the paper, in the prefactor in the equation under (i).

Definition at line 433 of file densemat.cpp.

◆ FNorm()

real_t mfem::DenseMatrix::FNorm ( ) const
inline

Compute the Frobenius norm of the matrix.

Definition at line 285 of file densemat.hpp.

◆ FNorm2()

real_t mfem::DenseMatrix::FNorm2 ( ) const
inline

Compute the square of the Frobenius norm of the matrix.

Definition at line 288 of file densemat.hpp.

◆ GetColumn() [1/3]

void mfem::DenseMatrix::GetColumn ( int c,
Vector & col ) const

Definition at line 1328 of file densemat.cpp.

◆ GetColumn() [2/3]

real_t * mfem::DenseMatrix::GetColumn ( int col)
inline

Definition at line 328 of file densemat.hpp.

◆ GetColumn() [3/3]

const real_t * mfem::DenseMatrix::GetColumn ( int col) const
inline

Definition at line 329 of file densemat.hpp.

◆ GetColumnReference()

void mfem::DenseMatrix::GetColumnReference ( int c,
Vector & col )
inline

Definition at line 331 of file densemat.hpp.

◆ GetData()

real_t * mfem::DenseMatrix::GetData ( ) const
inline

Returns the matrix data array. Warning: this method casts away constness.

Definition at line 126 of file densemat.hpp.

◆ GetDiag()

void mfem::DenseMatrix::GetDiag ( Vector & d) const

Returns the diagonal of the matrix.

Definition at line 1342 of file densemat.cpp.

◆ GetFromVector()

void mfem::DenseMatrix::GetFromVector ( int offset,
const Vector & v )

Get the matrix 'data' from the Vector 'v' at the given 'offset'.

Definition at line 2121 of file densemat.cpp.

◆ Getl1Diag()

void mfem::DenseMatrix::Getl1Diag ( Vector & l) const

Returns the l1 norm of the rows of the matrix v_i = sum_j |a_ij|.

Definition at line 1356 of file densemat.cpp.

◆ GetMemory() [1/2]

Memory< real_t > & mfem::DenseMatrix::GetMemory ( )
inline

Definition at line 128 of file densemat.hpp.

◆ GetMemory() [2/2]

const Memory< real_t > & mfem::DenseMatrix::GetMemory ( ) const
inline

Definition at line 130 of file densemat.hpp.

◆ GetRow()

void mfem::DenseMatrix::GetRow ( int r,
Vector & row ) const

Definition at line 1312 of file densemat.cpp.

◆ GetRowSums()

void mfem::DenseMatrix::GetRowSums ( Vector & l) const

Compute the row sums of the DenseMatrix.

Definition at line 1373 of file densemat.cpp.

◆ GetSubMatrix() [1/4]

void mfem::DenseMatrix::GetSubMatrix ( const Array< int > & idx,
DenseMatrix & A ) const

Get the square submatrix which corresponds to the given indices idx. Note: the A matrix will be resized to accommodate the data

Definition at line 1780 of file densemat.cpp.

◆ GetSubMatrix() [2/4]

void mfem::DenseMatrix::GetSubMatrix ( const Array< int > & idx_i,
const Array< int > & idx_j,
DenseMatrix & A ) const

Get the rectangular submatrix which corresponds to the given indices idx_i and idx_j. Note: the A matrix will be resized to accommodate the data

Definition at line 1801 of file densemat.cpp.

◆ GetSubMatrix() [3/4]

void mfem::DenseMatrix::GetSubMatrix ( int ibeg,
int iend,
DenseMatrix & A )

Get the square submatrix which corresponds to the range [ ibeg, iend ). Note: the A matrix will be resized to accommodate the data

Definition at line 1827 of file densemat.cpp.

◆ GetSubMatrix() [4/4]

void mfem::DenseMatrix::GetSubMatrix ( int ibeg,
int iend,
int jbeg,
int jend,
DenseMatrix & A )

Get the square submatrix which corresponds to the range i ∈ [ ibeg, iend ) and j ∈ [ jbeg, jend ) Note: the A matrix will be resized to accommodate the data

Definition at line 1851 of file densemat.cpp.

◆ GradToCurl()

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. In 2D this computes the scalar-valued curl of a 2D vector field

Definition at line 1475 of file densemat.cpp.

◆ GradToDiv()

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 1548 of file densemat.cpp.

◆ GradToVectorCurl2D()

void mfem::DenseMatrix::GradToVectorCurl2D ( DenseMatrix & curl)

Given a DShape matrix (from a scalar FE), stored in *this, returns the CurlShape matrix. This computes the vector-valued curl of a scalar field. this is N by 2 matrix and curl is N by 2 matrix as well.

Definition at line 1534 of file densemat.cpp.

◆ HostRead()

const real_t * mfem::DenseMatrix::HostRead ( ) const
inline

Shortcut for mfem::Read(GetMemory(), TotalSize(), false).

Definition at line 489 of file densemat.hpp.

◆ HostReadWrite()

real_t * mfem::DenseMatrix::HostReadWrite ( )
inline

Shortcut for mfem::ReadWrite(GetMemory(), TotalSize(), false).

Definition at line 501 of file densemat.hpp.

◆ HostWrite()

real_t * mfem::DenseMatrix::HostWrite ( )
inline

Shortcut for mfem::Write(GetMemory(), TotalSize(), false).

Definition at line 495 of file densemat.hpp.

◆ InnerProduct() [1/2]

real_t mfem::DenseMatrix::InnerProduct ( const real_t * x,
const real_t * y ) const

Compute y^t A x.

Definition at line 281 of file densemat.cpp.

◆ InnerProduct() [2/2]

real_t mfem::DenseMatrix::InnerProduct ( const Vector & x,
const Vector & y ) const
inline

Compute y^t A x.

Definition at line 216 of file densemat.hpp.

◆ Inverse()

MatrixInverse * mfem::DenseMatrix::Inverse ( ) const
overridevirtual

Returns a pointer to the inverse matrix.

Implements mfem::Matrix.

Definition at line 428 of file densemat.cpp.

◆ Invert()

void mfem::DenseMatrix::Invert ( )

Replaces the current matrix with its inverse.

Definition at line 674 of file densemat.cpp.

◆ InvLeftScaling()

void mfem::DenseMatrix::InvLeftScaling ( const Vector & s)

InvLeftScaling this = diag(1./s) * this.

Definition at line 312 of file densemat.cpp.

◆ InvRightScaling()

void mfem::DenseMatrix::InvRightScaling ( const Vector & s)

InvRightScaling: this = this * diag(1./s);.

Definition at line 340 of file densemat.cpp.

◆ InvSymmetricScaling()

void mfem::DenseMatrix::InvSymmetricScaling ( const Vector & s)

InvSymmetricScaling this = diag(sqrt(1./s)) * this * diag(sqrt(1./s))

Definition at line 382 of file densemat.cpp.

◆ LeftScaling()

void mfem::DenseMatrix::LeftScaling ( const Vector & s)

LeftScaling this = diag(s) * this.

Definition at line 299 of file densemat.cpp.

◆ Lump()

void mfem::DenseMatrix::Lump ( )

Definition at line 1461 of file densemat.cpp.

◆ MaxMaxNorm()

real_t mfem::DenseMatrix::MaxMaxNorm ( ) const

Compute the norm ||A|| = max_{ij} |A_{ij}|.

Definition at line 842 of file densemat.cpp.

◆ MemoryUsage()

std::size_t mfem::DenseMatrix::MemoryUsage ( ) const
inline

Definition at line 483 of file densemat.hpp.

◆ Mult() [1/4]

void mfem::DenseMatrix::Mult ( const real_t * x,
real_t * y ) const

Matrix vector multiplication.

Definition at line 108 of file densemat.cpp.

◆ Mult() [2/4]

void mfem::DenseMatrix::Mult ( const real_t * x,
Vector & y ) const

Matrix vector multiplication.

Definition at line 113 of file densemat.cpp.

◆ Mult() [3/4]

void mfem::DenseMatrix::Mult ( const Vector & x,
real_t * y ) const

Matrix vector multiplication.

Definition at line 120 of file densemat.cpp.

◆ Mult() [4/4]

void mfem::DenseMatrix::Mult ( const Vector & x,
Vector & y ) const
overridevirtual

Matrix vector multiplication.

Implements mfem::Operator.

Definition at line 127 of file densemat.cpp.

◆ MultTranspose() [1/4]

void mfem::DenseMatrix::MultTranspose ( const real_t * x,
real_t * y ) const

Multiply a vector with the transpose matrix.

Definition at line 158 of file densemat.cpp.

◆ MultTranspose() [2/4]

void mfem::DenseMatrix::MultTranspose ( const real_t * x,
Vector & y ) const

Multiply a vector with the transpose matrix.

Definition at line 163 of file densemat.cpp.

◆ MultTranspose() [3/4]

void mfem::DenseMatrix::MultTranspose ( const Vector & x,
real_t * y ) const

Multiply a vector with the transpose matrix.

Definition at line 170 of file densemat.cpp.

◆ MultTranspose() [4/4]

void mfem::DenseMatrix::MultTranspose ( const Vector & x,
Vector & y ) const
overridevirtual

Multiply a vector with the transpose matrix.

Reimplemented from mfem::Operator.

Definition at line 177 of file densemat.cpp.

◆ Neg()

void mfem::DenseMatrix::Neg ( )

(*this) = -(*this)

Definition at line 665 of file densemat.cpp.

◆ Norm2() [1/2]

void mfem::DenseMatrix::Norm2 ( real_t * v) const

Take the 2-norm of the columns of A and store in v.

Definition at line 829 of file densemat.cpp.

◆ Norm2() [2/2]

void mfem::DenseMatrix::Norm2 ( Vector & v) const
inline

Take the 2-norm of the columns of A and store in v.

Definition at line 275 of file densemat.hpp.

◆ operator()() [1/2]

real_t & mfem::DenseMatrix::operator() ( int i,
int j )
inline

Returns reference to a_{ij}.

Definition at line 1297 of file densemat.hpp.

◆ operator()() [2/2]

const real_t & mfem::DenseMatrix::operator() ( int i,
int j ) const
inline

Returns constant reference to a_{ij}.

Definition at line 1303 of file densemat.hpp.

◆ operator*()

real_t mfem::DenseMatrix::operator* ( const DenseMatrix & m) const

Matrix inner product: tr(A^t B)

Definition at line 143 of file densemat.cpp.

◆ operator*=()

DenseMatrix & mfem::DenseMatrix::operator*= ( real_t c)

Definition at line 655 of file densemat.cpp.

◆ operator+=() [1/2]

DenseMatrix & mfem::DenseMatrix::operator+= ( const DenseMatrix & m)

Definition at line 635 of file densemat.cpp.

◆ operator+=() [2/2]

DenseMatrix & mfem::DenseMatrix::operator+= ( const real_t * m)

Definition at line 629 of file densemat.cpp.

◆ operator-=()

DenseMatrix & mfem::DenseMatrix::operator-= ( const DenseMatrix & m)

Definition at line 642 of file densemat.cpp.

◆ operator=() [1/4]

DenseMatrix & mfem::DenseMatrix::operator= ( const DenseMatrix & )
default

Copy assignment (deep copy).

◆ operator=() [2/4]

DenseMatrix & mfem::DenseMatrix::operator= ( const real_t * d)

Copy the matrix entries from the given array.

Definition at line 619 of file densemat.cpp.

◆ operator=() [3/4]

DenseMatrix & mfem::DenseMatrix::operator= ( DenseMatrix && )
default

Move assignment.

◆ operator=() [4/4]

DenseMatrix & mfem::DenseMatrix::operator= ( real_t c)

Sets the matrix elements equal to constant c.

Definition at line 609 of file densemat.cpp.

◆ OwnsData()

bool mfem::DenseMatrix::OwnsData ( ) const
inline

Return the DenseMatrix data (host pointer) ownership flag.

Definition at line 133 of file densemat.hpp.

◆ Print()

void mfem::DenseMatrix::Print ( std::ostream & out = mfem::out,
int width_ = 4 ) const
overridevirtual

Prints matrix to stream out.

Reimplemented from mfem::Matrix.

Definition at line 2219 of file densemat.cpp.

◆ PrintMathematica()

void mfem::DenseMatrix::PrintMathematica ( std::ostream & out = mfem::out) const
virtual

Definition at line 2264 of file densemat.cpp.

◆ PrintMatlab()

void mfem::DenseMatrix::PrintMatlab ( std::ostream & out = mfem::out) const
overridevirtual

Prints operator in Matlab format.

Reimplemented from mfem::Operator.

Definition at line 2245 of file densemat.cpp.

◆ PrintT()

void mfem::DenseMatrix::PrintT ( std::ostream & out = mfem::out,
int width_ = 4 ) const
virtual

Prints the transpose matrix to stream out.

Definition at line 2293 of file densemat.cpp.

◆ Rank()

int mfem::DenseMatrix::Rank ( real_t tol) const

Definition at line 1250 of file densemat.cpp.

◆ Read()

const real_t * mfem::DenseMatrix::Read ( bool on_dev = true) const
inline

Shortcut for mfem::Read( GetMemory(), TotalSize(), on_dev).

Definition at line 486 of file densemat.hpp.

◆ ReadWrite()

real_t * mfem::DenseMatrix::ReadWrite ( bool on_dev = true)
inline

Shortcut for mfem::ReadWrite(GetMemory(), TotalSize(), on_dev).

Definition at line 498 of file densemat.hpp.

◆ Reset()

void mfem::DenseMatrix::Reset ( real_t * d,
int h,
int w )
inline

Change the data array and the size of the DenseMatrix.

The DenseMatrix does not assume ownership of the data array, i.e. it will not delete the new array d. This method will delete the current data array, if owned.

Definition at line 98 of file densemat.hpp.

◆ RightScaling()

void mfem::DenseMatrix::RightScaling ( const Vector & s)

RightScaling: this = this * diag(s);.

Definition at line 325 of file densemat.cpp.

◆ Set() [1/2]

void mfem::DenseMatrix::Set ( real_t alpha,
const DenseMatrix & A )
inline

Set the matrix to alpha * A.

Definition at line 242 of file densemat.hpp.

◆ Set() [2/2]

void mfem::DenseMatrix::Set ( real_t alpha,
const real_t * A )

Set the matrix to alpha * A, assuming that A has the same dimensions as the matrix and uses column-major layout.

Definition at line 580 of file densemat.cpp.

◆ SetCol() [1/3]

void mfem::DenseMatrix::SetCol ( int c,
const real_t * col )

Definition at line 2190 of file densemat.cpp.

◆ SetCol() [2/3]

void mfem::DenseMatrix::SetCol ( int c,
const Vector & col )

Definition at line 2199 of file densemat.cpp.

◆ SetCol() [3/3]

void mfem::DenseMatrix::SetCol ( int col,
real_t value )

Set all entries of a column to the specified value.

Definition at line 2167 of file densemat.cpp.

◆ SetRow() [1/3]

void mfem::DenseMatrix::SetRow ( int r,
const real_t * row )

Definition at line 2175 of file densemat.cpp.

◆ SetRow() [2/3]

void mfem::DenseMatrix::SetRow ( int r,
const Vector & row )

Definition at line 2184 of file densemat.cpp.

◆ SetRow() [3/3]

void mfem::DenseMatrix::SetRow ( int row,
real_t value )

Set all entries of a row to the specified value.

Definition at line 2159 of file densemat.cpp.

◆ SetSize() [1/2]

void mfem::DenseMatrix::SetSize ( int h,
int w )

Change the size of the DenseMatrix to h x w.

Definition at line 85 of file densemat.cpp.

◆ SetSize() [2/2]

void mfem::DenseMatrix::SetSize ( int s)
inline

Change the size of the DenseMatrix to s x s.

Definition at line 116 of file densemat.hpp.

◆ SetSubMatrix() [1/4]

void mfem::DenseMatrix::SetSubMatrix ( const Array< int > & idx,
const DenseMatrix & A )

Set (*this)(idx[i],idx[j]) = A(i,j)

Definition at line 1884 of file densemat.cpp.

◆ SetSubMatrix() [2/4]

void mfem::DenseMatrix::SetSubMatrix ( const Array< int > & idx_i,
const Array< int > & idx_j,
const DenseMatrix & A )

Set (*this)(idx_i[i],idx_j[j]) = A(i,j)

Definition at line 1913 of file densemat.cpp.

◆ SetSubMatrix() [3/4]

void mfem::DenseMatrix::SetSubMatrix ( int ibeg,
const DenseMatrix & A )

Set a submatrix of (this) to the given matrix A with row and column offset ibeg

Definition at line 1943 of file densemat.cpp.

◆ SetSubMatrix() [4/4]

void mfem::DenseMatrix::SetSubMatrix ( int ibeg,
int jbeg,
const DenseMatrix & A )

Set a submatrix of (this) to the given matrix A with row and column offset ibeg and jbeg respectively

Definition at line 1969 of file densemat.cpp.

◆ SingularValues()

void mfem::DenseMatrix::SingularValues ( Vector & sv) const

Definition at line 1210 of file densemat.cpp.

◆ Size()

int mfem::DenseMatrix::Size ( ) const
inline

For backward compatibility define Size to be synonym of Width()

Definition at line 110 of file densemat.hpp.

◆ SquareRootInverse()

void mfem::DenseMatrix::SquareRootInverse ( )

Replaces the current matrix with its square root inverse.

Definition at line 784 of file densemat.cpp.

◆ Swap()

void mfem::DenseMatrix::Swap ( DenseMatrix & other)

Definition at line 2333 of file densemat.cpp.

◆ SymmetricScaling()

void mfem::DenseMatrix::SymmetricScaling ( const Vector & s)

SymmetricScaling this = diag(sqrt(s)) * this * diag(sqrt(s))

Definition at line 354 of file densemat.cpp.

◆ Symmetrize()

void mfem::DenseMatrix::Symmetrize ( )

(*this) = 1/2 ((*this) + (*this)^t)

Definition at line 1450 of file densemat.cpp.

◆ TestInversion()

void mfem::DenseMatrix::TestInversion ( )

Invert and print the numerical conditioning of the inversion.

Definition at line 2319 of file densemat.cpp.

◆ Threshold()

void mfem::DenseMatrix::Threshold ( real_t eps)

Replace small entries, abs(a_ij) <= eps, with zero.

Definition at line 2205 of file densemat.cpp.

◆ TotalSize()

int mfem::DenseMatrix::TotalSize ( ) const
inline

Total size = width*height.

Definition at line 113 of file densemat.hpp.

◆ Trace()

real_t mfem::DenseMatrix::Trace ( ) const

Trace of a square matrix.

Definition at line 409 of file densemat.cpp.

◆ Transpose() [1/2]

void mfem::DenseMatrix::Transpose ( )

(*this) = (*this)^t

Definition at line 1417 of file densemat.cpp.

◆ Transpose() [2/2]

void mfem::DenseMatrix::Transpose ( const DenseMatrix & A)

(*this) = A^t

Definition at line 1439 of file densemat.cpp.

◆ UseExternalData()

void mfem::DenseMatrix::UseExternalData ( real_t * d,
int h,
int w )
inline

Change the data array and the size of the DenseMatrix.

The DenseMatrix does not assume ownership of the data array, i.e. it will not delete the data array d.

Definition at line 88 of file densemat.hpp.

◆ Weight()

real_t mfem::DenseMatrix::Weight ( ) const

Definition at line 553 of file densemat.cpp.

◆ Write()

real_t * mfem::DenseMatrix::Write ( bool on_dev = true)
inline

Shortcut for mfem::Write(GetMemory(), TotalSize(), on_dev).

Definition at line 492 of file densemat.hpp.

Friends And Related Symbol Documentation

◆ DenseMatrixInverse

friend class DenseMatrixInverse
friend

Definition at line 26 of file densemat.hpp.

◆ DenseTensor

friend class DenseTensor
friend

Definition at line 25 of file densemat.hpp.


The documentation for this class was generated from the following files: