Support equivalent and identical units in capgen and ccpp-prebuild

scripts/common.py

@@ -148,6 +148,21 @@ def isstring(s):
     """Return true if a variable is a string"""
     return isinstance(s, str)
 
+def insert_plus_sign_for_positive_exponents(string):
+    """Parse a string (a unit string) and insert plus (+) signs
+    for positive exponents where needed"""
+    # Break up the string by spaces
+    items = string.split()
+    # Identify units with positive exponents
+    # without a plus sign (m2 instead of m+2).
+    pattern = re.compile(r"([a-zA-Z]+)([0-9]+)")
+    for index, item in enumerate(items):
+        match = pattern.match(item)
+        if match:
+            items[index] = "+".join(match.groups())
+    # Recombine items to string
+    return " ".join(items)
+
 def string_to_python_identifier(string):
     """Replaces forbidden characters in strings with standard substitutions
     so that the result is a valid Python object (variable, function) name.

scripts/conversion_tools/unit_conversion.py

@@ -176,3 +176,23 @@ def W_m_minus_2__to__erg_cm_minus_2_s_minus_1():
 def erg_cm_minus_2_s_minus_1__to__W_m_minus_2():
     """Convert erg per square centimeter and second to Watt per square meter"""
     return '1.0E-3{kind}*{var}'
+
+####################
+# Equivalent units #
+####################
+
+def m_plus_2_s_minus_2__to__J_kg_minus_1():
+    """Equivalent units"""
+    return '{var}'
+
+def J_kg_minus_1__to__m_plus_2_s_minus_2():
+    """Equivalent units"""
+    return '{var}'
+
+def V_A__to__W():
+    """Equivalent units"""
+    return '{var}'
+
+def W__to__V_A():
+    """Equivalent units"""
+    return '{var}'

scripts/metadata_parser.py

@@ -10,6 +10,7 @@
 
 from common import encode_container, CCPP_STAGES
 from common import CCPP_ERROR_CODE_VARIABLE, CCPP_ERROR_MSG_VARIABLE
+from common import insert_plus_sign_for_positive_exponents
 from mkcap import Var
 
 sys.path.append(os.path.join(os.path.split(__file__)[0], 'fortran_tools'))
@@ -230,9 +231,14 @@ def read_new_metadata(filename, module_name, table_name, scheme_name = None, sub
                 #kind = new_var.get_prop_value('kind')
                 # *DH 20210812
 
+                # Workaround to support units with positive exponents with
+                # and without a plus (+) sign. Internally, we convert all
+                # units from capgen to the "+"-format (i.e. "m2 s-2" --> "m+2 s-2")
+                units = insert_plus_sign_for_positive_exponents(new_var.get_prop_value('units'))
+
                 var = Var(standard_name = standard_name,
                           long_name     = new_var.get_prop_value('long_name') + legacy_note,
-                          units         = new_var.get_prop_value('units'),
+                          units         = units,
                           local_name    = new_var.get_prop_value('local_name'),
                           type          = new_var.get_prop_value('type').lower(),
                           dimensions    = dimensions,

scripts/parse_tools/parse_checkers.py

@@ -16,7 +16,8 @@
 _NON_LEADING_ZERO_NUM          = "[1-9]\d*"
 _CHAR_WITH_UNDERSCORE          = "([a-zA-Z]+_[a-zA-Z]+)+"
 _NEGATIVE_NON_LEADING_ZERO_NUM = f"[-]{_NON_LEADING_ZERO_NUM}"
-_UNIT_EXPONENT                 = f"({_NEGATIVE_NON_LEADING_ZERO_NUM}|{_NON_LEADING_ZERO_NUM})"
+_POSITIVE_NON_LEADING_ZERO_NUM = f"[+]{_NON_LEADING_ZERO_NUM}"
+_UNIT_EXPONENT                 = f"({_NEGATIVE_NON_LEADING_ZERO_NUM}|{_POSITIVE_NON_LEADING_ZERO_NUM}|{_NON_LEADING_ZERO_NUM})"
 _UNIT_REGEX                    = f"[a-zA-Z]+{_UNIT_EXPONENT}?"
 _UNITS_REGEX                   = f"^({_CHAR_WITH_UNDERSCORE}|{_UNIT_REGEX}(\s{_UNIT_REGEX})*|{_UNITLESS_REGEX})$"
 _UNITS_RE                      = re.compile(_UNITS_REGEX)
@@ -29,6 +30,10 @@ def check_units(test_val, prop_dict, error):
     'm s-1'
     >>> check_units('kg m-3', None, True)
     'kg m-3'
+    >>> check_units('m2 s-2', None, True)
+    'm2 s-2'
+    >>> check_units('m+2 s-2', None, True)
+    'm+2 s-2'
     >>> check_units('1', None, True)
     '1'
     >>> check_units('', None, False)

scripts/var_props.py

@@ -852,6 +852,20 @@ class VarCompatObj:
                      "var_stdname", "real", "kind_phys", "km",['horizontal_dimension', 'vertical_layer_dimension'], "var2_lname", True, \
                      _DOCTEST_RUNENV).reverse_transform("var1_lname", "var2_lname", ('i','k'), ('i','nk-k+1'))
     'var1_lname(i,nk-k+1) = 1.0E+3_kind_phys*var2_lname(i,k)'
+
+    # Test that a 2-D var with equivalent units works and that it
+    # skips any unit transformations
+    >>> VarCompatObj("var_stdname", "real", "kind_phys", "m2 s-2",  ['horizontal_dimension'], "var1_lname", False, \
+                     "var_stdname", "real", "kind_phys", "J kg-1", ['horizontal_dimension'], "var2_lname", False, \
+                     _DOCTEST_RUNENV).forward_transform("var1_lname", "var2_lname", 'i', 'i')
+    'var1_lname(i) = var2_lname(i)'
+
+    # Test that a 2-D var with identical units works and that it
+    # skips any unit transformations
+    >>> VarCompatObj("var_stdname", "real", "kind_phys", "m2 s-2",  ['horizontal_dimension'], "var1_lname", False, \
+                     "var_stdname", "real", "kind_phys", "m+2 s-2", ['horizontal_dimension'], "var2_lname", False, \
+                     _DOCTEST_RUNENV).forward_transform("var1_lname", "var2_lname", 'i', 'i')
+    'var1_lname(i) = var2_lname(i)'
     """
 
     def __init__(self, var1_stdname, var1_type, var1_kind, var1_units,
@@ -960,20 +974,27 @@ def __init__(self, var1_stdname, var1_type, var1_kind, var1_units,
                 # A tendency variable's units should be "<var2_units> s-1"
                 tendency_split_units = var1_units.split('s-1')[0].strip()
                 if tendency_split_units != var2_units:
-                    # We don't currently support unit conversions for tendency variables
-                    emsg = f"\nMismatch tendency variable units '{var1_units}'"
-                    emsg += f" for variable '{var1_stdname}'."
-                    emsg += " No variable transforms supported for tendencies."
-                    emsg += f" Tendency units should be '{var2_units} s-1' to match state variable."
-                    self.__equiv = False
-                    self.__compat = False
-                    incompat_reason.append(emsg)
+                    # We don't currently support unit conversions for tendency variables,
+                    # but we can check if the units are identical or equivalent
+                    unit_transforms = self._get_unit_convstrs(tendency_split_units,
+                                                              var2_units)
+                    if not unit_transforms == (None, None):
+                        emsg = f"\nMismatch tendency variable units '{var1_units}'"
+                        emsg += f" for variable '{var1_stdname}'."
+                        emsg += " No variable transforms supported for tendencies."
+                        emsg += f" Tendency units should be '{var2_units} s-1' to match state variable."
+                        self.__equiv = False
+                        self.__compat = False
+                        incompat_reason.append(emsg)
                 # end if
             elif var1_units != var2_units:
                 # Try to find a set of unit conversions
-                self.__equiv = False
-                self.__unit_transforms = self._get_unit_convstrs(var1_units,
-                                                                 var2_units)
+                unit_transforms = self._get_unit_convstrs(var1_units,
+                                                          var2_units)
+                # Handle equivalent or identical units = (None, None)
+                if not unit_transforms == (None, None):
+                    self.__equiv = False
+                    self.__unit_transforms = unit_transforms
             # end if
         # end if
         if self.__compat:
@@ -1148,7 +1169,8 @@ def _get_kind_convstrs(self, var1_kind, var2_kind, run_env):
     def _get_unit_convstrs(self, var1_units, var2_units):
         """Attempt to retrieve the forward and reverse unit transformations
         for transforming a variable in <var1_units> to / from a variable in
-        <var2_units>.
+        <var2_units>. Return (None, None) if units are equivalent or identical
+        after parsing (this can happen when comparing m2 and m+2).
 
         # Initial setup
         >>> from parse_tools import init_log, set_log_to_null
@@ -1177,6 +1199,14 @@ def _get_unit_convstrs(self, var1_units, var2_units):
         ('1.0E+3{kind}*{var}', '1.0E-3{kind}*{var}')
         >>> _DOCTEST_VCOMPAT._get_unit_convstrs('C', 'K')
         ('{var}+273.15{kind}', '{var}-273.15{kind}')
+        >>> _DOCTEST_VCOMPAT._get_unit_convstrs('V A', 'W')
+        (None, None)
+        >>> _DOCTEST_VCOMPAT._get_unit_convstrs('m2 s-2', 'J kg-1')
+        (None, None)
+        >>> _DOCTEST_VCOMPAT._get_unit_convstrs('m+2 s-2', 'J kg-1')
+        (None, None)
+        >>> _DOCTEST_VCOMPAT._get_unit_convstrs('m+2 s-2', 'm2 s-2')
+        (None, None)
 
         # Try an invalid conversion
         >>> _DOCTEST_VCOMPAT._get_unit_convstrs('1', 'none') #doctest: +ELLIPSIS
@@ -1192,6 +1222,11 @@ def _get_unit_convstrs(self, var1_units, var2_units):
         """
         u1_str = self.units_to_string(var1_units, self.__v1_context)
         u2_str = self.units_to_string(var2_units, self.__v2_context)
+        # If u1_str and u2_str are identical, for example after parsing
+        # "m2 s-2" and "m+2 s-2", return (None, None) to signal that
+        # the units are in fact identical
+        if u1_str == u2_str:
+            return (None, None)
         unit_conv_str = "{0}__to__{1}".format(u1_str, u2_str)
         try:
             forward_transform = getattr(unit_conversion, unit_conv_str)()
@@ -1210,7 +1245,11 @@ def _get_unit_convstrs(self, var1_units, var2_units):
                                                self.__stdname,
                                                context=self.__v1_context))
         # end if
