[P]arallel [Hi]gh-order [Li]brary for [P]DEs  Latest
Parallel High-Order Library for PDEs through hp-adaptive Discontinuous Galerkin methods
convection_diffusion.h
1 #ifndef __CONVECTION_DIFFUSION__
2 #define __CONVECTION_DIFFUSION__
3 
4 #include <deal.II/base/tensor.h>
5 
6 #include "parameters/all_parameters.h"
7 #include "physics.h"
8 #include "parameters/parameters_manufactured_solution.h"
9 
10 namespace PHiLiP {
11 namespace Physics {
13 
28 template <int dim, int nstate, typename real>
29 class ConvectionDiffusion : public PhysicsBase <dim, nstate, real>
30 {
31 protected:
32  // For overloading the virtual functions defined in PhysicsBase
41 protected:
43  double linear_advection_velocity[3] = { 1.1, -atan(1)*4.0 / exp(1), exp(1)/(atan(1)*4.0) };
45  double diffusion_scaling_coeff = 0.1*atan(1)*4.0/exp(1);
46 public:
47  const bool hasConvection;
48 
49  const bool hasDiffusion;
52 
55  const Parameters::AllParameters *const parameters_input,
56  const bool convection = true,
57  const bool diffusion = true,
58  const dealii::Tensor<2,3,double> input_diffusion_tensor = Parameters::ManufacturedSolutionParam::get_default_diffusion_tensor(),
59  const dealii::Tensor<1,3,double> input_advection_vector = Parameters::ManufacturedSolutionParam::get_default_advection_vector(),
60  const double input_diffusion_coefficient = Parameters::ManufacturedSolutionParam::get_default_diffusion_coefficient(),
62  const Parameters::AllParameters::TestType parameters_test = Parameters::AllParameters::TestType::run_control,
63  const bool has_nonzero_physical_source = false) :
64  PhysicsBase<dim,nstate,real>(parameters_input, diffusion, has_nonzero_physical_source, input_diffusion_tensor, manufactured_solution_function),
65  linear_advection_velocity{input_advection_vector[0], input_advection_vector[1], input_advection_vector[2]},
66  diffusion_scaling_coeff(input_diffusion_coefficient),
67  hasConvection(convection),
68  hasDiffusion(diffusion),
69  test_type(parameters_test)
70  {};
71 
73  std::array<dealii::Tensor<1,dim,real>,nstate> convective_flux (const std::array<real,nstate> &solution) const;
74 
76  std::array<dealii::Tensor<1,dim,real>,nstate> convective_numerical_split_flux (
77  const std::array<real,nstate> &soln1,
78  const std::array<real,nstate> &soln2) const override;
79 
81  std::array<real,nstate> compute_entropy_variables (
82  const std::array<real,nstate> &conservative_soln) const;
83 
86  const std::array<real,nstate> &entropy_var) const;
87 
89  std::array<real,nstate> convective_eigenvalues (
90  const std::array<real,nstate> &/*solution*/,
91  const dealii::Tensor<1,dim,real> &/*normal*/) const;
92 
94  real max_convective_eigenvalue (const std::array<real,nstate> &soln) const;
95 
97  real max_viscous_eigenvalue (const std::array<real,nstate> &soln) const;
98 
99  // /// Diffusion matrix is identity
100  // std::array<dealii::Tensor<1,dim,real>,nstate> apply_diffusion_matrix (
101  // const std::array<real,nstate> &solution,
102  // const std::array<dealii::Tensor<1,dim,real>,nstate> &solution_grad) const;
103 
105  std::array<dealii::Tensor<1,dim,real>,nstate> dissipative_flux (
106  const std::array<real,nstate> &solution,
107  const std::array<dealii::Tensor<1,dim,real>,nstate> &solution_gradient,
108  const dealii::types::global_dof_index cell_index) const;
109 
111  std::array<dealii::Tensor<1,dim,real>,nstate> dissipative_flux (
112  const std::array<real,nstate> &solution,
113  const std::array<dealii::Tensor<1,dim,real>,nstate> &solution_gradient) const;
114 
116  std::array<real,nstate> source_term (
117  const dealii::Point<dim,real> &pos,
118  const std::array<real,nstate> &solution,
119  const real current_time,
120  const dealii::types::global_dof_index cell_index) const;
121 
123  std::array<real,nstate> source_term (
124  const dealii::Point<dim,real> &pos,
125  const std::array<real,nstate> &solution,
126  const real current_time) const;
127 
129 
132  void boundary_face_values (
133  const int /*boundary_type*/,
134  const dealii::Point<dim, real> &/*pos*/,
135  const dealii::Tensor<1,dim,real> &/*normal*/,
136  const std::array<real,nstate> &/*soln_int*/,
137  const std::array<dealii::Tensor<1,dim,real>,nstate> &/*soln_grad_int*/,
138  std::array<real,nstate> &/*soln_bc*/,
139  std::array<dealii::Tensor<1,dim,real>,nstate> &/*soln_grad_bc*/) const;
140 
141 protected:
143  dealii::Tensor<1,dim,real> advection_speed () const;
145  real diffusion_coefficient () const;
146 };
147 
148 
149 } // Physics namespace
150 } // PHiLiP namespace
151 
152 #endif
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 &#39;c&#39;.
std::array< real, nstate > compute_entropy_variables(const std::array< real, nstate > &conservative_soln) const
Computes the entropy variables.
ConvectionDiffusion(const Parameters::AllParameters *const parameters_input, const bool convection=true, const bool diffusion=true, const dealii::Tensor< 2, 3, double > input_diffusion_tensor=Parameters::ManufacturedSolutionParam::get_default_diffusion_tensor(), const dealii::Tensor< 1, 3, double > input_advection_vector=Parameters::ManufacturedSolutionParam::get_default_advection_vector(), const double input_diffusion_coefficient=Parameters::ManufacturedSolutionParam::get_default_diffusion_coefficient(), std::shared_ptr< ManufacturedSolutionFunction< dim, real > > manufactured_solution_function=nullptr, const Parameters::AllParameters::TestType parameters_test=Parameters::AllParameters::TestType::run_control, const bool has_nonzero_physical_source=false)
Constructor.
static dealii::Tensor< 1, 3, double > get_default_advection_vector()
gets the default advection vector
TestType
Possible integration tests to run.
Base class from which Advection, Diffusion, ConvectionDiffusion, and Euler is derived.
Definition: physics.h:34
Manufactured solution used for grid studies to check convergence orders.
Convection-diffusion with linear advective and diffusive term. Derived from PhysicsBase.
std::array< dealii::Tensor< 1, dim, real >, nstate > convective_numerical_split_flux(const std::array< real, nstate > &soln1, const std::array< real, nstate > &soln2) const override
Convective numerical split flux for split form.
Files for the baseline physics.
Definition: ADTypes.hpp:10
const Parameters::AllParameters::TestType test_type
Allows convection diffusion to distinguish between different unsteady test types. ...
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
If diffusion is present, assign Dirichlet boundary condition.
std::array< real, nstate > source_term(const dealii::Point< dim, real > &pos, const std::array< real, nstate > &solution, const real current_time, const dealii::types::global_dof_index cell_index) const
Source term is zero or depends on manufactured solution.
double diffusion_scaling_coeff
Diffusion scaling coefficient in front of the diffusion tensor.
const bool hasConvection
Turns ON/OFF convection term.
Main parameter class that contains the various other sub-parameter classes.
double linear_advection_velocity[3]
Linear advection velocity in x, y, and z directions.
real max_viscous_eigenvalue(const std::array< real, nstate > &soln) const
Maximum viscous eigenvalue.
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.
static dealii::Tensor< 2, 3, double > get_default_diffusion_tensor()
gets the default diffusion tensor
real max_convective_eigenvalue(const std::array< real, nstate > &soln) const
Maximum convective eigenvalue.
dealii::Tensor< 1, dim, real > advection_speed() const
Linear advection speed: c.
const bool has_nonzero_physical_source
Flag to signal that physical source term is non-zero.
Definition: physics.h:62
std::array< dealii::Tensor< 1, dim, real >, nstate > dissipative_flux(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
Dissipative flux: u.
real diffusion_coefficient() const
Diffusion coefficient.
static double get_default_diffusion_coefficient()
gets the default diffusion coefficient;
std::shared_ptr< ManufacturedSolutionFunction< dim, real > > manufactured_solution_function
Manufactured solution function.
Definition: physics.h:71
std::array< dealii::Tensor< 1, dim, real >, nstate > convective_flux(const std::array< real, nstate > &solution) const
Convective flux: .