libMesh
second_order_unsteady_solver_test.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/equation_systems.h>
8 #include <libmesh/mesh.h>
9 #include <libmesh/mesh_generation.h>
10 #include <libmesh/fem_system.h>
11 #include <libmesh/quadrature.h>
12 #include <libmesh/diff_solver.h>
13 #include <libmesh/newmark_solver.h>
14 #include <libmesh/euler_solver.h>
15 #include <libmesh/euler2_solver.h>
16 
17 #include "test_comm.h"
18 
20 
21 // THE CPPUNIT_TEST_SUITE_END macro expands to code that involves
22 // std::auto_ptr, which in turn produces -Wdeprecated-declarations
23 // warnings. These can be ignored in GCC as long as we wrap the
24 // offending code in appropriate pragmas. We can't get away with a
25 // single ignore_warnings.h inclusion at the beginning of this file,
26 // since the libmesh headers pull in a restore_warnings.h at some
27 // point. We also don't bother restoring warnings at the end of this
28 // file since it's not a header.
29 #include <libmesh/ignore_warnings.h>
30 
32 template<typename SystemBase>
33 class ConstantSecondOrderODE : public SystemBase
34 {
35 public:
36  ConstantSecondOrderODE(EquationSystems & es,
37  const std::string & name_in,
38  const unsigned int number_in)
39  : SystemBase(es, name_in, number_in)
40  {}
41 
42  virtual Number F( FEMContext & /*context*/, unsigned int /*qp*/ ) libmesh_override
43  { return -2.71; }
44 
45  virtual Number C( FEMContext & /*context*/, unsigned int /*qp*/ ) libmesh_override
46  { return 0.0; }
47 
48  virtual Number M( FEMContext & /*context*/, unsigned int /*qp*/ ) libmesh_override
49  { return 3.14; }
50 
51  virtual Number u( Real t ) libmesh_override
52  { return 2.71/3.14*0.5*t*t; }
53 };
54 
56 template<typename SystemBase>
57 class LinearTimeSecondOrderODE : public SystemBase
58 {
59 public:
60  LinearTimeSecondOrderODE(EquationSystems & es,
61  const std::string & name_in,
62  const unsigned int number_in)
63  : SystemBase(es, name_in, number_in)
64  {}
65 
66  virtual Number F( FEMContext & context, unsigned int /*qp*/ ) libmesh_override
67  { return -6.0*context.get_time()-2.0; }
68 
69  virtual Number C( FEMContext & /*context*/, unsigned int /*qp*/ ) libmesh_override
70  { return 0.0; }
71 
72  virtual Number M( FEMContext & /*context*/, unsigned int /*qp*/ ) libmesh_override
73  { return 1.0; }
74 
75  virtual Number u( Real t ) libmesh_override
76  { return t*t*t+t*t; }
77 };
78 
80 {
81 public:
83  : TimeSolverTestImplementation<NewmarkSolver>(),
84  _beta(0.25)
85  {}
86 
87 protected:
88 
89  virtual void aux_time_solver_init( NewmarkSolver & time_solver ) libmesh_override
90  { time_solver.set_beta(_beta);
91  time_solver.compute_initial_accel(); }
92 
93  void set_beta( Real beta )
94  { _beta = beta; }
95 
97 };
98 
99 class NewmarkSolverTest : public CppUnit::TestCase,
100  public NewmarkSolverTestBase
101 {
102 public:
103  CPPUNIT_TEST_SUITE( NewmarkSolverTest );
104 
105  CPPUNIT_TEST( testNewmarkSolverConstantSecondOrderODESecondOrderStyle );
106  CPPUNIT_TEST( testNewmarkSolverLinearTimeSecondOrderODESecondOrderStyle );
107  CPPUNIT_TEST( testNewmarkSolverConstantSecondOrderODEFirstOrderStyle );
108  CPPUNIT_TEST( testNewmarkSolverLinearTimeSecondOrderODEFirstOrderStyle );
109 
110  CPPUNIT_TEST_SUITE_END();
111 
112 public:
113 
115  {
116  this->run_test_with_exact_soln<ConstantSecondOrderODE<SecondOrderScalarSystemSecondOrderTimeSolverBase>>(0.5,10);
117  }
118 
120  {
121  // For \beta = 1/6, we have the "linear acceleration method" for which
122  // we should be able to exactly integrate linear (in time) acceleration
123  // functions.
124  this->set_beta(1.0/6.0);
125  this->run_test_with_exact_soln<LinearTimeSecondOrderODE<SecondOrderScalarSystemSecondOrderTimeSolverBase>>(0.5,10);
126  }
127 
129  {
130  this->run_test_with_exact_soln<ConstantSecondOrderODE<SecondOrderScalarSystemFirstOrderTimeSolverBase>>(0.5,10);
131  }
132 
134  {
135  // For \beta = 1/6, we have the "linear acceleration method" for which
136  // we should be able to exactly integrate linear (in time) acceleration
137  // functions.
138  this->set_beta(1.0/6.0);
139  this->run_test_with_exact_soln<LinearTimeSecondOrderODE<SecondOrderScalarSystemFirstOrderTimeSolverBase>>(0.5,10);
140  }
141 
142 };
143 
144 template<typename TimeSolverType>
145 class ThetaSolverTestBase : public TimeSolverTestImplementation<TimeSolverType>
146 {
147 public:
149  : TimeSolverTestImplementation<TimeSolverType>(),
150  _theta(1.0)
151  {}
152 
153 protected:
154 
155  virtual void aux_time_solver_init( TimeSolverType & time_solver )
156  { time_solver.theta = _theta; }
157 
158  void set_theta( Real theta )
159  { _theta = theta; }
160 
161  Real _theta;
162 };
163 
164 class EulerSolverSecondOrderTest : public CppUnit::TestCase,
165  public ThetaSolverTestBase<EulerSolver>
166 {
167 public:
168  CPPUNIT_TEST_SUITE( EulerSolverSecondOrderTest );
169 
170  CPPUNIT_TEST( testEulerSolverConstantSecondOrderODE );
171 
172  CPPUNIT_TEST_SUITE_END();
173 
174 public:
175 
177  {
178  this->set_theta(0.5);
179  this->run_test_with_exact_soln<ConstantSecondOrderODE<SecondOrderScalarSystemFirstOrderTimeSolverBase>>(0.5,10);
180  }
181 
182 };
183 
184 class Euler2SolverSecondOrderTest : public CppUnit::TestCase,
185  public ThetaSolverTestBase<Euler2Solver>
186 {
187 public:
188  CPPUNIT_TEST_SUITE( Euler2SolverSecondOrderTest );
189 
190  CPPUNIT_TEST( testEuler2SolverConstantSecondOrderODE );
191 
192  CPPUNIT_TEST_SUITE_END();
193 
194 public:
196  {
197  this->set_theta(0.5);
198  this->run_test_with_exact_soln<ConstantSecondOrderODE<SecondOrderScalarSystemFirstOrderTimeSolverBase>>(0.5,10);
199  }
200 };
201 
virtual Number M(FEMContext &, unsigned int) libmesh_override
Implements ODE: 3.14{u} = 2.71, u(0) = 0,.
virtual Number C(FEMContext &, unsigned int) libmesh_override
virtual Number F(FEMContext &, unsigned int) libmesh_override
virtual Number F(FEMContext &context, unsigned int) libmesh_override
virtual Number M(FEMContext &, unsigned int) libmesh_override
virtual Number u(Real t) libmesh_override
CPPUNIT_TEST_SUITE_REGISTRATION(NewmarkSolverTest)
virtual Number C(FEMContext &, unsigned int) libmesh_override
ConstantSecondOrderODE(EquationSystems &es, const std::string &name_in, const unsigned int number_in)
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
virtual void aux_time_solver_init(TimeSolverType &time_solver)
Implements ODE: 1.0{u} = 6.0*t+2.0, u(0) = 0,.
LinearTimeSecondOrderODE(EquationSystems &es, const std::string &name_in, const unsigned int number_in)
virtual void aux_time_solver_init(NewmarkSolver &time_solver) libmesh_override
virtual Number u(Real t) libmesh_override