v0.15.0
Loading...
Searching...
No Matches
Schur.hpp
Go to the documentation of this file.
1/**
2 * @file Schur.hpp
3 * @brief Assemble Schur complement
4 * @date 2023-02-02
5 *
6 * @copyright Copyright (c) 2023
7 *
8 * To create nested system of Schur complements, you push sequence of operator,
9 * to set up data on entities, and then assemble complements.
10 *
11 * \note Try separate floating points operations from book keeping. Also, align
12 * memory that blocks follow floating point operations.
13 *
14 */
15
16#ifndef __SCHUR_HPP__
17#define __SCHUR_HPP__
18
19namespace MoFEM {
20
21constexpr const char MoFEM_BLOCK_MAT[] = "mofem_block_mat";
22
23/**
24 * @brief Structure to register events for Schur block assembly and solver
25 */
27 static PetscLogEvent MOFEM_EVENT_schurMatSetValues;
28 static PetscLogEvent MOFEM_EVENT_opSchurAssembleEnd;
30 static PetscLogEvent MOFEM_EVENT_BlockStructureMult;
31 static PetscLogEvent MOFEM_EVENT_BlockStructureSolve;
32 static PetscLogEvent MOFEM_EVENT_AssembleSchurMat;
33 static PetscLogEvent MOFEM_EVENT_zeroRowsAndCols;
35};
36
37struct SchurElemMats;
40
42
44
45private:
46 using UserDataOperator::UserDataOperator;
47};
48
50
51/**
52 * @brief Construct a new Op Schur Assemble End object
53 *
54 * @param fields_name list of fields (can be empty)
55 * @param field_ents list of entities on which schur complement is applied
56 * @param schur_aos maps dofs indices from main problem to schur complement
57 * @param schur_mat schur matrix
58 * @param sym_schur true if schur (matrix) complement is symmetric
59 * @param symm_op true if block diagonal is symmetric
60 */
62createOpSchurAssembleEnd(std::vector<std::string> fields_name,
63 std::vector<boost::shared_ptr<Range>> field_ents,
66 bool sym_schur = false, bool symm_op = false);
67
68/**
69 * @brief Create a Op Schur Zero Rows And Cols object
70 *
71 * @param marker_ptr
72 * @return OpSchurAssembleBase*
73 */
75 boost::shared_ptr<std::vector<unsigned char>> marker_ptr, double diag_val);
76
77using SchurFieldPair = std::pair<std::string, std::string>;
78
79using SchurFEOpsFEandFields = std::vector<
80
81 std::pair<std::string, std::vector<SchurFieldPair>>
82
83 >;
84
85/**
86 * @brief Create a Mat Diag Blocks object
87 *
88 * @return Mat
89 */
90boost::shared_ptr<BlockStructure> createBlockMatStructure(
91
92 DM dm, //< dm
93 SchurFEOpsFEandFields schur_fe_op_vec //< block elements
94
95);
96
98 std::pair<SmartPetscObj<Mat>, boost::shared_ptr<BlockStructure>>;
99
100/**
101 * @brief Create a Schur Mat object
102 *
103 * @param dm
104 * @param data
105 * @return std::pair<SmartPetscObj<Mat>, boost::shared_ptr<BlockStructure>>
106 */
107SchurShellMatData createBlockMat(DM dm, boost::shared_ptr<BlockStructure> data);
108
109using NestSchurData = std::tuple<
110
111 std::array<SmartPetscObj<Mat>, 4>,
112 std::array<boost::shared_ptr<BlockStructure>, 4>,
113 boost::shared_ptr<BlockStructure>,
114 std::pair<SmartPetscObj<IS>, SmartPetscObj<IS>>
115
116 >;
117
118/**
119 * @brief Get the Schur Nest Mat Array object
120 *
121 * @param dms schur dm, and block dm
122 * @param block mat A data
123 * @param fields_name list of fields
124 * @param field_ents list of entities on which schur complement is applied
125 * @param add_preconditioner_block add block for preconditioner
126 * @return boost::shared_ptr<NestSchurData>
127 */
128boost::shared_ptr<NestSchurData> createSchurNestedMatrixStruture(
129 std::pair<SmartPetscObj<DM>, SmartPetscObj<DM>> dms,
130 boost::shared_ptr<BlockStructure> block_mat_data,
131
132 std::vector<std::string> fields_name, //< a00 fields
133 std::vector<boost::shared_ptr<Range>> field_ents, //< a00 ents
134 bool add_preconditioner_block = false
135
136);
137
138/**
139 * @brief Assemble Schur matrix
140 *
141 * @param m_field
142 * @param B
143 * @param S
144 * @param fields_name
145 * @param field_ents
146 * @param ao
147 * @return MoFEMErrorCode
148 */
150assembleBlockMatSchur(MoFEM::Interface &m_field, Mat B, Mat S,
151 std::vector<std::string> fields_name,
152 std::vector<boost::shared_ptr<Range>> field_ents,
154
155/**
156 * @brief Get the Block Storage object
157 *
158 * @return std::vector<double>&
159 */
160boost::shared_ptr<std::vector<double>> getBlockMatStorageMat(Mat B);
161
162/**
163 * @brief Get the Block Storage object
164 *
165 * @return std::vector<double>&
166 */
167boost::shared_ptr<std::vector<double>>
169
170/**
171 * @brief Switch preconditioner
172 *
173 * @param block_mat_data
174 * @return MoFEMErrorCode
175 */
177schurSwitchPreconditioner(boost::shared_ptr<BlockStructure> block_mat_data);
178
179/**
180 * @brief Save block matrix as a mesh
181 *
182 * @param block_mat_data
183 * @param filename
184 * @return MoFEMErrorCode
185 */
187schurSaveBlockMesh(boost::shared_ptr<BlockStructure> block_mat_data,
188 std::string filename);
189
190/**
191 * @brief Create a Mat Diag Blocks object
192 *
193 * \code {.cpp}
194 *
195 * auto [nested_mat, nested_data_ptr] = createSchurNestedMatrix(
196 *
197 * createSchurNestedMatrixStruture(
198 *
199 * {schur_dm, block_dm}, shell_data,
200 *
201 * {"TENSOR"}, {nullptr}
202 *
203 * )
204 *
205 * );
206 *
207 * \endcode
208 *
209 * @return Mat
210 */
211std::pair<SmartPetscObj<Mat>, boost::shared_ptr<NestSchurData>>
212createSchurNestedMatrix(boost::shared_ptr<NestSchurData> schur_net_data_ptr);
213
214template <>
215MoFEMErrorCode DMMoFEMSetNestSchurData(DM dm, boost::shared_ptr<NestSchurData>);
216
217/**
218 * @brief Set PC for A00 block
219 *
220 * @param pc
221 * @return MoFEMErrorCode
222 */
224
226
229 using MatSetValuesPtr = boost::function<MoFEMErrorCode(
230 Mat M, const EntitiesFieldData::EntData &row_data,
231 const EntitiesFieldData::EntData &col_data, const MatrixDouble &mat,
232 InsertMode iora)>;
233 static MatSetValuesPtr matSetValuesPtr; ///< backend assembly function
234 static MatSetValuesPtr
235 matSetValuesBlockPtr; ///< backend assembly block mat function
236 static MatSetValuesPtr
237 matSetValuesPreconditionedBlockPtr; ///< backend assembly block
238 ///< preconditioner mat function
239};
240
241template <>
244 const EntitiesFieldData::EntData &col_data,
245 const MatrixDouble &mat, InsertMode iora);
246
247template <>
250 const VectorDouble &nf, InsertMode iora) {
251 return VecSetValues<EssentialBcStorage>(V, data, nf, iora);
252}
253
254template <>
256 Mat M, const EntitiesFieldData::EntData &row_data,
257 const EntitiesFieldData::EntData &col_data, const MatrixDouble &mat,
258 InsertMode iora);
259
260template <>
262 Vec V, const EntitiesFieldData::EntData &data, const VectorDouble &nf,
263 InsertMode iora) {
264 return VecSetValues<SchurElemMats>(V, data, nf, iora);
265}
266
267/***
268 * @brief Specialization of MatSetValues for BlockStructure
269 */
270template <>
272MatSetValues<BlockStructure>(Mat M, const EntitiesFieldData::EntData &row_data,
273 const EntitiesFieldData::EntData &col_data,
274 const MatrixDouble &mat, InsertMode iora);
275
276/***
277 * @brief Specialisation of MatSetValues for AssemblyTypeSelector<BLOCK_MAT>
278 */
279template <>
281 Mat M, const EntitiesFieldData::EntData &row_data,
282 const EntitiesFieldData::EntData &col_data, const MatrixDouble &mat,
283 InsertMode iora) {
284 return MatSetValues<BlockStructure>(M, row_data, col_data, mat, iora);
285}
286
287/***
288 * @brief Specialisation of VecSetValues for AssemblyTypeSelector<BLOCK_MAT>
289 */
290template <>
292 Vec V, const EntitiesFieldData::EntData &data, const VectorDouble &nf,
293 InsertMode iora) {
294 return VecSetValues<EssentialBcStorage>(V, data, nf, iora);
295}
296
297/***
298 * @brief Specialization of MatSetValues for SchurElemMatsBlock
299 */
300template <>
303 const EntitiesFieldData::EntData &row_data,
304 const EntitiesFieldData::EntData &col_data,
305 const MatrixDouble &mat, InsertMode iora);
306
307/***
308 * @brief Specialization of VecSetValues for SchurElemMatsBlock
309 */
310template <>
313 const VectorDouble &nf, InsertMode iora) {
314 return VecSetValues<EssentialBcStorage>(V, data, nf, iora);
315}
316
317/***
318 * @brief Specialisation of MatSetValues for AssemblyTypeSelector<BLOCK_SCHUR>
319 */
320template <>
322 Mat M, const EntitiesFieldData::EntData &row_data,
323 const EntitiesFieldData::EntData &col_data, const MatrixDouble &mat,
324 InsertMode iora) {
325 return MatSetValues<SchurElemMatsBlock>(M, row_data, col_data, mat, iora);
326}
327
328/***
329 * @brief Specialisation of VecSetValues for AssemblyTypeSelector<BLOCK_SCHUR>
330 */
331template <>
333 Vec V, const EntitiesFieldData::EntData &data, const VectorDouble &nf,
334 InsertMode iora) {
335 return VecSetValues<EssentialBcStorage>(V, data, nf, iora);
336}
337
338/***
339 * @brief Specialization of MatSetValues for SchurElemMatsPreconditionedBlock
340 */
341template <>
343 Mat M, const EntitiesFieldData::EntData &row_data,
344 const EntitiesFieldData::EntData &col_data, const MatrixDouble &mat,
345 InsertMode iora);
346
347/***
348 * @brief Specialisation of MatSetValues for
349 * AssemblyTypeSelector<BLOCK_PRECONDITIONER_SCHUR>
350 */
351template <>
354 Mat M, const EntitiesFieldData::EntData &row_data,
355 const EntitiesFieldData::EntData &col_data, const MatrixDouble &mat,
356 InsertMode iora) {
357 return MatSetValues<SchurElemMatsPreconditionedBlock>(M, row_data, col_data,
358 mat, iora);
359}
360
361/***
362 * @brief Specialisation of VecSetValues for AssemblyTypeSelector<BLOCK_SCHUR>
363 */
364template <>
367 Vec V, const EntitiesFieldData::EntData &data, const VectorDouble &nf,
368 InsertMode iora) {
369 return ::VecSetValues(V, data.getIndices().size(),
370 &*data.getIndices().begin(), &*nf.begin(), iora);
371}
372
373/**
374 * @deprecated do not use
375 */
376DEPRECATED OpSchurAssembleBase *createOpSchurAssembleEnd(
377 std::vector<std::string> fields_name,
378 std::vector<boost::shared_ptr<Range>> field_ents,
379 std::vector<SmartPetscObj<AO>> sequence_of_aos,
380 std::vector<SmartPetscObj<Mat>> sequence_of_mats,
381 std::vector<bool> sym_schur, bool symm_op,
382 boost::shared_ptr<BlockStructure> diag_blocks = nullptr);
383
384/**
385 * @deprecated do not use
386 */
387DEPRECATED OpSchurAssembleBase *createOpSchurAssembleEnd(
388 std::vector<std::string> fields_name,
389 std::vector<boost::shared_ptr<Range>> field_ents,
390 std::vector<SmartPetscObj<AO>> sequence_of_aos,
391 std::vector<SmartPetscObj<Mat>> sequence_of_mats,
392 std::vector<bool> sym_schur, std::vector<double> diag_eps, bool symm_op,
393 boost::shared_ptr<BlockStructure> diag_blocks = nullptr);
394
395/** @deprecated do not use, use createSchurNestedMatrixStruture instead */
396DEPRECATED inline boost::shared_ptr<NestSchurData>
398 boost::shared_ptr<BlockStructure> block_mat_data,
399
400 std::vector<std::string> fields_name, //< a00 fields
401 std::vector<boost::shared_ptr<Range>> field_ents, //< a00 ents
402 bool add_preconditioner_block = false
403
404) {
405 return createSchurNestedMatrixStruture(dms, block_mat_data, fields_name,
406 field_ents, add_preconditioner_block);
407}
408
409} // namespace MoFEM
410
411#endif //__SCHUR_HPP__
#define DEPRECATED
Definition definitions.h:17
MoFEMErrorCode VecSetValues< EssentialBcStorage >(Vec V, const EntitiesFieldData::EntData &data, const double *ptr, InsertMode iora)
Set values to vector in operator.
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
UBlasMatrix< double > MatrixDouble
Definition Types.hpp:77
UBlasVector< double > VectorDouble
Definition Types.hpp:68
implementation of Data Operators for Forces and Sources
Definition Common.hpp:10
MoFEMErrorCode MatSetValues< SchurElemMats >(Mat M, const EntitiesFieldData::EntData &row_data, const EntitiesFieldData::EntData &col_data, const MatrixDouble &mat, InsertMode iora)
Definition Schur.cpp:2653
OpSchurAssembleBase * createOpSchurAssembleEnd(std::vector< std::string > fields_name, std::vector< boost::shared_ptr< Range > > field_ents, SmartPetscObj< AO > ao, SmartPetscObj< Mat > schur, bool sym_schur, bool symm_op)
Construct a new Op Schur Assemble End object.
Definition Schur.cpp:2585
MoFEMErrorCode MatSetValues< SchurElemMatsBlock >(Mat M, const EntitiesFieldData::EntData &row_data, const EntitiesFieldData::EntData &col_data, const MatrixDouble &mat, InsertMode iora)
Definition Schur.cpp:2729
MoFEMErrorCode schurSwitchPreconditioner(boost::shared_ptr< BlockStructure > block_mat_data)
Switch preconditioner.
Definition Schur.cpp:2769
MoFEMErrorCode setSchurA00MatSolvePC(SmartPetscObj< PC > pc)
Set PC for A00 block.
Definition Schur.cpp:2627
MoFEMErrorCode MatSetValues< SchurElemMatsPreconditionedBlock >(Mat M, const EntitiesFieldData::EntData &row_data, const EntitiesFieldData::EntData &col_data, const MatrixDouble &mat, InsertMode iora)
Definition Schur.cpp:2760
std::pair< SmartPetscObj< Mat >, boost::shared_ptr< BlockStructure > > SchurShellMatData
Definition Schur.hpp:97
std::pair< SmartPetscObj< Mat >, boost::shared_ptr< NestSchurData > > createSchurNestedMatrix(boost::shared_ptr< NestSchurData > schur_net_data_ptr)
Create a Mat Diag Blocks object.
Definition Schur.cpp:2552
MoFEMErrorCode MatSetValues< BlockStructure >(Mat M, const EntitiesFieldData::EntData &row_data, const EntitiesFieldData::EntData &col_data, const MatrixDouble &mat, InsertMode iora)
Definition Schur.cpp:2697
std::pair< std::string, std::string > SchurFieldPair
Definition Schur.hpp:77
MoFEMErrorCode setSchurMatSolvePC(SmartPetscObj< PC > pc)
Definition Schur.cpp:2645
boost::shared_ptr< std::vector< double > > getBlockMatPrconditionerStorageMat(Mat B)
Get the Block Storage object.
Definition Schur.cpp:2139
DEPRECATED boost::shared_ptr< NestSchurData > getNestSchurData(std::pair< SmartPetscObj< DM >, SmartPetscObj< DM > > dms, boost::shared_ptr< BlockStructure > block_mat_data, std::vector< std::string > fields_name, std::vector< boost::shared_ptr< Range > > field_ents, bool add_preconditioner_block=false)
Definition Schur.hpp:397
boost::shared_ptr< BlockStructure > createBlockMatStructure(DM dm, SchurFEOpsFEandFields schur_fe_op_vec)
Create a Mat Diag Blocks object.
Definition Schur.cpp:1082
boost::shared_ptr< NestSchurData > createSchurNestedMatrixStruture(std::pair< SmartPetscObj< DM >, SmartPetscObj< DM > > dms, boost::shared_ptr< BlockStructure > block_mat_data_ptr, std::vector< std::string > fields_names, std::vector< boost::shared_ptr< Range > > field_ents, bool add_preconditioner_block)
Get the Schur Nest Mat Array object.
Definition Schur.cpp:2343
MoFEMErrorCode MatSetValues(Mat M, const EntitiesFieldData::EntData &row_data, const EntitiesFieldData::EntData &col_data, const double *ptr, InsertMode iora)
Assemble PETSc matrix.
OpSchurAssembleBase * createOpSchurZeroRowsAndCols(boost::shared_ptr< std::vector< unsigned char > > marker_ptr, double diag_val)
Create a Op Schur Zero Rows And Cols object.
Definition Schur.cpp:2597
MoFEMErrorCode schurSaveBlockMesh(boost::shared_ptr< BlockStructure > block_mat_data, std::string filename)
Save block matrix as a mesh.
Definition Schur.cpp:2784
constexpr const char MoFEM_BLOCK_MAT[]
Definition Schur.hpp:21
std::vector< std::pair< std::string, std::vector< SchurFieldPair > > > SchurFEOpsFEandFields
Definition Schur.hpp:79
MoFEMErrorCode VecSetValues< SchurElemMatsBlock >(Vec V, const EntitiesFieldData::EntData &data, const VectorDouble &nf, InsertMode iora)
Definition Schur.hpp:312
MoFEMErrorCode DMMoFEMSetNestSchurData(DM dm, boost::shared_ptr< NestSchurData >)
Definition DMMoFEM.cpp:1554
MoFEMErrorCode VecSetValues< SchurElemMats >(Vec V, const EntitiesFieldData::EntData &data, const VectorDouble &nf, InsertMode iora)
Definition Schur.hpp:249
MoFEMErrorCode VecSetValues(Vec V, const EntitiesFieldData::EntData &data, const double *ptr, InsertMode iora)
Assemble PETSc vector.
SchurShellMatData createBlockMat(DM dm, boost::shared_ptr< BlockStructure > data)
Create a Schur Mat object.
Definition Schur.cpp:1452
MoFEMErrorCode assembleBlockMatSchur(MoFEM::Interface &m_field, Mat B, Mat S, std::vector< std::string > fields_name, std::vector< boost::shared_ptr< Range > > field_ents, SmartPetscObj< AO > ao)
Assemble Schur matrix.
Definition Schur.cpp:1817
OpSchurAssembleBase * createOpSchurAssembleBegin()
Definition Schur.cpp:2580
std::tuple< std::array< SmartPetscObj< Mat >, 4 >, std::array< boost::shared_ptr< BlockStructure >, 4 >, boost::shared_ptr< BlockStructure >, std::pair< SmartPetscObj< IS >, SmartPetscObj< IS > > > NestSchurData
Definition Schur.hpp:109
boost::shared_ptr< std::vector< double > > getBlockMatStorageMat(Mat B)
Get the Block Storage object.
Definition Schur.cpp:2132
Deprecated interface functions.
Data on single entity (This is passed as argument to DataOperator::doWork)
static MatSetValuesPtr matSetValuesPreconditionedBlockPtr
Definition Schur.hpp:237
static MatSetValuesPtr matSetValuesPtr
backend assembly function
Definition Schur.hpp:233
boost::function< MoFEMErrorCode( Mat M, const EntitiesFieldData::EntData &row_data, const EntitiesFieldData::EntData &col_data, const MatrixDouble &mat, InsertMode iora)> MatSetValuesPtr
Definition Schur.hpp:229
static MatSetValuesPtr matSetValuesBlockPtr
backend assembly block mat function
Definition Schur.hpp:235
Schur complement data storage.
Definition Schur.cpp:161
Structure to register events for Schur block assembly and solver.
Definition Schur.hpp:26
static PetscLogEvent MOFEM_EVENT_BlockStructureSetValues
Definition Schur.hpp:29
static PetscLogEvent MOFEM_EVENT_BlockStructureMult
Definition Schur.hpp:30
static PetscLogEvent MOFEM_EVENT_zeroRowsAndCols
Definition Schur.hpp:33
static PetscLogEvent MOFEM_EVENT_opSchurAssembleEnd
Definition Schur.hpp:28
static PetscLogEvent MOFEM_EVENT_BlockStructureSolve
Definition Schur.hpp:31
static PetscLogEvent MOFEM_EVENT_schurMatSetValues
Definition Schur.hpp:27
static PetscLogEvent MOFEM_EVENT_AssembleSchurMat
Definition Schur.hpp:32
intrusive_ptr for managing petsc objects