v0.13.2
Loading...
Searching...
No Matches
TsCtx.cpp
Go to the documentation of this file.
1
2
3namespace MoFEM {
4
5MoFEMErrorCode TsCtx::clearLoops() {
6 MoFEMFunctionBeginHot;
7 loopsIJacobian.clear();
8 loopsIFunction.clear();
9 loopsMonitor.clear();
10 loopsRHSJacobian.clear();
11 loopsRHSFunction.clear();
12 preProcessIJacobian.clear();
13 postProcessIJacobian.clear();
14 preProcessIFunction.clear();
15 postProcessIFunction.clear();
16 preProcessMonitor.clear();
17 postProcessMonitor.clear();
18 preProcessRHSJacobian.clear();
19 preProcessRHSFunction.clear();
20 postProcessRHSJacobian.clear();
21 postProcessRHSFunction.clear();
22 MoFEMFunctionReturnHot(0);
23}
24
25PetscErrorCode TsSetIFunction(TS ts, PetscReal t, Vec u, Vec u_t, Vec F,
26 void *ctx) {
27 MoFEMFunctionBegin;
28 TsCtx *ts_ctx = static_cast<TsCtx *>(ctx);
29 PetscLogEventBegin(ts_ctx->MOFEM_EVENT_TsCtxIFunction, 0, 0, 0, 0);
30 CHKERR VecGhostUpdateBegin(u, INSERT_VALUES, SCATTER_FORWARD);
31 CHKERR VecGhostUpdateEnd(u, INSERT_VALUES, SCATTER_FORWARD);
32 CHKERR VecGhostUpdateBegin(u_t, INSERT_VALUES, SCATTER_FORWARD);
33 CHKERR VecGhostUpdateEnd(u_t, INSERT_VALUES, SCATTER_FORWARD);
34 CHKERR ts_ctx->mField.getInterface<VecManager>()->setLocalGhostVector(
35 ts_ctx->problemName, COL, u, INSERT_VALUES, SCATTER_REVERSE);
36
37 auto zero_ghost_vec = [](Vec g) {
38 MoFEMFunctionBegin;
39 Vec l;
40 CHKERR VecGhostGetLocalForm(g, &l);
41 double *a;
42 CHKERR VecGetArray(l, &a);
43 int s;
44 CHKERR VecGetLocalSize(l, &s);
45 for (int i = 0; i != s; ++i)
46 a[i] = 0;
47 CHKERR VecRestoreArray(l, &a);
48 CHKERR VecGhostRestoreLocalForm(g, &l);
49 MoFEMFunctionReturn(0);
50 };
51 CHKERR zero_ghost_vec(F);
52
53 ts_ctx->vecAssembleSwitch = boost::movelib::make_unique<bool>(true);
54
55 int step;
56#if PETSC_VERSION_GE(3, 8, 0)
57 CHKERR TSGetStepNumber(ts, &step);
58#else
59 CHKERR TSGetTimeStepNumber(ts, &step);
60#endif
61
62 auto set = [&](auto &fe) {
63 fe.ts = ts;
64 fe.ts_u = u;
65 fe.ts_u_t = u_t;
66 fe.ts_F = F;
67 fe.ts_t = t;
68 fe.ts_step = step;
69 fe.ts_ctx = TSMethod::CTX_TSSETIFUNCTION;
70 fe.snes_ctx = SnesMethod::CTX_SNESSETFUNCTION;
71 fe.ksp_ctx = KspMethod::CTX_SETFUNCTION;
72 fe.data_ctx = PetscData::CtxSetF | PetscData::CtxSetX |
73 PetscData::CtxSetX_T | PetscData::CtxSetTime;
74 CHK_THROW_MESSAGE(TSGetTimeStep(ts, &fe.ts_dt), "get time step failed");
75 };
76
77 auto unset = [&](auto &fe) {
78 fe.ts_ctx = TSMethod::CTX_TSNONE;
79 fe.snes_ctx = SnesMethod::CTX_SNESNONE;
80 fe.ksp_ctx = KspMethod::CTX_KSPNONE;
81 fe.data_ctx = PetscData::CtxSetNone;
82 };
83
84 // preprocess
85 for (auto &bit : ts_ctx->preProcessIFunction) {
86 bit->vecAssembleSwitch = boost::move(ts_ctx->vecAssembleSwitch);
87 set(*bit);
88 CHKERR ts_ctx->mField.problem_basic_method_preProcess(ts_ctx->problemName,
89 *bit);
90 unset(*bit);
91 ts_ctx->vecAssembleSwitch = boost::move(bit->vecAssembleSwitch);
92 }
93
94 auto cache_ptr = boost::make_shared<CacheTuple>();
95 CHKERR ts_ctx->mField.cache_problem_entities(ts_ctx->problemName, cache_ptr);
96
97 // fe loops
98 for (auto &lit : ts_ctx->loopsIFunction) {
99 lit.second->vecAssembleSwitch = boost::move(ts_ctx->vecAssembleSwitch);
100 set(*lit.second);
101 CHKERR ts_ctx->mField.loop_finite_elements(ts_ctx->problemName, lit.first,
102 *(lit.second), nullptr,
103 ts_ctx->bH, cache_ptr);
104 unset(*lit.second);
105 ts_ctx->vecAssembleSwitch = boost::move(lit.second->vecAssembleSwitch);
106 }
107
108 // post process
109 for (auto &bit : ts_ctx->postProcessIFunction) {
110 bit->vecAssembleSwitch = boost::move(ts_ctx->vecAssembleSwitch);
111 set(*bit);
112 CHKERR ts_ctx->mField.problem_basic_method_postProcess(ts_ctx->problemName,
113 *bit);
114 unset(*bit);
115 ts_ctx->vecAssembleSwitch = boost::move(bit->vecAssembleSwitch);
116 }
117
118 if (*ts_ctx->vecAssembleSwitch) {
119 CHKERR VecGhostUpdateBegin(F, ADD_VALUES, SCATTER_REVERSE);
120 CHKERR VecGhostUpdateEnd(F, ADD_VALUES, SCATTER_REVERSE);
121 CHKERR VecAssemblyBegin(F);
122 CHKERR VecAssemblyEnd(F);
123 }
124
125 PetscLogEventEnd(ts_ctx->MOFEM_EVENT_TsCtxIFunction, 0, 0, 0, 0);
126 MoFEMFunctionReturn(0);
127}
128
129PetscErrorCode TsSetIJacobian(TS ts, PetscReal t, Vec u, Vec u_t, PetscReal a,
130 Mat A, Mat B, void *ctx) {
131 MoFEMFunctionBegin;
132
133 TsCtx *ts_ctx = static_cast<TsCtx *>(ctx);
134 PetscLogEventBegin(ts_ctx->MOFEM_EVENT_TsCtxIFunction, 0, 0, 0, 0);
135 CHKERR VecGhostUpdateBegin(u, INSERT_VALUES, SCATTER_FORWARD);
136 CHKERR VecGhostUpdateEnd(u, INSERT_VALUES, SCATTER_FORWARD);
137 CHKERR VecGhostUpdateBegin(u_t, INSERT_VALUES, SCATTER_FORWARD);
138 CHKERR VecGhostUpdateEnd(u_t, INSERT_VALUES, SCATTER_FORWARD);
139 CHKERR ts_ctx->mField.getInterface<VecManager>()->setLocalGhostVector(
140 ts_ctx->problemName, COL, u, INSERT_VALUES, SCATTER_REVERSE);
141 if (ts_ctx->zeroMatrix) {
142 CHKERR MatZeroEntries(B);
143 }
144 int step;
145#if PETSC_VERSION_GE(3, 8, 0)
146 CHKERR TSGetStepNumber(ts, &step);
147#else
148 CHKERR TSGetTimeStepNumber(ts, &step);
149#endif
150
151 ts_ctx->matAssembleSwitch =
152 boost::movelib::make_unique<bool>(ts_ctx->zeroMatrix);
153
154 auto set = [&](auto &fe) {
155 fe.ts = ts;
156 fe.ts_u = u;
157 fe.ts_u_t = u_t;
158 fe.ts_A = A;
159 fe.ts_B = B;
160 fe.ts_t = t;
161 fe.ts_a = a;
162 fe.ts_step = step;
163 fe.ts_ctx = TSMethod::CTX_TSSETIJACOBIAN;
164 fe.snes_ctx = SnesMethod::CTX_SNESSETJACOBIAN;
165 fe.ksp_ctx = KspMethod::CTX_OPERATORS;
166 fe.data_ctx = PetscData::CtxSetA | PetscData::CtxSetB | PetscData::CtxSetX |
167 PetscData::CtxSetX_T | PetscData::CtxSetTime;
168 CHK_THROW_MESSAGE(TSGetTimeStep(ts, &fe.ts_dt), "get time step failed");
169 };
170
171 auto unset = [&](auto &fe) {
172 fe.ts_ctx = TSMethod::CTX_TSNONE;
173 fe.snes_ctx = SnesMethod::CTX_SNESNONE;
174 fe.ksp_ctx = KspMethod::CTX_KSPNONE;
175 fe.data_ctx = PetscData::CtxSetNone;
176 };
177
178 // preproces
179 for (auto &bit : ts_ctx->preProcessIJacobian) {
180 bit->matAssembleSwitch = boost::move(ts_ctx->matAssembleSwitch);
181 set(*bit);
182 CHKERR ts_ctx->mField.problem_basic_method_preProcess(ts_ctx->problemName,
183 *bit);
184 unset(*bit);
185 ts_ctx->matAssembleSwitch = boost::move(bit->matAssembleSwitch);
186 }
187
188 auto cache_ptr = boost::make_shared<CacheTuple>();
189 CHKERR ts_ctx->mField.cache_problem_entities(ts_ctx->problemName, cache_ptr);
190
191 for (auto &lit : ts_ctx->loopsIJacobian) {
192 lit.second->matAssembleSwitch = boost::move(ts_ctx->matAssembleSwitch);
193 set(*lit.second);
194 CHKERR ts_ctx->mField.loop_finite_elements(ts_ctx->problemName, lit.first,
195 *(lit.second), nullptr,
196 ts_ctx->bH, cache_ptr);
197 unset(*lit.second);
198 ts_ctx->matAssembleSwitch = boost::move(lit.second->matAssembleSwitch);
199 }
200
201 // post process
202 for (auto &bit : ts_ctx->postProcessIJacobian) {
203 bit->matAssembleSwitch = boost::move(ts_ctx->matAssembleSwitch);
204 set(*bit);
205 CHKERR ts_ctx->mField.problem_basic_method_postProcess(ts_ctx->problemName,
206 *bit);
207 unset(*bit);
208 ts_ctx->matAssembleSwitch = boost::move(bit->matAssembleSwitch);
209 }
210
211 if (ts_ctx->matAssembleSwitch) {
212 CHKERR MatAssemblyBegin(B, MAT_FINAL_ASSEMBLY);
213 CHKERR MatAssemblyEnd(B, MAT_FINAL_ASSEMBLY);
214 }
215 PetscLogEventEnd(ts_ctx->MOFEM_EVENT_TsCtxIFunction, 0, 0, 0, 0);
216 MoFEMFunctionReturn(0);
217}
218
219PetscErrorCode TsMonitorSet(TS ts, PetscInt step, PetscReal t, Vec u,
220 void *ctx) {
221 MoFEMFunctionBegin;
222 TsCtx *ts_ctx = static_cast<TsCtx *>(ctx);
223 PetscLogEventBegin(ts_ctx->MOFEM_EVENT_TsCtxMonitor, 0, 0, 0, 0);
224 CHKERR VecGhostUpdateBegin(u, INSERT_VALUES, SCATTER_FORWARD);
225 CHKERR VecGhostUpdateEnd(u, INSERT_VALUES, SCATTER_FORWARD);
226 CHKERR ts_ctx->mField.getInterface<VecManager>()->setLocalGhostVector(
227 ts_ctx->problemName, COL, u, INSERT_VALUES, SCATTER_REVERSE);
228
229 auto set = [&](auto &fe) {
230 fe.ts = ts;
231 fe.ts_u = u;
232 fe.ts_t = t;
233 fe.ts_step = step;
234 fe.ts_F = PETSC_NULL;
235 fe.ts_ctx = TSMethod::CTX_TSTSMONITORSET;
236 fe.snes_ctx = SnesMethod::CTX_SNESNONE;
237 fe.ksp_ctx = KspMethod::CTX_KSPNONE;
238 fe.data_ctx = PetscData::CtxSetX | PetscData::CtxSetTime;
239 CHK_THROW_MESSAGE(TSGetTimeStep(ts, &fe.ts_dt), "get time step failed");
240 };
241
242 auto unset = [&](auto &fe) {
243 fe.ts_ctx = TSMethod::CTX_TSNONE;
244 fe.data_ctx = PetscData::CtxSetNone;
245 };
246
247 // preproces
248 for (auto &bit : ts_ctx->preProcessMonitor) {
249 set(*bit);
250 CHKERR ts_ctx->mField.problem_basic_method_preProcess(ts_ctx->problemName,
251 *bit);
252 unset(*bit);
253 }
254
255 auto cache_ptr = boost::make_shared<CacheTuple>();
256 CHKERR ts_ctx->mField.cache_problem_entities(ts_ctx->problemName, cache_ptr);
257
258 for (auto &lit : ts_ctx->loopsMonitor) {
259 set(*lit.second);
260 CHKERR ts_ctx->mField.loop_finite_elements(ts_ctx->problemName, lit.first,
261 *(lit.second), nullptr,
262 ts_ctx->bH, cache_ptr);
263 unset(*lit.second);
264 }
265
266 // post process
267 for (auto &bit : ts_ctx->postProcessMonitor) {
268 set(*bit);
269 CHKERR ts_ctx->mField.problem_basic_method_postProcess(ts_ctx->problemName,
270 *bit);
271 unset(*bit);
272 }
273
274 PetscLogEventEnd(ts_ctx->MOFEM_EVENT_TsCtxMonitor, 0, 0, 0, 0);
275 MoFEMFunctionReturn(0);
276}
277
278PetscErrorCode TsSetRHSFunction(TS ts, PetscReal t, Vec u, Vec F, void *ctx) {
279 MoFEMFunctionBegin;
280 TsCtx *ts_ctx = static_cast<TsCtx *>(ctx);
281 PetscLogEventBegin(ts_ctx->MOFEM_EVENT_TsCtxRHSFunction, 0, 0, 0, 0);
282 CHKERR VecGhostUpdateBegin(u, INSERT_VALUES, SCATTER_FORWARD);
283 CHKERR VecGhostUpdateEnd(u, INSERT_VALUES, SCATTER_FORWARD);
284 CHKERR ts_ctx->mField.getInterface<VecManager>()->setLocalGhostVector(
285 ts_ctx->problemName, COL, u, INSERT_VALUES, SCATTER_REVERSE);
286
287 auto zero_ghost_vec = [](Vec g) {
288 MoFEMFunctionBegin;
289 Vec l;
290 CHKERR VecGhostGetLocalForm(g, &l);
291 double *a;
292 CHKERR VecGetArray(l, &a);
293 int s;
294 CHKERR VecGetLocalSize(l, &s);
295 for (int i = 0; i != s; ++i)
296 a[i] = 0;
297 CHKERR VecRestoreArray(l, &a);
298 CHKERR VecGhostRestoreLocalForm(g, &l);
299 MoFEMFunctionReturn(0);
300 };
301 CHKERR zero_ghost_vec(F);
302
303 ts_ctx->vecAssembleSwitch = boost::movelib::make_unique<bool>(true);
304
305 int step;
306#if PETSC_VERSION_GE(3, 8, 0)
307 CHKERR TSGetStepNumber(ts, &step);
308#else
309 CHKERR TSGetTimeStepNumber(ts, &step);
310#endif
311
312 auto set = [&](auto &fe) {
313 fe.ts_u = u;
314 fe.ts_F = F;
315 fe.ts_t = t;
316 fe.ts = ts;
317 fe.ts_step = step;
318 fe.ts_ctx = TSMethod::CTX_TSSETRHSFUNCTION;
319 fe.snes_ctx = SnesMethod::CTX_SNESSETFUNCTION;
320 fe.ksp_ctx = KspMethod::CTX_SETFUNCTION;
321 fe.data_ctx =
322 PetscData::CtxSetF | PetscData::CtxSetX | PetscData::CtxSetTime;
323 CHK_THROW_MESSAGE(TSGetTimeStep(ts, &fe.ts_dt), "get time step failed");
324 };
325
326 auto unset = [&](auto &fe) {
327 fe.ts_ctx = TSMethod::CTX_TSNONE;
328 fe.snes_ctx = SnesMethod::CTX_SNESNONE;
329 fe.ksp_ctx = KspMethod::CTX_KSPNONE;
330 fe.data_ctx = PetscData::CtxSetNone;
331 };
332
333 for (auto &bit : ts_ctx->preProcessRHSJacobian) {
334 bit->vecAssembleSwitch = boost::move(ts_ctx->vecAssembleSwitch);
335 set(*bit);
336 CHKERR ts_ctx->mField.problem_basic_method_preProcess(ts_ctx->problemName,
337 *bit);
338 unset(*bit);
339 ts_ctx->vecAssembleSwitch = boost::move(bit->vecAssembleSwitch);
340 }
341
342 auto cache_ptr = boost::make_shared<CacheTuple>();
343 CHKERR ts_ctx->mField.cache_problem_entities(ts_ctx->problemName, cache_ptr);
344
345 // fe loops
346 for (auto &lit : ts_ctx->loopsRHSFunction) {
347 lit.second->vecAssembleSwitch = boost::move(ts_ctx->vecAssembleSwitch);
348 set(*lit.second);
349 CHKERR ts_ctx->mField.loop_finite_elements(ts_ctx->problemName, lit.first,
350 *(lit.second), nullptr,
351 ts_ctx->bH, cache_ptr);
352 unset(*lit.second);
353 ts_ctx->vecAssembleSwitch = boost::move(lit.second->vecAssembleSwitch);
354 }
355
356 // post process
357 for (auto &bit : ts_ctx->postProcessRHSJacobian) {
358 bit->vecAssembleSwitch = boost::move(ts_ctx->vecAssembleSwitch);
359 set(*bit);
360 CHKERR ts_ctx->mField.problem_basic_method_postProcess(ts_ctx->problemName,
361 *bit);
362 unset(*bit);
363 ts_ctx->vecAssembleSwitch = boost::move(bit->vecAssembleSwitch);
364 }
365
366 if (*ts_ctx->vecAssembleSwitch) {
367 CHKERR VecGhostUpdateBegin(F, ADD_VALUES, SCATTER_REVERSE);
368 CHKERR VecGhostUpdateEnd(F, ADD_VALUES, SCATTER_REVERSE);
369 CHKERR VecAssemblyBegin(F);
370 CHKERR VecAssemblyEnd(F);
371 }
372
373 PetscLogEventEnd(ts_ctx->MOFEM_EVENT_TsCtxRHSFunction, 0, 0, 0, 0);
374 MoFEMFunctionReturn(0);
375}
376
377PetscErrorCode TsSetRHSJacobian(TS ts, PetscReal t, Vec u, Mat A, Mat B,
378 void *ctx) {
379 MoFEMFunctionBegin;
380 TsCtx *ts_ctx = static_cast<TsCtx *>(ctx);
381 PetscLogEventBegin(ts_ctx->MOFEM_EVENT_TsCtxRHSJacobian, 0, 0, 0, 0);
382 CHKERR VecGhostUpdateBegin(u, INSERT_VALUES, SCATTER_FORWARD);
383 CHKERR VecGhostUpdateEnd(u, INSERT_VALUES, SCATTER_FORWARD);
384 CHKERR ts_ctx->mField.getInterface<VecManager>()->setLocalGhostVector(
385 ts_ctx->problemName, COL, u, INSERT_VALUES, SCATTER_REVERSE);
386
387 if (ts_ctx->zeroMatrix) {
388 CHKERR MatZeroEntries(B);
389 }
390
391 ts_ctx->matAssembleSwitch =
392 boost::movelib::make_unique<bool>(ts_ctx->zeroMatrix);
393
394 int step;
395#if PETSC_VERSION_GE(3, 8, 0)
396 CHKERR TSGetStepNumber(ts, &step);
397#else
398 CHKERR TSGetTimeStepNumber(ts, &step);
399#endif
400
401 auto set = [&](auto &fe) {
402 fe.ts_u = u;
403 fe.ts_A = A;
404 fe.ts_B = B;
405 fe.ts_t = t;
406 fe.ts_step = step;
407 fe.ts_ctx = TSMethod::CTX_TSSETRHSJACOBIAN;
408 fe.snes_ctx = SnesMethod::CTX_SNESSETJACOBIAN;
409 fe.ksp_ctx = KspMethod::CTX_OPERATORS;
410 fe.data_ctx = PetscData::CtxSetA | PetscData::CtxSetB | PetscData::CtxSetX |
411 PetscData::CtxSetTime;
412 fe.ts = ts;
413 CHK_THROW_MESSAGE(TSGetTimeStep(ts, &fe.ts_dt), "get time step failed");
414 };
415
416 auto unset = [&](auto &fe) {
417 fe.ts_ctx = TSMethod::CTX_TSNONE;
418 fe.snes_ctx = SnesMethod::CTX_SNESNONE;
419 fe.ksp_ctx = KspMethod::CTX_KSPNONE;
420 fe.data_ctx = PetscData::CtxSetNone;
421 };
422
423 // preprocess
424 for (auto &bit : ts_ctx->preProcessRHSJacobian) {
425 bit->matAssembleSwitch = boost::move(ts_ctx->matAssembleSwitch);
426 set(*bit);
427 CHKERR ts_ctx->mField.problem_basic_method_preProcess(ts_ctx->problemName,
428 *bit);
429 unset(*bit);
430 ts_ctx->matAssembleSwitch = boost::move(bit->matAssembleSwitch);
431 }
432
433 auto cache_ptr = boost::make_shared<CacheTuple>();
434 CHKERR ts_ctx->mField.cache_problem_entities(ts_ctx->problemName, cache_ptr);
435
436 // fe loops
437 for (auto &lit : ts_ctx->loopsRHSJacobian) {
438 lit.second->matAssembleSwitch = boost::move(ts_ctx->matAssembleSwitch);
439 set(*lit.second);
440 CHKERR ts_ctx->mField.loop_finite_elements(ts_ctx->problemName, lit.first,
441 *(lit.second), nullptr,
442 ts_ctx->bH, cache_ptr);
443 unset(*lit.second);
444 ts_ctx->matAssembleSwitch = boost::move(lit.second->matAssembleSwitch);
445 }
446
447 // post process
448 for (auto &bit : ts_ctx->postProcessRHSJacobian) {
449 bit->matAssembleSwitch = boost::move(ts_ctx->matAssembleSwitch);
450 set(*bit);
451 CHKERR ts_ctx->mField.problem_basic_method_postProcess(ts_ctx->problemName,
452 *bit);
453 unset(*bit);
454 ts_ctx->matAssembleSwitch = boost::move(bit->matAssembleSwitch);
455 }
456
457 if (ts_ctx->matAssembleSwitch) {
458 CHKERR MatAssemblyBegin(B, MAT_FINAL_ASSEMBLY);
459 CHKERR MatAssemblyEnd(B, MAT_FINAL_ASSEMBLY);
460 }
461
462 PetscLogEventEnd(ts_ctx->MOFEM_EVENT_TsCtxRHSJacobian, 0, 0, 0, 0);
463 MoFEMFunctionReturn(0);
464}
465
466PetscErrorCode TsSetI2Jacobian(TS ts, PetscReal t, Vec u, Vec u_t, Vec u_tt,
467 PetscReal a, PetscReal aa, Mat A, Mat B,
468 void *ctx) {
469 MoFEMFunctionBegin;
470
471 TsCtx *ts_ctx = static_cast<TsCtx *>(ctx);
472 PetscLogEventBegin(ts_ctx->MOFEM_EVENT_TsCtxI2Function, 0, 0, 0, 0);
473 CHKERR VecGhostUpdateBegin(u, INSERT_VALUES, SCATTER_FORWARD);
474 CHKERR VecGhostUpdateEnd(u, INSERT_VALUES, SCATTER_FORWARD);
475 CHKERR VecGhostUpdateBegin(u_t, INSERT_VALUES, SCATTER_FORWARD);
476 CHKERR VecGhostUpdateEnd(u_t, INSERT_VALUES, SCATTER_FORWARD);
477 CHKERR VecGhostUpdateBegin(u_tt, INSERT_VALUES, SCATTER_FORWARD);
478 CHKERR VecGhostUpdateEnd(u_tt, INSERT_VALUES, SCATTER_FORWARD);
479 CHKERR ts_ctx->mField.getInterface<VecManager>()->setLocalGhostVector(
480 ts_ctx->problemName, COL, u, INSERT_VALUES, SCATTER_REVERSE);
481 if (ts_ctx->zeroMatrix) {
482 CHKERR MatZeroEntries(B);
483 }
484 int step;
485#if PETSC_VERSION_GE(3, 8, 0)
486 CHKERR TSGetStepNumber(ts, &step);
487#else
488 CHKERR TSGetTimeStepNumber(ts, &step);
489#endif
490
491 ts_ctx->matAssembleSwitch =
492 boost::movelib::make_unique<bool>(ts_ctx->zeroMatrix);
493
494 auto set = [&](auto &fe) {
495 fe.ts_u = u;
496 fe.ts_u_t = u_t;
497 fe.ts_u_tt = u_tt;
498 fe.ts_A = A;
499 fe.ts_B = B;
500 fe.ts_t = t;
501 fe.ts_a = a;
502 fe.ts_aa = aa;
503 fe.ts_step = step;
504
505 fe.ts_ctx = TSMethod::CTX_TSSETIJACOBIAN;
506 fe.snes_ctx = SnesMethod::CTX_SNESSETJACOBIAN;
507 fe.ksp_ctx = KspMethod::CTX_OPERATORS;
508 fe.data_ctx = PetscData::CtxSetA | PetscData::CtxSetB | PetscData::CtxSetX |
509 PetscData::CtxSetX_T | PetscData::CtxSetX_TT |
510 PetscData::CtxSetTime;
511 fe.ts = ts;
512 CHK_THROW_MESSAGE(TSGetTimeStep(ts, &fe.ts_dt), "get time step failed");
513 };
514
515 auto unset = [&](auto &fe) {
516 fe.ts_ctx = TSMethod::CTX_TSNONE;
517 fe.snes_ctx = SnesMethod::CTX_SNESNONE;
518 fe.ksp_ctx = KspMethod::CTX_KSPNONE;
519 fe.data_ctx = PetscData::CtxSetNone;
520 };
521
522 // preproces
523 for (auto &bit : ts_ctx->preProcessIJacobian) {
524 bit->matAssembleSwitch = boost::move(ts_ctx->matAssembleSwitch);
525 set(*bit);
526 CHKERR ts_ctx->mField.problem_basic_method_preProcess(ts_ctx->problemName,
527 *bit);
528 unset(*bit);
529 ts_ctx->matAssembleSwitch = boost::move(bit->matAssembleSwitch);
530 }
531
532 auto cache_ptr = boost::make_shared<CacheTuple>();
533 CHKERR ts_ctx->mField.cache_problem_entities(ts_ctx->problemName, cache_ptr);
534
535 for (auto &lit : ts_ctx->loopsIJacobian) {
536 lit.second->matAssembleSwitch = boost::move(ts_ctx->matAssembleSwitch);
537 set(*lit.second);
538 CHKERR ts_ctx->mField.loop_finite_elements(ts_ctx->problemName, lit.first,
539 *(lit.second), nullptr,
540 ts_ctx->bH, cache_ptr);
541 unset(*lit.second);
542 ts_ctx->matAssembleSwitch = boost::move(lit.second->matAssembleSwitch);
543 }
544
545 // post process
546 for (auto &bit : ts_ctx->postProcessIJacobian) {
547 bit->matAssembleSwitch = boost::move(ts_ctx->matAssembleSwitch);
548 set(*bit);
549 CHKERR ts_ctx->mField.problem_basic_method_postProcess(ts_ctx->problemName,
550 *bit);
551 unset(*bit);
552 ts_ctx->matAssembleSwitch = boost::move(bit->matAssembleSwitch);
553 }
554
555 if (ts_ctx->matAssembleSwitch) {
556 CHKERR MatAssemblyBegin(B, MAT_FINAL_ASSEMBLY);
557 CHKERR MatAssemblyEnd(B, MAT_FINAL_ASSEMBLY);
558 }
559 PetscLogEventEnd(ts_ctx->MOFEM_EVENT_TsCtxI2Function, 0, 0, 0, 0);
560 MoFEMFunctionReturn(0);
561}
562
563PetscErrorCode TsSetI2Function(TS ts, PetscReal t, Vec u, Vec u_t, Vec u_tt,
564 Vec F, void *ctx) {
565 MoFEMFunctionBegin;
566 TsCtx *ts_ctx = static_cast<TsCtx *>(ctx);
567 PetscLogEventBegin(ts_ctx->MOFEM_EVENT_TsCtxIFunction, 0, 0, 0, 0);
568 CHKERR VecGhostUpdateBegin(u, INSERT_VALUES, SCATTER_FORWARD);
569 CHKERR VecGhostUpdateEnd(u, INSERT_VALUES, SCATTER_FORWARD);
570 CHKERR VecGhostUpdateBegin(u_t, INSERT_VALUES, SCATTER_FORWARD);
571 CHKERR VecGhostUpdateEnd(u_t, INSERT_VALUES, SCATTER_FORWARD);
572 CHKERR VecGhostUpdateBegin(u_tt, INSERT_VALUES, SCATTER_FORWARD);
573 CHKERR VecGhostUpdateEnd(u_tt, INSERT_VALUES, SCATTER_FORWARD);
574 CHKERR ts_ctx->mField.getInterface<VecManager>()->setLocalGhostVector(
575 ts_ctx->problemName, COL, u, INSERT_VALUES, SCATTER_REVERSE);
576
577 auto zero_ghost_vec = [](Vec g) {
578 MoFEMFunctionBegin;
579 Vec l;
580 CHKERR VecGhostGetLocalForm(g, &l);
581 double *a;
582 CHKERR VecGetArray(l, &a);
583 int s;
584 CHKERR VecGetLocalSize(l, &s);
585 for (int i = 0; i != s; ++i)
586 a[i] = 0;
587 CHKERR VecRestoreArray(l, &a);
588 CHKERR VecGhostRestoreLocalForm(g, &l);
589 MoFEMFunctionReturn(0);
590 };
591 CHKERR zero_ghost_vec(F);
592
593 ts_ctx->vecAssembleSwitch = boost::movelib::make_unique<bool>(true);
594
595 int step;
596#if PETSC_VERSION_GE(3, 8, 0)
597 CHKERR TSGetStepNumber(ts, &step);
598#else
599 CHKERR TSGetTimeStepNumber(ts, &step);
600#endif
601
602 auto set = [&](auto &fe) {
603 fe.ts_u = u;
604 fe.ts_u_t = u_t;
605 fe.ts_u_tt = u_tt;
606 fe.ts_F = F;
607 fe.ts_t = t;
608 fe.ts_step = step;
609 fe.ts_ctx = TSMethod::CTX_TSSETIFUNCTION;
610 fe.snes_ctx = SnesMethod::CTX_SNESSETFUNCTION;
611 fe.ksp_ctx = KspMethod::CTX_SETFUNCTION;
612 fe.data_ctx = PetscData::CtxSetF | PetscData::CtxSetX |
613 PetscData::CtxSetX_T | PetscData::CtxSetX_TT |
614 PetscData::CtxSetTime;
615 fe.ts = ts;
616 CHK_THROW_MESSAGE(TSGetTimeStep(ts, &fe.ts_dt), "get time step failed");
617 };
618
619 auto unset = [&](auto &fe) {
620 fe.ts_ctx = TSMethod::CTX_TSNONE;
621 fe.snes_ctx = SnesMethod::CTX_SNESNONE;
622 fe.ksp_ctx = KspMethod::CTX_KSPNONE;
623 fe.data_ctx = PetscData::CtxSetNone;
624 };
625
626 // preprocess
627 for (auto &bit : ts_ctx->preProcessIFunction) {
628 bit->vecAssembleSwitch = boost::move(ts_ctx->vecAssembleSwitch);
629 set(*bit);
630 CHKERR ts_ctx->mField.problem_basic_method_preProcess(ts_ctx->problemName,
631 *bit);
632 unset(*bit);
633 ts_ctx->vecAssembleSwitch = boost::move(bit->vecAssembleSwitch);
634 }
635
636 auto cache_ptr = boost::make_shared<CacheTuple>();
637 CHKERR ts_ctx->mField.cache_problem_entities(ts_ctx->problemName, cache_ptr);
638
639 // fe loops
640 for (auto &lit : ts_ctx->loopsIFunction) {
641 lit.second->vecAssembleSwitch = boost::move(ts_ctx->vecAssembleSwitch);
642 set(*lit.second);
643 CHKERR ts_ctx->mField.loop_finite_elements(ts_ctx->problemName, lit.first,
644 *(lit.second), nullptr,
645 ts_ctx->bH, cache_ptr);
646 unset(*lit.second);
647 ts_ctx->vecAssembleSwitch = boost::move(lit.second->vecAssembleSwitch);
648 }
649
650 // post process
651 for (auto &bit : ts_ctx->postProcessIFunction) {
652 bit->vecAssembleSwitch = boost::move(ts_ctx->vecAssembleSwitch);
653 set(*bit);
654 CHKERR ts_ctx->mField.problem_basic_method_postProcess(ts_ctx->problemName,
655 *bit);
656 unset(*bit);
657 ts_ctx->vecAssembleSwitch = boost::move(bit->vecAssembleSwitch);
658 }
659
660 if (*ts_ctx->vecAssembleSwitch) {
661 CHKERR VecGhostUpdateBegin(F, ADD_VALUES, SCATTER_REVERSE);
662 CHKERR VecGhostUpdateEnd(F, ADD_VALUES, SCATTER_REVERSE);
663 CHKERR VecAssemblyBegin(F);
664 CHKERR VecAssemblyEnd(F);
665 }
666
667 PetscLogEventEnd(ts_ctx->MOFEM_EVENT_TsCtxIFunction, 0, 0, 0, 0);
668 MoFEMFunctionReturn(0);
669}
670
671} // namespace MoFEM
a
constexpr double a
Definition: approx_sphere.cpp:30
COL
@ COL
Definition: definitions.h:123
CHK_THROW_MESSAGE
#define CHK_THROW_MESSAGE(err, msg)
Check and throw MoFEM exception.
Definition: definitions.h:595
MoFEMFunctionReturnHot
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:447
MoFEMFunctionBegin
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:346
MoFEMFunctionReturn
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:416
CHKERR
#define CHKERR
Inline error check.
Definition: definitions.h:535
MoFEMFunctionBeginHot
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:440
MoFEM::CoreInterface::problem_basic_method_postProcess
virtual MoFEMErrorCode problem_basic_method_postProcess(const Problem *problem_ptr, BasicMethod &method, int verb=DEFAULT_VERBOSITY)=0
Set data for BasicMethod.
MoFEM::CoreInterface::loop_finite_elements
virtual MoFEMErrorCode loop_finite_elements(const std::string problem_name, const std::string &fe_name, FEMethod &method, boost::shared_ptr< NumeredEntFiniteElement_multiIndex > fe_ptr=nullptr, MoFEMTypes bh=MF_EXIST, CacheTupleWeakPtr cache_ptr=CacheTupleSharedPtr(), int verb=DEFAULT_VERBOSITY)=0
Make a loop over finite elements.
bit
auto bit
set bit
Definition: hanging_node_approx.cpp:75
i
FTensor::Index< 'i', SPACE_DIM > i
Definition: hcurl_divergence_operator_2d.cpp:27
ts_ctx
MoFEM::TsCtx * ts_ctx
Definition: level_set.cpp:1876
l
FTensor::Index< 'l', 3 > l
Definition: matrix_function.cpp:21
MoFEM::Exceptions::MoFEMErrorCode
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:56
MoFEM
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:10
MoFEM::TsSetIJacobian
PetscErrorCode TsSetIJacobian(TS ts, PetscReal t, Vec u, Vec u_t, PetscReal a, Mat A, Mat B, void *ctx)
Set function evaluating jacobina in TS solver.
Definition: TsCtx.cpp:129
MoFEM::TsMonitorSet
PetscErrorCode TsMonitorSet(TS ts, PetscInt step, PetscReal t, Vec u, void *ctx)
Set monitor for TS solver.
Definition: TsCtx.cpp:219
MoFEM::TsSetI2Function
PetscErrorCode TsSetI2Function(TS ts, PetscReal t, Vec u, Vec u_t, Vec u_tt, Vec F, void *ctx)
Calculation the right hand side for second order PDE in time.
Definition: TsCtx.cpp:563
MoFEM::TsSetRHSFunction
PetscErrorCode TsSetRHSFunction(TS ts, PetscReal t, Vec u, Vec F, void *ctx)
TS solver function.
Definition: TsCtx.cpp:278
MoFEM::TsSetIFunction
PetscErrorCode TsSetIFunction(TS ts, PetscReal t, Vec u, Vec u_t, Vec F, void *ctx)
Set IFunction for TS solver.
Definition: TsCtx.cpp:25
MoFEM::TsSetI2Jacobian
PetscErrorCode TsSetI2Jacobian(TS ts, PetscReal t, Vec u, Vec u_t, Vec u_tt, PetscReal a, PetscReal aa, Mat A, Mat B, void *ctx)
Calculation Jaconian for second order PDE in time.
Definition: TsCtx.cpp:466
MoFEM::TsSetRHSJacobian
PetscErrorCode TsSetRHSJacobian(TS ts, PetscReal t, Vec u, Mat A, Mat B, void *ctx)
TS solver function.
Definition: TsCtx.cpp:377
A
constexpr AssemblyType A
Definition: operators_tests.cpp:30
t
constexpr double t
plate stiffness
Definition: plate.cpp:59
g
constexpr double g
Definition: shallow_wave.cpp:63
MoFEM::CoreInterface::cache_problem_entities
virtual MoFEMErrorCode cache_problem_entities(const std::string prb_name, CacheTupleWeakPtr cache_ptr)=0
Cache variables.
MoFEM::CoreInterface::problem_basic_method_preProcess
virtual MoFEMErrorCode problem_basic_method_preProcess(const Problem *problem_ptr, BasicMethod &method, int verb=DEFAULT_VERBOSITY)=0
Set data for BasicMethod.
MoFEM::PetscData::CtxSetA
static constexpr Switches CtxSetA
Definition: LoopMethods.hpp:37
MoFEM::PetscData::CtxSetX
static constexpr Switches CtxSetX
Definition: LoopMethods.hpp:39
MoFEM::PetscData::CtxSetX_TT
static constexpr Switches CtxSetX_TT
Definition: LoopMethods.hpp:41
MoFEM::PetscData::CtxSetNone
static constexpr Switches CtxSetNone
Definition: LoopMethods.hpp:35
MoFEM::PetscData::CtxSetF
static constexpr Switches CtxSetF
Definition: LoopMethods.hpp:36
MoFEM::PetscData::CtxSetX_T
static constexpr Switches CtxSetX_T
Definition: LoopMethods.hpp:40
MoFEM::PetscData::CtxSetB
static constexpr Switches CtxSetB
Definition: LoopMethods.hpp:38
MoFEM::PetscData::CtxSetTime
static constexpr Switches CtxSetTime
Definition: LoopMethods.hpp:42
MoFEM::TsCtx
Interface for Time Stepping (TS) solver.
Definition: TsCtx.hpp:15
MoFEM::TsCtx::loopsIJacobian
FEMethodsSequence loopsIJacobian
Definition: TsCtx.hpp:27
MoFEM::TsCtx::zeroMatrix
bool zeroMatrix
Definition: TsCtx.hpp:43
MoFEM::TsCtx::bH
MoFEMTypes bH
If set to MF_EXIST check if element exist.
Definition: TsCtx.hpp:21
MoFEM::TsCtx::MOFEM_EVENT_TsCtxRHSJacobian
PetscLogEvent MOFEM_EVENT_TsCtxRHSJacobian
Definition: TsCtx.hpp:219
MoFEM::TsCtx::loopsRHSFunction
FEMethodsSequence loopsRHSFunction
Definition: TsCtx.hpp:31
MoFEM::TsCtx::MOFEM_EVENT_TsCtxMonitor
PetscLogEvent MOFEM_EVENT_TsCtxMonitor
Definition: TsCtx.hpp:222
MoFEM::TsCtx::postProcessRHSFunction
BasicMethodsSequence postProcessRHSFunction
Definition: TsCtx.hpp:41
MoFEM::TsCtx::preProcessMonitor
BasicMethodsSequence preProcessMonitor
Definition: TsCtx.hpp:36
MoFEM::TsCtx::MOFEM_EVENT_TsCtxRHSFunction
PetscLogEvent MOFEM_EVENT_TsCtxRHSFunction
Definition: TsCtx.hpp:218
MoFEM::TsCtx::preProcessIJacobian
BasicMethodsSequence preProcessIJacobian
Definition: TsCtx.hpp:32
MoFEM::TsCtx::postProcessIFunction
BasicMethodsSequence postProcessIFunction
Definition: TsCtx.hpp:35
MoFEM::TsCtx::loopsIFunction
FEMethodsSequence loopsIFunction
Definition: TsCtx.hpp:28
MoFEM::TsCtx::loopsRHSJacobian
FEMethodsSequence loopsRHSJacobian
Definition: TsCtx.hpp:30
MoFEM::TsCtx::preProcessIFunction
BasicMethodsSequence preProcessIFunction
Definition: TsCtx.hpp:34
MoFEM::TsCtx::loopsMonitor
FEMethodsSequence loopsMonitor
Definition: TsCtx.hpp:29
MoFEM::TsCtx::postProcessMonitor
BasicMethodsSequence postProcessMonitor
Definition: TsCtx.hpp:37
MoFEM::TsCtx::preProcessRHSJacobian
BasicMethodsSequence preProcessRHSJacobian
Definition: TsCtx.hpp:38
MoFEM::TsCtx::vecAssembleSwitch
boost::movelib::unique_ptr< bool > vecAssembleSwitch
Definition: TsCtx.hpp:226
MoFEM::TsCtx::matAssembleSwitch
boost::movelib::unique_ptr< bool > matAssembleSwitch
Definition: TsCtx.hpp:227
MoFEM::TsCtx::clearLoops
MoFEMErrorCode clearLoops()
Clear loops.
Definition: TsCtx.cpp:5
MoFEM::TsCtx::problemName
std::string problemName
Definition: TsCtx.hpp:20
MoFEM::TsCtx::mField
MoFEM::Interface & mField
Definition: TsCtx.hpp:17
MoFEM::TsCtx::preProcessRHSFunction
BasicMethodsSequence preProcessRHSFunction
Definition: TsCtx.hpp:39
MoFEM::TsCtx::postProcessIJacobian
BasicMethodsSequence postProcessIJacobian
Definition: TsCtx.hpp:33
MoFEM::TsCtx::MOFEM_EVENT_TsCtxIFunction
PetscLogEvent MOFEM_EVENT_TsCtxIFunction
Definition: TsCtx.hpp:220
MoFEM::TsCtx::postProcessRHSJacobian
BasicMethodsSequence postProcessRHSJacobian
Definition: TsCtx.hpp:40
MoFEM::TsCtx::MOFEM_EVENT_TsCtxI2Function
PetscLogEvent MOFEM_EVENT_TsCtxI2Function
Definition: TsCtx.hpp:223
MoFEM::UnknownInterface::getInterface
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface refernce to pointer of interface.
Definition: UnknownInterface.hpp:93
MoFEM::VecManager
Vector manager is used to create vectors \mofem_vectors.
Definition: VecManager.hpp:23
Generated by Doxygen 1.9.5 and hosted at University of Glasgow