diff --git a/include/geom/cell_prism.h b/include/geom/cell_prism.h index fbe7688280..20bfa80560 100644 --- a/include/geom/cell_prism.h +++ b/include/geom/cell_prism.h @@ -211,7 +211,7 @@ class Prism : public Cell /** * Master element node locations */ - static const Real _master_points[18][3]; + static const Real _master_points[21][3]; /** * This maps the \f$ j^{th} \f$ Prism vertex to the (in this case) 3 diff --git a/include/geom/cell_pyramid.h b/include/geom/cell_pyramid.h index fe347d1953..7f8a9c4f39 100644 --- a/include/geom/cell_pyramid.h +++ b/include/geom/cell_pyramid.h @@ -211,7 +211,7 @@ class Pyramid : public Cell /** * Master element node locations */ - static const Real _master_points[14][3]; + static const Real _master_points[18][3]; #ifdef LIBMESH_ENABLE_AMR diff --git a/include/geom/face_tri.h b/include/geom/face_tri.h index bfc249a0cb..595ee04014 100644 --- a/include/geom/face_tri.h +++ b/include/geom/face_tri.h @@ -205,7 +205,7 @@ class Tri : public Face /** * Master element node locations */ - static const Real _master_points[6][3]; + static const Real _master_points[7][3]; /** * This maps the \f$ j^{th} \f$ node to the (in this case) 2 side diff --git a/src/fe/fe_abstract.C b/src/fe/fe_abstract.C index 20d435ae44..74c4ffe05e 100644 --- a/src/fe/fe_abstract.C +++ b/src/fe/fe_abstract.C @@ -584,10 +584,10 @@ void FEAbstract::get_refspace_nodes(const ElemType itemType, std::vector case PYRAMID18: { // triangle centers - nodes[14] = Point (-2/Real(3),0.,1/Real(3)); - nodes[15] = Point (0.,2/Real(3),1/Real(3)); - nodes[16] = Point (2/Real(3),0.,1/Real(3)); - nodes[17] = Point (0.,-2/Real(3),1/Real(3)); + nodes[14] = Point (0., -2/Real(3),1/Real(3)); + nodes[15] = Point (2/Real(3),0.,1/Real(3)); + nodes[16] = Point (0.,2/Real(3),1/Real(3)); + nodes[17] = Point (-2/Real(3),0.,1/Real(3)); libmesh_fallthrough(); } diff --git a/src/geom/cell_prism.C b/src/geom/cell_prism.C index 618d681051..180720f34c 100644 --- a/src/geom/cell_prism.C +++ b/src/geom/cell_prism.C @@ -34,7 +34,7 @@ const int Prism::num_sides; const int Prism::num_edges; const int Prism::num_children; -const Real Prism::_master_points[18][3] = +const Real Prism::_master_points[21][3] = { {0, 0, -1}, {1, 0, -1}, @@ -53,7 +53,10 @@ const Real Prism::_master_points[18][3] = {0, 0.5, 1}, {0.5, 0, 0}, {0.5, 0.5, 0}, - {0, 0.5, 0} + {0, 0.5, 0}, + {Real(1)/3, Real(1)/3, -1}, + {Real(1)/3, Real(1)/3, 1}, + {Real(1)/3, Real(1)/3, 0} }; const unsigned int Prism::edge_sides_map[9][2] = diff --git a/src/geom/cell_pyramid.C b/src/geom/cell_pyramid.C index 04bdf5298a..1169693a8d 100644 --- a/src/geom/cell_pyramid.C +++ b/src/geom/cell_pyramid.C @@ -33,7 +33,7 @@ const int Pyramid::num_sides; const int Pyramid::num_edges; const int Pyramid::num_children; -const Real Pyramid::_master_points[14][3] = +const Real Pyramid::_master_points[18][3] = { {-1, -1, 0}, {1, -1, 0}, @@ -44,11 +44,15 @@ const Real Pyramid::_master_points[14][3] = {1, 0, 0}, {0, 1, 0}, {-1, 0, 0}, - {0, -0.5, 0.5}, - {0.5, 0, 0.5}, - {0, 0.5, 0.5}, - {-0.5, 0, 0.5}, - {0, 0, 0} + {-0.5, -0.5, 0.5}, + {0.5, -0.5, 0.5}, + {0.5, 0.5, 0.5}, + {-0.5, 0.5, 0.5}, + {0, 0, 0}, + {0, -Real(2)/3, Real(1)/3}, + { Real(2)/3, 0, Real(1)/3}, + {0, Real(2)/3, Real(1)/3}, + {-Real(2)/3, 0, Real(1)/3} }; const unsigned int Pyramid::edge_sides_map[8][2] = diff --git a/src/geom/face_tri.C b/src/geom/face_tri.C index f63a640954..b6ebd6a879 100644 --- a/src/geom/face_tri.C +++ b/src/geom/face_tri.C @@ -35,14 +35,15 @@ const int Tri::num_children; // Note: we can omit initialization of the third entry of each row because // static variables are automatically zero-initialized. -const Real Tri::_master_points[6][3] = +const Real Tri::_master_points[7][3] = { {0, 0}, {1, 0}, {0, 1}, {0.5, 0}, {0.5, 0.5}, - {0, 0.5} + {0, 0.5}, + {Real(1)/3, Real(1)/3} }; const unsigned int Tri::adjacent_sides_map[/*num_vertices*/3][/*n_adjacent_sides*/2] = diff --git a/tests/fe/fe_test.h b/tests/fe/fe_test.h index 468d982d5f..d2a6e49836 100644 --- a/tests/fe/fe_test.h +++ b/tests/fe/fe_test.h @@ -25,6 +25,7 @@ #define FETEST \ CPPUNIT_TEST( testFEInterface ); \ + CPPUNIT_TEST( testRefspaceNodes ); \ CPPUNIT_TEST( testU ); \ CPPUNIT_TEST( testPartitionOfUnity ); \ CPPUNIT_TEST( testGradU ); \ @@ -668,8 +669,33 @@ class FETest : public FETestBase { #endif // LIBMESH_ENABLE_EXCEPTIONS } + void testRefspaceNodes() + { + LOG_UNIT_TEST; + + if (!this->_elem) + return; + + // We don't currently support reference infinite elements + if (this->_elem->infinite()) + return; + + // We'll probably never support reference polytopes + if (this->_elem->runtime_topology()) + return; + + std::vector nodes; + this->_fe->get_refspace_nodes(this->_elem->type(), nodes); + for (auto n : index_range(nodes)) + LIBMESH_ASSERT_REALVEC_EQUAL(nodes[n], + this->_elem->master_point(n), + TOLERANCE*TOLERANCE); + } + void testPartitionOfUnity() { + LOG_UNIT_TEST; + if (!this->_elem) return; diff --git a/tests/geom/elem_test.C b/tests/geom/elem_test.C index 41ced6f6b7..943e4f8554 100644 --- a/tests/geom/elem_test.C +++ b/tests/geom/elem_test.C @@ -50,6 +50,31 @@ public: } } + void test_ref_elem() + { + LOG_UNIT_TEST; + + for (const auto & elem : this->_mesh->active_local_element_ptr_range()) + { + // We don't currently support reference infinite elements + if (elem->infinite()) + continue; + + // We'll probably never support reference polytopes + if (elem->runtime_topology()) + continue; + + const Elem & ref_elem = *elem->reference_elem(); + + for (auto i : elem->node_index_range()) + { + LIBMESH_ASSERT_REALVEC_EQUAL(elem->master_point(i), + ref_elem.point(i), + TOLERANCE*TOLERANCE); + } + } + } + void test_quality() { LOG_UNIT_TEST; @@ -944,6 +969,7 @@ public: #define ELEMTEST \ CPPUNIT_TEST( test_bounding_box ); \ + CPPUNIT_TEST( test_ref_elem ); \ CPPUNIT_TEST( test_quality ); \ CPPUNIT_TEST( test_node_edge_map_consistency ); \ CPPUNIT_TEST( test_maps ); \