add RH_ice to moist env thermodynamic state

PySDM/backends/impl_numba/methods/physics_methods.py

@@ -80,24 +80,28 @@ def _a_w_ice_body(self):
         ff = self.formulae_flattened
 
         @numba.njit(**self.default_jit_flags)
-        def body(*, T_in, p_in, RH_in, water_vapour_mixing_ratio_in, a_w_ice_out):
+        def body(
+            *, T_in, p_in, RH_in, water_vapour_mixing_ratio_in, a_w_ice_out, RH_ice_out
+        ):
             for i in prange(T_in.shape[0]):  # pylint: disable=not-an-iterable
                 pvi = ff.saturation_vapour_pressure__pvs_ice(T_in[i])
                 pv = ff.state_variable_triplet__pv(
                     p_in[i], water_vapour_mixing_ratio_in[i]
                 )
                 pvs = pv / RH_in[i]
                 a_w_ice_out[i] = pvi / pvs
+                RH_ice_out[i] = pv / pvi
 
         return body
 
-    def a_w_ice(self, *, T, p, RH, water_vapour_mixing_ratio, a_w_ice):
+    def a_w_ice(self, *, T, p, RH, water_vapour_mixing_ratio, a_w_ice, RH_ice):
         self._a_w_ice_body(
             T_in=T.data,
             p_in=p.data,
             RH_in=RH.data,
             water_vapour_mixing_ratio_in=water_vapour_mixing_ratio.data,
             a_w_ice_out=a_w_ice.data,
+            RH_ice_out=RH_ice.data,
         )
 
     @cached_property

PySDM/environments/impl/moist.py

@@ -11,7 +11,7 @@ class Moist:
     def __init__(self, dt, mesh, variables, mixed_phase=False):
         variables += ["water_vapour_mixing_ratio", "thd", "T", "p", "RH"]
         if mixed_phase:
-            variables += ["a_w_ice"]
+            variables += ["a_w_ice", "RH_ice"]
         all_vars_unique = len(variables) == len(set(variables))
         assert all_vars_unique
 
@@ -84,6 +84,7 @@ def sync(self):
                 RH=target["RH"],
                 water_vapour_mixing_ratio=target["water_vapour_mixing_ratio"],
                 a_w_ice=target["a_w_ice"],
+                RH_ice=target["RH_ice"],
             )
         if "air density" in self.variables:
             self.particulator.backend.air_density(

tests/unit_tests/environments/test_moist.py

@@ -0,0 +1,43 @@
+""" check for ice-phase-related commons in environment codes """
+
+import pytest
+import numpy as np
+from PySDM.environments import Parcel
+from PySDM.physics import si
+from PySDM import Builder
+from PySDM.backends import GPU
+
+
+@pytest.mark.parametrize(
+    "env",
+    (
+        Parcel(
+            mixed_phase=True,
+            dt=np.nan,
+            mass_of_dry_air=np.nan,
+            p0=1000 * si.hPa,
+            initial_water_vapour_mixing_ratio=20 * si.g / si.kg,
+            T0=300 * si.K,
+            w=np.nan,
+        ),
+    ),
+)
+def test_ice_properties(backend_instance, env):
+    """checks ice-related values in recalculated thermodynamic state make sense"""
+    if isinstance(backend_instance, GPU):
+        pytest.skip("TODO #1495")
+
+    # arrange
+    builder = Builder(n_sd=0, backend=backend_instance, environment=env)
+
+    # act
+    thermo = {
+        key: builder.particulator.environment[key].to_ndarray()[0]
+        for key in ("RH", "RH_ice", "a_w_ice")
+    }
+
+    # assert
+    assert 1 > thermo["RH"] > thermo["RH_ice"] > 0
+    np.testing.assert_approx_equal(
+        thermo["a_w_ice"] * thermo["RH_ice"], thermo["RH"], significant=10
+    )