v0.9.0
LoopMethods.hpp
Go to the documentation of this file.
1 /** \file LoopMethods.hpp
2  * \brief MoFEM interface
3  *
4  * Data structures for making loops over finite elements and entities in the
5  * problem or MoFEM database.
6  *
7  */
8 
9 /*
10  * MoFEM is distributed in the hope that it will be useful, but WITHOUT
11  * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
12  * FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
13  * License for more details.
14  *
15  * You should have received a copy of the GNU Lesser General Public
16  * License along with MoFEM. If not, see <http://www.gnu.org/licenses/>
17  */
18 
19 #ifndef __LOOPMETHODS_HPP__
20 #define __LOOPMETHODS_HPP__
21 
22 namespace MoFEM {
23 
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; ///< snes solver
141  Vec snes_x; ///< state vector
142  Vec snes_f; ///< residual
143  Mat snes_A; ///< jacobian matrix
144  Mat snes_B; ///< preconditioner of jacobian matrix
145 };
146 
147 /**
148  * \brief data structure for TS (time stepping) context
149  * \ingroup mofem_loops
150  *
151  * Structure stores context data which are set in functions run by PETSc Time
152  * Stepping functions.
153  */
154 struct TSMethod : virtual public UnknownInterface {
155 
157  UnknownInterface **iface) const {
158  if (uuid == IDD_MOFEMTsMethod) {
159  *iface = const_cast<TSMethod *>(this);
161  }
162  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "unknown interface");
163  }
164 
165  enum TSContext {
172  };
173 
176  : ts_ctx(CTX_TSNONE), ts_u(PETSC_NULL), ts_u_t(PETSC_NULL),
177  ts_F(PETSC_NULL), ts_A(PETSC_NULL), ts_B(PETSC_NULL), ts_step(-1),
178  ts_a(0), ts_t(0) {}
179 
180  virtual ~TSMethod() {}
181 
182  /// \brief Set Ts context
183  MoFEMErrorCode setTsCtx(const TSContext &ctx);
184 
185  /// \brief Copy TS solver data
187 
188  /// \brief Set TS solver
189  MoFEMErrorCode setTs(TS _ts);
190 
191  TS ts; ///< time solver
192  Vec ts_u; ///< state vector
193  Vec ts_u_t; ///< time derivative of state vector
194  Vec ts_u_tt; ///< second time derivative of state vector
195  Vec ts_F; ///< residual vector
196 
197  Mat ts_A; ///< Jacobian of G(U) = F(t,U,W+v*U,W'+a*U), equivalent to dF/dU +
198  ///< v*dF/dU_t + a*dF/dU_tt
199  Mat ts_B; ///< Preconditioner for ts_A
200 
201  PetscInt ts_step; ///< time step
202  PetscReal ts_a; ///< shift for U_tt (see PETSc Time Solver)
203  PetscReal ts_v; ///< shift for U_t shift for U_t
204  PetscReal ts_t; ///< time
205 };
206 
207 /**
208  * \brief Data structure to exchange data between mofem and User Loop Methods.
209  * \ingroup mofem_loops
210  *
211  * It allows to exchange data between MoFEM and user functions. It stores
212  * information about multi-indices.
213  *
214  */
216 
218  UnknownInterface **iface) const {
219  if (uuid == IDD_MOFEMBasicMethod) {
220  *iface = const_cast<BasicMethod *>(this);
222  }
223  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "unknown interface");
224  }
225 
226  BasicMethod();
227  virtual ~BasicMethod(){};
228 
229  /**
230  * @brief number currently of processed method
231  */
232  int nInTheLoop;
233 
234  /**
235  * @brief local number oe methods to process
236  */
237  int loopSize;
238 
239  /** \brief get number of evaluated element in the loop
240  */
241  inline int getNinTheLoop() const { return nInTheLoop; }
242 
243  /** \brief get loop size
244  */
245  inline int getLoopSize() const { return loopSize; }
246 
247  int rAnk; ///< processor rank
248 
249  int sIze; ///< number of processors in communicator
250 
252  *refinedEntitiesPtr; ///< container of mofem dof entities
253 
255  *refinedFiniteElementsPtr; ///< container of mofem finite element entities
256 
257  const Problem *problemPtr; ///< raw pointer to problem
258 
259  const Field_multiIndex *fieldsPtr; ///< raw pointer to fields container
260 
262  *entitiesPtr; ///< raw pointer to container of field entities
263 
264  const DofEntity_multiIndex *dofsPtr; ///< raw pointer container of dofs
265 
267  *finiteElementsPtr; ///< raw pointer to container finite elements
268 
270  *finiteElementsEntitiesPtr; ///< raw pointer to container finite elements
271  ///< entities
272 
274  *adjacenciesPtr; ///< raw pointer to container to adjacencies between dofs
275  ///< and finite elements
276 
277  /**
278  * @brief Copy data from other base method to this base method
279  *
280  * @param basic
281  * @return MoFEMErrorCode
282  */
284 
285  /**
286  * @brief Hook function for pre-processing
287  */
288  boost::function<MoFEMErrorCode()> preProcessHook;
289 
290  /**
291  * @brief Hook function for post-processing
292  */
293  boost::function<MoFEMErrorCode()> postProcessHook;
294 
295  /**
296  * @brief Hook function for operator
297  */
298  boost::function<MoFEMErrorCode()> operatorHook;
299 
300  /** \brief function is run at the beginning of loop
301  *
302  * It is used to zeroing matrices and vectors, calculation of shape
303  * functions on reference element, preprocessing boundary conditions, etc.
304  */
305  virtual MoFEMErrorCode preProcess();
306 
307  /** \brief function is run for every finite element
308  *
309  * It is used to calculate element local matrices and assembly. It can be
310  * used for post-processing.
311  */
312  virtual MoFEMErrorCode operator()();
313 
314  /** \brief function is run at the end of loop
315  *
316  * It is used to assembly matrices and vectors, calculating global variables,
317  * f.e. total internal energy, ect.
318  *
319  * Iterating over dofs:
320  * Example1 iterating over dofs in row by name of the field
321  * for(_IT_GET_FEROW_BY_NAME_DOFS_FOR_LOOP_(this,"DISPLACEMENT",it)) { ... }
322  *
323  *
324  */
325  virtual MoFEMErrorCode postProcess();
326 
327  boost::movelib::unique_ptr<bool> vecAssembleSwitch;
328  boost::movelib::unique_ptr<bool> matAssembleSwitch;
329 };
330 
331 /**
332  * \brief structure for User Loop Methods on finite elements
333  * \ingroup mofem_loops
334  *
335  * It can be used to calculate stiffness matrices, residuals, load vectors etc.
336  * It is low level class however in some class users looking for speed and
337  * efficiency, can use it directly.
338  *
339  * This class is used with Interface::loop_finite_elements, where
340  * user overloaded operator FEMethod::operator() is executed for each element in
341  * the problem. Class have to additional methods which are overloaded by user,
342  * FEMethod::preProcess() and FEMethod::postProcess() executed at beginning and
343  * end of the loop over problem elements, respectively.
344  *
345  */
346 struct FEMethod : public BasicMethod {
347 
349  UnknownInterface **iface) const {
351  if (uuid == IDD_MOFEMFEMethod) {
352  *iface = const_cast<FEMethod *>(this);
354  }
355 
356  ierr = query_interface(uuid, iface);
357  CHKERRG(ierr);
359  }
360 
361  FEMethod();
362 
363  std::string feName; ///< Name of finite element
364 
365  boost::shared_ptr<const NumeredEntFiniteElement>
366  numeredEntFiniteElementPtr; ///< Pointer to finite element database
367  ///< structure
368  boost::shared_ptr<const FENumeredDofEntity_multiIndex>
369  rowPtr; ///< Pointer to finite element rows dofs view
370  boost::shared_ptr<const FENumeredDofEntity_multiIndex>
371  colPtr; ///< Pointer to finite element columns dofs view
372  boost::shared_ptr<const FEDofEntity_multiIndex>
373  dataPtr; ///< Pointer to finite element data dofs
374 
375  boost::shared_ptr<const FieldEntity_vector_view>
376  rowFieldEntsPtr; ///< Pointer to finite element field entities row view
377  boost::shared_ptr<const FieldEntity_vector_view>
378  colFieldEntsPtr; ///< Pointer to finite element field entities column view
379  boost::shared_ptr<const FieldEntity_multiIndex_spaceType_view>
380  dataFieldEntsPtr; ///< Pointer to finite element field entities data view
381 
382 /** \brief loop over all dofs which are on a particular FE row
383  * \ingroup mofem_loops
384  */
385 #define _IT_GET_FEROW_DOFS_FOR_LOOP_(FE, IT) \
386  auto IT = FE->rowPtr->begin(); \
387  IT != FE->rowPtr->end(); \
388  IT++
389 
390 /** \brief loop over all dofs which are on a particular FE column
391  * \ingroup mofem_loops
392  */
393 #define _IT_GET_FECOL_DOFS_FOR_LOOP_(FE, IT) \
394  auto IT = FE->colPtr->begin(); \
395  IT != FE->colPtr->end(); \
396  IT++
397 
398 /** \brief loop over all dofs which are on a particular FE data
399  * \ingroup mofem_loops
400  */
401 #define _IT_GET_FEDATA_DOFS_FOR_LOOP_(FE, IT) \
402  auto IT = FE->dataPtr->begin(); \
403  IT != FE->dataPtr->end(); \
404  IT++
405 
406  template <class MULTIINDEX>
407  typename MULTIINDEX::iterator
408  get_begin(const MULTIINDEX &index, const std::string &field_name,
409  const EntityType type, const int side_number) const {
410  return index.lower_bound(boost::make_tuple(field_name, type, side_number));
411  }
412 
413  template <class MULTIINDEX>
414  typename MULTIINDEX::iterator
415  get_end(const MULTIINDEX &index, const std::string &field_name,
416  const EntityType type, const int side_number) const {
417  return index.upper_bound(boost::make_tuple(field_name, type, side_number));
418  }
419 
420  /** \brief loop over all dofs which are on a particular FE row, field,
421  * entity type and canonical side number \ingroup mofem_loops
422  *
423  * \param FE finite elements
424  * \param Name field name
425  * \param Type moab entity type (MBVERTEX, MBEDGE etc)
426  * \param Side side canonical number
427  * \param IT the interator in use
428  */
429 #define _IT_GET_FEROW_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->rowPtr->get<Composite_mi_tag>(), NAME, TYPE, SIDE); \
434  IT != FE->get_end< \
435  FENumeredDofEntity_multiIndex::index<Composite_mi_tag>::type>( \
436  FE->rowPtr->get<Composite_mi_tag>(), NAME, TYPE, SIDE); \
437  IT++
438 
439 /** \brief loop over all dofs which are on a particular FE column, field, entity
440  * type and canonical side number \ingroup mofem_loops
441  */
442 #define _IT_GET_FECOL_BY_SIDE_DOFS_FOR_LOOP_(FE, NAME, TYPE, SIDE, IT) \
443  FENumeredDofEntity_multiIndex::index<Composite_mi_tag>::type::iterator IT = \
444  FE->get_begin< \
445  FENumeredDofEntity_multiIndex::index<Composite_mi_tag>::type>( \
446  FE->colPtr->get<Composite_mi_tag>(), NAME, TYPE, SIDE); \
447  IT != FE->get_end< \
448  FENumeredDofEntity_multiIndex::index<Composite_mi_tag>::type>( \
449  FE->colPtr->get<Composite_mi_tag>(), NAME, TYPE, SIDE); \
450  IT++
451 
452 /** \brief loop over all dofs which are on a particular FE data, field, entity
453  * type and canonical side number \ingroup mofem_loops
454  */
455 #define _IT_GET_FEDATA_BY_SIDE_DOFS_FOR_LOOP_(FE, NAME, TYPE, SIDE, IT) \
456  FEDofEntity_multiIndex::index<Composite_mi_tag>::type::iterator IT = \
457  FE->get_begin<FEDofEntity_multiIndex::index<Composite_mi_tag>::type>( \
458  FE->dataPtr->get<Composite_mi_tag>(), NAME, TYPE, SIDE); \
459  IT != FE->get_end<FEDofEntity_multiIndex::index<Composite_mi_tag>::type>( \
460  FE->dataPtr->get<Composite_mi_tag>(), NAME, TYPE, SIDE); \
461  IT++
462 
463  template <class MULTIINDEX>
464  typename MULTIINDEX::iterator get_begin(const MULTIINDEX &index,
465  const std::string &field_name,
466  const EntityType type) const {
467  return index.lower_bound(boost::make_tuple(field_name, type));
468  }
469  template <class MULTIINDEX>
470  typename MULTIINDEX::iterator get_end(const MULTIINDEX &index,
471  const std::string &field_name,
472  const EntityType type) const {
473  return index.upper_bound(boost::make_tuple(field_name, type));
474  }
475 
476 /** \brief loop over all dofs which are on a particular FE row, field and entity
477  * type \ingroup mofem_loops
478  */
479 #define _IT_GET_FEROW_BY_TYPE_DOFS_FOR_LOOP_(FE, NAME, TYPE, IT) \
480  auto IT = FE->get_begin<FENumeredDofEntityByNameAndType>( \
481  FE->rowPtr->get<Composite_Name_And_Type_mi_tag>(), NAME, TYPE); \
482  IT != FE->get_end<FENumeredDofEntityByNameAndType>( \
483  FE->rowPtr->get<Composite_Name_And_Type_mi_tag>(), NAME, TYPE); \
484  IT++
485 
486 /** \brief loop over all dofs which are on a particular FE column, field and
487  * entity type \ingroup mofem_loops
488  */
489 #define _IT_GET_FECOL_BY_TYPE_DOFS_FOR_LOOP_(FE, NAME, TYPE, IT) \
490  auto IT = FE->get_begin<FENumeredDofEntityByNameAndType>( \
491  FE->colPtr->get<Composite_Name_And_Type_mi_tag>(), NAME, TYPE); \
492  IT != FE->get_end<FENumeredDofEntityByNameAndType>( \
493  FE->colPtr->get<Composite_Name_And_Type_mi_tag>(), NAME, TYPE); \
494  IT++
495 
496 /** \brief loop over all dofs which are on a particular FE data, field and
497  * entity type \ingroup mofem_loops
498  */
499 #define _IT_GET_FEDATA_BY_TYPE_DOFS_FOR_LOOP_(FE, NAME, TYPE, IT) \
500  auto IT = FE->get_begin<FEDofEntityByNameAndType>( \
501  FE->dataPtr->get<Composite_Name_And_Type_mi_tag>(), NAME, TYPE); \
502  IT != FE->get_end<FEDofEntityByNameAndType>( \
503  FE->dataPtr->get<Composite_Name_And_Type_mi_tag>(), NAME, TYPE); \
504  IT++
505 
506  template <class MULTIINDEX>
507  typename MULTIINDEX::iterator get_begin(const MULTIINDEX &index,
508  const std::string &field_name) const {
509  return index.lower_bound(field_name);
510  }
511  template <class MULTIINDEX>
512  typename MULTIINDEX::iterator get_end(const MULTIINDEX &index,
513  const std::string &field_name) const {
514  return index.upper_bound(field_name);
515  }
516 
517 /** \brief loop over all dofs which are on a particular FE row and field
518  * \ingroup mofem_loops
519  */
520 #define _IT_GET_FEROW_BY_NAME_DOFS_FOR_LOOP_(FE, NAME, IT) \
521  auto IT = FE->get_begin<FENumeredDofEntityByFieldName>( \
522  FE->rowPtr->get<FieldName_mi_tag>(), NAME); \
523  IT != FE->get_end<FENumeredDofEntityByFieldName>( \
524  FE->rowPtr->get<FieldName_mi_tag>(), NAME); \
525  IT++
526 
527 /** \brief loop over all dofs which are on a particular FE column and field
528  * \ingroup mofem_loops
529  */
530 #define _IT_GET_FECOL_BY_NAME_DOFS_FOR_LOOP_(FE, NAME, IT) \
531  auto IT = FE->get_begin<FENumeredDofEntityByFieldName>( \
532  FE->colPtr->get<FieldName_mi_tag>(), NAME); \
533  IT != FE->get_end<FENumeredDofEntityByFieldName>( \
534  FE->colPtr->get<FieldName_mi_tag>(), NAME); \
535  IT++
536 
537 /** \brief loop over all dofs which are on a particular FE data and field
538  * \ingroup mofem_loops
539  */
540 #define _IT_GET_FEDATA_BY_NAME_DOFS_FOR_LOOP_(FE, NAME, IT) \
541  auto IT = FE->get_begin<FEDofEntityByFieldName>( \
542  FE->dataPtr->get<FieldName_mi_tag>(), NAME); \
543  IT != FE->get_end<FEDofEntityByFieldName>( \
544  FE->dataPtr->get<FieldName_mi_tag>(), NAME); \
545  IT++
546 
547  template <class MULTIINDEX>
548  typename MULTIINDEX::iterator get_begin(const MULTIINDEX &index,
549  const EntityHandle ent) const {
550  return index.lower_bound(ent);
551  }
552  template <class MULTIINDEX>
553  typename MULTIINDEX::iterator get_end(const MULTIINDEX &index,
554  const EntityHandle ent) const {
555  return index.upper_bound(ent);
556  }
557 
558 /** \brief loop over all dofs which are on a particular FE row and given element
559  * entity (handle from moab) \ingroup mofem_loops
560  */
561 #define _IT_GET_FEROW_DOFS_BY_ENT_FOR_LOOP_(FE, ENT, IT) \
562  auto IT = FE->get_begin<FENumeredDofEntityByEnt>( \
563  FE->rowPtr->get<Ent_mi_tag>(), ENT); \
564  IT != FE->get_end<FENumeredDofEntityByEnt>(FE->rowPtr->get<Ent_mi_tag>(), \
565  ENT); \
566  IT++
567 
568 /** \brief loop over all dofs which are on a particular FE column and given
569  * element entity (handle from moab) \ingroup mofem_loops
570  */
571 #define _IT_GET_FECOL_DOFS_BY_ENT_FOR_LOOP_(FE, ENT, IT) \
572  auto IT = FE->get_begin<FENumeredDofEntityByEnt>( \
573  FE->colPtr->get<Ent_mi_tag>(), ENT); \
574  IT != FE->get_end<FENumeredDofEntityByEnt>(FE->colPtr->get<Ent_mi_tag>(), \
575  ENT); \
576  IT++
577 
578 /** \brief loop over all dofs which are on a particular FE data and given
579  * element entity (handle from moab) \ingroup mofem_loops
580  */
581 #define _IT_GET_FEDATA_DOFS_BY_ENT_FOR_LOOP_(FE, ENT, IT) \
582  auto IT = FE->get_begin<FEDofEntity_multiIndex::index<Ent_mi_tag>::type>( \
583  FE->dataPtr->get<Ent_mi_tag>(), ENT); \
584  IT != FE->get_end<FEDofEntity_multiIndex::index<Ent_mi_tag>::type>( \
585  FE->dataPtr->get<Ent_mi_tag>(), ENT); \
586  IT++
587 
588  template <class MULTIINDEX>
589  typename MULTIINDEX::iterator get_begin(const MULTIINDEX &index,
590  const std::string &field_name,
591  const EntityHandle ent) const {
592  return index.lower_bound(boost::make_tuple(field_name, ent));
593  }
594  template <class MULTIINDEX>
595  typename MULTIINDEX::iterator get_end(const MULTIINDEX &index,
596  const std::string &field_name,
597  const EntityHandle ent) const {
598  return index.upper_bound(boost::make_tuple(field_name, ent));
599  }
600 
601 /** \brief loop over all dofs which are on a particular FE row, field and given
602  * element entity (handle from moab) \ingroup mofem_loops
603  */
604 #define _IT_GET_FEROW_DOFS_BY_NAME_AND_ENT_FOR_LOOP_(FE, NAME, ENT, IT) \
605  auto IT = FE->get_begin<FENumeredDofEntityByNameAndEnt>( \
606  FE->rowPtr->get<Composite_Name_And_Ent_mi_tag>(), NAME, ENT); \
607  IT != FE->get_end<FENumeredDofEntityByNameAndEnt>( \
608  FE->rowPtr->get<Composite_Name_And_Ent_mi_tag>(), NAME, ENT); \
609  IT++
610 
611 /** \brief loop over all dofs which are on a particular FE column, field and
612  * given element entity (handle from moab) \ingroup mofem_loops
613  */
614 #define _IT_GET_FECOL_DOFS_BY_NAME_AND_ENT_FOR_LOOP_(FE, NAME, ENT, IT) \
615  auto IT = FE->get_begin<FENumeredDofEntityByNameAndEnt>( \
616  FE->colPtr->get<Composite_Name_And_Ent_mi_tag>(), NAME, ENT); \
617  IT != FE->get_end<FENumeredDofEntityByNameAndEnt>( \
618  FE->colPtr->get<Composite_Name_And_Ent_mi_tag>(), NAME, ENT); \
619  IT++
620 
621 /** \brief loop over all dofs which are on a particular FE data, field and given
622  * element entity (handle from moab) \ingroup mofem_loops
623  */
624 #define _IT_GET_FEDATA_DOFS_BY_NAME_AND_ENT_FOR_LOOP_(FE, NAME, ENT, IT) \
625  auto IT = FE->get_begin<FEDofEntityByNameAndEnt>( \
626  FE->dataPtr->get<Composite_Name_And_Ent_mi_tag>(), NAME, ENT); \
627  IT != FE->get_end<FEDofEntityByNameAndEnt>( \
628  FE->dataPtr->get<Composite_Name_And_Ent_mi_tag>(), NAME, ENT); \
629  IT++
630 };
631 
632 /**
633  * \brief Data structure to exchange data between mofem and User Loop Methods on
634  * entities. \ingroup mofem_loops
635  *
636  * It allows to exchange data between MoFEM and user functions. It stores
637  * information about multi-indices.
638  */
639 struct EntityMethod : public BasicMethod {
640 
642  UnknownInterface **iface) const {
644  if (uuid == IDD_MOFEMEntityMethod) {
645  *iface = const_cast<EntityMethod *>(this);
647  }
648  CHKERR query_interface(uuid, iface);
650  }
651 
652  EntityMethod();
653 
654  boost::shared_ptr<Field> fieldPtr;
655  boost::shared_ptr<FieldEntity> entPtr;
656 };
657 
658 /**
659  * \brief Data structure to exchange data between mofem and User Loop Methods on
660  * entities. \ingroup mofem_loops
661  *
662  * It allows to exchange data between MoFEM and user functions. It stores
663  * information about multi-indices.
664  */
665 struct DofMethod : public BasicMethod {
666 
668  UnknownInterface **iface) const {
670  if (uuid == IDD_MOFEMDofMethod) {
671  *iface = const_cast<DofMethod *>(this);
673  }
674 
675  CHKERR query_interface(uuid, iface);
677  }
678 
679  DofMethod();
680 
681  boost::shared_ptr<Field> fieldPtr;
682  boost::shared_ptr<DofEntity> dofPtr;
683  boost::shared_ptr<NumeredDofEntity> dofNumeredPtr;
684 };
685 
686 /// \deprecated name changed use DofMethod insead EntMethod
688 
689 } // namespace MoFEM
690 
691 #endif // __LOOPMETHODS_HPP__
692 
693 /**
694  * \defgroup mofem_loops Loops
695  * \ingroup mofem
696  */
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.
Vec snes_f
residual
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
Vec ts_u_t
time derivative of state vector
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
PetscReal ts_t
time
MoFEMErrorCode query_interface(const MOFEMuuid &uuid, UnknownInterface **iface) const
Mat snes_A
jacobian matrix
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.
PetscReal ts_a
shift for U_tt (see PETSc Time Solver)
Vec ts_F
residual vector
int sIze
number of processors in communicator
PetscInt ts_step
time step
MoFEMErrorCode setTs(TS _ts)
Set TS solver.
Definition: LoopMethods.cpp:70
TS ts
time solver
multi_index_container< boost::shared_ptr< DofEntity >, indexed_by< ordered_unique< tag< Unique_mi_tag >, member< DofEntity, UId, &DofEntity::globalUId > >, ordered_non_unique< tag< Composite_Name_And_Ent_And_EntDofIdx_mi_tag >, composite_key< DofEntity, const_mem_fun< DofEntity::interface_type_Field, boost::string_ref, &DofEntity::getNameRef >, const_mem_fun< DofEntity, EntityHandle, &DofEntity::getEnt >, const_mem_fun< DofEntity, DofIdx, &DofEntity::getEntDofIdx > > >, ordered_non_unique< tag< Unique_Ent_mi_tag >, const_mem_fun< DofEntity, const UId &, &DofEntity::getEntGlobalUniqueId > >, ordered_non_unique< tag< FieldName_mi_tag >, const_mem_fun< DofEntity::interface_type_Field, boost::string_ref, &DofEntity::getNameRef > >, ordered_non_unique< tag< Ent_mi_tag >, const_mem_fun< DofEntity, EntityHandle, &DofEntity::getEnt > >, ordered_non_unique< tag< Composite_Name_And_Ent_mi_tag >, composite_key< DofEntity, const_mem_fun< DofEntity::interface_type_Field, boost::string_ref, &DofEntity::getNameRef >, const_mem_fun< DofEntity, EntityHandle, &DofEntity::getEnt > > >, ordered_non_unique< tag< Composite_Name_And_Type_mi_tag >, composite_key< DofEntity, const_mem_fun< DofEntity::interface_type_Field, boost::string_ref, &DofEntity::getNameRef >, const_mem_fun< DofEntity::interface_type_RefEntity, EntityType, &DofEntity::getEntType > > > > > DofEntity_multiIndex
MultiIndex container keeps DofEntity.
Vec ts_u_tt
second time derivative of state vector
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
Vec snes_x
state vector
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.
Vec ts_u
state vector
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
SNES snes
snes solver
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
Mat snes_B
preconditioner of jacobian matrix
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.
PetscReal ts_v
shift for U_t shift for U_t
MoFEMErrorCode query_interface(const MOFEMuuid &uuid, UnknownInterface **iface) const
MoFEMErrorCode query_interface(const MOFEMuuid &uuid, UnknownInterface **iface) const
Mat ts_B
Preconditioner for ts_A.
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