[P]arallel [Hi]gh-order [Li]brary for [P]DEs  Latest
Parallel High-Order Library for PDEs through hp-adaptive Discontinuous Galerkin methods
convective_numerical_flux.cpp
1 #include "ADTypes.hpp"
2 
3 #include "convective_numerical_flux.hpp"
4 
5 namespace PHiLiP {
6 namespace NumericalFlux {
7 
8 using AllParam = Parameters::AllParameters;
9 
10 // Protyping low level functions
11 template<int nstate, typename real_tensor>
12 std::array<real_tensor, nstate> array_average(
13  const std::array<real_tensor, nstate> &array1,
14  const std::array<real_tensor, nstate> &array2)
15 {
16  std::array<real_tensor,nstate> array_average;
17  for (int s=0; s<nstate; s++) {
18  array_average[s] = 0.5*(array1[s] + array2[s]);
19  }
20  return array_average;
21 }
22 
23 template <int dim, int nstate, typename real>
25  std::unique_ptr< BaselineNumericalFluxConvective<dim,nstate,real> > baseline_input,
26  std::unique_ptr< RiemannSolverDissipation<dim,nstate,real> > riemann_solver_dissipation_input)
27  : baseline(std::move(baseline_input))
28  , riemann_solver_dissipation(std::move(riemann_solver_dissipation_input))
29 { }
30 
31 template<int dim, int nstate, typename real>
32 std::array<real, nstate> NumericalFluxConvective<dim,nstate,real>
34  const std::array<real, nstate> &soln_int,
35  const std::array<real, nstate> &soln_ext,
36  const dealii::Tensor<1,dim,real> &normal_int) const
37 {
38  // baseline flux (without upwind dissipation)
39  const std::array<real, nstate> baseline_flux_dot_n
40  = this->baseline->evaluate_flux(soln_int, soln_ext, normal_int);
41 
42  // Riemann solver dissipation
43  const std::array<real, nstate> riemann_solver_dissipation_dot_n
44  = this->riemann_solver_dissipation->evaluate_riemann_solver_dissipation(soln_int, soln_ext, normal_int);
45 
46  // convective numerical flux: sum of baseline and Riemann solver dissipation term
47  std::array<real, nstate> numerical_flux_dot_n;
48  for (int s=0; s<nstate; s++) {
49  numerical_flux_dot_n[s] = baseline_flux_dot_n[s] + riemann_solver_dissipation_dot_n[s];
50  }
51  return numerical_flux_dot_n;
52 }
53 
54 template <int dim, int nstate, typename real>
56  std::shared_ptr<Physics::PhysicsBase<dim, nstate, real>> physics_input)
57  : NumericalFluxConvective<dim,nstate,real>(
58  std::make_unique< CentralBaselineNumericalFluxConvective<dim, nstate, real> > (physics_input),
59  std::make_unique< LaxFriedrichsRiemannSolverDissipation<dim, nstate, real> > (physics_input))
60 {}
61 
62 template <int dim, int nstate, typename real>
64  std::shared_ptr<Physics::PhysicsBase<dim, nstate, real>> physics_input)
65  : NumericalFluxConvective<dim,nstate,real>(
66  std::make_unique< CentralBaselineNumericalFluxConvective<dim, nstate, real> > (physics_input),
67  std::make_unique< RoePikeRiemannSolverDissipation<dim, nstate, real> > (physics_input))
68 {}
69 
70 template <int dim, int nstate, typename real>
72  std::shared_ptr<Physics::PhysicsBase<dim, nstate, real>> physics_input)
73  : NumericalFluxConvective<dim,nstate,real>(
74  std::make_unique< CentralBaselineNumericalFluxConvective<dim, nstate, real> > (physics_input),
75  std::make_unique< L2RoeRiemannSolverDissipation<dim, nstate, real> > (physics_input))
76 {}
77 
78 template <int dim, int nstate, typename real>
80  std::shared_ptr<Physics::PhysicsBase<dim, nstate, real>> physics_input)
81  : NumericalFluxConvective<dim,nstate,real>(
82  std::make_unique< CentralBaselineNumericalFluxConvective<dim, nstate, real> > (physics_input),
83  std::make_unique< ZeroRiemannSolverDissipation<dim, nstate, real> > ())
84 {}
85 
86 template <int dim, int nstate, typename real>
88  std::shared_ptr<Physics::PhysicsBase<dim, nstate, real>> physics_input)
89  : NumericalFluxConvective<dim,nstate,real>(
90  std::make_unique< EntropyConservingBaselineNumericalFluxConvective<dim, nstate, real> > (physics_input),
91  std::make_unique< ZeroRiemannSolverDissipation<dim, nstate, real> > ())
92 {}
93 
94 template <int dim, int nstate, typename real>
96  std::shared_ptr<Physics::PhysicsBase<dim, nstate, real>> physics_input)
97  : NumericalFluxConvective<dim,nstate,real>(
98  std::make_unique< EntropyConservingBaselineNumericalFluxConvective<dim, nstate, real> > (physics_input),
99  std::make_unique< LaxFriedrichsRiemannSolverDissipation<dim, nstate, real> > (physics_input))
100 {}
101 
102 template <int dim, int nstate, typename real>
104  std::shared_ptr<Physics::PhysicsBase<dim, nstate, real>> physics_input)
105  : NumericalFluxConvective<dim,nstate,real>(
106  std::make_unique< EntropyConservingBaselineNumericalFluxConvective<dim, nstate, real> > (physics_input),
107  std::make_unique< RoePikeRiemannSolverDissipation<dim, nstate, real> > (physics_input))
108 {}
109 
110 template <int dim, int nstate, typename real>
112  std::shared_ptr<Physics::PhysicsBase<dim, nstate, real>> physics_input)
113  : NumericalFluxConvective<dim,nstate,real>(
114  std::make_unique< EntropyConservingBaselineNumericalFluxConvective<dim, nstate, real> > (physics_input),
115  std::make_unique< L2RoeRiemannSolverDissipation<dim, nstate, real> > (physics_input))
116 {}
117 
118 template <int dim, int nstate, typename real>
120  const std::array<real, nstate> &soln_int,
121  const std::array<real, nstate> &soln_ext,
122  const dealii::Tensor<1,dim,real> &normal_int) const
123 {
124  using RealArrayVector = std::array<dealii::Tensor<1,dim,real>,nstate>;
125  RealArrayVector conv_phys_flux_int;
126  RealArrayVector conv_phys_flux_ext;
127 
128  conv_phys_flux_int = pde_physics->convective_flux (soln_int);
129  conv_phys_flux_ext = pde_physics->convective_flux (soln_ext);
130 
131  RealArrayVector flux_avg;
132  for (int s=0; s<nstate; s++) {
133  flux_avg[s] = 0.0;
134  for (int d=0; d<dim; ++d) {
135  flux_avg[s][d] = 0.5*(conv_phys_flux_int[s][d] + conv_phys_flux_ext[s][d]);
136  }
137  }
138 
139  std::array<real, nstate> numerical_flux_dot_n;
140  for (int s=0; s<nstate; s++) {
141  real flux_dot_n = 0.0;
142  for (int d=0; d<dim; ++d) {
143  flux_dot_n += flux_avg[s][d]*normal_int[d];
144  }
145  numerical_flux_dot_n[s] = flux_dot_n;
146  }
147  return numerical_flux_dot_n;
148 }
149 
150 template <int dim, int nstate, typename real>
152  const std::array<real, nstate> &soln_int,
153  const std::array<real, nstate> &soln_ext,
154  const dealii::Tensor<1,dim,real> &normal_int) const
155 {
156  using RealArrayVector = std::array<dealii::Tensor<1,dim,real>,nstate>;
157  RealArrayVector conv_phys_split_flux;
158 
159  conv_phys_split_flux = pde_physics->convective_numerical_split_flux (soln_int,soln_ext);
160 
161  // Scalar dissipation
162  std::array<real, nstate> numerical_flux_dot_n;
163  for (int s=0; s<nstate; s++) {
164  real flux_dot_n = 0.0;
165  for (int d=0; d<dim; ++d) {
166  flux_dot_n += conv_phys_split_flux[s][d] * normal_int[d];
167  }
168  numerical_flux_dot_n[s] = flux_dot_n;
169  }
170  return numerical_flux_dot_n;
171 }
172 
173 template<int dim, int nstate, typename real>
174 std::array<real, nstate> ZeroRiemannSolverDissipation<dim,nstate,real>
176  const std::array<real, nstate> &/*soln_int*/,
177  const std::array<real, nstate> &/*soln_ext*/,
178  const dealii::Tensor<1,dim,real> &/*normal_int*/) const
179 {
180  // zero upwind dissipation
181  std::array<real, nstate> numerical_flux_dot_n;
182  numerical_flux_dot_n.fill(0.0);
183  return numerical_flux_dot_n;
184 }
185 
186 template<int dim, int nstate, typename real>
189  const std::array<real, nstate> &soln_int,
190  const std::array<real, nstate> &soln_ext,
191  const dealii::Tensor<1,dim,real> &normal_int) const
192 {
193  const real conv_max_eig_int = pde_physics->max_convective_normal_eigenvalue(soln_int,normal_int);
194  const real conv_max_eig_ext = pde_physics->max_convective_normal_eigenvalue(soln_ext,normal_int);
195  // Replaced the std::max with an if-statement for the AD to work properly.
196  //const real conv_max_eig = std::max(conv_max_eig_int, conv_max_eig_ext);
197  real conv_max_eig;
198  if (conv_max_eig_int > conv_max_eig_ext) {
199  conv_max_eig = conv_max_eig_int;
200  } else {
201  conv_max_eig = conv_max_eig_ext;
202  }
203  //conv_max_eig = std::max(conv_max_eig_int, conv_max_eig_ext);
204  // Scalar dissipation
205  std::array<real, nstate> numerical_flux_dot_n;
206  for (int s=0; s<nstate; s++) {
207  numerical_flux_dot_n[s] = - 0.5 * conv_max_eig * (soln_ext[s]-soln_int[s]);
208  }
209 
210  return numerical_flux_dot_n;
211 }
212 
213 template <int dim, int nstate, typename real>
216  const std::array<real, 3> &eig_L,
217  const std::array<real, 3> &eig_R,
218  std::array<real, 3> &eig_RoeAvg,
219  const real /*vel2_ravg*/,
220  const real /*sound_ravg*/) const
221 {
222  // Harten's entropy fix
223  // -- See Blazek 2015, p.103-105
224  for(int e=0;e<3;e++) {
225  const real eps = std::max(abs(eig_RoeAvg[e] - eig_L[e]), abs(eig_R[e] - eig_RoeAvg[e]));
226  if(eig_RoeAvg[e] < eps) {
227  eig_RoeAvg[e] = 0.5*(eig_RoeAvg[e] * eig_RoeAvg[e]/eps + eps);
228  }
229  }
230 }
231 
232 template <int dim, int nstate, typename real>
235  const std::array<real, nstate> &/*soln_int*/,
236  const std::array<real, nstate> &/*soln_ext*/,
237  const std::array<real, 3> &/*eig_L*/,
238  const std::array<real, 3> &/*eig_R*/,
239  real &/*dV_normal*/,
240  dealii::Tensor<1,dim,real> &/*dV_tangent*/
241  ) const
242 {
243  // No additional modifications for the Roe-Pike scheme
244 }
245 
246 template <int dim, int nstate, typename real>
249  const std::array<real, 3> &eig_L,
250  const std::array<real, 3> &eig_R,
251  int &ssw_LEFT,
252  int &ssw_RIGHT) const
253 {
254  // Shock indicator of Wada & Liou (1994 Flux) -- Eq.(39)
255  // -- See also p.74 of Osswald et al. (2016 L2Roe)
256 
257  ssw_LEFT=0; ssw_RIGHT=0; // initialize
258 
259  // ssw_L: i=L --> j=R
260  if((eig_L[0]>0.0 && eig_R[0]<0.0) || (eig_L[2]>0.0 && eig_R[2]<0.0)) {
261  ssw_LEFT = 1;
262  }
263 
264  // ssw_R: i=R --> j=L
265  if((eig_R[0]>0.0 && eig_L[0]<0.0) || (eig_R[2]>0.0 && eig_L[2]<0.0)) {
266  ssw_RIGHT = 1;
267  }
268 }
269 
270 template <int dim, int nstate, typename real>
273  const std::array<real, 3> &eig_L,
274  const std::array<real, 3> &eig_R,
275  std::array<real, 3> &eig_RoeAvg,
276  const real vel2_ravg,
277  const real sound_ravg) const
278 {
279  // Van Leer et al. (1989 Sonic) entropy fix for acoustic waves
280  // -- p.74 of Osswald et al. (2016 L2Roe)
281  for(int e=0;e<3;e++) {
282  if(e!=1) {
283  // const real deig = std::max((eig_R[e]-eig_L[e]), 0.0);
284  const real deig = std::max(static_cast<real>(eig_R[e] - eig_L[e]), static_cast<real>(0.0));
285  if(eig_RoeAvg[e] < 2.0*deig) {
286  eig_RoeAvg[e] = 0.25*(eig_RoeAvg[e] * eig_RoeAvg[e]/deig) + deig;
287  }
288  }
289  }
290 
291  // Entropy fix of Liou (2000 Mass)
292  // -- p.74 of Osswald et al. (2016 L2Roe)
293  int ssw_L, ssw_R;
294  evaluate_shock_indicator(eig_L,eig_R,ssw_L,ssw_R);
295  if(ssw_L!=0 || ssw_R!=0) {
296  eig_RoeAvg[1] = std::max(sound_ravg, static_cast<real>(sqrt(vel2_ravg)));
297  }
298 }
299 
300 template <int dim, int nstate, typename real>
303  const std::array<real, nstate> &soln_int,
304  const std::array<real, nstate> &soln_ext,
305  const std::array<real, 3> &eig_L,
306  const std::array<real, 3> &eig_R,
307  real &dV_normal,
308  dealii::Tensor<1,dim,real> &dV_tangent) const
309 {
310  const real mach_number_L = this->euler_physics->compute_mach_number(soln_int);
311  const real mach_number_R = this->euler_physics->compute_mach_number(soln_ext);
312 
313  // Osswald's two modifications to Roe-Pike scheme --> L2Roe
314  // - Blending factor (variable 'z' in reference)
315  const real blending_factor = std::min(static_cast<real>(1.0), std::max(mach_number_L,mach_number_R));
316  // - Scale jump in (1) normal and (2) tangential velocities
317  int ssw_L, ssw_R;
318  evaluate_shock_indicator(eig_L,eig_R,ssw_L,ssw_R);
319  if(ssw_L==0 && ssw_R==0)
320  {
321  dV_normal *= blending_factor;
322  for (int d=0;d<dim;d++)
323  {
324  dV_tangent[d] *= blending_factor;
325  }
326  }
327 }
328 
329 template <int dim, int nstate, typename real>
332  const std::array<real, nstate> &soln_int,
333  const std::array<real, nstate> &soln_ext,
334  const dealii::Tensor<1,dim,real> &normal_int) const
335 {
336  // See Blazek 2015, p.103-105
337  // -- Note: Modified calculation of alpha_{3,4} to use
338  // dVt (jump in tangential velocities);
339  // expressions are equivalent
340 
341  // Blazek 2015
342  // p. 103-105
343  // Note: This is in fact the Roe-Pike method of Roe & Pike (1984 - Efficient)
344  const std::array<real,nstate> prim_soln_int = euler_physics->convert_conservative_to_primitive(soln_int);
345  const std::array<real,nstate> prim_soln_ext = euler_physics->convert_conservative_to_primitive(soln_ext);
346  // Left cell
347  const real density_L = prim_soln_int[0];
348  const dealii::Tensor< 1,dim,real > velocities_L = euler_physics->extract_velocities_from_primitive(prim_soln_int);
349  const real pressure_L = prim_soln_int[nstate-1];
350 
351  //const real normal_vel_L = velocities_L*normal_int;
352  real normal_vel_L = 0.0;
353  for (int d=0; d<dim; ++d) {
354  normal_vel_L+= velocities_L[d]*normal_int[d];
355  }
356  const real specific_enthalpy_L = euler_physics->compute_specific_enthalpy(soln_int, pressure_L);
357 
358  // Right cell
359  const real density_R = prim_soln_ext[0];
360  const dealii::Tensor< 1,dim,real > velocities_R = euler_physics->extract_velocities_from_primitive(prim_soln_ext);
361  const real pressure_R = prim_soln_ext[nstate-1];
362 
363  //const real normal_vel_R = velocities_R*normal_int;
364  real normal_vel_R = 0.0;
365  for (int d=0; d<dim; ++d) {
366  normal_vel_R+= velocities_R[d]*normal_int[d];
367  }
368  const real specific_enthalpy_R = euler_physics->compute_specific_enthalpy(soln_ext, pressure_R);
369 
370  // Roe-averaged states
371  const real r = sqrt(density_R/density_L);
372  const real rp1 = r+1.0;
373 
374  const real density_ravg = r*density_L;
375  //const dealii::Tensor< 1,dim,real > velocities_ravg = (r*velocities_R + velocities_L) / rp1;
376  dealii::Tensor< 1,dim,real > velocities_ravg;
377  for (int d=0; d<dim; ++d) {
378  velocities_ravg[d] = (r*velocities_R[d] + velocities_L[d]) / rp1;
379  }
380  const real specific_total_enthalpy_ravg = (r*specific_enthalpy_R + specific_enthalpy_L) / rp1;
381 
382  const real vel2_ravg = euler_physics->compute_velocity_squared (velocities_ravg);
383  //const real normal_vel_ravg = velocities_ravg*normal_int;
384  real normal_vel_ravg = 0.0;
385  for (int d=0; d<dim; ++d) {
386  normal_vel_ravg += velocities_ravg[d]*normal_int[d];
387  }
388 
389  const real sound2_ravg = euler_physics->gamm1*(specific_total_enthalpy_ravg-0.5*vel2_ravg);
390  real sound_ravg = 1e10;
391  if (sound2_ravg > 0.0) {
392  sound_ravg = sqrt(sound2_ravg);
393  }
394 
395  // Compute eigenvalues
396  std::array<real, 3> eig_ravg;
397  eig_ravg[0] = abs(normal_vel_ravg-sound_ravg);
398  eig_ravg[1] = abs(normal_vel_ravg);
399  eig_ravg[2] = abs(normal_vel_ravg+sound_ravg);
400 
401  const real sound_L = euler_physics->compute_sound(density_L, pressure_L);
402  std::array<real, 3> eig_L;
403  eig_L[0] = abs(normal_vel_L-sound_L);
404  eig_L[1] = abs(normal_vel_L);
405  eig_L[2] = abs(normal_vel_L+sound_L);
406 
407  const real sound_R = euler_physics->compute_sound(density_R, pressure_R);
408  std::array<real, 3> eig_R;
409  eig_R[0] = abs(normal_vel_R-sound_R);
410  eig_R[1] = abs(normal_vel_R);
411  eig_R[2] = abs(normal_vel_R+sound_R);
412 
413  // Jumps in pressure and density
414  const real dp = pressure_R - pressure_L;
415  const real drho = density_R - density_L;
416 
417  // Jump in normal velocity
418  real dVn = normal_vel_R-normal_vel_L;
419 
420  // Jumps in tangential velocities
421  dealii::Tensor<1,dim,real> dVt;
422  for (int d=0;d<dim;d++) {
423  dVt[d] = (velocities_R[d] - velocities_L[d]) - dVn*normal_int[d];
424  }
425 
426  // Evaluate entropy fix on wave speeds
427  evaluate_entropy_fix (eig_L, eig_R, eig_ravg, vel2_ravg, sound_ravg);
428 
429  // Evaluate additional modifications to the Roe-Pike scheme (if applicable)
430  evaluate_additional_modifications (soln_int, soln_ext, eig_L, eig_R, dVn, dVt);
431 
432  // Product of eigenvalues and wave strengths
433  real coeff[4];
434  coeff[0] = eig_ravg[0]*(dp-density_ravg*sound_ravg*dVn)/(2.0*sound2_ravg);
435  coeff[1] = eig_ravg[1]*(drho - dp/sound2_ravg);
436  coeff[2] = eig_ravg[1]*density_ravg;
437  coeff[3] = eig_ravg[2]*(dp+density_ravg*sound_ravg*dVn)/(2.0*sound2_ravg);
438 
439  // Evaluate |A_Roe| * (W_R - W_L)
440  std::array<real,nstate> AdW;
441 
442  // Vn-c (i=1)
443  AdW[0] = coeff[0] * 1.0;
444  for (int d=0;d<dim;d++) {
445  AdW[1+d] = coeff[0] * (velocities_ravg[d] - sound_ravg * normal_int[d]);
446  }
447  AdW[nstate-1] = coeff[0] * (specific_total_enthalpy_ravg - sound_ravg*normal_vel_ravg);
448 
449  // Vn (i=2)
450  AdW[0] += coeff[1] * 1.0;
451  for (int d=0;d<dim;d++) {
452  AdW[1+d] += coeff[1] * velocities_ravg[d];
453  }
454  AdW[nstate-1] += coeff[1] * vel2_ravg * 0.5;
455 
456  // (i=3,4)
457  AdW[0] += coeff[2] * 0.0;
458  real dVt_dot_vel_ravg = 0.0;
459  for (int d=0;d<dim;d++) {
460  AdW[1+d] += coeff[2]*dVt[d];
461  dVt_dot_vel_ravg += velocities_ravg[d]*dVt[d];
462  }
463  AdW[nstate-1] += coeff[2]*dVt_dot_vel_ravg;
464 
465  // Vn+c (i=5)
466  AdW[0] += coeff[3] * 1.0;
467  for (int d=0;d<dim;d++) {
468  AdW[1+d] += coeff[3] * (velocities_ravg[d] + sound_ravg * normal_int[d]);
469  }
470  AdW[nstate-1] += coeff[3] * (specific_total_enthalpy_ravg + sound_ravg*normal_vel_ravg);
471 
472  std::array<real, nstate> numerical_flux_dot_n;
473  for (int s=0; s<nstate; s++) {
474  numerical_flux_dot_n[s] = - 0.5 * AdW[s];
475  }
476 
477  return numerical_flux_dot_n;
478 }
479 
480 // Instantiation
499 
512 
543 
554 
555 template class Central<PHILIP_DIM, 1, double>;
556 template class Central<PHILIP_DIM, 2, double>;
557 template class Central<PHILIP_DIM, 3, double>;
558 template class Central<PHILIP_DIM, 4, double>;
559 template class Central<PHILIP_DIM, 5, double>;
560 //template class Central<PHILIP_DIM, 6, double>;
561 template class Central<PHILIP_DIM, 1, FadType >;
562 template class Central<PHILIP_DIM, 2, FadType >;
563 template class Central<PHILIP_DIM, 3, FadType >;
564 template class Central<PHILIP_DIM, 4, FadType >;
565 template class Central<PHILIP_DIM, 5, FadType >;
566 //template class Central<PHILIP_DIM, 6, FadType >;
567 template class Central<PHILIP_DIM, 1, RadType >;
568 template class Central<PHILIP_DIM, 2, RadType >;
569 template class Central<PHILIP_DIM, 3, RadType >;
570 template class Central<PHILIP_DIM, 4, RadType >;
571 template class Central<PHILIP_DIM, 5, RadType >;
572 //template class Central<PHILIP_DIM, 6, RadType >;
578 //template class Central<PHILIP_DIM, 6, FadFadType >;
584 //template class Central<PHILIP_DIM, 6, RadFadType >;
585 
591 //template class EntropyConserving<PHILIP_DIM, 6, double>;
597 //template class EntropyConserving<PHILIP_DIM, 6, FadType >;
603 //template class EntropyConserving<PHILIP_DIM, 6, RadType >;
609 //template class EntropyConserving<PHILIP_DIM, 6, FadFadType >;
615 //template class EntropyConserving<PHILIP_DIM, 6, RadFadType >;
616 
622 //template class EntropyConservingWithLaxFriedrichsDissipation<PHILIP_DIM, 6, double>;
628 //template class EntropyConservingWithLaxFriedrichsDissipation<PHILIP_DIM, 6, FadType >;
634 //template class EntropyConservingWithLaxFriedrichsDissipation<PHILIP_DIM, 6, RadType >;
640 //template class EntropyConservingWithLaxFriedrichsDissipation<PHILIP_DIM, 6, FadFadType >;
646 //template class EntropyConservingWithLaxFriedrichsDissipation<PHILIP_DIM, 6, RadFadType >;
647 
653 
659 
660 //template class EntropyConservingWithRoeDissipation<PHILIP_DIM, PHILIP_DIM+3, double>;
661 //template class EntropyConservingWithRoeDissipation<PHILIP_DIM, PHILIP_DIM+3, FadType >;
662 //template class EntropyConservingWithRoeDissipation<PHILIP_DIM, PHILIP_DIM+3, RadType >;
663 //template class EntropyConservingWithRoeDissipation<PHILIP_DIM, PHILIP_DIM+3, FadFadType >;
664 //template class EntropyConservingWithRoeDissipation<PHILIP_DIM, PHILIP_DIM+3, RadFadType >;
665 //
666 //template class EntropyConservingWithL2RoeDissipation<PHILIP_DIM, PHILIP_DIM+3, double>;
667 //template class EntropyConservingWithL2RoeDissipation<PHILIP_DIM, PHILIP_DIM+3, FadType >;
668 //template class EntropyConservingWithL2RoeDissipation<PHILIP_DIM, PHILIP_DIM+3, RadType >;
669 //template class EntropyConservingWithL2RoeDissipation<PHILIP_DIM, PHILIP_DIM+3, FadFadType >;
670 //template class EntropyConservingWithL2RoeDissipation<PHILIP_DIM, PHILIP_DIM+3, RadFadType >;
671 
702 
733 
739 //template class EntropyConservingBaselineNumericalFluxConvective<PHILIP_DIM, 6, double>;
745 //template class EntropyConservingBaselineNumericalFluxConvective<PHILIP_DIM, 6, FadType >;
751 //template class EntropyConservingBaselineNumericalFluxConvective<PHILIP_DIM, 6, RadType >;
757 //template class EntropyConservingBaselineNumericalFluxConvective<PHILIP_DIM, 6, FadFadType >;
763 //template class EntropyConservingBaselineNumericalFluxConvective<PHILIP_DIM, 6, RadFadType >;
764 
795 
801 //template class ZeroRiemannSolverDissipation<PHILIP_DIM, 6, double>;
807 //template class ZeroRiemannSolverDissipation<PHILIP_DIM, 6, FadType >;
813 //template class ZeroRiemannSolverDissipation<PHILIP_DIM, 6, RadType >;
819 //template class ZeroRiemannSolverDissipation<PHILIP_DIM, 6, FadFadType >;
825 //template class ZeroRiemannSolverDissipation<PHILIP_DIM, 6, RadFadType >;
826 
857 
863 
869 
875 
876 } // NumericalFlux namespace
877 } // PHiLiP namespace
std::unique_ptr< BaselineNumericalFluxConvective< dim, nstate, real > > baseline
Baseline convective numerical flux object.
EntropyConservingWithRoeDissipation(std::shared_ptr< Physics::PhysicsBase< dim, nstate, real >> physics_input)
Constructor.
std::unique_ptr< RiemannSolverDissipation< dim, nstate, real > > riemann_solver_dissipation
Upwind convective numerical flux object.
RoePike(std::shared_ptr< Physics::PhysicsBase< dim, nstate, real >> physics_input)
Constructor.
Central numerical flux. Derived from BaselineNumericalFluxConvective.
Base class from which Advection, Diffusion, ConvectionDiffusion, and Euler is derived.
Definition: physics.h:34
void evaluate_entropy_fix(const std::array< real, 3 > &eig_L, const std::array< real, 3 > &eig_R, std::array< real, 3 > &eig_RoeAvg, const real vel2_ravg, const real sound_ravg) const
Entropy conserving numerical flux with Roe dissipation. Derived from NumericalFluxConvective.
std::array< real, nstate > evaluate_riemann_solver_dissipation(const std::array< real, nstate > &soln_int, const std::array< real, nstate > &soln_ext, const dealii::Tensor< 1, dim, real > &normal1) const
Returns zeros.
Files for the baseline physics.
Definition: ADTypes.hpp:10
EntropyConservingWithL2RoeDissipation(std::shared_ptr< Physics::PhysicsBase< dim, nstate, real >> physics_input)
Constructor.
EntropyConservingWithLaxFriedrichsDissipation(std::shared_ptr< Physics::PhysicsBase< dim, nstate, real >> physics_input)
Constructor.
Roe-Pike numerical flux. Derived from NumericalFluxConvective.
std::array< real, nstate > evaluate_riemann_solver_dissipation(const std::array< real, nstate > &soln_int, const std::array< real, nstate > &soln_ext, const dealii::Tensor< 1, dim, real > &normal1) const
Entropy conserving numerical flux with L2Roe dissipation. Derived from NumericalFluxConvective.
Base class of Roe (Roe-Pike) flux with entropy fix. Derived from RiemannSolverDissipation.
Base class of Riemann solver dissipation (i.e. upwind-term) for numerical flux associated with convec...
void evaluate_additional_modifications(const std::array< real, nstate > &soln_int, const std::array< real, nstate > &soln_ext, const std::array< real, 3 > &eig_L, const std::array< real, 3 > &eig_R, real &dV_normal, dealii::Tensor< 1, dim, real > &dV_tangent) const
Entropy conserving numerical flux. Derived from NumericalFluxConvective.
void evaluate_entropy_fix(const std::array< real, 3 > &eig_L, const std::array< real, 3 > &eig_R, std::array< real, 3 > &eig_RoeAvg, const real vel2_ravg, const real sound_ravg) const
NumericalFluxConvective(std::unique_ptr< BaselineNumericalFluxConvective< dim, nstate, real > > baseline_input, std::unique_ptr< RiemannSolverDissipation< dim, nstate, real > > riemann_solver_dissipation_input)
Constructor.
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.
LaxFriedrichs(std::shared_ptr< Physics::PhysicsBase< dim, nstate, real >> physics_input)
Constructor.
Lax-Friedrichs Riemann solver dissipation. Derived from RiemannSolverDissipation. ...
RoePike flux with entropy fix. Derived from RoeBase.
std::array< real, nstate > evaluate_riemann_solver_dissipation(const std::array< real, nstate > &soln_int, const std::array< real, nstate > &soln_ext, const dealii::Tensor< 1, dim, real > &normal1) const
Entropy conserving numerical flux with Lax Friedrichs dissipation. Derived from NumericalFluxConvecti...
Zero Riemann solver dissipation. Derived from RiemannSolverDissipation.
void evaluate_additional_modifications(const std::array< real, nstate > &soln_int, const std::array< real, nstate > &soln_ext, const std::array< real, 3 > &eig_L, const std::array< real, 3 > &eig_R, real &dV_normal, dealii::Tensor< 1, dim, real > &dV_tangent) const
Empty function. No additional modifications for the Roe-Pike scheme.
EntropyConserving(std::shared_ptr< Physics::PhysicsBase< dim, nstate, real >> physics_input)
Constructor.
Central numerical flux. Derived from NumericalFluxConvective.
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
Returns the convective numerical flux at an interface.
void evaluate_shock_indicator(const std::array< real, 3 > &eig_L, const std::array< real, 3 > &eig_R, int &ssw_LEFT, int &ssw_RIGHT) const
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
Returns the convective numerical flux at an interface.
Lax-Friedrichs numerical flux. Derived from NumericalFluxConvective.
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
Returns the convective numerical flux at an interface.
Central(std::shared_ptr< Physics::PhysicsBase< dim, nstate, real >> physics_input)
Constructor.
L2Roe(std::shared_ptr< Physics::PhysicsBase< dim, nstate, real >> physics_input)
Constructor.
Entropy Conserving Numerical Flux. Derived from BaselineNumericalFluxConvective.