v0.9.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 
28 
29 /**
30  * \brief Matrix manager is used to build and partition problems
31  * \ingroup mofem_mat_interface
32  *
33  */
35 
37  UnknownInterface **iface) const;
38 
40  MatrixManager(const MoFEM::Core &core);
41 
42  /**
43  * \brief Destructor
44  */
46 
47  /**
48  * @brief Creates a MPI AIJ matrix using arrays that contain in standard CSR
49  * format the local rows.
50  *
51  * <a
52  href=https://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/Mat/MatCreateMPIAIJWithArrays.html>See
53  PETSc for details</a>
54 
55  * @tparam Tag
56  * @param name
57  * @param Aij
58  * @param verb
59  * @return MoFEMErrorCode
60  */
61  template <class Tag>
62  MoFEMErrorCode createMPIAIJWithArrays(const std::string name, Mat *Aij,
63  int verb = QUIET) {
64  static_assert(!std::is_same<Tag, Tag>::value, "not implemented");
65  return 0;
66  }
67 
68  /**
69  * @brief Creates a MPI AIJ matrix using arrays that contain in standard CSR
70  * format the local rows.
71  *
72  * <a
73  href=https://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/Mat/MatCreateMPIAIJWithArrays.html>See
74  PETSc for details</a>
75 
76  * @tparam Tag
77  * @param name
78  * @param Aij (SmartPetscObj)
79  * @param verb
80  * @return MoFEMErrorCode
81  */
82  template <class Tag>
83  MoFEMErrorCode createMPIAIJWithArrays(const std::string name,
84  SmartPetscObj<Mat> &aij_ptr,
85  int verb = QUIET) {
87  Mat aij;
88  CHKERR createMPIAIJWithArrays<Tag>(name, &aij, verb);
89  aij_ptr.reset(aij, false);
91  }
92 
93  /**
94  * @brief Creates a sparse matrix representing an adjacency list.
95  *
96  * The matrix
97  * does not have numerical values associated with it, but is intended for
98  * ordering (to reduce bandwidth etc) and partitioning.
99  *
100  * <a
101  * href=https://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/Mat/MatCreateMPIAdj.html>See
102  * PETSc for details</a>
103  *
104  * \note This matrix object does not support most matrix operations, include
105  * MatSetValues(). You must NOT free the ii, values and jj arrays yourself.
106  * PETSc will free them when the matrix is destroyed; you must allocate them
107  * with PetscMalloc(). If you call from Fortran you need not create the arrays
108  * with PetscMalloc(). Should not include the matrix diagonals.
109  *
110  * @tparam Tag
111  * @param name
112  * @param Adj
113  * @param verb
114  * @return MoFEMErrorCode
115  */
116  template <class Tag>
117  MoFEMErrorCode createMPIAdjWithArrays(const std::string name, Mat *Adj,
118  int verb = QUIET) {
119  static_assert(!std::is_same<Tag, Tag>::value, "not implemented");
120  return MOFEM_NOT_IMPLEMENTED;
121  }
122 
123  /**
124  * @brief Create sequencial matrix
125  *
126  * Creates a sparse matrix in AIJ (compressed row) format (the default
127  * parallel PETSc format). For good matrix assembly performance the user
128  * should preallocate the matrix storage by setting the parameter nz (or the
129  * array nnz). By setting these parameters accurately, performance during
130  * matrix assembly can be increased by more than a factor of 50.
131  *
132  * <a
133  * href=https://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/Mat/MatCreateSeqAIJ.html>See
134  * PETSc for details</a>
135  *
136  * @tparam Tag
137  * @param name
138  * @param Aij
139  * @param i
140  * @param j
141  * @param v
142  * @param verb
143  * @return MoFEMErrorCode
144  */
145  template <class Tag>
146  MoFEMErrorCode createSeqAIJWithArrays(const std::string name, Mat *Aij,
147  int verb = QUIET) {
148  static_assert(!std::is_same<Tag, Tag>::value, "not implemented");
149  return 0;
150  }
151 
152  /**
153  * @brief Create sequencial matrix
154  *
155  * Creates a sparse matrix in AIJ (compressed row) format (the default
156  * parallel PETSc format). For good matrix assembly performance the user
157  * should preallocate the matrix storage by setting the parameter nz (or the
158  * array nnz). By setting these parameters accurately, performance during
159  * matrix assembly can be increased by more than a factor of 50.
160  *
161  * <a
162  * href=https://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/Mat/MatCreateSeqAIJ.html>See
163  * PETSc for details</a>
164  *
165  * @tparam Tag
166  * @param name
167  * @param Aij (SmartPetscObj)
168  * @param i
169  * @param j
170  * @param v
171  * @param verb
172  * @return MoFEMErrorCode
173  */
174  template <class Tag>
175  MoFEMErrorCode createSeqAIJWithArrays(const std::string name,
176  SmartPetscObj<Mat> &aij_ptr,
177  int verb = QUIET) {
179  Mat aij;
180  CHKERR createSeqAIJWithArrays<Tag>(name, &aij, verb);
181  aij_ptr.reset(aij, false);
183  }
184 
185  /** \brief check if matrix fill in correspond to finite element indices
186 
187  This is used to check consistency of code. If problem is notices with
188  additional non-zero elements in matrix, this function can help detect
189  problem. Should be used as a part of atom tests
190 
191  * @tparam Tag
192  * @param problem_name
193  * @param row print info at particular row
194  * @param col print info at particular col
195  * @return MoFEMErrorCode
196 
197  */
198  template <class Tag>
200  checkMPIAIJWithArraysMatrixFillIn(const std::string problem_name,
201  int row_print, int col_print,
202  int verb = QUIET) {
203  static_assert(!std::is_same<Tag, Tag>::value, "not implemented");
204  return 0;
205  }
206 
207 private:
212 };
213 
214 template <>
215 MoFEMErrorCode MatrixManager::createMPIAIJWithArrays<PetscGlobalIdx_mi_tag>(
216  const std::string name, Mat *Aij, int verb);
217 
218 template <>
220 MatrixManager::createMPIAdjWithArrays<Idx_mi_tag>(const std::string name,
221  Mat *Adj, int verb);
222 
223 template <>
224 MoFEMErrorCode MatrixManager::createSeqAIJWithArrays<PetscLocalIdx_mi_tag>(
225  const std::string name, Mat *Aij, int verb);
226 
227 template <>
229 MatrixManager::checkMPIAIJWithArraysMatrixFillIn<PetscGlobalIdx_mi_tag>(
230  const std::string problem_name, int row_print, int col_print, int verb);
231 
232 } // namespace MoFEM
233 
234 #endif // __MATMANAGER_HPP__a
235 
236 /**
237  * \defgroup mofem_mat_interface Matrix Manager
238  * \brief Creating and managing matrices
239  *
240  * \ingroup mofem
241  */
MoFEM interface unique ID.
base class for all interface classes
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:477
PetscLogEvent MOFEM_EVENT_checkMPIAIJWithArraysMatrixFillIn
Core (interface) class.
Definition: Core.hpp:50
MoFEM interface.
MoFEMErrorCode createMPIAdjWithArrays(const std::string name, Mat *Adj, int verb=QUIET)
Creates a sparse matrix representing an adjacency list.
PetscLogEvent MOFEM_EVENT_createSeqAIJWithArrays
std::bitset< BITINTERFACEUID_SIZE > BitIntefaceId
Definition: Types.hpp:56
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:21
~MatrixManager()
Destructor.
MoFEMErrorCode createSeqAIJWithArrays(const std::string name, Mat *Aij, int verb=QUIET)
Create sequencial matrix.
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:66
PetscLogEvent MOFEM_EVENT_createMPIAdjWithArrays
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
MatrixManager(const MoFEM::Core &core)
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.
#define CHKERR
Inline error check.
Definition: definitions.h:596
MoFEMErrorCode createSeqAIJWithArrays(const std::string name, SmartPetscObj< Mat > &aij_ptr, int verb=QUIET)
Create sequencial matrix.
static const MOFEMuuid IDD_MOFEMMatrixManager
PetscLogEvent MOFEM_EVENT_createMPIAIJWithArrays
MoFEMErrorCode query_interface(const MOFEMuuid &uuid, UnknownInterface **iface) const
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:407
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.
Matrix manager is used to build and partition problems.