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