[P]arallel [Hi]gh-order [Li]brary for [P]DEs  Latest
Parallel High-Order Library for PDEs through hp-adaptive Discontinuous Galerkin methods
euler.h
1 #ifndef __EULER__
2 #define __EULER__
3 
4 #include <deal.II/base/tensor.h>
5 #include "physics.h"
6 #include "parameters/all_parameters.h"
7 #include "parameters/parameters_manufactured_solution.h"
8 
9 namespace PHiLiP {
10 namespace Physics {
11 
13 
77 template <int dim, int nstate, typename real>
78 class Euler : public PhysicsBase <dim, nstate, real>
79 {
80 protected:
81  // For overloading the virtual functions defined in PhysicsBase
89 public:
92  Euler (
93  const Parameters::AllParameters *const parameters_input,
94  const double ref_length,
95  const double gamma_gas,
96  const double mach_inf,
97  const double angle_of_attack,
98  const double side_slip_angle,
100  const two_point_num_flux_enum two_point_num_flux_type = two_point_num_flux_enum::KG,
101  const bool has_nonzero_diffusion = false,
102  const bool has_nonzero_physical_source = false);
103 
104  const double ref_length;
105  const double gam;
106  const double gamm1;
107 
111  const double density_inf;
112 
113  const double mach_inf;
114  const double mach_inf_sqr;
115 
118  const double angle_of_attack;
120 
122  const double side_slip_angle;
123 
124 
125  const double sound_inf;
126  const double pressure_inf;
127  const double entropy_inf;
128  const two_point_num_flux_enum two_point_num_flux_type;
131 
132  //const double internal_energy_inf;
134 
136  dealii::Tensor<1,dim,double> velocities_inf; // should be const
137 
138 
139  // dealii::Tensor<1,dim,double> compute_velocities_inf() const;
140 
141  // std::array<real,nstate> manufactured_solution (const dealii::Point<dim,double> &pos) const;
142 
144  std::array<dealii::Tensor<1,dim,real>,nstate> convective_flux (
145  const std::array<real,nstate> &conservative_soln) const override;
146 
148  std::array<real,nstate> convective_normal_flux (const std::array<real,nstate> &conservative_soln, const dealii::Tensor<1,dim,real> &normal) const;
149 
151  dealii::Tensor<2,nstate,real> convective_flux_directional_jacobian (
152  const std::array<real,nstate> &conservative_soln,
153  const dealii::Tensor<1,dim,real> &normal) const;
154 
156  std::array<real,nstate> convective_eigenvalues (
157  const std::array<real,nstate> &/*conservative_soln*/,
158  const dealii::Tensor<1,dim,real> &/*normal*/) const override;
159 
161  real max_convective_eigenvalue (const std::array<real,nstate> &soln) const override;
162 
164 
166  const std::array<real,nstate> &soln,
167  const dealii::Tensor<1,dim,real> &normal) const override;
168 
170  real max_viscous_eigenvalue (const std::array<real,nstate> &soln) const override;
171 
173  std::array<dealii::Tensor<1,dim,real>,nstate> dissipative_flux (
174  const std::array<real,nstate> &conservative_soln,
175  const std::array<dealii::Tensor<1,dim,real>,nstate> &solution_gradient,
176  const dealii::types::global_dof_index cell_index) const;
177 
179  virtual std::array<dealii::Tensor<1,dim,real>,nstate> dissipative_flux (
180  const std::array<real,nstate> &conservative_soln,
181  const std::array<dealii::Tensor<1,dim,real>,nstate> &solution_gradient) const;
182 
184  std::array<real,nstate> source_term (
185  const dealii::Point<dim,real> &pos,
186  const std::array<real,nstate> &conservative_soln,
187  const real current_time,
188  const dealii::types::global_dof_index cell_index) const;
189 
191  virtual std::array<real,nstate> source_term (
192  const dealii::Point<dim,real> &pos,
193  const std::array<real,nstate> &conservative_soln,
194  const real current_time) const;
195 
197  std::array<real,nstate> convective_source_term (
198  const dealii::Point<dim,real> &pos) const;
199 
200 protected:
202 
204  template<typename real2>
205  bool check_positive_quantity(real2 &quantity, const std::string qty_name) const;
206 
207 public:
212  template<typename real2>
213  std::array<real2,nstate> convert_conservative_to_primitive ( const std::array<real2,nstate> &conservative_soln ) const;
214 
219  std::array<real,nstate> convert_primitive_to_conservative ( const std::array<real,nstate> &primitive_soln ) const;
220 
222  template<typename real2>
223  real2 compute_pressure ( const std::array<real2,nstate> &conservative_soln ) const;
224 
226  template<typename real2>
227  real2 compute_entropy (const real2 density, const real2 pressure) const;
228 
230  real compute_specific_enthalpy ( const std::array<real,nstate> &conservative_soln, const real pressure) const;
231 
233  real compute_numerical_entropy_function(const std::array<real,nstate> &conservative_soln) const;
234 
236  real compute_sound ( const std::array<real,nstate> &conservative_soln ) const;
238  real compute_sound ( const real density, const real pressure ) const;
239 
241  template<typename real2>
242  dealii::Tensor<1,dim,real2> compute_velocities ( const std::array<real2,nstate> &conservative_soln ) const;
244  template<typename real2>
245  real2 compute_velocity_squared ( const dealii::Tensor<1,dim,real2> &velocities ) const;
246 
248  template<typename real2>
249  dealii::Tensor<1,dim,real2> extract_velocities_from_primitive ( const std::array<real2,nstate> &primitive_soln ) const;
251 
254  real compute_total_energy ( const std::array<real,nstate> &primitive_soln ) const;
255 
257  real compute_kinetic_energy_from_primitive_solution ( const std::array<real,nstate> &primitive_soln ) const;
258 
260  real compute_kinetic_energy_from_conservative_solution ( const std::array<real,nstate> &conservative_soln ) const;
261 
263 
271  real compute_entropy_measure ( const std::array<real,nstate> &conservative_soln ) const;
272 
274  real compute_entropy_measure ( const real density, const real pressure ) const;
275 
277  real compute_mach_number ( const std::array<real,nstate> &conservative_soln ) const;
278 
280 
281  template<typename real2>
282  real2 compute_temperature ( const std::array<real2,nstate> &primitive_soln ) const;
283 
285 
286  real compute_density_from_pressure_temperature ( const real pressure, const real temperature ) const;
287 
289 
290  real compute_temperature_from_density_pressure ( const real density, const real pressure ) const;
291 
293 
294  real compute_pressure_from_density_temperature ( const real density, const real temperature ) const;
295 
297  std::array<dealii::Tensor<1,dim,real>,nstate> convective_numerical_split_flux (
298  const std::array<real,nstate> &conservative_soln1,
299  const std::array<real,nstate> &conservative_soln2) const override;
300 
304  std::array<real,nstate> compute_entropy_variables (
305  const std::array<real,nstate> &conservative_soln) const;
306 
310  const std::array<real,nstate> &entropy_var) const;
311 
313  std::array<real,nstate> compute_kinetic_energy_variables (
314  const std::array<real,nstate> &conservative_soln) const;
315 
317 
320  const std::array<real,nstate> &conservative_soln1,
321  const std::array<real,nstate> &convervative_soln2) const;
322 
324 
327  const std::array<real,nstate> &conservative_soln1,
328  const std::array<real,nstate> &convervative_soln2) const;
329 
331 
333  dealii::Tensor<1,dim,real> compute_mean_velocities(
334  const std::array<real,nstate> &conservative_soln1,
335  const std::array<real,nstate> &convervative_soln2) const;
336 
338 
341  const std::array<real,nstate> &conservative_soln1,
342  const std::array<real,nstate> &convervative_soln2) const;
343 
345  void boundary_face_values (
346  const int /*boundary_type*/,
347  const dealii::Point<dim, real> &/*pos*/,
348  const dealii::Tensor<1,dim,real> &/*normal*/,
349  const std::array<real,nstate> &/*soln_int*/,
350  const std::array<dealii::Tensor<1,dim,real>,nstate> &/*soln_grad_int*/,
351  std::array<real,nstate> &/*soln_bc*/,
352  std::array<dealii::Tensor<1,dim,real>,nstate> &/*soln_grad_bc*/) const;
353 
355  virtual dealii::Vector<double> post_compute_derived_quantities_vector (
356  const dealii::Vector<double> &uh,
357  const std::vector<dealii::Tensor<1,dim> > &duh,
358  const std::vector<dealii::Tensor<2,dim> > &dduh,
359  const dealii::Tensor<1,dim> &normals,
360  const dealii::Point<dim> &evaluation_points) const;
361 
363  virtual std::vector<std::string> post_get_names () const;
364 
366  virtual std::vector<dealii::DataComponentInterpretation::DataComponentInterpretation> post_get_data_component_interpretation () const;
367 
369  virtual dealii::UpdateFlags post_get_needed_update_flags () const;
370 
371 protected:
378  void boundary_slip_wall (
379  const dealii::Tensor<1,dim,real> &normal_int,
380  const std::array<real,nstate> &soln_int,
381  const std::array<dealii::Tensor<1,dim,real>,nstate> &soln_grad_int,
382  std::array<real,nstate> &soln_bc,
383  std::array<dealii::Tensor<1,dim,real>,nstate> &soln_grad_bc) const;
384 
386  virtual void boundary_wall (
387  const dealii::Tensor<1,dim,real> &normal_int,
388  const std::array<real,nstate> &soln_int,
389  const std::array<dealii::Tensor<1,dim,real>,nstate> &soln_grad_int,
390  std::array<real,nstate> &soln_bc,
391  std::array<dealii::Tensor<1,dim,real>,nstate> &soln_grad_bc) const;
392 
394  virtual void boundary_manufactured_solution (
395  const dealii::Point<dim, real> &pos,
396  const dealii::Tensor<1,dim,real> &normal_int,
397  const std::array<real,nstate> &soln_int,
398  const std::array<dealii::Tensor<1,dim,real>,nstate> &soln_grad_int,
399  std::array<real,nstate> &soln_bc,
400  std::array<dealii::Tensor<1,dim,real>,nstate> &soln_grad_bc) const;
401 
405  const real total_inlet_pressure,
406  const real back_pressure,
407  const std::array<real,nstate> &soln_int,
408  std::array<real,nstate> &soln_bc) const;
409 
412  void boundary_inflow (
413  const real total_inlet_pressure,
414  const real total_inlet_temperature,
415  const dealii::Tensor<1,dim,real> &normal_int,
416  const std::array<real,nstate> &soln_int,
417  std::array<real,nstate> &soln_bc) const;
418 
421  void boundary_riemann (
422  const dealii::Tensor<1,dim,real> &normal_int,
423  const std::array<real,nstate> &soln_int,
424  std::array<real,nstate> &soln_bc) const;
425 
427  void boundary_farfield (
428  std::array<real,nstate> &soln_bc) const;
429 
431  std::array<real,nstate> get_manufactured_solution_value(
432  const dealii::Point<dim,real> &pos) const;
433 
435  std::array<dealii::Tensor<1,dim,real>,nstate> get_manufactured_solution_gradient(
436  const dealii::Point<dim,real> &pos) const;
437 
440  std::array<dealii::Tensor<1,dim,real>,nstate> convective_numerical_split_flux_kennedy_gruber (
441  const std::array<real,nstate> &conservative_soln1,
442  const std::array<real,nstate> &conservative_soln2) const;
443 
446  const std::array<real,nstate> &primitive_soln) const;
447 
449  real compute_ismail_roe_logarithmic_mean(const real val1, const real val2) const;
450 
453  std::array<dealii::Tensor<1,dim,real>,nstate> convective_numerical_split_flux_ismail_roe (
454  const std::array<real,nstate> &conservative_soln1,
455  const std::array<real,nstate> &conservative_soln2) const;
456 
458  std::array<dealii::Tensor<1,dim,real>,nstate> convective_numerical_split_flux_chandrashekar (
459  const std::array<real,nstate> &conservative_soln1,
460  const std::array<real,nstate> &conservative_soln2) const;
461 
463  std::array<dealii::Tensor<1,dim,real>,nstate> convective_numerical_split_flux_ranocha (
464  const std::array<real,nstate> &conservative_soln1,
465  const std::array<real,nstate> &conservative_soln2) const;
466 };
467 
468 } // Physics namespace
469 } // PHiLiP namespace
470 
471 #endif
std::array< real, nstate > convective_eigenvalues(const std::array< real, nstate > &, const dealii::Tensor< 1, dim, real > &) const override
Spectral radius of convective term Jacobian is &#39;c&#39;.
Definition: euler.cpp:886
real compute_ismail_roe_logarithmic_mean(const real val1, const real val2) const
Compute Ismail-Roe logarithmic mean.
Definition: euler.cpp:501
std::array< dealii::Tensor< 1, dim, real >, nstate > get_manufactured_solution_gradient(const dealii::Point< dim, real > &pos) const
Get manufactured solution gradient.
Definition: euler.cpp:107
Euler(const Parameters::AllParameters *const parameters_input, const double ref_length, const double gamma_gas, const double mach_inf, const double angle_of_attack, const double side_slip_angle, std::shared_ptr< ManufacturedSolutionFunction< dim, real > > manufactured_solution_function=nullptr, const two_point_num_flux_enum two_point_num_flux_type=two_point_num_flux_enum::KG, const bool has_nonzero_diffusion=false, const bool has_nonzero_physical_source=false)
Constructor.
Definition: euler.cpp:13
const double mach_inf
Farfield Mach number.
Definition: euler.h:113
real max_viscous_eigenvalue(const std::array< real, nstate > &soln) const override
Maximum viscous eigenvalue.
Definition: euler.cpp:937
std::array< real, nstate > compute_conservative_variables_from_entropy_variables(const std::array< real, nstate > &entropy_var) const
Definition: euler.cpp:704
const double angle_of_attack
Angle of attack.
Definition: euler.h:118
const bool has_nonzero_diffusion
Flag to signal that diffusion term is non-zero.
Definition: physics.h:59
Base class from which Advection, Diffusion, ConvectionDiffusion, and Euler is derived.
Definition: physics.h:34
Manufactured solution used for grid studies to check convergence orders.
const double gam
Constant heat capacity ratio of fluid.
Definition: euler.h:105
const double density_inf
Definition: euler.h:111
void boundary_farfield(std::array< real, nstate > &soln_bc) const
Simple farfield boundary conditions based on freestream values.
Definition: euler.cpp:1307
std::array< dealii::Tensor< 1, dim, real >, nstate > convective_numerical_split_flux(const std::array< real, nstate > &conservative_soln1, const std::array< real, nstate > &conservative_soln2) const override
Evaluates convective flux based on the chosen split form.
Definition: euler.cpp:440
virtual void boundary_manufactured_solution(const dealii::Point< dim, real > &pos, const dealii::Tensor< 1, dim, real > &normal_int, const std::array< real, nstate > &soln_int, const std::array< dealii::Tensor< 1, dim, real >, nstate > &soln_grad_int, std::array< real, nstate > &soln_bc, std::array< dealii::Tensor< 1, dim, real >, nstate > &soln_grad_bc) const
Evaluate the manufactured solution boundary conditions.
Definition: euler.cpp:1103
std::array< dealii::Tensor< 1, dim, real >, nstate > dissipative_flux(const std::array< real, nstate > &conservative_soln, const std::array< dealii::Tensor< 1, dim, real >, nstate > &solution_gradient, const dealii::types::global_dof_index cell_index) const
Dissipative flux: 0.
Definition: euler.cpp:944
const double gamm1
Constant heat capacity ratio (Gamma-1.0) used often.
Definition: euler.h:106
Files for the baseline physics.
Definition: ADTypes.hpp:10
real compute_numerical_entropy_function(const std::array< real, nstate > &conservative_soln) const
Compute numerical entropy function -rho s.
Definition: euler.cpp:319
void boundary_inflow(const real total_inlet_pressure, const real total_inlet_temperature, const dealii::Tensor< 1, dim, real > &normal_int, const std::array< real, nstate > &soln_int, std::array< real, nstate > &soln_bc) const
Definition: euler.cpp:1200
const double side_slip_angle
Sideslip angle.
Definition: euler.h:122
real2 compute_temperature(const std::array< real2, nstate > &primitive_soln) const
Given primitive variables, returns NON-DIMENSIONALIZED temperature using free-stream non-dimensionali...
Definition: euler.cpp:334
double temperature_inf
Non-dimensionalized temperature* at infinity. Should equal 1/density*(inf)
Definition: euler.h:129
void boundary_slip_wall(const dealii::Tensor< 1, dim, real > &normal_int, const std::array< real, nstate > &soln_int, const std::array< dealii::Tensor< 1, dim, real >, nstate > &soln_grad_int, std::array< real, nstate > &soln_bc, std::array< dealii::Tensor< 1, dim, real >, nstate > &soln_grad_bc) const
Definition: euler.cpp:1042
std::array< real, nstate > convective_normal_flux(const std::array< real, nstate > &conservative_soln, const dealii::Tensor< 1, dim, real > &normal) const
Convective normal flux: .
Definition: euler.cpp:812
std::array< dealii::Tensor< 1, dim, real >, nstate > convective_numerical_split_flux_kennedy_gruber(const std::array< real, nstate > &conservative_soln1, const std::array< real, nstate > &conservative_soln2) const
Definition: euler.cpp:459
void boundary_pressure_outflow(const real total_inlet_pressure, const real back_pressure, const std::array< real, nstate > &soln_int, std::array< real, nstate > &soln_bc) const
Definition: euler.cpp:1160
std::array< real, nstate > convective_source_term(const dealii::Point< dim, real > &pos) const
Convective flux contribution to the source term.
Definition: euler.cpp:123
std::array< real, nstate > source_term(const dealii::Point< dim, real > &pos, const std::array< real, nstate > &conservative_soln, const real current_time, const dealii::types::global_dof_index cell_index) const
Source term is zero or depends on manufactured solution.
Definition: euler.cpp:70
Main parameter class that contains the various other sub-parameter classes.
real compute_mach_number(const std::array< real, nstate > &conservative_soln) const
Given conservative variables, returns Mach number.
Definition: euler.cpp:427
std::array< dealii::Tensor< 1, dim, real >, nstate > convective_numerical_split_flux_ismail_roe(const std::array< real, nstate > &conservative_soln1, const std::array< real, nstate > &conservative_soln2) const
Definition: euler.cpp:522
const double sound_inf
Non-dimensionalized sound* at infinity.
Definition: euler.h:125
real2 compute_velocity_squared(const dealii::Tensor< 1, dim, real2 > &velocities) const
Given the velocity vector , returns the dot-product .
Definition: euler.cpp:236
real compute_mean_specific_total_energy(const std::array< real, nstate > &conservative_soln1, const std::array< real, nstate > &convervative_soln2) const
Mean specific total energy given two sets of conservative solutions.
Definition: euler.cpp:778
bool check_positive_quantity(real2 &quantity, const std::string qty_name) const
Check positive quantity and modify it according to handle_non_physical_result()
Definition: euler.cpp:154
real compute_total_energy(const std::array< real, nstate > &primitive_soln) const
Given primitive variables, returns total energy.
Definition: euler.cpp:258
virtual dealii::Vector< double > post_compute_derived_quantities_vector(const dealii::Vector< double > &uh, const std::vector< dealii::Tensor< 1, dim > > &duh, const std::vector< dealii::Tensor< 2, dim > > &dduh, const dealii::Tensor< 1, dim > &normals, const dealii::Point< dim > &evaluation_points) const
For post processing purposes (update comment later)
Definition: euler.cpp:1373
TwoPointNumericalFlux
Two point numerical flux type for split form.
virtual void boundary_wall(const dealii::Tensor< 1, dim, real > &normal_int, const std::array< real, nstate > &soln_int, const std::array< dealii::Tensor< 1, dim, real >, nstate > &soln_grad_int, std::array< real, nstate > &soln_bc, std::array< dealii::Tensor< 1, dim, real >, nstate > &soln_grad_bc) const
Wall boundary condition.
Definition: euler.cpp:1090
real compute_pressure_from_density_temperature(const real density, const real temperature) const
Given density and temperature, returns NON-DIMENSIONALIZED pressure using free-stream non-dimensional...
Definition: euler.cpp:360
void boundary_riemann(const dealii::Tensor< 1, dim, real > &normal_int, const std::array< real, nstate > &soln_int, std::array< real, nstate > &soln_bc) const
Definition: euler.cpp:968
std::array< dealii::Tensor< 1, dim, real >, nstate > convective_numerical_split_flux_ranocha(const std::array< real, nstate > &conservative_soln1, const std::array< real, nstate > &conservative_soln2) const
Ranocha pressure equilibrium preserving, entropy and energy conserving flux.
Definition: euler.cpp:631
const double entropy_inf
Entropy measure at infinity.
Definition: euler.h:127
real compute_mean_pressure(const std::array< real, nstate > &conservative_soln1, const std::array< real, nstate > &convervative_soln2) const
Mean pressure given two sets of conservative solutions.
Definition: euler.cpp:754
void boundary_face_values(const int, const dealii::Point< dim, real > &, const dealii::Tensor< 1, dim, real > &, const std::array< real, nstate > &, const std::array< dealii::Tensor< 1, dim, real >, nstate > &, std::array< real, nstate > &, std::array< dealii::Tensor< 1, dim, real >, nstate > &) const
Boundary condition handler.
Definition: euler.cpp:1324
virtual dealii::UpdateFlags post_get_needed_update_flags() const
For post processing purposes (update comment later)
Definition: euler.cpp:1478
virtual std::vector< dealii::DataComponentInterpretation::DataComponentInterpretation > post_get_data_component_interpretation() const
For post processing purposes, sets the interpretation of each computed quantity as either scalar or v...
Definition: euler.cpp:1428
Euler equations. Derived from PhysicsBase.
Definition: euler.h:78
real compute_sound(const std::array< real, nstate > &conservative_soln) const
Evaluate speed of sound from conservative variables.
Definition: euler.cpp:407
std::array< dealii::Tensor< 1, dim, real >, nstate > convective_numerical_split_flux_chandrashekar(const std::array< real, nstate > &conservative_soln1, const std::array< real, nstate > &conservative_soln2) const
Chandrashekar entropy conserving flux.
Definition: euler.cpp:582
real compute_temperature_from_density_pressure(const real density, const real pressure) const
Given density and pressure, returns NON-DIMENSIONALIZED temperature using free-stream non-dimensional...
Definition: euler.cpp:352
dealii::Tensor< 1, dim, real2 > extract_velocities_from_primitive(const std::array< real2, nstate > &primitive_soln) const
Given primitive variables, returns velocities.
Definition: euler.cpp:249
std::array< real, nstate > compute_kinetic_energy_variables(const std::array< real, nstate > &conservative_soln) const
Computes the kinetic energy variables.
Definition: euler.cpp:728
real max_convective_eigenvalue(const std::array< real, nstate > &soln) const override
Maximum convective eigenvalue.
Definition: euler.cpp:906
std::array< real, nstate > compute_ismail_roe_parameter_vector_from_primitive(const std::array< real, nstate > &primitive_soln) const
Compute Ismail-Roe parameter vector from primitive solution.
Definition: euler.cpp:486
real compute_density_from_pressure_temperature(const real pressure, const real temperature) const
Given pressure and temperature, returns NON-DIMENSIONALIZED density using free-stream non-dimensional...
Definition: euler.cpp:344
const bool has_nonzero_physical_source
Flag to signal that physical source term is non-zero.
Definition: physics.h:62
real max_convective_normal_eigenvalue(const std::array< real, nstate > &soln, const dealii::Tensor< 1, dim, real > &normal) const override
Maximum convective normal eigenvalue (used in Lax-Friedrichs)
Definition: euler.cpp:921
double dynamic_pressure_inf
Non-dimensionalized dynamic pressure* at infinity.
Definition: euler.h:130
std::array< dealii::Tensor< 1, dim, real >, nstate > convective_flux(const std::array< real, nstate > &conservative_soln) const override
Convective flux: .
Definition: euler.cpp:787
real2 compute_entropy(const real2 density, const real2 pressure) const
Evaluate physical entropy = log(p ^{-}) from pressure and density.
Definition: euler.cpp:388
real2 compute_pressure(const std::array< real2, nstate > &conservative_soln) const
Evaluate pressure from conservative variables.
Definition: euler.cpp:369
real compute_specific_enthalpy(const std::array< real, nstate > &conservative_soln, const real pressure) const
Evaluate pressure from conservative variables.
Definition: euler.cpp:309
dealii::Tensor< 2, nstate, real > convective_flux_directional_jacobian(const std::array< real, nstate > &conservative_soln, const dealii::Tensor< 1, dim, real > &normal) const
Convective flux Jacobian: .
Definition: euler.cpp:839
std::array< real, nstate > get_manufactured_solution_value(const dealii::Point< dim, real > &pos) const
Get manufactured solution value.
Definition: euler.cpp:92
real compute_kinetic_energy_from_conservative_solution(const std::array< real, nstate > &conservative_soln) const
Given conservative variables, returns kinetic energy.
Definition: euler.cpp:279
virtual std::vector< std::string > post_get_names() const
For post processing purposes, sets the base names (with no prefix or suffix) of the computed quantiti...
Definition: euler.cpp:1456
std::array< real2, nstate > convert_conservative_to_primitive(const std::array< real2, nstate > &conservative_soln) const
Definition: euler.cpp:178
dealii::Tensor< 1, dim, real2 > compute_velocities(const std::array< real2, nstate > &conservative_soln) const
Evaluate velocities from conservative variables.
Definition: euler.cpp:225
real compute_mean_density(const std::array< real, nstate > &conservative_soln1, const std::array< real, nstate > &convervative_soln2) const
Mean density given two sets of conservative solutions.
Definition: euler.cpp:746
real compute_kinetic_energy_from_primitive_solution(const std::array< real, nstate > &primitive_soln) const
Given primitive variables, returns kinetic energy.
Definition: euler.cpp:268
const two_point_num_flux_enum two_point_num_flux_type
Two point numerical flux type (for split form)
Definition: euler.h:128
const double ref_length
Reference length.
Definition: euler.h:104
std::array< real, nstate > convert_primitive_to_conservative(const std::array< real, nstate > &primitive_soln) const
Definition: euler.cpp:199
std::shared_ptr< ManufacturedSolutionFunction< dim, real > > manufactured_solution_function
Manufactured solution function.
Definition: physics.h:71
const double pressure_inf
Non-dimensionalized pressure* at infinity.
Definition: euler.h:126
const double mach_inf_sqr
Definition: euler.h:114
dealii::Tensor< 1, dim, double > velocities_inf
Non-dimensionalized Velocity vector at farfield.
Definition: euler.h:136
dealii::Tensor< 1, dim, real > compute_mean_velocities(const std::array< real, nstate > &conservative_soln1, const std::array< real, nstate > &convervative_soln2) const
Mean velocities given two sets of conservative solutions.
Definition: euler.cpp:764
std::array< real, nstate > compute_entropy_variables(const std::array< real, nstate > &conservative_soln) const
Definition: euler.cpp:682
real compute_entropy_measure(const std::array< real, nstate > &conservative_soln) const
Evaluate entropy from conservative variables.
Definition: euler.cpp:288