libMesh
Public Member Functions | Protected Member Functions | Protected Attributes | List of all members
MeshPerElemTest< elem_type > Class Template Reference

#include <mesh_elem_test.h>

Inheritance diagram for MeshPerElemTest< elem_type >:
[legend]

Public Member Functions

void setUp ()
 

Protected Member Functions

bool meshes_equal_enough (Mesh &other_mesh, bool double_precision)
 

Protected Attributes

std::unique_ptr< Mesh_mesh
 
std::string libmesh_suite_name
 

Detailed Description

template<ElemType elem_type>
class MeshPerElemTest< elem_type >

Definition at line 11 of file mesh_elem_test.h.

Member Function Documentation

◆ meshes_equal_enough()

template<ElemType elem_type>
bool MeshPerElemTest< elem_type >::meshes_equal_enough ( Mesh other_mesh,
bool  double_precision 
)
inlineprotected

Definition at line 15 of file mesh_elem_test.h.

References libMesh::MeshBase::elem_ref(), libMesh::DofObject::id(), libMesh::libmesh_ignore(), libMesh::make_range(), libMesh::DistributedMesh::max_elem_id(), libMesh::DistributedMesh::max_node_id(), libMesh::DistributedMesh::node_ptr(), libMesh::Elem::node_ref_range(), libMesh::TypeVector< T >::norm(), libMesh::DofObject::processor_id(), libMesh::DistributedMesh::renumber_elem(), libMesh::DistributedMesh::renumber_node(), libMesh::DistributedMesh::set_next_unique_id(), libMesh::DofObject::set_unique_id(), and libMesh::MeshBase::sub_point_locator().

