v0.10.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)
 
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 209 of file TsCtx.cpp.

210  {
212  TsCtx *ts_ctx = static_cast<TsCtx *>(ctx);
213  PetscLogEventBegin(ts_ctx->MOFEM_EVENT_TsCtxMonitor, 0, 0, 0, 0);
214  CHKERR VecGhostUpdateBegin(u, INSERT_VALUES, SCATTER_FORWARD);
215  CHKERR VecGhostUpdateEnd(u, INSERT_VALUES, SCATTER_FORWARD);
216  CHKERR ts_ctx->mField.getInterface<VecManager>()->setLocalGhostVector(
217  ts_ctx->problemName, COL, u, INSERT_VALUES, SCATTER_REVERSE);
218 
219  auto set = [&](auto &fe) {
220  fe.ts = ts;
221  fe.ts_u = u;
222  fe.ts_t = t;
223  fe.ts_step = step;
224  fe.ts_F = PETSC_NULL;
225  fe.ts_ctx = TSMethod::CTX_TSTSMONITORSET;
226  fe.snes_ctx = SnesMethod::CTX_SNESNONE;
227  fe.ksp_ctx = KspMethod::CTX_KSPNONE;
229  };
230 
231  auto unset = [&](auto &fe) {
232  fe.ts_ctx = TSMethod::CTX_TSNONE;
233  fe.data_ctx = PetscData::CtxSetNone;
234  };
235 
236  auto cache_ptr = boost::make_shared<CacheTuple>();
237  CHKERR ts_ctx->mField.cache_problem_entities(ts_ctx->problemName, cache_ptr);
238 
239  // preproces
240  for (auto &bit : ts_ctx->preProcess_Monitor) {
241  set(*bit);
242  CHKERR ts_ctx->mField.problem_basic_method_preProcess(ts_ctx->problemName,
243  *bit);
244  unset(*bit);
245  }
246 
247  for (auto &lit : ts_ctx->loops_to_do_Monitor) {
248  set(*lit.second);
249  CHKERR ts_ctx->mField.loop_finite_elements(ts_ctx->problemName, lit.first,
250  *(lit.second), nullptr,
251  ts_ctx->bH, cache_ptr);
252  unset(*lit.second);
253  }
254 
255  // post process
256  for (auto &bit : ts_ctx->postProcess_Monitor) {
257  set(*bit);
258  CHKERR ts_ctx->mField.problem_basic_method_postProcess(ts_ctx->problemName,
259  *bit);
260  unset(*bit);
261  }
262 
263  PetscLogEventEnd(ts_ctx->MOFEM_EVENT_TsCtxMonitor, 0, 0, 0, 0);
265 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:485
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:604
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:415
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
t
u
u_t
u_tt
F
ctx
Returns
PetscErrorCode

Definition at line 549 of file TsCtx.cpp.

550  {
552  TsCtx *ts_ctx = static_cast<TsCtx *>(ctx);
553  PetscLogEventBegin(ts_ctx->MOFEM_EVENT_TsCtxIFunction, 0, 0, 0, 0);
554  CHKERR VecGhostUpdateBegin(u, INSERT_VALUES, SCATTER_FORWARD);
555  CHKERR VecGhostUpdateEnd(u, INSERT_VALUES, SCATTER_FORWARD);
556  CHKERR VecGhostUpdateBegin(u_t, INSERT_VALUES, SCATTER_FORWARD);
557  CHKERR VecGhostUpdateEnd(u_t, INSERT_VALUES, SCATTER_FORWARD);
558  CHKERR VecGhostUpdateBegin(u_tt, INSERT_VALUES, SCATTER_FORWARD);
559  CHKERR VecGhostUpdateEnd(u_tt, INSERT_VALUES, SCATTER_FORWARD);
560  CHKERR ts_ctx->mField.getInterface<VecManager>()->setLocalGhostVector(
561  ts_ctx->problemName, COL, u, INSERT_VALUES, SCATTER_REVERSE);
562 
563  auto zero_ghost_vec = [](Vec g) {
565  Vec l;
566  CHKERR VecGhostGetLocalForm(g, &l);
567  double *a;
568  CHKERR VecGetArray(l, &a);
569  int s;
570  CHKERR VecGetLocalSize(l, &s);
571  for (int i = 0; i != s; ++i)
572  a[i] = 0;
573  CHKERR VecRestoreArray(l, &a);
574  CHKERR VecGhostRestoreLocalForm(g, &l);
576  };
577  CHKERR zero_ghost_vec(F);
578 
579  ts_ctx->vecAssembleSwitch = boost::movelib::make_unique<bool>(true);
580 
581  int step;
582 #if PETSC_VERSION_GE(3, 8, 0)
583  CHKERR TSGetStepNumber(ts, &step);
584 #else
585  CHKERR TSGetTimeStepNumber(ts, &step);
586 #endif
587 
588  auto set = [&](auto &fe) {
589  fe.ts_u = u;
590  fe.ts_u_t = u_t;
591  fe.ts_u_tt = u_tt;
592  fe.ts_F = F;
593  fe.ts_t = t;
594  fe.ts_step = step;
595  fe.ts_ctx = TSMethod::CTX_TSSETIFUNCTION;
596  fe.snes_ctx = SnesMethod::CTX_SNESSETFUNCTION;
597  fe.ksp_ctx = KspMethod::CTX_SETFUNCTION;
598  fe.data_ctx = PetscData::CtxSetF | PetscData::CtxSetX |
601  fe.ts = ts;
602  };
603 
604  auto unset = [&](auto &fe) {
605  fe.ts_ctx = TSMethod::CTX_TSNONE;
606  fe.snes_ctx = SnesMethod::CTX_SNESNONE;
607  fe.ksp_ctx = KspMethod::CTX_KSPNONE;
608  fe.data_ctx = PetscData::CtxSetNone;
609  };
610 
611  auto cache_ptr = boost::make_shared<CacheTuple>();
612  CHKERR ts_ctx->mField.cache_problem_entities(ts_ctx->problemName, cache_ptr);
613 
614  // preprocess
615  for (auto &bit : ts_ctx->preProcess_IFunction) {
616  bit->vecAssembleSwitch = boost::move(ts_ctx->vecAssembleSwitch);
617  set(*bit);
618  CHKERR ts_ctx->mField.problem_basic_method_preProcess(ts_ctx->problemName,
619  *bit);
620  unset(*bit);
621  ts_ctx->vecAssembleSwitch = boost::move(bit->vecAssembleSwitch);
622  }
623 
624  // fe loops
625  for (auto &lit : ts_ctx->loops_to_do_IFunction) {
626  lit.second->vecAssembleSwitch = boost::move(ts_ctx->vecAssembleSwitch);
627  set(*lit.second);
628  CHKERR ts_ctx->mField.loop_finite_elements(ts_ctx->problemName, lit.first,
629  *(lit.second), nullptr,
630  ts_ctx->bH, cache_ptr);
631  unset(*lit.second);
632  ts_ctx->vecAssembleSwitch = boost::move(lit.second->vecAssembleSwitch);
633  }
634 
635  // post process
636  for (auto &bit : ts_ctx->postProcess_IFunction) {
637  bit->vecAssembleSwitch = boost::move(ts_ctx->vecAssembleSwitch);
638  set(*bit);
639  CHKERR ts_ctx->mField.problem_basic_method_postProcess(ts_ctx->problemName,
640  *bit);
641  unset(*bit);
642  ts_ctx->vecAssembleSwitch = boost::move(bit->vecAssembleSwitch);
643  }
644 
645  if (*ts_ctx->vecAssembleSwitch) {
646  CHKERR VecGhostUpdateBegin(F, ADD_VALUES, SCATTER_REVERSE);
647  CHKERR VecGhostUpdateEnd(F, ADD_VALUES, SCATTER_REVERSE);
648  CHKERR VecAssemblyBegin(F);
649  CHKERR VecAssemblyEnd(F);
650  }
651 
652  PetscLogEventEnd(ts_ctx->MOFEM_EVENT_TsCtxIFunction, 0, 0, 0, 0);
654 }
static constexpr Switches CtxSetX_T
Definition: LoopMethods.hpp:66
FTensor::Index< 'i', 3 > i
FTensor::Index< 'l', 3 > l
const FTensor::Tensor2< T, Dim, Dim > Vec
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:485
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:604
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:415
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.

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

455  {
457 
458  TsCtx *ts_ctx = static_cast<TsCtx *>(ctx);
459  PetscLogEventBegin(ts_ctx->MOFEM_EVENT_TsCtxI2Function, 0, 0, 0, 0);
460  CHKERR VecGhostUpdateBegin(u, INSERT_VALUES, SCATTER_FORWARD);
461  CHKERR VecGhostUpdateEnd(u, INSERT_VALUES, SCATTER_FORWARD);
462  CHKERR VecGhostUpdateBegin(u_t, INSERT_VALUES, SCATTER_FORWARD);
463  CHKERR VecGhostUpdateEnd(u_t, INSERT_VALUES, SCATTER_FORWARD);
464  CHKERR VecGhostUpdateBegin(u_tt, INSERT_VALUES, SCATTER_FORWARD);
465  CHKERR VecGhostUpdateEnd(u_tt, INSERT_VALUES, SCATTER_FORWARD);
466  CHKERR ts_ctx->mField.getInterface<VecManager>()->setLocalGhostVector(
467  ts_ctx->problemName, COL, u, INSERT_VALUES, SCATTER_REVERSE);
468  if (ts_ctx->zeroMatrix) {
469  CHKERR MatZeroEntries(B);
470  }
471  int step;
472 #if PETSC_VERSION_GE(3, 8, 0)
473  CHKERR TSGetStepNumber(ts, &step);
474 #else
475  CHKERR TSGetTimeStepNumber(ts, &step);
476 #endif
477 
478  ts_ctx->matAssembleSwitch =
479  boost::movelib::make_unique<bool>(ts_ctx->zeroMatrix);
480 
481  auto set = [&](auto &fe) {
482  fe.ts_u = u;
483  fe.ts_u_t = u_t;
484  fe.ts_u_tt = u_tt;
485  fe.ts_A = A;
486  fe.ts_B = B;
487  fe.ts_t = t;
488  fe.ts_a = a;
489  fe.ts_aa = aa;
490  fe.ts_step = step;
491 
492  fe.ts_ctx = TSMethod::CTX_TSSETIJACOBIAN;
493  fe.snes_ctx = SnesMethod::CTX_SNESSETJACOBIAN;
494  fe.ksp_ctx = KspMethod::CTX_OPERATORS;
498  fe.ts = ts;
499  };
500 
501  auto unset = [&](auto &fe) {
502  fe.ts_ctx = TSMethod::CTX_TSNONE;
503  fe.snes_ctx = SnesMethod::CTX_SNESNONE;
504  fe.ksp_ctx = KspMethod::CTX_KSPNONE;
505  fe.data_ctx = PetscData::CtxSetNone;
506  };
507 
508  auto cache_ptr = boost::make_shared<CacheTuple>();
509  CHKERR ts_ctx->mField.cache_problem_entities(ts_ctx->problemName, cache_ptr);
510 
511  // preproces
512  for (auto &bit : ts_ctx->preProcess_IJacobian) {
513  bit->matAssembleSwitch = boost::move(ts_ctx->matAssembleSwitch);
514  set(*bit);
515  CHKERR ts_ctx->mField.problem_basic_method_preProcess(ts_ctx->problemName,
516  *bit);
517  unset(*bit);
518  ts_ctx->matAssembleSwitch = boost::move(bit->matAssembleSwitch);
519  }
520 
521  for (auto &lit : ts_ctx->loops_to_do_IJacobian) {
522  lit.second->matAssembleSwitch = boost::move(ts_ctx->matAssembleSwitch);
523  set(*lit.second);
524  CHKERR ts_ctx->mField.loop_finite_elements(ts_ctx->problemName, lit.first,
525  *(lit.second), nullptr,
526  ts_ctx->bH, cache_ptr);
527  unset(*lit.second);
528  ts_ctx->matAssembleSwitch = boost::move(lit.second->matAssembleSwitch);
529  }
530 
531  // post process
532  for (auto &bit : ts_ctx->postProcess_IJacobian) {
533  bit->matAssembleSwitch = boost::move(ts_ctx->matAssembleSwitch);
534  set(*bit);
535  CHKERR ts_ctx->mField.problem_basic_method_postProcess(ts_ctx->problemName,
536  *bit);
537  unset(*bit);
538  ts_ctx->matAssembleSwitch = boost::move(bit->matAssembleSwitch);
539  }
540 
541  if (ts_ctx->matAssembleSwitch) {
542  CHKERR MatAssemblyBegin(B, MAT_FINAL_ASSEMBLY);
543  CHKERR MatAssemblyEnd(B, MAT_FINAL_ASSEMBLY);
544  }
545  PetscLogEventEnd(ts_ctx->MOFEM_EVENT_TsCtxI2Function, 0, 0, 0, 0);
547 }
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:485
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:604
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:415
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  auto cache_ptr = boost::make_shared<CacheTuple>();
76  CHKERR ts_ctx->mField.cache_problem_entities(ts_ctx->problemName, cache_ptr);
77 
78  // preprocess
79  for (auto &bit : ts_ctx->preProcess_IFunction) {
80  bit->vecAssembleSwitch = boost::move(ts_ctx->vecAssembleSwitch);
81  set(*bit);
82  CHKERR ts_ctx->mField.problem_basic_method_preProcess(ts_ctx->problemName,
83  *bit);
84  unset(*bit);
85  ts_ctx->vecAssembleSwitch = boost::move(bit->vecAssembleSwitch);
86  }
87 
88  // fe loops
89  for (auto &lit : ts_ctx->loops_to_do_IFunction) {
90  lit.second->vecAssembleSwitch = boost::move(ts_ctx->vecAssembleSwitch);
91  set(*lit.second);
92  CHKERR ts_ctx->mField.loop_finite_elements(ts_ctx->problemName, lit.first,
93  *(lit.second), nullptr,
94  ts_ctx->bH, cache_ptr);
95  unset(*lit.second);
96  ts_ctx->vecAssembleSwitch = boost::move(lit.second->vecAssembleSwitch);
97  }
98 
99  // post process
100  for (auto &bit : ts_ctx->postProcess_IFunction) {
101  bit->vecAssembleSwitch = boost::move(ts_ctx->vecAssembleSwitch);
102  set(*bit);
103  CHKERR ts_ctx->mField.problem_basic_method_postProcess(ts_ctx->problemName,
104  *bit);
105  unset(*bit);
106  ts_ctx->vecAssembleSwitch = boost::move(bit->vecAssembleSwitch);
107  }
108 
109  if (*ts_ctx->vecAssembleSwitch) {
110  CHKERR VecGhostUpdateBegin(F, ADD_VALUES, SCATTER_REVERSE);
111  CHKERR VecGhostUpdateEnd(F, ADD_VALUES, SCATTER_REVERSE);
112  CHKERR VecAssemblyBegin(F);
113  CHKERR VecAssemblyEnd(F);
114  }
115 
116  PetscLogEventEnd(ts_ctx->MOFEM_EVENT_TsCtxIFunction, 0, 0, 0, 0);
118 }
static constexpr Switches CtxSetX_T
Definition: LoopMethods.hpp:66
FTensor::Index< 'i', 3 > i
FTensor::Index< 'l', 3 > l
const FTensor::Tensor2< T, Dim, Dim > Vec
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:485
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:604
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:415
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 120 of file TsCtx.cpp.

