Merge pull request #1506 from CliMA/js/phasenonequil

use PhaseNonEquil instead of PhaseEquil

Author: juliasloan25

NEWS.md

@@ -6,6 +6,9 @@ ClimaCoupler.jl Release Notes
 
 ### ClimaCoupler features
 
+#### Replace `TD.PhaseEquil_ρTq` with `TD.PhaseNonEquil_ρTq` PR[#1506](https://github.com/CliMA/ClimaCoupler.jl/pull/1506)
+This should be more physically correct.
+
 #### Use individual surface model temperatures to compute turbulent fluxes PR[#1498](https://github.com/CliMA/ClimaCoupler.jl/pull/1498)
 We use a partitioned approach to computing turbulent fluxes, so we should be using
 the individual surface temperature (and humidity, air density, etc) for each component

src/FluxCalculator.jl

@@ -218,7 +218,8 @@ All fields should be on the exchange grid.
 - `thermo_params`: [TD.Parameters.ThermodynamicsParameters] the thermodynamic parameters.
 """
 function compute_surface_humidity!(q_sfc, T_atmos, q_atmos, ρ_atmos, T_sfc, thermo_params)
-    thermo_state_atmos = TD.PhaseEquil_ρTq.(thermo_params, ρ_atmos, T_atmos, q_atmos)
+    thermo_state_atmos =
+        TD.PhaseNonEquil_ρTq.(thermo_params, ρ_atmos, T_atmos, TD.PhasePartition.(q_atmos))
     ρ_sfc = FluxCalculator.extrapolate_ρ_to_sfc.(thermo_params, thermo_state_atmos, T_sfc)
 
     T_freeze = TDP.T_freeze(thermo_params)
@@ -279,7 +280,12 @@ function compute_surface_fluxes!(
 
     # construct the atmospheric thermo state
     thermo_state_atmos =
-        TD.PhaseEquil_ρTq.(thermo_params, csf.ρ_atmos, csf.T_atmos, csf.q_atmos)
+        TD.PhaseNonEquil_ρTq.(
+            thermo_params,
+            csf.ρ_atmos,
+            csf.T_atmos,
+            TD.PhasePartition.(csf.q_atmos),
+        )
 
     # compute surface humidity from the surface temperature, surface density, and phase
     Interfacer.get_field!(csf.scalar_temp1, sim, Val(:surface_temperature))
@@ -299,7 +305,8 @@ function compute_surface_fluxes!(
 
     # construct the surface thermo state
     # after this we can reuse `scalar_temp1` and `scalar_temp2` again
-    thermo_state_sfc = TD.PhaseEquil_ρTq.(thermo_params, ρ_sfc, T_sfc, q_sfc)
+    thermo_state_sfc =
+        TD.PhaseNonEquil_ρTq.(thermo_params, ρ_sfc, T_sfc, TD.PhasePartition.(q_sfc))
 
     # get area fraction (min = 0, max = 1)
     Interfacer.get_field!(csf.scalar_temp1, sim, Val(:area_fraction))

test/flux_calculator_tests.jl

@@ -191,11 +191,11 @@ for FT in (Float32, Float64)
                 Interfacer.get_field(atmos_sim, Val(:v_int)),
             )
         thermo_state_atmos =
-            TD.PhaseEquil_ρTq.(
+            TD.PhaseNonEquil_ρTq.(
                 thermo_params,
                 atmos_sim.integrator.ρ,
                 atmos_sim.integrator.T,
-                atmos_sim.integrator.q,
+                TD.PhasePartition.(atmos_sim.integrator.q),
             )
 
         # Get surface properties
@@ -217,7 +217,8 @@ for FT in (Float32, Float64)
         )
         ρ_sfc =
             FluxCalculator.extrapolate_ρ_to_sfc.(thermo_params, thermo_state_atmos, T_sfc)
-        thermo_state_sfc = TD.PhaseEquil_ρTq.(thermo_params, ρ_sfc, T_sfc, q_sfc)
+        thermo_state_sfc =
+            TD.PhaseNonEquil_ρTq.(thermo_params, ρ_sfc, T_sfc, TD.PhasePartition.(q_sfc))
 
         # Use SurfaceFluxes.jl to compute the expected fluxes
         inputs = @. SF.ValuesOnly(