12 namespace EquationOfStateLibrary
29 this->add_parameter(
"A", capA,
"The A constant");
32 this->add_parameter(
"B", capB,
"The B constant");
35 this->add_parameter(
"R1", R1,
"The R1 constant");
38 this->add_parameter(
"R2", R2,
"The R2 constant");
41 this->add_parameter(
"omega", omega,
"The Gruneisen coefficient");
44 this->add_parameter(
"rho_0", rho_0,
"The reference density");
47 this->add_parameter(
"q_0", q_0,
"The specific internal energy offset");
51 "c_v", cv_,
"The specific heat capacity at constant volume");
62 double pressure(
double rho,
double e)
const final
65 const auto ratio = rho / rho_0;
67 const auto first_term =
68 capA * (1. - omega / R1 * ratio) * std::exp(-R1 * 1. / ratio);
69 const auto second_term =
70 capB * (1. - omega / R2 * ratio) * std::exp(-R2 * 1. / ratio);
72 return first_term + second_term + omega * rho * (e + q_0);
85 const auto ratio = rho / rho_0;
87 const auto first_term =
88 capA * (1. - omega / R1 * ratio) * std::exp(-R1 * 1. / ratio);
89 const auto second_term =
90 capB * (1. - omega / R2 * ratio) * std::exp(-R2 * 1. / ratio);
92 return (p - first_term - second_term) / (rho * omega);
105 const auto ratio = rho / rho_0;
107 const auto first_term = capA / R1 * std::exp(-R1 * 1. / ratio);
108 const auto second_term = capB / R2 * std::exp(-R2 * 1. / ratio);
110 return (e + q_0 - 1. / rho_0 * (first_term + second_term)) / cv_;
120 const auto t1 = omega * rho / (R1 * rho_0);
121 const auto factor1 = omega * (1. - t1) * (1. + 1. / t1) - t1;
122 const auto first_term =
123 capA / rho * factor1 * std::exp(-1. / t1 / omega);
125 const auto t2 = omega * rho / (R2 * rho_0);
126 const auto factor2 = omega * (1. - t2) * (1. + 1. / t2) - t2;
127 const auto second_term =
128 capB / rho * factor2 * std::exp(-1. / t2 / omega);
130 const auto third_term = omega * (omega + 1.) * e;
132 return std::sqrt(first_term + second_term + third_term);
double pressure(double rho, double e) const final
JonesWilkinsLee(const std::string &subsection)
double specific_internal_energy(double rho, double p) const final
double speed_of_sound(double rho, double e) const final
double temperature(double rho, double e) const final