libMesh
Public Member Functions | Protected Member Functions | Static Protected Member Functions | Protected Attributes | Private Attributes | List of all members
libMesh::DenseSubMatrix< T > Class Template Reference

Defines a dense submatrix for use in Finite Element-type computations. More...

#include <dense_submatrix.h>

Inheritance diagram for libMesh::DenseSubMatrix< T >:
[legend]

Public Member Functions

 DenseSubMatrix (DenseMatrix< T > &new_parent, const unsigned int ioff=0, const unsigned int joff=0, const unsigned int m=0, const unsigned int n=0)
 Constructor. More...
 
 DenseSubMatrix (DenseSubMatrix &&)=default
 The 5 special functions can be defaulted for this class, as it does not manage any memory itself. More...
 
 DenseSubMatrix (const DenseSubMatrix &)=default
 
DenseSubMatrixoperator= (const DenseSubMatrix &)=default
 
DenseSubMatrixoperator= (DenseSubMatrix &&)=default
 
virtual ~DenseSubMatrix ()=default
 
DenseMatrix< T > & parent ()
 
virtual void zero () override final
 Set every element in the matrix to 0. More...
 
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...
 
virtual void right_multiply (const DenseMatrixBase< T > &M3) override final
 Performs the operation: (*this) <- (*this) * M3. More...
 
void reposition (const unsigned int ioff, const unsigned int joff, const unsigned int new_m, const unsigned int new_n)
 Changes the location of the submatrix in the parent matrix. More...
 
unsigned int i_off () const
 
unsigned int j_off () const
 
void condense (const unsigned int i, const unsigned int j, const T val, DenseSubVector< T > &rhs)
 Condense-out the (i,j) entry of the matrix, forcing it to take on the value val. More...
 
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...
 

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 Attributes

DenseMatrix< T > & _parent_matrix
 The parent matrix that contains this submatrix. More...
 
unsigned int _i_off
 The row offset into the parent matrix. More...
 
unsigned int _j_off
 The column offset into the parent matrix. More...
 

Detailed Description

template<typename T>
class libMesh::DenseSubMatrix< T >

Defines a dense submatrix for use in Finite Element-type computations.

Useful for storing element stiffness matrices before summation into a global matrix, particularly when you have systems of equations. All overridden virtual functions are documented in dense_matrix_base.h.

Author
Benjamin S. Kirk
Date
2003

Definition at line 45 of file dense_submatrix.h.

Constructor & Destructor Documentation

◆ DenseSubMatrix() [1/3]

template<typename T >
libMesh::DenseSubMatrix< T >::DenseSubMatrix ( DenseMatrix< T > &  new_parent,
const unsigned int  ioff = 0,
const unsigned int  joff = 0,
const unsigned int  m = 0,
const unsigned int  n = 0 
)
inline

Constructor.

Creates a dense submatrix of the matrix parent. The submatrix has dimensions \((m \times n)\), and the \((0,0)\) entry of the submatrix is located at the \((ioff,joff)\) location in the parent matrix.

Definition at line 159 of file dense_submatrix.h.

References libMesh::DenseSubMatrix< T >::reposition().

163  :
164  DenseMatrixBase<T>(new_m,new_n),
165  _parent_matrix(new_parent)
166 {
167  this->reposition (ioff, joff, new_m, new_n);
168 }
void reposition(const unsigned int ioff, const unsigned int joff, const unsigned int new_m, const unsigned int new_n)
Changes the location of the submatrix in the parent matrix.
DenseMatrix< T > & _parent_matrix
The parent matrix that contains this submatrix.

◆ DenseSubMatrix() [2/3]

template<typename T>
libMesh::DenseSubMatrix< T >::DenseSubMatrix ( DenseSubMatrix< T > &&  )
default

The 5 special functions can be defaulted for this class, as it does not manage any memory itself.

◆ DenseSubMatrix() [3/3]

template<typename T>
libMesh::DenseSubMatrix< T >::DenseSubMatrix ( const DenseSubMatrix< T > &  )
default

◆ ~DenseSubMatrix()

template<typename T>
virtual libMesh::DenseSubMatrix< T >::~DenseSubMatrix ( )
virtualdefault

