v0.14.0
Loading...
Searching...
No Matches
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.

Examples
Remodeling.cpp, level_set.cpp, and nonlinear_dynamics.cpp.

Definition at line 15 of file TsCtx.hpp.

Member Typedef Documentation

◆ BasicMethodsSequence

Definition at line 25 of file TsCtx.hpp.

◆ FEMethodsSequence

Definition at line 24 of file TsCtx.hpp.

◆ PairNameFEMethodPtr

Definition at line 23 of file TsCtx.hpp.

Constructor & Destructor Documentation

◆ TsCtx()

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

Definition at line 45 of file TsCtx.hpp.

46 : mField(m_field), moab(m_field.get_moab()), problemName(problem_name),
47 bH(MF_EXIST), zeroMatrix(true) {
48 PetscLogEventRegister("LoopTsIFunction", 0, &MOFEM_EVENT_TsCtxIFunction);
49 PetscLogEventRegister("LoopTsIJacobian", 0, &MOFEM_EVENT_TsCtxIJacobian);
50 PetscLogEventRegister("LoopTsRHSFunction", 0,
52 PetscLogEventRegister("LoopTsRHSJacobian", 0,
54 PetscLogEventRegister("LoopTsMonitor", 0, &MOFEM_EVENT_TsCtxMonitor);
55 PetscLogEventRegister("LoopTsI2Function", 0, &MOFEM_EVENT_TsCtxI2Function);
56 PetscLogEventRegister("LoopTsI2Jacobian", 0, &MOFEM_EVENT_TsCtxI2Jacobian);
57 }
@ MF_EXIST
Definition: definitions.h:100
virtual moab::Interface & get_moab()=0
PetscLogEvent MOFEM_EVENT_TsCtxI2Jacobian
Definition: TsCtx.hpp:224
bool zeroMatrix
Definition: TsCtx.hpp:43
MoFEMTypes bH
If set to MF_EXIST check if element exist.
Definition: TsCtx.hpp:21
moab::Interface & moab
Definition: TsCtx.hpp:18
PetscLogEvent MOFEM_EVENT_TsCtxRHSJacobian
Definition: TsCtx.hpp:219
PetscLogEvent MOFEM_EVENT_TsCtxMonitor
Definition: TsCtx.hpp:222
PetscLogEvent MOFEM_EVENT_TsCtxRHSFunction
Definition: TsCtx.hpp:218
PetscLogEvent MOFEM_EVENT_TsCtxIJacobian
Definition: TsCtx.hpp:221
std::string problemName
Definition: TsCtx.hpp:20
MoFEM::Interface & mField
Definition: TsCtx.hpp:17
PetscLogEvent MOFEM_EVENT_TsCtxIFunction
Definition: TsCtx.hpp:220
PetscLogEvent MOFEM_EVENT_TsCtxI2Function
Definition: TsCtx.hpp:223

◆ ~TsCtx()

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

Member Function Documentation

◆ clearLoops()

MoFEMErrorCode MoFEM::TsCtx::clearLoops ( )

Clear loops.

Returns
MoFEMErrorCode

Definition at line 5 of file TsCtx.cpp.

5 {
7 loopsIJacobian.clear();
8 loopsIFunction.clear();
9 loopsMonitor.clear();
10 loopsRHSJacobian.clear();
11 loopsRHSFunction.clear();
12 preProcessIJacobian.clear();
14 preProcessIFunction.clear();
16 preProcessMonitor.clear();
17 postProcessMonitor.clear();
23}
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:447
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:440
FEMethodsSequence loopsIJacobian
Definition: TsCtx.hpp:27
FEMethodsSequence loopsRHSFunction
Definition: TsCtx.hpp:31
BasicMethodsSequence postProcessRHSFunction
Definition: TsCtx.hpp:41
BasicMethodsSequence preProcessMonitor
Definition: TsCtx.hpp:36
BasicMethodsSequence preProcessIJacobian
Definition: TsCtx.hpp:32
BasicMethodsSequence postProcessIFunction
Definition: TsCtx.hpp:35
FEMethodsSequence loopsIFunction
Definition: TsCtx.hpp:28
FEMethodsSequence loopsRHSJacobian
Definition: TsCtx.hpp:30
BasicMethodsSequence preProcessIFunction
Definition: TsCtx.hpp:34
FEMethodsSequence loopsMonitor
Definition: TsCtx.hpp:29
BasicMethodsSequence postProcessMonitor
Definition: TsCtx.hpp:37
BasicMethodsSequence preProcessRHSJacobian
Definition: TsCtx.hpp:38
BasicMethodsSequence preProcessRHSFunction
Definition: TsCtx.hpp:39
BasicMethodsSequence postProcessIJacobian
Definition: TsCtx.hpp:33
BasicMethodsSequence postProcessRHSJacobian
Definition: TsCtx.hpp:40

◆ getLoopsIFunction()

FEMethodsSequence & MoFEM::TsCtx::getLoopsIFunction ( )
inline

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&
Examples
nonlinear_dynamics.cpp.

Definition at line 69 of file TsCtx.hpp.

69{ return loopsIFunction; }

◆ getLoopsIJacobian()

FEMethodsSequence & MoFEM::TsCtx::getLoopsIJacobian ( )
inline

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

89{ return loopsIJacobian; }

◆ getLoopsMonitor()

FEMethodsSequence & MoFEM::TsCtx::getLoopsMonitor ( )
inline

Get the loops to do Monitor object.

It is sequence used to monitor solution of time solver.

Returns
FEMethodsSequence&
Examples
EshelbianPlasticity.cpp, and nonlinear_dynamics.cpp.

Definition at line 108 of file TsCtx.hpp.

108{ return loopsMonitor; }

◆ getLoopsRHSFunction()

FEMethodsSequence & MoFEM::TsCtx::getLoopsRHSFunction ( )
inline

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

79{ return loopsRHSFunction; }

◆ getLoopsRHSJacobian()

FEMethodsSequence & MoFEM::TsCtx::getLoopsRHSJacobian ( )
inline

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

99{ return loopsRHSJacobian; }

◆ getPostProcessIFunction()

BasicMethodsSequence & MoFEM::TsCtx::getPostProcessIFunction ( )
inline

Get the postProcess to do IFunction object.

Returns
BasicMethodsSequence&
Examples
nonlinear_dynamics.cpp.

Definition at line 122 of file TsCtx.hpp.

122 {
124 }

◆ getPostProcessIJacobian()

BasicMethodsSequence & MoFEM::TsCtx::getPostProcessIJacobian ( )
inline

Get the postProcess to do IJacobian object.

Returns
BasicMethodsSequence&
Examples
nonlinear_dynamics.cpp.

Definition at line 138 of file TsCtx.hpp.

138 {
140 }

◆ getPostProcessMonitor()

BasicMethodsSequence & MoFEM::TsCtx::getPostProcessMonitor ( )
inline

Get the postProcess to do Monitor object.

Returns
BasicMethodsSequence&
Examples
Remodeling.cpp, and UnsaturatedFlow.hpp.

Definition at line 154 of file TsCtx.hpp.

154{ return postProcessMonitor; }

◆ getPostProcessRHSFunction()

BasicMethodsSequence & MoFEM::TsCtx::getPostProcessRHSFunction ( )
inline

Get the postProcess to do RHSFunction object.

Returns
BasicMethodsSequence&

Definition at line 188 of file TsCtx.hpp.

188 {
190 }

◆ getPostProcessRHSJacobian()

BasicMethodsSequence & MoFEM::TsCtx::getPostProcessRHSJacobian ( )
inline

Get the postProcess to do RHSJacobian object.

Returns
BasicMethodsSequence&

Definition at line 170 of file TsCtx.hpp.

170 {
172 }

◆ getPreProcessIFunction()

BasicMethodsSequence & MoFEM::TsCtx::getPreProcessIFunction ( )
inline

Get the preProcess to do IFunction object.

Returns
BasicMethodsSequence&
Examples
nonlinear_dynamics.cpp.

Definition at line 115 of file TsCtx.hpp.

115{ return preProcessIFunction; }

◆ getPreProcessIJacobian()

BasicMethodsSequence & MoFEM::TsCtx::getPreProcessIJacobian ( )
inline

Get the preProcess to do IJacobian object.

Returns
BasicMethodsSequence&
Examples
nonlinear_dynamics.cpp.

Definition at line 131 of file TsCtx.hpp.

131{ return preProcessIJacobian; }

◆ getPreProcessMonitor()

BasicMethodsSequence & MoFEM::TsCtx::getPreProcessMonitor ( )
inline

Get the preProcess to do Monitor object.

Returns
BasicMethodsSequence&

Definition at line 147 of file TsCtx.hpp.

147{ return preProcessMonitor; }

◆ getPreProcessRHSFunction()

BasicMethodsSequence & MoFEM::TsCtx::getPreProcessRHSFunction ( )
inline

Get the preProcess to do RHSFunction object.

Returns
BasicMethodsSequence&

Definition at line 179 of file TsCtx.hpp.

179 {
181 }

◆ getPreProcessRHSJacobian()

BasicMethodsSequence & MoFEM::TsCtx::getPreProcessRHSJacobian ( )
inline

Get the preProcess to do RHSJacobian object.

Returns
BasicMethodsSequence&

Definition at line 161 of file TsCtx.hpp.

161 {
163 }

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

222 {
224 TsCtx *ts_ctx = static_cast<TsCtx *>(ctx);
225 PetscLogEventBegin(ts_ctx->MOFEM_EVENT_TsCtxMonitor, 0, 0, 0, 0);
226 CHKERR VecGhostUpdateBegin(u, INSERT_VALUES, SCATTER_FORWARD);
227 CHKERR VecGhostUpdateEnd(u, INSERT_VALUES, SCATTER_FORWARD);
228 CHKERR ts_ctx->mField.getInterface<VecManager>()->setLocalGhostVector(
229 ts_ctx->problemName, COL, u, INSERT_VALUES, SCATTER_REVERSE);
230
231 auto cache_ptr = boost::make_shared<CacheTuple>();
233
234 auto set = [&](auto &fe) {
235 fe.ts = ts;
236 fe.ts_u = u;
237 fe.ts_t = t;
238 fe.ts_step = step;
239 fe.ts_F = PETSC_NULL;
241 fe.snes_ctx = SnesMethod::CTX_SNESNONE;
242 fe.ksp_ctx = KspMethod::CTX_KSPNONE;
244 fe.cacheWeakPtr = cache_ptr;
245 CHK_THROW_MESSAGE(TSGetTimeStep(ts, &fe.ts_dt), "get time step failed");
246 };
247
248 auto unset = [&](auto &fe) {
249 fe.ts_ctx = TSMethod::CTX_TSNONE;
250 fe.data_ctx = PetscData::CtxSetNone;
251 };
252
253 // preprocess
254 for (auto &bit : ts_ctx->preProcessMonitor) {
255 set(*bit);
257 *bit);
258 unset(*bit);
259 }
260
261 for (auto &lit : ts_ctx->loopsMonitor) {
262 set(*lit.second);
264 *(lit.second), nullptr,
265 ts_ctx->bH, cache_ptr);
266 unset(*lit.second);
267 }
268
269 // post process
270 for (auto &bit : ts_ctx->postProcessMonitor) {
271 set(*bit);
273 *bit);
274 unset(*bit);
275 }
276
277 PetscLogEventEnd(ts_ctx->MOFEM_EVENT_TsCtxMonitor, 0, 0, 0, 0);
279}
@ COL
Definition: definitions.h:123
#define CHK_THROW_MESSAGE(err, msg)
Check and throw MoFEM exception.
Definition: definitions.h:595
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:346
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:416
#define CHKERR
Inline error check.
Definition: definitions.h:535
virtual MoFEMErrorCode problem_basic_method_postProcess(const Problem *problem_ptr, BasicMethod &method, int verb=DEFAULT_VERBOSITY)=0
Set data for BasicMethod.
virtual MoFEMErrorCode loop_finite_elements(const std::string problem_name, const std::string &fe_name, FEMethod &method, boost::shared_ptr< NumeredEntFiniteElement_multiIndex > fe_ptr=nullptr, MoFEMTypes bh=MF_EXIST, CacheTupleWeakPtr cache_ptr=CacheTupleSharedPtr(), int verb=DEFAULT_VERBOSITY)=0
Make a loop over finite elements.
auto bit
set bit
MoFEM::TsCtx * ts_ctx
Definition: level_set.cpp:1884
constexpr double t
plate stiffness
Definition: plate.cpp:59
virtual MoFEMErrorCode cache_problem_entities(const std::string prb_name, CacheTupleWeakPtr cache_ptr)=0
Cache variables.
virtual MoFEMErrorCode problem_basic_method_preProcess(const Problem *problem_ptr, BasicMethod &method, int verb=DEFAULT_VERBOSITY)=0
Set data for BasicMethod.
static constexpr Switches CtxSetX
Definition: LoopMethods.hpp:39
static constexpr Switches CtxSetNone
Definition: LoopMethods.hpp:35
static constexpr Switches CtxSetTime
Definition: LoopMethods.hpp:42
TsCtx(MoFEM::Interface &m_field, const std::string &problem_name)
Definition: TsCtx.hpp:45
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface refernce to pointer of interface.

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

