Generate reference mesh on single element.
Each element is subdivided on smaller elements, i.e. a reference mesh on single element is created. Nodes of such reference mesh are used as integration points at which field values are calculated and to each node a "moab" tag is attached to store those values.
316 auto get_nb_of_ref_levels_from_options = [
this] {
319 PetscBool flg = PETSC_TRUE;
321 "-my_max_post_proc_ref_level", &max_level, &flg);
329 auto generate_for_hex = [&]() {
335 auto create_reference_element = [&moab_ref]() {
337 constexpr
double base_coords[] = {0, 0, 0, 1, 0, 0, 1, 1, 0, 0, 1, 0,
338 0, 0, 1, 1, 0, 1, 1, 1, 1, 0, 1, 1};
340 for (
int nn = 0; nn < 8; nn++)
341 CHKERR moab_ref.create_vertex(&base_coords[3 * nn], nodes[nn]);
343 CHKERR moab_ref.create_element(MBHEX, nodes, 8, hex);
347 auto add_ho_nodes = [&]() {
350 CHKERR moab_ref.get_entities_by_type(0, MBHEX, hexes,
true);
352 CHKERR moab_ref.create_meshset(MESHSET_SET, meshset);
353 CHKERR moab_ref.add_entities(meshset, hexes);
354 CHKERR moab_ref.convert_entities(meshset,
true,
true,
true);
355 CHKERR moab_ref.delete_entities(&meshset, 1);
359 auto set_gauss_pts = [&](std::map<EntityHandle, int> &little_map) {
362 CHKERR moab_ref.get_entities_by_type(0, MBHEX, hexes,
true);
364 CHKERR moab_ref.get_connectivity(hexes, hexes_nodes,
false);
366 gauss_pts.resize(hexes_nodes.size(), 4,
false);
368 for (
auto node : hexes_nodes) {
369 CHKERR moab_ref.get_coords(&node, 1, &gauss_pts(gg, 0));
370 little_map[node] = gg;
373 gauss_pts = trans(gauss_pts);
377 auto set_ref_hexes = [&](std::map<EntityHandle, int> &little_map) {
380 CHKERR moab_ref.get_entities_by_type(0, MBHEX, hexes,
true);
383 for (
auto hex : hexes) {
386 CHKERR moab_ref.get_connectivity(hex, conn, num_nodes,
false);
387 if (ref_hexes.size2() != num_nodes) {
388 ref_hexes.resize(hexes.size(), num_nodes);
390 for (
int nn = 0; nn != num_nodes; ++nn) {
391 ref_hexes(hh, nn) = little_map[conn[nn]];
398 auto set_shape_functions = [&]() {
402 const auto nb_gauss_pts = gauss_pts.size2();
403 shape_functions.resize(nb_gauss_pts, 8);
404 for (
int gg = 0; gg != nb_gauss_pts; ++gg) {
405 const double ksi = gauss_pts(0, gg);
406 const double zeta = gauss_pts(1, gg);
407 const double eta = gauss_pts(2, gg);
424 CHKERR create_reference_element();
427 std::map<EntityHandle, int> little_map;
428 CHKERR set_gauss_pts(little_map);
429 CHKERR set_ref_hexes(little_map);
430 CHKERR set_shape_functions();
435 auto generate_for_tet = [&]() {
438 const int max_level = get_nb_of_ref_levels_from_options();
443 auto create_reference_element = [&moab_ref]() {
445 constexpr
double base_coords[] = {0, 0, 0, 1, 0, 0, 0, 1, 0, 0, 0, 1};
447 for (
int nn = 0; nn < 4; nn++) {
449 moab_ref.create_vertex(&base_coords[3 * nn], nodes[nn]);
452 CHKERR moab_ref.create_element(MBTET, nodes, 4, tet);
459 auto refine_ref_tetrahedron = [
this, &m_field_ref, max_level]() {
464 m_field_ref.
getInterface<BitRefManager>()->setBitRefLevelByDim(
466 for (
int ll = 0; ll != max_level; ++ll) {
468 "Refine Level %d", ll);
471 ->getEntitiesByTypeAndRefLevel(
475 ->getEntitiesByTypeAndRefLevel(
478 MeshRefinement *m_ref;
480 CHKERR m_ref->addVerticesInTheMiddleOfEdges(
487 auto get_ref_gauss_pts_and_shape_functions = [
this, max_level, &moab_ref,
490 for (
int ll = 0; ll != max_level + 1; ++ll) {
494 ->getEntitiesByTypeAndRefLevel(
498 CHKERR moab_ref.create_meshset(MESHSET_SET, meshset);
499 CHKERR moab_ref.add_entities(meshset, tets);
500 CHKERR moab_ref.convert_entities(meshset,
true,
false,
false);
501 CHKERR moab_ref.delete_entities(&meshset, 1);
504 CHKERR moab_ref.get_connectivity(tets, elem_nodes,
false);
507 gauss_pts.resize(elem_nodes.size(), 4,
false);
508 std::map<EntityHandle, int> little_map;
509 Range::iterator nit = elem_nodes.begin();
510 for (
int gg = 0; nit != elem_nodes.end(); nit++, gg++) {
511 CHKERR moab_ref.get_coords(&*nit, 1, &gauss_pts(gg, 0));
512 little_map[*nit] = gg;
514 gauss_pts = trans(gauss_pts);
517 Range::iterator tit = tets.begin();
518 for (
int tt = 0; tit != tets.end(); ++tit, ++tt) {
521 CHKERR moab_ref.get_connectivity(*tit, conn, num_nodes,
false);
523 ref_tets.resize(tets.size(), num_nodes);
525 for (
int nn = 0; nn != num_nodes; ++nn) {
526 ref_tets(tt, nn) = little_map[conn[nn]];
531 shape_functions.resize(elem_nodes.size(), 4);
533 &gauss_pts(1, 0), &gauss_pts(2, 0),
543 CHKERR create_reference_element();
544 CHKERR refine_ref_tetrahedron();
545 CHKERR get_ref_gauss_pts_and_shape_functions();
550 CHKERR generate_for_hex();
551 CHKERR generate_for_tet();
#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.
PetscErrorCode PetscOptionsGetInt(PetscOptions *, const char pre[], const char name[], PetscInt *ivalue, PetscBool *set)
DeprecatedCoreInterface Interface
Deprecated interface functions.
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface refernce to pointer of interface.
std::vector< MatrixDouble > levelGaussPtsOnRefMeshTets
std::vector< ublas::matrix< int > > levelRefTets
std::vector< MatrixDouble > levelShapeFunctionsHexes
std::vector< MatrixDouble > levelGaussPtsOnRefMeshHexes
std::vector< MatrixDouble > levelShapeFunctionsTets
std::vector< ublas::matrix< int > > levelRefHexes