Member Function Documentation

◆ add()

template<typename T >
template<typename T2 , typename T3 >
boostcopy::enable_if_c< ScalarTraits< T2 >::value, void >::type libMesh::DenseMatrixBase< T >::add ( const T2  factor,
const DenseMatrixBase< T3 > &  mat 
)
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().

197 {
198  libmesh_assert_equal_to (this->m(), mat.m());
199  libmesh_assert_equal_to (this->n(), mat.n());
200 
201  for (auto j : make_range(this->n()))
202  for (auto i : make_range(this->m()))
203  this->el(i,j) += factor*mat.el(i,j);
204 }
unsigned int m() const
virtual T el(const unsigned int i, const unsigned int j) const =0
IntRange< T > make_range(T beg, T end)
The 2-parameter make_range() helper function returns an IntRange<T> when both input parameters are of...
Definition: int_range.h:134
unsigned int n() const

◆ condense() [1/2]

template<typename T>
void libMesh::DenseSubMatrix< T >::condense ( const unsigned int  i,
const unsigned int  j,
const T  val,
DenseSubVector< T > &  rhs 
)
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 126 of file dense_submatrix.h.

References libMesh::DenseSubMatrix< T >::i_off(), libMesh::DenseSubMatrix< T >::j_off(), libMesh::DenseSubVector< T >::parent(), and libMesh::DenseSubMatrix< T >::parent().

130  {
131  this->parent().condense(this->i_off()+i,
132  this->j_off()+j,
133  val, rhs.parent());
134  }
unsigned int j_off() const
DenseMatrix< T > & parent()
unsigned int i_off() const

◆ condense() [2/2]

template<typename T >
void libMesh::DenseMatrixBase< T >::condense ( const unsigned int  i,
const unsigned int  j,
const T  val,
DenseVectorBase< T > &  rhs 
)
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().

77  {
78  libmesh_assert_equal_to (this->_m, rhs.size());
79  libmesh_assert_equal_to (iv, jv);
80 
81 
82  // move the known value into the RHS
83  // and zero the column
84  for (auto i : make_range(this->m()))
85  {
86  rhs.el(i) -= this->el(i,jv)*val;
87  this->el(i,jv) = 0.;
88  }
89 
90  // zero the row
91  for (auto j : make_range(this->n()))
92  this->el(iv,j) = 0.;
93 
94  this->el(iv,jv) = 1.;
95  rhs.el(iv) = val;
96 
97  }
unsigned int m() const
virtual T el(const unsigned int i, const unsigned int j) const =0
IntRange< T > make_range(T beg, T end)
The 2-parameter make_range() helper function returns an IntRange<T> when both input parameters are of...
Definition: int_range.h:134
unsigned int n() const
unsigned int _m
The row dimension.

◆ diagonal()

template<typename T >
DenseVector< T > libMesh::DenseMatrixBase< T >::diagonal ( ) const
inherited

Return the matrix diagonal.

Definition at line 37 of file dense_matrix_base_impl.h.

Referenced by libMesh::DiagonalMatrix< T >::add_matrix().

38  {
39  DenseVector<T> ret(_m);
40  for (decltype(_m) i = 0; i < _m; ++i)
41  ret(i) = el(i, i);
42  return ret;
43  }
virtual T el(const unsigned int i, const unsigned int j) const =0
unsigned int _m
The row dimension.

◆ el() [1/2]

template<typename T>
virtual T libMesh::DenseSubMatrix< T >::el ( const unsigned int  i,
const unsigned int  j 
) const
inlinefinaloverridevirtual
Returns
The (i,j) element of the matrix. Since internal data representations may differ, you must redefine this function.

Implements libMesh::DenseMatrixBase< T >.

Definition at line 90 of file dense_submatrix.h.

92  { return (*this)(i,j); }

◆ el() [2/2]

