Renaming pyfunc to kernels

docs/getting_started/explanation_concepts.md

@@ -181,7 +181,7 @@ dt = np.timedelta64(5, "m")
 runtime = np.timedelta64(1, "D")
 
 # Run the simulation
-pset.execute(pyfunc=kernels, dt=dt, runtime=runtime)
+pset.execute(kernels=kernels, dt=dt, runtime=runtime)
 ```
 
 ### Output

docs/user_guide/examples/tutorial_interaction.ipynb

@@ -139,7 +139,7 @@
     "]\n",
     "\n",
     "pset.execute(\n",
-    "    pyfunc=kernels,\n",
+    "    kernels=kernels,\n",
     "    runtime=np.timedelta64(60, \"s\"),\n",
     "    dt=np.timedelta64(1, \"s\"),\n",
     "    output_file=output_file,\n",
@@ -331,7 +331,7 @@
     "]\n",
     "\n",
     "pset.execute(\n",
-    "    pyfunc=kernels,\n",
+    "    kernels=kernels,\n",
     "    runtime=np.timedelta64(60, \"s\"),\n",
     "    dt=np.timedelta64(1, \"s\"),\n",
     "    output_file=output_file,\n",

docs/user_guide/examples_v3/tutorial_stommel_uxarray.ipynb

@@ -334,7 +334,7 @@
     "        pset.execute(\n",
     "            endtime=endtime,\n",
     "            dt=timedelta(seconds=60),\n",
-    "            pyfunc=AdvectionEE,\n",
+    "            kernels=AdvectionEE,\n",
     "            verbose_progress=False,\n",
     "        )\n",
     "    except FieldOutOfBoundError:\n",

docs/user_guide/v4-migration.md

@@ -17,6 +17,7 @@ Version 4 of Parcels is unreleased at the moment. The information in this migrat
 - The `InteractionKernel` class has been removed. Since normal Kernels now have access to _all_ particles, particle-particle interaction can be performed within normal Kernels.
 - Users need to explicitly use `convert_z_to_sigma_croco` in sampling kernels (such as the `AdvectionRK4_3D_CROCO` or `SampleOMegaCroco` kernels) when working with CROCO data, as the automatic conversion from depth to sigma grids under the hood has been removed.
 - We added a new AdvectionRK2 Kernel. The AdvectionRK4 kernel is still available, but RK2 is now the recommended default advection scheme as it is faster while the accuracy is comparable for most applications. See also the Choosing an integration method tutorial.
+- Functions shouldn't be converted to Kernels before adding to a pset.execute() call. Instead, simply pass the function(s) as a list to pset.execute(), which will convert them to Kernels internally.
 
 ## FieldSet
 

src/parcels/_core/kernel.py

@@ -49,8 +49,8 @@ class Kernel:
         FieldSet object providing the field information (possibly None)
     ptype :
         PType object for the kernel particle
-    pyfunc :
-        (aggregated) Kernel function
+    kernels :
+        list of Kernel functions
 
     Notes
     -----
@@ -60,32 +60,35 @@ class Kernel:
 
     def __init__(
         self,
+        kernels: list[types.FunctionType],
         fieldset,
         ptype,
-        pyfuncs: list[types.FunctionType],
     ):
-        for f in pyfuncs:
+        if not isinstance(kernels, list):
+            kernels = [kernels]
+
+        for f in kernels:
             if not isinstance(f, types.FunctionType):
-                raise TypeError(f"Argument pyfunc should be a function or list of functions. Got {type(f)}")
+                raise TypeError(f"Argument `kernels` should be a function or list of functions. Got {type(f)}")
             assert_same_function_signature(f, ref=AdvectionRK4, context="Kernel")
 
-        if len(pyfuncs) == 0:
-            raise ValueError("List of `pyfuncs` should have at least one function.")
+        if len(kernels) == 0:
+            raise ValueError("List of `kernels` should have at least one function.")
 
         self._fieldset = fieldset
         self._ptype = ptype
 
         self._positionupdate_kernel_added = False
 
-        for f in pyfuncs:
+        for f in kernels:
             self.check_fieldsets_in_kernels(f)
 
-        self._pyfuncs: list[Callable] = pyfuncs
+        self._kernels: list[Callable] = kernels
 
     @property  #! Ported from v3. To be removed in v4? (/find another way to name kernels in output file)
     def funcname(self):
         ret = ""
-        for f in self._pyfuncs:
+        for f in self._kernels:
             ret += f.__name__
         return ret
 
@@ -123,21 +126,21 @@ def PositionUpdate(particles, fieldset):  # pragma: no cover
                 # Update dt in case it's increased in RK45 kernel
                 particles.dt = particles.next_dt
 
-        self._pyfuncs = (PositionUpdate + self)._pyfuncs
+        self._kernels = (PositionUpdate + self)._kernels
 
-    def check_fieldsets_in_kernels(self, pyfunc):  # TODO v4: this can go into another method? assert_is_compatible()?
+    def check_fieldsets_in_kernels(self, kernel):  # TODO v4: this can go into another method? assert_is_compatible()?
         """
         Checks the integrity of the fieldset with the kernels.
 
-        This function is to be called from the derived class when setting up the 'pyfunc'.
+        This function is to be called from the derived class when setting up the 'kernel'.
         """
         if self.fieldset is not None:
-            if pyfunc is AdvectionAnalytical:
+            if kernel is AdvectionAnalytical:
                 if self._fieldset.U.interp_method != "cgrid_velocity":
                     raise NotImplementedError("Analytical Advection only works with C-grids")
                 if self._fieldset.U.grid._gtype not in [GridType.CurvilinearZGrid, GridType.RectilinearZGrid]:
                     raise NotImplementedError("Analytical Advection only works with Z-grids in the vertical")
-            elif pyfunc is AdvectionRK45:
+            elif kernel is AdvectionRK45:
                 if "next_dt" not in [v.name for v in self.ptype.variables]:
                     raise ValueError('ParticleClass requires a "next_dt" for AdvectionRK45 Kernel.')
                 if not hasattr(self.fieldset, "RK45_tol"):
@@ -174,48 +177,21 @@ def merge(self, kernel):
         assert self.ptype == kernel.ptype, "Cannot merge kernels with different particle types"
 
         return type(self)(
+            self._kernels + kernel._kernels,
             self.fieldset,
             self.ptype,
-            pyfuncs=self._pyfuncs + kernel._pyfuncs,
         )
 
     def __add__(self, kernel):
         if isinstance(kernel, types.FunctionType):
-            kernel = type(self)(self.fieldset, self.ptype, pyfuncs=[kernel])
+            kernel = type(self)([kernel], self.fieldset, self.ptype)
         return self.merge(kernel)
 
     def __radd__(self, kernel):
         if isinstance(kernel, types.FunctionType):
-            kernel = type(self)(self.fieldset, self.ptype, pyfuncs=[kernel])
+            kernel = type(self)([kernel], self.fieldset, self.ptype)
         return kernel.merge(self)
 
-    @classmethod
-    def from_list(cls, fieldset, ptype, pyfunc_list):
-        """Create a combined kernel from a list of functions.
-
-        Takes a list of functions, converts them to kernels, and joins them
-        together.
-
-        Parameters
-        ----------
-        fieldset : parcels.Fieldset
-            FieldSet object providing the field information (possibly None)
-        ptype :
-            PType object for the kernel particle
-        pyfunc_list : list of functions
-            List of functions to be combined into a single kernel.
-        *args :
-            Additional arguments passed to first kernel during construction.
-        **kwargs :
-            Additional keyword arguments passed to first kernel during construction.
-        """
-        if not isinstance(pyfunc_list, list):
-            raise TypeError(f"Argument `pyfunc_list` should be a list of functions. Got {type(pyfunc_list)}")
-        if not all([isinstance(f, types.FunctionType) for f in pyfunc_list]):
-            raise ValueError("Argument `pyfunc_list` should be a list of functions.")
-
-        return cls(fieldset, ptype, pyfunc_list)
-
     def execute(self, pset, endtime, dt):
         """Execute this Kernel over a ParticleSet for several timesteps.
 
@@ -248,7 +224,7 @@ def execute(self, pset, endtime, dt):
                 pset.dt = np.minimum(np.maximum(pset.dt, -time_to_endtime), 0)
 
             # run kernels for all particles that need to be evaluated
-            for f in self._pyfuncs:
+            for f in self._kernels:
                 f(pset[evaluate_particles], self._fieldset)
 
                 # check for particles that have to be repeated

src/parcels/_core/particleset.py

@@ -290,29 +290,6 @@ def from_particlefile(cls, fieldset, pclass, filename, restart=True, restarttime
             "ParticleSet.from_particlefile is not yet implemented in v4."
         )  # TODO implement this when ParticleFile is implemented in v4
 
-    def Kernel(self, pyfunc):
-        """Wrapper method to convert a `pyfunc` into a :class:`parcels.kernel.Kernel` object.
-
-        Conversion is based on `fieldset` and `ptype` of the ParticleSet.
-
-        Parameters
-        ----------
-        pyfunc : function or list of functions
-            Python function to convert into kernel. If a list of functions is provided,
-            the functions will be converted to kernels and combined into a single kernel.
-        """
-        if isinstance(pyfunc, list):
-            return Kernel.from_list(
-                self.fieldset,
-                self._ptype,
-                pyfunc,
-            )
-        return Kernel(
-            self.fieldset,
-            self._ptype,
-            pyfuncs=[pyfunc],
-        )
-
     def data_indices(self, variable_name, compare_values, invert=False):
         """Get the indices of all particles where the value of `variable_name` equals (one of) `compare_values`.
 
@@ -376,7 +353,7 @@ def set_variable_write_status(self, var, write_status):
 
     def execute(
         self,
-        pyfunc,
+        kernels,
         dt: datetime.timedelta | np.timedelta64 | float,
         endtime: np.timedelta64 | np.datetime64 | None = None,
         runtime: datetime.timedelta | np.timedelta64 | float | None = None,
@@ -390,10 +367,9 @@ def execute(
 
         Parameters
         ----------
-        pyfunc :
-            Kernel function to execute. This can be the name of a
+        kernels :
+            List of Kernel functions to execute. This can be the name of a
             defined Python function or a :class:`parcels.kernel.Kernel` object.
-            Kernels can be concatenated using the + operator.
         dt (np.timedelta64 or float):
             Timestep interval (as a np.timedelta64 object of float in seconds) to be passed to the kernel.
             Use a negative value for a backward-in-time simulation.
@@ -417,10 +393,9 @@ def execute(
         if len(self) == 0:
             return
 
-        if not isinstance(pyfunc, Kernel):
-            pyfunc = self.Kernel(pyfunc)
-
-        self._kernel = pyfunc
+        if not isinstance(kernels, Kernel):
+            kernels = Kernel(kernels, self.fieldset, self._ptype)
+        self._kernel = kernels
 
         if output_file is not None:
             output_file.set_metadata(self.fieldset.gridset[0]._mesh)

tests-v3/test_kernel_language.py

@@ -18,7 +18,7 @@
 def expr_kernel(name, pset, expr):
     pycode = (f"def {name}(particle, fieldset, time):\n"
               f"    particle.p = {expr}")  # fmt: skip
-    return Kernel(pset.fieldset, pset.particledata.ptype, pyfunc=None, funccode=pycode, funcname=name)
+    return Kernel(kernels=None, fieldset=pset.fieldset, ptype=pset._ptype, funccode=pycode, funcname=name)
 
 
 @pytest.fixture

tests/test_diffusion.py

@@ -34,7 +34,7 @@ def test_fieldKh_Brownian(mesh):
 
     np.random.seed(1234)
     pset = ParticleSet(fieldset=fieldset, lon=np.zeros(npart), lat=np.zeros(npart))
-    pset.execute(pset.Kernel(DiffusionUniformKh), runtime=runtime, dt=np.timedelta64(1, "h"))
+    pset.execute(DiffusionUniformKh, runtime=runtime, dt=np.timedelta64(1, "h"))
 
     expected_std_lon = np.sqrt(2 * kh_zonal * mesh_conversion**2 * timedelta_to_float(runtime))
     expected_std_lat = np.sqrt(2 * kh_meridional * mesh_conversion**2 * timedelta_to_float(runtime))
@@ -70,7 +70,7 @@ def test_fieldKh_SpatiallyVaryingDiffusion(mesh, kernel):
 
     np.random.seed(1636)
     pset = ParticleSet(fieldset=fieldset, lon=np.zeros(npart), lat=np.zeros(npart))
-    pset.execute(pset.Kernel(kernel), runtime=np.timedelta64(3, "h"), dt=np.timedelta64(1, "h"))
+    pset.execute(kernel, runtime=np.timedelta64(3, "h"), dt=np.timedelta64(1, "h"))
 
     tol = 2000 * mesh_conversion  # effectively 2000 m errors (because of low numbers of particles)
     assert np.allclose(np.mean(pset.lon), 0, atol=tol)

tests/test_kernel.py

@@ -35,22 +35,22 @@ def ErrorKernel(particles, fieldset):  # pragma: no cover
 
 
 def test_kernel_init(fieldset):
-    Kernel(fieldset, ptype=Particle, pyfuncs=[AdvectionRK4])
+    Kernel(kernels=[AdvectionRK4], fieldset=fieldset, ptype=Particle)
 
 
 def test_kernel_merging(fieldset):
-    k1 = Kernel(fieldset, ptype=Particle, pyfuncs=[AdvectionRK4])
-    k2 = Kernel(fieldset, ptype=Particle, pyfuncs=[MoveEast, MoveNorth])
+    k1 = Kernel(kernels=[AdvectionRK4], fieldset=fieldset, ptype=Particle)
+    k2 = Kernel(kernels=[MoveEast, MoveNorth], fieldset=fieldset, ptype=Particle)
 
     merged_kernel = k1 + k2
     assert merged_kernel.funcname == "AdvectionRK4MoveEastMoveNorth"
-    assert len(merged_kernel._pyfuncs) == 3
-    assert merged_kernel._pyfuncs == [AdvectionRK4, MoveEast, MoveNorth]
+    assert len(merged_kernel._kernels) == 3
+    assert merged_kernel._kernels == [AdvectionRK4, MoveEast, MoveNorth]
 
     merged_kernel = k2 + k1
     assert merged_kernel.funcname == "MoveEastMoveNorthAdvectionRK4"
-    assert len(merged_kernel._pyfuncs) == 3
-    assert merged_kernel._pyfuncs == [MoveEast, MoveNorth, AdvectionRK4]
+    assert len(merged_kernel._kernels) == 3
+    assert merged_kernel._kernels == [MoveEast, MoveNorth, AdvectionRK4]
 
 
 def test_kernel_from_list(fieldset):
@@ -61,8 +61,8 @@ def test_kernel_from_list(fieldset):
     mixed functions and kernel objects.
     """
     pset = ParticleSet(fieldset, lon=[0.5], lat=[0.5])
-    kernels_single = pset.Kernel([AdvectionRK4])
-    kernels_functions = pset.Kernel([AdvectionRK4, MoveEast, MoveNorth])
+    kernels_single = Kernel(kernels=[AdvectionRK4], fieldset=fieldset, ptype=pset._ptype)
+    kernels_functions = Kernel(kernels=[AdvectionRK4, MoveEast, MoveNorth], fieldset=fieldset, ptype=pset._ptype)
 
     # Check if the kernels were combined correctly
     assert kernels_single.funcname == "AdvectionRK4"
@@ -77,14 +77,18 @@ def test_kernel_from_list_error_checking(fieldset):
     """
     pset = ParticleSet(fieldset, lon=[0.5], lat=[0.5])
 
-    with pytest.raises(ValueError, match="List of `pyfuncs` should have at least one function."):
-        pset.Kernel([])
+    with pytest.raises(ValueError, match="List of `kernels` should have at least one function."):
+        Kernel(kernels=[], fieldset=fieldset, ptype=pset._ptype)
 
-    with pytest.raises(ValueError, match="Argument `pyfunc_list` should be a list of functions."):
-        pset.Kernel([AdvectionRK4, "something else"])
+    with pytest.raises(TypeError, match=r"Argument `kernels` should be a function or list of functions.*"):
+        Kernel(kernels=[AdvectionRK4, "something else"], fieldset=fieldset, ptype=pset._ptype)
 
-    with pytest.raises(ValueError, match="Argument `pyfunc_list` should be a list of functions."):
-        kernels_mixed = pset.Kernel([pset.Kernel(AdvectionRK4), MoveEast, MoveNorth])
+    with pytest.raises(TypeError, match=r"Argument `kernels` should be a function or list of functions.*"):
+        kernels_mixed = Kernel(
+            kernels=[Kernel(kernels=[AdvectionRK4], fieldset=fieldset, ptype=pset._ptype), MoveEast, MoveNorth],
+            fieldset=fieldset,
+            ptype=pset._ptype,
+        )
         assert kernels_mixed.funcname == "AdvectionRK4MoveEastMoveNorth"
 
 
@@ -93,7 +97,7 @@ def test_RK45Kernel_error_no_next_dt(fieldset):
     pset = ParticleSet(fieldset, lon=[0.5], lat=[0.5])
 
     with pytest.raises(ValueError, match='ParticleClass requires a "next_dt" for AdvectionRK45 Kernel.'):
-        pset.Kernel(AdvectionRK45)
+        Kernel(kernels=AdvectionRK45, fieldset=fieldset, ptype=pset._ptype)
 
 
 def test_kernel_signature(fieldset):
@@ -114,23 +118,23 @@ def kernel_switched_args(fieldset, particle):
     def kernel_with_forced_kwarg(particles, *, fieldset=0):
         pass
 
-    pset.Kernel(good_kernel)
+    Kernel(kernels=good_kernel, fieldset=fieldset, ptype=pset._ptype)
 
     with pytest.raises(ValueError, match="Kernel function must have 2 parameters, got 3"):
-        pset.Kernel(version_3_kernel)
+        Kernel(kernels=version_3_kernel, fieldset=fieldset, ptype=pset._ptype)
 
     with pytest.raises(
         ValueError, match="Parameter 'particle' has incorrect name. Expected 'particles', got 'particle'"
     ):
-        pset.Kernel(version_3_kernel_without_time)
+        Kernel(kernels=version_3_kernel_without_time, fieldset=fieldset, ptype=pset._ptype)
 
     with pytest.raises(
         ValueError, match="Parameter 'fieldset' has incorrect name. Expected 'particles', got 'fieldset'"
     ):
-        pset.Kernel(kernel_switched_args)
+        Kernel(kernels=kernel_switched_args, fieldset=fieldset, ptype=pset._ptype)
 
     with pytest.raises(
         ValueError,
         match="Parameter 'fieldset' has incorrect parameter kind. Expected POSITIONAL_OR_KEYWORD, got KEYWORD_ONLY",
     ):
-        pset.Kernel(kernel_with_forced_kwarg)
+        Kernel(kernels=kernel_with_forced_kwarg, fieldset=fieldset, ptype=pset._ptype)