[P]arallel [Hi]gh-order [Li]brary for [P]DEs  Latest
Parallel High-Order Library for PDEs through hp-adaptive Discontinuous Galerkin methods
PHiLiP::Physics::ReynoldsAveragedNavierStokesBase< dim, nstate, real > Class Template Referenceabstract

Reynolds-Averaged Navier-Stokes (RANS) equations. Derived from Navier-Stokes for modifying the stress tensor and heat flux, which is derived from PhysicsBase. More...

#include <reynolds_averaged_navier_stokes.h>

Inheritance diagram for PHiLiP::Physics::ReynoldsAveragedNavierStokesBase< dim, nstate, real >:
Collaboration diagram for PHiLiP::Physics::ReynoldsAveragedNavierStokesBase< dim, nstate, real >:

Public Types

using thermal_boundary_condition_enum = Parameters::NavierStokesParam::ThermalBoundaryCondition
 
using two_point_num_flux_enum = Parameters::AllParameters::TwoPointNumericalFlux
 

Public Member Functions

 ReynoldsAveragedNavierStokesBase (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, const double prandtl_number, const double reynolds_number_inf, const bool use_constant_viscosity, const double constant_viscosity, const double turbulent_prandtl_number, const double temperature_inf=273.15, const double isothermal_wall_temperature=1.0, const thermal_boundary_condition_enum thermal_boundary_condition_type=thermal_boundary_condition_enum::adiabatic, 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)
 Constructor.
 
std::array< dealii::Tensor< 1, dim, real >, nstate > convective_flux (const std::array< real, nstate > &conservative_soln) const
 Additional convective flux of RANS + convective flux of turbulence model.
 
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
 Additional viscous flux of RANS + viscous flux of turbulence model.
 
std::array< real, nstate > convective_eigenvalues (const std::array< real, nstate > &, const dealii::Tensor< 1, dim, real > &) const
 Convective eigenvalues of the additional models' PDEs. More...
 
real max_convective_eigenvalue (const std::array< real, nstate > &soln) const
 Maximum convective eigenvalue of the additional models' PDEs. More...
 
real max_convective_normal_eigenvalue (const std::array< real, nstate > &soln, const dealii::Tensor< 1, dim, real > &normal) const
 Maximum convective normal eigenvalue (used in Lax-Friedrichs) of the additional models' PDEs. More...
 
std::array< real, nstate > source_term (const dealii::Point< dim, real > &pos, const std::array< real, nstate > &conservative_solution, const real current_time, const dealii::types::global_dof_index cell_index) const
 Source term for manufactured solution functions.
 
std::array< real, nstate > convective_dissipative_source_term (const dealii::Point< dim, real > &pos, const std::array< real, nstate > &conservative_solution, const dealii::types::global_dof_index cell_index) const
 Convective and dissipative source term for manufactured solution functions.
 
std::array< real, nstate > physical_source_term (const dealii::Point< dim, real > &pos, const std::array< real, nstate > &conservative_solution, const std::array< dealii::Tensor< 1, dim, real >, nstate > &solution_gradient, const dealii::types::global_dof_index cell_index) const override
 Physical source term.
 
virtual dealii::Tensor< 2, dim, real > compute_Reynolds_stress_tensor (const std::array< real, nstate_navier_stokes > &primitive_soln_rans, const std::array< dealii::Tensor< 1, dim, real >, nstate_navier_stokes > &primitive_soln_gradient_rans, const std::array< real, nstate_turbulence_model > &primitive_soln_turbulence_model) const =0
 Nondimensionalized Reynolds stress tensor, (tau^reynolds)*.
 
virtual dealii::Tensor< 1, dim, real > compute_Reynolds_heat_flux (const std::array< real, nstate_navier_stokes > &primitive_soln_rans, const std::array< dealii::Tensor< 1, dim, real >, nstate_navier_stokes > &primitive_soln_gradient_rans, const std::array< real, nstate_turbulence_model > &primitive_soln_turbulence_model) const =0
 Nondimensionalized Reynolds heat flux, (q^reynolds)*.
 
virtual dealii::Tensor< 2, dim, FadTypecompute_Reynolds_stress_tensor_fad (const std::array< FadType, nstate_navier_stokes > &primitive_soln_rans, const std::array< dealii::Tensor< 1, dim, FadType >, nstate_navier_stokes > &primitive_soln_gradient_rans, const std::array< FadType, nstate_turbulence_model > &primitive_soln_turbulence_model) const =0
 Nondimensionalized Reynolds stress tensor, (tau^reynolds)* (Automatic Differentiation Type: FadType)
 