template<typename T>
virtual T& libMesh::DenseSubMatrix< T >::el ( const unsigned int  i,
const unsigned int  j 
)
inlinefinaloverridevirtual
Returns
The (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 94 of file dense_submatrix.h.

96  { return (*this)(i,j); }

◆ i_off()

template<typename T>
unsigned int libMesh::DenseSubMatrix< T >::i_off ( ) const
inline
Returns
The row offset into the parent matrix.

Definition at line 113 of file dense_submatrix.h.

References libMesh::DenseSubMatrix< T >::_i_off.

Referenced by libMesh::DenseSubMatrix< T >::condense().

113 { return _i_off; }
unsigned int _i_off
The row offset into the parent matrix.

◆ j_off()

template<typename T>
unsigned int libMesh::DenseSubMatrix< T >::j_off ( ) const
inline
Returns
The column offset into the parent matrix.

Definition at line 118 of file dense_submatrix.h.

References libMesh::DenseSubMatrix< T >::_j_off.

Referenced by libMesh::DenseSubMatrix< T >::condense().

118 { return _j_off; }
unsigned int _j_off
The column offset into the parent matrix.

◆ left_multiply()

template<typename T >
void libMesh::DenseSubMatrix< T >::left_multiply ( const DenseMatrixBase< T > &  M2)
finaloverridevirtual

Performs the operation: (*this) <- M2 * (*this)

Implements libMesh::DenseMatrixBase< T >.

Definition at line 31 of file dense_submatrix.C.

32 {
33  // (*this) <- M2 * M3
34  // Where:
35  // (*this) = (m x n),
36  // M2 = (m x p),
37  // M3 = (p x n)
38 
39  // M3 is a simply a copy of *this
40  DenseSubMatrix<T> M3(*this);
41 
42  // Call the multiply function in the base class
43  this->multiply(*this, M2, M3);
44 }
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)...

◆ m()

template<typename T>
unsigned int libMesh::DenseMatrixBase< T >::m ( ) const
inlineinherited
Returns
The row-dimension of the matrix.

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< libMesh::Number >::add_block_matrix(), libMesh::SparseMatrix< ValOut >::add_block_matrix(), libMesh::EigenSparseMatrix< T >::add_matrix(), libMesh::DiagonalMatrix< T >::add_matrix(), libMesh::LaspackMatrix< T >::add_matrix(), libMesh::EpetraMatrix< T >::add_matrix(), libMesh::PetscMatrix< libMesh::Number >::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().

104 { return _m; }
unsigned int _m
The row dimension.

◆ multiply()

template<typename T >
void libMesh::DenseMatrixBase< T >::multiply ( DenseMatrixBase< T > &  M1,
const DenseMatrixBase< T > &  M2,
const DenseMatrixBase< T > &  M3 
)
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().

49  {
50  // Assertions to make sure we have been
51  // passed matrices of the correct dimension.
52  libmesh_assert_equal_to (M1.m(), M2.m());
53  libmesh_assert_equal_to (M1.n(), M3.n());
54  libmesh_assert_equal_to (M2.n(), M3.m());
55 
56  const unsigned int m_s = M2.m();
57  const unsigned int p_s = M2.n();
58  const unsigned int n_s = M1.n();
59 
60  // Do it this way because there is a
61  // decent chance (at least for constraint matrices)
62  // that M3(k,j) = 0. when right-multiplying.
63  for (unsigned int k=0; k<p_s; k++)
64  for (unsigned int j=0; j<n_s; j++)
65  if (M3.el(k,j) != 0.)
66  for (unsigned int i=0; i<m_s; i++)
67  M1.el(i,j) += M2.el(i,k) * M3.el(k,j);
68  }

◆ n()

template<typename T>
unsigned int libMesh::DenseMatrixBase< T >::n ( ) const
inlineinherited
Returns
The column-dimension of the matrix.

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< libMesh::Number >::add_block_matrix(), libMesh::SparseMatrix< ValOut >::add_block_matrix(), libMesh::EigenSparseMatrix< T >::add_matrix(), libMesh::DiagonalMatrix< T >::add_matrix(), libMesh::LaspackMatrix< T >::add_matrix(), libMesh::EpetraMatrix< T >::add_matrix(), libMesh::PetscMatrix< libMesh::Number >::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().

109 { return _n; }
unsigned int _n
The column dimension.

◆ operator()() [1/2]

template<typename T >
T libMesh::DenseSubMatrix< T >::operator() ( const unsigned int  i,
const unsigned int  j 
) const
inline
Returns
The (i,j) element of the submatrix.

