4 #include <deal.II/base/tensor.h> 6 #include "parameters/all_parameters.h" 7 #include "parameters/parameters_manufactured_solution.h" 77 template <
int dim,
int nstate,
typename real>
95 const double gamma_gas,
145 const std::array<real,nstate> &conservative_soln)
const override;
148 std::array<real,nstate>
convective_normal_flux (
const std::array<real,nstate> &conservative_soln,
const dealii::Tensor<1,dim,real> &normal)
const;
152 const std::array<real,nstate> &conservative_soln,
153 const dealii::Tensor<1,dim,real> &normal)
const;
157 const std::array<real,nstate> &,
158 const dealii::Tensor<1,dim,real> &)
const override;
166 const std::array<real,nstate> &soln,
167 const dealii::Tensor<1,dim,real> &normal)
const override;
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;
180 const std::array<real,nstate> &conservative_soln,
181 const std::array<dealii::Tensor<1,dim,real>,nstate> &solution_gradient)
const;
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;
192 const dealii::Point<dim,real> &pos,
193 const std::array<real,nstate> &conservative_soln,
194 const real current_time)
const;
198 const dealii::Point<dim,real> &pos)
const;
204 template<
typename real2>
212 template<
typename real2>
222 template<
typename real2>
223 real2
compute_pressure (
const std::array<real2,nstate> &conservative_soln )
const;
226 template<
typename real2>
227 real2
compute_entropy (
const real2 density,
const real2 pressure)
const;
236 real
compute_sound (
const std::array<real,nstate> &conservative_soln )
const;
238 real
compute_sound (
const real density,
const real pressure )
const;
241 template<
typename real2>
242 dealii::Tensor<1,dim,real2>
compute_velocities (
const std::array<real2,nstate> &conservative_soln )
const;
244 template<
typename real2>
248 template<
typename real2>
281 template<
typename real2>
298 const std::array<real,nstate> &conservative_soln1,
299 const std::array<real,nstate> &conservative_soln2)
const override;
305 const std::array<real,nstate> &conservative_soln)
const;
310 const std::array<real,nstate> &entropy_var)
const;
314 const std::array<real,nstate> &conservative_soln)
const;
320 const std::array<real,nstate> &conservative_soln1,
321 const std::array<real,nstate> &convervative_soln2)
const;
327 const std::array<real,nstate> &conservative_soln1,
328 const std::array<real,nstate> &convervative_soln2)
const;
334 const std::array<real,nstate> &conservative_soln1,
335 const std::array<real,nstate> &convervative_soln2)
const;
341 const std::array<real,nstate> &conservative_soln1,
342 const std::array<real,nstate> &convervative_soln2)
const;
347 const dealii::Point<dim, real> &,
348 const dealii::Tensor<1,dim,real> &,
349 const std::array<real,nstate> &,
350 const std::array<dealii::Tensor<1,dim,real>,nstate> &,
351 std::array<real,nstate> &,
352 std::array<dealii::Tensor<1,dim,real>,nstate> &)
const;
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;
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;
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;
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;
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;
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;
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;
428 std::array<real,nstate> &soln_bc)
const;
432 const dealii::Point<dim,real> &pos)
const;
436 const dealii::Point<dim,real> &pos)
const;
441 const std::array<real,nstate> &conservative_soln1,
442 const std::array<real,nstate> &conservative_soln2)
const;
446 const std::array<real,nstate> &primitive_soln)
const;
454 const std::array<real,nstate> &conservative_soln1,
455 const std::array<real,nstate> &conservative_soln2)
const;
459 const std::array<real,nstate> &conservative_soln1,
460 const std::array<real,nstate> &conservative_soln2)
const;
464 const std::array<real,nstate> &conservative_soln1,
465 const std::array<real,nstate> &conservative_soln2)
const;
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 'c'.
real compute_ismail_roe_logarithmic_mean(const real val1, const real val2) const
Compute Ismail-Roe logarithmic mean.
std::array< dealii::Tensor< 1, dim, real >, nstate > get_manufactured_solution_gradient(const dealii::Point< dim, real > &pos) const
Get manufactured solution gradient.
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.
const double mach_inf
Farfield Mach number.
real max_viscous_eigenvalue(const std::array< real, nstate > &soln) const override
Maximum viscous eigenvalue.
std::array< real, nstate > compute_conservative_variables_from_entropy_variables(const std::array< real, nstate > &entropy_var) const
const double angle_of_attack
Angle of attack.
const bool has_nonzero_diffusion
Flag to signal that diffusion term is non-zero.
Base class from which Advection, Diffusion, ConvectionDiffusion, and Euler is derived.
Manufactured solution used for grid studies to check convergence orders.
const double gam
Constant heat capacity ratio of fluid.
void boundary_farfield(std::array< real, nstate > &soln_bc) const
Simple farfield boundary conditions based on freestream values.
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.
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.
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.
const double gamm1
Constant heat capacity ratio (Gamma-1.0) used often.
Files for the baseline physics.
real compute_numerical_entropy_function(const std::array< real, nstate > &conservative_soln) const
Compute numerical entropy function -rho s.
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
const double side_slip_angle
Sideslip angle.
real2 compute_temperature(const std::array< real2, nstate > &primitive_soln) const
Given primitive variables, returns NON-DIMENSIONALIZED temperature using free-stream non-dimensionali...
double temperature_inf
Non-dimensionalized temperature* at infinity. Should equal 1/density*(inf)
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
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: .
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
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
std::array< real, nstate > convective_source_term(const dealii::Point< dim, real > &pos) const
Convective flux contribution to the source term.
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.
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.
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
const double sound_inf
Non-dimensionalized sound* at infinity.
real2 compute_velocity_squared(const dealii::Tensor< 1, dim, real2 > &velocities) const
Given the velocity vector , returns the dot-product .
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.
bool check_positive_quantity(real2 &quantity, const std::string qty_name) const
Check positive quantity and modify it according to handle_non_physical_result()
real compute_total_energy(const std::array< real, nstate > &primitive_soln) const
Given primitive variables, returns total energy.
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)
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.
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...
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
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.
const double entropy_inf
Entropy measure at infinity.
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.
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.
virtual dealii::UpdateFlags post_get_needed_update_flags() const
For post processing purposes (update comment later)
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...
Euler equations. Derived from PhysicsBase.
real compute_sound(const std::array< real, nstate > &conservative_soln) const
Evaluate speed of sound from conservative variables.
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.
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...
dealii::Tensor< 1, dim, real2 > extract_velocities_from_primitive(const std::array< real2, nstate > &primitive_soln) const
Given primitive variables, returns velocities.
std::array< real, nstate > compute_kinetic_energy_variables(const std::array< real, nstate > &conservative_soln) const
Computes the kinetic energy variables.
real max_convective_eigenvalue(const std::array< real, nstate > &soln) const override
Maximum convective eigenvalue.
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.
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...
const bool has_nonzero_physical_source
Flag to signal that physical source term is non-zero.
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)
double dynamic_pressure_inf
Non-dimensionalized dynamic pressure* at infinity.
std::array< dealii::Tensor< 1, dim, real >, nstate > convective_flux(const std::array< real, nstate > &conservative_soln) const override
Convective flux: .
real2 compute_entropy(const real2 density, const real2 pressure) const
Evaluate physical entropy = log(p ^{-}) from pressure and density.
real2 compute_pressure(const std::array< real2, nstate > &conservative_soln) const
Evaluate pressure from conservative variables.
real compute_specific_enthalpy(const std::array< real, nstate > &conservative_soln, const real pressure) const
Evaluate pressure from conservative variables.
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: .
std::array< real, nstate > get_manufactured_solution_value(const dealii::Point< dim, real > &pos) const
Get manufactured solution value.
real compute_kinetic_energy_from_conservative_solution(const std::array< real, nstate > &conservative_soln) const
Given conservative variables, returns kinetic energy.
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...
std::array< real2, nstate > convert_conservative_to_primitive(const std::array< real2, nstate > &conservative_soln) const
dealii::Tensor< 1, dim, real2 > compute_velocities(const std::array< real2, nstate > &conservative_soln) const
Evaluate velocities from conservative variables.
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.
real compute_kinetic_energy_from_primitive_solution(const std::array< real, nstate > &primitive_soln) const
Given primitive variables, returns kinetic energy.
const two_point_num_flux_enum two_point_num_flux_type
Two point numerical flux type (for split form)
const double ref_length
Reference length.
std::array< real, nstate > convert_primitive_to_conservative(const std::array< real, nstate > &primitive_soln) const
std::shared_ptr< ManufacturedSolutionFunction< dim, real > > manufactured_solution_function
Manufactured solution function.
const double pressure_inf
Non-dimensionalized pressure* at infinity.
const double mach_inf_sqr
dealii::Tensor< 1, dim, double > velocities_inf
Non-dimensionalized Velocity vector at farfield.
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.
std::array< real, nstate > compute_entropy_variables(const std::array< real, nstate > &conservative_soln) const
real compute_entropy_measure(const std::array< real, nstate > &conservative_soln) const
Evaluate entropy from conservative variables.