Python では、2 次元の中央有限差分スキーム (辞書編集) 係数を開発しようとしています。私のコードは
import numpy as np
import scipy as sp
x = np.linspace(0, 9, num=10)
nx = len(x)
hx = (x[-1] - x[0])/(nx - 1)
one_x = np.ones((nx, 1))
x_coeffs = np.hstack((one_x, one_x))
x_coeffs = np.insert(x_coeffs, 1, -2, axis=1)
# So far so good, then carry on
x_diff = spdiags(x_coeffs.transpose(), np.asarray([0, 1, 2]), nx-2, nx)
nx_id = identity(nx)
Dx = np.kron(nx_id, x_diff.transpose())
Dx = 1/hx*Dx
ダミーの例で試してみましょう。
x = np.linspace(0, 4, num=5)
f = np.power(x, 2)
dfdx = Dx * f
しかし、このダミーの例の結果は です(array([ -4., -6., -8., -10.])
。何がうまくいかないかについて何か考えはありますか?