diff --git a/components/eamxx/src/diagnostics/potential_temperature.cpp b/components/eamxx/src/diagnostics/potential_temperature.cpp index 602f3ffc218..988bb555224 100644 --- a/components/eamxx/src/diagnostics/potential_temperature.cpp +++ b/components/eamxx/src/diagnostics/potential_temperature.cpp @@ -50,6 +50,7 @@ void PotentialTemperatureDiagnostic::set_grids(const std::shared_ptr("T_mid", scalar3d_layout_mid, K, grid_name, ps); add_field("p_mid", scalar3d_layout_mid, Pa, grid_name, ps); // Only needed for LiqPotentialTemperature, but put it here for ease + // TODO: only request it if it is needed add_field("qc", scalar3d_layout_mid, Q, grid_name, ps); // Construct and allocate the diagnostic field diff --git a/components/eamxx/src/diagnostics/tests/potential_temperature_test.cpp b/components/eamxx/src/diagnostics/tests/potential_temperature_test.cpp index d82f85fa8ce..1f584ec7e71 100644 --- a/components/eamxx/src/diagnostics/tests/potential_temperature_test.cpp +++ b/components/eamxx/src/diagnostics/tests/potential_temperature_test.cpp @@ -91,7 +91,6 @@ void run(std::mt19937_64& engine, int int_ptype) auto diag = diag_factory.create("PotentialTemperature",comm,params); diag->set_grids(gm); - // Set the required fields for the diagnostic. std::map input_fields; for (const auto& req : diag->get_required_field_requests()) { diff --git a/components/eamxx/src/physics/shoc/eamxx_shoc_process_interface.hpp b/components/eamxx/src/physics/shoc/eamxx_shoc_process_interface.hpp index 8d9feb496b8..d3f5b05f88d 100644 --- a/components/eamxx/src/physics/shoc/eamxx_shoc_process_interface.hpp +++ b/components/eamxx/src/physics/shoc/eamxx_shoc_process_interface.hpp @@ -75,7 +75,6 @@ class SHOCMacrophysics : public scream::AtmosphereProcess const int i = team.league_rank(); const Real zvir = C::ZVIR; - const Real latvap = C::LatVap; const Real cpair = C::Cpair; const Real ggr = C::gravit; const Real inv_ggr = 1/ggr; @@ -107,12 +106,9 @@ class SHOCMacrophysics : public scream::AtmosphereProcess qw(i,k) = qv(i,k) + qc(i,k); // Temperature - // NOTE: vtheta (thv) is inconsistent with one in HOMME - // SEE gh issue #2813 and gh pr #2812 + // NOTE: theta_v (thv) is intentionally different from one in HOMME const auto theta_zt = PF::calculate_theta_from_T(T_mid(i,k),p_mid(i,k)); - // TODO? Can use PF::calculate_thetal_from_theta instead? (added in gh pr #2812) - thlm(i,k) = theta_zt-(theta_zt/T_mid(i,k))*(latvap/cpair)*qc(i,k); - // TODO? Can use PF::calculate_virtual_temperature instead? (also be consistent with HOMME) + thlm(i,k) = PF::calculate_thetal_from_theta(theta_zt,T_mid(i,k),qc(i,k)); thv(i,k) = theta_zt*(1 + zvir*qv(i,k) - qc(i,k)); // Vertical layer thickness diff --git a/components/eamxx/src/share/util/scream_common_physics_functions.hpp b/components/eamxx/src/share/util/scream_common_physics_functions.hpp index b926a7b7bc8..87ae6604d99 100644 --- a/components/eamxx/src/share/util/scream_common_physics_functions.hpp +++ b/components/eamxx/src/share/util/scream_common_physics_functions.hpp @@ -133,14 +133,8 @@ struct PhysicsFunctions static ScalarT calculate_temperature_from_virtual_temperature(const ScalarT& T_virtual, const ScalarT& qv); //-----------------------------------------------------------------------------------------------// - // FIXME: this comment is different from what's actually happening in the code: - // FIXME: T_virtual = temperature * ( one + c1*qv ) where c1 = - one + one/ep_2 - // FIXME: The code seems to be correct, and the comment is wrong. - // FIXME: ep_2 = 0.622, so the comment produces vpt = pt * (0.75) for qv = 0.5 - // FIXME: but the code produces vpt = pt * (1.30) for qv = 0.5 - // XREFS: pr #1059, issue #2812 // Compute virtual temperature - // T_virtual = temperature * (qv+ep_2)/(qv+1) + // T_virtual = temperature * (1+(-1+1/ep_2)*qv) // where // ep_2 is ratio of molecular mass of water to the molecular mass of dry air // temperature is the atmospheric temperature. Units in [K]. @@ -411,12 +405,12 @@ struct PhysicsFunctions const InputProviderP& pressure, const view_1d& theta); - template + template KOKKOS_INLINE_FUNCTION static void calculate_thetal_from_theta(const MemberType& team, - const InputProviderT& theta, + const InputProviderTheta& theta, const InputProviderT& temperature, - const InputProviderP& qc, + const InputProviderQ& qc, const view_1d& thetal); template KOKKOS_INLINE_FUNCTION diff --git a/components/eamxx/src/share/util/scream_common_physics_functions_impl.hpp b/components/eamxx/src/share/util/scream_common_physics_functions_impl.hpp index 769622a30ad..f41e389bf67 100644 --- a/components/eamxx/src/share/util/scream_common_physics_functions_impl.hpp +++ b/components/eamxx/src/share/util/scream_common_physics_functions_impl.hpp @@ -123,10 +123,8 @@ KOKKOS_INLINE_FUNCTION ScalarT PhysicsFunctions::calculate_thetal_from_theta(const ScalarT& theta, const ScalarT& temperature, const ScalarT& qc) { using C = scream::physics::Constants; - static constexpr auto cpair = C::Cpair; - static constexpr auto latvap = C::LatVap; - return theta - (theta / temperature) * (latvap/cpair) * qc; + return theta - (theta / temperature) * (C::LatVap/C::Cpair) * qc; } template @@ -144,13 +142,13 @@ void PhysicsFunctions::calculate_theta_from_T(const MemberType& team, } template -template +template KOKKOS_INLINE_FUNCTION void PhysicsFunctions::calculate_thetal_from_theta(const MemberType& team, - const InputProviderT& theta, - const InputProviderT& temperature, - const InputProviderP& qc, - const view_1d& thetal) + const InputProviderTheta& theta, + const InputProviderT& temperature, + const InputProviderQ& qc, + const view_1d& thetal) { Kokkos::parallel_for(Kokkos::TeamVectorRange(team,thetal.extent(0)), [&] (const int k) {