libMesh
dense_matrix_base.C
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 // C++ includes
19 #include <iomanip> // for std::setw()
20 
21 // Local Includes
22 #include "libmesh/dense_matrix_base.h"
23 #include "libmesh/dense_vector_base.h"
24 
25 namespace libMesh
26 {
27 
28 
29 
30 template<typename T>
32  const DenseMatrixBase<T> & M2,
33  const DenseMatrixBase<T> & M3)
34 {
35  // Assertions to make sure we have been
36  // passed matrices of the correct dimension.
37  libmesh_assert_equal_to (M1.m(), M2.m());
38  libmesh_assert_equal_to (M1.n(), M3.n());
39  libmesh_assert_equal_to (M2.n(), M3.m());
40 
41  const unsigned int m_s = M2.m();
42  const unsigned int p_s = M2.n();
43  const unsigned int n_s = M1.n();
44 
45  // Do it this way because there is a
46  // decent chance (at least for constraint matrices)
47  // that M3(k,j) = 0. when right-multiplying.
48  for (unsigned int k=0; k<p_s; k++)
49  for (unsigned int j=0; j<n_s; j++)
50  if (M3.el(k,j) != 0.)
51  for (unsigned int i=0; i<m_s; i++)
52  M1.el(i,j) += M2.el(i,k) * M3.el(k,j);
53 }
54 
55 
56 
57 template<typename T>
58 void DenseMatrixBase<T>::condense(const unsigned int iv,
59  const unsigned int jv,
60  const T val,
61  DenseVectorBase<T> & rhs)
62 {
63  libmesh_assert_equal_to (this->_m, rhs.size());
64  libmesh_assert_equal_to (iv, jv);
65 
66 
67  // move the known value into the RHS
68  // and zero the column
69  for (unsigned int i=0; i<this->m(); i++)
70  {
71  rhs.el(i) -= this->el(i,jv)*val;
72  this->el(i,jv) = 0.;
73  }
74 
75  // zero the row
76  for (unsigned int j=0; j<this->n(); j++)
77  this->el(iv,j) = 0.;
78 
79  this->el(iv,jv) = 1.;
80  rhs.el(iv) = val;
81 
82 }
83 
84 
85 template<typename T>
86 void DenseMatrixBase<T>::print_scientific (std::ostream & os, unsigned precision) const
87 {
88 #ifndef LIBMESH_BROKEN_IOSTREAM
89 
90  // save the initial format flags
91  std::ios_base::fmtflags os_flags = os.flags();
92 
93  // Print the matrix entries.
94  for (unsigned int i=0; i<this->m(); i++)
95  {
96  for (unsigned int j=0; j<this->n(); j++)
97  os << std::setw(15)
98  << std::scientific
99  << std::setprecision(precision)
100  << this->el(i,j) << " ";
101 
102  os << std::endl;
103  }
104 
105  // reset the original format flags
106  os.flags(os_flags);
107 
108 #else
109 
110  // Print the matrix entries.
111  for (unsigned int i=0; i<this->m(); i++)
112  {
113  for (unsigned int j=0; j<this->n(); j++)
114  os << std::setprecision(precision)
115  << this->el(i,j)
116  << " ";
117 
118  os << std::endl;
119  }
120 
121 
122 #endif
123 }
124 
125 
126 
127 template<typename T>
128 void DenseMatrixBase<T>::print (std::ostream & os) const
129 {
130  for (unsigned int i=0; i<this->m(); i++)
131  {
132  for (unsigned int j=0; j<this->n(); j++)
133  os << std::setw(8)
134  << this->el(i,j) << " ";
135 
136  os << std::endl;
137  }
138 
139  return;
140 }
141 
142 
143 
144 
145 
146 
147 
148 
149 
150 //--------------------------------------------------------------
151 // Explicit instantiations
152 template class DenseMatrixBase<Real>;
153 
154 #ifdef LIBMESH_USE_COMPLEX_NUMBERS
155 template class DenseMatrixBase<Complex>;
156 #endif
157 
158 } // namespace libMesh
void print_scientific(std::ostream &os, unsigned precision=8) const
Prints the matrix entries with more decimal places in scientific notation.
unsigned int n() const
The libMesh namespace provides an interface to certain functionality in the library.
unsigned int m() const
virtual T el(const unsigned int i) const =0
virtual unsigned int size() const =0
virtual T el(const unsigned int i, const unsigned int j) const =0
Defines an abstract dense vector base class for use in Finite Element-type computations.
Definition: dof_map.h: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)...
Defines an abstract dense matrix base class for use in Finite Element-type computations.
void print(std::ostream &os=libMesh::out) const
Pretty-print the matrix, by default to libMesh::out.
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.