v0.13.0
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)
 
virtual ~TsCtx ()=default
 
FEMethodsSequencegetLoopsIFunction ()
 Get the loops to do IFunction object. More...
 
FEMethodsSequencegetLoopsRHSFunction ()
 Get the loops to do RHSFunction object. More...
 
FEMethodsSequencegetLoopsIJacobian ()
 Get the loops to do IJacobian object. More...
 
FEMethodsSequencegetLoopsRHSJacobian ()
 Get the loops to do RHSJacobian object. More...
 
FEMethodsSequencegetLoopsMonitor ()
 Get the loops to do Monitor object. More...
 
BasicMethodsSequencegetPreProcessIFunction ()
 Get the preProcess to do IFunction object. More...
 
BasicMethodsSequencegetPostProcessIFunction ()
 Get the postProcess to do IFunction object. More...
 
BasicMethodsSequencegetPreProcessIJacobian ()
 Get the preProcess to do IJacobian object. More...
 
BasicMethodsSequencegetPostProcessIJacobian ()
 Get the postProcess to do IJacobian object. More...
 
BasicMethodsSequencegetPreProcessMonitor ()
 Get the preProcess to do Monitor object. More...
 
BasicMethodsSequencegetPostProcessMonitor ()
 Get the postProcess to do Monitor object. More...
 
BasicMethodsSequencegetPreProcessRHSJacobian ()
 Get the preProcess to do RHSJacobian object. More...
 
BasicMethodsSequencegetPostProcessRHSJacobian ()
 Get the postProcess to do RHSJacobian object. More...
 
BasicMethodsSequencegetPreProcessRHSFunction ()
 Get the preProcess to do RHSFunction object. More...
 
BasicMethodsSequencegetPostProcessRHSFunction ()
 Get the postProcess to do RHSFunction object. More...
 
MoFEMErrorCode clearLoops ()
 Clear loops. More...
 

Public Attributes

MoFEM::InterfacemField
 
moab::Interface & moab
 
std::string problemName
 
MoFEMTypes bH
 If set to MF_EXIST check if element exist. More...
 
FEMethodsSequence loopsIJacobian
 
FEMethodsSequence loopsIFunction
 
FEMethodsSequence loopsMonitor
 
FEMethodsSequence loopsRHSJacobian
 
FEMethodsSequence loopsRHSFunction
 
BasicMethodsSequence preProcessIJacobian
 
BasicMethodsSequence postProcessIJacobian
 
BasicMethodsSequence preProcessIFunction
 
BasicMethodsSequence postProcessIFunction
 
BasicMethodsSequence preProcessMonitor
 
BasicMethodsSequence postProcessMonitor
 
BasicMethodsSequence preProcessRHSJacobian
 
BasicMethodsSequence preProcessRHSFunction
 
BasicMethodsSequence postProcessRHSJacobian
 
BasicMethodsSequence postProcessRHSFunction
 
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.

Definition at line 27 of file TsCtx.hpp.

Member Typedef Documentation

◆ BasicMethodsSequence

Definition at line 37 of file TsCtx.hpp.

◆ FEMethodsSequence

Examples
nonlinear_dynamics.cpp.

Definition at line 36 of file TsCtx.hpp.

◆ PairNameFEMethodPtr

Examples
EshelbianPlasticity.cpp, and nonlinear_dynamics.cpp.

Definition at line 35 of file TsCtx.hpp.

Constructor & Destructor Documentation

◆ TsCtx()

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

Definition at line 57 of file TsCtx.hpp.