567 {
569 TsCtx *ts_ctx = static_cast<TsCtx *>(ctx);
570 PetscLogEventBegin(ts_ctx->MOFEM_EVENT_TsCtxIFunction, 0, 0, 0, 0);
571 CHKERR VecGhostUpdateBegin(u, INSERT_VALUES, SCATTER_FORWARD);
572 CHKERR VecGhostUpdateEnd(u, INSERT_VALUES, SCATTER_FORWARD);
573 CHKERR VecGhostUpdateBegin(u_t, INSERT_VALUES, SCATTER_FORWARD);
574 CHKERR VecGhostUpdateEnd(u_t, INSERT_VALUES, SCATTER_FORWARD);
575 CHKERR VecGhostUpdateBegin(u_tt, INSERT_VALUES, SCATTER_FORWARD);
576 CHKERR VecGhostUpdateEnd(u_tt, INSERT_VALUES, SCATTER_FORWARD);
577 CHKERR ts_ctx->mField.getInterface<VecManager>()->setLocalGhostVector(
578 ts_ctx->problemName, COL, u, INSERT_VALUES, SCATTER_REVERSE);
579
580 auto zero_ghost_vec = [](Vec g) {
582 Vec l;
583 CHKERR VecGhostGetLocalForm(g, &l);
584 double *a;
585 CHKERR VecGetArray(l, &a);
586 int s;
587 CHKERR VecGetLocalSize(l, &s);
588 for (int i = 0; i != s; ++i)
589 a[i] = 0;
590 CHKERR VecRestoreArray(l, &a);
591 CHKERR VecGhostRestoreLocalForm(g, &l);
593 };
594 CHKERR zero_ghost_vec(F);
595
596 ts_ctx->vecAssembleSwitch = boost::movelib::make_unique<bool>(true);
597 auto cache_ptr = boost::make_shared<CacheTuple>();
599
600 int step;
601#if PETSC_VERSION_GE(3, 8, 0)
602 CHKERR TSGetStepNumber(ts, &step);
603#else
604 CHKERR TSGetTimeStepNumber(ts, &step);
605#endif
606
607 auto set = [&](auto &fe) {
608 fe.ts_u = u;
609 fe.ts_u_t = u_t;
610 fe.ts_u_tt = u_tt;
611 fe.ts_F = F;
612 fe.ts_t = t;
613 fe.ts_step = step;
616 fe.ksp_ctx = KspMethod::CTX_SETFUNCTION;
617 fe.data_ctx = PetscData::CtxSetF | PetscData::CtxSetX |
620 fe.ts = ts;
621 fe.cacheWeakPtr = cache_ptr;
622 CHK_THROW_MESSAGE(TSGetTimeStep(ts, &fe.ts_dt), "get time step failed");
623 };
624
625 auto unset = [&](auto &fe) {
626 fe.ts_ctx = TSMethod::CTX_TSNONE;
627 fe.snes_ctx = SnesMethod::CTX_SNESNONE;
628 fe.ksp_ctx = KspMethod::CTX_KSPNONE;
629 fe.data_ctx = PetscData::CtxSetNone;
630 };
631
632 // preprocess
633 for (auto &bit : ts_ctx->preProcessIFunction) {
634 bit->vecAssembleSwitch = boost::move(ts_ctx->vecAssembleSwitch);
635 set(*bit);
637 *bit);
638 unset(*bit);
639 ts_ctx->vecAssembleSwitch = boost::move(bit->vecAssembleSwitch);
640 }
641
642 // fe loops
643 for (auto &lit : ts_ctx->loopsIFunction) {
644 lit.second->vecAssembleSwitch = boost::move(ts_ctx->vecAssembleSwitch);
645 set(*lit.second);
647 *(lit.second), nullptr,
648 ts_ctx->bH, cache_ptr);
649 unset(*lit.second);
650 ts_ctx->vecAssembleSwitch = boost::move(lit.second->vecAssembleSwitch);
651 }
652
653 // post process
654 for (auto &bit : ts_ctx->postProcessIFunction) {
655 bit->vecAssembleSwitch = boost::move(ts_ctx->vecAssembleSwitch);
656 set(*bit);
658 *bit);
659 unset(*bit);
660 ts_ctx->vecAssembleSwitch = boost::move(bit->vecAssembleSwitch);
661 }
662
664 CHKERR VecGhostUpdateBegin(F, ADD_VALUES, SCATTER_REVERSE);
665 CHKERR VecGhostUpdateEnd(F, ADD_VALUES, SCATTER_REVERSE);
666 CHKERR VecAssemblyBegin(F);
667 CHKERR VecAssemblyEnd(F);
668 }
669
670 PetscLogEventEnd(ts_ctx->MOFEM_EVENT_TsCtxIFunction, 0, 0, 0, 0);
672}
constexpr double a
@ F
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:41
static constexpr Switches CtxSetF
Definition: LoopMethods.hpp:36
static constexpr Switches CtxSetX_T
Definition: LoopMethods.hpp:40
boost::movelib::unique_ptr< bool > vecAssembleSwitch
Definition: TsCtx.hpp:226

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

