如何表示二阶数值微分的截断误差是一阶数值微分的平方

Question

我试图证明二阶数值微分的截断误差确实给了我们比一阶数值微分更高的精度（考虑到机器 error/round-off 误差 eps()）

这是我在 Julia 中的代码：

function first_order_numerical_D(f)
  function df(x)
    h = sqrt(eps(x))        
    (f(x+h) - f(x))/h
  end
  df
end

function second_order_numerical_D(f)
  function df(x)
    h = sqrt(eps(x))          
    (f(x+h) - f(x-h))/(2.0*h)
  end
  df
end

function analytical_diff_exp(x)
    return exp(x)
end
function analytical_diff_sin(x)
    return cos(x)
end
function analytical_diff_cos(x)
    return -sin(x)
end
function analytical_diff_sqrt(x)
    return 1/(2.0*sqrt(x))
end

function first_order_error_exp(x)
    return first_order_numerical_D(exp)(x) - analytical_diff_exp(x) 
end
function first_order_error_sin(x)
    return first_order_numerical_D(sin)(x) - analytical_diff_sin(x) 
end
function first_order_error_cos(x)
    return first_order_numerical_D(cos)(x) - analytical_diff_cos(x) 
end
function first_order_error_sqrt(x)
    return first_order_numerical_D(sqrt)(x) - analytical_diff_sqrt(x) 
end

function second_order_error_exp(x)
    return second_order_numerical_D(exp)(x) - analytical_diff_exp(x)
end
function second_order_error_sin(x)
    return second_order_numerical_D(sin)(x) - analytical_diff_sin(x)
end
function second_order_error_cos(x)
    return second_order_numerical_D(cos)(x) - analytical_diff_cos(x) 
end
function second_order_error_sqrt(x)
    return second_order_numerical_D(sqrt)(x) - analytical_diff_sqrt(x) 
end

function round_off_err_exp(x)
    return 2.0*sqrt(eps(x))*exp(x)
end
function round_off_err_sin(x)
    return 2.0*sqrt(eps(x))*sin(x)
end
function round_off_err_cos(x)
    return 2.0*sqrt(eps(x))*cos(x)
end
function round_off_err_sqrt(x)
    return 2.0*sqrt(eps(x))*sqrt(x)
end

function first_order_truncation_err_exp(x)
    return abs(first_order_error_exp(x)+round_off_err_exp(x))
end
function first_order_truncation_err_sin(x)
    return abs(first_order_error_sin(x)+round_off_err_sin(x))
end
function first_order_truncation_err_cos(x)
    return abs(first_order_error_cos(x)+round_off_err_cos(x))
end
function first_order_truncation_err_sqrt(x)
    return abs(first_order_error_sqrt(x)+round_off_err_sqrt(x))
end

function second_order_truncation_err_exp(x)
    return abs(second_order_error_exp(x)+0.5*round_off_err_exp(x))
end
function second_order_truncation_err_sin(x)
    return abs(second_order_error_sin(x)+0.5*round_off_err_sin(x))
end
function second_order_truncation_err_cos(x)
    return abs(second_order_error_cos(x)+0.5*round_off_err_cos(x))
end
function second_order_truncation_err_sqrt(x)
    return abs(second_order_error_sqrt(x)+0.5*round_off_err_sqrt(x))
end

如果我减去，这应该给我正确的截断误差（这里我使用加法，因为实际的泰勒展开表明舍入误差和截断误差前面都有一个负号） round_off_err_f项。

分析derivation/proof见： https://www.uio.no/studier/emner/matnat/math/MAT-INF1100/h10/kompendiet/kap11.pdf http://www2.math.umd.edu/~dlevy/classes/amsc466/lecture-notes/differentiation-chap.pdf

但是结果显示：

first_order_truncation_err_exp(0.5), first_order_truncation_err_sin(0.5), first_order_truncation_err_cos(0.5), first_order_truncation_err_sqrt(0.5)

(4.6783240139052204e-8, 1.2990419187857229e-8, 2.8342226290287478e-9, 4.364449135429996e-9)

second_order_truncation_err_exp(0.5), second_order_truncation_err_sin(0.5), second_order_truncation_err_cos(0.5), second_order_truncation_err_sqrt(0.5)

(1.8874426561390482e-8, 7.938850300905947e-9, 4.1240999200086055e-9, 7.45058059692383e-9)

其中：

eps(0.5)=1.1102230246251565e-16

second_order_truncation_err_f()应该是1e-16的顺序而不是1e-8，不知道为什么不行

Answer 1

那是因为你计算的是舍入误差。即：

julia> round_off_err_sqrt(0.5)
1.490116119384766e-8

并且为了查看您的 "second order" 导数之间的差异（通常我看到术语中心差异）您需要选择更大的步长。在文献中通常看到h = cbrt(eps())。

function second_order_numerical_D(f)
    function df(x)
        h = cbrt(eps(x))          
        (f(x+h/2) - f(x-h/2))/h
    end
    df
end
julia> second_order_error_exp(0.5)
1.308131380994837e-11

如何表示二阶数值微分的截断误差是一阶数值微分的平方

How to show truncation error of second-order numerical differentiation is square of the first-order numerical differentiation

numeric

numerical-computing

numerical-analysis

julia