[P]arallel [Hi]gh-order [Li]brary for [P]DEs  Latest
Parallel High-Order Library for PDEs through hp-adaptive Discontinuous Galerkin methods
PHiLiP::FlowSolver::PeriodicTurbulence< dim, nstate > Class Template Reference
Inheritance diagram for PHiLiP::FlowSolver::PeriodicTurbulence< dim, nstate >:
Collaboration diagram for PHiLiP::FlowSolver::PeriodicTurbulence< dim, nstate >:

Public Member Functions

 PeriodicTurbulence (const Parameters::AllParameters *const parameters_input)
 Constructor. More...
 
void compute_and_update_integrated_quantities (DGBase< dim, double > &dg)
 
double get_integrated_kinetic_energy () const
 
double get_integrated_enstrophy () const
 
double get_vorticity_based_dissipation_rate () const
 
double get_pressure_dilatation_based_dissipation_rate () const
 
double get_deviatoric_strain_rate_tensor_based_dissipation_rate () const
 
double get_strain_rate_tensor_based_dissipation_rate () const
 
void output_velocity_field (std::shared_ptr< DGBase< dim, double >> dg, const unsigned int output_file_index, const double current_time) const
 Output the velocity field to file.
 
double compute_current_integrated_numerical_entropy (const std::shared_ptr< DGBase< dim, double >> dg) const
 Calculate numerical entropy by matrix-vector product.
 
void update_numerical_entropy (const double FR_entropy_contribution_RRK_solver, const unsigned int current_iteration, const std::shared_ptr< DGBase< dim, double >> dg)
 Update numerical entropy variables.
 
double get_numerical_entropy (const std::shared_ptr< DGBase< dim, double >>) const
 Retrieves cumulative_numerical_entropy_change_FRcorrected.
 
- Public Member Functions inherited from PHiLiP::FlowSolver::PeriodicCubeFlow< dim, nstate >
 PeriodicCubeFlow (const Parameters::AllParameters *const parameters_input)
 Constructor.
 
std::shared_ptr< Triangulation > generate_grid () const override
 Function to generate the grid.
 
- Public Member Functions inherited from PHiLiP::FlowSolver::CubeFlow_UniformGrid< dim, nstate >
 CubeFlow_UniformGrid (const Parameters::AllParameters *const parameters_input)
 
double get_adaptive_time_step (std::shared_ptr< DGBase< dim, double >> dg) const override
 Function to compute the adaptive time step.
 
double get_adaptive_time_step_initial (std::shared_ptr< DGBase< dim, double >> dg) override
 Function to compute the initial adaptive time step.
 
void update_maximum_local_wave_speed (DGBase< dim, double > &dg)
 Updates the maximum local wave speed.
 
- Public Member Functions inherited from PHiLiP::FlowSolver::FlowSolverCaseBase< dim, nstate >
 FlowSolverCaseBase (const Parameters::AllParameters *const parameters_input)
 Constructor.
 
virtual ~FlowSolverCaseBase ()=default
 Destructor.
 
void display_flow_solver_setup (std::shared_ptr< DGBase< dim, double >> dg) const
 Displays the flow setup parameters.
 
virtual void set_higher_order_grid (std::shared_ptr< DGBase< dim, double >> dg) const
 Set higher order grid.
 
virtual void compute_unsteady_data_and_write_to_table (const unsigned int current_iteration, const double current_time, const std::shared_ptr< DGBase< dim, double >> dg, const std::shared_ptr< dealii::TableHandler > unsteady_data_table)
 Virtual function to write unsteady snapshot data to table.
 
virtual void steady_state_postprocessing (std::shared_ptr< DGBase< dim, double >> dg) const
 Virtual function for postprocessing when solving for steady state.
 
void set_time_step (const double time_step_input)
 Setter for time step.
 

Protected Types

enum  IntegratedQuantitiesEnum {
  kinetic_energy, enstrophy, pressure_dilatation, deviatoric_strain_rate_tensor_magnitude_sqr,
  strain_rate_tensor_magnitude_sqr
}
 List of possible integrated quantities over the domain.
 

Protected Member Functions

void display_additional_flow_case_specific_parameters () const override
 Display additional more specific flow case parameters.
 
double get_constant_time_step (std::shared_ptr< DGBase< dim, double >> dg) const override
 Function to compute the constant time step.
 
void compute_unsteady_data_and_write_to_table (const std::shared_ptr< ODE::ODESolverBase< dim, double >> ode_solver, const std::shared_ptr< DGBase< dim, double >> dg, const std::shared_ptr< dealii::TableHandler > unsteady_data_table) override
 Compute the desired unsteady data and write it to a table.
 