471 {
473
474 TsCtx *ts_ctx = static_cast<TsCtx *>(ctx);
475 PetscLogEventBegin(ts_ctx->MOFEM_EVENT_TsCtxI2Function, 0, 0, 0, 0);
476 CHKERR VecGhostUpdateBegin(u, INSERT_VALUES, SCATTER_FORWARD);
477 CHKERR VecGhostUpdateEnd(u, INSERT_VALUES, SCATTER_FORWARD);
478 CHKERR VecGhostUpdateBegin(u_t, INSERT_VALUES, SCATTER_FORWARD);
479 CHKERR VecGhostUpdateEnd(u_t, INSERT_VALUES, SCATTER_FORWARD);
480 CHKERR VecGhostUpdateBegin(u_tt, INSERT_VALUES, SCATTER_FORWARD);
481 CHKERR VecGhostUpdateEnd(u_tt, INSERT_VALUES, SCATTER_FORWARD);
482 CHKERR ts_ctx->mField.getInterface<VecManager>()->setLocalGhostVector(
483 ts_ctx->problemName, COL, u, INSERT_VALUES, SCATTER_REVERSE);
484 if (ts_ctx->zeroMatrix) {
485 CHKERR MatZeroEntries(B);
486 }
487 int step;
488#if PETSC_VERSION_GE(3, 8, 0)
489 CHKERR TSGetStepNumber(ts, &step);
490#else
491 CHKERR TSGetTimeStepNumber(ts, &step);
492#endif
493
495 boost::movelib::make_unique<bool>(ts_ctx->zeroMatrix);
496 auto cache_ptr = boost::make_shared<CacheTuple>();
498
499 auto set = [&](auto &fe) {
500 fe.ts_u = u;
501 fe.ts_u_t = u_t;
502 fe.ts_u_tt = u_tt;
503 fe.ts_A = A;
504 fe.ts_B = B;
505 fe.ts_t = t;
506 fe.ts_a = a;
507 fe.ts_aa = aa;
508 fe.ts_step = step;
509
512 fe.ksp_ctx = KspMethod::CTX_OPERATORS;
516 fe.ts = ts;
517 fe.cacheWeakPtr = cache_ptr;
518 CHK_THROW_MESSAGE(TSGetTimeStep(ts, &fe.ts_dt), "get time step failed");
519 };
520
521 auto unset = [&](auto &fe) {
522 fe.ts_ctx = TSMethod::CTX_TSNONE;
523 fe.snes_ctx = SnesMethod::CTX_SNESNONE;
524 fe.ksp_ctx = KspMethod::CTX_KSPNONE;
525 fe.data_ctx = PetscData::CtxSetNone;
526 };
527
528 // preprocess
529 for (auto &bit : ts_ctx->preProcessIJacobian) {
530 bit->matAssembleSwitch = boost::move(ts_ctx->matAssembleSwitch);
531 set(*bit);
533 *bit);
534 unset(*bit);
535 ts_ctx->matAssembleSwitch = boost::move(bit->matAssembleSwitch);
536 }
537
538 for (auto &lit : ts_ctx->loopsIJacobian) {
539 lit.second->matAssembleSwitch = boost::move(ts_ctx->matAssembleSwitch);
540 set(*lit.second);
542 *(lit.second), nullptr,
543 ts_ctx->bH, cache_ptr);
544 unset(*lit.second);
545 ts_ctx->matAssembleSwitch = boost::move(lit.second->matAssembleSwitch);
546 }
547
548 // post process
549 for (auto &bit : ts_ctx->postProcessIJacobian) {
550 bit->matAssembleSwitch = boost::move(ts_ctx->matAssembleSwitch);
551 set(*bit);
553 *bit);
554 unset(*bit);
555 ts_ctx->matAssembleSwitch = boost::move(bit->matAssembleSwitch);
556 }
557
558 if (*(ts_ctx->matAssembleSwitch)) {
559 CHKERR MatAssemblyBegin(B, MAT_FINAL_ASSEMBLY);
560 CHKERR MatAssemblyEnd(B, MAT_FINAL_ASSEMBLY);
561 }
562 PetscLogEventEnd(ts_ctx->MOFEM_EVENT_TsCtxI2Function, 0, 0, 0, 0);
564}
constexpr AssemblyType A
static constexpr Switches CtxSetA
Definition: LoopMethods.hpp:37
static constexpr Switches CtxSetB
Definition: LoopMethods.hpp:38
boost::movelib::unique_ptr< bool > matAssembleSwitch
Definition: TsCtx.hpp:227

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

