Alsvinn  0.5.3
The fast FVM simulator with UQ support
HLL3.hpp
Go to the documentation of this file.
1 /* Copyright (c) 2018 ETH Zurich, Kjetil Olsen Lye
2  * This program is free software: you can redistribute it and/or modify
3  * it under the terms of the GNU General Public License as published by
4  * the Free Software Foundation, either version 3 of the License, or
5  * (at your option) any later version.
6  *
7  * This program is distributed in the hope that it will be useful,
8  * but WITHOUT ANY WARRANTY; without even the implied warranty of
9  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
10  * GNU General Public License for more details.
11  *
12  * You should have received a copy of the GNU General Public License
13  * along with this program. If not, see <http://www.gnu.org/licenses/>.
14  */
15 
16 #pragma once
20 #include <algorithm>
21 #include <numeric>
23 
24 namespace alsfvm {
25 namespace numflux {
26 namespace euler {
27 
36 template<int nsd>
37 class HLL3 {
38 public:
42  static const std::string name;
43 
44  typedef typename Types<nsd>::rvec rvec;
45 
52  template<int direction>
53  __device__ __host__ inline static real computeFlux(const
57 
58  static_assert(direction < 3, "We only support three dimensions.");
59  real speedLeft;
60  real speedRight;
61  real cs;
62  equation::euler::ConservedVariables<nsd> fluxLeft, fluxRight;
63 
64  computeHLLSpeeds<direction>(eq, left, right, speedLeft, speedRight, cs);
65  eq.template computePointFlux<direction>(left, fluxLeft);
66  eq.template computePointFlux<direction>(right, fluxRight);
67 
68  const real pressureLeft = left.p;
69  const real pressureRight = right.p;
70 
71  const real udl = left.u[direction] - speedLeft;
72  const real udr = right.u[direction] - speedRight;
73 
74  real aa = udr * right.rho - udl * left.rho;
75  real sm = (right.m[direction] * udr - left.m[direction] * udl + pressureRight -
76  pressureLeft) / aa;
77 
78  rvec us = ( fluxRight.m - fluxLeft.m - speedRight * right.m + speedLeft *
79  left.m ) / aa;
80  us[direction] = sm;
81 
82  if (speedLeft == 0) {
84  } else if (speedLeft > 0) {
85  F = fluxLeft;
86  } else if (speedRight < 0) {
87  F = fluxRight;
88  } else if (sm >= 0) {
90  middle.rho = left.rho * udl / (sm - speedLeft);
91  middle.m = middle.rho * us;
92  const real p = pressureLeft + left.rho * (left.u[direction] - sm) * udl;
93  middle.E = (udl * left.E + pressureLeft * left.u[direction] - p * sm) /
94  (sm - speedLeft);
95  F = fluxLeft + speedLeft * (middle - left.conserved());
96  } else {
98  middle.rho = right.rho * udr / (sm - speedRight);
99  middle.m = middle.rho * us;
100  real p = pressureRight + right.rho * (right.u[direction] - sm) * udr;
101  middle.E = (udr * right.E + pressureRight * right.u[direction] - p * sm) /
102  (sm - speedRight);
103  F = fluxRight + speedRight * (middle - right.conserved());
104  }
105 
106  return fmax(fabs(speedLeft), fabs(speedRight));
107  }
108 
119  template<int direction>
121  const equation::euler::Euler<nsd>& eq,
123  const equation::euler::AllVariables<nsd>& right, real& speedLeft,
124  real& speedRight, real& cs) {
125 
126  const real waveLeft = sqrt(left.rho);
127  const real waveRight = sqrt(right.rho);
128 
129  const real rho = (left.rho + right.rho) / 2;
130  const rvec u = (waveLeft * left.u + waveRight * right.u) /
131  (waveLeft + waveRight);
132  const real p = (waveLeft * left.p + waveRight * right.p) /
133  (waveLeft + waveRight);
134 
135  // calculate extended fast speed cf at the left boundary
136  const real cfLeft = sqrt( eq.getGamma() * left.p / left.rho );
137 
138  // calculate extended fast speed cf at the right boundary
139  const real cfRight = sqrt( eq.getGamma() * right.p / right.rho );
140 
141  const real correct = 0.5 * fmax(real(0),
142  left.u[direction] - right.u[direction]);
143 
144  cs = sqrt(eq.getGamma() * p / rho);
145  speedLeft = fmin(left.u[direction] + correct - cfLeft, u[direction] - cs);
146  speedRight = fmax(right.u[direction] - correct + cfRight, u[direction] + cs);
147 
148  }
149 };
150 } // namespace alsfvm
151 } // namespace numflux
152 } // namespace euler
153 
Definition: types.hpp:104
__device__ __host__ real getGamma() const
Definition: Euler.hpp:299
Definition: HLL3.hpp:37
Types< nsd >::rvec rvec
Definition: HLL3.hpp:44
#define __host__
Definition: types.hpp:46
Definition: ConservedVariables.hpp:30
double real
Definition: types.hpp:65
__device__ static __host__ void computeHLLSpeeds(const equation::euler::Euler< nsd > &eq, const equation::euler::AllVariables< nsd > &left, const equation::euler::AllVariables< nsd > &right, real &speedLeft, real &speedRight, real &cs)
Definition: HLL3.hpp:120
int rho
Definition: sodshocktube.py:3
rvec m
Definition: ConservedVariables.hpp:93
rvec u
Definition: ExtraVariables.hpp:59
#define static_assert(x, y)
Definition: types.hpp:52
static const std::string name
name is "hll3"
Definition: HLL3.hpp:42
#define __device__
Definition: types.hpp:45
Definition: Euler.hpp:36
real p
Definition: ExtraVariables.hpp:58
real rho
Definition: ConservedVariables.hpp:92
Definition: AllVariables.hpp:27
Various utility functions to implement the tecno flux.
Definition: types.hpp:30
float p
Definition: sodshocktube.py:5
__device__ static __host__ real computeFlux(const equation::euler::Euler< nsd > &eq, const equation::euler::AllVariables< nsd > &left, const equation::euler::AllVariables< nsd > &right, equation::euler::ConservedVariables< nsd > &F)
Definition: HLL3.hpp:53
real E
Definition: ConservedVariables.hpp:94