libMesh
boundary_info.C
Go to the documentation of this file.
1 // Ignore unused parameter warnings coming from cppunit headers
2 #include <libmesh/ignore_warnings.h>
3 #include <cppunit/extensions/HelperMacros.h>
4 #include <cppunit/TestCase.h>
5 #include <libmesh/restore_warnings.h>
6 
7 #include <libmesh/mesh.h>
8 #include <libmesh/mesh_generation.h>
9 #include <libmesh/boundary_info.h>
10 #include <libmesh/elem.h>
11 #include <libmesh/face_quad4_shell.h>
12 #include <libmesh/equation_systems.h>
13 #include <libmesh/zero_function.h>
14 #include <libmesh/dirichlet_boundaries.h>
15 #include <libmesh/dof_map.h>
16 
17 #include "test_comm.h"
18 
19 // THE CPPUNIT_TEST_SUITE_END macro expands to code that involves
20 // std::auto_ptr, which in turn produces -Wdeprecated-declarations
21 // warnings. These can be ignored in GCC as long as we wrap the
22 // offending code in appropriate pragmas. We can't get away with a
23 // single ignore_warnings.h inclusion at the beginning of this file,
24 // since the libmesh headers pull in a restore_warnings.h at some
25 // point. We also don't bother restoring warnings at the end of this
26 // file since it's not a header.
27 #include <libmesh/ignore_warnings.h>
28 
29 using namespace libMesh;
30 
31 class BoundaryInfoTest : public CppUnit::TestCase {
35 public:
36  CPPUNIT_TEST_SUITE( BoundaryInfoTest );
37 
38  CPPUNIT_TEST( testMesh );
39  CPPUNIT_TEST( testEdgeBoundaryConditions );
40  CPPUNIT_TEST( testShellFaceConstraints );
41 
42  CPPUNIT_TEST_SUITE_END();
43 
44 protected:
45 
46 public:
47  void setUp()
48  {
49  }
50 
51  void tearDown()
52  {
53  }
54 
55  void testMesh()
56  {
58 
60  2, 2,
61  0., 1.,
62  0., 1.,
63  QUAD4);
64 
65  BoundaryInfo & bi = mesh.get_boundary_info();
66 
67  // Side lists should be cleared and refilled by each call
68  std::vector<dof_id_type> element_id_list;
69  std::vector<unsigned short int> side_list;
70  std::vector<boundary_id_type> bc_id_list;
71 
72  // build_square adds boundary_ids 0,1,2,3 for the bottom, right,
73  // top, and left sides, respectively.
74 
75  // On a ReplicatedMesh, we should see all 4 ids on each processor
76  if (mesh.is_serial())
77  CPPUNIT_ASSERT_EQUAL(static_cast<std::size_t>(4), bi.n_boundary_ids());
78 
79  // On any mesh, we should see each id on *some* processor
80  {
81  const std::set<boundary_id_type> & bc_ids = bi.get_boundary_ids();
82  for (unsigned int i = 0 ; i != 4; ++i)
83  {
84  bool has_bcid = bc_ids.count(i);
85  mesh.comm().max(has_bcid);
86  CPPUNIT_ASSERT(has_bcid);
87  }
88  }
89 
90  // Build the side list
91  bi.build_side_list (element_id_list, side_list, bc_id_list);
92 
93  // Check that there are exactly 8 sides in the BoundaryInfo for a
94  // replicated mesh
95  if (mesh.is_serial())
96  CPPUNIT_ASSERT_EQUAL(static_cast<std::size_t>(8), element_id_list.size());
97 
98  // Let's test that we can remove them successfully.
99  bi.remove_id(0);
100 
101  // Check that there are now only 3 boundary ids total on the Mesh.
102  if (mesh.is_serial())
103  CPPUNIT_ASSERT_EQUAL(static_cast<std::size_t>(3), bi.n_boundary_ids());
104 
105  {
106  const std::set<boundary_id_type> & bc_ids = bi.get_boundary_ids();
107  CPPUNIT_ASSERT(!bc_ids.count(0));
108  for (unsigned int i = 1 ; i != 4; ++i)
109  {
110  bool has_bcid = bc_ids.count(i);
111  mesh.comm().max(has_bcid);
112  CPPUNIT_ASSERT(has_bcid);
113  }
114  }
115 
116  // Build the side list again
117  bi.build_side_list (element_id_list, side_list, bc_id_list);
118 
119  // Check that there are now exactly 6 sides left in the
120  // BoundaryInfo on a replicated mesh
121  if (mesh.is_serial())
122  CPPUNIT_ASSERT_EQUAL(static_cast<std::size_t>(6), element_id_list.size());
123 
124  // Check that the removed ID is really removed
125  CPPUNIT_ASSERT(std::find(bc_id_list.begin(), bc_id_list.end(), 0) == bc_id_list.end());
126 
127  // Remove the same id again, make sure nothing changes.
128  bi.remove_id(0);
129  if (mesh.is_serial())
130  CPPUNIT_ASSERT_EQUAL(static_cast<std::size_t>(3), bi.n_boundary_ids());
131 
132  // Remove the remaining IDs, verify that we have no sides left and
133  // that we can safely reuse the same vectors in the
134  // build_side_list() call.
135  bi.remove_id(1);
136  bi.remove_id(2);
137  bi.remove_id(3);
138  bi.build_side_list (element_id_list, side_list, bc_id_list);
139 
140  CPPUNIT_ASSERT_EQUAL(static_cast<std::size_t>(0), bi.n_boundary_ids());
141  CPPUNIT_ASSERT_EQUAL(static_cast<std::size_t>(0), element_id_list.size());
142  }
143 
145  {
146  const unsigned int n_elem = 5;
147  const std::string mesh_filename = "cube_mesh.xda";
148 
149  {
152  n_elem, n_elem, n_elem,
153  0., 1.,
154  0., 1.,
155  0., 1.,
156  HEX8);
157 
158  BoundaryInfo & bi = mesh.get_boundary_info();
159 
160  // build_cube does not add any edge boundary IDs
161  CPPUNIT_ASSERT_EQUAL(static_cast<std::size_t>(0), bi.n_edge_conds());
162 
163  // Let's now add some edge boundary IDs.
164  // We loop over all elements (not just local elements) so that
165  // all processors know about the boundary IDs
166  const boundary_id_type BOUNDARY_ID_MAX_X = 2;
167  const boundary_id_type BOUNDARY_ID_MIN_Y = 1;
168  const boundary_id_type EDGE_BOUNDARY_ID = 20;
169 
171  const MeshBase::const_element_iterator end_el = mesh.elements_end();
172  for ( ; el != end_el; ++el)
173  {
174  const Elem * elem = *el;
175 
176  unsigned int side_max_x = 0, side_min_y = 0;
177  bool found_side_max_x = false, found_side_min_y = false;
178 
179  for (unsigned int side=0; side<elem->n_sides(); side++)
180  {
181  if (mesh.get_boundary_info().has_boundary_id(elem, side, BOUNDARY_ID_MAX_X))
182  {
183  side_max_x = side;
184  found_side_max_x = true;
185  }
186 
187  if (mesh.get_boundary_info().has_boundary_id(elem, side, BOUNDARY_ID_MIN_Y))
188  {
189  side_min_y = side;
190  found_side_min_y = true;
191  }
192  }
193 
194  // If elem has sides on boundaries
195  // BOUNDARY_ID_MAX_X and BOUNDARY_ID_MIN_Y
196  // then let's set an edge boundary condition
197  if (found_side_max_x && found_side_min_y)
198  for (unsigned int e=0; e<elem->n_edges(); e++)
199  if (elem->is_edge_on_side(e, side_max_x) &&
200  elem->is_edge_on_side(e, side_min_y))
201  bi.add_edge(elem, e, EDGE_BOUNDARY_ID);
202  }
203 
204  // Check that we have the expected number of edge boundary IDs after
205  // updating bi
206  CPPUNIT_ASSERT_EQUAL(static_cast<std::size_t>(n_elem), bi.n_edge_conds());
207 
208  mesh.write(mesh_filename);
209  }
210 
212  mesh.read(mesh_filename);
213 
214  // Check that writing and reading preserves the edge boundary IDs
215  BoundaryInfo & bi = mesh.get_boundary_info();
216  CPPUNIT_ASSERT_EQUAL(static_cast<std::size_t>(n_elem), bi.n_edge_conds());
217  }
218 
220  {
221  // Make a simple two element mesh that we can use to test constraints
222  Mesh mesh(*TestCommWorld, 3);
223 
224  // (0,1) (1,1)
225  // x---------------x
226  // | |
227  // | |
228  // | |
229  // | |
230  // | |
231  // x---------------x
232  // (0,0) (1,0)
233  // | |
234  // | |
235  // | |
236  // | |
237  // x---------------x
238  // (0,-1) (1,-1)
239 
240  mesh.add_point( Point(0.0,-1.0), 4 );
241  mesh.add_point( Point(1.0,-1.0), 5 );
242  mesh.add_point( Point(1.0, 0.0), 1 );
243  mesh.add_point( Point(1.0, 1.0), 2 );
244  mesh.add_point( Point(0.0, 1.0), 3 );
245  mesh.add_point( Point(0.0, 0.0), 0 );
246 
247  Elem* elem_top = mesh.add_elem( new QuadShell4 );
248  elem_top->set_node(0) = mesh.node_ptr(0);
249  elem_top->set_node(1) = mesh.node_ptr(1);
250  elem_top->set_node(2) = mesh.node_ptr(2);
251  elem_top->set_node(3) = mesh.node_ptr(3);
252 
253  Elem* elem_bottom = mesh.add_elem( new QuadShell4 );
254  elem_bottom->set_node(0) = mesh.node_ptr(4);
255  elem_bottom->set_node(1) = mesh.node_ptr(5);
256  elem_bottom->set_node(2) = mesh.node_ptr(1);
257  elem_bottom->set_node(3) = mesh.node_ptr(0);
258 
259  BoundaryInfo & bi = mesh.get_boundary_info();
260  bi.add_shellface(elem_top, 0, 10);
261  bi.add_shellface(elem_bottom, 1, 20);
262 
263  mesh.prepare_for_use(false /*skip_renumber*/);
264 
265  CPPUNIT_ASSERT_EQUAL(static_cast<std::size_t>(2), bi.n_shellface_conds());
266 
267  EquationSystems es(mesh);
268  System & system = es.add_system<System> ("SimpleSystem");
269  system.add_variable("u", FIRST);
270 
271  // Add a Dirichlet constraint to check that we impose constraints
272  // correctly on shell faces.
273  std::vector<unsigned int> variables;
274  variables.push_back(0);
275  std::set<boundary_id_type> shellface_ids;
276  shellface_ids.insert(20);
277  ZeroFunction<> zf;
278  DirichletBoundary dirichlet_bc(shellface_ids,
279  variables,
280  &zf);
281  system.get_dof_map().add_dirichlet_boundary(dirichlet_bc);
282  es.init();
283 
284  // Find elem_bottom again if we have it (it may have been deleted
285  // in a DistributedMesh or renumbered in theory)
286  elem_bottom = NULL;
287  for (unsigned int e = 0; e != mesh.max_elem_id(); ++e)
288  {
289  Elem *elem = mesh.query_elem_ptr(e);
290  if (elem && elem->point(3) == Point(0,0))
291  elem_bottom = elem;
292  }
293 
294  // We expect to have a dof constraint on all four dofs of
295  // elem_bottom
296  CPPUNIT_ASSERT_EQUAL(static_cast<std::size_t>(4), static_cast<std::size_t>(system.n_constrained_dofs()));
297 
298  // But we may only know the details of that
299  // constraint on the processor which owns elem_bottom.
300  if (elem_bottom &&
301  elem_bottom->processor_id() == mesh.processor_id())
302  {
303  std::vector<dof_id_type> dof_indices;
304  system.get_dof_map().dof_indices(elem_bottom, dof_indices);
305  CPPUNIT_ASSERT_EQUAL(static_cast<std::size_t>(4), dof_indices.size());
306 
307  for(unsigned int i=0; i<dof_indices.size(); i++)
308  {
309  dof_id_type dof_id = dof_indices[i];
310  CPPUNIT_ASSERT( system.get_dof_map().is_constrained_dof(dof_id) );
311  }
312  }
313  }
314 
315 };
316 
dof_id_type n_constrained_dofs() const
Definition: system.C:155
const BoundaryInfo & get_boundary_info() const
The information about boundary ids on the mesh.
Definition: mesh_base.h:117
This is the EquationSystems class.
virtual void read(const std::string &name, void *mesh_data=libmesh_nullptr, bool skip_renumber_nodes_and_elements=false, bool skip_find_neighbors=false) libmesh_override
Reads the file specified by name.
virtual Node *& set_node(const unsigned int i)
Definition: elem.h:1941
ConstFunction that simply returns 0.
Definition: zero_function.h:35
void testShellFaceConstraints()
virtual bool is_edge_on_side(const unsigned int e, const unsigned int s) const =0
void max(T &r) const
Take a local variable and replace it with the maximum of it&#39;s values on all processors.
dof_id_type n_elem(const MeshBase::const_element_iterator &begin, const MeshBase::const_element_iterator &end)
Count up the number of elements of a specific type (as defined by an iterator range).
Definition: mesh_tools.C:656
libMesh::Parallel::Communicator * TestCommWorld
Definition: driver.C:28
void remove_id(boundary_id_type id)
Removes all entities (nodes, sides, edges, shellfaces) with boundary id id from their respective cont...
The definition of the const_element_iterator struct.
Definition: mesh_base.h:1494
virtual unsigned int n_edges() const =0
unsigned int add_variable(const std::string &var, const FEType &type, const std::set< subdomain_id_type > *const active_subdomains=libmesh_nullptr)
Adds the variable var to the list of variables for this system.
Definition: system.C:1101
CPPUNIT_TEST_SUITE_REGISTRATION(BoundaryInfoTest)
unsigned short int side
Definition: xdr_io.C:49
std::size_t n_edge_conds() const
This is the base class from which all geometric element types are derived.
Definition: elem.h:89
MeshBase & mesh
This class allows one to associate Dirichlet boundary values with a given set of mesh boundary ids an...
void build_square(UnstructuredMesh &mesh, const unsigned int nx, const unsigned int ny, const Real xmin=0., const Real xmax=1., const Real ymin=0., const Real ymax=1., const ElemType type=INVALID_ELEM, const bool gauss_lobatto_grid=false)
A specialized build_cube() for 2D meshes.
The libMesh namespace provides an interface to certain functionality in the library.
dof_id_type dof_id
Definition: xdr_io.C:48
void build_side_list(std::vector< dof_id_type > &element_id_list, std::vector< unsigned short int > &side_list, std::vector< boundary_id_type > &bc_id_list) const
Creates a list of element numbers, sides, and ids for those sides.
const std::set< boundary_id_type > & get_boundary_ids() const
std::size_t n_shellface_conds() const
void testEdgeBoundaryConditions()
bool is_constrained_dof(const dof_id_type dof) const
Definition: dof_map.h:1734
const DofMap & get_dof_map() const
Definition: system.h:2030
int8_t boundary_id_type
Definition: id_types.h:51
virtual const Elem * query_elem_ptr(const dof_id_type i) const libmesh_override
This is the base class for classes which contain information related to any physical process that mig...
Definition: system.h:76
virtual void write(const std::string &name) libmesh_override
Write the file specified by name.
The BoundaryInfo class contains information relevant to boundary conditions including storing faces...
Definition: boundary_info.h:56
virtual bool is_serial() const libmesh_override
virtual const Node * node_ptr(const dof_id_type i) const libmesh_override
void prepare_for_use(const bool skip_renumber_nodes_and_elements=false, const bool skip_find_neighbors=false)
Prepare a newly created (or read) mesh for use.
Definition: mesh_base.C:174
virtual unsigned int n_sides() const =0
QuadShell4 is almost identical to Quad4.
virtual Elem * add_elem(Elem *e) libmesh_override
Add elem e to the end of the element array.
const Point & point(const unsigned int i) const
Definition: elem.h:1809
std::size_t n_boundary_ids() const
const Parallel::Communicator & comm() const
bool has_boundary_id(const Node *const node, const boundary_id_type id) const
void add_shellface(const dof_id_type elem, const unsigned short int shellface, const boundary_id_type id)
Add shell face shellface of element number elem with boundary id id to the boundary information data ...
void add_dirichlet_boundary(const DirichletBoundary &dirichlet_boundary)
Adds a copy of the specified Dirichlet boundary to the system.
The Mesh class is a thin wrapper, around the ReplicatedMesh class by default.
Definition: mesh.h:50
virtual element_iterator elements_begin() libmesh_override
Elem iterator accessor functions.
A Point defines a location in LIBMESH_DIM dimensional Real space.
Definition: point.h:38
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.
virtual Node * add_point(const Point &p, const dof_id_type id=DofObject::invalid_id, const processor_id_type proc_id=DofObject::invalid_processor_id) libmesh_override
functions for adding /deleting nodes elements.
virtual dof_id_type max_elem_id() const libmesh_override
uint8_t dof_id_type
Definition: id_types.h:64
processor_id_type processor_id() const
processor_id_type processor_id() const
Definition: dof_object.h:694
void add_edge(const dof_id_type elem, const unsigned short int edge, const boundary_id_type id)
Add edge edge of element number elem with boundary id id to the boundary information data structure...
void dof_indices(const Elem *const elem, std::vector< dof_id_type > &di) const
Fills the vector di with the global degree of freedom indices for the element.
Definition: dof_map.C:1917
virtual element_iterator elements_end() libmesh_override