Re-ordering of Kernel actions to vectorially sum particle displacements

[erikvansebille]

(note: this PR supersedes #1388, which included attempts to write using parquet and sqlite. But the performance of these approaches was poor, so this PR will keep zarr writing.)

Until now, the way that Kernel-concatenation worked in Parcels meant that the execution of multiple Kernels was generally not commutable: for example first advection a particle and then letting it sink would lead to a different result than first sinking and then advecting. That is because particle positions were updated within each kernel, instead of the vectorial summing of displacements that would be more appropriate (and would ensure commutability).

More problematically, while the particle positions would be updated within a kernel, the time was generally only updated at the end. This meant that Field-sampling could give unexpected results (e.g. when sampling was done at the new location but the old time).

To solve these issues, in this PR I propose that individual Kernels compute displacements (particle_dlon, particle_dlat and particle_ddepth) and that these displacements are only added to the particle locations when all kernels have been executed for a timestep dt.

This means that the workflow for each particle p and each time t becomes (here for particle_dlon, same for particle_dlat and particle_ddepth):

1. Set particle_dlon = 0
2. Set particle.lon = particle.lon_nextloop
3. Execute all Kernels, including sampling Kernels, and let the Kernels update particle_dlon, when required
4. Update time to t += dt and particle locations to particle.lon_nextloop = particle.lon + particle_dlon
5. If particle p is written this timestep, write particle.lon

(note, the concept for this loop was updated after #1422)

Note that this workflow means that particles will be written at t=0 (and there is no need anymore for a separate initial values sampling step!), but not anymore at t=runtime (the last tilmestep), unless we also change that the particle loop is executed one time more at the end (tbd).

Also, since particle locations in all Kernels will now be vectorially summed, we don't need the SummedFields class anymore, simplifying the codebase a bit.

Issues/questions that need to be explored or implemented:

• a. Move the writing step to after the Kernel execution but before the location and time update
• b. Investigate whether it's wise to perform the loop standard one extra time to also write output at t=runtime (see comment above)
• c. Check if the particles_backups are stil needed now that particles are not moved in a Kernel
• d. Confirm that RK45 still works
• e. Check how to deal with if-statement in Kernels based on 'current' locations (see [OUTDATED; DO NOT MERGE] Vectorial summing of particle displacement at end of kernel loop #1388 (comment))
• f. Explore how to deal with Error Kernels
• g. Remove the SummedFields class
• h. Fix the writing in AnalyticalAdvection
• i. Fix interactive particles
• j. Update this new workflow in documentation and tutorials
• k. Check performance in production simulations

Since this is a major overhaul that can change results of a pet.execute() (so is not per se backward-compatible), I suggest to make this PR (when merged) part of v3.0.0