Line data Source code
1 : // The libMesh Finite Element Library.
2 : // Copyright (C) 2002-2025 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_SUBMATRIX_H
21 : #define LIBMESH_DENSE_SUBMATRIX_H
22 :
23 : // Local Includes
24 : #include "libmesh/libmesh_common.h"
25 : #include "libmesh/dense_matrix.h"
26 : #include "libmesh/dense_subvector.h"
27 : #include "libmesh/int_range.h"
28 :
29 : // C++ includes
30 :
31 : namespace libMesh
32 : {
33 :
34 : /**
35 : * Defines a dense submatrix for use in Finite Element-type
36 : * computations. Useful for storing element stiffness matrices before
37 : * summation into a global matrix, particularly when you have systems
38 : * of equations. All overridden virtual functions are documented in
39 : * dense_matrix_base.h.
40 : *
41 : * \author Benjamin S. Kirk
42 : * \date 2003
43 : */
44 : template<typename T>
45 : class DenseSubMatrix : public DenseMatrixBase<T>
46 : {
47 : public:
48 :
49 : /**
50 : * Constructor. Creates a dense submatrix of the matrix
51 : * \p parent. The submatrix has dimensions \f$(m \times n)\f$,
52 : * and the \f$(0,0)\f$ entry of the submatrix is located
53 : * at the \f$(ioff,joff)\f$ location in the parent matrix.
54 : */
55 : DenseSubMatrix(DenseMatrix<T> & new_parent,
56 : const unsigned int ioff=0,
57 : const unsigned int joff=0,
58 : const unsigned int m=0,
59 : const unsigned int n=0);
60 :
61 : /**
62 : * The 5 special functions can be defaulted for this class, as it
63 : * does not manage any memory itself.
64 : */
65 0 : DenseSubMatrix (DenseSubMatrix &&) = default;
66 0 : DenseSubMatrix (const DenseSubMatrix &) = default;
67 : DenseSubMatrix & operator= (const DenseSubMatrix &) = default;
68 : DenseSubMatrix & operator= (DenseSubMatrix &&) = default;
69 1216855 : virtual ~DenseSubMatrix() = default;
70 :
71 : /**
72 : * \returns A reference to the parent matrix.
73 : */
74 0 : DenseMatrix<T> & parent () { return _parent_matrix; }
75 :
76 : virtual void zero() override final;
77 :
78 : /**
79 : * \returns The \p (i,j) element of the submatrix.
80 : */
81 : T operator() (const unsigned int i,
82 : const unsigned int j) const;
83 :
84 : /**
85 : * \returns The \p (i,j) element of the submatrix as a writable reference.
86 : */
87 : T & operator() (const unsigned int i,
88 : const unsigned int j);
89 :
90 0 : virtual T el(const unsigned int i,
91 : const unsigned int j) const override final
92 0 : { return (*this)(i,j); }
93 :
94 0 : virtual T & el(const unsigned int i,
95 : const unsigned int j) override final
96 0 : { return (*this)(i,j); }
97 :
98 : virtual void left_multiply (const DenseMatrixBase<T> & M2) override final;
99 :
100 : virtual void right_multiply (const DenseMatrixBase<T> & M3) override final;
101 :
102 : /**
103 : * Changes the location of the submatrix in the parent matrix.
104 : */
105 : void reposition(const unsigned int ioff,
106 : const unsigned int joff,
107 : const unsigned int new_m,
108 : const unsigned int new_n);
109 :
110 : /**
111 : * \returns The row offset into the parent matrix.
112 : */
113 4128033455 : unsigned int i_off() const { return _i_off; }
114 :
115 : /**
116 : * \returns The column offset into the parent matrix.
117 : */
118 4128033455 : unsigned int j_off() const { return _j_off; }
119 :
120 : /**
121 : * Condense-out the \p (i,j) entry of the matrix, forcing
122 : * it to take on the value \p val. This is useful in numerical
123 : * simulations for applying boundary conditions. Preserves the
124 : * symmetry of the matrix.
125 : */
126 0 : void condense(const unsigned int i,
127 : const unsigned int j,
128 : const T val,
129 : DenseSubVector<T> & rhs)
130 : {
131 0 : this->parent().condense(this->i_off()+i,
132 0 : this->j_off()+j,
133 : val, rhs.parent());
134 0 : }
135 :
136 : private:
137 :
138 : /**
139 : * The parent matrix that contains this submatrix.
140 : */
141 : DenseMatrix<T> & _parent_matrix;
142 :
143 : /**
144 : * The row offset into the parent matrix.
145 : */
146 : unsigned int _i_off;
147 :
148 : /**
149 : * The column offset into the parent matrix.
150 : */
151 : unsigned int _j_off;
152 : };
153 :
154 :
155 : // --------------------------------------------------
156 : // Constructor
157 : template<typename T>
158 : inline
159 1215239 : DenseSubMatrix<T>::DenseSubMatrix(DenseMatrix<T> & new_parent,
160 : const unsigned int ioff,
161 : const unsigned int joff,
162 : const unsigned int new_m,
163 : const unsigned int new_n) :
164 : DenseMatrixBase<T>(new_m,new_n),
165 1215239 : _parent_matrix(new_parent)
166 : {
167 51997 : this->reposition (ioff, joff, new_m, new_n);
168 51997 : }
169 :
170 :
171 : template<typename T>
172 : inline
173 7656972 : void DenseSubMatrix<T>::reposition(const unsigned int ioff,
174 : const unsigned int joff,
175 : const unsigned int new_m,
176 : const unsigned int new_n)
177 : {
178 72876221 : _i_off = ioff;
179 72871421 : _j_off = joff;
180 72823976 : this->_m = new_m;
181 32039584 : this->_n = new_n;
182 :
183 : // Make sure we still fit in the parent matrix.
184 7656972 : libmesh_assert_less_equal ((this->i_off() + this->m()), _parent_matrix.m());
185 7656972 : libmesh_assert_less_equal ((this->j_off() + this->n()), _parent_matrix.n());
186 48387056 : }
187 :
188 :
189 :
190 : template<typename T>
191 : inline
192 0 : void DenseSubMatrix<T>::zero()
193 : {
194 0 : for (auto i : make_range(this->m()))
195 0 : for (auto j : make_range(this->n()))
196 0 : _parent_matrix(i + this->i_off(),
197 0 : j + this->j_off()) = 0.;
198 0 : }
199 :
200 :
201 :
202 : template<typename T>
203 : inline
204 0 : T DenseSubMatrix<T>::operator () (const unsigned int i,
205 : const unsigned int j) const
206 : {
207 0 : libmesh_assert_less (i, this->m());
208 0 : libmesh_assert_less (j, this->n());
209 0 : libmesh_assert_less (i + this->i_off(), _parent_matrix.m());
210 0 : libmesh_assert_less (j + this->j_off(), _parent_matrix.n());
211 :
212 0 : return _parent_matrix (i + this->i_off(),
213 0 : j + this->j_off());
214 : }
215 :
216 :
217 : template<typename T>
218 : inline
219 0 : T & DenseSubMatrix<T>::operator () (const unsigned int i,
220 : const unsigned int j)
221 : {
222 0 : libmesh_assert_less (i, this->m());
223 0 : libmesh_assert_less (j, this->n());
224 0 : libmesh_assert_less (i + this->i_off(), _parent_matrix.m());
225 0 : libmesh_assert_less (j + this->j_off(), _parent_matrix.n());
226 :
227 5367885624 : return _parent_matrix (i + this->i_off(),
228 229049802 : j + this->j_off());
229 : }
230 :
231 :
232 : } // namespace libMesh
233 :
234 :
235 : #endif // LIBMESH_DENSE_SUBMATRIX_H
|