26 {
28 TsCtx *ts_ctx = static_cast<TsCtx *>(ctx);
29 PetscLogEventBegin(ts_ctx->MOFEM_EVENT_TsCtxIFunction, 0, 0, 0, 0);
30 CHKERR VecGhostUpdateBegin(u, INSERT_VALUES, SCATTER_FORWARD);
31 CHKERR VecGhostUpdateEnd(u, INSERT_VALUES, SCATTER_FORWARD);
32 CHKERR VecGhostUpdateBegin(u_t, INSERT_VALUES, SCATTER_FORWARD);
33 CHKERR VecGhostUpdateEnd(u_t, INSERT_VALUES, SCATTER_FORWARD);
34 CHKERR ts_ctx->mField.getInterface<VecManager>()->setLocalGhostVector(
35 ts_ctx->problemName, COL, u, INSERT_VALUES, SCATTER_REVERSE);
36
37 auto zero_ghost_vec = [](Vec g) {
39 Vec l;
40 CHKERR VecGhostGetLocalForm(g, &l);
41 double *a;
42 CHKERR VecGetArray(l, &a);
43 int s;
44 CHKERR VecGetLocalSize(l, &s);
45 for (int i = 0; i != s; ++i)
46 a[i] = 0;
47 CHKERR VecRestoreArray(l, &a);
48 CHKERR VecGhostRestoreLocalForm(g, &l);
50 };
51 CHKERR zero_ghost_vec(F);
52
53 ts_ctx->vecAssembleSwitch = boost::movelib::make_unique<bool>(true);
54
55 int step;
56#if PETSC_VERSION_GE(3, 8, 0)
57 CHKERR TSGetStepNumber(ts, &step);
58#else
59 CHKERR TSGetTimeStepNumber(ts, &step);
60#endif
61
62 auto cache_ptr = boost::make_shared<CacheTuple>();
64
65 auto set = [&](auto &fe) {
66 fe.ts = ts;
67 fe.ts_u = u;
68 fe.ts_u_t = u_t;
69 fe.ts_F = F;
70 fe.ts_t = t;
71 fe.ts_step = step;
74 fe.ksp_ctx = KspMethod::CTX_SETFUNCTION;
77 fe.cacheWeakPtr = cache_ptr;
78 CHK_THROW_MESSAGE(TSGetTimeStep(ts, &fe.ts_dt), "get time step failed");
79 };
80
81 auto unset = [&](auto &fe) {
82 fe.ts_ctx = TSMethod::CTX_TSNONE;
83 fe.snes_ctx = SnesMethod::CTX_SNESNONE;
84 fe.ksp_ctx = KspMethod::CTX_KSPNONE;
85 fe.data_ctx = PetscData::CtxSetNone;
86 };
87
88 // preprocess
89 for (auto &bit : ts_ctx->preProcessIFunction) {
90 bit->vecAssembleSwitch = boost::move(ts_ctx->vecAssembleSwitch);
91 set(*bit);
93 *bit);
94 unset(*bit);
95 ts_ctx->vecAssembleSwitch = boost::move(bit->vecAssembleSwitch);
96 }
97
98
99 // fe loops
100 for (auto &lit : ts_ctx->loopsIFunction) {
101 lit.second->vecAssembleSwitch = boost::move(ts_ctx->vecAssembleSwitch);
102 set(*lit.second);
104 *(lit.second), nullptr,
105 ts_ctx->bH, cache_ptr);
106 unset(*lit.second);
107 ts_ctx->vecAssembleSwitch = boost::move(lit.second->vecAssembleSwitch);
108 }
109
110 // post process
111 for (auto &bit : ts_ctx->postProcessIFunction) {
112 bit->vecAssembleSwitch = boost::move(ts_ctx->vecAssembleSwitch);
113 set(*bit);
115 *bit);
116 unset(*bit);
117 ts_ctx->vecAssembleSwitch = boost::move(bit->vecAssembleSwitch);
118 }
119
121 CHKERR VecGhostUpdateBegin(F, ADD_VALUES, SCATTER_REVERSE);
122 CHKERR VecGhostUpdateEnd(F, ADD_VALUES, SCATTER_REVERSE);
123 CHKERR VecAssemblyBegin(F);
124 CHKERR VecAssemblyEnd(F);
125 }
126
127 PetscLogEventEnd(ts_ctx->MOFEM_EVENT_TsCtxIFunction, 0, 0, 0, 0);
129}

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

