v0.9.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 
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;
76  Vec x_t;
77  Vec x_tt;
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 
112  Vec &ksp_f;
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 {
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_tt (see PETSc Time Solver)
199  PetscReal ts_v; ///< shift for U_t shift for U_t
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  /**
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  boost::movelib::unique_ptr<bool> vecAssembleSwitch;
338  boost::movelib::unique_ptr<bool> matAssembleSwitch;
339 };
340 
341 /**
342  * \brief structure for User Loop Methods on finite elements
343  * \ingroup mofem_loops
344  *
345  * It can be used to calculate stiffness matrices, residuals, load vectors etc.
346  * It is low level class however in some class users looking for speed and
347  * efficiency, can use it directly.
348  *
349  * This class is used with Interface::loop_finite_elements, where
350  * user overloaded operator FEMethod::operator() is executed for each element in
351  * the problem. Class have to additional methods which are overloaded by user,
352  * FEMethod::preProcess() and FEMethod::postProcess() executed at beginning and
353  * end of the loop over problem elements, respectively.
354  *
355  */
356 struct FEMethod : public BasicMethod {
357 
359  UnknownInterface **iface) const {
361  if (uuid == IDD_MOFEMFEMethod) {
362  *iface = const_cast<FEMethod *>(this);
364  }
365 
366  ierr = query_interface(uuid, iface);
367  CHKERRG(ierr);
369  }
370 
371  FEMethod();
372 
373  std::string feName; ///< Name of finite element
374 
375  boost::shared_ptr<const NumeredEntFiniteElement>
376  numeredEntFiniteElementPtr; ///< Pointer to finite element database
377  ///< structure
378  boost::shared_ptr<const FENumeredDofEntity_multiIndex>
379  rowPtr; ///< Pointer to finite element rows dofs view
380  boost::shared_ptr<const FENumeredDofEntity_multiIndex>
381  colPtr; ///< Pointer to finite element columns dofs view
382  boost::shared_ptr<const FEDofEntity_multiIndex>
383  dataPtr; ///< Pointer to finite element data dofs
384 
385  boost::shared_ptr<const FieldEntity_vector_view>
386  rowFieldEntsPtr; ///< Pointer to finite element field entities row view
387  boost::shared_ptr<const FieldEntity_vector_view>
388  colFieldEntsPtr; ///< Pointer to finite element field entities column view
389  boost::shared_ptr<const FieldEntity_multiIndex_spaceType_view>
390  dataFieldEntsPtr; ///< Pointer to finite element field entities data view
391 
392  /// \brief Get number of DOFs on element
393  MoFEMErrorCode getNumberOfNodes(int &num_nodes) const;
394 
395  inline EntityHandle getFEEntityHandle() const;
396 
397  MoFEMErrorCode getNodeData(const std::string field_name, VectorDouble &data,
398  const bool reset_dofs = true);
399 
400  template <class MULTIINDEX>
401  typename MULTIINDEX::iterator get_begin(const MULTIINDEX &index,
402  const std::string &field_name,
403  const EntityType type) const {
404  return index.lower_bound(boost::make_tuple(field_name, type));
405  }
406  template <class MULTIINDEX>
407  typename MULTIINDEX::iterator get_end(const MULTIINDEX &index,
408  const std::string &field_name,
409  const EntityType type) const {
410  return index.upper_bound(boost::make_tuple(field_name, type));
411  }
412 
413 /** \brief loop over all dofs which are on a particular FE row, field and entity
414  * type \ingroup mofem_loops
415  */
416 #define _IT_GET_FEROW_BY_TYPE_DOFS_FOR_LOOP_(FE, NAME, TYPE, IT) \
417  auto IT = FE->get_begin<FENumeredDofEntityByNameAndType>( \
418  FE->rowPtr->get<Composite_Name_And_Type_mi_tag>(), NAME, TYPE); \
419  IT != FE->get_end<FENumeredDofEntityByNameAndType>( \
420  FE->rowPtr->get<Composite_Name_And_Type_mi_tag>(), NAME, TYPE); \
421  IT++
422 
423 /** \brief loop over all dofs which are on a particular FE column, field and
424  * entity type \ingroup mofem_loops
425  */
426 #define _IT_GET_FECOL_BY_TYPE_DOFS_FOR_LOOP_(FE, NAME, TYPE, IT) \
427  auto IT = FE->get_begin<FENumeredDofEntityByNameAndType>( \
428  FE->colPtr->get<Composite_Name_And_Type_mi_tag>(), NAME, TYPE); \
429  IT != FE->get_end<FENumeredDofEntityByNameAndType>( \
430  FE->colPtr->get<Composite_Name_And_Type_mi_tag>(), NAME, TYPE); \
431  IT++
432 
433 /** \brief loop over all dofs which are on a particular FE data, field and
434  * entity type \ingroup mofem_loops
435  */
436 #define _IT_GET_FEDATA_BY_TYPE_DOFS_FOR_LOOP_(FE, NAME, TYPE, IT) \
437  auto IT = FE->get_begin<FEDofEntityByNameAndType>( \
438  FE->dataPtr->get<Composite_Name_And_Type_mi_tag>(), NAME, TYPE); \
439  IT != FE->get_end<FEDofEntityByNameAndType>( \
440  FE->dataPtr->get<Composite_Name_And_Type_mi_tag>(), NAME, TYPE); \
441  IT++
442 
443  template <class MULTIINDEX>
444  typename MULTIINDEX::iterator get_begin(const MULTIINDEX &index,
445  const std::string &field_name) const {
446  return index.lower_bound(field_name);
447  }
448  template <class MULTIINDEX>
449  typename MULTIINDEX::iterator get_end(const MULTIINDEX &index,
450  const std::string &field_name) const {
451  return index.upper_bound(field_name);
452  }
453 
454 /** \brief loop over all dofs which are on a particular FE row and field
455  * \ingroup mofem_loops
456  */
457 #define _IT_GET_FEROW_BY_NAME_DOFS_FOR_LOOP_(FE, NAME, IT) \
458  auto IT = FE->get_begin<FENumeredDofEntityByFieldName>( \
459  FE->rowPtr->get<FieldName_mi_tag>(), NAME); \
460  IT != FE->get_end<FENumeredDofEntityByFieldName>( \
461  FE->rowPtr->get<FieldName_mi_tag>(), NAME); \
462  IT++
463 
464 /** \brief loop over all dofs which are on a particular FE column and field
465  * \ingroup mofem_loops
466  */
467 #define _IT_GET_FECOL_BY_NAME_DOFS_FOR_LOOP_(FE, NAME, IT) \
468  auto IT = FE->get_begin<FENumeredDofEntityByFieldName>( \
469  FE->colPtr->get<FieldName_mi_tag>(), NAME); \
470  IT != FE->get_end<FENumeredDofEntityByFieldName>( \
471  FE->colPtr->get<FieldName_mi_tag>(), NAME); \
472  IT++
473 
474 /** \brief loop over all dofs which are on a particular FE data and field
475  * \ingroup mofem_loops
476  */
477 #define _IT_GET_FEDATA_BY_NAME_DOFS_FOR_LOOP_(FE, NAME, IT) \
478  auto IT = FE->get_begin<FEDofEntityByFieldName>( \
479  FE->dataPtr->get<FieldName_mi_tag>(), NAME); \
480  IT != FE->get_end<FEDofEntityByFieldName>( \
481  FE->dataPtr->get<FieldName_mi_tag>(), NAME); \
482  IT++
483 };
484 
486  return numeredEntFiniteElementPtr->getEnt();
487 }
488 
489 /**
490  * \brief Data structure to exchange data between mofem and User Loop Methods on
491  * entities. \ingroup mofem_loops
492  *
493  * It allows to exchange data between MoFEM and user functions. It stores
494  * information about multi-indices.
495  */
496 struct EntityMethod : public BasicMethod {
497 
499  UnknownInterface **iface) const {
501  if (uuid == IDD_MOFEMEntityMethod) {
502  *iface = const_cast<EntityMethod *>(this);
504  }
505  CHKERR query_interface(uuid, iface);
507  }
508 
509  EntityMethod();
510 
511  boost::shared_ptr<Field> fieldPtr;
512  boost::shared_ptr<FieldEntity> entPtr;
513 };
514 
515 /**
516  * \brief Data structure to exchange data between mofem and User Loop Methods on
517  * entities. \ingroup mofem_loops
518  *
519  * It allows to exchange data between MoFEM and user functions. It stores
520  * information about multi-indices.
521  */
522 struct DofMethod : public BasicMethod {
523 
525  UnknownInterface **iface) const {
527  if (uuid == IDD_MOFEMDofMethod) {
528  *iface = const_cast<DofMethod *>(this);
530  }
531 
532  CHKERR query_interface(uuid, iface);
534  }
535 
536  DofMethod();
537 
538  boost::shared_ptr<Field> fieldPtr;
539  boost::shared_ptr<DofEntity> dofPtr;
540  boost::shared_ptr<NumeredDofEntity> dofNumeredPtr;
541 };
542 
543 /// \deprecated name changed use DofMethod insead EntMethod
545 
546 } // namespace MoFEM
547 
548 #endif // __LOOPMETHODS_HPP__
549 
550 /**
551  * \defgroup mofem_loops Loops
552  * \ingroup mofem
553  */
MoFEMErrorCode copyBasicMethod(const BasicMethod &basic)
Copy data from other base method to this base method.
boost::shared_ptr< DofEntity > dofPtr
structure for User Loop Methods on finite elementsIt can be used to calculate stiffness matrices,...
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.
boost::shared_ptr< Field > fieldPtr
KSP ksp
KSP solver.
MULTIINDEX::iterator get_end(const MULTIINDEX &index, const std::string &field_name) const
multi_index_container< boost::shared_ptr< RefElement >, indexed_by< ordered_unique< tag< Ent_mi_tag >, const_mem_fun< RefElement::interface_type_RefEntity, EntityHandle, &RefElement::getRefEnt > >, ordered_non_unique< tag< EntType_mi_tag >, const_mem_fun< RefElement::interface_type_RefEntity, EntityType, &RefElement::getEntType > > > > RefElement_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.
const Field_multiIndex * fieldsPtr
raw pointer to fields container
MoFEMErrorCode query_interface(const MOFEMuuid &uuid, UnknownInterface **iface) const
Definition: LoopMethods.cpp:46
MoFEMErrorCode query_interface(const MOFEMuuid &uuid, UnknownInterface **iface) const
int getLoopSize() const
get loop size
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.
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
Data structure to exchange data between mofem and User Loop Methods on entities.It allows to exchange...
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:507
MULTIINDEX::iterator get_begin(const MULTIINDEX &index, const std::string &field_name) const
const FieldEntity_multiIndex * entitiesPtr
raw pointer to container of field entities
Mat & ts_B
Preconditioner for ts_A.
base class for all interface classes
multi_index_container< boost::shared_ptr< FieldEntity >, indexed_by< ordered_unique< tag< Unique_mi_tag >, member< FieldEntity, UId, &FieldEntity::globalUId > >, ordered_non_unique< tag< FieldName_mi_tag >, const_mem_fun< FieldEntity::interface_type_Field, boost::string_ref, &FieldEntity::getNameRef > >, ordered_non_unique< tag< Ent_mi_tag >, const_mem_fun< FieldEntity, EntityHandle, &FieldEntity::getEnt > >, ordered_non_unique< tag< Composite_Name_And_Ent_mi_tag >, composite_key< FieldEntity, const_mem_fun< FieldEntity::interface_type_Field, boost::string_ref, &FieldEntity::getNameRef >, const_mem_fun< FieldEntity, EntityHandle, &FieldEntity::getEnt > > > > > FieldEntity_multiIndex
MultiIndex container keeps FieldEntity.
boost::shared_ptr< const FieldEntity_vector_view > colFieldEntsPtr
Pointer to finite element field entities column view.
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:483
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:550
std::string feName
Name of finite element.
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.
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...
boost::shared_ptr< Field > fieldPtr
int getNinTheLoop() const
get number of evaluated element in the loop
int rAnk
processor rank
boost::shared_ptr< const FieldEntity_vector_view > rowFieldEntsPtr
Pointer to finite element field entities row view.
static const MOFEMuuid IDD_MOFEMBasicMethod
Definition: LoopMethods.hpp:31
MoFEMErrorCode query_interface(const MOFEMuuid &uuid, UnknownInterface **iface) const
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:514
Data structure to exchange data between mofem and User Loop Methods on entities.It allows to exchange...
keeps basic data about problemThis is low level structure with information about problem,...
virtual MoFEMErrorCode preProcess()
function is run at the beginning of loop
KSPContext ksp_ctx
Context.
boost::shared_ptr< const FEDofEntity_multiIndex > dataPtr
Pointer to finite element data dofs.
std::bitset< BITINTERFACEUID_SIZE > BitIntefaceId
Definition: Types.hpp:55
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 query_interface(const MOFEMuuid &uuid, UnknownInterface **iface) const
boost::shared_ptr< const NumeredEntFiniteElement > numeredEntFiniteElementPtr
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
boost::movelib::unique_ptr< bool > matAssembleSwitch
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
multi_index_container< boost::shared_ptr< DofEntity >, indexed_by< ordered_unique< tag< Unique_mi_tag >, member< DofEntity, UId, &DofEntity::globalUId > >, ordered_non_unique< tag< Composite_Name_And_Ent_And_EntDofIdx_mi_tag >, composite_key< DofEntity, const_mem_fun< DofEntity::interface_type_Field, boost::string_ref, &DofEntity::getNameRef >, const_mem_fun< DofEntity, EntityHandle, &DofEntity::getEnt >, const_mem_fun< DofEntity, DofIdx, &DofEntity::getEntDofIdx > > >, ordered_non_unique< tag< Unique_Ent_mi_tag >, const_mem_fun< DofEntity, const UId &, &DofEntity::getEntGlobalUniqueId > >, ordered_non_unique< tag< FieldName_mi_tag >, const_mem_fun< DofEntity::interface_type_Field, boost::string_ref, &DofEntity::getNameRef > >, ordered_non_unique< tag< Ent_mi_tag >, const_mem_fun< DofEntity, EntityHandle, &DofEntity::getEnt > >, ordered_non_unique< tag< Composite_Name_And_Ent_mi_tag >, composite_key< DofEntity, const_mem_fun< DofEntity::interface_type_Field, boost::string_ref, &DofEntity::getNameRef >, const_mem_fun< DofEntity, EntityHandle, &DofEntity::getEnt > > >, ordered_non_unique< tag< Composite_Name_And_Type_mi_tag >, composite_key< DofEntity, const_mem_fun< DofEntity::interface_type_Field, boost::string_ref, &DofEntity::getNameRef >, const_mem_fun< DofEntity::interface_type_RefEntity, EntityType, &DofEntity::getEntType > > > > > DofEntity_multiIndex
MultiIndex container keeps DofEntity.
DEPRECATED MoFEMErrorCode setTsCtx(TSContext ctx)
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:67
data structure for ksp (linear solver) contextStruture stores context data which are set in functions...
Definition: LoopMethods.hpp:88
MoFEMErrorCode getNodeData(const std::string field_name, VectorDouble &data, const bool reset_dofs=true)
Vec & snes_x
state vector
EntityHandle getFEEntityHandle() const
static constexpr Switches CtxSetF
Definition: LoopMethods.hpp:62
KSPContext
pass information about context of KSP/DM for with finite element is computed
Definition: LoopMethods.hpp:96
virtual ~SnesMethod()=default
DEPRECATED typedef DofMethod EntMethod
DEPRECATED MoFEMErrorCode setSnesCtx(SNESContext ctx)
data structure for TS (time stepping) contextStructure stores context data which are set in functions...
static constexpr Switches CtxSetX_TT
Definition: LoopMethods.hpp:67
#define CHKERR
Inline error check.
Definition: definitions.h:602
multi_index_container< boost::shared_ptr< EntFiniteElement >, indexed_by< ordered_unique< tag< Unique_mi_tag >, member< EntFiniteElement, UId, &EntFiniteElement::globalUId > >, ordered_non_unique< tag< Ent_mi_tag >, const_mem_fun< EntFiniteElement, 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< BitFEId_mi_tag >, const_mem_fun< EntFiniteElement::interface_type_FiniteElement, BitFEId, &EntFiniteElement::getId >, LtBit< BitFEId > >, ordered_non_unique< tag< EntType_mi_tag >, const_mem_fun< EntFiniteElement::interface_type_RefEntity, EntityType, &EntFiniteElement::getEntType > >, 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, EntityHandle, &EntFiniteElement::getEnt > > > > > EntFiniteElement_multiIndex
MultiIndex container for EntFiniteElement.
int nInTheLoop
number currently of processed method
static const MOFEMuuid IDD_MOFEMEntityMethod
Definition: LoopMethods.hpp:34
static constexpr Switches CtxSetA
Definition: LoopMethods.hpp:63
Vec & ts_F
residual vector
boost::shared_ptr< NumeredDofEntity > dofNumeredPtr
SNES snes
snes solver
boost::shared_ptr< FieldEntity > entPtr
virtual ~TSMethod()=default
multi_index_container< boost::shared_ptr< RefEntity >, indexed_by< ordered_unique< tag< Ent_mi_tag >, member< RefEntity::BasicEntity, EntityHandle, &RefEntity::ent > >, ordered_non_unique< tag< Ent_Ent_mi_tag >, const_mem_fun< RefEntity, EntityHandle, &RefEntity::getParentEnt > >, ordered_non_unique< tag< EntType_mi_tag >, const_mem_fun< RefEntity::BasicEntity, EntityType, &RefEntity::getEntType > >, ordered_non_unique< tag< ParentEntType_mi_tag >, const_mem_fun< RefEntity, EntityType, &RefEntity::getParentEntType > >, ordered_non_unique< tag< Composite_EntType_and_ParentEntType_mi_tag >, composite_key< RefEntity, const_mem_fun< RefEntity::BasicEntity, 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::BasicEntity, EntityType, &RefEntity::getEntType >, const_mem_fun< RefEntity, EntityHandle, &RefEntity::getParentEnt > > > > > RefEntity_multiIndex
#define DEPRECATED
Definition: definitions.h:29
MoFEMErrorCode copySnes(const SnesMethod &snes)
Copy snes data.
Definition: LoopMethods.cpp:84
virtual ~KspMethod()=default
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
boost::shared_ptr< const FieldEntity_multiIndex_spaceType_view > dataFieldEntsPtr
Pointer to finite element field entities data view.
PetscReal ts_v
shift for U_t shift for U_t
Mat & snes_A
jacobian matrix
MoFEMErrorCode query_interface(const MOFEMuuid &uuid, UnknownInterface **iface) const
MoFEMErrorCode query_interface(const MOFEMuuid &uuid, UnknownInterface **iface) const
Definition: LoopMethods.cpp:96
MULTIINDEX::iterator get_end(const MULTIINDEX &index, const std::string &field_name, const EntityType type) const
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:413
static const MOFEMuuid IDD_MOFEMTsMethod
Definition: LoopMethods.hpp:30
TSContext ts_ctx
ublas::vector< double, DoubleAllocator > VectorDouble
Definition: Types.hpp:73
static const MOFEMuuid IDD_MOFEMKspMethod
Definition: LoopMethods.hpp:26
static const MOFEMuuid IDD_MOFEMDofMethod
Definition: LoopMethods.hpp:36
Vec & ts_u_tt
second time derivative of state vector
boost::shared_ptr< const FENumeredDofEntity_multiIndex > colPtr
Pointer to finite element columns dofs view.
boost::function< MoFEMErrorCode()> preProcessHook
Hook function for pre-processing.
virtual ~PetscData()=default
static constexpr Switches CtxSetX
Definition: LoopMethods.hpp:65
static const MOFEMuuid IDD_MOFEMFEMethod
Definition: LoopMethods.hpp:33
Data structure to exchange data between mofem and User Loop Methods.It allows to exchange data betwee...
virtual ~BasicMethod()=default
boost::movelib::unique_ptr< bool > vecAssembleSwitch
boost::shared_ptr< const FENumeredDofEntity_multiIndex > rowPtr
Pointer to finite element rows dofs view.
MoFEMErrorCode getNumberOfNodes(int &num_nodes) const
Get number of DOFs on element.
MULTIINDEX::iterator get_begin(const MULTIINDEX &index, const std::string &field_name, const EntityType type) const
static constexpr Switches CtxSetNone
Definition: LoopMethods.hpp:61
static const MOFEMuuid IDD_MOFEMSnesMethod
Definition: LoopMethods.hpp:28