121  {
123 
124  TsCtx *ts_ctx = static_cast<TsCtx *>(ctx);
125  PetscLogEventBegin(ts_ctx->MOFEM_EVENT_TsCtxIFunction, 0, 0, 0, 0);
126  CHKERR VecGhostUpdateBegin(u, INSERT_VALUES, SCATTER_FORWARD);
127  CHKERR VecGhostUpdateEnd(u, INSERT_VALUES, SCATTER_FORWARD);
128  CHKERR VecGhostUpdateBegin(u_t, INSERT_VALUES, SCATTER_FORWARD);
129  CHKERR VecGhostUpdateEnd(u_t, INSERT_VALUES, SCATTER_FORWARD);
130  CHKERR ts_ctx->mField.getInterface<VecManager>()->setLocalGhostVector(
131  ts_ctx->problemName, COL, u, INSERT_VALUES, SCATTER_REVERSE);
132  if (ts_ctx->zeroMatrix) {
133  CHKERR MatZeroEntries(B);
134  }
135  int step;
136 #if PETSC_VERSION_GE(3, 8, 0)
137  CHKERR TSGetStepNumber(ts, &step);
138 #else
139  CHKERR TSGetTimeStepNumber(ts, &step);
140 #endif
141 
142  ts_ctx->matAssembleSwitch =
143  boost::movelib::make_unique<bool>(ts_ctx->zeroMatrix);
144 
145  auto set = [&](auto &fe) {
146  fe.ts = ts;
147  fe.ts_u = u;
148  fe.ts_u_t = u_t;
149  fe.ts_A = A;
150  fe.ts_B = B;
151  fe.ts_t = t;
152  fe.ts_a = a;
153  fe.ts_step = step;
154  fe.ts_ctx = TSMethod::CTX_TSSETIJACOBIAN;
155  fe.snes_ctx = SnesMethod::CTX_SNESSETJACOBIAN;
156  fe.ksp_ctx = KspMethod::CTX_OPERATORS;
159  };
160 
161  auto unset = [&](auto &fe) {
162  fe.ts_ctx = TSMethod::CTX_TSNONE;
163  fe.snes_ctx = SnesMethod::CTX_SNESNONE;
164  fe.ksp_ctx = KspMethod::CTX_KSPNONE;
165  fe.data_ctx = PetscData::CtxSetNone;
166  };
167 
168  auto cache_ptr = boost::make_shared<CacheTuple>();
169  CHKERR ts_ctx->mField.cache_problem_entities(ts_ctx->problemName, cache_ptr);
170 
171  // preproces
172  for (auto &bit : ts_ctx->preProcess_IJacobian) {
173  bit->matAssembleSwitch = boost::move(ts_ctx->matAssembleSwitch);
174  set(*bit);
175  CHKERR ts_ctx->mField.problem_basic_method_preProcess(ts_ctx->problemName,
176  *bit);
177  unset(*bit);
178  ts_ctx->matAssembleSwitch = boost::move(bit->matAssembleSwitch);
179  }
180 
181  for (auto &lit : ts_ctx->loops_to_do_IJacobian) {
182  lit.second->matAssembleSwitch = boost::move(ts_ctx->matAssembleSwitch);
183  set(*lit.second);
184  CHKERR ts_ctx->mField.loop_finite_elements(ts_ctx->problemName, lit.first,
185  *(lit.second), nullptr,
186  ts_ctx->bH, cache_ptr);
187  unset(*lit.second);
188  ts_ctx->matAssembleSwitch = boost::move(lit.second->matAssembleSwitch);
189  }
190 
191  // post process
192  for (auto &bit : ts_ctx->postProcess_IJacobian) {
193  bit->matAssembleSwitch = boost::move(ts_ctx->matAssembleSwitch);
194  set(*bit);
195  CHKERR ts_ctx->mField.problem_basic_method_postProcess(ts_ctx->problemName,
196  *bit);
197  unset(*bit);
198  ts_ctx->matAssembleSwitch = boost::move(bit->matAssembleSwitch);
199  }
200 
201  if (ts_ctx->matAssembleSwitch) {
202  CHKERR MatAssemblyBegin(B, MAT_FINAL_ASSEMBLY);
203  CHKERR MatAssemblyEnd(B, MAT_FINAL_ASSEMBLY);
204  }
205  PetscLogEventEnd(ts_ctx->MOFEM_EVENT_TsCtxIFunction, 0, 0, 0, 0);
207 }
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:485
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:604
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:415
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 267 of file TsCtx.cpp.