16  {
17  // We'll need to fix up processor_id() and unique_id() values
18  // before we can operator== these meshes. But worse: our gold
19  // meshes might have been numbered differently to our generated
20  // meshes. Some of our generated mesh options practically
21  // *require* renumbering (e.g. after interior HEX20 nodes are
22  // deleted, ExodusII still wants to see a contiguous numbering),
23  // but ReplicatedMesh and DistributedMesh renumber differently.
24  //
25  // So, let's renumber too.
26 
27  MeshSerializer serialthis(*this->_mesh);
28  MeshSerializer serialother(other_mesh);
29 
30  const dof_id_type max_elem_id = this->_mesh->max_elem_id();
31  const dof_id_type max_node_id = this->_mesh->max_node_id();
32 
33  CPPUNIT_ASSERT_EQUAL(max_elem_id, other_mesh.max_elem_id());
34  CPPUNIT_ASSERT_EQUAL(max_node_id, other_mesh.max_node_id());
35 
36  auto locator = other_mesh.sub_point_locator();
37 
38  for (Elem * e1 : this->_mesh->element_ptr_range())
39  {
40  const Elem * e2c = (*locator)(e1->vertex_average());
41  CPPUNIT_ASSERT(e2c);
42  Elem & e2 = other_mesh.elem_ref(e2c->id());
43  e1->processor_id() = 0;
44  e2.processor_id() = 0;
45 
46  const dof_id_type e1_id = e1->id();
47  const dof_id_type e2_id = e2.id();
48  // Do a swap if necessary, using a free temporary id
49  if (e1_id != e2_id)
50  {
51  other_mesh.renumber_elem(e1_id, max_elem_id);
52  other_mesh.renumber_elem(e2_id, e1_id);
53  other_mesh.renumber_elem(max_elem_id, e2_id);
54  }
55 
56 #ifdef LIBMESH_ENABLE_UNIQUE_ID
57  e2.set_unique_id(e1->unique_id());
58 #endif
59  }
60 
61  for (Node * n1 : this->_mesh->node_ptr_range())
62  {
63  const Elem * e1c = (*locator)(*n1);
64  Node * n2 = nullptr;
65  for (const Node & n : e1c->node_ref_range())
66  {
67 #if defined(LIBMESH_DEFAULT_QUADRUPLE_PRECISION) || defined(LIBMESH_DEFAULT_TRIPLE_PRECISION)
68  if (double_precision)
69  {
70  const Point diff = Point(*n1)-Point(n);
71 
72  // We may be testing against ExodusII input, and if
73  // we're in triple or quadruple precision that means our
74  // lovely higher-precision node coordinates got
75  // truncated to double to be written. We need to adjust
76  // ours or they won't satisfy operator== later.
77 
78  // We're *also* testing against gold files that were
79  // calculated at double precision, so just casting a
80  // higher precision calculation to double won't give the
81  // exact same result, we have to account for error.
82  if (diff.norm() < 1e-15)
83  for (auto d : make_range(LIBMESH_DIM))
84  (*n1)(d) = double(n(d));
85  }
86 #else
87  libmesh_ignore(double_precision);
88 #endif
89  if (Point(*n1) == Point(n))
90  n2 = other_mesh.node_ptr(n.id());
91  }
92  CPPUNIT_ASSERT(n2);
93  n1->processor_id() = 0;
94  n2->processor_id() = 0;
95 
96  const dof_id_type n1_id = n1->id();
97  const dof_id_type n2_id = n2->id();
98  // Do a swap if necessary, using a free temporary id
99  if (n1_id != n2_id)
100  {
101  other_mesh.renumber_node(n1_id,max_node_id);
102  other_mesh.renumber_node(n2_id,n1_id);
103  other_mesh.renumber_node(max_node_id, n2_id);
104  }
105 
106 #ifdef LIBMESH_ENABLE_UNIQUE_ID
107  n2->set_unique_id(n1->unique_id());
108 #endif
109  }
110 
111 #ifdef LIBMESH_ENABLE_UNIQUE_ID
112  other_mesh.set_next_unique_id(this->_mesh->parallel_max_unique_id());
113  this->_mesh->set_next_unique_id(this->_mesh->parallel_max_unique_id());
114 #endif
115 
116  return *this->_mesh == other_mesh;
117  }
A Node is like a Point, but with more information.
Definition: node.h:52
auto norm() const
Definition: type_vector.h:908
std::unique_ptr< PointLocatorBase > sub_point_locator() const
Definition: mesh_base.C:1826
virtual void renumber_elem(dof_id_type old_id, dof_id_type new_id) override final
Changes the id of element old_id, both by changing elem(old_id)->id() and by moving elem(old_id) in t...
std::unique_ptr< Mesh > _mesh
Definition: elem_test.h:26
This is the base class from which all geometric element types are derived.
Definition: elem.h:94
void libmesh_ignore(const Args &...)
dof_id_type id() const
Definition: dof_object.h:819
virtual void set_next_unique_id(unique_id_type id) override
Sets the next available unique id to be used.
SimpleRange< NodeRefIter > node_ref_range()
Returns a range with all nodes of an element, usable in range-based for loops.
Definition: elem.h:2679
virtual void renumber_node(dof_id_type old_id, dof_id_type new_id) override final
Changes the id of node old_id, both by changing node(old_id)->id() and by moving node(old_id) in the ...
void set_unique_id(unique_id_type new_id)
Sets the unique_id for this DofObject.
Definition: dof_object.h:848
Temporarily serialize a DistributedMesh for non-distributed-mesh capable code paths.
virtual dof_id_type max_node_id() const override final
virtual const Elem & elem_ref(const dof_id_type i) const
Definition: mesh_base.h:778
IntRange< T > make_range(T beg, T end)
The 2-parameter make_range() helper function returns an IntRange<T> when both input parameters are of...
Definition: int_range.h:176
virtual dof_id_type max_elem_id() const override final
processor_id_type processor_id() const
Definition: dof_object.h:881
A Point defines a location in LIBMESH_DIM dimensional Real space.
Definition: point.h:39
uint8_t dof_id_type
Definition: id_types.h:67
virtual const Node * node_ptr(const dof_id_type i) const override final

◆ setUp()

template<ElemType elem_type>
void PerElemTest< elem_type >::setUp ( )
inlineinherited

Definition at line 30 of file elem_test.h.

References libMesh::Elem::build(), libMesh::MeshTools::Generation::build_cube(), libMesh::C0POLYGON, libMesh::C0POLYHEDRON, dim, libMesh::index_range(), libMesh::INFHEX16, libMesh::INFHEX18, libMesh::INFHEX8, libMesh::INFPRISM12, libMesh::INFPRISM6, libMesh::INFQUAD4, libMesh::INFQUAD6, libMesh::make_range(), libMesh::Real, and libMesh::Elem::set_node().

