diff --git a/docs/user_guide/examples/tutorial_Argofloats.ipynb b/docs/user_guide/examples/tutorial_Argofloats.ipynb index d183577d5..0a37193ce 100644 --- a/docs/user_guide/examples/tutorial_Argofloats.ipynb +++ b/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", diff --git a/docs/user_guide/examples/tutorial_croco_3D.ipynb b/docs/user_guide/examples/tutorial_croco_3D.ipynb index 89491ca26..a84b0cc74 100644 --- a/docs/user_guide/examples/tutorial_croco_3D.ipynb +++ b/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", diff --git a/docs/user_guide/examples/tutorial_diffusion.ipynb b/docs/user_guide/examples/tutorial_diffusion.ipynb index fdc0ef83f..6e42b4ce9 100644 --- a/docs/user_guide/examples/tutorial_diffusion.ipynb +++ b/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)" ] }, { diff --git a/docs/user_guide/v4-migration.md b/docs/user_guide/v4-migration.md index a4026b7bf..7db6451cc 100644 --- a/docs/user_guide/v4-migration.md +++ b/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 diff --git a/src/parcels/_core/fieldset.py b/src/parcels/_core/fieldset.py index 86cf7e0be..fdb47f8e9 100644 --- a/src/parcels/_core/fieldset.py +++ b/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 = [] diff --git a/src/parcels/_core/kernel.py b/src/parcels/_core/kernel.py index 90bdd9038..0f6493e44 100644 --- a/src/parcels/_core/kernel.py +++ b/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)): diff --git a/tests/test_advection.py b/tests/test_advection.py index a1dbaa849..c5d6a9ebf 100644 --- a/tests/test_advection.py +++ b/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] diff --git a/tests/test_diffusion.py b/tests/test_diffusion.py index d9bfb19e0..873bfea22 100644 --- a/tests/test_diffusion.py +++ b/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) diff --git a/tests/test_fieldset.py b/tests/test_fieldset.py index c75a3e3da..8056d4617 100644 --- a/tests/test_fieldset.py +++ b/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) diff --git a/tests/test_particlefile.py b/tests/test_particlefile.py index 1f8ed550b..8d559eb45 100755 --- a/tests/test_particlefile.py +++ b/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")] diff --git a/tests/test_sigmagrids.py b/tests/test_sigmagrids.py index f87a81c29..de437c8fb 100644 --- a/tests/test_sigmagrids.py +++ b/tests/test_sigmagrids.py @@ -31,7 +31,7 @@ def test_conversion_3DCROCO(): ds_fset = parcels.convert.croco_to_sgrid(fields=fields, coords=ds_fields) fieldset = parcels.FieldSet.from_sgrid_conventions(ds_fset) - fieldset.constants["hc"] = ds_fields.hc.item() + fieldset.add_constant("hc", ds_fields.hc.item()) s_xroms = np.array([-1.0, -0.9, -0.8, -0.7, -0.6, -0.5, -0.4, -0.3, -0.2, -0.1, 0.0], dtype=np.float32) z_xroms = np.array([-1.26000000e02, -1.10585846e02, -9.60985413e01, -8.24131317e01, -6.94126511e01, -5.69870148e01, -4.50318756e01, -3.34476166e01, -2.21383114e01, -1.10107975e01, 2.62768921e-02,], dtype=np.float32,) # fmt: skip @@ -63,7 +63,7 @@ def test_advection_3DCROCO(): ds_fset = parcels.convert.croco_to_sgrid(fields=fields, coords=ds_fields) fieldset = parcels.FieldSet.from_sgrid_conventions(ds_fset) - fieldset.constants["hc"] = ds_fields.hc.item() + fieldset.add_constant("hc", ds_fields.hc.item()) runtime = 10_000 X, Z = np.meshgrid([40e3, 80e3, 120e3], [-10, -130])