267  {
269  TsCtx *ts_ctx = static_cast<TsCtx *>(ctx);
270  PetscLogEventBegin(ts_ctx->MOFEM_EVENT_TsCtxRHSFunction, 0, 0, 0, 0);
271  CHKERR VecGhostUpdateBegin(u, INSERT_VALUES, SCATTER_FORWARD);
272  CHKERR VecGhostUpdateEnd(u, INSERT_VALUES, SCATTER_FORWARD);
273  CHKERR ts_ctx->mField.getInterface<VecManager>()->setLocalGhostVector(
274  ts_ctx->problemName, COL, u, INSERT_VALUES, SCATTER_REVERSE);
275 
276  auto zero_ghost_vec = [](Vec g) {
278  Vec l;
279  CHKERR VecGhostGetLocalForm(g, &l);
280  double *a;
281  CHKERR VecGetArray(l, &a);
282  int s;
283  CHKERR VecGetLocalSize(l, &s);
284  for (int i = 0; i != s; ++i)
285  a[i] = 0;
286  CHKERR VecRestoreArray(l, &a);
287  CHKERR VecGhostRestoreLocalForm(g, &l);
289  };
290  CHKERR zero_ghost_vec(F);
291 
292  ts_ctx->vecAssembleSwitch = boost::movelib::make_unique<bool>(true);
293 
294  int step;
295 #if PETSC_VERSION_GE(3, 8, 0)
296  CHKERR TSGetStepNumber(ts, &step);
297 #else
298  CHKERR TSGetTimeStepNumber(ts, &step);
299 #endif
300 
301  auto set = [&](auto &fe) {
302  fe.ts_u = u;
303  fe.ts_F = F;
304  fe.ts_t = t;
305  fe.ts = ts;
306  fe.ts_step = step;
307  fe.ts_ctx = TSMethod::CTX_TSSETRHSFUNCTION;
308  fe.snes_ctx = SnesMethod::CTX_SNESSETFUNCTION;
309  fe.ksp_ctx = KspMethod::CTX_SETFUNCTION;
310  fe.data_ctx =
312  };
313 
314  auto unset = [&](auto &fe) {
315  fe.ts_ctx = TSMethod::CTX_TSNONE;
316  fe.snes_ctx = SnesMethod::CTX_SNESNONE;
317  fe.ksp_ctx = KspMethod::CTX_KSPNONE;
318  fe.data_ctx = PetscData::CtxSetNone;
319  };
320 
321  auto cache_ptr = boost::make_shared<CacheTuple>();
322  CHKERR ts_ctx->mField.cache_problem_entities(ts_ctx->problemName, cache_ptr);
323 
324  for (auto &bit : ts_ctx->preProcess_RHSJacobian) {
325  bit->vecAssembleSwitch = boost::move(ts_ctx->vecAssembleSwitch);
326  set(*bit);
327  CHKERR ts_ctx->mField.problem_basic_method_preProcess(ts_ctx->problemName,
328  *bit);
329  unset(*bit);
330  ts_ctx->vecAssembleSwitch = boost::move(bit->vecAssembleSwitch);
331  }
332 
333  // fe loops
334  for (auto &lit : ts_ctx->loops_to_do_RHSFunction) {
335  lit.second->vecAssembleSwitch = boost::move(ts_ctx->vecAssembleSwitch);
336  set(*lit.second);
337  CHKERR ts_ctx->mField.loop_finite_elements(ts_ctx->problemName, lit.first,
338  *(lit.second), nullptr,
339  ts_ctx->bH, cache_ptr);
340  unset(*lit.second);
341  ts_ctx->vecAssembleSwitch = boost::move(lit.second->vecAssembleSwitch);
342  }
343 
344  // post process
345  for (auto &bit : ts_ctx->postProcess_RHSJacobian) {
346  bit->vecAssembleSwitch = boost::move(ts_ctx->vecAssembleSwitch);
347  set(*bit);
348  CHKERR ts_ctx->mField.problem_basic_method_postProcess(ts_ctx->problemName,
349  *bit);
350  unset(*bit);
351  ts_ctx->vecAssembleSwitch = boost::move(bit->vecAssembleSwitch);
352  }
353 
354  if (*ts_ctx->vecAssembleSwitch) {
355  CHKERR VecGhostUpdateBegin(F, ADD_VALUES, SCATTER_REVERSE);
356  CHKERR VecGhostUpdateEnd(F, ADD_VALUES, SCATTER_REVERSE);
357  CHKERR VecAssemblyBegin(F);
358  CHKERR VecAssemblyEnd(F);
359  }
360 
361  PetscLogEventEnd(ts_ctx->MOFEM_EVENT_TsCtxRHSFunction, 0, 0, 0, 0);
363 }
FTensor::Index< 'i', 3 > i
FTensor::Index< 'l', 3 > l
const FTensor::Tensor2< T, Dim, Dim > Vec
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:485
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:604
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:415
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 365 of file TsCtx.cpp.

