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  boost::movelib::unique_ptr<bool> vecAssembleSwitch;
318  boost::movelib::unique_ptr<bool> matAssembleSwitch;
319 
320 };
321 
322 /**
323  * \brief structure for User Loop Methods on finite elements
324  * \ingroup mofem_loops
325  *
326  * It can be used to calculate stiffness matrices, residuals, load vectors etc.
327  * It is low level class however in some class users looking for speed and
328  * efficiency, can use it directly.
329  *
330  * This class is used with Interface::loop_finite_elements, where
331  * user overloaded operator FEMethod::operator() is executed for each element in
332  * the problem. Class have to additional methods which are overloaded by user,
333  * FEMethod::preProcess() and FEMethod::postProcess() executed at beginning and
334  * end of the loop over problem elements, respectively.
335  *
336  */
337 struct FEMethod : public BasicMethod {
338 
340  UnknownInterface **iface) const {
342  if (uuid == IDD_MOFEMFEMethod) {
343  *iface = const_cast<FEMethod *>(this);
345  }
346 
347  ierr = query_interface(uuid, iface);
348  CHKERRG(ierr);
350  }
351 
352  FEMethod();
353 
354  std::string feName; ///< Name of finite element
355 
356  boost::shared_ptr<const NumeredEntFiniteElement>
357  numeredEntFiniteElementPtr; ///< Pointer to finite element database
358  ///< structure
359  boost::shared_ptr<const FENumeredDofEntity_multiIndex>
360  rowPtr; ///< Pointer to finite element rows dofs view
361  boost::shared_ptr<const FENumeredDofEntity_multiIndex>
362  colPtr; ///< Pointer to finite element columns dofs view
363  boost::shared_ptr<const FEDofEntity_multiIndex>
364  dataPtr; ///< Pointer to finite element data dofs
365 
366  boost::shared_ptr<const FieldEntity_vector_view>
367  rowFieldEntsPtr; ///< Pointer to finite element field entities row view
368  boost::shared_ptr<const FieldEntity_vector_view>
369  colFieldEntsPtr; ///< Pointer to finite element field entities column view
370  boost::shared_ptr<const FieldEntity_multiIndex_spaceType_view>
371  dataFieldEntsPtr; ///< Pointer to finite element field entities data view
372 
373 /** \brief loop over all dofs which are on a particular FE row
374  * \ingroup mofem_loops
375  */
376 #define _IT_GET_FEROW_DOFS_FOR_LOOP_(FE, IT) \
377  auto IT = FE->rowPtr->begin(); \
378  IT != FE->rowPtr->end(); \
379  IT++
380 
381 /** \brief loop over all dofs which are on a particular FE column
382  * \ingroup mofem_loops
383  */
384 #define _IT_GET_FECOL_DOFS_FOR_LOOP_(FE, IT) \
385  auto IT = FE->colPtr->begin(); \
386  IT != FE->colPtr->end(); \
387  IT++
388 
389 /** \brief loop over all dofs which are on a particular FE data
390  * \ingroup mofem_loops
391  */
392 #define _IT_GET_FEDATA_DOFS_FOR_LOOP_(FE, IT) \
393  auto IT = FE->dataPtr->begin(); \
394  IT != FE->dataPtr->end(); \
395  IT++
396 
397  template <class MULTIINDEX>
398  typename MULTIINDEX::iterator
399  get_begin(const MULTIINDEX &index, const std::string &field_name,
400  const EntityType type, const int side_number) const {
401  return index.lower_bound(boost::make_tuple(field_name, type, side_number));
402  }
403 
404  template <class MULTIINDEX>
405  typename MULTIINDEX::iterator
406  get_end(const MULTIINDEX &index, const std::string &field_name,
407  const EntityType type, const int side_number) const {
408  return index.upper_bound(boost::make_tuple(field_name, type, side_number));
409  }
410 
411  /** \brief loop over all dofs which are on a particular FE row, field,
412  * entity type and canonical side number \ingroup mofem_loops
413  *
414  * \param FE finite elements
415  * \param Name field name
416  * \param Type moab entity type (MBVERTEX, MBEDGE etc)
417  * \param Side side canonical number
418  * \param IT the interator in use
419  */
420 #define _IT_GET_FEROW_BY_SIDE_DOFS_FOR_LOOP_(FE, NAME, TYPE, SIDE, IT) \
421  FENumeredDofEntity_multiIndex::index<Composite_mi_tag>::type::iterator IT = \
422  FE->get_begin< \
423  FENumeredDofEntity_multiIndex::index<Composite_mi_tag>::type>( \
424  FE->rowPtr->get<Composite_mi_tag>(), NAME, TYPE, SIDE); \
425  IT != FE->get_end< \
426  FENumeredDofEntity_multiIndex::index<Composite_mi_tag>::type>( \
427  FE->rowPtr->get<Composite_mi_tag>(), NAME, TYPE, SIDE); \
428  IT++
429 
430 /** \brief loop over all dofs which are on a particular FE column, field, entity
431  * type and canonical side number \ingroup mofem_loops
432  */
433 #define _IT_GET_FECOL_BY_SIDE_DOFS_FOR_LOOP_(FE, NAME, TYPE, SIDE, IT) \
434  FENumeredDofEntity_multiIndex::index<Composite_mi_tag>::type::iterator IT = \
435  FE->get_begin< \
436  FENumeredDofEntity_multiIndex::index<Composite_mi_tag>::type>( \
437  FE->colPtr->get<Composite_mi_tag>(), NAME, TYPE, SIDE); \
438  IT != FE->get_end< \
439  FENumeredDofEntity_multiIndex::index<Composite_mi_tag>::type>( \
440  FE->colPtr->get<Composite_mi_tag>(), NAME, TYPE, SIDE); \
441  IT++
442 
443 /** \brief loop over all dofs which are on a particular FE data, field, entity
444  * type and canonical side number \ingroup mofem_loops
445  */
446 #define _IT_GET_FEDATA_BY_SIDE_DOFS_FOR_LOOP_(FE, NAME, TYPE, SIDE, IT) \
447  FEDofEntity_multiIndex::index<Composite_mi_tag>::type::iterator IT = \
448  FE->get_begin<FEDofEntity_multiIndex::index<Composite_mi_tag>::type>( \
449  FE->dataPtr->get<Composite_mi_tag>(), NAME, TYPE, SIDE); \
450  IT != FE->get_end<FEDofEntity_multiIndex::index<Composite_mi_tag>::type>( \
451  FE->dataPtr->get<Composite_mi_tag>(), NAME, TYPE, SIDE); \
452  IT++
453 
454  template <class MULTIINDEX>
455  typename MULTIINDEX::iterator get_begin(const MULTIINDEX &index,
456  const std::string &field_name,
457  const EntityType type) const {
458  return index.lower_bound(boost::make_tuple(field_name, type));
459  }
460  template <class MULTIINDEX>
461  typename MULTIINDEX::iterator get_end(const MULTIINDEX &index,
462  const std::string &field_name,
463  const EntityType type) const {
464  return index.upper_bound(boost::make_tuple(field_name, type));
465  }
466 
467 /** \brief loop over all dofs which are on a particular FE row, field and entity
468  * type \ingroup mofem_loops
469  */
470 #define _IT_GET_FEROW_BY_TYPE_DOFS_FOR_LOOP_(FE, NAME, TYPE, IT) \
471  auto IT = FE->get_begin<FENumeredDofEntityByNameAndType>( \
472  FE->rowPtr->get<Composite_Name_And_Type_mi_tag>(), NAME, TYPE); \
473  IT != FE->get_end<FENumeredDofEntityByNameAndType>( \
474  FE->rowPtr->get<Composite_Name_And_Type_mi_tag>(), NAME, TYPE); \
475  IT++
476 
477 /** \brief loop over all dofs which are on a particular FE column, field and
478  * entity type \ingroup mofem_loops
479  */
480 #define _IT_GET_FECOL_BY_TYPE_DOFS_FOR_LOOP_(FE, NAME, TYPE, IT) \
481  auto IT = FE->get_begin<FENumeredDofEntityByNameAndType>( \
482  FE->colPtr->get<Composite_Name_And_Type_mi_tag>(), NAME, TYPE); \
483  IT != FE->get_end<FENumeredDofEntityByNameAndType>( \
484  FE->colPtr->get<Composite_Name_And_Type_mi_tag>(), NAME, TYPE); \
485  IT++
486 
487 /** \brief loop over all dofs which are on a particular FE data, field and
488  * entity type \ingroup mofem_loops
489  */
490 #define _IT_GET_FEDATA_BY_TYPE_DOFS_FOR_LOOP_(FE, NAME, TYPE, IT) \
491  auto IT = FE->get_begin<FEDofEntityByNameAndType>( \
492  FE->dataPtr->get<Composite_Name_And_Type_mi_tag>(), NAME, TYPE); \
493  IT != FE->get_end<FEDofEntityByNameAndType>( \
494  FE->dataPtr->get<Composite_Name_And_Type_mi_tag>(), NAME, TYPE); \
495  IT++
496 
497  template <class MULTIINDEX>
498  typename MULTIINDEX::iterator get_begin(const MULTIINDEX &index,
499  const std::string &field_name) const {
500  return index.lower_bound(field_name);
501  }
502  template <class MULTIINDEX>
503  typename MULTIINDEX::iterator get_end(const MULTIINDEX &index,
504  const std::string &field_name) const {
505  return index.upper_bound(field_name);
506  }
507 
508 /** \brief loop over all dofs which are on a particular FE row and field
509  * \ingroup mofem_loops
510  */
511 #define _IT_GET_FEROW_BY_NAME_DOFS_FOR_LOOP_(FE, NAME, IT) \
512  auto IT = FE->get_begin<FENumeredDofEntityByFieldName>( \
513  FE->rowPtr->get<FieldName_mi_tag>(), NAME); \
514  IT != FE->get_end<FENumeredDofEntityByFieldName>( \
515  FE->rowPtr->get<FieldName_mi_tag>(), NAME); \
516  IT++
517 
518 /** \brief loop over all dofs which are on a particular FE column and field
519  * \ingroup mofem_loops
520  */
521 #define _IT_GET_FECOL_BY_NAME_DOFS_FOR_LOOP_(FE, NAME, IT) \
522  auto IT = FE->get_begin<FENumeredDofEntityByFieldName>( \
523  FE->colPtr->get<FieldName_mi_tag>(), NAME); \
524  IT != FE->get_end<FENumeredDofEntityByFieldName>( \
525  FE->colPtr->get<FieldName_mi_tag>(), NAME); \
526  IT++
527 
528 /** \brief loop over all dofs which are on a particular FE data and field
529  * \ingroup mofem_loops
530  */
531 #define _IT_GET_FEDATA_BY_NAME_DOFS_FOR_LOOP_(FE, NAME, IT) \
532  auto IT = FE->get_begin<FEDofEntityByFieldName>( \
533  FE->dataPtr->get<FieldName_mi_tag>(), NAME); \
534  IT != FE->get_end<FEDofEntityByFieldName>( \
535  FE->dataPtr->get<FieldName_mi_tag>(), NAME); \
536  IT++
537 
538  template <class MULTIINDEX>
539  typename MULTIINDEX::iterator get_begin(const MULTIINDEX &index,
540  const EntityHandle ent) const {
541  return index.lower_bound(ent);
542  }
543  template <class MULTIINDEX>
544  typename MULTIINDEX::iterator get_end(const MULTIINDEX &index,
545  const EntityHandle ent) const {
546  return index.upper_bound(ent);
547  }
548 
549 /** \brief loop over all dofs which are on a particular FE row and given element
550  * entity (handle from moab) \ingroup mofem_loops
551  */
552 #define _IT_GET_FEROW_DOFS_BY_ENT_FOR_LOOP_(FE, ENT, IT) \
553  auto IT = FE->get_begin<FENumeredDofEntityByEnt>( \
554  FE->rowPtr->get<Ent_mi_tag>(), ENT); \
555  IT != FE->get_end<FENumeredDofEntityByEnt>(FE->rowPtr->get<Ent_mi_tag>(), \
556  ENT); \
557  IT++
558 
559 /** \brief loop over all dofs which are on a particular FE column and given
560  * element entity (handle from moab) \ingroup mofem_loops
561  */
562 #define _IT_GET_FECOL_DOFS_BY_ENT_FOR_LOOP_(FE, ENT, IT) \
563  auto IT = FE->get_begin<FENumeredDofEntityByEnt>( \
564  FE->colPtr->get<Ent_mi_tag>(), ENT); \
565  IT != FE->get_end<FENumeredDofEntityByEnt>(FE->colPtr->get<Ent_mi_tag>(), \
566  ENT); \
567  IT++
568 
569 /** \brief loop over all dofs which are on a particular FE data and given
570  * element entity (handle from moab) \ingroup mofem_loops
571  */
572 #define _IT_GET_FEDATA_DOFS_BY_ENT_FOR_LOOP_(FE, ENT, IT) \
573  auto IT = FE->get_begin<FEDofEntity_multiIndex::index<Ent_mi_tag>::type>( \
574  FE->dataPtr->get<Ent_mi_tag>(), ENT); \
575  IT != FE->get_end<FEDofEntity_multiIndex::index<Ent_mi_tag>::type>( \
576  FE->dataPtr->get<Ent_mi_tag>(), ENT); \
577  IT++
578 
579  template <class MULTIINDEX>
580  typename MULTIINDEX::iterator get_begin(const MULTIINDEX &index,
581  const std::string &field_name,
582  const EntityHandle ent) const {
583  return index.lower_bound(boost::make_tuple(field_name, ent));
584  }
585  template <class MULTIINDEX>
586  typename MULTIINDEX::iterator get_end(const MULTIINDEX &index,
587  const std::string &field_name,
588  const EntityHandle ent) const {
589  return index.upper_bound(boost::make_tuple(field_name, ent));
590  }
591 
592 /** \brief loop over all dofs which are on a particular FE row, field and given
593  * element entity (handle from moab) \ingroup mofem_loops
594  */
595 #define _IT_GET_FEROW_DOFS_BY_NAME_AND_ENT_FOR_LOOP_(FE, NAME, ENT, IT) \
596  auto IT = FE->get_begin<FENumeredDofEntityByNameAndEnt>( \
597  FE->rowPtr->get<Composite_Name_And_Ent_mi_tag>(), NAME, ENT); \
598  IT != FE->get_end<FENumeredDofEntityByNameAndEnt>( \
599  FE->rowPtr->get<Composite_Name_And_Ent_mi_tag>(), NAME, ENT); \
600  IT++
601 
602 /** \brief loop over all dofs which are on a particular FE column, field and
603  * given element entity (handle from moab) \ingroup mofem_loops
604  */
605 #define _IT_GET_FECOL_DOFS_BY_NAME_AND_ENT_FOR_LOOP_(FE, NAME, ENT, IT) \
606  auto IT = FE->get_begin<FENumeredDofEntityByNameAndEnt>( \
607  FE->colPtr->get<Composite_Name_And_Ent_mi_tag>(), NAME, ENT); \
608  IT != FE->get_end<FENumeredDofEntityByNameAndEnt>( \
609  FE->colPtr->get<Composite_Name_And_Ent_mi_tag>(), NAME, ENT); \
610  IT++
611 
612 /** \brief loop over all dofs which are on a particular FE data, field and given
613  * element entity (handle from moab) \ingroup mofem_loops
614  */
615 #define _IT_GET_FEDATA_DOFS_BY_NAME_AND_ENT_FOR_LOOP_(FE, NAME, ENT, IT) \
616  auto IT = FE->get_begin<FEDofEntityByNameAndEnt>( \
617  FE->dataPtr->get<Composite_Name_And_Ent_mi_tag>(), NAME, ENT); \
618  IT != FE->get_end<FEDofEntityByNameAndEnt>( \
619  FE->dataPtr->get<Composite_Name_And_Ent_mi_tag>(), NAME, ENT); \
620  IT++
621 };
622 
623 /**
624  * \brief Data structure to exchange data between mofem and User Loop Methods on
625  * entities. \ingroup mofem_loops
626  *
627  * It allows to exchange data between MoFEM and user functions. It stores
628  * information about multi-indices.
629  */
630 struct EntityMethod : public BasicMethod {
631 
633  UnknownInterface **iface) const {
635  if (uuid == IDD_MOFEMEntityMethod) {
636  *iface = const_cast<EntityMethod *>(this);
638  }
639  CHKERR query_interface(uuid, iface);
641  }
642 
643  EntityMethod();
644 
645  boost::shared_ptr<Field> fieldPtr;
646  boost::shared_ptr<FieldEntity> entPtr;
647 };
648 
649 /**
650  * \brief Data structure to exchange data between mofem and User Loop Methods on
651  * entities. \ingroup mofem_loops
652  *
653  * It allows to exchange data between MoFEM and user functions. It stores
654  * information about multi-indices.
655  */
656 struct DofMethod : public BasicMethod {
657 
659  UnknownInterface **iface) const {
661  if (uuid == IDD_MOFEMDofMethod) {
662  *iface = const_cast<DofMethod *>(this);
664  }
665 
666  CHKERR query_interface(uuid, iface);
668  }
669 
670  DofMethod();
671 
672  boost::shared_ptr<Field> fieldPtr;
673  boost::shared_ptr<DofEntity> dofPtr;
674  boost::shared_ptr<NumeredDofEntity> dofNumeredPtr;
675 };
676 
677 /// \deprecated name changed use DofMethod insead EntMethod
679 
680 } // namespace MoFEM
681 
682 #endif // __LOOPMETHODS_HPP__
683 
684 /**
685  * \defgroup mofem_loops Loops
686  * \ingroup mofem
687  */
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:501
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:477
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:544
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:508
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
boost::movelib::unique_ptr< bool > matAssembleSwitch
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:596
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:407
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::movelib::unique_ptr< bool > vecAssembleSwitch
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