v0.14.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 #ifndef __MATMANAGER_HPP__
10 #define __MATMANAGER_HPP__
11 
12 #include "UnknownInterface.hpp"
13 
14 namespace MoFEM {
15 
16 /**
17  * \brief Matrix manager is used to build and partition problems
18  * \ingroup mofem_mat_interface
19  *
20  */
22 
23  MoFEMErrorCode query_interface(boost::typeindex::type_index type_index,
24  UnknownInterface **iface) const;
25 
27  MatrixManager(const MoFEM::Core &core);
28 
29  /**
30  * \brief Destructor
31  */
32  virtual ~MatrixManager() = default;
33 
34  /**
35  * @brief Creates a MPI AIJ matrix using arrays that contain in standard CSR
36  * format the local rows.
37  *
38  * <a
39  href=https://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/Mat/MatCreateMPIAIJWithArrays.html>See
40  PETSc for details</a>
41 
42  * @tparam Tag
43  * @param name
44  * @param Aij
45  * @param verb
46  * @return MoFEMErrorCode
47  */
48  template <class Tag>
49  MoFEMErrorCode createMPIAIJWithArrays(const std::string name, Mat *Aij,
50  int verb = QUIET) {
51  static_assert(!std::is_same<Tag, Tag>::value, "not implemented");
52  return 0;
53  }
54 
55  /** \copydoc MoFEM::MatrixManager::createMPIAIJWithArrays
56  */
57  template <class Tag>
58  MoFEMErrorCode createMPIAIJWithArrays(const std::string name,
59  SmartPetscObj<Mat> &aij_ptr,
60  int verb = QUIET) {
62  Mat aij;
63  CHKERR createMPIAIJWithArrays<Tag>(name, &aij, verb);
64  aij_ptr.reset(aij, false);
66  }
67 
68  template <class Tag>
70  Mat *Aij, int verb = QUIET) {
71  static_assert(!std::is_same<Tag, Tag>::value, "not implemented");
72  return 0;
73  }
74 
75  /** \copydoc MoFEM::MatrixManager::createMPIAIJCUSPARSEWithArrays
76  */
77  template <class Tag>
79  SmartPetscObj<Mat> &aij_ptr,
80  int verb = QUIET) {
82  Mat aij;
83  CHKERR createMPIAIJCUSPARSEWithArrays<Tag>(name, &aij, verb);
84  aij_ptr.reset(aij, false);
86  }
87 
88  template <class Tag>
90  Mat *Aij, int verb = QUIET) {
91  static_assert(!std::is_same<Tag, Tag>::value, "not implemented");
92  return 0;
93  }
94 
95  /** \copydoc MoFEM::MatrixManager::createSeqAIJCUSPARSEWithArrays
96  */
97  template <class Tag>
99  SmartPetscObj<Mat> &aij_ptr,
100  int verb = QUIET) {
102  Mat aij;
103  CHKERR createSeqAIJCUSPARSEWithArrays<Tag>(name, &aij, verb);
104  aij_ptr.reset(aij, false);
106  }
107 
108  /**
109  * @brief Creates a MPI AIJ matrix using arrays that contain in standard CSR
110  * format the local rows.
111  *
112  * <a
113  href=https://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/Mat/MATMPIAIJ.html>See
114  PETSc for details</a>
115 
116  * @tparam Tag
117  * @param name
118  * @param Aij
119  * @param verb
120  * @return MoFEMErrorCode
121  */
122  template <class Tag>
123  MoFEMErrorCode createMPIAIJ(const std::string name, Mat *Aij,
124  int verb = QUIET) {
125  static_assert(!std::is_same<Tag, Tag>::value, "not implemented");
126  return 0;
127  }
128 
129  /** \copydoc MoFEM::MatrixManager::createMPIAIJ
130  */
131  template <class Tag>
132  MoFEMErrorCode createMPIAIJ(const std::string name,
133  SmartPetscObj<Mat> &aij_ptr,
134  int verb = QUIET) {
136  Mat aij;
137  CHKERR createMPIAIJ<Tag>(name, &aij, verb);
138  aij_ptr.reset(aij, false);
140  }
141 
142  /**
143  * @brief Creates a sparse matrix representing an adjacency list.
144  *
145  * The matrix
146  * does not have numerical values associated with it, but is intended for
147  * ordering (to reduce bandwidth etc) and partitioning.
148  *
149  * <a
150  * href=https://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/Mat/MatCreateMPIAdj.html>See
151  * PETSc for details</a>
152  *
153  * \note This matrix object does not support most matrix operations, include
154  * MatSetValues(). You must NOT free the ii, values and jj arrays yourself.
155  * PETSc will free them when the matrix is destroyed; you must allocate them
156  * with PetscMalloc(). If you call from Fortran you need not create the arrays
157  * with PetscMalloc(). Should not include the matrix diagonals.
158  *
159  * @tparam Tag
160  * @param name
161  * @param Adj
162  * @param verb
163  * @return MoFEMErrorCode
164  */
165  template <class Tag>
166  MoFEMErrorCode createMPIAdjWithArrays(const std::string name, Mat *Adj,
167  int verb = QUIET) {
168  static_assert(!std::is_same<Tag, Tag>::value, "not implemented");
169  return MOFEM_NOT_IMPLEMENTED;
170  }
171 
172  /**
173  * @brief Create sequencial matrix
174  *
175  * Creates a sparse matrix in AIJ (compressed row) format (the default
176  * parallel PETSc format). For good matrix assembly performance the user
177  * should preallocate the matrix storage by setting the parameter nz (or the
178  * array nnz). By setting these parameters accurately, performance during
179  * matrix assembly can be increased by more than a factor of 50.
180  *
181  * <a
182  * href=https://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/Mat/MatCreateSeqAIJ.html>See
183  * PETSc for details</a>
184  *
185  * @tparam Tag
186  * @param name
187  * @param Aij
188  * @param i
189  * @param j
190  * @param v
191  * @param verb
192  * @return MoFEMErrorCode
193  */
194  template <class Tag>
195  MoFEMErrorCode createSeqAIJWithArrays(const std::string name, Mat *Aij,
196  int verb = QUIET) {
197  static_assert(!std::is_same<Tag, Tag>::value, "not implemented");
198  return 0;
199  }
200 
201  /**
202  * @brief Create sequencial matrix
203  *
204  * Creates a sparse matrix in AIJ (compressed row) format (the default
205  * parallel PETSc format). For good matrix assembly performance the user
206  * should preallocate the matrix storage by setting the parameter nz (or the
207  * array nnz). By setting these parameters accurately, performance during
208  * matrix assembly can be increased by more than a factor of 50.
209  *
210  * <a
211  * href=https://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/Mat/MatCreateSeqAIJ.html>See
212  * PETSc for details</a>
213  *
214  * @tparam Tag
215  * @param name
216  * @param Aij (SmartPetscObj)
217  * @param i
218  * @param j
219  * @param v
220  * @param verb
221  * @return MoFEMErrorCode
222  */
223  template <class Tag>
224  MoFEMErrorCode createSeqAIJWithArrays(const std::string name,
225  SmartPetscObj<Mat> &aij_ptr,
226  int verb = QUIET) {
228  Mat aij;
229  CHKERR createSeqAIJWithArrays<Tag>(name, &aij, verb);
230  aij_ptr.reset(aij, false);
232  }
233 
234  /** \brief check if matrix fill in correspond to finite element indices
235 
236  This is used to check consistency of code. If problem is notices with
237  additional non-zero elements in matrix, this function can help detect
238  problem. Should be used as a part of atom tests
239 
240  * @param problem_name
241  * @param row print info at particular row
242  * @param col print info at particular col
243  * @return MoFEMErrorCode
244 
245  */
246  MoFEMErrorCode checkMatrixFillIn(const std::string problem_name,
247  int row_print, int col_print, Mat A,
248  int verb = QUIET);
249 
250  /** \brief check if matrix fill in correspond to finite element indices
251 
252  This is used to check consistency of code. If problem is notices with
253  additional non-zero elements in matrix, this function can help detect
254  problem. Should be used as a part of atom tests
255 
256  * @tparam Tag
257  * @param problem_name
258  * @param row print info at particular row
259  * @param col print info at particular col
260  * @return MoFEMErrorCode
261 
262  */
263  template <class Tag>
265  checkMPIAIJWithArraysMatrixFillIn(const std::string problem_name,
266  int row_print, int col_print,
267  int verb = QUIET) {
268  static_assert(!std::is_same<Tag, Tag>::value, "not implemented");
269  return 0;
270  }
271 
272  /** \brief check if matrix fill in correspond to finite element indices
273 
274  This is used to check consistency of code. If problem is notices with
275  additional non-zero elements in matrix, this function can help detect
276  problem. Should be used as a part of atom tests
277 
278  * @tparam Tag
279  * @param problem_name
280  * @param row print info at particular row
281  * @param col print info at particular col
282  * @return MoFEMErrorCode
283 
284  */
285  template <class Tag>
286  MoFEMErrorCode checkMPIAIJMatrixFillIn(const std::string problem_name,
287  int row_print, int col_print,
288  int verb = QUIET) {
289  static_assert(!std::is_same<Tag, Tag>::value, "not implemented");
290  return 0;
291  }
292 
293  /**
294  * @brief Create a Hybrid MPIAij object
295  *
296  * It assumes that L2 field is on skeleton, and dim+1 is bridge adjacency.
297  *
298  * @tparam Tag
299  * @param problem_name
300  * @param aij_ptr
301  * @param verb
302  * @return MoFEMErrorCode
303  */
304  template <class Tag>
305  MoFEMErrorCode createHybridL2MPIAIJ(const std::string problem_name,
306  SmartPetscObj<Mat> &aij_ptr,
307  int verb = QUIET) {
308  static_assert(!std::is_same<Tag, Tag>::value, "not implemented");
309  return 0;
310  }
311 
312 private:
320 };
321 
322 template <>
323 MoFEMErrorCode MatrixManager::createMPIAIJWithArrays<PetscGlobalIdx_mi_tag>(
324  const std::string name, Mat *Aij, int verb);
325 
326 template <>
328 MatrixManager::createMPIAIJCUSPARSEWithArrays<PetscGlobalIdx_mi_tag>(
329  const std::string name, Mat *Aij, int verb);
330 
331 template <>
333 MatrixManager::createMPIAIJ<PetscGlobalIdx_mi_tag>(const std::string name,
334  Mat *Aij, int verb);
335 
336 template <>
338 MatrixManager::createMPIAdjWithArrays<Idx_mi_tag>(const std::string name,
339  Mat *Adj, int verb);
340 
341 template <>
342 MoFEMErrorCode MatrixManager::createSeqAIJWithArrays<PetscLocalIdx_mi_tag>(
343  const std::string name, Mat *Aij, int verb);
344 
345 template <>
347 MatrixManager::createSeqAIJCUSPARSEWithArrays<PetscLocalIdx_mi_tag>(
348  const std::string name, Mat *Aij, int verb);
349 
350 template <>
352 MatrixManager::checkMPIAIJWithArraysMatrixFillIn<PetscGlobalIdx_mi_tag>(
353  const std::string problem_name, int row_print, int col_print, int verb);
354 
355 template <>
356 MoFEMErrorCode MatrixManager::checkMPIAIJMatrixFillIn<PetscGlobalIdx_mi_tag>(
357  const std::string problem_name, int row_print, int col_print, int verb);
358 
359 template <>
360 MoFEMErrorCode MatrixManager::createHybridL2MPIAIJ<PetscGlobalIdx_mi_tag>(
361  const std::string problem_name, SmartPetscObj<Mat> &aij_ptr, int verb);
362 
363 } // namespace MoFEM
364 
365 #endif // __MATMANAGER_HPP__a
366 
367 /**
368  * \defgroup mofem_mat_interface Matrix Manager
369  * \brief Creating and managing matrices
370  *
371  * \ingroup mofem
372  */
MoFEM::MatrixManager::MOFEM_EVENT_createMPIAIJCUSPARSEWithArrays
PetscLogEvent MOFEM_EVENT_createMPIAIJCUSPARSEWithArrays
Definition: MatrixManager.hpp:317
MoFEM::CoreTmp< 0 >
Core (interface) class.
Definition: Core.hpp:82
MoFEM::MatrixManager::createHybridL2MPIAIJ
MoFEMErrorCode createHybridL2MPIAIJ(const std::string problem_name, SmartPetscObj< Mat > &aij_ptr, int verb=QUIET)
Create a Hybrid MPIAij object.
Definition: MatrixManager.hpp:305
MoFEM::Exceptions::MoFEMErrorCode
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:56
MoFEM::MatrixManager::createMPIAIJ
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.
Definition: MatrixManager.hpp:132
MoFEM::MatrixManager::~MatrixManager
virtual ~MatrixManager()=default
Destructor.
A
constexpr AssemblyType A
Definition: operators_tests.cpp:30
MoFEM::MatrixManager::checkMPIAIJMatrixFillIn
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
Definition: MatrixManager.hpp:286
MoFEM::MatrixManager::MOFEM_EVENT_createSeqAIJCUSPARSEWithArrays
PetscLogEvent MOFEM_EVENT_createSeqAIJCUSPARSEWithArrays
Definition: MatrixManager.hpp:318
MoFEM::MatrixManager::MOFEM_EVENT_checkMatrixFillIn
PetscLogEvent MOFEM_EVENT_checkMatrixFillIn
Definition: MatrixManager.hpp:319
MoFEM::MatrixManager::query_interface
MoFEMErrorCode query_interface(boost::typeindex::type_index type_index, UnknownInterface **iface) const
Definition: MatrixManager.cpp:600
MoFEM::MatrixManager::createSeqAIJCUSPARSEWithArrays
MoFEMErrorCode createSeqAIJCUSPARSEWithArrays(const std::string name, Mat *Aij, int verb=QUIET)
Definition: MatrixManager.hpp:89
MoFEM::MatrixManager::MOFEM_EVENT_createMPIAIJWithArrays
PetscLogEvent MOFEM_EVENT_createMPIAIJWithArrays
Definition: MatrixManager.hpp:314
CHKERR
#define CHKERR
Inline error check.
Definition: definitions.h:548
MoFEM
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:10
MoFEM::MatrixManager::createMPIAIJCUSPARSEWithArrays
MoFEMErrorCode createMPIAIJCUSPARSEWithArrays(const std::string name, Mat *Aij, int verb=QUIET)
Definition: MatrixManager.hpp:69
MoFEM::MatrixManager::MOFEM_EVENT_createMPIAdjWithArrays
PetscLogEvent MOFEM_EVENT_createMPIAdjWithArrays
Definition: MatrixManager.hpp:315
MoFEM::MatrixManager::createSeqAIJCUSPARSEWithArrays
MoFEMErrorCode createSeqAIJCUSPARSEWithArrays(const std::string name, SmartPetscObj< Mat > &aij_ptr, int verb=QUIET)
Definition: MatrixManager.hpp:98
MoFEM::MatrixManager::MatrixManager
MatrixManager(const MoFEM::Core &core)
Definition: MatrixManager.cpp:606
UnknownInterface.hpp
MoFEM interface.
MoFEM::MatrixManager::MOFEM_EVENT_createMPIAIJ
PetscLogEvent MOFEM_EVENT_createMPIAIJ
Definition: MatrixManager.hpp:313
MoFEM::MatrixManager
Matrix manager is used to build and partition problems.
Definition: MatrixManager.hpp:21
MoFEM::MatrixManager::createMPIAIJWithArrays
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.
Definition: MatrixManager.hpp:49
MoFEM::UnknownInterface
base class for all interface classes
Definition: UnknownInterface.hpp:34
MoFEM::MatrixManager::createMPIAdjWithArrays
MoFEMErrorCode createMPIAdjWithArrays(const std::string name, Mat *Adj, int verb=QUIET)
Creates a sparse matrix representing an adjacency list.
Definition: MatrixManager.hpp:166
MoFEM::MatrixManager::createSeqAIJWithArrays
MoFEMErrorCode createSeqAIJWithArrays(const std::string name, Mat *Aij, int verb=QUIET)
Create sequencial matrix.
Definition: MatrixManager.hpp:195
MoFEM::MatrixManager::createSeqAIJWithArrays
MoFEMErrorCode createSeqAIJWithArrays(const std::string name, SmartPetscObj< Mat > &aij_ptr, int verb=QUIET)
Create sequencial matrix.
Definition: MatrixManager.hpp:224
MoFEM::MatrixManager::checkMPIAIJWithArraysMatrixFillIn
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
Definition: MatrixManager.hpp:265
MoFEM::MatrixManager::createMPIAIJWithArrays
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.
Definition: MatrixManager.hpp:58
MoFEM::MatrixManager::createMPIAIJ
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.
Definition: MatrixManager.hpp:123
QUIET
@ QUIET
Definition: definitions.h:221
MoFEM::SmartPetscObj< Mat >
MoFEMFunctionReturn
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:429
MoFEM::MatrixManager::cOre
MoFEM::Core & cOre
Definition: MatrixManager.hpp:26
MOFEM_NOT_IMPLEMENTED
@ MOFEM_NOT_IMPLEMENTED
Definition: definitions.h:32
MoFEMFunctionBegin
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:359
MoFEM::MatrixManager::checkMatrixFillIn
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
Definition: MatrixManager.cpp:899
MoFEM::MatrixManager::MOFEM_EVENT_createSeqAIJWithArrays
PetscLogEvent MOFEM_EVENT_createSeqAIJWithArrays
Definition: MatrixManager.hpp:316
MoFEM::MatrixManager::createMPIAIJCUSPARSEWithArrays
MoFEMErrorCode createMPIAIJCUSPARSEWithArrays(const std::string name, SmartPetscObj< Mat > &aij_ptr, int verb=QUIET)
Definition: MatrixManager.hpp:78