Line data Source code
1 : #ifndef SOLVER_LOCAL_RIEMANN_H
2 : #define SOLVER_LOCAL_RIEMANN_H
3 :
4 : #include "Set/Base.H"
5 : #include "Set/Set.H"
6 : #include "Model/Gas/Gas.H"
7 : #include <memory>
8 :
9 : namespace Solver
10 : {
11 : namespace Local
12 : {
13 : namespace Riemann
14 : {
15 : struct State {
16 : Set::Scalar rho = NAN;
17 : Set::Scalar M_normal = NAN;
18 : Set::Scalar M_tangent = NAN;
19 : Set::Scalar E = NAN;
20 : // Construtor for convenience
21 : State() { rho = 0.0; M_normal = 0.0, M_tangent = 0.0; E = 0.0;}
22 279746024 : State(Set::Scalar a_rho, Set::Scalar a_M_normal, Set::Scalar a_M_tangent, Set::Scalar a_E)
23 279746024 : : rho(a_rho), M_normal(a_M_normal), M_tangent(a_M_tangent), E(a_E) {}
24 158822400 : State(Set::Patch<const Set::Scalar> density_mf, Set::Patch<const Set::Scalar> momentum_mf, Set::Patch<const Set::Scalar> energy_mf, int i, int j, int k, int direction)
25 158822400 : {
26 158822400 : rho = density_mf(i,j,k);
27 158822400 : if (direction == 0)
28 : {
29 79411200 : M_normal = momentum_mf(i,j,k,0);
30 79411200 : M_tangent = momentum_mf(i,j,k,1);
31 : }
32 79411200 : else if (direction == 1)
33 : {
34 79411200 : M_normal = momentum_mf(i,j,k,1);
35 79411200 : M_tangent = momentum_mf(i,j,k,0);
36 : }
37 : else
38 : {
39 0 : Util::Abort(INFO, "Not supported yet");
40 : }
41 158822400 : E = energy_mf(i,j,k);
42 158822400 : }
43 :
44 :
45 0 : friend std::ostream &operator<<(std::ostream &os, const State &state)
46 : {
47 0 : os << "rho=" << state.rho << ", ";
48 0 : os << "Mn=" << state.M_normal << ", ";
49 0 : os << "Mt=" << state.M_tangent << ", ";
50 0 : os << "E=" << state.E << ", ";
51 : // do stuf
52 0 : return os;
53 : }
54 : void operator += (const State &a)
55 : {
56 : rho += a.rho; M_normal += a.M_normal; M_tangent += a.M_tangent; E += a.E;
57 : };
58 : void operator -= (const State &a)
59 : {
60 : rho -= a.rho; M_normal -= a.M_normal; M_tangent -= a.M_tangent; E -= a.E;
61 : };
62 8401992 : void operator *= (const Set::Scalar alpha)
63 : {
64 8401992 : rho *= alpha; M_normal *= alpha; M_tangent *= alpha; E *= alpha;
65 8401992 : };
66 : void operator /= (const Set::Scalar alpha)
67 : {
68 : rho /= alpha; M_normal /= alpha; M_tangent /= alpha; E /= alpha;
69 : };
70 : friend State operator + (const State &a, const State &b)
71 : {
72 : return State(a.rho+b.rho, a.M_normal+b.M_normal, a.M_tangent+b.M_tangent, a.E+b.E);
73 : }
74 95717648 : friend State operator - (const State &a, const State &b)
75 : {
76 95717648 : return State(a.rho-b.rho, a.M_normal-b.M_normal, a.M_tangent-b.M_tangent, a.E-b.E);
77 : }
78 87813192 : friend State operator * (const Set::Scalar alpha, const State &b)
79 : {
80 87813192 : return State(b.rho*alpha, b.M_normal*alpha, b.M_tangent*alpha, b.E*alpha);
81 : }
82 : friend State operator * (const State &b, const Set::Scalar alpha)
83 : {
84 : return State(b.rho*alpha, b.M_normal*alpha, b.M_tangent*alpha, b.E*alpha);
85 : }
86 79411200 : friend State operator / (const State &b, const Set::Scalar alpha)
87 : {
88 79411200 : return State(b.rho/alpha, b.M_normal/alpha, b.M_tangent/alpha, b.E/alpha);
89 : }
90 : };
91 :
92 : struct Flux {
93 : Set::Scalar mass;
94 : Set::Scalar momentum_normal;
95 : Set::Scalar momentum_tangent;
96 : Set::Scalar energy;
97 88985600 : Flux() : mass(0.0), momentum_normal(0.0), momentum_tangent(0.0), energy(0.0) {}
98 133857528 : Flux(Set::Scalar a_mass, Set::Scalar a_momentum_normal, Set::Scalar a_momentum_tangent, Set::Scalar a_energy) :
99 133857528 : mass(a_mass), momentum_normal(a_momentum_normal),
100 133857528 : momentum_tangent(a_momentum_tangent), energy(a_energy) {}
101 :
102 16306448 : Flux(const State &in) :
103 16306448 : mass(in.rho),
104 16306448 : momentum_normal(in.M_normal),
105 16306448 : momentum_tangent(in.M_tangent),
106 16306448 : energy(in.E) {}
107 :
108 : friend std::ostream &operator<<(std::ostream &os, const Flux &flux)
109 : {
110 : os << "mass=" << flux.mass << ", ";
111 : os << "Mn=" << flux.momentum_normal << ", ";
112 : os << "Mt=" << flux.momentum_tangent << ", ";
113 : os << "E=" << flux.energy << ", ";
114 : // do stuff
115 : return os;
116 : }
117 : void operator += (const Flux &a)
118 : {
119 : mass += a.mass;
120 : momentum_normal += a.momentum_normal; momentum_tangent += a.momentum_tangent;
121 : energy += a.energy;
122 : }
123 : void operator -= (const Flux &a)
124 : {
125 : mass -= a.mass;
126 : momentum_normal -= a.momentum_normal; momentum_tangent -= a.momentum_tangent;
127 : energy -= a.energy;
128 : }
129 : void operator *= (const Set::Scalar alpha)
130 : {
131 : mass *= alpha;
132 : momentum_normal *= alpha;
133 : momentum_tangent *= alpha;
134 : energy *= alpha;
135 : }
136 : void operator /= (const Set::Scalar alpha)
137 : {
138 : mass /= alpha;
139 : momentum_normal /= alpha;
140 : momentum_tangent /= alpha;
141 : energy /= alpha;
142 : }
143 16306448 : friend Flux operator + (const Flux &a, const Flux &b)
144 : {
145 16306448 : return Flux(a.mass+b.mass, a.momentum_normal+b.momentum_normal,
146 16306448 : a.momentum_tangent+b.momentum_tangent, a.energy+b.energy);
147 : }
148 7904456 : friend Flux operator - (const Flux &a, const Flux &b)
149 : {
150 7904456 : return Flux(a.mass-b.mass, a.momentum_normal-b.momentum_normal,
151 7904456 : a.momentum_tangent-b.momentum_tangent, a.energy-b.energy);
152 : }
153 52940800 : friend Flux operator * (const Flux &a, const Set::Scalar beta)
154 : {
155 52940800 : return Flux(a.mass*beta, a.momentum_normal*beta,
156 52940800 : a.momentum_tangent*beta, a.energy*beta);
157 : }
158 23713368 : friend Flux operator * (const Set::Scalar beta,const Flux &a)
159 : {
160 23713368 : return Flux(a.mass*beta, a.momentum_normal*beta,
161 23713368 : a.momentum_tangent*beta, a.energy*beta);
162 : }
163 7904456 : friend Flux operator / (const Flux &a, const Set::Scalar beta)
164 : {
165 7904456 : return Flux(a.mass/beta, a.momentum_normal/beta,
166 7904456 : a.momentum_tangent/beta, a.energy/beta);
167 : }
168 :
169 : };
170 :
171 :
172 : class Riemann
173 : {
174 : public:
175 : virtual Flux Solve(State lo, State hi, Model::Gas::Gas& gas, Set::Patch<const Set::Scalar>& X, int i, int j, int k, int side, Set::Scalar small) = 0;
176 8 : virtual ~Riemann() = default;
177 : };
178 :
179 :
180 : }
181 : }
182 : }
183 :
184 :
185 :
186 :
187 : #endif
|
