v0.8.23
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 href=https://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/Mat/MatCreateMPIAIJWithArrays.html>See PETSc for details</a>
52 
53  * @tparam Tag
54  * @param name
55  * @param Aij
56  * @param verb
57  * @return MoFEMErrorCode
58  */
59  template <class Tag>
60  MoFEMErrorCode createMPIAIJWithArrays(const std::string name, Mat *Aij,
61  int verb = QUIET) {
62  static_assert(!std::is_same<Tag, Tag>::value, "not implemented");
63  return 0;
64  }
65 
66  /**
67  * @brief Creates a sparse matrix representing an adjacency list.
68  *
69  * The matrix
70  * does not have numerical values associated with it, but is intended for
71  * ordering (to reduce bandwidth etc) and partitioning.
72  *
73  * <a href=https://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/Mat/MatCreateMPIAdj.html>See PETSc for details</a>
74  *
75  * \note This matrix object does not support most matrix operations, include
76  * MatSetValues(). You must NOT free the ii, values and jj arrays yourself.
77  * PETSc will free them when the matrix is destroyed; you must allocate them
78  * with PetscMalloc(). If you call from Fortran you need not create the arrays
79  * with PetscMalloc(). Should not include the matrix diagonals.
80  *
81  * @tparam Tag
82  * @param name
83  * @param Adj
84  * @param verb
85  * @return MoFEMErrorCode
86  */
87  template <class Tag>
88  MoFEMErrorCode createMPIAdjWithArrays(const std::string name, Mat *Adj,
89  int verb = QUIET) {
90  static_assert(!std::is_same<Tag, Tag>::value, "not implemented");
91  }
92 
93  /**
94  * @brief Create sequencial matrix
95  *
96  * Creates a sparse matrix in AIJ (compressed row) format (the default
97  * parallel PETSc format). For good matrix assembly performance the user
98  * should preallocate the matrix storage by setting the parameter nz (or the
99  * array nnz). By setting these parameters accurately, performance during
100  * matrix assembly can be increased by more than a factor of 50.
101  *
102  * <a href=https://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/Mat/MatCreateSeqAIJ.html>See PETSc for details</a>
103  *
104  * @tparam Tag
105  * @param name
106  * @param Aij
107  * @param i
108  * @param j
109  * @param v
110  * @param verb
111  * @return MoFEMErrorCode
112  */
113  template <class Tag>
114  MoFEMErrorCode createSeqAIJWithArrays(const std::string name, Mat *Aij,
115  int verb = QUIET) {
116  static_assert(!std::is_same<Tag, Tag>::value, "not implemented");
117  return 0;
118  }
119 
120  /** \brief check if matrix fill in correspond to finite element indices
121 
122  This is used to check consistency of code. If problem is notices with
123  additional non-zero elements in matrix, this function can help detect
124  problem. Should be used as a part of atom tests
125 
126  * @tparam Tag
127  * @param problem_name
128  * @param row print info at particular row
129  * @param col print info at particular col
130  * @return MoFEMErrorCode
131 
132  */
133  template <class Tag>
135  checkMPIAIJWithArraysMatrixFillIn(const std::string problem_name,
136  int row_print, int col_print,
137  int verb = QUIET) {
138  static_assert(!std::is_same<Tag, Tag>::value, "not implemented");
139  return 0;
140  }
141 
142 private:
147 };
148 
149 template <>
150 MoFEMErrorCode MatrixManager::createMPIAIJWithArrays<PetscGlobalIdx_mi_tag>(
151  const std::string name, Mat *Aij, int verb);
152 
153 template <>
155 MatrixManager::createMPIAdjWithArrays<Idx_mi_tag>(const std::string name,
156  Mat *Adj, int verb);
157 
158 template <>
159 MoFEMErrorCode MatrixManager::createSeqAIJWithArrays<PetscLocalIdx_mi_tag>(
160  const std::string name, Mat *Aij, int verb);
161 
162 template <>
164 MatrixManager::checkMPIAIJWithArraysMatrixFillIn<PetscGlobalIdx_mi_tag>(
165  const std::string problem_name, int row_print, int col_print, int verb);
166 
167 } // namespace MoFEM
168 
169 #endif // __MATMANAGER_HPP__a
170 
171 /**
172  * \defgroup mofem_mat_interface Matrix Manager
173  * \brief Creating and managing matrices
174  *
175  * \ingroup mofem
176  */
MoFEM interface unique ID.
base class for all interface classes
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.
static const MOFEMuuid IDD_MOFEMMatrixManager
PetscLogEvent MOFEM_EVENT_createMPIAIJWithArrays
MoFEMErrorCode query_interface(const MOFEMuuid &uuid, UnknownInterface **iface) const
Matrix manager is used to build and partition problems.