generating the dual mesh

[Chilipp]

hey @ChrisBarker-NOAA and @jay-hennen (and anyone else who's interested)!

I am currently discussing with one of my colleagues how to best visualize variables that are defined on the nodes and edges of a mesh (see psyplot/psy-maps#32). I already have a visualization for the faces (see https://psyplot.github.io/examples/maps/example_ugrid.html for instance), but nodes and edges are not yet implemented.

We think the best way would be to visualize these variables on the dual mesh.

For a face_node_connectivity like this

this looks like the following magenta lines for variables on the nodes

and the following magenta lines for variables on the the edges.

Currently we know of no library to compute this mesh. The ESMpy library defines a convert_to_dual parameter, but this is only implemented for SCRIP.

I'd be happy to implement this and could think of a method for the UGrid class with something like

    def build_dual_mesh(self, location: Union[Literal["node"], Literal["edge"]]) -> UGrid:

Would this be an interesting enhancement for gridded?

I do not have a clear plan of the algorithm yet, but I'd probably implement it in cython, because it would require quite a lot of lookups and looping. I'd be happy to come up with a PR if this is of interest for you as the developers of gridded.

ping @platipodium

[jay-hennen]

What sort of visualization are you looking for? Similar to the heatmaps in the psplot ugrid example you linked?

Also are we working exclusively with triangle grids here? In all cases this makes grid handling and visualization easier, which is why I bring it up.

In my opinion the best way to visualize polygon data on the vertices is simple Gouraud shading (linear interpolation) as seen in this basic example from matplotlib. This assumes of course that smoothing out the value inconsistency between nodes is what you want. The converse (maintaining the discontinuities) is very challenging indeed, and would probably warrant the enhancement you're suggesting here.

I never quite understood how edge data would get visualized in the triangle case. It sort of makes sense to me in the staggered quad grid case, when you can basically interpolate across the two opposing sides of a quad to get something coherent. What are you trying to visualize here?

[Chilipp]

thanks for the quick response @jay-hennen!

Also are we working exclusively with triangle grids here? In all cases this makes grid handling and visualization easier, which is why I bring it up.

nope, it's important that the algorithm works for 2D flexible mesh topologies, too. This is also what the proposed function build_dual_mesh function would return. For the dual mesh of the edges, for instance, you'd have triangles at the domain boundary and quadrilaterals in between (see the image above).

In my opinion the best way to visualize polygon data on the vertices is simple Gouraud shading (linear interpolation) as seen in this basic example from matplotlib.

well, I tend to disagree here. For model development, it is often useful to see the exact value in a specific grid cell, without any interpolation.

What are you trying to visualize here?

I must admit that I am also quite new to this thematic (which is why it is not yet implemented in my visualization framework). Nonetheless, a classic example I guess would be a flow from one grid cell into another. left to right would be positive, right to left would be negative. There are of course many different ways to visualize this (e.g. via a quiver plot), but the most generic one without any manipulation on the data or the need for adjustments (e.g. arrowsize) would be the visualization on a mesh as in the last figure above.

[jay-hennen]

Nonetheless, a classic example I guess would be a flow from one grid cell into another. left to right would be positive, right to left would be negative.

Oh I see the quads in the edges case now, but I don't know what visualization you would use with those quads that makes any sense for visualizing flow. What's left-to-right and right-to-left look like in an unstructured grid? I have a feeling it's pretty model-dependent how this gets represented too.

For model development, it is often useful to see the exact value in a specific grid cell, without any interpolation.

Hold on. This makes sense when talking about face data, but how does it make sense when talking about vertex data? If I wanted a visualization of node data to help with actual model development (not a presentation, grant proposal, etc) (and the 'in-betweens' values didn't matter), I would probably just draw the grid, and then put dots with labels on the vertices and edges with their values, and maybe toss in some color/size scaling. It's ugly, but much more informative.

Sorry man. While I can see definite computational use cases for this sort of 'dual mesh' representation, I'm still struggling to see how it helps visualize anything. Perhaps our use case is too simplistic and it's blocking my thinking. We often use gridded to represent an actual physical quantity (say, water Temperature) at specific times at specific locations in X/Y/Z space. In these sort of cases the visualizations are usually fairly obvious.

[Chilipp]

While I can see definite computational use cases for this sort of 'dual mesh' representation, I'm still struggling to see how it helps visualize anything.

ha, that's fine 😅 the purpose of this issue is not to talk about visualization. That's just my personal motivation for the moment and that would also happen outside of the gridded framework. My question is just, whether you would be interested if I implement the dual mesh generation as a PR to the gridded project, or whether I should rather make my own thing (which I'd prefer not to do).

[jay-hennen]

My question is just, whether you would be interested if I implement the dual mesh generation as a PR to the gridded project,

Oh well, I don't see why not. Go for it.

[ChrisBarker-NOAA]

I am also somewhat confused about what the point of this "dual mesh" is. Is the idea that you want to fill in a complete cell with one color corresponding to the value at a given node?

Nonetheless: the dual of a triangular mesh is the Voronoi Polygons (well, if the triangles are Delaunay -- and if not, you would still build it the same way).

That is -- a node of a triangular mesh can be shared by any number of triangles -- so it's "dual" cell is an polygon, not necessarily a triangle -- and in fact very rarely a triangle -- in your example, you got quads, but they will often be more than four sided.

Voronoi Polygons are build into at least some of the triangularization codes -- but if you want to implement it in gridded -- sure! It's pretty straightforward, if a bit computationally expensive.

I'm also confused about what the purpose of a "dual" to the edges is, as data "on the edge" usually represent fluxes. But I think is is well defined -- so go for that too :-)

[Chilipp]

thanks for the feedback @jay-hennen and @ChrisBarker-NOAA! I'll come back to this within the next week probably.

It's pretty straightforward, if a bit computationally expensive.

true, the basic theory is quite straight-forward. I think there are three challenges:

1. the computational efficiency
2. the domain boundary (for the nodes, for instance, you need to take midpoint of the edge if the edge does belong to one triangle only)
3. the ordering to make sure the nodes of the dual are in a counterclockwise order.

but let's see when I am getting a bit more concrete next week