v0.13.1
TsCtx.cpp
Go to the documentation of this file.
1/* This file is part of MoFEM.
2 * MoFEM is free software: you can redistribute it and/or modify it under
3 * the terms of the GNU Lesser General Public License as published by the
4 * Free Software Foundation, either version 3 of the License, or (at your
5 * option) any later version.
6 *
7 * MoFEM is distributed in the hope that it will be useful, but WITHOUT
8 * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
9 * FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
10 * License for more details.
11 *
12 * You should have received a copy of the GNU Lesser General Public
13 * License along with MoFEM. If not, see <http://www.gnu.org/licenses/>. */
14
15namespace MoFEM {
16
19 loopsIJacobian.clear();
20 loopsIFunction.clear();
21 loopsMonitor.clear();
22 loopsRHSJacobian.clear();
23 loopsRHSFunction.clear();
24 preProcessIJacobian.clear();
26 preProcessIFunction.clear();
28 preProcessMonitor.clear();
29 postProcessMonitor.clear();
35}
36
37PetscErrorCode TsSetIFunction(TS ts, PetscReal t, Vec u, Vec u_t, Vec F,
38 void *ctx) {
40 TsCtx *ts_ctx = static_cast<TsCtx *>(ctx);
41 PetscLogEventBegin(ts_ctx->MOFEM_EVENT_TsCtxIFunction, 0, 0, 0, 0);
42 CHKERR VecGhostUpdateBegin(u, INSERT_VALUES, SCATTER_FORWARD);
43 CHKERR VecGhostUpdateEnd(u, INSERT_VALUES, SCATTER_FORWARD);
44 CHKERR VecGhostUpdateBegin(u_t, INSERT_VALUES, SCATTER_FORWARD);
45 CHKERR VecGhostUpdateEnd(u_t, INSERT_VALUES, SCATTER_FORWARD);
46 CHKERR ts_ctx->mField.getInterface<VecManager>()->setLocalGhostVector(
47 ts_ctx->problemName, COL, u, INSERT_VALUES, SCATTER_REVERSE);
48
49 auto zero_ghost_vec = [](Vec g) {
51 Vec l;
52 CHKERR VecGhostGetLocalForm(g, &l);
53 double *a;
54 CHKERR VecGetArray(l, &a);
55 int s;
56 CHKERR VecGetLocalSize(l, &s);
57 for (int i = 0; i != s; ++i)
58 a[i] = 0;
59 CHKERR VecRestoreArray(l, &a);
60 CHKERR VecGhostRestoreLocalForm(g, &l);
62 };
63 CHKERR zero_ghost_vec(F);
64
65 ts_ctx->vecAssembleSwitch = boost::movelib::make_unique<bool>(true);
66
67 int step;
68#if PETSC_VERSION_GE(3, 8, 0)
69 CHKERR TSGetStepNumber(ts, &step);
70#else
71 CHKERR TSGetTimeStepNumber(ts, &step);
72#endif
73
74 auto set = [&](auto &fe) {
75 fe.ts = ts;
76 fe.ts_u = u;
77 fe.ts_u_t = u_t;
78 fe.ts_F = F;
79 fe.ts_t = t;
80 fe.ts_step = step;
83 fe.ksp_ctx = KspMethod::CTX_SETFUNCTION;
86 };
87
88 auto unset = [&](auto &fe) {
89 fe.ts_ctx = TSMethod::CTX_TSNONE;
90 fe.snes_ctx = SnesMethod::CTX_SNESNONE;
91 fe.ksp_ctx = KspMethod::CTX_KSPNONE;
92 fe.data_ctx = PetscData::CtxSetNone;
93 };
94
95 // preprocess
96 for (auto &bit : ts_ctx->preProcessIFunction) {
97 bit->vecAssembleSwitch = boost::move(ts_ctx->vecAssembleSwitch);
98 set(*bit);
100 *bit);
101 unset(*bit);
102 ts_ctx->vecAssembleSwitch = boost::move(bit->vecAssembleSwitch);
103 }
104
105 auto cache_ptr = boost::make_shared<CacheTuple>();
106 CHKERR ts_ctx->mField.cache_problem_entities(ts_ctx->problemName, cache_ptr);
107
108 // fe loops
109 for (auto &lit : ts_ctx->loopsIFunction) {
110 lit.second->vecAssembleSwitch = boost::move(ts_ctx->vecAssembleSwitch);
111 set(*lit.second);
112 CHKERR ts_ctx->mField.loop_finite_elements(ts_ctx->problemName, lit.first,
113 *(lit.second), nullptr,
114 ts_ctx->bH, cache_ptr);
115 unset(*lit.second);
116 ts_ctx->vecAssembleSwitch = boost::move(lit.second->vecAssembleSwitch);
117 }
118
119 // post process
120 for (auto &bit : ts_ctx->postProcessIFunction) {
121 bit->vecAssembleSwitch = boost::move(ts_ctx->vecAssembleSwitch);
122 set(*bit);
124 *bit);
125 unset(*bit);
126 ts_ctx->vecAssembleSwitch = boost::move(bit->vecAssembleSwitch);
127 }
128
129 if (*ts_ctx->vecAssembleSwitch) {
130 CHKERR VecGhostUpdateBegin(F, ADD_VALUES, SCATTER_REVERSE);
131 CHKERR VecGhostUpdateEnd(F, ADD_VALUES, SCATTER_REVERSE);
132 CHKERR VecAssemblyBegin(F);
133 CHKERR VecAssemblyEnd(F);
134 }
135
136 PetscLogEventEnd(ts_ctx->MOFEM_EVENT_TsCtxIFunction, 0, 0, 0, 0);
138}
139
140PetscErrorCode TsSetIJacobian(TS ts, PetscReal t, Vec u, Vec u_t, PetscReal a,
141 Mat A, Mat B, void *ctx) {
143
144 TsCtx *ts_ctx = static_cast<TsCtx *>(ctx);
145 PetscLogEventBegin(ts_ctx->MOFEM_EVENT_TsCtxIFunction, 0, 0, 0, 0);
146 CHKERR VecGhostUpdateBegin(u, INSERT_VALUES, SCATTER_FORWARD);
147 CHKERR VecGhostUpdateEnd(u, INSERT_VALUES, SCATTER_FORWARD);
148 CHKERR VecGhostUpdateBegin(u_t, INSERT_VALUES, SCATTER_FORWARD);
149 CHKERR VecGhostUpdateEnd(u_t, INSERT_VALUES, SCATTER_FORWARD);
150 CHKERR ts_ctx->mField.getInterface<VecManager>()->setLocalGhostVector(
151 ts_ctx->problemName, COL, u, INSERT_VALUES, SCATTER_REVERSE);
152 if (ts_ctx->zeroMatrix) {
153 CHKERR MatZeroEntries(B);
154 }
155 int step;
156#if PETSC_VERSION_GE(3, 8, 0)
157 CHKERR TSGetStepNumber(ts, &step);
158#else
159 CHKERR TSGetTimeStepNumber(ts, &step);
160#endif
161
162 ts_ctx->matAssembleSwitch =
163 boost::movelib::make_unique<bool>(ts_ctx->zeroMatrix);
164
165 auto set = [&](auto &fe) {
166 fe.ts = ts;
167 fe.ts_u = u;
168 fe.ts_u_t = u_t;
169 fe.ts_A = A;
170 fe.ts_B = B;
171 fe.ts_t = t;
172 fe.ts_a = a;
173 fe.ts_step = step;
176 fe.ksp_ctx = KspMethod::CTX_OPERATORS;
179 };
180
181 auto unset = [&](auto &fe) {
182 fe.ts_ctx = TSMethod::CTX_TSNONE;
183 fe.snes_ctx = SnesMethod::CTX_SNESNONE;
184 fe.ksp_ctx = KspMethod::CTX_KSPNONE;
185 fe.data_ctx = PetscData::CtxSetNone;
186 };
187
188 // preproces
189 for (auto &bit : ts_ctx->preProcessIJacobian) {
190 bit->matAssembleSwitch = boost::move(ts_ctx->matAssembleSwitch);
191 set(*bit);
193 *bit);
194 unset(*bit);
195 ts_ctx->matAssembleSwitch = boost::move(bit->matAssembleSwitch);
196 }
197
198 auto cache_ptr = boost::make_shared<CacheTuple>();
199 CHKERR ts_ctx->mField.cache_problem_entities(ts_ctx->problemName, cache_ptr);
200
201 for (auto &lit : ts_ctx->loopsIJacobian) {
202 lit.second->matAssembleSwitch = boost::move(ts_ctx->matAssembleSwitch);
203 set(*lit.second);
204 CHKERR ts_ctx->mField.loop_finite_elements(ts_ctx->problemName, lit.first,
205 *(lit.second), nullptr,
206 ts_ctx->bH, cache_ptr);
207 unset(*lit.second);
208 ts_ctx->matAssembleSwitch = boost::move(lit.second->matAssembleSwitch);
209 }
210
211 // post process
212 for (auto &bit : ts_ctx->postProcessIJacobian) {
213 bit->matAssembleSwitch = boost::move(ts_ctx->matAssembleSwitch);
214 set(*bit);
216 *bit);
217 unset(*bit);
218 ts_ctx->matAssembleSwitch = boost::move(bit->matAssembleSwitch);
219 }
220
221 if (ts_ctx->matAssembleSwitch) {
222 CHKERR MatAssemblyBegin(B, MAT_FINAL_ASSEMBLY);
223 CHKERR MatAssemblyEnd(B, MAT_FINAL_ASSEMBLY);
224 }
225 PetscLogEventEnd(ts_ctx->MOFEM_EVENT_TsCtxIFunction, 0, 0, 0, 0);
227}
228
229PetscErrorCode TsMonitorSet(TS ts, PetscInt step, PetscReal t, Vec u,
230 void *ctx) {
232 TsCtx *ts_ctx = static_cast<TsCtx *>(ctx);
233 PetscLogEventBegin(ts_ctx->MOFEM_EVENT_TsCtxMonitor, 0, 0, 0, 0);
234 CHKERR VecGhostUpdateBegin(u, INSERT_VALUES, SCATTER_FORWARD);
235 CHKERR VecGhostUpdateEnd(u, INSERT_VALUES, SCATTER_FORWARD);
236 CHKERR ts_ctx->mField.getInterface<VecManager>()->setLocalGhostVector(
237 ts_ctx->problemName, COL, u, INSERT_VALUES, SCATTER_REVERSE);
238
239 auto set = [&](auto &fe) {
240 fe.ts = ts;
241 fe.ts_u = u;
242 fe.ts_t = t;
243 fe.ts_step = step;
244 fe.ts_F = PETSC_NULL;
246 fe.snes_ctx = SnesMethod::CTX_SNESNONE;
247 fe.ksp_ctx = KspMethod::CTX_KSPNONE;
249 };
250
251 auto unset = [&](auto &fe) {
252 fe.ts_ctx = TSMethod::CTX_TSNONE;
253 fe.data_ctx = PetscData::CtxSetNone;
254 };
255
256 // preproces
257 for (auto &bit : ts_ctx->preProcessMonitor) {
258 set(*bit);
260 *bit);
261 unset(*bit);
262 }
263
264 auto cache_ptr = boost::make_shared<CacheTuple>();
265 CHKERR ts_ctx->mField.cache_problem_entities(ts_ctx->problemName, cache_ptr);
266
267 for (auto &lit : ts_ctx->loopsMonitor) {
268 set(*lit.second);
269 CHKERR ts_ctx->mField.loop_finite_elements(ts_ctx->problemName, lit.first,
270 *(lit.second), nullptr,
271 ts_ctx->bH, cache_ptr);
272 unset(*lit.second);
273 }
274
275 // post process
276 for (auto &bit : ts_ctx->postProcessMonitor) {
277 set(*bit);
279 *bit);
280 unset(*bit);
281 }
282
283 PetscLogEventEnd(ts_ctx->MOFEM_EVENT_TsCtxMonitor, 0, 0, 0, 0);
285}
286
287PetscErrorCode TsSetRHSFunction(TS ts, PetscReal t, Vec u, Vec F, void *ctx) {
289 TsCtx *ts_ctx = static_cast<TsCtx *>(ctx);
290 PetscLogEventBegin(ts_ctx->MOFEM_EVENT_TsCtxRHSFunction, 0, 0, 0, 0);
291 CHKERR VecGhostUpdateBegin(u, INSERT_VALUES, SCATTER_FORWARD);
292 CHKERR VecGhostUpdateEnd(u, INSERT_VALUES, SCATTER_FORWARD);
293 CHKERR ts_ctx->mField.getInterface<VecManager>()->setLocalGhostVector(
294 ts_ctx->problemName, COL, u, INSERT_VALUES, SCATTER_REVERSE);
295
296 auto zero_ghost_vec = [](Vec g) {
298 Vec l;
299 CHKERR VecGhostGetLocalForm(g, &l);
300 double *a;
301 CHKERR VecGetArray(l, &a);
302 int s;
303 CHKERR VecGetLocalSize(l, &s);
304 for (int i = 0; i != s; ++i)
305 a[i] = 0;
306 CHKERR VecRestoreArray(l, &a);
307 CHKERR VecGhostRestoreLocalForm(g, &l);
309 };
310 CHKERR zero_ghost_vec(F);
311
312 ts_ctx->vecAssembleSwitch = boost::movelib::make_unique<bool>(true);
313
314 int step;
315#if PETSC_VERSION_GE(3, 8, 0)
316 CHKERR TSGetStepNumber(ts, &step);
317#else
318 CHKERR TSGetTimeStepNumber(ts, &step);
319#endif
320
321 auto set = [&](auto &fe) {
322 fe.ts_u = u;
323 fe.ts_F = F;
324 fe.ts_t = t;
325 fe.ts = ts;
326 fe.ts_step = step;
329 fe.ksp_ctx = KspMethod::CTX_SETFUNCTION;
330 fe.data_ctx =
332 };
333
334 auto unset = [&](auto &fe) {
335 fe.ts_ctx = TSMethod::CTX_TSNONE;
336 fe.snes_ctx = SnesMethod::CTX_SNESNONE;
337 fe.ksp_ctx = KspMethod::CTX_KSPNONE;
338 fe.data_ctx = PetscData::CtxSetNone;
339 };
340
341 for (auto &bit : ts_ctx->preProcessRHSJacobian) {
342 bit->vecAssembleSwitch = boost::move(ts_ctx->vecAssembleSwitch);
343 set(*bit);
345 *bit);
346 unset(*bit);
347 ts_ctx->vecAssembleSwitch = boost::move(bit->vecAssembleSwitch);
348 }
349
350 auto cache_ptr = boost::make_shared<CacheTuple>();
351 CHKERR ts_ctx->mField.cache_problem_entities(ts_ctx->problemName, cache_ptr);
352
353 // fe loops
354 for (auto &lit : ts_ctx->loopsRHSFunction) {
355 lit.second->vecAssembleSwitch = boost::move(ts_ctx->vecAssembleSwitch);
356 set(*lit.second);
357 CHKERR ts_ctx->mField.loop_finite_elements(ts_ctx->problemName, lit.first,
358 *(lit.second), nullptr,
359 ts_ctx->bH, cache_ptr);
360 unset(*lit.second);
361 ts_ctx->vecAssembleSwitch = boost::move(lit.second->vecAssembleSwitch);
362 }
363
364 // post process
365 for (auto &bit : ts_ctx->postProcessRHSJacobian) {
366 bit->vecAssembleSwitch = boost::move(ts_ctx->vecAssembleSwitch);
367 set(*bit);
369 *bit);
370 unset(*bit);
371 ts_ctx->vecAssembleSwitch = boost::move(bit->vecAssembleSwitch);
372 }
373
374 if (*ts_ctx->vecAssembleSwitch) {
375 CHKERR VecGhostUpdateBegin(F, ADD_VALUES, SCATTER_REVERSE);
376 CHKERR VecGhostUpdateEnd(F, ADD_VALUES, SCATTER_REVERSE);
377 CHKERR VecAssemblyBegin(F);
378 CHKERR VecAssemblyEnd(F);
379 }
380
381 PetscLogEventEnd(ts_ctx->MOFEM_EVENT_TsCtxRHSFunction, 0, 0, 0, 0);
383}
384
385PetscErrorCode TsSetRHSJacobian(TS ts, PetscReal t, Vec u, Mat A, Mat B,
386 void *ctx) {
388 TsCtx *ts_ctx = static_cast<TsCtx *>(ctx);
389 PetscLogEventBegin(ts_ctx->MOFEM_EVENT_TsCtxRHSJacobian, 0, 0, 0, 0);
390 CHKERR VecGhostUpdateBegin(u, INSERT_VALUES, SCATTER_FORWARD);
391 CHKERR VecGhostUpdateEnd(u, INSERT_VALUES, SCATTER_FORWARD);
392 CHKERR ts_ctx->mField.getInterface<VecManager>()->setLocalGhostVector(
393 ts_ctx->problemName, COL, u, INSERT_VALUES, SCATTER_REVERSE);
394
395 if (ts_ctx->zeroMatrix) {
396 CHKERR MatZeroEntries(B);
397 }
398
399 ts_ctx->matAssembleSwitch =
400 boost::movelib::make_unique<bool>(ts_ctx->zeroMatrix);
401
402 int step;
403#if PETSC_VERSION_GE(3, 8, 0)
404 CHKERR TSGetStepNumber(ts, &step);
405#else
406 CHKERR TSGetTimeStepNumber(ts, &step);
407#endif
408
409 auto set = [&](auto &fe) {
410 fe.ts_u = u;
411 fe.ts_A = A;
412 fe.ts_B = B;
413 fe.ts_t = t;
414 fe.ts_step = step;
417 fe.ksp_ctx = KspMethod::CTX_OPERATORS;
420 fe.ts = ts;
421 };
422
423 auto unset = [&](auto &fe) {
424 fe.ts_ctx = TSMethod::CTX_TSNONE;
425 fe.snes_ctx = SnesMethod::CTX_SNESNONE;
426 fe.ksp_ctx = KspMethod::CTX_KSPNONE;
427 fe.data_ctx = PetscData::CtxSetNone;
428 };
429
430 // preprocess
431 for (auto &bit : ts_ctx->preProcessRHSJacobian) {
432 bit->matAssembleSwitch = boost::move(ts_ctx->matAssembleSwitch);
433 set(*bit);
435 *bit);
436 unset(*bit);
437 ts_ctx->matAssembleSwitch = boost::move(bit->matAssembleSwitch);
438 }
439
440 auto cache_ptr = boost::make_shared<CacheTuple>();
441 CHKERR ts_ctx->mField.cache_problem_entities(ts_ctx->problemName, cache_ptr);
442
443 // fe loops
444 for (auto &lit : ts_ctx->loopsRHSJacobian) {
445 lit.second->matAssembleSwitch = boost::move(ts_ctx->matAssembleSwitch);
446 set(*lit.second);
447 CHKERR ts_ctx->mField.loop_finite_elements(ts_ctx->problemName, lit.first,
448 *(lit.second), nullptr,
449 ts_ctx->bH, cache_ptr);
450 unset(*lit.second);
451 ts_ctx->matAssembleSwitch = boost::move(lit.second->matAssembleSwitch);
452 }
453
454 // post process
455 for (auto &bit : ts_ctx->postProcessRHSJacobian) {
456 bit->matAssembleSwitch = boost::move(ts_ctx->matAssembleSwitch);
457 set(*bit);
459 *bit);
460 unset(*bit);
461 ts_ctx->matAssembleSwitch = boost::move(bit->matAssembleSwitch);
462 }
463
464 if (ts_ctx->matAssembleSwitch) {
465 CHKERR MatAssemblyBegin(B, MAT_FINAL_ASSEMBLY);
466 CHKERR MatAssemblyEnd(B, MAT_FINAL_ASSEMBLY);
467 }
468
469 PetscLogEventEnd(ts_ctx->MOFEM_EVENT_TsCtxRHSJacobian, 0, 0, 0, 0);
471}
472
473PetscErrorCode TsSetI2Jacobian(TS ts, PetscReal t, Vec u, Vec u_t, Vec u_tt,
474 PetscReal a, PetscReal aa, Mat A, Mat B,
475 void *ctx) {
477
478 TsCtx *ts_ctx = static_cast<TsCtx *>(ctx);
479 PetscLogEventBegin(ts_ctx->MOFEM_EVENT_TsCtxI2Function, 0, 0, 0, 0);
480 CHKERR VecGhostUpdateBegin(u, INSERT_VALUES, SCATTER_FORWARD);
481 CHKERR VecGhostUpdateEnd(u, INSERT_VALUES, SCATTER_FORWARD);
482 CHKERR VecGhostUpdateBegin(u_t, INSERT_VALUES, SCATTER_FORWARD);
483 CHKERR VecGhostUpdateEnd(u_t, INSERT_VALUES, SCATTER_FORWARD);
484 CHKERR VecGhostUpdateBegin(u_tt, INSERT_VALUES, SCATTER_FORWARD);
485 CHKERR VecGhostUpdateEnd(u_tt, INSERT_VALUES, SCATTER_FORWARD);
486 CHKERR ts_ctx->mField.getInterface<VecManager>()->setLocalGhostVector(
487 ts_ctx->problemName, COL, u, INSERT_VALUES, SCATTER_REVERSE);
488 if (ts_ctx->zeroMatrix) {
489 CHKERR MatZeroEntries(B);
490 }
491 int step;
492#if PETSC_VERSION_GE(3, 8, 0)
493 CHKERR TSGetStepNumber(ts, &step);
494#else
495 CHKERR TSGetTimeStepNumber(ts, &step);
496#endif
497
498 ts_ctx->matAssembleSwitch =
499 boost::movelib::make_unique<bool>(ts_ctx->zeroMatrix);
500
501 auto set = [&](auto &fe) {
502 fe.ts_u = u;
503 fe.ts_u_t = u_t;
504 fe.ts_u_tt = u_tt;
505 fe.ts_A = A;
506 fe.ts_B = B;
507 fe.ts_t = t;
508 fe.ts_a = a;
509 fe.ts_aa = aa;
510 fe.ts_step = step;
511
514 fe.ksp_ctx = KspMethod::CTX_OPERATORS;
518 fe.ts = ts;
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 // preproces
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 auto cache_ptr = boost::make_shared<CacheTuple>();
539 CHKERR ts_ctx->mField.cache_problem_entities(ts_ctx->problemName, cache_ptr);
540
541 for (auto &lit : ts_ctx->loopsIJacobian) {
542 lit.second->matAssembleSwitch = boost::move(ts_ctx->matAssembleSwitch);
543 set(*lit.second);
544 CHKERR ts_ctx->mField.loop_finite_elements(ts_ctx->problemName, lit.first,
545 *(lit.second), nullptr,
546 ts_ctx->bH, cache_ptr);
547 unset(*lit.second);
548 ts_ctx->matAssembleSwitch = boost::move(lit.second->matAssembleSwitch);
549 }
550
551 // post process
552 for (auto &bit : ts_ctx->postProcessIJacobian) {
553 bit->matAssembleSwitch = boost::move(ts_ctx->matAssembleSwitch);
554 set(*bit);
556 *bit);
557 unset(*bit);
558 ts_ctx->matAssembleSwitch = boost::move(bit->matAssembleSwitch);
559 }
560
561 if (ts_ctx->matAssembleSwitch) {
562 CHKERR MatAssemblyBegin(B, MAT_FINAL_ASSEMBLY);
563 CHKERR MatAssemblyEnd(B, MAT_FINAL_ASSEMBLY);
564 }
565 PetscLogEventEnd(ts_ctx->MOFEM_EVENT_TsCtxI2Function, 0, 0, 0, 0);
567}
568
569PetscErrorCode TsSetI2Function(TS ts, PetscReal t, Vec u, Vec u_t, Vec u_tt,
570 Vec F, void *ctx) {
572 TsCtx *ts_ctx = static_cast<TsCtx *>(ctx);
573 PetscLogEventBegin(ts_ctx->MOFEM_EVENT_TsCtxIFunction, 0, 0, 0, 0);
574 CHKERR VecGhostUpdateBegin(u, INSERT_VALUES, SCATTER_FORWARD);
575 CHKERR VecGhostUpdateEnd(u, INSERT_VALUES, SCATTER_FORWARD);
576 CHKERR VecGhostUpdateBegin(u_t, INSERT_VALUES, SCATTER_FORWARD);
577 CHKERR VecGhostUpdateEnd(u_t, INSERT_VALUES, SCATTER_FORWARD);
578 CHKERR VecGhostUpdateBegin(u_tt, INSERT_VALUES, SCATTER_FORWARD);
579 CHKERR VecGhostUpdateEnd(u_tt, INSERT_VALUES, SCATTER_FORWARD);
580 CHKERR ts_ctx->mField.getInterface<VecManager>()->setLocalGhostVector(
581 ts_ctx->problemName, COL, u, INSERT_VALUES, SCATTER_REVERSE);
582
583 auto zero_ghost_vec = [](Vec g) {
585 Vec l;
586 CHKERR VecGhostGetLocalForm(g, &l);
587 double *a;
588 CHKERR VecGetArray(l, &a);
589 int s;
590 CHKERR VecGetLocalSize(l, &s);
591 for (int i = 0; i != s; ++i)
592 a[i] = 0;
593 CHKERR VecRestoreArray(l, &a);
594 CHKERR VecGhostRestoreLocalForm(g, &l);
596 };
597 CHKERR zero_ghost_vec(F);
598
599 ts_ctx->vecAssembleSwitch = boost::movelib::make_unique<bool>(true);
600
601 int step;
602#if PETSC_VERSION_GE(3, 8, 0)
603 CHKERR TSGetStepNumber(ts, &step);
604#else
605 CHKERR TSGetTimeStepNumber(ts, &step);
606#endif
607
608 auto set = [&](auto &fe) {
609 fe.ts_u = u;
610 fe.ts_u_t = u_t;
611 fe.ts_u_tt = u_tt;
612 fe.ts_F = F;
613 fe.ts_t = t;
614 fe.ts_step = step;
617 fe.ksp_ctx = KspMethod::CTX_SETFUNCTION;
618 fe.data_ctx = PetscData::CtxSetF | PetscData::CtxSetX |
621 fe.ts = ts;
622 };
623
624 auto unset = [&](auto &fe) {
625 fe.ts_ctx = TSMethod::CTX_TSNONE;
626 fe.snes_ctx = SnesMethod::CTX_SNESNONE;
627 fe.ksp_ctx = KspMethod::CTX_KSPNONE;
628 fe.data_ctx = PetscData::CtxSetNone;
629 };
630
631 // preprocess
632 for (auto &bit : ts_ctx->preProcessIFunction) {
633 bit->vecAssembleSwitch = boost::move(ts_ctx->vecAssembleSwitch);
634 set(*bit);
636 *bit);
637 unset(*bit);
638 ts_ctx->vecAssembleSwitch = boost::move(bit->vecAssembleSwitch);
639 }
640
641 auto cache_ptr = boost::make_shared<CacheTuple>();
642 CHKERR ts_ctx->mField.cache_problem_entities(ts_ctx->problemName, cache_ptr);
643
644 // fe loops
645 for (auto &lit : ts_ctx->loopsIFunction) {
646 lit.second->vecAssembleSwitch = boost::move(ts_ctx->vecAssembleSwitch);
647 set(*lit.second);
648 CHKERR ts_ctx->mField.loop_finite_elements(ts_ctx->problemName, lit.first,
649 *(lit.second), nullptr,
650 ts_ctx->bH, cache_ptr);
651 unset(*lit.second);
652 ts_ctx->vecAssembleSwitch = boost::move(lit.second->vecAssembleSwitch);
653 }
654
655 // post process
656 for (auto &bit : ts_ctx->postProcessIFunction) {
657 bit->vecAssembleSwitch = boost::move(ts_ctx->vecAssembleSwitch);
658 set(*bit);
660 *bit);
661 unset(*bit);
662 ts_ctx->vecAssembleSwitch = boost::move(bit->vecAssembleSwitch);
663 }
664
665 if (*ts_ctx->vecAssembleSwitch) {
666 CHKERR VecGhostUpdateBegin(F, ADD_VALUES, SCATTER_REVERSE);
667 CHKERR VecGhostUpdateEnd(F, ADD_VALUES, SCATTER_REVERSE);
668 CHKERR VecAssemblyBegin(F);
669 CHKERR VecAssemblyEnd(F);
670 }
671
672 PetscLogEventEnd(ts_ctx->MOFEM_EVENT_TsCtxIFunction, 0, 0, 0, 0);
674}
675
676} // namespace MoFEM
constexpr double a
@ COL
Definition: definitions.h:136
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:460
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:359
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:429
#define CHKERR
Inline error check.
Definition: definitions.h:548
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:453
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
FTensor::Index< 'l', 3 > l
const FTensor::Tensor2< T, Dim, Dim > Vec
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:67
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:21
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:140
PetscErrorCode TsMonitorSet(TS ts, PetscInt step, PetscReal t, Vec u, void *ctx)
Set monitor for TS solver.
Definition: TsCtx.cpp:229
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:569
PetscErrorCode TsSetRHSFunction(TS ts, PetscReal t, Vec u, Vec F, void *ctx)
TS solver function.
Definition: TsCtx.cpp:287
PetscErrorCode TsSetIFunction(TS ts, PetscReal t, Vec u, Vec u_t, Vec F, void *ctx)
Set IFunction for TS solver.
Definition: TsCtx.cpp:37
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:473
PetscErrorCode TsSetRHSJacobian(TS ts, PetscReal t, Vec u, Mat A, Mat B, void *ctx)
TS solver function.
Definition: TsCtx.cpp:385
double A
constexpr double t
plate stiffness
Definition: plate.cpp:76
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:47
static constexpr Switches CtxSetX
Definition: LoopMethods.hpp:49
static constexpr Switches CtxSetX_TT
Definition: LoopMethods.hpp:51
static constexpr Switches CtxSetNone
Definition: LoopMethods.hpp:45
static constexpr Switches CtxSetF
Definition: LoopMethods.hpp:46
static constexpr Switches CtxSetX_T
Definition: LoopMethods.hpp:50
static constexpr Switches CtxSetB
Definition: LoopMethods.hpp:48
static constexpr Switches CtxSetTime
Definition: LoopMethods.hpp:52
Interface for Time Stepping (TS) solver.
Definition: TsCtx.hpp:27
FEMethodsSequence loopsIJacobian
Definition: TsCtx.hpp:39
bool zeroMatrix
Definition: TsCtx.hpp:55
MoFEMTypes bH
If set to MF_EXIST check if element exist.
Definition: TsCtx.hpp:33
PetscLogEvent MOFEM_EVENT_TsCtxRHSJacobian
Definition: TsCtx.hpp:231
FEMethodsSequence loopsRHSFunction
Definition: TsCtx.hpp:43
PetscLogEvent MOFEM_EVENT_TsCtxMonitor
Definition: TsCtx.hpp:234
BasicMethodsSequence postProcessRHSFunction
Definition: TsCtx.hpp:53
BasicMethodsSequence preProcessMonitor
Definition: TsCtx.hpp:48
PetscLogEvent MOFEM_EVENT_TsCtxRHSFunction
Definition: TsCtx.hpp:230
BasicMethodsSequence preProcessIJacobian
Definition: TsCtx.hpp:44
BasicMethodsSequence postProcessIFunction
Definition: TsCtx.hpp:47
FEMethodsSequence loopsIFunction
Definition: TsCtx.hpp:40
FEMethodsSequence loopsRHSJacobian
Definition: TsCtx.hpp:42
BasicMethodsSequence preProcessIFunction
Definition: TsCtx.hpp:46
FEMethodsSequence loopsMonitor
Definition: TsCtx.hpp:41
BasicMethodsSequence postProcessMonitor
Definition: TsCtx.hpp:49
BasicMethodsSequence preProcessRHSJacobian
Definition: TsCtx.hpp:50
boost::movelib::unique_ptr< bool > vecAssembleSwitch
Definition: TsCtx.hpp:238
boost::movelib::unique_ptr< bool > matAssembleSwitch
Definition: TsCtx.hpp:239
MoFEMErrorCode clearLoops()
Clear loops.
Definition: TsCtx.cpp:17
std::string problemName
Definition: TsCtx.hpp:32
MoFEM::Interface & mField
Definition: TsCtx.hpp:29
BasicMethodsSequence preProcessRHSFunction
Definition: TsCtx.hpp:51
BasicMethodsSequence postProcessIJacobian
Definition: TsCtx.hpp:45
PetscLogEvent MOFEM_EVENT_TsCtxIFunction
Definition: TsCtx.hpp:232
BasicMethodsSequence postProcessRHSJacobian
Definition: TsCtx.hpp:52
PetscLogEvent MOFEM_EVENT_TsCtxI2Function
Definition: TsCtx.hpp:235
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface refernce to pointer of interface.
Vector manager is used to create vectors \mofem_vectors.
Definition: VecManager.hpp:33