MAST
Multidisciplinary-design Adaptation and Sensitivity Toolkit (MAST)
check_fluid_jacobian.cpp
Go to the documentation of this file.
1 /*
2  * MAST: Multidisciplinary-design Adaptation and Sensitivity Toolkit
3  * Copyright (C) 2013-2020 Manav Bhatia and MAST authors
4  *
5  * This library is free software; you can redistribute it and/or
6  * modify it under the terms of the GNU Lesser General Public
7  * License as published by the Free Software Foundation; either
8  * version 2.1 of the License, or (at your option) any later version.
9  *
10  * This library is distributed in the hope that it will be useful,
11  * but WITHOUT ANY WARRANTY; without even the implied warranty of
12  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
13  * Lesser General Public License for more details.
14  *
15  * You should have received a copy of the GNU Lesser General Public
16  * License along with this library; if not, write to the Free Software
17  * Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA
18  */
19 
20 
21 #define BOOST_TEST_DYN_LINK
22 #define BOOST_TEST_MODULE MAST_TESTS
23 #include <boost/test/unit_test.hpp>
24 #include <boost/version.hpp>
25 
26 // C++ includes
27 #include <memory>
28 
29 // MAST includes
30 #include "base/mast_data_types.h"
31 #include "base/parameter.h"
32 
33 // libMesh includes
34 #include "libmesh/libmesh.h"
35 
36 
37 libMesh::LibMeshInit *_libmesh_init = nullptr;
38 const Real _frac = 1.e-4;
39 const Real _delta = 1.e-4;
40 const Real _tol = 1.e-5;
41 
42 struct GlobalTestFixture {
43 
45 
46  // create the libMeshInit function
48  new libMesh::LibMeshInit(boost::unit_test::framework::master_test_suite().argc,
49  boost::unit_test::framework::master_test_suite().argv);
50  }
51 
53 
54  delete _libmesh_init;
55  }
56 
57 };
58 
59 
60 #if BOOST_VERSION > 106100
61 BOOST_TEST_GLOBAL_FIXTURE( GlobalTestFixture );
62 #else
64 #endif
65 
66 
67 // Test includes
69 #include "base/test_comparisons.h"
70 
71 
72 BOOST_FIXTURE_TEST_SUITE(ConservativeFluidElemJacobians, BuildFluidElem)
73 
74 
75 // sensitivity wrt rho is all zero. Hence, the first one tested is wrt u1
76 BOOST_AUTO_TEST_CASE(FarFieldJacobian) {
77 
78  struct Check {
79  bool jac_xdot;
80  Real frac;
81  Real delta;
82  Real tol;
83  unsigned int side;
84  BuildFluidElem* e;
85  void compute(bool jac, RealVectorX& f, RealMatrixX& j) {
86  e->_fluid_elem->far_field_surface_residual(jac, f, j, side,
87  *e->_far_field_bc);
88  }
89  };
90 
91  Check val;
92  val.jac_xdot = false;
93  val.frac = _frac;
94  val.delta = _delta;
95  val.tol = _tol;
96  val.e = this;
97  this->init(false);
98 
99  // check for each side.
100  for (val.side=0; val.side<_fluid_elem->elem().get_reference_elem().n_sides(); val.side++)
101  BOOST_CHECK(check_jacobian(val));
102 
103 }
104 
105 
106 
107 BOOST_AUTO_TEST_CASE(FarFieldResSens) {
108 
109  struct Check {
110  bool jac_xdot;
111  Real frac;
112  Real delta;
113  Real tol;
114  unsigned int side;
115  BuildFluidElem* e;
116  void compute(bool jac, RealVectorX& f, RealMatrixX& j) {
117  e->_fluid_elem->far_field_surface_residual(jac, f, j, side,
118  *e->_far_field_bc);
119  }
120  };
121 
122  Check val;
123  val.jac_xdot = false;
124  val.frac = _frac;
125  val.delta = _delta;
126  val.tol = _tol;
127  val.e = this;
128  this->init(false);
129 
130  Real
131  rho0 = this->_flight_cond->gas_property.rho;
132 
134  sol = RealVectorX::Zero(this->_sys->n_dofs()),
135  res_sens = RealVectorX::Zero(this->_sys->n_dofs()),
136  res_up = RealVectorX::Zero(this->_sys->n_dofs()),
137  res_lo = RealVectorX::Zero(this->_sys->n_dofs());
138 
140  jac;
141 
142  MAST::Parameter p("dummy", 0.);
143 
144  this->init_solution_for_elem(sol);
145  this->_fluid_elem->set_solution(sol);
146 
147  // check for each side.
148  for (val.side=0; val.side<_fluid_elem->elem().get_reference_elem().n_sides(); val.side++) {
149 
150  res_sens.setZero();
151  res_up.setZero();
152  res_lo.setZero();
153 
154  this->_flight_cond->gas_property.rho = rho0;
155  this->_flight_cond->init();
156  this->_fluid_elem->far_field_surface_residual_sensitivity(p, false, res_sens, jac, val.side, *this->_far_field_bc);
157 
158 
159 
160  this->_flight_cond->gas_property.rho = rho0 * (1. + _frac);
161  this->_flight_cond->init();
162  this->_fluid_elem->far_field_surface_residual(false, res_up, jac, val.side, *this->_far_field_bc);
163 
164  this->_flight_cond->gas_property.rho = rho0 * (1. - _frac);
165  this->_flight_cond->init();
166  this->_fluid_elem->far_field_surface_residual(false, res_lo, jac, val.side, *this->_far_field_bc);
167 
168  res_up -= res_lo;
169  res_up /= (2.*_frac*rho0);
170 
171  BOOST_CHECK(MAST::compare_vector(res_up, res_sens, _tol));
172  }
173 }
174 
175 
176 BOOST_AUTO_TEST_CASE(SlipWallJacobian) {
177 
178  struct Check {
179  bool jac_xdot;
180  Real frac;
181  Real delta;
182  Real tol;
183  unsigned int side;
184  BuildFluidElem* e;
185  void compute(bool jac, RealVectorX& f, RealMatrixX& j) {
186  e->_fluid_elem->slip_wall_surface_residual(jac, f, j, side,
187  *e->_slip_wall_bc);
188  }
189  };
190 
191  Check val;
192  val.jac_xdot = false;
193  val.frac = _frac;
194  val.delta = _delta;
195  val.tol = _tol;
196  val.e = this;
197  this->init(false);
198 
199  // check for each side.
200  for (val.side=0; val.side<_fluid_elem->elem().get_reference_elem().n_sides(); val.side++)
201  BOOST_CHECK(check_jacobian(val));
202 
203 }
204 
205 
206 
207 BOOST_AUTO_TEST_CASE(NoSlipWallJacobian) {
208 
209  struct Check {
210  bool jac_xdot;
211  Real frac;
212  Real delta;
213  Real tol;
214  unsigned int side;
215  BuildFluidElem* e;
216  void compute(bool jac, RealVectorX& f, RealMatrixX& j) {
217  e->_fluid_elem->noslip_wall_surface_residual(jac, f, j, side,
218  *e->_noslip_wall_bc);
219  }
220  };
221 
222  Check val;
223  val.jac_xdot = false;
224  val.frac = _frac;
225  val.delta = _delta;
226  val.tol = _tol;
227  val.e = this;
228  this->_delta = 0.e-4;
229  this->init(true);
230 
231  // check for each side.
232  for (val.side=0; val.side<_fluid_elem->elem().get_reference_elem().n_sides(); val.side++)
233  BOOST_CHECK(check_jacobian(val));
234 
235 }
236 
237 
238 BOOST_AUTO_TEST_CASE(InternalResidualJacobianInviscid) {
239 
240  struct Check {
241  bool jac_xdot;
242  Real frac;
243  Real delta;
244  Real tol;
245  BuildFluidElem* e;
246  void compute(bool jac, RealVectorX& f, RealMatrixX& j) {
247  e->_fluid_elem->internal_residual(jac, f, j);
248  }
249  };
250 
251  Check val;
252  val.jac_xdot = false;
253  val.frac = _frac;
254  val.delta = _delta;
255  // a smaller tolerance is required for the internal resisudal since
256  // an exact Jacobian is not computed for the stabilization terms
257  val.tol = 1.e-2;
258  val.e = this;
259  this->init(false);
260 
261  BOOST_CHECK(check_jacobian(val));
262 
263 }
264 
265 
266 BOOST_AUTO_TEST_CASE(InternalResidualJacobianViscous) {
267 
268  struct Check {
269  bool jac_xdot;
270  Real frac;
271  Real delta;
272  Real tol;
273  BuildFluidElem* e;
274  void compute(bool jac, RealVectorX& f, RealMatrixX& j) {
275  e->_fluid_elem->internal_residual(jac, f, j);
276  }
277  };
278 
279  Check val;
280  val.jac_xdot = false;
281  val.frac = _frac;
282  val.delta = _delta;
283  // a smaller tolerance is required for the internal resisudal since
284  // an exact Jacobian is not computed for the stabilization terms
285  val.tol = 1.e-2;
286  val.e = this;
287  this->_delta = 1.e-4;
288  this->init(true);
289 
290  BOOST_CHECK(check_jacobian(val));
291 
292 }
293 
294 
295 BOOST_AUTO_TEST_CASE(VelocityResidualXJacobian) {
296 
297  struct Check {
298  bool jac_xdot;
299  Real frac;
300  Real delta;
301  Real tol;
302  BuildFluidElem* e;
303  void compute(bool jac, RealVectorX& f, RealMatrixX& j) {
304  RealMatrixX jj = j; // dummy for jac wrt x_dot
305  e->_fluid_elem->velocity_residual(jac, f, jj, j);
306  }
307  };
308 
309  Check val;
310  val.jac_xdot = false;
311  val.frac = _frac;
312  val.delta = _delta;
313  val.tol = _tol;
314  val.e = this;
315  this->init(false);
316 
317  //BOOST_CHECK(check_jacobian(val));
318 
319 }
320 
321 
322 BOOST_AUTO_TEST_CASE(VelocityResidualXdotJacobian) {
323 
324  struct Check {
325  bool jac_xdot;
326  Real frac;
327  Real delta;
328  Real tol;
329  BuildFluidElem* e;
330  void compute(bool jac, RealVectorX& f, RealMatrixX& j) {
331  RealMatrixX jj = j; // dummy for jac wrt x
332  e->_fluid_elem->velocity_residual(jac, f, j, jj);
333  }
334  };
335 
336  Check val;
337  val.jac_xdot = true;
338  val.frac = _frac;
339  val.delta = _delta;
340  val.tol = _tol;
341  val.e = this;
342  this->init(false);
343 
344  BOOST_CHECK(check_jacobian(val));
345 
346 }
347 
348 
349 BOOST_AUTO_TEST_SUITE_END()
libMesh::LibMeshInit * _libmesh_init
bool compare_vector(const RealVectorX &v0, const RealVectorX &v, const Real tol)
This is a scalar function whose value can be changed and one that can be used as a design variable in...
Definition: parameter.h:35
MAST::BoundaryConditionBase * _far_field_bc
BOOST_GLOBAL_FIXTURE(GlobalTestFixture)
virtual bool velocity_residual(bool request_jacobian, RealVectorX &f, RealMatrixX &jac_xdot, RealMatrixX &jac)
inertial force contribution to system residual
BOOST_FIXTURE_TEST_SUITE(PanelSmallDisturbanceFrequencyDomain2D, MAST::PanelInviscidSmallDisturbanceFrequencyDomain2DAnalysis) BOOST_AUTO_TEST_CASE(FreqDomainSensitivityWrtOmega)
libMesh::Real Real
virtual bool slip_wall_surface_residual(bool request_jacobian, RealVectorX &f, RealMatrixX &jac, const unsigned int s, MAST::BoundaryConditionBase &bc)
Matrix< Real, Dynamic, Dynamic > RealMatrixX
BOOST_AUTO_TEST_CASE(FarFieldJacobian)
MAST::ConservativeFluidElementBase * _fluid_elem
const Real _frac
const Real _tol
Matrix< Real, Dynamic, 1 > RealVectorX
virtual bool noslip_wall_surface_residual(bool request_jacobian, RealVectorX &f, RealMatrixX &jac, const unsigned int s, MAST::BoundaryConditionBase &bc)
virtual bool internal_residual(bool request_jacobian, RealVectorX &f, RealMatrixX &jac)
internal force contribution to system residual
MAST::BoundaryConditionBase * _noslip_wall_bc
virtual bool far_field_surface_residual(bool request_jacobian, RealVectorX &f, RealMatrixX &jac, const unsigned int s, MAST::BoundaryConditionBase &bc)
MAST::BoundaryConditionBase * _slip_wall_bc
const Real _delta