v0.15.0
Loading...
Searching...
No Matches
KspCtx.cpp
Go to the documentation of this file.
1/**
2 * @brief Context for PETSc KSP, i.e. linear solver
3 * @author Anonymous contributor under MIT license
4 *
5 */
6
7namespace MoFEM {
8
9struct KspCtx::KspCtx::KspCtxImpl {
10
11 KspCtxImpl(MoFEM::Interface &m_field, std::string problem_name);
12
13 virtual ~KspCtxImpl() = default;
14
15 /**
16 * @return return reference to vector with FEMethod to calculate matrix
17 */
19
20 /**
21 * @return return vector to vector with FEMethod to vector
22 */
24
25 /**
26 * The sequence of BasicMethod is executed before residual is calculated. It
27 * can be used to setup data structures, e.g. zero global variable which is
28 * integrated in domain, e.g. for calculation of strain energy.
29 *
30 * @return reference to BasicMethod for preprocessing
31 */
33
34 /**
35 * The sequence of BasicMethod is executed after residual is calculated. It
36 * can be used to setup data structures, e.g. aggregate data from processors
37 * or to apply essential boundary conditions.
38 *
39 * @return reference to BasicMethod for postprocessing
40 */
42
43 /**
44 * @return reference to BasicMethod for preprocessing
45 */
47
48 /**
49 * The sequence of BasicMethod is executed after tangent matrix is calculated.
50 * It can be used to setup data structures, e.g. aggregate data from
51 * processors or to apply essential boundary conditions.
52 *
53 * @return reference to BasicMethod for postprocessing
54 */
58
59 friend PetscErrorCode KspRhs(KSP ksp, Vec f, void *ctx);
60 friend PetscErrorCode KspMat(KSP ksp, Mat A, Mat B, void *ctx);
61
62 /**
63 * @brief Clear loops
64 *
65 * @return MoFEMErrorCode
66 */
67 MoFEMErrorCode clearLoops();
68
69protected:
70 PetscLogEvent MOFEM_EVENT_KspRhs;
71 PetscLogEvent MOFEM_EVENT_KspMat;
72
73 boost::movelib::unique_ptr<bool> vecAssembleSwitch;
74 boost::movelib::unique_ptr<bool> matAssembleSwitch;
75
77 moab::Interface &moab;
78
79 std::string problemName; ///< Problem name
80 MoFEMTypes bH; ///< If set to MF_EXIST check if element exist
81
82 FEMethodsSequence loopsOperator; ///< Sequence of finite elements instances
83 ///< assembling tangent matrix
84 FEMethodsSequence loopsRhs; ///< Sequence of finite elements
85 ///< instances assembling residual vector
86 BasicMethodsSequence preProcessOperator; ///< Sequence of methods run before
87 ///< tangent matrix is assembled
88 BasicMethodsSequence postProcessOperator; ///< Sequence of methods run after
89 ///< tangent matrix is assembled
90 BasicMethodsSequence preProcessRhs; ///< Sequence of methods run before
91 ///< residual is assembled
92 BasicMethodsSequence postProcessRhs; ///< Sequence of methods run after
93 ///< residual is assembled
94};
95
96KspCtx::KspCtx(MoFEM::Interface &m_field, std::string problem_name)
97 : kspCtxImpl(boost::movelib::make_unique<KspCtx::KspCtxImpl>(
98 m_field, problem_name)) {}
99
100KspCtx::~KspCtx() = default;
101
103 return kspCtxImpl->getSetOperators();
104}
105
107 return kspCtxImpl->getComputeRhs();
108}
109
111 return kspCtxImpl->getPreProcComputeRhs();
112}
113
115 return kspCtxImpl->getPostProcComputeRhs();
116}
117
119 return kspCtxImpl->getPreProcSetOperators();
120}
121
123 return kspCtxImpl->getPostProcSetOperators();
124}
125
127
128KspCtx::KspCtxImpl::KspCtxImpl(MoFEM::Interface &m_field,
129 std::string problem_name)
130 : mField(m_field), moab(m_field.get_moab()), problemName(problem_name),
131 bH(MF_EXIST) {
132 PetscLogEventRegister("LoopKSPRhs", 0, &MOFEM_EVENT_KspRhs);
133 PetscLogEventRegister("LoopKSPMat", 0, &MOFEM_EVENT_KspMat);
134}
135
136MoFEMErrorCode KspCtx::KspCtxImpl::clearLoops() {
138 loopsOperator.clear();
139 loopsRhs.clear();
140 preProcessOperator.clear();
141 postProcessOperator.clear();
142 preProcessRhs.clear();
144}
145
146PetscErrorCode KspRhs(KSP ksp, Vec f, void *ctx) {
148 KspCtx::KspCtxImpl *ksp_ctx = static_cast<KspCtx *>(ctx)->kspCtxImpl.get();
149 PetscLogEventBegin(ksp_ctx->MOFEM_EVENT_KspRhs, 0, 0, 0, 0);
150
151 ksp_ctx->vecAssembleSwitch = boost::movelib::make_unique<bool>(true);
152
153 auto cache_ptr = boost::make_shared<CacheTuple>();
154 CHKERR ksp_ctx->mField.cache_problem_entities(ksp_ctx->problemName,
155 cache_ptr);
156
157 auto set = [&](auto &fe) {
158 fe.ksp = ksp;
159 fe.ksp_ctx = KspMethod::CTX_SETFUNCTION;
160 fe.data_ctx = PetscData::CtxSetF;
161 fe.ksp_f = f;
162 fe.cacheWeakPtr = cache_ptr;
163 };
164
165 auto unset = [&](auto &fe) {
166 fe.ksp_ctx = KspMethod::CTX_KSPNONE;
167 fe.data_ctx = PetscData::CtxSetNone;
168 };
169
170 // pre-process
171 for (auto &bit : ksp_ctx->preProcessRhs) {
172 bit->vecAssembleSwitch = boost::move(ksp_ctx->vecAssembleSwitch);
173 set(*bit);
174 CHKERR ksp_ctx->mField.problem_basic_method_preProcess(ksp_ctx->problemName,
175 *bit);
176 unset(*bit);
177 ksp_ctx->vecAssembleSwitch = boost::move(bit->vecAssembleSwitch);
178 }
179
180 // operators
181 for (auto &lit : ksp_ctx->loopsRhs) {
182 lit.second->vecAssembleSwitch = boost::move(ksp_ctx->vecAssembleSwitch);
183 set(*lit.second);
184 CHKERR ksp_ctx->mField.loop_finite_elements(ksp_ctx->problemName, lit.first,
185 *(lit.second), nullptr,
186 ksp_ctx->bH, cache_ptr);
187 unset(*lit.second);
188 ksp_ctx->vecAssembleSwitch = boost::move(lit.second->vecAssembleSwitch);
189 }
190
191 // post-process
192 for (auto &bit : ksp_ctx->postProcessRhs) {
193 bit->vecAssembleSwitch = boost::move(ksp_ctx->vecAssembleSwitch);
194 set(*bit);
195 CHKERR ksp_ctx->mField.problem_basic_method_postProcess(
196 ksp_ctx->problemName, *bit);
197 unset(*bit);
198 ksp_ctx->vecAssembleSwitch = boost::move(bit->vecAssembleSwitch);
199 }
200
201 if (*ksp_ctx->vecAssembleSwitch) {
202 CHKERR VecGhostUpdateBegin(f, ADD_VALUES, SCATTER_REVERSE);
203 CHKERR VecGhostUpdateEnd(f, ADD_VALUES, SCATTER_REVERSE);
204 CHKERR VecAssemblyBegin(f);
205 CHKERR VecAssemblyEnd(f);
206 }
207 PetscLogEventEnd(ksp_ctx->MOFEM_EVENT_KspRhs, 0, 0, 0, 0);
209}
210PetscErrorCode KspMat(KSP ksp, Mat A, Mat B, void *ctx) {
211 // PetscValidHeaderSpecific(ksp,KSP_CLASSID,1);
213 KspCtx::KspCtxImpl *ksp_ctx = static_cast<KspCtx *>(ctx)->kspCtxImpl.get();
214 PetscLogEventBegin(ksp_ctx->MOFEM_EVENT_KspMat, 0, 0, 0, 0);
215
216 ksp_ctx->matAssembleSwitch = boost::movelib::make_unique<bool>(true);
217 auto cache_ptr = boost::make_shared<CacheTuple>();
218 CHKERR ksp_ctx->mField.cache_problem_entities(ksp_ctx->problemName,
219 cache_ptr);
220
221 auto set = [&](auto &fe) {
222 fe.ksp = ksp;
223 fe.ksp_A = A;
224 fe.ksp_B = B;
225 fe.ksp_ctx = KspMethod::CTX_OPERATORS;
227 fe.cacheWeakPtr = cache_ptr;
228 };
229
230 auto unset = [&](auto &fe) {
231 fe.ksp_ctx = KspMethod::CTX_KSPNONE;
232 fe.data_ctx = PetscData::CtxSetNone;
233 };
234
235 auto ent_data_cache = boost::make_shared<std::vector<EntityCacheDofs>>();
236 auto ent_row_cache =
237 boost::make_shared<std::vector<EntityCacheNumeredDofs>>();
238 auto ent_col_cache =
239 boost::make_shared<std::vector<EntityCacheNumeredDofs>>();
240
241 // pre-procsess
242 for (auto &bit : ksp_ctx->preProcessOperator) {
243 bit->matAssembleSwitch = boost::move(ksp_ctx->matAssembleSwitch);
244 set(*bit);
245 CHKERR ksp_ctx->mField.problem_basic_method_preProcess(ksp_ctx->problemName,
246 *bit);
247 unset(*bit);
248 ksp_ctx->matAssembleSwitch = boost::move(bit->matAssembleSwitch);
249 }
250
251 // operators
252 for (auto &lit : ksp_ctx->loopsOperator) {
253 lit.second->matAssembleSwitch = boost::move(ksp_ctx->matAssembleSwitch);
254 set(*lit.second);
255 CHKERR ksp_ctx->mField.loop_finite_elements(ksp_ctx->problemName, lit.first,
256 *(lit.second), nullptr,
257 ksp_ctx->bH, cache_ptr);
258 unset(*lit.second);
259 ksp_ctx->matAssembleSwitch = boost::move(lit.second->matAssembleSwitch);
260 }
261
262 // post-process
263 for (auto &bit : ksp_ctx->postProcessOperator) {
264 bit->matAssembleSwitch = boost::move(ksp_ctx->matAssembleSwitch);
265 set(*bit);
266 CHKERR ksp_ctx->mField.problem_basic_method_postProcess(
267 ksp_ctx->problemName, *bit);
268 unset(*bit);
269 ksp_ctx->matAssembleSwitch = boost::move(bit->matAssembleSwitch);
270 }
271
272 if (ksp_ctx->matAssembleSwitch) {
273 CHKERR MatAssemblyBegin(B, MAT_FINAL_ASSEMBLY);
274 CHKERR MatAssemblyEnd(B, MAT_FINAL_ASSEMBLY);
275 }
276 PetscLogEventEnd(ksp_ctx->MOFEM_EVENT_KspMat, 0, 0, 0, 0);
278}
279
280} // namespace MoFEM
MoFEMTypes
Those types control how functions respond on arguments, f.e. error handling.
@ MF_EXIST
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
#define CHKERR
Inline error check.
auto bit
set bit
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
implementation of Data Operators for Forces and Sources
Definition Common.hpp:10
std::deque< BasicMethodPtr > BasicMethodsSequence
Definition AuxPETSc.hpp:57
std::deque< PairNameFEMethodPtr > FEMethodsSequence
Definition AuxPETSc.hpp:56
constexpr AssemblyType A
Deprecated interface functions.
MoFEM::Interface & mField
Definition KspCtx.cpp:76
BasicMethodsSequence preProcessOperator
Definition KspCtx.cpp:86
boost::movelib::unique_ptr< bool > matAssembleSwitch
Definition KspCtx.cpp:74
BasicMethodsSequence & getPreProcSetOperators()
Definition KspCtx.cpp:46
BasicMethodsSequence postProcessOperator
Definition KspCtx.cpp:88
std::string problemName
Problem name.
Definition KspCtx.cpp:79
FEMethodsSequence & getSetOperators()
Definition KspCtx.cpp:18
MoFEMTypes bH
If set to MF_EXIST check if element exist.
Definition KspCtx.cpp:80
FEMethodsSequence loopsOperator
Definition KspCtx.cpp:82
BasicMethodsSequence & getPostProcComputeRhs()
Definition KspCtx.cpp:41
FEMethodsSequence & getComputeRhs()
Definition KspCtx.cpp:23
BasicMethodsSequence postProcessRhs
Definition KspCtx.cpp:92
boost::movelib::unique_ptr< bool > vecAssembleSwitch
Definition KspCtx.cpp:73
BasicMethodsSequence & getPostProcSetOperators()
Definition KspCtx.cpp:55
BasicMethodsSequence preProcessRhs
Definition KspCtx.cpp:90
BasicMethodsSequence & getPreProcComputeRhs()
Definition KspCtx.cpp:32
Interface for linear (KSP) solver.
Definition KspCtx.hpp:14
boost::movelib::unique_ptr< KspCtxImpl > kspCtxImpl
Definition KspCtx.hpp:74
BasicMethodsSequence & getPostProcComputeRhs()
Definition KspCtx.cpp:114
BasicMethodsSequence & getPreProcComputeRhs()
Definition KspCtx.cpp:110
BasicMethodsSequence & getPreProcSetOperators()
Definition KspCtx.cpp:118
FEMethodsSequence & getSetOperators()
Definition KspCtx.cpp:102
friend PetscErrorCode KspRhs(KSP ksp, Vec f, void *ctx)
Run over elements in the lists.
Definition KspCtx.cpp:146
FEMethodsSequence & getComputeRhs()
Definition KspCtx.cpp:106
MoFEMErrorCode clearLoops()
Clear loops.
Definition KspCtx.cpp:126
virtual ~KspCtx()
friend PetscErrorCode KspMat(KSP ksp, Mat A, Mat B, void *ctx)
Run over elements in the list.
Definition KspCtx.cpp:210
BasicMethodsSequence & getPostProcSetOperators()
Definition KspCtx.cpp:122
@ CTX_KSPNONE
No specific KSP context.
@ CTX_SETFUNCTION
Setting up the linear system function.
@ CTX_OPERATORS
Setting up linear operators.
static constexpr Switches CtxSetA
Jacobian matrix switch.
static constexpr Switches CtxSetNone
No data switch.
static constexpr Switches CtxSetF
Residual vector switch.
static constexpr Switches CtxSetB
Preconditioner matrix switch.