- Protected Member Functions inherited from PHiLiP::FlowSolver::PeriodicCubeFlow< dim, nstate >
void display_grid_parameters () const
 Display grid parameters.
 
- Protected Member Functions inherited from PHiLiP::FlowSolver::FlowSolverCaseBase< dim, nstate >
void add_value_to_data_table (const double value, const std::string value_string, const std::shared_ptr< dealii::TableHandler > data_table) const
 Add a value to a given data table with scientific format.
 
double get_time_step () const
 Getter for time step.
 

Protected Attributes

const std::string unsteady_data_table_filename_with_extension
 Filename (with extension) for the unsteady data table.
 
const unsigned int number_of_times_to_output_velocity_field
 Number of times to output the velocity field.
 
const bool output_velocity_field_at_fixed_times
 Flag for outputting velocity field at fixed times.
 
const bool output_vorticity_magnitude_field_in_addition_to_velocity
 Flag for outputting vorticity magnitude field in addition to velocity field at fixed times.
 
const std::string output_flow_field_files_directory_name
 Directory for writting flow field files.
 
const bool output_solution_at_exact_fixed_times
 Flag for outputting the solution at exact fixed times by decreasing the time step on the fly.
 
std::shared_ptr< Physics::NavierStokes< dim, dim+2, double > > navier_stokes_physics
 Pointer to Navier-Stokes physics object for computing things on the fly.
 
bool is_taylor_green_vortex = false
 Identifies if taylor green vortex case; initialized as false.
 
bool is_decaying_homogeneous_isotropic_turbulence = false
 Identified if DHIT case; initialized as false.
 
bool is_viscous_flow = true
 Identifies if viscous flow; initialized as true.
 
bool do_calculate_numerical_entropy = false
 Identifies if numerical entropy should be calculated; initialized as false.
 
std::array< double, NUMBER_OF_INTEGRATED_QUANTITIESintegrated_quantities
 Array for storing the integrated quantities; done for computational efficiency.
 
double maximum_local_wave_speed
 Maximum local wave speed (i.e. convective eigenvalue)
 
double previous_numerical_entropy = 0
 Numerical entropy at previous timestep.
 
double cumulative_numerical_entropy_change_FRcorrected = 0
 Cumulative change in numerical entropy.
 
double initial_numerical_entropy_abs = 0
 Numerical entropy at initial time.
 
dealii::Table< 1, double > output_velocity_field_times
 Times at which to output the velocity field.
 
unsigned int index_of_current_desired_time_to_output_velocity_field
 Index of current desired time to output velocity field.
 
std::string flow_field_quantity_filename_prefix
 Flow field quantity filename prefix.
 
std::shared_ptr< dealii::TableHandler > exact_output_times_of_velocity_field_files_table
 Data table storing the exact output times for the velocity field files.
 
- Protected Attributes inherited from PHiLiP::FlowSolver::PeriodicCubeFlow< dim, nstate >
const int number_of_cells_per_direction
 Number of cells per direction for the grid.
 
const double domain_left
 Domain left-boundary value for generating the grid.
 
const double domain_right
 Domain right-boundary value for generating the grid.
 
const double domain_size
 Domain size (length in 1D, area in 2D, and volume in 3D)
 
- Protected Attributes inherited from PHiLiP::FlowSolver::CubeFlow_UniformGrid< dim, nstate >
double maximum_local_wave_speed
 Maximum local wave speed (i.e. convective eigenvalue)
 
std::shared_ptr< Physics::PhysicsBase< dim, nstate, double > > pde_physics
 Pointer to Physics object for computing things on the fly.
 
- Protected Attributes inherited from PHiLiP::FlowSolver::FlowSolverCaseBase< dim, nstate >
const Parameters::AllParameters all_param
 All parameters.
 
const MPI_Comm mpi_communicator
 MPI communicator.
 
const int mpi_rank
 MPI rank.
 
const int n_mpi
 Number of MPI processes.
 
dealii::ConditionalOStream pcout
 ConditionalOStream. More...
 

Static Private Attributes

static const int NUMBER_OF_INTEGRATED_QUANTITIES = 5
 

Additional Inherited Members

- Public Attributes inherited from PHiLiP::FlowSolver::FlowSolverCaseBase< dim, nstate >
std::shared_ptr< InitialConditionFunction< dim, nstate, double > > initial_condition_function
 Initial condition function.
 

Detailed Description

template<int dim, int nstate>
class PHiLiP::FlowSolver::PeriodicTurbulence< dim, nstate >

Definition at line 14 of file periodic_turbulence.h.