-        return (forward_transform, reverse_transform)
+        # For equivalent units, return (None, None)
+        if forward_transform == '{var}' and reverse_transform == '{var}':
+            return (None, None)
+        else:
+            return (forward_transform, reverse_transform)
 
     def _get_dim_transforms(self, var1_dims, var2_dims):
         """Attempt to find forward and reverse permutations for transforming a
@@ -1355,8 +1394,17 @@ def units_to_string(self, units, context=None):
         """Replace variable unit description with string that is a legal
         Python identifier.
         If the resulting string is a Python keyword, raise an exception."""
-        # Replace each whitespace with an underscore
-        string = units.replace(" ","_")
+        # Start with breaking up the string by spaces
+        items = units.split()
+        # Identify units with positive exponents
+        # without a plus sign (m2 instead of m+2).
+        pattern = re.compile(r"([a-zA-Z]+)([0-9]+)")
+        for index, item in enumerate(items):
+            match = pattern.match(item)
+            if match:
+                items[index] = "+".join(match.groups())
+        # Combine list into string using underscores
+        string = "_".join(items)
         # Replace each minus sign with '_minus_'
         string = string.replace("-","_minus_")
         # Replace each plus sign with '_plus_'
@@ -1458,7 +1506,7 @@ def has_unit_transforms(self):
            and <var> arguments to produce code to transform one variable into
            the correct units of the other.
         """
-        return self.__unit_transforms is not None
+        return self.__unit_transforms is not None and self.__unit_transforms[0]
 
     def __bool__(self):
         """Return True if this object describes two Var objects which are

test/unit_tests/test_var_transforms.py

@@ -428,27 +428,47 @@ def test_compatible_tendency_variable(self):
         self.assertFalse(compat.has_dim_transforms)
         self.assertFalse(compat.has_unit_transforms)
 
+    def test_compatible_tendency_variable_equivalent_units(self):
+        """Test that a given tendency variable is compatible with
+        its corresponding state variable"""
+        real_array1 = self._new_var('real_stdname1', 'V A',
+                                    ['horizontal_dimension',
+                                     'vertical_layer_dimension'],
+                                     'real', vkind='kind_phys')
+        real_array2 = self._new_var('tendency_of_real_stdname1', 'W s-1',
+                                    ['horizontal_dimension',
+                                     'vertical_layer_dimension'],
+                                     'real', vkind='kind_phys')
+        compat = real_array2.compatible(real_array1, self.__run_env, is_tend=True)
+        self.assertIsInstance(compat, VarCompatObj,
+                              msg=self.__inst_emsg.format(type(compat)))
+        self.assertTrue(compat)
+        self.assertTrue(compat.compat)
+        self.assertEqual(compat.incompat_reason, '')
+        self.assertFalse(compat.has_kind_transforms)
+        self.assertFalse(compat.has_dim_transforms)
+        self.assertFalse(compat.has_unit_transforms)
+
     def test_incompatible_tendency_variable(self):
         """Test that the correct error is returned when a given tendency
         variable has inconsistent units vs the state variable"""
-        real_array1 = self._new_var('real_stdname1', 'C',
+        real_array1 = self._new_var('real_stdname1', 'm',
                                     ['horizontal_dimension',
                                      'vertical_layer_dimension'],
                                      'real', vkind='kind_phys')
-        real_array2 = self._new_var('tendency_of_real_stdname1', 'C kg s-1',
+        real_array2 = self._new_var('tendency_of_real_stdname1', 'cm s-1',
                                     ['horizontal_dimension',
                                      'vertical_layer_dimension'],
                                      'real', vkind='kind_phys')
         compat = real_array2.compatible(real_array1, self.__run_env, is_tend=True)
         self.assertIsInstance(compat, VarCompatObj,
                               msg=self.__inst_emsg.format(type(compat)))
-        #Verify correct error message returned
-        emsg = "\nMismatch tendency variable units 'C kg s-1' for variable 'tendency_of_real_stdname1'. No variable transforms supported for tendencies. Tendency units should be 'C s-1' to match state variable."
+        # Verify correct error message returned
+        emsg = "\nMismatch tendency variable units 'cm s-1' for variable 'tendency_of_real_stdname1'. No variable transforms supported for tendencies. Tendency units should be 'm s-1' to match state variable."
         self.assertEqual(compat.incompat_reason, emsg)
         self.assertFalse(compat.has_kind_transforms)
         self.assertFalse(compat.has_dim_transforms)
         self.assertFalse(compat.has_unit_transforms)
-        #Verify correct error message returned
 
 
 if __name__ == "__main__":

test/var_compatibility_test/effr_calc.F90

@@ -37,7 +37,8 @@ end subroutine effr_calc_init
    !!
    subroutine effr_calc_run(ncol, nlev, effrr_in, effrg_in, ncg_in, nci_out, &
                             effrl_inout, effri_out, effrs_inout, ncl_out,    &
-                            has_graupel, scalar_var, errmsg, errflg)
+                            has_graupel, scalar_var, tke_inout, tke2_inout,  &
+                            errmsg, errflg)
 
       integer,            intent(in)    :: ncol
       integer,            intent(in)    :: nlev
@@ -53,6 +54,9 @@ subroutine effr_calc_run(ncol, nlev, effrr_in, effrg_in, ncg_in, nci_out, &
       character(len=512), intent(out)   :: errmsg
       integer,            intent(out)   :: errflg
       real(kind_phys),    intent(out),optional   :: ncl_out(:,:)
+      real(kind_phys),    intent(inout) :: tke_inout
+      real(kind_phys),    intent(inout) :: tke2_inout
+
       !----------------------------------------------------------------
 
       real(kind_phys), parameter :: re_qc_min = 2.5 ! microns

test/var_compatibility_test/effr_calc.meta

@@ -130,6 +130,22 @@
   type = real
   kind = kind_phys
   intent = inout
+[ tke_inout ]
+  standard_name = turbulent_kinetic_energy
+  long_name = turbulent_kinetic_energy
+  units = m2 s-2
+  dimensions = ()
+  type = real
+  kind = kind_phys
+  intent = inout
+[ tke2_inout ]
+  standard_name = turbulent_kinetic_energy2
+  long_name = turbulent_kinetic_energy2
+  units = m+2 s-2
+  dimensions = ()
+  type = real
+  kind = kind_phys
+  intent = inout
 [ errmsg ]
   standard_name = ccpp_error_message
   long_name = Error message for error handling in CCPP

test/var_compatibility_test/run_test

@@ -152,9 +152,10 @@ required_vars_var_compatibility="${required_vars_var_compatibility},scalar_varia
 required_vars_var_compatibility="${required_vars_var_compatibility},scalar_variable_for_testing_c"
 required_vars_var_compatibility="${required_vars_var_compatibility},scheme_order_in_suite"
 required_vars_var_compatibility="${required_vars_var_compatibility},shortwave_radiation_fluxes"
+required_vars_var_compatibility="${required_vars_var_compatibility},turbulent_kinetic_energy"
+required_vars_var_compatibility="${required_vars_var_compatibility},turbulent_kinetic_energy2"
 required_vars_var_compatibility="${required_vars_var_compatibility},vertical_layer_dimension"
 input_vars_var_compatibility="cloud_graupel_number_concentration"
-#input_vars_var_compatibility="${input_vars_var_compatibility},cloud_ice_number_concentration"
 input_vars_var_compatibility="${input_vars_var_compatibility},effective_radius_of_stratiform_cloud_graupel"
 input_vars_var_compatibility="${input_vars_var_compatibility},effective_radius_of_stratiform_cloud_liquid_water_particle"
 input_vars_var_compatibility="${input_vars_var_compatibility},effective_radius_of_stratiform_cloud_rain_particle"
@@ -172,6 +173,8 @@ input_vars_var_compatibility="${input_vars_var_compatibility},scalar_variable_fo
 input_vars_var_compatibility="${input_vars_var_compatibility},scalar_variable_for_testing_c"
 input_vars_var_compatibility="${input_vars_var_compatibility},scheme_order_in_suite"
 input_vars_var_compatibility="${input_vars_var_compatibility},shortwave_radiation_fluxes"
+input_vars_var_compatibility="${input_vars_var_compatibility},turbulent_kinetic_energy"
+input_vars_var_compatibility="${input_vars_var_compatibility},turbulent_kinetic_energy2"
 input_vars_var_compatibility="${input_vars_var_compatibility},vertical_layer_dimension"
 output_vars_var_compatibility="ccpp_error_code,ccpp_error_message"
 output_vars_var_compatibility="${output_vars_var_compatibility},cloud_ice_number_concentration"
@@ -183,6 +186,8 @@ output_vars_var_compatibility="${output_vars_var_compatibility},longwave_radiati
 output_vars_var_compatibility="${output_vars_var_compatibility},scalar_variable_for_testing"
 output_vars_var_compatibility="${output_vars_var_compatibility},scheme_order_in_suite"
 output_vars_var_compatibility="${output_vars_var_compatibility},shortwave_radiation_fluxes"
+output_vars_var_compatibility="${output_vars_var_compatibility},turbulent_kinetic_energy"
+output_vars_var_compatibility="${output_vars_var_compatibility},turbulent_kinetic_energy2"
 
 ##
 ## Run a database report and check the return string