v0.8.23
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 
36 
37 /**
38  * \brief data structure for ksp (linear solver) context
39  * \ingroup mofem_loops
40  *
41  * Struture stores context data which are set in functions run by PETSc SNES
42  * functions.
43  *
44  */
45 struct KspMethod : virtual public UnknownInterface {
46 
48  UnknownInterface **iface) const {
50  if (uuid == IDD_MOFEMKspMethod) {
51  *iface = const_cast<KspMethod *>(this);
53  }
54  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "unknown interface");
55  }
56 
57  /**
58  * \brief pass information about context of KSP/DM for with finite element is
59  * computed
60  */
62 
64  : ksp_ctx(CTX_KSPNONE), ksp(PETSC_NULL), ksp_f(PETSC_NULL),
65  ksp_A(PETSC_NULL), ksp_B(PETSC_NULL) {}
66 
67  virtual ~KspMethod() {}
68 
69  /**
70  * \brief set operator type
71  * @param ctx Context, CTX_SETFUNCTION, CTX_OPERATORS, CTX_KSPNONE
72  * @return error code
73  */
75 
76  /**
77  * \brief set solver
78  * @param ksp solver
79  * @return error code
80  */
82 
83  /**
84  * \brief copy data form another method
85  * @param ksp ksp method
86  * @return error code
87  */
89 
90  KSPContext ksp_ctx; ///< Context
91  KSP ksp; ///< KSP solver
92  Vec ksp_f; ///< the right hand side vector
93  Mat ksp_A; ///< matrix
94  Mat ksp_B; ///< preconditioner matrix
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 public UnknownInterface {
106 
108  UnknownInterface **iface) const {
109  if (uuid == IDD_MOFEMSnesMethod) {
110  *iface = const_cast<SnesMethod *>(this);
112  }
113  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "unknown interface");
114  }
115 
117 
120  : snes_ctx(CTX_SNESNONE), snes_x(PETSC_NULL), snes_f(PETSC_NULL),
121  snes_A(PETSC_NULL), snes_B(PETSC_NULL) {}
122 
123  virtual ~SnesMethod() {}
124 
125  /**
126  * \brief Set SNES context
127  */
129 
130  /**
131  * \brief Set SNES instance
132  */
134 
135  /**
136  * \brief Copy snes data
137  */
139 
140  SNES snes;
143 };
144 
145 /**
146  * \brief data structure for TS (time stepping) context
147  * \ingroup mofem_loops
148  *
149  * Structure stores context data which are set in functions run by PETSc Time
150  * Stepping functions.
151  */
152 struct TSMethod : virtual public UnknownInterface {
153 
155  UnknownInterface **iface) const {
156  if (uuid == IDD_MOFEMTsMethod) {
157  *iface = const_cast<TSMethod *>(this);
159  }
160  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "unknown interface");
161  }
162 
163  enum TSContext {
170  };
171 
174  : ts_ctx(CTX_TSNONE), ts_u(PETSC_NULL), ts_u_t(PETSC_NULL),
175  ts_F(PETSC_NULL), ts_A(PETSC_NULL), ts_B(PETSC_NULL), ts_step(-1),
176  ts_a(0), ts_t(0) {}
177 
178  virtual ~TSMethod() {}
179 
180  /// \brief Set Ts context
181  MoFEMErrorCode setTsCtx(const TSContext &ctx);
182 
183  /// \brief Copy TS solver data
185 
186  /// \brief Set TS solver
187  MoFEMErrorCode setTs(TS _ts);
188 
189  TS ts;
190  Vec ts_u, ts_u_t, ts_F;
191  Mat ts_A, ts_B;
192 
193  PetscInt ts_step;
194  PetscReal ts_a, ts_t;
195 };
196 
197 /**
198  * \brief Data structure to exchange data between mofem and User Loop Methods.
199  * \ingroup mofem_loops
200  *
201  * It allows to exchange data between MoFEM and user functions. It stores
202  * information about multi-indices.
203  *
204  */
206 
208  UnknownInterface **iface) const {
209  if (uuid == IDD_MOFEMBasicMethod) {
210  *iface = const_cast<BasicMethod *>(this);
212  }
213  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "unknown interface");
214  }
215 
216  BasicMethod();
217  virtual ~BasicMethod(){};
218 
219  /**
220  * @brief number currently of processed method
221  */
222  int nInTheLoop;
223 
224  /**
225  * @brief local number oe methods to process
226  */
227  int loopSize;
228 
229  /** \brief get number of evaluated element in the loop
230  */
231  inline int getNinTheLoop() const { return nInTheLoop; }
232 
233  /** \brief get loop size
234  */
235  inline int getLoopSize() const { return loopSize; }
236 
237  int rAnk; ///< processor rank
238 
239  int sIze; ///< number of processors in communicator
240 
242  *refinedEntitiesPtr; ///< container of mofem dof entities
243 
245  *refinedFiniteElementsPtr; ///< container of mofem finite element entities
246 
247  const Problem *problemPtr; ///< raw pointer to problem
248 
249  const Field_multiIndex *fieldsPtr; ///< raw pointer to fields container
250 
252  *entitiesPtr; ///< raw pointer to container of field entities
253 
254  const DofEntity_multiIndex *dofsPtr; ///< raw pointer container of dofs
255 
257  *finiteElementsPtr; ///< raw pointer to container finite elements
258 
260  *finiteElementsEntitiesPtr; ///< raw pointer to container finite elements
261  ///< entities
262 
264  *adjacenciesPtr; ///< raw pointer to container to adjacencies between dofs
265  ///< and finite elements
266 
267  /**
268  * @brief Copy data from other base method to this base method
269  *
270  * @param basic
271  * @return MoFEMErrorCode
272  */
274 
275  /**
276  * @brief Hook function for pre-processing
277  */
278  boost::function<MoFEMErrorCode()> preProcessHook;
279 
280  /**
281  * @brief Hook function for post-processing
282  */
283  boost::function<MoFEMErrorCode()> postProcessHook;
284 
285  /**
286  * @brief Hook function for operator
287  */
288  boost::function<MoFEMErrorCode()> operatorHook;
289 
290  /** \brief function is run at the beginning of loop
291  *
292  * It is used to zeroing matrices and vectors, calculation of shape
293  * functions on reference element, preprocessing boundary conditions, etc.
294  */
295  virtual MoFEMErrorCode preProcess();
296 
297  /** \brief function is run for every finite element
298  *
299  * It is used to calculate element local matrices and assembly. It can be
300  * used for post-processing.
301  */
302  virtual MoFEMErrorCode operator()();
303 
304  /** \brief function is run at the end of loop
305  *
306  * It is used to assembly matrices and vectors, calculating global variables,
307  * f.e. total internal energy, ect.
308  *
309  * Iterating over dofs:
310  * Example1 iterating over dofs in row by name of the field
311  * for(_IT_GET_FEROW_BY_NAME_DOFS_FOR_LOOP_(this,"DISPLACEMENT",it)) { ... }
312  *
313  *
314  */
315  virtual MoFEMErrorCode postProcess();
316 };
317 
318 /**
319  * \brief structure for User Loop Methods on finite elements
320  * \ingroup mofem_loops
321  *
322  * It can be used to calculate stiffness matrices, residuals, load vectors etc.
323  * It is low level class however in some class users looking for speed and
324  * efficiency, can use it directly.
325  *
326  * This class is used with Interface::loop_finite_elements, where
327  * user overloaded operator FEMethod::operator() is executed for each element in
328  * the problem. Class have to additional methods which are overloaded by user,
329  * FEMethod::preProcess() and FEMethod::postProcess() executed at beginning and
330  * end of the loop over problem elements, respectively.
331  *
332  */
333 struct FEMethod : public BasicMethod {
334 
336  UnknownInterface **iface) const {
338  if (uuid == IDD_MOFEMFEMethod) {
339  *iface = const_cast<FEMethod *>(this);
341  }
342 
343  ierr = query_interface(uuid, iface);
344  CHKERRG(ierr);
346  }
347 
348  FEMethod();
349 
350  std::string feName; ///< Name of finite element
351 
352  boost::shared_ptr<const NumeredEntFiniteElement>
353  numeredEntFiniteElementPtr; ///< Pointer to finite element database
354  ///< structure
355  boost::shared_ptr<const FENumeredDofEntity_multiIndex>
356  rowPtr; ///< Pointer to finite element rows dofs view
357  boost::shared_ptr<const FENumeredDofEntity_multiIndex>
358  colPtr; ///< Pointer to finite element columns dofs view
359  boost::shared_ptr<const FEDofEntity_multiIndex>
360  dataPtr; ///< Pointer to finite element data dofs
361 
362  boost::shared_ptr<const FieldEntity_vector_view>
363  rowFieldEntsPtr; ///< Pointer to finite element field entities row view
364  boost::shared_ptr<const FieldEntity_vector_view>
365  colFieldEntsPtr; ///< Pointer to finite element field entities column view
366  boost::shared_ptr<const FieldEntity_multiIndex_spaceType_view>
367  dataFieldEntsPtr; ///< Pointer to finite element field entities data view
368 
369 /** \brief loop over all dofs which are on a particular FE row
370  * \ingroup mofem_loops
371  */
372 #define _IT_GET_FEROW_DOFS_FOR_LOOP_(FE, IT) \
373  auto IT = FE->rowPtr->begin(); \
374  IT != FE->rowPtr->end(); \
375  IT++
376 
377 /** \brief loop over all dofs which are on a particular FE column
378  * \ingroup mofem_loops
379  */
380 #define _IT_GET_FECOL_DOFS_FOR_LOOP_(FE, IT) \
381  auto IT = FE->colPtr->begin(); \
382  IT != FE->colPtr->end(); \
383  IT++
384 
385 /** \brief loop over all dofs which are on a particular FE data
386  * \ingroup mofem_loops
387  */
388 #define _IT_GET_FEDATA_DOFS_FOR_LOOP_(FE, IT) \
389  auto IT = FE->dataPtr->begin(); \
390  IT != FE->dataPtr->end(); \
391  IT++
392 
393  template <class MULTIINDEX>
394  typename MULTIINDEX::iterator
395  get_begin(const MULTIINDEX &index, const std::string &field_name,
396  const EntityType type, const int side_number) const {
397  return index.lower_bound(boost::make_tuple(field_name, type, side_number));
398  }
399 
400  template <class MULTIINDEX>
401  typename MULTIINDEX::iterator
402  get_end(const MULTIINDEX &index, const std::string &field_name,
403  const EntityType type, const int side_number) const {
404  return index.upper_bound(boost::make_tuple(field_name, type, side_number));
405  }
406 
407  /** \brief loop over all dofs which are on a particular FE row, field,
408  * entity type and canonical side number \ingroup mofem_loops
409  *
410  * \param FE finite elements
411  * \param Name field name
412  * \param Type moab entity type (MBVERTEX, MBEDGE etc)
413  * \param Side side canonical number
414  * \param IT the interator in use
415  */
416 #define _IT_GET_FEROW_BY_SIDE_DOFS_FOR_LOOP_(FE, NAME, TYPE, SIDE, IT) \
417  FENumeredDofEntity_multiIndex::index<Composite_mi_tag>::type::iterator IT = \
418  FE->get_begin< \
419  FENumeredDofEntity_multiIndex::index<Composite_mi_tag>::type>( \
420  FE->rowPtr->get<Composite_mi_tag>(), NAME, TYPE, SIDE); \
421  IT != FE->get_end< \
422  FENumeredDofEntity_multiIndex::index<Composite_mi_tag>::type>( \
423  FE->rowPtr->get<Composite_mi_tag>(), NAME, TYPE, SIDE); \
424  IT++
425 
426 /** \brief loop over all dofs which are on a particular FE column, field, entity
427  * type and canonical side number \ingroup mofem_loops
428  */
429 #define _IT_GET_FECOL_BY_SIDE_DOFS_FOR_LOOP_(FE, NAME, TYPE, SIDE, IT) \
430  FENumeredDofEntity_multiIndex::index<Composite_mi_tag>::type::iterator IT = \
431  FE->get_begin< \
432  FENumeredDofEntity_multiIndex::index<Composite_mi_tag>::type>( \
433  FE->colPtr->get<Composite_mi_tag>(), NAME, TYPE, SIDE); \
434  IT != FE->get_end< \
435  FENumeredDofEntity_multiIndex::index<Composite_mi_tag>::type>( \
436  FE->colPtr->get<Composite_mi_tag>(), NAME, TYPE, SIDE); \
437  IT++
438 
439 /** \brief loop over all dofs which are on a particular FE data, field, entity
440  * type and canonical side number \ingroup mofem_loops
441  */
442 #define _IT_GET_FEDATA_BY_SIDE_DOFS_FOR_LOOP_(FE, NAME, TYPE, SIDE, IT) \
443  FEDofEntity_multiIndex::index<Composite_mi_tag>::type::iterator IT = \
444  FE->get_begin<FEDofEntity_multiIndex::index<Composite_mi_tag>::type>( \
445  FE->dataPtr->get<Composite_mi_tag>(), NAME, TYPE, SIDE); \
446  IT != FE->get_end<FEDofEntity_multiIndex::index<Composite_mi_tag>::type>( \
447  FE->dataPtr->get<Composite_mi_tag>(), NAME, TYPE, SIDE); \
448  IT++
449 
450  template <class MULTIINDEX>
451  typename MULTIINDEX::iterator get_begin(const MULTIINDEX &index,
452  const std::string &field_name,
453  const EntityType type) const {
454  return index.lower_bound(boost::make_tuple(field_name, type));
455  }
456  template <class MULTIINDEX>
457  typename MULTIINDEX::iterator get_end(const MULTIINDEX &index,
458  const std::string &field_name,
459  const EntityType type) const {
460  return index.upper_bound(boost::make_tuple(field_name, type));
461  }
462 
463 /** \brief loop over all dofs which are on a particular FE row, field and entity
464  * type \ingroup mofem_loops
465  */
466 #define _IT_GET_FEROW_BY_TYPE_DOFS_FOR_LOOP_(FE, NAME, TYPE, IT) \
467  auto IT = FE->get_begin<FENumeredDofEntityByNameAndType>( \
468  FE->rowPtr->get<Composite_Name_And_Type_mi_tag>(), NAME, TYPE); \
469  IT != FE->get_end<FENumeredDofEntityByNameAndType>( \
470  FE->rowPtr->get<Composite_Name_And_Type_mi_tag>(), NAME, TYPE); \
471  IT++
472 
473 /** \brief loop over all dofs which are on a particular FE column, field and
474  * entity type \ingroup mofem_loops
475  */
476 #define _IT_GET_FECOL_BY_TYPE_DOFS_FOR_LOOP_(FE, NAME, TYPE, IT) \
477  auto IT = FE->get_begin<FENumeredDofEntityByNameAndType>( \
478  FE->colPtr->get<Composite_Name_And_Type_mi_tag>(), NAME, TYPE); \
479  IT != FE->get_end<FENumeredDofEntityByNameAndType>( \
480  FE->colPtr->get<Composite_Name_And_Type_mi_tag>(), NAME, TYPE); \
481  IT++
482 
483 /** \brief loop over all dofs which are on a particular FE data, field and
484  * entity type \ingroup mofem_loops
485  */
486 #define _IT_GET_FEDATA_BY_TYPE_DOFS_FOR_LOOP_(FE, NAME, TYPE, IT) \
487  auto IT = FE->get_begin<FEDofEntityByNameAndType>( \
488  FE->dataPtr->get<Composite_Name_And_Type_mi_tag>(), NAME, TYPE); \
489  IT != FE->get_end<FEDofEntityByNameAndType>( \
490  FE->dataPtr->get<Composite_Name_And_Type_mi_tag>(), NAME, TYPE); \
491  IT++
492 
493  template <class MULTIINDEX>
494  typename MULTIINDEX::iterator get_begin(const MULTIINDEX &index,
495  const std::string &field_name) const {
496  return index.lower_bound(field_name);
497  }
498  template <class MULTIINDEX>
499  typename MULTIINDEX::iterator get_end(const MULTIINDEX &index,
500  const std::string &field_name) const {
501  return index.upper_bound(field_name);
502  }
503 
504 /** \brief loop over all dofs which are on a particular FE row and field
505  * \ingroup mofem_loops
506  */
507 #define _IT_GET_FEROW_BY_NAME_DOFS_FOR_LOOP_(FE, NAME, IT) \
508  auto IT = FE->get_begin<FENumeredDofEntityByFieldName>( \
509  FE->rowPtr->get<FieldName_mi_tag>(), NAME); \
510  IT != FE->get_end<FENumeredDofEntityByFieldName>( \
511  FE->rowPtr->get<FieldName_mi_tag>(), NAME); \
512  IT++
513 
514 /** \brief loop over all dofs which are on a particular FE column and field
515  * \ingroup mofem_loops
516  */
517 #define _IT_GET_FECOL_BY_NAME_DOFS_FOR_LOOP_(FE, NAME, IT) \
518  auto IT = FE->get_begin<FENumeredDofEntityByFieldName>( \
519  FE->colPtr->get<FieldName_mi_tag>(), NAME); \
520  IT != FE->get_end<FENumeredDofEntityByFieldName>( \
521  FE->colPtr->get<FieldName_mi_tag>(), NAME); \
522  IT++
523 
524 /** \brief loop over all dofs which are on a particular FE data and field
525  * \ingroup mofem_loops
526  */
527 #define _IT_GET_FEDATA_BY_NAME_DOFS_FOR_LOOP_(FE, NAME, IT) \
528  auto IT = FE->get_begin<FEDofEntityByFieldName>( \
529  FE->dataPtr->get<FieldName_mi_tag>(), NAME); \
530  IT != FE->get_end<FEDofEntityByFieldName>( \
531  FE->dataPtr->get<FieldName_mi_tag>(), NAME); \
532  IT++
533 
534  template <class MULTIINDEX>
535  typename MULTIINDEX::iterator get_begin(const MULTIINDEX &index,
536  const EntityHandle ent) const {
537  return index.lower_bound(ent);
538  }
539  template <class MULTIINDEX>
540  typename MULTIINDEX::iterator get_end(const MULTIINDEX &index,
541  const EntityHandle ent) const {
542  return index.upper_bound(ent);
543  }
544 
545 /** \brief loop over all dofs which are on a particular FE row and given element
546  * entity (handle from moab) \ingroup mofem_loops
547  */
548 #define _IT_GET_FEROW_DOFS_BY_ENT_FOR_LOOP_(FE, ENT, IT) \
549  auto IT = FE->get_begin<FENumeredDofEntityByEnt>( \
550  FE->rowPtr->get<Ent_mi_tag>(), ENT); \
551  IT != FE->get_end<FENumeredDofEntityByEnt>(FE->rowPtr->get<Ent_mi_tag>(), \
552  ENT); \
553  IT++
554 
555 /** \brief loop over all dofs which are on a particular FE column and given
556  * element entity (handle from moab) \ingroup mofem_loops
557  */
558 #define _IT_GET_FECOL_DOFS_BY_ENT_FOR_LOOP_(FE, ENT, IT) \
559  auto IT = FE->get_begin<FENumeredDofEntityByEnt>( \
560  FE->colPtr->get<Ent_mi_tag>(), ENT); \
561  IT != FE->get_end<FENumeredDofEntityByEnt>(FE->colPtr->get<Ent_mi_tag>(), \
562  ENT); \
563  IT++
564 
565 /** \brief loop over all dofs which are on a particular FE data and given
566  * element entity (handle from moab) \ingroup mofem_loops
567  */
568 #define _IT_GET_FEDATA_DOFS_BY_ENT_FOR_LOOP_(FE, ENT, IT) \
569  auto IT = FE->get_begin<FEDofEntity_multiIndex::index<Ent_mi_tag>::type>( \
570  FE->dataPtr->get<Ent_mi_tag>(), ENT); \
571  IT != FE->get_end<FEDofEntity_multiIndex::index<Ent_mi_tag>::type>( \
572  FE->dataPtr->get<Ent_mi_tag>(), ENT); \
573  IT++
574 
575  template <class MULTIINDEX>
576  typename MULTIINDEX::iterator get_begin(const MULTIINDEX &index,
577  const std::string &field_name,
578  const EntityHandle ent) const {
579  return index.lower_bound(boost::make_tuple(field_name, ent));
580  }
581  template <class MULTIINDEX>
582  typename MULTIINDEX::iterator get_end(const MULTIINDEX &index,
583  const std::string &field_name,
584  const EntityHandle ent) const {
585  return index.upper_bound(boost::make_tuple(field_name, ent));
586  }
587 
588 /** \brief loop over all dofs which are on a particular FE row, field and given
589  * element entity (handle from moab) \ingroup mofem_loops
590  */
591 #define _IT_GET_FEROW_DOFS_BY_NAME_AND_ENT_FOR_LOOP_(FE, NAME, ENT, IT) \
592  auto IT = FE->get_begin<FENumeredDofEntityByNameAndEnt>( \
593  FE->rowPtr->get<Composite_Name_And_Ent_mi_tag>(), NAME, ENT); \
594  IT != FE->get_end<FENumeredDofEntityByNameAndEnt>( \
595  FE->rowPtr->get<Composite_Name_And_Ent_mi_tag>(), NAME, ENT); \
596  IT++
597 
598 /** \brief loop over all dofs which are on a particular FE column, field and
599  * given element entity (handle from moab) \ingroup mofem_loops
600  */
601 #define _IT_GET_FECOL_DOFS_BY_NAME_AND_ENT_FOR_LOOP_(FE, NAME, ENT, IT) \
602  auto IT = FE->get_begin<FENumeredDofEntityByNameAndEnt>( \
603  FE->colPtr->get<Composite_Name_And_Ent_mi_tag>(), NAME, ENT); \
604  IT != FE->get_end<FENumeredDofEntityByNameAndEnt>( \
605  FE->colPtr->get<Composite_Name_And_Ent_mi_tag>(), NAME, ENT); \
606  IT++
607 
608 /** \brief loop over all dofs which are on a particular FE data, field and given
609  * element entity (handle from moab) \ingroup mofem_loops
610  */
611 #define _IT_GET_FEDATA_DOFS_BY_NAME_AND_ENT_FOR_LOOP_(FE, NAME, ENT, IT) \
612  auto IT = FE->get_begin<FEDofEntityByNameAndEnt>( \
613  FE->dataPtr->get<Composite_Name_And_Ent_mi_tag>(), NAME, ENT); \
614  IT != FE->get_end<FEDofEntityByNameAndEnt>( \
615  FE->dataPtr->get<Composite_Name_And_Ent_mi_tag>(), NAME, ENT); \
616  IT++
617 };
618 
619 /**
620  * \brief Data structure to exchange data between mofem and User Loop Methods on
621  * entities. \ingroup mofem_loops
622  *
623  * It allows to exchange data between MoFEM and user functions. It stores
624  * information about multi-indices.
625  */
626 struct EntityMethod : public BasicMethod {
627 
629  UnknownInterface **iface) const {
631  if (uuid == IDD_MOFEMEntityMethod) {
632  *iface = const_cast<EntityMethod *>(this);
634  }
635  CHKERR query_interface(uuid, iface);
637  }
638 
639  EntityMethod();
640 
641  boost::shared_ptr<Field> fieldPtr;
642  boost::shared_ptr<FieldEntity> entPtr;
643 };
644 
645 /**
646  * \brief Data structure to exchange data between mofem and User Loop Methods on
647  * entities. \ingroup mofem_loops
648  *
649  * It allows to exchange data between MoFEM and user functions. It stores
650  * information about multi-indices.
651  */
652 struct DofMethod : public BasicMethod {
653 
655  UnknownInterface **iface) const {
657  if (uuid == IDD_MOFEMDofMethod) {
658  *iface = const_cast<DofMethod *>(this);
660  }
661 
662  CHKERR query_interface(uuid, iface);
664  }
665 
666  DofMethod();
667 
668  boost::shared_ptr<Field> fieldPtr;
669  boost::shared_ptr<DofEntity> dofPtr;
670  boost::shared_ptr<NumeredDofEntity> dofNumeredPtr;
671 };
672 
673 /// \deprecated name changed use DofMethod insead EntMethod
675 
676 } // namespace MoFEM
677 
678 #endif // __LOOPMETHODS_HPP__
679 
680 /**
681  * \defgroup mofem_loops Loops
682  * \ingroup mofem
683  */
MULTIINDEX::iterator get_end(const MULTIINDEX &index, const EntityHandle ent) const
MoFEMErrorCode copyBasicMethod(const BasicMethod &basic)
Copy data from other base method to this base method.
Definition: LoopMethods.cpp:99
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
MoFEM interface unique ID.
boost::shared_ptr< Field > fieldPtr
KSP ksp
KSP solver.
Definition: LoopMethods.hpp:91
MULTIINDEX::iterator get_begin(const MULTIINDEX &index, const std::string &field_name, const EntityType type, const int side_number) const
Mat ksp_B
preconditioner matrix
Definition: LoopMethods.hpp:94
MULTIINDEX::iterator get_end(const MULTIINDEX &index, const std::string &field_name) const
MULTIINDEX::iterator get_begin(const MULTIINDEX &index, const std::string &field_name, const EntityHandle ent) 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.hpp:47
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.
const FiniteElement_multiIndex * finiteElementsPtr
raw pointer to container finite elements
MoFEMErrorCode setSnes(SNES snes)
Set SNES instance.
Definition: LoopMethods.cpp:48
int loopSize
local number oe methods to process
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:500
MULTIINDEX::iterator get_begin(const MULTIINDEX &index, const std::string &field_name) const
const FieldEntity_multiIndex * entitiesPtr
raw pointer to container of field entities
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:476
virtual ~BasicMethod()
MULTIINDEX::iterator get_end(const MULTIINDEX &index, const std::string &field_name, const EntityType type, const int side_number) const
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:543
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
virtual ~SnesMethod()
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:29
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:507
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.
Definition: LoopMethods.hpp:90
boost::shared_ptr< const FEDofEntity_multiIndex > dataPtr
Pointer to finite element data dofs.
std::bitset< BITINTERFACEUID_SIZE > BitIntefaceId
Definition: Types.hpp:56
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:21
MULTIINDEX::iterator get_end(const MULTIINDEX &index, const std::string &field_name, const EntityHandle ent) const
boost::function< MoFEMErrorCode()> postProcessHook
Hook function for post-processing.
const Problem * problemPtr
raw pointer to problem
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:32
virtual MoFEMErrorCode operator()()
function is run for every finite element
virtual MoFEMErrorCode postProcess()
function is run at the end of loop
const RefEntity_multiIndex * refinedEntitiesPtr
container of mofem dof entities
MoFEMErrorCode copyTs(const TSMethod &ts)
Copy TS solver data.
Definition: LoopMethods.cpp:75
MULTIINDEX::iterator get_begin(const MULTIINDEX &index, const EntityHandle ent) const
MoFEMErrorCode setSnesCtx(const SNESContext &ctx)
Set SNES context.
Definition: LoopMethods.cpp:43
boost::function< MoFEMErrorCode()> operatorHook
Hook function for operator.
int sIze
number of processors in communicator
MoFEMErrorCode setTs(TS _ts)
Set TS solver.
Definition: LoopMethods.cpp:70
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.
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:66
data structure for ksp (linear solver) contextStruture stores context data which are set in functions...
Definition: LoopMethods.hpp:45
MoFEMErrorCode setKsp(KSP ksp)
set solver
Definition: LoopMethods.cpp:27
KSPContext
pass information about context of KSP/DM for with finite element is computed
Definition: LoopMethods.hpp:61
DEPRECATED typedef DofMethod EntMethod
data structure for TS (time stepping) contextStructure stores context data which are set in functions...
virtual ~KspMethod()
Definition: LoopMethods.hpp:67
#define CHKERR
Inline error check.
Definition: definitions.h:595
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
MoFEMErrorCode setKspCtx(const KSPContext &ctx)
set operator type
Definition: LoopMethods.cpp:22
static const MOFEMuuid IDD_MOFEMEntityMethod
Definition: LoopMethods.hpp:32
boost::shared_ptr< NumeredDofEntity > dofNumeredPtr
MoFEMErrorCode setTsCtx(const TSContext &ctx)
Set Ts context.
Definition: LoopMethods.cpp:65
boost::shared_ptr< FieldEntity > entPtr
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
virtual ~TSMethod()
MoFEMErrorCode copySnes(const SnesMethod &snes)
Copy snes data.
Definition: LoopMethods.cpp:53
MoFEMErrorCode query_interface(const MOFEMuuid &uuid, UnknownInterface **iface) const
SNESContext snes_ctx
const EntFiniteElement_multiIndex * finiteElementsEntitiesPtr
boost::shared_ptr< const FieldEntity_multiIndex_spaceType_view > dataFieldEntsPtr
Pointer to finite element field entities data view.
MoFEMErrorCode query_interface(const MOFEMuuid &uuid, UnknownInterface **iface) const
MoFEMErrorCode query_interface(const MOFEMuuid &uuid, UnknownInterface **iface) const
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:406
Vec ksp_f
the right hand side vector
Definition: LoopMethods.hpp:92
static const MOFEMuuid IDD_MOFEMTsMethod
Definition: LoopMethods.hpp:28
TSContext ts_ctx
static const MOFEMuuid IDD_MOFEMKspMethod
Definition: LoopMethods.hpp:24
static const MOFEMuuid IDD_MOFEMDofMethod
Definition: LoopMethods.hpp:34
boost::shared_ptr< const FENumeredDofEntity_multiIndex > colPtr
Pointer to finite element columns dofs view.
boost::function< MoFEMErrorCode()> preProcessHook
Hook function for pre-processing.
static const MOFEMuuid IDD_MOFEMFEMethod
Definition: LoopMethods.hpp:31
Data structure to exchange data between mofem and User Loop Methods.It allows to exchange data betwee...
boost::shared_ptr< const FENumeredDofEntity_multiIndex > rowPtr
Pointer to finite element rows dofs view.
MULTIINDEX::iterator get_begin(const MULTIINDEX &index, const std::string &field_name, const EntityType type) const
static const MOFEMuuid IDD_MOFEMSnesMethod
Definition: LoopMethods.hpp:26