31  {
32  const Real minpos = 1.5, maxpos = 5.5;
33  const unsigned int N = 2;
34 
35  _mesh = std::make_unique<Mesh>(*TestCommWorld);
36 
37  std::unique_ptr<Elem> test_elem;
38 
39  if (elem_type != C0POLYHEDRON)
40  test_elem = Elem::build(elem_type);
41 
42 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
43 #if LIBMESH_DIM > 1
44  if (test_elem.get() && test_elem->infinite())
45  {
46  Elem * elem = _mesh->add_elem(std::move(test_elem));
47 
48  const auto add_point =
49  [this, elem](const unsigned int i,
50  const Real x,
51  const Real y,
52  const Real
53 #if LIBMESH_DIM == 3
54  z
55 #endif
56  )
57  {
58 #if LIBMESH_DIM == 2
59  auto node = _mesh->add_point(Point(x, y), i);
60 #else
61  auto node = _mesh->add_point(Point(x, y, z), i);
62 #endif
63  elem->set_node(i, node);
64  };
65 
66  const Real halfpos = (minpos + maxpos) / 2.;
67 
68  if (elem_type == INFQUAD4 || elem_type == INFQUAD6 ||
69  elem_type == INFHEX8 || elem_type == INFHEX16 || elem_type == INFHEX18)
70  {
71  const bool is_quad = (elem_type == INFQUAD4 || elem_type == INFQUAD6);
72 
73  add_point(0, minpos, minpos, minpos);
74  add_point(1, maxpos, minpos, minpos);
75  add_point(2+is_quad, maxpos, maxpos, minpos);
76  add_point(3-is_quad, minpos, maxpos, minpos);
77 
78  if (elem_type == INFQUAD6)
79  {
80  add_point(4, halfpos, minpos, minpos);
81  add_point(5, halfpos, maxpos, minpos);
82  }
83  }
84  if (elem_type == INFHEX8 || elem_type == INFHEX16 || elem_type == INFHEX18)
85  {
86  add_point(4, minpos, minpos, maxpos);
87  add_point(5, maxpos, minpos, maxpos);
88  add_point(6, maxpos, maxpos, maxpos);
89  add_point(7, minpos, maxpos, maxpos);
90 
91  if (elem_type == INFHEX16 || elem_type == INFHEX18)
92  {
93  add_point(8, halfpos, minpos, minpos);
94  add_point(9, maxpos, halfpos, minpos);
95  add_point(10, halfpos, maxpos, minpos);
96  add_point(11, minpos, halfpos, minpos);
97  add_point(12, halfpos, minpos, maxpos);
98  add_point(13, maxpos, halfpos, maxpos);
99  add_point(14, halfpos, maxpos, maxpos);
100  add_point(15, minpos, halfpos, maxpos);
101  }
102  if (elem_type == INFHEX18)
103  {
104  add_point(16, halfpos, halfpos, minpos);
105  add_point(17, halfpos, halfpos, maxpos);
106  }
107  }
108  if (elem_type == INFPRISM6 || elem_type == INFPRISM12)
109  {
110  add_point(0, minpos, minpos, minpos);
111  add_point(1, maxpos, minpos, minpos);
112  add_point(2, halfpos, maxpos, minpos);
113  add_point(3, minpos, minpos, maxpos);
114  add_point(4, maxpos, minpos, maxpos);
115  add_point(5, halfpos, maxpos, maxpos);
116 
117  if (elem_type == INFPRISM12)
118  {
119  add_point(6, halfpos, minpos, minpos);
120  add_point(7, (halfpos + maxpos) / 2., halfpos, minpos);
121  add_point(8, (halfpos + minpos) / 2., halfpos, minpos);
122  add_point(9, halfpos, minpos, maxpos);
123  add_point(10, (halfpos + maxpos) / 2., halfpos, maxpos);
124  add_point(11, (halfpos + minpos) / 2., halfpos, maxpos);
125  }
126  }
127 
128  _mesh->prepare_for_use();
129  }
130  else
131 #endif // LIBMESH_DIM > 1
132 #endif // LIBMESH_ENABLE_INFINITE_ELEMENTS
133  if (elem_type == C0POLYGON)
134  {
135  // We're not going to implement build_square for e.g.
136  // pentagons any time soon.
137  //
138  // We should probably implement build_dual_mesh, though, and
139  // run it on a perturbed or triangle input so we don't just
140  // get quads in the output...
141 
142  _mesh->add_point(Point(0, 0), 0);
143  _mesh->add_point(Point(1, 0), 1);
144  _mesh->add_point(Point(1.5, 0.5), 2);
145  _mesh->add_point(Point(1, 1), 3);
146  _mesh->add_point(Point(0, 1), 4);
147 
148  std::unique_ptr<Elem> polygon = std::make_unique<C0Polygon>(5);
149  for (auto i : make_range(5))
150  polygon->set_node(i, _mesh->node_ptr(i));
151  polygon->set_id() = 0;
152 
153  _mesh->add_elem(std::move(polygon));
154  _mesh->prepare_for_use();
155  }
156  else if (elem_type == C0POLYHEDRON)
157  {
158  // There's even less point in having a build_cube for general
159  // polyhedra, so again we'll hand-make one for testing and
160  // we'll plan on making a build_dual_mesh for the future.
161 
162 #if 0
163  // Try a truncated cube for relatively simple verification.
164 
165  // We have some differing orientations on the top sides, to
166  // test handling of that.
167  // Lower octagon points, counterclockwise
168  _mesh->add_point(Point(1/Real(3), 0, 0), 0);
169  _mesh->add_point(Point(2/Real(3), 0, 0), 1);
170  _mesh->add_point(Point(1, 1/Real(3), 0), 2);
171  _mesh->add_point(Point(1, 2/Real(3), 0), 3);
172  _mesh->add_point(Point(2/Real(3), 1, 0), 4);
173  _mesh->add_point(Point(1/Real(3), 1, 0), 5);
174  _mesh->add_point(Point(0, 2/Real(3), 0), 6);
175  _mesh->add_point(Point(0, 1/Real(3), 0), 7);
176 
177  // Lower-middle points, counterclockwise
178  _mesh->add_point(Point(0, 0, 1/Real(3)), 8);
179  _mesh->add_point(Point(1, 0, 1/Real(3)), 9);
180  _mesh->add_point(Point(1, 1, 1/Real(3)), 10);
181  _mesh->add_point(Point(0, 1, 1/Real(3)), 11);
182 
183  // Upper-middle points, counterclockwise
184  _mesh->add_point(Point(0, 0, 2/Real(3)), 12);
185  _mesh->add_point(Point(1, 0, 2/Real(3)), 13);
186  _mesh->add_point(Point(1, 1, 2/Real(3)), 14);
187  _mesh->add_point(Point(0, 1, 2/Real(3)), 15);
188 
189  // Upper octagon points, counterclockwise
190  _mesh->add_point(Point(1/Real(3), 0, 1), 16);
191  _mesh->add_point(Point(2/Real(3), 0, 1), 17);
192  _mesh->add_point(Point(1, 1/Real(3), 1), 18);
193  _mesh->add_point(Point(1, 2/Real(3), 1), 19);
194  _mesh->add_point(Point(2/Real(3), 1, 1), 20);
195  _mesh->add_point(Point(1/Real(3), 1, 1), 21);
196  _mesh->add_point(Point(0, 2/Real(3), 1), 22);
197  _mesh->add_point(Point(0, 1/Real(3), 1), 23);
198 
199  const std::vector<std::vector<unsigned int>> nodes_on_side =
200  { {0, 1, 2, 3, 4, 5, 6, 7}, // min z
201  {0, 1, 9, 13, 17, 16, 12, 8}, // min y
202  {2, 3, 10, 14, 19, 18, 13, 9}, // max x
203  {4, 5, 11, 15, 21, 20, 14, 10}, // max y
204  {6, 7, 8, 12, 23, 22, 15, 11}, // min x
205  {16, 17, 18, 19, 20, 21, 22, 23}, // max z
206  {7, 0, 8}, // max nothing
207  {1, 2, 9}, // max x
208  {3, 4, 10}, // max xy
209  {5, 6, 11}, // max y
210  {23, 16, 12}, // max z
211  {17, 18, 13}, // max xz
212  {19, 20, 14}, // max xyz
213  {21, 22, 15} }; // max yz
214 #endif
215 
216 #if 1
217  // Or try a plain cube, for simpler debugging
218  _mesh->add_point(Point(0, 0, 0), 0);
219  _mesh->add_point(Point(1, 0, 0), 1);
220  _mesh->add_point(Point(1, 1, 0), 2);
221  _mesh->add_point(Point(0, 1, 0), 3);
222  _mesh->add_point(Point(0, 0, 1), 4);
223  _mesh->add_point(Point(1, 0, 1), 5);
224  _mesh->add_point(Point(1, 1, 1), 6);
225  _mesh->add_point(Point(0, 1, 1), 7);
226 
227  // With some combinations of face triangulations, even a
228  // simple cube has no tetrahedralization that doesn't have
229  // either interior discontinuities or a 0-volume pancake tet!
230  //
231  // The initial "natural" way to orient our sides is commented
232  // out here; permutations that fix diagonals for us are used
233  // instead.
234  const std::vector<std::vector<unsigned int>> nodes_on_side =
235  { {0, 1, 2, 3}, // min z
236  {0, 1, 5, 4}, // min y
237  // {1, 2, 6, 5}, // max x - bad
238  {2, 6, 5, 1}, // max x
239  {2, 3, 7, 6}, // max y
240  // {3, 0, 4, 7}, // min x - bad
241  {0, 4, 7, 3}, // min x
242  // {4, 5, 6, 7} }; // max z - bad
243  {5, 6, 7, 4} }; // max z
244 #endif
245 
246  // Build all the sides.
247  std::vector<std::shared_ptr<Polygon>> sides(nodes_on_side.size());
248 
249  for (auto s : index_range(nodes_on_side))
250  {
251  const auto & nodes_on_s = nodes_on_side[s];
252  sides[s] = std::make_shared<C0Polygon>(nodes_on_s.size());
253  for (auto i : index_range(nodes_on_s))
254  sides[s]->set_node(i, _mesh->node_ptr(nodes_on_s[i]));
255  }
256 
257  std::unique_ptr<libMesh::Node> mid_elem_node;
258  std::unique_ptr<Elem> polyhedron = std::make_unique<C0Polyhedron>(sides, mid_elem_node);
259  _mesh->add_elem(std::move(polyhedron));
260  if (mid_elem_node)
261  _mesh->add_node(std::move(mid_elem_node));
262  _mesh->prepare_for_use();
263  }
264  else
265  {
266  const unsigned int dim = test_elem->dim();
267  const unsigned int use_x = dim > 0;
268  const unsigned int use_y = dim > 1;
269  const unsigned int use_z = dim > 2;
270 
272  N*use_x, N*use_y, N*use_z,
273  minpos, maxpos,
274  minpos, use_y*maxpos,
275  minpos, use_z*maxpos,
276  elem_type);
277  }
278 
279  // Use non-default subdomain ids so we can properly test their
280  // preservation.
281  //
282  // Use long subdomain names so we can test that we're not
283  // truncating as badly as we used to in ExodusII.
284  for (const auto & elem :
285  this->_mesh->element_ptr_range())
286  {
287  const subdomain_id_type sbdid = 10 + (elem->id() % 10);
288  elem->subdomain_id() = sbdid;
289  std::ostringstream sbdname;
290  sbdname <<
291  "a_very_long_subdomain_name_for_the_subdomain_with_number_" << sbdid;
292  this->_mesh->subdomain_name(sbdid) = sbdname.str();
293  }
294 
295  // Make sure our mesh's cache knows about them all for later
296  // test comparisons
297  this->_mesh->cache_elem_data();
298 
299  // We may have updated only a portion of subdomain names if we're
300  // on a distributed mesh, but every processor should know about
301  // every name.
302  this->_mesh->sync_subdomain_name_map();
303  }
virtual Node *& set_node(const unsigned int i)
Definition: elem.h:2564
std::unique_ptr< Mesh > _mesh
Definition: elem_test.h:26
unsigned int dim
This is the base class from which all geometric element types are derived.
Definition: elem.h:94
dof_id_type id() const
Definition: dof_object.h:819
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
subdomain_id_type subdomain_id() const
Definition: elem.h:2588
IntRange< T > make_range(T beg, T end)
The 2-parameter make_range() helper function returns an IntRange<T> when both input parameters are of...
Definition: int_range.h:176
A Point defines a location in LIBMESH_DIM dimensional Real space.
Definition: point.h:39
auto index_range(const T &sizable)
Helper function that returns an IntRange<std::size_t> representing all the indices of the passed-in v...
Definition: int_range.h:153
void build_cube(UnstructuredMesh &mesh, const unsigned int nx=0, const unsigned int ny=0, const unsigned int nz=0, const Real xmin=0., const Real xmax=1., const Real ymin=0., const Real ymax=1., const Real zmin=0., const Real zmax=1., const ElemType type=INVALID_ELEM, const bool gauss_lobatto_grid=false)
Builds a (elements) cube.

Member Data Documentation

◆ _mesh

template<ElemType elem_type>
std::unique_ptr<Mesh> PerElemTest< elem_type >::_mesh
protectedinherited

Definition at line 26 of file elem_test.h.

◆ libmesh_suite_name

template<ElemType elem_type>
std::string PerElemTest< elem_type >::libmesh_suite_name
protectedinherited

Definition at line 27 of file elem_test.h.


The documentation for this class was generated from the following file: