v0.9.1
LoopMethods.cpp
Go to the documentation of this file.
1 /** \file LoopMethods.cpp
2 \brief methods for managing loops
3 */
4 
5 /* MoFEM is free software: you can redistribute it and/or modify it under
6  * the terms of the GNU Lesser General Public License as published by the
7  * Free Software Foundation, either version 3 of the License, or (at your
8  * option) any later version.
9  *
10  * MoFEM is distributed in the hope that it will be useful, but WITHOUT
11  * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
12  * FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
13  * License for more details.
14  *
15  * You should have received a copy of the GNU Lesser General Public
16  * License along with MoFEM. If not, see <http://www.gnu.org/licenses/>
17  */
18 
19 namespace MoFEM {
20 
29 
30 // PetscData
32  UnknownInterface **iface) const {
34  if (uuid == IDD_MOFEMPetscDataMethod) {
35  *iface = const_cast<PetscData *>(this);
37  }
38  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "unknown interface");
39  }
40 
42  : f(PETSC_NULL), A(PETSC_NULL), B(PETSC_NULL), x(PETSC_NULL),
43  x_t(PETSC_NULL), x_tt(PETSC_NULL) {}
44 
45 // KSP
47  UnknownInterface **iface) const {
49  if (uuid == IDD_MOFEMKspMethod) {
50  *iface = const_cast<KspMethod *>(this);
52  }
53  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "unknown interface");
54 };
55 
57  : ksp_ctx(CTX_KSPNONE), ksp(PETSC_NULL), ksp_f(PetscData::f),
58  ksp_A(PetscData::A), ksp_B(PetscData::B) {}
59 
62  this->ksp_ctx = ksp.ksp_ctx;
63  this->ksp = ksp.ksp;
64  this->ksp_f = ksp.ksp_f;
65  this->ksp_A = ksp.ksp_A;
66  this->ksp_B = ksp.ksp_B;
68 }
69 
70 // SNES
72  UnknownInterface **iface) const {
73  if (uuid == IDD_MOFEMSnesMethod) {
74  *iface = const_cast<SnesMethod *>(this);
76  }
77  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "unknown interface");
78  }
79 
81  : snes_ctx(CTX_SNESNONE), snes_x(PetscData::x), snes_f(PetscData::f),
82  snes_A(PetscData::A), snes_B(PetscData::B) {}
83 
86  this->snes_ctx = snes.snes_ctx;
87  this->snes = snes.snes;
88  this->snes_x = snes.snes_x;
89  this->snes_f = snes.snes_f;
90  this->snes_A = snes.snes_A;
91  this->snes_B = snes.snes_B;
93 }
94 
95 // TS
97  UnknownInterface **iface) const {
98  if (uuid == IDD_MOFEMTsMethod) {
99  *iface = const_cast<TSMethod *>(this);
101  }
102  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "unknown interface");
103  }
104 
106  : ts_ctx(CTX_TSNONE), ts_step(-1), ts_a(0), ts_t(0), ts_u(PetscData::x),
107  ts_u_t(PetscData::x_t), ts_u_tt(PetscData::x_tt), ts_F(PetscData::f),
108  ts_A(PetscData::A), ts_B(PetscData::B) {}
109 
112  this->ts_ctx = ts.ts_ctx;
113  this->ts = ts.ts;
114  this->ts_u = ts.ts_u;
115  this->ts_u_t = ts.ts_u_t;
116  this->ts_u_tt = ts.ts_u_tt;
117  this->ts_F = ts.ts_F;
118  this->ts_A = ts.ts_A;
119  this->ts_B = ts.ts_B;
120  this->ts_step = ts.ts_step;
121  this->ts_a = ts.ts_a;
122  this->ts_t = ts.ts_t;
124 }
125 
126 // BasicMethod
128  : PetscData(), KspMethod(), SnesMethod(), TSMethod(), nInTheLoop(0),
129  loopSize(0), rAnk(-1), sIze(-1), refinedEntitiesPtr(nullptr),
130  refinedFiniteElementsPtr(nullptr), problemPtr(nullptr),
131  fieldsPtr(nullptr), entitiesPtr(nullptr), dofsPtr(nullptr),
132  finiteElementsPtr(nullptr), finiteElementsEntitiesPtr(nullptr),
133  adjacenciesPtr(nullptr) {}
134 
137 
138  this->nInTheLoop = basic.nInTheLoop;
139  this->loopSize = basic.loopSize;
140  this->rAnk = basic.rAnk;
141  this->sIze = basic.sIze;
144  this->problemPtr = basic.problemPtr;
145  this->fieldsPtr = basic.fieldsPtr;
146  this->entitiesPtr = basic.entitiesPtr;
147  this->dofsPtr = basic.dofsPtr;
148  this->finiteElementsPtr = basic.finiteElementsPtr;
150  this->adjacenciesPtr = basic.adjacenciesPtr;
151 
153 }
154 
157  if (preProcessHook) {
158  ierr = preProcessHook();
159  CHKERRG(ierr);
160  } else {
161  SETERRQ(PETSC_COMM_SELF, MOFEM_OPERATION_UNSUCCESSFUL,
162  "should be implemented by user in derived class (preProcess)");
163  }
165 }
168  if (postProcessHook) {
169  ierr = postProcessHook();
170  CHKERRG(ierr);
171  } else {
172  SETERRQ(PETSC_COMM_SELF, MOFEM_OPERATION_UNSUCCESSFUL,
173  "should be implemented by user in derived class (postProcess)");
174  }
176 }
179  if (operatorHook) {
180  ierr = operatorHook();
181  CHKERRG(ierr);
182  } else {
183  SETERRQ(PETSC_COMM_SELF, MOFEM_OPERATION_UNSUCCESSFUL,
184  "should be implemented by user in derived class (operator)");
185  }
187 }
188 
189 // FEMethod
191 
192 // Entity method
194 
195 // DofMethod
197 
198 } // namespace MoFEM
MoFEMErrorCode copyBasicMethod(const BasicMethod &basic)
Copy data from other base method to this base method.
const DofEntity_multiIndex * dofsPtr
raw pointer container of dofs
const FieldEntityEntFiniteElementAdjacencyMap_multiIndex * adjacenciesPtr
static constexpr Switches CtxSetX_T
Definition: LoopMethods.hpp:66
MoFEM interface unique ID.
KSP ksp
KSP solver.
const Field_multiIndex * fieldsPtr
raw pointer to fields container
MoFEMErrorCode query_interface(const MOFEMuuid &uuid, UnknownInterface **iface) const
Definition: LoopMethods.cpp:46
static const MOFEMuuid IDD_MOFEMPetscDataMethod
Definition: LoopMethods.hpp:24
const FiniteElement_multiIndex * finiteElementsPtr
raw pointer to container finite elements
static constexpr Switches CtxSetB
Definition: LoopMethods.hpp:64
int loopSize
local number oe methods to process
std::bitset< 8 > Switches
Definition: LoopMethods.hpp:59
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:506
const FieldEntity_multiIndex * entitiesPtr
raw pointer to container of field entities
Mat & ts_B
Preconditioner for ts_A.
base class for all interface classes
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:549
const RefElement_multiIndex * refinedFiniteElementsPtr
container of mofem finite element entities
data structure for snes (nonlinear solver) contextStructure stores context data which are set in func...
int rAnk
processor rank
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:513
virtual MoFEMErrorCode preProcess()
function is run at the beginning of loop
KSPContext ksp_ctx
Context.
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:21
boost::function< MoFEMErrorCode()> postProcessHook
Hook function for post-processing.
const Problem * problemPtr
raw pointer to problem
PetscReal ts_t
time
MoFEMErrorCode query_interface(const MOFEMuuid &uuid, UnknownInterface **iface) const
Definition: LoopMethods.cpp:31
MoFEMErrorCode copyKsp(const KspMethod &ksp)
copy data form another method
Definition: LoopMethods.cpp:60
Vec & snes_f
residual
virtual MoFEMErrorCode operator()()
function is run for every finite element
virtual MoFEMErrorCode postProcess()
function is run at the end of loop
const RefEntity_multiIndex * refinedEntitiesPtr
container of mofem dof entities
MoFEMErrorCode copyTs(const TSMethod &ts)
Copy TS solver data.
static constexpr Switches CtxSetTime
Definition: LoopMethods.hpp:68
Vec & ts_u_t
time derivative of state vector
boost::function< MoFEMErrorCode()> operatorHook
Hook function for operator.
PetscReal ts_a
shift for U_tt (see PETSc Time Solver)
Mat & snes_B
preconditioner of jacobian matrix
int sIze
number of processors in communicator
PetscInt ts_step
time step
TS ts
time solver
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:66
data structure for ksp (linear solver) contextStruture stores context data which are set in functions...
Definition: LoopMethods.hpp:88
Vec & snes_x
state vector
static constexpr Switches CtxSetF
Definition: LoopMethods.hpp:62
data structure for TS (time stepping) contextStructure stores context data which are set in functions...
static constexpr Switches CtxSetX_TT
Definition: LoopMethods.hpp:67
int nInTheLoop
number currently of processed method
static constexpr Switches CtxSetA
Definition: LoopMethods.hpp:63
Vec & ts_F
residual vector
SNES snes
snes solver
MoFEMErrorCode copySnes(const SnesMethod &snes)
Copy snes data.
Definition: LoopMethods.cpp:84
MoFEMErrorCode query_interface(const MOFEMuuid &uuid, UnknownInterface **iface) const
Definition: LoopMethods.cpp:71
Vec & ts_u
state vector
SNESContext snes_ctx
const EntFiniteElement_multiIndex * finiteElementsEntitiesPtr
Mat & snes_A
jacobian matrix
MoFEMErrorCode query_interface(const MOFEMuuid &uuid, UnknownInterface **iface) const
Definition: LoopMethods.cpp:96
static const MOFEMuuid IDD_MOFEMTsMethod
Definition: LoopMethods.hpp:30
TSContext ts_ctx
static const MOFEMuuid IDD_MOFEMKspMethod
Definition: LoopMethods.hpp:26
Vec & ts_u_tt
second time derivative of state vector
boost::function< MoFEMErrorCode()> preProcessHook
Hook function for pre-processing.
static constexpr Switches CtxSetX
Definition: LoopMethods.hpp:65
Data structure to exchange data between mofem and User Loop Methods.It allows to exchange data betwee...
static constexpr Switches CtxSetNone
Definition: LoopMethods.hpp:61
static const MOFEMuuid IDD_MOFEMSnesMethod
Definition: LoopMethods.hpp:28