58  : mField(m_field), moab(m_field.get_moab()), problemName(problem_name),
59  bH(MF_EXIST), zeroMatrix(true) {
60  PetscLogEventRegister("LoopTsIFunction", 0, &MOFEM_EVENT_TsCtxIFunction);
61  PetscLogEventRegister("LoopTsIJacobian", 0, &MOFEM_EVENT_TsCtxIJacobian);
62  PetscLogEventRegister("LoopTsRHSFunction", 0,
64  PetscLogEventRegister("LoopTsRHSJacobian", 0,
66  PetscLogEventRegister("LoopTsMonitor", 0, &MOFEM_EVENT_TsCtxMonitor);
67  PetscLogEventRegister("LoopTsI2Function", 0, &MOFEM_EVENT_TsCtxI2Function);
68  PetscLogEventRegister("LoopTsI2Jacobian", 0, &MOFEM_EVENT_TsCtxI2Jacobian);
69  }
@ MF_EXIST
Definition: definitions.h:113
virtual moab::Interface & get_moab()=0
PetscLogEvent MOFEM_EVENT_TsCtxI2Jacobian
Definition: TsCtx.hpp:236
bool zeroMatrix
Definition: TsCtx.hpp:55
MoFEMTypes bH
If set to MF_EXIST check if element exist.
Definition: TsCtx.hpp:33
moab::Interface & moab
Definition: TsCtx.hpp:30
PetscLogEvent MOFEM_EVENT_TsCtxRHSJacobian
Definition: TsCtx.hpp:231
PetscLogEvent MOFEM_EVENT_TsCtxMonitor
Definition: TsCtx.hpp:234
PetscLogEvent MOFEM_EVENT_TsCtxRHSFunction
Definition: TsCtx.hpp:230
PetscLogEvent MOFEM_EVENT_TsCtxIJacobian
Definition: TsCtx.hpp:233
std::string problemName
Definition: TsCtx.hpp:32
MoFEM::Interface & mField
Definition: TsCtx.hpp:29
PetscLogEvent MOFEM_EVENT_TsCtxIFunction
Definition: TsCtx.hpp:232
PetscLogEvent MOFEM_EVENT_TsCtxI2Function
Definition: TsCtx.hpp:235

◆ ~TsCtx()

virtual MoFEM::TsCtx::~TsCtx ( )
virtualdefault

Member Function Documentation

◆ clearLoops()

MoFEMErrorCode MoFEM::TsCtx::clearLoops ( )

Clear loops.

Returns
MoFEMErrorCode

Definition at line 17 of file TsCtx.cpp.

17  {
19  loopsIJacobian.clear();
20  loopsIFunction.clear();
21  loopsMonitor.clear();
22  loopsRHSJacobian.clear();
23  loopsRHSFunction.clear();
24  preProcessIJacobian.clear();
25  postProcessIJacobian.clear();
26  preProcessIFunction.clear();
27  postProcessIFunction.clear();
28  preProcessMonitor.clear();
29  postProcessMonitor.clear();
30  preProcessRHSJacobian.clear();
31  preProcessRHSFunction.clear();
32  postProcessRHSJacobian.clear();
33  postProcessRHSFunction.clear();
35 }
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:460
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:453
FEMethodsSequence loopsIJacobian
Definition: TsCtx.hpp:39
FEMethodsSequence loopsRHSFunction
Definition: TsCtx.hpp:43
BasicMethodsSequence postProcessRHSFunction
Definition: TsCtx.hpp:53
BasicMethodsSequence preProcessMonitor
Definition: TsCtx.hpp:48
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
BasicMethodsSequence preProcessRHSFunction
Definition: TsCtx.hpp:51
BasicMethodsSequence postProcessIJacobian
Definition: TsCtx.hpp:45
BasicMethodsSequence postProcessRHSJacobian
Definition: TsCtx.hpp:52

◆ getLoopsIFunction()

FEMethodsSequence& MoFEM::TsCtx::getLoopsIFunction ( )

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 81 of file TsCtx.hpp.

81 { return loopsIFunction; }

◆ getLoopsIJacobian()

FEMethodsSequence& MoFEM::TsCtx::getLoopsIJacobian ( )

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 101 of file TsCtx.hpp.

101 { return loopsIJacobian; }

◆ getLoopsMonitor()

FEMethodsSequence& MoFEM::TsCtx::getLoopsMonitor ( )

Get the loops to do Monitor object.

It is sequence used to monitor solution of time solver.

Returns
FEMethodsSequence&

Definition at line 120 of file TsCtx.hpp.

120 { return loopsMonitor; }

◆ getLoopsRHSFunction()

FEMethodsSequence& MoFEM::TsCtx::getLoopsRHSFunction ( )

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 91 of file TsCtx.hpp.

91 { return loopsRHSFunction; }

◆ getLoopsRHSJacobian()

FEMethodsSequence& MoFEM::TsCtx::getLoopsRHSJacobian ( )

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 111 of file TsCtx.hpp.

111 { return loopsRHSJacobian; }

◆ getPostProcessIFunction()

BasicMethodsSequence& MoFEM::TsCtx::getPostProcessIFunction ( )

Get the postProcess to do IFunction object.

Returns
BasicMethodsSequence&

Definition at line 134 of file TsCtx.hpp.

134  {
135  return postProcessIFunction;
136  }

◆ getPostProcessIJacobian()

BasicMethodsSequence& MoFEM::TsCtx::getPostProcessIJacobian ( )

Get the postProcess to do IJacobian object.

Returns
BasicMethodsSequence&

Definition at line 150 of file TsCtx.hpp.

150  {
151  return postProcessIJacobian;
152  }

◆ getPostProcessMonitor()

BasicMethodsSequence& MoFEM::TsCtx::getPostProcessMonitor ( )

Get the postProcess to do Monitor object.

Returns
BasicMethodsSequence&

Definition at line 166 of file TsCtx.hpp.

166 { return postProcessMonitor; }

◆ getPostProcessRHSFunction()

BasicMethodsSequence& MoFEM::TsCtx::getPostProcessRHSFunction ( )

Get the postProcess to do RHSFunction object.

Returns
BasicMethodsSequence&

Definition at line 200 of file TsCtx.hpp.

200  {
201  return postProcessRHSFunction;
202  }

◆ getPostProcessRHSJacobian()

BasicMethodsSequence& MoFEM::TsCtx::getPostProcessRHSJacobian ( )

Get the postProcess to do RHSJacobian object.

Returns
BasicMethodsSequence&

Definition at line 182 of file TsCtx.hpp.

182  {
183  return postProcessRHSJacobian;
184  }

◆ getPreProcessIFunction()

BasicMethodsSequence& MoFEM::TsCtx::getPreProcessIFunction ( )

Get the preProcess to do IFunction object.

Returns
BasicMethodsSequence&

Definition at line 127 of file TsCtx.hpp.

127 { return preProcessIFunction; }

◆ getPreProcessIJacobian()

BasicMethodsSequence& MoFEM::TsCtx::getPreProcessIJacobian ( )

Get the preProcess to do IJacobian object.

Returns
BasicMethodsSequence&

Definition at line 143 of file TsCtx.hpp.

143 { return preProcessIJacobian; }

◆ getPreProcessMonitor()

BasicMethodsSequence& MoFEM::TsCtx::getPreProcessMonitor ( )

Get the preProcess to do Monitor object.

Returns
BasicMethodsSequence&

Definition at line 159 of file TsCtx.hpp.

159 { return preProcessMonitor; }

◆ getPreProcessRHSFunction()

BasicMethodsSequence& MoFEM::TsCtx::getPreProcessRHSFunction ( )

Get the preProcess to do RHSFunction object.

Returns
BasicMethodsSequence&

Definition at line 191 of file TsCtx.hpp.

191  {
192  return preProcessRHSFunction;
193  }

◆ getPreProcessRHSJacobian()

BasicMethodsSequence& MoFEM::TsCtx::getPreProcessRHSJacobian ( )

Get the preProcess to do RHSJacobian object.

Returns
BasicMethodsSequence&

Definition at line 173 of file TsCtx.hpp.

173  {
174  return preProcessRHSJacobian;
175  }

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 229 of file TsCtx.cpp.

230  {
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;
245  fe.ts_ctx = TSMethod::CTX_TSTSMONITORSET;
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);
259  CHKERR ts_ctx->mField.problem_basic_method_preProcess(ts_ctx->problemName,
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);
278  CHKERR ts_ctx->mField.problem_basic_method_postProcess(ts_ctx->problemName,
279  *bit);
280  unset(*bit);
281  }
282 
283  PetscLogEventEnd(ts_ctx->MOFEM_EVENT_TsCtxMonitor, 0, 0, 0, 0);
285 }
@ COL
Definition: definitions.h:136
#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
auto bit
set bit
constexpr double t
plate stiffness
Definition: plate.cpp:76
static constexpr Switches CtxSetX
Definition: LoopMethods.hpp:49
static constexpr Switches CtxSetNone
Definition: LoopMethods.hpp:45
static constexpr Switches CtxSetTime
Definition: LoopMethods.hpp:52
TsCtx(MoFEM::Interface &m_field, const std::string &problem_name)
Definition: TsCtx.hpp:57

