https://mooseframework.inl.gov
Functions
NumericsTest.C File Reference

Go to the source code of this file.

Functions

 TEST (NumericsTest, test_absoluteFuzzyEqualVectors_True)
 
 TEST (NumericsTest, test_absoluteFuzzyEqualVectors_False)
 
 TEST (NumericsTest, test_areParallelVectors_True)
 
 TEST (NumericsTest, test_areParallelVectors_False)
 
 TEST (NumericsTest, test_haveSameDirection_True)
 
 TEST (NumericsTest, test_haveSameDirection_False)
 
 TEST (NumericsTest, Reynolds)
 
 TEST (NumericsTest, Prandtl)
 
 TEST (NumericsTest, Peclet)
 
 TEST (NumericsTest, Grashof)
 
 TEST (NumericsTest, Laplace)
 
 TEST (NumericsTest, viscosityNumber)
 
 TEST (NumericsTest, wallHeatTransferCoefficient)
 
 TEST (NumericsTest, Dean)
 
 TEST (NumericsTest, vel)
 
 TEST (NumericsTest, rho)
 
 TEST (NumericsTest, specific_volume)
 
 TEST (NumericsTest, internal_energy)
 
 TEST (NumericsTest, total_energy)
 
 TEST (NumericsTest, specific_enthalpy)
 
 TEST (NumericsTest, xlets)
 

Function Documentation

◆ TEST() [1/21]

TEST ( NumericsTest  ,
test_absoluteFuzzyEqualVectors_True   
)

Definition at line 14 of file NumericsTest.C.

15 {
16  const RealVectorValue a(1.0, 2.0, 3.0);
17  const RealVectorValue b(1.0 + 1.0e-15, 2.0, 3.0);
18  EXPECT_TRUE(THM::absoluteFuzzyEqualVectors(a, b));
19 }
bool absoluteFuzzyEqualVectors(const RealVectorValue &a, const RealVectorValue &b, const Real &tol=libMesh::TOLERANCE *libMesh::TOLERANCE)
Tests if two real-valued vectors are equal within some absolute tolerance.
Definition: Numerics.C:18

◆ TEST() [2/21]

TEST ( NumericsTest  ,
test_absoluteFuzzyEqualVectors_False   
)

Definition at line 21 of file NumericsTest.C.

22 {
23  const RealVectorValue a(1.0, 2.0, 3.0);
24  const RealVectorValue b(1.1, 2.0, 3.0);
25  EXPECT_FALSE(THM::absoluteFuzzyEqualVectors(a, b));
26 }
bool absoluteFuzzyEqualVectors(const RealVectorValue &a, const RealVectorValue &b, const Real &tol=libMesh::TOLERANCE *libMesh::TOLERANCE)
Tests if two real-valued vectors are equal within some absolute tolerance.
Definition: Numerics.C:18

◆ TEST() [3/21]

TEST ( NumericsTest  ,
test_areParallelVectors_True   
)

Definition at line 28 of file NumericsTest.C.

29 {
30  const RealVectorValue a(1.0, 2.0, 3.0);
31  const RealVectorValue b(-2.0, -4.0, -6.0);
32  EXPECT_TRUE(THM::areParallelVectors(a, b));
33 }
bool areParallelVectors(const RealVectorValue &a, const RealVectorValue &b, const Real &tol=libMesh::TOLERANCE *libMesh::TOLERANCE)
Tests if two real-valued vectors are parallel within some absolute tolerance.
Definition: Numerics.C:26

◆ TEST() [4/21]

TEST ( NumericsTest  ,
test_areParallelVectors_False   
)

Definition at line 35 of file NumericsTest.C.

36 {
37  const RealVectorValue a(1.0, 2.0, 3.0);
38  const RealVectorValue b(2.0, 3.0, 1.0);
39  EXPECT_FALSE(THM::areParallelVectors(a, b));
40 }
bool areParallelVectors(const RealVectorValue &a, const RealVectorValue &b, const Real &tol=libMesh::TOLERANCE *libMesh::TOLERANCE)
Tests if two real-valued vectors are parallel within some absolute tolerance.
Definition: Numerics.C:26

