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 
57 
62  PetscLogEvent MOFEM_EVENT_TsCtxMonitor;
63 
64  bool zeroMatrix;
65  TsCtx(MoFEM::Interface &m_field, const std::string &problem_name)
66  : mField(m_field), moab(m_field.get_moab()), problemName(problem_name),
67  bH(MF_EXIST), zeroMatrix(true) {
68  PetscLogEventRegister("LoopTsIFunction", 0, &MOFEM_EVENT_TsCtxIFunction);
69  PetscLogEventRegister("LoopTsIJacobian", 0, &MOFEM_EVENT_TsCtxIJacobian);
70  PetscLogEventRegister("LoopTsRHSFunction", 0,
72  PetscLogEventRegister("LoopTsRHSJacobian", 0,
74  PetscLogEventRegister("LoopTsMonitor", 0, &MOFEM_EVENT_TsCtxMonitor);
75  }
76 
77  /**
78  * @brief Get the loops to do IFunction object
79  *
80  * It is sequence of finite elements used to evaluate the right hand side of
81  * implicit time solver.
82  *
83  * @return FEMethodsSequence&
84  */
86  return loops_to_do_IFunction;
87  }
88 
89  /**
90  * @brief Get the loops to do IJacobian object
91  *
92  * It is sequence of finite elements used to evalite the left hand sie of
93  * implimcit time solver.
94  *
95  * @return FEMethodsSequence&
96  */
98  return loops_to_do_IJacobian;
99  }
100 
101  /**
102  * @brief Get the loops to do Monitor object
103  *
104  * It is sequence used to monitor solution of time solver.
105  *
106  * @return FEMethodsSequence&
107  */
109 
110  /**
111  * @brief Get the preProcess to do IFunction object
112  *
113  * @return BasicMethodsSequence&
114  */
116  return preProcess_IFunction;
117  }
118 
119  /**
120  * @brief Get the postProcess to do IFunction object
121  *
122  * @return BasicMethodsSequence&
123  */
125  return postProcess_IFunction;
126  }
127 
128  /**
129  * @brief Get the preProcess to do IJacobian object
130  *
131  * @return BasicMethodsSequence&
132  */
134  return preProcess_IJacobian;
135  }
136 
137  /**
138  * @brief Get the postProcess to do IJacobian object
139  *
140  * @return BasicMethodsSequence&
141  */
143  return postProcess_IJacobian;
144  }
145 
146  /**
147  * @brief Get the preProcess to do Monitor object
148  *
149  * @return BasicMethodsSequence&
150  */
152  return preProcess_Monitor;
153  }
154 
155  /**
156  * @brief Get the postProcess to do Monitor object
157  *
158  * @return BasicMethodsSequence&
159  */
161  return postProcess_Monitor;
162  }
163 
164  friend PetscErrorCode TsSetIFunction(TS ts, PetscReal t, Vec u, Vec u_t,
165  Vec F, void *ctx);
166  friend PetscErrorCode TsSetIJacobian(TS ts, PetscReal t, Vec u, Vec U_t,
167  PetscReal a, Mat A, Mat B, void *ctx);
168  friend PetscErrorCode TsMonitorSet(TS ts, PetscInt step, PetscReal t, Vec u,
169  void *ctx);
170 };
171 
172 /**
173  * @brief Set IFunction for TS solver
174  *
175  * <a href=https://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/TS/TSSetIFunction.html>See petsc for details</a>
176  *
177  * @param ts
178  * @param t
179  * @param u
180  * @param u_t
181  * @param F
182  * @param ctx
183  * @return PetscErrorCode
184  */
185 PetscErrorCode TsSetIFunction(TS ts, PetscReal t, Vec u, Vec u_t, Vec F,
186  void *ctx);
187 
188 /**
189  * @brief Set function evaluating jacobina in TS solver
190  *
191  * <a href=https://www.mcs.anl.gov/petsc/petsc-3.1/docs/manualpages/TS/TSSetIJacobian.html>See PETSc for details</a>
192  *
193  * @param ts
194  * @param t
195  * @param u
196  * @param u_t
197  * @param a
198  * @param A
199  * @param B
200  * @param ctx
201  * @return PetscErrorCode
202  */
203 PetscErrorCode TsSetIJacobian(TS ts, PetscReal t, Vec u, Vec u_t, PetscReal a,
204  Mat A, Mat B, void *ctx);
205 
206 /**
207  * @brief Set monitor for TS solver
208  *
209  * <a href=https://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/TS/TSMonitorSet.html>See PETSc for details</a>
210  *
211  * @param ts
212  * @param step
213  * @param t
214  * @param u
215  * @param ctx
216  * @return PetscErrorCode
217  */
218 PetscErrorCode TsMonitorSet(TS ts, PetscInt step, PetscReal t, Vec u,
219  void *ctx);
220 
221 DEPRECATED inline PetscErrorCode f_TSSetIFunction(TS ts, PetscReal t, Vec u,
222  Vec u_t, Vec F, void *ctx) {
223  return TsSetIFunction(ts, t, u, u_t, F, ctx);
224 }
225 
226 DEPRECATED inline PetscErrorCode f_TSSetIJacobian(TS ts, PetscReal t, Vec u,
227  Vec u_t, PetscReal a, Mat A,
228  Mat B, void *ctx) {
229  return TsSetIJacobian(ts, t, u, u_t, a, A, B, ctx);
230 }
231 
232 DEPRECATED inline PetscErrorCode f_TSMonitorSet(TS ts, PetscInt step,
233  PetscReal t, Vec u, void *ctx) {
234  return TsMonitorSet(ts, step, t, u, ctx);
235 }
236 
237 } // namespace MoFEM
238 
239 #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:133
PetscLogEvent MOFEM_EVENT_TsCtxIFunction
Definition: TsCtx.hpp:60
bool zeroMatrix
Definition: TsCtx.hpp:64
BasicMethodsSequence & get_preProcess_to_do_Monitor()
Get the preProcess to do Monitor object.
Definition: TsCtx.hpp:151
DEPRECATED typedef MoFEM::FEMethodsSequence loops_to_do_type
Definition: TsCtx.hpp:39
PetscLogEvent MOFEM_EVENT_TsCtxIJacobian
Definition: TsCtx.hpp:61
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:160
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:51
PetscErrorCode TsMonitorSet(TS ts, PetscInt step, PetscReal t, Vec u, void *ctx)
Set monitor for TS solver.
Definition: TsCtx.cpp:175
std::string problemName
Definition: TsCtx.hpp:32
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:103
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)
Definition: TsCtx.hpp:221
PetscLogEvent MOFEM_EVENT_TsCtxRHSFunction
Definition: TsCtx.hpp:58
PetscErrorCode TsSetIFunction(TS ts, PetscReal t, Vec u, Vec u_t, Vec F, void *ctx)
Set IFunction for TS solver.
Definition: TsCtx.cpp:23
TsCtx(MoFEM::Interface &m_field, const std::string &problem_name)
Definition: TsCtx.hpp:65
std::vector< BasicMethodPtr > BasicMethodsSequence
Definition: AuxPETSc.hpp:65
BasicMethodsSequence postProcess_IFunction
Definition: TsCtx.hpp:54
BasicMethodsSequence preProcess_IFunction
Definition: TsCtx.hpp:53
BasicMethodsSequence & get_preProcess_to_do_IFunction()
Get the preProcess to do IFunction object.
Definition: TsCtx.hpp:115
BasicMethodsSequence postProcess_IJacobian
Definition: TsCtx.hpp:52
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:142
FEMethodsSequence & get_loops_to_do_Monitor()
Get the loops to do Monitor object.
Definition: TsCtx.hpp:108
std::vector< PairNameFEMethodPtr > FEMethodsSequence
Definition: AuxPETSc.hpp:64
PetscLogEvent MOFEM_EVENT_TsCtxMonitor
Definition: TsCtx.hpp:62
friend PetscErrorCode TsSetIFunction(TS ts, PetscReal t, Vec u, Vec u_t, Vec F, void *ctx)
Set IFunction for TS solver.
Definition: TsCtx.cpp:23
#define DEPRECATED
Definition: definitions.h:29
BasicMethodsSequence & get_postProcess_to_do_IFunction()
Get the postProcess to do IFunction object.
Definition: TsCtx.hpp:124
DEPRECATED typedef MoFEM::BasicMethodsSequence basic_method_to_do
Definition: TsCtx.hpp:42
FEMethodsSequence loops_to_do_IFunction
Definition: TsCtx.hpp:49
MoFEM::Interface & mField
Definition: TsCtx.hpp:29
BasicMethodsSequence preProcess_Monitor
Definition: TsCtx.hpp:55
DEPRECATED PetscErrorCode f_TSSetIJacobian(TS ts, PetscReal t, Vec u, Vec u_t, PetscReal a, Mat A, Mat B, void *ctx)
Definition: TsCtx.hpp:226
MoFEM::FEMethodsSequence FEMethodsSequence
Definition: TsCtx.hpp:45
DEPRECATED PetscErrorCode f_TSMonitorSet(TS ts, PetscInt step, PetscReal t, Vec u, void *ctx)
Definition: TsCtx.hpp:232
FEMethodsSequence & get_loops_to_do_IJacobian()
Get the loops to do IJacobian object.
Definition: TsCtx.hpp:97
PetscLogEvent MOFEM_EVENT_TsCtxRHSJacobian
Definition: TsCtx.hpp:59
FEMethodsSequence & get_loops_to_do_IFunction()
Get the loops to do IFunction object.
Definition: TsCtx.hpp:85
MoFEMTypes
Those types control how functions respond on arguments, f.e. error handling.
Definition: definitions.h:182
BasicMethodsSequence postProcess_Monitor
Definition: TsCtx.hpp:56
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:103
friend PetscErrorCode TsMonitorSet(TS ts, PetscInt step, PetscReal t, Vec u, void *ctx)
Set monitor for TS solver.
Definition: TsCtx.cpp:175