v0.14.0
Loading...
Searching...
No Matches
SnesCtx.hpp
Go to the documentation of this file.
1/** \file SnesCtx.hpp
2 * \brief Context for PETSc SNES, i.e. nonlinear solver
3 */
4
5#ifndef __SNESCTX_HPP__
6#define __SNESCTX_HPP__
7
8namespace MoFEM {
9
10/** \brief Interface for nonlinear (SNES) solver
11 * \ingroup mofem_petsc_solvers
12 */
13struct SnesCtx {
14
15 MoFEM::Interface &mField; ///< database Interface
16 moab::Interface &moab; ///< moab Interface
17
18 std::string problemName; ///< problem name
20 bH; ///< If set to MF_EXIST check if element exist, default MF_EXIST
21 bool zeroPreCondMatrixB; ///< If true zero matrix, otherwise user need to do
22 ///< it, default true
23 MatAssemblyType typeOfAssembly; ///< type of assembly at the end
24 bool vErify; ///< If true verify vector
25
29
30 FEMethodsSequence loops_to_do_Mat; ///< Sequence of finite elements instances
31 ///< assembling tangent matrix
32 FEMethodsSequence loops_to_do_Rhs; ///< Sequence of finite elements instances
33 ///< assembling residual vector
34 BasicMethodsSequence preProcess_Mat; ///< Sequence of methods run before
35 ///< tangent matrix is assembled
36 BasicMethodsSequence postProcess_Mat; ///< Sequence of methods run after
37 ///< tangent matrix is assembled
39 preProcess_Rhs; ///< Sequence of methods run before residual is assembled
41 postProcess_Rhs; ///< Sequence of methods run after residual is assembled
42
43 SnesCtx(Interface &m_field, const std::string &problem_name)
44 : mField(m_field), moab(m_field.get_moab()), problemName(problem_name),
46 typeOfAssembly(MAT_FINAL_ASSEMBLY), vErify(false) {
47 PetscLogEventRegister("LoopSNESRhs", 0, &MOFEM_EVENT_SnesRhs);
48 PetscLogEventRegister("LoopSNESMat", 0, &MOFEM_EVENT_SnesMat);
49 if (!LogManager::checkIfChannelExist("SNES_WORLD")) {
50 auto core_log = logging::core::get();
51 core_log->add_sink(
53 LogManager::setLog("SNES_WORLD");
54 MOFEM_LOG_TAG("SNES_WORLD", "SNES");
55 }
56 }
57
58 virtual ~SnesCtx() = default;
59
60 /**
61 * @return return reference to vector with FEMethod to calculate tangent
62 * matrix
63 */
65
66 /**
67 * @return return vector to vector with FEMethod to calculate residual
68 */
70
71 /**
72 * The sequence of BasicMethod is executed before residual is calculated. It
73 * can be used to setup data structures, e.g. zero global variable which is
74 * integrated in domain, e.g. for calculation of strain energy.
75 *
76 * @return reference to BasicMethod for preprocessing
77 */
79
80 /**
81 * The sequence of BasicMethod is executed after residual is calculated. It
82 * can be used to setup data structures, e.g. aggregate data from processors
83 * or to apply essential boundary conditions.
84 *
85 * @return reference to BasicMethod for postprocessing
86 */
88
89 /**
90 * @return reference to BasicMethod for preprocessing
91 */
93
94 /**
95 * The sequence of BasicMethod is executed after tangent matrix is calculated.
96 * It can be used to setup data structures, e.g. aggregate data from
97 * processors or to apply essential boundary conditions.
98 *
99 * @return reference to BasicMethod for postprocessing
100 */
102
103 /**
104 * \brief Copy sequences from other SNES contex
105 * @param snes_ctx SNES contex from which Sequence is copied from
106 * @return error code
107 */
108 MoFEMErrorCode copyLoops(const SnesCtx &snes_ctx);
109
110 /**
111 * @brief Clear loops
112 *
113 * @return MoFEMErrorCode
114 */
116
117 friend PetscErrorCode SnesRhs(SNES snes, Vec x, Vec f, void *ctx);
118 friend PetscErrorCode SnesMat(SNES snes, Vec x, Mat A, Mat B, void *ctx);
119
121 MatAssemblyType type);
122 friend MoFEMErrorCode SnesMoFEMSetBehavior(SNES snes, MoFEMTypes bh);
123
124 /** @deprecated use getSetOperators */
126 return getSetOperators();
127 }
128
129 /** @deprecated use getComputeRhs */
131 return getComputeRhs();
132 }
133
134 /** @deprecated use getPreProcComputeRhs */
136 return getPreProcComputeRhs();
137 }
138
139 /** @deprecated use getPostProcComputeRhs */
141 return getPostProcComputeRhs();
142 }
143
144 /** @deprecated use getPreProcSetOperators */
146 return getPreProcSetOperators();
147 }
148
149 /** @deprecated use getPostProcSetOperators */
152 }
153
154private:
155 boost::movelib::unique_ptr<bool> vecAssembleSwitch;
156 boost::movelib::unique_ptr<bool> matAssembleSwitch;
157 PetscLogEvent MOFEM_EVENT_SnesRhs; ///< Log events to assemble residual
158 PetscLogEvent MOFEM_EVENT_SnesMat; ///< Log events to assemble tangent matrix
159};
160
161/**
162 * \brief This is MoFEM implementation for the right hand side (residual vector)
163 * evaluation in SNES solver
164 *
165 * For more information pleas look to PETSc manual, i.e. SNESSetFunction
166 * <http://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/SNES/SNESSetFunction.html>
167 *
168 * @param snes SNES solver
169 * @param x Solution vector at current iteration
170 * @param f The right hand side vector
171 * @param ctx Pointer to context i.e. SnesCtx
172 * @return Error code
173 */
174PetscErrorCode SnesRhs(SNES snes, Vec x, Vec f, void *ctx);
175
176/**
177 * \brief This is MoFEM implementation for the left hand side (tangent matrix)
178 * evaluation in SNES solver
179 *
180 * For more information pleas look to PETSc manual, i.e. SNESSetJacobian
181 * <http://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/SNES/SNESSetJacobian.html#SNESSetJacobian>
182 *
183 * @param snes SNES solver
184 * @param x Solution vector at current iteration
185 * @param A Tangent matrix
186 * @param B Preconditioner tangent matrix
187 * @param ctx Pointer to context i.e. SnesCtx
188 * @return Error code
189 */
190PetscErrorCode SnesMat(SNES snes, Vec x, Mat A, Mat B, void *ctx);
191
192/**
193 * \brief Set assembly type at the end of SnesMat
194 *
195 * \note Note that tangent matrix need have to have final assembly, you would
196 * use flush assembly in special case that you call SnesMat form other function
197 * set to SNESSetJacobian
198 *
199 * @param snes
200 * @param type type of assembly, either MAT_FLUSH_ASSEMBLY or
201 * MAT_FINAL_ASSEMBLY
202 * @return error code
203 */
204MoFEMErrorCode SnesMoFEMSetAssemblyType(SNES snes, MatAssemblyType type);
205
206/**
207 * \brief Set behavior if finite element in sequence does not exist
208 * @param snes
209 * @param bh If set to MF_EXIST check if element exist, default MF_EXIST.
210 * Otherwise set MF_ZERO
211 * @return error code
212 */
214
215/**
216 * @brief Sens monitor printing residual field by field
217 *
218 */
219MoFEMErrorCode MoFEMSNESMonitorFields(SNES snes, PetscInt its, PetscReal fgnorm,
220 SnesCtx *snes_ctx);
221
222} // namespace MoFEM
223
224#endif // __SNESCTX_HPP__
MoFEMTypes
Those types control how functions respond on arguments, f.e. error handling.
Definition: definitions.h:97
@ MF_EXIST
Definition: definitions.h:100
#define DEPRECATED
Definition: definitions.h:17
static LoggerType & setLog(const std::string channel)
Set ans resset chanel logger.
Definition: LogManager.cpp:389
#define MOFEM_LOG_TAG(channel, tag)
Tag channel.
Definition: LogManager.hpp:339
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:56
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:136
MoFEMErrorCode MoFEMSNESMonitorFields(SNES snes, PetscInt its, PetscReal fgnorm, SnesCtx *snes_ctx)
Sens monitor printing residual field by field.
Definition: SnesCtx.cpp:226
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:27
MoFEMErrorCode SnesMoFEMSetAssemblyType(SNES snes, MatAssemblyType type)
Set assembly type at the end of SnesMat.
Definition: SnesCtx.cpp:209
MoFEMErrorCode SnesMoFEMSetBehavior(SNES snes, MoFEMTypes bh)
Set behavior if finite element in sequence does not exist.
Definition: SnesCtx.cpp:218
std::deque< BasicMethodPtr > BasicMethodsSequence
Definition: AuxPETSc.hpp:57
std::deque< PairNameFEMethodPtr > FEMethodsSequence
Definition: AuxPETSc.hpp:56
constexpr AssemblyType A
Deprecated interface functions.
static boost::shared_ptr< SinkType > createSink(boost::shared_ptr< std::ostream > stream_ptr, std::string comm_filter)
Create a sink object.
Definition: LogManager.cpp:298
static boost::shared_ptr< std::ostream > getStrmWorld()
Get the strm world object.
Definition: LogManager.cpp:344
static bool checkIfChannelExist(const std::string channel)
Check if channel exist.
Definition: LogManager.cpp:399
Interface for nonlinear (SNES) solver.
Definition: SnesCtx.hpp:13
DEPRECATED FEMethodsSequence & get_loops_to_do_Mat()
Definition: SnesCtx.hpp:125
PetscLogEvent MOFEM_EVENT_SnesMat
Log events to assemble tangent matrix.
Definition: SnesCtx.hpp:158
BasicMethodsSequence preProcess_Rhs
Sequence of methods run before residual is assembled.
Definition: SnesCtx.hpp:39
BasicMethodsSequence postProcess_Mat
Definition: SnesCtx.hpp:36
FEMethodsSequence & getSetOperators()
Definition: SnesCtx.hpp:64
friend MoFEMErrorCode SnesMoFEMSetBehavior(SNES snes, MoFEMTypes bh)
Set behavior if finite element in sequence does not exist.
Definition: SnesCtx.cpp:218
SnesCtx(Interface &m_field, const std::string &problem_name)
Definition: SnesCtx.hpp:43
MoFEMTypes bH
If set to MF_EXIST check if element exist, default MF_EXIST.
Definition: SnesCtx.hpp:20
BasicMethodsSequence postProcess_Rhs
Sequence of methods run after residual is assembled.
Definition: SnesCtx.hpp:41
moab::Interface & moab
moab Interface
Definition: SnesCtx.hpp:16
MoFEMErrorCode copyLoops(const SnesCtx &snes_ctx)
Copy sequences from other SNES contex.
Definition: SnesCtx.cpp:5
DEPRECATED BasicMethodsSequence & get_preProcess_to_do_Mat()
Definition: SnesCtx.hpp:145
MatAssemblyType typeOfAssembly
type of assembly at the end
Definition: SnesCtx.hpp:23
BasicMethodsSequence & getPreProcSetOperators()
Definition: SnesCtx.hpp:92
bool zeroPreCondMatrixB
Definition: SnesCtx.hpp:21
friend MoFEMErrorCode SNESMoFEMSetAssmblyType(SNES snes, MatAssemblyType type)
MoFEMErrorCode clearLoops()
Clear loops.
Definition: SnesCtx.cpp:16
DEPRECATED BasicMethodsSequence & get_preProcess_to_do_Rhs()
Definition: SnesCtx.hpp:135
BasicMethodsSequence preProcess_Mat
Definition: SnesCtx.hpp:34
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:27
virtual ~SnesCtx()=default
MoFEM::PairNameFEMethodPtr PairNameFEMethodPtr
Definition: SnesCtx.hpp:26
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:136
bool vErify
If true verify vector.
Definition: SnesCtx.hpp:24
MoFEM::FEMethodsSequence FEMethodsSequence
Definition: SnesCtx.hpp:27
BasicMethodsSequence & getPostProcComputeRhs()
Definition: SnesCtx.hpp:87
FEMethodsSequence loops_to_do_Mat
Definition: SnesCtx.hpp:30
BasicMethodsSequence & getPostProcSetOperators()
Definition: SnesCtx.hpp:101
BasicMethodsSequence & getPreProcComputeRhs()
Definition: SnesCtx.hpp:78
std::string problemName
problem name
Definition: SnesCtx.hpp:18
FEMethodsSequence loops_to_do_Rhs
Definition: SnesCtx.hpp:32
FEMethodsSequence & getComputeRhs()
Definition: SnesCtx.hpp:69
DEPRECATED BasicMethodsSequence & get_postProcess_to_do_Rhs()
Definition: SnesCtx.hpp:140
boost::movelib::unique_ptr< bool > matAssembleSwitch
Definition: SnesCtx.hpp:156
DEPRECATED FEMethodsSequence & get_loops_to_do_Rhs()
Definition: SnesCtx.hpp:130
DEPRECATED BasicMethodsSequence & get_postProcess_to_do_Mat()
Definition: SnesCtx.hpp:150
MoFEM::Interface & mField
database Interface
Definition: SnesCtx.hpp:15
boost::movelib::unique_ptr< bool > vecAssembleSwitch
Definition: SnesCtx.hpp:155
MoFEM::BasicMethodsSequence BasicMethodsSequence
Definition: SnesCtx.hpp:28
PetscLogEvent MOFEM_EVENT_SnesRhs
Log events to assemble residual.
Definition: SnesCtx.hpp:157