Ariane, a Lagrangian analysis application for ocean circulation models
developed by Bruno Blanke and Nicolas Grima
Calculation of trajectories
In this fifth part of the course,
we come to the actual calculation of trajectories.
A classical approach for the calculation of trajectories
from a known initial position and a known velocity field
is to interpolate finely in space and time
the velocity field on the current position of the particle,
then to perform an infinitesimal displacement
with this interpolated velocity,
and finally to repeat the interpolation on the new position,
with a velocity field updated in time.
We can work with a temporal scheme called Runge-Kutta,
in order to solve the integration dx/dt = U,
knowing that the spatial interpolation
can be performed
by including a neighbourhood more or less close
to the point on which the velocity must be estimated.
An important parameter of this approach
is the time step ∆t
on which the Runge-Kutta scheme is constructed.
The smaller the ∆t chosen,
the more accurate the results will be,
with of course the penalty of a higher computational cost.
The interpolation in space of the velocity field
also requires special attention,
especially when working near a coastline or topography
on which the velocity values are no longer defined.
The approach used for Ariane is very different.
The spatial and temporal interpolations of the velocity field
are coarse, but robust.
The advantage of this approach is a fast,
reliable, time-reversible scheme,
which preserves the non-divergence of the flow,
even at scales finer than the model mesh.
The velocity field is not interpolated in time.
The archived instants of the model are used sequentially,
over the duration of their validity.
When needed, the velocity field is updated to its new value
and is kept constant until the next necessary update.
This simplicity leads to disadvantages
linked to a stepwise integration,
such as temporal discontinuities of the velocity locally,
with infinite accelerations, thus not physical.
However, nothing prevents the user
from providing a velocity field with a better time resolution
by sampling the calculations of the ocean model more finely
to limit the height of the steps.
The spatial interpolation is limited
to a linear and independent formulation
for each of the 3 transport components in each individual grid cell.
In this interpolation,
the velocity, or transport,
does not depend on time
since we are working on piecewise constant velocities.
The chosen interpolation allows us to linearly link
the known velocity components on 2 opposite faces of a grid cell.
Again, this coarse interpolation
introduces discontinuities, this time spatial, of the velocity field.
For example, at the transition between 2 grid cells,
only the transverse component of the velocity field is continuous,
while the 2 other components are calculated independently
(in a discontinuous way)
upstream and downstream of the grid cell change.
The great advantage of the method
is that the interpolated velocity field at the sub-grid cell scale
is non-divergent,
like the velocity field defined at the grid cell scale in the ocean model.
This property is very interesting
when using trajectory calculations for water mass tracking.
Indeed, in Lagrangian mode, we respect
one of the fundamental properties present
in the equations solved by the ocean model.
The other main advantage
is that the trajectory equation is solved analytically.
From an input position in a grid cell of the model,
simple relations determine
the output position on one of the other faces of the grid cell.
Some illustrations allow this calculation to be visualised.
It is presented here in 2 dimensions,
ignoring the vertical component of motion.
The transport field is known on the west, east, south and north faces of a model grid cell.
Linear interpolation of each component of the transport
allows the velocity to be determined
analytically and virtually on the entire interior of the grid cell.
In 2 dimensions, and for a velocity field kept constant piecewise in time,
the trajectories are the streamlines of the flow.
Their formulation is simple to write.
In practice, and in 3D, we have in each grid cell 3 independent equations
to continuously define the velocity components U, V and W inside that grid cell.
From a given initial position on the edge or in the grid cell,
it is easy to determine the virtual travel time
to leave the grid cell in each of the 3 directions.
The minimum of the 3 travel times
(in i, j and k; or x, y and z, if we refer to geographical coordinates)
gives the actual time needed to cross the grid cell,
and is used to precisely diagnose the exit position.
Let's go back to the previous example,
with the choice of an initial position
in the middle of the western edge of the grid cell.
For the chosen initial position, we can see
that the virtual travel time to the eastern edge of the grid cell
(log(3)/2)
is greater than the travel time to the northern edge
(log(2)/2).
It is thus by the northern edge that the particle will escape,
and the analytical relations give its precise position
on this edge between two grid cells.
Ariane is therefore characterised by a significant bias
in relation to spatial and temporal interpolation of the velocity field
on the current position of a particle.
The advantages of the formalism chosen for Ariane are:
1) Calculations performed as close as possible to the model grid cell
(which produced the velocity field used);
- A possible adaptation for online calculations,
i.e., performed in the core of the ocean model
as it computes the 3D velocity field.
- A simple and inexpensive analytical resolution of individual particle trajectories.
- The absence of any adjustable parameter
(the trajectories are unique, intimately linked to the definition of the velocity field).
- The possibility of backward calculations,
virtually superimposable on direct calculations.
To do this, we only need to order in reverse chronological order
the archived velocity fields and to multiply them by -1.
2) A continuity relation (non-divergence relation)
preserved at the local scale of the particles,
which is a great advantage from the point of view of tracking the properties of water masses.
3) A formalism directly linked to the ocean model grid,
and calculations performed directly on transports (rather than on velocities),
which introduces normalisation factors into the core of the calculations made by Ariane.
It is essential to keep in mind the drawbacks of using this scheme.
Ariane works on the basis of stepwise knowledge
of the archived time series of a velocity field.
It is up to the user to ensure that this archive
includes all the relevant scales of variability that are to be analysed.
And Ariane works at the sub-grid cell scale on a velocity field interpolated in a rather crude way,
leading to spatial discontinuities of certain velocity components,
and thus to non-physical accelerations.
To present this trade-off between simplicity and robustness in a way favourable to Ariane,
here is a comparison taken from a DESS internship.
The comparison concerns the spatial interpolation of the velocity field,
leaving aside its temporal variations.
The objective is to determine the trajectories corresponding to Stommel's solution
for the circulation that is established in a rectangular basin
under the effect of a surface wind characteristic of mid-latitudes.
The exact solution is expressed as a stream function,
whose contours are the virtual trajectories of particles in the fluid.
By construction, these contours are closed curves,
with a fast and narrow northward motion in the western part of the basin,
and a slower and more diffuse southward motion
for the return current in the eastern part of the basin.
With a 16-point Lagrangian interpolation of the velocity field
and a Runge-Kutta integration of order 4,
the trajectories are close to the theoretical solution,
but they hardly loop back on themselves.
Spatial interpolation of the velocity field in the western part of the basin,
where velocities and their gradients are strongest,
and near the coast,
leads to a progressive loss of information and to spiral trajectories.
This defect can be corrected, without being completely eliminated,
by refining the time step used in the Runge-Kutta algorithm,
at the price of more expensive calculations.
The solution computed by Ariane integrates the conservation of mass
i.e., the respect of the non-divergence of the flow,
down to the motion of a point particle.
Consequently, the trajectories formally follow closed contours,
showing very good agreement with the theoretical solution,
for a moderate computational cost.