◆ TEST() [5/21]

TEST ( NumericsTest  ,
test_haveSameDirection_True   
)

Definition at line 42 of file NumericsTest.C.

43 {
44  const RealVectorValue a(1.0, 2.0, 3.0);
45  const RealVectorValue b(2.0, 4.0, 6.0);
46  EXPECT_TRUE(THM::haveSameDirection(a, b));
47 }
bool haveSameDirection(const RealVectorValue &a, const RealVectorValue &b, const Real &tol=libMesh::TOLERANCE *libMesh::TOLERANCE)
Tests if two real-valued vectors are in the same direction.
Definition: Numerics.C:33

◆ TEST() [6/21]

TEST ( NumericsTest  ,
test_haveSameDirection_False   
)

Definition at line 49 of file NumericsTest.C.

50 {
51  const RealVectorValue a(1.0, 2.0, 3.0);
52  const RealVectorValue b(-1.0, 2.0, 3.0);
53  EXPECT_FALSE(THM::haveSameDirection(a, b));
54 }
bool haveSameDirection(const RealVectorValue &a, const RealVectorValue &b, const Real &tol=libMesh::TOLERANCE *libMesh::TOLERANCE)
Tests if two real-valued vectors are in the same direction.
Definition: Numerics.C:33

◆ TEST() [7/21]

TEST ( NumericsTest  ,
Reynolds   
)

Definition at line 56 of file NumericsTest.C.

57 {
58  ABS_TEST(THM::Reynolds(0.1, 999., 0.5, 2.0e-2, 0.9), 1.11, 1e-13);
59  ABS_TEST(THM::Reynolds(0.1, 999., -0.5, 2.0e-2, 0.9), 1.11, 1e-13);
60 }
auto Reynolds(const T1 &volume_fraction, const T2 &rho, const T3 &vel, const T4 &D_h, const T5 &mu)
Compute Reynolds number.
Definition: Numerics.h:118

◆ TEST() [8/21]

TEST ( NumericsTest  ,
Prandtl   
)

Definition at line 62 of file NumericsTest.C.

62 { ABS_TEST(THM::Prandtl(10., 0.1, 2.), 0.5, 1e-13); }
auto Prandtl(const T1 &cp, const T2 &mu, const T3 &k)
Compute Prandtl number.
Definition: Numerics.h:133

◆ TEST() [9/21]

TEST ( NumericsTest  ,
Peclet   
)

Definition at line 64 of file NumericsTest.C.

65 {
66  ABS_TEST(THM::Peclet(0.1, 999., 999., 0.5, 2.0e-2, 0.9), 1108.89, 1e-13);
67 }
auto Peclet(const T1 &volume_fraction, const T2 &cp, const T3 &rho, const T4 &vel, const T5 &D_h, const T6 &k)
Compute Peclet number.
Definition: Numerics.h:153

◆ TEST() [10/21]

TEST ( NumericsTest  ,
Grashof   
)

Definition at line 69 of file NumericsTest.C.

70 {
71  ABS_TEST(THM::Grashof(0.1, 1., 2.0e-2, 999., 0.05, 9.81), 3.1329247392e3, 1e-13);
72 }
auto Grashof(const T1 &beta, const T2 &dT, const T3 &D_h, const T4 &rho_liquid, const T5 &mu_liquid, const Real &gravity_magnitude)
Compute Grashof number.
Definition: Numerics.h:178

◆ TEST() [11/21]

TEST ( NumericsTest  ,
Laplace   
)

Definition at line 74 of file NumericsTest.C.

74 { ABS_TEST(THM::Laplace(0.001, 1., 9.81), 0.010096375546923, 1e-13); }
auto Laplace(const T1 &surf_tension, const T2 &delta_rho, const Real &gravity_magnitude)
Compute Laplace number (or coefficient)
Definition: Numerics.h:200

◆ TEST() [12/21]

TEST ( NumericsTest  ,
viscosityNumber   
)

Definition at line 76 of file NumericsTest.C.

77 {
78  ABS_TEST(THM::viscosityNumber(0.05, 0.02, 999., 2., 9.81), 0.062602188259, 1e-13);
79 }
auto viscosityNumber(const T1 &viscosity, const T2 &surf_tension, const T3 &rho_k, const T4 &delta_rho, const Real &gravity_magnitude)
Compute viscosity number (or coefficient)
Definition: Numerics.h:218

◆ TEST() [13/21]

TEST ( NumericsTest  ,
wallHeatTransferCoefficient   
)

Definition at line 81 of file NumericsTest.C.

82 {
83  ABS_TEST(THM::wallHeatTransferCoefficient(2., 6., 3.), 4, 1e-13);
84 }
auto wallHeatTransferCoefficient(const T1 &Nu, const T2 &k, const T3 &D_h)
Compute wall heat transfer coefficient.

◆ TEST() [14/21]

TEST ( NumericsTest  ,
Dean   
)

Definition at line 86 of file NumericsTest.C.

86 { ABS_TEST(THM::Dean(1., 4.), 2, 1e-15); }
auto Dean(const T1 &Re, const T2 &doD)
Compute Dean number.
Definition: Numerics.h:239

◆ TEST() [15/21]

TEST ( NumericsTest  ,
vel   
)

Definition at line 88 of file NumericsTest.C.

89 {
90  const Real arhoA = 1;
91  const Real arhouA = 2;
92 
95  ABS_TEST(vel, 2., 1e-15);
96  ABS_TEST(dvel_darhoA, -2, 1e-15);
97  ABS_TEST(dvel_darhouA, 1., 1e-15);
98 
99  ABS_TEST(THM::dvel_darhoA(arhoA, arhouA), -2., 1e-15);
100  ABS_TEST(THM::dvel_darhouA(arhoA), 1., 1e-15);
101 }
void vel_from_arhoA_arhouA(Real arhoA, Real arhouA, Real &vel, Real &dvel_darhoA, Real &dvel_darhouA)
Computes velocity and its derivatives from alpha*rho*A and alpha*rho*u*A.
Definition: Numerics.C:58
Real dvel_darhouA(Real arhoA)
Derivative of velocity w.r.t.
Definition: Numerics.C:78
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
Real dvel_darhoA(Real arhoA, Real arhouA)
Derivative of velocity w.r.t.
Definition: Numerics.C:72

◆ TEST() [16/21]

TEST ( NumericsTest  ,
rho   
)

Definition at line 103 of file NumericsTest.C.

104 {
105  const Real arhoA = 2;
106  const Real alpha = 0.1;
107  const Real A = 1;
108 
109  Real rho, drho_darhoA, drho_dalpha;
110  THM::rho_from_arhoA_alpha_A(arhoA, alpha, A, rho, drho_darhoA, drho_dalpha);
111  ABS_TEST(rho, 20., 1e-15);
112  ABS_TEST(drho_darhoA, 10., 1e-15);
113  ABS_TEST(drho_dalpha, -200., 1e-13);
114 }
void rho_from_arhoA_alpha_A(Real arhoA, Real alpha, Real A, Real &rho, Real &drho_darhoA, Real &drho_dalpha)
Computes density and its derivatives from alpha*rho*A, alpha, and area.
Definition: Numerics.C:84
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
static const std::string alpha
Definition: NS.h:134

◆ TEST() [17/21]

TEST ( NumericsTest  ,
specific_volume   
)

Definition at line 116 of file NumericsTest.C.

117 {
118  const Real arhoA = 2;
119  const Real A = 1;
120 
121  Real v, dv_drhoA;
122  THM::v_from_rhoA_A(arhoA, A, v, dv_drhoA);
123  ABS_TEST(v, 0.5, 1e-15);
124  ABS_TEST(dv_drhoA, -0.25, 1e-15);
125 
126  const Real alpha = 0.1;
127  Real dv_darhoA, dv_dalpha;
128  THM::v_from_arhoA_alpha_A(arhoA, alpha, A, v, dv_darhoA, dv_dalpha);
129  ABS_TEST(v, 0.05, 1e-15);
130  ABS_TEST(dv_darhoA, -0.025, 1e-15);
131  ABS_TEST(dv_dalpha, 0.5, 1e-13);
132 
133  const Real rho = 20;
134  Real dv_drho;
135  THM::v_from_rho(rho, v, dv_drho);
136  ABS_TEST(v, 0.05, 1e-15);
137  ABS_TEST(dv_darhoA, -0.025, 1e-15);
138 
139  ABS_TEST(THM::dv_dalpha_liquid(A, arhoA, true), 0.5, 1e-15);
140  ABS_TEST(THM::dv_dalpha_liquid(A, arhoA, false), -0.5, 1e-15);
141 
142  ABS_TEST(THM::dv_darhoA(A, arhoA), -0.25, 1e-15);
143 }
void v_from_rhoA_A(Real rhoA, Real A, Real &v, Real &dv_drhoA)
Computes specific volume and its derivatives from rho*A, and area.
Definition: Numerics.C:100
Real dv_darhoA(Real area, Real arhoA)
Derivative of specific volume wrt density equation solution variable.
Definition: Numerics.C:141
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
static const std::string v
Definition: NS.h:84
static const std::string alpha
Definition: NS.h:134
void v_from_rho(Real rho, Real &v, Real &dv_drho)
Computes specific volume and its derivative with respect to density.
Definition: Numerics.C:127
Real dv_dalpha_liquid(Real area, Real arhoA, bool is_liquid)
Derivative of specific volume wrt alpha_liquid.
Definition: Numerics.C:134
void v_from_arhoA_alpha_A(Real arhoA, Real alpha, Real A, Real &v, Real &dv_darhoA, Real &dv_dalpha)
Computes specific volume and its derivatives from alpha*rho*A, volume fraction, and area...
Definition: Numerics.C:113

◆ TEST() [18/21]

TEST ( NumericsTest  ,
internal_energy   
)

Definition at line 145 of file NumericsTest.C.

146 {
147  const Real arhoA = 2;
148  const Real arhouA = 4;
149  const Real arhoEA = 8;
150 
152  THM::e_from_arhoA_arhouA_arhoEA(arhoA, arhouA, arhoEA, e, de_darhoA, de_darhouA, de_darhoEA);
153  ABS_TEST(e, 2., 1e-15);
154  ABS_TEST(de_darhoA, 0., 1e-15);
155  ABS_TEST(de_darhouA, -1., 1e-15);
156  ABS_TEST(de_darhoEA, 0.5, 1e-15);
157 
158  const Real E = 10;
159  const Real vel = 2;
160  Real de_dE, de_dvel;
161  THM::e_from_E_vel(E, vel, e, de_dE, de_dvel);
162  ABS_TEST(e, 8., 1e-15);
163  ABS_TEST(de_dE, 1., 1e-15);
164  ABS_TEST(de_dvel, -2., 1e-15);
165 
166  ABS_TEST(THM::de_darhoA(arhoA, arhouA, arhoEA), 0., 1e-15);
167  ABS_TEST(THM::de_darhouA(arhoA, arhouA), -1., 1e-15);
168  ABS_TEST(THM::de_darhoEA(arhoA), 0.5, 1e-15);
169 }
void e_from_arhoA_arhouA_arhoEA(Real arhoA, Real arhouA, Real arhoEA, Real &e, Real &de_darhoA, Real &de_darhouA, Real &de_darhoEA)
Computes specific internal energy and its derivatives from alpha*rho*A, alpha*rho*u*A, and alpha*rho*E*A.
Definition: Numerics.C:147
Real de_darhoA(Real arhoA, Real arhouA, Real arhoEA)
Derivative of specific internal energy wrt density of the phase (rhoA or arhoA)
Definition: Numerics.C:181
Real de_darhoEA(Real arhoA)
Derivative of specific internal energy wrt total energy of the phase (rhoEA or arhoEA) ...
Definition: Numerics.C:193
void e_from_E_vel(Real E, Real vel, Real &e, Real &de_dE, Real &de_dvel)
Computes specific internal energy and its derivatives from specific total energy and velocity...
Definition: Numerics.C:168
Real de_darhouA(Real arhoA, Real arhouA)
Derivative of specific internal energy wrt momentum of the phase (rhouA or arhouA) ...
Definition: Numerics.C:187
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real

