libMesh
dense_matrix.h
Go to the documentation of this file.
1 // The libMesh Finite Element Library.
2 // Copyright (C) 2002-2017 Benjamin S. Kirk, John W. Peterson, Roy H. Stogner
3 
4 // This library is free software; you can redistribute it and/or
5 // modify it under the terms of the GNU Lesser General Public
6 // License as published by the Free Software Foundation; either
7 // version 2.1 of the License, or (at your option) any later version.
8 
9 // This library is distributed in the hope that it will be useful,
10 // but WITHOUT ANY WARRANTY; without even the implied warranty of
11 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
12 // Lesser General Public License for more details.
13 
14 // You should have received a copy of the GNU Lesser General Public
15 // License along with this library; if not, write to the Free Software
16 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
17 
18 
19 
20 #ifndef LIBMESH_DENSE_MATRIX_H
21 #define LIBMESH_DENSE_MATRIX_H
22 
23 // Local Includes
24 #include "libmesh/libmesh_common.h"
25 #include "libmesh/dense_matrix_base.h"
26 
27 // For the definition of PetscBLASInt.
28 #if (LIBMESH_HAVE_PETSC)
29 # include "libmesh/petsc_macro.h"
30 # include <petscblaslapack.h>
31 #endif
32 
33 // C++ includes
34 #include <vector>
35 #include <algorithm>
36 
37 namespace libMesh
38 {
39 
40 // Forward Declarations
41 template <typename T> class DenseVector;
42 
53 template<typename T>
54 class DenseMatrix : public DenseMatrixBase<T>
55 {
56 public:
57 
61  DenseMatrix(const unsigned int new_m=0,
62  const unsigned int new_n=0);
63 
67  virtual ~DenseMatrix() {}
68 
69 
70  virtual void zero() libmesh_override;
71 
75  T operator() (const unsigned int i,
76  const unsigned int j) const;
77 
81  T & operator() (const unsigned int i,
82  const unsigned int j);
83 
84  virtual T el(const unsigned int i,
85  const unsigned int j) const libmesh_override
86  { return (*this)(i,j); }
87 
88  virtual T & el(const unsigned int i,
89  const unsigned int j) libmesh_override
90  { return (*this)(i,j); }
91 
92  virtual void left_multiply (const DenseMatrixBase<T> & M2) libmesh_override;
93 
97  template <typename T2>
98  void left_multiply (const DenseMatrixBase<T2> & M2);
99 
100  virtual void right_multiply (const DenseMatrixBase<T> & M2) libmesh_override;
101 
105  template <typename T2>
106  void right_multiply (const DenseMatrixBase<T2> & M2);
107 
112  void vector_mult (DenseVector<T> & dest,
113  const DenseVector<T> & arg) const;
114 
120  template <typename T2>
122  const DenseVector<T2> & arg) const;
123 
129  const DenseVector<T> & arg) const;
130 
136  template <typename T2>
138  const DenseVector<T2> & arg) const;
139 
144  void vector_mult_add (DenseVector<T> & dest,
145  const T factor,
146  const DenseVector<T> & arg) const;
147 
153  template <typename T2, typename T3>
154  void vector_mult_add (DenseVector<typename CompareTypes<T, typename CompareTypes<T2,T3>::supertype>::supertype> & dest,
155  const T2 factor,
156  const DenseVector<T3> & arg) const;
157 
161  void get_principal_submatrix (unsigned int sub_m, unsigned int sub_n, DenseMatrix<T> & dest) const;
162 
166  void get_principal_submatrix (unsigned int sub_m, DenseMatrix<T> & dest) const;
167 
173  DenseMatrix<T> & operator = (const DenseMatrix<T> & other_matrix);
174 
184  template <typename T2>
185  DenseMatrix<T> & operator = (const DenseMatrix<T2> & other_matrix);
186 
190  void swap(DenseMatrix<T> & other_matrix);
191 
196  void resize(const unsigned int new_m,
197  const unsigned int new_n);
198 
202  void scale (const T factor);
203 
207  void scale_column (const unsigned int col, const T factor);
208 
214  DenseMatrix<T> & operator *= (const T factor);
215 
221  template<typename T2, typename T3>
222  typename boostcopy::enable_if_c<
223  ScalarTraits<T2>::value, void >::type add (const T2 factor,
224  const DenseMatrix<T3> & mat);
225 
229  bool operator== (const DenseMatrix<T> & mat) const;
230 
234  bool operator!= (const DenseMatrix<T> & mat) const;
235 
242 
249 
254  Real min () const;
255 
260  Real max () const;
261 
270  Real l1_norm () const;
271 
280  Real linfty_norm () const;
281 
286 
291  template <typename T2>
292  void left_multiply_transpose (const DenseMatrix<T2> & A);
293 
294 
298  void right_multiply_transpose (const DenseMatrix<T> & A);
299 
304  template <typename T2>
305  void right_multiply_transpose (const DenseMatrix<T2> & A);
306 
310  T transpose (const unsigned int i,
311  const unsigned int j) const;
312 
316  void get_transpose(DenseMatrix<T> & dest) const;
317 
325  std::vector<T> & get_values() { return _val; }
326 
330  const std::vector<T> & get_values() const { return _val; }
331 
338  void condense(const unsigned int i,
339  const unsigned int j,
340  const T val,
341  DenseVector<T> & rhs)
342  { DenseMatrixBase<T>::condense (i, j, val, rhs); }
343 
349  void lu_solve (const DenseVector<T> & b,
350  DenseVector<T> & x);
351 
363  template <typename T2>
364  void cholesky_solve(const DenseVector<T2> & b,
365  DenseVector<T2> & x);
366 
375  void svd(DenseVector<Real> & sigma);
376 
388  void svd(DenseVector<Real> & sigma,
390  DenseMatrix<Number> & VT);
391 
409  void svd_solve(const DenseVector<T> & rhs,
410  DenseVector<T> & x,
411  Real rcond=std::numeric_limits<Real>::epsilon()) const;
412 
420  void evd(DenseVector<T> & lambda_real,
421  DenseVector<T> & lambda_imag);
422 
441  void evd_left(DenseVector<T> & lambda_real,
442  DenseVector<T> & lambda_imag,
443  DenseMatrix<T> & VL);
444 
463  void evd_right(DenseVector<T> & lambda_real,
464  DenseVector<T> & lambda_imag,
465  DenseMatrix<T> & VR);
466 
477  void evd_left_and_right(DenseVector<T> & lambda_real,
478  DenseVector<T> & lambda_imag,
479  DenseMatrix<T> & VL,
480  DenseMatrix<T> & VR);
481 
489  T det();
490 
502  // void inverse();
503 
510 
511 private:
512 
516  std::vector<T> _val;
517 
523  void _lu_decompose ();
524 
530  void _lu_back_substitute (const DenseVector<T> & b,
531  DenseVector<T> & x) const;
532 
540  void _cholesky_decompose();
541 
549  template <typename T2>
551  DenseVector<T2> & x) const;
552 
561 
567 
577  };
578 
586  void _multiply_blas(const DenseMatrixBase<T> & other,
587  _BLAS_Multiply_Flag flag);
588 
598  void _lu_decompose_lapack();
599 
605  void _svd_lapack(DenseVector<Real> & sigma);
606 
612  void _svd_lapack(DenseVector<Real> & sigma,
614  DenseMatrix<Number> & VT);
615 
619  void _svd_solve_lapack(const DenseVector<T> & rhs,
620  DenseVector<T> & x,
621  Real rcond) const;
622 
627  void _svd_helper (char JOBU,
628  char JOBVT,
629  std::vector<Real> & sigma_val,
630  std::vector<Number> & U_val,
631  std::vector<Number> & VT_val);
632 
641  void _evd_lapack(DenseVector<T> & lambda_real,
642  DenseVector<T> & lambda_imag,
645 
651 #if (LIBMESH_HAVE_PETSC && LIBMESH_USE_REAL_NUMBERS)
652  typedef PetscBLASInt pivot_index_t;
653 #else
654  typedef int pivot_index_t;
655 #endif
656  std::vector<pivot_index_t> _pivots;
657 
668  DenseVector<T> & x);
669 
682  void _matvec_blas(T alpha, T beta,
683  DenseVector<T> & dest,
684  const DenseVector<T> & arg,
685  bool trans=false) const;
686 };
687 
688 
689 
690 
691 
692 // ------------------------------------------------------------
696 namespace DenseMatrices
697 {
698 
704 
713 
714 }
715 
716 
717 
718 using namespace DenseMatrices;
719 
720 
721 
722 
723 
724 // ------------------------------------------------------------
725 // Dense Matrix member functions
726 template<typename T>
727 inline
728 DenseMatrix<T>::DenseMatrix(const unsigned int new_m,
729  const unsigned int new_n) :
730  DenseMatrixBase<T>(new_m,new_n),
731 #if defined(LIBMESH_HAVE_PETSC) && defined(LIBMESH_USE_REAL_NUMBERS) && defined(LIBMESH_DEFAULT_DOUBLE_PRECISION)
732  use_blas_lapack(true),
733 #else
734  use_blas_lapack(false),
735 #endif
736  _val(),
738 {
739  this->resize(new_m,new_n);
740 }
741 
742 
743 
744 template<typename T>
745 inline
747 {
748  std::swap(this->_m, other_matrix._m);
749  std::swap(this->_n, other_matrix._n);
750  _val.swap(other_matrix._val);
752  _decomposition_type = other_matrix._decomposition_type;
753  other_matrix._decomposition_type = _temp;
754 }
755 
756 
757 template <typename T>
758 template <typename T2>
759 inline
762 {
763  unsigned int mat_m = mat.m(), mat_n = mat.n();
764  this->resize(mat_m, mat_n);
765  for (unsigned int i=0; i<mat_m; i++)
766  for (unsigned int j=0; j<mat_n; j++)
767  (*this)(i,j) = mat(i,j);
768 
769  return *this;
770 }
771 
772 
773 
774 template<typename T>
775 inline
776 void DenseMatrix<T>::resize(const unsigned int new_m,
777  const unsigned int new_n)
778 {
779  _val.resize(new_m*new_n);
780 
781  this->_m = new_m;
782  this->_n = new_n;
783 
784  // zero and set decomposition_type to NONE
785  this->zero();
786 }
787 
788 
789 
790 template<typename T>
791 inline
793 {
795 
796  std::fill (_val.begin(), _val.end(), static_cast<T>(0));
797 }
798 
799 
800 
801 template<typename T>
802 inline
804 {
805  this->_m = other_matrix._m;
806  this->_n = other_matrix._n;
807 
808  _val = other_matrix._val;
809  _decomposition_type = other_matrix._decomposition_type;
810 
811  return *this;
812 }
813 
814 
815 
816 template<typename T>
817 inline
818 T DenseMatrix<T>::operator () (const unsigned int i,
819  const unsigned int j) const
820 {
821  libmesh_assert_less (i*j, _val.size());
822  libmesh_assert_less (i, this->_m);
823  libmesh_assert_less (j, this->_n);
824 
825 
826  // return _val[(i) + (this->_m)*(j)]; // col-major
827  return _val[(i)*(this->_n) + (j)]; // row-major
828 }
829 
830 
831 
832 template<typename T>
833 inline
834 T & DenseMatrix<T>::operator () (const unsigned int i,
835  const unsigned int j)
836 {
837  libmesh_assert_less (i*j, _val.size());
838  libmesh_assert_less (i, this->_m);
839  libmesh_assert_less (j, this->_n);
840 
841  //return _val[(i) + (this->_m)*(j)]; // col-major
842  return _val[(i)*(this->_n) + (j)]; // row-major
843 }
844 
845 
846 
847 
848 
849 template<typename T>
850 inline
851 void DenseMatrix<T>::scale (const T factor)
852 {
853  for (std::size_t i=0; i<_val.size(); i++)
854  _val[i] *= factor;
855 }
856 
857 
858 template<typename T>
859 inline
860 void DenseMatrix<T>::scale_column (const unsigned int col, const T factor)
861 {
862  for (unsigned int i=0; i<this->m(); i++)
863  (*this)(i, col) *= factor;
864 }
865 
866 
867 
868 template<typename T>
869 inline
871 {
872  this->scale(factor);
873  return *this;
874 }
875 
876 
877 
878 template<typename T>
879 template<typename T2, typename T3>
880 inline
881 typename boostcopy::enable_if_c<
882  ScalarTraits<T2>::value, void >::type
883 DenseMatrix<T>::add (const T2 factor,
884  const DenseMatrix<T3> & mat)
885 {
886  libmesh_assert_equal_to (this->m(), mat.m());
887  libmesh_assert_equal_to (this->n(), mat.n());
888 
889  for (unsigned int i=0; i<this->m(); i++)
890  for (unsigned int j=0; j<this->n(); j++)
891  (*this)(i,j) += factor * mat(i,j);
892 }
893 
894 
895 
896 template<typename T>
897 inline
899 {
900  for (std::size_t i=0; i<_val.size(); i++)
901  if (_val[i] != mat._val[i])
902  return false;
903 
904  return true;
905 }
906 
907 
908 
909 template<typename T>
910 inline
912 {
913  for (std::size_t i=0; i<_val.size(); i++)
914  if (_val[i] != mat._val[i])
915  return true;
916 
917  return false;
918 }
919 
920 
921 
922 template<typename T>
923 inline
925 {
926  for (std::size_t i=0; i<_val.size(); i++)
927  _val[i] += mat._val[i];
928 
929  return *this;
930 }
931 
932 
933 
934 template<typename T>
935 inline
937 {
938  for (std::size_t i=0; i<_val.size(); i++)
939  _val[i] -= mat._val[i];
940 
941  return *this;
942 }
943 
944 
945 
946 template<typename T>
947 inline
949 {
950  libmesh_assert (this->_m);
951  libmesh_assert (this->_n);
952  Real my_min = libmesh_real((*this)(0,0));
953 
954  for (unsigned int i=0; i!=this->_m; i++)
955  {
956  for (unsigned int j=0; j!=this->_n; j++)
957  {
958  Real current = libmesh_real((*this)(i,j));
959  my_min = (my_min < current? my_min : current);
960  }
961  }
962  return my_min;
963 }
964 
965 
966 
967 template<typename T>
968 inline
970 {
971  libmesh_assert (this->_m);
972  libmesh_assert (this->_n);
973  Real my_max = libmesh_real((*this)(0,0));
974 
975  for (unsigned int i=0; i!=this->_m; i++)
976  {
977  for (unsigned int j=0; j!=this->_n; j++)
978  {
979  Real current = libmesh_real((*this)(i,j));
980  my_max = (my_max > current? my_max : current);
981  }
982  }
983  return my_max;
984 }
985 
986 
987 
988 template<typename T>
989 inline
991 {
992  libmesh_assert (this->_m);
993  libmesh_assert (this->_n);
994 
995  Real columnsum = 0.;
996  for (unsigned int i=0; i!=this->_m; i++)
997  {
998  columnsum += std::abs((*this)(i,0));
999  }
1000  Real my_max = columnsum;
1001  for (unsigned int j=1; j!=this->_n; j++)
1002  {
1003  columnsum = 0.;
1004  for (unsigned int i=0; i!=this->_m; i++)
1005  {
1006  columnsum += std::abs((*this)(i,j));
1007  }
1008  my_max = (my_max > columnsum? my_max : columnsum);
1009  }
1010  return my_max;
1011 }
1012 
1013 
1014 
1015 template<typename T>
1016 inline
1018 {
1019  libmesh_assert (this->_m);
1020  libmesh_assert (this->_n);
1021 
1022  Real rowsum = 0.;
1023  for (unsigned int j=0; j!=this->_n; j++)
1024  {
1025  rowsum += std::abs((*this)(0,j));
1026  }
1027  Real my_max = rowsum;
1028  for (unsigned int i=1; i!=this->_m; i++)
1029  {
1030  rowsum = 0.;
1031  for (unsigned int j=0; j!=this->_n; j++)
1032  {
1033  rowsum += std::abs((*this)(i,j));
1034  }
1035  my_max = (my_max > rowsum? my_max : rowsum);
1036  }
1037  return my_max;
1038 }
1039 
1040 
1041 
1042 template<typename T>
1043 inline
1044 T DenseMatrix<T>::transpose (const unsigned int i,
1045  const unsigned int j) const
1046 {
1047  // Implement in terms of operator()
1048  return (*this)(j,i);
1049 }
1050 
1051 
1052 
1053 
1054 
1055 // template<typename T>
1056 // inline
1057 // void DenseMatrix<T>::condense(const unsigned int iv,
1058 // const unsigned int jv,
1059 // const T val,
1060 // DenseVector<T> & rhs)
1061 // {
1062 // libmesh_assert_equal_to (this->_m, rhs.size());
1063 // libmesh_assert_equal_to (iv, jv);
1064 
1065 
1066 // // move the known value into the RHS
1067 // // and zero the column
1068 // for (unsigned int i=0; i<this->m(); i++)
1069 // {
1070 // rhs(i) -= ((*this)(i,jv))*val;
1071 // (*this)(i,jv) = 0.;
1072 // }
1073 
1074 // // zero the row
1075 // for (unsigned int j=0; j<this->n(); j++)
1076 // (*this)(iv,j) = 0.;
1077 
1078 // (*this)(iv,jv) = 1.;
1079 // rhs(iv) = val;
1080 
1081 // }
1082 
1083 
1084 
1085 
1086 } // namespace libMesh
1087 
1088 
1089 
1090 
1091 #endif // LIBMESH_DENSE_MATRIX_H
T libmesh_real(T a)
bool operator!=(const DenseMatrix< T > &mat) const
Definition: dense_matrix.h:911
double abs(double a)
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...
Definition: dense_matrix.C:989
void _lu_decompose_lapack()
Computes an LU factorization of the matrix using the Lapack routine "getrf".
void _lu_back_substitute(const DenseVector< T > &b, DenseVector< T > &x) const
Solves the system Ax=b through back substitution.
Definition: dense_matrix.C:648
void swap(DenseMatrix< T > &other_matrix)
STL-like swap method.
Definition: dense_matrix.h:746
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.
Definition: dense_matrix.h:338
unsigned int n() const
DenseMatrix< Complex > ComplexDenseMatrix
This typedef may be either a real-only matrix, or a truly complex matrix, depending on how Number was...
Definition: dense_matrix.h:712
Real max() const
Definition: dense_matrix.h:969
void svd_solve(const DenseVector< T > &rhs, DenseVector< T > &x, Real rcond=std::numeric_limits< Real >::epsilon()) const
Solve the system of equations for in the least-squares sense.
Definition: dense_matrix.C:796
virtual void zero() libmesh_override
Set every element in the matrix to 0.
Definition: dense_matrix.h:792
_BLAS_Multiply_Flag
Enumeration used to determine the behavior of the _multiply_blas function.
Definition: dense_matrix.h:572
const std::vector< T > & get_values() const
Definition: dense_matrix.h:330
DenseMatrix< T > & operator+=(const DenseMatrix< T > &mat)
Adds mat to this matrix.
Definition: dense_matrix.h:924
DecompositionType _decomposition_type
This flag keeps track of which type of decomposition has been performed on the matrix.
Definition: dense_matrix.h:566
void vector_mult_transpose(DenseVector< T > &dest, const DenseVector< T > &arg) const
Performs the matrix-vector multiplication, dest := (*this)^T * arg.
Definition: dense_matrix.C:438
bool operator==(const DenseMatrix< T > &mat) const
Definition: dense_matrix.h:898
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...
Definition: dense_matrix.C:816
Real l1_norm() const
Definition: dense_matrix.h:990
const class libmesh_nullptr_t libmesh_nullptr
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...
Definition: dense_matrix.C:838
void _lu_back_substitute_lapack(const DenseVector< T > &b, DenseVector< T > &x)
Companion function to _lu_decompose_lapack().
void get_transpose(DenseMatrix< T > &dest) const
Put the tranposed matrix into dest.
Definition: dense_matrix.C:576
DenseMatrix< T > & operator-=(const DenseMatrix< T > &mat)
Subtracts mat from this matrix.
Definition: dense_matrix.h:936
void _svd_lapack(DenseVector< Real > &sigma)
Computes an SVD of the matrix using the Lapack routine "getsvd".
The libMesh namespace provides an interface to certain functionality in the library.
bool use_blas_lapack
Computes the inverse of the dense matrix (assuming it is invertible) by first computing the LU decomp...
Definition: dense_matrix.h:509
DenseMatrix(const unsigned int new_m=0, const unsigned int new_n=0)
Constructor.
Definition: dense_matrix.h:728
libmesh_assert(j)
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.
Real linfty_norm() const
void scale_column(const unsigned int col, const T factor)
Multiplies every element in the column col matrix by factor.
Definition: dense_matrix.h:860
virtual ~DenseMatrix()
Destructor.
Definition: dense_matrix.h:67
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. ...
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.
PetscErrorCode Vec x
virtual void right_multiply(const DenseMatrixBase< T > &M2) libmesh_override
Performs the operation: (*this) <- (*this) * M3.
Definition: dense_matrix.C:210
unsigned int m() const
void svd(DenseVector< Real > &sigma)
Compute the singular value decomposition of the matrix.
Definition: dense_matrix.C:777
void _lu_decompose()
Form the LU decomposition of the matrix.
Definition: dense_matrix.C:699
Real min() const
Definition: dense_matrix.h:948
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. ...
Definition: dense_matrix.C:508
T transpose(const unsigned int i, const unsigned int j) const
boostcopy::enable_if_c< ScalarTraits< T2 >::value, void >::type add(const T2 factor, const DenseMatrix< T3 > &mat)
Adds factor times mat to this matrix.
Definition: dense_matrix.h:883
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.
Definition: dense_matrix.C:553
DenseMatrix< Real > RealDenseMatrix
Convenient definition of a real-only dense matrix.
Definition: dense_matrix.h:703
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...
Definition: dense_matrix.C:827
virtual T el(const unsigned int i, const unsigned int j) const libmesh_override
Definition: dense_matrix.h:84
void _cholesky_decompose()
Decomposes a symmetric positive definite matrix into a product of two lower triangular matrices accor...
Definition: dense_matrix.C:943
std::vector< T > & get_values()
Definition: dense_matrix.h:325
virtual void left_multiply(const DenseMatrixBase< T > &M2) libmesh_override
Performs the operation: (*this) <- M2 * (*this)
Definition: dense_matrix.C:37
T operator()(const unsigned int i, const unsigned int j) const
Definition: dense_matrix.h:818
unsigned int _n
The column dimension.
void _svd_solve_lapack(const DenseVector< T > &rhs, DenseVector< T > &x, Real rcond) const
Called by svd_solve(rhs).
void _evd_lapack(DenseVector< T > &lambda_real, DenseVector< T > &lambda_imag, DenseMatrix< T > *VL=libmesh_nullptr, DenseMatrix< T > *VR=libmesh_nullptr)
Computes the eigenvalues of the matrix using the Lapack routine "DGEEV".
void right_multiply_transpose(const DenseMatrix< T > &A)
Right multiplies by the transpose of the matrix A.
Definition: dense_matrix.C:257
void scale(const T factor)
Multiplies every element in the matrix by factor.
Definition: dense_matrix.h:851
void lu_solve(const DenseVector< T > &b, DenseVector< T > &x)
Solve the system Ax=b given the input vector b.
Definition: dense_matrix.C:589
DenseMatrix< T > & operator=(const DenseMatrix< T > &other_matrix)
Assignment operator.
Definition: dense_matrix.h:803
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
void evd(DenseVector< T > &lambda_real, DenseVector< T > &lambda_imag)
Compute the eigenvalues (both real and imaginary parts) of a general matrix.
Definition: dense_matrix.C:806
void left_multiply_transpose(const DenseMatrix< T > &A)
Left multiplies by the transpose of the matrix A.
Definition: dense_matrix.C:85
static PetscErrorCode Mat * A
void swap(Iterator &lhs, Iterator &rhs)
swap, used to implement op=
virtual T & el(const unsigned int i, const unsigned int j) libmesh_override
Definition: dense_matrix.h:88
void resize(const unsigned int new_m, const unsigned int new_n)
Resize the matrix.
Definition: dense_matrix.h:776
std::vector< T > _val
The actual data values, stored as a 1D array.
Definition: dense_matrix.h:516
Defines a dense vector for use in Finite Element-type computations.
DenseMatrix< T > & operator*=(const T factor)
Multiplies every element in the matrix by factor.
Definition: dense_matrix.h:870
Defines an abstract dense matrix base class for use in Finite Element-type computations.
if(!eq) SETERRQ2(((PetscObject) dm) -> comm, PETSC_ERR_ARG_WRONG,"DM of type %s, not of type %s",((PetscObject) dm) ->type, DMLIBMESH)
void vector_mult(DenseVector< T > &dest, const DenseVector< T > &arg) const
Performs the matrix-vector multiplication, dest := (*this) * arg.
Definition: dense_matrix.C:382
DecompositionType
The decomposition schemes above change the entries of the matrix A.
Definition: dense_matrix.h:560
Defines a dense matrix for use in Finite Element-type computations.
Definition: dof_map.h:64
std::vector< pivot_index_t > _pivots
Definition: dense_matrix.h:656
void cholesky_solve(const DenseVector< T2 > &b, DenseVector< T2 > &x)
For symmetric positive definite (SPD) matrices.
Definition: dense_matrix.C:911
PetscBLASInt pivot_index_t
Array used to store pivot indices.
Definition: dense_matrix.h:652
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.
unsigned int _m
The row dimension.