libMesh
|
Defines a dense matrix for use in Finite Element-type computations. More...
#include <dof_map.h>
Classes | |
struct | UseBlasLapack |
Helper structure for determining whether to use blas_lapack. More... | |
Public Member Functions | |
DenseMatrix (const unsigned int new_m=0, const unsigned int new_n=0) | |
Constructor. More... | |
template<typename T2 > | |
DenseMatrix (unsigned int nrow, unsigned int ncol, std::initializer_list< T2 > init_list) | |
Constructor taking the number of rows, columns, and an initializer_list, which must be of length nrow * ncol, of row-major values to initialize the DenseMatrix with. More... | |
DenseMatrix (DenseMatrix &&)=default | |
The 5 special functions can be defaulted for this class, as it does not manage any memory itself. More... | |
DenseMatrix (const DenseMatrix &)=default | |
DenseMatrix & | operator= (const DenseMatrix &)=default |
DenseMatrix & | operator= (DenseMatrix &&)=default |
virtual | ~DenseMatrix ()=default |
virtual void | zero () override final |
Sets all elements of the matrix to 0 and resets any decomposition flag which may have been previously set. More... | |
DenseMatrix | sub_matrix (unsigned int row_id, unsigned int row_size, unsigned int col_id, unsigned int col_size) const |
Get submatrix with the smallest row and column indices and the submatrix size. More... | |
T | operator() (const unsigned int i, const unsigned int j) const |
T & | operator() (const unsigned int i, const unsigned int j) |
virtual T | el (const unsigned int i, const unsigned int j) const override final |
virtual T & | el (const unsigned int i, const unsigned int j) override final |
virtual void | left_multiply (const DenseMatrixBase< T > &M2) override final |
Performs the operation: (*this) <- M2 * (*this) More... | |
template<typename T2 > | |
void | left_multiply (const DenseMatrixBase< T2 > &M2) |
Left multiplies by the matrix M2 of different type. More... | |
virtual void | right_multiply (const DenseMatrixBase< T > &M2) override final |
Performs the operation: (*this) <- (*this) * M3. More... | |
template<typename T2 > | |
void | right_multiply (const DenseMatrixBase< T2 > &M2) |
Right multiplies by the matrix M2 of different type. More... | |
void | vector_mult (DenseVector< T > &dest, const DenseVector< T > &arg) const |
Performs the matrix-vector multiplication, dest := (*this) * arg . More... | |
template<typename T2 > | |
void | vector_mult (DenseVector< typename CompareTypes< T, T2 >::supertype > &dest, const DenseVector< T2 > &arg) const |
Performs the matrix-vector multiplication, dest := (*this) * arg on mixed types. More... | |
void | vector_mult_transpose (DenseVector< T > &dest, const DenseVector< T > &arg) const |
Performs the matrix-vector multiplication, dest := (*this)^T * arg . More... | |
template<typename T2 > | |
void | vector_mult_transpose (DenseVector< typename CompareTypes< T, T2 >::supertype > &dest, const DenseVector< T2 > &arg) const |
Performs the matrix-vector multiplication, dest := (*this)^T * arg . More... | |
void | vector_mult_add (DenseVector< T > &dest, const T factor, const DenseVector< T > &arg) const |
Performs the scaled matrix-vector multiplication, dest += factor * (*this) * arg . More... | |
template<typename T2 , typename T3 > | |
void | vector_mult_add (DenseVector< typename CompareTypes< T, typename CompareTypes< T2, T3 >::supertype >::supertype > &dest, const T2 factor, const DenseVector< T3 > &arg) const |
Performs the scaled matrix-vector multiplication, dest += factor * (*this) * arg . More... | |
void | get_principal_submatrix (unsigned int sub_m, unsigned int sub_n, DenseMatrix< T > &dest) const |
Put the sub_m x sub_n principal submatrix into dest . More... | |
void | get_principal_submatrix (unsigned int sub_m, DenseMatrix< T > &dest) const |
Put the sub_m x sub_m principal submatrix into dest . More... | |
void | outer_product (const DenseVector< T > &a, const DenseVector< T > &b) |
Computes the outer (dyadic) product of two vectors and stores in (*this). More... | |
template<typename T2 > | |
DenseMatrix< T > & | operator= (const DenseMatrix< T2 > &other_matrix) |
Assignment-from-other-matrix-type operator. More... | |
void | swap (DenseMatrix< T > &other_matrix) |
STL-like swap method. More... | |
void | resize (const unsigned int new_m, const unsigned int new_n) |
Resizes the matrix to the specified size and calls zero(). More... | |
void | scale (const T factor) |
Multiplies every element in the matrix by factor . More... | |
void | scale_column (const unsigned int col, const T factor) |
Multiplies every element in the column col matrix by factor . More... | |
DenseMatrix< T > & | operator*= (const T factor) |
Multiplies every element in the matrix by factor . More... | |
template<typename T2 , typename T3 > | |
boostcopy::enable_if_c< ScalarTraits< T2 >::value, void >::type | add (const T2 factor, const DenseMatrix< T3 > &mat) |
Adds factor times mat to this matrix. More... | |
bool | operator== (const DenseMatrix< T > &mat) const |
bool | operator!= (const DenseMatrix< T > &mat) const |
DenseMatrix< T > & | operator+= (const DenseMatrix< T > &mat) |
Adds mat to this matrix. More... | |
DenseMatrix< T > & | operator-= (const DenseMatrix< T > &mat) |
Subtracts mat from this matrix. More... | |
auto | min () const -> decltype(libmesh_real(T(0))) |
auto | max () const -> decltype(libmesh_real(T(0))) |
auto | l1_norm () const -> decltype(std::abs(T(0))) |
auto | linfty_norm () const -> decltype(std::abs(T(0))) |
void | left_multiply_transpose (const DenseMatrix< T > &A) |
Left multiplies by the transpose of the matrix A . More... | |
template<typename T2 > | |
void | left_multiply_transpose (const DenseMatrix< T2 > &A) |
Left multiplies by the transpose of the matrix A which contains a different numerical type. More... | |
void | right_multiply_transpose (const DenseMatrix< T > &A) |
Right multiplies by the transpose of the matrix A . More... | |
template<typename T2 > | |
void | right_multiply_transpose (const DenseMatrix< T2 > &A) |
Right multiplies by the transpose of the matrix A which contains a different numerical type. More... | |
T | transpose (const unsigned int i, const unsigned int j) const |
void | get_transpose (DenseMatrix< T > &dest) const |
Put the tranposed matrix into dest . More... | |
std::vector< T > & | get_values () |
const std::vector< T > & | get_values () const |
void | condense (const unsigned int i, const unsigned int j, const T val, DenseVector< T > &rhs) |
Condense-out the (i,j) entry of the matrix, forcing it to take on the value val . More... | |
void | lu_solve (const DenseVector< T > &b, DenseVector< T > &x) |
Solve the system Ax=b given the input vector b. More... | |
template<typename T2 > | |
void | cholesky_solve (const DenseVector< T2 > &b, DenseVector< T2 > &x) |
For symmetric positive definite (SPD) matrices. More... | |
void | svd (DenseVector< Real > &sigma) |
Compute the singular value decomposition of the matrix. More... | |
void | svd (DenseVector< Real > &sigma, DenseMatrix< Number > &U, DenseMatrix< Number > &VT) |
Compute the "reduced" singular value decomposition of the matrix. More... | |
void | svd_solve (const DenseVector< T > &rhs, DenseVector< T > &x, Real rcond=std::numeric_limits< Real >::epsilon()) const |
Solve the system of equations \( A x = rhs \) for \( x \) in the least-squares sense. More... | |
void | evd (DenseVector< T > &lambda_real, DenseVector< T > &lambda_imag) |
Compute the eigenvalues (both real and imaginary parts) of a general matrix. More... | |
void | evd_left (DenseVector< T > &lambda_real, DenseVector< T > &lambda_imag, DenseMatrix< T > &VL) |
Compute the eigenvalues (both real and imaginary parts) and left eigenvectors of a general matrix, \( A \). More... | |
void | evd_right (DenseVector< T > &lambda_real, DenseVector< T > &lambda_imag, DenseMatrix< T > &VR) |
Compute the eigenvalues (both real and imaginary parts) and right eigenvectors of a general matrix, \( A \). More... | |
void | evd_left_and_right (DenseVector< T > &lambda_real, DenseVector< T > &lambda_imag, DenseMatrix< T > &VL, DenseMatrix< T > &VR) |
Compute the eigenvalues (both real and imaginary parts) as well as the left and right eigenvectors of a general matrix. More... | |
T | det () |
unsigned int | m () const |
unsigned int | n () const |
void | print (std::ostream &os=libMesh::out) const |
Pretty-print the matrix, by default to libMesh::out . More... | |
void | print_scientific (std::ostream &os, unsigned precision=8) const |
Prints the matrix entries with more decimal places in scientific notation. More... | |
template<typename T2 , typename T3 > | |
boostcopy::enable_if_c< ScalarTraits< T2 >::value, void >::type | add (const T2 factor, const DenseMatrixBase< T3 > &mat) |
Adds factor to every element in the matrix. More... | |
DenseVector< T > | diagonal () const |
Return the matrix diagonal. More... | |
Public Attributes | |
bool | use_blas_lapack |
Computes the inverse of the dense matrix (assuming it is invertible) by first computing the LU decomposition and then performing multiple back substitution steps. More... | |
Protected Member Functions | |
void | condense (const unsigned int i, const unsigned int j, const T val, DenseVectorBase< T > &rhs) |
Condense-out the (i,j) entry of the matrix, forcing it to take on the value val . More... | |
Static Protected Member Functions | |
static void | multiply (DenseMatrixBase< T > &M1, const DenseMatrixBase< T > &M2, const DenseMatrixBase< T > &M3) |
Helper function - Performs the computation M1 = M2 * M3 where: M1 = (m x n) M2 = (m x p) M3 = (p x n) More... | |
Protected Attributes | |
unsigned int | _m |
The row dimension. More... | |
unsigned int | _n |
The column dimension. More... | |
Private Types | |
enum | DecompositionType { LU =0, CHOLESKY =1, LU_BLAS_LAPACK, NONE } |
The decomposition schemes above change the entries of the matrix A. More... | |
enum | _BLAS_Multiply_Flag { LEFT_MULTIPLY = 0, RIGHT_MULTIPLY, LEFT_MULTIPLY_TRANSPOSE, RIGHT_MULTIPLY_TRANSPOSE } |
Enumeration used to determine the behavior of the _multiply_blas function. More... | |
typedef PetscBLASInt | pivot_index_t |
Array used to store pivot indices. More... | |
typedef int | pivot_index_t |
Private Member Functions | |
void | _lu_decompose () |
Form the LU decomposition of the matrix. More... | |
void | _lu_back_substitute (const DenseVector< T > &b, DenseVector< T > &x) const |
Solves the system Ax=b through back substitution. More... | |
void | _cholesky_decompose () |
Decomposes a symmetric positive definite matrix into a product of two lower triangular matrices according to A = LL^T. More... | |
template<typename T2 > | |
void | _cholesky_back_substitute (const DenseVector< T2 > &b, DenseVector< T2 > &x) const |
Solves the equation Ax=b for the unknown value x and rhs b based on the Cholesky factorization of A. More... | |
void | _multiply_blas (const DenseMatrixBase< T > &other, _BLAS_Multiply_Flag flag) |
The _multiply_blas function computes A <- op(A) * op(B) using BLAS gemm function. More... | |
void | _lu_decompose_lapack () |
Computes an LU factorization of the matrix using the Lapack routine "getrf". More... | |
void | _svd_lapack (DenseVector< Real > &sigma) |
Computes an SVD of the matrix using the Lapack routine "getsvd". More... | |
void | _svd_lapack (DenseVector< Real > &sigma, DenseMatrix< Number > &U, DenseMatrix< Number > &VT) |
Computes a "reduced" SVD of the matrix using the Lapack routine "getsvd". More... | |
void | _svd_solve_lapack (const DenseVector< T > &rhs, DenseVector< T > &x, Real rcond) const |
Called by svd_solve(rhs). More... | |
void | _svd_helper (char JOBU, char JOBVT, std::vector< Real > &sigma_val, std::vector< Number > &U_val, std::vector< Number > &VT_val) |
Helper function that actually performs the SVD. More... | |
void | _evd_lapack (DenseVector< T > &lambda_real, DenseVector< T > &lambda_imag, DenseMatrix< T > *VL=nullptr, DenseMatrix< T > *VR=nullptr) |
Computes the eigenvalues of the matrix using the Lapack routine "DGEEV". More... | |
void | _lu_back_substitute_lapack (const DenseVector< T > &b, DenseVector< T > &x) |
Companion function to _lu_decompose_lapack(). More... | |
void | _matvec_blas (T alpha, T beta, DenseVector< T > &dest, const DenseVector< T > &arg, bool trans=false) const |
Uses the BLAS GEMV function (through PETSc) to compute. More... | |
template<typename T2 > | |
void | _left_multiply_transpose (const DenseMatrix< T2 > &A) |
Left multiplies by the transpose of the matrix A which may contain a different numerical type. More... | |
template<typename T2 > | |
void | _right_multiply_transpose (const DenseMatrix< T2 > &A) |
Right multiplies by the transpose of the matrix A which may contain a different numerical type. More... | |
Private Attributes | |
std::vector< T > | _val |
The actual data values, stored as a 1D array. More... | |
DecompositionType | _decomposition_type |
This flag keeps track of which type of decomposition has been performed on the matrix. More... | |
std::vector< pivot_index_t > | _pivots |
Defines a dense matrix for use in Finite Element-type computations.
Useful for storing element stiffness matrices before summation into a global matrix. All overridden virtual functions are documented in dense_matrix_base.h.
|
private |
Array used to store pivot indices.
May be used by whatever factorization is currently active, clients of the class should not rely on it for any reason.
Definition at line 736 of file dense_matrix.h.
|
private |
Definition at line 738 of file dense_matrix.h.
|
private |
Enumeration used to determine the behavior of the _multiply_blas function.
Enumerator | |
---|---|
LEFT_MULTIPLY | |
RIGHT_MULTIPLY | |
LEFT_MULTIPLY_TRANSPOSE | |
RIGHT_MULTIPLY_TRANSPOSE |
Definition at line 656 of file dense_matrix.h.
|
private |
The decomposition schemes above change the entries of the matrix A.
It is therefore an error to call A.lu_solve() and subsequently call A.cholesky_solve() since the result will probably not match any desired outcome. This typedef keeps track of which decomposition has been called for this matrix.
Enumerator | |
---|---|
LU | |
CHOLESKY | |
LU_BLAS_LAPACK | |
NONE |
Definition at line 644 of file dense_matrix.h.
|
inline |
Constructor.
Creates a dense matrix of dimension m
by n
.
Definition at line 836 of file dense_matrix.h.
libMesh::DenseMatrix< T >::DenseMatrix | ( | unsigned int | nrow, |
unsigned int | ncol, | ||
std::initializer_list< T2 > | init_list | ||
) |
Constructor taking the number of rows, columns, and an initializer_list, which must be of length nrow * ncol, of row-major values to initialize the DenseMatrix with.
Definition at line 848 of file dense_matrix.h.
|
default |
The 5 special functions can be defaulted for this class, as it does not manage any memory itself.
|
default |
|
virtualdefault |
|
private |
Solves the equation Ax=b for the unknown value x and rhs b based on the Cholesky factorization of A.
Definition at line 919 of file dense_matrix_impl.h.
|
private |
Decomposes a symmetric positive definite matrix into a product of two lower triangular matrices according to A = LL^T.
Definition at line 873 of file dense_matrix_impl.h.
|
private |
Computes the eigenvalues of the matrix using the Lapack routine "DGEEV".
If VR and/or VL are not nullptr, then the matrix of right and/or left eigenvectors is also computed and returned by this function.
[ Implementation in dense_matrix_blas_lapack.C ]
Definition at line 695 of file dense_matrix_blas_lapack.C.
|
private |
Left multiplies by the transpose of the matrix A
which may contain a different numerical type.
Definition at line 112 of file dense_matrix_impl.h.
|
private |
Solves the system Ax=b through back substitution.
This function is private since it is only called as part of the implementation of the lu_solve(...) function.
Definition at line 580 of file dense_matrix_impl.h.
|
private |
Companion function to _lu_decompose_lapack().
Do not use directly, called through the public lu_solve() interface. This function is logically const in that it does not modify the matrix, but since we are just calling LAPACK routines, it's less const_cast hassle to just declare the function non-const. [ Implementation in dense_matrix_blas_lapack.C ]
Definition at line 901 of file dense_matrix_blas_lapack.C.
|
private |
Form the LU decomposition of the matrix.
This function is private since it is only called as part of the implementation of the lu_solve(...) function.
Definition at line 631 of file dense_matrix_impl.h.
|
private |
Computes an LU factorization of the matrix using the Lapack routine "getrf".
This routine should only be used by the "use_blas_lapack" branch of the lu_solve() function. After the call to this function, the matrix is replaced by its factorized version, and the DecompositionType is set to LU_BLAS_LAPACK. [ Implementation in dense_matrix_blas_lapack.C ]
Definition at line 210 of file dense_matrix_blas_lapack.C.
|
private |
Uses the BLAS GEMV function (through PETSc) to compute.
dest := alpha*A*arg + beta*dest
where alpha and beta are scalars, A is this matrix, and arg and dest are input vectors of appropriate size. If trans is true, the transpose matvec is computed instead. By default, trans==false.
[ Implementation in dense_matrix_blas_lapack.C ]
Definition at line 990 of file dense_matrix_blas_lapack.C.
|
private |
The _multiply_blas function computes A <- op(A) * op(B) using BLAS gemm function.
Used in the right_multiply(), left_multiply(), right_multiply_transpose(), and left_multiply_transpose() routines. [ Implementation in dense_matrix_blas_lapack.C ]
Definition at line 37 of file dense_matrix_blas_lapack.C.
|
private |
Right multiplies by the transpose of the matrix A
which may contain a different numerical type.
Definition at line 239 of file dense_matrix_impl.h.
|
private |
Helper function that actually performs the SVD.
[ Implementation in dense_matrix_blas_lapack.C ]
Definition at line 384 of file dense_matrix_blas_lapack.C.
|
private |
Computes an SVD of the matrix using the Lapack routine "getsvd".
[ Implementation in dense_matrix_blas_lapack.C ]
Definition at line 277 of file dense_matrix_blas_lapack.C.
|
private |
Computes a "reduced" SVD of the matrix using the Lapack routine "getsvd".
[ Implementation in dense_matrix_blas_lapack.C ]
Definition at line 318 of file dense_matrix_blas_lapack.C.
|
private |
Called by svd_solve(rhs).
Definition at line 529 of file dense_matrix_blas_lapack.C.
|
inlineinherited |
Adds factor
to every element in the matrix.
This should only work if T += T2 * T3 is valid C++ and if T2 is scalar. Return type is void
Definition at line 195 of file dense_matrix_base.h.
References libMesh::DenseMatrixBase< T >::el(), libMesh::DenseMatrixBase< T >::m(), libMesh::make_range(), and libMesh::DenseMatrixBase< T >::n().
|
inline |
Adds factor
times mat
to this matrix.
Definition at line 1018 of file dense_matrix.h.
Referenced by assemble_wave(), libMesh::FEMSystem::assembly(), libMesh::TransientRBEvaluation::rb_solve(), and libMesh::RBEvaluation::rb_solve().
template LIBMESH_EXPORT void libMesh::DenseMatrix< T >::cholesky_solve | ( | const DenseVector< T2 > & | b, |
DenseVector< T2 > & | x | ||
) |
For symmetric positive definite (SPD) matrices.
A Cholesky factorization of A such that A = L L^T is about twice as fast as a standard LU factorization. Therefore you can use this method if you know a-priori that the matrix is SPD. If the matrix is not SPD, an error is generated. One nice property of Cholesky decompositions is that they do not require pivoting for stability.
Important note: once you call cholesky_solve(), you must not modify the entries of the matrix via calls to operator(i,j) and call cholesky_solve() again without first calling either zero() or resize(), otherwise the code will skip computing the decomposition of the matrix and go directly to the back substitution step. This is done on purpose for efficiency, so that the same decomposition can be used with multiple right-hand sides, but it does also make it possible to "shoot yourself in the foot", so be careful!
Definition at line 841 of file dense_matrix_impl.h.
Referenced by libMesh::FEGenericBase< FEOutputType< T >::type >::coarsened_dof_values(), libMesh::FEGenericBase< FEOutputType< T >::type >::compute_periodic_constraints(), libMesh::FEGenericBase< FEOutputType< T >::type >::compute_proj_constraints(), LinearElasticity::compute_stresses(), libMesh::GenericProjector< FFunctor, GFunctor, FValue, ProjectionAction >::SubProjector::construct_projection(), libMesh::BoundaryProjectSolution::operator()(), and libMesh::HPCoarsenTest::select_refinement().
|
protectedinherited |
Condense-out the (i,j) entry of the matrix, forcing it to take on the value
val
.
This is useful in numerical simulations for applying boundary conditions. Preserves the symmetry of the matrix.
Definition at line 73 of file dense_matrix_base_impl.h.
References libMesh::DenseVectorBase< T >::el(), libMesh::make_range(), and libMesh::DenseVectorBase< T >::size().
Referenced by libMesh::DenseMatrix< Real >::condense().
|
inline |
Condense-out the (i,j) entry of the matrix, forcing it to take on the value
val
.
This is useful in numerical simulations for applying boundary conditions. Preserves the symmetry of the matrix.
Definition at line 395 of file dense_matrix.h.
T libMesh::DenseMatrix< T >::det | ( | ) |
Definition at line 780 of file dense_matrix_impl.h.
|
inherited |
Return the matrix diagonal.
Definition at line 37 of file dense_matrix_base_impl.h.
Referenced by libMesh::DiagonalMatrix< T >::add_matrix().
|
inlinefinaloverridevirtual |
(i,j) element of the matrix. Since internal data representations may differ, you must redefine this function. Implements libMesh::DenseMatrixBase< T >.
Definition at line 126 of file dense_matrix.h.
Referenced by libMesh::RBConstruction::add_scaled_matrix_and_vector(), libMesh::TransientRBConstruction::enrich_RB_space(), libMesh::RBEIMConstruction::train_eim_approximation_with_POD(), and libMesh::RBConstruction::train_reduced_basis_with_POD().
|
inlinefinaloverridevirtual |
(i,j) element of the matrix as a writable reference. Since internal data representations may differ, you must redefine this function. Implements libMesh::DenseMatrixBase< T >.
Definition at line 130 of file dense_matrix.h.
void libMesh::DenseMatrix< T >::evd | ( | DenseVector< T > & | lambda_real, |
DenseVector< T > & | lambda_imag | ||
) |
Compute the eigenvalues (both real and imaginary parts) of a general matrix.
Warning: the contents of *this
are overwritten by this function!
The implementation requires the LAPACKgeev_ function which is wrapped by PETSc.
Definition at line 736 of file dense_matrix_impl.h.
void libMesh::DenseMatrix< T >::evd_left | ( | DenseVector< T > & | lambda_real, |
DenseVector< T > & | lambda_imag, | ||
DenseMatrix< T > & | VL | ||
) |
Compute the eigenvalues (both real and imaginary parts) and left eigenvectors of a general matrix, \( A \).
Warning: the contents of *this
are overwritten by this function!
The left eigenvector \( u_j \) of \( A \) satisfies: \( u_j^H A = lambda_j u_j^H \) where \( u_j^H \) denotes the conjugate-transpose of \( u_j \).
If the j-th and (j+1)-st eigenvalues form a complex conjugate pair, then the j-th and (j+1)-st columns of VL "share" their real-valued storage in the following way: u_j = VL(:,j) + i*VL(:,j+1) and u_{j+1} = VL(:,j) - i*VL(:,j+1).
The implementation requires the LAPACKgeev_ routine which is provided by PETSc.
Definition at line 746 of file dense_matrix_impl.h.
void libMesh::DenseMatrix< T >::evd_left_and_right | ( | DenseVector< T > & | lambda_real, |
DenseVector< T > & | lambda_imag, | ||
DenseMatrix< T > & | VL, | ||
DenseMatrix< T > & | VR | ||
) |
Compute the eigenvalues (both real and imaginary parts) as well as the left and right eigenvectors of a general matrix.
Warning: the contents of *this
are overwritten by this function!
See the documentation of the evd_left()
and evd_right()
functions for more information. The implementation requires the LAPACKgeev_ routine which is provided by PETSc.
Definition at line 768 of file dense_matrix_impl.h.
Referenced by DenseMatrixTest::testEVD_helper().
void libMesh::DenseMatrix< T >::evd_right | ( | DenseVector< T > & | lambda_real, |
DenseVector< T > & | lambda_imag, | ||
DenseMatrix< T > & | VR | ||
) |
Compute the eigenvalues (both real and imaginary parts) and right eigenvectors of a general matrix, \( A \).
Warning: the contents of *this
are overwritten by this function!
The right eigenvector \( v_j \) of \( A \) satisfies: \( A v_j = lambda_j v_j \) where \( lambda_j \) is its corresponding eigenvalue.
The implementation requires the LAPACKgeev_ routine which is provided by PETSc.
Definition at line 757 of file dense_matrix_impl.h.
void libMesh::DenseMatrix< T >::get_principal_submatrix | ( | unsigned int | sub_m, |
unsigned int | sub_n, | ||
DenseMatrix< T > & | dest | ||
) | const |
Put the sub_m
x sub_n
principal submatrix into dest
.
Definition at line 470 of file dense_matrix_impl.h.
Referenced by libMesh::RBEIMConstruction::compute_max_eim_error(), libMesh::RBEIMEvaluation::rb_eim_solve(), libMesh::RBEIMEvaluation::rb_eim_solves(), and libMesh::TransientRBConstruction::update_RB_initial_condition_all_N().
void libMesh::DenseMatrix< T >::get_principal_submatrix | ( | unsigned int | sub_m, |
DenseMatrix< T > & | dest | ||
) | const |
Put the sub_m
x sub_m
principal submatrix into dest
.
Definition at line 500 of file dense_matrix_impl.h.
void libMesh::DenseMatrix< T >::get_transpose | ( | DenseMatrix< T > & | dest | ) | const |
Put the tranposed matrix into dest
.
Definition at line 508 of file dense_matrix_impl.h.
Referenced by libMesh::RBConstruction::add_scaled_matrix_and_vector().
|
inline |
This should be used with caution (i.e. one should not change the size of the vector, etc.) but is useful for interoperating with low level BLAS routines which expect a simple array.
Definition at line 382 of file dense_matrix.h.
Referenced by libMesh::DenseMatrix< Real >::_evd_lapack(), libMesh::DenseMatrix< Real >::_matvec_blas(), libMesh::DenseMatrix< Real >::_multiply_blas(), libMesh::DenseMatrix< Real >::_svd_lapack(), libMesh::DenseMatrix< Real >::_svd_solve_lapack(), libMesh::PetscMatrix< T >::add_block_matrix(), libMesh::EpetraMatrix< T >::add_matrix(), libMesh::PetscMatrix< T >::add_matrix(), and OverlappingCouplingGhostingTest::run_sparsity_pattern_test().
|
inline |
Definition at line 387 of file dense_matrix.h.
|
inline |
\( |M|_1 = max_{all columns j} \sum_{all rows i} |M_ij| \),
This is the natural matrix norm that is compatible to the l1-norm for vectors, i.e. \( |Mv|_1 \leq |M|_1 |v|_1 \).
Definition at line 1125 of file dense_matrix.h.
Referenced by libMesh::FEMSystem::assembly().
|
finaloverridevirtual |
Performs the operation: (*this) <- M2 * (*this)
Implements libMesh::DenseMatrixBase< T >.
Definition at line 45 of file dense_matrix_impl.h.
void libMesh::DenseMatrix< T >::left_multiply | ( | const DenseMatrixBase< T2 > & | M2 | ) |
Left multiplies by the matrix M2
of different type.
Definition at line 72 of file dense_matrix_impl.h.
void libMesh::DenseMatrix< T >::left_multiply_transpose | ( | const DenseMatrix< T > & | A | ) |
Left multiplies by the transpose of the matrix A
.
Definition at line 93 of file dense_matrix_impl.h.
Referenced by libMesh::DofMap::heterogeneously_constrain_element_jacobian_and_residual().
void libMesh::DenseMatrix< T >::left_multiply_transpose | ( | const DenseMatrix< T2 > & | A | ) |
Left multiplies by the transpose of the matrix A
which contains a different numerical type.
Definition at line 104 of file dense_matrix_impl.h.
|
inline |
\( |M|_\infty = max_{all rows i} \sum_{all columns j} |M_ij| \),
This is the natural matrix norm that is compatible to the linfty-norm of vectors, i.e. \( |Mv|_\infty \leq |M|_\infty |v|_\infty \).
Definition at line 1152 of file dense_matrix.h.
void libMesh::DenseMatrix< T >::lu_solve | ( | const DenseVector< T > & | b, |
DenseVector< T > & | x | ||
) |
Solve the system Ax=b given the input vector b.
Partial pivoting is performed by default in order to keep the algorithm stable to the effects of round-off error.
Important note: once you call lu_solve(), you must not modify the entries of the matrix via calls to operator(i,j) and call lu_solve() again without first calling either zero() or resize(), otherwise the code will skip computing the decomposition of the matrix and go directly to the back substitution step. This is done on purpose for efficiency, so that the same LU decomposition can be used with multiple right-hand sides, but it does also make it possible to "shoot yourself in the foot", so be careful!
Definition at line 521 of file dense_matrix_impl.h.
Referenced by compute_enriched_soln(), libMesh::RBEIMConstruction::compute_max_eim_error(), fe_assembly(), libMesh::WeightedPatchRecoveryErrorEstimator::EstimateError::operator()(), libMesh::PatchRecoveryErrorEstimator::EstimateError::operator()(), libMesh::RBEIMEvaluation::rb_eim_solve(), libMesh::RBEIMEvaluation::rb_eim_solves(), libMesh::TransientRBEvaluation::rb_solve(), libMesh::RBEvaluation::rb_solve(), libMesh::TransientRBEvaluation::rb_solve_again(), libMesh::TransientRBConstruction::set_error_temporal_data(), and libMesh::TransientRBConstruction::update_RB_initial_condition_all_N().
|
inlineinherited |
Definition at line 104 of file dense_matrix_base.h.
References libMesh::DenseMatrixBase< T >::_m.
Referenced by libMesh::DenseMatrix< Real >::_left_multiply_transpose(), libMesh::DenseMatrix< Real >::_multiply_blas(), libMesh::DenseMatrix< Real >::_right_multiply_transpose(), libMesh::DenseMatrix< Real >::_svd_solve_lapack(), libMesh::DenseMatrixBase< T >::add(), libMesh::PetscMatrix< T >::add_block_matrix(), libMesh::SparseMatrix< ValOut >::add_block_matrix(), libMesh::StaticCondensation::add_matrix(), libMesh::EigenSparseMatrix< T >::add_matrix(), libMesh::LaspackMatrix< T >::add_matrix(), libMesh::DiagonalMatrix< T >::add_matrix(), libMesh::EpetraMatrix< T >::add_matrix(), libMesh::PetscMatrix< T >::add_matrix(), libMesh::RBConstruction::add_scaled_matrix_and_vector(), libMesh::DofMap::build_constraint_matrix(), libMesh::DofMap::build_constraint_matrix_and_vector(), libMesh::DofMap::extract_local_vector(), libMesh::DenseMatrix< Real >::get_transpose(), libMesh::DofMap::heterogeneously_constrain_element_jacobian_and_residual(), libMesh::DofMap::heterogeneously_constrain_element_residual(), libMesh::DenseMatrix< Real >::left_multiply(), libMesh::DofMap::max_constraint_error(), libMesh::DenseMatrixBase< T >::multiply(), libMesh::WeightedPatchRecoveryErrorEstimator::EstimateError::operator()(), libMesh::PatchRecoveryErrorEstimator::EstimateError::operator()(), libMesh::DenseMatrix< Real >::right_multiply(), libMesh::RBEIMEvaluation::set_interpolation_matrix_entry(), DualShapeTest::testEdge2Lagrange(), DenseMatrixTest::testEVD_helper(), DenseMatrixTest::testSVD(), and MetaPhysicL::RawType< libMesh::DenseMatrix< T > >::value().
|
inline |
Definition at line 1104 of file dense_matrix.h.
|
inline |
Definition at line 1083 of file dense_matrix.h.
|
staticprotectedinherited |
Helper function - Performs the computation M1 = M2 * M3 where: M1 = (m x n) M2 = (m x p) M3 = (p x n)
Definition at line 46 of file dense_matrix_base_impl.h.
References libMesh::DenseMatrixBase< T >::el(), libMesh::DenseMatrixBase< T >::m(), and libMesh::DenseMatrixBase< T >::n().
|
inlineinherited |
Definition at line 109 of file dense_matrix_base.h.
References libMesh::DenseMatrixBase< T >::_n.
Referenced by libMesh::DenseMatrix< Real >::_left_multiply_transpose(), libMesh::DenseMatrix< Real >::_multiply_blas(), libMesh::DenseMatrix< Real >::_right_multiply_transpose(), libMesh::DenseMatrix< Real >::_svd_solve_lapack(), libMesh::DenseMatrixBase< T >::add(), libMesh::PetscMatrix< T >::add_block_matrix(), libMesh::SparseMatrix< ValOut >::add_block_matrix(), libMesh::StaticCondensation::add_matrix(), libMesh::EigenSparseMatrix< T >::add_matrix(), libMesh::LaspackMatrix< T >::add_matrix(), libMesh::DiagonalMatrix< T >::add_matrix(), libMesh::EpetraMatrix< T >::add_matrix(), libMesh::PetscMatrix< T >::add_matrix(), libMesh::RBConstruction::add_scaled_matrix_and_vector(), libMesh::DofMap::build_constraint_matrix(), libMesh::DofMap::build_constraint_matrix_and_vector(), libMesh::DofMap::constrain_element_residual(), libMesh::DofMap::extract_local_vector(), libMesh::DenseMatrix< Real >::get_transpose(), libMesh::DofMap::heterogeneously_constrain_element_jacobian_and_residual(), libMesh::DofMap::heterogeneously_constrain_element_residual(), libMesh::DenseMatrix< Real >::left_multiply(), libMesh::DofMap::max_constraint_error(), libMesh::DenseMatrixBase< T >::multiply(), libMesh::WeightedPatchRecoveryErrorEstimator::EstimateError::operator()(), libMesh::PatchRecoveryErrorEstimator::EstimateError::operator()(), libMesh::DenseMatrix< Real >::right_multiply(), libMesh::RBEIMEvaluation::set_interpolation_matrix_entry(), DualShapeTest::testEdge2Lagrange(), DenseMatrixTest::testSVD(), and MetaPhysicL::RawType< libMesh::DenseMatrix< T > >::value().
|
inline |
true
if mat
is not exactly equal to this matrix, false otherwise. Definition at line 1046 of file dense_matrix.h.
|
inline |
(i,j) element of the matrix. Definition at line 953 of file dense_matrix.h.
|
inline |
(i,j) element of the matrix as a writable reference. Definition at line 969 of file dense_matrix.h.
|
inline |
Multiplies every element in the matrix by factor
.
Definition at line 1005 of file dense_matrix.h.
|
inline |
Adds mat
to this matrix.
Definition at line 1059 of file dense_matrix.h.
|
inline |
Subtracts mat
from this matrix.
Definition at line 1071 of file dense_matrix.h.
|
default |
|
default |
|
inline |
Assignment-from-other-matrix-type operator.
Copies the dense matrix of type T2 into the present matrix. This is useful for copying real matrices into complex ones for further operations.
Definition at line 880 of file dense_matrix.h.
|
inline |
true
if mat
is exactly equal to this matrix, false
otherwise. Definition at line 1033 of file dense_matrix.h.
void libMesh::DenseMatrix< T >::outer_product | ( | const DenseVector< T > & | a, |
const DenseVector< T > & | b | ||
) |
Computes the outer (dyadic) product of two vectors and stores in (*this).
The outer product of two real-valued vectors \(\mathbf{a}\) and \(\mathbf{b}\) is
\[ (\mathbf{a}\mathbf{b}^T)_{i,j} = \mathbf{a}_i \mathbf{b}_j . \]
The outer product of two complex-valued vectors \(\mathbf{a}\) and \(\mathbf{b}\) is
\[ (\mathbf{a}\mathbf{b}^H)_{i,j} = \mathbf{a}_i \mathbf{b}^*_j , \]
where \(H\) denotes the conjugate transpose of the vector and \(*\) denotes the complex conjugate.
[in] | a | Vector whose entries correspond to rows in the product matrix. |
[in] | b | Vector whose entries correspond to columns in the product matrix. |
Definition at line 485 of file dense_matrix_impl.h.
Referenced by DenseMatrixTest::testOuterProduct().
|
inherited |
Pretty-print the matrix, by default to libMesh::out
.
Definition at line 125 of file dense_matrix_base_impl.h.
References libMesh::make_range().
|
inherited |
Prints the matrix entries with more decimal places in scientific notation.
Definition at line 101 of file dense_matrix_base_impl.h.
References libMesh::make_range().
|
inline |
Resizes the matrix to the specified size and calls zero().
Will never free memory, but may allocate more. Note: when the matrix is zero()'d, any decomposition (LU, Cholesky, etc.) is also cleared, forcing a new decomposition to be computed the next time e.g. lu_solve() is called.
Definition at line 895 of file dense_matrix.h.
Referenced by libMesh::DenseMatrix< Real >::_evd_lapack(), libMesh::DenseMatrix< Real >::_lu_decompose(), libMesh::DenseMatrix< Real >::_svd_lapack(), libMesh::HPCoarsenTest::add_projection(), alternative_fe_assembly(), assemble(), LinearElasticity::assemble(), assemble_1D(), AssembleOptimization::assemble_A_and_F(), assemble_biharmonic(), assemble_cd(), assemble_divgrad(), assemble_elasticity(), assemble_ellipticdg(), assemble_func(), assemble_graddiv(), assemble_helmholtz(), libMesh::ClawSystem::assemble_jump_coupling_matrix(), assemble_laplace(), assemble_mass(), libMesh::ClawSystem::assemble_mass_matrix(), assemble_matrices(), assemble_matrix_and_rhs(), assemble_poisson(), assemble_SchroedingerEquation(), assemble_shell(), assemble_stokes(), assemble_wave(), libMesh::DofMap::build_constraint_matrix(), libMesh::DofMap::build_constraint_matrix_and_vector(), libMesh::TransientRBEvaluation::cache_online_residual_terms(), NonlinearNeoHookeCurrentConfig::calculate_tangent(), libMesh::RBEIMConstruction::clear(), libMesh::RBEIMEvaluation::clear(), libMesh::StaticCondensation::close(), libMesh::FEGenericBase< FEOutputType< T >::type >::coarsened_dof_values(), compute_enriched_soln(), compute_jacobian(), libMesh::FEGenericBase< FEOutputType< T >::type >::compute_periodic_constraints(), libMesh::FEGenericBase< FEOutputType< T >::type >::compute_proj_constraints(), compute_stresses(), fe_assembly(), form_matrixA(), NonlinearNeoHookeCurrentConfig::get_linearized_stiffness(), libMesh::DenseMatrix< Real >::get_principal_submatrix(), NonlinearNeoHookeCurrentConfig::get_residual(), libMesh::DenseMatrix< Real >::get_transpose(), libMesh::FESubdivision::init_subdivision_matrix(), libMesh::HDGProblem::jacobian(), LaplaceYoung::jacobian(), LargeDeformationElasticity::jacobian(), libMesh::DGFEMContext::neighbor_side_fe_reinit(), libMesh::BoundaryProjectSolution::operator()(), periodic_bc_test_poisson(), libMesh::FEMContext::pre_fe_reinit(), libMesh::TransientRBEvaluation::rb_solve(), libMesh::RBEIMConstruction::reinit_eim_projection_matrix(), Biharmonic::JR::residual_and_jacobian(), LinearElasticityWithContact::residual_and_jacobian(), libMesh::DenseMatrix< Real >::resize(), libMesh::TransientRBEvaluation::resize_data_structures(), libMesh::RBEvaluation::resize_data_structures(), libMesh::RBEIMEvaluation::resize_data_structures(), OverlappingCouplingGhostingTest::run_sparsity_pattern_test(), libMesh::HPCoarsenTest::select_refinement(), libMesh::StaticCondensation::StaticCondensation(), libMesh::RBEIMConstruction::train_eim_approximation_with_greedy(), and libMesh::RBEIMConstruction::train_eim_approximation_with_POD().
|
finaloverridevirtual |
Performs the operation: (*this) <- (*this) * M3.
Implements libMesh::DenseMatrixBase< T >.
Definition at line 171 of file dense_matrix_impl.h.
Referenced by libMesh::DofMap::build_constraint_matrix(), libMesh::DofMap::build_constraint_matrix_and_vector(), NonlinearNeoHookeCurrentConfig::get_linearized_stiffness(), libMesh::DofMap::heterogeneously_constrain_element_jacobian_and_residual(), and libMesh::FESubdivision::init_shape_functions().
void libMesh::DenseMatrix< T >::right_multiply | ( | const DenseMatrixBase< T2 > & | M2 | ) |
Right multiplies by the matrix M2
of different type.
Definition at line 197 of file dense_matrix_impl.h.
void libMesh::DenseMatrix< T >::right_multiply_transpose | ( | const DenseMatrix< T > & | A | ) |
Right multiplies by the transpose of the matrix A
.
Definition at line 218 of file dense_matrix_impl.h.
Referenced by NonlinearNeoHookeCurrentConfig::get_linearized_stiffness().
void libMesh::DenseMatrix< T >::right_multiply_transpose | ( | const DenseMatrix< T2 > & | A | ) |
Right multiplies by the transpose of the matrix A
which contains a different numerical type.
Definition at line 230 of file dense_matrix_impl.h.
|
inline |
Multiplies every element in the matrix by factor
.
Definition at line 986 of file dense_matrix.h.
Referenced by compute_stresses(), and SolidSystem::element_time_derivative().
|
inline |
Multiplies every element in the column col
matrix by factor
.
Definition at line 995 of file dense_matrix.h.
|
inline |
Get submatrix with the smallest row and column indices and the submatrix size.
Definition at line 922 of file dense_matrix.h.
Referenced by DenseMatrixTest::testSubMatrix().
void libMesh::DenseMatrix< T >::svd | ( | DenseVector< Real > & | sigma | ) |
Compute the singular value decomposition of the matrix.
On exit, sigma holds all of the singular values (in descending order).
The implementation uses PETSc's interface to BLAS/LAPACK. If this is not available, this function throws an error.
Definition at line 707 of file dense_matrix_impl.h.
Referenced by libMesh::TransientRBConstruction::enrich_RB_space(), DenseMatrixTest::testComplexSVD(), DenseMatrixTest::testSVD(), libMesh::RBEIMConstruction::train_eim_approximation_with_POD(), and libMesh::RBConstruction::train_reduced_basis_with_POD().
void libMesh::DenseMatrix< T >::svd | ( | DenseVector< Real > & | sigma, |
DenseMatrix< Number > & | U, | ||
DenseMatrix< Number > & | VT | ||
) |
Compute the "reduced" singular value decomposition of the matrix.
On exit, sigma holds all of the singular values (in descending order), U holds the left singular vectors, and VT holds the transpose of the right singular vectors. In the reduced SVD, U has min(m,n) columns and VT has min(m,n) rows. (In the "full" SVD, U and VT would be square.)
The implementation uses PETSc's interface to BLAS/LAPACK. If this is not available, this function throws an error.
Definition at line 715 of file dense_matrix_impl.h.
void libMesh::DenseMatrix< T >::svd_solve | ( | const DenseVector< T > & | rhs, |
DenseVector< T > & | x, | ||
Real | rcond = std::numeric_limits<Real>::epsilon() |
||
) | const |
Solve the system of equations \( A x = rhs \) for \( x \) in the least-squares sense.
\( A \) may be non-square and/or rank-deficient. You can control which singular values are treated as zero by changing the "rcond" parameter. Singular values S(i) for which S(i) <= rcond*S(1) are treated as zero for purposes of the solve. Passing a negative number for rcond forces a "machine precision" value to be used instead.
This function is marked const, since due to various implementation details, we do not need to modify the contents of A in order to compute the SVD (a copy is made internally instead).
Requires PETSc >= 3.1 since this was the first version to provide the LAPACKgelss_ wrapper.
Definition at line 726 of file dense_matrix_impl.h.
Referenced by fe_assembly().
|
inline |
STL-like swap method.
Definition at line 865 of file dense_matrix.h.
Referenced by libMesh::EulerSolver::_general_residual(), libMesh::Euler2Solver::_general_residual(), and libMesh::NewmarkSolver::_general_residual().
|
inline |
(i,j) element of the transposed matrix. Definition at line 1179 of file dense_matrix.h.
Referenced by libMesh::DenseMatrix< Real >::_left_multiply_transpose().
void libMesh::DenseMatrix< T >::vector_mult | ( | DenseVector< T > & | dest, |
const DenseVector< T > & | arg | ||
) | const |
Performs the matrix-vector multiplication, dest
:= (*this) * arg
.
Definition at line 299 of file dense_matrix_impl.h.
Referenced by libMesh::FEMap::compute_single_point_map(), NonlinearNeoHookeCurrentConfig::get_residual(), libMesh::TransientRBEvaluation::rb_solve(), and libMesh::TransientRBEvaluation::rb_solve_again().
void libMesh::DenseMatrix< T >::vector_mult | ( | DenseVector< typename CompareTypes< T, T2 >::supertype > & | dest, |
const DenseVector< T2 > & | arg | ||
) | const |
Performs the matrix-vector multiplication, dest
:= (*this) * arg
on mixed types.
Definition at line 330 of file dense_matrix_impl.h.
void libMesh::DenseMatrix< T >::vector_mult_add | ( | DenseVector< T > & | dest, |
const T | factor, | ||
const DenseVector< T > & | arg | ||
) | const |
Performs the scaled matrix-vector multiplication, dest
+= factor
* (*this) * arg
.
Definition at line 425 of file dense_matrix_impl.h.
Referenced by libMesh::DofMap::build_constraint_matrix_and_vector().
void libMesh::DenseMatrix< T >::vector_mult_add | ( | DenseVector< typename CompareTypes< T, typename CompareTypes< T2, T3 >::supertype >::supertype > & | dest, |
const T2 | factor, | ||
const DenseVector< T3 > & | arg | ||
) | const |
Performs the scaled matrix-vector multiplication, dest
+= factor
* (*this) * arg
.
on mixed types
Definition at line 450 of file dense_matrix_impl.h.
void libMesh::DenseMatrix< T >::vector_mult_transpose | ( | DenseVector< T > & | dest, |
const DenseVector< T > & | arg | ||
) | const |
Performs the matrix-vector multiplication, dest
:= (*this)^T * arg
.
Definition at line 355 of file dense_matrix_impl.h.
Referenced by libMesh::DofMap::constrain_element_residual(), libMesh::DofMap::heterogeneously_constrain_element_jacobian_and_residual(), and libMesh::DofMap::heterogeneously_constrain_element_residual().
void libMesh::DenseMatrix< T >::vector_mult_transpose | ( | DenseVector< typename CompareTypes< T, T2 >::supertype > & | dest, |
const DenseVector< T2 > & | arg | ||
) | const |
Performs the matrix-vector multiplication, dest
:= (*this)^T * arg
.
on mixed types
Definition at line 394 of file dense_matrix_impl.h.
|
inlinefinaloverridevirtual |
Sets all elements of the matrix to 0 and resets any decomposition flag which may have been previously set.
This allows e.g. a new LU decomposition to be computed while reusing the same storage.
Implements libMesh::DenseMatrixBase< T >.
Definition at line 911 of file dense_matrix.h.
Referenced by libMesh::HPCoarsenTest::add_projection(), libMesh::FEMSystem::assembly(), libMesh::FEGenericBase< FEOutputType< T >::type >::coarsened_dof_values(), libMesh::BoundaryProjectSolution::operator()(), libMesh::TransientRBEvaluation::rb_solve(), libMesh::RBEvaluation::rb_solve(), libMesh::HPCoarsenTest::select_refinement(), and DiagonalMatrixTest::testNumerics().
|
private |
This flag keeps track of which type of decomposition has been performed on the matrix.
Definition at line 650 of file dense_matrix.h.
|
protectedinherited |
The row dimension.
Definition at line 176 of file dense_matrix_base.h.
Referenced by libMesh::DenseMatrixBase< T >::m().
|
protectedinherited |
The column dimension.
Definition at line 181 of file dense_matrix_base.h.
Referenced by libMesh::DenseMatrixBase< T >::n().
|
private |
Definition at line 740 of file dense_matrix.h.
|
private |
The actual data values, stored as a 1D array.
Definition at line 600 of file dense_matrix.h.
Referenced by libMesh::DenseMatrix< Real >::get_values().
bool libMesh::DenseMatrix< T >::use_blas_lapack |
Computes the inverse of the dense matrix (assuming it is invertible) by first computing the LU decomposition and then performing multiple back substitution steps.
Follows the algorithm from Numerical Recipes in C that is available on the web.
This routine is commented out since it is not really a memory- or computationally- efficient implementation. Also, you typically don't need the actual inverse for anything, and can use something like lu_solve() instead. Run-time selectable option to turn on/off BLAS support. This was primarily used for testing purposes, and could be removed...
Definition at line 585 of file dense_matrix.h.