Addressing field staggering

[pletzer]

Problem

The CF conventions have dealt with Arakawa fields for a long time, which are staggered. For instance, temperature might be on cell nodes, velocity on cell faces, vorticity on cell edges and density attached to cells. The traditional way was to provide additional coordinates for these fields, corresponding to their staggering. However, this does not accurately represent the fields. For instance, a face velocity is not located at the middle of a face middle; rather there is a mean velocity attached to that face. Likewise, a density is not at the cell centre but instead there is a mean density across the cell.

This distinction could easily be dismissed as unimportant in the case of a uniform lat-lon grid where everyone agrees where the middle point of a face is, but becomes more complicated in the case of the cubed-sphere, which has very distorted faces when projected onto a lat-lon coordinate system.

This issue can be addressed by using "extensive" fields instead of "intensive" fields. Extensive fields represent integrals of density, velocity, etc. Thus, with an extensive field we're saying there is a total mass inside each cell, a total flux across each face, etc. Extensive fields have many advantages over intensive fields: (1) there is no need for a separate masking array, (2) they don't suffer from pole singularities to the same extent and (3) they are critical for conserving mass, energy, fluxes, etc.

The problem is that we don't have any way to specify whether a field is integrated (extensive) or defined on nodes (intensive). If I'm not mistaken, neither CF nor UGRID have any attribute to distinguish extensive from intensive fields.

Possible solution

While it would be best to extend the CF and UGRID conventions, the presence/absence of the "coordinates" attribute could be used to distinguish intensive/extensive fields. Extensive fields would not have the "coordinates" attribute. Intensive fields are required to have the coordinates attribute.

[zklaus]

It's worth noting that this seems to be done in the SGRID convention, which was indeed inspired by UGRID and that both conventions are strengthening ties with CF, albeit at the typical CF speed, see cf-convention/cf-conventions#153.

[zklaus]

Thanks, @pletzer for this discussion. It's certainly interesting and important to see how we can develop support for newer grids and perhaps features of existing grids that have received suboptimal support so far.

I am not sure how much of this discussion should be in this regridding interface project and how much should perhaps be moved upstream, for example to the Iris "mother" project for support in the infrastructure and data model, i.e. the cube, or even up to the CF conventions. But let's see if we can sort out a couple of points already here.

First, I think the presence or absence of the coordinates attribute can not be used for this purpose because it is already important for too many other things.

Second, cell_methods are by no means limited to time. It's perfectly idiomatic CF to use those for other dimensions. Perhaps that's all that is needed?

Third, I do recognize the importance of intensive vs extensive fields and find that this generally should be considered when deciding on a regridding method, though up to now that is considered expert knowledge that the user has to bring to the table.
I am, however, a bit fuzzy on the general advantages that you bring up. Could you elaborate on why one wouldn't need a mask? It still seems important to be able to distinguish, say, ocean points from land points. Points (2) and (3) are certainly true, but I find it hard to count them as general advantages. To me, they seem to be simply inherent properties. But any model will need both intensive and extensive quantities, the extensive ones for the reasons you mention, the intensive ones for other reasons, e.g. local accuracy.

[pletzer]

Thanks @zklaus. That's a good point, maybe cell_methods might be sufficient to indicate that a field is extensive!

Only intensive quantities require masking. Extensive quantities, because the field is integrated over edges, faces or cells, don't need masking. If a cell is outside the ocean domain (ie on land) then the water content in this cell is just zero. If the cell is partially valid then the total amount of stuff in that cell will be lower than that of a neighbouring cell -- there is no need to define a masking array when using extensive field as it is contained by the integral.

Masking arrays become a pain when some fields are on edges, while others are on faces and yet other on cells. To rely on intensive fields, one would then need to carry around arrays denoting how much each edge is valid, how much each face is valid and how much each cell is valid. ESMF wants us to pass the fractional cell arrays (srcfracfield and dstfracfield) when regridding. However, these quantities are just multiplied to the field values times the cell volumes, essentially they are converting the intensive fields to extensive fields.

[zklaus]

I'm still a bit unconvinced about the masking. Sure, you can just put zero and be happy, but "zero" and "masked" is not necessarily the same thing. There is, for example, a difference between an ocean point that happens to have no sea-ice and a land point.

But regardless, perhaps we can turn the discussion now to practical implications for this project. Suppose we make it so that intensive and extensive fields can be clearly identified. What should iris-esmf-regrid do about it?

[pletzer]

Thanks @zklaus.

(1) Masking is required for nodal fields but only for nodal fields. Masking can be dispensed with with extensive edge, face and cell fields. (As an aside, intensive and extensive are the same for a nodal field since there is nothing to integrate over.)

(2) Assuming there is some attribute that tells us if a field is intensive or extensive and whether it is nodal, edge, face or cell. Let's say we have a Regridder instance (rgd), it knows about the src and tgt meshes/grids. Here are some examples:
(a) tgt_intensive_nodal = rgd.regrid(src_intensive_nodal, src_mask=None, tgt_mask=None) # bilinear, src/tgt masks are optional
(b) tgt_intensive_edge = rgd.regrid(src_intensive_edge, src_invalid_frac=None, tgt_invalid_frac=None) # edge regridding with fraction of invalid edges
... similar for face data...
(c) tgt_intensive_cell = rgd.regrid(src_intensive_cell, src_invalid_frac=None, tgt_invalid_frac=None) # area weighted regridding with fraction of invalid cells

(One could also set the masking in the GridInfo or MeshInfo?)

Now for the extensive fields
(d) tgt_extensive_edge = rgd.regrid(src_extensive_edge) # mimetic from MINT
(e) tgt_extensive_face = rgd.regrid(src_extensive_face) # mimetic from MINT
(f) tgt_extensive_cell = rgd.regrid(src_extensive_cell) # area weighted from ESMF

We might also need some method to go back and forth between intensive and extensive fields. Since the extensive fields regridding interface is much cleaner we could ask the user to convert their field to extensive prior to regrid.

[zklaus]

Thanks, @pletzer.

(1) I feel like sea-ice as mentioned above can be an example of an extensive field that would also like to have a mask. But I also think this issue is a bit tangential. At the end of the day, it seems very simple to me: If the user gives us a field with a mask (keep in mind that the mask is part of the cube), the regridding has to take the mask into account and produce a masked field as a result, if they don't, we have no mask to work with and need to rely on the field as a whole being valid.

(2) So if I understand correctly, you are saying that the regridder should choose some defaults based on characteristics of the field, including whether it's extensive or intensive, whether it lives on edges, faces, cells (or vertices?), but this might conceivably be extended to other characteristics such as the relative resolution of source and target.

That does make sense to me, but I think it should always remain possible for the user to override these defaults and specify their own choices, and of course, it presupposes that the option is available in the first place.

If this understanding is correct, I suppose the question is can be split up into (a) designing the API, (b) implementing the necessary options, (c) put in place a nice framework to devise defaults from the fields, (d) choose and document the default choices.