v0.9.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 petsc_context_struture
26  */
27 struct TsCtx {
28 
30  moab::Interface &moab;
31 
32  std::string problemName;
33  MoFEMTypes bH; ///< If set to MF_EXIST check if element exist
34 
35  /// \deprecated use PairNameFEMethodPtr
37 
38  /// \deprecated use FEMethodsSequence
40 
41  /// \deprecated use BasicMethodsSequence
43 
47 
63 
64  bool zeroMatrix;
65 
66  TsCtx(MoFEM::Interface &m_field, const std::string &problem_name)
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  }
79 
80  /**
81  * @brief Get the loops to do IFunction object
82  *
83  * It is sequence of finite elements used to evaluate the right hand side of
84  * implicit time solver.
85  *
86  * @return FEMethodsSequence&
87  */
89  return loops_to_do_IFunction;
90  }
91 
92  /**
93  * @brief Get the loops to do RHSFunction object
94  *
95  * It is sequence of finite elements used to evaluate the right hand side of
96  * implicit time solver.
97  *
98  * @return FEMethodsSequence&
99  */
102  }
103 
104  /**
105  * @brief Get the loops to do IJacobian object
106  *
107  * It is sequence of finite elements used to evalite the left hand sie of
108  * implimcit time solver.
109  *
110  * @return FEMethodsSequence&
111  */
113  return loops_to_do_IJacobian;
114  }
115 
116  /**
117  * @brief Get the loops to do RHSJacobian object
118  *
119  * It is sequence of finite elements used to evalite the left hand sie of
120  * implimcit time solver.
121  *
122  * @return FEMethodsSequence&
123  */
126  }
127 
128  /**
129  * @brief Get the loops to do Monitor object
130  *
131  * It is sequence used to monitor solution of time solver.
132  *
133  * @return FEMethodsSequence&
134  */
136 
137  /**
138  * @brief Get the preProcess to do IFunction object
139  *
140  * @return BasicMethodsSequence&
141  */
143  return preProcess_IFunction;
144  }
145 
146  /**
147  * @brief Get the postProcess to do IFunction object
148  *
149  * @return BasicMethodsSequence&
150  */
152  return postProcess_IFunction;
153  }
154 
155  /**
156  * @brief Get the preProcess to do IJacobian object
157  *
158  * @return BasicMethodsSequence&
159  */
161  return preProcess_IJacobian;
162  }
163 
164  /**
165  * @brief Get the postProcess to do IJacobian object
166  *
167  * @return BasicMethodsSequence&
168  */
170  return postProcess_IJacobian;
171  }
172 
173  /**
174  * @brief Get the preProcess to do Monitor object
175  *
176  * @return BasicMethodsSequence&
177  */
179  return preProcess_Monitor;
180  }
181 
182  /**
183  * @brief Get the postProcess to do Monitor object
184  *
185  * @return BasicMethodsSequence&
186  */
188  return postProcess_Monitor;
189  }
190 
191  /**
192  * @brief Get the preProcess to do RHSJacobian object
193  *
194  * @return BasicMethodsSequence&
195  */
197  return preProcess_RHSJacobian;
198  }
199 
200  /**
201  * @brief Get the postProcess to do RHSJacobian object
202  *
203  * @return BasicMethodsSequence&
204  */
207  }
208 
209  /**
210  * @brief Get the preProcess to do RHSFunction object
211  *
212  * @return BasicMethodsSequence&
213  */
215  return preProcess_RHSFunction;
216  }
217 
218  /**
219  * @brief Get the postProcess to do RHSFunction object
220  *
221  * @return BasicMethodsSequence&
222  */
225  }
226 
227  friend PetscErrorCode TsSetIFunction(TS ts, PetscReal t, Vec u, Vec u_t,
228  Vec F, void *ctx);
229  friend PetscErrorCode TsSetIJacobian(TS ts, PetscReal t, Vec u, Vec U_t,
230  PetscReal a, Mat A, Mat B, void *ctx);
231  friend PetscErrorCode TsMonitorSet(TS ts, PetscInt step, PetscReal t, Vec u,
232  void *ctx);
233  friend PetscErrorCode TsSetRHSFunction(TS ts, PetscReal t, Vec u, Vec F,
234  void *ctx);
235  friend PetscErrorCode TsSetRHSJacobian(TS ts, PetscReal t, Vec u, Mat A,
236  Mat B, void *ctx);
237 
238  friend PetscErrorCode TsSetI2Function(TS ts, PetscReal t, Vec U, Vec U_t,
239  Vec U_tt, Vec F, void *ctx);
240 
241  friend PetscErrorCode TsSetI2Jacobian(TS ts, PetscReal t, Vec U, Vec U_t,
242  Vec U_tt, PetscReal v, PetscReal a,
243  Mat J, Mat P, void *ctx);
244 
245 private:
253 
254  boost::movelib::unique_ptr<bool> vecAssembleSwitch;
255  boost::movelib::unique_ptr<bool> matAssembleSwitch;
256 };
257 
258 /**
259  * @brief Set IFunction for TS solver
260  *
261  * <a
262  * href=https://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/TS/TSSetIFunction.html>See
263  * petsc for details</a>
264  *
265  * @param ts
266  * @param t
267  * @param u
268  * @param u_t
269  * @param F
270  * @param ctx
271  * @return PetscErrorCode
272  */
273 PetscErrorCode TsSetIFunction(TS ts, PetscReal t, Vec u, Vec u_t, Vec F,
274  void *ctx);
275 
276 /**
277  * @brief Set function evaluating jacobina in TS solver
278  *
279  * <a
280  * href=https://www.mcs.anl.gov/petsc/petsc-3.1/docs/manualpages/TS/TSSetIJacobian.html>See
281  * PETSc for details</a>
282  *
283  * @param ts
284  * @param t
285  * @param u
286  * @param u_t
287  * @param a
288  * @param A
289  * @param B
290  * @param ctx
291  * @return PetscErrorCode
292  */
293 PetscErrorCode TsSetIJacobian(TS ts, PetscReal t, Vec u, Vec u_t, PetscReal a,
294  Mat A, Mat B, void *ctx);
295 
296 /**
297  * @brief Set monitor for TS solver
298  *
299  * <a
300  * href=https://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/TS/TSMonitorSet.html>See
301  * PETSc for details</a>
302  *
303  * @param ts
304  * @param step
305  * @param t
306  * @param u
307  * @param ctx
308  * @return PetscErrorCode
309  */
310 PetscErrorCode TsMonitorSet(TS ts, PetscInt step, PetscReal t, Vec u,
311  void *ctx);
312 
313 /// \deprecate Do not use, change to TsSetIFunction
314 DEPRECATED inline PetscErrorCode f_TSSetIFunction(TS ts, PetscReal t, Vec u,
315  Vec u_t, Vec F, void *ctx) {
316  return TsSetIFunction(ts, t, u, u_t, F, ctx);
317 }
318 
319 /// \deprecated Do not use, change to TsSetIJacobian
320 DEPRECATED inline PetscErrorCode f_TSSetIJacobian(TS ts, PetscReal t, Vec u,
321  Vec u_t, PetscReal a, Mat A,
322  Mat B, void *ctx) {
323  return TsSetIJacobian(ts, t, u, u_t, a, A, B, ctx);
324 }
325 
326 /// \deprecated Do not use, change to TsMonitorSet
327 DEPRECATED inline PetscErrorCode f_TSMonitorSet(TS ts, PetscInt step,
328  PetscReal t, Vec u, void *ctx) {
329  return TsMonitorSet(ts, step, t, u, ctx);
330 }
331 
332 /**
333  * @brief TS solver function
334  *
335  * <a
336  * href=https://www.mcs.anl.gov/petsc/petsc-3.11/docs/manualpages/TS/TSSetRHSFunction.html#TSSetRHSFunction>See
337  * PETSc for details</a>
338  *
339  * @param ts
340  * @param t
341  * @param u
342  * @param F
343  * @param ctx
344  * @return PetscErrorCode
345  */
346 PetscErrorCode TsSetRHSFunction(TS ts, PetscReal t, Vec u, Vec F, void *ctx);
347 
348 /**
349  * @brief TS solver function
350  *
351  * <a
352  * href=https://www.mcs.anl.gov/petsc/petsc-3.11/docs/manualpages/TS/TSSetRHSJacobian.html#TSSetRHSJacobian>See
353  * PETSc for details</a>
354  *
355  * @param ts
356  * @param t
357  * @param u
358  * @param A
359  * @param B
360  * @param ctx
361  * @return PetscErrorCode
362  */
363 PetscErrorCode TsSetRHSJacobian(TS ts, PetscReal t, Vec u, Mat A, Mat B,
364  void *ctx);
365 
366 /**
367  * @brief Calculation Jaconian for second order PDE in time
368  *
369  * <a
370  * href=https://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/TS/TSSetI2Jacobian.html>PETSc
371  * for details</a>
372  *
373  * @param ts
374  * @param J
375  * @param P
376  * @param jac
377  * @param ctx
378  * @return PetscErrorCode
379  */
380 PetscErrorCode TsSetI2Jacobian(TS ts, PetscReal t, Vec u, Vec u_t, Vec u_tt,
381  PetscReal v, PetscReal a, Mat J, Mat P,
382  void *ctx);
383 
384 /**
385  * @brief Calculation the right hand side for second order PDE in time
386  *
387  * <a
388  * href=https://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/TS/TSSetI2Function.html>PETSc
389  * for details</a>
390  *
391  * @param ts
392  * @param F
393  * @param fun
394  * @param ctx
395  * @return PetscErrorCode
396  */
397 PetscErrorCode TsSetI2Function(TS ts, PetscReal t, Vec u, Vec u_t, Vec u_tt,
398  Vec F, void *ctx);
399 
400 } // namespace MoFEM
401 
402 #endif // __TSCTX_HPP__
PetscErrorCode TsSetI2Jacobian(TS ts, PetscReal t, Vec u, Vec u_t, Vec u_tt, PetscReal v, PetscReal a, Mat A, Mat B, void *ctx)
Calculation Jaconian for second order PDE in time.
Definition: TsCtx.cpp:401
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:401
FEMethodsSequence loops_to_do_Monitor
Definition: TsCtx.hpp:50
BasicMethodsSequence & get_preProcess_to_do_IJacobian()
Get the preProcess to do IJacobian object.
Definition: TsCtx.hpp:160
PetscLogEvent MOFEM_EVENT_TsCtxIFunction
Definition: TsCtx.hpp:248
bool zeroMatrix
Definition: TsCtx.hpp:64
BasicMethodsSequence preProcess_RHSJacobian
Definition: TsCtx.hpp:59
friend PetscErrorCode TsSetRHSJacobian(TS ts, PetscReal t, Vec u, Mat A, Mat B, void *ctx)
TS solver function.
Definition: TsCtx.cpp:325
BasicMethodsSequence & get_preProcess_to_do_Monitor()
Get the preProcess to do Monitor object.
Definition: TsCtx.hpp:178
DEPRECATED typedef MoFEM::FEMethodsSequence loops_to_do_type
Definition: TsCtx.hpp:39
BasicMethodsSequence postProcess_RHSFunction
Definition: TsCtx.hpp:62
PetscLogEvent MOFEM_EVENT_TsCtxIJacobian
Definition: TsCtx.hpp:249
DEPRECATED typedef MoFEM::PairNameFEMethodPtr loop_pair_type
Definition: TsCtx.hpp:36
BasicMethodsSequence & get_postProcess_to_do_Monitor()
Get the postProcess to do Monitor object.
Definition: TsCtx.hpp:187
Interface for Time Stepping (TS) solver.
Definition: TsCtx.hpp:27
MoFEMTypes bH
If set to MF_EXIST check if element exist.
Definition: TsCtx.hpp:33
BasicMethodsSequence preProcess_IJacobian
Definition: TsCtx.hpp:53
PetscErrorCode TsSetRHSFunction(TS ts, PetscReal t, Vec u, Vec F, void *ctx)
TS solver function.
Definition: TsCtx.cpp:241
BasicMethodsSequence & get_postProcess_to_do_RHSFunction()
Get the postProcess to do RHSFunction object.
Definition: TsCtx.hpp:223
PetscErrorCode TsMonitorSet(TS ts, PetscInt step, PetscReal t, Vec u, void *ctx)
Set monitor for TS solver.
Definition: TsCtx.cpp:190
std::string problemName
Definition: TsCtx.hpp:32
boost::movelib::unique_ptr< bool > vecAssembleSwitch
Definition: TsCtx.hpp:254
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:108
moab::Interface & moab
Definition: TsCtx.hpp:30
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:21
MoFEM::PairNameFEMethodPtr PairNameFEMethodPtr
Definition: TsCtx.hpp:44
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:314
PetscLogEvent MOFEM_EVENT_TsCtxRHSFunction
Definition: TsCtx.hpp:246
BasicMethodsSequence & get_preProcess_to_do_RHSFunction()
Get the preProcess to do RHSFunction object.
Definition: TsCtx.hpp:214
PetscErrorCode TsSetIFunction(TS ts, PetscReal t, Vec u, Vec u_t, Vec F, void *ctx)
Set IFunction for TS solver.
Definition: TsCtx.cpp:17
TsCtx(MoFEM::Interface &m_field, const std::string &problem_name)
Definition: TsCtx.hpp:66
std::vector< BasicMethodPtr > BasicMethodsSequence
Definition: AuxPETSc.hpp:112
FEMethodsSequence loops_to_do_RHSJacobian
Definition: TsCtx.hpp:51
BasicMethodsSequence postProcess_IFunction
Definition: TsCtx.hpp:56
BasicMethodsSequence preProcess_IFunction
Definition: TsCtx.hpp:55
PetscLogEvent MOFEM_EVENT_TsCtxI2Jacobian
Definition: TsCtx.hpp:252
BasicMethodsSequence & get_preProcess_to_do_IFunction()
Get the preProcess to do IFunction object.
Definition: TsCtx.hpp:142
BasicMethodsSequence postProcess_IJacobian
Definition: TsCtx.hpp:54
BasicMethodsSequence preProcess_RHSFunction
Definition: TsCtx.hpp:60
FEMethodsSequence loops_to_do_IJacobian
Definition: TsCtx.hpp:48
MoFEM::BasicMethodsSequence BasicMethodsSequence
Definition: TsCtx.hpp:46
BasicMethodsSequence & get_postProcess_to_do_IJacobian()
Get the postProcess to do IJacobian object.
Definition: TsCtx.hpp:169
friend PetscErrorCode TsSetRHSFunction(TS ts, PetscReal t, Vec u, Vec F, void *ctx)
TS solver function.
Definition: TsCtx.cpp:241
FEMethodsSequence & get_loops_to_do_Monitor()
Get the loops to do Monitor object.
Definition: TsCtx.hpp:135
std::vector< PairNameFEMethodPtr > FEMethodsSequence
Definition: AuxPETSc.hpp:111
PetscLogEvent MOFEM_EVENT_TsCtxI2Function
Definition: TsCtx.hpp:251
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:492
PetscLogEvent MOFEM_EVENT_TsCtxMonitor
Definition: TsCtx.hpp:250
boost::movelib::unique_ptr< bool > matAssembleSwitch
Definition: TsCtx.hpp:255
BasicMethodsSequence & get_preProcess_to_do_RHSJacobian()
Get the preProcess to do RHSJacobian object.
Definition: TsCtx.hpp:196
friend PetscErrorCode TsSetIFunction(TS ts, PetscReal t, Vec u, Vec u_t, Vec F, void *ctx)
Set IFunction for TS solver.
Definition: TsCtx.cpp:17
#define DEPRECATED
Definition: definitions.h:29
BasicMethodsSequence & get_postProcess_to_do_IFunction()
Get the postProcess to do IFunction object.
Definition: TsCtx.hpp:151
DEPRECATED typedef MoFEM::BasicMethodsSequence basic_method_to_do
Definition: TsCtx.hpp:42
PetscErrorCode TsSetRHSJacobian(TS ts, PetscReal t, Vec u, Mat A, Mat B, void *ctx)
TS solver function.
Definition: TsCtx.cpp:325
FEMethodsSequence loops_to_do_IFunction
Definition: TsCtx.hpp:49
BasicMethodsSequence & get_postProcess_to_do_RHSJacobian()
Get the postProcess to do RHSJacobian object.
Definition: TsCtx.hpp:205
MoFEM::Interface & mField
Definition: TsCtx.hpp:29
BasicMethodsSequence preProcess_Monitor
Definition: TsCtx.hpp:57
DEPRECATED PetscErrorCode f_TSSetIJacobian(TS ts, PetscReal t, Vec u, Vec u_t, PetscReal a, Mat A, Mat B, void *ctx)
Definition: TsCtx.hpp:320
MoFEM::FEMethodsSequence FEMethodsSequence
Definition: TsCtx.hpp:45
DEPRECATED PetscErrorCode f_TSMonitorSet(TS ts, PetscInt step, PetscReal t, Vec u, void *ctx)
Definition: TsCtx.hpp:327
FEMethodsSequence loops_to_do_RHSFunction
Definition: TsCtx.hpp:52
FEMethodsSequence & get_loops_to_do_IJacobian()
Get the loops to do IJacobian object.
Definition: TsCtx.hpp:112
PetscLogEvent MOFEM_EVENT_TsCtxRHSJacobian
Definition: TsCtx.hpp:247
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:492
FEMethodsSequence & get_loops_to_do_RHSJacobian()
Get the loops to do RHSJacobian object.
Definition: TsCtx.hpp:124
FEMethodsSequence & get_loops_to_do_IFunction()
Get the loops to do IFunction object.
Definition: TsCtx.hpp:88
MoFEMTypes
Those types control how functions respond on arguments, f.e. error handling.
Definition: definitions.h:183
FEMethodsSequence & get_loops_to_do_RHSFunction()
Get the loops to do RHSFunction object.
Definition: TsCtx.hpp:100
BasicMethodsSequence postProcess_Monitor
Definition: TsCtx.hpp:58
BasicMethodsSequence postProcess_RHSJacobian
Definition: TsCtx.hpp:61
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:108
friend PetscErrorCode TsMonitorSet(TS ts, PetscInt step, PetscReal t, Vec u, void *ctx)
Set monitor for TS solver.
Definition: TsCtx.cpp:190