Tripolar B-grid issue in the North Pole (MOM5 velocity fields)

[eavellashaw]

Hi !

I am currently trying to advect particles (3D advection) from the North-Atlantic to the Arctic with MOM5 velocity fields (outputs from GFDL CM2.5 and CM2.6 models runs). My grid is curvilinear and tripolar, here's a quick snapshot of my geolons and geolats on a T-grid:

My velocity fields have the following dimensions:

u (st_ocean, yu_ocean, xu_ocean)
v (st_ocean, yu_ocean, xu_ocean)
w (sw_ocean, yt_ocean, xt_ocean)

In the beginning I thought about using xu_ocean and yu_ocean, but the tripolarity of my grid is then not accounted for (I do not include it here but it leads to erroneous results). So I then switched to using geolon_c and geolat_c (basically the same as in the figures above but on the U grid).

This again lead to weird results, so I was wondering how Parcels was handling this kind of tripolarity, and if I had to change something in its code due to the fact that the velocities are defined on the 1-D vectors coordinates (y and x) and not on the 2-D ones (geolat and geolon)?

This is a kind of follow-up question to #860, but the question was not raised as the tripolar grid does not influence the coordinates in the Southern Ocean. The script I am using is available here.

Thanks in advance for your help!

[erikvansebille]

Thanks for raising this @eavellashaw. I'm surprised that this doesn't work; so far it has always worked on C-grids. See e.g. the homepage animation tutorial, where particles move happily through the Arctic in curvilinear NEMO simulation (see also the Curvilinear Grid tutorial.

Can you show what is wrong with your particles?

[eavellashaw]

Thanks so much for your fast reply @erikvansebille! A fellow student used Parcels in a NEMO configuration and it worked, that is why I am confused too ...

Here are two figures I made showing the advection of the particles, colors corresponding to the depth the particles were at at time = 0. The first one is the result when run with xu_ocean and yu_ocean, and the second one when run with geolon_c and geolat_c.

The particles are generally supposed to circulate along the slope (black arrow), not go straight toward the North Pole (1st figure). In the second one, it is also acting weirdly as the particles seem to be "blocked", but they are not that deep and the bathymetry should not be a problem at their depths.

I was thinking that when interpolating, Parcels may get the wrong indices due to the fact that the dimensions of the velocity fields are x and y, not geolon and geolat. I am fairly new to all of this so I might be completely wrong...

If there is any kind of unclarity don't hesitate to let me know!

[erikvansebille]

Dear @eavellashaw, I must say that the stop at exactly 80N in the right figure looks conspicuous. But I have no idea why this would go wrong. As you say, there is no problem with NEMO's C-Grid, so it might have to do something with MOM5's B-grid?

I'm sorry, but I can't help you much more. I don't have access to B-grid models (we don't use them in my team) and thus also no time to dig much deeper into this. Perhaps other MOM5/Parcels users might have solutions?

[eavellashaw]

I think it definitely has something to do with MOM5's B-grid, I am now reaching out to people having worked with MOM5 B-grid velocities (typically the ones mentioned in #860), I'll let you know if a solution pops out!

Thanks for your help!

[eavellashaw]

Hi @erikvansebille, as I have a bit more insight on the situation, I have a follow-up question.

So, I ran two same simulations but creating the fieldset with from_mom5() in the first one and with from_nemo() in the second one. I chose to do so to be sure that the error lied elsewhere than in the curvilinearity and tripolarity of my grid.

I don't know if you are familiar with the surface circulation in the Arctic, but the results in the second figure are way more convincing than in the first one.

My question is the following: Do you know what is the main difference when indexing the latitude and longitude coordinates of the particles between from_nemo() and from_mom5() when performing the interpolation? I think the error might come from parcels looking at the wrong coordinates.

As said earlier, here are the dimensions of my velocity fields:

u (st_ocean, yu_ocean, xu_ocean)
v (st_ocean, yu_ocean, xu_ocean)
w (sw_ocean, yt_ocean, xt_ocean)

These are 1-D arrays and when indexing, parcels might be searching for indices in these arrays and not in geolat_c and geolon_c.

This is of course an hypothesis, tell me if it makes sense and if you know where I should look into the parcels code for testing it! I already looked into the code but I am not sure where to look specifically and what to change.

Cheers,

[erikvansebille]

Thanks for diving into this, @eavellashaw. I agree with you that the second plot (with from_nemo()) looks more realistic.

The biggest difference between from_mom5 and from_nemo is that the first expects an Arakawa B-grid, whereas from_nemo interprets the data as on a C-grid.

So if your data is indeed on a C-grid, it makes much more sense to use from_nemo (or from_mitgcm, which has slightly different conventions on where the grid corners are)

Anyways, the specific calls for the MOM5 interpolation are at every instance where MOM5 is used in the code here https://github.com/OceanParcels/parcels/blob/dad04e30fbc4b64f687b861fa55b791d64ba65df/parcels/include/parcels.h#L437-L551 and also https://github.com/OceanParcels/parcels/blob/dad04e30fbc4b64f687b861fa55b791d64ba65df/parcels/include/index_search.h#L101-L105

But note that the code above is the C/JIT version, which may be complicated to read. If you prefer to look at the scripy version, that is here https://github.com/OceanParcels/parcels/blob/dad04e30fbc4b64f687b861fa55b791d64ba65df/parcels/field.py#L1076-L1080 and https://github.com/OceanParcels/parcels/blob/dad04e30fbc4b64f687b861fa55b791d64ba65df/parcels/field.py#L698-L703

[eavellashaw]

Thank you for this fast reply, I'll dive into the code and let you know if I make an interesting discovery!

I am working on an Arakawa B-grid with mom5 velocity fields (CM2-O models suite), that is why I used from_mom5(). But I thought about using from_nemo() as the biggest difference between these two are the interpolation methods and the locations of the points in the vertical. This was to check if something was wrong with my datasets, or if the problem might lie in Parcels, handling in the wrong way the coordinates of my velocity fields.

Anyway, I'll let you know!

[eavellashaw]

This issue has been solved in #1509 !