-
-
Save spencerkclark/bae91200c3d2c1b59adcd831b1594a2e to your computer and use it in GitHub Desktop.
Hi @kuchaale,
Sorry I just saw this now (for some reason I wasn't notified by your tagging of me in your post). Many thanks for trying things out!
For right now I would say the boundary behavior of xdiff
you experienced is expected, because I've only allowed for the mode='wrap'
option, which is enabled by default, and the interior values look good.
But it's clear you've uncovered a bug in xgradient
. I think I've fixed things in the updated gist:
In [1]: test = xr.DataArray([30., 33., 39., 36., 41., 37.], dims=['x'])
Out [1]: <xarray.DataArray (x: 6)>
array([-1. , 2.25, 0.75, 0.5 , 0.25, -2.75])
Dimensions without coordinates: x
In [2]: xgradient(test, 'x', spacing=2., accuracy=2)
Out [2]: <xarray.DataArray (x: 6)>
array([ 0.75, 2.25, 0.75, 0.5 , 0.25, -4.25])
Dimensions without coordinates: x
Note that my result using xgradient
is different than NCL's on the boundaries, because I'm doing second-order accurate differencing everywhere. NCL uses second-order centered differencing on the interior and first-order forward or backward differencing on the boundaries.
In the updated gist I also added some logic in xgradient
that checks to see if the underlying data is stored in a dask array; only then will it chunk the left and right components. (This way you don't need to call .compute()
at the end anymore).
Hi @spencerkclark,
I like your approach very much. However, when I tried to produce similar results to ncl example 1 using center_finite_diff_n (see below), I receive different results especially on the boundary points. See below:
t = np.arange(6); x = [30.,33.,39.,36.,41.,37.]; test = xr.DataArray(x, coords=[t], dims=['t'])
Using xdiff
xdiff(test, 't', accuracy=2, method='centered', spacing=2.)
produced:
`<xarray.DataArray (t: 6)>
array([-1. , 2.25, 0.75, 0.5 , 0.25, -2.75])
Coordinates:
Using xgradient
xgradient(test, 't', accuracy=2, spacing=2.).compute()
produced:
`<xarray.DataArray (t: 6)>
array([ 0.75, 5.25, 0.75, 0.5 , 4.5 , -4.25])
Coordinates:
Am I doing something wrong? Thx for your answer,
Ales