Add custom gradient function

[dopplershift]

It seems np.gradient doesn't support arrays for dx, etc. For a case where we want to compute the geostrophic wind on a lat/lon grid (like GFS), we have a varying dx. Time for our own I guess.

[deeplycloudy]

A while ago I wrote a 4th order accurate gradient and third order accurate second derivative (along each axis) function. From a quick glance at it, it looks like it might adapt trivially to arrays of dx. I put it in a gist, which you're free to use.

https://gist.github.com/deeplycloudy/1b9fa46d5290314d9be02a5156b48741

[sgdecker]

I believe this issue explains why the geostrophic winds in the example at https://github.com/Unidata/notebook-gallery/blob/master/notebooks/Geostrophic_Wind_Example.ipynb are close to but not exactly parallel to the height contours.

[dopplershift]

Move up in priority since numpy 1.12 has broken the ability to pass arrays for dx/dy. We just need to move to our own, likely based on @deeplycloudy 's implementation. (cc @kgoebber )

[deeplycloudy]

It's mentioned in the comment at the top of the gist, but my implementation should drop in easily, since it was based on the version in numpy. It's worth noting that the implementation extends easily to vectorizing other derivative stencils of various orders of accuracy.

[dopplershift]

I agree, but I can't help but think of this:

[deeplycloudy]

I believe that's called Mathematica

[laufers]

Or someone’s master thesis…..

…

On Apr 24, at 5:43 PM, deeplycloudy ***@***.***> wrote: I believe that's called Mathematica — You are receiving this because you are subscribed to this thread. Reply to this email directly, view it on GitHub <#174 (comment)>, or mute the thread <https://github.com/notifications/unsubscribe-auth/AAzdAjKQYXVzqU5jZd8_wyWdG_AK4Iiaks5rzSWGgaJpZM4JCjj9>.

[lesserwhirls]

...a few of us resemble that remark...

[jrleeman]

@tjwixtrom just hit this in working up some new examples. Will be nice to fix by 0.6 if possible. Probably a decent amount of example/doc updating to go with it.

[dopplershift]

Not sure, but numpy/numpy#9408 (for 1.13.2) may restore the variable differences.

[eric-wieser]

@dopplershift: Passing anything other than 1D arrays was never supported, and computed incorrect results. 1.13.2 just raises an error in those cases instead

[dopplershift]

Good to know, thx.

[sgdecker]

The following reference (especially the first three equations) might have some useful info for computing the gradient (and other quantities that involve the del operator like divergence) properly on non-uniform grids:
http://www1.arl.noaa.gov/documents/reports/cmapf.pdf
Looks like the routines still exist at https://www.arl.noaa.gov/HYSPLIT_cmapf.php as well.

[dopplershift]

Thanks @sgdecker ! We'll have a look.

[dopplershift]

Here's an updated link for the NOAA stuff: https://www.arl.noaa.gov/hysplit/cmapf-mapping-routines/

I think we'll defer on the issues with projected coordinate systems in the short term, but I should get in a PR in a bit that tackles the non-uniform issue.

[dopplershift]

Re-opening to remind us to tackle the coordinate system issues.

[dopplershift]

We can continue to talk about the map projection issues in #893.