34 {
36
38
39 moab::Core core_ref;
40 moab::Interface &moab_ref = core_ref;
41
42 auto create_reference_element = [&moab_ref]() {
44 constexpr double base_coords[] = {0, 0, 0, 1, 0, 0, 0, 1, 0, 0, 0, 1};
46 for (int nn = 0; nn < 4; nn++) {
48 moab_ref.create_vertex(&base_coords[3 * nn], nodes[nn]);
49 }
51 CHKERR moab_ref.create_element(MBTET, nodes, 4, tet);
53 };
54
57
58 auto refine_ref_tetrahedron = [this, &m_field_ref, max_level]() {
60
61
63 m_field_ref.
getInterface<BitRefManager>()->setBitRefLevelByDim(
65 for (int ll = 0; ll != max_level; ++ll) {
67 ll);
70 ->getEntitiesByTypeAndRefLevel(
BitRefLevel().set(ll),
74 ->getEntitiesByTypeAndRefLevel(
BitRefLevel().set(ll),
76
77 MeshRefinement *m_ref;
79 CHKERR m_ref->addVerticesInTheMiddleOfEdges(edges,
82 }
84 };
85
86 auto get_ref_gauss_pts_and_shape_functions = [this, max_level, &moab_ref,
87 &m_field_ref]() {
89 for (int ll = 0; ll != max_level + 1; ++ll) {
92 m_field_ref.
getInterface<BitRefManager>()->getEntitiesByTypeAndRefLevel(
96 CHKERR moab_ref.create_meshset(MESHSET_SET, meshset);
97 CHKERR moab_ref.add_entities(meshset, tets);
98 CHKERR moab_ref.convert_entities(meshset,
true,
false,
false);
99 CHKERR moab_ref.delete_entities(&meshset, 1);
100 }
102 CHKERR moab_ref.get_connectivity(tets, elem_nodes,
false);
103
105 gauss_pts.resize(elem_nodes.size(), 4, false);
106 std::map<EntityHandle, int> little_map;
107 Range::iterator nit = elem_nodes.begin();
108 for (int gg = 0; nit != elem_nodes.end(); nit++, gg++) {
109 CHKERR moab_ref.get_coords(&*nit, 1, &gauss_pts(gg, 0));
110 little_map[*nit] = gg;
111 }
112 gauss_pts = trans(gauss_pts);
113
115 Range::iterator tit = tets.begin();
116 for (int tt = 0; tit != tets.end(); ++tit, ++tt) {
118 int num_nodes;
119 CHKERR moab_ref.get_connectivity(*tit, conn, num_nodes,
false);
120 if (tt == 0) {
121 ref_tets.resize(tets.size(), num_nodes);
122 }
123 for (int nn = 0; nn != num_nodes; ++nn) {
124 ref_tets(tt, nn) = little_map[conn[nn]];
125 }
126 }
127
129 shape_functions.resize(elem_nodes.size(), 4);
131 &gauss_pts(1, 0), &gauss_pts(2, 0), elem_nodes.size());
132 }
134 };
135
139
140 CHKERR create_reference_element();
141 CHKERR refine_ref_tetrahedron();
142 CHKERR get_ref_gauss_pts_and_shape_functions();
143
145}
#define MOFEM_TAG_AND_LOG_C(channel, severity, tag, format,...)
Tag and log in channel.
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
#define CHKERR
Inline error check.
std::bitset< BITREFLEVEL_SIZE > BitRefLevel
Bit structure attached to each entity identifying to what mesh entity is attached.
Deprecated interface functions.
std::vector< ublas::matrix< int > > levelRef
std::vector< MatrixDouble > levelShapeFunctions
std::vector< MatrixDouble > levelGaussPtsOnRefMesh
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface refernce to pointer of interface.