Definition at line 204 of file dense_submatrix.h.

206 {
207  libmesh_assert_less (i, this->m());
208  libmesh_assert_less (j, this->n());
209  libmesh_assert_less (i + this->i_off(), _parent_matrix.m());
210  libmesh_assert_less (j + this->j_off(), _parent_matrix.n());
211 
212  return _parent_matrix (i + this->i_off(),
213  j + this->j_off());
214 }
unsigned int j_off() const
unsigned int m() const
unsigned int n() const
DenseMatrix< T > & _parent_matrix
The parent matrix that contains this submatrix.
unsigned int i_off() const

◆ operator()() [2/2]

template<typename T >
T & libMesh::DenseSubMatrix< T >::operator() ( const unsigned int  i,
const unsigned int  j 
)
inline
Returns
The (i,j) element of the submatrix as a writable reference.

Definition at line 219 of file dense_submatrix.h.

221 {
222  libmesh_assert_less (i, this->m());
223  libmesh_assert_less (j, this->n());
224  libmesh_assert_less (i + this->i_off(), _parent_matrix.m());
225  libmesh_assert_less (j + this->j_off(), _parent_matrix.n());
226 
227  return _parent_matrix (i + this->i_off(),
228  j + this->j_off());
229 }
unsigned int j_off() const
unsigned int m() const
unsigned int n() const
DenseMatrix< T > & _parent_matrix
The parent matrix that contains this submatrix.
unsigned int i_off() const

◆ operator=() [1/2]

template<typename T>
DenseSubMatrix& libMesh::DenseSubMatrix< T >::operator= ( const DenseSubMatrix< T > &  )
default

◆ operator=() [2/2]

template<typename T>
DenseSubMatrix& libMesh::DenseSubMatrix< T >::operator= ( DenseSubMatrix< T > &&  )
default

◆ parent()

template<typename T>
DenseMatrix<T>& libMesh::DenseSubMatrix< T >::parent ( )
inline
Returns
A reference to the parent matrix.

Definition at line 74 of file dense_submatrix.h.

References libMesh::DenseSubMatrix< T >::_parent_matrix.

Referenced by libMesh::DenseSubMatrix< T >::condense().

74 { return _parent_matrix; }
DenseMatrix< T > & _parent_matrix
The parent matrix that contains this submatrix.

◆ print()

template<typename T >
void libMesh::DenseMatrixBase< T >::print ( std::ostream &  os = libMesh::out) const
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().

126  {
127  for (auto i : make_range(this->m()))
128  {
129  for (auto j : make_range(this->n()))
130  os << std::setw(8)
131  << this->el(i,j) << " ";
132 
133  os << std::endl;
134  }
135 
136  return;
137  }
unsigned int m() const
virtual T el(const unsigned int i, const unsigned int j) const =0
IntRange< T > make_range(T beg, T end)
The 2-parameter make_range() helper function returns an IntRange<T> when both input parameters are of...
Definition: int_range.h:134
unsigned int n() const

◆ print_scientific()

template<typename T >
void libMesh::DenseMatrixBase< T >::print_scientific ( std::ostream &  os,
unsigned  precision = 8 
) const
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().

102  {
103  // save the initial format flags
104  std::ios_base::fmtflags os_flags = os.flags();
105 
106  // Print the matrix entries.
107  for (auto i : make_range(this->m()))
108  {
109  for (auto j : make_range(this->n()))
110  os << std::setw(15)
111  << std::scientific
112  << std::setprecision(precision)
113  << this->el(i,j) << " ";
114 
115  os << std::endl;
116  }
117 
118  // reset the original format flags
119  os.flags(os_flags);
120  }
unsigned int m() const
virtual T el(const unsigned int i, const unsigned int j) const =0
IntRange< T > make_range(T beg, T end)
The 2-parameter make_range() helper function returns an IntRange<T> when both input parameters are of...
Definition: int_range.h:134
unsigned int n() const

◆ reposition()

template<typename T >
void libMesh::DenseSubMatrix< T >::reposition ( const unsigned int  ioff,
const unsigned int  joff,
const unsigned int  new_m,
const unsigned int  new_n 
)
inline

