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