Bugfix array equality

docs/src/whatsnew/3.12.rst

@@ -215,6 +215,12 @@ v3.12.2 (09 May 2025)
 #. `@trexfeathers`_ refactored Iris loading and saving to make it compatible
    with Dask version ``2025.4.0`` and above. (:pull:`6451`)
 
+#. `@trexfeathers`_ and `@ukmo-ccbunney`_ adapted array comparison in response
+   to NumPy v1.25 deprecating comparison of un-broadcastable arrays. It is
+   hoped that users will see no difference in behaviour, but please get in touch
+   if you notice anything. See `NumPy v1.25 expired deprecations`_ and
+   `numpy#22707`_ for more. (:pull:`6665`)
+
 
 📚 Documentation
 ================
@@ -271,3 +277,5 @@ v3.12.2 (09 May 2025)
 
 .. _SPEC 0: https://scientific-python.org/specs/spec-0000/
 .. _Running setuptools commands: https://setuptools.pypa.io/en/latest/deprecated/commands.html
+.. _NumPy v1.25 expired deprecations: https://numpy.org/doc/stable/release/1.25.0-notes.html#expired-deprecations
+.. _numpy#22707: https://github.com/numpy/numpy/pull/22707

lib/iris/_concatenate.py

@@ -16,6 +16,7 @@
 from xxhash import xxh3_64
 
 from iris._lazy_data import concatenate as concatenate_arrays
+from iris.common.metadata import hexdigest
 import iris.coords
 from iris.coords import AncillaryVariable, AuxCoord, CellMeasure, DimCoord
 import iris.cube
@@ -786,7 +787,7 @@ def _coordinate_differences(self, other, attr, reason="metadata"):
             diff_names = []
             for self_key, self_value in self_dict.items():
                 other_value = other_dict[self_key]
