maroba / findiff

Python package for numerical derivatives and partial differential equations in any number of dimensions.
MIT License
409 stars 59 forks source link

multi dimensional coordinates #44

Open markcollier opened 3 years ago

markcollier commented 3 years ago

Hi, is it possible to generate derivatives of multi dimensional coordinates, for example curvilinear (2d) latitude longitude grids? I have tried the following example which does not work:

lat = np.array(( ((50,60), (55,65), (60,70)) ),dtype=float) lon = np.array(( ((0,180), (0,180), (0,180)) ),dtype=float)

print('lat.shape=',lat.shape) print('lon.shape=',lon.shape)

d_dlon = FinDiff(1, lon, acc=2)

print(d_dlon)

time = np.array((1,3,5,7), dtype=float)

data = np.array(( (((3,5), (6,2), (4,3)), ((3,5), (6,2), (4,3)), ((3,5), (6,2), (4,3)), ((3,5), (6,2), (4,3))) ),dtype=float)

print('data.shape=',data.shape)

ddata_dlon = d_dlon(data)

print('ddata_dlon=',ddata_dlon)

maroba commented 2 years ago

Do you mean something like the following?

import numpy as np

# latitude coordinates:
thetas = np.linspace(-np.pi/2, np.pi/2, 101) 

# longitude coordinates:
phis = np.linspace(0, 2*np.pi, 101) 

# set up 2D latitude/longitude grid:
Thetas, Phis = np.meshgrid(thetas, phis, indexing='ij')

# then define some 2D array:
f = np.sin(Phis)**2 * np.sin(Thetas)**2

print(
 "theta = %f, phi = %f, f(theta, phi) = %f" % (thetas[17], phis[21], f[17, 21])
)

# partial derivative with respect to Theta:
delta_theta = thetas[1] - thetas[0]
d_dtheta = FinDiff(0, delta_theta, 1)

# apply it:
f_theta = d_dtheta(f)