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