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

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

x        f(x)
=============
0.0      -6.0
0.1      -5.89483
0.3      -5.65014
0.6      -5.17788
1.0      -4.28172
1.1      -3.99583
=============

using Printf

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



function newton_err(x::Array{Float64,1}, da::Float64 , xa::Float64)
    result = da
    for i=n-1:-1:1
        result = result*(xa - x[i])
    end
    return result
end


x= [0.0 , 0.1 , 0.3 ,0.6 , 1.0  , 1.1]
f= [  [-6.0     ,0.0 , 0.0   ,0.0  ,0.0 ,0.0 ],
      [-5.89483 ,0.0 , 0.0   ,0.0  ,0.0 ,0.0 ],
      [-5.65014 ,0.0 , 0.0   ,0.0  ,0.0 ,0.0 ],
      [-5.17788 ,0.0 , 0.0   ,0.0  ,0.0 ,0.0 ],
      [-4.28172 ,0.0 , 0.0   ,0.0  ,0.0 ,0.0 ],
      [-3.99583 ,0.0 , 0.0   ,0.0  ,0.0 ,0.0 ]]

n=length(x)
println("  Divided Difference Table: ")
println("=============================")
for j=2:n
    for i=1:n-j+1
        f[i][j]=(f[i+1][j-1]-f[i][j-1])/(x[i+j-1]-x[i])
    end
end 

print("i\tx(i)\t\tf(i)\t\tf(i,i+1)\tf(i,i+1.i+2),  ......................\n")
for i=1:n
    s=@sprintf("%d\t%8.5f",i,x[i])
    print(s)
    for j=1:n-i+1
        s=@sprintf("\t%8.5f",f[i][j])
        print(s)
    end
    println()
 
end

d=[0.0 for i=1:n]
println(d)
for i=1:n
    d[i]=f[1][i]
end 

s=@sprintf("牛頓前向差除表")
println(s,d)

xa=0.2
p = newtonform(x,d,xa)
s=@sprintf("%0.5f",p)
print("\nNewton內插法的差除表")
println("P(",xa,")=",s)

da=d[n]
p = newton_err(x,da,xa)
s=@sprintf("%0.8f",abs(p))
print("\nNewton內插法的誤差")
println("P(",xa,")=", s)

輸出畫面

  Divided Difference Table: 
=============================
i x(i) f(i) f(i,i+1) f(i,i+1.i+2),  ......................
1 0.00000 -6.00000 1.05170 0.57250 0.21500 0.06302 0.01416
2 0.10000 -5.89483 1.22345 0.70150 0.27802 0.07859
3 0.30000 -5.65014 1.57420 0.95171 0.35661
4 0.60000 -5.17788 2.24040 1.23700
5 1.00000 -4.28172 2.85890
6 1.10000 -3.99583
[0.0, 0.0, 0.0, 0.0, 0.0, 0.0]
牛頓前向差除表[-6.0, 1.0517, 0.5725, 0.215, 0.0630159, 0.0141595]

Newton內插法的差除表P(0.2)=-5.77860

Newton內插法的誤差P(0.2)=0.00000906