www.mooseframework.org
Functions
MatrixTools Namespace Reference

Functions

void inverse (const std::vector< std::vector< Real >> &m, std::vector< std::vector< Real >> &m_inv)
 Inverse the dense square matrix m using LAPACK routines. More...
 
void inverse (std::vector< PetscScalar > &A, unsigned int n)
 Inverts the dense "matrix" A using LAPACK routines. More...
 

Function Documentation

void MatrixTools::inverse ( const std::vector< std::vector< Real >> &  m,
std::vector< std::vector< Real >> &  m_inv 
)

Inverse the dense square matrix m using LAPACK routines.

If you need to invert a matrix "in place", make the two arguments the same @ param m The matrix to invert @ param[out] m_inv The inverse of m, which must be of the same size as m. @ return if zero then the inversion was successful. Otherwise m was not square, contained illegal entries or was singular

Definition at line 25 of file MatrixTools.C.

Referenced by RankFourTensor::invSymm().

26 {
27  unsigned int n = m.size();
28 
29  // check the matrix m exists and is square
30  if (n == 0)
31  throw MooseException("Input matrix empty during matrix inversion.");
32  if (n != m_inv.size() || n != m[0].size() || n != m_inv[0].size())
33  throw MooseException("Input and output matrix are not same size square matrices.");
34 
35  // build the vectorial representation
36  std::vector<PetscScalar> A;
37  for (const auto & rowvec : m)
38  for (const auto & matrix_entry : rowvec)
39  A.push_back(matrix_entry);
40 
41  inverse(A, n);
42 
43  // build the inverse
44  unsigned int i = 0;
45  for (auto & rowvec : m_inv)
46  for (auto & inv_entry : rowvec)
47  inv_entry = A[i++];
48 }
void inverse(const std::vector< std::vector< Real >> &m, std::vector< std::vector< Real >> &m_inv)
Inverse the dense square matrix m using LAPACK routines.
Definition: MatrixTools.C:25
static PetscErrorCode Mat * A
PetscInt m
Provides a way for users to bail out of the current solve.
PetscInt n
void MatrixTools::inverse ( std::vector< PetscScalar > &  A,
unsigned int  n 
)

Inverts the dense "matrix" A using LAPACK routines.

Parameters
Aupon input this is a row vector representing a square matrix of size sqrt(n)*sqrt(n). Upon output it is the inverse (as a row-vector)
nsize of the vector A
Returns
if zero then inversion was successful. Otherwise A contained illegal entries or was singular

Definition at line 51 of file MatrixTools.C.

52 {
53  mooseAssert(n >= 1, "MatrixTools::inverse - n (leading dimension) needs to be positive");
54  mooseAssert(n <= std::numeric_limits<int>::max(),
55  "MatrixTools::inverse - n (leading dimension) too large");
56 
57  std::vector<PetscBLASInt> ipiv(n);
58  std::vector<PetscScalar> buffer(n * 64);
59 
60  // Following does a LU decomposition of "square matrix A"
61  // upon return "A = P*L*U" if return_value == 0
62  // Here I use quotes because A is actually an array of length n^2, not a matrix of size n-by-n
63  int return_value;
64  LAPACKgetrf_(reinterpret_cast<int *>(&n),
65  reinterpret_cast<int *>(&n),
66  &A[0],
67  reinterpret_cast<int *>(&n),
68  &ipiv[0],
69  &return_value);
70 
71  if (return_value != 0)
72  throw MooseException(
73  return_value < 0
74  ? "Argument " + Moose::stringify(-return_value) +
75  " was invalid during LU factorization in MatrixTools::inverse."
76  : "Matrix on-diagonal entry " + Moose::stringify(return_value) +
77  " was exactly zero during LU factorization in MatrixTools::inverse.");
78 
79  // get the inverse of A
80  int buffer_size = buffer.size();
81 #if PETSC_VERSION_LESS_THAN(3, 5, 0)
82  FORTRAN_CALL(dgetri)
83  (reinterpret_cast<int *>(&n),
84  &A[0],
85  reinterpret_cast<int *>(&n),
86  &ipiv[0],
87  &buffer[0],
88  &buffer_size,
89  &return_value);
90 #else
91  LAPACKgetri_(reinterpret_cast<int *>(&n),
92  &A[0],
93  reinterpret_cast<int *>(&n),
94  &ipiv[0],
95  &buffer[0],
96  &buffer_size,
97  &return_value);
98 #endif
99 
100  if (return_value != 0)
101  throw MooseException(return_value < 0
102  ? "Argument " + Moose::stringify(-return_value) +
103  " was invalid during invert in MatrixTools::inverse."
104  : "Matrix on-diagonal entry " + Moose::stringify(return_value) +
105  " was exactly zero during invert in MatrixTools::inverse.");
106 }
void FORTRAN_CALL() dgetri(...)
static PetscErrorCode Mat * A
std::string stringify(const T &t)
conversion to string
Definition: Conversion.h:66
Provides a way for users to bail out of the current solve.
PetscInt n