libMesh
contains_point.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/libmesh.h>
8 #include <libmesh/elem.h>
9 
10 #include "test_comm.h"
11 
12 // THE CPPUNIT_TEST_SUITE_END macro expands to code that involves
13 // std::auto_ptr, which in turn produces -Wdeprecated-declarations
14 // warnings. These can be ignored in GCC as long as we wrap the
15 // offending code in appropriate pragmas. We can't get away with a
16 // single ignore_warnings.h inclusion at the beginning of this file,
17 // since the libmesh headers pull in a restore_warnings.h at some
18 // point. We also don't bother restoring warnings at the end of this
19 // file since it's not a header.
20 #include <libmesh/ignore_warnings.h>
21 
22 using namespace libMesh;
23 
24 class ContainsPointTest : public CppUnit::TestCase
25 {
31 public:
32  CPPUNIT_TEST_SUITE( ContainsPointTest );
33 
34  CPPUNIT_TEST( testContainsPointTri3 );
35 
36  CPPUNIT_TEST( testContainsPointTet4 );
37 
38  CPPUNIT_TEST_SUITE_END();
39 
40 public:
41  void setUp() {}
42 
43  void tearDown() {}
44 
45  // TRI3 test
47  {
48  Point a(3,1,2), b(1,2,3), c(2,3,1);
49  containsPointTri3Helper(a, b, c, a);
50 
51  // Ben's 1st "small triangle" test case.
52  containsPointTri3Helper(a/5000., b/5000., c/5000., c/5000.);
53 
54  // Ben's 2nd "small triangle" test case.
55  containsPointTri3Helper(Point(0.000808805, 0.0047284, 0.),
56  Point(0.0011453, 0.00472831, 0.),
57  Point(0.000982249, 0.00508037, 0.),
58  Point(0.001, 0.005, 0.));
59  }
60 
61 
62 
63  // TET4 test
65  {
66  // Construct unit Tet.
67  {
68  Node zero (0., 0., 0., 0);
69  Node one (1., 0., 0., 1);
70  Node two (0., 1., 0., 2);
71  Node three (0., 0., 1., 3);
72 
74  elem->set_node(0) = &zero;
75  elem->set_node(1) = &one;
76  elem->set_node(2) = &two;
77  elem->set_node(3) = &three;
78 
79  // The centroid must be inside the element.
80  CPPUNIT_ASSERT (elem->contains_point(elem->centroid()));
81 
82  // The vertices must be contained in the element.
83  CPPUNIT_ASSERT (elem->contains_point(zero));
84  CPPUNIT_ASSERT (elem->contains_point(one));
85  CPPUNIT_ASSERT (elem->contains_point(two));
86  CPPUNIT_ASSERT (elem->contains_point(three));
87 
88  // Make sure that outside points are not contained.
89  CPPUNIT_ASSERT (!elem->contains_point(Point(.34, .34, .34)));
90  CPPUNIT_ASSERT (!elem->contains_point(Point(.33, .33, -.1)));
91  CPPUNIT_ASSERT (!elem->contains_point(Point(0., -.1, .5)));
92  }
93 
94 
95  // Construct a nearly degenerate (sliver) tet. A unit tet with
96  // nodes "one" and "two" moved to within a distance of epsilon
97  // from the origin.
98  {
99  Real epsilon = 1.e-4;
100 
101  Node zero (0., 0., 0., 0);
102  Node one (epsilon, 0., 0., 1);
103  Node two (0., epsilon, 0., 2);
104  Node three (0., 0., 1., 3);
105 
107  elem->set_node(0) = &zero;
108  elem->set_node(1) = &one;
109  elem->set_node(2) = &two;
110  elem->set_node(3) = &three;
111 
112  // The centroid must be inside the element.
113  CPPUNIT_ASSERT (elem->contains_point(elem->centroid()));
114 
115  // The vertices must be contained in the element.
116  CPPUNIT_ASSERT (elem->contains_point(zero));
117  CPPUNIT_ASSERT (elem->contains_point(one));
118  CPPUNIT_ASSERT (elem->contains_point(two));
119  CPPUNIT_ASSERT (elem->contains_point(three));
120 
121  // Verify that a point which is on a mid-edge is contained in the element.
122  CPPUNIT_ASSERT (elem->contains_point(Point(epsilon/2, 0, 0.5)));
123 
124  // Make sure that "just barely" outside points are outside.
125  CPPUNIT_ASSERT (!elem->contains_point(Point(epsilon, epsilon, epsilon/2)));
126  CPPUNIT_ASSERT (!elem->contains_point(Point(epsilon/10, epsilon/10, 1.0)));
127  CPPUNIT_ASSERT (!elem->contains_point(Point(epsilon/2, -epsilon/10, 0.5)));
128  }
129  }
130 
131 protected:
132  // The first 3 arguments are the points of the triangle, the last argument is a
133  // custom point that should be inside the Tri.
134  void containsPointTri3Helper(Point a_in, Point b_in, Point c_in, Point p)
135  {
136  // vertices
137  Node a(a_in, 0);
138  Node b(b_in, 1);
139  Node c(c_in, 2);
140 
141  // Build the test Elem
143 
144  elem->set_node(0) = &a;
145  elem->set_node(1) = &b;
146  elem->set_node(2) = &c;
147 
148  // helper vectors to span the tri and point out of plane
149  Point va(a-c);
150  Point vb(b-c);
151  Point oop(va.cross(vb));
152  Point oop_unit = oop.unit();
153 
154  // The triangle must contain its own centroid
155  CPPUNIT_ASSERT (elem->contains_point(elem->centroid()));
156 
157  // triangle should contain all its vertices
158  CPPUNIT_ASSERT (elem->contains_point(a));
159  CPPUNIT_ASSERT (elem->contains_point(b));
160  CPPUNIT_ASSERT (elem->contains_point(c));
161 
162  // out of plane from the centroid, should *not* be in the element
163  CPPUNIT_ASSERT (!elem->contains_point(elem->centroid() + std::sqrt(TOLERANCE) * oop_unit));
164 
165  // These points should be outside the triangle
166  CPPUNIT_ASSERT (!elem->contains_point(a + va * TOLERANCE * 10));
167  CPPUNIT_ASSERT (!elem->contains_point(b + vb * TOLERANCE * 10));
168  CPPUNIT_ASSERT (!elem->contains_point(c - (va + vb) * TOLERANCE * 10));
169 
170  // Test the custom point
171  CPPUNIT_ASSERT (elem->contains_point(p));
172  }
173 };
174 
175 
void containsPointTri3Helper(Point a_in, Point b_in, Point c_in, Point p)
CPPUNIT_TEST_SUITE_REGISTRATION(ContainsPointTest)
A Node is like a Point, but with more information.
Definition: node.h:52
static UniquePtr< Elem > build(const ElemType type, Elem *p=libmesh_nullptr)
Definition: elem.C:238
static const Real TOLERANCE
The libMesh namespace provides an interface to certain functionality in the library.
const Number zero
.
Definition: libmesh.h:178
std::unique_ptr< T > UniquePtr
Definition: auto_ptr.h:46
void testContainsPointTet4()
TypeVector< typename CompareTypes< T, T2 >::supertype > cross(const TypeVector< T2 > &v) const
Definition: type_vector.h:879
void testContainsPointTri3()
TypeVector< T > unit() const
Definition: type_vector.C:37
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
A Point defines a location in LIBMESH_DIM dimensional Real space.
Definition: point.h:38