{
"cells": [
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# Getting started with Parcels: general structure"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"The flexibility of Parcels allows for a wide range of applications and to build complex simulations. In order to help structuring your code, this tutorial describes the structure that a Parcels script uses.\n",
"\n",
"Code that uses Parcels is generally build up from four different components:\n",
"1. [**FieldSet**](#1.-FieldSet). Load and set up the fields. These can be velocity fields that are used to advect the particles, but it can also be e.g. temperature.\n",
"2. [**ParticleSet**](#2.-ParticleSet). Define the type of particles. Also additional `Variables` can be added to the particles (e.g. temperature, to keep track of the temperature that particles experience).\n",
"3. [**Kernels**](#3.-Kernels). Define and compile kernels. Kernels perform some specific operation on the particles every time step (e.g. interpolate the temperature from the temperature field to the particle location).\n",
"4. [**Execution and output**](#4.-Execution-and-Output). Execute the simulation and write and store the output in a NetCDF file.\n",
"\n",
"We discuss each component in more detail below."
]
},
{
"cell_type": "code",
"execution_count": 1,
"metadata": {},
"outputs": [
{
"data": {
"image/svg+xml": [
""
],
"text/plain": [
""
]
},
"execution_count": 1,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"from IPython.display import SVG\n",
"# Illustration of Parcels code structure\n",
"SVG(filename='parcels_user_diagram.svg') "
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# 1. FieldSet"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Parcels provides a framework to simulate the movement of particles **within an existing flow field environment**. To start a parcels simulation we must define this environment with the [**`FieldSet`** class](https://oceanparcels.org/gh-pages/html/#module-parcels.fieldset). The minimal requirements for this Fieldset are that it must contain the `'U'` and `'V'` fields: the 2D hydrodynamic data that will move the particles in a horizontal direction. Additionally, it can contain e.g. a temperature or vertical flow field."
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"A fieldset can be loaded with [**`FieldSet.from_netcdf`**](https://oceanparcels.org/gh-pages/html/#parcels.fieldset.FieldSet.from_netcdf), if the model output with the fields is written in NetCDF files. This function requires `filenames`, `variables` and `dimensions`. "
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"In this example, we only load the `'U'` and `'V'` fields, which represent the zonal and meridional flow velocity. First, `fname` points to the location of the model output. "
]
},
{
"cell_type": "code",
"execution_count": 2,
"metadata": {},
"outputs": [],
"source": [
"fname = 'GlobCurrent_example_data/*.nc'\n",
"filenames = {'U': fname, 'V': fname}"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Second, we have to specify the `'U'` and `'V'` variable names, as given by the NetCDF files."
]
},
{
"cell_type": "code",
"execution_count": 3,
"metadata": {},
"outputs": [],
"source": [
"variables = {'U': 'eastward_eulerian_current_velocity', 'V': 'northward_eulerian_current_velocity'}"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Third, we specify the names of the variable dimensions, as given by the NetCDF files."
]
},
{
"cell_type": "code",
"execution_count": 4,
"metadata": {},
"outputs": [],
"source": [
"# In the GlobCurrent data the dimensions are also called 'lon', 'lat' and 'time\n",
"dimensions = {'U': {'lat': 'lat', 'lon': 'lon', 'time': 'time'},\n",
" 'V': {'lat': 'lat', 'lon': 'lon', 'time': 'time'}}"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Finally, we load the fieldset."
]
},
{
"cell_type": "code",
"execution_count": 5,
"metadata": {},
"outputs": [],
"source": [
"from parcels import FieldSet\n",
"fieldset = FieldSet.from_netcdf(filenames, variables, dimensions)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"#### For more advanced tutorials on creating `FieldSets`:\n",
"\n",
"* [Implement **periodic boundaries**](https://nbviewer.jupyter.org/github/OceanParcels/parcels/blob/master/parcels/examples/tutorial_periodic_boundaries.ipynb)\n",
"* [How to **interpolate** field data for different fields](https://nbviewer.jupyter.org/github/OceanParcels/parcels/blob/master/parcels/examples/tutorial_interpolation.ipynb)\n",
"* [**Converting units** in the field data](https://nbviewer.jupyter.org/github/OceanParcels/parcels/blob/master/parcels/examples/tutorial_unitconverters.ipynb)\n",
"* [Working around incompatible **time coordinates**](https://nbviewer.jupyter.org/github/OceanParcels/parcels/blob/master/parcels/examples/tutorial_timestamps.ipynb)\n",
"\n",
"#### If you are working with field data on different grids:\n",
"* [Grid **indexing** on different grids](https://nbviewer.jupyter.org/github/OceanParcels/parcels/blob/master/parcels/examples/documentation_indexing.ipynb)\n",
"* [Load field data from **Curvilinear grids**](https://nbviewer.jupyter.org/github/OceanParcels/parcels/blob/master/parcels/examples/tutorial_nemo_curvilinear.ipynb)\n",
"* [Load field data from **3D C-grids**](https://nbviewer.jupyter.org/github/OceanParcels/parcels/blob/master/parcels/examples/tutorial_nemo_3D.ipynb)\n",
"\n",
"#### If you want to combine different velocity fields:\n",
"* [Add different velocity fields in a **SummedField**](https://nbviewer.jupyter.org/github/OceanParcels/parcels/blob/master/parcels/examples/tutorial_SummedFields.ipynb)\n",
"* [Nest velocity fields of different regions or resolutions in a **NestedField**](https://nbviewer.jupyter.org/github/OceanParcels/parcels/blob/master/parcels/examples/tutorial_NestedFields.ipynb)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# 2. ParticleSet"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Once the environment has a `FieldSet` object, you can start defining your particles in a [**`ParticleSet`** object](https://oceanparcels.org/gh-pages/html/#module-parcels.particleset). This object requires:\n",
"1. The [**`FieldSet`**](https://oceanparcels.org/gh-pages/html/#module-parcels.fieldset) object in which the particles will be released.\n",
"2. The type of [**`Particle`**](https://oceanparcels.org/gh-pages/html/#module-parcels.particle): Either a JITParticle or ScipyParticle.\n",
"3. Initial conditions for each [**`Variable`**](https://oceanparcels.org/gh-pages/html/#parcels.particle.Variable) defined in the `Particle`, most notably the release locations in `lon` and `lat`."
]
},