◆ 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
t
u
u_t
u_tt
F
ctx
Returns
PetscErrorCode

Definition at line 569 of file TsCtx.cpp.

570  {
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;
615  fe.ts_ctx = TSMethod::CTX_TSSETIFUNCTION;
616  fe.snes_ctx = SnesMethod::CTX_SNESSETFUNCTION;
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);
635  CHKERR ts_ctx->mField.problem_basic_method_preProcess(ts_ctx->problemName,
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);
659  CHKERR ts_ctx->mField.problem_basic_method_postProcess(ts_ctx->problemName,
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 }
constexpr double a
FTensor::Index< 'i', SPACE_DIM > i
FTensor::Index< 'l', 3 > l
const FTensor::Tensor2< T, Dim, Dim > Vec
constexpr double g
static constexpr Switches CtxSetX_TT
Definition: LoopMethods.hpp:51
static constexpr Switches CtxSetF
Definition: LoopMethods.hpp:46
static constexpr Switches CtxSetX_T
Definition: LoopMethods.hpp:50

◆ 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.

See PETSc for details

Parameters
ts
ttime at step/stage being solved
ustate vectora
u_ttime derivative of state vector
u_ttsecond time derivative of state vector
ashift for u_t
aashift for u_tt
AJacobian of G(U) = F(t,U,W+v*U,W'+a*U), equivalent to dF/dU + v*dF/dU_t + a*dF/dU_tt
Bpreconditioning matrix for J, may be same as J
ctxTsCtx context for matrix evaluation routine
Returns
PetscErrorCode

Definition at line 473 of file TsCtx.cpp.

475  {
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 
512  fe.ts_ctx = TSMethod::CTX_TSSETIJACOBIAN;
513  fe.snes_ctx = SnesMethod::CTX_SNESSETJACOBIAN;
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);
532  CHKERR ts_ctx->mField.problem_basic_method_preProcess(ts_ctx->problemName,
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);
555  CHKERR ts_ctx->mField.problem_basic_method_postProcess(ts_ctx->problemName,
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 }
double A
static constexpr Switches CtxSetA
Definition: LoopMethods.hpp:47
static constexpr Switches CtxSetB
Definition: LoopMethods.hpp:48

◆ 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 37 of file TsCtx.cpp.

38  {
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;
81  fe.ts_ctx = TSMethod::CTX_TSSETIFUNCTION;
82  fe.snes_ctx = SnesMethod::CTX_SNESSETFUNCTION;
83  fe.ksp_ctx = KspMethod::CTX_SETFUNCTION;
84  fe.data_ctx = PetscData::CtxSetF | PetscData::CtxSetX |
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);
99  CHKERR ts_ctx->mField.problem_basic_method_preProcess(ts_ctx->problemName,
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);
123  CHKERR ts_ctx->mField.problem_basic_method_postProcess(ts_ctx->problemName,
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 }

◆ 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 140 of file TsCtx.cpp.

141  {
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;
174  fe.ts_ctx = TSMethod::CTX_TSSETIJACOBIAN;
175  fe.snes_ctx = SnesMethod::CTX_SNESSETJACOBIAN;
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);
192  CHKERR ts_ctx->mField.problem_basic_method_preProcess(ts_ctx->problemName,
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);
215  CHKERR ts_ctx->mField.problem_basic_method_postProcess(ts_ctx->problemName,
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 }

◆ 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 287 of file TsCtx.cpp.

287  {
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;
327  fe.ts_ctx = TSMethod::CTX_TSSETRHSFUNCTION;
328  fe.snes_ctx = SnesMethod::CTX_SNESSETFUNCTION;
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);
344  CHKERR ts_ctx->mField.problem_basic_method_preProcess(ts_ctx->problemName,
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);
368  CHKERR ts_ctx->mField.problem_basic_method_postProcess(ts_ctx->problemName,
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 }

◆ 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 385 of file TsCtx.cpp.

386  {
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;
415  fe.ts_ctx = TSMethod::CTX_TSSETRHSJACOBIAN;
416  fe.snes_ctx = SnesMethod::CTX_SNESSETJACOBIAN;
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);
434  CHKERR ts_ctx->mField.problem_basic_method_preProcess(ts_ctx->problemName,
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);
458  CHKERR ts_ctx->mField.problem_basic_method_postProcess(ts_ctx->problemName,
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 }

Member Data Documentation

◆ bH

MoFEMTypes MoFEM::TsCtx::bH

If set to MF_EXIST check if element exist.

Definition at line 33 of file TsCtx.hpp.

◆ loopsIFunction

FEMethodsSequence MoFEM::TsCtx::loopsIFunction

Definition at line 40 of file TsCtx.hpp.

◆ loopsIJacobian

FEMethodsSequence MoFEM::TsCtx::loopsIJacobian

Definition at line 39 of file TsCtx.hpp.

◆ loopsMonitor

FEMethodsSequence MoFEM::TsCtx::loopsMonitor

Definition at line 41 of file TsCtx.hpp.

◆ loopsRHSFunction

FEMethodsSequence MoFEM::TsCtx::loopsRHSFunction

Definition at line 43 of file TsCtx.hpp.

◆ loopsRHSJacobian

FEMethodsSequence MoFEM::TsCtx::loopsRHSJacobian

Definition at line 42 of file TsCtx.hpp.

◆ matAssembleSwitch

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

Definition at line 239 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 235 of file TsCtx.hpp.

◆ MOFEM_EVENT_TsCtxI2Jacobian

PetscLogEvent MoFEM::TsCtx::MOFEM_EVENT_TsCtxI2Jacobian
private

Definition at line 236 of file TsCtx.hpp.

◆ MOFEM_EVENT_TsCtxIFunction

PetscLogEvent MoFEM::TsCtx::MOFEM_EVENT_TsCtxIFunction
private

Definition at line 232 of file TsCtx.hpp.

◆ MOFEM_EVENT_TsCtxIJacobian

PetscLogEvent MoFEM::TsCtx::MOFEM_EVENT_TsCtxIJacobian
private

Definition at line 233 of file TsCtx.hpp.

◆ MOFEM_EVENT_TsCtxMonitor

PetscLogEvent MoFEM::TsCtx::MOFEM_EVENT_TsCtxMonitor
private

Definition at line 234 of file TsCtx.hpp.

◆ MOFEM_EVENT_TsCtxRHSFunction

PetscLogEvent MoFEM::TsCtx::MOFEM_EVENT_TsCtxRHSFunction
private

Definition at line 230 of file TsCtx.hpp.

◆ MOFEM_EVENT_TsCtxRHSJacobian

PetscLogEvent MoFEM::TsCtx::MOFEM_EVENT_TsCtxRHSJacobian
private

Definition at line 231 of file TsCtx.hpp.

◆ postProcessIFunction

BasicMethodsSequence MoFEM::TsCtx::postProcessIFunction

Definition at line 47 of file TsCtx.hpp.

◆ postProcessIJacobian

BasicMethodsSequence MoFEM::TsCtx::postProcessIJacobian

Definition at line 45 of file TsCtx.hpp.

◆ postProcessMonitor

BasicMethodsSequence MoFEM::TsCtx::postProcessMonitor

Definition at line 49 of file TsCtx.hpp.

◆ postProcessRHSFunction

BasicMethodsSequence MoFEM::TsCtx::postProcessRHSFunction

Definition at line 53 of file TsCtx.hpp.

◆ postProcessRHSJacobian

BasicMethodsSequence MoFEM::TsCtx::postProcessRHSJacobian

Definition at line 52 of file TsCtx.hpp.

◆ preProcessIFunction

BasicMethodsSequence MoFEM::TsCtx::preProcessIFunction

Definition at line 46 of file TsCtx.hpp.

◆ preProcessIJacobian

BasicMethodsSequence MoFEM::TsCtx::preProcessIJacobian

Definition at line 44 of file TsCtx.hpp.

◆ preProcessMonitor

BasicMethodsSequence MoFEM::TsCtx::preProcessMonitor

Definition at line 48 of file TsCtx.hpp.

◆ preProcessRHSFunction

BasicMethodsSequence MoFEM::TsCtx::preProcessRHSFunction

Definition at line 51 of file TsCtx.hpp.

◆ preProcessRHSJacobian

BasicMethodsSequence MoFEM::TsCtx::preProcessRHSJacobian

Definition at line 50 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 238 of file TsCtx.hpp.

◆ zeroMatrix

bool MoFEM::TsCtx::zeroMatrix

Definition at line 55 of file TsCtx.hpp.


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