-                if self_value != other_value:
+                if hexdigest(self_value) != hexdigest(other_value):
                     diff_names.append(self_key)
             result = (
                 " " + reason,

lib/iris/_constraints.py

@@ -531,6 +531,8 @@ def __init__(self, **attributes):
         super().__init__(cube_func=self._cube_func)
 
     def __eq__(self, other):
+        # Note: equality means that NumPy arrays are not supported for
+        #  AttributeConstraints (get the truth ambiguity error).
         eq = (
             isinstance(other, AttributeConstraint)
             and self._attributes == other._attributes
@@ -553,6 +555,8 @@ def _cube_func(self, cube):
                         match = False
                         break
                 else:
+                    # Note: equality means that NumPy arrays are not supported
+                    #  for AttributeConstraints (get the truth ambiguity error).
                     if cube_attr != value:
                         match = False
                         break

lib/iris/coords.py

@@ -29,6 +29,7 @@
     DimCoordMetadata,
     metadata_manager_factory,
 )
+from iris.common.metadata import hexdigest
 import iris.exceptions
 import iris.time
 import iris.util
@@ -760,7 +761,7 @@ def is_compatible(self, other, ignore=None):
                     ignore = (ignore,)
                 common_keys = common_keys.difference(ignore)
             for key in common_keys:
-                if np.any(self.attributes[key] != other.attributes[key]):
+                if hexdigest(self.attributes[key]) != hexdigest(other.attributes[key]):
                     compatible = False
                     break
 

lib/iris/cube.py

@@ -40,7 +40,7 @@
 import iris.aux_factory
 from iris.aux_factory import AuxCoordFactory
 from iris.common import CFVariableMixin, CubeMetadata, metadata_manager_factory
-from iris.common.metadata import CoordMetadata, metadata_filter
+from iris.common.metadata import CoordMetadata, hexdigest, metadata_filter
 from iris.common.mixin import LimitedAttributeDict
 import iris.coord_systems
 import iris.coords
@@ -1435,7 +1435,7 @@ def is_compatible(
                     ignore = (ignore,)
                 common_keys = common_keys.difference(ignore)
             for key in common_keys:
-                if np.any(self.attributes[key] != other.attributes[key]):
+                if hexdigest(self.attributes[key]) != hexdigest(other.attributes[key]):
                     compatible = False
                     break
 

lib/iris/fileformats/_structured_array_identification.py

@@ -111,7 +111,9 @@ def __eq__(self, other):
 
         result = NotImplemented
         if stride is not None or arr is not None:
-            result = stride == self.stride and np.all(self.unique_ordered_values == arr)
+            result = stride == self.stride and np.array_equal(
+                self.unique_ordered_values, arr
+            )
         return result
 
     def __ne__(self, other):
@@ -284,7 +286,7 @@ def from_array(cls, arr):
             # Do one last sanity check - does the array we've just described
             # actually compute the correct array?
             constructed_array = structure.construct_array(arr.size)
-            if not np.all(constructed_array == arr):
+            if not np.array_equal(constructed_array, arr):
                 structure = None
 
         return structure

lib/iris/fileformats/netcdf/saver.py

@@ -41,6 +41,7 @@
     OceanSigmaFactory,
     OceanSigmaZFactory,
 )
+from iris.common import hexdigest
 import iris.config
 import iris.coord_systems
 import iris.coords
@@ -2590,12 +2591,7 @@ def save(
         def attr_values_equal(val1, val2):
             # An equality test which also works when some values are numpy arrays (!)
             # As done in :meth:`iris.common.mixin.LimitedAttributeDict.__eq__`.
-            match = val1 == val2
-            try:
-                match = bool(match)
-            except ValueError:
-                match = match.all()
-            return match
+            return hexdigest(val1) == hexdigest(val2)
 
         cube0 = cubes[0]
         invalid_globals = set(
@@ -2682,7 +2678,7 @@ def attr_values_equal(val1, val2):
             common_keys.intersection_update(keys)
             different_value_keys = []
             for key in common_keys:
-                if np.any(attributes[key] != cube.attributes[key]):
+                if hexdigest(attributes[key]) != hexdigest(cube.attributes[key]):
                     different_value_keys.append(key)
             common_keys.difference_update(different_value_keys)
             local_keys.update(different_value_keys)

lib/iris/fileformats/pp_load_rules.py

@@ -139,8 +139,15 @@ def _convert_vertical_coords(
         )
         coords_and_dims.append((coord, dim))
 
+    # Common calc for Depth
+    try:
+        svd_lev_eq = brsvd1 == brlev
+    except ValueError:
+        # In case of broadcasting errors.
+        svd_lev_eq = False
+
     # Depth - unbound.
-    if (len(lbcode) != 5) and (lbvc == 2) and np.all(brsvd1 == brlev):
+    if (len(lbcode) != 5) and (lbvc == 2) and np.all(svd_lev_eq):
         coord = _dim_or_aux(
             blev,
             standard_name="depth",
@@ -150,7 +157,7 @@ def _convert_vertical_coords(
         coords_and_dims.append((coord, dim))
 
     # Depth - bound.
-    if (len(lbcode) != 5) and (lbvc == 2) and np.all(brsvd1 != brlev):
+    if (len(lbcode) != 5) and (lbvc == 2) and np.all(~svd_lev_eq):
         coord = _dim_or_aux(
             blev,
             standard_name="depth",
@@ -164,10 +171,10 @@ def _convert_vertical_coords(
     if (
         (len(lbcode) != 5)
         and (lbvc == 2)
-        and (np.any(brsvd1 == brlev) and np.any(brsvd1 != brlev))
+        and (np.any(svd_lev_eq) and np.any(~svd_lev_eq))
     ):
-        lower = np.where(brsvd1 == brlev, blev, brsvd1)
-        upper = np.where(brsvd1 == brlev, blev, brlev)
+        lower = np.where(svd_lev_eq, blev, brsvd1)
+        upper = np.where(svd_lev_eq, blev, brlev)
         coord = _dim_or_aux(
             blev,
             standard_name="depth",
@@ -189,7 +196,7 @@ def _convert_vertical_coords(
                 units="1",
             )
             coords_and_dims.append((coord, dim))
-        elif np.any(brsvd1 != brlev):
+        elif np.any(~svd_lev_eq):
             # UM populates metadata CORRECTLY,
             # so treat it as the expected (bounded) soil depth.
             coord = _dim_or_aux(

lib/iris/tests/unit/fileformats/pp_load_rules/test__convert_vertical_coords.py

@@ -8,6 +8,7 @@
 """
 
 import numpy as np
+import pytest
 
 from iris.aux_factory import HybridHeightFactory, HybridPressureFactory
 from iris.coords import AuxCoord, DimCoord
@@ -280,6 +281,29 @@ def test_unbounded__vector_no_depth(self):
             dim=0,
         )
 
+    def test_unbounded_incompatible_vectors(self):
+        # Confirm this is not vulnerable to the non-broadcastable error.
+        lblev = [1, 2, 3]
+        blev = [10, 20, 30]
+        brsvd1 = [5, 15, 25, 35]
+        brlev = [5, 15, 25]
+        avoided_error = "operands could not be broadcast together"
+        try:
+            self._check_depth(
+                _lbcode(1),
+                lblev=lblev,
+                blev=blev,
+                brsvd1=brsvd1,
+                brlev=brlev,
+                expect_bounds=False,
+                dim=1,
+            )
+        except ValueError as err:
+            if avoided_error in str(err):
+                message = f'Test failed to avoid specified error: "{err}"'
+                pytest.fail(message)
+            pass
+
     def test_bounded(self):
         self._check_depth(
             _lbcode(1), lblev=23.0, brlev=22.5, brsvd1=23.5, expect_bounds=True
@@ -300,6 +324,29 @@ def test_bounded__vector(self):
             dim=1,
         )
 
+    def test_bounded_incompatible_vectors(self):
+        # Confirm this is not vulnerable to the non-broadcastable error.
+        lblev = [1, 2, 3]
+        blev = [10, 20, 30]
+        brsvd1 = [5, 15, 25, 35]
+        brlev = [15, 25, 35]
+        avoided_error = "operands could not be broadcast together"
+        try:
+            self._check_depth(
+                _lbcode(1),
+                lblev=lblev,
+                blev=blev,
+                brsvd1=brsvd1,
+                brlev=brlev,
+                expect_bounds=True,
+                dim=1,
+            )
+        except ValueError as err:
+            if avoided_error in str(err):
+                message = f'Test failed to avoid specified error: "{err}"'
+                pytest.fail(message)
+            pass
+
     def test_cross_section(self):
         self._check_depth(_lbcode(ix=1, iy=2), lblev=23.0, expect_match=False)
 
@@ -360,6 +407,37 @@ def test_normal__vector(self):
         lblev = np.arange(10)
         self._check_soil_level(_lbcode(0), lblev=lblev, dim=0)
 
+    def test_normal_incompatible_vectors(self):
+        # Confirm this is not vulnerable to the non-broadcastable error.
+        lbvc = 6
+        stash = STASH(1, 1, 1)
+        lbcode = _lbcode(0)
+        lblev = np.arange(10)
+        brsvd1 = [1] * len(lblev)
+        brlev = brsvd1 + [1]
+        blev, bhlev, bhrlev, brsvd2 = None, None, None, None
+
+        avoided_error = "operands could not be broadcast together"
+        try:
+            _ = _convert_vertical_coords(
+                lbcode=lbcode,
+                lbvc=lbvc,
+                blev=blev,
+                lblev=lblev,
+                stash=stash,
+                bhlev=bhlev,
+                bhrlev=bhrlev,
+                brsvd1=brsvd1,
+                brsvd2=brsvd2,
+                brlev=brlev,
+                dim=0,
+            )
+        except ValueError as err:
+            if avoided_error in str(err):
+                message = f'Test failed to avoid specified error: "{err}"'
+                pytest.fail(message)
+            pass
+
     def test_cross_section(self):
         self._check_soil_level(_lbcode(ix=1, iy=2), expect_match=False)
 

lib/iris/tests/unit/fileformats/structured_array_identification/test_ArrayStructure.py

@@ -93,6 +93,13 @@ def test_irregular_3d(self):
         a[0, 0, 0] = 5
         assert self.struct_from_arr(a) is None
 
+    def test_irregular_1d(self):
+        # Note this tests a unique scenario where NumPy raises a broadcasting
+        #  error (having previously allowed elementwise comparison in earlier
+        #  versions).
+        a = np.array([0, 0, 0, 0, 1, 1, 1, 1, 2])
+        assert self.struct_from_arr(a) is None
+
     def test_repeated_3d(self):
         sub = np.array([-1, 3, 1, 2])
         a = construct_nd(sub, 2, (3, 2, 4))
@@ -121,6 +128,11 @@ def test_multi_dim_array(self):
         with pytest.raises(ValueError):
             ArrayStructure.from_array(np.arange(12).reshape(3, 4))
 
+    def test_eq_incompatible_shapes(self):
+        struct1 = ArrayStructure(1, np.array([1, 2]))
+        struct2 = ArrayStructure(1, np.array([1, 2, 3]))
+        assert struct1 != struct2
+
 
 class TestNdarrayAndDimsCases:
     """Defines the test functionality for nd_array_and_dims. This class

lib/iris/tests/unit/util/test_describe_diff.py

@@ -58,3 +58,10 @@ def test_different_array_attributes(self):
                 "incompatible_array_attrs.str.txt",
             ],
         )
+
+    def test_incompatible_array_attributes(self):
+        # test incompatible array attribute
+        self.cube_a.attributes["test_array"] = np.array([1, 2, 3])
+        self.cube_b.attributes["test_array"] = np.array([1, 2])
+        with pytest.raises(ValueError, match="Error comparing test_array attributes"):
+            _ = self._compare_result(self.cube_a, self.cube_b)

lib/iris/tests/unit/util/test_equalise_attributes.py

@@ -129,6 +129,18 @@ def test_array_same(self):
         cubes = [self.cube_a1b5v1, self.cube_a1b6v1]
         self._test(cubes, {"a": 1, "v": self.v1}, [{"b": 5}, {"b": 6}])
 
+    @pytest.fixture
+    def make_a1b6v2_incompatible(self):
+        v_array = self.cube_a1b6v2.attributes["v"]
+        self.cube_a1b6v2.attributes["v"] = np.repeat(v_array, 2)
+
+    def test_array_incompatible(self, make_a1b6v2_incompatible):
+        cubes = [self.cube_a1b5v1, self.cube_a1b6v2]
+        with pytest.raises(
+            ValueError, match=r"Error comparing \('local', 'v'\) attributes"
+        ):
+            self._test(cubes, {}, [])
+
     @_shared_utils.skip_data
     def test_complex_nonecommon(self):
         # Example with cell methods and factories, but no common attributes.

lib/iris/util.py

@@ -232,7 +232,12 @@ def describe_diff(cube_a, cube_b, output_file=None):
     else:
         common_keys = set(cube_a.attributes).intersection(cube_b.attributes)
         for key in common_keys:
-            if np.any(cube_a.attributes[key] != cube_b.attributes[key]):
+            try:
+                eq = np.any(cube_a.attributes[key] == cube_b.attributes[key])
+            except ValueError as err:
+                message = f"Error comparing {key} attributes: {err}"
+                raise ValueError(message)
+            if not eq:
                 output_file.write(
                     '"%s" cube_a attribute value "%s" is not '
                     "compatible with cube_b "
@@ -427,7 +432,12 @@ def _masked_array_equal(
             else:
                 ignore |= nanmask
 
-        eqs = ma.getdata(array1) == ma.getdata(array2)
+        try:
+            eqs = ma.getdata(array1) == ma.getdata(array2)
+        except ValueError:
+            # In case of broadcasting errors.
+            eqs = np.zeros(array1.shape, dtype=bool)
+
         if ignore is not None:
             eqs = np.where(ignore, True, eqs)
 
@@ -2106,11 +2116,15 @@ def equalise_attributes(cubes):
     for attrs in cube_attrs[1:]:
         cube_keys = list(attrs.keys())
         keys_to_remove.update(cube_keys)
-        common_keys = [
-            key
-            for key in common_keys
-            if (key in cube_keys and np.all(attrs[key] == cube_attrs[0][key]))
-        ]
+
+        def eq(key):
+            try:
+                return np.all(attrs[key] == cube_attrs[0][key])
+            except ValueError as err:
+                message = f"Error comparing {key} attributes: {err}"
+                raise ValueError(message)
+
+        common_keys = [key for key in common_keys if (key in cube_keys and eq(key))]
     keys_to_remove.difference_update(common_keys)
 
     # Convert back from the resulting 'paired' keys set, extracting just the