132 {
134
135 TsCtx *ts_ctx = static_cast<TsCtx *>(ctx);
136 PetscLogEventBegin(ts_ctx->MOFEM_EVENT_TsCtxIFunction, 0, 0, 0, 0);
137 CHKERR VecGhostUpdateBegin(u, INSERT_VALUES, SCATTER_FORWARD);
138 CHKERR VecGhostUpdateEnd(u, INSERT_VALUES, SCATTER_FORWARD);
139 CHKERR VecGhostUpdateBegin(u_t, INSERT_VALUES, SCATTER_FORWARD);
140 CHKERR VecGhostUpdateEnd(u_t, INSERT_VALUES, SCATTER_FORWARD);
141 CHKERR ts_ctx->mField.getInterface<VecManager>()->setLocalGhostVector(
142 ts_ctx->problemName, COL, u, INSERT_VALUES, SCATTER_REVERSE);
143 if (ts_ctx->zeroMatrix) {
144 CHKERR MatZeroEntries(B);
145 }
146 int step;
147#if PETSC_VERSION_GE(3, 8, 0)
148 CHKERR TSGetStepNumber(ts, &step);
149#else
150 CHKERR TSGetTimeStepNumber(ts, &step);
151#endif
152
154 boost::movelib::make_unique<bool>(ts_ctx->zeroMatrix);
155 auto cache_ptr = boost::make_shared<CacheTuple>();
157
158 auto set = [&](auto &fe) {
159 fe.ts = ts;
160 fe.ts_u = u;
161 fe.ts_u_t = u_t;
162 fe.ts_A = A;
163 fe.ts_B = B;
164 fe.ts_t = t;
165 fe.ts_a = a;
166 fe.ts_step = step;
169 fe.ksp_ctx = KspMethod::CTX_OPERATORS;
172 fe.cacheWeakPtr = cache_ptr;
173 CHK_THROW_MESSAGE(TSGetTimeStep(ts, &fe.ts_dt), "get time step failed");
174 };
175
176 auto unset = [&](auto &fe) {
177 fe.ts_ctx = TSMethod::CTX_TSNONE;
178 fe.snes_ctx = SnesMethod::CTX_SNESNONE;
179 fe.ksp_ctx = KspMethod::CTX_KSPNONE;
180 fe.data_ctx = PetscData::CtxSetNone;
181 };
182
183 // preproces
184 for (auto &bit : ts_ctx->preProcessIJacobian) {
185 bit->matAssembleSwitch = boost::move(ts_ctx->matAssembleSwitch);
186 set(*bit);
188 *bit);
189 unset(*bit);
190 ts_ctx->matAssembleSwitch = boost::move(bit->matAssembleSwitch);
191 }
192
193 for (auto &lit : ts_ctx->loopsIJacobian) {
194 lit.second->matAssembleSwitch = boost::move(ts_ctx->matAssembleSwitch);
195 set(*lit.second);
197 *(lit.second), nullptr,
198 ts_ctx->bH, cache_ptr);
199 unset(*lit.second);
200 ts_ctx->matAssembleSwitch = boost::move(lit.second->matAssembleSwitch);
201 }
202
203 // post process
204 for (auto &bit : ts_ctx->postProcessIJacobian) {
205 bit->matAssembleSwitch = boost::move(ts_ctx->matAssembleSwitch);
206 set(*bit);
208 *bit);
209 unset(*bit);
210 ts_ctx->matAssembleSwitch = boost::move(bit->matAssembleSwitch);
211 }
212
213 if (*(ts_ctx->matAssembleSwitch)) {
214 CHKERR MatAssemblyBegin(B, MAT_FINAL_ASSEMBLY);
215 CHKERR MatAssemblyEnd(B, MAT_FINAL_ASSEMBLY);
216 }
217 PetscLogEventEnd(ts_ctx->MOFEM_EVENT_TsCtxIFunction, 0, 0, 0, 0);
219}

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

