v0.13.0
MatrixManager.hpp
Go to the documentation of this file.
1 /** \file MatrixManager.hpp
2  * \brief Interface for creating matrices and managing matrices
3  * \ingroup mofem_mat_interface
4  *
5  * Creating and managing matrices
6  *
7  */
8 
9 /*
10  * MoFEM is distributed in the hope that it will be useful, but WITHOUT
11  * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
12  * FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
13  * License for more details.
14  *
15  * You should have received a copy of the GNU Lesser General Public
16  * License along with MoFEM. If not, see <http://www.gnu.org/licenses/>
17  */
18 
19 #ifndef __MATMANAGER_HPP__
20 #define __MATMANAGER_HPP__
21 
22 #include "UnknownInterface.hpp"
23 
24 namespace MoFEM {
25 
26 /**
27  * \brief Matrix manager is used to build and partition problems
28  * \ingroup mofem_mat_interface
29  *
30  */
32 
33  MoFEMErrorCode query_interface(boost::typeindex::type_index type_index,
34  UnknownInterface **iface) const;
35 
37  MatrixManager(const MoFEM::Core &core);
38 
39  /**
40  * \brief Destructor
41  */
42  virtual ~MatrixManager() = default;
43 
44  /**
45  * @brief Creates a MPI AIJ matrix using arrays that contain in standard CSR
46  * format the local rows.
47  *
48  * <a
49  href=https://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/Mat/MatCreateMPIAIJWithArrays.html>See
50  PETSc for details</a>
51 
52  * @tparam Tag
53  * @param name
54  * @param Aij
55  * @param verb
56  * @return MoFEMErrorCode
57  */
58  template <class Tag>
59  MoFEMErrorCode createMPIAIJWithArrays(const std::string name, Mat *Aij,
60  int verb = QUIET) {
61  static_assert(!std::is_same<Tag, Tag>::value, "not implemented");
62  return 0;
63  }
64 
65  /** \copydoc MoFEM::MatrixManager::createMPIAIJWithArrays
66  */
67  template <class Tag>
68  MoFEMErrorCode createMPIAIJWithArrays(const std::string name,
69  SmartPetscObj<Mat> &aij_ptr,
70  int verb = QUIET) {
72  Mat aij;
73  CHKERR createMPIAIJWithArrays<Tag>(name, &aij, verb);
74  aij_ptr.reset(aij, false);
76  }
77 
78  /**
79  * @brief Creates a MPI AIJ matrix using arrays that contain in standard CSR
80  * format the local rows.
81  *
82  * <a
83  href=https://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/Mat/MATMPIAIJ.html>See
84  PETSc for details</a>
85 
86  * @tparam Tag
87  * @param name
88  * @param Aij
89  * @param verb
90  * @return MoFEMErrorCode
91  */
92  template <class Tag>
93  MoFEMErrorCode createMPIAIJ(const std::string name, Mat *Aij,
94  int verb = QUIET) {
95  static_assert(!std::is_same<Tag, Tag>::value, "not implemented");
96  return 0;
97  }
98 
99  /** \copydoc MoFEM::MatrixManager::createMPIAIJ
100  */
101  template <class Tag>
102  MoFEMErrorCode createMPIAIJ(const std::string name,
103  SmartPetscObj<Mat> &aij_ptr,
104  int verb = QUIET) {
106  Mat aij;
107  CHKERR createMPIAIJ<Tag>(name, &aij, verb);
108  aij_ptr.reset(aij, false);
110  }
111 
112  /**
113  * @brief Creates a sparse matrix representing an adjacency list.
114  *
115  * The matrix
116  * does not have numerical values associated with it, but is intended for
117  * ordering (to reduce bandwidth etc) and partitioning.
118  *
119  * <a
120  * href=https://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/Mat/MatCreateMPIAdj.html>See
121  * PETSc for details</a>
122  *
123  * \note This matrix object does not support most matrix operations, include
124  * MatSetValues(). You must NOT free the ii, values and jj arrays yourself.
125  * PETSc will free them when the matrix is destroyed; you must allocate them
126  * with PetscMalloc(). If you call from Fortran you need not create the arrays
127  * with PetscMalloc(). Should not include the matrix diagonals.
128  *
129  * @tparam Tag
130  * @param name
131  * @param Adj
132  * @param verb
133  * @return MoFEMErrorCode
134  */
135  template <class Tag>
136  MoFEMErrorCode createMPIAdjWithArrays(const std::string name, Mat *Adj,
137  int verb = QUIET) {
138  static_assert(!std::is_same<Tag, Tag>::value, "not implemented");
139  return MOFEM_NOT_IMPLEMENTED;
140  }
141 
142  /**
143  * @brief Create sequencial matrix
144  *
145  * Creates a sparse matrix in AIJ (compressed row) format (the default
146  * parallel PETSc format). For good matrix assembly performance the user
147  * should preallocate the matrix storage by setting the parameter nz (or the
148  * array nnz). By setting these parameters accurately, performance during
149  * matrix assembly can be increased by more than a factor of 50.
150  *
151  * <a
152  * href=https://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/Mat/MatCreateSeqAIJ.html>See
153  * PETSc for details</a>
154  *
155  * @tparam Tag
156  * @param name
157  * @param Aij
158  * @param i
159  * @param j
160  * @param v
161  * @param verb
162  * @return MoFEMErrorCode
163  */
164  template <class Tag>
165  MoFEMErrorCode createSeqAIJWithArrays(const std::string name, Mat *Aij,
166  int verb = QUIET) {
167  static_assert(!std::is_same<Tag, Tag>::value, "not implemented");
168  return 0;
169  }
170 
171  /**
172  * @brief Create sequencial matrix
173  *
174  * Creates a sparse matrix in AIJ (compressed row) format (the default
175  * parallel PETSc format). For good matrix assembly performance the user
176  * should preallocate the matrix storage by setting the parameter nz (or the
177  * array nnz). By setting these parameters accurately, performance during
178  * matrix assembly can be increased by more than a factor of 50.
179  *
180  * <a
181  * href=https://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/Mat/MatCreateSeqAIJ.html>See
182  * PETSc for details</a>
183  *
184  * @tparam Tag
185  * @param name
186  * @param Aij (SmartPetscObj)
187  * @param i
188  * @param j
189  * @param v
190  * @param verb
191  * @return MoFEMErrorCode
192  */
193  template <class Tag>
194  MoFEMErrorCode createSeqAIJWithArrays(const std::string name,
195  SmartPetscObj<Mat> &aij_ptr,
196  int verb = QUIET) {
198  Mat aij;
199  CHKERR createSeqAIJWithArrays<Tag>(name, &aij, verb);
200  aij_ptr.reset(aij, false);
202  }
203 
204  /** \brief check if matrix fill in correspond to finite element indices
205 
206  This is used to check consistency of code. If problem is notices with
207  additional non-zero elements in matrix, this function can help detect
208  problem. Should be used as a part of atom tests
209 
210  * @param problem_name
211  * @param row print info at particular row
212  * @param col print info at particular col
213  * @return MoFEMErrorCode
214 
215  */
216  MoFEMErrorCode checkMatrixFillIn(const std::string problem_name,
217  int row_print, int col_print, Mat A,
218  int verb = QUIET);
219 
220  /** \brief check if matrix fill in correspond to finite element indices
221 
222  This is used to check consistency of code. If problem is notices with
223  additional non-zero elements in matrix, this function can help detect
224  problem. Should be used as a part of atom tests
225 
226  * @tparam Tag
227  * @param problem_name
228  * @param row print info at particular row
229  * @param col print info at particular col
230  * @return MoFEMErrorCode
231 
232  */
233  template <class Tag>
235  checkMPIAIJWithArraysMatrixFillIn(const std::string problem_name,
236  int row_print, int col_print,
237  int verb = QUIET) {
238  static_assert(!std::is_same<Tag, Tag>::value, "not implemented");
239  return 0;
240  }
241 
242  /** \brief check if matrix fill in correspond to finite element indices
243 
244  This is used to check consistency of code. If problem is notices with
245  additional non-zero elements in matrix, this function can help detect
246  problem. Should be used as a part of atom tests
247 
248  * @tparam Tag
249  * @param problem_name
250  * @param row print info at particular row
251  * @param col print info at particular col
252  * @return MoFEMErrorCode
253 
254  */
255  template <class Tag>
256  MoFEMErrorCode checkMPIAIJMatrixFillIn(const std::string problem_name,
257  int row_print, int col_print,
258  int verb = QUIET) {
259  static_assert(!std::is_same<Tag, Tag>::value, "not implemented");
260  return 0;
261  }
262 
263 private:
269 };
270 
271 template <>
272 MoFEMErrorCode MatrixManager::createMPIAIJWithArrays<PetscGlobalIdx_mi_tag>(
273  const std::string name, Mat *Aij, int verb);
274 
275 template <>
277 MatrixManager::createMPIAIJ<PetscGlobalIdx_mi_tag>(const std::string name,
278  Mat *Aij, int verb);
279 
280 template <>
282 MatrixManager::createMPIAdjWithArrays<Idx_mi_tag>(const std::string name,
283  Mat *Adj, int verb);
284 
285 template <>
286 MoFEMErrorCode MatrixManager::createSeqAIJWithArrays<PetscLocalIdx_mi_tag>(
287  const std::string name, Mat *Aij, int verb);
288 
289 template <>
291 MatrixManager::checkMPIAIJWithArraysMatrixFillIn<PetscGlobalIdx_mi_tag>(
292  const std::string problem_name, int row_print, int col_print, int verb);
293 
294 template <>
295 MoFEMErrorCode MatrixManager::checkMPIAIJMatrixFillIn<PetscGlobalIdx_mi_tag>(
296  const std::string problem_name, int row_print, int col_print, int verb);
297 
298 } // namespace MoFEM
299 
300 #endif // __MATMANAGER_HPP__a
301 
302 /**
303  * \defgroup mofem_mat_interface Matrix Manager
304  * \brief Creating and managing matrices
305  *
306  * \ingroup mofem
307  */
MoFEM interface.
@ QUIET
Definition: definitions.h:221
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:359
@ MOFEM_NOT_IMPLEMENTED
Definition: definitions.h:45
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:429
#define CHKERR
Inline error check.
Definition: definitions.h:548
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:67
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:21
double A
Core (interface) class.
Definition: Core.hpp:92
Matrix manager is used to build and partition problems.
MoFEMErrorCode query_interface(boost::typeindex::type_index type_index, UnknownInterface **iface) const
MoFEMErrorCode createMPIAdjWithArrays(const std::string name, Mat *Adj, int verb=QUIET)
Creates a sparse matrix representing an adjacency list.
MoFEMErrorCode createMPIAIJ(const std::string name, Mat *Aij, int verb=QUIET)
Creates a MPI AIJ matrix using arrays that contain in standard CSR format the local rows.
PetscLogEvent MOFEM_EVENT_createMPIAdjWithArrays
virtual ~MatrixManager()=default
Destructor.
PetscLogEvent MOFEM_EVENT_checkMatrixFillIn
PetscLogEvent MOFEM_EVENT_createMPIAIJ
MoFEMErrorCode createMPIAIJWithArrays(const std::string name, Mat *Aij, int verb=QUIET)
Creates a MPI AIJ matrix using arrays that contain in standard CSR format the local rows.
PetscLogEvent MOFEM_EVENT_createMPIAIJWithArrays
MoFEMErrorCode checkMatrixFillIn(const std::string problem_name, int row_print, int col_print, Mat A, int verb=QUIET)
check if matrix fill in correspond to finite element indices
MoFEMErrorCode createMPIAIJ(const std::string name, SmartPetscObj< Mat > &aij_ptr, int verb=QUIET)
Creates a MPI AIJ matrix using arrays that contain in standard CSR format the local rows.
PetscLogEvent MOFEM_EVENT_createSeqAIJWithArrays
MoFEMErrorCode createSeqAIJWithArrays(const std::string name, Mat *Aij, int verb=QUIET)
Create sequencial matrix.
MatrixManager(const MoFEM::Core &core)
MoFEMErrorCode checkMPIAIJMatrixFillIn(const std::string problem_name, int row_print, int col_print, int verb=QUIET)
check if matrix fill in correspond to finite element indices
MoFEMErrorCode createMPIAIJWithArrays(const std::string name, SmartPetscObj< Mat > &aij_ptr, int verb=QUIET)
Creates a MPI AIJ matrix using arrays that contain in standard CSR format the local rows.
MoFEMErrorCode createSeqAIJWithArrays(const std::string name, SmartPetscObj< Mat > &aij_ptr, int verb=QUIET)
Create sequencial matrix.
MoFEMErrorCode checkMPIAIJWithArraysMatrixFillIn(const std::string problem_name, int row_print, int col_print, int verb=QUIET)
check if matrix fill in correspond to finite element indices
base class for all interface classes