Julia語言 例題1-8 已知5點的座標值 使用牛頓內插法求差除表？

Julia語言 例題1-8 已知5點的座標值 使用牛頓內插法求差除表？
x = [0.0, 1.0  , 2.0  , 3.0  , 4.0 ]
f = [0.0  , -3.0 , 0.0 , 15.0, 48.0]
並計算P(1.5)與P(2.5)之值=?

程式

using Printf

function divdif(x::Array{Float64,1},f::Array{Float64,1})
    #
    # Returns the vector of divided differences for the
    # Newton form of the interpolating polynomial.
    #
    # ON ENTRY:
    # x abscisses, given as a column vector;
    # f ordinates, given as a column vector.
    #
    # ON RETURN:
    # d divided differences
    #
    #
    n = length(x)
    d = deepcopy(f)
    for i=2:n
        for j=1:i-1
            d[i] = (d[j] - d[i])/(x[j] - x[i])
        end
    end
    return d
end

function newtonform(x::Array{Float64,1},
    d::Array{Float64,1},
    xx::Float64)
    #
    # Evaluates the Newton form of the
    # interpolating polynomial, with abscisses
    # in x and divided differences in d at xx.
    #
    n = length(d)
    result = d[n]
    for i=n-1:-1:1
        result = result*(xx - x[i]) + d[i]
    end
    return result
end

function newton(x::Array{Float64,1},f::Array{Float64,1},xx::Float64)
    #
    # implements the interpolation algorithm of Newton
    #
    # ON ENTRY :
    # x abscisses, given as a column vector;
    # f ordinates, given as a column vector;
    # xx point where to evaluate the interpolating
    # polynomial through (x[i],f[i]).
    #
    # ON RETURN :
    # d divided differences, computed from and f;
    # p value of the interpolating polynomial at xx.
    #
    # EXAMPLE :
    divided = divdif(x,f)
    result = newtonform(x,divided,xx)
    return divided, result
end


x = [0.0, 1. , 2. , 3. , 4. ]
f = [0. , -3.0 , 0. , 15.0, 48.0]
xx = 2.5
d, p = newton(x,f,xx)
print("Newton內插法的差除表")
println(d)
s = @sprintf("Newton內插法的Pn(%0.2f)值= %0.5f " , float(xx) , float(p) )
println(s)

xx = 1.5
d, p = newton(x,f,xx)
print("Newton內插法的差除表")
println(d)
s = @sprintf("Newton內插法的Pn(%0.2f)值= %0.5f " , float(xx) , float(p) )
println(s)

輸出結果

Newton內插法的差除表[0.0, -3.0, 3.0, 1.0, -0.0]
Newton內插法的Pn(2.50)值= 5.62500 
Newton內插法的差除表[0.0, -3.0, 3.0, 1.0, -0.0]
Newton內插法的Pn(1.50)值= -2.62500