281 {
283 TsCtx *ts_ctx = static_cast<TsCtx *>(ctx);
284 PetscLogEventBegin(ts_ctx->MOFEM_EVENT_TsCtxRHSFunction, 0, 0, 0, 0);
285 CHKERR VecGhostUpdateBegin(u, INSERT_VALUES, SCATTER_FORWARD);
286 CHKERR VecGhostUpdateEnd(u, INSERT_VALUES, SCATTER_FORWARD);
287 CHKERR ts_ctx->mField.getInterface<VecManager>()->setLocalGhostVector(
288 ts_ctx->problemName, COL, u, INSERT_VALUES, SCATTER_REVERSE);
289
290 auto zero_ghost_vec = [](Vec g) {
292 Vec l;
293 CHKERR VecGhostGetLocalForm(g, &l);
294 double *a;
295 CHKERR VecGetArray(l, &a);
296 int s;
297 CHKERR VecGetLocalSize(l, &s);
298 for (int i = 0; i != s; ++i)
299 a[i] = 0;
300 CHKERR VecRestoreArray(l, &a);
301 CHKERR VecGhostRestoreLocalForm(g, &l);
303 };
304 CHKERR zero_ghost_vec(F);
305
306 ts_ctx->vecAssembleSwitch = boost::movelib::make_unique<bool>(true);
307 auto cache_ptr = boost::make_shared<CacheTuple>();
309
310 int step;
311#if PETSC_VERSION_GE(3, 8, 0)
312 CHKERR TSGetStepNumber(ts, &step);
313#else
314 CHKERR TSGetTimeStepNumber(ts, &step);
315#endif
316
317 auto set = [&](auto &fe) {
318 fe.ts_u = u;
319 fe.ts_F = F;
320 fe.ts_t = t;
321 fe.ts = ts;
322 fe.ts_step = step;
325 fe.ksp_ctx = KspMethod::CTX_SETFUNCTION;
326 fe.data_ctx =
328 fe.cacheWeakPtr = cache_ptr;
329 CHK_THROW_MESSAGE(TSGetTimeStep(ts, &fe.ts_dt), "get time step failed");
330 };
331
332 auto unset = [&](auto &fe) {
333 fe.ts_ctx = TSMethod::CTX_TSNONE;
334 fe.snes_ctx = SnesMethod::CTX_SNESNONE;
335 fe.ksp_ctx = KspMethod::CTX_KSPNONE;
336 fe.data_ctx = PetscData::CtxSetNone;
337 };
338
339 for (auto &bit : ts_ctx->preProcessRHSFunction) {
340 bit->vecAssembleSwitch = boost::move(ts_ctx->vecAssembleSwitch);
341 set(*bit);
343 *bit);
344 unset(*bit);
345 ts_ctx->vecAssembleSwitch = boost::move(bit->vecAssembleSwitch);
346 }
347
348 // fe loops
349 for (auto &lit : ts_ctx->loopsRHSFunction) {
350 lit.second->vecAssembleSwitch = boost::move(ts_ctx->vecAssembleSwitch);
351 set(*lit.second);
353 *(lit.second), nullptr,
354 ts_ctx->bH, cache_ptr);
355 unset(*lit.second);
356 ts_ctx->vecAssembleSwitch = boost::move(lit.second->vecAssembleSwitch);
357 }
358
359 // post process
360 for (auto &bit : ts_ctx->postProcessRHSFunction) {
361 bit->vecAssembleSwitch = boost::move(ts_ctx->vecAssembleSwitch);
362 set(*bit);
364 *bit);
365 unset(*bit);
366 ts_ctx->vecAssembleSwitch = boost::move(bit->vecAssembleSwitch);
367 }
368
370 CHKERR VecGhostUpdateBegin(F, ADD_VALUES, SCATTER_REVERSE);
371 CHKERR VecGhostUpdateEnd(F, ADD_VALUES, SCATTER_REVERSE);
372 CHKERR VecAssemblyBegin(F);
373 CHKERR VecAssemblyEnd(F);
374 }
375
376 PetscLogEventEnd(ts_ctx->MOFEM_EVENT_TsCtxRHSFunction, 0, 0, 0, 0);
378}

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