virtual dealii::Tensor< 1, dim, FadTypecompute_Reynolds_heat_flux_fad (const std::array< FadType, nstate_navier_stokes > &primitive_soln_rans, const std::array< dealii::Tensor< 1, dim, FadType >, nstate_navier_stokes > &primitive_soln_gradient_rans, const std::array< FadType, nstate_turbulence_model > &primitive_soln_turbulence_model) const =0
 Nondimensionalized Reynolds heat flux, (q^reynolds)* (Automatic Differentiation Type: FadType)
 
virtual std::array< real, nstate_turbulence_modelcompute_effective_viscosity_turbulence_model (const std::array< real, nstate_navier_stokes > &primitive_soln_rans, const std::array< real, nstate_turbulence_model > &primitive_soln_turbulence_model) const =0
 Nondimensionalized effective (total) viscosities for the turbulence model.
 
virtual std::array< FadType, nstate_turbulence_modelcompute_effective_viscosity_turbulence_model_fad (const std::array< FadType, nstate_navier_stokes > &primitive_soln_rans, const std::array< FadType, nstate_turbulence_model > &primitive_soln_turbulence_model) const =0
 Nondimensionalized effective (total) viscosities for the turbulence model (Automatic Differentiation Type: FadType)
 
virtual std::array< real, nstate > compute_production_dissipation_cross_term (const dealii::Point< dim, real > &pos, const std::array< real, nstate > &conservative_solution, const std::array< dealii::Tensor< 1, dim, real >, nstate > &solution_gradient) const =0
 Physical source term (production, dissipation source terms and source term with cross derivatives) in the turbulence model.
 
- Public Member Functions inherited from PHiLiP::Physics::ModelBase< dim, nstate, real >
 ModelBase (std::shared_ptr< ManufacturedSolutionFunction< dim, real > > manufactured_solution_function_input=nullptr)
 Constructor.
 
virtual ~ModelBase ()=default
 Virtual destructor required for abstract classes.
 
