libMesh
second_order_unsteady_solver_test.C
Go to the documentation of this file.
1 #include <libmesh/equation_systems.h>
2 #include <libmesh/mesh.h>
3 #include <libmesh/mesh_generation.h>
4 #include <libmesh/fem_system.h>
5 #include <libmesh/quadrature.h>
6 #include <libmesh/diff_solver.h>
7 #include <libmesh/newmark_solver.h>
8 #include <libmesh/euler_solver.h>
9 #include <libmesh/euler2_solver.h>
10 
12 
13 
15 template<typename SystemBase>
16 class ConstantSecondOrderODE : public SystemBase
17 {
18 public:
19  ConstantSecondOrderODE(EquationSystems & es,
20  const std::string & name_in,
21  const unsigned int number_in)
22  : SystemBase(es, name_in, number_in)
23  {}
24 
25  virtual Number F( FEMContext & /*context*/, unsigned int /*qp*/ ) override
26  { return -Real(271)/100; }
27 
28  virtual Number C( FEMContext & /*context*/, unsigned int /*qp*/ ) override
29  { return 0.0; }
30 
31  virtual Number M( FEMContext & /*context*/, unsigned int /*qp*/ ) override
32  { return Real(314)/100; }
33 
34  virtual Number u( Real t ) override
35  { return Real(271)/Real(314)*0.5*t*t; }
36 };
37 
39 template<typename SystemBase>
40 class LinearTimeSecondOrderODE : public SystemBase
41 {
42 public:
43  LinearTimeSecondOrderODE(EquationSystems & es,
44  const std::string & name_in,
45  const unsigned int number_in)
46  : SystemBase(es, name_in, number_in)
47  {}
48 
49  virtual Number F( FEMContext & context, unsigned int /*qp*/ ) override
50  { return -6.0*context.get_time()-2.0; }
51 
52  virtual Number C( FEMContext & /*context*/, unsigned int /*qp*/ ) override
53  { return 0.0; }
54 
55  virtual Number M( FEMContext & /*context*/, unsigned int /*qp*/ ) override
56  { return 1.0; }
57 
58  virtual Number u( Real t ) override
59  { return t*t*t+t*t; }
60 };
61 
63 {
64 public:
66  : TimeSolverTestImplementation<NewmarkSolver>(),
67  _beta(0.25)
68  {}
69 
70 protected:
71 
72  virtual void aux_time_solver_init( NewmarkSolver & time_solver ) override
73  { time_solver.set_beta(_beta);
74  time_solver.compute_initial_accel(); }
75 
76  void set_beta( Real beta )
77  { _beta = beta; }
78 
80 };
81 
82 class NewmarkSolverTest : public CppUnit::TestCase,
84 {
85 public:
87 
88 #ifdef LIBMESH_HAVE_SOLVER
93 #endif
94 
96 
97 public:
98 
100  {
101  LOG_UNIT_TEST;
102 
103  this->run_test_with_exact_soln<ConstantSecondOrderODE<SecondOrderScalarSystemSecondOrderTimeSolverBase>>(0.5,10);
104  }
105 
107  {
108  LOG_UNIT_TEST;
109 
110  // For \beta = 1/6, we have the "linear acceleration method" for which
111  // we should be able to exactly integrate linear (in time) acceleration
112  // functions.
113  this->set_beta(Real(1)/Real(6));
114  this->run_test_with_exact_soln<LinearTimeSecondOrderODE<SecondOrderScalarSystemSecondOrderTimeSolverBase>>(0.5,10);
115  }
116 
118  {
119  LOG_UNIT_TEST;
120 
121  this->run_test_with_exact_soln<ConstantSecondOrderODE<SecondOrderScalarSystemFirstOrderTimeSolverBase>>(0.5,10);
122  }
123 
125  {
126  LOG_UNIT_TEST;
127 
128  // For \beta = 1/6, we have the "linear acceleration method" for which
129  // we should be able to exactly integrate linear (in time) acceleration
130  // functions.
131  this->set_beta(Real(1)/Real(6));
132  this->run_test_with_exact_soln<LinearTimeSecondOrderODE<SecondOrderScalarSystemFirstOrderTimeSolverBase>>(0.5,10);
133  }
134 
135 };
136 
137 template<typename TimeSolverType>
138 class ThetaSolverTestBase : public TimeSolverTestImplementation<TimeSolverType>
139 {
140 public:
142  : TimeSolverTestImplementation<TimeSolverType>(),
143  _theta(1.0)
144  {}
145 
146 protected:
147 
148  virtual void aux_time_solver_init( TimeSolverType & time_solver )
149  { time_solver.theta = _theta; }
150 
151  void set_theta( Real theta )
152  { _theta = theta; }
153 
154  Real _theta;
155 };
156 
157 class EulerSolverSecondOrderTest : public CppUnit::TestCase,
158  public ThetaSolverTestBase<EulerSolver>
159 {
160 public:
162 
163 #ifdef LIBMESH_HAVE_SOLVER
165 #endif
166 
168 
169 public:
170 
172  {
173  LOG_UNIT_TEST;
174 
175  this->set_theta(0.5);
176  this->run_test_with_exact_soln<ConstantSecondOrderODE<SecondOrderScalarSystemFirstOrderTimeSolverBase>>(0.5,10);
177  }
178 
179 };
180 
181 class Euler2SolverSecondOrderTest : public CppUnit::TestCase,
182  public ThetaSolverTestBase<Euler2Solver>
183 {
184 public:
186 
187 #ifdef LIBMESH_HAVE_SOLVER
189 #endif
190 
192 
193 public:
195  {
196  LOG_UNIT_TEST;
197 
198  this->set_theta(0.5);
199  this->run_test_with_exact_soln<ConstantSecondOrderODE<SecondOrderScalarSystemFirstOrderTimeSolverBase>>(0.5,10);
200  }
201 };
202 
virtual void aux_time_solver_init(NewmarkSolver &time_solver) override
Implements ODE: 3.14{u} = 2.71, u(0) = 0,.
virtual Number u(Real t) override
CPPUNIT_TEST(testEulerSolverConstantSecondOrderODE)
virtual Number F(FEMContext &context, unsigned int) override
LIBMESH_CPPUNIT_TEST_SUITE(NewmarkSolverTest)
virtual Number F(FEMContext &, unsigned int) override
virtual Number C(FEMContext &, unsigned int) override
CPPUNIT_TEST_SUITE_REGISTRATION(NewmarkSolverTest)
ConstantSecondOrderODE(EquationSystems &es, const std::string &name_in, const unsigned int number_in)
CPPUNIT_TEST(testNewmarkSolverConstantSecondOrderODESecondOrderStyle)
virtual Number M(FEMContext &, unsigned int) override
LIBMESH_CPPUNIT_TEST_SUITE(EulerSolverSecondOrderTest)
LIBMESH_CPPUNIT_TEST_SUITE(Euler2SolverSecondOrderTest)
CPPUNIT_TEST(testEuler2SolverConstantSecondOrderODE)
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
virtual Number C(FEMContext &, unsigned int) override
virtual Number u(Real t) override
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 Number M(FEMContext &, unsigned int) override