Line data Source code
1 : #ifndef SOLVER_LOCAL_RIEMANN_H
2 : #define SOLVER_LOCAL_RIEMANN_H
3 :
4 :
5 : namespace Solver
6 : {
7 : namespace Local
8 : {
9 : namespace Riemann
10 : {
11 : struct State {
12 : Set::Scalar rho = NAN;
13 : Set::Scalar M_normal = NAN;
14 : Set::Scalar M_tangent = NAN;
15 : Set::Scalar E = NAN;
16 : // Construtor for convenience
17 : State() { rho = 0.0; M_normal = 0.0, M_tangent = 0.0; E = 0.0;}
18 112896024 : State(Set::Scalar a_rho, Set::Scalar a_M_normal, Set::Scalar a_M_tangent, Set::Scalar a_E)
19 112896024 : : rho(a_rho), M_normal(a_M_normal), M_tangent(a_M_tangent), E(a_E) {}
20 75264000 : 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)
21 75264000 : {
22 75264000 : rho = density_mf(i,j,k);
23 75264000 : if (direction == 0)
24 : {
25 37632000 : M_normal = momentum_mf(i,j,k,0);
26 37632000 : M_tangent = momentum_mf(i,j,k,1);
27 : }
28 37632000 : else if (direction == 1)
29 : {
30 37632000 : M_normal = momentum_mf(i,j,k,1);
31 37632000 : M_tangent = momentum_mf(i,j,k,0);
32 : }
33 : else
34 : {
35 0 : Util::Abort(INFO, "Not supported yet");
36 : }
37 75264000 : E = energy_mf(i,j,k);
38 75264000 : }
39 :
40 :
41 0 : friend std::ostream &operator<<(std::ostream &os, const State &state)
42 : {
43 0 : os << "rho=" << state.rho << ", ";
44 0 : os << "Mn=" << state.M_normal << ", ";
45 0 : os << "Mt=" << state.M_tangent << ", ";
46 0 : os << "E=" << state.E << ", ";
47 : // do stuf
48 0 : return os;
49 : }
50 : void operator += (const State &a)
51 : {
52 : rho += a.rho; M_normal += a.M_normal; M_tangent += a.M_tangent; E += a.E;
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 Set::Scalar alpha)
59 : {
60 : rho *= alpha; M_normal *= alpha; M_tangent *= alpha; E *= alpha;
61 : };
62 : void operator /= (const Set::Scalar alpha)
63 : {
64 : rho /= alpha; M_normal /= alpha; M_tangent /= alpha; E /= alpha;
65 : };
66 : friend State operator + (const State &a, const State &b)
67 : {
68 : return State(a.rho+b.rho, a.M_normal+b.M_normal, a.M_tangent+b.M_tangent, a.E+b.E);
69 : }
70 37632000 : friend State operator - (const State &a, const State &b)
71 : {
72 37632000 : return State(a.rho-b.rho, a.M_normal-b.M_normal, a.M_tangent-b.M_tangent, a.E-b.E);
73 : }
74 37632000 : friend State operator * (const Set::Scalar alpha, const State &b)
75 : {
76 37632000 : return State(b.rho*alpha, b.M_normal*alpha, b.M_tangent*alpha, b.E*alpha);
77 : }
78 : friend State operator * (const State &b, const Set::Scalar alpha)
79 : {
80 : return State(b.rho*alpha, b.M_normal*alpha, b.M_tangent*alpha, b.E*alpha);
81 : }
82 37632000 : friend State operator / (const State &b, const Set::Scalar alpha)
83 : {
84 37632000 : return State(b.rho/alpha, b.M_normal/alpha, b.M_tangent/alpha, b.E/alpha);
85 : }
86 : };
87 :
88 : struct Flux {
89 : Set::Scalar mass;
90 : Set::Scalar momentum_normal;
91 : Set::Scalar momentum_tangent;
92 : Set::Scalar energy;
93 50176016 : Flux() : mass(0.0), momentum_normal(0.0), momentum_tangent(0.0), energy(0.0) {}
94 25088000 : Flux(Set::Scalar a_mass, Set::Scalar a_momentum_normal, Set::Scalar a_momentum_tangent, Set::Scalar a_energy) :
95 25088000 : mass(a_mass), momentum_normal(a_momentum_normal),
96 25088000 : momentum_tangent(a_momentum_tangent), energy(a_energy) {}
97 0 : friend std::ostream &operator<<(std::ostream &os, const Flux &flux)
98 : {
99 0 : os << "mass=" << flux.mass << ", ";
100 0 : os << "Mn=" << flux.momentum_normal << ", ";
101 0 : os << "Mt=" << flux.momentum_tangent << ", ";
102 0 : os << "E=" << flux.energy << ", ";
103 : // do stuff
104 0 : return os;
105 : }
106 : void operator += (const Flux &a)
107 : {
108 : mass += a.mass;
109 : momentum_normal += a.momentum_normal; momentum_tangent += a.momentum_tangent;
110 : energy += a.energy;
111 : }
112 : void operator -= (const Flux &a)
113 : {
114 : mass -= a.mass;
115 : momentum_normal -= a.momentum_normal; momentum_tangent -= a.momentum_tangent;
116 : energy -= a.energy;
117 : }
118 : void operator *= (const Set::Scalar alpha)
119 : {
120 : mass *= alpha;
121 : momentum_normal *= alpha;
122 : momentum_tangent *= alpha;
123 : energy *= alpha;
124 : }
125 : void operator /= (const Set::Scalar alpha)
126 : {
127 : mass /= alpha;
128 : momentum_normal /= alpha;
129 : momentum_tangent /= alpha;
130 : energy /= alpha;
131 : }
132 : friend Flux operator + (const Flux &a, const Flux &b)
133 : {
134 : return Flux(a.mass+b.mass, a.momentum_normal+b.momentum_normal,
135 : a.momentum_tangent+b.momentum_tangent, a.energy+b.energy);
136 : }
137 : friend Flux operator - (const Flux &a, const Flux &b)
138 : {
139 : return Flux(a.mass-b.mass, a.momentum_normal-b.momentum_normal,
140 : a.momentum_tangent-b.momentum_tangent, a.energy-b.energy);
141 : }
142 25088000 : friend Flux operator * (const Flux &a, const Set::Scalar beta)
143 : {
144 25088000 : return Flux(a.mass*beta, a.momentum_normal*beta,
145 25088000 : a.momentum_tangent*beta, a.energy*beta);
146 : }
147 : friend Flux operator * (const Set::Scalar beta,const Flux &a)
148 : {
149 : return Flux(a.mass*beta, a.momentum_normal*beta,
150 : a.momentum_tangent*beta, a.energy*beta);
151 : }
152 : friend Flux operator / (const Flux &a, const Set::Scalar beta)
153 : {
154 : return Flux(a.mass/beta, a.momentum_normal/beta,
155 : a.momentum_tangent/beta, a.energy/beta);
156 : }
157 :
158 : };
159 :
160 :
161 :
162 : }
163 : }
164 : }
165 :
166 :
167 :
168 :
169 : #endif
|
