v0.15.0
Loading...
Searching...
No Matches
SnesCtx.hpp
Go to the documentation of this file.
1/**
2 * \file SnesCtx.hpp
3 * \brief Context for PETSc SNES, i.e. nonlinear solver
4 * \author Anonymous contributor under MIT license
5 */
6
7#ifndef __SNESCTX_HPP__
8#define __SNESCTX_HPP__
9
10namespace MoFEM {
11
12/** \brief Interface for nonlinear (SNES) solver
13 * \ingroup mofem_petsc_solvers
14 */
15struct SnesCtx {
16
20 using HookFunction = boost::function<MoFEMErrorCode(
21 SNES snes, Vec x, Vec F, boost::shared_ptr<CacheTuple> cache_ptr,
22 void *ctx)>;
23
24 SnesCtx(Interface &m_field, std::string problem_name);
25
26 virtual ~SnesCtx();
27
28 /**
29 * @return return reference to vector with FEMethod to calculate tangent
30 * matrix
31 */
33
34 /**
35 * @return return vector to vector with FEMethod to calculate residual
36 */
38
39 /**
40 * The sequence of BasicMethod is executed before residual is calculated. It
41 * can be used to setup data structures, e.g. zero global variable which is
42 * integrated in domain, e.g. for calculation of strain energy.
43 *
44 * @return reference to BasicMethod for preprocessing
45 */
47
48 /**
49 * The sequence of BasicMethod is executed after residual is calculated. It
50 * can be used to setup data structures, e.g. aggregate data from processors
51 * or to apply essential boundary conditions.
52 *
53 * @return reference to BasicMethod for postprocessing
54 */
56
57 /**
58 * @return reference to BasicMethod for preprocessing
59 */
61
62 /**
63 * The sequence of BasicMethod is executed after tangent matrix is calculated.
64 * It can be used to setup data structures, e.g. aggregate data from
65 * processors or to apply essential boundary conditions.
66 *
67 * @return reference to BasicMethod for postprocessing
68 */
70
71 /**
72 * @brief Get the finite element pipeline for FunctionFn, that calculate
73 * tangent of function load for arc length control
74 *
75 * @return FEMethodsSequence&
76 */
78
79 /**
80 * @brief Get the BasicMethod sequence for preprocessing of FunctionFn
81 *
82 * @return BasicMethodsSequence&
83 */
85
86 /**
87 * @brief Get the BasicMethod sequence for postprocessing of FunctionFn
88 *
89 * @return BasicMethodsSequence&
90 */
92
93 /**
94 * @brief Get the Load Tangent Hook function
95 *
96 * Using hook you can debug or post-process while load vector tangent is
97 * evaluated
98 *
99 * @return HookFunction&
100 */
102
103 /**
104 * @brief Get the Right Hand Side Hook function
105 *
106 * Using hook you can debug or post-process while right hand side vector is
107 * evaluated
108 *
109 * @return HookFunction&
110 */
112
113 /**
114 * \brief Copy sequences from other SNES context
115 * @param snes_ctx SNES context from which Sequence is copied from
116 * @return error code
117 */
118 MoFEMErrorCode copyLoops(const SnesCtx &snes_ctx);
119
120 /**
121 * @brief Clear loops
122 *
123 * @return MoFEMErrorCode
124 */
126
127 friend PetscErrorCode SnesRhs(SNES snes, Vec x, Vec f, void *ctx);
128 friend PetscErrorCode SnesMat(SNES snes, Vec x, Mat A, Mat B, void *ctx);
129
131 MatAssemblyType type);
132 friend MoFEMErrorCode SnesMoFEMSetBehavior(SNES snes, MoFEMTypes bh);
133
135 MatAssemblyType type);
136
137 friend MoFEMErrorCode MoFEMSNESMonitorFields(SNES snes, PetscInt its,
138 PetscReal fgnorm,
139 SnesCtx *snes_ctx);
140
141 friend MoFEMErrorCode MoFEMSNESMonitorEnergy(SNES snes, PetscInt its,
142 PetscReal fgnorm,
143 SnesCtx *snes_ctx);
144
145 friend PetscErrorCode SnesLoadTangent(SNES snes, Vec u, Vec F, void *ctx);
146
147 struct SnesCtxImpl;
148
149protected:
150 boost::movelib::unique_ptr<SnesCtxImpl> snesCtxImpl;
151};
152
153/**
154 * \brief This is MoFEM implementation for the right hand side (residual vector)
155 * evaluation in SNES solver
156 *
157 * For more information pleas look to PETSc manual, i.e. SNESSetFunction
158 * <http://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/SNES/SNESSetFunction.html>
159 *
160 * @param snes SNES solver
161 * @param x Solution vector at current iteration
162 * @param f The right hand side vector
163 * @param ctx Pointer to context i.e. SnesCtx
164 * @return Error code
165 */
166PetscErrorCode SnesRhs(SNES snes, Vec x, Vec f, void *ctx);
167
168/**
169 * \brief This is MoFEM implementation for the left hand side (tangent matrix)
170 * evaluation in SNES solver
171 *
172 * For more information pleas look to PETSc manual, i.e. SNESSetJacobian
173 * <http://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/SNES/SNESSetJacobian.html#SNESSetJacobian>
174 *
175 * @param snes SNES solver
176 * @param x Solution vector at current iteration
177 * @param A Tangent matrix
178 * @param B Preconditioner tangent matrix
179 * @param ctx Pointer to context i.e. SnesCtx
180 * @return Error code
181 */
182PetscErrorCode SnesMat(SNES snes, Vec x, Mat A, Mat B, void *ctx);
183
184/**
185 * \brief This function calls finite element pipeline to compute tangent of of
186 * load vector in arc length control
187 *
188 * For more information pleas look to PETSc manual, i.e. SNESSetFunction
189 * <https://petsc.org/main/manualpages/SNES/SNESSetFunction/>
190 *
191 * @param snes SNES solver
192 * @param u Solution vector at current iteration
193 * @param F The right hand side vector
194 * @param ctx Pointer to context i.e. SnesCtx
195 * @return Error code
196 */
197PetscErrorCode SnesLoadTangent(SNES snes, Vec u, Vec F, void *ctx);
198
199/**
200 * \brief Set assembly type at the end of SnesMat
201 *
202 * \note Note that tangent matrix need have to have final assembly, you would
203 * use flush assembly in special case that you call SnesMat form other function
204 * set to SNESSetJacobian
205 *
206 * @param snes
207 * @param type type of assembly, either MAT_FLUSH_ASSEMBLY or
208 * MAT_FINAL_ASSEMBLY
209 * @return error code
210 */
211MoFEMErrorCode SnesMoFEMSetAssemblyType(SNES snes, MatAssemblyType type);
212
213/**
214 * \brief Set behavior if finite element in sequence does not exist
215 * @param snes
216 * @param bh If set to MF_EXIST check if element exist, default MF_EXIST.
217 * Otherwise set MF_ZERO
218 * @return error code
219 */
221
222/**
223 * @brief Sens monitor printing residual field by field
224 *
225 */
226MoFEMErrorCode MoFEMSNESMonitorFields(SNES snes, PetscInt its, PetscReal fgnorm,
227 SnesCtx *snes_ctx);
228
229/**
230 * @brief Sens monitor printing residual field by field
231 *
232 */
233MoFEMErrorCode MoFEMSNESMonitorEnergy(SNES snes, PetscInt its, PetscReal fgnorm,
234 SnesCtx *snes_ctx);
235
236} // namespace MoFEM
237
238#endif // __SNESCTX_HPP__
MoFEMTypes
Those types control how functions respond on arguments, f.e. error handling.
@ F
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
implementation of Data Operators for Forces and Sources
Definition Common.hpp:10
PetscErrorCode SnesMat(SNES snes, Vec x, Mat A, Mat B, void *ctx)
This is MoFEM implementation for the left hand side (tangent matrix) evaluation in SNES solver.
Definition SnesCtx.cpp:483
MoFEMErrorCode MoFEMSNESMonitorFields(SNES snes, PetscInt its, PetscReal fgnorm, SnesCtx *ctx)
Sens monitor printing residual field by field.
Definition SnesCtx.cpp:592
PetscErrorCode SnesRhs(SNES snes, Vec x, Vec f, void *ctx)
This is MoFEM implementation for the right hand side (residual vector) evaluation in SNES solver.
Definition SnesCtx.cpp:223
MoFEMErrorCode SnesMoFEMSetAssemblyType(SNES snes, MatAssemblyType type)
Set assembly type at the end of SnesMat.
Definition SnesCtx.cpp:564
MoFEMErrorCode SnesMoFEMSetBehavior(SNES snes, MoFEMTypes bh)
Set behavior if finite element in sequence does not exist.
Definition SnesCtx.cpp:578
std::deque< BasicMethodPtr > BasicMethodsSequence
Definition AuxPETSc.hpp:57
PetscErrorCode SnesLoadTangent(SNES snes, Vec x, Vec f, void *ctx)
This function calls finite element pipeline to compute tangent of of load vector in arc length contro...
Definition SnesCtx.cpp:353
std::deque< PairNameFEMethodPtr > FEMethodsSequence
Definition AuxPETSc.hpp:56
MoFEMErrorCode MoFEMSNESMonitorEnergy(SNES snes, PetscInt its, PetscReal fgnorm, SnesCtx *ctx)
Sens monitor printing residual field by field.
Definition SnesCtx.cpp:648
Deprecated interface functions.
Interface for nonlinear (SNES) solver.
Definition SnesCtx.hpp:15
MoFEM::BasicMethodsSequence BasicMethodsSequence
Definition SnesCtx.hpp:19
friend PetscErrorCode SnesLoadTangent(SNES snes, Vec u, Vec F, void *ctx)
This function calls finite element pipeline to compute tangent of of load vector in arc length contro...
Definition SnesCtx.cpp:353
FEMethodsSequence & getSetOperators()
Definition SnesCtx.cpp:131
friend MoFEMErrorCode SnesMoFEMSetBehavior(SNES snes, MoFEMTypes bh)
Set behavior if finite element in sequence does not exist.
Definition SnesCtx.cpp:578
friend MoFEMErrorCode MoFEMSNESMonitorEnergy(SNES snes, PetscInt its, PetscReal fgnorm, SnesCtx *snes_ctx)
Sens monitor printing residual field by field.
Definition SnesCtx.cpp:648
SnesCtx(Interface &m_field, std::string problem_name)
Definition SnesCtx.cpp:125
BasicMethodsSequence & getPreProcLoadTangent()
Get the BasicMethod sequence for preprocessing of FunctionFn.
Definition SnesCtx.cpp:159
BasicMethodsSequence & getPreProcSetOperators()
Definition SnesCtx.cpp:147
friend MoFEMErrorCode SNESMoFEMSetAssmblyType(SNES snes, MatAssemblyType type)
MoFEM::FEMethodsSequence FEMethodsSequence
Definition SnesCtx.hpp:18
BasicMethodsSequence & getPostProcLoadTangent()
Get the BasicMethod sequence for postprocessing of FunctionFn.
Definition SnesCtx.cpp:163
MoFEMErrorCode copyLoops(const SnesCtx &snes_ctx)
Copy sequences from other SNES context.
Definition SnesCtx.cpp:167
boost::function< MoFEMErrorCode( SNES snes, Vec x, Vec F, boost::shared_ptr< CacheTuple > cache_ptr, void *ctx)> HookFunction
Definition SnesCtx.hpp:20
friend PetscErrorCode SnesRhs(SNES snes, Vec x, Vec f, void *ctx)
This is MoFEM implementation for the right hand side (residual vector) evaluation in SNES solver.
Definition SnesCtx.cpp:223
friend PetscErrorCode SnesMat(SNES snes, Vec x, Mat A, Mat B, void *ctx)
This is MoFEM implementation for the left hand side (tangent matrix) evaluation in SNES solver.
Definition SnesCtx.cpp:483
BasicMethodsSequence & getPostProcComputeRhs()
Definition SnesCtx.cpp:143
HookFunction & getLoadTangentHook()
Get the Load Tangent Hook function.
Definition SnesCtx.cpp:173
BasicMethodsSequence & getPostProcSetOperators()
Definition SnesCtx.cpp:151
boost::movelib::unique_ptr< SnesCtxImpl > snesCtxImpl
Definition SnesCtx.hpp:150
BasicMethodsSequence & getPreProcComputeRhs()
Definition SnesCtx.cpp:139
FEMethodsSequence & getLoadTangent()
Get the finite element pipeline for FunctionFn, that calculate tangent of function load for arc lengt...
Definition SnesCtx.cpp:155
FEMethodsSequence & getComputeRhs()
Definition SnesCtx.cpp:135
HookFunction & getRhsHook()
Get the Right Hand Side Hook function.
Definition SnesCtx.cpp:177
virtual ~SnesCtx()
friend MoFEMErrorCode MoFEMSNESMonitorFields(SNES snes, PetscInt its, PetscReal fgnorm, SnesCtx *snes_ctx)
Sens monitor printing residual field by field.
Definition SnesCtx.cpp:592
MoFEMErrorCode clearLoops()
Clear loops.
Definition SnesCtx.cpp:171
friend MoFEMErrorCode SnesMoFEMSetAssemblyType(SNES snes, MatAssemblyType type)
Set assembly type at the end of SnesMat.
Definition SnesCtx.cpp:564