381 {
383 TsCtx *ts_ctx = static_cast<TsCtx *>(ctx);
384 PetscLogEventBegin(ts_ctx->MOFEM_EVENT_TsCtxRHSJacobian, 0, 0, 0, 0);
385 CHKERR VecGhostUpdateBegin(u, INSERT_VALUES, SCATTER_FORWARD);
386 CHKERR VecGhostUpdateEnd(u, INSERT_VALUES, SCATTER_FORWARD);
387 CHKERR ts_ctx->mField.getInterface<VecManager>()->setLocalGhostVector(
388 ts_ctx->problemName, COL, u, INSERT_VALUES, SCATTER_REVERSE);
389
390 if (ts_ctx->zeroMatrix) {
391 CHKERR MatZeroEntries(B);
392 }
393
395 boost::movelib::make_unique<bool>(ts_ctx->zeroMatrix);
396 auto cache_ptr = boost::make_shared<CacheTuple>();
398
399 int step;
400#if PETSC_VERSION_GE(3, 8, 0)
401 CHKERR TSGetStepNumber(ts, &step);
402#else
403 CHKERR TSGetTimeStepNumber(ts, &step);
404#endif
405
406 auto set = [&](auto &fe) {
407 fe.ts_u = u;
408 fe.ts_A = A;
409 fe.ts_B = B;
410 fe.ts_t = t;
411 fe.ts_step = step;
414 fe.ksp_ctx = KspMethod::CTX_OPERATORS;
417 fe.ts = ts;
418 fe.cacheWeakPtr = cache_ptr;
419 CHK_THROW_MESSAGE(TSGetTimeStep(ts, &fe.ts_dt), "get time step failed");
420 };
421
422 auto unset = [&](auto &fe) {
423 fe.ts_ctx = TSMethod::CTX_TSNONE;
424 fe.snes_ctx = SnesMethod::CTX_SNESNONE;
425 fe.ksp_ctx = KspMethod::CTX_KSPNONE;
426 fe.data_ctx = PetscData::CtxSetNone;
427 };
428
429 // preprocess
430 for (auto &bit : ts_ctx->preProcessRHSJacobian) {
431 bit->matAssembleSwitch = boost::move(ts_ctx->matAssembleSwitch);
432 set(*bit);
434 *bit);
435 unset(*bit);
436 ts_ctx->matAssembleSwitch = boost::move(bit->matAssembleSwitch);
437 }
438
439 // fe loops
440 for (auto &lit : ts_ctx->loopsRHSJacobian) {
441 lit.second->matAssembleSwitch = boost::move(ts_ctx->matAssembleSwitch);
442 set(*lit.second);
444 *(lit.second), nullptr,
445 ts_ctx->bH, cache_ptr);
446 unset(*lit.second);
447 ts_ctx->matAssembleSwitch = boost::move(lit.second->matAssembleSwitch);
448 }
449
450 // post process
451 for (auto &bit : ts_ctx->postProcessRHSJacobian) {
452 bit->matAssembleSwitch = boost::move(ts_ctx->matAssembleSwitch);
453 set(*bit);
455 *bit);
456 unset(*bit);
457 ts_ctx->matAssembleSwitch = boost::move(bit->matAssembleSwitch);
458 }
459
460 if (*(ts_ctx->matAssembleSwitch)) {
461 CHKERR MatAssemblyBegin(B, MAT_FINAL_ASSEMBLY);
462 CHKERR MatAssemblyEnd(B, MAT_FINAL_ASSEMBLY);
463 }
464
465 PetscLogEventEnd(ts_ctx->MOFEM_EVENT_TsCtxRHSJacobian, 0, 0, 0, 0);
467}

Member Data Documentation

◆ bH

MoFEMTypes MoFEM::TsCtx::bH

If set to MF_EXIST check if element exist.

Definition at line 21 of file TsCtx.hpp.

◆ loopsIFunction

FEMethodsSequence MoFEM::TsCtx::loopsIFunction

Definition at line 28 of file TsCtx.hpp.

◆ loopsIJacobian

FEMethodsSequence MoFEM::TsCtx::loopsIJacobian

Definition at line 27 of file TsCtx.hpp.

◆ loopsMonitor

FEMethodsSequence MoFEM::TsCtx::loopsMonitor

Definition at line 29 of file TsCtx.hpp.

◆ loopsRHSFunction

FEMethodsSequence MoFEM::TsCtx::loopsRHSFunction

Definition at line 31 of file TsCtx.hpp.

◆ loopsRHSJacobian

FEMethodsSequence MoFEM::TsCtx::loopsRHSJacobian

Definition at line 30 of file TsCtx.hpp.

◆ matAssembleSwitch

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

Definition at line 227 of file TsCtx.hpp.

◆ mField

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

Definition at line 17 of file TsCtx.hpp.

◆ moab

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

Definition at line 18 of file TsCtx.hpp.

◆ MOFEM_EVENT_TsCtxI2Function

PetscLogEvent MoFEM::TsCtx::MOFEM_EVENT_TsCtxI2Function
private

Definition at line 223 of file TsCtx.hpp.

◆ MOFEM_EVENT_TsCtxI2Jacobian

PetscLogEvent MoFEM::TsCtx::MOFEM_EVENT_TsCtxI2Jacobian
private

Definition at line 224 of file TsCtx.hpp.

◆ MOFEM_EVENT_TsCtxIFunction

PetscLogEvent MoFEM::TsCtx::MOFEM_EVENT_TsCtxIFunction
private

Definition at line 220 of file TsCtx.hpp.

◆ MOFEM_EVENT_TsCtxIJacobian

PetscLogEvent MoFEM::TsCtx::MOFEM_EVENT_TsCtxIJacobian
private

Definition at line 221 of file TsCtx.hpp.

◆ MOFEM_EVENT_TsCtxMonitor

PetscLogEvent MoFEM::TsCtx::MOFEM_EVENT_TsCtxMonitor
private

Definition at line 222 of file TsCtx.hpp.

◆ MOFEM_EVENT_TsCtxRHSFunction

PetscLogEvent MoFEM::TsCtx::MOFEM_EVENT_TsCtxRHSFunction
private

Definition at line 218 of file TsCtx.hpp.

◆ MOFEM_EVENT_TsCtxRHSJacobian

PetscLogEvent MoFEM::TsCtx::MOFEM_EVENT_TsCtxRHSJacobian
private

Definition at line 219 of file TsCtx.hpp.

◆ postProcessIFunction

BasicMethodsSequence MoFEM::TsCtx::postProcessIFunction

Definition at line 35 of file TsCtx.hpp.

◆ postProcessIJacobian

BasicMethodsSequence MoFEM::TsCtx::postProcessIJacobian

Definition at line 33 of file TsCtx.hpp.

◆ postProcessMonitor

BasicMethodsSequence MoFEM::TsCtx::postProcessMonitor

Definition at line 37 of file TsCtx.hpp.

◆ postProcessRHSFunction

BasicMethodsSequence MoFEM::TsCtx::postProcessRHSFunction

Definition at line 41 of file TsCtx.hpp.

◆ postProcessRHSJacobian

BasicMethodsSequence MoFEM::TsCtx::postProcessRHSJacobian

Definition at line 40 of file TsCtx.hpp.

◆ preProcessIFunction

BasicMethodsSequence MoFEM::TsCtx::preProcessIFunction

Definition at line 34 of file TsCtx.hpp.

◆ preProcessIJacobian

BasicMethodsSequence MoFEM::TsCtx::preProcessIJacobian

Definition at line 32 of file TsCtx.hpp.

◆ preProcessMonitor

BasicMethodsSequence MoFEM::TsCtx::preProcessMonitor

Definition at line 36 of file TsCtx.hpp.

◆ preProcessRHSFunction

BasicMethodsSequence MoFEM::TsCtx::preProcessRHSFunction

Definition at line 39 of file TsCtx.hpp.

◆ preProcessRHSJacobian

BasicMethodsSequence MoFEM::TsCtx::preProcessRHSJacobian

Definition at line 38 of file TsCtx.hpp.

◆ problemName

std::string MoFEM::TsCtx::problemName

Definition at line 20 of file TsCtx.hpp.

◆ vecAssembleSwitch

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

Definition at line 226 of file TsCtx.hpp.

◆ zeroMatrix

bool MoFEM::TsCtx::zeroMatrix

Definition at line 43 of file TsCtx.hpp.


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