v0.14.0
Loading...
Searching...
No Matches
ProblemsMultiIndices.hpp
Go to the documentation of this file.
1/** \file ProblemsMultiIndices.hpp
2 * \brief Multi-index containers, data structures for problems and other
3 * low-level functions
4 */
5
6#ifndef __PROBLEMSMULTIINDICES_HPP__
7#define __PROBLEMSMULTIINDICES_HPP__
8
9namespace MoFEM {
10
11struct Problem;
12
13/**
14 * Data structure created when composite problem is created
15 */
17
18 std::vector<const Problem *> rowProblemsAdd;
19 std::vector<const Problem *> colProblemsAdd;
20
21 std::vector<SmartPetscObj<IS>> rowIs;
22 std::vector<SmartPetscObj<IS>> colIs;
23
24 /**
25 * @brief Get the col sub dm IS
26 *
27 * @param is sub problem IS
28 * @param pp problem number
29 * @return MoFEMErrorCode
30 */
31 inline MoFEMErrorCode getRowIs(IS *is, const unsigned int pp) const;
32
33 /**
34 * @brief Get the Col sub dm IS object
35 *
36 * @param is sub problem IS
37 * @param pp problem number
38 * @return MoFEMErrorCode
39 */
40 inline MoFEMErrorCode getColIs(IS *is, const unsigned int pp) const;
41
42 virtual ~ComposedProblemsData() = default;
43};
44
45/** \brief keeps basic data about problem
46 * \ingroup problems_multi_indices
47 *
48 * This is low level structure with information about problem, what elements
49 * compose problem and what DOFs are on rows and columns.
50 *
51 * \todo fix names following name convention
52 *
53 */
54struct Problem {
55
56 EntityHandle meshset; ///< Problem meshset (on tag of this meshset all data
57 ///< related to problem are stored)
58 BitProblemId *tagId; ///< Unique problem ID
59 const char *tagName; ///< Problem name
60 int tagNameSize; ///< Size of problem name
61 BitFEId *tagBitFEId; ///< IDs of finite elements in problem
62 BitRefLevel *tagBitRefLevel; ///< BitRef level of finite elements in problem
63 BitRefLevel *tagBitRefLevelMask; ///< BItRefMask of elements in problem
64
65 mutable DofIdx nbDofsRow; ///< Global number of DOFs in row
66 mutable DofIdx nbDofsCol; ///< Global number of DOFs in col
67 mutable DofIdx nbLocDofsRow; ///< Local number of DOFs in row
68 mutable DofIdx nbLocDofsCol; ///< Local number of DOFs in colIs
69 mutable DofIdx nbGhostDofsRow; ///< Number of ghost DOFs in row
70 mutable DofIdx nbGhostDofsCol; ///< Number of ghost DOFs in col
71
72 mutable boost::shared_ptr<NumeredDofEntity_multiIndex>
73 numeredRowDofsPtr; ///< store DOFs on rows for this problem
74 mutable boost::shared_ptr<NumeredDofEntity_multiIndex>
75 numeredColDofsPtr; ///< store DOFs on columns for this problem
76 mutable boost::shared_ptr<NumeredEntFiniteElement_multiIndex>
77 numeredFiniteElementsPtr; ///< store finite elements
78
79 /**
80 * \brief get access to numeredRowDofsPtr storing DOFs on rows
81 */
82 inline auto &getNumeredRowDofsPtr() const { return numeredRowDofsPtr; }
83
84 /**
85 * \brief get access to numeredColDofsPtr storing DOFs on cols
86 */
87 inline auto &getNumeredColDofsPtr() const { return numeredColDofsPtr; }
88
89 /**
90 * \brief get access to reference for multi-index storing finite elements
91 */
92 inline const auto &getNumeredFiniteElementsPtr() const {
94 }
95
96 /**
97 * @brief Erase elements by entities
98 *
99 * @param entities
100 * @return MoFEMErrorCode
101 */
102 MoFEMErrorCode eraseElements(Range entities) const;
103
104 /**
105 * \brief Subproblem problem data
106 */
107 struct SubProblemData;
108
109 /**
110 * Pointer to data structure. This pointer has allocated data only for
111 * sub problems.
112 */
113 mutable boost::shared_ptr<SubProblemData> subProblemData;
114
115 /**
116 * \brief Get main problem of sub-problem is
117 * @return sub problem data structure
118 */
119 inline boost::shared_ptr<SubProblemData> &getSubData() const {
120 return subProblemData;
121 }
122
123 /**
124 * Pointer to data structure from which this problem is composed
125 */
126 mutable boost::shared_ptr<ComposedProblemsData> composedProblemsData;
127
128 /**
129 * \brief Het composed problems data structure
130 */
131 inline auto &getComposedProblemsData() const { return composedProblemsData; }
132
133 /**
134 * \brief get DOFs from problem
135 *
136 * Note that \e ent_dof_idx is not coefficient number, is local number of DOFs
137 * on the entity. The coefficient number and local index of DOFs or entity are
138 * the same on vertices and H1 approximation.
139 *
140 * @param field_bit_number field name field_bit_number = (use
141 * m_field.get_field_bit_number(field_name);
142 * @param ent entity handle
143 * @param ent_dof_idx index of DOFs on entity
144 * @param row_or_col ROW or COL
145 * @param dof_ptr shared pointer to DOFs if found
146 * @return error code
147 */
149 const int field_bit_number, const EntityHandle ent, const int ent_dof_idx,
150 const RowColData row_or_col,
151 boost::shared_ptr<NumeredDofEntity> &dof_ptr) const;
152
153/**
154 * \brief use with loops to iterate row DOFs
155 * \ingroup problems_multi_indices
156 *
157 * \code
158 * for(_IT_NUMEREDDOF_ROW_FOR_LOOP_(PROBLEMPTR,IT)) {
159 * ...
160 * }
161 * \endcode
162 *
163 */
164#define _IT_NUMEREDDOF_ROW_FOR_LOOP_(PROBLEMPTR, IT) \
165 NumeredDofEntity_multiIndex::iterator IT = \
166 PROBLEMPTR->getNumeredRowDofsBegin(); \
167 IT != PROBLEMPTR->getNumeredRowDofsEnd(); \
168 IT++
169
170/**
171 * use with loops to iterate col DOFs
172 * \ingroup problems_multi_indices
173 *
174 * \code
175 * for(_IT_NUMEREDDOF_COL_FOR_LOOP_(PROBLEMPTR,IT)) {
176 * ...
177 * }
178 * \endcode
179 *
180 */
181#define _IT_NUMEREDDOF_COL_FOR_LOOP_(PROBLEMPTR, IT) \
182 NumeredDofEntity_multiIndex::iterator IT = \
183 PROBLEMPTR->getNumeredColDofsBegin(); \
184 IT != PROBLEMPTR->getNumeredColDofsEnd(); \
185 IT++
186
187 /// get begin iterator for numeredRowDofsPtr (insted you can use
188 /// #_IT_NUMEREDDOFMOFEMENTITY_ROW_FOR_LOOP_ for loops)
189 NumeredDofEntity_multiIndex::iterator getNumeredRowDofsBegin() const {
190 return numeredRowDofsPtr->begin();
191 }
192
193 /// get end iterator for numeredRowDofsPtr (insted you can use
194 /// #_IT_NUMEREDDOFMOFEMENTITY_ROW_FOR_LOOP_ for loops)
195 NumeredDofEntity_multiIndex::iterator getNumeredRowDofsEnd() const {
196 return numeredRowDofsPtr->end();
197 }
198
199 /// get begin iterator for numeredColDofsPtr (insted you can use
200 /// #_IT_NUMEREDDOFMOFEMENTITY_COL_FOR_LOOP_ for loops)
201 NumeredDofEntity_multiIndex::iterator getNumeredColDofsBegin() const {
202 return numeredColDofsPtr->begin();
203 }
204
205 /// get end iterator for numeredColDofsPtr (insted you can use
206 /// #_IT_NUMEREDDOFMOFEMENTITY_COL_FOR_LOOP_ for loops)
207 NumeredDofEntity_multiIndex::iterator getNumeredColDofsEnd() const {
208 return numeredColDofsPtr->end();
209 }
210
211/**
212 * \brief use with loops to iterate row DOFs
213 * \ingroup problems_multi_indices
214 *
215 * \code
216 * for(_IT_NUMEREDDOF_ROW_BY_LOCIDX_FOR_LOOP_(PROBLEMPTR,IT)) {
217 * ...
218 * }
219 * \endcode
220 *
221 */
222#define _IT_NUMEREDDOF_ROW_BY_LOCIDX_FOR_LOOP_(PROBLEMPTR, IT) \
223 NumeredDofEntityByLocalIdx::iterator IT = \
224 PROBLEMPTR->getNumeredRowDofsByLocIdxBegin(0); \
225 IT != PROBLEMPTR->getNumeredRowDofsByLocIdxEnd( \
226 PROBLEMPTR->getNbLocalDofsRow() - 1); \
227 IT++
228
229/**
230 * \brief use with loops to iterate col DOFs
231 *
232 * \code
233 * for(_IT_NUMEREDDOF_COL_BY_LOCIDX_FOR_LOOP_(PROBLEMPTR,IT)) {
234 * ...
235 * }
236 * \endcode
237 *
238 */
239#define _IT_NUMEREDDOF_COL_BY_LOCIDX_FOR_LOOP_(PROBLEMPTR, IT) \
240 NumeredDofEntityByUId::iterator IT = \
241 PROBLEMPTR->getNumeredColDofsByLocIdxBegin(0); \
242 IT != PROBLEMPTR->getNumeredColDofsByLocIdxEnd( \
243 PROBLEMPTR->getNbLocalDofsRow() - 1); \
244 IT++
245
246 /// get begin iterator for numeredRowDofsPtr (insted you can use
247 /// #_IT_NUMEREDDOF_ROW_FOR_LOOP_ for loops)
248 auto getNumeredRowDofsByLocIdxBegin(const DofIdx locidx) const {
249 return numeredRowDofsPtr->get<PetscLocalIdx_mi_tag>().lower_bound(locidx);
250 }
251
252 /// get end iterator for numeredRowDofsPtr (insted you can use
253 /// #_IT_NUMEREDDOF_ROW_FOR_LOOP_ for loops)
254 auto getNumeredRowDofsByLocIdxEnd(const DofIdx locidx) const {
255 return numeredRowDofsPtr->get<PetscLocalIdx_mi_tag>().upper_bound(locidx);
256 }
257
258 /// get begin iterator for numeredColDofsPtr (insted you can use
259 /// #_IT_NUMEREDDOF_COL_FOR_LOOP_ for loops)
260 auto getNumeredColDofsByLocIdxBegin(const DofIdx locidx) const {
261 return numeredColDofsPtr->get<PetscLocalIdx_mi_tag>().lower_bound(locidx);
262 }
263
264 /// get end iterator for numeredColDofsPtr (insted you can use
265 /// #_IT_NUMEREDDOF_COL_FOR_LOOP_ for loops)
266 auto getNumeredColDofsByLocIdxEnd(const DofIdx locidx) const {
267 return numeredColDofsPtr->get<PetscLocalIdx_mi_tag>().upper_bound(locidx);
268 }
269
270/**
271 * \brief use with loops to iterate row DOFs
272 * \ingroup problems_multi_indices
273 *
274 * \code
275 * for(_IT_NUMEREDDOF_BY_ENT_ROW_FOR_LOOP_(PROBLEMPTR,ENT,IT)) {
276 * ...
277 * }
278 * \endcode
279 *
280 */
281#define _IT_NUMEREDDOF_ROW_BY_ENT_FOR_LOOP_(PROBLEMPTR, ENT, IT) \
282 auto IT = PROBLEMPTR->getNumeredRowDofsByEntBegin(ENT); \
283 IT != PROBLEMPTR->getNumeredRowDofsByEntEnd(ENT); \
284 IT++
285
286/**
287 * \brief use with loops to iterate col DOFs
288 * \ingroup problems_multi_indices
289 *
290 * \code
291 * for(_IT_NUMEREDDOF_COL_BY_ENT_FOR_LOOP_(PROBLEMPTR,ENT,IT)) {
292 * ...
293 * }
294 * \endcode
295 *
296 */
297#define _IT_NUMEREDDOF_COL_BY_ENT_FOR_LOOP_(PROBLEMPTR, ENT, IT) \
298 auto IT = PROBLEMPTR->getNumeredColDofsByEntBegin(ENT); \
299 IT != PROBLEMPTR->getNumeredColDofsByEntEnd(ENT); \
300 IT++
301
302 /// get begin iterator for numeredRowDofsPtr (insted you can use
303 /// #_IT_NUMEREDDOF_ROW_BY_ENT_FOR_LOOP_ for loops)
305 return numeredRowDofsPtr->get<Ent_mi_tag>().lower_bound(ent);
306 }
307
308 /// get end iterator for numeredRowDofsPtr (insted you can use
309 /// #_IT_NUMEREDDOF_ROW_BY_ENT_FOR_LOOP_ for loops)
311 return numeredRowDofsPtr->get<Ent_mi_tag>().upper_bound(ent);
312 }
313
314 /// get begin iterator for numeredColDofsPtr (insted you can use
315 /// #_IT_NUMEREDDOF_COL_BY_ENT_FOR_LOOP_ for loops)
317 return numeredColDofsPtr->get<Ent_mi_tag>().lower_bound(ent);
318 }
319
320 /// get end iterator for numeredColDofsPtr (insted you can use
321 /// #_IT_NUMEREDDOF_COL_BY_ENT_FOR_LOOP_ for loops)
323 return numeredColDofsPtr->get<Ent_mi_tag>().upper_bound(ent);
324 }
325
326/**
327 * use with loops to iterate row DOFs
328 * \ingroup problems_multi_indices
329 *
330 * \code
331 * for(_IT_NUMEREDDOF_BY_NAME_ROW_FOR_LOOP_(PROBLEMPTR,m_field.get_field_bit_number(FIELD_BIT_NUMBER),IT))
332 * {
333 * ...
334 * }
335 * \endcode
336 *
337 */
338#define _IT_NUMEREDDOF_ROW_BY_BITNUMBER_FOR_LOOP_(PROBLEMPTR, \
339 FIELD_BIT_NUMBER, IT) \
340 auto IT = PROBLEMPTR->numeredRowDofsPtr->lower_bound( \
341 FieldEntity::getLoBitNumberUId(FIELD_BIT_NUMBER)); \
342 IT != PROBLEMPTR->numeredRowDofsPtr->upper_bound( \
343 FieldEntity::getHiBitNumberUId(FIELD_BIT_NUMBER)); \
344 IT++
345
346/**
347 * \brief use with loops to iterate col DOFs
348 * \ingroup problems_multi_indices
349 *
350 * \code
351 * for(_IT_NUMEREDDOF_COL_BY_BITNUMBER_FOR_LOOP_(PROBLEMPTR,m_field.get_field_bit_number(FIELD_BIT_NUMBER),IT))
352 * {
353 * ...
354 * }
355 * \endcode
356 *
357 */
358#define _IT_NUMEREDDOF_COL_BY_BITNUMBER_FOR_LOOP_(PROBLEMPTR, \
359 FIELD_BIT_NUMBER, IT) \
360 auto IT = PROBLEMPTR->numeredColDofsPtr->lower_bound( \
361 FieldEntity::getLoBitNumberUId(FIELD_BIT_NUMBER)); \
362 IT != PROBLEMPTR->numeredColDofsPtr->upper_bound( \
363 FieldEntity::getHiBitNumberUId(FIELD_BIT_NUMBER)); \
364 IT++
365
366 Problem(Interface &moab, const EntityHandle meshset);
367
368 virtual ~Problem() = default;
369
370 inline BitProblemId getId() const { return *((BitProblemId *)tagId); }
371
372 inline auto getName() const {
373 return std::string((char *)tagName, tagNameSize);
374 }
375
376 inline DofIdx getNbDofsRow() const { return nbDofsRow; }
377 inline DofIdx getNbDofsCol() const { return nbDofsCol; }
378 inline DofIdx getNbLocalDofsRow() const { return nbLocDofsRow; }
379 inline DofIdx getNbLocalDofsCol() const { return nbLocDofsCol; }
380 inline DofIdx getNbGhostDofsRow() const { return nbGhostDofsRow; }
381 inline DofIdx getNbGhostDofsCol() const { return nbGhostDofsCol; }
382
383 inline BitRefLevel getBitRefLevel() const { return *tagBitRefLevel; }
385
386 // DEPRECATED inline BitRefLevel getMaskBitRefLevel() const {
387 // return *tagBitRefLevelMask;
388 // }
389
390 /**
391 * @brief Get the Row Dofs By Petsc Global Dof Idx object
392 *
393 * @param idx
394 * @param dof_ptr
395 * @param bh
396 * @return MoFEMErrorCode
397 */
399 const NumeredDofEntity **dof_ptr,
400 MoFEMTypes bh = MF_EXIST) const;
401
402 /**
403 * @brief Get the Col Dofs By Petsc Global Dof Idx object
404 *
405 * @param idx
406 * @param dof_ptr
407 * @param bh
408 * @return MoFEMErrorCode
409 */
411 const NumeredDofEntity **dof_ptr,
412 MoFEMTypes bh = MF_EXIST) const;
413
414 /**
415 * @brief Get the Row Dofs By Petsc Global Dof Idx object
416 *
417 * @param idx
418 * @return boost::weak_ptr<NumeredDofEntity>
419 */
420 boost::weak_ptr<NumeredDofEntity>
422
423 /**
424 * @brief Get the Col Dofs By Petsc Global Dof Idx object
425 *
426 * @param idx
427 * @return boost::weak_ptr<NumeredDofEntity>
428 */
429 boost::weak_ptr<NumeredDofEntity>
431
432 BitFEId getBitFEId() const;
433
434 friend std::ostream &operator<<(std::ostream &os, const Problem &e);
435
436 /**
437 * \brief Get number of finite elements by name on processors
438 *
439 * Size retuned IS is equal to size of processors.
440 *
441 * What PetscLayout see,
442 <http://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/IS/index.html>
443 *
444 * Example of usage for layout
445 * \code
446
447 PetscInt rstart, rend;
448 ierr = PetscLayoutGetRange(layout, &rstart, &rend); CHKERRG(ierr);
449 int global_size;
450 ierr = PetscLayoutGetSize(layout,&global_size); CHKERRG(ierr);
451
452 \endcode
453 *
454 * @param comm Communicator
455 * @param name Finite element name
456 * @param layout Get number of elements on each processor
457 * @return error code
458 */
460 const std::string name,
461 PetscLayout *layout) const;
462
463 /**
464 * \brief Get number of finite elements on processors
465 *
466 * Size retuned IS is equal to size of processors.
467 *
468 * What PetscLayout see,
469 <http://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/IS/index.html>
470 *
471 * Example of usage for layout
472 * \code
473
474 PetscInt rstart, rend;
475 ierr = PetscLayoutGetRange(layout, &rstart, &rend); CHKERRG(ierr);
476 int global_size;
477 ierr = PetscLayoutGetSize(layout,&global_size); CHKERRG(ierr);
478
479 \endcode
480 *
481 * @param comm Communicator
482 * @param layout Get number of elements on each processor
483 * @return error code
484 */
486 PetscLayout *layout) const;
487
488 typedef multi_index_container<boost::weak_ptr<std::vector<NumeredDofEntity>>,
489 indexed_by<sequenced<>>>
491
492 /**
493 * \brief Get reference to sequence data numbered dof container
494 *
495 * In sequence data container data are physically stored. The purpose of this
496 * is to allocate NumeredDofEntity data in bulk, having only one allocation
497 instead
498 * each time entity is inserted. That makes code efficient.
499 *
500 * The vector in sequence is destroyed if last entity inside that vector is
501 * destroyed. All MoFEM::NumeredDofEntity have aliased shared_ptr which points
502 to the vector.
503
504 * @return MoFEM::Problem::SequenceDofContainer
505 */
506 inline auto &getRowDofsSequence() const { return sequenceRowDofContainer; }
507
508 /**
509 * \brief Get reference to sequence data numbered dof container
510 *
511 * In sequence data container data are physically stored. The purpose of this
512 * is to allocate NumeredDofEntity data in bulk, having only one allocation
513 instead
514 * each time entity is inserted. That makes code efficient.
515 *
516 * The vector in sequence is destroyed if last entity inside that vector is
517 * destroyed. All MoFEM::NumeredDofEntity have aliased shared_ptr which points
518 to the vector.
519
520 * @return MoFEM::Problem::SequenceDofContainer
521 */
522 inline auto &getColDofsSequence() const { return sequenceColDofContainer; }
523
524 using EmptyFieldBlocks = std::pair<BitFieldId, BitFieldId>;
525
526 /**
527 * @brief Get the empty field blocks
528 *
529 * Emtpy field blocks is a pair contains IDs of the fields for which matrix
530 * has zero entries.
531 *
532 * @return EmptyFieldBlocks&
533 */
535 return emptyFieldBlocks;
536 }
537
538 /**
539 * @brief Add fields to the empty field blocks
540 *
541 * Emtpy field blocks is a pair contains IDs of the fields for which matrix
542 * has zero entries.
543 *
544 * @param add_fields
545 * @return EmptyFieldBlocks&
546 */
547 inline EmptyFieldBlocks &
549 emptyFieldBlocks.first |= add_fields.first;
550 emptyFieldBlocks.second |= add_fields.second;
551 return emptyFieldBlocks;
552 }
553
554private:
555 // Keep vector of DoFS on entity
556 mutable boost::shared_ptr<SequenceDofContainer> sequenceRowDofContainer;
557 mutable boost::shared_ptr<SequenceDofContainer> sequenceColDofContainer;
558
560};
561
563
564/**
565 * \brief Subproblem problem data
566 */
568
570 rowIs; ///< indices of main problem of which sub problem is this
572 colIs; ///< indices of main problem of which sub problem is this
574 rowMap; ///< mapping form main problem indices to sub-problem indices
576 colMap; ///< mapping form main problem indices to sub-problem indices
577
578 inline auto getSmartRowIs() { return rowIs; }
579 inline auto getSmartColIs() { return colIs; }
580 inline auto getSmartRowMap() { return rowMap; }
581 inline auto getSmartColMap() { return colMap; }
582
583 /**
584 * get row Is for sub problem
585 * @param is create is
586 * @return error code
587 */
588 inline MoFEMErrorCode getRowIs(IS *is) {
590 *is = rowIs;
591 PetscObjectReference((PetscObject)(*is));
593 }
594
595 /**
596 * get col Is for sub problem
597 * @param is create is
598 * @return error code
599 */
600 inline MoFEMErrorCode getColIs(IS *is) {
602 *is = colIs;
603 PetscObjectReference((PetscObject)(*is));
605 };
606
607 /**
608 * get row AO mapping for sub problem
609 * @param ao get mapping
610 * @return error code
611 */
612 inline MoFEMErrorCode getRowMap(AO *ao) {
614 *ao = rowMap;
615 PetscObjectReference((PetscObject)(*ao));
617 }
618
619 /**
620 * get col AO mapping for sub problem
621 * @param ao get mapping
622 * @return error code
623 */
624 inline MoFEMErrorCode getColMap(AO *ao) {
626 *ao = colMap;
627 PetscObjectReference((PetscObject)(*ao));
629 }
630
631 SubProblemData() = default;
632 virtual ~SubProblemData() = default;
633};
634
636
637/**
638 * @relates multi_index_container
639 * \brief MultiIndex for entities for Problem
640 * \ingroup fe_multi_indices
641 */
642typedef multi_index_container<
643 Problem,
644 indexed_by<
645 ordered_unique<tag<Meshset_mi_tag>,
646 member<Problem, EntityHandle, &Problem::meshset>>,
647 hashed_unique<tag<BitProblemId_mi_tag>,
648 const_mem_fun<Problem, BitProblemId, &Problem::getId>,
650 hashed_unique<tag<Problem_mi_tag>,
651 const_mem_fun<Problem, std::string, &Problem::getName>>>>
653
654/** \brief add ref level to problem
655 * \ingroup problems_multi_indices
656 */
660 void operator()(Problem &p) { *(p.tagBitRefLevel) |= bit; };
661};
662
663/** \brief set ref level to problem
664 * \ingroup problems_multi_indices
665 */
669 void operator()(Problem &p) { *(p.tagBitRefLevel) = bit; };
670};
671
672/** \brief set prof dof bit ref mask
673 * \ingroup problems_multi_indices
674 */
678 void operator()(Problem &p) { *(p.tagBitRefLevelMask) = bit; };
679};
680
681/** \brief add finite element to problem
682 * \ingroup problems_multi_indices
683 */
687 void operator()(Problem &p);
688};
689
690/** \brief set prof dof bit ref mask
691 * \ingroup problems_multi_indices
692 */
696 void operator()(Problem &p) { *(p.tagBitRefLevelMask) |= bit; };
697};
698
699/** \brief remove finite element from problem
700 * \ingroup problems_multi_indices
701 */
705 void operator()(Problem &p);
706};
707
708/** \brief zero nb. of DOFs in row
709 * \ingroup problems_multi_indices
710 */
712 void operator()(Problem &e);
713};
714
715/** \brief zero nb. of DOFs in col
716 * \ingroup problems_multi_indices
717 */
719 void operator()(Problem &e);
720};
721
722/** \brief clear problem finite elements
723 * \ingroup problems_multi_indices
724 */
726 void operator()(Problem &e);
727};
728
729/**
730 * \brief Clear sub-problem data structure
731 */
733 void operator()(Problem &e) { e.subProblemData.reset(); }
734};
735
736/**
737 * \brief Clear composed problem data structure
738 */
741};
742
743inline MoFEMErrorCode
744ComposedProblemsData::getRowIs(IS *is, const unsigned int pp) const {
746 if (pp <= rowIs.size()) {
747 SETERRQ1(PETSC_COMM_WORLD, MOFEM_INVALID_DATA, "Exceed size of array pp<%d",
748 rowIs.size());
749 }
750 *is = rowIs[pp].get();
751 PetscObjectReference((PetscObject)(*is));
753}
754
755inline MoFEMErrorCode
756ComposedProblemsData::getColIs(IS *is, const unsigned int pp) const {
758 if (pp <= colIs.size()) {
759 SETERRQ1(PETSC_COMM_WORLD, MOFEM_INVALID_DATA, "Exceed size of array pp<%d",
760 colIs.size());
761 }
762 *is = colIs[pp].get();
763 PetscObjectReference((PetscObject)(*is));
765}
766
767} // namespace MoFEM
768
769#endif //__PROBLEMSMULTIINDICES_HPP__
770
771/**
772 * \defgroup problems_multi_indices Problems structures and multi-indices
773 * \ingroup mofem
774 ******************************************************************************/
static Index< 'p', 3 > p
RowColData
RowColData.
Definition: definitions.h:123
MoFEMTypes
Those types control how functions respond on arguments, f.e. error handling.
Definition: definitions.h:97
@ MF_EXIST
Definition: definitions.h:100
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:447
@ MOFEM_INVALID_DATA
Definition: definitions.h:36
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:440
multi_index_container< Problem, indexed_by< ordered_unique< tag< Meshset_mi_tag >, member< Problem, EntityHandle, &Problem::meshset > >, hashed_unique< tag< BitProblemId_mi_tag >, const_mem_fun< Problem, BitProblemId, &Problem::getId >, HashBit< BitProblemId >, EqBit< BitProblemId > >, hashed_unique< tag< Problem_mi_tag >, const_mem_fun< Problem, std::string, &Problem::getName > > > > Problem_multiIndex
MultiIndex for entities for Problem.
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:56
std::bitset< BITFEID_SIZE > BitFEId
Finite element Id.
Definition: Types.hpp:43
int DofIdx
Index of DOF.
Definition: Types.hpp:18
std::bitset< BITPROBLEMID_SIZE > BitProblemId
Problem Id.
Definition: Types.hpp:44
std::bitset< BITREFLEVEL_SIZE > BitRefLevel
Bit structure attached to each entity identifying to what mesh entity is attached.
Definition: Types.hpp:40
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:10
Problem::EmptyFieldBlocks EmptyFieldBlocks
DeprecatedCoreInterface Interface
Definition: Interface.hpp:1975
std::vector< const Problem * > rowProblemsAdd
std::vector< SmartPetscObj< IS > > colIs
virtual ~ComposedProblemsData()=default
MoFEMErrorCode getRowIs(IS *is, const unsigned int pp) const
Get the col sub dm IS.
MoFEMErrorCode getColIs(IS *is, const unsigned int pp) const
Get the Col sub dm IS object.
std::vector< SmartPetscObj< IS > > rowIs
std::vector< const Problem * > colProblemsAdd
keeps information about indexed dofs for the problem
SmartPetscObj< AO > rowMap
mapping form main problem indices to sub-problem indices
virtual ~SubProblemData()=default
SmartPetscObj< IS > rowIs
indices of main problem of which sub problem is this
SmartPetscObj< AO > colMap
mapping form main problem indices to sub-problem indices
SmartPetscObj< IS > colIs
indices of main problem of which sub problem is this
ProblemChangeRefLevelBitAdd(const BitRefLevel _bit)
ProblemChangeRefLevelBitDofMaskAdd(const BitRefLevel _bit)
ProblemChangeRefLevelBitDofMaskSet(const BitRefLevel _bit)
ProblemChangeRefLevelBitSet(const BitRefLevel _bit)
Clear composed problem data structure.
Clear sub-problem data structure.
remove finite element from problem
keeps basic data about problem
MoFEMErrorCode getNumberOfElementsByPart(MPI_Comm comm, PetscLayout *layout) const
Get number of finite elements on processors.
auto getNumeredRowDofsByEntBegin(const EntityHandle ent) const
NumeredDofEntity_multiIndex::iterator getNumeredRowDofsBegin() const
auto getNumeredRowDofsByLocIdxBegin(const DofIdx locidx) const
auto getNumeredColDofsByLocIdxEnd(const DofIdx locidx) const
DofIdx getNbLocalDofsRow() const
MoFEMErrorCode getNumberOfElementsByNameAndPart(MPI_Comm comm, const std::string name, PetscLayout *layout) const
Get number of finite elements by name on processors.
DofIdx getNbDofsRow() const
const char * tagName
Problem name.
BitRefLevel getBitRefLevel() const
boost::shared_ptr< SequenceDofContainer > sequenceColDofContainer
BitProblemId * tagId
Unique problem ID.
auto & getColDofsSequence() const
Get reference to sequence data numbered dof container.
std::pair< BitFieldId, BitFieldId > EmptyFieldBlocks
MoFEMErrorCode getDofByNameEntAndEntDofIdx(const int field_bit_number, const EntityHandle ent, const int ent_dof_idx, const RowColData row_or_col, boost::shared_ptr< NumeredDofEntity > &dof_ptr) const
get DOFs from problem
virtual ~Problem()=default
auto getNumeredColDofsByEntEnd(const EntityHandle ent) const
BitRefLevel getBitRefLevelMask() const
auto getNumeredRowDofsByEntEnd(const EntityHandle ent) const
EmptyFieldBlocks & getEmptyFieldBlocks() const
Get the empty field blocks.
boost::shared_ptr< SubProblemData > & getSubData() const
Get main problem of sub-problem is.
DofIdx nbGhostDofsCol
Number of ghost DOFs in col.
boost::shared_ptr< NumeredDofEntity_multiIndex > numeredRowDofsPtr
store DOFs on rows for this problem
EmptyFieldBlocks & addFieldToEmptyFieldBlocks(const EmptyFieldBlocks add_fields) const
Add fields to the empty field blocks.
boost::shared_ptr< NumeredEntFiniteElement_multiIndex > numeredFiniteElementsPtr
store finite elements
NumeredDofEntity_multiIndex::iterator getNumeredRowDofsEnd() const
auto getNumeredColDofsByLocIdxBegin(const DofIdx locidx) const
friend std::ostream & operator<<(std::ostream &os, const Problem &e)
boost::shared_ptr< SequenceDofContainer > sequenceRowDofContainer
MoFEMErrorCode getRowDofsByPetscGlobalDofIdx(DofIdx idx, const NumeredDofEntity **dof_ptr, MoFEMTypes bh=MF_EXIST) const
Get the Row Dofs By Petsc Global Dof Idx object.
boost::shared_ptr< ComposedProblemsData > composedProblemsData
MoFEMErrorCode eraseElements(Range entities) const
Erase elements by entities.
const auto & getNumeredFiniteElementsPtr() const
get access to reference for multi-index storing finite elements
DofIdx getNbGhostDofsCol() const
EmptyFieldBlocks emptyFieldBlocks
auto getNumeredColDofsByEntBegin(const EntityHandle ent) const
DofIdx nbLocDofsRow
Local number of DOFs in row.
BitFEId getBitFEId() const
BitProblemId getId() const
DofIdx nbDofsCol
Global number of DOFs in col.
auto & getComposedProblemsData() const
Het composed problems data structure.
DofIdx getNbDofsCol() const
DofIdx getNbLocalDofsCol() const
DofIdx nbGhostDofsRow
Number of ghost DOFs in row.
auto & getNumeredColDofsPtr() const
get access to numeredColDofsPtr storing DOFs on cols
boost::shared_ptr< SubProblemData > subProblemData
DofIdx nbLocDofsCol
Local number of DOFs in colIs.
int tagNameSize
Size of problem name.
MoFEMErrorCode getColDofsByPetscGlobalDofIdx(DofIdx idx, const NumeredDofEntity **dof_ptr, MoFEMTypes bh=MF_EXIST) const
Get the Col Dofs By Petsc Global Dof Idx object.
BitFEId * tagBitFEId
IDs of finite elements in problem.
BitRefLevel * tagBitRefLevelMask
BItRefMask of elements in problem.
NumeredDofEntity_multiIndex::iterator getNumeredColDofsEnd() const
DofIdx nbDofsRow
Global number of DOFs in row.
DofIdx getNbGhostDofsRow() const
multi_index_container< boost::weak_ptr< std::vector< NumeredDofEntity > >, indexed_by< sequenced<> > > SequenceDofContainer
boost::shared_ptr< NumeredDofEntity_multiIndex > numeredColDofsPtr
store DOFs on columns for this problem
BitRefLevel * tagBitRefLevel
BitRef level of finite elements in problem.
NumeredDofEntity_multiIndex::iterator getNumeredColDofsBegin() const
auto getNumeredRowDofsByLocIdxEnd(const DofIdx locidx) const
auto & getNumeredRowDofsPtr() const
get access to numeredRowDofsPtr storing DOFs on rows
auto & getRowDofsSequence() const
Get reference to sequence data numbered dof container.
intrusive_ptr for managing petsc objects