Revert "Remove add constant method" (#2503)

Author: VeckoTheGecko

docs/user_guide/examples/tutorial_Argofloats.ipynb

@@ -131,7 +131,7 @@
     "# Convert to SGRID-compliant dataset and create FieldSet\n",
     "ds_fset = parcels.convert.copernicusmarine_to_sgrid(fields=fields)\n",
     "fieldset = parcels.FieldSet.from_sgrid_conventions(ds_fset)\n",
-    "fieldset.constants[\"mindepth\"] = 1.0\n",
+    "fieldset.add_constant(\"mindepth\", 1.0)\n",
     "\n",
     "# Define a new Particle type including extra Variables\n",
     "ArgoParticle = parcels.Particle.add_variable(\n",

docs/user_guide/examples/tutorial_croco_3D.ipynb

@@ -91,7 +91,7 @@
     "fieldset = parcels.FieldSet.from_sgrid_conventions(ds_fset)\n",
     "\n",
     "# Add the critical depth (`hc`) as a constant to the fieldset\n",
-    "fieldset.constants[\"hc\"] = ds_fields.hc.item()"
+    "fieldset.add_constant(\"hc\", ds_fields.hc.item())"
    ]
   },
   {
@@ -204,7 +204,7 @@
    "source": [
     "fieldset_noW = parcels.FieldSet.from_sgrid_conventions(ds_fset)\n",
     "fieldset_noW.W.data[:] = 0.0\n",
-    "fieldset_noW.constants[\"hc\"] = ds_fields.hc.item()\n",
+    "fieldset_noW.add_constant(\"hc\", ds_fields.hc.item())\n",
     "\n",
     "X, Z = np.meshgrid(\n",
     "    [40e3, 80e3, 120e3],\n",

docs/user_guide/examples/tutorial_diffusion.ipynb

@@ -90,7 +90,7 @@
     "\n",
     "Just like velocities, diffusivities are passed to Parcels in the form of `Field` objects. When using `DiffusionUniformKh`, they should be added to the `FieldSet` object as constant fields, e.g. `fieldset.add_constant_field(\"Kh_zonal\", 1, mesh=\"flat\")`.\n",
     "\n",
-    "To make a central difference approximation for computing the gradient in diffusivity, a resolution for this approximation `dres` is needed: _Parcels_ approximates the gradients in diffusivities by using their values at the particle's location ± `dres` (in both $x$ and $y$). A value of `dres` must be specified and added to the FieldSet by the user (e.g. `fieldset.constants[\"dres\"] = 0.01`). Currently, it is unclear what the best value of `dres` is. From experience, the size of `dres` should be smaller than the spatial resolution of the data, but within reasonable limits of machine precision to avoid numerical errors. We are working on a method to compute gradients differently so that specifying `dres` is not necessary anymore.\n",
+    "To make a central difference approximation for computing the gradient in diffusivity, a resolution for this approximation `dres` is needed: _Parcels_ approximates the gradients in diffusivities by using their values at the particle's location ± `dres` (in both $x$ and $y$). A value of `dres` must be specified and added to the FieldSet by the user (e.g. `fieldset.add_constant(\"dres\", 0.01)`). Currently, it is unclear what the best value of `dres` is. From experience, the size of `dres` should be smaller than the spatial resolution of the data, but within reasonable limits of machine precision to avoid numerical errors. We are working on a method to compute gradients differently so that specifying `dres` is not necessary anymore.\n",
     "\n",
     "## Example: Impermeable Diffusivity Profile\n",
     "\n",
@@ -206,7 +206,7 @@
    "source": [
     "fieldset = parcels.FieldSet.from_sgrid_conventions(ds, mesh=\"flat\")\n",
     "fieldset.add_constant_field(\"Kh_zonal\", 1, mesh=\"flat\")\n",
-    "fieldset.constants[\"dres\"] = 0.00005"
+    "fieldset.add_constant(\"dres\", 0.00005)"
    ]
   },
   {
@@ -537,7 +537,7 @@
     ")\n",
     "fieldset.add_field(cell_areas_field)\n",
     "\n",
-    "fieldset.constants[\"Cs\"] = 0.1"
+    "fieldset.add_constant(\"Cs\", 0.1)"
    ]
   },
   {

docs/user_guide/v4-migration.md

@@ -21,7 +21,6 @@ Version 4 of Parcels is unreleased at the moment. The information in this migrat
 ## FieldSet
 
 - `interp_method` has to be an Interpolation function, instead of a string.
-- `add_constant` has been removed as a method. Instead, use direct dictionary assignment (`fieldset.constants[key] = value`)
 
 ## Particle
 

src/parcels/_core/fieldset.py

@@ -12,6 +12,7 @@
 
 from parcels._core.field import Field, VectorField
 from parcels._core.utils import sgrid
+from parcels._core.utils.string import _assert_str_and_python_varname
 from parcels._core.utils.time import get_datetime_type_calendar
 from parcels._core.utils.time import is_compatible as datetime_is_compatible
 from parcels._core.uxgrid import UxGrid
@@ -152,6 +153,25 @@ def add_constant_field(self, name: str, value, mesh: Mesh = "spherical"):
         grid = XGrid(xgrid, mesh=mesh)
         self.add_field(Field(name, ds[name], grid, interp_method=XConstantField))
 
+    def add_constant(self, name, value):
+        """Add a constant to the FieldSet.
+
+        Parameters
+        ----------
+        name : str
+            Name of the constant
+        value :
+            Value of the constant
+
+        """
+        _assert_str_and_python_varname(name)
+
+        if name in self.constants:
+            raise ValueError(f"FieldSet already has a constant with name '{name}'")
+        if not isinstance(value, (float, np.floating, int, np.integer)):
+            raise ValueError(f"FieldSet constants have to be of type float or int, got a {type(value)}")
+        self.constants[name] = value
+
     @property
     def gridset(self) -> list[BaseGrid]:
         grids = []

src/parcels/_core/kernel.py

@@ -142,29 +142,29 @@ def check_fieldsets_in_kernels(self, pyfunc):  # TODO v4: this can go into anoth
                     raise ValueError('ParticleClass requires a "next_dt" for AdvectionRK45 Kernel.')
                 if not hasattr(self.fieldset, "RK45_tol"):
                     warnings.warn(
-                        "Setting RK45 tolerance to 10 m. Set fieldset.constants['RK45_tol'] = [distance] to change.",
+                        "Setting RK45 tolerance to 10 m. Use fieldset.add_constant('RK45_tol', [distance]) to change.",
                         KernelWarning,
                         stacklevel=2,
                     )
-                    self.fieldset.constants["RK45_tol"] = 10
+                    self.fieldset.add_constant("RK45_tol", 10)
                 if self.fieldset.U.grid._mesh == "spherical":
                     self.fieldset.RK45_tol /= (
                         1852 * 60
                     )  # TODO does not account for zonal variation in meter -> degree conversion
                 if not hasattr(self.fieldset, "RK45_min_dt"):
                     warnings.warn(
-                        "Setting RK45 minimum timestep to 1 s. Set fieldset.constants['RK45_min_dt'] = [timestep] to change.",
+                        "Setting RK45 minimum timestep to 1 s. Use fieldset.add_constant('RK45_min_dt', [timestep]) to change.",
                         KernelWarning,
                         stacklevel=2,
                     )
-                    self.fieldset.constants["RK45_min_dt"] = 1
+                    self.fieldset.add_constant("RK45_min_dt", 1)
                 if not hasattr(self.fieldset, "RK45_max_dt"):
                     warnings.warn(
-                        "Setting RK45 maximum timestep to 1 day. Set fieldset.constants['RK45_max_dt'] = [timestep] to change.",
+                        "Setting RK45 maximum timestep to 1 day. Use fieldset.add_constant('RK45_max_dt', [timestep]) to change.",
                         KernelWarning,
                         stacklevel=2,
                     )
-                    self.fieldset.constants["RK45_max_dt"] = 60 * 60 * 24
+                    self.fieldset.add_constant("RK45_max_dt", 60 * 60 * 24)
 
     def merge(self, kernel):
         if not isinstance(kernel, type(self)):

tests/test_advection.py

@@ -297,14 +297,14 @@ def test_moving_eddy(kernel, rtol):
         ds["Kh"] = (["time", "depth", "YG", "XG"], np.full(ds["U"].shape, 0))
         fieldset.add_field(Field("Kh", ds["Kh"], grid, interp_method=XLinear), "Kh_zonal")
         fieldset.add_field(Field("Kh", ds["Kh"], grid, interp_method=XLinear), "Kh_meridional")
-        fieldset.constants["dres"] = 0.1
+        fieldset.add_constant("dres", 0.1)
 
     start_lon, start_lat, start_z = 12000, 12500, 12500
     dt = np.timedelta64(30, "m")
     endtime = np.timedelta64(1, "h")
 
     if kernel == AdvectionRK45:
-        fieldset.constants["RK45_tol"] = rtol
+        fieldset.add_constant("RK45_tol", rtol)
 
     pset = ParticleSet(
         fieldset, pclass=DEFAULT_PARTICLES[kernel], lon=start_lon, lat=start_lat, z=start_z, time=np.timedelta64(0, "s")
@@ -346,8 +346,8 @@ def test_decaying_moving_eddy(kernel, rtol):
     endtime = np.timedelta64(23, "h")
 
     if kernel == AdvectionRK45:
-        fieldset.constants["RK45_tol"] = rtol
-        fieldset.constants["RK45_min_dt"] = 10 * 60
+        fieldset.add_constant("RK45_tol", rtol)
+        fieldset.add_constant("RK45_min_dt", 10 * 60)
 
     pset = ParticleSet(
         fieldset, pclass=DEFAULT_PARTICLES[kernel], lon=start_lon, lat=start_lat, time=np.timedelta64(0, "s")
@@ -403,7 +403,7 @@ def test_stommelgyre_fieldset(kernel, rtol, grid_type):
     )
 
     if kernel == AdvectionRK45:
-        fieldset.constants["RK45_tol"] = rtol
+        fieldset.add_constant("RK45_tol", rtol)
 
     def UpdateP(particles, fieldset):  # pragma: no cover
         particles.p = fieldset.P[particles.time, particles.z, particles.lat, particles.lon]
@@ -443,7 +443,7 @@ def test_peninsula_fieldset(kernel, rtol, grid_type):
     )
 
     if kernel == AdvectionRK45:
-        fieldset.constants["RK45_tol"] = rtol
+        fieldset.add_constant("RK45_tol", rtol)
 
     def UpdateP(particles, fieldset):  # pragma: no cover
         particles.p = fieldset.P[particles.time, particles.z, particles.lat, particles.lon]

tests/test_diffusion.py

@@ -64,7 +64,7 @@ def test_fieldKh_SpatiallyVaryingDiffusion(mesh, kernel):
     ds["Kh_zonal"] = (["time", "depth", "YG", "XG"], np.full((2, 1, ydim, xdim), Kh))
     ds["Kh_meridional"] = (["time", "depth", "YG", "XG"], np.full((2, 1, ydim, xdim), Kh))
     fieldset = FieldSet.from_sgrid_conventions(ds, mesh=mesh)
-    fieldset.constants["dres"] = float(ds["lon"][1] - ds["lon"][0])
+    fieldset.add_constant("dres", float(ds["lon"][1] - ds["lon"][0]))
 
     npart = 10000
 
@@ -84,7 +84,7 @@ def test_randomexponential(lambd):
     npart = 1000
 
     # Rate parameter for random.expovariate
-    fieldset.constants["lambd"] = lambd
+    fieldset.add_constant("lambd", lambd)
 
     # Set random seed
     np.random.seed(1234)

tests/test_fieldset.py

@@ -36,6 +36,22 @@ def test_fieldset_init_wrong_types():
         FieldSet([1.0, 2.0, 3.0])
 
 
+def test_fieldset_add_constant(fieldset):
+    fieldset.add_constant("test_constant", 1.0)
+    assert fieldset.test_constant == 1.0
+
+
+def test_fieldset_add_constant_int_name(fieldset):
+    with pytest.raises(TypeError, match="Expected a string for variable name, got int instead."):
+        fieldset.add_constant(123, 1.0)
+
+
+@pytest.mark.parametrize("name", ["a b", "123", "while"])
+def test_fieldset_add_constant_invalid_name(fieldset, name):
+    with pytest.raises(ValueError, match=r"Received invalid Python variable name.*"):
+        fieldset.add_constant(name, 1.0)
+
+
 def test_fieldset_add_constant_field(fieldset):
     fieldset.add_constant_field("test_constant_field", 1.0)
 

tests/test_particlefile.py

@@ -189,7 +189,7 @@ def test_variable_written_once():
 def test_pset_repeated_release_delayed_adding_deleting(fieldset, tmp_zarrfile, dt, maxvar):
     """Tests that if particles are released and deleted based on age that resulting output file is correct."""
     npart = 10
-    fieldset.constants["maxvar"] = maxvar
+    fieldset.add_constant("maxvar", maxvar)
 
     MyParticle = Particle.add_variable(
         [Variable("sample_var", initial=0.0), Variable("v_once", dtype=np.float64, initial=0.0, to_write="once")]