v0.8.23
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
 
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...
 

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  }
PetscLogEvent MOFEM_EVENT_TsCtxIFunction
Definition: TsCtx.hpp:242
bool zeroMatrix
Definition: TsCtx.hpp:64
virtual moab::Interface & get_moab()=0
PetscLogEvent MOFEM_EVENT_TsCtxIJacobian
Definition: TsCtx.hpp:243
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:240
PetscLogEvent MOFEM_EVENT_TsCtxMonitor
Definition: TsCtx.hpp:244
MoFEM::Interface & mField
Definition: TsCtx.hpp:29
PetscLogEvent MOFEM_EVENT_TsCtxRHSJacobian
Definition: TsCtx.hpp:241

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

87  {
88  return loops_to_do_IFunction;
89  }
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 111 of file TsCtx.hpp.

111  {
112  return loops_to_do_IJacobian;
113  }
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 99 of file TsCtx.hpp.

99  {
101  }
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 190 of file TsCtx.cpp.

191  {
193  TsCtx *ts_ctx = static_cast<TsCtx *>(ctx);
194  PetscLogEventBegin(ts_ctx->MOFEM_EVENT_TsCtxMonitor, 0, 0, 0, 0);
195  CHKERR VecGhostUpdateBegin(u, INSERT_VALUES, SCATTER_FORWARD);
196  CHKERR VecGhostUpdateEnd(u, INSERT_VALUES, SCATTER_FORWARD);
197  CHKERR ts_ctx->mField.getInterface<VecManager>()->setLocalGhostVector(
198  ts_ctx->problemName, COL, u, INSERT_VALUES, SCATTER_REVERSE);
199 
200  // preproces
201  for (auto &bit : ts_ctx->preProcess_Monitor) {
202  bit->ts_u = u;
203  bit->ts_t = t;
204  bit->ts_step = step;
205  bit->ts_F = PETSC_NULL;
206  CHKERR bit->setTsCtx(TSMethod::CTX_TSTSMONITORSET);
207  CHKERR bit->setTs(ts);
208  CHKERR ts_ctx->mField.problem_basic_method_preProcess(ts_ctx->problemName,
209  *bit);
210  CHKERR bit->setTsCtx(TSMethod::CTX_TSNONE);
211  }
212 
213  for (auto &lit : ts_ctx->loops_to_do_Monitor) {
214  lit.second->ts_u = u;
215  lit.second->ts_t = t;
216  lit.second->ts_step = step;
217  lit.second->ts_F = PETSC_NULL;
218  CHKERR lit.second->setTsCtx(TSMethod::CTX_TSTSMONITORSET);
219  CHKERR lit.second->setTs(ts);
220  CHKERR ts_ctx->mField.loop_finite_elements(
221  ts_ctx->problemName, lit.first, *(lit.second), nullptr, ts_ctx->bH);
222  CHKERR lit.second->setTsCtx(TSMethod::CTX_TSNONE);
223  }
224 
225  // post process
226  for (auto &bit : ts_ctx->postProcess_Monitor) {
227  bit->ts_u = u;
228  bit->ts_t = t;
229  bit->ts_step = step;
230  bit->ts_F = PETSC_NULL;
231  CHKERR bit->setTsCtx(TSMethod::CTX_TSTSMONITORSET);
232  CHKERR bit->setTs(ts);
233  CHKERR ts_ctx->mField.problem_basic_method_postProcess(ts_ctx->problemName,
234  *bit);
235  CHKERR bit->setTsCtx(TSMethod::CTX_TSNONE);
236  }
237  PetscLogEventEnd(ts_ctx->MOFEM_EVENT_TsCtxMonitor, 0, 0, 0, 0);
239 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:476
TsCtx(MoFEM::Interface &m_field, const std::string &problem_name)
Definition: TsCtx.hpp:66
#define CHKERR
Inline error check.
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:406

◆ 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  CHKERR TSGetTimeStepNumber(ts, &step);
49  // preprocess
50  for (auto &bit : ts_ctx->preProcess_IFunction) {
51  bit->vecAssembleSwitch = boost::move(ts_ctx->vecAssembleSwitch);
52  bit->ts_u = u;
53  bit->ts_u_t = u_t;
54  bit->ts_F = F;
55  bit->ts_t = t;
56  bit->ts_step = step;
57  CHKERR bit->setTsCtx(TSMethod::CTX_TSSETIFUNCTION);
58  CHKERR bit->setTs(ts);
59  CHKERR ts_ctx->mField.problem_basic_method_preProcess(ts_ctx->problemName,
60  *bit);
61  CHKERR bit->setTsCtx(TSMethod::CTX_TSNONE);
62  ts_ctx->vecAssembleSwitch = boost::move(bit->vecAssembleSwitch);
63  }
64 
65  // fe loops
66  for (auto &lit : ts_ctx->loops_to_do_IFunction) {
67  lit.second->vecAssembleSwitch = boost::move(ts_ctx->vecAssembleSwitch);
68  lit.second->ts_u = u;
69  lit.second->ts_u_t = u_t;
70  lit.second->ts_F = F;
71  lit.second->ts_t = t;
72  lit.second->ts_step = step;
73  CHKERR lit.second->setTsCtx(TSMethod::CTX_TSSETIFUNCTION);
74  CHKERR lit.second->setTs(ts);
75  CHKERR ts_ctx->mField.loop_finite_elements(
76  ts_ctx->problemName, lit.first, *(lit.second), nullptr, ts_ctx->bH);
77  CHKERR lit.second->setTsCtx(TSMethod::CTX_TSNONE);
78  ts_ctx->vecAssembleSwitch = boost::move(lit.second->vecAssembleSwitch);
79  }
80 
81  // post process
82  for (auto &bit : ts_ctx->postProcess_IFunction) {
83  bit->vecAssembleSwitch = boost::move(ts_ctx->vecAssembleSwitch);
84  bit->ts_u = u;
85  bit->ts_u_t = u_t;
86  bit->ts_F = F;
87  bit->ts_t = t;
88  bit->ts_step = step;
89  CHKERR bit->setTsCtx(TSMethod::CTX_TSSETIFUNCTION);
90  CHKERR bit->setTs(ts);
91  CHKERR ts_ctx->mField.problem_basic_method_postProcess(ts_ctx->problemName,
92  *bit);
93  CHKERR bit->setTsCtx(TSMethod::CTX_TSNONE);
94  ts_ctx->vecAssembleSwitch = boost::move(bit->vecAssembleSwitch);
95  }
96 
97  if (*ts_ctx->vecAssembleSwitch) {
98  CHKERR VecGhostUpdateBegin(F, ADD_VALUES, SCATTER_REVERSE);
99  CHKERR VecGhostUpdateEnd(F, ADD_VALUES, SCATTER_REVERSE);
100  CHKERR VecAssemblyBegin(F);
101  CHKERR VecAssemblyEnd(F);
102  }
103 
104  PetscLogEventEnd(ts_ctx->MOFEM_EVENT_TsCtxIFunction, 0, 0, 0, 0);
106 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:476
TsCtx(MoFEM::Interface &m_field, const std::string &problem_name)
Definition: TsCtx.hpp:66
#define CHKERR
Inline error check.
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:406

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

109  {
111 
112  TsCtx *ts_ctx = static_cast<TsCtx *>(ctx);
113  PetscLogEventBegin(ts_ctx->MOFEM_EVENT_TsCtxIFunction, 0, 0, 0, 0);
114  CHKERR VecGhostUpdateBegin(u, INSERT_VALUES, SCATTER_FORWARD);
115  CHKERR VecGhostUpdateEnd(u, INSERT_VALUES, SCATTER_FORWARD);
116  CHKERR VecGhostUpdateBegin(u_t, INSERT_VALUES, SCATTER_FORWARD);
117  CHKERR VecGhostUpdateEnd(u_t, INSERT_VALUES, SCATTER_FORWARD);
118  CHKERR ts_ctx->mField.getInterface<VecManager>()->setLocalGhostVector(
119  ts_ctx->problemName, COL, u, INSERT_VALUES, SCATTER_REVERSE);
120  if (ts_ctx->zeroMatrix) {
121  CHKERR MatZeroEntries(B);
122  }
123  int step;
124  CHKERR TSGetTimeStepNumber(ts, &step);
125 
126  ts_ctx->matAssembleSwitch =
127  boost::movelib::make_unique<bool>(ts_ctx->zeroMatrix);
128 
129  // preproces
130  for (auto &bit : ts_ctx->preProcess_IJacobian) {
131  bit->matAssembleSwitch = boost::move(ts_ctx->matAssembleSwitch);
132  bit->ts_u = u;
133  bit->ts_u_t = u_t;
134  bit->ts_A = A;
135  bit->ts_B = B;
136  bit->ts_t = t;
137  bit->ts_a = a;
138  bit->ts_step = step;
139  CHKERR bit->setTsCtx(TSMethod::CTX_TSSETIJACOBIAN);
140  CHKERR bit->setTs(ts);
141  CHKERR ts_ctx->mField.problem_basic_method_preProcess(ts_ctx->problemName,
142  *bit);
143  CHKERR bit->setTsCtx(TSMethod::CTX_TSNONE);
144  ts_ctx->matAssembleSwitch = boost::move(bit->matAssembleSwitch);
145  }
146 
147  for (auto &lit : ts_ctx->loops_to_do_IJacobian) {
148  lit.second->matAssembleSwitch = boost::move(ts_ctx->matAssembleSwitch);
149  lit.second->ts_u = u;
150  lit.second->ts_u_t = u_t;
151  lit.second->ts_A = A;
152  lit.second->ts_B = B;
153  lit.second->ts_t = t;
154  lit.second->ts_a = a;
155  lit.second->ts_step = step;
156  CHKERR lit.second->setTsCtx(TSMethod::CTX_TSSETIJACOBIAN);
157  CHKERR lit.second->setTs(ts);
158  CHKERR ts_ctx->mField.loop_finite_elements(
159  ts_ctx->problemName, lit.first, *(lit.second), nullptr, ts_ctx->bH);
160  CHKERR lit.second->setTsCtx(TSMethod::CTX_TSNONE);
161  ts_ctx->matAssembleSwitch = boost::move(lit.second->matAssembleSwitch);
162  }
163 
164  // post process
165  for (auto &bit : ts_ctx->postProcess_IJacobian) {
166  bit->matAssembleSwitch = boost::move(ts_ctx->matAssembleSwitch);
167  bit->ts_u = u;
168  bit->ts_u_t = u_t;
169  bit->ts_A = A;
170  bit->ts_B = B;
171  bit->ts_t = t;
172  bit->ts_a = a;
173  bit->ts_step = step;
174  CHKERR bit->setTsCtx(TSMethod::CTX_TSSETIJACOBIAN);
175  CHKERR bit->setTs(ts);
176  CHKERR ts_ctx->mField.problem_basic_method_postProcess(ts_ctx->problemName,
177  *bit);
178  CHKERR bit->setTsCtx(TSMethod::CTX_TSNONE);
179  ts_ctx->matAssembleSwitch = boost::move(bit->matAssembleSwitch);
180  }
181 
182  if (ts_ctx->matAssembleSwitch) {
183  CHKERR MatAssemblyBegin(B, MAT_FINAL_ASSEMBLY);
184  CHKERR MatAssemblyEnd(B, MAT_FINAL_ASSEMBLY);
185  }
186  PetscLogEventEnd(ts_ctx->MOFEM_EVENT_TsCtxIFunction, 0, 0, 0, 0);
188 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:476
TsCtx(MoFEM::Interface &m_field, const std::string &problem_name)
Definition: TsCtx.hpp:66
#define CHKERR
Inline error check.
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:406

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

241  {
243  TsCtx *ts_ctx = static_cast<TsCtx *>(ctx);
244  PetscLogEventBegin(ts_ctx->MOFEM_EVENT_TsCtxRHSFunction, 0, 0, 0, 0);
245  CHKERR VecGhostUpdateBegin(u, INSERT_VALUES, SCATTER_FORWARD);
246  CHKERR VecGhostUpdateEnd(u, INSERT_VALUES, SCATTER_FORWARD);
247  CHKERR ts_ctx->mField.getInterface<VecManager>()->setLocalGhostVector(
248  ts_ctx->problemName, COL, u, INSERT_VALUES, SCATTER_REVERSE);
249 
250  auto zero_ghost_vec = [](Vec g) {
252  Vec l;
253  CHKERR VecGhostGetLocalForm(g, &l);
254  double *a;
255  CHKERR VecGetArray(l, &a);
256  int s;
257  CHKERR VecGetLocalSize(l, &s);
258  for (int i = 0; i != s; ++i)
259  a[i] = 0;
260  CHKERR VecRestoreArray(l, &a);
261  CHKERR VecGhostRestoreLocalForm(g, &l);
263  };
264  CHKERR zero_ghost_vec(F);
265 
266  ts_ctx->vecAssembleSwitch = boost::movelib::make_unique<bool>(true);
267 
268  int step;
269  CHKERR TSGetTimeStepNumber(ts, &step);
270  for (auto &bit : ts_ctx->preProcess_RHSJacobian) {
271  bit->vecAssembleSwitch = boost::move(ts_ctx->vecAssembleSwitch);
272  bit->ts_u = u;
273  bit->ts_F = F;
274  bit->ts_t = t;
275  bit->ts_step = step;
276  CHKERR bit->setTsCtx(TSMethod::CTX_TSSETRHSFUNCTION);
277  CHKERR bit->setTs(ts);
278  CHKERR ts_ctx->mField.problem_basic_method_preProcess(ts_ctx->problemName,
279  *bit);
280  CHKERR bit->setTsCtx(TSMethod::CTX_TSNONE);
281  ts_ctx->vecAssembleSwitch = boost::move(bit->vecAssembleSwitch);
282  }
283 
284  // fe loops
285  for (auto &lit : ts_ctx->loops_to_do_RHSFunction) {
286  lit.second->vecAssembleSwitch = boost::move(ts_ctx->vecAssembleSwitch);
287  lit.second->ts_u = u;
288  lit.second->ts_F = F;
289  lit.second->ts_t = t;
290  lit.second->ts_step = step;
291  CHKERR lit.second->setTsCtx(TSMethod::CTX_TSSETRHSFUNCTION);
292  CHKERR lit.second->setTs(ts);
293  CHKERR ts_ctx->mField.loop_finite_elements(
294  ts_ctx->problemName, lit.first, *(lit.second), nullptr, ts_ctx->bH);
295  CHKERR lit.second->setTsCtx(TSMethod::CTX_TSNONE);
296  ts_ctx->vecAssembleSwitch = boost::move(lit.second->vecAssembleSwitch);
297  }
298 
299  // post process
300  for (auto &bit : ts_ctx->postProcess_RHSJacobian) {
301  bit->vecAssembleSwitch = boost::move(ts_ctx->vecAssembleSwitch);
302  bit->ts_u = u;
303  bit->ts_F = F;
304  bit->ts_t = t;
305  bit->ts_step = step;
306  CHKERR bit->setTsCtx(TSMethod::CTX_TSSETRHSFUNCTION);
307  CHKERR bit->setTs(ts);
308  CHKERR ts_ctx->mField.problem_basic_method_postProcess(ts_ctx->problemName,
309  *bit);
310  CHKERR bit->setTsCtx(TSMethod::CTX_TSNONE);
311  ts_ctx->vecAssembleSwitch = boost::move(bit->vecAssembleSwitch);
312  }
313 
314  if (*ts_ctx->vecAssembleSwitch) {
315  CHKERR VecGhostUpdateBegin(F, ADD_VALUES, SCATTER_REVERSE);
316  CHKERR VecGhostUpdateEnd(F, ADD_VALUES, SCATTER_REVERSE);
317  CHKERR VecAssemblyBegin(F);
318  CHKERR VecAssemblyEnd(F);
319  }
320 
321  PetscLogEventEnd(ts_ctx->MOFEM_EVENT_TsCtxRHSFunction, 0, 0, 0, 0);
323 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:476
TsCtx(MoFEM::Interface &m_field, const std::string &problem_name)
Definition: TsCtx.hpp:66
#define CHKERR
Inline error check.
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:406

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

326  {
328  TsCtx *ts_ctx = static_cast<TsCtx *>(ctx);
329  PetscLogEventBegin(ts_ctx->MOFEM_EVENT_TsCtxRHSJacobian, 0, 0, 0, 0);
330  CHKERR VecGhostUpdateBegin(u, INSERT_VALUES, SCATTER_FORWARD);
331  CHKERR VecGhostUpdateEnd(u, INSERT_VALUES, SCATTER_FORWARD);
332  CHKERR ts_ctx->mField.getInterface<VecManager>()->setLocalGhostVector(
333  ts_ctx->problemName, COL, u, INSERT_VALUES, SCATTER_REVERSE);
334 
335  if (ts_ctx->zeroMatrix) {
336  CHKERR MatZeroEntries(B);
337  }
338 
339  ts_ctx->matAssembleSwitch =
340  boost::movelib::make_unique<bool>(ts_ctx->zeroMatrix);
341 
342  int step;
343  CHKERR TSGetTimeStepNumber(ts, &step);
344  // preprocess
345  for (auto &bit : ts_ctx->preProcess_RHSJacobian) {
346  bit->matAssembleSwitch = boost::move(ts_ctx->matAssembleSwitch);
347  bit->ts_u = u;
348  bit->ts_A = A;
349  bit->ts_B = B;
350  bit->ts_t = t;
351  bit->ts_step = step;
352  CHKERR bit->setTsCtx(TSMethod::CTX_TSSETRHSJACOBIAN);
353  CHKERR bit->setTs(ts);
354  CHKERR ts_ctx->mField.problem_basic_method_preProcess(ts_ctx->problemName,
355  *bit);
356  CHKERR bit->setTsCtx(TSMethod::CTX_TSNONE);
357  ts_ctx->matAssembleSwitch = boost::move(bit->matAssembleSwitch);
358  }
359 
360  // fe loops
361  for (auto &lit : ts_ctx->loops_to_do_RHSJacobian) {
362  lit.second->matAssembleSwitch = boost::move(ts_ctx->matAssembleSwitch);
363  lit.second->ts_u = u;
364  lit.second->ts_A = A;
365  lit.second->ts_B = B;
366  lit.second->ts_t = t;
367  lit.second->ts_step = step;
368  CHKERR lit.second->setTsCtx(TSMethod::CTX_TSSETRHSJACOBIAN);
369  CHKERR lit.second->setTs(ts);
370  CHKERR ts_ctx->mField.loop_finite_elements(
371  ts_ctx->problemName, lit.first, *(lit.second), nullptr, ts_ctx->bH);
372  CHKERR lit.second->setTsCtx(TSMethod::CTX_TSNONE);
373  ts_ctx->matAssembleSwitch = boost::move(lit.second->matAssembleSwitch);
374  }
375 
376  // post process
377  for (auto &bit : ts_ctx->postProcess_RHSJacobian) {
378  bit->matAssembleSwitch = boost::move(ts_ctx->matAssembleSwitch);
379  bit->ts_u = u;
380  bit->ts_A = A;
381  bit->ts_B = B;
382  bit->ts_t = t;
383  bit->ts_step = step;
384  CHKERR bit->setTsCtx(TSMethod::CTX_TSSETRHSJACOBIAN);
385  CHKERR bit->setTs(ts);
386  CHKERR ts_ctx->mField.problem_basic_method_postProcess(ts_ctx->problemName,
387  *bit);
388  CHKERR bit->setTsCtx(TSMethod::CTX_TSNONE);
389  ts_ctx->matAssembleSwitch = boost::move(bit->matAssembleSwitch);
390  }
391 
392  if (ts_ctx->matAssembleSwitch) {
393  CHKERR MatAssemblyBegin(B, MAT_FINAL_ASSEMBLY);
394  CHKERR MatAssemblyEnd(B, MAT_FINAL_ASSEMBLY);
395  }
396 
397  PetscLogEventEnd(ts_ctx->MOFEM_EVENT_TsCtxRHSJacobian, 0, 0, 0, 0);
399 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:476
TsCtx(MoFEM::Interface &m_field, const std::string &problem_name)
Definition: TsCtx.hpp:66
#define CHKERR
Inline error check.
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:406

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 247 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_TsCtxIFunction

PetscLogEvent MoFEM::TsCtx::MOFEM_EVENT_TsCtxIFunction
private

Definition at line 242 of file TsCtx.hpp.

◆ MOFEM_EVENT_TsCtxIJacobian

PetscLogEvent MoFEM::TsCtx::MOFEM_EVENT_TsCtxIJacobian
private

Definition at line 243 of file TsCtx.hpp.

◆ MOFEM_EVENT_TsCtxMonitor

PetscLogEvent MoFEM::TsCtx::MOFEM_EVENT_TsCtxMonitor
private

Definition at line 244 of file TsCtx.hpp.

◆ MOFEM_EVENT_TsCtxRHSFunction

PetscLogEvent MoFEM::TsCtx::MOFEM_EVENT_TsCtxRHSFunction
private

Definition at line 240 of file TsCtx.hpp.

◆ MOFEM_EVENT_TsCtxRHSJacobian

PetscLogEvent MoFEM::TsCtx::MOFEM_EVENT_TsCtxRHSJacobian
private

Definition at line 241 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 246 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: