diff --git a/EOS/breakout/actual_eos.H b/EOS/breakout/actual_eos.H index 0674931872..ac0ba24bca 100644 --- a/EOS/breakout/actual_eos.H +++ b/EOS/breakout/actual_eos.H @@ -8,8 +8,6 @@ #include #include -using namespace eos_rp; - const std::string eos_name = "breakout"; inline @@ -17,8 +15,8 @@ void actual_eos_init () { // constant ratio of specific heats - if (eos_gamma > 0.e0_rt) { - gamma_const = eos_gamma; + if (eos_rp::eos_gamma > 0.e0_rt) { + gamma_const = eos_rp::eos_gamma; } else { gamma_const = 5.0_rt / 3.0_rt; } diff --git a/EOS/gamma_law/actual_eos.H b/EOS/gamma_law/actual_eos.H index 41eba41423..174a870ef6 100644 --- a/EOS/gamma_law/actual_eos.H +++ b/EOS/gamma_law/actual_eos.H @@ -6,8 +6,6 @@ #include #include -using namespace eos_rp; - // This is a constant gamma equation of state, using an ideal gas. // // The gas may either be completely ionized or completely neutral. @@ -22,7 +20,7 @@ inline void actual_eos_init() { // constant ratio of specific heats - if (eos_gamma <= 0.0) { + if (eos_rp::eos_gamma <= 0.0) { amrex::Error("gamma_const cannot be < 0"); } @@ -56,7 +54,7 @@ void actual_eos (I input, T& state) const amrex::Real m_nucleon = C::m_u; if constexpr (has_xn::value) { - if (eos_assume_neutral) { + if (eos_rp::eos_assume_neutral) { state.mu = state.abar; } else { amrex::Real sum = 0.0; @@ -88,7 +86,7 @@ void actual_eos (I input, T& state) // h = e + p/rho = (p/rho)*[1 + 1/(gamma-1)] = (p/rho)*gamma/(gamma-1) if constexpr (has_enthalpy::value) { - state.T = (state.h * state.mu * m_nucleon / C::k_B)*(eos_gamma - 1.0)/eos_gamma; + state.T = (state.h * state.mu * m_nucleon / C::k_B)*(eos_rp::eos_gamma - 1.0)/eos_rp::eos_gamma; } break; @@ -127,7 +125,7 @@ void actual_eos (I input, T& state) // e = k T / [(mu m_nucleon)*(gamma-1)] if constexpr (has_energy::value) { - state.T = state.e * state.mu * m_nucleon * (eos_gamma - 1.0) / C::k_B; + state.T = state.e * state.mu * m_nucleon * (eos_rp::eos_gamma - 1.0) / C::k_B; } break; @@ -162,7 +160,7 @@ void actual_eos (I input, T& state) // Solve for temperature and density if constexpr (has_pressure::value && has_enthalpy::value) { - state.rho = state.p / state.h * eos_gamma / (eos_gamma - 1.0); + state.rho = state.p / state.h * eos_rp::eos_gamma / (eos_rp::eos_gamma - 1.0); state.T = state.p * state.mu * m_nucleon / (C::k_B * state.rho); } @@ -199,7 +197,7 @@ void actual_eos (I input, T& state) // Compute the pressure simply from the ideal gas law, and the // specific internal energy using the gamma-law EOS relation. amrex::Real pressure = state.rho * state.T * C::k_B / (state.mu * m_nucleon); - amrex::Real energy = pressure / (eos_gamma - 1.0) * rhoinv; + amrex::Real energy = pressure / (eos_rp::eos_gamma - 1.0) * rhoinv; if constexpr (has_pressure::value) { state.p = pressure; } @@ -256,17 +254,17 @@ void actual_eos (I input, T& state) state.cv = state.dedT; if constexpr (has_pressure::value) { - state.cp = eos_gamma * state.cv; + state.cp = eos_rp::eos_gamma * state.cv; - state.gam1 = eos_gamma; + state.gam1 = eos_rp::eos_gamma; state.dpdr_e = state.dpdr - state.dpdT * state.dedr * (1.0 / state.dedT); state.dpde = state.dpdT * (1.0 / state.dedT); // sound speed - state.cs = std::sqrt(eos_gamma * state.p * rhoinv); + state.cs = std::sqrt(eos_rp::eos_gamma * state.p * rhoinv); if constexpr (has_G::value) { - state.G = 0.5 * (1.0 + eos_gamma); + state.G = 0.5 * (1.0 + eos_rp::eos_gamma); } } } @@ -278,7 +276,7 @@ void actual_eos (I input, T& state) state.dedA = - state.e * (1.0 / state.abar); } - if (eos_assume_neutral) { + if (eos_rp::eos_assume_neutral) { if constexpr (has_dpdZ::value) { state.dpdZ = 0.0; } diff --git a/EOS/multigamma/actual_eos.H b/EOS/multigamma/actual_eos.H index e673474f7c..7c8189b0d7 100644 --- a/EOS/multigamma/actual_eos.H +++ b/EOS/multigamma/actual_eos.H @@ -14,8 +14,6 @@ #include #include -using namespace eos_rp; - const std::string eos_name = "multigamma"; inline @@ -26,24 +24,24 @@ void actual_eos_init () // that can override the default gammas for a few named species. for (int n = 0; n < NumSpec; ++n) { - gammas[n] = eos_gamma_default; + gammas[n] = eos_rp::eos_gamma_default; } int idx; - idx = network_spec_index(species_a_name); + idx = network_spec_index(eos_rp::species_a_name); if (idx >= 0) { - gammas[idx] = species_a_gamma; + gammas[idx] = eos_rp::species_a_gamma; } - idx = network_spec_index(species_b_name); + idx = network_spec_index(eos_rp::species_b_name); if (idx >= 0) { - gammas[idx] = species_b_gamma; + gammas[idx] = eos_rp::species_b_gamma; } - idx = network_spec_index(species_c_name); + idx = network_spec_index(eos_rp::species_c_name); if (idx >= 0) { - gammas[idx] = species_c_gamma; + gammas[idx] = eos_rp::species_c_gamma; } } diff --git a/EOS/polytrope/actual_eos.H b/EOS/polytrope/actual_eos.H index 3f0108a8e5..639fa72a0e 100644 --- a/EOS/polytrope/actual_eos.H +++ b/EOS/polytrope/actual_eos.H @@ -28,8 +28,6 @@ #include #include -using namespace eos_rp; - const std::string eos_name = "polytrope"; inline @@ -41,9 +39,9 @@ void actual_eos_init () // 1: Non-relativistic, fully degenerate electron gas // 2: Relativistic, fully degenerate electron gas - if (polytrope_type > 0) { - mu_e = polytrope_mu_e; - polytrope = polytrope_type; + if (eos_rp::polytrope_type > 0) { + mu_e = eos_rp::polytrope_mu_e; + polytrope = eos_rp::polytrope_type; if (polytrope == 1) { gamma_const = 5.0_rt / 3.0_rt; @@ -59,9 +57,9 @@ void actual_eos_init () amrex::Error("EOS: Polytrope type currently not defined"); } } - else if (polytrope_gamma > 0.0_rt && polytrope_K > 0.0_rt) { - gamma_const = polytrope_gamma; - K_const = polytrope_K; + else if (eos_rp::polytrope_gamma > 0.0_rt && eos_rp::polytrope_K > 0.0_rt) { + gamma_const = eos_rp::polytrope_gamma; + K_const = eos_rp::polytrope_K; mu_e = 2.0_rt; // This will not be used } else { diff --git a/EOS/primordial_chem/actual_eos.H b/EOS/primordial_chem/actual_eos.H index 3d50fa76b9..9d743432ab 100644 --- a/EOS/primordial_chem/actual_eos.H +++ b/EOS/primordial_chem/actual_eos.H @@ -12,8 +12,6 @@ #include #include -using namespace eos_rp; - const std::string eos_name = "multigamma"; inline @@ -23,94 +21,94 @@ void actual_eos_init () // Set the gammas & masses for the species for (int n = 0; n < NumSpec; ++n) { - gammas[n] = eos_gamma_default; + gammas[n] = eos_rp::eos_gamma_default; spmasses[n] = 1.67353251819e-24; } int idx; - idx = network_spec_index(species_a_name); + idx = network_spec_index(eos_rp::species_a_name); if (idx >= 0) { - gammas[idx] = species_a_gamma; - spmasses[idx] = species_a_mass; + gammas[idx] = eos_rp::species_a_gamma; + spmasses[idx] = eos_rp::species_a_mass; } - idx = network_spec_index(species_b_name); + idx = network_spec_index(eos_rp::species_b_name); if (idx >= 0) { - gammas[idx] = species_b_gamma; - spmasses[idx] = species_b_mass; + gammas[idx] = eos_rp::species_b_gamma; + spmasses[idx] = eos_rp::species_b_mass; } - idx = network_spec_index(species_c_name); + idx = network_spec_index(eos_rp::species_c_name); if (idx >= 0) { - gammas[idx] = species_c_gamma; - spmasses[idx] = species_c_mass; + gammas[idx] = eos_rp::species_c_gamma; + spmasses[idx] = eos_rp::species_c_mass; } - idx = network_spec_index(species_d_name); + idx = network_spec_index(eos_rp::species_d_name); if (idx >= 0) { - gammas[idx] = species_d_gamma; - spmasses[idx] = species_d_mass; + gammas[idx] = eos_rp::species_d_gamma; + spmasses[idx] = eos_rp::species_d_mass; } - idx = network_spec_index(species_e_name); + idx = network_spec_index(eos_rp::species_e_name); if (idx >= 0) { - gammas[idx] = species_e_gamma; - spmasses[idx] = species_e_mass; + gammas[idx] = eos_rp::species_e_gamma; + spmasses[idx] = eos_rp::species_e_mass; } - idx = network_spec_index(species_f_name); + idx = network_spec_index(eos_rp::species_f_name); if (idx >= 0) { - gammas[idx] = species_f_gamma; - spmasses[idx] = species_f_mass; + gammas[idx] = eos_rp::species_f_gamma; + spmasses[idx] = eos_rp::species_f_mass; } - idx = network_spec_index(species_g_name); + idx = network_spec_index(eos_rp::species_g_name); if (idx >= 0) { - gammas[idx] = species_g_gamma; - spmasses[idx] = species_g_mass; + gammas[idx] = eos_rp::species_g_gamma; + spmasses[idx] = eos_rp::species_g_mass; } - idx = network_spec_index(species_h_name); + idx = network_spec_index(eos_rp::species_h_name); if (idx >= 0) { - gammas[idx] = species_h_gamma; - spmasses[idx] = species_h_mass; + gammas[idx] = eos_rp::species_h_gamma; + spmasses[idx] = eos_rp::species_h_mass; } - idx = network_spec_index(species_i_name); + idx = network_spec_index(eos_rp::species_i_name); if (idx >= 0) { - gammas[idx] = species_i_gamma; - spmasses[idx] = species_i_mass; + gammas[idx] = eos_rp::species_i_gamma; + spmasses[idx] = eos_rp::species_i_mass; } - idx = network_spec_index(species_j_name); + idx = network_spec_index(eos_rp::species_j_name); if (idx >= 0) { - gammas[idx] = species_j_gamma; - spmasses[idx] = species_j_mass; + gammas[idx] = eos_rp::species_j_gamma; + spmasses[idx] = eos_rp::species_j_mass; } - idx = network_spec_index(species_k_name); + idx = network_spec_index(eos_rp::species_k_name); if (idx >= 0) { - gammas[idx] = species_k_gamma; - spmasses[idx] = species_k_mass; + gammas[idx] = eos_rp::species_k_gamma; + spmasses[idx] = eos_rp::species_k_mass; } - idx = network_spec_index(species_l_name); + idx = network_spec_index(eos_rp::species_l_name); if (idx >= 0) { - gammas[idx] = species_l_gamma; - spmasses[idx] = species_l_mass; + gammas[idx] = eos_rp::species_l_gamma; + spmasses[idx] = eos_rp::species_l_mass; } - idx = network_spec_index(species_m_name); + idx = network_spec_index(eos_rp::species_m_name); if (idx >= 0) { - gammas[idx] = species_m_gamma; - spmasses[idx] = species_m_mass; + gammas[idx] = eos_rp::species_m_gamma; + spmasses[idx] = eos_rp::species_m_mass; } - idx = network_spec_index(species_n_name); + idx = network_spec_index(eos_rp::species_n_name); if (idx >= 0) { - gammas[idx] = species_n_gamma; - spmasses[idx] = species_n_mass; + gammas[idx] = eos_rp::species_n_gamma; + spmasses[idx] = eos_rp::species_n_mass; } } diff --git a/EOS/rad_power_law/actual_eos.H b/EOS/rad_power_law/actual_eos.H index 5a1fabbca6..4a3ed428b2 100644 --- a/EOS/rad_power_law/actual_eos.H +++ b/EOS/rad_power_law/actual_eos.H @@ -18,19 +18,17 @@ #include #include -using namespace eos_rp; - const std::string eos_name = "rad_power_law"; inline void actual_eos_init () { - if (eos_const_c_v <= 0.e0_rt) { + if (eos_rp::eos_const_c_v <= 0.e0_rt) { amrex::Error("eos_const_c_v must be > 0"); } - if (eos_c_v_exp_n == 1.0e0_rt) { + if (eos_rp::eos_c_v_exp_n == 1.0e0_rt) { amrex::Error("eos_c_v_exp_n == 1 is unsupported"); } @@ -77,8 +75,8 @@ void actual_eos (I input, T& state) case eos_input_rt: if constexpr (has_energy::value) { - state.cv = eos_const_c_v * std::pow(state.rho, eos_c_v_exp_m) * std::pow(state.T, -eos_c_v_exp_n); - state.e = eos_const_c_v * std::pow(state.rho, eos_c_v_exp_m) * std::pow(state.T, 1 - eos_c_v_exp_n) / (1 - eos_c_v_exp_n); + state.cv = eos_rp::eos_const_c_v * std::pow(state.rho, eos_rp::eos_c_v_exp_m) * std::pow(state.T, -eos_rp::eos_c_v_exp_n); + state.e = eos_rp::eos_const_c_v * std::pow(state.rho, eos_rp::eos_c_v_exp_m) * std::pow(state.T, 1 - eos_rp::eos_c_v_exp_n) / (1 - eos_rp::eos_c_v_exp_n); } break; @@ -86,8 +84,8 @@ void actual_eos (I input, T& state) case eos_input_re: if constexpr (has_energy::value) { - state.T = std::pow((1 - eos_c_v_exp_n) * state.e * std::pow(state.rho, -eos_c_v_exp_m) / eos_const_c_v, 1.0_rt / (1.0_rt - eos_c_v_exp_n)); - state.cv = eos_const_c_v * std::pow(state.rho, eos_c_v_exp_m) * std::pow(state.T, -eos_c_v_exp_n); + state.T = std::pow((1 - eos_rp::eos_c_v_exp_n) * state.e * std::pow(state.rho, -eos_rp::eos_c_v_exp_m) / eos_rp::eos_const_c_v, 1.0_rt / (1.0_rt - eos_rp::eos_c_v_exp_n)); + state.cv = eos_rp::eos_const_c_v * std::pow(state.rho, eos_rp::eos_c_v_exp_m) * std::pow(state.T, -eos_rp::eos_c_v_exp_n); } break; diff --git a/EOS/tillotson/actual_eos.H b/EOS/tillotson/actual_eos.H index c9df21832a..a7805d4d15 100644 --- a/EOS/tillotson/actual_eos.H +++ b/EOS/tillotson/actual_eos.H @@ -7,8 +7,6 @@ #include #include -using namespace eos_rp; - const std::string eos_name = "tillotson"; inline @@ -19,7 +17,7 @@ void actual_eos_init() { // (Note that rho_0 is not the minimum density but rather just // a reference density.) - EOSData::mine = eos_e_0; + EOSData::mine = eos_rp::eos_e_0; } @@ -56,55 +54,55 @@ void actual_eos (I input, T& state) // that we are consistent with eos_input_re. if (input == eos_input_re) { - state.T = amrex::max((state.e - eos_e_0) / eos_c_v, EOSData::mintemp); + state.T = amrex::max((state.e - eos_rp::eos_e_0) / eos_rp::eos_c_v, EOSData::mintemp); } else if (input == eos_input_rt) { - state.e = amrex::max(eos_e_0 + eos_c_v * state.T, EOSData::mine); + state.e = amrex::max(eos_rp::eos_e_0 + eos_rp::eos_c_v * state.T, EOSData::mine); } - amrex::Real eta = state.rho / eos_rho_0; + amrex::Real eta = state.rho / eos_rp::eos_rho_0; amrex::Real mu = eta - 1.0_rt; - amrex::Real omega = state.e / (eos_e_0 * eta * eta) + 1.0_rt; + amrex::Real omega = state.e / (eos_rp::eos_e_0 * eta * eta) + 1.0_rt; amrex::Real z = (1.0_rt / eta - 1.0_rt); - amrex::Real g_c = eos_la + eos_lb / omega; - amrex::Real g_e = eos_la + eos_lb / omega * std::exp(-eos_beta * z * z); + amrex::Real g_c = eos_rp::eos_la + eos_rp::eos_lb / omega; + amrex::Real g_e = eos_rp::eos_la + eos_rp::eos_lb / omega * std::exp(-eos_rp::eos_beta * z * z); - amrex::Real P_c = g_c * state.e * state.rho + eos_A * mu + eos_B * mu * mu; - amrex::Real P_e = g_e * state.e * state.rho + eos_A * mu * std::exp(-eos_alpha * z - eos_beta * z * z); + amrex::Real P_c = g_c * state.e * state.rho + eos_rp::eos_A * mu + eos_rp::eos_B * mu * mu; + amrex::Real P_e = g_e * state.e * state.rho + eos_rp::eos_A * mu * std::exp(-eos_rp::eos_alpha * z - eos_rp::eos_beta * z * z); // Floor the pressure since the above expressions can be negative P_c = amrex::max(P_c, EOSData::minp); P_e = amrex::max(P_e, EOSData::minp); - amrex::Real cs_squared_c = (g_c + 2.0_rt * eos_lb * (omega - 1.0_rt) / (omega * omega)) * state.e + - (eos_A + 2.0_rt * eos_B * mu) / eos_rho_0; - amrex::Real cs_squared_e = (g_e + (2.0_rt * eos_lb / omega) * std::exp(-eos_beta * z * z) * - ((omega - 1.0_rt) / omega + eos_beta * z / eta)) * state.e + - (eos_A / eos_rho_0) * std::exp(-eos_alpha * z - eos_beta * z * z) * - (1.0_rt + mu / (eta * eta) * (eos_alpha + 2.0_rt * eos_beta * z)); + amrex::Real cs_squared_c = (g_c + 2.0_rt * eos_rp::eos_lb * (omega - 1.0_rt) / (omega * omega)) * state.e + + (eos_rp::eos_A + 2.0_rt * eos_rp::eos_B * mu) / eos_rp::eos_rho_0; + amrex::Real cs_squared_e = (g_e + (2.0_rt * eos_rp::eos_lb / omega) * std::exp(-eos_rp::eos_beta * z * z) * + ((omega - 1.0_rt) / omega + eos_rp::eos_beta * z / eta)) * state.e + + (eos_rp::eos_A / eos_rp::eos_rho_0) * std::exp(-eos_rp::eos_alpha * z - eos_rp::eos_beta * z * z) * + (1.0_rt + mu / (eta * eta) * (eos_rp::eos_alpha + 2.0_rt * eos_rp::eos_beta * z)); // Floor the sound speed as well (and make sure we get the floor dimensionally right) cs_squared_c = amrex::max(cs_squared_c, EOSData::minp / state.rho); cs_squared_e = amrex::max(cs_squared_e, EOSData::minp / state.rho); - amrex::Real dPde_c = eos_la * state.rho + eos_lb * state.rho / (omega * omega); - amrex::Real dPde_e = eos_la * state.rho + eos_lb * state.rho / (omega * omega) * std::exp(-eos_beta * z * z); + amrex::Real dPde_c = eos_rp::eos_la * state.rho + eos_rp::eos_lb * state.rho / (omega * omega); + amrex::Real dPde_e = eos_rp::eos_la * state.rho + eos_rp::eos_lb * state.rho / (omega * omega) * std::exp(-eos_rp::eos_beta * z * z); // At this point we should have valid (rho, e); compute P = P(rho, e). amrex::Real pres, dpdr_e, dpde; - if (state.rho >= eos_rho_0 || (state.rho < eos_rho_0 && state.e < eos_e_s)) { + if (state.rho >= eos_rp::eos_rho_0 || (state.rho < eos_rp::eos_rho_0 && state.e < eos_rp::eos_e_s)) { // Region I, II (compressed state) pres = P_c; dpdr_e = cs_squared_c; dpde = dPde_c; } - else if (state.rho < eos_rho_0 && state.e > eos_e_s_prime) { + else if (state.rho < eos_rp::eos_rho_0 && state.e > eos_rp::eos_e_s_prime) { // Region IV (expanded state) pres = P_e; @@ -114,9 +112,9 @@ void actual_eos (I input, T& state) else { // Region III (interpolate between compressed and expanded state) - amrex::Real denom = (eos_e_s_prime - eos_e_s); - amrex::Real xi_c = (state.e - eos_e_s) / denom; - amrex::Real xi_e = (eos_e_s_prime - state.e) / denom; + amrex::Real denom = (eos_rp::eos_e_s_prime - eos_rp::eos_e_s); + amrex::Real xi_c = (state.e - eos_rp::eos_e_s) / denom; + amrex::Real xi_e = (eos_rp::eos_e_s_prime - state.e) / denom; pres = P_c * xi_c + P_e * xi_e; dpdr_e = cs_squared_c * xi_c + cs_squared_e * xi_e; @@ -137,7 +135,7 @@ void actual_eos (I input, T& state) state.dpdr = state.dpdr_e; // Use chain rule to get dp/dT from de/dT - state.dpdT = state.dpde * eos_c_v; + state.dpdT = state.dpde * eos_rp::eos_c_v; } // Sound speed @@ -154,15 +152,15 @@ void actual_eos (I input, T& state) // Fill in heat capacity (c_v from provided data, c_p from c_v) - state.cv = eos_c_v; + state.cv = eos_rp::eos_c_v; if constexpr (has_pressure::value) { - state.cp = state.gam1 * eos_c_v; + state.cp = state.gam1 * eos_rp::eos_c_v; } // de/dT == specific heat at constant volume - state.dedT = eos_c_v; + state.dedT = eos_rp::eos_c_v; // de/drho is not relevant for this EOS