matlab - 尝试在 MATLAB 中实现差异公式

Question

我正在尝试实现差异公式

f'(x) ≈ [ f(x+h) - f(x) ] / h

对 x=1 和 h=10^-k 使用 MATLAB，其中 k=0,...,16。此外，我想绘制错误。

下面是我的代码。我看到错误在 3 左右，我认为它太大了。它应该接近于 0。

syms f(x)
f(x) = tan(x);
df = diff(f,x);
x = 1;
for k = 0:16
    h = 10^-k;
    finitediff = double((f(x+h)-f(x))/h);
    err = double(abs(finitediff-df(x)));
end

Answer 1

您的代码没有任何问题，有限差分公式运行良好，您得到的错误在于计算以下数值的项目:

• 当a和b都非常非常小的时候计算a/b产生的错误。

  <

• 当 a 和 b 非常非常接近时计算 a-b MATLAB 将给出 0。

这是当k从1变为15时的结果:

感谢@CrisLuengo 富有洞察力的评论!

这表明 err 立即下降到几乎为零，并在 h 变为 1e-9 时再次上升(情况 1 发生在这之后) .最后，当h变为1e-14时，df变为0(这里出现情况2)。

我在你的代码中添加了几行代码来展示这一点:

clc;
clear;
format long
syms f(x)
f(x) = tan(x);
h=1;
df = diff(f,x);
double(df(1));
x=1;


range=1:15;
[finitediff,err]=deal(zeros(size(range)));
for k=range
    h=10^-k;
    finitediff(k)=double((f(x+h)-f(x))/h);
    err(k)=double(abs(finitediff(k)-df(1)));
end

figure(1)

subplot(1,2,1)
hold on
plot(err)
plot(err,'bx')
set(gca,'yscale','log')
title('err')

subplot(1,2,2)
hold on
ezplot(df)
axis([0.5 1.5 0 5])
plot(ones(size(range)),finitediff,'rx','MarkerSize',7)
for ii=range
  text(1,finitediff(ii),['h=1e-' num2str(ii)])
end