v0.9.1
Public Types | Public Member Functions | Public Attributes | Private Attributes | Friends | List of all members
MoFEM::TsCtx Struct Reference

Interface for Time Stepping (TS) solver. More...

#include <src/petsc/TsCtx.hpp>

Collaboration diagram for MoFEM::TsCtx:
[legend]

Public Types

typedef MoFEM::PairNameFEMethodPtr PairNameFEMethodPtr
 
typedef MoFEM::FEMethodsSequence FEMethodsSequence
 
typedef MoFEM::BasicMethodsSequence BasicMethodsSequence
 

Public Member Functions

 TsCtx (MoFEM::Interface &m_field, const std::string &problem_name)
 
FEMethodsSequenceget_loops_to_do_IFunction ()
 Get the loops to do IFunction object. More...
 
FEMethodsSequenceget_loops_to_do_RHSFunction ()
 Get the loops to do RHSFunction object. More...
 
FEMethodsSequenceget_loops_to_do_IJacobian ()
 Get the loops to do IJacobian object. More...
 
FEMethodsSequenceget_loops_to_do_RHSJacobian ()
 Get the loops to do RHSJacobian object. More...
 
FEMethodsSequenceget_loops_to_do_Monitor ()
 Get the loops to do Monitor object. More...
 
BasicMethodsSequenceget_preProcess_to_do_IFunction ()
 Get the preProcess to do IFunction object. More...
 
BasicMethodsSequenceget_postProcess_to_do_IFunction ()
 Get the postProcess to do IFunction object. More...
 
BasicMethodsSequenceget_preProcess_to_do_IJacobian ()
 Get the preProcess to do IJacobian object. More...
 
BasicMethodsSequenceget_postProcess_to_do_IJacobian ()
 Get the postProcess to do IJacobian object. More...
 
BasicMethodsSequenceget_preProcess_to_do_Monitor ()
 Get the preProcess to do Monitor object. More...
 
BasicMethodsSequenceget_postProcess_to_do_Monitor ()
 Get the postProcess to do Monitor object. More...
 
BasicMethodsSequenceget_preProcess_to_do_RHSJacobian ()
 Get the preProcess to do RHSJacobian object. More...
 
BasicMethodsSequenceget_postProcess_to_do_RHSJacobian ()
 Get the postProcess to do RHSJacobian object. More...
 
BasicMethodsSequenceget_preProcess_to_do_RHSFunction ()
 Get the preProcess to do RHSFunction object. More...
 
BasicMethodsSequenceget_postProcess_to_do_RHSFunction ()
 Get the postProcess to do RHSFunction object. More...
 

Public Attributes

DEPRECATED typedef MoFEM::PairNameFEMethodPtr loop_pair_type
 
DEPRECATED typedef MoFEM::FEMethodsSequence loops_to_do_type
 
DEPRECATED typedef MoFEM::BasicMethodsSequence basic_method_to_do
 
MoFEM::InterfacemField
 
moab::Interface & moab
 
std::string problemName
 
MoFEMTypes bH
 If set to MF_EXIST check if element exist. More...
 
FEMethodsSequence loops_to_do_IJacobian
 
FEMethodsSequence loops_to_do_IFunction
 
FEMethodsSequence loops_to_do_Monitor
 
FEMethodsSequence loops_to_do_RHSJacobian
 
FEMethodsSequence loops_to_do_RHSFunction
 
BasicMethodsSequence preProcess_IJacobian
 
BasicMethodsSequence postProcess_IJacobian
 
BasicMethodsSequence preProcess_IFunction
 
BasicMethodsSequence postProcess_IFunction
 
BasicMethodsSequence preProcess_Monitor
 
BasicMethodsSequence postProcess_Monitor
 
BasicMethodsSequence preProcess_RHSJacobian
 
BasicMethodsSequence preProcess_RHSFunction
 
BasicMethodsSequence postProcess_RHSJacobian
 
BasicMethodsSequence postProcess_RHSFunction
 
bool zeroMatrix
 

Private Attributes

PetscLogEvent MOFEM_EVENT_TsCtxRHSFunction
 
PetscLogEvent MOFEM_EVENT_TsCtxRHSJacobian
 
PetscLogEvent MOFEM_EVENT_TsCtxIFunction
 
PetscLogEvent MOFEM_EVENT_TsCtxIJacobian
 
PetscLogEvent MOFEM_EVENT_TsCtxMonitor
 
PetscLogEvent MOFEM_EVENT_TsCtxI2Function
 
PetscLogEvent MOFEM_EVENT_TsCtxI2Jacobian
 
boost::movelib::unique_ptr< boolvecAssembleSwitch
 
boost::movelib::unique_ptr< boolmatAssembleSwitch
 

Friends

PetscErrorCode TsSetIFunction (TS ts, PetscReal t, Vec u, Vec u_t, Vec F, void *ctx)
 Set IFunction for TS solver. More...
 
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. More...
 
PetscErrorCode TsMonitorSet (TS ts, PetscInt step, PetscReal t, Vec u, void *ctx)
 Set monitor for TS solver. More...
 
PetscErrorCode TsSetRHSFunction (TS ts, PetscReal t, Vec u, Vec F, void *ctx)
 TS solver function. More...
 
PetscErrorCode TsSetRHSJacobian (TS ts, PetscReal t, Vec u, Mat A, Mat B, void *ctx)
 TS solver function. More...
 
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. More...
 
PetscErrorCode TsSetI2Jacobian (TS ts, PetscReal t, Vec U, Vec U_t, Vec U_tt, PetscReal v, PetscReal a, Mat J, Mat P, void *ctx)
 Calculation Jaconian for second order PDE in time. More...
 

Detailed Description

Interface for Time Stepping (TS) solver.

Examples
Remodeling.cpp.

Definition at line 27 of file TsCtx.hpp.

Member Typedef Documentation

◆ BasicMethodsSequence

Definition at line 46 of file TsCtx.hpp.

◆ FEMethodsSequence

Definition at line 45 of file TsCtx.hpp.

◆ PairNameFEMethodPtr

Definition at line 44 of file TsCtx.hpp.

Constructor & Destructor Documentation

◆ TsCtx()

MoFEM::TsCtx::TsCtx ( MoFEM::Interface m_field,
const std::string &  problem_name 
)

Definition at line 66 of file TsCtx.hpp.