Changes the location of the submatrix in the parent matrix.

Definition at line 173 of file dense_submatrix.h.

Referenced by assemble_elasticity(), assemble_laplace(), assemble_shell(), assemble_stokes(), and libMesh::DenseSubMatrix< T >::DenseSubMatrix().

177 {
178  _i_off = ioff;
179  _j_off = joff;
180  this->_m = new_m;
181  this->_n = new_n;
182 
183  // Make sure we still fit in the parent matrix.
184  libmesh_assert_less_equal ((this->i_off() + this->m()), _parent_matrix.m());
185  libmesh_assert_less_equal ((this->j_off() + this->n()), _parent_matrix.n());
186 }
unsigned int j_off() const
unsigned int m() const
unsigned int _i_off
The row offset into the parent matrix.
unsigned int _j_off
The column offset into the parent matrix.
unsigned int _n
The column dimension.
unsigned int n() const
unsigned int _m
The row dimension.
DenseMatrix< T > & _parent_matrix
The parent matrix that contains this submatrix.
unsigned int i_off() const

◆ right_multiply()

template<typename T >
void libMesh::DenseSubMatrix< T >::right_multiply ( const DenseMatrixBase< T > &  M3)
finaloverridevirtual

Performs the operation: (*this) <- (*this) * M3.

Implements libMesh::DenseMatrixBase< T >.

Definition at line 49 of file dense_submatrix.C.

50 {
51  // (*this) <- M2 * M3
52  // Where:
53  // (*this) = (m x n),
54  // M2 = (m x p),
55  // M3 = (p x n)
56 
57  // M2 is simply a copy of *this
58  DenseSubMatrix<T> M2(*this);
59 
60  // Call the multiply function in the base class
61  this->multiply(*this, M2, M3);
62 }
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)...

◆ zero()

template<typename T >
void libMesh::DenseSubMatrix< T >::zero ( )
inlinefinaloverridevirtual

Set every element in the matrix to 0.

You must redefine what you mean by zeroing the matrix since it depends on how your values are stored.

Implements libMesh::DenseMatrixBase< T >.

Definition at line 192 of file dense_submatrix.h.

References libMesh::make_range().

193 {
194  for (auto i : make_range(this->m()))
195  for (auto j : make_range(this->n()))
196  _parent_matrix(i + this->i_off(),
197  j + this->j_off()) = 0.;
198 }
unsigned int j_off() const
unsigned int m() const
IntRange< T > make_range(T beg, T end)
The 2-parameter make_range() helper function returns an IntRange<T> when both input parameters are of...
Definition: int_range.h:134
unsigned int n() const
DenseMatrix< T > & _parent_matrix
The parent matrix that contains this submatrix.
unsigned int i_off() const

Member Data Documentation

◆ _i_off

template<typename T>
unsigned int libMesh::DenseSubMatrix< T >::_i_off
private

The row offset into the parent matrix.

Definition at line 146 of file dense_submatrix.h.

Referenced by libMesh::DenseSubMatrix< T >::i_off().

◆ _j_off

template<typename T>
unsigned int libMesh::DenseSubMatrix< T >::_j_off
private

The column offset into the parent matrix.

Definition at line 151 of file dense_submatrix.h.

Referenced by libMesh::DenseSubMatrix< T >::j_off().

◆ _m

template<typename T>
unsigned int libMesh::DenseMatrixBase< T >::_m
protectedinherited

The row dimension.

Definition at line 176 of file dense_matrix_base.h.

Referenced by libMesh::DenseMatrixBase< T >::m().

◆ _n

template<typename T>
unsigned int libMesh::DenseMatrixBase< T >::_n
protectedinherited

The column dimension.

Definition at line 181 of file dense_matrix_base.h.

Referenced by libMesh::DenseMatrixBase< T >::n().

◆ _parent_matrix

template<typename T>
DenseMatrix<T>& libMesh::DenseSubMatrix< T >::_parent_matrix
private

The parent matrix that contains this submatrix.

Definition at line 141 of file dense_submatrix.h.

Referenced by libMesh::DenseSubMatrix< T >::parent().


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