v0.8.23
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  }
77 
78 
79  /**
80  * @brief Get the loops to do IFunction object
81  *
82  * It is sequence of finite elements used to evaluate the right hand side of
83  * implicit time solver.
84  *
85  * @return FEMethodsSequence&
86  */
88  return loops_to_do_IFunction;
89  }
90 
91  /**
92  * @brief Get the loops to do RHSFunction object
93  *
94  * It is sequence of finite elements used to evaluate the right hand side of
95  * implicit time solver.
96  *
97  * @return FEMethodsSequence&
98  */
101  }
102 
103  /**
104  * @brief Get the loops to do IJacobian 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  return loops_to_do_IJacobian;
113  }
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 private:
239 
245 
246  boost::movelib::unique_ptr<bool> vecAssembleSwitch;
247  boost::movelib::unique_ptr<bool> matAssembleSwitch;
248 
249 };
250 
251 /**
252  * @brief Set IFunction for TS solver
253  *
254  * <a
255  * href=https://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/TS/TSSetIFunction.html>See
256  * petsc for details</a>
257  *
258  * @param ts
259  * @param t
260  * @param u
261  * @param u_t
262  * @param F
263  * @param ctx
264  * @return PetscErrorCode
265  */
266 PetscErrorCode TsSetIFunction(TS ts, PetscReal t, Vec u, Vec u_t, Vec F,
267  void *ctx);
268 
269 /**
270  * @brief Set function evaluating jacobina in TS solver
271  *
272  * <a
273  * href=https://www.mcs.anl.gov/petsc/petsc-3.1/docs/manualpages/TS/TSSetIJacobian.html>See
274  * PETSc for details</a>
275  *
276  * @param ts
277  * @param t
278  * @param u
279  * @param u_t
280  * @param a
281  * @param A
282  * @param B
283  * @param ctx
284  * @return PetscErrorCode
285  */
286 PetscErrorCode TsSetIJacobian(TS ts, PetscReal t, Vec u, Vec u_t, PetscReal a,
287  Mat A, Mat B, void *ctx);
288 
289 /**
290  * @brief Set monitor for TS solver
291  *
292  * <a
293  * href=https://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/TS/TSMonitorSet.html>See
294  * PETSc for details</a>
295  *
296  * @param ts
297  * @param step
298  * @param t
299  * @param u
300  * @param ctx
301  * @return PetscErrorCode
302  */
303 PetscErrorCode TsMonitorSet(TS ts, PetscInt step, PetscReal t, Vec u,
304  void *ctx);
305 
306 /// \deprecate Do not use, change to TsSetIFunction
307 DEPRECATED inline PetscErrorCode f_TSSetIFunction(TS ts, PetscReal t, Vec u,
308  Vec u_t, Vec F, void *ctx) {
309  return TsSetIFunction(ts, t, u, u_t, F, ctx);
310 }
311 
312 /// \deprecated Do not use, change to TsSetIJacobian
313 DEPRECATED inline PetscErrorCode f_TSSetIJacobian(TS ts, PetscReal t, Vec u,
314  Vec u_t, PetscReal a, Mat A,
315  Mat B, void *ctx) {
316  return TsSetIJacobian(ts, t, u, u_t, a, A, B, ctx);
317 }
318 
319 /// \deprecated Do not use, change to TsMonitorSet
320 DEPRECATED inline PetscErrorCode f_TSMonitorSet(TS ts, PetscInt step,
321  PetscReal t, Vec u, void *ctx) {
322  return TsMonitorSet(ts, step, t, u, ctx);
323 }
324 
325 /**
326  * @brief TS solver function
327  *
328  * <a
329  * href=https://www.mcs.anl.gov/petsc/petsc-3.11/docs/manualpages/TS/TSSetRHSFunction.html#TSSetRHSFunction>See
330  * PETSc for details</a>
331  *
332  * @param ts
333  * @param t
334  * @param u
335  * @param F
336  * @param ctx
337  * @return PetscErrorCode
338  */
339 PetscErrorCode TsSetRHSFunction(TS ts, PetscReal t, Vec u, Vec F, void *ctx);
340 
341 /**
342  * @brief TS solver function
343  *
344  * <a
345  * href=https://www.mcs.anl.gov/petsc/petsc-3.11/docs/manualpages/TS/TSSetRHSJacobian.html#TSSetRHSJacobian>See
346  * PETSc for details</a>
347  *
348  * @param ts
349  * @param t
350  * @param u
351  * @param A
352  * @param B
353  * @param ctx
354  * @return PetscErrorCode
355  */
356 PetscErrorCode TsSetRHSJacobian(TS ts, PetscReal t, Vec u, Mat A, Mat B,
357  void *ctx);
358 
359 } // namespace MoFEM
360 
361 #endif // __TSCTX_HPP__
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:242
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:243
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:246
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:307
PetscLogEvent MOFEM_EVENT_TsCtxRHSFunction
Definition: TsCtx.hpp:240
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:108
FEMethodsSequence loops_to_do_RHSJacobian
Definition: TsCtx.hpp:51
BasicMethodsSequence postProcess_IFunction
Definition: TsCtx.hpp:56
BasicMethodsSequence preProcess_IFunction
Definition: TsCtx.hpp:55
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:107
PetscLogEvent MOFEM_EVENT_TsCtxMonitor
Definition: TsCtx.hpp:244
boost::movelib::unique_ptr< bool > matAssembleSwitch
Definition: TsCtx.hpp:247
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:313
MoFEM::FEMethodsSequence FEMethodsSequence
Definition: TsCtx.hpp:45
DEPRECATED PetscErrorCode f_TSMonitorSet(TS ts, PetscInt step, PetscReal t, Vec u, void *ctx)
Definition: TsCtx.hpp:320
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:111
PetscLogEvent MOFEM_EVENT_TsCtxRHSJacobian
Definition: TsCtx.hpp:241
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:87
MoFEMTypes
Those types control how functions respond on arguments, f.e. error handling.
Definition: definitions.h:182
FEMethodsSequence & get_loops_to_do_RHSFunction()
Get the loops to do RHSFunction object.
Definition: TsCtx.hpp:99
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