155 <<
"Refinement for hexes is not implemented";
161 auto create_reference_element = [&moab_ref]() {
163 constexpr
double base_coords[] = {0, 0, 0, 1, 0, 0, 1, 1, 0, 0, 1, 0,
164 0, 0, 1, 1, 0, 1, 1, 1, 1, 0, 1, 1};
166 for (
int nn = 0; nn < 8; nn++)
167 CHKERR moab_ref.create_vertex(&base_coords[3 * nn], nodes[nn]);
169 CHKERR moab_ref.create_element(MBHEX, nodes, 8, hex);
173 auto add_ho_nodes = [&]() {
176 CHKERR moab_ref.get_entities_by_type(0, MBHEX, hexes,
true);
178 CHKERR moab_ref.create_meshset(MESHSET_SET, meshset);
179 CHKERR moab_ref.add_entities(meshset, hexes);
180 CHKERR moab_ref.convert_entities(meshset,
true,
true,
true);
181 CHKERR moab_ref.delete_entities(&meshset, 1);
185 auto set_gauss_pts = [&](std::map<EntityHandle, int> &little_map) {
188 CHKERR moab_ref.get_entities_by_type(0, MBHEX, hexes,
true);
190 CHKERR moab_ref.get_connectivity(hexes, hexes_nodes,
false);
192 gauss_pts.resize(hexes_nodes.size(), 4,
false);
194 for (
auto node : hexes_nodes) {
195 CHKERR moab_ref.get_coords(&node, 1, &gauss_pts(gg, 0));
196 little_map[node] = gg;
199 gauss_pts = trans(gauss_pts);
203 auto set_ref_hexes = [&](std::map<EntityHandle, int> &little_map) {
206 CHKERR moab_ref.get_entities_by_type(0, MBHEX, hexes,
true);
209 for (
auto hex : hexes) {
212 CHKERR moab_ref.get_connectivity(hex, conn, num_nodes,
false);
213 if (ref_hexes.size2() != num_nodes) {
214 ref_hexes.resize(hexes.size(), num_nodes);
216 for (
int nn = 0; nn != num_nodes; ++nn) {
217 ref_hexes(hh, nn) = little_map[conn[nn]];
224 auto set_shape_functions = [&]() {
228 const auto nb_gauss_pts = gauss_pts.size2();
229 shape_functions.resize(nb_gauss_pts, 8);
230 for (
int gg = 0; gg != nb_gauss_pts; ++gg) {
231 const double ksi = gauss_pts(0, gg);
232 const double zeta = gauss_pts(1, gg);
233 const double eta = gauss_pts(2, gg);
250 CHKERR create_reference_element();
253 std::map<EntityHandle, int> little_map;
254 CHKERR set_gauss_pts(little_map);
255 CHKERR set_ref_hexes(little_map);
256 CHKERR set_shape_functions();