
{
"cell_type": "code",
"execution_count": 6,
"metadata": {},
"outputs": [],
"source": [
"from parcels import Variable, JITParticle, ParticleSet\n",
"\n",
"# Define a new particle class\n",
"class AgeParticle(JITParticle): # It is a JIT particle\n",
" age = Variable('age',\n",
" initial=0) # Variable 'age' is added with initial value 0.\n",
"\n",
"pset = ParticleSet(fieldset=fieldset, # the fields that the particleset uses\n",
" pclass=AgeParticle, # define the type of particle\n",
" lon=29, # release longitude\n",
" lat=-33) # release latitude"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"#### For more advanced tutorials on how to setup your `ParticleSet`:\n",
"* [**Releasing particles** at different times](https://nbviewer.jupyter.org/github/OceanParcels/parcels/blob/master/parcels/examples/tutorial_delaystart.ipynb)\n",
"* [The difference between **JITParticles and ScipyParticles**](https://nbviewer.jupyter.org/github/OceanParcels/parcels/blob/master/parcels/examples/tutorial_jit_vs_scipy.ipynb)\n",
"\n",
"For more information on how to implement **`Particle` types** with specific behaviour, see the [section on writing your own kernels](#For-more-advanced-tutorials-on-writing-custom-kernels-that-work-on-custom-particles:)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# 3. Kernels"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Kernels are little snippets of code, which are applied to every particle in the `ParticleSet`, for every time step during a simulation.\n",
"Basic kernels are [included in Parcels](https://oceanparcels.org/gh-pages/html/#module-parcels.kernels.advection), among which several types of advection kernels. `AdvectionRK4` is the main advection kernel for two-dimensional advection, which is also used in this example.\n",
"\n",
"One can also write custom kernels, to add certain types of 'behaviour' to the particles."
]
},