1 #ifndef __CONVECTIVE_NUMERICAL_FLUX__ 2 #define __CONVECTIVE_NUMERICAL_FLUX__ 4 #include <deal.II/base/tensor.h> 5 #include "physics/physics.h" 6 #include "physics/euler.h" 9 namespace NumericalFlux {
11 using AllParam = Parameters::AllParameters;
14 template<
int dim,
int nstate,
typename real>
23 const std::array<real, nstate> &soln_int,
24 const std::array<real, nstate> &soln_ext,
25 const dealii::Tensor<1,dim,real> &normal1)
const = 0;
29 template<
int dim,
int nstate,
typename real>
35 : pde_physics(physics_input) {};
39 const std::array<real, nstate> &soln_int,
40 const std::array<real, nstate> &soln_ext,
41 const dealii::Tensor<1,dim,real> &normal1)
const;
45 const std::shared_ptr < Physics::PhysicsBase<dim, nstate, real> >
pde_physics;
49 template<
int dim,
int nstate,
typename real>
55 : pde_physics(physics_input) {};
59 const std::array<real, nstate> &soln_int,
60 const std::array<real, nstate> &soln_ext,
61 const dealii::Tensor<1,dim,real> &normal1)
const;
65 const std::shared_ptr < Physics::PhysicsBase<dim, nstate, real> >
pde_physics;
69 template<
int dim,
int nstate,
typename real>
77 virtual std::array<real, nstate> evaluate_riemann_solver_dissipation (
78 const std::array<real, nstate> &soln_int,
79 const std::array<real, nstate> &soln_ext,
80 const dealii::Tensor<1,dim,real> &normal1)
const = 0;
84 template<
int dim,
int nstate,
typename real>
89 std::array<real, nstate> evaluate_riemann_solver_dissipation (
90 const std::array<real, nstate> &soln_int,
91 const std::array<real, nstate> &soln_ext,
92 const dealii::Tensor<1,dim,real> &normal1)
const;
96 template<
int dim,
int nstate,
typename real>
102 : pde_physics(physics_input) {};
110 std::array<real, nstate> evaluate_riemann_solver_dissipation (
111 const std::array<real, nstate> &soln_int,
112 const std::array<real, nstate> &soln_ext,
113 const dealii::Tensor<1,dim,real> &normal1)
const;
117 const std::shared_ptr < Physics::PhysicsBase<dim, nstate, real> >
pde_physics;
121 template<
int dim,
int nstate,
typename real>
131 : euler_physics(
std::dynamic_pointer_cast<Physics::Euler<dim,nstate,real>>(physics_input)) {};
134 virtual void evaluate_entropy_fix (
135 const std::array<real, 3> &eig_L,
136 const std::array<real, 3> &eig_R,
137 std::array<real, 3> &eig_RoeAvg,
138 const real vel2_ravg,
139 const real sound_ravg)
const = 0;
142 virtual void evaluate_additional_modifications (
143 const std::array<real, nstate> &soln_int,
144 const std::array<real, nstate> &soln_ext,
145 const std::array<real, 3> &eig_L,
146 const std::array<real, 3> &eig_R,
148 dealii::Tensor<1,dim,real> &dV_tangent)
const = 0;
155 std::array<real, nstate> evaluate_riemann_solver_dissipation (
156 const std::array<real, nstate> &soln_int,
157 const std::array<real, nstate> &soln_ext,
158 const dealii::Tensor<1,dim,real> &normal1)
const;
162 template<
int dim,
int nstate,
typename real>
172 void evaluate_entropy_fix(
173 const std::array<real, 3> &eig_L,
174 const std::array<real, 3> &eig_R,
175 std::array<real, 3> &eig_RoeAvg,
176 const real vel2_ravg,
177 const real sound_ravg)
const;
180 void evaluate_additional_modifications(
181 const std::array<real, nstate> &soln_int,
182 const std::array<real, nstate> &soln_ext,
183 const std::array<real, 3> &eig_L,
184 const std::array<real, 3> &eig_R,
186 dealii::Tensor<1,dim,real> &dV_tangent)
const;
191 template<
int dim,
int nstate,
typename real>
202 void evaluate_entropy_fix(
203 const std::array<real, 3> &eig_L,
204 const std::array<real, 3> &eig_R,
205 std::array<real, 3> &eig_RoeAvg,
206 const real vel2_ravg,
207 const real sound_ravg)
const;
211 void evaluate_additional_modifications(
212 const std::array<real, nstate> &soln_int,
213 const std::array<real, nstate> &soln_ext,
214 const std::array<real, 3> &eig_L,
215 const std::array<real, 3> &eig_R,
217 dealii::Tensor<1,dim,real> &dV_tangent)
const;
222 void evaluate_shock_indicator(
223 const std::array<real, 3> &eig_L,
224 const std::array<real, 3> &eig_R,
226 int &ssw_RIGHT)
const;
230 template<
int dim,
int nstate,
typename real>
241 std::unique_ptr< BaselineNumericalFluxConvective<dim,nstate,real> >
baseline;
249 const std::array<real, nstate> &soln_int,
250 const std::array<real, nstate> &soln_ext,
251 const dealii::Tensor<1,dim,real> &normal1)
const;
255 template<
int dim,
int nstate,
typename real>
265 template<
int dim,
int nstate,
typename real>
274 template<
int dim,
int nstate,
typename real>
283 template<
int dim,
int nstate,
typename real>
292 template<
int dim,
int nstate,
typename real>
301 template<
int dim,
int nstate,
typename real>
310 template<
int dim,
int nstate,
typename real>
319 template<
int dim,
int nstate,
typename real>
std::unique_ptr< BaselineNumericalFluxConvective< dim, nstate, real > > baseline
Baseline convective numerical flux object.
const std::shared_ptr< Physics::PhysicsBase< dim, nstate, real > > pde_physics
Numerical flux requires physics to evaluate split form convective flux.
std::unique_ptr< RiemannSolverDissipation< dim, nstate, real > > riemann_solver_dissipation
Upwind convective numerical flux object.
virtual ~BaselineNumericalFluxConvective()=default
Base class destructor required for abstract classes.
Central numerical flux. Derived from BaselineNumericalFluxConvective.
Base class from which Advection, Diffusion, ConvectionDiffusion, and Euler is derived.
Entropy conserving numerical flux with Roe dissipation. Derived from NumericalFluxConvective.
Files for the baseline physics.
EntropyConservingBaselineNumericalFluxConvective(std::shared_ptr< Physics::PhysicsBase< dim, nstate, real >> physics_input)
Constructor.
CentralBaselineNumericalFluxConvective(std::shared_ptr< Physics::PhysicsBase< dim, nstate, real >> physics_input)
Constructor.
Roe-Pike numerical flux. Derived from NumericalFluxConvective.
Entropy conserving numerical flux with L2Roe dissipation. Derived from NumericalFluxConvective.
Base class of Roe (Roe-Pike) flux with entropy fix. Derived from RiemannSolverDissipation.
L2RoeRiemannSolverDissipation(std::shared_ptr< Physics::PhysicsBase< dim, nstate, real >> physics_input)
Constructor.
Base class of Riemann solver dissipation (i.e. upwind-term) for numerical flux associated with convec...
Entropy conserving numerical flux. Derived from NumericalFluxConvective.
virtual std::array< real, nstate > evaluate_flux(const std::array< real, nstate > &soln_int, const std::array< real, nstate > &soln_ext, const dealii::Tensor< 1, dim, real > &normal1) const =0
Returns the convective numerical flux at an interface.
Base class of baseline numerical flux (without upwind term) associated with convection.
Base class of numerical flux associated with convection.
L2Roe numerical flux. Derived from NumericalFluxConvective.
const std::shared_ptr< Physics::PhysicsBase< dim, nstate, real > > pde_physics
Numerical flux requires physics to evaluate convective flux.
Lax-Friedrichs Riemann solver dissipation. Derived from RiemannSolverDissipation. ...
RoePike flux with entropy fix. Derived from RoeBase.
Entropy conserving numerical flux with Lax Friedrichs dissipation. Derived from NumericalFluxConvecti...
Zero Riemann solver dissipation. Derived from RiemannSolverDissipation.
LaxFriedrichsRiemannSolverDissipation(std::shared_ptr< Physics::PhysicsBase< dim, nstate, real >> physics_input)
Constructor.
Central numerical flux. Derived from NumericalFluxConvective.
const std::shared_ptr< Physics::PhysicsBase< dim, nstate, real > > pde_physics
Numerical flux requires physics to evaluate convective eigenvalues.
const std::shared_ptr< Physics::Euler< dim, nstate, real > > euler_physics
Numerical flux requires physics to evaluate convective eigenvalues.
Lax-Friedrichs numerical flux. Derived from NumericalFluxConvective.
RoeBaseRiemannSolverDissipation(std::shared_ptr< Physics::PhysicsBase< dim, nstate, real >> physics_input)
Constructor.
RoePikeRiemannSolverDissipation(std::shared_ptr< Physics::PhysicsBase< dim, nstate, real >> physics_input)
Constructor.
Entropy Conserving Numerical Flux. Derived from BaselineNumericalFluxConvective.