v0.13.2
Loading...
Searching...
No Matches
MeshProjectionDataOperators.hpp
Go to the documentation of this file.
1/** \file MeshProjectionDataOperators.hpp
2 * \brief Mesh projection operators
3
4*/
5
6#ifndef __MESH_PROJECTION_DATA_OPERATORS_HPP__
7#define __MESH_PROJECTION_DATA_OPERATORS_HPP__
8
9namespace MoFEM {
10
11/**
12 * @brief Operator to execute finite element instance on parent element.
13 * This operator is typically used to project field from parent to child, or
14 * vice versa. It enables to evaluate filed data of parent entity on chile
15 * entity integration points.
16 */
18
19 /**
20 * @brief Construct a new Op Run Parent object
21 *
22 * @note Finite element instance usually has to be class which has overloaded
23 * method from projecting integration points from child tp parent.
24 *
25 * @note Typically parent_ele_ptr and bit_this_mask is the same instance
26 *
27 * @param parent_ele_ptr finite element instance executed on parent entity
28 * @param bit_parent bit of parent entity
29 * @param bit_parent_mask mask of parent entity
30 * @param this_ele_ptr "this" element instance
31 * @param bit_this bit of entity on which "this" finite element is executed
32 * @param bit_this_mask mask of entity on which "this" finite element instance
33 * is executed
34 * @param verb verbosity level
35 * @param sev logging severity level
36 */
37 OpRunParent(boost::shared_ptr<ForcesAndSourcesCore> parent_ele_ptr,
38 BitRefLevel bit_parent, BitRefLevel bit_parent_mask,
39 boost::shared_ptr<ForcesAndSourcesCore> this_ele_ptr,
40 BitRefLevel bit_this, BitRefLevel bit_this_mask, int verb = QUIET,
41 Sev sev = Sev::noisy);
42
43 MoFEMErrorCode doWork(int side, EntityType type,
45
46private:
47 boost::shared_ptr<ForcesAndSourcesCore> parentElePtr;
48 boost::shared_ptr<ForcesAndSourcesCore> thisElePtr;
55};
56
57/**
58 * @brief Operator to project base functions from parent entity to child
59 *
60 * This operator project base functions, field data (i.e. indices, field values
61 * of dofs, etc.), into parent element. Operator can be called as a hierarchy to
62 * get access to information on lower refinement levels.
63 *
64 */
66
67 /**
68 * @brief Construct a new Op Add Parent Ent Data object
69 *
70 * @param field_name field name DOFs projected from parent
71 * @param op_parent_type type of user data operator
72 * @param parent_ele_ptr parent finite element instance
73 * @param bit_child bit of child entity
74 * @param bit_child_mask bit mask of child
75 * @param bit_parent_ent bit of parent entity
76 * @param bit_parent_ent_mask bit mask of parent
77 * @param verb verbosity level
78 * @param sev severity level for logging
79 */
80 OpAddParentEntData(std::string field_name, OpType op_parent_type,
81 boost::shared_ptr<ForcesAndSourcesCore> parent_ele_ptr,
82 BitRefLevel bit_child, BitRefLevel bit_child_mask,
83 BitRefLevel bit_parent_ent,
84 BitRefLevel bit_parent_ent_mask, int verb = QUIET,
85 Sev sev = Sev::noisy);
86
87 /**
88 * @brief Construct a new Op Add Parent Ent Data object
89 *
90 * @param space field space
91 * @param op_parent_type type of user data operator
92 * @param parent_ele_ptr parent finite element instance
93 * @param bit_child bit of child entity
94 * @param bit_child_mask bit mask of child
95 * @param bit_parent_ent bit of parent entity
96 * @param bit_parent_ent_mask bit mask of parent
97 * @param verb verbosity level
98 * @param sev severity level for logging
99 */
100 OpAddParentEntData(FieldSpace space, OpType op_parent_type,
101 boost::shared_ptr<ForcesAndSourcesCore> parent_ele_ptr,
102 BitRefLevel bit_child, BitRefLevel bit_child_mask,
103 BitRefLevel bit_parent_ent,
104 BitRefLevel bit_parent_ent_mask, int verb = QUIET,
105 Sev sev = Sev::noisy);
106
108 const bool error_if_no_base = false);
109
110private:
111 std::string fieldName;
114 boost::shared_ptr<ForcesAndSourcesCore> parentElePtr;
121
122 boost::ptr_deque<EntitiesFieldData::EntData> poolEntsVector;
123};
124
125/**
126 * @brief Create adjacency to parent elements.
127 *
128 * That class is used during entity finite element construction.
129 *
130 * @tparam DIM dimension of parent element
131 */
132template <int DIM> struct ParentFiniteElementAdjacencyFunction {
133
135 BitRefLevel bit_parent_mask,
136 BitRefLevel bit_ent,
137 BitRefLevel bit_ent_mask)
138 : bitParent(bit_parent), bitParentMask(bit_parent_mask), bitEnt(bit_ent),
139 bitEntMask(bit_ent_mask) {}
140
141 /**
142 * @brief Function setting adjacencies to DOFs of parent element
143 *
144 * @note elements form child, see dofs from parent, so DOFs located on
145 * adjacencies of parent entity has adjacent to dofs of child.
146 *
147 * @tparam DIM dimension of the element entity
148 * @param moab
149 * @param field
150 * @param fe
151 * @param adjacency
152 * @return MoFEMErrorCode
153 */
154 MoFEMErrorCode operator()(moab::Interface &moab, const Field &field,
155 const EntFiniteElement &fe,
156 std::vector<EntityHandle> &adjacency) {
158
159 static_assert(DIM >= 0 && DIM <= 3, "DIM is out of scope");
160
161 adjacency.clear();
162
163 if (field.getSpace() != NOFIELD) {
164
165 auto basic_entity_data_ptr = fe.getBasicDataPtr();
166 auto th_parent_handle = basic_entity_data_ptr->th_RefParentHandle;
167 auto th_bit_level = basic_entity_data_ptr->th_RefBitLevel;
168
169 std::vector<EntityHandle> parents;
170 parents.reserve(BITREFLEVEL_SIZE);
171
172 CHKERR getParent(fe.getEnt(), parents, moab, th_parent_handle,
173 th_bit_level);
174
175 CHKERR getParentsAdjacencies(field, moab, parents, adjacency);
176 }
177
178 adjTmp.clear();
179 CHKERR getDefaultAdjacencies(moab, field, fe, adjTmp);
180 adjacency.insert(adjacency.end(), adjTmp.begin(), adjTmp.end());
181
182 std::sort(adjacency.begin(), adjacency.end());
183 auto it = std::unique(adjacency.begin(), adjacency.end());
184 adjacency.resize(std::distance(adjacency.begin(), it));
185
186 for (auto e : adjacency) {
187 auto &side_table = fe.getSideNumberTable();
188 if (side_table.find(e) == side_table.end())
189 const_cast<SideNumber_multiIndex &>(side_table)
190 .insert(boost::shared_ptr<SideNumber>(new SideNumber(e, -1, 0, 0)));
191 }
192
194 }
195
196protected:
197 MoFEMErrorCode getParent(EntityHandle fe, std::vector<EntityHandle> &parents,
198 moab::Interface &moab, Tag th_parent_handle,
199 Tag th_bit_level) {
201
202 auto check = [](auto &b, auto &m, auto &bit) {
203 return ((bit & b).any()) && ((bit & m) == bit);
204 };
205
206 BitRefLevel bit_fe;
207 CHKERR moab.tag_get_data(th_bit_level, &fe, 1, &bit_fe);
208 if (check(bitEnt, bitEntMask, bit_fe)) {
209
210 using GetParent = boost::function<MoFEMErrorCode(
211 EntityHandle fe, std::vector<EntityHandle> & parents)>;
212 /**
213 * @brief this function os called recursively, until all stack of parents
214 * is found.
215 *
216 */
217 GetParent get_parent = [&](EntityHandle fe,
218 std::vector<EntityHandle> &parents) {
220 EntityHandle fe_parent;
221
222 CHKERR moab.tag_get_data(th_parent_handle, &fe, 1, &fe_parent);
223 auto parent_type = type_from_handle(fe_parent);
224 auto back_type = type_from_handle(fe);
225 BitRefLevel bit_parent;
226 CHKERR moab.tag_get_data(th_bit_level, &fe_parent, 1, &bit_parent);
227 if (check(bitParent, bitParentMask, bit_parent)) {
228 if ((fe_parent != 0) && (fe_parent != fe) &&
229 (parent_type == back_type)) {
230 parents.push_back(fe_parent);
231 CHKERR get_parent(parents.back(), parents);
232 }
233 }
235 };
236
237 CHKERR get_parent(fe, parents);
238 }
240 }
241
243 moab::Interface &moab,
244 std::vector<EntityHandle> &parents,
245 std::vector<EntityHandle> &adjacency) {
247
248 switch (field.getSpace()) {
249 case H1:
250 for (auto fe_ent : parents)
251 CHKERR moab.get_connectivity(&*parents.begin(), parents.size(),
252 adjacency, true);
253 case HCURL:
254 if constexpr (DIM >= 1)
255 for (auto fe_ent : parents)
256 CHKERR moab.get_adjacencies(&fe_ent, 1, 1, false, adjacency,
257 moab::Interface::UNION);
258 case HDIV:
259 if constexpr (DIM == 2)
260 for (auto fe_ent : parents)
261 CHKERR moab.get_adjacencies(&fe_ent, 1, 2, false, adjacency,
262 moab::Interface::UNION);
263 case L2:
264 for (auto fe_ent : parents)
265 adjacency.push_back(fe_ent);
266 break;
267 default:
268 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
269 "this field is not implemented for face finite element");
270 }
271
272 if (adjacency.size()) {
273 std::sort(adjacency.begin(), adjacency.end());
274 auto it = std::unique(adjacency.begin(), adjacency.end());
275 adjacency.resize(std::distance(adjacency.begin(), it));
276 }
277
279 }
280
282 const Field &field,
283 const EntFiniteElement &fe,
284 std::vector<EntityHandle> &adjacency) {
286 if constexpr (DIM == 3)
287 CHKERR DefaultElementAdjacency::defaultVolume(moab, field, fe, adjacency);
288 if constexpr (DIM == 2)
289 CHKERR DefaultElementAdjacency::defaultFace(moab, field, fe, adjacency);
290 else if constexpr (DIM == 1)
291 CHKERR DefaultElementAdjacency::defaultEdge(moab, field, fe, adjacency);
292 else if constexpr (DIM == 0)
293 CHKERR DefaultElementAdjacency::defaultVertex(moab, field, fe, adjacency);
295 };
296
301 std::vector<EntityHandle> adjTmp;
302};
303
304/**
305 * @brief Create adjacency to parent skeleton elements.
306 *
307 * That class is used during entity finite element construction.
308 *
309 * @tparam DIM dimension of parent element
310 */
311template <int DIM>
314
316 BitRefLevel bit_parent_mask,
317 BitRefLevel bit_ent,
318 BitRefLevel bit_ent_mask)
319 : ParentFiniteElementAdjacencyFunction<DIM>(bit_parent, bit_parent_mask,
320 bit_ent, bit_ent_mask) {}
321
322 /**
323 * @brief Function setting adjacencies to DOFs of parent element
324 *
325 * @note elements form child, see dofs from parent, so DOFs located on
326 * adjacencies of parent entity has adjacent to dofs of child.
327 *
328 * @tparam DIM dimension of the element entity
329 * @param moab
330 * @param field
331 * @param fe
332 * @param adjacency
333 * @return MoFEMErrorCode
334 */
335 MoFEMErrorCode operator()(moab::Interface &moab, const Field &field,
336 const EntFiniteElement &fe,
337 std::vector<EntityHandle> &adjacency) {
339
340 adjacency.clear();
341 CHKERR this->getDefaultAdjacencies(moab, field, fe, adjacency);
342
343 if (field.getSpace() != NOFIELD) {
344
345 const auto fe_ent = fe.getEnt();
346 brideAdjacencyEdge.clear();
347 CHKERR moab.get_adjacencies(&fe_ent, 1, DIM + 1, false,
349
350 std::vector<EntityHandle> parents;
351
352 if (this->bitParent.any()) {
353 auto basic_entity_data_ptr = fe.getBasicDataPtr();
354 auto th_parent_handle = basic_entity_data_ptr->th_RefParentHandle;
355 auto th_bit_level = basic_entity_data_ptr->th_RefBitLevel;
356 parents.reserve(BITREFLEVEL_SIZE);
357 for (auto bridge_fe : brideAdjacencyEdge) {
358 CHKERR this->getParent(bridge_fe, parents, moab, th_parent_handle,
359 th_bit_level);
360 };
361 parents.insert(parents.end(), brideAdjacencyEdge.begin(),
362 brideAdjacencyEdge.end());
363 } else {
364 parents.swap(brideAdjacencyEdge);
365 }
366
367 CHKERR this->getParentsAdjacencies(field, moab, parents, adjacency);
368
369 std::sort(adjacency.begin(), adjacency.end());
370 auto it = std::unique(adjacency.begin(), adjacency.end());
371 adjacency.resize(std::distance(adjacency.begin(), it));
372
373 for (auto e : adjacency) {
374 auto &side_table = fe.getSideNumberTable();
375 if (side_table.find(e) == side_table.end())
376 const_cast<SideNumber_multiIndex &>(side_table)
377 .insert(
378 boost::shared_ptr<SideNumber>(new SideNumber(e, -1, 0, 0)));
379 }
380 }
381
383 }
384
385protected:
386 std::vector<EntityHandle> brideAdjacencyEdge;
387};
388
389} // namespace MoFEM
390
391#endif //__MESH_PROJECTION_DATA_OPERATORS_HPP__
@ QUIET
Definition: definitions.h:208
#define BITREFLEVEL_SIZE
max number of refinements
Definition: definitions.h:219
FieldSpace
approximation spaces
Definition: definitions.h:82
@ L2
field with C-1 continuity
Definition: definitions.h:88
@ NOFIELD
scalar or vector of scalars describe (no true field)
Definition: definitions.h:84
@ H1
continuous field
Definition: definitions.h:85
@ HCURL
field with continuous tangents
Definition: definitions.h:86
@ HDIV
field with continuous normal traction
Definition: definitions.h:87
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:346
@ MOFEM_NOT_IMPLEMENTED
Definition: definitions.h:32
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:416
#define CHKERR
Inline error check.
Definition: definitions.h:535
FTensor::Index< 'm', SPACE_DIM > m
multi_index_container< boost::shared_ptr< SideNumber >, indexed_by< ordered_unique< member< SideNumber, EntityHandle, &SideNumber::ent > >, ordered_non_unique< composite_key< SideNumber, const_mem_fun< SideNumber, EntityType, &SideNumber::getEntType >, member< SideNumber, signed char, &SideNumber::side_number > > > > > SideNumber_multiIndex
SideNumber_multiIndex for SideNumber.
SeverityLevel
Severity levels.
Definition: LogManager.hpp:33
auto bit
set bit
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:56
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
auto type_from_handle(const EntityHandle h)
get type from entity handle
Definition: Templates.hpp:1594
constexpr auto field_name
static MoFEMErrorCode defaultVertex(Interface &moab, const Field &field, const EntFiniteElement &fe, std::vector< EntityHandle > &adjacency)
static MoFEMErrorCode defaultVolume(Interface &moab, const Field &field, const EntFiniteElement &fe, std::vector< EntityHandle > &adjacency)
static MoFEMErrorCode defaultEdge(Interface &moab, const Field &field, const EntFiniteElement &fe, std::vector< EntityHandle > &adjacency)
static MoFEMErrorCode defaultFace(Interface &moab, const Field &field, const EntFiniteElement &fe, std::vector< EntityHandle > &adjacency)
Finite element data for entity.
Data on single entity (This is passed as argument to DataOperator::doWork)
data structure for finite element entity
Provide data structure for (tensor) field approximation.
FieldSpace getSpace() const
Get field approximation space.
OpType
Controls loop over entities on element.
Operator to project base functions from parent entity to child.
MoFEMErrorCode opRhs(EntitiesFieldData &data, const bool error_if_no_base=false)
boost::ptr_deque< EntitiesFieldData::EntData > poolEntsVector
boost::shared_ptr< ForcesAndSourcesCore > parentElePtr
Operator to execute finite element instance on parent element. This operator is typically used to pro...
boost::shared_ptr< ForcesAndSourcesCore > parentElePtr
boost::shared_ptr< ForcesAndSourcesCore > thisElePtr
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
Operator for linear form, usually to calculate values on right hand side.
MoFEMErrorCode operator()(moab::Interface &moab, const Field &field, const EntFiniteElement &fe, std::vector< EntityHandle > &adjacency)
Function setting adjacencies to DOFs of parent element.
MoFEMErrorCode getParent(EntityHandle fe, std::vector< EntityHandle > &parents, moab::Interface &moab, Tag th_parent_handle, Tag th_bit_level)
MoFEMErrorCode getParentsAdjacencies(const Field &field, moab::Interface &moab, std::vector< EntityHandle > &parents, std::vector< EntityHandle > &adjacency)
ParentFiniteElementAdjacencyFunction(BitRefLevel bit_parent, BitRefLevel bit_parent_mask, BitRefLevel bit_ent, BitRefLevel bit_ent_mask)
MoFEMErrorCode getDefaultAdjacencies(moab::Interface &moab, const Field &field, const EntFiniteElement &fe, std::vector< EntityHandle > &adjacency)
ParentFiniteElementAdjacencyFunctionSkeleton(BitRefLevel bit_parent, BitRefLevel bit_parent_mask, BitRefLevel bit_ent, BitRefLevel bit_ent_mask)
MoFEMErrorCode operator()(moab::Interface &moab, const Field &field, const EntFiniteElement &fe, std::vector< EntityHandle > &adjacency)
Function setting adjacencies to DOFs of parent element.
keeps information about side number for the finite element
SideNumber_multiIndex & getSideNumberTable() const
const boost::shared_ptr< BasicEntityData > getBasicDataPtr() const