v0.9.0
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 
15 namespace MoFEM {
16 
17 PetscErrorCode TsSetIFunction(TS ts, PetscReal t, Vec u, Vec u_t, Vec F,
18  void *ctx) {
20  TsCtx *ts_ctx = static_cast<TsCtx *>(ctx);
21  PetscLogEventBegin(ts_ctx->MOFEM_EVENT_TsCtxIFunction, 0, 0, 0, 0);
22  CHKERR VecGhostUpdateBegin(u, INSERT_VALUES, SCATTER_FORWARD);
23  CHKERR VecGhostUpdateEnd(u, INSERT_VALUES, SCATTER_FORWARD);
24  CHKERR VecGhostUpdateBegin(u_t, INSERT_VALUES, SCATTER_FORWARD);
25  CHKERR VecGhostUpdateEnd(u_t, INSERT_VALUES, SCATTER_FORWARD);
26  CHKERR ts_ctx->mField.getInterface<VecManager>()->setLocalGhostVector(
27  ts_ctx->problemName, COL, u, INSERT_VALUES, SCATTER_REVERSE);
28 
29  auto zero_ghost_vec = [](Vec g) {
31  Vec l;
32  CHKERR VecGhostGetLocalForm(g, &l);
33  double *a;
34  CHKERR VecGetArray(l, &a);
35  int s;
36  CHKERR VecGetLocalSize(l, &s);
37  for (int i = 0; i != s; ++i)
38  a[i] = 0;
39  CHKERR VecRestoreArray(l, &a);
40  CHKERR VecGhostRestoreLocalForm(g, &l);
42  };
43  CHKERR zero_ghost_vec(F);
44 
45  ts_ctx->vecAssembleSwitch = boost::movelib::make_unique<bool>(true);
46 
47  int step;
48  CHKERR TSGetTimeStepNumber(ts, &step);
49  // preprocess
50  for (auto &bit : ts_ctx->preProcess_IFunction) {
51  bit->vecAssembleSwitch = boost::move(ts_ctx->vecAssembleSwitch);
52  bit->ts_u = u;
53  bit->ts_u_t = u_t;
54  bit->ts_F = F;
55  bit->ts_t = t;
56  bit->ts_step = step;
57  CHKERR bit->setTsCtx(TSMethod::CTX_TSSETIFUNCTION);
58  CHKERR bit->setTs(ts);
60  *bit);
61  CHKERR bit->setTsCtx(TSMethod::CTX_TSNONE);
62  ts_ctx->vecAssembleSwitch = boost::move(bit->vecAssembleSwitch);
63  }
64 
65  // fe loops
66  for (auto &lit : ts_ctx->loops_to_do_IFunction) {
67  lit.second->vecAssembleSwitch = boost::move(ts_ctx->vecAssembleSwitch);
68  lit.second->ts_u = u;
69  lit.second->ts_u_t = u_t;
70  lit.second->ts_F = F;
71  lit.second->ts_t = t;
72  lit.second->ts_step = step;
73  CHKERR lit.second->setTsCtx(TSMethod::CTX_TSSETIFUNCTION);
74  CHKERR lit.second->setTs(ts);
76  ts_ctx->problemName, lit.first, *(lit.second), nullptr, ts_ctx->bH);
77  CHKERR lit.second->setTsCtx(TSMethod::CTX_TSNONE);
78  ts_ctx->vecAssembleSwitch = boost::move(lit.second->vecAssembleSwitch);
79  }
80 
81  // post process
82  for (auto &bit : ts_ctx->postProcess_IFunction) {
83  bit->vecAssembleSwitch = boost::move(ts_ctx->vecAssembleSwitch);
84  bit->ts_u = u;
85  bit->ts_u_t = u_t;
86  bit->ts_F = F;
87  bit->ts_t = t;
88  bit->ts_step = step;
89  CHKERR bit->setTsCtx(TSMethod::CTX_TSSETIFUNCTION);
90  CHKERR bit->setTs(ts);
92  *bit);
93  CHKERR bit->setTsCtx(TSMethod::CTX_TSNONE);
94  ts_ctx->vecAssembleSwitch = boost::move(bit->vecAssembleSwitch);
95  }
96 
97  if (*ts_ctx->vecAssembleSwitch) {
98  CHKERR VecGhostUpdateBegin(F, ADD_VALUES, SCATTER_REVERSE);
99  CHKERR VecGhostUpdateEnd(F, ADD_VALUES, SCATTER_REVERSE);
100  CHKERR VecAssemblyBegin(F);
101  CHKERR VecAssemblyEnd(F);
102  }
103 
104  PetscLogEventEnd(ts_ctx->MOFEM_EVENT_TsCtxIFunction, 0, 0, 0, 0);
106 }
107 
108 PetscErrorCode TsSetIJacobian(TS ts, PetscReal t, Vec u, Vec u_t, PetscReal a,
109  Mat A, Mat B, void *ctx) {
111 
112  TsCtx *ts_ctx = static_cast<TsCtx *>(ctx);
113  PetscLogEventBegin(ts_ctx->MOFEM_EVENT_TsCtxIFunction, 0, 0, 0, 0);
114  CHKERR VecGhostUpdateBegin(u, INSERT_VALUES, SCATTER_FORWARD);
115  CHKERR VecGhostUpdateEnd(u, INSERT_VALUES, SCATTER_FORWARD);
116  CHKERR VecGhostUpdateBegin(u_t, INSERT_VALUES, SCATTER_FORWARD);
117  CHKERR VecGhostUpdateEnd(u_t, INSERT_VALUES, SCATTER_FORWARD);
118  CHKERR ts_ctx->mField.getInterface<VecManager>()->setLocalGhostVector(
119  ts_ctx->problemName, COL, u, INSERT_VALUES, SCATTER_REVERSE);
120  if (ts_ctx->zeroMatrix) {
121  CHKERR MatZeroEntries(B);
122  }
123  int step;
124  CHKERR TSGetTimeStepNumber(ts, &step);
125 
126  ts_ctx->matAssembleSwitch =
127  boost::movelib::make_unique<bool>(ts_ctx->zeroMatrix);
128 
129  // preproces
130  for (auto &bit : ts_ctx->preProcess_IJacobian) {
131  bit->matAssembleSwitch = boost::move(ts_ctx->matAssembleSwitch);
132  bit->ts_u = u;
133  bit->ts_u_t = u_t;
134  bit->ts_A = A;
135  bit->ts_B = B;
136  bit->ts_t = t;
137  bit->ts_a = a;
138  bit->ts_step = step;
139  CHKERR bit->setTsCtx(TSMethod::CTX_TSSETIJACOBIAN);
140  CHKERR bit->setTs(ts);
142  *bit);
143  CHKERR bit->setTsCtx(TSMethod::CTX_TSNONE);
144  ts_ctx->matAssembleSwitch = boost::move(bit->matAssembleSwitch);
145  }
146 
147  for (auto &lit : ts_ctx->loops_to_do_IJacobian) {
148  lit.second->matAssembleSwitch = boost::move(ts_ctx->matAssembleSwitch);
149  lit.second->ts_u = u;
150  lit.second->ts_u_t = u_t;
151  lit.second->ts_A = A;
152  lit.second->ts_B = B;
153  lit.second->ts_t = t;
154  lit.second->ts_a = a;
155  lit.second->ts_step = step;
156  CHKERR lit.second->setTsCtx(TSMethod::CTX_TSSETIJACOBIAN);
157  CHKERR lit.second->setTs(ts);
159  ts_ctx->problemName, lit.first, *(lit.second), nullptr, ts_ctx->bH);
160  CHKERR lit.second->setTsCtx(TSMethod::CTX_TSNONE);
161  ts_ctx->matAssembleSwitch = boost::move(lit.second->matAssembleSwitch);
162  }
163 
164  // post process
165  for (auto &bit : ts_ctx->postProcess_IJacobian) {
166  bit->matAssembleSwitch = boost::move(ts_ctx->matAssembleSwitch);
167  bit->ts_u = u;
168  bit->ts_u_t = u_t;
169  bit->ts_A = A;
170  bit->ts_B = B;
171  bit->ts_t = t;
172  bit->ts_a = a;
173  bit->ts_step = step;
174  CHKERR bit->setTsCtx(TSMethod::CTX_TSSETIJACOBIAN);
175  CHKERR bit->setTs(ts);
177  *bit);
178  CHKERR bit->setTsCtx(TSMethod::CTX_TSNONE);
179  ts_ctx->matAssembleSwitch = boost::move(bit->matAssembleSwitch);
180  }
181 
182  if (ts_ctx->matAssembleSwitch) {
183  CHKERR MatAssemblyBegin(B, MAT_FINAL_ASSEMBLY);
184  CHKERR MatAssemblyEnd(B, MAT_FINAL_ASSEMBLY);
185  }
186  PetscLogEventEnd(ts_ctx->MOFEM_EVENT_TsCtxIFunction, 0, 0, 0, 0);
188 }
189 
190 PetscErrorCode TsMonitorSet(TS ts, PetscInt step, PetscReal t, Vec u,
191  void *ctx) {
193  TsCtx *ts_ctx = static_cast<TsCtx *>(ctx);
194  PetscLogEventBegin(ts_ctx->MOFEM_EVENT_TsCtxMonitor, 0, 0, 0, 0);
195  CHKERR VecGhostUpdateBegin(u, INSERT_VALUES, SCATTER_FORWARD);
196  CHKERR VecGhostUpdateEnd(u, INSERT_VALUES, SCATTER_FORWARD);
197  CHKERR ts_ctx->mField.getInterface<VecManager>()->setLocalGhostVector(
198  ts_ctx->problemName, COL, u, INSERT_VALUES, SCATTER_REVERSE);
199 
200  // preproces
201  for (auto &bit : ts_ctx->preProcess_Monitor) {
202  bit->ts_u = u;
203  bit->ts_t = t;
204  bit->ts_step = step;
205  bit->ts_F = PETSC_NULL;
206  CHKERR bit->setTsCtx(TSMethod::CTX_TSTSMONITORSET);
207  CHKERR bit->setTs(ts);
209  *bit);
210  CHKERR bit->setTsCtx(TSMethod::CTX_TSNONE);
211  }
212 
213  for (auto &lit : ts_ctx->loops_to_do_Monitor) {
214  lit.second->ts_u = u;
215  lit.second->ts_t = t;
216  lit.second->ts_step = step;
217  lit.second->ts_F = PETSC_NULL;
218  CHKERR lit.second->setTsCtx(TSMethod::CTX_TSTSMONITORSET);
219  CHKERR lit.second->setTs(ts);
221  ts_ctx->problemName, lit.first, *(lit.second), nullptr, ts_ctx->bH);
222  CHKERR lit.second->setTsCtx(TSMethod::CTX_TSNONE);
223  }
224 
225  // post process
226  for (auto &bit : ts_ctx->postProcess_Monitor) {
227  bit->ts_u = u;
228  bit->ts_t = t;
229  bit->ts_step = step;
230  bit->ts_F = PETSC_NULL;
231  CHKERR bit->setTsCtx(TSMethod::CTX_TSTSMONITORSET);
232  CHKERR bit->setTs(ts);
234  *bit);
235  CHKERR bit->setTsCtx(TSMethod::CTX_TSNONE);
236  }
237  PetscLogEventEnd(ts_ctx->MOFEM_EVENT_TsCtxMonitor, 0, 0, 0, 0);
239 }
240 
241 PetscErrorCode TsSetRHSFunction(TS ts, PetscReal t, Vec u, Vec F, void *ctx) {
243  TsCtx *ts_ctx = static_cast<TsCtx *>(ctx);
244  PetscLogEventBegin(ts_ctx->MOFEM_EVENT_TsCtxRHSFunction, 0, 0, 0, 0);
245  CHKERR VecGhostUpdateBegin(u, INSERT_VALUES, SCATTER_FORWARD);
246  CHKERR VecGhostUpdateEnd(u, INSERT_VALUES, SCATTER_FORWARD);
247  CHKERR ts_ctx->mField.getInterface<VecManager>()->setLocalGhostVector(
248  ts_ctx->problemName, COL, u, INSERT_VALUES, SCATTER_REVERSE);
249 
250  auto zero_ghost_vec = [](Vec g) {
252  Vec l;
253  CHKERR VecGhostGetLocalForm(g, &l);
254  double *a;
255  CHKERR VecGetArray(l, &a);
256  int s;
257  CHKERR VecGetLocalSize(l, &s);
258  for (int i = 0; i != s; ++i)
259  a[i] = 0;
260  CHKERR VecRestoreArray(l, &a);
261  CHKERR VecGhostRestoreLocalForm(g, &l);
263  };
264  CHKERR zero_ghost_vec(F);
265 
266  ts_ctx->vecAssembleSwitch = boost::movelib::make_unique<bool>(true);
267 
268  int step;
269  CHKERR TSGetTimeStepNumber(ts, &step);
270  for (auto &bit : ts_ctx->preProcess_RHSJacobian) {
271  bit->vecAssembleSwitch = boost::move(ts_ctx->vecAssembleSwitch);
272  bit->ts_u = u;
273  bit->ts_F = F;
274  bit->ts_t = t;
275  bit->ts_step = step;
276  CHKERR bit->setTsCtx(TSMethod::CTX_TSSETRHSFUNCTION);
277  CHKERR bit->setTs(ts);
279  *bit);
280  CHKERR bit->setTsCtx(TSMethod::CTX_TSNONE);
281  ts_ctx->vecAssembleSwitch = boost::move(bit->vecAssembleSwitch);
282  }
283 
284  // fe loops
285  for (auto &lit : ts_ctx->loops_to_do_RHSFunction) {
286  lit.second->vecAssembleSwitch = boost::move(ts_ctx->vecAssembleSwitch);
287  lit.second->ts_u = u;
288  lit.second->ts_F = F;
289  lit.second->ts_t = t;
290  lit.second->ts_step = step;
291  CHKERR lit.second->setTsCtx(TSMethod::CTX_TSSETRHSFUNCTION);
292  CHKERR lit.second->setTs(ts);
294  ts_ctx->problemName, lit.first, *(lit.second), nullptr, ts_ctx->bH);
295  CHKERR lit.second->setTsCtx(TSMethod::CTX_TSNONE);
296  ts_ctx->vecAssembleSwitch = boost::move(lit.second->vecAssembleSwitch);
297  }
298 
299  // post process
300  for (auto &bit : ts_ctx->postProcess_RHSJacobian) {
301  bit->vecAssembleSwitch = boost::move(ts_ctx->vecAssembleSwitch);
302  bit->ts_u = u;
303  bit->ts_F = F;
304  bit->ts_t = t;
305  bit->ts_step = step;
306  CHKERR bit->setTsCtx(TSMethod::CTX_TSSETRHSFUNCTION);
307  CHKERR bit->setTs(ts);
309  *bit);
310  CHKERR bit->setTsCtx(TSMethod::CTX_TSNONE);
311  ts_ctx->vecAssembleSwitch = boost::move(bit->vecAssembleSwitch);
312  }
313 
314  if (*ts_ctx->vecAssembleSwitch) {
315  CHKERR VecGhostUpdateBegin(F, ADD_VALUES, SCATTER_REVERSE);
316  CHKERR VecGhostUpdateEnd(F, ADD_VALUES, SCATTER_REVERSE);
317  CHKERR VecAssemblyBegin(F);
318  CHKERR VecAssemblyEnd(F);
319  }
320 
321  PetscLogEventEnd(ts_ctx->MOFEM_EVENT_TsCtxRHSFunction, 0, 0, 0, 0);
323 }
324 
325 PetscErrorCode TsSetRHSJacobian(TS ts, PetscReal t, Vec u, Mat A, Mat B,
326  void *ctx) {
328  TsCtx *ts_ctx = static_cast<TsCtx *>(ctx);
329  PetscLogEventBegin(ts_ctx->MOFEM_EVENT_TsCtxRHSJacobian, 0, 0, 0, 0);
330  CHKERR VecGhostUpdateBegin(u, INSERT_VALUES, SCATTER_FORWARD);
331  CHKERR VecGhostUpdateEnd(u, INSERT_VALUES, SCATTER_FORWARD);
332  CHKERR ts_ctx->mField.getInterface<VecManager>()->setLocalGhostVector(
333  ts_ctx->problemName, COL, u, INSERT_VALUES, SCATTER_REVERSE);
334 
335  if (ts_ctx->zeroMatrix) {
336  CHKERR MatZeroEntries(B);
337  }
338 
339  ts_ctx->matAssembleSwitch =
340  boost::movelib::make_unique<bool>(ts_ctx->zeroMatrix);
341 
342  int step;
343  CHKERR TSGetTimeStepNumber(ts, &step);
344  // preprocess
345  for (auto &bit : ts_ctx->preProcess_RHSJacobian) {
346  bit->matAssembleSwitch = boost::move(ts_ctx->matAssembleSwitch);
347  bit->ts_u = u;
348  bit->ts_A = A;
349  bit->ts_B = B;
350  bit->ts_t = t;
351  bit->ts_step = step;
352  CHKERR bit->setTsCtx(TSMethod::CTX_TSSETRHSJACOBIAN);
353  CHKERR bit->setTs(ts);
355  *bit);
356  CHKERR bit->setTsCtx(TSMethod::CTX_TSNONE);
357  ts_ctx->matAssembleSwitch = boost::move(bit->matAssembleSwitch);
358  }
359 
360  // fe loops
361  for (auto &lit : ts_ctx->loops_to_do_RHSJacobian) {
362  lit.second->matAssembleSwitch = boost::move(ts_ctx->matAssembleSwitch);
363  lit.second->ts_u = u;
364  lit.second->ts_A = A;
365  lit.second->ts_B = B;
366  lit.second->ts_t = t;
367  lit.second->ts_step = step;
368  CHKERR lit.second->setTsCtx(TSMethod::CTX_TSSETRHSJACOBIAN);
369  CHKERR lit.second->setTs(ts);
371  ts_ctx->problemName, lit.first, *(lit.second), nullptr, ts_ctx->bH);
372  CHKERR lit.second->setTsCtx(TSMethod::CTX_TSNONE);
373  ts_ctx->matAssembleSwitch = boost::move(lit.second->matAssembleSwitch);
374  }
375 
376  // post process
377  for (auto &bit : ts_ctx->postProcess_RHSJacobian) {
378  bit->matAssembleSwitch = boost::move(ts_ctx->matAssembleSwitch);
379  bit->ts_u = u;
380  bit->ts_A = A;
381  bit->ts_B = B;
382  bit->ts_t = t;
383  bit->ts_step = step;
384  CHKERR bit->setTsCtx(TSMethod::CTX_TSSETRHSJACOBIAN);
385  CHKERR bit->setTs(ts);
387  *bit);
388  CHKERR bit->setTsCtx(TSMethod::CTX_TSNONE);
389  ts_ctx->matAssembleSwitch = boost::move(bit->matAssembleSwitch);
390  }
391 
392  if (ts_ctx->matAssembleSwitch) {
393  CHKERR MatAssemblyBegin(B, MAT_FINAL_ASSEMBLY);
394  CHKERR MatAssemblyEnd(B, MAT_FINAL_ASSEMBLY);
395  }
396 
397  PetscLogEventEnd(ts_ctx->MOFEM_EVENT_TsCtxRHSJacobian, 0, 0, 0, 0);
399 }
400 
401 PetscErrorCode TsSetI2Jacobian(TS ts, PetscReal t, Vec u, Vec u_t, Vec u_tt,
402  PetscReal v, PetscReal a, Mat A, Mat B,
403  void *ctx) {
405 
406  TsCtx *ts_ctx = static_cast<TsCtx *>(ctx);
407  PetscLogEventBegin(ts_ctx->MOFEM_EVENT_TsCtxI2Function, 0, 0, 0, 0);
408  CHKERR VecGhostUpdateBegin(u, INSERT_VALUES, SCATTER_FORWARD);
409  CHKERR VecGhostUpdateEnd(u, INSERT_VALUES, SCATTER_FORWARD);
410  CHKERR VecGhostUpdateBegin(u_t, INSERT_VALUES, SCATTER_FORWARD);
411  CHKERR VecGhostUpdateEnd(u_t, INSERT_VALUES, SCATTER_FORWARD);
412  CHKERR VecGhostUpdateBegin(u_tt, INSERT_VALUES, SCATTER_FORWARD);
413  CHKERR VecGhostUpdateEnd(u_tt, INSERT_VALUES, SCATTER_FORWARD);
414  CHKERR ts_ctx->mField.getInterface<VecManager>()->setLocalGhostVector(
415  ts_ctx->problemName, COL, u, INSERT_VALUES, SCATTER_REVERSE);
416  if (ts_ctx->zeroMatrix) {
417  CHKERR MatZeroEntries(B);
418  }
419  int step;
420  CHKERR TSGetTimeStepNumber(ts, &step);
421 
422  ts_ctx->matAssembleSwitch =
423  boost::movelib::make_unique<bool>(ts_ctx->zeroMatrix);
424 
425  // preproces
426  for (auto &bit : ts_ctx->preProcess_IJacobian) {
427  bit->matAssembleSwitch = boost::move(ts_ctx->matAssembleSwitch);
428  bit->ts_u = u;
429  bit->ts_u_t = u_t;
430  bit->ts_u_tt = u_tt;
431  bit->ts_A = A;
432  bit->ts_B = B;
433  bit->ts_t = t;
434  bit->ts_a = a;
435  bit->ts_step = step;
436 
437  CHKERR bit->setTsCtx(TSMethod::CTX_TSSETIJACOBIAN);
438  CHKERR bit->setTs(ts);
440  *bit);
441  CHKERR bit->setTsCtx(TSMethod::CTX_TSNONE);
442  ts_ctx->matAssembleSwitch = boost::move(bit->matAssembleSwitch);
443  }
444 
445  for (auto &lit : ts_ctx->loops_to_do_IJacobian) {
446  lit.second->matAssembleSwitch = boost::move(ts_ctx->matAssembleSwitch);
447  lit.second->ts_u = u;
448  lit.second->ts_u_t = u_t;
449  lit.second->ts_u_tt = u_tt;
450  lit.second->ts_A = A;
451  lit.second->ts_B = B;
452  lit.second->ts_t = t;
453  lit.second->ts_v = v;
454  lit.second->ts_a = a;
455  lit.second->ts_step = step;
456  CHKERR lit.second->setTsCtx(TSMethod::CTX_TSSETIJACOBIAN);
457  CHKERR lit.second->setTs(ts);
459  ts_ctx->problemName, lit.first, *(lit.second), nullptr, ts_ctx->bH);
460  CHKERR lit.second->setTsCtx(TSMethod::CTX_TSNONE);
461  ts_ctx->matAssembleSwitch = boost::move(lit.second->matAssembleSwitch);
462  }
463 
464  // post process
465  for (auto &bit : ts_ctx->postProcess_IJacobian) {
466  bit->matAssembleSwitch = boost::move(ts_ctx->matAssembleSwitch);
467  bit->ts_u = u;
468  bit->ts_u_t = u_t;
469  bit->ts_u_tt = u_tt;
470  bit->ts_A = A;
471  bit->ts_B = B;
472  bit->ts_t = t;
473  bit->ts_v = v;
474  bit->ts_a = a;
475  bit->ts_step = step;
476  CHKERR bit->setTsCtx(TSMethod::CTX_TSSETIJACOBIAN);
477  CHKERR bit->setTs(ts);
479  *bit);
480  CHKERR bit->setTsCtx(TSMethod::CTX_TSNONE);
481  ts_ctx->matAssembleSwitch = boost::move(bit->matAssembleSwitch);
482  }
483 
484  if (ts_ctx->matAssembleSwitch) {
485  CHKERR MatAssemblyBegin(B, MAT_FINAL_ASSEMBLY);
486  CHKERR MatAssemblyEnd(B, MAT_FINAL_ASSEMBLY);
487  }
488  PetscLogEventEnd(ts_ctx->MOFEM_EVENT_TsCtxI2Function, 0, 0, 0, 0);
490 }
491 
492 PetscErrorCode TsSetI2Function(TS ts, PetscReal t, Vec u, Vec u_t, Vec u_tt,
493  Vec F, void *ctx) {
495  TsCtx *ts_ctx = static_cast<TsCtx *>(ctx);
496  PetscLogEventBegin(ts_ctx->MOFEM_EVENT_TsCtxIFunction, 0, 0, 0, 0);
497  CHKERR VecGhostUpdateBegin(u, INSERT_VALUES, SCATTER_FORWARD);
498  CHKERR VecGhostUpdateEnd(u, INSERT_VALUES, SCATTER_FORWARD);
499  CHKERR VecGhostUpdateBegin(u_t, INSERT_VALUES, SCATTER_FORWARD);
500  CHKERR VecGhostUpdateEnd(u_t, INSERT_VALUES, SCATTER_FORWARD);
501  CHKERR VecGhostUpdateBegin(u_tt, INSERT_VALUES, SCATTER_FORWARD);
502  CHKERR VecGhostUpdateEnd(u_tt, INSERT_VALUES, SCATTER_FORWARD);
503  CHKERR ts_ctx->mField.getInterface<VecManager>()->setLocalGhostVector(
504  ts_ctx->problemName, COL, u, INSERT_VALUES, SCATTER_REVERSE);
505 
506  auto zero_ghost_vec = [](Vec g) {
508  Vec l;
509  CHKERR VecGhostGetLocalForm(g, &l);
510  double *a;
511  CHKERR VecGetArray(l, &a);
512  int s;
513  CHKERR VecGetLocalSize(l, &s);
514  for (int i = 0; i != s; ++i)
515  a[i] = 0;
516  CHKERR VecRestoreArray(l, &a);
517  CHKERR VecGhostRestoreLocalForm(g, &l);
519  };
520  CHKERR zero_ghost_vec(F);
521 
522  ts_ctx->vecAssembleSwitch = boost::movelib::make_unique<bool>(true);
523 
524  int step;
525  CHKERR TSGetTimeStepNumber(ts, &step);
526  // preprocess
527  for (auto &bit : ts_ctx->preProcess_IFunction) {
528  bit->vecAssembleSwitch = boost::move(ts_ctx->vecAssembleSwitch);
529  bit->ts_u = u;
530  bit->ts_u_t = u_t;
531  bit->ts_u_tt = u_tt;
532  bit->ts_F = F;
533  bit->ts_t = t;
534  bit->ts_step = step;
535  CHKERR bit->setTsCtx(TSMethod::CTX_TSSETIFUNCTION);
536  CHKERR bit->setTs(ts);
538  *bit);
539  CHKERR bit->setTsCtx(TSMethod::CTX_TSNONE);
540  ts_ctx->vecAssembleSwitch = boost::move(bit->vecAssembleSwitch);
541  }
542 
543  // fe loops
544  for (auto &lit : ts_ctx->loops_to_do_IFunction) {
545  lit.second->vecAssembleSwitch = boost::move(ts_ctx->vecAssembleSwitch);
546  lit.second->ts_u = u;
547  lit.second->ts_u_t = u_t;
548  lit.second->ts_u_tt = u_tt;
549  lit.second->ts_F = F;
550  lit.second->ts_t = t;
551  lit.second->ts_step = step;
552  CHKERR lit.second->setTsCtx(TSMethod::CTX_TSSETIFUNCTION);
553  CHKERR lit.second->setTs(ts);
555  ts_ctx->problemName, lit.first, *(lit.second), nullptr, ts_ctx->bH);
556  CHKERR lit.second->setTsCtx(TSMethod::CTX_TSNONE);
557  ts_ctx->vecAssembleSwitch = boost::move(lit.second->vecAssembleSwitch);
558  }
559 
560  // post process
561  for (auto &bit : ts_ctx->postProcess_IFunction) {
562  bit->vecAssembleSwitch = boost::move(ts_ctx->vecAssembleSwitch);
563  bit->ts_u = u;
564  bit->ts_u_t = u_t;
565  bit->ts_u_tt = u_tt;
566  bit->ts_F = F;
567  bit->ts_t = t;
568  bit->ts_step = step;
569  CHKERR bit->setTsCtx(TSMethod::CTX_TSSETIFUNCTION);
570  CHKERR bit->setTs(ts);
572  *bit);
573  CHKERR bit->setTsCtx(TSMethod::CTX_TSNONE);
574  ts_ctx->vecAssembleSwitch = boost::move(bit->vecAssembleSwitch);
575  }
576 
577  if (*ts_ctx->vecAssembleSwitch) {
578  CHKERR VecGhostUpdateBegin(F, ADD_VALUES, SCATTER_REVERSE);
579  CHKERR VecGhostUpdateEnd(F, ADD_VALUES, SCATTER_REVERSE);
580  CHKERR VecAssemblyBegin(F);
581  CHKERR VecAssemblyEnd(F);
582  }
583 
584  PetscLogEventEnd(ts_ctx->MOFEM_EVENT_TsCtxIFunction, 0, 0, 0, 0);
586 }
587 
588 } // namespace MoFEM
PetscErrorCode TsSetI2Jacobian(TS ts, PetscReal t, Vec u, Vec u_t, Vec u_tt, PetscReal v, PetscReal a, Mat A, Mat B, void *ctx)
Calculation Jaconian for second order PDE in time.
Definition: TsCtx.cpp:401
FEMethodsSequence loops_to_do_Monitor
Definition: TsCtx.hpp:50
PetscLogEvent MOFEM_EVENT_TsCtxIFunction
Definition: TsCtx.hpp:248
bool zeroMatrix
Definition: TsCtx.hpp:64
BasicMethodsSequence preProcess_RHSJacobian
Definition: TsCtx.hpp:59
Interface for Time Stepping (TS) solver.
Definition: TsCtx.hpp:27
MoFEMTypes bH
If set to MF_EXIST check if element exist.
Definition: TsCtx.hpp:33
BasicMethodsSequence preProcess_IJacobian
Definition: TsCtx.hpp:53
PetscErrorCode TsSetRHSFunction(TS ts, PetscReal t, Vec u, Vec F, void *ctx)
TS solver function.
Definition: TsCtx.cpp:241
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:477
PetscErrorCode TsMonitorSet(TS ts, PetscInt step, PetscReal t, Vec u, void *ctx)
Set monitor for TS solver.
Definition: TsCtx.cpp:190
std::string problemName
Definition: TsCtx.hpp:32
boost::movelib::unique_ptr< bool > vecAssembleSwitch
Definition: TsCtx.hpp:254
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:21
PetscLogEvent MOFEM_EVENT_TsCtxRHSFunction
Definition: TsCtx.hpp:246
MoFEMErrorCode getInterface(const MOFEMuuid &uuid, IFACE *&iface) const
Get interface by uuid and return reference to pointer of interface.
PetscErrorCode TsSetIFunction(TS ts, PetscReal t, Vec u, Vec u_t, Vec F, void *ctx)
Set IFunction for TS solver.
Definition: TsCtx.cpp:17
FEMethodsSequence loops_to_do_RHSJacobian
Definition: TsCtx.hpp:51
BasicMethodsSequence postProcess_IFunction
Definition: TsCtx.hpp:56
BasicMethodsSequence preProcess_IFunction
Definition: TsCtx.hpp:55
BasicMethodsSequence postProcess_IJacobian
Definition: TsCtx.hpp:54
Vector manager is used to create vectors \mofem_vectors.
Definition: VecManager.hpp:36
FEMethodsSequence loops_to_do_IJacobian
Definition: TsCtx.hpp:48
virtual MoFEMErrorCode problem_basic_method_preProcess(const Problem *problem_ptr, BasicMethod &method, int verb=DEFAULT_VERBOSITY)=0
Set data for BasicMethod.
#define CHKERR
Inline error check.
Definition: definitions.h:596
PetscLogEvent MOFEM_EVENT_TsCtxI2Function
Definition: TsCtx.hpp:251
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:492
PetscLogEvent MOFEM_EVENT_TsCtxMonitor
Definition: TsCtx.hpp:250
boost::movelib::unique_ptr< bool > matAssembleSwitch
Definition: TsCtx.hpp:255
PetscErrorCode TsSetRHSJacobian(TS ts, PetscReal t, Vec u, Mat A, Mat B, void *ctx)
TS solver function.
Definition: TsCtx.cpp:325
FEMethodsSequence loops_to_do_IFunction
Definition: TsCtx.hpp:49
MoFEM::Interface & mField
Definition: TsCtx.hpp:29
BasicMethodsSequence preProcess_Monitor
Definition: TsCtx.hpp:57
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:407
FEMethodsSequence loops_to_do_RHSFunction
Definition: TsCtx.hpp:52
PetscLogEvent MOFEM_EVENT_TsCtxRHSJacobian
Definition: TsCtx.hpp:247
virtual MoFEMErrorCode problem_basic_method_postProcess(const Problem *problem_ptr, BasicMethod &method, int verb=DEFAULT_VERBOSITY)=0
Set data for BasicMethodThis function set data about problem, adjacencies and other multi-indices in ...
BasicMethodsSequence postProcess_Monitor
Definition: TsCtx.hpp:58
BasicMethodsSequence postProcess_RHSJacobian
Definition: TsCtx.hpp:61
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:108
DEPRECATED MoFEMErrorCode loop_finite_elements(const Problem *problem_ptr, const std::string &fe_name, FEMethod &method, int lower_rank, int upper_rank, MoFEMTypes bh, int verb=DEFAULT_VERBOSITY)