Line data Source code
1 : //* This file is part of the MOOSE framework
2 : //* https://mooseframework.inl.gov
3 : //*
4 : //* All rights reserved, see COPYRIGHT for full restrictions
5 : //* https://github.com/idaholab/moose/blob/master/COPYRIGHT
6 : //*
7 : //* Licensed under LGPL 2.1, please see LICENSE for details
8 : //* https://www.gnu.org/licenses/lgpl-2.1.html
9 :
10 : #include "Numerics.h"
11 : #include "MooseUtils.h"
12 : #include "ADReal.h"
13 :
14 : namespace THM
15 : {
16 :
17 : bool
18 30110 : absoluteFuzzyEqualVectors(const RealVectorValue & a, const RealVectorValue & b, const Real & tol)
19 : {
20 30106 : return MooseUtils::absoluteFuzzyEqual(a(0), b(0), tol) &&
21 60213 : MooseUtils::absoluteFuzzyEqual(a(1), b(1), tol) &&
22 30110 : MooseUtils::absoluteFuzzyEqual(a(2), b(2), tol);
23 : }
24 :
25 : bool
26 30108 : areParallelVectors(const RealVectorValue & a, const RealVectorValue & b, const Real & tol)
27 : {
28 30108 : const RealVectorValue c = a.cross(b);
29 30108 : return absoluteFuzzyEqualVectors(c, RealVectorValue(0, 0, 0), tol);
30 : }
31 :
32 : bool
33 2 : haveSameDirection(const RealVectorValue & a, const RealVectorValue & b, const Real & tol)
34 : {
35 2 : return areParallelVectors(a, b, tol) && a * b > 0;
36 : }
37 :
38 : Real
39 0 : applyQuotientRule(const Real & num, const Real & den, const Real & dnum_dy, const Real & dden_dy)
40 : {
41 0 : return (dnum_dy * den - num * dden_dy) / (den * den);
42 : }
43 :
44 : DenseVector<Real>
45 0 : applyQuotientRule(const Real & num,
46 : const Real & den,
47 : const DenseVector<Real> & dnum_dy,
48 : const DenseVector<Real> & dden_dy)
49 : {
50 : DenseVector<Real> d_dy = dnum_dy;
51 0 : d_dy *= 1.0 / den;
52 0 : d_dy.add(-num / std::pow(den, 2), dden_dy);
53 :
54 0 : return d_dy;
55 : }
56 :
57 : void
58 991 : vel_from_arhoA_arhouA(Real arhoA, Real arhouA, Real & vel, Real & dvel_darhoA, Real & dvel_darhouA)
59 : {
60 991 : vel = arhouA / arhoA;
61 991 : dvel_darhoA = -arhouA / (arhoA * arhoA);
62 991 : dvel_darhouA = 1.0 / arhoA;
63 991 : }
64 :
65 : ADReal
66 0 : vel_from_arhoA_arhouA(ADReal arhoA, ADReal arhouA)
67 : {
68 0 : return arhouA / arhoA;
69 : }
70 :
71 : Real
72 1 : dvel_darhoA(Real arhoA, Real arhouA)
73 : {
74 1 : return -arhouA / (arhoA * arhoA);
75 : }
76 :
77 : Real
78 1 : dvel_darhouA(Real arhoA)
79 : {
80 1 : return 1.0 / arhoA;
81 : }
82 :
83 : void
84 1 : rho_from_arhoA_alpha_A(
85 : Real arhoA, Real alpha, Real A, Real & rho, Real & drho_darhoA, Real & drho_dalpha)
86 : {
87 1 : rho = arhoA / (alpha * A);
88 :
89 1 : drho_darhoA = 1.0 / (alpha * A);
90 1 : drho_dalpha = -arhoA / (alpha * alpha * A);
91 1 : }
92 :
93 : ADReal
94 0 : rho_from_arhoA_alpha_A(ADReal arhoA, ADReal alpha, ADReal A)
95 : {
96 0 : return arhoA / (alpha * A);
97 : }
98 :
99 : void
100 145567 : v_from_rhoA_A(Real rhoA, Real A, Real & v, Real & dv_drhoA)
101 : {
102 145567 : v = A / rhoA;
103 145567 : dv_drhoA = -A / (rhoA * rhoA);
104 145567 : }
105 :
106 : ADReal
107 26122 : v_from_rhoA_A(ADReal rhoA, ADReal A)
108 : {
109 26122 : return A / rhoA;
110 : }
111 :
112 : void
113 1 : v_from_arhoA_alpha_A(Real arhoA, Real alpha, Real A, Real & v, Real & dv_darhoA, Real & dv_dalpha)
114 : {
115 1 : v = (alpha * A) / arhoA;
116 1 : dv_darhoA = -(alpha * A) / (arhoA * arhoA);
117 1 : dv_dalpha = A / arhoA;
118 1 : }
119 :
120 : ADReal
121 0 : v_from_arhoA_alpha_A(ADReal arhoA, ADReal alpha, ADReal A)
122 : {
123 0 : return (alpha * A) / arhoA;
124 : }
125 :
126 : void
127 1 : v_from_rho(Real rho, Real & v, Real & dv_drho)
128 : {
129 1 : v = 1.0 / rho;
130 1 : dv_drho = -1.0 / (rho * rho);
131 1 : }
132 :
133 : Real
134 2 : dv_dalpha_liquid(Real area, Real arhoA, bool is_liquid)
135 : {
136 2 : const Real sign = is_liquid ? 1.0 : -1.0;
137 2 : return sign * (area / arhoA);
138 : }
139 :
140 : Real
141 991 : dv_darhoA(Real area, Real arhoA)
142 : {
143 991 : return -area / arhoA / arhoA;
144 : }
145 :
146 : void
147 1 : e_from_arhoA_arhouA_arhoEA(Real arhoA,
148 : Real arhouA,
149 : Real arhoEA,
150 : Real & e,
151 : Real & de_darhoA_val,
152 : Real & de_darhouA_val,
153 : Real & de_darhoEA_val)
154 : {
155 1 : e = arhoEA / arhoA - 0.5 * arhouA * arhouA / (arhoA * arhoA);
156 1 : de_darhoA_val = de_darhoA(arhoA, arhouA, arhoEA);
157 1 : de_darhouA_val = de_darhouA(arhoA, arhouA);
158 1 : de_darhoEA_val = de_darhoEA(arhoA);
159 1 : }
160 :
161 : ADReal
162 10679 : e_from_arhoA_arhouA_arhoEA(ADReal arhoA, ADReal arhouA, ADReal arhoEA)
163 : {
164 10679 : return arhoEA / arhoA - 0.5 * arhouA * arhouA / (arhoA * arhoA);
165 : }
166 :
167 : void
168 1 : e_from_E_vel(Real E, Real vel, Real & e, Real & de_dE, Real & de_dvel)
169 : {
170 1 : e = E - 0.5 * vel * vel;
171 1 : de_dE = 1;
172 1 : de_dvel = -vel;
173 1 : }
174 : ADReal
175 0 : e_from_E_vel(ADReal E, ADReal vel)
176 : {
177 0 : return E - 0.5 * vel * vel;
178 : }
179 :
180 : Real
181 992 : de_darhoA(Real arhoA, Real arhouA, Real arhoEA)
182 : {
183 992 : return (-arhoEA / arhoA / arhoA + arhouA * arhouA / arhoA / arhoA / arhoA);
184 : }
185 :
186 : Real
187 992 : de_darhouA(Real arhoA, Real arhouA)
188 : {
189 992 : return (-arhouA / arhoA / arhoA);
190 : }
191 :
192 : Real
193 992 : de_darhoEA(Real arhoA)
194 : {
195 992 : return (1 / arhoA);
196 : }
197 :
198 : void
199 1 : E_from_arhoA_arhoEA(Real arhoA, Real arhoEA, Real & E, Real & dE_darhoA, Real & dE_darhoEA)
200 : {
201 1 : E = arhoEA / arhoA;
202 1 : dE_darhoA = -arhoEA / (arhoA * arhoA);
203 1 : dE_darhoEA = 1.0 / arhoA;
204 1 : }
205 : ADReal
206 0 : E_from_arhoA_arhoEA(ADReal arhoA, ADReal arhoEA)
207 : {
208 0 : return arhoEA / arhoA;
209 : }
210 :
211 : void
212 1 : E_from_e_vel(Real e, Real vel, Real & E, Real & dE_de, Real & dE_dvel)
213 : {
214 1 : E = e + 0.5 * vel * vel;
215 1 : dE_de = 1.0;
216 1 : dE_dvel = vel;
217 1 : }
218 :
219 : void
220 1 : h_from_e_p_rho(Real e, Real p, Real rho, Real & h, Real & dh_de, Real & dh_dp, Real & dh_drho)
221 : {
222 1 : h = e + p / rho;
223 1 : dh_de = 1.0;
224 1 : dh_dp = 1.0 / rho;
225 1 : dh_drho = -p / (rho * rho);
226 1 : }
227 :
228 : ADReal
229 10563 : h_from_e_p_rho(ADReal e, ADReal p, ADReal rho)
230 : {
231 10563 : return e + p / rho;
232 : }
233 :
234 : bool
235 65132 : isInlet(Real vel, Real normal)
236 : {
237 65132 : return (vel * normal) < 0;
238 : }
239 :
240 : bool
241 16361 : isInlet(ADReal vel, Real normal)
242 : {
243 16361 : return (MetaPhysicL::raw_value(vel) * normal) < 0;
244 : }
245 :
246 : bool
247 4 : isOutlet(Real vel, Real normal)
248 : {
249 4 : return (vel * normal) >= 0;
250 : }
251 :
252 : bool
253 21108 : isOutlet(ADReal vel, Real normal)
254 : {
255 21108 : return (MetaPhysicL::raw_value(vel) * normal) >= 0;
256 : }
257 :
258 : } // namespace THM
|