v0.14.0
TsCtx.hpp
Go to the documentation of this file.
1 /** \file TsCtx.hpp
2  * \brief Context for PETSc Time Stepping
3  */
4 
5 #ifndef __TSCTX_HPP__
6 #define __TSCTX_HPP__
7 
8 #include <petsc/private/tsimpl.h>
9 
10 #define TSADAPTMOFEM "mofem"
11 
12 namespace MoFEM {
13 
14 /** \brief Interface for Time Stepping (TS) solver
15  * \ingroup mofem_petsc_solvers
16  */
17 struct TsCtx {
18 
21 
22  std::string problemName;
23  MoFEMTypes bH; ///< If set to MF_EXIST check if element exist
24 
28 
44 
45  boost::function<MoFEMErrorCode(TS ts, PetscReal t, Vec u, Vec u_t, Vec u_tt,
46  Vec F, void *ctx)>
48 
49  bool zeroMatrix;
50 
51  TsCtx(MoFEM::Interface &m_field, const std::string &problem_name);
52 
53  virtual ~TsCtx() = default;
54 
55  /**
56  * @brief Get the loops to do IFunction object
57  *
58  * It is sequence of finite elements used to evaluate the right hand side of
59  * implicit time solver.
60  *
61  * @return FEMethodsSequence&
62  */
64 
65  /**
66  * @brief Get the loops to do RHSFunction object
67  *
68  * It is sequence of finite elements used to evaluate the right hand side of
69  * implicit time solver.
70  *
71  * @return FEMethodsSequence&
72  */
74 
75  /**
76  * @brief Get the loops to do IJacobian object
77  *
78  * It is sequence of finite elements used to evalite the left hand sie of
79  * implicit time solver.
80  *
81  * @return FEMethodsSequence&
82  */
84 
85  /**
86  * @brief Get the loops to do RHSJacobian object
87  *
88  * It is sequence of finite elements used to evalite the left hand sie of
89  * implicit time solver.
90  *
91  * @return FEMethodsSequence&
92  */
94 
95  /**
96  * @brief Get the loops to do Monitor object
97  *
98  * It is sequence used to monitor solution of time solver.
99  *
100  * @return FEMethodsSequence&
101  */
103 
104  /**
105  * @brief Get the preProcess to do IFunction object
106  *
107  * @return BasicMethodsSequence&
108  */
110 
111  /**
112  * @brief Get the postProcess to do IFunction object
113  *
114  * @return BasicMethodsSequence&
115  */
117  return postProcessIFunction;
118  }
119 
120  /**
121  * @brief Get the preProcess to do IJacobian object
122  *
123  * @return BasicMethodsSequence&
124  */
126 
127  /**
128  * @brief Get the postProcess to do IJacobian object
129  *
130  * @return BasicMethodsSequence&
131  */
133  return postProcessIJacobian;
134  }
135 
136  /**
137  * @brief Get the preProcess to do Monitor object
138  *
139  * @return BasicMethodsSequence&
140  */
142 
143  /**
144  * @brief Get the postProcess to do Monitor object
145  *
146  * @return BasicMethodsSequence&
147  */
149 
150  /**
151  * @brief Get the preProcess to do RHSJacobian object
152  *
153  * @return BasicMethodsSequence&
154  */
156  return preProcessRHSJacobian;
157  }
158 
159  /**
160  * @brief Get the postProcess to do RHSJacobian object
161  *
162  * @return BasicMethodsSequence&
163  */
165  return postProcessRHSJacobian;
166  }
167 
168  /**
169  * @brief Get the preProcess to do RHSFunction object
170  *
171  * @return BasicMethodsSequence&
172  */
174  return preProcessRHSFunction;
175  }
176 
177  /**
178  * @brief Get the postProcess to do RHSFunction object
179  *
180  * @return BasicMethodsSequence&
181  */
183  return postProcessRHSFunction;
184  }
185 
186  /**
187  * @brief Clear loops
188  *
189  * @return MoFEMErrorCode
190  */
192 
193  friend PetscErrorCode TsSetIFunction(TS ts, PetscReal t, Vec u, Vec u_t,
194  Vec F, void *ctx);
195  friend PetscErrorCode TsSetIJacobian(TS ts, PetscReal t, Vec u, Vec U_t,
196  PetscReal a, Mat A, Mat B, void *ctx);
197  friend PetscErrorCode TsMonitorSet(TS ts, PetscInt step, PetscReal t, Vec u,
198  void *ctx);
199  friend PetscErrorCode TsSetRHSFunction(TS ts, PetscReal t, Vec u, Vec F,
200  void *ctx);
201  friend PetscErrorCode TsSetRHSJacobian(TS ts, PetscReal t, Vec u, Mat A,
202  Mat B, void *ctx);
203 
204  friend PetscErrorCode TsSetI2Function(TS ts, PetscReal t, Vec U, Vec U_t,
205  Vec U_tt, Vec F, void *ctx);
206 
207  friend PetscErrorCode TsSetI2Jacobian(TS ts, PetscReal t, Vec U, Vec U_t,
208  Vec U_tt, PetscReal v, PetscReal a,
209  Mat J, Mat P, void *ctx);
210 
211 private:
219 
220  boost::movelib::unique_ptr<bool> vecAssembleSwitch;
221  boost::movelib::unique_ptr<bool> matAssembleSwitch;
222 };
223 
224 /**
225  * @brief Set IFunction for TS solver
226  *
227  * <a
228  * href=https://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/TS/TSSetIFunction.html>See
229  * petsc for details</a>
230  *
231  * @param ts
232  * @param t
233  * @param u
234  * @param u_t
235  * @param F
236  * @param ctx
237  * @return PetscErrorCode
238  */
239 PetscErrorCode TsSetIFunction(TS ts, PetscReal t, Vec u, Vec u_t, Vec F,
240  void *ctx);
241 
242 /**
243  * @brief Set function evaluating jacobian in TS solver
244  *
245  * <a
246  * href=https://www.mcs.anl.gov/petsc/petsc-3.1/docs/manualpages/TS/TSSetIJacobian.html>See
247  * PETSc for details</a>
248  *
249  * @param ts
250  * @param t
251  * @param u
252  * @param u_t
253  * @param a
254  * @param A
255  * @param B
256  * @param ctx
257  * @return PetscErrorCode
258  */
259 PetscErrorCode TsSetIJacobian(TS ts, PetscReal t, Vec u, Vec u_t, PetscReal a,
260  Mat A, Mat B, void *ctx);
261 
262 /**
263  * @brief Set monitor for TS solver
264  *
265  * <a
266  * href=https://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/TS/TSMonitorSet.html>See
267  * PETSc for details</a>
268  *
269  * @param ts
270  * @param step
271  * @param t
272  * @param u
273  * @param ctx
274  * @return PetscErrorCode
275  */
276 PetscErrorCode TsMonitorSet(TS ts, PetscInt step, PetscReal t, Vec u,
277  void *ctx);
278 
279 /**
280  * @brief TS solver function
281  *
282  * <a
283  * href=https://www.mcs.anl.gov/petsc/petsc-3.11/docs/manualpages/TS/TSSetRHSFunction.html#TSSetRHSFunction>See
284  * PETSc for details</a>
285  *
286  * @param ts
287  * @param t
288  * @param u
289  * @param F
290  * @param ctx
291  * @return PetscErrorCode
292  */
293 PetscErrorCode TsSetRHSFunction(TS ts, PetscReal t, Vec u, Vec F, void *ctx);
294 
295 /**
296  * @brief TS solver function
297  *
298  * <a
299  * href=https://www.mcs.anl.gov/petsc/petsc-3.11/docs/manualpages/TS/TSSetRHSJacobian.html#TSSetRHSJacobian>See
300  * PETSc for details</a>
301  *
302  * @param ts
303  * @param t
304  * @param u
305  * @param A
306  * @param B
307  * @param ctx
308  * @return PetscErrorCode
309  */
310 PetscErrorCode TsSetRHSJacobian(TS ts, PetscReal t, Vec u, Mat A, Mat B,
311  void *ctx);
312 
313 /**
314  * @brief Calculation Jacobian for second order PDE in time
315  *
316  * <a
317  * href=https://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/TS/TSSetI2Jacobian.html>See
318  * PETSc for details</a>
319  *
320  * @param ts
321  * @param t time at step/stage being solved
322  * @param u state vectora
323  * @param u_t time derivative of state vector
324  * @param u_tt second time derivative of state vector
325  * @param a shift for u_t
326  * @param aa shift for u_tt
327  * @param A Jacobian of G(U) = F(t,U,W+v*U,W'+a*U), equivalent to dF/dU +
328  * v*dF/dU_t + a*dF/dU_tt
329  * @param B preconditioning matrix for J, may be same as J
330  * @param ctx TsCtx context for matrix evaluation routine
331  * @return PetscErrorCode
332  */
333 PetscErrorCode TsSetI2Jacobian(TS ts, PetscReal t, Vec u, Vec u_t, Vec u_tt,
334  PetscReal a, PetscReal aa, Mat A, Mat B,
335  void *ctx);
336 
337 /**
338  * @brief Calculation the right hand side for second order PDE in time
339  *
340  * <a
341  * href=https://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/TS/TSSetI2Function.html>PETSc
342  * for details</a>
343  *
344  * @param ts
345  * @param t
346  * @param u
347  * @param u_t
348  * @param u_tt
349  * @param F
350  * @param ctx
351  * @return PetscErrorCode
352  */
353 PetscErrorCode TsSetI2Function(TS ts, PetscReal t, Vec u, Vec u_t, Vec u_tt,
354  Vec F, void *ctx);
355 
356 static PetscErrorCode TSAdaptChooseMoFEM(TSAdapt adapt, TS ts, PetscReal h,
357  PetscInt *next_sc, PetscReal *next_h,
358  PetscBool *accept, PetscReal *wlte,
359  PetscReal *wltea, PetscReal *wlter);
360 static PetscErrorCode TSAdaptResetMoFEM(TSAdapt adapt);
361 
362 static PetscErrorCode TSAdaptDestroyMoFEM(TSAdapt adapt);
363 
364 /**
365  * @brief Craete MOFEM adapt
366  *
367  * \code
368  * CHKERR TSAdaptRegister(TSADAPTMOFEM, TSAdaptCreateMoFEM);
369  * TSAdapt adapt;
370  * CHKERR TSGetAdapt(solver, &adapt);
371  * CHKERR TSAdaptSetType(adapt, TSADAPTMOFEM);
372  * \endcode
373  *
374  * @param adapt
375  * @return PetscErrorCode
376  */
377 PetscErrorCode TSAdaptCreateMoFEM(TSAdapt adapt);
378 
379 } // namespace MoFEM
380 
381 #endif // __TSCTX_HPP__
MoFEM::TsCtx::preProcessRHSJacobian
BasicMethodsSequence preProcessRHSJacobian
Definition: TsCtx.hpp:40
MoFEM::TsCtx::TsSetRHSFunction
friend PetscErrorCode TsSetRHSFunction(TS ts, PetscReal t, Vec u, Vec F, void *ctx)
TS solver function.
Definition: TsCtx.cpp:327
MoFEM::TsCtx::getLoopsRHSFunction
FEMethodsSequence & getLoopsRHSFunction()
Get the loops to do RHSFunction object.
Definition: TsCtx.hpp:73
MoFEM::TsCtx::clearLoops
MoFEMErrorCode clearLoops()
Clear loops.
Definition: TsCtx.cpp:36
MoFEM::TsCtx::getLoopsRHSJacobian
FEMethodsSequence & getLoopsRHSJacobian()
Get the loops to do RHSJacobian object.
Definition: TsCtx.hpp:93
MoFEM::TsSetIFunction
PetscErrorCode TsSetIFunction(TS ts, PetscReal t, Vec u, Vec u_t, Vec F, void *ctx)
Set IFunction for TS solver.
Definition: TsCtx.cpp:56
MoFEM::TsCtx::MOFEM_EVENT_TsCtxIJacobian
PetscLogEvent MOFEM_EVENT_TsCtxIJacobian
Definition: TsCtx.hpp:215
MoFEM::TSAdaptDestroyMoFEM
PetscErrorCode TSAdaptDestroyMoFEM(TSAdapt adapt)
Definition: TsCtx.cpp:823
MoFEM::TsCtx::getLoopsMonitor
FEMethodsSequence & getLoopsMonitor()
Get the loops to do Monitor object.
Definition: TsCtx.hpp:102
MoFEM::TsCtx::moab
moab::Interface & moab
Definition: TsCtx.hpp:20
MoFEM::TsCtx::getPostProcessMonitor
BasicMethodsSequence & getPostProcessMonitor()
Get the postProcess to do Monitor object.
Definition: TsCtx.hpp:148
MoFEM::TsCtx::MOFEM_EVENT_TsCtxI2Function
PetscLogEvent MOFEM_EVENT_TsCtxI2Function
Definition: TsCtx.hpp:217
MoFEM::TsCtx::getLoopsIFunction
FEMethodsSequence & getLoopsIFunction()
Get the loops to do IFunction object.
Definition: TsCtx.hpp:63
MoFEM::TsCtx::MOFEM_EVENT_TsCtxIFunction
PetscLogEvent MOFEM_EVENT_TsCtxIFunction
Definition: TsCtx.hpp:214
MoFEM::Exceptions::MoFEMErrorCode
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:56
MoFEM::FEMethodsSequence
std::deque< PairNameFEMethodPtr > FEMethodsSequence
Definition: AuxPETSc.hpp:56
MoFEM::TsCtx::MOFEM_EVENT_TsCtxMonitor
PetscLogEvent MOFEM_EVENT_TsCtxMonitor
Definition: TsCtx.hpp:216
MoFEM::TsCtx::getPostProcessRHSFunction
BasicMethodsSequence & getPostProcessRHSFunction()
Get the postProcess to do RHSFunction object.
Definition: TsCtx.hpp:182
A
constexpr AssemblyType A
Definition: operators_tests.cpp:30
MoFEM::TsCtx::getPostProcessIJacobian
BasicMethodsSequence & getPostProcessIJacobian()
Get the postProcess to do IJacobian object.
Definition: TsCtx.hpp:132
MoFEM::TsSetI2Jacobian
PetscErrorCode TsSetI2Jacobian(TS ts, PetscReal t, Vec u, Vec u_t, Vec u_tt, PetscReal a, PetscReal aa, Mat A, Mat B, void *ctx)
Calculation Jacobian for second order PDE in time.
Definition: TsCtx.cpp:519
MoFEM::TsCtx::vecAssembleSwitch
boost::movelib::unique_ptr< bool > vecAssembleSwitch
Definition: TsCtx.hpp:220
J
FTensor::Index< 'J', DIM1 > J
Definition: level_set.cpp:30
MoFEM::TsCtx::FEMethodsSequence
MoFEM::FEMethodsSequence FEMethodsSequence
Definition: TsCtx.hpp:26
MoFEM::TsCtx::loopsIJacobian
FEMethodsSequence loopsIJacobian
Definition: TsCtx.hpp:29
MoFEM::DeprecatedCoreInterface
Deprecated interface functions.
Definition: DeprecatedCoreInterface.hpp:16
MoFEM::TsCtx::bH
MoFEMTypes bH
If set to MF_EXIST check if element exist.
Definition: TsCtx.hpp:23
MoFEM::TsCtx::problemName
std::string problemName
Definition: TsCtx.hpp:22
MoFEM::TsCtx::loopsRHSJacobian
FEMethodsSequence loopsRHSJacobian
Definition: TsCtx.hpp:32
MoFEM::Interface
DeprecatedCoreInterface Interface
Definition: Interface.hpp:2010
MoFEM::TsMonitorSet
PetscErrorCode TsMonitorSet(TS ts, PetscInt step, PetscReal t, Vec u, void *ctx)
Set monitor for TS solver.
Definition: TsCtx.cpp:263
MoFEM::TSAdaptResetMoFEM
PetscErrorCode TSAdaptResetMoFEM(TSAdapt adapt)
Definition: TsCtx.cpp:818
MoFEM::TsCtx::TsSetRHSJacobian
friend PetscErrorCode TsSetRHSJacobian(TS ts, PetscReal t, Vec u, Mat A, Mat B, void *ctx)
TS solver function.
Definition: TsCtx.cpp:430
MoFEM::TsCtx::postProcessIJacobian
BasicMethodsSequence postProcessIJacobian
Definition: TsCtx.hpp:35
MoFEM::TsCtx::getPostProcessRHSJacobian
BasicMethodsSequence & getPostProcessRHSJacobian()
Get the postProcess to do RHSJacobian object.
Definition: TsCtx.hpp:164
MoFEM
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:10
MoFEM::TsCtx::TsSetIFunction
friend PetscErrorCode TsSetIFunction(TS ts, PetscReal t, Vec u, Vec u_t, Vec F, void *ctx)
Set IFunction for TS solver.
Definition: TsCtx.cpp:56
a
constexpr double a
Definition: approx_sphere.cpp:30
MoFEM::TsCtx::mField
MoFEM::Interface & mField
Definition: TsCtx.hpp:19
MoFEM::TsCtx::getPreProcessIJacobian
BasicMethodsSequence & getPreProcessIJacobian()
Get the preProcess to do IJacobian object.
Definition: TsCtx.hpp:125
MoFEM::TsCtx::BasicMethodsSequence
MoFEM::BasicMethodsSequence BasicMethodsSequence
Definition: TsCtx.hpp:27
MoFEM::TsCtx::postProcessMonitor
BasicMethodsSequence postProcessMonitor
Definition: TsCtx.hpp:39
MoFEM::TsCtx::TsMonitorSet
friend PetscErrorCode TsMonitorSet(TS ts, PetscInt step, PetscReal t, Vec u, void *ctx)
Set monitor for TS solver.
Definition: TsCtx.cpp:263
MoFEM::TsCtx::~TsCtx
virtual ~TsCtx()=default
MoFEM::TsCtx::loopsMonitor
FEMethodsSequence loopsMonitor
Definition: TsCtx.hpp:31
MoFEM::BasicMethodsSequence
std::deque< BasicMethodPtr > BasicMethodsSequence
Definition: AuxPETSc.hpp:57
MoFEM::PairNameFEMethodPtr
Definition: AuxPETSc.hpp:12
EshelbianPlasticity::P
@ P
Definition: EshelbianContact.cpp:201
h
double h
Definition: photon_diffusion.cpp:60
MoFEM::TsCtx::zeroMatrix
bool zeroMatrix
Definition: TsCtx.hpp:49
MoFEM::TsCtx::getPreProcessMonitor
BasicMethodsSequence & getPreProcessMonitor()
Get the preProcess to do Monitor object.
Definition: TsCtx.hpp:141
MoFEM::TsCtx
Interface for Time Stepping (TS) solver.
Definition: TsCtx.hpp:17
MoFEM::TsCtx::getPostProcessIFunction
BasicMethodsSequence & getPostProcessIFunction()
Get the postProcess to do IFunction object.
Definition: TsCtx.hpp:116
MoFEM::TsCtx::matAssembleSwitch
boost::movelib::unique_ptr< bool > matAssembleSwitch
Definition: TsCtx.hpp:221
MoFEM::TsCtx::preProcessMonitor
BasicMethodsSequence preProcessMonitor
Definition: TsCtx.hpp:38
MoFEM::TsCtx::postProcessIFunction
BasicMethodsSequence postProcessIFunction
Definition: TsCtx.hpp:37
MoFEM::TsCtx::MOFEM_EVENT_TsCtxRHSFunction
PetscLogEvent MOFEM_EVENT_TsCtxRHSFunction
Definition: TsCtx.hpp:212
MoFEM::TsCtx::getPreProcessRHSJacobian
BasicMethodsSequence & getPreProcessRHSJacobian()
Get the preProcess to do RHSJacobian object.
Definition: TsCtx.hpp:155
MoFEM::TsCtx::preProcessIFunction
BasicMethodsSequence preProcessIFunction
Definition: TsCtx.hpp:36
MoFEM::TsCtx::loopsIFunction
FEMethodsSequence loopsIFunction
Definition: TsCtx.hpp:30
MoFEM::TsSetIJacobian
PetscErrorCode TsSetIJacobian(TS ts, PetscReal t, Vec u, Vec u_t, PetscReal a, Mat A, Mat B, void *ctx)
Set function evaluating jacobian in TS solver.
Definition: TsCtx.cpp:169
MoFEM::TsCtx::preProcessIJacobian
BasicMethodsSequence preProcessIJacobian
Definition: TsCtx.hpp:34
MoFEM::TsCtx::PairNameFEMethodPtr
MoFEM::PairNameFEMethodPtr PairNameFEMethodPtr
Definition: TsCtx.hpp:25
t
constexpr double t
plate stiffness
Definition: plate.cpp:58
MoFEM::TsCtx::loopsRHSFunction
FEMethodsSequence loopsRHSFunction
Definition: TsCtx.hpp:33
MoFEM::TsCtx::MOFEM_EVENT_TsCtxRHSJacobian
PetscLogEvent MOFEM_EVENT_TsCtxRHSJacobian
Definition: TsCtx.hpp:213
MoFEM::TsCtx::getLoopsIJacobian
FEMethodsSequence & getLoopsIJacobian()
Get the loops to do IJacobian object.
Definition: TsCtx.hpp:83
MoFEM::TSAdaptChooseMoFEM
PetscErrorCode TSAdaptChooseMoFEM(TSAdapt adapt, TS ts, PetscReal h, PetscInt *next_sc, PetscReal *next_h, PetscBool *accept, PetscReal *wlte, PetscReal *wltea, PetscReal *wlter)
Definition: TsCtx.cpp:769
MoFEM::TsCtx::TsSetIJacobian
friend PetscErrorCode TsSetIJacobian(TS ts, PetscReal t, Vec u, Vec U_t, PetscReal a, Mat A, Mat B, void *ctx)
Set function evaluating jacobian in TS solver.
Definition: TsCtx.cpp:169
MoFEM::TsCtx::getPreProcessIFunction
BasicMethodsSequence & getPreProcessIFunction()
Get the preProcess to do IFunction object.
Definition: TsCtx.hpp:109
MoFEM::TsCtx::TsSetI2Jacobian
friend 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 Jacobian for second order PDE in time.
Definition: TsCtx.cpp:519
MoFEM::TsCtx::postProcessRHSFunction
BasicMethodsSequence postProcessRHSFunction
Definition: TsCtx.hpp:43
v
const double v
phase velocity of light in medium (cm/ns)
Definition: initial_diffusion.cpp:40
MoFEM::TsCtx::MOFEM_EVENT_TsCtxI2Jacobian
PetscLogEvent MOFEM_EVENT_TsCtxI2Jacobian
Definition: TsCtx.hpp:218
EigenMatrix::Vec
const FTensor::Tensor2< T, Dim, Dim > Vec
Definition: MatrixFunction.hpp:64
MoFEMTypes
MoFEMTypes
Those types control how functions respond on arguments, f.e. error handling.
Definition: definitions.h:110
MoFEM::TsCtx::postProcessRHSJacobian
BasicMethodsSequence postProcessRHSJacobian
Definition: TsCtx.hpp:42
MoFEM::TsSetRHSJacobian
PetscErrorCode TsSetRHSJacobian(TS ts, PetscReal t, Vec u, Mat A, Mat B, void *ctx)
TS solver function.
Definition: TsCtx.cpp:430
MoFEM::TsCtx::TsSetI2Function
friend 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.
Definition: TsCtx.cpp:620
MoFEM::TsSetRHSFunction
PetscErrorCode TsSetRHSFunction(TS ts, PetscReal t, Vec u, Vec F, void *ctx)
TS solver function.
Definition: TsCtx.cpp:327
MoFEM::TsSetI2Function
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.
Definition: TsCtx.cpp:620
MoFEM::TsCtx::TsCtx
TsCtx(MoFEM::Interface &m_field, const std::string &problem_name)
Definition: TsCtx.cpp:5
MoFEM::TsCtx::getPreProcessRHSFunction
BasicMethodsSequence & getPreProcessRHSFunction()
Get the preProcess to do RHSFunction object.
Definition: TsCtx.hpp:173
MoFEM::TsCtx::preProcessRHSFunction
BasicMethodsSequence preProcessRHSFunction
Definition: TsCtx.hpp:41
EshelbianPlasticity::U
@ U
Definition: EshelbianContact.cpp:201
MoFEM::TsCtx::tsDebugHook
boost::function< MoFEMErrorCode(TS ts, PetscReal t, Vec u, Vec u_t, Vec u_tt, Vec F, void *ctx)> tsDebugHook
Definition: TsCtx.hpp:47
F
@ F
Definition: free_surface.cpp:396
MoFEM::TSAdaptCreateMoFEM
PetscErrorCode TSAdaptCreateMoFEM(TSAdapt adapt)
Craete MOFEM adapt.
Definition: TsCtx.cpp:829