366  {
368  TsCtx *ts_ctx = static_cast<TsCtx *>(ctx);
369  PetscLogEventBegin(ts_ctx->MOFEM_EVENT_TsCtxRHSJacobian, 0, 0, 0, 0);
370  CHKERR VecGhostUpdateBegin(u, INSERT_VALUES, SCATTER_FORWARD);
371  CHKERR VecGhostUpdateEnd(u, INSERT_VALUES, SCATTER_FORWARD);
372  CHKERR ts_ctx->mField.getInterface<VecManager>()->setLocalGhostVector(
373  ts_ctx->problemName, COL, u, INSERT_VALUES, SCATTER_REVERSE);
374 
375  if (ts_ctx->zeroMatrix) {
376  CHKERR MatZeroEntries(B);
377  }
378 
379  ts_ctx->matAssembleSwitch =
380  boost::movelib::make_unique<bool>(ts_ctx->zeroMatrix);
381 
382  int step;
383 #if PETSC_VERSION_GE(3, 8, 0)
384  CHKERR TSGetStepNumber(ts, &step);
385 #else
386  CHKERR TSGetTimeStepNumber(ts, &step);
387 #endif
388 
389  auto set = [&](auto &fe) {
390  fe.ts_u = u;
391  fe.ts_A = A;
392  fe.ts_B = B;
393  fe.ts_t = t;
394  fe.ts_step = step;
395  fe.ts_ctx = TSMethod::CTX_TSSETRHSJACOBIAN;
396  fe.snes_ctx = SnesMethod::CTX_SNESSETJACOBIAN;
397  fe.ksp_ctx = KspMethod::CTX_OPERATORS;
400  fe.ts = ts;
401  };
402 
403  auto unset = [&](auto &fe) {
404  fe.ts_ctx = TSMethod::CTX_TSNONE;
405  fe.snes_ctx = SnesMethod::CTX_SNESNONE;
406  fe.ksp_ctx = KspMethod::CTX_KSPNONE;
407  fe.data_ctx = PetscData::CtxSetNone;
408  };
409 
410  auto cache_ptr = boost::make_shared<CacheTuple>();
411  CHKERR ts_ctx->mField.cache_problem_entities(ts_ctx->problemName, cache_ptr);
412 
413  // preprocess
414  for (auto &bit : ts_ctx->preProcess_RHSJacobian) {
415  bit->matAssembleSwitch = boost::move(ts_ctx->matAssembleSwitch);
416  set(*bit);
417  CHKERR ts_ctx->mField.problem_basic_method_preProcess(ts_ctx->problemName,
418  *bit);
419  unset(*bit);
420  ts_ctx->matAssembleSwitch = boost::move(bit->matAssembleSwitch);
421  }
422 
423  // fe loops
424  for (auto &lit : ts_ctx->loops_to_do_RHSJacobian) {
425  lit.second->matAssembleSwitch = boost::move(ts_ctx->matAssembleSwitch);
426  set(*lit.second);
427  CHKERR ts_ctx->mField.loop_finite_elements(ts_ctx->problemName, lit.first,
428  *(lit.second), nullptr,
429  ts_ctx->bH, cache_ptr);
430  unset(*lit.second);
431  ts_ctx->matAssembleSwitch = boost::move(lit.second->matAssembleSwitch);
432  }
433 
434  // post process
435  for (auto &bit : ts_ctx->postProcess_RHSJacobian) {
436  bit->matAssembleSwitch = boost::move(ts_ctx->matAssembleSwitch);
437  set(*bit);
438  CHKERR ts_ctx->mField.problem_basic_method_postProcess(ts_ctx->problemName,
439  *bit);
440  unset(*bit);
441  ts_ctx->matAssembleSwitch = boost::move(bit->matAssembleSwitch);
442  }
443 
444  if (ts_ctx->matAssembleSwitch) {
445  CHKERR MatAssemblyBegin(B, MAT_FINAL_ASSEMBLY);
446  CHKERR MatAssemblyEnd(B, MAT_FINAL_ASSEMBLY);
447  }
448 
449  PetscLogEventEnd(ts_ctx->MOFEM_EVENT_TsCtxRHSJacobian, 0, 0, 0, 0);
451 }
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:485
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:604
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:415
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: