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