cosima-recipes icon indicating copy to clipboard operation
cosima-recipes copied to clipboard

xgcm vorticity formula

Open aekiss opened this issue 3 years ago • 6 comments

@navidcy should https://github.com/COSIMA/cosima-recipes/blob/master/DocumentedExamples/RelativeVorticity.ipynb use dxt instead of dxu in cell 18?

ζ_xgcm = ( grid.interp( grid.diff(ds.v, 'X') / grid.interp(ds.dxu, 'X'), 'Y', boundary='extend')
          - grid.interp( grid.diff(ds.u, 'Y', boundary='extend') / grid.interp(ds.dyt, 'X'), 'X') )

dxu is inconsistent with the subsequent text.

aekiss avatar Oct 28 '21 04:10 aekiss

What do you mean "inconsistent with the subsequent text"?

The objective in this cell is to reproduce ζ_mom5. If you can do that with dyu then go for it and make a PR.

I'm not suggesting that what I did is definitely the only way forward. I don't even recall my thinking. But at least by plotting the difference you can see that (the values at least) of ζ_xgcm match those of ζ_mom5. The plot shows almost-zero (see the figure below). However, ζ_xgcm lives on t-grid while ζ_mom5 lives on u-grid.

navidcy avatar Oct 28 '21 06:10 navidcy

I haven't looked at this in great detail. It just seems odd to have dxu and dyt in the same formula, especially as the "simpler formula" given afterwards uses dxtand dyt:

ζ_xgcm = grid.interp(grid.diff(ds.v, 'X'), 'Y', boundary='extend')/ds.dxt - grid.interp(grid.diff(ds.u, 'Y', boundary='extend'), 'X')/ds.dyt

aekiss avatar Oct 28 '21 07:10 aekiss

It's odd. If you can help sort it out in a cleaner way it'll be great.

But the objective of the first (seemingly weird) method was to replicate was MOM5 is actually doing. I couldn't figure out a way to do that using dyt! And also I couldn't figure how to do this computation and end up on the u-grid (like where ζ_mom5 lives). That is, the way the notebook is now ζ_mom5 and ζ_xgcm have the same values but live on different grids! If you can figure out how to replicate the values of ζ_mom5 and end up on the u-grid let me know. However, what (I think I understood that) you suggested doesn't even give same values!

ζ_xgcm.values[200:205, 100:104]

array([[ 1.1507019e-06,  1.0019297e-06,  5.4093033e-07,  2.8089801e-07],
       [ 9.9383203e-07,  9.8258420e-07,  1.0025971e-06,  9.3004383e-07],
       [ 1.9963846e-07,  5.3464299e-07,  7.5498042e-07,  4.4450223e-07],
       [-1.6323114e-07, -1.6445910e-08,  2.2685566e-07, -1.9235793e-07],
       [-6.1921597e-08,  2.8606735e-08,  3.5336069e-07,  2.6428950e-08]],
      dtype=float32)

ζ_aekiss = (grid.interp( grid.diff(ds.v, 'X') / grid.interp(ds.dxu, 'X'), 'X') 
             -  grid.interp( grid.diff(ds.u, 'Y', boundary='extend') / grid.interp(ds.dyu, 'Y', boundary='extend'), 'Y', boundary='extend'))
             - 
ζ_aekiss.values[200:205, 100:104]

array([[ 1.0400078e-06,  8.1025718e-07,  6.4869403e-07,  9.7753821e-07],
       [ 7.0336637e-07,  8.3177400e-07,  8.4166714e-07,  7.4260828e-07],
       [ 1.4705354e-07,  4.5017367e-07,  3.7383137e-07, -1.6373960e-07],
       [-4.7659199e-08,  1.5902492e-07,  1.2845433e-07, -4.5021255e-07],
       [-3.4233143e-08,  2.7488650e-07,  3.5775838e-07, -1.7128500e-08]],
      dtype=float32)

ζ_mom5.values[200:205, 100:104]

array([[ 1.1507033e-06,  1.0019312e-06,  5.4093209e-07,  2.8090085e-07],
       [ 9.9383305e-07,  9.8258579e-07,  1.0025997e-06,  9.3004758e-07],
       [ 1.9963758e-07,  5.3464328e-07,  7.5498144e-07,  4.4450331e-07],
       [-1.6323250e-07, -1.6447160e-08,  2.2685441e-07, -1.9235949e-07],
       [-6.1921909e-08,  2.8606166e-08,  3.5335984e-07,  2.6427927e-08]],
      dtype=float32)

But if replicating the MOM5 diagnostic is not part of the requirements then you should just use the "simpler"

ζ_simpler = grid.interp(grid.diff(ds.v, 'X'), 'Y', boundary='extend')/ds.dxt - grid.interp(grid.diff(ds.u, 'Y', boundary='extend'), 'X')/ds.dyt

which lives on the t-grid, or the even simpler

ζ_evensimpler = ds.v / ds.dxu - ds.u / ds.dyu

which lives on the u-grid.

Bottom line is that it all depends what you wanna compute. Maybe I have missed your point entirely? Sorry if I have done that... But changing dxt->dxu just for the sake of symmetry ignoring

navidcy avatar Oct 29 '21 00:10 navidcy

I think you have a point here. Or at least there is something fishy.... because the ζ_xgcm lives on the cell centers while ζ_mom5 lives at the corners.

navidcy avatar Jan 23 '23 06:01 navidcy

Question: the vorticity_z diagnostic from MOM5 in which location is it output??

navidcy avatar Jan 23 '23 06:01 navidcy

Is this still relevant @aekiss ?

adele-morrison avatar Apr 16 '24 05:04 adele-morrison