v0.13.2
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 */
12
13#ifndef __SCHUR_HPP__
14#define __SCHUR_HPP__
15
16namespace MoFEM {
17
18/**
19 * @brief Clear Schur complement internal data
20 *
21 */
23
26
27 MoFEMErrorCode doWork(int side, EntityType type,
29};
30
31/**
32 * @brief Assemble Schur complement (Implementation)
33 *
34 */
36
37 using MatSetValuesRaw = boost::function<MoFEMErrorCode(
38 Mat mat, PetscInt m, const PetscInt idxm[], PetscInt n,
39 const PetscInt idxn[], const PetscScalar v[], InsertMode addv)>;
40
43
44 /**
45 * @brief Construct a new Op Schur Assemble End object
46 *
47 * @param fields_name list of fields
48 * @param field_ents list of entities on which schur complement is applied (can be empty)
49 * @param sequence_of_aos list of maps from base problem to Schur complement matrix
50 * @param sequence_of_mats list of Schur complement matrices
51 * @param sym_schur true if Schur complement is symmetric
52 * @param symm_op true if block diagonal is symmetric
53 */
54 OpSchurAssembleEndImpl(std::vector<std::string> fields_name,
55 std::vector<boost::shared_ptr<Range>> field_ents,
56 std::vector<SmartPetscObj<AO>> sequence_of_aos,
57 std::vector<SmartPetscObj<Mat>> sequence_of_mats,
58 std::vector<bool> sym_schur, bool symm_op = true);
59
60 /**
61 * @brief Construct a new Op Schur Assemble End object
62 *
63 * @param fields_name list of fields
64 * @param field_ents list of entities on which schur complement is applied (can be empty)
65 * @param sequence_of_aos list of maps from base problem to Schur complement matrix
66 * @param sequence_of_mats list of Schur complement matrices
67 * @param sym_schur true if Schur complement is symmetric
68 * @param diag_eps add epsilon on diagonal of inverted matrix
69 * @param symm_op true if block diagonal is symmetric
70 */
71 OpSchurAssembleEndImpl(std::vector<std::string> fields_name,
72 std::vector<boost::shared_ptr<Range>> field_ents,
73 std::vector<SmartPetscObj<AO>> sequence_of_aos,
74 std::vector<SmartPetscObj<Mat>> sequence_of_mats,
75 std::vector<bool> sym_schur,
76 std::vector<double> diag_eps, bool symm_op = true);
77
78protected:
79
80 template <typename I>
83
84 std::vector<std::string> fieldsName;
85 std::vector<boost::shared_ptr<Range>> fieldEnts;
86 std::vector<SmartPetscObj<AO>> sequenceOfAOs;
87 std::vector<SmartPetscObj<Mat>> sequenceOfMats;
88 std::vector<bool> symSchur;
89 std::vector<double> diagEps;
90
95};
96
97struct SCHUR_DSYSV; ///< SY symmetric
98struct SCHUR_DGESV; ///< GE general (i.e., nonsymmetric, in some cases
99 ///< rectangular)
100
101/**
102 * @brief Assemble Schur complement
103 *
104 */
105template <typename I> struct OpSchurAssembleEnd;
106
107template <>
110 MoFEMErrorCode doWork(int side, EntityType type,
112};
113
114template <>
117 MoFEMErrorCode doWork(int side, EntityType type,
119};
120
121/**
122 * @brief Schur complement data storage
123 *
124 */
125struct SchurL2Mats : public boost::enable_shared_from_this<SchurL2Mats> {
126
127 SchurL2Mats(const size_t idx, const UId uid_row, const UId uid_col);
128 virtual ~SchurL2Mats() = default;
129
130 const UId uidRow;
131 const UId uidCol;
132
133 inline auto &getMat() const { return locMats[iDX]; }
134 inline auto &getRowInd() const { return rowIndices[iDX]; }
135 inline auto &getColInd() const { return colIndices[iDX]; }
136
137 using MatSetValuesPtr = boost::function<MoFEMErrorCode(
138 Mat M, const EntitiesFieldData::EntData &row_data,
139 const EntitiesFieldData::EntData &col_data, const MatrixDouble &mat,
140 InsertMode iora)>;
141
142 static MatSetValuesPtr matSetValuesPtr; ///< backend assembly function
143
144 static MoFEMErrorCode MatSetValues(Mat M,
145 const EntitiesFieldData::EntData &row_data,
146 const EntitiesFieldData::EntData &col_data,
147 const MatrixDouble &mat, InsertMode iora);
148
149private:
150 const size_t iDX;
151
154
155 struct idx_mi_tag {};
156 struct uid_mi_tag {};
157 struct row_mi_tag {};
158 struct col_mi_tag {};
159
160 using SchurL2Storage = multi_index_container<
162 indexed_by<
163
164 ordered_unique<
165 tag<uid_mi_tag>,
166 composite_key<
168
169 member<SchurL2Mats, const UId, &SchurL2Mats::uidRow>,
170 member<SchurL2Mats, const UId, &SchurL2Mats::uidCol>
171
172 >>,
173
174 ordered_non_unique<tag<row_mi_tag>, member<SchurL2Mats, const UId,
176
177 ordered_non_unique<tag<col_mi_tag>, member<SchurL2Mats, const UId,
179
180 >>;
181
182 static boost::ptr_vector<MatrixDouble> locMats;
183 static boost::ptr_vector<VectorInt> rowIndices;
184 static boost::ptr_vector<VectorInt> colIndices;
186};
187
188template <>
191 const EntitiesFieldData::EntData &col_data,
192 const MatrixDouble &mat, InsertMode iora) {
193 return SchurL2Mats::MatSetValues(M, row_data, col_data, mat, iora);
194}
195
196template <>
199 const VectorDouble &nf, InsertMode iora) {
200 return VecSetValues<EssentialBcStorage>(V, data, nf, iora);
201}
202
203template <>
204inline MoFEMErrorCode MatSetValues<AssemblyTypeSelector<SCHUR>>(
205 Mat M, const EntitiesFieldData::EntData &row_data,
206 const EntitiesFieldData::EntData &col_data, const MatrixDouble &mat,
207 InsertMode iora) {
208 return MatSetValues<SchurL2Mats>(M, row_data, col_data, mat, iora);
209}
210
211template <>
212inline MoFEMErrorCode VecSetValues<AssemblyTypeSelector<SCHUR>>(
213 Vec V, const EntitiesFieldData::EntData &data, const VectorDouble &nf,
214 InsertMode iora) {
215 return VecSetValues<SchurL2Mats>(V, data, nf, iora);
216}
217
218} // namespace MoFEM
219
220#endif //__SCHUR_HPP__
static Index< 'M', 3 > M
ForcesAndSourcesCore::UserDataOperator UserDataOperator
@ NOSPACE
Definition: definitions.h:83
FTensor::Index< 'n', SPACE_DIM > n
FTensor::Index< 'm', SPACE_DIM > m
MoFEMErrorCode VecSetValues< EssentialBcStorage >(Vec V, const EntitiesFieldData::EntData &data, const double *ptr, InsertMode iora)
Set values to vector in operator.
const double v
phase velocity of light in medium (cm/ns)
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:56
uint128_t UId
Unique Id.
Definition: Types.hpp:31
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(Mat M, const EntitiesFieldData::EntData &row_data, const EntitiesFieldData::EntData &col_data, const double *ptr, InsertMode iora)
Assemble PETSc matrix.
MoFEMErrorCode VecSetValues< SchurL2Mats >(Vec V, const EntitiesFieldData::EntData &data, const VectorDouble &nf, InsertMode iora)
Definition: Schur.hpp:198
MoFEMErrorCode MatSetValues< SchurL2Mats >(Mat M, const EntitiesFieldData::EntData &row_data, const EntitiesFieldData::EntData &col_data, const MatrixDouble &mat, InsertMode iora)
Definition: Schur.hpp:190
Data on single entity (This is passed as argument to DataOperator::doWork)
@ OPSPACE
operator do Work is execute on space data
structure to get information form mofem into EntitiesFieldData
Clear Schur complement internal data.
Definition: Schur.hpp:22
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
Operator for linear form, usually to calculate values on right hand side.
Definition: Schur.cpp:116
Assemble Schur complement.
Definition: Schur.hpp:105
Assemble Schur complement (Implementation)
Definition: Schur.hpp:35
MatrixDouble offMatInvDiagOffMat
Definition: Schur.hpp:93
std::vector< boost::shared_ptr< Range > > fieldEnts
Definition: Schur.hpp:85
boost::function< MoFEMErrorCode(Mat mat, PetscInt m, const PetscInt idxm[], PetscInt n, const PetscInt idxn[], const PetscScalar v[], InsertMode addv)> MatSetValuesRaw
Definition: Schur.hpp:39
MatSetValuesRaw matSetValuesSchurRaw
Definition: Schur.hpp:42
MoFEMErrorCode doWorkImpl(int side, EntityType type, EntitiesFieldData::EntData &data)
Definition: Schur.cpp:159
OpSchurAssembleEndImpl(std::vector< std::string > fields_name, std::vector< boost::shared_ptr< Range > > field_ents, std::vector< SmartPetscObj< AO > > sequence_of_aos, std::vector< SmartPetscObj< Mat > > sequence_of_mats, std::vector< bool > sym_schur, bool symm_op=true)
Construct a new Op Schur Assemble End object.
Definition: Schur.cpp:146
MatSetValuesRaw matSetValuesRaw
Definition: Schur.hpp:41
std::vector< bool > symSchur
Definition: Schur.hpp:88
std::vector< std::string > fieldsName
Definition: Schur.hpp:84
MatrixDouble invDiagOffMat
Definition: Schur.hpp:92
std::vector< SmartPetscObj< Mat > > sequenceOfMats
Definition: Schur.hpp:87
std::vector< SmartPetscObj< AO > > sequenceOfAOs
Definition: Schur.hpp:86
MatrixDouble transOffMatInvDiagOffMat
Definition: Schur.hpp:94
std::vector< double > diagEps
Definition: Schur.hpp:89
Schur complement data storage.
Definition: Schur.hpp:125
virtual ~SchurL2Mats()=default
const UId uidCol
Definition: Schur.hpp:131
const UId uidRow
Definition: Schur.hpp:130
static MatSetValuesPtr matSetValuesPtr
backend assembly function
Definition: Schur.hpp:142
static boost::ptr_vector< VectorInt > rowIndices
Definition: Schur.hpp:183
friend OpSchurAssembleEndImpl
Definition: Schur.hpp:153
static SchurL2Storage schurL2Storage
Definition: Schur.hpp:185
multi_index_container< SchurL2Mats, indexed_by< ordered_unique< tag< uid_mi_tag >, composite_key< SchurL2Mats, member< SchurL2Mats, const UId, &SchurL2Mats::uidRow >, member< SchurL2Mats, const UId, &SchurL2Mats::uidCol > > >, ordered_non_unique< tag< row_mi_tag >, member< SchurL2Mats, const UId, &SchurL2Mats::uidRow > >, ordered_non_unique< tag< col_mi_tag >, member< SchurL2Mats, const UId, &SchurL2Mats::uidCol > > > > SchurL2Storage
Definition: Schur.hpp:180
static boost::ptr_vector< VectorInt > colIndices
Definition: Schur.hpp:184
const size_t iDX
Definition: Schur.hpp:150
boost::function< MoFEMErrorCode(Mat M, const EntitiesFieldData::EntData &row_data, const EntitiesFieldData::EntData &col_data, const MatrixDouble &mat, InsertMode iora)> MatSetValuesPtr
Definition: Schur.hpp:140
auto & getMat() const
Definition: Schur.hpp:133
static MoFEMErrorCode MatSetValues(Mat M, const EntitiesFieldData::EntData &row_data, const EntitiesFieldData::EntData &col_data, const MatrixDouble &mat, InsertMode iora)
Definition: Schur.cpp:32
friend OpSchurAssembleBegin
Definition: Schur.hpp:152
static boost::ptr_vector< MatrixDouble > locMats
Definition: Schur.hpp:182
auto & getRowInd() const
Definition: Schur.hpp:134
auto & getColInd() const
Definition: Schur.hpp:135
intrusive_ptr for managing petsc objects