正弦函数分析导数和近似导数的比较

Question

我试着比较正弦函数的解析导数和近似导数，用后向、前向、中心、3阶后差法，它们之间有很大的差别。这一定是错误的。你能给我一个提示，我是在哪里犯的错误吗？先谢谢您了!

%matplotlib notebook
import numpy as np
import matplotlib.pyplot as plt

from math import pi
from numpy import *



# spatial domain
xmin = 0
xmax = 1
n = 50 # num of grid points

# x grid of n points
x = np.linspace(xmin, xmax, n+1);

k=2


def f1(x):
    return np.sin(2*pi*k*x)
def f2(x):
    return np.pi*k*2*cos(2*pi*k*x)   #analytical derivative of f1 function

f = zeros(n+1,dtype=float)        # array to store values of f
fd=zeros(n+1,dtype=float)         # array to store values of fd
dfrwd = zeros(n+1,dtype=float)    # array to store values of calulated derivative with forward difference
dbwd = zeros(n+1,dtype=float)     # array to store values of calulated derivative with backward difference
dzd = zeros(n+1,dtype=float)      # array to store values of calulated derivative with central difference
ddo = zeros(n+1,dtype=float)      # array to store values of calulated derivative with 3.order backward difference

for i in range(0,n): # adds values to arrays for x and f(x)

    f[i] = f1(x[i])
    fd[i] = f2(x[i])
step=f[2]-f[1]

# periodic boundary conditions


dfrwd[n] = (f[0]-f[n])/step

dbwd[0] = (f[0]-f[n])/step

dzd[0]=(f[1]-f[n])/(2*step)
dzd[n]=(f[0]-f[n-1])/(2*step)

ddo[n] = (2*f[0]+3*f[n]-6*f[n-1]+f[n-2])/(6*step)
ddo[1] = (2*f[2]+3*f[1]-6*f[0]+f[n])/(6*step)
ddo[0] = (2*f[1]+3*f[0]-6*f[n]+f[n-1])/(6*step)


for i in range(0,n-1): # add values to array for derivative with forward difference
    dfrwd[i] = (f[i+1]-f[i])/step
for i in range(1,n): # add values to array for derivative with backward difference
    dbwd[i] = (f[i]-f[i-1])/step
for i in range(1,n-1): # add values to array for derivative with central difference
    dzd[i] = (f[i+1]-f[i-1])/(2*step)
for i in range(2,n-1): # add values to array for derivative with 3.order backward difference
    ddo[i] = (2*f[i+1]+3*f[i]-6*f[i-1]+f[i-2])/(6*step)




plt.plot(x,fd)
plt.plot(x,dzd)

Answer 1

你可能想计算 step 作为 x[1]-x[0].

---

清理你的导入。您导入并使用 numpy 作为 np 但也是所有内容的总导入。np.sin 在 f1 但只是 cos 在 f2. 你使用 pi 从数学和 np.pi 在同一个公式中。这一切并没有错，但却让人困惑，给调试带来了难度。

你可以使用 numpy 向量操作。f = f1(x); fd = f2(x). 同样，数值导数也可以通过指数移动来计算。或者你也可以将它们实现为 dfrwd = (f1(x+step)-f1(x))/step 等。效率稍低，但读起来比较干净。