void boundary_face_values (const int boundary_type, const dealii::Point< dim, real > &pos, const dealii::Tensor< 1, dim, real > &normal, 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
 Boundary condition handler.
 
virtual dealii::Vector< double > post_compute_derived_quantities_vector (const dealii::Vector< double > &uh, const std::vector< dealii::Tensor< 1, dim > > &, const std::vector< dealii::Tensor< 2, dim > > &, const dealii::Tensor< 1, dim > &, const dealii::Point< dim > &) const
 Returns current vector solution to be used by PhysicsPostprocessor to output current solution. More...
 
virtual std::vector< dealii::DataComponentInterpretation::DataComponentInterpretation > post_get_data_component_interpretation () const
 Returns DataComponentInterpretation of the solution to be used by PhysicsPostprocessor to output current solution. More...
 
virtual std::vector< std::string > post_get_names () const
 Returns names of the solution to be used by PhysicsPostprocessor to output current solution. More...
 

Public Attributes

const double turbulent_prandtl_number
 Turbulent Prandtl number.
 
std::unique_ptr< NavierStokes< dim, nstate_navier_stokes, real > > navier_stokes_physics
 Pointer to Navier-Stokes physics object.
 
- Public Attributes inherited from PHiLiP::Physics::ModelBase< dim, nstate, real >
std::shared_ptr< ManufacturedSolutionFunction< dim, real > > manufactured_solution_function
 Manufactured solution function.
 
dealii::LinearAlgebra::distributed::Vector< int > cellwise_poly_degree
 Cellwise polynomial degree.
 
dealii::LinearAlgebra::distributed::Vector< double > cellwise_volume
 

Static Public Attributes

static const int nstate_navier_stokes = dim+2
 Number of PDEs for RANS equations.
 
static const int nstate_turbulence_model = nstate-(dim+2)
 Number of PDEs for RANS turbulence model.
 

Protected Member Functions

template<typename real2 >
real2 get_vector_magnitude_sqr (const dealii::Tensor< 1, 3, real2 > &vector) const
 Returns the square of the magnitude of the vector.
 
template<typename real2 >
real2 get_tensor_magnitude_sqr (const dealii::Tensor< 2, dim, real2 > &tensor) const
 Returns the square of the magnitude of the tensor (i.e. the double dot product of a tensor with itself)
 
template<typename real2 >
std::array< dealii::Tensor< 1, dim, real2 >, nstate > dissipative_flux_templated (const std::array< real2, nstate > &conservative_soln, const std::array< dealii::Tensor< 1, dim, real2 >, nstate > &solution_gradient, const dealii::types::global_dof_index cell_index) const
 Templated additional dissipative (i.e. viscous) flux.
 
template<typename real2 >
std::array< dealii::Tensor< 1, dim, real2 >, nstate > convective_flux_templated (const std::array< real2, nstate > &conservative_soln) const
 Templated additional convective flux.
 
template<typename real2 >
std::array< real2, dim+2 > extract_rans_conservative_solution (const std::array< real2, nstate > &conservative_soln) const
 Returns the conservative solutions of Reynolds-averaged Navier-Stokes equations (without additional RANS turbulence model)
 
template<typename real2 >
std::array< dealii::Tensor< 1, dim, real2 >, dim+2 > extract_rans_solution_gradient (const std::array< dealii::Tensor< 1, dim, real2 >, nstate > &solution_gradient) const
 Returns the conservative solutions gradient of Reynolds-averaged Navier-Stokes equations (without additional RANS turbulence model)
 
template<typename real2 >
std::array< dealii::Tensor< 1, dim, real2 >, nstate-(dim+2)> dissipative_flux_turbulence_model (const std::array< real2, nstate_navier_stokes > &primitive_soln_rans, const std::array< real2, nstate_turbulence_model > &primitive_soln_turbulence_model, const std::array< dealii::Tensor< 1, dim, real2 >, nstate_turbulence_model > &primitive_solution_gradient_turbulence_model) const
 Templated Additional viscous flux of RANS + viscous flux of turbulence model.
 
template<typename real2 >
std::array< real2, nstate-(dim+2)> convert_conservative_to_primitive_turbulence_model (const std::array< real2, nstate > &conservative_soln) const
 
template<typename real2 >
std::array< dealii::Tensor< 1, dim, real2 >, nstate-(dim+2)> convert_conservative_gradient_to_primitive_gradient_turbulence_model (const std::array< real2, nstate > &conservative_soln, const std::array< dealii::Tensor< 1, dim, real2 >, nstate > &solution_gradient) const
 
std::array< real, nstate-(dim+2)> compute_mean_turbulence_property (const std::array< real, nstate > &conservative_soln1, const std::array< real, nstate > &conservative_soln2) const
 Mean turbulence properties given two sets of conservative solutions. More...
 
dealii::Tensor< 2, nstate, real > convective_flux_directional_jacobian (const std::array< real, nstate > &conservative_soln, const dealii::Tensor< 1, dim, real > &normal) const
 
dealii::Tensor< 2, nstate, real > dissipative_flux_directional_jacobian (const std::array< real, nstate > &conservative_soln, const std::array< dealii::Tensor< 1, dim, real >, nstate > &solution_gradient, const dealii::Tensor< 1, dim, real > &normal, const dealii::types::global_dof_index cell_index) const
 
dealii::Tensor< 2, nstate, real > dissipative_flux_directional_jacobian_wrt_gradient_component (const std::array< real, nstate > &conservative_soln, const std::array< dealii::Tensor< 1, dim, real >, nstate > &solution_gradient, const dealii::Tensor< 1, dim, real > &normal, const int d_gradient, const dealii::types::global_dof_index cell_index) const
 
std::array< real, nstate > get_manufactured_solution_value (const dealii::Point< dim, real > &pos) const
 Get manufactured solution value.
 
std::array< dealii::Tensor< 1, dim, real >, nstate > get_manufactured_solution_gradient (const dealii::Point< dim, real > &pos) const
 Get manufactured solution value.
 
std::array< real, nstate > convective_source_term_computed_from_manufactured_solution (const dealii::Point< dim, real > &pos) const
 
std::array< real, nstate > dissipative_source_term_computed_from_manufactured_solution (const dealii::Point< dim, real > &pos, const dealii::types::global_dof_index cell_index) const
 
std::array< real, nstate > physical_source_term_computed_from_manufactured_solution (const dealii::Point< dim, real > &pos, const dealii::types::global_dof_index cell_index) const
 
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 override
 Evaluate the manufactured solution boundary conditions.
 
- Protected Member Functions inherited from PHiLiP::Physics::ModelBase< dim, nstate, real >
virtual void boundary_wall (std::array< real, nstate > &soln_bc, std::array< dealii::Tensor< 1, dim, real >, nstate > &soln_grad_bc) const
 Wall boundary condition.
 
virtual void boundary_outflow (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
 Outflow Boundary Condition.
 
virtual void boundary_inflow (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
 Inflow boundary conditions.
 
virtual void boundary_farfield (std::array< real, nstate > &soln_bc) const
 Farfield boundary conditions based on freestream values.
 
virtual 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
 
virtual 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
 

Additional Inherited Members

- Protected Attributes inherited from PHiLiP::Physics::ModelBase< dim, nstate, real >
const MPI_Comm mpi_communicator
 MPI communicator.
 
dealii::ConditionalOStream pcout
 Parallel std::cout that only outputs on mpi_rank==0.
 

Detailed Description

template<int dim, int nstate, typename real>
class PHiLiP::Physics::ReynoldsAveragedNavierStokesBase< dim, nstate, real >

Reynolds-Averaged Navier-Stokes (RANS) equations. Derived from Navier-Stokes for modifying the stress tensor and heat flux, which is derived from PhysicsBase.

Definition at line 13 of file reynolds_averaged_navier_stokes.h.

Member Function Documentation

◆ compute_mean_turbulence_property()

template<int dim, int nstate, typename real >
std::array< real, nstate-(dim+2)> PHiLiP::Physics::ReynoldsAveragedNavierStokesBase< dim, nstate, real >::compute_mean_turbulence_property ( const std::array< real, nstate > &  conservative_soln1,
const std::array< real, nstate > &  conservative_soln2 
) const
protected

Mean turbulence properties given two sets of conservative solutions.

Used in the implementation of the split form.

Definition at line 281 of file reynolds_averaged_navier_stokes.cpp.

◆ convective_eigenvalues()

template<int dim, int nstate, typename real >
std::array< real, nstate > PHiLiP::Physics::ReynoldsAveragedNavierStokesBase< dim, nstate, real >::convective_eigenvalues ( const std::array< real, nstate > &  conservative_soln,
const dealii::Tensor< 1, dim, real > &  normal 
) const
virtual

Convective eigenvalues of the additional models' PDEs.

For RANS model, all entries associated with RANS are assigned to be zero all entries associated with turbulence model are assigned to be the corresponding eigenvalues

Implements PHiLiP::Physics::ModelBase< dim, nstate, real >.

Definition at line 298 of file reynolds_averaged_navier_stokes.cpp.

◆ convective_flux_directional_jacobian()

template<int dim, int nstate, typename real >
dealii::Tensor< 2, nstate, real > PHiLiP::Physics::ReynoldsAveragedNavierStokesBase< dim, nstate, real >::convective_flux_directional_jacobian ( const std::array< real, nstate > &  conservative_soln,
const dealii::Tensor< 1, dim, real > &  normal 
) const
protected

Convective flux Jacobian Note: Only used for computing the manufactured solution source term; computed using automatic differentiation

Definition at line 419 of file reynolds_averaged_navier_stokes.cpp.

◆ convective_source_term_computed_from_manufactured_solution()

template<int dim, int nstate, typename real >
std::array< real, nstate > PHiLiP::Physics::ReynoldsAveragedNavierStokesBase< dim, nstate, real >::convective_source_term_computed_from_manufactured_solution ( const dealii::Point< dim, real > &  pos) const
protected

Convective flux contribution to the source term Note: Only used for computing the manufactured solution source term; computed using automatic differentiation

Definition at line 567 of file reynolds_averaged_navier_stokes.cpp.

◆ convert_conservative_gradient_to_primitive_gradient_turbulence_model()

template<int dim, int nstate, typename real >
template<typename real2 >
template std::array< dealii::Tensor< 1, PHILIP_DIM, FadType >, 1 > PHiLiP::Physics::ReynoldsAveragedNavierStokesBase< dim, nstate, real >::convert_conservative_gradient_to_primitive_gradient_turbulence_model< FadType > ( const std::array< real2, nstate > &  conservative_soln,
const std::array< dealii::Tensor< 1, dim, real2 >, nstate > &  solution_gradient 
) const
protected

Given conservative variable gradients of turbulence model Return primitive variable gradients of turbulence model

Definition at line 263 of file reynolds_averaged_navier_stokes.cpp.

◆ convert_conservative_to_primitive_turbulence_model()

template<int dim, int nstate, typename real >
template<typename real2 >
template std::array< FadType, 1 > PHiLiP::Physics::ReynoldsAveragedNavierStokesBase< dim, nstate, real >::convert_conservative_to_primitive_turbulence_model< FadType > ( const std::array< real2, nstate > &  conservative_soln) const
protected

Given conservative variables of turbulence model Return primitive variables of turbulence model

Definition at line 249 of file reynolds_averaged_navier_stokes.cpp.

◆ dissipative_flux_directional_jacobian()

template<int dim, int nstate, typename real >
dealii::Tensor< 2, nstate, real > PHiLiP::Physics::ReynoldsAveragedNavierStokesBase< dim, nstate, real >::dissipative_flux_directional_jacobian ( const std::array< real, nstate > &  conservative_soln,
const std::array< dealii::Tensor< 1, dim, real >, nstate > &  solution_gradient,
const dealii::Tensor< 1, dim, real > &  normal,
const dealii::types::global_dof_index  cell_index 
) const
protected

Dissipative flux Jacobian Note: Only used for computing the manufactured solution source term; computed using automatic differentiation

Definition at line 452 of file reynolds_averaged_navier_stokes.cpp.

◆ dissipative_flux_directional_jacobian_wrt_gradient_component()

template<int dim, int nstate, typename real >
dealii::Tensor< 2, nstate, real > PHiLiP::Physics::ReynoldsAveragedNavierStokesBase< dim, nstate, real >::dissipative_flux_directional_jacobian_wrt_gradient_component ( const std::array< real, nstate > &  conservative_soln,
const std::array< dealii::Tensor< 1, dim, real >, nstate > &  solution_gradient,
const dealii::Tensor< 1, dim, real > &  normal,
const int  d_gradient,
const dealii::types::global_dof_index  cell_index 
) const
protected

Dissipative flux Jacobian wrt gradient component Note: Only used for computing the manufactured solution source term; computed using automatic differentiation

Definition at line 491 of file reynolds_averaged_navier_stokes.cpp.

◆ dissipative_source_term_computed_from_manufactured_solution()

template<int dim, int nstate, typename real >
std::array< real, nstate > PHiLiP::Physics::ReynoldsAveragedNavierStokesBase< dim, nstate, real >::dissipative_source_term_computed_from_manufactured_solution ( const dealii::Point< dim, real > &  pos,
const dealii::types::global_dof_index  cell_index 
) const
protected

Dissipative flux contribution to the source term Note: Only used for computing the manufactured solution source term; computed using automatic differentiation

This function is same as the one in NavierStokes; TO DO: Reduce this code repetition in the future by moving it elsewhere

Definition at line 601 of file reynolds_averaged_navier_stokes.cpp.

◆ max_convective_eigenvalue()

template<int dim, int nstate, typename real >
real PHiLiP::Physics::ReynoldsAveragedNavierStokesBase< dim, nstate, real >::max_convective_eigenvalue ( const std::array< real, nstate > &  soln) const
virtual

Maximum convective eigenvalue of the additional models' PDEs.

For RANS model, this value is assigned to be the maximum eigenvalue of the turbulence model In most of the RANS models, the maximum eigenvalue of the convective flux is the flow velocity

Implements PHiLiP::Physics::ModelBase< dim, nstate, real >.

Definition at line 316 of file reynolds_averaged_navier_stokes.cpp.

◆ max_convective_normal_eigenvalue()

template<int dim, int nstate, typename real >
real PHiLiP::Physics::ReynoldsAveragedNavierStokesBase< dim, nstate, real >::max_convective_normal_eigenvalue ( const std::array< real, nstate > &  soln,
const dealii::Tensor< 1, dim, real > &  normal 
) const
virtual

Maximum convective normal eigenvalue (used in Lax-Friedrichs) of the additional models' PDEs.

For RANS model, this value is assigned to be the maximum eigenvalue of the turbulence model In most of the RANS models, the maximum normal eigenvalue of the convective flux is the flow normal velocity

Implements PHiLiP::Physics::ModelBase< dim, nstate, real >.

Definition at line 331 of file reynolds_averaged_navier_stokes.cpp.

◆ physical_source_term_computed_from_manufactured_solution()

template<int dim, int nstate, typename real >
std::array< real, nstate > PHiLiP::Physics::ReynoldsAveragedNavierStokesBase< dim, nstate, real >::physical_source_term_computed_from_manufactured_solution ( const dealii::Point< dim, real > &  pos,
const dealii::types::global_dof_index  cell_index 
) const
protected

Physical source contribution to the source term Note: Only used for computing the manufactured solution source term;

Definition at line 673 of file reynolds_averaged_navier_stokes.cpp.


The documentation for this class was generated from the following files: