Add fix for spherical uxgrid mesh type

Author: fluidnumerics-joe

parcels/spatialhash.py

@@ -138,12 +138,12 @@ def __init__(
                 _xbound, _ybound, _zbound = _latlon_rad_to_xyz(np.deg2rad(lat), np.deg2rad(lon))
 
                 # Compute bounding box of each face
-                self._xlow = np.min(_xbound, axis=-1)
-                self._xhigh = np.max(_xbound, axis=-1)
-                self._ylow = np.min(_ybound, axis=-1)
-                self._yhigh = np.max(_ybound, axis=-1)
-                self._zlow = np.min(_zbound, axis=-1)
-                self._zhigh = np.max(_zbound, axis=-1)
+                self._xlow = np.atleast_2d(np.min(_xbound, axis=-1))
+                self._xhigh = np.atleast_2d(np.max(_xbound, axis=-1))
+                self._ylow = np.atleast_2d(np.min(_ybound, axis=-1))
+                self._yhigh = np.atleast_2d(np.max(_ybound, axis=-1))
+                self._zlow = np.atleast_2d(np.min(_zbound, axis=-1))
+                self._zhigh = np.atleast_2d(np.max(_zbound, axis=-1))
 
             else:
                 # Boundaries of the hash grid are the bounding box of the source grid

tests/v4/test_field.py

@@ -185,8 +185,8 @@ def test_field_constant_in_time():
 
     # Assert that the field can be evaluated at any time, and returns the same value
     time = np.datetime64("2000-01-01T00:00:00")
-    P1 = P.eval(time=time, z=10.0, y=30.0, x=30.0, applyConversion=False)
-    P2 = P.eval(time=time + np.timedelta64(1, "D"), z=10.0, y=30.0, x=30.0, applyConversion=False)
+    P1 = P.eval(time=time, z=[10.0], y=[30.0], x=[30.0], applyConversion=False)
+    P2 = P.eval(time=time + np.timedelta64(1, "D"), z=[10.0], y=[30.0], x=[30.0], applyConversion=False)
     assert np.isclose(P1, P2)
 
 

tests/v4/test_uxarray_fieldset.py

@@ -119,10 +119,10 @@ def test_fesom2_square_delaunay_uniform_z_coordinate_eval():
     P = Field(name="p", data=ds.p, grid=UxGrid(ds.uxgrid, z=ds.coords["nz"]), interp_method=UXPiecewiseLinearNode)
     fieldset = FieldSet([UVW, P, UVW.U, UVW.V, UVW.W])
 
-    assert fieldset.U.eval(time=ds.time[0].values, z=1.0, y=30.0, x=30.0, applyConversion=False) == 1.0
-    assert fieldset.V.eval(time=ds.time[0].values, z=1.0, y=30.0, x=30.0, applyConversion=False) == 1.0
-    assert fieldset.W.eval(time=ds.time[0].values, z=1.0, y=30.0, x=30.0, applyConversion=False) == 0.0
-    assert fieldset.p.eval(time=ds.time[0].values, z=1.0, y=30.0, x=30.0, applyConversion=False) == 1.0
+    assert fieldset.U.eval(time=ds.time[0].values, z=[1.0], y=[30.0], x=[30.0], applyConversion=False) == 1.0
+    assert fieldset.V.eval(time=ds.time[0].values, z=[1.0], y=[30.0], x=[30.0], applyConversion=False) == 1.0
+    assert fieldset.W.eval(time=ds.time[0].values, z=[1.0], y=[30.0], x=[30.0], applyConversion=False) == 0.0
+    assert fieldset.p.eval(time=ds.time[0].values, z=[1.0], y=[30.0], x=[30.0], applyConversion=False) == 1.0
 
 
 def test_fesom2_square_delaunay_antimeridian_eval():
@@ -132,10 +132,19 @@ def test_fesom2_square_delaunay_antimeridian_eval():
     Since the underlying data is constant, we can check that the values are as expected.
     """
     ds = datasets_unstructured["fesom2_square_delaunay_antimeridian"]
-    P = Field(name="p", data=ds.p, grid=UxGrid(ds.uxgrid, z=ds.coords["nz"]), interp_method=UXPiecewiseLinearNode)
+    P = Field(
+        name="p",
+        data=ds.p,
+        grid=UxGrid(ds.uxgrid, z=ds.coords["nz"], mesh="spherical"),
+        interp_method=UXPiecewiseLinearNode,
+    )
     fieldset = FieldSet([P])
 
-    assert np.isclose(fieldset.p.eval(time=ds.time[0].values, z=1.0, y=30.0, x=-170.0, applyConversion=False), 1.0)
-    assert np.isclose(fieldset.p.eval(time=ds.time[0].values, z=1.0, y=30.0, x=-180.0, applyConversion=False), 1.0)
-    assert np.isclose(fieldset.p.eval(time=ds.time[0].values, z=1.0, y=30.0, x=180.0, applyConversion=False), 1.0)
-    assert np.isclose(fieldset.p.eval(time=ds.time[0].values, z=1.0, y=30.0, x=170.0, applyConversion=False), 1.0)
+    assert np.isclose(
+        fieldset.p.eval(time=ds.time[0].values, z=[1.0], y=[30.0], x=[-170.0], applyConversion=False), 1.0
+    )
+    assert np.isclose(
+        fieldset.p.eval(time=ds.time[0].values, z=[1.0], y=[30.0], x=[-180.0], applyConversion=False), 1.0
+    )
+    assert np.isclose(fieldset.p.eval(time=ds.time[0].values, z=[1.0], y=[30.0], x=[180.0], applyConversion=False), 1.0)
+    assert np.isclose(fieldset.p.eval(time=ds.time[0].values, z=[1.0], y=[30.0], x=[170.0], applyConversion=False), 1.0)