0

在 python 中,我正在尝试开发二维中心有限差分方案(字典)系数。我的代码是

    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.]). 任何想法可能会出错?

4

0 回答 0