
{
"cell_type": "code",
"execution_count": 7,
"metadata": {},
"outputs": [],
"source": [
"# load basic advection kernels\n",
"from parcels import AdvectionRK4\n",
"\n",
"# Create a custom kernel which displaces each particle southward\n",
"def NorthVel(particle, fieldset, time):\n",
" if time > 10*86400 and time < 10.2*86400:\n",
" vvel = -1e-4\n",
" particle.lat += vvel * particle.dt\n",
" \n",
"# Create a custom kernel which keeps track of the particle age (minutes)\n",
"def Age(particle, fieldset, time):\n",
" particle.age += particle.dt / 3600\n",
"\n",
"# define all kernels to be executed on particles\n",
"kernels = pset.Kernel(Age) + pset.Kernel(NorthVel) + AdvectionRK4"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"\n",
"Some key limitations exist to the Kernels that everyone who wants to write their own should be aware of:\n",
"* Every Kernel must be a function with the following (and only those) arguments: `(particle, fieldset, time)`\n",
"* In order to run successfully in JIT mode, Kernel definitions can only contain the following types of commands:\n",
" * Basic arithmetical operators (`+`, `-`, `*`, `/`, `**`) and assignments (`=`).\n",
" * Basic logical operators (`<`, `==`, `!=`, `>`, `&`, `|`). Note that you can use a statement like `particle.lon != particle.lon` to check if `particle.lon` is NaN (since `math.nan != math.nan`).\n",
" * `if` and `while` loops, as well as `break` statements. Note that `for`-loops are not supported in JIT mode.\n",
" * Interpolation of a `Field` from the `FieldSet` at a `[time, depth, lat, lon]` point, using square brackets notation.\n",
" For example, to interpolate the zonal velocity (U) field at the particle location, use the following statement:\n",
" `value = fieldset.U[time, particle.depth, particle.lat, particle.lon]`\n",
" * Functions from the `math` standard library and from the custom `ParcelsRandom` library at [`parcels.rng`](https://oceanparcels.org/gh-pages/html/#module-parcels.rng)\n",
" * Simple `print` statements, such as:\n",
" * `print(\"Some print\")`\n",
" * `print(particle.lon)`\n",
" * `print(\"particle id: %d\" % particle.id)`\n",
" * `print(\"lon: %f, lat: %f\" % (particle.lon, particle.lat))`\n",
" \n",
" Although note that these `print` statements are not shown in Jupyter notebooks in JIT mode, see [this long-standing Issue](https://github.com/OceanParcels/parcels/issues/369).\n",
" * Local variables can be used in Kernels, and these variables will be accessible in all concatenated Kernels. Note that these local variables are not shared between particles, and also not between time steps.\n",
" * Note that one has to be careful with writing kernels for vector fields on Curvilinear grids. While Parcels automatically rotates the U and V field when necessary, this is not the case for for example wind data. In that case, a custom rotation function will have to be written."
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"#### For more advanced tutorials on writing custom kernels that work on custom particles:\n",
"* [Sample other fields like temperature](https://nbviewer.jupyter.org/github/OceanParcels/parcels/blob/master/parcels/examples/tutorial_sampling.ipynb).\n",
"* [Mimic the behavior of ARGO floats](https://nbviewer.jupyter.org/github/OceanParcels/parcels/blob/master/parcels/examples/tutorial_Argofloats.ipynb).\n",
"* [Adding diffusion to approximate subgrid-scale processes and unresolved physics](https://nbviewer.jupyter.org/github/OceanParcels/parcels/blob/master/parcels/examples/tutorial_diffusion.ipynb).\n",
"* [Converting between units in m/s and degree/s](https://nbviewer.jupyter.org/github/OceanParcels/parcels/blob/master/parcels/examples/tutorial_unitconverters.ipynb)."
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# 4. Execution and output"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"The final part executes the simulation, given the `ParticleSet`, `FieldSet` and `Kernels`, that have been defined in the previous steps. If you like to store the particle data generated in the simulation, you define the [**`ParticleFile`**](https://oceanparcels.org/gh-pages/html/#parcels.particlefile.ParticleFile) to which the output of the kernel execution will be written. Then, on the [`ParticleSet`](#ParticleSet) you have defined, you can use the method [**`ParticleSet.execute()`**](https://oceanparcels.org/gh-pages/html/#parcels.particleset.ParticleSet.execute) which requires the following arguments:\n",
"1. The kernels to be executed.\n",
"2. The `runtime` defining how long the execution loop runs. Alternatively, you may define the `endtime` at which the execution loop stops.\n",
"3. The timestep `dt` at which to execute the kernels.\n",
"4. (Optional) The `ParticleFile` object to write the output to."
]
},
{
"cell_type": "code",
"execution_count": 8,
"metadata": {},
"outputs": [
{
"name": "stderr",
"output_type": "stream",
"text": [
"INFO: Compiled ArrayAgeParticleAgeNorthVelAdvectionRK4 ==> /var/folders/s5/gxtkk3c12yqgd7hkt1b_s0vr0000gq/T/parcels-503/lib0e6e33cfac25b37419e8147f88e31705_0.so\n"
]
}
],
"source": [
"output_file = pset.ParticleFile(name=\"GCParticles.nc\", outputdt=3600) # the file name and the time step of the outputs\n",
"\n",
"pset.execute(kernels, # the kernel (which defines how particles move)\n",
" runtime=86400*24, # the total length of the run in seconds\n",
" dt=300, # the timestep of the kernel in seconds\n",
" output_file=output_file)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"While executing the `ParticleSet`, parcels stores the data in **npy** files in an output folder. To take all the data and store them in a netcdf file, you can use [**`ParticleFile.export()`**](https://oceanparcels.org/gh-pages/html/#parcels.particlefile.ParticleFile.export) if you want to keep the folder with npy files; or [**`ParticleFile.close()`**](https://oceanparcels.org/gh-pages/html/#parcels.particlefile.ParticleFile.close) if you only want to keep the netcdf file:"
]
},
{
"cell_type": "code",
"execution_count": 9,
"metadata": {},
"outputs": [],
"source": [
"output_file.export()\n",
"output_file.close()"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"After running your simulation, you probably want to analyse the output. Although there is some simple plotting functionality built into Parcels, we **recommend you write your own code** to analyse your specific output and you can probably separate the analysis from the simulation. \n",
"#### For more tutorials on the parcels output:\n",
"* [How the output is structured and how to start your own analysis](https://nbviewer.jupyter.org/github/OceanParcels/parcels/blob/master/parcels/examples/tutorial_output.ipynb)"
]
}
],
"metadata": {
"celltoolbar": "Metagegevens bewerken",
"kernelspec": {
"display_name": "Python 3 (ipykernel)",
"language": "python",
"name": "python3"
},
"language_info": {
"codemirror_mode": {
"name": "ipython",
"version": 3
},
"file_extension": ".py",
"mimetype": "text/x-python",
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython3",
"version": "3.8.13"
}
},
"nbformat": 4,
"nbformat_minor": 4
}