v0.13.1
LoopMethods.hpp
Go to the documentation of this file.
1/** \file LoopMethods.hpp
2 * \brief MoFEM interface
3 *
4 * Data structures for making loops over finite elements and entities in the
5 * problem or MoFEM database.
6 *
7 */
8
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#ifndef __LOOPMETHODS_HPP__
20#define __LOOPMETHODS_HPP__
21
22namespace MoFEM {
23struct PetscData : public UnknownInterface {
24
25 MoFEMErrorCode query_interface(boost::typeindex::type_index type_index,
26 UnknownInterface **iface) const;
27
28 PetscData();
29
30 virtual ~PetscData() = default;
31
34 CTX_SET_F = 1 << 0,
35 CTX_SET_A = 1 << 1,
36 CTX_SET_B = 1 << 2,
37 CTX_SET_X = 1 << 3,
38 CTX_SET_X_T = 1 << 4,
39 CTX_SET_X_TT = 1 << 6,
40 CTX_SET_TIME = 1 << 7
41 };
42
43 using Switches = std::bitset<8>;
44
53
55
56 MoFEMErrorCode copyPetscData(const PetscData &petsc_data);
57
59 Mat A;
60 Mat B;
64};
65
66/**
67 * \brief data structure for ksp (linear solver) context
68 * \ingroup mofem_loops
69 *
70 * Struture stores context data which are set in functions run by PETSc SNES
71 * functions.
72 *
73 */
74struct KspMethod : virtual public PetscData {
75
76 MoFEMErrorCode query_interface(boost::typeindex::type_index type_index,
77 UnknownInterface **iface) const;
78
79 /**
80 * \brief pass information about context of KSP/DM for with finite element is
81 * computed
82 */
84
85 KspMethod();
86
87 virtual ~KspMethod() = default;
88
89 /**
90 * \brief copy data form another method
91 * @param ksp ksp method
92 * @return error code
93 */
95
96 KSPContext ksp_ctx; ///< Context
97 KSP ksp; ///< KSP solver
98
100 Mat &ksp_A;
101 Mat &ksp_B;
102};
103
104/**
105 * \brief data structure for snes (nonlinear solver) context
106 * \ingroup mofem_loops
107 *
108 * Structure stores context data which are set in functions run by PETSc SNES
109 * functions.
110 *
111 */
112struct SnesMethod : virtual protected PetscData {
113
114 MoFEMErrorCode query_interface(boost::typeindex::type_index type_index,
115 UnknownInterface **iface) const;
116
118
119 SnesMethod();
120
121 virtual ~SnesMethod() = default;
122
123 /**
124 * \brief Copy snes data
125 */
127
129
130 SNES snes; ///< snes solver
131 Vec &snes_x; ///< state vector
132 Vec &snes_f; ///< residual
133 Mat &snes_A; ///< jacobian matrix
134 Mat &snes_B; ///< preconditioner of jacobian matrix
135};
136
137/**
138 * \brief data structure for TS (time stepping) context
139 * \ingroup mofem_loops
140 *
141 * Structure stores context data which are set in functions run by PETSc Time
142 * Stepping functions.
143 */
144struct TSMethod : virtual protected PetscData {
145
146 MoFEMErrorCode query_interface(boost::typeindex::type_index type_index,
147 UnknownInterface **iface) const;
148
156 };
157
158 TSMethod();
159
160 virtual ~TSMethod() = default;
161
162 /// \brief Copy TS solver data
164
165 TS ts; ///< time solver
166
168
169 PetscInt ts_step; ///< time step
170 PetscReal ts_a; ///< shift for U_t (see PETSc Time Solver)
171 PetscReal ts_aa; ///< shift for U_tt shift for U_tt
172 PetscReal ts_t; ///< time
173
174 Vec &ts_u; ///< state vector
175 Vec &ts_u_t; ///< time derivative of state vector
176 Vec &ts_u_tt; ///< second time derivative of state vector
177 Vec &ts_F; ///< residual vector
178
179 Mat &ts_A; ///< Jacobian of G(U) = F(t,U,W+v*U,W'+a*U), equivalent to dF/dU +
180 ///< v*dF/dU_t + a*dF/dU_tt
181 Mat &ts_B; ///< Preconditioner for ts_A
182};
183
184/**
185 * \brief Data structure to exchange data between mofem and User Loop Methods.
186 * \ingroup mofem_loops
187 *
188 * It allows to exchange data between MoFEM and user functions. It stores
189 * information about multi-indices.
190 *
191 */
193
194 MoFEMErrorCode query_interface(boost::typeindex::type_index type_index,
195 UnknownInterface **iface) const {
197 *iface = const_cast<BasicMethod *>(this);
199 }
200
201 BasicMethod();
202 virtual ~BasicMethod() = default;
203
204 /**
205 * @brief number currently of processed method
206 */
208
209 /**
210 * @brief local number oe methods to process
211 */
213
214 /** \brief get number of evaluated element in the loop
215 */
216 inline int getNinTheLoop() const { return nInTheLoop; }
217
218 /** \brief get loop size
219 */
220 inline int getLoopSize() const { return loopSize; }
221
222 /**
223 * @brief Llo and hi processor rank of iterated entities
224 *
225 */
226 std::pair<int, int> loHiFERank;
227
228 /**
229 * @brief Get lo and hi processor rank of iterated entities
230 *
231 * @return raturn std::pair<int, int> loHiFERank
232 */
233 inline auto getLoHiFERank() const { return loHiFERank; }
234
235 /**
236 * @brief Get upper rank in loop for iterating elements
237 *
238 * @return loHiFERank.first
239 */
240 inline auto getLoFERank() const { return loHiFERank.first; }
241
242 /**
243 * @brief Get upper rank in loop for iterating elements
244 *
245 * @return loHiFERank.first
246 */
247 inline auto getHiFERank() const { return loHiFERank.second; }
248
249 int rAnk; ///< processor rank
250
251 int sIze; ///< number of processors in communicator
252
254 *refinedEntitiesPtr; ///< container of mofem dof entities
255
257 *refinedFiniteElementsPtr; ///< container of mofem finite element entities
258
259 const Problem *problemPtr; ///< raw pointer to problem
260
261 const Field_multiIndex *fieldsPtr; ///< raw pointer to fields container
262
264 *entitiesPtr; ///< raw pointer to container of field entities
265
266 const DofEntity_multiIndex *dofsPtr; ///< raw pointer container of dofs
267
269 *finiteElementsPtr; ///< raw pointer to container finite elements
270
272 *finiteElementsEntitiesPtr; ///< raw pointer to container finite elements
273 ///< entities
274
276 *adjacenciesPtr; ///< raw pointer to container to adjacencies between dofs
277 ///< and finite elements
278
279 inline unsigned int getFieldBitNumber(std::string field_name) const {
280 if (fieldsPtr) {
281 auto field_it = fieldsPtr->get<FieldName_mi_tag>().find(field_name);
282 if (field_it != fieldsPtr->get<FieldName_mi_tag>().end())
283 return (*field_it)->getBitNumber();
284 else
285 return BITFEID_SIZE;
286 } else {
287 THROW_MESSAGE("Pointer to fields multi-index is not set");
288 return BITFEID_SIZE;
289 }
290 }
291
292 /**
293 * @brief Copy data from other base method to this base method
294 *
295 * @param basic
296 * @return MoFEMErrorCode
297 */
299
300 /**
301 * @brief Hook function for pre-processing
302 */
303 boost::function<MoFEMErrorCode()> preProcessHook;
304
305 /**
306 * @brief Hook function for post-processing
307 */
308 boost::function<MoFEMErrorCode()> postProcessHook;
309
310 /**
311 * @brief Hook function for operator
312 */
313 boost::function<MoFEMErrorCode()> operatorHook;
314
315 /** \brief function is run at the beginning of loop
316 *
317 * It is used to zeroing matrices and vectors, calculation of shape
318 * functions on reference element, preprocessing boundary conditions, etc.
319 */
320 virtual MoFEMErrorCode preProcess();
321
322 /** \brief function is run for every finite element
323 *
324 * It is used to calculate element local matrices and assembly. It can be
325 * used for post-processing.
326 */
327 virtual MoFEMErrorCode operator()();
328
329 /** \brief function is run at the end of loop
330 *
331 * It is used to assembly matrices and vectors, calculating global variables,
332 * f.e. total internal energy, ect.
333 *
334 * Iterating over dofs:
335 * Example1 iterating over dofs in row by name of the field
336 * for(_IT_GET_FEROW_BY_NAME_DOFS_FOR_LOOP_(this,"DISPLACEMENT",it)) { ... }
337 *
338 *
339 */
340 virtual MoFEMErrorCode postProcess();
341
342 /**
343 * @brief Get the cache weak ptr object
344 *
345 * \note This store problem information on entities about DOFs. Each problem
346 * store different information. If you iterate over finite elements in
347 * preprocessor of TS solve element, us TS cache in the loop. Otherwise you
348 * will create undetermined behaviour or segmentation error. This is necessary
349 * compromise over bug resilience for memory saving and preformans.
350 *
351 * @return boost::weak_ptr<CacheTuple>
352 */
353 inline boost::weak_ptr<CacheTuple> getCacheWeakPtr() const {
354 return cacheWeakPtr;
355 }
356
357 boost::movelib::unique_ptr<bool> vecAssembleSwitch;
358 boost::movelib::unique_ptr<bool> matAssembleSwitch;
359
360 boost::weak_ptr<CacheTuple> cacheWeakPtr; // cache pointer entity data
361};
362
363/**
364 * \brief structure for User Loop Methods on finite elements
365 * \ingroup mofem_loops
366 *
367 * It can be used to calculate stiffness matrices, residuals, load vectors etc.
368 * It is low level class however in some class users looking for speed and
369 * efficiency, can use it directly.
370 *
371 * This class is used with Interface::loop_finite_elements, where
372 * user overloaded operator FEMethod::operator() is executed for each element in
373 * the problem. Class have to additional methods which are overloaded by user,
374 * FEMethod::preProcess() and FEMethod::postProcess() executed at beginning and
375 * end of the loop over problem elements, respectively.
376 *
377 */
378struct FEMethod : public BasicMethod {
379
380 MoFEMErrorCode query_interface(boost::typeindex::type_index type_index,
381 UnknownInterface **iface) const {
383 *iface = const_cast<FEMethod *>(this);
385 }
386
387 FEMethod() = default;
388
389 std::string feName; ///< Name of finite element
390
391 boost::shared_ptr<const NumeredEntFiniteElement>
392 numeredEntFiniteElementPtr; ///< Pointer to finite element database
393 ///< structure
394
395 /**
396 * @brief get finite element name
397 *
398 * @return std::string
399 */
400 inline auto getFEName() const { return feName; }
401
402 /**
403 * @brief Tet if element to skip element
404 *
405 * If is set and return false elemnent us skiped in
406 * MoFEM::Core::loop_finite_elements
407 *
408 * \note That functionality is used to run elements on particular bit levels
409 *
410 */
411 boost::function<bool(FEMethod *fe_method_ptr)> exeTestHook;
412
413 inline auto getDataDofsPtr() const {
414 return numeredEntFiniteElementPtr->getDataDofsPtr();
415 };
416
417 inline auto getDataVectorDofsPtr() const {
418 return numeredEntFiniteElementPtr->getDataVectorDofsPtr();
419 };
420
422 return numeredEntFiniteElementPtr->getDataFieldEnts();
423 }
424
425 inline boost::shared_ptr<FieldEntity_vector_view> &
427 return const_cast<NumeredEntFiniteElement *>(
430 }
431
432 inline auto &getRowFieldEnts() const {
433 return numeredEntFiniteElementPtr->getRowFieldEnts();
434 };
435
436 inline auto &getRowFieldEntsPtr() const {
437 return numeredEntFiniteElementPtr->getRowFieldEntsPtr();
438 };
439
440 inline auto &getColFieldEnts() const {
441 return numeredEntFiniteElementPtr->getColFieldEnts();
442 };
443
444 inline auto &getColFieldEntsPtr() const {
445 return numeredEntFiniteElementPtr->getColFieldEntsPtr();
446 };
447
448 inline auto getRowDofsPtr() const {
449 return numeredEntFiniteElementPtr->getRowDofsPtr();
450 };
451
452 inline auto getColDofsPtr() const {
453 return numeredEntFiniteElementPtr->getColDofsPtr();
454 };
455
456 inline auto getNumberOfNodes() const;
457
458 inline EntityHandle getFEEntityHandle() const;
459
460 MoFEMErrorCode getNodeData(const std::string field_name, VectorDouble &data,
461 const bool reset_dofs = true);
462
463};
464
465inline auto FEMethod::getNumberOfNodes() const {
466 return moab::CN::VerticesPerEntity(numeredEntFiniteElementPtr->getEntType());
467};
468
469inline EntityHandle FEMethod::getFEEntityHandle() const {
470 return numeredEntFiniteElementPtr->getEnt();
471}
472
473/**
474 * \brief Data structure to exchange data between mofem and User Loop Methods on
475 * entities. \ingroup mofem_loops
476 *
477 * It allows to exchange data between MoFEM and user functions. It stores
478 * information about multi-indices.
479 */
480struct EntityMethod : public BasicMethod {
481
482 MoFEMErrorCode query_interface(boost::typeindex::type_index type_index,
483 UnknownInterface **iface) const {
485 *iface = const_cast<EntityMethod *>(this);
487 }
488
489 EntityMethod() = default;
490
491 boost::shared_ptr<Field> fieldPtr;
492 boost::shared_ptr<FieldEntity> entPtr;
493};
494
495/**
496 * \brief Data structure to exchange data between mofem and User Loop Methods on
497 * entities. \ingroup mofem_loops
498 *
499 * It allows to exchange data between MoFEM and user functions. It stores
500 * information about multi-indices.
501 */
502struct DofMethod : public BasicMethod {
503
504 MoFEMErrorCode query_interface(boost::typeindex::type_index type_index,
505 UnknownInterface **iface) const {
507 *iface = const_cast<DofMethod *>(this);
509 }
510
511 DofMethod() = default;
512
513 boost::shared_ptr<Field> fieldPtr;
514 boost::shared_ptr<DofEntity> dofPtr;
515 boost::shared_ptr<NumeredDofEntity> dofNumeredPtr;
516};
517
518/// \deprecated name changed use DofMethod insead EntMethod
520
521} // namespace MoFEM
522
523#endif // __LOOPMETHODS_HPP__
524
525/**
526 * \defgroup mofem_loops Loops
527 * \ingroup mofem
528 */
multi_index_container< FieldEntityEntFiniteElementAdjacencyMap, indexed_by< ordered_unique< tag< Composite_Unique_mi_tag >, composite_key< FieldEntityEntFiniteElementAdjacencyMap, const_mem_fun< FieldEntityEntFiniteElementAdjacencyMap, UId, &FieldEntityEntFiniteElementAdjacencyMap::getEntUniqueId >, const_mem_fun< FieldEntityEntFiniteElementAdjacencyMap, UId, &FieldEntityEntFiniteElementAdjacencyMap::getFeUniqueId > > >, ordered_non_unique< tag< Unique_mi_tag >, const_mem_fun< FieldEntityEntFiniteElementAdjacencyMap, UId, &FieldEntityEntFiniteElementAdjacencyMap::getEntUniqueId > >, ordered_non_unique< tag< FE_Unique_mi_tag >, const_mem_fun< FieldEntityEntFiniteElementAdjacencyMap, UId, &FieldEntityEntFiniteElementAdjacencyMap::getFeUniqueId > >, ordered_non_unique< tag< FEEnt_mi_tag >, const_mem_fun< FieldEntityEntFiniteElementAdjacencyMap, EntityHandle, &FieldEntityEntFiniteElementAdjacencyMap::getFeHandle > >, ordered_non_unique< tag< Ent_mi_tag >, const_mem_fun< FieldEntityEntFiniteElementAdjacencyMap, EntityHandle, &FieldEntityEntFiniteElementAdjacencyMap::getEntHandle > > > > FieldEntityEntFiniteElementAdjacencyMap_multiIndex
MultiIndex container keeps Adjacencies Element and dof entities adjacencies and vice versa.
multi_index_container< boost::shared_ptr< Field >, indexed_by< hashed_unique< tag< BitFieldId_mi_tag >, const_mem_fun< Field, const BitFieldId &, &Field::getId >, HashBit< BitFieldId >, EqBit< BitFieldId > >, ordered_unique< tag< Meshset_mi_tag >, member< Field, EntityHandle, &Field::meshSet > >, ordered_unique< tag< FieldName_mi_tag >, const_mem_fun< Field, boost::string_ref, &Field::getNameRef > >, ordered_non_unique< tag< BitFieldId_space_mi_tag >, const_mem_fun< Field, FieldSpace, &Field::getSpace > > > > Field_multiIndex
Field_multiIndex for Field.
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:460
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:359
#define DEPRECATED
Definition: definitions.h:30
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:429
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:453
#define BITFEID_SIZE
max number of finite elements
Definition: definitions.h:234
#define THROW_MESSAGE(msg)
Throw MoFEM exception.
Definition: definitions.h:574
multi_index_container< boost::shared_ptr< DofEntity >, indexed_by< ordered_unique< tag< Unique_mi_tag >, const_mem_fun< DofEntity, UId, &DofEntity::getLocalUniqueId > >, ordered_non_unique< tag< Ent_mi_tag >, const_mem_fun< DofEntity, EntityHandle, &DofEntity::getEnt > > > > DofEntity_multiIndex
MultiIndex container keeps DofEntity.
multi_index_container< boost::shared_ptr< RefEntity >, indexed_by< ordered_unique< tag< Ent_mi_tag >, const_mem_fun< RefEntity, EntityHandle, &RefEntity::getEnt > >, ordered_non_unique< tag< Ent_Ent_mi_tag >, const_mem_fun< RefEntity, EntityHandle, &RefEntity::getParentEnt > >, ordered_non_unique< tag< Composite_EntType_and_ParentEntType_mi_tag >, composite_key< RefEntity, const_mem_fun< RefEntity, EntityType, &RefEntity::getEntType >, const_mem_fun< RefEntity, EntityType, &RefEntity::getParentEntType > > >, ordered_non_unique< tag< Composite_ParentEnt_And_EntType_mi_tag >, composite_key< RefEntity, const_mem_fun< RefEntity, EntityType, &RefEntity::getEntType >, const_mem_fun< RefEntity, EntityHandle, &RefEntity::getParentEnt > > > > > RefEntity_multiIndex
multi_index_container< boost::shared_ptr< RefElement >, indexed_by< ordered_unique< tag< Ent_mi_tag >, const_mem_fun< RefElement::interface_type_RefEntity, EntityHandle, &RefElement::getEnt > > > > RefElement_multiIndex
multi_index_container< boost::shared_ptr< EntFiniteElement >, indexed_by< ordered_unique< tag< Unique_mi_tag >, const_mem_fun< EntFiniteElement, UId, &EntFiniteElement::getLocalUniqueId > >, ordered_non_unique< tag< Ent_mi_tag >, const_mem_fun< EntFiniteElement::interface_type_RefEntity, EntityHandle, &EntFiniteElement::getEnt > >, ordered_non_unique< tag< FiniteElement_name_mi_tag >, const_mem_fun< EntFiniteElement::interface_type_FiniteElement, boost::string_ref, &EntFiniteElement::getNameRef > >, ordered_non_unique< tag< Composite_Name_And_Ent_mi_tag >, composite_key< EntFiniteElement, const_mem_fun< EntFiniteElement::interface_type_FiniteElement, boost::string_ref, &EntFiniteElement::getNameRef >, const_mem_fun< EntFiniteElement::interface_type_RefEntity, EntityHandle, &EntFiniteElement::getEnt > > > > > EntFiniteElement_multiIndex
MultiIndex container for EntFiniteElement.
multi_index_container< boost::shared_ptr< FiniteElement >, indexed_by< hashed_unique< tag< FiniteElement_Meshset_mi_tag >, member< FiniteElement, EntityHandle, &FiniteElement::meshset > >, hashed_unique< tag< BitFEId_mi_tag >, const_mem_fun< FiniteElement, BitFEId, &FiniteElement::getId >, HashBit< BitFEId >, EqBit< BitFEId > >, ordered_unique< tag< FiniteElement_name_mi_tag >, const_mem_fun< FiniteElement, boost::string_ref, &FiniteElement::getNameRef > > > > FiniteElement_multiIndex
MultiIndex for entities for FiniteElement.
const FTensor::Tensor2< T, Dim, Dim > Vec
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:67
UBlasVector< double > VectorDouble
Definition: Types.hpp:79
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:21
DEPRECATED typedef DofMethod EntMethod
std::vector< boost::weak_ptr< FieldEntity > > FieldEntity_vector_view
Field_multiIndex::index< FieldName_mi_tag >::type::iterator field_it
multi_index_container< boost::shared_ptr< FieldEntity >, indexed_by< ordered_unique< tag< Unique_mi_tag >, member< FieldEntity, UId, &FieldEntity::localUId > >, ordered_non_unique< tag< Ent_mi_tag >, const_mem_fun< FieldEntity::interface_type_RefEntity, EntityHandle, &FieldEntity::getEnt > > > > FieldEntity_multiIndex
constexpr auto field_name
Data structure to exchange data between mofem and User Loop Methods.
int loopSize
local number oe methods to process
boost::movelib::unique_ptr< bool > matAssembleSwitch
const FieldEntity_multiIndex * entitiesPtr
raw pointer to container of field entities
unsigned int getFieldBitNumber(std::string field_name) const
virtual ~BasicMethod()=default
MoFEMErrorCode copyBasicMethod(const BasicMethod &basic)
Copy data from other base method to this base method.
auto getHiFERank() const
Get upper rank in loop for iterating elements.
const FiniteElement_multiIndex * finiteElementsPtr
raw pointer to container finite elements
const RefElement_multiIndex * refinedFiniteElementsPtr
container of mofem finite element entities
int getLoopSize() const
get loop size
std::pair< int, int > loHiFERank
Llo and hi processor rank of iterated entities.
boost::function< MoFEMErrorCode()> preProcessHook
Hook function for pre-processing.
const DofEntity_multiIndex * dofsPtr
raw pointer container of dofs
int rAnk
processor rank
boost::weak_ptr< CacheTuple > getCacheWeakPtr() const
Get the cache weak ptr object.
MoFEMErrorCode query_interface(boost::typeindex::type_index type_index, UnknownInterface **iface) const
const Field_multiIndex * fieldsPtr
raw pointer to fields container
int nInTheLoop
number currently of processed method
auto getLoHiFERank() const
Get lo and hi processor rank of iterated entities.
const RefEntity_multiIndex * refinedEntitiesPtr
container of mofem dof entities
boost::movelib::unique_ptr< bool > vecAssembleSwitch
virtual MoFEMErrorCode postProcess()
function is run at the end of loop
const FieldEntityEntFiniteElementAdjacencyMap_multiIndex * adjacenciesPtr
virtual MoFEMErrorCode operator()()
function is run for every finite element
virtual MoFEMErrorCode preProcess()
function is run at the beginning of loop
const EntFiniteElement_multiIndex * finiteElementsEntitiesPtr
const Problem * problemPtr
raw pointer to problem
auto getLoFERank() const
Get upper rank in loop for iterating elements.
boost::function< MoFEMErrorCode()> operatorHook
Hook function for operator.
int sIze
number of processors in communicator
boost::weak_ptr< CacheTuple > cacheWeakPtr
int getNinTheLoop() const
get number of evaluated element in the loop
boost::function< MoFEMErrorCode()> postProcessHook
Hook function for post-processing.
Data structure to exchange data between mofem and User Loop Methods on entities.
boost::shared_ptr< DofEntity > dofPtr
DofMethod()=default
MoFEMErrorCode query_interface(boost::typeindex::type_index type_index, UnknownInterface **iface) const
boost::shared_ptr< NumeredDofEntity > dofNumeredPtr
boost::shared_ptr< Field > fieldPtr
Data structure to exchange data between mofem and User Loop Methods on entities.
EntityMethod()=default
boost::shared_ptr< Field > fieldPtr
MoFEMErrorCode query_interface(boost::typeindex::type_index type_index, UnknownInterface **iface) const
boost::shared_ptr< FieldEntity > entPtr
structure for User Loop Methods on finite elements
std::string feName
Name of finite element.
auto & getColFieldEntsPtr() const
FEMethod()=default
auto & getRowFieldEnts() const
auto getFEName() const
get finite element name
MoFEMErrorCode getNodeData(const std::string field_name, VectorDouble &data, const bool reset_dofs=true)
auto getDataDofsPtr() const
const FieldEntity_vector_view & getDataFieldEnts() const
auto getDataVectorDofsPtr() const
boost::shared_ptr< FieldEntity_vector_view > & getDataFieldEntsPtr() const
auto getColDofsPtr() const
EntityHandle getFEEntityHandle() const
MoFEMErrorCode query_interface(boost::typeindex::type_index type_index, UnknownInterface **iface) const
auto & getRowFieldEntsPtr() const
auto getNumberOfNodes() const
auto getRowDofsPtr() const
boost::shared_ptr< const NumeredEntFiniteElement > numeredEntFiniteElementPtr
auto & getColFieldEnts() const
boost::function< bool(FEMethod *fe_method_ptr)> exeTestHook
Tet if element to skip element.
MultiIndex Tag for field name.
data structure for ksp (linear solver) context
Definition: LoopMethods.hpp:74
MoFEMErrorCode copyKsp(const KspMethod &ksp)
copy data form another method
Definition: LoopMethods.cpp:68
KSPContext ksp_ctx
Context.
Definition: LoopMethods.hpp:96
virtual ~KspMethod()=default
KSPContext
pass information about context of KSP/DM for with finite element is computed
Definition: LoopMethods.hpp:83
KSP ksp
KSP solver.
Definition: LoopMethods.hpp:97
MoFEMErrorCode query_interface(boost::typeindex::type_index type_index, UnknownInterface **iface) const
Definition: LoopMethods.cpp:58
Partitioned (Indexed) Finite Element in Problem.
static constexpr Switches CtxSetA
Definition: LoopMethods.hpp:47
static constexpr Switches CtxSetX
Definition: LoopMethods.hpp:49
static constexpr Switches CtxSetX_TT
Definition: LoopMethods.hpp:51
static constexpr Switches CtxSetNone
Definition: LoopMethods.hpp:45
static constexpr Switches CtxSetF
Definition: LoopMethods.hpp:46
std::bitset< 8 > Switches
Definition: LoopMethods.hpp:43
virtual ~PetscData()=default
static constexpr Switches CtxSetX_T
Definition: LoopMethods.hpp:50
MoFEMErrorCode query_interface(boost::typeindex::type_index type_index, UnknownInterface **iface) const
Definition: LoopMethods.cpp:34
static constexpr Switches CtxSetB
Definition: LoopMethods.hpp:48
static constexpr Switches CtxSetTime
Definition: LoopMethods.hpp:52
MoFEMErrorCode copyPetscData(const PetscData &petsc_data)
Definition: LoopMethods.cpp:44
keeps basic data about problem
data structure for snes (nonlinear solver) context
MoFEMErrorCode copySnes(const SnesMethod &snes)
Copy snes data.
Definition: LoopMethods.cpp:88
Vec & snes_f
residual
MoFEMErrorCode query_interface(boost::typeindex::type_index type_index, UnknownInterface **iface) const
Definition: LoopMethods.cpp:78
Vec & snes_x
state vector
virtual ~SnesMethod()=default
Mat & snes_B
preconditioner of jacobian matrix
SNESContext snes_ctx
Mat & snes_A
jacobian matrix
SNES snes
snes solver
data structure for TS (time stepping) context
TS ts
time solver
PetscReal ts_t
time
Vec & ts_F
residual vector
MoFEMErrorCode query_interface(boost::typeindex::type_index type_index, UnknownInterface **iface) const
Definition: LoopMethods.cpp:98
Vec & ts_u_tt
second time derivative of state vector
TSContext ts_ctx
PetscReal ts_a
shift for U_t (see PETSc Time Solver)
Vec & ts_u_t
time derivative of state vector
PetscInt ts_step
time step
Vec & ts_u
state vector
virtual ~TSMethod()=default
PetscReal ts_aa
shift for U_tt shift for U_tt
MoFEMErrorCode copyTs(const TSMethod &ts)
Copy TS solver data.
Mat & ts_B
Preconditioner for ts_A.
base class for all interface classes