MAST
fluid_elem_base.h
Go to the documentation of this file.
1 /*
2  * MAST: Multidisciplinary-design Adaptation and Sensitivity Toolkit
3  * Copyright (C) 2013-2019 Manav Bhatia
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 #ifndef __mast__fluid_elem_base__
21 #define __mast__fluid_elem_base__
22 
23 // C++ includes
24 #include <ostream>
25 #include <map>
26 
27 // MAST includes
28 #include "base/mast_data_types.h"
30 //#include "BoundaryConditions/boundary_surface_motion.h"
31 
32 // libMesh includes
33 #include "libmesh/fe_base.h"
34 
35 
36 
37 namespace MAST {
38 
39  // Forward declerations
40  class FlightCondition;
41  class PrimitiveSolution;
42  template <typename ValType> class SmallPerturbationPrimitiveSolution;
43  class FEBase;
44 
54  };
55 
56 
66  };
67 
68 
73  class FluidElemBase {
74 
75  public:
76 
80  FluidElemBase(const unsigned int dimension,
81  const MAST::FlightCondition& f);
82 
83 
84  virtual ~FluidElemBase();
85 
86 
87  void get_infinity_vars( RealVectorX& vars_inf ) const;
88 
94  //MAST::SurfaceMotionBase* surface_motion;
95 
97 
98  unsigned int dim;
99 
100  bool if_viscous() const {
101  return _if_viscous;
102  }
103 
104  void calculate_dxidX (const unsigned int qp,
105  const MAST::FEBase& fe,
106  RealMatrixX& dxi_dX,
107  RealMatrixX& dX_dxi);
108 
109 
110  void
111  update_solution_at_quadrature_point(const unsigned int qp,
112  const MAST::FEBase& fe,
113  const RealVectorX& elem_solution,
114  RealVectorX& conservative_sol,
115  MAST::PrimitiveSolution& primitive_sol,
117  std::vector<MAST::FEMOperatorMatrix>& dB_mat);
118 
119 
120  void
121  calculate_advection_flux(const unsigned int calculate_dim,
122  const MAST::PrimitiveSolution& sol,
123  RealVectorX& flux);
124 
125  void
126  calculate_diffusion_flux(const unsigned int calculate_dim,
127  const MAST::PrimitiveSolution& sol,
128  const RealMatrixX& stress_tensor,
129  const RealVectorX& temp_gradient,
130  RealVectorX& flux);
131 
135  void
137  const std::vector<MAST::FEMOperatorMatrix>& dB_mat,
138  const RealMatrixX& dprim_dcons,
139  const MAST::PrimitiveSolution& psol,
140  RealMatrixX& stress_tensor,
141  RealVectorX& temp_gradient);
142 
143 
147  void
149  const RealVectorX& elem_sol_sens,
150  const std::vector<MAST::FEMOperatorMatrix>& dB_mat,
151  const RealMatrixX& dprim_dcons,
152  const MAST::PrimitiveSolution& psol,
154  RealMatrixX& stress_tensor_sens,
155  RealVectorX& temp_gradient_sens);
156 
157  void
159  RealMatrixX& dcons_dprim,
160  RealMatrixX& dprim_dcons);
161 
162  void
163  calculate_advection_flux_jacobian(const unsigned int calculate_dim,
164  const MAST::PrimitiveSolution& sol,
165  RealMatrixX& mat);
166 
167  void
168  calculate_advection_flux_jacobian_rho_derivative(const unsigned int calculate_dim,
169  const MAST::PrimitiveSolution& sol,
170  RealMatrixX& mat);
171 
172 
173  void
174  calculate_advection_flux_jacobian_u1_derivative(const unsigned int calculate_dim,
175  const MAST::PrimitiveSolution& sol,
176  RealMatrixX& mat);
177 
178  void
179  calculate_advection_flux_jacobian_u2_derivative(const unsigned int calculate_dim,
180  const MAST::PrimitiveSolution& sol,
181  RealMatrixX& mat);
182 
183 
184  void
185  calculate_advection_flux_jacobian_u3_derivative(const unsigned int calculate_dim,
186  const MAST::PrimitiveSolution& sol,
187  RealMatrixX& mat);
188 
189 
190  void
191  calculate_advection_flux_jacobian_T_derivative(const unsigned int calculate_dim,
192  const MAST::PrimitiveSolution& sol,
193  RealMatrixX& mat);
194 
195 
196  template <typename ValType>
197  void
198  calculate_advection_flux_jacobian_vec_mult_second_derivative(const unsigned int calculate_dim,
199  const MAST::PrimitiveSolution &sol,
200  const ValType& lin_sol,
201  RealMatrixX &mat);
202 
203  template <typename ValType>
204  void
206  const MAST::PrimitiveSolution &sol,
207  const ValType& lin_sol,
208  RealMatrixX &mat);
209 
210 
211  void calculate_diffusion_flux_jacobian(const unsigned int flux_dim,
212  const unsigned int deriv_dim,
213  const MAST::PrimitiveSolution& sol,
214  RealMatrixX& mat);
215 
216 
217  void
218  calculate_diffusion_flux_jacobian_primitive_vars (const unsigned int flux_dim,
219  const unsigned int deriv_dim,
220  const RealVectorX& uvec,
221  const bool zero_kth,
222  const MAST::PrimitiveSolution& sol,
223  RealMatrixX& mat);
224 
225 
227  (const unsigned int calculate_dim,
228  const MAST::PrimitiveSolution& sol,
229  std::vector<RealMatrixX >& mat);
230 
231 
233  (const unsigned int calculate_dim,
234  const unsigned int primitive_var,
235  const MAST::PrimitiveSolution& sol,
236  RealMatrixX& mat);
237 
238 
240  (const MAST::PrimitiveSolution& sol,
241  const libMesh::Point& normal,
242  RealVectorX& eig_vals,
243  RealMatrixX& l_eig_mat,
244  RealMatrixX& l_eig_mat_inv_tr);
245 
246 
248  (const MAST::PrimitiveSolution& sol,
249  const libMesh::Point& normal,
250  RealMatrixX& eig_vals,
251  RealMatrixX& l_eig_mat,
252  RealMatrixX& l_eig_mat_inv_tr);
253 
254 
256  (const MAST::PrimitiveSolution& sol,
257  const libMesh::Point& normal,
258  RealMatrixX& eig_vals,
259  RealMatrixX& l_eig_mat,
260  RealMatrixX& l_eig_mat_inv_tr);
261 
262 
264  (const MAST::PrimitiveSolution& sol,
265  const libMesh::Point& normal,
266  RealMatrixX& eig_vals,
267  RealMatrixX& l_eig_mat,
268  RealMatrixX& l_eig_mat_inv_tr);
269 
270 
272  (const MAST::PrimitiveSolution& sol,
273  const libMesh::Point& normal,
274  RealMatrixX& eig_vals,
275  RealMatrixX& l_eig_mat,
276  RealMatrixX& l_eig_mat_inv_tr);
277 
278 
280  (const MAST::PrimitiveSolution& sol,
281  const libMesh::Point& normal,
282  RealMatrixX& eig_vals,
283  RealMatrixX& l_eig_mat,
284  RealMatrixX& l_eig_mat_inv_tr);
285 
286 
287  template <typename ValType>
288  void
290  const libMesh::Point &normal,
291  const ValType& lin_sol,
292  RealMatrixX &mat);
293 
294 
295  template <typename ValType>
296  void
298  const libMesh::Point &normal,
299  const ValType& lin_sol,
300  RealMatrixX &mat);
301 
302 
304  RealMatrixX& dUdV,
305  RealMatrixX& dVdU);
306 
307 
308 
309  void
311  RealVectorX& dpdX);
312 
313 
314  void
316  const Real ui_ni,
317  const libMesh::Point& nvec,
318  const RealVectorX& dnvec,
319  RealMatrixX& mat);
320 
321 
322  void
324  (const MAST::PrimitiveSolution& sol,
325  const Real ui_ni,
326  const libMesh::Point& nvec,
327  const RealVectorX& dnvec,
328  RealMatrixX& mat);
329 
330  void
332  (const MAST::PrimitiveSolution& sol,
333  const Real ui_ni,
334  const libMesh::Point& nvec,
335  const RealVectorX& dnvec,
336  RealMatrixX& mat);
337 
338  void
340  (const MAST::PrimitiveSolution& sol,
341  const Real ui_ni,
342  const libMesh::Point& nvec,
343  const RealVectorX& dnvec,
344  RealMatrixX& mat);
345 
346  void
348  (const MAST::PrimitiveSolution& sol,
349  const Real ui_ni,
350  const libMesh::Point& nvec,
351  const RealVectorX& dnvec,
352  RealMatrixX& mat);
353 
354  void
356  (const MAST::PrimitiveSolution& sol,
357  const Real ui_ni,
358  const libMesh::Point& nvec,
359  const RealVectorX& dnvec,
360  RealMatrixX& mat);
361 
362 
363  template <typename ValType>
364  void
366  (const MAST::PrimitiveSolution &sol,
367  const Real ui_ni,
368  const libMesh::Point &nvec,
369  const RealVectorX &dnvec,
370  const ValType& lin_sol,
371  RealMatrixX &mat);
372 
373  template <typename ValType>
374  void
376  (const MAST::PrimitiveSolution &sol,
377  const Real ui_ni,
378  const libMesh::Point &nvec,
379  const RealVectorX &dnvec,
380  const ValType& lin_sol,
381  RealMatrixX &mat);
382 
383 
384 
385  bool calculate_barth_tau_matrix(const unsigned int qp,
386  const MAST::FEBase& fe,
387  const MAST::PrimitiveSolution& sol,
388  RealMatrixX& tau,
389  std::vector<RealMatrixX >& tau_sens);
390 
391  bool calculate_aliabadi_tau_matrix(const unsigned int qp,
392  const MAST::FEBase& fe,
393  const MAST::PrimitiveSolution& sol,
394  RealMatrixX& tau,
395  std::vector<RealMatrixX >& tau_sens);
396 
398  (const unsigned int qp,
399  const MAST::FEBase& fe,
400  const MAST::PrimitiveSolution& sol,
401  const RealVectorX& elem_solution,
402  const std::vector<MAST::FEMOperatorMatrix>& dB_mat,
403  const RealMatrixX& Ai_Bi_advection,
404  RealVectorX& discontinuity_val);
405 
406 
408  (const unsigned int qp,
409  const MAST::FEBase& fe,
410  const MAST::PrimitiveSolution& sol,
411  const RealVectorX& elem_solution,
412  const std::vector<MAST::FEMOperatorMatrix>& dB_mat,
413  const RealMatrixX& Ai_Bi_advection,
414  RealVectorX& discontinuity_val);
415 
416 
417  template <typename ValType>
419  (const unsigned int qp,
420  const MAST::FEBase& fe,
421  const MAST::PrimitiveSolution& sol,
423  const RealVectorX& elem_solution,
424  const std::vector<MAST::FEMOperatorMatrix>& dB_mat,
425  const RealMatrixX& Ai_Bi_advection,
426  RealVectorX& discontinuity_val);
427 
428 
430  (const unsigned int qp,
431  const MAST::FEBase& fe,
432  const RealVectorX& elem_solution,
433  const MAST::PrimitiveSolution& sol,
434  const MAST::FEMOperatorMatrix& B_mat,
435  const std::vector<MAST::FEMOperatorMatrix>& dB_mat,
436  const std::vector<RealMatrixX >& Ai_advection,
437  const RealMatrixX& Ai_Bi_advection,
438  const std::vector<std::vector<RealMatrixX > >& Ai_sens,
439  RealMatrixX& LS_operator,
440  RealMatrixX& LS_sens);
441 
442  protected:
443 
444  std::vector<MAST::FluidPrimitiveVars> _active_primitive_vars;
445 
446  std::vector<FluidConservativeVars> _active_conservative_vars;
447 
449 
451 
453  };
454 
455 
456 }
457 
458 #endif /* defined(__mast__fluid_elem_base__) */
void calculate_advection_flux_jacobian_vec_adjoint_mult_second_derivative_for_moving_solid_wall_boundary(const MAST::PrimitiveSolution &sol, const Real ui_ni, const libMesh::Point &nvec, const RealVectorX &dnvec, const ValType &lin_sol, RealMatrixX &mat)
Class defines basic operations and calculation of the small disturbance primitive variables...
void calculate_advection_flux(const unsigned int calculate_dim, const MAST::PrimitiveSolution &sol, RealVectorX &flux)
bool if_viscous() const
void calculate_advection_flux_jacobian_rho_derivative(const unsigned int calculate_dim, const MAST::PrimitiveSolution &sol, RealMatrixX &mat)
void calculate_advection_flux_jacobian_for_moving_solid_wall_boundary(const MAST::PrimitiveSolution &sol, const Real ui_ni, const libMesh::Point &nvec, const RealVectorX &dnvec, RealMatrixX &mat)
void calculate_advection_flux_jacobian_u3_derivative_for_moving_solid_wall_boundary(const MAST::PrimitiveSolution &sol, const Real ui_ni, const libMesh::Point &nvec, const RealVectorX &dnvec, RealMatrixX &mat)
void calculate_advection_left_eigenvector_and_inverse_T_derivative_for_normal(const MAST::PrimitiveSolution &sol, const libMesh::Point &normal, RealMatrixX &eig_vals, RealMatrixX &l_eig_mat, RealMatrixX &l_eig_mat_inv_tr)
void calculate_advection_flux_jacobian_vec_adjoint_mult_second_derivative(const unsigned int calculate_dim, const MAST::PrimitiveSolution &sol, const ValType &lin_sol, RealMatrixX &mat)
FluidElemBase(const unsigned int dimension, const MAST::FlightCondition &f)
Constructor.
This class provides the necessary functions to evaluate the flux vectors and their Jacobians for both...
void calculate_advection_flux_jacobian_T_derivative(const unsigned int calculate_dim, const MAST::PrimitiveSolution &sol, RealMatrixX &mat)
void calculate_diffusion_tensors_sensitivity(const RealVectorX &elem_sol, const RealVectorX &elem_sol_sens, const std::vector< MAST::FEMOperatorMatrix > &dB_mat, const RealMatrixX &dprim_dcons, const MAST::PrimitiveSolution &psol, const MAST::SmallPerturbationPrimitiveSolution< Real > &dsol, RealMatrixX &stress_tensor_sens, RealVectorX &temp_gradient_sens)
calculates and returns the stress tensor in stress_tensor.
Class defines the conversion and some basic operations on primitive fluid variables used in calculati...
void calculate_conservative_variable_jacobian(const MAST::PrimitiveSolution &sol, RealMatrixX &dcons_dprim, RealMatrixX &dprim_dcons)
void calculate_advection_flux_jacobian_sensitivity_for_primitive_variable(const unsigned int calculate_dim, const unsigned int primitive_var, const MAST::PrimitiveSolution &sol, RealMatrixX &mat)
void update_solution_at_quadrature_point(const unsigned int qp, const MAST::FEBase &fe, const RealVectorX &elem_solution, RealVectorX &conservative_sol, MAST::PrimitiveSolution &primitive_sol, MAST::FEMOperatorMatrix &B_mat, std::vector< MAST::FEMOperatorMatrix > &dB_mat)
bool calculate_aliabadi_tau_matrix(const unsigned int qp, const MAST::FEBase &fe, const MAST::PrimitiveSolution &sol, RealMatrixX &tau, std::vector< RealMatrixX > &tau_sens)
void calculate_advection_outflow_flux_jacobian_vec_mult_second_derivative(const MAST::PrimitiveSolution &sol, const libMesh::Point &normal, const ValType &lin_sol, RealMatrixX &mat)
void calculate_advection_left_eigenvector_and_inverse_u1_derivative_for_normal(const MAST::PrimitiveSolution &sol, const libMesh::Point &normal, RealMatrixX &eig_vals, RealMatrixX &l_eig_mat, RealMatrixX &l_eig_mat_inv_tr)
void calculate_advection_flux_jacobian_sensitivity_for_conservative_variable(const unsigned int calculate_dim, const MAST::PrimitiveSolution &sol, std::vector< RealMatrixX > &mat)
libMesh::Real Real
void calculate_diffusion_flux(const unsigned int calculate_dim, const MAST::PrimitiveSolution &sol, const RealMatrixX &stress_tensor, const RealVectorX &temp_gradient, RealVectorX &flux)
void calculate_advection_flux_jacobian_vec_mult_second_derivative(const unsigned int calculate_dim, const MAST::PrimitiveSolution &sol, const ValType &lin_sol, RealMatrixX &mat)
void calculate_diffusion_tensors(const RealVectorX &elem_sol, const std::vector< MAST::FEMOperatorMatrix > &dB_mat, const RealMatrixX &dprim_dcons, const MAST::PrimitiveSolution &psol, RealMatrixX &stress_tensor, RealVectorX &temp_gradient)
calculates and returns the stress tensor in stress_tensor.
void calculate_diffusion_flux_jacobian(const unsigned int flux_dim, const unsigned int deriv_dim, const MAST::PrimitiveSolution &sol, RealMatrixX &mat)
std::vector< FluidConservativeVars > _active_conservative_vars
const MAST::FlightCondition * flight_condition
This defines the surface motion for use with the nonlinear fluid solver.
void calculate_advection_outflow_flux_jacobian_vec_adjoint_mult_second_derivative(const MAST::PrimitiveSolution &sol, const libMesh::Point &normal, const ValType &lin_sol, RealMatrixX &mat)
void calculate_small_disturbance_aliabadi_discontinuity_operator(const unsigned int qp, const MAST::FEBase &fe, const MAST::PrimitiveSolution &sol, const MAST::SmallPerturbationPrimitiveSolution< ValType > &dsol, const RealVectorX &elem_solution, const std::vector< MAST::FEMOperatorMatrix > &dB_mat, const RealMatrixX &Ai_Bi_advection, RealVectorX &discontinuity_val)
void calculate_entropy_variable_jacobian(const MAST::PrimitiveSolution &sol, RealMatrixX &dUdV, RealMatrixX &dVdU)
FluidConservativeVars
enumeration of the conservative fluid variables
void calculate_advection_flux_jacobian(const unsigned int calculate_dim, const MAST::PrimitiveSolution &sol, RealMatrixX &mat)
void calculate_advection_flux_jacobian_u2_derivative_for_moving_solid_wall_boundary(const MAST::PrimitiveSolution &sol, const Real ui_ni, const libMesh::Point &nvec, const RealVectorX &dnvec, RealMatrixX &mat)
Matrix< Real, Dynamic, Dynamic > RealMatrixX
void calculate_advection_flux_jacobian_u2_derivative(const unsigned int calculate_dim, const MAST::PrimitiveSolution &sol, RealMatrixX &mat)
void calculate_advection_left_eigenvector_and_inverse_u3_derivative_for_normal(const MAST::PrimitiveSolution &sol, const libMesh::Point &normal, RealMatrixX &eig_vals, RealMatrixX &l_eig_mat, RealMatrixX &l_eig_mat_inv_tr)
void calculate_aliabadi_discontinuity_operator(const unsigned int qp, const MAST::FEBase &fe, const MAST::PrimitiveSolution &sol, const RealVectorX &elem_solution, const std::vector< MAST::FEMOperatorMatrix > &dB_mat, const RealMatrixX &Ai_Bi_advection, RealVectorX &discontinuity_val)
void calculate_advection_left_eigenvector_and_inverse_rho_derivative_for_normal(const MAST::PrimitiveSolution &sol, const libMesh::Point &normal, RealMatrixX &eig_vals, RealMatrixX &l_eig_mat, RealMatrixX &l_eig_mat_inv_tr)
Matrix< Real, Dynamic, 1 > RealVectorX
void calculate_differential_operator_matrix(const unsigned int qp, const MAST::FEBase &fe, const RealVectorX &elem_solution, const MAST::PrimitiveSolution &sol, const MAST::FEMOperatorMatrix &B_mat, const std::vector< MAST::FEMOperatorMatrix > &dB_mat, const std::vector< RealMatrixX > &Ai_advection, const RealMatrixX &Ai_Bi_advection, const std::vector< std::vector< RealMatrixX > > &Ai_sens, RealMatrixX &LS_operator, RealMatrixX &LS_sens)
void calculate_advection_flux_jacobian_vec_mult_second_derivative_for_moving_solid_wall_boundary(const MAST::PrimitiveSolution &sol, const Real ui_ni, const libMesh::Point &nvec, const RealVectorX &dnvec, const ValType &lin_sol, RealMatrixX &mat)
void calculate_advection_left_eigenvector_and_inverse_u2_derivative_for_normal(const MAST::PrimitiveSolution &sol, const libMesh::Point &normal, RealMatrixX &eig_vals, RealMatrixX &l_eig_mat, RealMatrixX &l_eig_mat_inv_tr)
void calculate_advection_flux_jacobian_T_derivative_for_moving_solid_wall_boundary(const MAST::PrimitiveSolution &sol, const Real ui_ni, const libMesh::Point &nvec, const RealVectorX &dnvec, RealMatrixX &mat)
void calculate_pressure_derivative_wrt_conservative_variables(const MAST::PrimitiveSolution &sol, RealVectorX &dpdX)
FluidPrimitiveVars
enumeration of the primitive fluid variables
void calculate_hartmann_discontinuity_operator(const unsigned int qp, const MAST::FEBase &fe, const MAST::PrimitiveSolution &sol, const RealVectorX &elem_solution, const std::vector< MAST::FEMOperatorMatrix > &dB_mat, const RealMatrixX &Ai_Bi_advection, RealVectorX &discontinuity_val)
void calculate_advection_flux_jacobian_u1_derivative_for_moving_solid_wall_boundary(const MAST::PrimitiveSolution &sol, const Real ui_ni, const libMesh::Point &nvec, const RealVectorX &dnvec, RealMatrixX &mat)
void calculate_dxidX(const unsigned int qp, const MAST::FEBase &fe, RealMatrixX &dxi_dX, RealMatrixX &dX_dxi)
void calculate_advection_flux_jacobian_u1_derivative(const unsigned int calculate_dim, const MAST::PrimitiveSolution &sol, RealMatrixX &mat)
void calculate_advection_flux_jacobian_u3_derivative(const unsigned int calculate_dim, const MAST::PrimitiveSolution &sol, RealMatrixX &mat)
void calculate_diffusion_flux_jacobian_primitive_vars(const unsigned int flux_dim, const unsigned int deriv_dim, const RealVectorX &uvec, const bool zero_kth, const MAST::PrimitiveSolution &sol, RealMatrixX &mat)
void calculate_advection_flux_jacobian_rho_derivative_for_moving_solid_wall_boundary(const MAST::PrimitiveSolution &sol, const Real ui_ni, const libMesh::Point &nvec, const RealVectorX &dnvec, RealMatrixX &mat)
bool calculate_barth_tau_matrix(const unsigned int qp, const MAST::FEBase &fe, const MAST::PrimitiveSolution &sol, RealMatrixX &tau, std::vector< RealMatrixX > &tau_sens)
void calculate_advection_left_eigenvector_and_inverse_for_normal(const MAST::PrimitiveSolution &sol, const libMesh::Point &normal, RealVectorX &eig_vals, RealMatrixX &l_eig_mat, RealMatrixX &l_eig_mat_inv_tr)
std::vector< MAST::FluidPrimitiveVars > _active_primitive_vars
void get_infinity_vars(RealVectorX &vars_inf) const