◆ TEST() [19/21]

TEST ( NumericsTest  ,
total_energy   
)

Definition at line 171 of file NumericsTest.C.

172 {
173  const Real arhoA = 2.;
174  const Real arhoEA = 8.;
175 
176  Real E, dE_darhoA, dE_darhoEA;
177  THM::E_from_arhoA_arhoEA(arhoA, arhoEA, E, dE_darhoA, dE_darhoEA);
178  ABS_TEST(E, 4., 1e-15);
179  ABS_TEST(dE_darhoA, -2., 1e-15);
180  ABS_TEST(dE_darhoEA, 0.5, 1e-15);
181 
182  const Real e = 4.;
183  const Real vel = 2.;
184 
185  Real dE_de, dE_dvel;
186  THM::E_from_e_vel(e, vel, E, dE_de, dE_dvel);
187 
188  ABS_TEST(E, 6., 1e-15);
189  ABS_TEST(dE_de, 1., 1e-15);
190  ABS_TEST(dE_dvel, 2., 1e-15);
191 }
void E_from_arhoA_arhoEA(Real arhoA, Real arhoEA, Real &E, Real &dE_darhoA, Real &dE_darhoEA)
Computes specific total energy and its derivatives from alpha*rho*A and alpha*rho*E*A.
Definition: Numerics.C:199
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
void E_from_e_vel(Real e, Real vel, Real &E, Real &dE_de, Real &dE_dvel)
Computes specific total energy and its derivatives from specific internal energy and velocity...
Definition: Numerics.C:212

◆ TEST() [20/21]

TEST ( NumericsTest  ,
specific_enthalpy   
)

Definition at line 193 of file NumericsTest.C.

194 {
195  const Real e = 6;
196  const Real p = 4;
197  const Real rho = 2;
198 
199  Real h, dh_de, dh_dp, dh_drho;
200  THM::h_from_e_p_rho(e, p, rho, h, dh_de, dh_dp, dh_drho);
201  ABS_TEST(h, 8., 1e-15);
202  ABS_TEST(dh_de, 1., 1e-15);
203  ABS_TEST(dh_dp, 0.5, 1e-15);
204  ABS_TEST(dh_drho, -1., 1e-15);
205 }
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
void h_from_e_p_rho(Real e, Real p, Real rho, Real &h, Real &dh_de, Real &dh_dp, Real &dh_drho)
Computes specific enthalpy and its derivatives from specific internal energy, pressure, and density.
Definition: Numerics.C:220

◆ TEST() [21/21]

TEST ( NumericsTest  ,
xlets   
)

Definition at line 207 of file NumericsTest.C.

208 {
209  EXPECT_FALSE(THM::isInlet(2, 1));
210  EXPECT_TRUE(THM::isInlet(2, -1));
211  EXPECT_TRUE(THM::isInlet(-2, 1));
212  EXPECT_FALSE(THM::isInlet(-2, -1));
213 
214  EXPECT_TRUE(THM::isOutlet(2, 1));
215  EXPECT_FALSE(THM::isOutlet(2, -1));
216  EXPECT_FALSE(THM::isOutlet(-2, 1));
217  EXPECT_TRUE(THM::isOutlet(-2, -1));
218 }
bool isOutlet(Real vel, Real normal)
Determine if outlet boundary condition should be applied.
Definition: Numerics.C:247
bool isInlet(Real vel, Real normal)
Determine if inlet boundary condition should be applied.
Definition: Numerics.C:235