Fix combining Fields with and without time in one FieldSet (#2741)

* Adding unit test to show combining fields with and without time doesn't currently work

* Support for time_interval=None for individual fields

* Remove time from kernelloop explanation

Author: erikvansebille

docs/user_guide/examples/explanation_kernelloop.md

@@ -59,14 +59,14 @@ import parcels.tutorial
 ds_fields = parcels.tutorial.open_dataset("CopernicusMarine_data_for_Argo_tutorial/data")
 
 # Create an idealised wind field and add it to the dataset
-tdim, ydim, xdim = (len(ds_fields.time),len(ds_fields.latitude), len(ds_fields.longitude))
+ydim, xdim = len(ds_fields.latitude), len(ds_fields.longitude)
 ds_fields["UWind"] = xr.DataArray(
-    data=0.5 * np.ones((tdim, ydim, xdim)) * np.sin(ds_fields.latitude.values - ds_fields.latitude.values.mean())[None, :, None],
-    coords=[ds_fields.time, ds_fields.latitude, ds_fields.longitude])
+    data=0.5 * np.ones((ydim, xdim)) * np.sin(ds_fields.latitude.values - ds_fields.latitude.values.mean())[:, None],
+    coords=[ds_fields.latitude, ds_fields.longitude])
 
 ds_fields["VWind"] = xr.DataArray(
-    data=np.zeros((tdim, ydim, xdim)),
-    coords=[ds_fields.time, ds_fields.latitude, ds_fields.longitude])
+    data=np.zeros((ydim, xdim)),
+    coords=[ds_fields.latitude, ds_fields.longitude])
 
 fields = {
     "U": ds_fields["uo"],

src/parcels/_core/field.py

@@ -109,6 +109,9 @@ def grid(self):  # TODO PR: Remove in favour of referencing model grid directly
 
     @property
     def time_interval(self):  # TODO PR: Remove in favour of referencing model time_interval directly
+        if "time" not in self.model.data[self.name].dims:
+            # This field does not have a time-dimension
+            return None
         return self.model.time_interval
 
     def __repr__(self):

src/parcels/_core/index_search.py

@@ -77,8 +77,8 @@ def _search_time_index(field: Field, time: np.ndarray):
     if field.time_interval is None:
         return {
             "T": {
-                "index": np.zeros(shape=time.shape, dtype=np.int32),
-                "bcoord": np.zeros(shape=time.shape, dtype=np.float32),
+                "index": np.zeros(shape=np.atleast_1d(time).shape, dtype=np.int32),
+                "bcoord": np.zeros(shape=np.atleast_1d(time).shape, dtype=np.float32),
             }
         }
 

tests/test_fieldset.py

@@ -6,6 +6,7 @@
 import numpy as np
 import pandas as pd
 import pytest
+import xarray as xr
 
 import parcels.tutorial
 from parcels import Field, ParticleFile, ParticleSet, XGrid, convert, open_raw_zarr
@@ -360,6 +361,17 @@ def test_fieldset_add():
     assert set(fields_before) == set(fset.fields.keys())
 
 
+def test_vectorfields_without_time():
+    """Test that vector fields without a time dimension can be evaluated."""
+    ds1 = datasets_structured["ds_2d_left"][["U_A_grid", "V_A_grid", "grid"]].rename({"U_A_grid": "U", "V_A_grid": "V"})
+    ds2 = ds1.isel(time=0).drop_vars("time").rename({"U": "U_const", "V": "V_const"})
+    ds = xr.merge([ds1, ds2])
+
+    fset = FieldSet.from_sgrid_conventions(ds, mesh="flat", vector_fields={"UV_const": ("U_const", "V_const")})
+    fset.UV_const.eval(t=0, z=0, y=0, x=0)
+    fset.U_const.eval(t=0, z=0, y=0, x=0)
+
+
 def test_fieldset_add_error_on_duplicate_context_values():
     """Test that adding FieldSets with overlapping context value names raises a ValueError."""
     ds1 = datasets_structured["ds_2d_left"][["U_A_grid", "grid"]].rename({"U_A_grid": "U1"})