67  : mField(m_field), moab(m_field.get_moab()), problemName(problem_name),
68  bH(MF_EXIST), zeroMatrix(true) {
69  PetscLogEventRegister("LoopTsIFunction", 0, &MOFEM_EVENT_TsCtxIFunction);
70  PetscLogEventRegister("LoopTsIJacobian", 0, &MOFEM_EVENT_TsCtxIJacobian);
71  PetscLogEventRegister("LoopTsRHSFunction", 0,
73  PetscLogEventRegister("LoopTsRHSJacobian", 0,
75  PetscLogEventRegister("LoopTsMonitor", 0, &MOFEM_EVENT_TsCtxMonitor);
76  PetscLogEventRegister("LoopTsI2Function", 0, &MOFEM_EVENT_TsCtxI2Function);
77  PetscLogEventRegister("LoopTsI2Jacobian", 0, &MOFEM_EVENT_TsCtxI2Jacobian);
78  }
PetscLogEvent MOFEM_EVENT_TsCtxIFunction
Definition: TsCtx.hpp:248
bool zeroMatrix
Definition: TsCtx.hpp:64
virtual moab::Interface & get_moab()=0
PetscLogEvent MOFEM_EVENT_TsCtxIJacobian
Definition: TsCtx.hpp:249
MoFEMTypes bH
If set to MF_EXIST check if element exist.
Definition: TsCtx.hpp:33
std::string problemName
Definition: TsCtx.hpp:32
moab::Interface & moab
Definition: TsCtx.hpp:30
PetscLogEvent MOFEM_EVENT_TsCtxRHSFunction
Definition: TsCtx.hpp:246
PetscLogEvent MOFEM_EVENT_TsCtxI2Jacobian
Definition: TsCtx.hpp:252
PetscLogEvent MOFEM_EVENT_TsCtxI2Function
Definition: TsCtx.hpp:251
PetscLogEvent MOFEM_EVENT_TsCtxMonitor
Definition: TsCtx.hpp:250
MoFEM::Interface & mField
Definition: TsCtx.hpp:29
PetscLogEvent MOFEM_EVENT_TsCtxRHSJacobian
Definition: TsCtx.hpp:247

Member Function Documentation

◆ get_loops_to_do_IFunction()

FEMethodsSequence& MoFEM::TsCtx::get_loops_to_do_IFunction ( )

Get the loops to do IFunction object.

It is sequence of finite elements used to evaluate the right hand side of implicit time solver.

Returns
FEMethodsSequence&

Definition at line 88 of file TsCtx.hpp.

88  {
89  return loops_to_do_IFunction;
90  }
FEMethodsSequence loops_to_do_IFunction
Definition: TsCtx.hpp:49

◆ get_loops_to_do_IJacobian()

FEMethodsSequence& MoFEM::TsCtx::get_loops_to_do_IJacobian ( )

Get the loops to do IJacobian object.

It is sequence of finite elements used to evalite the left hand sie of implimcit time solver.

Returns
FEMethodsSequence&

Definition at line 112 of file TsCtx.hpp.

112  {
113  return loops_to_do_IJacobian;
114  }
FEMethodsSequence loops_to_do_IJacobian
Definition: TsCtx.hpp:48

◆ get_loops_to_do_Monitor()

FEMethodsSequence& MoFEM::TsCtx::get_loops_to_do_Monitor ( )

Get the loops to do Monitor object.

It is sequence used to monitor solution of time solver.

Returns
FEMethodsSequence&

Definition at line 135 of file TsCtx.hpp.

135 { return loops_to_do_Monitor; }
FEMethodsSequence loops_to_do_Monitor
Definition: TsCtx.hpp:50

◆ get_loops_to_do_RHSFunction()

FEMethodsSequence& MoFEM::TsCtx::get_loops_to_do_RHSFunction ( )

Get the loops to do RHSFunction object.

It is sequence of finite elements used to evaluate the right hand side of implicit time solver.

Returns
FEMethodsSequence&

Definition at line 100 of file TsCtx.hpp.

100  {
102  }
FEMethodsSequence loops_to_do_RHSFunction
Definition: TsCtx.hpp:52

◆ get_loops_to_do_RHSJacobian()

FEMethodsSequence& MoFEM::TsCtx::get_loops_to_do_RHSJacobian ( )

Get the loops to do RHSJacobian object.

It is sequence of finite elements used to evalite the left hand sie of implimcit time solver.

Returns
FEMethodsSequence&

Definition at line 124 of file TsCtx.hpp.

124  {
126  }
FEMethodsSequence loops_to_do_RHSJacobian
Definition: TsCtx.hpp:51

◆ get_postProcess_to_do_IFunction()

BasicMethodsSequence& MoFEM::TsCtx::get_postProcess_to_do_IFunction ( )

Get the postProcess to do IFunction object.

Returns
BasicMethodsSequence&

Definition at line 151 of file TsCtx.hpp.

151  {
152  return postProcess_IFunction;
153  }
BasicMethodsSequence postProcess_IFunction
Definition: TsCtx.hpp:56

◆ get_postProcess_to_do_IJacobian()

BasicMethodsSequence& MoFEM::TsCtx::get_postProcess_to_do_IJacobian ( )

Get the postProcess to do IJacobian object.

Returns
BasicMethodsSequence&

Definition at line 169 of file TsCtx.hpp.

169  {
170  return postProcess_IJacobian;
171  }
BasicMethodsSequence postProcess_IJacobian
Definition: TsCtx.hpp:54

◆ get_postProcess_to_do_Monitor()

BasicMethodsSequence& MoFEM::TsCtx::get_postProcess_to_do_Monitor ( )

Get the postProcess to do Monitor object.

Returns
BasicMethodsSequence&
Examples
Remodeling.cpp.

Definition at line 187 of file TsCtx.hpp.

187  {
188  return postProcess_Monitor;
189  }
BasicMethodsSequence postProcess_Monitor
Definition: TsCtx.hpp:58

◆ get_postProcess_to_do_RHSFunction()

BasicMethodsSequence& MoFEM::TsCtx::get_postProcess_to_do_RHSFunction ( )

Get the postProcess to do RHSFunction object.

Returns
BasicMethodsSequence&

Definition at line 223 of file TsCtx.hpp.

223  {
225  }
BasicMethodsSequence postProcess_RHSFunction
Definition: TsCtx.hpp:62

◆ get_postProcess_to_do_RHSJacobian()

BasicMethodsSequence& MoFEM::TsCtx::get_postProcess_to_do_RHSJacobian ( )

Get the postProcess to do RHSJacobian object.

Returns
BasicMethodsSequence&

Definition at line 205 of file TsCtx.hpp.

205  {
207  }
BasicMethodsSequence postProcess_RHSJacobian
Definition: TsCtx.hpp:61

◆ get_preProcess_to_do_IFunction()

BasicMethodsSequence& MoFEM::TsCtx::get_preProcess_to_do_IFunction ( )

Get the preProcess to do IFunction object.

Returns
BasicMethodsSequence&

Definition at line 142 of file TsCtx.hpp.

142  {
143  return preProcess_IFunction;
144  }
BasicMethodsSequence preProcess_IFunction
Definition: TsCtx.hpp:55

◆ get_preProcess_to_do_IJacobian()

BasicMethodsSequence& MoFEM::TsCtx::get_preProcess_to_do_IJacobian ( )

Get the preProcess to do IJacobian object.

Returns
BasicMethodsSequence&

Definition at line 160 of file TsCtx.hpp.

160  {
161  return preProcess_IJacobian;
162  }
BasicMethodsSequence preProcess_IJacobian
Definition: TsCtx.hpp:53

◆ get_preProcess_to_do_Monitor()

BasicMethodsSequence& MoFEM::TsCtx::get_preProcess_to_do_Monitor ( )

Get the preProcess to do Monitor object.

Returns
BasicMethodsSequence&

Definition at line 178 of file TsCtx.hpp.

178  {
179  return preProcess_Monitor;
180  }
BasicMethodsSequence preProcess_Monitor
Definition: TsCtx.hpp:57

◆ get_preProcess_to_do_RHSFunction()

BasicMethodsSequence& MoFEM::TsCtx::get_preProcess_to_do_RHSFunction ( )

Get the preProcess to do RHSFunction object.

Returns
BasicMethodsSequence&

Definition at line 214 of file TsCtx.hpp.

214  {
215  return preProcess_RHSFunction;
216  }
BasicMethodsSequence preProcess_RHSFunction
Definition: TsCtx.hpp:60

◆ get_preProcess_to_do_RHSJacobian()

BasicMethodsSequence& MoFEM::TsCtx::get_preProcess_to_do_RHSJacobian ( )

Get the preProcess to do RHSJacobian object.

Returns
BasicMethodsSequence&

Definition at line 196 of file TsCtx.hpp.

196  {
197  return preProcess_RHSJacobian;
198  }
BasicMethodsSequence preProcess_RHSJacobian
Definition: TsCtx.hpp:59

Friends And Related Function Documentation

◆ TsMonitorSet

PetscErrorCode TsMonitorSet ( TS  ts,
PetscInt  step,
PetscReal  t,
Vec  u,
void *  ctx 
)
friend

Set monitor for TS solver.

See PETSc for details

Parameters
ts
step
t
u
ctx
Returns
PetscErrorCode

Definition at line 202 of file TsCtx.cpp.

203  {
205  TsCtx *ts_ctx = static_cast<TsCtx *>(ctx);
206  PetscLogEventBegin(ts_ctx->MOFEM_EVENT_TsCtxMonitor, 0, 0, 0, 0);
207  CHKERR VecGhostUpdateBegin(u, INSERT_VALUES, SCATTER_FORWARD);
208  CHKERR VecGhostUpdateEnd(u, INSERT_VALUES, SCATTER_FORWARD);
209  CHKERR ts_ctx->mField.getInterface<VecManager>()->setLocalGhostVector(
210  ts_ctx->problemName, COL, u, INSERT_VALUES, SCATTER_REVERSE);
211 
212  auto set = [&](auto &fe) {
213  fe.ts = ts;
214  fe.ts_u = u;
215  fe.ts_t = t;
216  fe.ts_step = step;
217  fe.ts_F = PETSC_NULL;
218  fe.ts_ctx = TSMethod::CTX_TSTSMONITORSET;
219  fe.snes_ctx = SnesMethod::CTX_SNESNONE;
220  fe.ksp_ctx = KspMethod::CTX_KSPNONE;
222  };
223 
224  auto unset = [&](auto &fe) {
225  fe.ts_ctx = TSMethod::CTX_TSNONE;
226  fe.data_ctx = PetscData::CtxSetNone;
227  };
228 
229  // preproces
230  for (auto &bit : ts_ctx->preProcess_Monitor) {
231  set(*bit);
232  CHKERR ts_ctx->mField.problem_basic_method_preProcess(ts_ctx->problemName,
233  *bit);
234  unset(*bit);
235  }
236 
237  for (auto &lit : ts_ctx->loops_to_do_Monitor) {
238  set(*lit.second);
239  CHKERR ts_ctx->mField.loop_finite_elements(
240  ts_ctx->problemName, lit.first, *(lit.second), nullptr, ts_ctx->bH);
241  unset(*lit.second);
242  }
243 
244  // post process
245  for (auto &bit : ts_ctx->postProcess_Monitor) {
246  set(*bit);
247  CHKERR ts_ctx->mField.problem_basic_method_postProcess(ts_ctx->problemName,
248  *bit);
249  unset(*bit);
250  }
251 
252  PetscLogEventEnd(ts_ctx->MOFEM_EVENT_TsCtxMonitor, 0, 0, 0, 0);
254 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:482
static constexpr Switches CtxSetTime
Definition: LoopMethods.hpp:68
TsCtx(MoFEM::Interface &m_field, const std::string &problem_name)
Definition: TsCtx.hpp:66
#define CHKERR
Inline error check.
Definition: definitions.h:601
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:412
static constexpr Switches CtxSetX
Definition: LoopMethods.hpp:65
static constexpr Switches CtxSetNone
Definition: LoopMethods.hpp:61

◆ TsSetI2Function

PetscErrorCode TsSetI2Function ( TS  ts,
PetscReal  t,
Vec  U,
Vec  U_t,
Vec  U_tt,
Vec  F,
void *  ctx 
)
friend

Calculation the right hand side for second order PDE in time.

PETSc for details

Parameters
ts
F
fun
ctx
Returns
PetscErrorCode

Definition at line 525 of file TsCtx.cpp.

526  {
528  TsCtx *ts_ctx = static_cast<TsCtx *>(ctx);
529  PetscLogEventBegin(ts_ctx->MOFEM_EVENT_TsCtxIFunction, 0, 0, 0, 0);
530  CHKERR VecGhostUpdateBegin(u, INSERT_VALUES, SCATTER_FORWARD);
531  CHKERR VecGhostUpdateEnd(u, INSERT_VALUES, SCATTER_FORWARD);
532  CHKERR VecGhostUpdateBegin(u_t, INSERT_VALUES, SCATTER_FORWARD);
533  CHKERR VecGhostUpdateEnd(u_t, INSERT_VALUES, SCATTER_FORWARD);
534  CHKERR VecGhostUpdateBegin(u_tt, INSERT_VALUES, SCATTER_FORWARD);
535  CHKERR VecGhostUpdateEnd(u_tt, INSERT_VALUES, SCATTER_FORWARD);
536  CHKERR ts_ctx->mField.getInterface<VecManager>()->setLocalGhostVector(
537  ts_ctx->problemName, COL, u, INSERT_VALUES, SCATTER_REVERSE);
538 
539  auto zero_ghost_vec = [](Vec g) {
541  Vec l;
542  CHKERR VecGhostGetLocalForm(g, &l);
543  double *a;
544  CHKERR VecGetArray(l, &a);
545  int s;
546  CHKERR VecGetLocalSize(l, &s);
547  for (int i = 0; i != s; ++i)
548  a[i] = 0;
549  CHKERR VecRestoreArray(l, &a);
550  CHKERR VecGhostRestoreLocalForm(g, &l);
552  };
553  CHKERR zero_ghost_vec(F);
554 
555  ts_ctx->vecAssembleSwitch = boost::movelib::make_unique<bool>(true);
556 
557  int step;
558 #if PETSC_VERSION_GE(3, 8, 0)
559  CHKERR TSGetStepNumber(ts, &step);
560 #else
561  CHKERR TSGetTimeStepNumber(ts, &step);
562 #endif
563 
564  auto set = [&](auto &fe) {
565  fe.ts_u = u;
566  fe.ts_u_t = u_t;
567  fe.ts_u_tt = u_tt;
568  fe.ts_F = F;
569  fe.ts_t = t;
570  fe.ts_step = step;
571  fe.ts_ctx = TSMethod::CTX_TSSETIFUNCTION;
572  fe.snes_ctx = SnesMethod::CTX_SNESSETFUNCTION;
573  fe.ksp_ctx = KspMethod::CTX_SETFUNCTION;
574  fe.data_ctx = PetscData::CtxSetF | PetscData::CtxSetX |
577  fe.ts = ts;
578  };
579 
580  auto unset = [&](auto &fe) {
581  fe.ts_ctx = TSMethod::CTX_TSNONE;
582  fe.snes_ctx = SnesMethod::CTX_SNESNONE;
583  fe.ksp_ctx = KspMethod::CTX_KSPNONE;
584  fe.data_ctx = PetscData::CtxSetNone;
585  };
586 
587  // preprocess
588  for (auto &bit : ts_ctx->preProcess_IFunction) {
589  bit->vecAssembleSwitch = boost::move(ts_ctx->vecAssembleSwitch);
590  set(*bit);
591  CHKERR ts_ctx->mField.problem_basic_method_preProcess(ts_ctx->problemName,
592  *bit);
593  unset(*bit);
594  ts_ctx->vecAssembleSwitch = boost::move(bit->vecAssembleSwitch);
595  }
596 
597  // fe loops
598  for (auto &lit : ts_ctx->loops_to_do_IFunction) {
599  lit.second->vecAssembleSwitch = boost::move(ts_ctx->vecAssembleSwitch);
600  set(*lit.second);
601  CHKERR ts_ctx->mField.loop_finite_elements(
602  ts_ctx->problemName, lit.first, *(lit.second), nullptr, ts_ctx->bH);
603  unset(*lit.second);
604  ts_ctx->vecAssembleSwitch = boost::move(lit.second->vecAssembleSwitch);
605  }
606 
607  // post process
608  for (auto &bit : ts_ctx->postProcess_IFunction) {
609  bit->vecAssembleSwitch = boost::move(ts_ctx->vecAssembleSwitch);
610  set(*bit);
611  CHKERR ts_ctx->mField.problem_basic_method_postProcess(ts_ctx->problemName,
612  *bit);
613  unset(*bit);
614  ts_ctx->vecAssembleSwitch = boost::move(bit->vecAssembleSwitch);
615  }
616 
617  if (*ts_ctx->vecAssembleSwitch) {
618  CHKERR VecGhostUpdateBegin(F, ADD_VALUES, SCATTER_REVERSE);
619  CHKERR VecGhostUpdateEnd(F, ADD_VALUES, SCATTER_REVERSE);
620  CHKERR VecAssemblyBegin(F);
621  CHKERR VecAssemblyEnd(F);
622  }
623 
624  PetscLogEventEnd(ts_ctx->MOFEM_EVENT_TsCtxIFunction, 0, 0, 0, 0);
626 }
static constexpr Switches CtxSetX_T
Definition: LoopMethods.hpp:66
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:482
static constexpr Switches CtxSetTime
Definition: LoopMethods.hpp:68
TsCtx(MoFEM::Interface &m_field, const std::string &problem_name)
Definition: TsCtx.hpp:66
static constexpr Switches CtxSetF
Definition: LoopMethods.hpp:62
static constexpr Switches CtxSetX_TT
Definition: LoopMethods.hpp:67
#define CHKERR
Inline error check.
Definition: definitions.h:601
FTensor::Index< 'i', 2 > i
[Common data]
Definition: ContactOps.hpp:24
FTensor::Index< 'l', 2 > l
Definition: ElasticOps.hpp:29
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:412
static constexpr Switches CtxSetX
Definition: LoopMethods.hpp:65
static constexpr Switches CtxSetNone
Definition: LoopMethods.hpp:61

◆ TsSetI2Jacobian

PetscErrorCode TsSetI2Jacobian ( TS  ts,
PetscReal  t,
Vec  U,
Vec  U_t,
Vec  U_tt,
PetscReal  v,
PetscReal  a,
Mat  J,
Mat  P,
void *  ctx 
)
friend

Calculation Jaconian for second order PDE in time.

PETSc for details

Parameters
ts
J
P
jac
ctx
Returns
PetscErrorCode

Definition at line 434 of file TsCtx.cpp.

436  {
438 
439  TsCtx *ts_ctx = static_cast<TsCtx *>(ctx);
440  PetscLogEventBegin(ts_ctx->MOFEM_EVENT_TsCtxI2Function, 0, 0, 0, 0);
441  CHKERR VecGhostUpdateBegin(u, INSERT_VALUES, SCATTER_FORWARD);
442  CHKERR VecGhostUpdateEnd(u, INSERT_VALUES, SCATTER_FORWARD);
443  CHKERR VecGhostUpdateBegin(u_t, INSERT_VALUES, SCATTER_FORWARD);
444  CHKERR VecGhostUpdateEnd(u_t, INSERT_VALUES, SCATTER_FORWARD);
445  CHKERR VecGhostUpdateBegin(u_tt, INSERT_VALUES, SCATTER_FORWARD);
446  CHKERR VecGhostUpdateEnd(u_tt, INSERT_VALUES, SCATTER_FORWARD);
447  CHKERR ts_ctx->mField.getInterface<VecManager>()->setLocalGhostVector(
448  ts_ctx->problemName, COL, u, INSERT_VALUES, SCATTER_REVERSE);
449  if (ts_ctx->zeroMatrix) {
450  CHKERR MatZeroEntries(B);
451  }
452  int step;
453 #if PETSC_VERSION_GE(3, 8, 0)
454  CHKERR TSGetStepNumber(ts, &step);
455 #else
456  CHKERR TSGetTimeStepNumber(ts, &step);
457 #endif
458 
459  ts_ctx->matAssembleSwitch =
460  boost::movelib::make_unique<bool>(ts_ctx->zeroMatrix);
461 
462  auto set = [&](auto &fe) {
463  fe.ts_u = u;
464  fe.ts_u_t = u_t;
465  fe.ts_u_tt = u_tt;
466  fe.ts_A = A;
467  fe.ts_B = B;
468  fe.ts_t = t;
469  fe.ts_a = a;
470  fe.ts_step = step;
471 
472  fe.ts_ctx = TSMethod::CTX_TSSETIJACOBIAN;
473  fe.snes_ctx = SnesMethod::CTX_SNESSETJACOBIAN;
474  fe.ksp_ctx = KspMethod::CTX_OPERATORS;
475  fe.data_ctx = PetscData::CtxSetA | PetscData::CtxSetB |
478  fe.ts = ts;
479  };
480 
481  auto unset = [&](auto &fe) {
482  fe.ts_ctx = TSMethod::CTX_TSNONE;
483  fe.snes_ctx = SnesMethod::CTX_SNESNONE;
484  fe.ksp_ctx = KspMethod::CTX_KSPNONE;
485  fe.data_ctx = PetscData::CtxSetNone;
486  };
487 
488  // preproces
489  for (auto &bit : ts_ctx->preProcess_IJacobian) {
490  bit->matAssembleSwitch = boost::move(ts_ctx->matAssembleSwitch);
491  set(*bit);
492  CHKERR ts_ctx->mField.problem_basic_method_preProcess(ts_ctx->problemName,
493  *bit);
494  unset(*bit);
495  ts_ctx->matAssembleSwitch = boost::move(bit->matAssembleSwitch);
496  }
497 
498  for (auto &lit : ts_ctx->loops_to_do_IJacobian) {
499  lit.second->matAssembleSwitch = boost::move(ts_ctx->matAssembleSwitch);
500  set(*lit.second);
501  CHKERR ts_ctx->mField.loop_finite_elements(
502  ts_ctx->problemName, lit.first, *(lit.second), nullptr, ts_ctx->bH);
503  unset(*lit.second);
504  ts_ctx->matAssembleSwitch = boost::move(lit.second->matAssembleSwitch);
505  }
506 
507  // post process
508  for (auto &bit : ts_ctx->postProcess_IJacobian) {
509  bit->matAssembleSwitch = boost::move(ts_ctx->matAssembleSwitch);
510  set(*bit);
511  CHKERR ts_ctx->mField.problem_basic_method_postProcess(ts_ctx->problemName,
512  *bit);
513  unset(*bit);
514  ts_ctx->matAssembleSwitch = boost::move(bit->matAssembleSwitch);
515  }
516 
517  if (ts_ctx->matAssembleSwitch) {
518  CHKERR MatAssemblyBegin(B, MAT_FINAL_ASSEMBLY);
519  CHKERR MatAssemblyEnd(B, MAT_FINAL_ASSEMBLY);
520  }
521  PetscLogEventEnd(ts_ctx->MOFEM_EVENT_TsCtxI2Function, 0, 0, 0, 0);
523 }
static constexpr Switches CtxSetX_T
Definition: LoopMethods.hpp:66
static constexpr Switches CtxSetB
Definition: LoopMethods.hpp:64
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:482
static constexpr Switches CtxSetTime
Definition: LoopMethods.hpp:68
TsCtx(MoFEM::Interface &m_field, const std::string &problem_name)
Definition: TsCtx.hpp:66
static constexpr Switches CtxSetX_TT
Definition: LoopMethods.hpp:67
#define CHKERR
Inline error check.
Definition: definitions.h:601
static constexpr Switches CtxSetA
Definition: LoopMethods.hpp:63
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:412
static constexpr Switches CtxSetX
Definition: LoopMethods.hpp:65
static constexpr Switches CtxSetNone
Definition: LoopMethods.hpp:61

◆ TsSetIFunction

PetscErrorCode TsSetIFunction ( TS  ts,
PetscReal  t,
Vec  u,
Vec  u_t,
Vec  F,
void *  ctx 
)
friend

Set IFunction for TS solver.

See petsc for details

Parameters
ts
t
u
u_t
F
ctx
Returns
PetscErrorCode

Definition at line 17 of file TsCtx.cpp.

18  {
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 #if PETSC_VERSION_GE(3, 8, 0)
49  CHKERR TSGetStepNumber(ts, &step);
50 #else
51  CHKERR TSGetTimeStepNumber(ts, &step);
52 #endif
53 
54  auto set = [&](auto &fe) {
55  fe.ts = ts;
56  fe.ts_u = u;
57  fe.ts_u_t = u_t;
58  fe.ts_F = F;
59  fe.ts_t = t;
60  fe.ts_step = step;
61  fe.ts_ctx = TSMethod::CTX_TSSETIFUNCTION;
62  fe.snes_ctx = SnesMethod::CTX_SNESSETFUNCTION;
63  fe.ksp_ctx = KspMethod::CTX_SETFUNCTION;
64  fe.data_ctx = PetscData::CtxSetF | PetscData::CtxSetX |
66  };
67 
68  auto unset = [&](auto &fe) {
69  fe.ts_ctx = TSMethod::CTX_TSNONE;
70  fe.snes_ctx = SnesMethod::CTX_SNESNONE;
71  fe.ksp_ctx = KspMethod::CTX_KSPNONE;
72  fe.data_ctx = PetscData::CtxSetNone;
73  };
74 
75  // preprocess
76  for (auto &bit : ts_ctx->preProcess_IFunction) {
77  bit->vecAssembleSwitch = boost::move(ts_ctx->vecAssembleSwitch);
78  set(*bit);
79  CHKERR ts_ctx->mField.problem_basic_method_preProcess(ts_ctx->problemName,
80  *bit);
81  unset(*bit);
82  ts_ctx->vecAssembleSwitch = boost::move(bit->vecAssembleSwitch);
83  }
84 
85  // fe loops
86  for (auto &lit : ts_ctx->loops_to_do_IFunction) {
87  lit.second->vecAssembleSwitch = boost::move(ts_ctx->vecAssembleSwitch);
88  set(*lit.second);
89  CHKERR ts_ctx->mField.loop_finite_elements(
90  ts_ctx->problemName, lit.first, *(lit.second), nullptr, ts_ctx->bH);
91  unset(*lit.second);
92  ts_ctx->vecAssembleSwitch = boost::move(lit.second->vecAssembleSwitch);
93  }
94 
95  // post process
96  for (auto &bit : ts_ctx->postProcess_IFunction) {
97  bit->vecAssembleSwitch = boost::move(ts_ctx->vecAssembleSwitch);
98  set(*bit);
99  CHKERR ts_ctx->mField.problem_basic_method_postProcess(ts_ctx->problemName,
100  *bit);
101  unset(*bit);
102  ts_ctx->vecAssembleSwitch = boost::move(bit->vecAssembleSwitch);
103  }
104 
105  if (*ts_ctx->vecAssembleSwitch) {
106  CHKERR VecGhostUpdateBegin(F, ADD_VALUES, SCATTER_REVERSE);
107  CHKERR VecGhostUpdateEnd(F, ADD_VALUES, SCATTER_REVERSE);
108  CHKERR VecAssemblyBegin(F);
109  CHKERR VecAssemblyEnd(F);
110  }
111 
112  PetscLogEventEnd(ts_ctx->MOFEM_EVENT_TsCtxIFunction, 0, 0, 0, 0);
114 }
static constexpr Switches CtxSetX_T
Definition: LoopMethods.hpp:66
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:482
static constexpr Switches CtxSetTime
Definition: LoopMethods.hpp:68
TsCtx(MoFEM::Interface &m_field, const std::string &problem_name)
Definition: TsCtx.hpp:66
static constexpr Switches CtxSetF
Definition: LoopMethods.hpp:62
#define CHKERR
Inline error check.
Definition: definitions.h:601
FTensor::Index< 'i', 2 > i
[Common data]
Definition: ContactOps.hpp:24
FTensor::Index< 'l', 2 > l
Definition: ElasticOps.hpp:29
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:412
static constexpr Switches CtxSetX
Definition: LoopMethods.hpp:65
static constexpr Switches CtxSetNone
Definition: LoopMethods.hpp:61

◆ TsSetIJacobian

PetscErrorCode TsSetIJacobian ( TS  ts,
PetscReal  t,
Vec  u,
Vec  U_t,
PetscReal  a,
Mat  A,
Mat  B,
void *  ctx 
)
friend

Set function evaluating jacobina in TS solver.

See PETSc for details

Parameters
ts
t
u
u_t
a
A
B
ctx
Returns
PetscErrorCode

Definition at line 116 of file TsCtx.cpp.

117  {
119 
120  TsCtx *ts_ctx = static_cast<TsCtx *>(ctx);
121  PetscLogEventBegin(ts_ctx->MOFEM_EVENT_TsCtxIFunction, 0, 0, 0, 0);
122  CHKERR VecGhostUpdateBegin(u, INSERT_VALUES, SCATTER_FORWARD);
123  CHKERR VecGhostUpdateEnd(u, INSERT_VALUES, SCATTER_FORWARD);
124  CHKERR VecGhostUpdateBegin(u_t, INSERT_VALUES, SCATTER_FORWARD);
125  CHKERR VecGhostUpdateEnd(u_t, INSERT_VALUES, SCATTER_FORWARD);
126  CHKERR ts_ctx->mField.getInterface<VecManager>()->setLocalGhostVector(
127  ts_ctx->problemName, COL, u, INSERT_VALUES, SCATTER_REVERSE);
128  if (ts_ctx->zeroMatrix) {
129  CHKERR MatZeroEntries(B);
130  }
131  int step;
132 #if PETSC_VERSION_GE(3, 8, 0)
133  CHKERR TSGetStepNumber(ts, &step);
134 #else
135  CHKERR TSGetTimeStepNumber(ts, &step);
136 #endif
137 
138  ts_ctx->matAssembleSwitch =
139  boost::movelib::make_unique<bool>(ts_ctx->zeroMatrix);
140 
141  auto set = [&](auto &fe) {
142  fe.ts = ts;
143  fe.ts_u = u;
144  fe.ts_u_t = u_t;
145  fe.ts_A = A;
146  fe.ts_B = B;
147  fe.ts_t = t;
148  fe.ts_a = a;
149  fe.ts_step = step;
150  fe.ts_ctx = TSMethod::CTX_TSSETIJACOBIAN;
151  fe.snes_ctx = SnesMethod::CTX_SNESSETJACOBIAN;
152  fe.ksp_ctx = KspMethod::CTX_OPERATORS;
153  fe.data_ctx = PetscData::CtxSetA | PetscData::CtxSetB |
156  };
157 
158  auto unset = [&](auto &fe) {
159  fe.ts_ctx = TSMethod::CTX_TSNONE;
160  fe.snes_ctx = SnesMethod::CTX_SNESNONE;
161  fe.ksp_ctx = KspMethod::CTX_KSPNONE;
162  fe.data_ctx = PetscData::CtxSetNone;
163  };
164 
165  // preproces
166  for (auto &bit : ts_ctx->preProcess_IJacobian) {
167  bit->matAssembleSwitch = boost::move(ts_ctx->matAssembleSwitch);
168  set(*bit);
169  CHKERR ts_ctx->mField.problem_basic_method_preProcess(ts_ctx->problemName,
170  *bit);
171  unset(*bit);
172  ts_ctx->matAssembleSwitch = boost::move(bit->matAssembleSwitch);
173  }
174 
175  for (auto &lit : ts_ctx->loops_to_do_IJacobian) {
176  lit.second->matAssembleSwitch = boost::move(ts_ctx->matAssembleSwitch);
177  set(*lit.second);
178  CHKERR ts_ctx->mField.loop_finite_elements(
179  ts_ctx->problemName, lit.first, *(lit.second), nullptr, ts_ctx->bH);
180  unset(*lit.second);
181  ts_ctx->matAssembleSwitch = boost::move(lit.second->matAssembleSwitch);
182  }
183 
184  // post process
185  for (auto &bit : ts_ctx->postProcess_IJacobian) {
186  bit->matAssembleSwitch = boost::move(ts_ctx->matAssembleSwitch);
187  set(*bit);
188  CHKERR ts_ctx->mField.problem_basic_method_postProcess(ts_ctx->problemName,
189  *bit);
190  unset(*bit);
191  ts_ctx->matAssembleSwitch = boost::move(bit->matAssembleSwitch);
192  }
193 
194  if (ts_ctx->matAssembleSwitch) {
195  CHKERR MatAssemblyBegin(B, MAT_FINAL_ASSEMBLY);
196  CHKERR MatAssemblyEnd(B, MAT_FINAL_ASSEMBLY);
197  }
198  PetscLogEventEnd(ts_ctx->MOFEM_EVENT_TsCtxIFunction, 0, 0, 0, 0);
200 }
static constexpr Switches CtxSetX_T
Definition: LoopMethods.hpp:66
static constexpr Switches CtxSetB
Definition: LoopMethods.hpp:64
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:482
static constexpr Switches CtxSetTime
Definition: LoopMethods.hpp:68
TsCtx(MoFEM::Interface &m_field, const std::string &problem_name)
Definition: TsCtx.hpp:66
#define CHKERR
Inline error check.
Definition: definitions.h:601
static constexpr Switches CtxSetA
Definition: LoopMethods.hpp:63
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:412
static constexpr Switches CtxSetX
Definition: LoopMethods.hpp:65
static constexpr Switches CtxSetNone
Definition: LoopMethods.hpp:61

◆ TsSetRHSFunction

PetscErrorCode TsSetRHSFunction ( TS  ts,
PetscReal  t,
Vec  u,
Vec  F,
void *  ctx 
)
friend

TS solver function.

See PETSc for details

Parameters
ts
t
u
F
ctx
Returns
PetscErrorCode

Definition at line 256 of file TsCtx.cpp.

256  {
258  TsCtx *ts_ctx = static_cast<TsCtx *>(ctx);
259  PetscLogEventBegin(ts_ctx->MOFEM_EVENT_TsCtxRHSFunction, 0, 0, 0, 0);
260  CHKERR VecGhostUpdateBegin(u, INSERT_VALUES, SCATTER_FORWARD);
261  CHKERR VecGhostUpdateEnd(u, INSERT_VALUES, SCATTER_FORWARD);
262  CHKERR ts_ctx->mField.getInterface<VecManager>()->setLocalGhostVector(
263  ts_ctx->problemName, COL, u, INSERT_VALUES, SCATTER_REVERSE);
264 
265  auto zero_ghost_vec = [](Vec g) {
267  Vec l;
268  CHKERR VecGhostGetLocalForm(g, &l);
269  double *a;
270  CHKERR VecGetArray(l, &a);
271  int s;
272  CHKERR VecGetLocalSize(l, &s);
273  for (int i = 0; i != s; ++i)
274  a[i] = 0;
275  CHKERR VecRestoreArray(l, &a);
276  CHKERR VecGhostRestoreLocalForm(g, &l);
278  };
279  CHKERR zero_ghost_vec(F);
280 
281  ts_ctx->vecAssembleSwitch = boost::movelib::make_unique<bool>(true);
282 
283  int step;
284 #if PETSC_VERSION_GE(3, 8, 0)
285  CHKERR TSGetStepNumber(ts, &step);
286 #else
287  CHKERR TSGetTimeStepNumber(ts, &step);
288 #endif
289 
290  auto set = [&](auto &fe) {
291  fe.ts_u = u;
292  fe.ts_F = F;
293  fe.ts_t = t;
294  fe.ts = ts;
295  fe.ts_step = step;
296  fe.ts_ctx = TSMethod::CTX_TSSETRHSFUNCTION;
297  fe.snes_ctx = SnesMethod::CTX_SNESSETFUNCTION;
298  fe.ksp_ctx = KspMethod::CTX_SETFUNCTION;
299  fe.data_ctx =
301  };
302 
303  auto unset = [&](auto &fe) {
304  fe.ts_ctx = TSMethod::CTX_TSNONE;
305  fe.snes_ctx = SnesMethod::CTX_SNESNONE;
306  fe.ksp_ctx = KspMethod::CTX_KSPNONE;
307  fe.data_ctx = PetscData::CtxSetNone;
308  };
309 
310  for (auto &bit : ts_ctx->preProcess_RHSJacobian) {
311  bit->vecAssembleSwitch = boost::move(ts_ctx->vecAssembleSwitch);
312  set(*bit);
313  CHKERR ts_ctx->mField.problem_basic_method_preProcess(ts_ctx->problemName,
314  *bit);
315  unset(*bit);
316  ts_ctx->vecAssembleSwitch = boost::move(bit->vecAssembleSwitch);
317  }
318 
319  // fe loops
320  for (auto &lit : ts_ctx->loops_to_do_RHSFunction) {
321  lit.second->vecAssembleSwitch = boost::move(ts_ctx->vecAssembleSwitch);
322  set(*lit.second);
323  CHKERR ts_ctx->mField.loop_finite_elements(
324  ts_ctx->problemName, lit.first, *(lit.second), nullptr, ts_ctx->bH);
325  unset(*lit.second);
326  ts_ctx->vecAssembleSwitch = boost::move(lit.second->vecAssembleSwitch);
327  }
328 
329  // post process
330  for (auto &bit : ts_ctx->postProcess_RHSJacobian) {
331  bit->vecAssembleSwitch = boost::move(ts_ctx->vecAssembleSwitch);
332  set(*bit);
333  CHKERR ts_ctx->mField.problem_basic_method_postProcess(ts_ctx->problemName,
334  *bit);
335  unset(*bit);
336  ts_ctx->vecAssembleSwitch = boost::move(bit->vecAssembleSwitch);
337  }
338 
339  if (*ts_ctx->vecAssembleSwitch) {
340  CHKERR VecGhostUpdateBegin(F, ADD_VALUES, SCATTER_REVERSE);
341  CHKERR VecGhostUpdateEnd(F, ADD_VALUES, SCATTER_REVERSE);
342  CHKERR VecAssemblyBegin(F);
343  CHKERR VecAssemblyEnd(F);
344  }
345 
346  PetscLogEventEnd(ts_ctx->MOFEM_EVENT_TsCtxRHSFunction, 0, 0, 0, 0);
348 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:482
static constexpr Switches CtxSetTime
Definition: LoopMethods.hpp:68
TsCtx(MoFEM::Interface &m_field, const std::string &problem_name)
Definition: TsCtx.hpp:66
static constexpr Switches CtxSetF
Definition: LoopMethods.hpp:62
#define CHKERR
Inline error check.
Definition: definitions.h:601
FTensor::Index< 'i', 2 > i
[Common data]
Definition: ContactOps.hpp:24
FTensor::Index< 'l', 2 > l
Definition: ElasticOps.hpp:29
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:412
static constexpr Switches CtxSetX
Definition: LoopMethods.hpp:65
static constexpr Switches CtxSetNone
Definition: LoopMethods.hpp:61

◆ TsSetRHSJacobian

PetscErrorCode TsSetRHSJacobian ( TS  ts,
PetscReal  t,
Vec  u,
Mat  A,
Mat  B,
void *  ctx 
)
friend

TS solver function.

See PETSc for details

Parameters
ts
t
u
A
B
ctx
Returns
PetscErrorCode

Definition at line 350 of file TsCtx.cpp.

351  {
353  TsCtx *ts_ctx = static_cast<TsCtx *>(ctx);
354  PetscLogEventBegin(ts_ctx->MOFEM_EVENT_TsCtxRHSJacobian, 0, 0, 0, 0);
355  CHKERR VecGhostUpdateBegin(u, INSERT_VALUES, SCATTER_FORWARD);
356  CHKERR VecGhostUpdateEnd(u, INSERT_VALUES, SCATTER_FORWARD);
357  CHKERR ts_ctx->mField.getInterface<VecManager>()->setLocalGhostVector(
358  ts_ctx->problemName, COL, u, INSERT_VALUES, SCATTER_REVERSE);
359 
360  if (ts_ctx->zeroMatrix) {
361  CHKERR MatZeroEntries(B);
362  }
363 
364  ts_ctx->matAssembleSwitch =
365  boost::movelib::make_unique<bool>(ts_ctx->zeroMatrix);
366 
367  int step;
368 #if PETSC_VERSION_GE(3, 8, 0)
369  CHKERR TSGetStepNumber(ts, &step);
370 #else
371  CHKERR TSGetTimeStepNumber(ts, &step);
372 #endif
373 
374  auto set = [&](auto &fe) {
375  fe.ts_u = u;
376  fe.ts_A = A;
377  fe.ts_B = B;
378  fe.ts_t = t;
379  fe.ts_step = step;
380  fe.ts_ctx = TSMethod::CTX_TSSETRHSJACOBIAN;
381  fe.snes_ctx = SnesMethod::CTX_SNESSETJACOBIAN;
382  fe.ksp_ctx = KspMethod::CTX_OPERATORS;
383  fe.data_ctx = PetscData::CtxSetA | PetscData::CtxSetB |
385  fe.ts = ts;
386  };
387 
388  auto unset = [&](auto &fe) {
389  fe.ts_ctx = TSMethod::CTX_TSNONE;
390  fe.snes_ctx = SnesMethod::CTX_SNESNONE;
391  fe.ksp_ctx = KspMethod::CTX_KSPNONE;
392  fe.data_ctx = PetscData::CtxSetNone;
393  };
394 
395  // preprocess
396  for (auto &bit : ts_ctx->preProcess_RHSJacobian) {
397  bit->matAssembleSwitch = boost::move(ts_ctx->matAssembleSwitch);
398  set(*bit);
399  CHKERR ts_ctx->mField.problem_basic_method_preProcess(ts_ctx->problemName,
400  *bit);
401  unset(*bit);
402  ts_ctx->matAssembleSwitch = boost::move(bit->matAssembleSwitch);
403  }
404 
405  // fe loops
406  for (auto &lit : ts_ctx->loops_to_do_RHSJacobian) {
407  lit.second->matAssembleSwitch = boost::move(ts_ctx->matAssembleSwitch);
408  set(*lit.second);
409  CHKERR ts_ctx->mField.loop_finite_elements(
410  ts_ctx->problemName, lit.first, *(lit.second), nullptr, ts_ctx->bH);
411  unset(*lit.second);
412  ts_ctx->matAssembleSwitch = boost::move(lit.second->matAssembleSwitch);
413  }
414 
415  // post process
416  for (auto &bit : ts_ctx->postProcess_RHSJacobian) {
417  bit->matAssembleSwitch = boost::move(ts_ctx->matAssembleSwitch);
418  set(*bit);
419  CHKERR ts_ctx->mField.problem_basic_method_postProcess(ts_ctx->problemName,
420  *bit);
421  unset(*bit);
422  ts_ctx->matAssembleSwitch = boost::move(bit->matAssembleSwitch);
423  }
424 
425  if (ts_ctx->matAssembleSwitch) {
426  CHKERR MatAssemblyBegin(B, MAT_FINAL_ASSEMBLY);
427  CHKERR MatAssemblyEnd(B, MAT_FINAL_ASSEMBLY);
428  }
429 
430  PetscLogEventEnd(ts_ctx->MOFEM_EVENT_TsCtxRHSJacobian, 0, 0, 0, 0);
432 }
static constexpr Switches CtxSetB
Definition: LoopMethods.hpp:64
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:482
static constexpr Switches CtxSetTime
Definition: LoopMethods.hpp:68
TsCtx(MoFEM::Interface &m_field, const std::string &problem_name)
Definition: TsCtx.hpp:66
#define CHKERR
Inline error check.
Definition: definitions.h:601
static constexpr Switches CtxSetA
Definition: LoopMethods.hpp:63
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:412
static constexpr Switches CtxSetX
Definition: LoopMethods.hpp:65
static constexpr Switches CtxSetNone
Definition: LoopMethods.hpp:61

Member Data Documentation

◆ basic_method_to_do

DEPRECATED typedef MoFEM::BasicMethodsSequence MoFEM::TsCtx::basic_method_to_do
Deprecated:
use BasicMethodsSequence

Definition at line 42 of file TsCtx.hpp.

◆ bH

MoFEMTypes MoFEM::TsCtx::bH

If set to MF_EXIST check if element exist.

Definition at line 33 of file TsCtx.hpp.

◆ loop_pair_type

DEPRECATED typedef MoFEM::PairNameFEMethodPtr MoFEM::TsCtx::loop_pair_type
Deprecated:
use PairNameFEMethodPtr

Definition at line 36 of file TsCtx.hpp.

◆ loops_to_do_IFunction

FEMethodsSequence MoFEM::TsCtx::loops_to_do_IFunction

Definition at line 49 of file TsCtx.hpp.

◆ loops_to_do_IJacobian

FEMethodsSequence MoFEM::TsCtx::loops_to_do_IJacobian

Definition at line 48 of file TsCtx.hpp.

◆ loops_to_do_Monitor

FEMethodsSequence MoFEM::TsCtx::loops_to_do_Monitor

Definition at line 50 of file TsCtx.hpp.

◆ loops_to_do_RHSFunction

FEMethodsSequence MoFEM::TsCtx::loops_to_do_RHSFunction

Definition at line 52 of file TsCtx.hpp.

◆ loops_to_do_RHSJacobian

FEMethodsSequence MoFEM::TsCtx::loops_to_do_RHSJacobian

Definition at line 51 of file TsCtx.hpp.

◆ loops_to_do_type

DEPRECATED typedef MoFEM::FEMethodsSequence MoFEM::TsCtx::loops_to_do_type
Deprecated:
use FEMethodsSequence

Definition at line 39 of file TsCtx.hpp.

◆ matAssembleSwitch

boost::movelib::unique_ptr<bool> MoFEM::TsCtx::matAssembleSwitch
private

Definition at line 255 of file TsCtx.hpp.

◆ mField

MoFEM::Interface& MoFEM::TsCtx::mField

Definition at line 29 of file TsCtx.hpp.

◆ moab

moab::Interface& MoFEM::TsCtx::moab

Definition at line 30 of file TsCtx.hpp.

◆ MOFEM_EVENT_TsCtxI2Function

PetscLogEvent MoFEM::TsCtx::MOFEM_EVENT_TsCtxI2Function
private

Definition at line 251 of file TsCtx.hpp.

◆ MOFEM_EVENT_TsCtxI2Jacobian

PetscLogEvent MoFEM::TsCtx::MOFEM_EVENT_TsCtxI2Jacobian
private

Definition at line 252 of file TsCtx.hpp.

◆ MOFEM_EVENT_TsCtxIFunction

PetscLogEvent MoFEM::TsCtx::MOFEM_EVENT_TsCtxIFunction
private

Definition at line 248 of file TsCtx.hpp.

◆ MOFEM_EVENT_TsCtxIJacobian

PetscLogEvent MoFEM::TsCtx::MOFEM_EVENT_TsCtxIJacobian
private

Definition at line 249 of file TsCtx.hpp.

◆ MOFEM_EVENT_TsCtxMonitor

PetscLogEvent MoFEM::TsCtx::MOFEM_EVENT_TsCtxMonitor
private

Definition at line 250 of file TsCtx.hpp.

◆ MOFEM_EVENT_TsCtxRHSFunction

PetscLogEvent MoFEM::TsCtx::MOFEM_EVENT_TsCtxRHSFunction
private

Definition at line 246 of file TsCtx.hpp.

◆ MOFEM_EVENT_TsCtxRHSJacobian

PetscLogEvent MoFEM::TsCtx::MOFEM_EVENT_TsCtxRHSJacobian
private

Definition at line 247 of file TsCtx.hpp.

◆ postProcess_IFunction

BasicMethodsSequence MoFEM::TsCtx::postProcess_IFunction

Definition at line 56 of file TsCtx.hpp.

◆ postProcess_IJacobian

BasicMethodsSequence MoFEM::TsCtx::postProcess_IJacobian

Definition at line 54 of file TsCtx.hpp.

◆ postProcess_Monitor

BasicMethodsSequence MoFEM::TsCtx::postProcess_Monitor

Definition at line 58 of file TsCtx.hpp.

◆ postProcess_RHSFunction

BasicMethodsSequence MoFEM::TsCtx::postProcess_RHSFunction

Definition at line 62 of file TsCtx.hpp.

◆ postProcess_RHSJacobian

BasicMethodsSequence MoFEM::TsCtx::postProcess_RHSJacobian

Definition at line 61 of file TsCtx.hpp.

◆ preProcess_IFunction

BasicMethodsSequence MoFEM::TsCtx::preProcess_IFunction

Definition at line 55 of file TsCtx.hpp.

◆ preProcess_IJacobian

BasicMethodsSequence MoFEM::TsCtx::preProcess_IJacobian

Definition at line 53 of file TsCtx.hpp.

◆ preProcess_Monitor

BasicMethodsSequence MoFEM::TsCtx::preProcess_Monitor

Definition at line 57 of file TsCtx.hpp.

◆ preProcess_RHSFunction

BasicMethodsSequence MoFEM::TsCtx::preProcess_RHSFunction

Definition at line 60 of file TsCtx.hpp.

◆ preProcess_RHSJacobian

BasicMethodsSequence MoFEM::TsCtx::preProcess_RHSJacobian

Definition at line 59 of file TsCtx.hpp.

◆ problemName

std::string MoFEM::TsCtx::problemName

Definition at line 32 of file TsCtx.hpp.

◆ vecAssembleSwitch

boost::movelib::unique_ptr<bool> MoFEM::TsCtx::vecAssembleSwitch
private

Definition at line 254 of file TsCtx.hpp.

◆ zeroMatrix

bool MoFEM::TsCtx::zeroMatrix

Definition at line 64 of file TsCtx.hpp.


The documentation for this struct was generated from the following file: