v0.10.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 /*
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 
22 namespace MoFEM {
23 
38 
39 struct PetscData : public UnknownInterface {
40 
42  UnknownInterface **iface) const;
43 
44  PetscData();
45 
46  virtual ~PetscData() = default;
47 
48  enum DataContext {
50  CTX_SET_F = 1 << 0,
51  CTX_SET_A = 1 << 1,
52  CTX_SET_B = 1 << 2,
53  CTX_SET_X = 1 << 3,
54  CTX_SET_X_T = 1 << 4,
55  CTX_SET_X_TT = 1 << 6,
56  CTX_SET_TIME = 1 << 7
57  };
58 
59  using Switches = std::bitset<8>;
60 
69 
71 
72  Vec f;
73  Mat A;
74  Mat B;
75  Vec x;
78 };
79 
80 /**
81  * \brief data structure for ksp (linear solver) context
82  * \ingroup mofem_loops
83  *
84  * Struture stores context data which are set in functions run by PETSc SNES
85  * functions.
86  *
87  */
88 struct KspMethod : virtual public PetscData {
89 
91  UnknownInterface **iface) const;
92  /**
93  * \brief pass information about context of KSP/DM for with finite element is
94  * computed
95  */
97 
98  KspMethod();
99 
100  virtual ~KspMethod() = default;
101 
102  /**
103  * \brief copy data form another method
104  * @param ksp ksp method
105  * @return error code
106  */
108 
109  KSPContext ksp_ctx; ///< Context
110  KSP ksp; ///< KSP solver
111 
113  Mat &ksp_A;
114  Mat &ksp_B;
115 };
116 
117 /**
118  * \brief data structure for snes (nonlinear solver) context
119  * \ingroup mofem_loops
120  *
121  * Structure stores context data which are set in functions run by PETSc SNES
122  * functions.
123  *
124  */
125 struct SnesMethod : virtual protected PetscData {
126 
128  UnknownInterface **iface) const;
129 
131 
132  SnesMethod();
133 
134  virtual ~SnesMethod() = default;
135 
136  /**
137  * \brief Copy snes data
138  */
140 
142 
143  /**
144  * @deprecated Avoid using values by hand.
145  */
147 
148  SNES snes; ///< snes solver
149  Vec &snes_x; ///< state vector
150  Vec &snes_f; ///< residual
151  Mat &snes_A; ///< jacobian matrix
152  Mat &snes_B; ///< preconditioner of jacobian matrix
153 };
154 
156  snes_ctx = ctx;
157  return 0;
158 }
159 
160 /**
161  * \brief data structure for TS (time stepping) context
162  * \ingroup mofem_loops
163  *
164  * Structure stores context data which are set in functions run by PETSc Time
165  * Stepping functions.
166  */
167 struct TSMethod : virtual protected PetscData {
168 
170  UnknownInterface **iface) const;
171 
172  enum TSContext {
178  CTX_TSNONE
179  };
180 
181  TSMethod();
182 
183  virtual ~TSMethod() = default;
184 
185  /// \brief Copy TS solver data
187 
188  TS ts; ///< time solver
189 
191 
192  /**
193  * @deprecated Avoid using values by hand.
194  */
196 
197  PetscInt ts_step; ///< time step
198  PetscReal ts_a; ///< shift for U_t (see PETSc Time Solver)
199  PetscReal ts_aa; ///< shift for U_tt shift for U_tt
200  PetscReal ts_t; ///< time
201 
202  Vec &ts_u; ///< state vector
203  Vec &ts_u_t; ///< time derivative of state vector
204  Vec &ts_u_tt; ///< second time derivative of state vector
205  Vec &ts_F; ///< residual vector
206 
207  Mat &ts_A; ///< Jacobian of G(U) = F(t,U,W+v*U,W'+a*U), equivalent to dF/dU +
208  ///< v*dF/dU_t + a*dF/dU_tt
209  Mat &ts_B; ///< Preconditioner for ts_A
210 };
211 
213  ts_ctx = ctx;
214  return 0;
215 }
216 
217 /**
218  * \brief Data structure to exchange data between mofem and User Loop Methods.
219  * \ingroup mofem_loops
220  *
221  * It allows to exchange data between MoFEM and user functions. It stores
222  * information about multi-indices.
223  *
224  */
226 
228  UnknownInterface **iface) const {
229  if (uuid == IDD_MOFEMBasicMethod) {
230  *iface = const_cast<BasicMethod *>(this);
232  }
233  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "unknown interface");
234  }
235 
236  BasicMethod();
237  virtual ~BasicMethod() = default;
238 
239  /**
240  * @brief number currently of processed method
241  */
243 
244  /**
245  * @brief local number oe methods to process
246  */
247  int loopSize;
248 
249  /** \brief get number of evaluated element in the loop
250  */
251  inline int getNinTheLoop() const { return nInTheLoop; }
252 
253  /** \brief get loop size
254  */
255  inline int getLoopSize() const { return loopSize; }
256 
257  int rAnk; ///< processor rank
258 
259  int sIze; ///< number of processors in communicator
260 
262  *refinedEntitiesPtr; ///< container of mofem dof entities
263 
265  *refinedFiniteElementsPtr; ///< container of mofem finite element entities
266 
267  const Problem *problemPtr; ///< raw pointer to problem
268 
269  const Field_multiIndex *fieldsPtr; ///< raw pointer to fields container
270 
272  *entitiesPtr; ///< raw pointer to container of field entities
273 
274  const DofEntity_multiIndex *dofsPtr; ///< raw pointer container of dofs
275 
277  *finiteElementsPtr; ///< raw pointer to container finite elements
278 
280  *finiteElementsEntitiesPtr; ///< raw pointer to container finite elements
281  ///< entities
282 
284  *adjacenciesPtr; ///< raw pointer to container to adjacencies between dofs
285  ///< and finite elements
286 
287  inline unsigned int getFieldBitNumber(std::string field_name) const {
288  if (fieldsPtr) {
289  auto field_it = fieldsPtr->get<FieldName_mi_tag>().find(field_name);
290  if (field_it != fieldsPtr->get<FieldName_mi_tag>().end())
291  return (*field_it)->getBitNumber();
292  else
293  return BITFEID_SIZE;
294  } else {
295  THROW_MESSAGE("Pointer to fields multi-index is not set");
296  return BITFEID_SIZE;
297  }
298  }
299 
300  /**
301  * @brief Copy data from other base method to this base method
302  *
303  * @param basic
304  * @return MoFEMErrorCode
305  */
307 
308  /**
309  * @brief Hook function for pre-processing
310  */
311  boost::function<MoFEMErrorCode()> preProcessHook;
312 
313  /**
314  * @brief Hook function for post-processing
315  */
316  boost::function<MoFEMErrorCode()> postProcessHook;
317 
318  /**
319  * @brief Hook function for operator
320  */
321  boost::function<MoFEMErrorCode()> operatorHook;
322 
323  /** \brief function is run at the beginning of loop
324  *
325  * It is used to zeroing matrices and vectors, calculation of shape
326  * functions on reference element, preprocessing boundary conditions, etc.
327  */
328  virtual MoFEMErrorCode preProcess();
329 
330  /** \brief function is run for every finite element
331  *
332  * It is used to calculate element local matrices and assembly. It can be
333  * used for post-processing.
334  */
335  virtual MoFEMErrorCode operator()();
336 
337  /** \brief function is run at the end of loop
338  *
339  * It is used to assembly matrices and vectors, calculating global variables,
340  * f.e. total internal energy, ect.
341  *
342  * Iterating over dofs:
343  * Example1 iterating over dofs in row by name of the field
344  * for(_IT_GET_FEROW_BY_NAME_DOFS_FOR_LOOP_(this,"DISPLACEMENT",it)) { ... }
345  *
346  *
347  */
348  virtual MoFEMErrorCode postProcess();
349 
350  boost::movelib::unique_ptr<bool> vecAssembleSwitch;
351  boost::movelib::unique_ptr<bool> matAssembleSwitch;
352 };
353 
354 /**
355  * \brief structure for User Loop Methods on finite elements
356  * \ingroup mofem_loops
357  *
358  * It can be used to calculate stiffness matrices, residuals, load vectors etc.
359  * It is low level class however in some class users looking for speed and
360  * efficiency, can use it directly.
361  *
362  * This class is used with Interface::loop_finite_elements, where
363  * user overloaded operator FEMethod::operator() is executed for each element in
364  * the problem. Class have to additional methods which are overloaded by user,
365  * FEMethod::preProcess() and FEMethod::postProcess() executed at beginning and
366  * end of the loop over problem elements, respectively.
367  *
368  */
369 struct FEMethod : public BasicMethod {
370 
372  UnknownInterface **iface) const {
374  if (uuid == IDD_MOFEMFEMethod) {
375  *iface = const_cast<FEMethod *>(this);
377  }
378 
379  ierr = query_interface(uuid, iface);
380  CHKERRG(ierr);
382  }
383 
384  FEMethod();
385 
386  std::string feName; ///< Name of finite element
387 
388  boost::shared_ptr<const NumeredEntFiniteElement>
389  numeredEntFiniteElementPtr; ///< Pointer to finite element database
390  ///< structure
391 
392  inline auto getDataDofsPtr() const {
393  return numeredEntFiniteElementPtr->getDataDofsPtr();
394  };
395 
396  inline auto getDataVectorDofsPtr() const {
397  return numeredEntFiniteElementPtr->getDataVectorDofsPtr();
398  };
399 
401  return numeredEntFiniteElementPtr->getDataFieldEnts();
402  }
403 
404  inline boost::shared_ptr<FieldEntity_vector_view> &
406  return const_cast<NumeredEntFiniteElement *>(
409  }
410 
411  inline auto &getRowFieldEnts() const {
412  return numeredEntFiniteElementPtr->getRowFieldEnts();
413  };
414 
415  inline auto &getRowFieldEntsPtr() const {
416  return numeredEntFiniteElementPtr->getRowFieldEntsPtr();
417  };
418 
419  inline auto &getColFieldEnts() const {
420  return numeredEntFiniteElementPtr->getColFieldEnts();
421  };
422 
423  inline auto &getColFieldEntsPtr() const {
424  return numeredEntFiniteElementPtr->getColFieldEntsPtr();
425  };
426 
427  inline auto getRowDofsPtr() const {
428  return numeredEntFiniteElementPtr->getRowDofsPtr();
429  };
430 
431  inline auto getColDofsPtr() const {
432  return numeredEntFiniteElementPtr->getColDofsPtr();
433  };
434 
435  /// \brief Get number of DOFs on element
436  MoFEMErrorCode getNumberOfNodes(int &num_nodes) const;
437 
438  inline EntityHandle getFEEntityHandle() const;
439 
440  MoFEMErrorCode getNodeData(const std::string field_name, VectorDouble &data,
441  const bool reset_dofs = true);
442 
443  template <class MULTIINDEX>
444  typename MULTIINDEX::iterator get_begin(const MULTIINDEX &index,
445  const std::string &field_name,
446  const EntityType type) const {
447  return index.lower_bound(boost::make_tuple(field_name, type));
448  }
449 
450  template <class MULTIINDEX>
451  typename MULTIINDEX::iterator get_end(const MULTIINDEX &index,
452  const std::string &field_name,
453  const EntityType type) const {
454  return index.upper_bound(boost::make_tuple(field_name, type));
455  }
456 
457 /** \brief loop over all dofs which are on a particular FE data, field and
458  * entity type \ingroup mofem_loops
459  */
460 #define _IT_GET_FEDATA_BY_TYPE_DOFS_FOR_LOOP_(FE, NAME, TYPE, IT) \
461  auto IT = FE->get_begin<FEDofEntityByNameAndEnt>( \
462  FE->dataPtr->get<FEDofEntityByNameAndEnt>(), NAME, \
463  get_id_for_min_type(TYPE)); \
464  IT != FE->get_end<FEDofEntityByNameAndEnt>( \
465  FE->dataPtr->get<FEDofEntityByNameAndEnt>(), NAME, \
466  get_id_for_max_type(TYPE)); \
467  IT++
468 
469  template <class MULTIINDEX>
470  typename MULTIINDEX::iterator get_begin(const MULTIINDEX &index,
471  const std::string &field_name) const {
472  return index.lower_bound(field_name);
473  }
474  template <class MULTIINDEX>
475  typename MULTIINDEX::iterator get_end(const MULTIINDEX &index,
476  const std::string &field_name) const {
477  return index.upper_bound(field_name);
478  }
479 
480 /** \brief loop over all dofs which are on a particular FE row and field
481  * \ingroup mofem_loops
482  */
483 #define _IT_GET_FEROW_BY_NAME_DOFS_FOR_LOOP_(FE, NAME, IT) \
484  auto IT = FE->get_begin<FENumeredDofEntityByUId>( \
485  FE->getRowDofs().get<Unique_mi_tag>(), \
486  FieldEntity::getLoBitNumberUId(FE->getFieldBitNumber(NAME))); \
487  IT != FE->get_end<FENumeredDofEntityByUId>( \
488  FE->getRowDofs()->get<Unique_mi_tag>(), \
489  FieldEntity::getHiBitNumberUId(FE->getFieldBitNumber(NAME))); \
490  IT++
491 
492 /** \brief loop over all dofs which are on a particular FE column and field
493  * \ingroup mofem_loops
494  */
495 #define _IT_GET_FECOL_BY_NAME_DOFS_FOR_LOOP_(FE, NAME, IT) \
496  auto IT = FE->get_begin<FENumeredDofEntityByUId>( \
497  FE->getColDofs().get<Unique_mi_tag>(), \
498  FieldEntity::getLoBitNumberUId(FE->getFieldBitNumber(NAME))); \
499  IT != FE->get_end<FENumeredDofEntityByUId>( \
500  FE->getColDofs()->get<Unique_mi_tag>(), \
501  FieldEntity::getHiBitNumberUId(FE->getFieldBitNumber(NAME))); \
502  IT++
503 
504 /** \brief loop over all dofs which are on a particular FE data and field
505  * \ingroup mofem_loops
506  */
507 #define _IT_GET_FEDATA_BY_NAME_DOFS_FOR_LOOP_(FE, NAME, IT) \
508  auto IT = FE->get_begin<FEDofEntityByUId>( \
509  FE->dataPtr->get<Unique_mi_tag>(), \
510  FieldEntity::getLoBitNumberUId(FE->getFieldBitNumber(NAME))); \
511  IT != FE->get_end<FEDofEntityByUId>( \
512  FE->dataPtr->get<Unique_mi_tag>(), \
513  FieldEntity::getHiBitNumberUId(FE->getFieldBitNumber(NAME))); \
514  IT++
515 };
516 
517 inline EntityHandle FEMethod::getFEEntityHandle() const {
518  return numeredEntFiniteElementPtr->getEnt();
519 }
520 
521 /**
522  * \brief Data structure to exchange data between mofem and User Loop Methods on
523  * entities. \ingroup mofem_loops
524  *
525  * It allows to exchange data between MoFEM and user functions. It stores
526  * information about multi-indices.
527  */
528 struct EntityMethod : public BasicMethod {
529 
531  UnknownInterface **iface) const {
533  if (uuid == IDD_MOFEMEntityMethod) {
534  *iface = const_cast<EntityMethod *>(this);
536  }
537  CHKERR query_interface(uuid, iface);
539  }
540 
541  EntityMethod();
542 
543  boost::shared_ptr<Field> fieldPtr;
544  boost::shared_ptr<FieldEntity> entPtr;
545 };
546 
547 /**
548  * \brief Data structure to exchange data between mofem and User Loop Methods on
549  * entities. \ingroup mofem_loops
550  *
551  * It allows to exchange data between MoFEM and user functions. It stores
552  * information about multi-indices.
553  */
554 struct DofMethod : public BasicMethod {
555 
557  UnknownInterface **iface) const {
559  if (uuid == IDD_MOFEMDofMethod) {
560  *iface = const_cast<DofMethod *>(this);
562  }
563 
564  CHKERR query_interface(uuid, iface);
566  }
567 
568  DofMethod();
569 
570  boost::shared_ptr<Field> fieldPtr;
571  boost::shared_ptr<DofEntity> dofPtr;
572  boost::shared_ptr<NumeredDofEntity> dofNumeredPtr;
573 };
574 
575 /// \deprecated name changed use DofMethod insead EntMethod
577 
578 } // namespace MoFEM
579 
580 #endif // __LOOPMETHODS_HPP__
581 
582 /**
583  * \defgroup mofem_loops Loops
584  * \ingroup mofem
585  */
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:516
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:415
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:552
@ MOFEM_DATA_INCONSISTENCY
Definition: definitions.h:123
#define DEPRECATED
Definition: definitions.h:29
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:485
#define CHKERR
Inline error check.
Definition: definitions.h:604
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:509
#define BITFEID_SIZE
max number of finite elements
Definition: definitions.h:288
@ DOF_METHOD
Definition: definitions.h:83
@ ENTITY_METHOD
Definition: definitions.h:82
@ KSP_METHOD
Definition: definitions.h:77
@ FE_METHOD
Definition: definitions.h:81
@ SNES_METHOD
Definition: definitions.h:78
@ TS_METHOD
Definition: definitions.h:79
@ PETSC_DATA_METHOD
Definition: definitions.h:76
@ BASIC_METHOD
Definition: definitions.h:80
#define THROW_MESSAGE(a)
Throw MoFEM exception.
Definition: definitions.h:628
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< 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
MultiIndex container keeps FieldEntity.
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
static MoFEMErrorCodeGeneric< PetscErrorCode > ierr
Definition: Exceptions.hpp:87
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:67
std::bitset< BITINTERFACEUID_SIZE > BitIntefaceId
Definition: Types.hpp:56
ublas::vector< double, DoubleAllocator > VectorDouble
Definition: Types.hpp:74
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:21
DEPRECATED typedef DofMethod EntMethod
static const MOFEMuuid IDD_MOFEMSnesMethod
Definition: LoopMethods.hpp:28
static const MOFEMuuid IDD_MOFEMKspMethod
Definition: LoopMethods.hpp:26
static const MOFEMuuid IDD_MOFEMDofMethod
Definition: LoopMethods.hpp:36
static const MOFEMuuid IDD_MOFEMTsMethod
Definition: LoopMethods.hpp:30
static const MOFEMuuid IDD_MOFEMPetscDataMethod
Definition: LoopMethods.hpp:24
static const MOFEMuuid IDD_MOFEMFEMethod
Definition: LoopMethods.hpp:33
std::vector< boost::weak_ptr< FieldEntity > > FieldEntity_vector_view
static const MOFEMuuid IDD_MOFEMBasicMethod
Definition: LoopMethods.hpp:31
static const MOFEMuuid IDD_MOFEMEntityMethod
Definition: LoopMethods.hpp:34
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.
MoFEMErrorCode query_interface(const MOFEMuuid &uuid, UnknownInterface **iface) const
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
boost::function< MoFEMErrorCode()> preProcessHook
Hook function for pre-processing.
const DofEntity_multiIndex * dofsPtr
raw pointer container of dofs
int rAnk
processor rank
const Field_multiIndex * fieldsPtr
raw pointer to fields container
int nInTheLoop
number currently of processed method
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
boost::function< MoFEMErrorCode()> operatorHook
Hook function for operator.
int sIze
number of processors in communicator
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
MoFEMErrorCode query_interface(const MOFEMuuid &uuid, 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.
MoFEMErrorCode query_interface(const MOFEMuuid &uuid, UnknownInterface **iface) const
boost::shared_ptr< Field > fieldPtr
boost::shared_ptr< FieldEntity > entPtr
structure for User Loop Methods on finite elements
std::string feName
Name of finite element.
MoFEMErrorCode query_interface(const MOFEMuuid &uuid, UnknownInterface **iface) const
MULTIINDEX::iterator get_begin(const MULTIINDEX &index, const std::string &field_name, const EntityType type) const
MULTIINDEX::iterator get_begin(const MULTIINDEX &index, const std::string &field_name) const
MoFEMErrorCode getNodeData(const std::string field_name, VectorDouble &data, const bool reset_dofs=true)
auto getDataDofsPtr() const
auto & getColFieldEntsPtr() const
MoFEMErrorCode getNumberOfNodes(int &num_nodes) const
Get number of DOFs on element.
auto getDataVectorDofsPtr() const
MULTIINDEX::iterator get_end(const MULTIINDEX &index, const std::string &field_name, const EntityType type) const
auto & getColFieldEnts() const
auto getColDofsPtr() const
EntityHandle getFEEntityHandle() const
MULTIINDEX::iterator get_end(const MULTIINDEX &index, const std::string &field_name) const
const FieldEntity_vector_view & getDataFieldEnts() const
boost::shared_ptr< FieldEntity_vector_view > & getDataFieldEntsPtr() const
auto & getRowFieldEntsPtr() const
auto getRowDofsPtr() const
boost::shared_ptr< const NumeredEntFiniteElement > numeredEntFiniteElementPtr
auto & getRowFieldEnts() const
MultiIndex Tag for field name.
data structure for ksp (linear solver) context
Definition: LoopMethods.hpp:88
MoFEMErrorCode copyKsp(const KspMethod &ksp)
copy data form another method
Definition: LoopMethods.cpp:62
KSPContext ksp_ctx
Context.
virtual ~KspMethod()=default
KSPContext
pass information about context of KSP/DM for with finite element is computed
Definition: LoopMethods.hpp:96
MoFEMErrorCode query_interface(const MOFEMuuid &uuid, UnknownInterface **iface) const
Definition: LoopMethods.cpp:48
KSP ksp
KSP solver.
MoFEM interface unique ID.
Partitioned (Indexed) Finite Element in Problem.
static constexpr Switches CtxSetA
Definition: LoopMethods.hpp:63
static constexpr Switches CtxSetX
Definition: LoopMethods.hpp:65
static constexpr Switches CtxSetX_TT
Definition: LoopMethods.hpp:67
static constexpr Switches CtxSetNone
Definition: LoopMethods.hpp:61
static constexpr Switches CtxSetF
Definition: LoopMethods.hpp:62
std::bitset< 8 > Switches
Definition: LoopMethods.hpp:59
virtual ~PetscData()=default
static constexpr Switches CtxSetX_T
Definition: LoopMethods.hpp:66
static constexpr Switches CtxSetB
Definition: LoopMethods.hpp:64
static constexpr Switches CtxSetTime
Definition: LoopMethods.hpp:68
MoFEMErrorCode query_interface(const MOFEMuuid &uuid, UnknownInterface **iface) const
Definition: LoopMethods.cpp:33
keeps basic data about problem
data structure for snes (nonlinear solver) context
MoFEMErrorCode copySnes(const SnesMethod &snes)
Copy snes data.
Definition: LoopMethods.cpp:86
Vec & snes_f
residual
DEPRECATED MoFEMErrorCode setSnesCtx(SNESContext ctx)
Vec & snes_x
state vector
virtual ~SnesMethod()=default
Mat & snes_B
preconditioner of jacobian matrix
SNESContext snes_ctx
Mat & snes_A
jacobian matrix
MoFEMErrorCode query_interface(const MOFEMuuid &uuid, UnknownInterface **iface) const
Definition: LoopMethods.cpp:73
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(const MOFEMuuid &uuid, 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.
DEPRECATED MoFEMErrorCode setTsCtx(TSContext ctx)
Mat & ts_B
Preconditioner for ts_A.
base class for all interface classes