v0.13.0
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 /* 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 __SNESCTX_HPP__
20 #define __SNESCTX_HPP__
21 
22 namespace MoFEM {
23 
24 /** \brief Interface for nonlinear (SNES) solver
25  * \ingroup mofem_petsc_solvers
26  */
27 struct SnesCtx {
28 
29  MoFEM::Interface &mField; ///< database Interface
30  moab::Interface &moab; ///< moab Interface
31 
32  std::string problemName; ///< problem name
34  bH; ///< If set to MF_EXIST check if element exist, default MF_EXIST
35  bool zeroPreCondMatrixB; ///< If true zero matrix, otherwise user need to do
36  ///< it, default true
37  MatAssemblyType typeOfAssembly; ///< type of assembly at the end
38  bool vErify; ///< If true verify vector
39 
43 
44  FEMethodsSequence loops_to_do_Mat; ///< Sequence of finite elements instances
45  ///< assembling tangent matrix
46  FEMethodsSequence loops_to_do_Rhs; ///< Sequence of finite elements instances
47  ///< assembling residual vector
48  BasicMethodsSequence preProcess_Mat; ///< Sequence of methods run before
49  ///< tangent matrix is assembled
50  BasicMethodsSequence postProcess_Mat; ///< Sequence of methods run after
51  ///< tangent matrix is assembled
53  preProcess_Rhs; ///< Sequence of methods run before residual is assembled
55  postProcess_Rhs; ///< Sequence of methods run after residual is assembled
56 
57 
58  SnesCtx(Interface &m_field, const std::string &problem_name)
59  : mField(m_field), moab(m_field.get_moab()), problemName(problem_name),
61  typeOfAssembly(MAT_FINAL_ASSEMBLY), vErify(false) {
62  PetscLogEventRegister("LoopSNESRhs", 0, &MOFEM_EVENT_SnesRhs);
63  PetscLogEventRegister("LoopSNESMat", 0, &MOFEM_EVENT_SnesMat);
64  }
65 
66  virtual ~SnesCtx() = default;
67 
68  /**
69  * @return return reference to vector with FEMethod to calculate tangent
70  * matrix
71  */
73 
74  /**
75  * @return return vector to vector with FEMethod to calculate residual
76  */
78 
79  /**
80  * The sequence of BasicMethod is executed before residual is calculated. It
81  * can be used to setup data structures, e.g. zero global variable which is
82  * integrated in domain, e.g. for calculation of strain energy.
83  *
84  * @return reference to BasicMethod for preprocessing
85  */
87 
88  /**
89  * The sequence of BasicMethod is executed after residual is calculated. It
90  * can be used to setup data structures, e.g. aggregate data from processors
91  * or to apply essential boundary conditions.
92  *
93  * @return reference to BasicMethod for postprocessing
94  */
96 
97  /**
98  * @return reference to BasicMethod for preprocessing
99  */
101 
102  /**
103  * The sequence of BasicMethod is executed after tangent matrix is calculated.
104  * It can be used to setup data structures, e.g. aggregate data from
105  * processors or to apply essential boundary conditions.
106  *
107  * @return reference to BasicMethod for postprocessing
108  */
110 
111  /**
112  * \brief Copy sequences from other SNES contex
113  * @param snes_ctx SNES contex from which Sequence is copied from
114  * @return error code
115  */
116  MoFEMErrorCode copyLoops(const SnesCtx &snes_ctx);
117 
118  /**
119  * @brief Clear loops
120  *
121  * @return MoFEMErrorCode
122  */
124 
125  friend PetscErrorCode SnesRhs(SNES snes, Vec x, Vec f, void *ctx);
126  friend PetscErrorCode SnesMat(SNES snes, Vec x, Mat A, Mat B, void *ctx);
127 
129  MatAssemblyType type);
130  friend MoFEMErrorCode SnesMoFEMSetBehavior(SNES snes, MoFEMTypes bh);
131 
132 private:
133 
134  boost::movelib::unique_ptr<bool> vecAssembleSwitch;
135  boost::movelib::unique_ptr<bool> matAssembleSwitch;
136  PetscLogEvent MOFEM_EVENT_SnesRhs; ///< Log events to assemble residual
137  PetscLogEvent MOFEM_EVENT_SnesMat; ///< Log events to assemble tangent matrix
138 
139 };
140 
141 /**
142  * \brief This is MoFEM implementation for the right hand side (residual vector)
143  * evaluation in SNES solver
144  *
145  * For more information pleas look to PETSc manual, i.e. SNESSetFunction
146  * <http://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/SNES/SNESSetFunction.html>
147  *
148  * @param snes SNES solver
149  * @param x Solution vector at current iteration
150  * @param f The right hand side vector
151  * @param ctx Pointer to context i.e. SnesCtx
152  * @return Error code
153  */
154 PetscErrorCode SnesRhs(SNES snes, Vec x, Vec f, void *ctx);
155 
156 /**
157  * \brief This is MoFEM implementation for the left hand side (tangent matrix)
158  * evaluation in SNES solver
159  *
160  * For more information pleas look to PETSc manual, i.e. SNESSetJacobian
161  * <http://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/SNES/SNESSetJacobian.html#SNESSetJacobian>
162  *
163  * @param snes SNES solver
164  * @param x Solution vector at current iteration
165  * @param A Tangent matrix
166  * @param B Preconditioner tangent matrix
167  * @param ctx Pointer to context i.e. SnesCtx
168  * @return Error code
169  */
170 PetscErrorCode SnesMat(SNES snes, Vec x, Mat A, Mat B, void *ctx);
171 
172 /**
173  * \brief Set assembly type at the end of SnesMat
174  *
175  * \note Note that tangent matrix need have to have final assembly, you would
176  * use flush assembly in special case that you call SnesMat form other function
177  * set to SNESSetJacobian
178  *
179  * @param snes
180  * @param type type of assembly, either MAT_FLUSH_ASSEMBLY or
181  * MAT_FINAL_ASSEMBLY
182  * @return error code
183  */
184 MoFEMErrorCode SnesMoFEMSetAssemblyType(SNES snes, MatAssemblyType type);
185 
186 /**
187  * \brief Set behavior if finite element in sequence does not exist
188  * @param snes
189  * @param bh If set to MF_EXIST check if element exist, default MF_EXIST.
190  * Otherwise set MF_ZERO
191  * @return error code
192  */
194 
195 } // namespace MoFEM
196 
197 #endif // __SNESCTX_HPP__
MoFEMTypes
Those types control how functions respond on arguments, f.e. error handling.
Definition: definitions.h:110
@ MF_EXIST
Definition: definitions.h:113
const FTensor::Tensor2< T, Dim, Dim > Vec
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:67
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:21
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:148
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:39
MoFEMErrorCode SnesMoFEMSetAssemblyType(SNES snes, MatAssemblyType type)
Set assembly type at the end of SnesMat.
Definition: SnesCtx.cpp:219
MoFEMErrorCode SnesMoFEMSetBehavior(SNES snes, MoFEMTypes bh)
Set behavior if finite element in sequence does not exist.
Definition: SnesCtx.cpp:228
std::deque< BasicMethodPtr > BasicMethodsSequence
Definition: AuxPETSc.hpp:65
std::deque< PairNameFEMethodPtr > FEMethodsSequence
Definition: AuxPETSc.hpp:64
DeprecatedCoreInterface Interface
Definition: Interface.hpp:1965
double A
Deprecated interface functions.
Interface for nonlinear (SNES) solver.
Definition: SnesCtx.hpp:27
BasicMethodsSequence & get_postProcess_to_do_Rhs()
Definition: SnesCtx.hpp:95
BasicMethodsSequence & get_preProcess_to_do_Mat()
Definition: SnesCtx.hpp:100
PetscLogEvent MOFEM_EVENT_SnesMat
Log events to assemble tangent matrix.
Definition: SnesCtx.hpp:137
BasicMethodsSequence preProcess_Rhs
Sequence of methods run before residual is assembled.
Definition: SnesCtx.hpp:53
BasicMethodsSequence postProcess_Mat
Definition: SnesCtx.hpp:50
friend MoFEMErrorCode SnesMoFEMSetBehavior(SNES snes, MoFEMTypes bh)
Set behavior if finite element in sequence does not exist.
Definition: SnesCtx.cpp:228
SnesCtx(Interface &m_field, const std::string &problem_name)
Definition: SnesCtx.hpp:58
MoFEMTypes bH
If set to MF_EXIST check if element exist, default MF_EXIST.
Definition: SnesCtx.hpp:34
BasicMethodsSequence postProcess_Rhs
Sequence of methods run after residual is assembled.
Definition: SnesCtx.hpp:55
moab::Interface & moab
moab Interface
Definition: SnesCtx.hpp:30
MoFEMErrorCode copyLoops(const SnesCtx &snes_ctx)
Copy sequences from other SNES contex.
Definition: SnesCtx.cpp:17
BasicMethodsSequence & get_preProcess_to_do_Rhs()
Definition: SnesCtx.hpp:86
MatAssemblyType typeOfAssembly
type of assembly at the end
Definition: SnesCtx.hpp:37
bool zeroPreCondMatrixB
Definition: SnesCtx.hpp:35
friend MoFEMErrorCode SNESMoFEMSetAssmblyType(SNES snes, MatAssemblyType type)
MoFEMErrorCode clearLoops()
Clear loops.
Definition: SnesCtx.cpp:28
BasicMethodsSequence preProcess_Mat
Definition: SnesCtx.hpp:48
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:39
virtual ~SnesCtx()=default
MoFEM::PairNameFEMethodPtr PairNameFEMethodPtr
Definition: SnesCtx.hpp:40
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:148
bool vErify
If true verify vector.
Definition: SnesCtx.hpp:38
MoFEM::FEMethodsSequence FEMethodsSequence
Definition: SnesCtx.hpp:41
FEMethodsSequence loops_to_do_Mat
Definition: SnesCtx.hpp:44
std::string problemName
problem name
Definition: SnesCtx.hpp:32
FEMethodsSequence loops_to_do_Rhs
Definition: SnesCtx.hpp:46
boost::movelib::unique_ptr< bool > matAssembleSwitch
Definition: SnesCtx.hpp:135
FEMethodsSequence & get_loops_to_do_Rhs()
Definition: SnesCtx.hpp:77
FEMethodsSequence & get_loops_to_do_Mat()
Definition: SnesCtx.hpp:72
BasicMethodsSequence & get_postProcess_to_do_Mat()
Definition: SnesCtx.hpp:109
MoFEM::Interface & mField
database Interface
Definition: SnesCtx.hpp:29
boost::movelib::unique_ptr< bool > vecAssembleSwitch
Definition: SnesCtx.hpp:134
MoFEM::BasicMethodsSequence BasicMethodsSequence
Definition: SnesCtx.hpp:42
PetscLogEvent MOFEM_EVENT_SnesRhs
Log events to assemble residual.
Definition: SnesCtx.hpp:136