[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::MHD< dim, nstate, real > Class Template Reference

Magnetohydrodynamics (MHD) equations. Derived from PhysicsBase. More...

#include <mhd.h>

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

Public Member Functions

 MHD (const Parameters::AllParameters *const parameters_input, const double gamma_gas, const dealii::Tensor< 2, 3, double > input_diffusion_tensor=Parameters::ManufacturedSolutionParam::get_default_diffusion_tensor(), std::shared_ptr< ManufacturedSolutionFunction< dim, real > > manufactured_solution_function=nullptr, const bool has_nonzero_diffusion=false, const bool has_nonzero_physical_source=false)
 Constructor.
 
std::array< dealii::Tensor< 1, dim, real >, nstate > convective_flux (const std::array< real, nstate > &conservative_soln) const
 Convective flux: \( \mathbf{F}_{conv} \).
 
std::array< real, nstate > convective_normal_flux (const std::array< real, nstate > &conservative_soln, const dealii::Tensor< 1, dim, real > &normal) const
 Convective flux: \( \mathbf{F}_{conv} \hat{n} \).
 
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: \( \frac{\partial \mathbf{F}_{conv}}{\partial w} \cdot \mathbf{n} \).
 
std::array< real, nstate > convective_eigenvalues (const std::array< real, nstate > &, const dealii::Tensor< 1, dim, real > &) const
 Spectral radius of convective term Jacobian is 'c'.
 
real max_convective_eigenvalue (const std::array< real, nstate > &soln) const
 Maximum convective eigenvalue.
 
real max_viscous_eigenvalue (const std::array< real, nstate > &soln) const
 Maximum viscous eigenvalue.
 
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.
 
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
 (function overload) Dissipative flux: 0
 
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.
 
std::array< real, nstate > source_term (const dealii::Point< dim, real > &pos, const std::array< real, nstate > &conservative_soln, const real current_time) const
 (function overload) Source term is zero or depends on manufactured solution
 
std::array< real, nstate > convert_conservative_to_primitive (const std::array< real, nstate > &conservative_soln) const
 
std::array< real, nstate > convert_primitive_to_conservative (const std::array< real, nstate > &primitive_soln) const
 
real compute_pressure (const std::array< real, nstate > &conservative_soln) const
 Evaluate pressure from conservative variables.
 
real compute_magnetic_energy (const std::array< real, nstate > &conservative_soln) const
 Evaluate Magnetic Energy.
 
real compute_pressure_from_enthalpy (const std::array< real, 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.
 
real compute_sound (const std::array< real, nstate > &conservative_soln) const
 Evaluate speed of sound from conservative variables.
 
real compute_sound (const real density, const real pressure) const
 Evaluate speed of sound from density and pressure.
 
dealii::Tensor< 1, dim, real > compute_velocities (const std::array< real, nstate > &conservative_soln) const
 Evaluate velocities from conservative variables.
 
real compute_velocity_squared (const dealii::Tensor< 1, dim, real > &velocities) const
 Given the velocity vector \( \mathbf{u} \), returns the dot-product \( \mathbf{u} \cdot \mathbf{u} \).
 
dealii::Tensor< 1, dim, real > extract_velocities_from_primitive (const std::array< real, nstate > &primitive_soln) const
 Given primitive variables, returns velocities.
 
real compute_total_energy (const std::array< real, nstate > &primitive_soln) const
 Given primitive variables, returns total energy.
 
real compute_entropy_measure (const std::array< real, nstate > &conservative_soln) const
 Evaluate entropy from conservative variables. More...
 
real compute_mach_number (const std::array< real, nstate > &conservative_soln) const
 Given conservative variables, returns Mach number.
 
real compute_dimensional_temperature (const std::array< real, nstate > &primitive_soln) const
 Given primitive variables, returns DIMENSIONALIZED temperature using the equation of state.
 
real compute_temperature (const std::array< real, nstate > &primitive_soln) const
 Given primitive variables, returns NON-DIMENSIONALIZED temperature using free-stream non-dimensionalization. More...
 
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-dimensionalization. More...
 
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-dimensionalization. More...
 
std::array< real, nstate > compute_entropy_variables (const std::array< real, nstate > &conservative_soln) const
 Computes the entropy variables.
 
std::array< real, nstate > compute_conservative_variables_from_entropy_variables (const std::array< real, nstate > &entropy_var) const
 Computes the conservative variables from the entropy 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. More...
 
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. More...
 
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. More...
 
real compute_mean_specific_energy (const std::array< real, nstate > &conservative_soln1, const std::array< real, nstate > &convervative_soln2) const
 Mean specific energy given two sets of conservative solutions. More...
 
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
 Evaluates boundary values and gradients on the other side of the face.
 
- Public Member Functions inherited from PHiLiP::Physics::PhysicsBase< dim, nstate, real >
 PhysicsBase (const Parameters::AllParameters *const parameters_input, const bool has_nonzero_diffusion_input, const bool has_nonzero_physical_source_input, const dealii::Tensor< 2, 3, double > input_diffusion_tensor=Parameters::ManufacturedSolutionParam::get_default_diffusion_tensor(), std::shared_ptr< ManufacturedSolutionFunction< dim, real > > manufactured_solution_function_input=nullptr)
 Default constructor that will set the constants.
 
 PhysicsBase (const Parameters::AllParameters *const parameters_input, const bool has_nonzero_diffusion_input, const bool has_nonzero_physical_source_input, std::shared_ptr< ManufacturedSolutionFunction< dim, real > > manufactured_solution_function_input=nullptr)
 Constructor that will call default constructor.
 
virtual ~PhysicsBase ()=default
 Virtual destructor required for abstract classes.
 
virtual 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
 Convective Numerical Split Flux for split form.
 
virtual 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)
 
virtual std::array< real, nstate > physical_source_term (const dealii::Point< dim, real > &pos, const std::array< real, nstate > &solution, const std::array< dealii::Tensor< 1, dim, real >, nstate > &solution_gradient, const dealii::types::global_dof_index cell_index) const
 Physical source term that does require differentiation.
 
virtual std::array< real, nstate > artificial_source_term (const real viscosity_coefficient, const dealii::Point< dim, real > &pos, const std::array< real, nstate > &solution) const
 Artificial source term that does not require differentiation stemming from artificial dissipation.
 
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 dealii::Vector< double > post_compute_derived_quantities_scalar (const double &uh, const dealii::Tensor< 1, dim > &, const dealii::Tensor< 2, dim > &, const dealii::Tensor< 1, dim > &, const dealii::Point< dim > &) const
 Returns current scalar 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...
 
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 dealii::UpdateFlags post_get_needed_update_flags () const
 Returns required update flags of the solution to be used by PhysicsPostprocessor to output current solution. More...
 
template<typename real2 >
real2 handle_non_physical_result (const std::string message="") const
 Function to handle nonphysical results. More...
 

Public Attributes

const double gam
 Constant heat capacity ratio of air.
 
const double gamm1
 Gamma-1.0 used often.
 
- Public Attributes inherited from PHiLiP::Physics::PhysicsBase< dim, nstate, real >
const bool has_nonzero_diffusion
 Flag to signal that diffusion term is non-zero.
 
const bool has_nonzero_physical_source
 Flag to signal that physical source term is non-zero.
 
const Parameters::AllParameters *const all_parameters
 Pointer to parameters object.
 
const NonPhysicalBehaviorEnum non_physical_behavior_type
 Determines type of nonphysical behavior.
 
std::shared_ptr< ManufacturedSolutionFunction< dim, real > > manufactured_solution_function
 Manufactured solution function.
 
const double BIG_NUMBER = 1e100
 BIG_NUMBER which is returned in place of NaN according to handle_non_physical_result() More...
 

Additional Inherited Members

- Public Types inherited from PHiLiP::Physics::PhysicsBase< dim, nstate, real >
using NonPhysicalBehaviorEnum = Parameters::AllParameters::NonPhysicalBehaviorEnum
 
- Protected Attributes inherited from PHiLiP::Physics::PhysicsBase< dim, nstate, real >
dealii::ConditionalOStream pcout
 ConditionalOStream. More...
 
dealii::Tensor< 2, dim, double > diffusion_tensor
 Anisotropic diffusion matrix. More...
 

Detailed Description

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

Magnetohydrodynamics (MHD) equations. Derived from PhysicsBase.

Only 2D and 3D State variable and convective fluxes given by

\[ \mathbf{w} = \begin{bmatrix} \rho \\ \rho v_1 \\ \rho v_2 \\ \rho v_3 \\ \rho E \end{bmatrix} , \qquad \mathbf{F}_{conv} = \begin{bmatrix} \mathbf{f}^x_{conv}, \mathbf{f}^y_{conv}, \mathbf{f}^z_{conv} \end{bmatrix} = \begin{bmatrix} \begin{bmatrix} \rho v_1 \\ \rho v_1 v_1 + p \\ \rho v_1 v_2 \\ \rho v_1 v_3 \\ v_1 (\rho e+p) \end{bmatrix} , \begin{bmatrix} \rho v_2 \\ \rho v_1 v_2 \\ \rho v_2 v_2 + p \\ \rho v_2 v_3 \\ v_2 (\rho e+p) \end{bmatrix} , \begin{bmatrix} \rho v_3 \\ \rho v_1 v_3 \\ \rho v_2 v_3 \\ \rho v_3 v_3 + p \\ v_3 (\rho e+p) \end{bmatrix} \end{bmatrix} \]

where, \( E \) is the specific total energy and \( e \) is the specific internal energy, related by

\[ E = e + |V|^2 / 2 \]

For a calorically perfect gas

\[ p=(\gamma -1)(\rho e-\frac{1}{2}\rho \|\mathbf{v}\|) \]

Dissipative flux \( \mathbf{F}_{diss} = \mathbf{0} \)

Source term \( s(\mathbf{x}) \)

Equation:

\[ \boldsymbol{\nabla} \cdot ( \mathbf{F}_{conv}( w ) + \mathbf{F}_{diss}( w, \boldsymbol{\nabla}(w) ) = s(\mathbf{x}) \]

Still need to provide functions to un-non-dimensionalize the variables. Like, given density_inf

Definition at line 77 of file mhd.h.

Member Function Documentation

◆ compute_density_from_pressure_temperature()

template<int dim, int nstate, typename real >
real PHiLiP::Physics::MHD< dim, nstate, real >::compute_density_from_pressure_temperature ( const real  pressure,
const real  temperature 
) const
inline

Given pressure and temperature, returns NON-DIMENSIONALIZED density using free-stream non-dimensionalization.

See the book I do like CFD, sec 4.14.2

mach_inf_sqr

Definition at line 164 of file mhd.cpp.

◆ compute_entropy_measure()

template<int dim, int nstate, typename real >
real PHiLiP::Physics::MHD< dim, nstate, real >::compute_entropy_measure ( const std::array< real, nstate > &  conservative_soln) const
inline

Evaluate entropy from conservative variables.

Note that it is not the actual entropy since it's missing some constants. Used to check entropy convergence See discussion in https://physics.stackexchange.com/questions/116779/entropy-is-constant-how-to-express-this-equation-in-terms-of-pressure-and-densi?answertab=votes#tab-top

Definition at line 123 of file mhd.cpp.

◆ compute_mean_density()

template<int dim, int nstate, typename real >
real PHiLiP::Physics::MHD< dim, nstate, real >::compute_mean_density ( const std::array< real, nstate > &  conservative_soln1,
const std::array< real, nstate > &  convervative_soln2 
) const
inline

Mean density given two sets of conservative solutions.

Used in the implementation of the split form.

Definition at line 261 of file mhd.cpp.

◆ compute_mean_pressure()

template<int dim, int nstate, typename real >
real PHiLiP::Physics::MHD< dim, nstate, real >::compute_mean_pressure ( const std::array< real, nstate > &  conservative_soln1,
const std::array< real, nstate > &  convervative_soln2 
) const
inline

Mean pressure given two sets of conservative solutions.

Used in the implementation of the split form.

Definition at line 269 of file mhd.cpp.

◆ compute_mean_specific_energy()

template<int dim, int nstate, typename real >
real PHiLiP::Physics::MHD< dim, nstate, real >::compute_mean_specific_energy ( const std::array< real, nstate > &  conservative_soln1,
const std::array< real, nstate > &  convervative_soln2 
) const
inline

Mean specific energy given two sets of conservative solutions.

Used in the implementation of the split form.

Definition at line 294 of file mhd.cpp.

◆ compute_mean_velocities()

template<int dim, int nstate, typename real >
dealii::Tensor< 1, dim, real > PHiLiP::Physics::MHD< dim, nstate, real >::compute_mean_velocities ( const std::array< real, nstate > &  conservative_soln1,
const std::array< real, nstate > &  convervative_soln2 
) const
inline

Mean velocities given two sets of conservative solutions.

Used in the implementation of the split form.

Definition at line 279 of file mhd.cpp.

◆ compute_temperature()

template<int dim, int nstate, typename real >
real PHiLiP::Physics::MHD< dim, nstate, real >::compute_temperature ( const std::array< real, nstate > &  primitive_soln) const
inline

Given primitive variables, returns NON-DIMENSIONALIZED temperature using free-stream non-dimensionalization.

See the book I do like CFD, sec 4.14.2

mach_inf_sqr

Definition at line 155 of file mhd.cpp.

◆ compute_temperature_from_density_pressure()

template<int dim, int nstate, typename real >
real PHiLiP::Physics::MHD< dim, nstate, real >::compute_temperature_from_density_pressure ( const real  density,
const real  pressure 
) const
inline

Given density and pressure, returns NON-DIMENSIONALIZED temperature using free-stream non-dimensionalization.

See the book I do like CFD, sec 4.14.2

Definition at line 171 of file mhd.cpp.

◆ convert_conservative_to_primitive()

template<int dim, int nstate, typename real >
std::array< real, nstate > PHiLiP::Physics::MHD< dim, nstate, real >::convert_conservative_to_primitive ( const std::array< real, nstate > &  conservative_soln) const
inline

Given conservative variables [density, [momentum], total energy], returns primitive variables [density, [velocities], pressure].

Opposite of convert_primitive_to_conservative

Definition at line 42 of file mhd.cpp.

◆ convert_primitive_to_conservative()

template<int dim, int nstate, typename real >
std::array< real, nstate > PHiLiP::Physics::MHD< dim, nstate, real >::convert_primitive_to_conservative ( const std::array< real, nstate > &  primitive_soln) const
inline

Given primitive variables [density, [velocities], pressure], returns conservative variables [density, [momentum], total energy].

Opposite of convert_primitive_to_conservative

Definition at line 61 of file mhd.cpp.


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