v0.8.19
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 /**
147  * \brief data structure for TS (time stepping) context
148  * \ingroup mofem_loops
149  *
150  * Structure stores context data which are set in functions run by PETSc Time
151  * Stepping functions.
152  */
153 struct TSMethod : virtual public UnknownInterface {
154 
156  UnknownInterface **iface) const {
157  if (uuid == IDD_MOFEMTsMethod) {
158  *iface = const_cast<TSMethod *>(this);
160  }
161  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "unknown interface");
162  }
163 
164  enum TSContext {
171  };
172 
175  : ts_ctx(CTX_TSNONE), ts_u(PETSC_NULL), ts_u_t(PETSC_NULL),
176  ts_F(PETSC_NULL), ts_A(PETSC_NULL), ts_B(PETSC_NULL), ts_step(-1),
177  ts_a(0), ts_t(0) {}
178 
179  virtual ~TSMethod() {}
180 
181  /// \brief Set Ts context
182  MoFEMErrorCode setTsCtx(const TSContext &ctx);
183 
184  /// \brief Copy TS solver data
186 
187  /// \brief Set TS solver
188  MoFEMErrorCode setTs(TS _ts);
189 
190  TS ts;
191  Vec ts_u, ts_u_t, ts_F;
192  Mat ts_A, ts_B;
193 
194  PetscInt ts_step;
195  PetscReal ts_a, ts_t;
196 
197 };
198 
199 /**
200  * \brief Data structure to exchange data between mofem and User Loop Methods.
201  * \ingroup mofem_loops
202  *
203  * It allows to exchange data between MoFEM and user functions. It stores
204  * information about multi-indices.
205  *
206  */
208 
210  UnknownInterface **iface) const {
211  if (uuid == IDD_MOFEMBasicMethod) {
212  *iface = const_cast<BasicMethod *>(this);
214  }
215  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "unknown interface");
216  }
217 
218  BasicMethod();
219  virtual ~BasicMethod(){};
220 
221  int nInTheLoop;
222  int loopSize;
223 
224  /** \brief get number of evaluated element in the loop
225  */
226  inline int getNinTheLoop() const { return nInTheLoop; }
227 
228  /** \brief get loop size
229  */
230  inline int getLoopSize() const { return loopSize; }
231 
232  int rAnk, sIze;
242 
244 
245  /**
246  * Hook function for pre-processing
247  */
248  boost::function<MoFEMErrorCode()> preProcessHook;
249 
250  /**
251  * Hook function for post-processing
252  */
253  boost::function<MoFEMErrorCode()> postProcessHook;
254 
255  /**
256  * Hook function for operator
257  */
258  boost::function<MoFEMErrorCode()> operatorHook;
259 
260  /** \brief function is run at the beginning of loop
261  *
262  * It is used to zeroing matrices and vectors, calculation of shape
263  * functions on reference element, preprocessing boundary conditions, etc.
264  */
265  virtual MoFEMErrorCode preProcess();
266 
267  /** \brief function is run for every finite element
268  *
269  * It is used to calculate element local matrices and assembly. It can be
270  * used for post-processing.
271  */
272  virtual MoFEMErrorCode operator()();
273 
274  /** \brief function is run at the end of loop
275  *
276  * It is used to assembly matrices and vectors, calculating global variables,
277  * f.e. total internal energy, ect.
278  *
279  * Iterating over dofs:
280  * Example1 iterating over dofs in row by name of the field
281  * for(_IT_GET_FEROW_BY_NAME_DOFS_FOR_LOOP_(this,"DISPLACEMENT",it)) { ... }
282  *
283  *
284  */
285  virtual MoFEMErrorCode postProcess();
286 
287 private:
288  void iNit();
289 };
290 
291 /**
292  * \brief structure for User Loop Methods on finite elements
293  * \ingroup mofem_loops
294  *
295  * It can be used to calculate stiffness matrices, residuals, load vectors etc.
296  * It is low level class however in some class users looking for speed and
297  * efficiency, can use it directly.
298  *
299  * This class is used with Interface::loop_finite_elements, where
300  * user overloaded operator FEMethod::operator() is executed for each element in
301  * the problem. Class have to additional methods which are overloaded by user,
302  * FEMethod::preProcess() and FEMethod::postProcess() executed at beginning and
303  * end of the loop over problem elements, respectively.
304  *
305  */
306 struct FEMethod : public BasicMethod {
307 
309  UnknownInterface **iface) const {
311  if (uuid == IDD_MOFEMFEMethod) {
312  *iface = const_cast<FEMethod *>(this);
314  }
315 
316  ierr = query_interface(uuid, iface);
317  CHKERRG(ierr);
319  }
320 
321  FEMethod();
322 
323  std::string feName; ///< Name of finite element
324 
325  boost::shared_ptr<const NumeredEntFiniteElement>
326  numeredEntFiniteElementPtr; ///< Pointer to finite element database
327  ///< structure
328  boost::shared_ptr<const FENumeredDofEntity_multiIndex>
329  rowPtr; ///< Pointer to finite element rows dofs view
330  boost::shared_ptr<const FENumeredDofEntity_multiIndex>
331  colPtr; ///< Pointer to finite element columns dofs view
332  boost::shared_ptr<const FEDofEntity_multiIndex>
333  dataPtr; ///< Pointer to finite element data dofs
334 
335  boost::shared_ptr<const FieldEntity_vector_view>
336  rowFieldEntsPtr; ///< Pointer to finite element field entities row view
337  boost::shared_ptr<const FieldEntity_vector_view>
338  colFieldEntsPtr; ///< Pointer to finite element field entities column view
339  boost::shared_ptr<const FieldEntity_multiIndex_spaceType_view>
340  dataFieldEntsPtr; ///< Pointer to finite element field entities data view
341 
342 /** \brief loop over all dofs which are on a particular FE row
343  * \ingroup mofem_loops
344  */
345 #define _IT_GET_FEROW_DOFS_FOR_LOOP_(FE, IT) \
346  auto IT = FE->rowPtr->begin(); \
347  IT != FE->rowPtr->end(); \
348  IT++
349 
350 /** \brief loop over all dofs which are on a particular FE column
351  * \ingroup mofem_loops
352  */
353 #define _IT_GET_FECOL_DOFS_FOR_LOOP_(FE, IT) \
354  auto IT = FE->colPtr->begin(); \
355  IT != FE->colPtr->end(); \
356  IT++
357 
358 /** \brief loop over all dofs which are on a particular FE data
359  * \ingroup mofem_loops
360  */
361 #define _IT_GET_FEDATA_DOFS_FOR_LOOP_(FE, IT) \
362  auto IT = FE->dataPtr->begin(); \
363  IT != FE->dataPtr->end(); \
364  IT++
365 
366  template <class MULTIINDEX>
367  typename MULTIINDEX::iterator
368  get_begin(const MULTIINDEX &index, const std::string &field_name,
369  const EntityType type, const int side_number) const {
370  return index.lower_bound(boost::make_tuple(field_name, type, side_number));
371  }
372 
373  template <class MULTIINDEX>
374  typename MULTIINDEX::iterator
375  get_end(const MULTIINDEX &index, const std::string &field_name,
376  const EntityType type, const int side_number) const {
377  return index.upper_bound(boost::make_tuple(field_name, type, side_number));
378  }
379 
380  /** \brief loop over all dofs which are on a particular FE row, field,
381  * entity type and canonical side number \ingroup mofem_loops
382  *
383  * \param FE finite elements
384  * \param Name field name
385  * \param Type moab entity type (MBVERTEX, MBEDGE etc)
386  * \param Side side canonical number
387  * \param IT the interator in use
388  */
389 #define _IT_GET_FEROW_BY_SIDE_DOFS_FOR_LOOP_(FE, NAME, TYPE, SIDE, IT) \
390  FENumeredDofEntity_multiIndex::index<Composite_mi_tag>::type::iterator IT = \
391  FE->get_begin< \
392  FENumeredDofEntity_multiIndex::index<Composite_mi_tag>::type>( \
393  FE->rowPtr->get<Composite_mi_tag>(), NAME, TYPE, SIDE); \
394  IT != FE->get_end< \
395  FENumeredDofEntity_multiIndex::index<Composite_mi_tag>::type>( \
396  FE->rowPtr->get<Composite_mi_tag>(), NAME, TYPE, SIDE); \
397  IT++
398 
399 /** \brief loop over all dofs which are on a particular FE column, field, entity
400  * type and canonical side number \ingroup mofem_loops
401  */
402 #define _IT_GET_FECOL_BY_SIDE_DOFS_FOR_LOOP_(FE, NAME, TYPE, SIDE, IT) \
403  FENumeredDofEntity_multiIndex::index<Composite_mi_tag>::type::iterator IT = \
404  FE->get_begin< \
405  FENumeredDofEntity_multiIndex::index<Composite_mi_tag>::type>( \
406  FE->colPtr->get<Composite_mi_tag>(), NAME, TYPE, SIDE); \
407  IT != FE->get_end< \
408  FENumeredDofEntity_multiIndex::index<Composite_mi_tag>::type>( \
409  FE->colPtr->get<Composite_mi_tag>(), NAME, TYPE, SIDE); \
410  IT++
411 
412 /** \brief loop over all dofs which are on a particular FE data, field, entity
413  * type and canonical side number \ingroup mofem_loops
414  */
415 #define _IT_GET_FEDATA_BY_SIDE_DOFS_FOR_LOOP_(FE, NAME, TYPE, SIDE, IT) \
416  FEDofEntity_multiIndex::index<Composite_mi_tag>::type::iterator IT = \
417  FE->get_begin<FEDofEntity_multiIndex::index<Composite_mi_tag>::type>( \
418  FE->dataPtr->get<Composite_mi_tag>(), NAME, TYPE, SIDE); \
419  IT != FE->get_end<FEDofEntity_multiIndex::index<Composite_mi_tag>::type>( \
420  FE->dataPtr->get<Composite_mi_tag>(), NAME, TYPE, SIDE); \
421  IT++
422 
423  template <class MULTIINDEX>
424  typename MULTIINDEX::iterator get_begin(const MULTIINDEX &index,
425  const std::string &field_name,
426  const EntityType type) const {
427  return index.lower_bound(boost::make_tuple(field_name, type));
428  }
429  template <class MULTIINDEX>
430  typename MULTIINDEX::iterator get_end(const MULTIINDEX &index,
431  const std::string &field_name,
432  const EntityType type) const {
433  return index.upper_bound(boost::make_tuple(field_name, type));
434  }
435 
436 /** \brief loop over all dofs which are on a particular FE row, field and entity
437  * type \ingroup mofem_loops
438  */
439 #define _IT_GET_FEROW_BY_TYPE_DOFS_FOR_LOOP_(FE, NAME, TYPE, IT) \
440  auto IT = FE->get_begin<FENumeredDofEntityByNameAndType>( \
441  FE->rowPtr->get<Composite_Name_And_Type_mi_tag>(), NAME, TYPE); \
442  IT != FE->get_end<FENumeredDofEntityByNameAndType>( \
443  FE->rowPtr->get<Composite_Name_And_Type_mi_tag>(), NAME, TYPE); \
444  IT++
445 
446 /** \brief loop over all dofs which are on a particular FE column, field and
447  * entity type \ingroup mofem_loops
448  */
449 #define _IT_GET_FECOL_BY_TYPE_DOFS_FOR_LOOP_(FE, NAME, TYPE, IT) \
450  auto IT = FE->get_begin<FENumeredDofEntityByNameAndType>( \
451  FE->colPtr->get<Composite_Name_And_Type_mi_tag>(), NAME, TYPE); \
452  IT != FE->get_end<FENumeredDofEntityByNameAndType>( \
453  FE->colPtr->get<Composite_Name_And_Type_mi_tag>(), NAME, TYPE); \
454  IT++
455 
456 /** \brief loop over all dofs which are on a particular FE data, field and
457  * entity type \ingroup mofem_loops
458  */
459 #define _IT_GET_FEDATA_BY_TYPE_DOFS_FOR_LOOP_(FE, NAME, TYPE, IT) \
460  auto IT = FE->get_begin<FEDofEntityByNameAndType>( \
461  FE->dataPtr->get<Composite_Name_And_Type_mi_tag>(), NAME, TYPE); \
462  IT != FE->get_end<FEDofEntityByNameAndType>( \
463  FE->dataPtr->get<Composite_Name_And_Type_mi_tag>(), NAME, TYPE); \
464  IT++
465 
466  template <class MULTIINDEX>
467  typename MULTIINDEX::iterator get_begin(const MULTIINDEX &index,
468  const std::string &field_name) const {
469  return index.lower_bound(field_name);
470  }
471  template <class MULTIINDEX>
472  typename MULTIINDEX::iterator get_end(const MULTIINDEX &index,
473  const std::string &field_name) const {
474  return index.upper_bound(field_name);
475  }
476 
477 /** \brief loop over all dofs which are on a particular FE row and field
478  * \ingroup mofem_loops
479  */
480 #define _IT_GET_FEROW_BY_NAME_DOFS_FOR_LOOP_(FE, NAME, IT) \
481  auto IT = FE->get_begin<FENumeredDofEntityByFieldName>( \
482  FE->rowPtr->get<FieldName_mi_tag>(), NAME); \
483  IT != FE->get_end<FENumeredDofEntityByFieldName>( \
484  FE->rowPtr->get<FieldName_mi_tag>(), NAME); \
485  IT++
486 
487 /** \brief loop over all dofs which are on a particular FE column and field
488  * \ingroup mofem_loops
489  */
490 #define _IT_GET_FECOL_BY_NAME_DOFS_FOR_LOOP_(FE, NAME, IT) \
491  auto IT = FE->get_begin<FENumeredDofEntityByFieldName>( \
492  FE->colPtr->get<FieldName_mi_tag>(), NAME); \
493  IT != FE->get_end<FENumeredDofEntityByFieldName>( \
494  FE->colPtr->get<FieldName_mi_tag>(), NAME); \
495  IT++
496 
497 /** \brief loop over all dofs which are on a particular FE data and field
498  * \ingroup mofem_loops
499  */
500 #define _IT_GET_FEDATA_BY_NAME_DOFS_FOR_LOOP_(FE, NAME, IT) \
501  auto IT = FE->get_begin<FEDofEntityByFieldName>( \
502  FE->dataPtr->get<FieldName_mi_tag>(), NAME); \
503  IT != FE->get_end<FEDofEntityByFieldName>( \
504  FE->dataPtr->get<FieldName_mi_tag>(), NAME); \
505  IT++
506 
507  template <class MULTIINDEX>
508  typename MULTIINDEX::iterator get_begin(const MULTIINDEX &index,
509  const EntityHandle ent) const {
510  return index.lower_bound(ent);
511  }
512  template <class MULTIINDEX>
513  typename MULTIINDEX::iterator get_end(const MULTIINDEX &index,
514  const EntityHandle ent) const {
515  return index.upper_bound(ent);
516  }
517 
518 /** \brief loop over all dofs which are on a particular FE row and given element
519  * entity (handle from moab) \ingroup mofem_loops
520  */
521 #define _IT_GET_FEROW_DOFS_BY_ENT_FOR_LOOP_(FE, ENT, IT) \
522  auto IT = FE->get_begin<FENumeredDofEntityByEnt>( \
523  FE->rowPtr->get<Ent_mi_tag>(), ENT); \
524  IT != FE->get_end<FENumeredDofEntityByEnt>(FE->rowPtr->get<Ent_mi_tag>(), \
525  ENT); \
526  IT++
527 
528 /** \brief loop over all dofs which are on a particular FE column and given
529  * element entity (handle from moab) \ingroup mofem_loops
530  */
531 #define _IT_GET_FECOL_DOFS_BY_ENT_FOR_LOOP_(FE, ENT, IT) \
532  auto IT = FE->get_begin<FENumeredDofEntityByEnt>( \
533  FE->colPtr->get<Ent_mi_tag>(), ENT); \
534  IT != FE->get_end<FENumeredDofEntityByEnt>(FE->colPtr->get<Ent_mi_tag>(), \
535  ENT); \
536  IT++
537 
538 /** \brief loop over all dofs which are on a particular FE data and given
539  * element entity (handle from moab) \ingroup mofem_loops
540  */
541 #define _IT_GET_FEDATA_DOFS_BY_ENT_FOR_LOOP_(FE, ENT, IT) \
542  auto IT = FE->get_begin<FEDofEntity_multiIndex::index<Ent_mi_tag>::type>( \
543  FE->dataPtr->get<Ent_mi_tag>(), ENT); \
544  IT != FE->get_end<FEDofEntity_multiIndex::index<Ent_mi_tag>::type>( \
545  FE->dataPtr->get<Ent_mi_tag>(), ENT); \
546  IT++
547 
548  template <class MULTIINDEX>
549  typename MULTIINDEX::iterator get_begin(const MULTIINDEX &index,
550  const std::string &field_name,
551  const EntityHandle ent) const {
552  return index.lower_bound(boost::make_tuple(field_name, ent));
553  }
554  template <class MULTIINDEX>
555  typename MULTIINDEX::iterator get_end(const MULTIINDEX &index,
556  const std::string &field_name,
557  const EntityHandle ent) const {
558  return index.upper_bound(boost::make_tuple(field_name, ent));
559  }
560 
561 /** \brief loop over all dofs which are on a particular FE row, field and given
562  * element entity (handle from moab) \ingroup mofem_loops
563  */
564 #define _IT_GET_FEROW_DOFS_BY_NAME_AND_ENT_FOR_LOOP_(FE, NAME, ENT, IT) \
565  auto IT = FE->get_begin<FENumeredDofEntityByNameAndEnt>( \
566  FE->rowPtr->get<Composite_Name_And_Ent_mi_tag>(), NAME, ENT); \
567  IT != FE->get_end<FENumeredDofEntityByNameAndEnt>( \
568  FE->rowPtr->get<Composite_Name_And_Ent_mi_tag>(), NAME, ENT); \
569  IT++
570 
571 /** \brief loop over all dofs which are on a particular FE column, field and
572  * given element entity (handle from moab) \ingroup mofem_loops
573  */
574 #define _IT_GET_FECOL_DOFS_BY_NAME_AND_ENT_FOR_LOOP_(FE, NAME, ENT, IT) \
575  auto IT = FE->get_begin<FENumeredDofEntityByNameAndEnt>( \
576  FE->colPtr->get<Composite_Name_And_Ent_mi_tag>(), NAME, ENT); \
577  IT != FE->get_end<FENumeredDofEntityByNameAndEnt>( \
578  FE->colPtr->get<Composite_Name_And_Ent_mi_tag>(), NAME, ENT); \
579  IT++
580 
581 /** \brief loop over all dofs which are on a particular FE data, field and given
582  * element entity (handle from moab) \ingroup mofem_loops
583  */
584 #define _IT_GET_FEDATA_DOFS_BY_NAME_AND_ENT_FOR_LOOP_(FE, NAME, ENT, IT) \
585  auto IT = FE->get_begin<FEDofEntityByNameAndEnt>( \
586  FE->dataPtr->get<Composite_Name_And_Ent_mi_tag>(), NAME, ENT); \
587  IT != FE->get_end<FEDofEntityByNameAndEnt>( \
588  FE->dataPtr->get<Composite_Name_And_Ent_mi_tag>(), NAME, ENT); \
589  IT++
590 };
591 
592 /**
593  * \brief Data structure to exchange data between mofem and User Loop Methods on
594  * entities. \ingroup mofem_loops
595  *
596  * \todo Add implementation to loop over entities in the problem
597  *
598  * It allows to exchange data between MoFEM and user functions. It stores
599  * information about multi-indices.
600  */
601 struct EntityMethod : public BasicMethod {
602 
604  UnknownInterface **iface) const {
606  if (uuid == IDD_MOFEMEntityMethod) {
607  *iface = const_cast<EntityMethod *>(this);
609  }
610  CHKERR query_interface(uuid, iface);
612  }
613 
614  EntityMethod();
615 
616  boost::shared_ptr<Field> fieldPtr;
617  boost::shared_ptr<FieldEntity> entPtr;
618 
619 };
620 
621 
622 /**
623  * \brief Data structure to exchange data between mofem and User Loop Methods on
624  * entities. \ingroup mofem_loops
625  *
626  * It allows to exchange data between MoFEM and user functions. It stores
627  * information about multi-indices.
628  */
629 struct DofMethod : public BasicMethod {
630 
632  UnknownInterface **iface) const {
634  if (uuid == IDD_MOFEMDofMethod) {
635  *iface = const_cast<DofMethod *>(this);
637  }
638 
639  CHKERR query_interface(uuid, iface);
641  }
642 
643  DofMethod();
644 
645  boost::shared_ptr<Field> fieldPtr;
646  boost::shared_ptr<DofEntity> dofPtr;
647  boost::shared_ptr<NumeredDofEntity> dofNumeredPtr;
648 };
649 
650 /// \deprecated name changed use DofMethod insead EntMethod
652 
653 } // namespace MoFEM
654 
655 #endif // __LOOPMETHODS_HPP__
656 
657 /***************************************************************************/ /**
658  * \defgroup mofem_loops Loops
659  * \ingroup mofem
660  ******************************************************************************/
MULTIINDEX::iterator get_end(const MULTIINDEX &index, const EntityHandle ent) const
boost::shared_ptr< DofEntity > dofPtr
structure for User Loop Methods on finite elementsIt can be used to calculate stiffness matrices...
const DofEntity_multiIndex * dofsPtr
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
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Common.hpp:60
virtual MoFEMErrorCode postProcess()
function is run at the end of loop
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
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
virtual MoFEMErrorCode operator()()
function is run for every finite element
MoFEMErrorCode setSnes(SNES snes)
Set SNES instance.
Definition: LoopMethods.cpp:48
Data structure to exchange data between mofem and User Loop Methods on entities.
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:499
MULTIINDEX::iterator get_begin(const MULTIINDEX &index, const std::string &field_name) const
const FieldEntity_multiIndex * entitiesPtr
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:475
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:542
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
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
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:506
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, what elements compose problem and what DOFs are on rows and columns.
KSPContext ksp_ctx
Context.
Definition: LoopMethods.hpp:90
boost::shared_ptr< const FEDofEntity_multiIndex > dataPtr
Pointer to finite element data dofs.
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
const Problem * problemPtr
virtual MoFEMErrorCode preProcess()
function is run at the beginning of loop
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
static MoFEMErrorCodeGeneric< PetscErrorCode > ierr
Definition: Common.hpp:80
const RefEntity_multiIndex * refinedEntitiesPtr
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
MoFEMErrorCode copyBasicMethod(const BasicMethod &basic)
Definition: LoopMethods.cpp:99
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.
Mat ksp_A
matrix
Definition: LoopMethods.hpp:93
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:594
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.
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.
std::bitset< BITINTERFACEUID_SIZE > BitIntefaceId
Definition: Common.hpp:152
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:405
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
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