Constructor & Destructor Documentation

◆ PeriodicTurbulence()

template<int dim, int nstate>
PHiLiP::FlowSolver::PeriodicTurbulence< dim, nstate >::PeriodicTurbulence ( const Parameters::AllParameters *const  parameters_input)
explicit

Constructor.

For outputting velocity field

Definition at line 25 of file periodic_turbulence.cpp.

Member Function Documentation

◆ compute_and_update_integrated_quantities()

template<int dim, int nstate>
void PHiLiP::FlowSolver::PeriodicTurbulence< dim, nstate >::compute_and_update_integrated_quantities ( DGBase< dim, double > &  dg)

Computes the integrated quantities over the domain simultaneously and updates the array storing them Note: For efficiency, this also simultaneously updates the local maximum wave speed

Definition at line 324 of file periodic_turbulence.cpp.

◆ get_deviatoric_strain_rate_tensor_based_dissipation_rate()

template<int dim, int nstate>
double PHiLiP::FlowSolver::PeriodicTurbulence< dim, nstate >::get_deviatoric_strain_rate_tensor_based_dissipation_rate ( ) const

Gets non-dimensional theoretical deviatoric strain-rate tensor based dissipation rate – Reference: Cox, Christopher, et al. "Accuracy, stability, and performance comparison between the spectral difference and flux reconstruction schemes." Computers & Fluids 221 (2021): 104922.

Definition at line 523 of file periodic_turbulence.cpp.

◆ get_integrated_enstrophy()

template<int dim, int nstate>
double PHiLiP::FlowSolver::PeriodicTurbulence< dim, nstate >::get_integrated_enstrophy ( ) const

Gets the nondimensional integrated enstrophy given a DG object from dg->solution – Reference: Cox, Christopher, et al. "Accuracy, stability, and performance comparison between the spectral difference and flux reconstruction schemes." Computers & Fluids 221 (2021): 104922.

Definition at line 499 of file periodic_turbulence.cpp.

◆ get_integrated_kinetic_energy()

template<int dim, int nstate>
double PHiLiP::FlowSolver::PeriodicTurbulence< dim, nstate >::get_integrated_kinetic_energy ( ) const

Gets the nondimensional integrated kinetic energy given a DG object from dg->solution – Reference: Cox, Christopher, et al. "Accuracy, stability, and performance comparison between the spectral difference and flux reconstruction schemes." Computers & Fluids 221 (2021): 104922.

Definition at line 486 of file periodic_turbulence.cpp.

◆ get_pressure_dilatation_based_dissipation_rate()

template<int dim, int nstate>
double PHiLiP::FlowSolver::PeriodicTurbulence< dim, nstate >::get_pressure_dilatation_based_dissipation_rate ( ) const

Evaluate non-dimensional theoretical pressure-dilatation dissipation rate – Reference: Cox, Christopher, et al. "Accuracy, stability, and performance comparison between the spectral difference and flux reconstruction schemes." Computers & Fluids 221 (2021): 104922.

Definition at line 516 of file periodic_turbulence.cpp.

◆ get_strain_rate_tensor_based_dissipation_rate()

template<int dim, int nstate>
double PHiLiP::FlowSolver::PeriodicTurbulence< dim, nstate >::get_strain_rate_tensor_based_dissipation_rate ( ) const

Gets non-dimensional theoretical strain-rate tensor based dissipation rate from integrated strain-rate tensor magnitude squared. – Reference: Navah, Farshad, et al. "A High-Order Variational Multiscale Approach to Turbulence for Compact Nodal Schemes."

Definition at line 535 of file periodic_turbulence.cpp.

◆ get_vorticity_based_dissipation_rate()

template<int dim, int nstate>
double PHiLiP::FlowSolver::PeriodicTurbulence< dim, nstate >::get_vorticity_based_dissipation_rate ( ) const

Gets non-dimensional theoretical vorticity tensor based dissipation rate Note: For incompressible flows or when dilatation effects are negligible – Reference: Cox, Christopher, et al. "Accuracy, stability, and performance comparison between the spectral difference and flux reconstruction schemes." Computers & Fluids 221 (2021): 104922.

Definition at line 505 of file periodic_turbulence.cpp.

Member Data Documentation

◆ NUMBER_OF_INTEGRATED_QUANTITIES

template<int dim, int nstate>
const int PHiLiP::FlowSolver::PeriodicTurbulence< dim, nstate >::NUMBER_OF_INTEGRATED_QUANTITIES = 5
staticprivate

Number of different computed quantities Corresponds to the number of items in IntegratedQuantitiesEnum

Definition at line 19 of file periodic_turbulence.h.


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