Как получить специальную производную от интерполированной функции - PullRequest
Новая
       58

Как получить специальную производную от интерполированной функции

0 голосов
/ 31 января 2019

Я создал файл h5 для простого куба, затем прочитал его с помощью python и, наконец, использую функцию RegularGridInterpolator для интерполяции.У меня все отлично работает.Но я хочу знать, как я могу изменить свой код, чтобы я мог получить производную от этой интерполированной функции?Для вашей доброй информации, я привел здесь мой код:

код для создания файла h5 
import numpy as np
import h5py

def f(x,y,z):
   return 2 * x**3 + 3 * y**2 - z

x = np.linspace(-1, 1, 2)
y = np.linspace(-1, 1, 2)
z = np.linspace(-1, 1, 2)
mesh_data = f(*np.meshgrid(x, y, z, indexing='ij', sparse=True))

h5file = h5py.File('cube.h5', 'w')
h5file.create_dataset('/x', data=x)
h5file.create_dataset('/y', data=y)
h5file.create_dataset('/z', data=z)
h5file.create_dataset('/mesh_data', data=mesh_data)

h5file.close()

код для чтения файла h5 и интерполяции 

import numpy as np   
import h5py
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D
from scipy.interpolate import RegularGridInterpolator
from mpl_toolkits.mplot3d.art3d import Poly3DCollection, Line3DCollection
f = h5py.File('cube.h5', 'r')  
list(f.keys())
dset = f[u'mesh_data']
dset.shape
dset.value.shape
dset[0:2,0:2,0:2]
x = np.linspace(-1, 1, 2)
y = np.linspace(-1, 1, 2)
z = np.linspace(-1, 1, 2)
my_interpolating_function = RegularGridInterpolator((x, y, z), dset.value, method='nearest')
pts = np.array([0.2, 0.9, 0.6]) 
my_interpolating_function(pts)

значение интерполяции равно 4,0.

1 Ответ

0 голосов
/ 01 февраля 2019

Я не уверен, что понимаю, что вы ищете.Вот 1D пример, иллюстрирующий важные моменты, которые следует учитывать при численной оценке производной функции: 

import numpy as np
import matplotlib.pyplot as plt
from scipy.interpolate import interp1d

def f(x, y, z):
    return 2 * x**3 + 3 * y**2 - z

x_fine = np.linspace(-1, 1, 50)  # used for the plots

# Coarse sampling, only two points:
x_coarse = np.linspace(-1, 1, 2)

# Interpolation
interpolator_coarse = interp1d(x_coarse, f(x_coarse, 0, 0), kind='linear')

plt.plot(x_fine, f(x_fine, 0, 0), label='analytical')
plt.plot(x_coarse, f(x_coarse, 0, 0), 'ok', label='coarse sampling')

plt.plot(x_fine, interpolator_coarse(x_fine), '--r', label='interpolation based on the sampling')

plt.xlabel('x'); plt.ylabel('f(x, 0, 0)');
plt.legend();

график:

graph

Значение «истинной» производной при x = 0 равно нулю (плоский наклон).Однако, если производная вычисляется на основе данных выборки (при x = -1 и x = 1), независимо от любого вида интерполяции, оценочное значение будет равно 2 ...

Количество выборокточки должны быть увеличены: 

import numpy as np
import matplotlib.pyplot as plt
from scipy.interpolate import interp1d

def f(x, y, z):
    return 2 * x**3 + 3 * y**2 - z

x_fine = np.linspace(-1, 1, 50)  # used for the plots

# Coarse sampling:
x_coarse = np.linspace(-1, 1, 4)

# Interpolation
interpolator_coarse = interp1d(x_coarse, f(x_coarse, 0, 0), kind='linear')
interpolator_cubic = interp1d(x_coarse, f(x_coarse, 0, 0), kind='cubic')

plt.plot(x_fine, f(x_fine, 0, 0), 'k', label='analytical')
plt.plot(x_coarse, f(x_coarse, 0, 0), 'ok', label='coarse sampling')

plt.plot(x_fine, interpolator_coarse(x_fine), '--r', label='linear interpolation')
plt.plot(x_fine, interpolator_cubic(x_fine), '--b', label='cubic interpolation')

plt.xlabel('x'); plt.ylabel('f(x, 0, 0)');
plt.legend();

with 4 points

Наклон при x = 0 теперь ближе к нулю.Следующая часть проблемы заключается в оценке производной от выбранных данных, см., Например, Числовая_дифференциация .

...