libMesh
parsed_fem_function_test.C
Go to the documentation of this file.
1 // libmesh includes
2 #include "libmesh/elem.h"
3 #include "libmesh/equation_systems.h"
4 #include "libmesh/mesh.h"
5 #include "libmesh/mesh_generation.h"
6 #include "libmesh/numeric_vector.h"
7 #include "libmesh/parsed_fem_function.h"
8 #include "libmesh/system.h"
9 #include <libmesh/auto_ptr.h> // libmesh_make_unique
10 
11 #ifdef LIBMESH_HAVE_FPARSER
12 
13 // test includes
14 #include "test_comm.h"
15 #include "libmesh_cppunit.h"
16 
17 // C++ includes
18 #include <memory>
19 
20 
21 using namespace libMesh;
22 
23 class ParsedFEMFunctionTest : public CppUnit::TestCase
24 {
25 public:
26  void setUp() {
27 #if LIBMESH_DIM > 2
28  mesh = libmesh_make_unique<Mesh>(*TestCommWorld);
29  MeshTools::Generation::build_cube(*mesh, 1, 1, 1);
30  es = libmesh_make_unique<EquationSystems>(*mesh);
31  sys = &(es->add_system<System> ("SimpleSystem"));
32  sys->add_variable("x2");
33  sys->add_variable("x3");
34  sys->add_variable("c05");
35  sys->add_variable("y4");
36  sys->add_variable("xy");
37  sys->add_variable("yz");
38  sys->add_variable("xyz");
39 
40  es->init();
41 
42  NumericVector<Number> & sol = *sys->solution;
43  Elem *elem = mesh->query_elem_ptr(0);
44 
45  if (elem && elem->processor_id() == TestCommWorld->rank())
46  {
47  // Set x2 = 2*x
48  sol.set(elem->node_ref(1).dof_number(0,0,0), 2);
49  sol.set(elem->node_ref(2).dof_number(0,0,0), 2);
50  sol.set(elem->node_ref(5).dof_number(0,0,0), 2);
51  sol.set(elem->node_ref(6).dof_number(0,0,0), 2);
52 
53  // Set x3 = 3*x
54  sol.set(elem->node_ref(1).dof_number(0,1,0), 3);
55  sol.set(elem->node_ref(2).dof_number(0,1,0), 3);
56  sol.set(elem->node_ref(5).dof_number(0,1,0), 3);
57  sol.set(elem->node_ref(6).dof_number(0,1,0), 3);
58 
59  // Set c05 = 0.5
60  sol.set(elem->node_ref(0).dof_number(0,2,0), 0.5);
61  sol.set(elem->node_ref(1).dof_number(0,2,0), 0.5);
62  sol.set(elem->node_ref(2).dof_number(0,2,0), 0.5);
63  sol.set(elem->node_ref(3).dof_number(0,2,0), 0.5);
64  sol.set(elem->node_ref(4).dof_number(0,2,0), 0.5);
65  sol.set(elem->node_ref(5).dof_number(0,2,0), 0.5);
66  sol.set(elem->node_ref(6).dof_number(0,2,0), 0.5);
67  sol.set(elem->node_ref(7).dof_number(0,2,0), 0.5);
68 
69  // Set y4 = 4*y
70  sol.set(elem->node_ref(2).dof_number(0,3,0), 4);
71  sol.set(elem->node_ref(3).dof_number(0,3,0), 4);
72  sol.set(elem->node_ref(6).dof_number(0,3,0), 4);
73  sol.set(elem->node_ref(7).dof_number(0,3,0), 4);
74 
75  // Set xy = x*y
76  sol.set(elem->node_ref(2).dof_number(0,4,0), 1);
77  sol.set(elem->node_ref(6).dof_number(0,4,0), 1);
78 
79  // Set yz = y*z
80  sol.set(elem->node_ref(6).dof_number(0,5,0), 1);
81  sol.set(elem->node_ref(7).dof_number(0,5,0), 1);
82 
83  // Set xyz = x*y*z
84  sol.set(elem->node_ref(6).dof_number(0,6,0), 1);
85  }
86 
87  sol.close();
88  sys->update();
89 
90  c = libmesh_make_unique<FEMContext>(*sys);
91  s = libmesh_make_unique<FEMContext>(*sys);
92  if (elem && elem->processor_id() == TestCommWorld->rank())
93  {
94  c->pre_fe_reinit(*sys, elem);
95  c->elem_fe_reinit();
96  s->pre_fe_reinit(*sys, elem);
97  s->side = 3;
98  s->side_fe_reinit();
99  }
100 #endif
101  }
102 
103  void tearDown() {
104  c.reset();
105  s.reset();
106  es.reset();
107  mesh.reset();
108  }
109 
110  CPPUNIT_TEST_SUITE(ParsedFEMFunctionTest);
111 
112 #if LIBMESH_DIM > 2
113  CPPUNIT_TEST(testValues);
114  CPPUNIT_TEST(testGradients);
115 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
116  CPPUNIT_TEST(testHessians);
117 #endif
118  CPPUNIT_TEST(testInlineGetter);
119  CPPUNIT_TEST(testInlineSetter);
120  CPPUNIT_TEST(testNormals);
121 #endif
122 
123  CPPUNIT_TEST_SUITE_END();
124 
125 
126 private:
127  std::unique_ptr<UnstructuredMesh> mesh;
128  std::unique_ptr<EquationSystems> es;
129  System * sys;
130  std::unique_ptr<FEMContext> c, s;
131 
132  void testValues()
133  {
134  if (c->has_elem() &&
135  c->get_elem().processor_id() == TestCommWorld->rank())
136  {
137  ParsedFEMFunction<Number> x2(*sys, "x2");
138 
139  // Test that copy constructor works
140  ParsedFEMFunction<Number> x2_copy(x2);
141 
142  LIBMESH_ASSERT_FP_EQUAL
143  (libmesh_real(x2_copy(*c,Point(0.5,0.5,0.5))), 1.0, TOLERANCE*TOLERANCE);
144 
145  ParsedFEMFunction<Number> xy8(*sys, "x2*y4");
146 
147  // Test that move constructor works
148  ParsedFEMFunction<Number> xy8_stolen(std::move(xy8));
149 
150  LIBMESH_ASSERT_FP_EQUAL
151  (libmesh_real(xy8_stolen(*c,Point(0.5,0.5,0.5))), 2.0, TOLERANCE*TOLERANCE);
152  }
153  }
154 
155 
156  void testGradients()
157  {
158  if (c->has_elem() &&
159  c->get_elem().processor_id() == TestCommWorld->rank())
160  {
161  ParsedFEMFunction<Number> c2(*sys, "grad_x_x2");
162 
163  // Test that copy/move assignment fails to compile. Note:
164  // ParsedFEMFunction is neither move-assignable nor
165  // copy-assignable because it contains a const reference.
166  // ParsedFEMFunction<Number> c2_assigned(*sys, "grad_y_xyz");
167  // c2_assigned = c2;
168 
169  LIBMESH_ASSERT_FP_EQUAL
170  (libmesh_real(c2(*c,Point(0.35,0.45,0.55))), 2.0, TOLERANCE*TOLERANCE);
171 
172  ParsedFEMFunction<Number> xz(*sys, "grad_y_xyz");
173 
174  LIBMESH_ASSERT_FP_EQUAL
175  (libmesh_real(xz(*c,Point(0.25,0.35,0.75))), 0.1875, TOLERANCE*TOLERANCE);
176 
177  ParsedFEMFunction<Number> xyz(*sys, "grad_y_xyz*grad_x_xy");
178 
179  LIBMESH_ASSERT_FP_EQUAL
180  (libmesh_real(xyz(*c,Point(0.25,0.5,0.75))), 0.09375, TOLERANCE*TOLERANCE);
181  }
182  }
183 
184 
185  void testHessians()
186  {
187  if (c->has_elem() &&
188  c->get_elem().processor_id() == TestCommWorld->rank())
189  {
190  ParsedFEMFunction<Number> c1(*sys, "hess_xy_xy");
191 
192  LIBMESH_ASSERT_FP_EQUAL
193  (libmesh_real(c1(*c,Point(0.35,0.45,0.55))), 1.0, TOLERANCE*TOLERANCE);
194 
195  ParsedFEMFunction<Number> x(*sys, "hess_yz_xyz");
196 
197  LIBMESH_ASSERT_FP_EQUAL
198  (libmesh_real(x(*c,Point(0.25,0.35,0.55))), 0.25, TOLERANCE*TOLERANCE);
199 
200  ParsedFEMFunction<Number> xz(*sys, "hess_yz_xyz*grad_y_yz");
201 
202  LIBMESH_ASSERT_FP_EQUAL
203  (libmesh_real(xz(*c,Point(0.25,0.4,0.75))), 0.1875, TOLERANCE*TOLERANCE);
204  }
205  }
206 
207  void testInlineGetter()
208  {
209  if (c->has_elem() &&
210  c->get_elem().processor_id() == TestCommWorld->rank())
211  {
212  ParsedFEMFunction<Number> ax2(*sys, "a:=4.5;a*x2");
213 
214  LIBMESH_ASSERT_FP_EQUAL
215  (libmesh_real(ax2(*c,Point(0.25,0.25,0.25))), 2.25, TOLERANCE*TOLERANCE);
216 
217  LIBMESH_ASSERT_FP_EQUAL
218  (libmesh_real(ax2.get_inline_value("a")), 4.5, TOLERANCE*TOLERANCE);
219 
220  ParsedFEMFunction<Number> cxy8
221  (*sys, "a := 4 ; b := a/2+1; c:=b-a+3.5; c*x2*y4");
222 
223  LIBMESH_ASSERT_FP_EQUAL
224  (libmesh_real(cxy8(*c,Point(0.5,0.5,0.5))), 5.0, TOLERANCE*TOLERANCE);
225 
226  LIBMESH_ASSERT_FP_EQUAL
227  (libmesh_real(cxy8.get_inline_value("b")), 3.0, TOLERANCE*TOLERANCE);
228 
229  LIBMESH_ASSERT_FP_EQUAL
230  (libmesh_real(cxy8.get_inline_value("c")), 2.5, TOLERANCE*TOLERANCE);
231  }
232  }
233 
234  void testInlineSetter()
235  {
236  if (c->has_elem() &&
237  c->get_elem().processor_id() == TestCommWorld->rank())
238  {
239  ParsedFEMFunction<Number> ax2(*sys, "a:=4.5;a*x2");
240  ax2.set_inline_value("a", 2.5);
241 
242  LIBMESH_ASSERT_FP_EQUAL
243  (libmesh_real(ax2(*c,Point(0.25,0.25,0.25))), 1.25, TOLERANCE*TOLERANCE);
244 
245  LIBMESH_ASSERT_FP_EQUAL
246  (libmesh_real(ax2.get_inline_value("a")), 2.5, TOLERANCE*TOLERANCE);
247 
248  ParsedFEMFunction<Number> cxy8
249  (*sys, "a := 4 ; b := a/2+1; c:=b-a+3.5; c*x2*y4");
250 
251  cxy8.set_inline_value("a", 2);
252 
253  LIBMESH_ASSERT_FP_EQUAL
254  (libmesh_real(cxy8(*c,Point(0.5,0.5,0.5))), 7.0, TOLERANCE*TOLERANCE);
255 
256  LIBMESH_ASSERT_FP_EQUAL
257  (libmesh_real(cxy8.get_inline_value("b")), 2.0, TOLERANCE*TOLERANCE);
258 
259  LIBMESH_ASSERT_FP_EQUAL
260  (libmesh_real(cxy8.get_inline_value("c")), 3.5, TOLERANCE*TOLERANCE);
261 
262  }
263  }
264 
265  void testNormals()
266  {
267  if (s->has_elem() &&
268  s->get_elem().processor_id() == TestCommWorld->rank())
269  {
270  ParsedFEMFunction<Number> nx(*sys, "n_x");
271 
272  ParsedFEMFunction<Number> ny(*sys, "n_y");
273 
274  ParsedFEMFunction<Number> nz(*sys, "n_z");
275 
276  const std::vector<Point> & xyz = s->get_side_fe(0)->get_xyz();
277 
278  // On side 3 of a hex the normal direction is +y
279  for (std::size_t qp=0; qp != xyz.size(); ++qp)
280  {
281  LIBMESH_ASSERT_FP_EQUAL
282  (libmesh_real(nx(*s,xyz[qp])), 0.0, TOLERANCE*TOLERANCE);
283  LIBMESH_ASSERT_FP_EQUAL
284  (libmesh_real(ny(*s,xyz[qp])), 1.0, TOLERANCE*TOLERANCE);
285  LIBMESH_ASSERT_FP_EQUAL
286  (libmesh_real(nz(*s,xyz[qp])), 0.0, TOLERANCE*TOLERANCE);
287  }
288  }
289  }
290 
291 };
292 
293 CPPUNIT_TEST_SUITE_REGISTRATION(ParsedFEMFunctionTest);
294 
295 #endif // #ifdef LIBMESH_HAVE_FPARSER
libMesh::System
Manages consistently variables, degrees of freedom, and coefficient vectors.
Definition: system.h:100
ParsedFEMFunctionTest::es
std::unique_ptr< EquationSystems > es
Definition: parsed_fem_function_test.C:128
libMesh::ParsedFEMFunction::get_inline_value
Output get_inline_value(const std::string &inline_var_name) const
Definition: parsed_fem_function.h:445
libMesh::libmesh_real
T libmesh_real(T a)
Definition: libmesh_common.h:166
libMesh
The libMesh namespace provides an interface to certain functionality in the library.
Definition: factoryfunction.C:55
libMesh::NumericVector::close
virtual void close()=0
Calls the NumericVector's internal assembly routines, ensuring that the values are consistent across ...
ParsedFEMFunctionTest::s
std::unique_ptr< FEMContext > s
Definition: parsed_fem_function_test.C:130
libMesh::MeshTools::Generation::build_cube
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.
Definition: mesh_generation.C:298
ParsedFEMFunctionTest::testNormals
void testNormals()
Definition: parsed_fem_function_test.C:265
libMesh::TOLERANCE
static const Real TOLERANCE
Definition: libmesh_common.h:128
mesh
MeshBase & mesh
Definition: mesh_communication.C:1257
libMesh::DofObject::processor_id
processor_id_type processor_id() const
Definition: dof_object.h:829
ParsedFEMFunctionTest
Definition: parsed_fem_function_test.C:23
libMesh::MeshBase::query_elem_ptr
virtual const Elem * query_elem_ptr(const dof_id_type i) const =0
libMesh::NumericVector< Number >
libMesh::DofObject::dof_number
dof_id_type dof_number(const unsigned int s, const unsigned int var, const unsigned int comp) const
Definition: dof_object.h:956
ParsedFEMFunctionTest::testInlineSetter
void testInlineSetter()
Definition: parsed_fem_function_test.C:234
ParsedFEMFunctionTest::testHessians
void testHessians()
Definition: parsed_fem_function_test.C:185
libMesh::ParsedFEMFunction::set_inline_value
void set_inline_value(const std::string &inline_var_name, Output newval)
Changes the value of an inline variable.
Definition: parsed_fem_function.h:518
ParsedFEMFunctionTest::mesh
std::unique_ptr< UnstructuredMesh > mesh
Definition: parsed_fem_function_test.C:127
TestCommWorld
libMesh::Parallel::Communicator * TestCommWorld
Definition: driver.C:111
libMesh::Point
A Point defines a location in LIBMESH_DIM dimensional Real space.
Definition: point.h:38
ParsedFEMFunctionTest::setUp
void setUp()
Definition: parsed_fem_function_test.C:26
ParsedFEMFunctionTest::testGradients
void testGradients()
Definition: parsed_fem_function_test.C:156
CPPUNIT_TEST_SUITE_REGISTRATION
CPPUNIT_TEST_SUITE_REGISTRATION(ParsedFEMFunctionTest)
ParsedFEMFunctionTest::sys
System * sys
Definition: parsed_fem_function_test.C:129
libMesh::NumericVector::set
virtual void set(const numeric_index_type i, const T value)=0
Sets v(i) = value.
libmesh_cppunit.h
libMesh::ParsedFEMFunction
ParsedFEMFunction provides support for FParser-based parsed functions in FEMSystem.
Definition: parsed_fem_function.h:57
libMesh::Elem
This is the base class from which all geometric element types are derived.
Definition: elem.h:100
ParsedFEMFunctionTest::testValues
void testValues()
Definition: parsed_fem_function_test.C:132
ParsedFEMFunctionTest::testInlineGetter
void testInlineGetter()
Definition: parsed_fem_function_test.C:207
libMesh::Elem::node_ref
const Node & node_ref(const unsigned int i) const
Definition: elem.h:2031
test_comm.h
ParsedFEMFunctionTest::tearDown
void tearDown()
Definition: parsed_fem_function_test.C:103
Generated on Sat Jan 25 2020 12:07:06 for libMesh by   doxygen 1.8.16