writing out model grid location

[JamiePringle]

Dear all --

I am trying to figure out how to write out the location of particles as grid-indices. I am using NEMO on the c-grid, and would prefer the output with respect to the rho/T/S grids. Is this data available to the JITParticle class? What are the names of the relevant variables that are accessible within JITParticle class if any?

I would prefer to have this written out by making it an explicit variable with the JITParticle class instead of post-processing, because going from location to grid-cell in a curved grid is surprisingly subtle and easy to screw up. 

In the end, the goal of this effort will be just to write out grid locations for use in calculating connectivity. This requires high frequency output in order to catch a particle moving briefly through a cell, and ultimately I will use this grid-cell information to reduce the amount of IO needed. But for now, I just want to save the grid-indices along with the other data.

Thank you
Jamie Pringle
University of New Hampshire

[erikvansebille]

Good question, @JamiePringle! I fully agree that this would be a useful feature to have. I developed a PR (#1147), which works for SciPy, but not JIT. Two comments

• We don't (yet?) have a way in JIT to sample particle.xi, let alone store it to a file. How to sample a particle in JIT mode (i.e. make the line below work in JIT too) be something that @CKehl has an idea for?https://github.com/OceanParcels/parcels/blob/140329373dc79933cc4e1925113f7e4377a4eaf4/tests/test_particle_file.py#L183 might
• Note that when using multiple grids, there is not one particle.xi, but a particle.xi for each grid. So users will anyways need to be careful for which grid the particle.xi are outputted.

[CKehl]

I totally agree with Erik here - apart from the reasons he mentions, there are straight-out 4 more reasons why I cannot support this request as a general feature.
I also suggest going with Erik's proposal - create a variable in your Particle definition yourself to store the index.
If you think about 'which index to use for accessing xi[...]' ? These grids are are added successively in the order they are defined when constructing the grid (i.e. fieldset) - hence, typically, U is index 0, V is index 1, W is index 2, and then the further fields you're using.

[erikvansebille]

hence, typically, U is index 0, V is index 1, W is index 2, and then the further fields you're using.

Well, this is not true for A-grids and C-grids, where U, V and W share one grid. This makes it even more tricky to assess which Field is associated with which grid...

[JamiePringle]

Thanks -- keeping track of which grid is not a problem. I will try this out. I am unclear of one thing in the answer -- does particle.xi only exist in the scipy particles, and not the JIT particles? In the JIT particle is there no visibility into the underlying indices? I am happy to create new variables to save in the JITParticle.

More generally, is there a list somewhere of what variables can be accessed as part of the scipy particles and JIT particles?

Jamie

[erikvansebille]

No, it's surprisingly tricky to access the xi in JIT, and it depends for example on whether you run in AoS or SoA ParticleSet mode.

Furthermore, this requires adapting the codecompiler.py, which is where Kernels are converted to C so they can be compiled. In order to make this work, we will need some code to parse particle.xi directly. But working with codecompiler.py is only for the brave ;-)

I won't have time to dig into this in the coming month, but perhaps at some point in May I could see if I can develop this feature. So am sorry; but not anytime soon..