Julia語言 例題4-4利用辛普森 理則 （Simpson's Rule)與梯形法 計算定積分

Julia語言 例題4-4利用辛普森 理則 （Simpson's Rule) 與梯形法
#        (a) 計算 exp( 1/x)  在[1, 2]的定積分
#        (b) 計算 ( 1/x)       在[1, 2]的定積分

#========================================================
/* ex4-4.jl based on Simpson's Rule to compute
 * definite integral with domain [a,b] and
 * n even-grid. n must be even.
 */
========================================================#
using Printf

function F1(x::Float64) #// 欲微分函數
    return  (exp(1/x))
end


function F2(x::Float64) #// 欲微分函數
    return  (1/x)
end

a=1.0
b=2.0
n=4
m=n/2
h=(b-a)/n
sum1=0.0
sum2=0.0
for i=1:2*m-1
    x=a+i*h;
    if(i%2==0)
        sum2=sum2+F1(x)
        s=@sprintf("i=%2d ,x=%0.3f --- F(x)=%0.6f",i,x,F1(x))
        println(s)
    else
        sum1=sum1+F1(x)
        s=@sprintf("i=%2d ,x=%0.3f --- F(x)=%0.6f",i,x,F1(x))
        println(s)
    end   
end
sn=  (h/3.0)*(F1(a)+F1(b)+2.0*sum2+4.0*sum1)
print("辛普森積分法  ")
s=@sprintf("S%d=%lf\n",n,sn)
println(s)

println("---------------------------------")
n=10
a=1.0
b=2.0
n=4
h=(b-a)/n
x=a
result=0.0
for i=1:n-1
    x=x+h
    s=@sprintf("i=%2d , x=%0.2f ---- f(x)=%0.6f",i,x,F1(abs(x)))
    println(s)
    result=result+F1(abs(x))
end

tn=(h/2.0)*(F1(abs(a))+F1(abs(b))+2.0*result)

print("梯形積分法  ")
s=@sprintf("T%d=%0.7lf\n",n,tn)
println(s)


println("---------------------------------")
n=10
m=n/2
h=(b-a)/n
sum1=0.0
sum2=0.0
for i=1:2*m-1
    x=a+i*h;
    if(i%2==0)
        sum2=sum2+F2(x)
        s=@sprintf("i=%2d ,x=%0.3f --- F(x)=%0.6f",i,x,F2(x))
        println(s)
    else
        sum1=sum1+F2(x)
        s=@sprintf("i=%2d ,x=%0.3f --- F(x)=%0.6f",i,x,F2(x))
        println(s)
    end   
end
sn= (h/3.0)*(F2(a)+F2(b)+2.0*sum2+4.0*sum1)
print("辛普森積分法  ")
s=@sprintf("S%d=%lf\n",n,sn)
println(s)


println("---------------------------------")
n=10
a=1.0
b=2.0
n=10
h=(b-a)/n
x=a
result=0.0
for i=1:n-1
    x=x+h
    s=@sprintf("i=%2d , x=%0.2f ---- f(x)=%0.6f",i,x,F2(abs(x)))
    println(s)
    result=result+F2(abs(x))
end

tn=(h/2.0)*(F2(abs(a))+F2(abs(b))+2.0*result)

print("梯形積分法  ")
s=@sprintf("T%d=%0.7lf\n",n,tn)
println(s)

輸出畫面

i= 1 ,x=1.250 --- F(x)=2.225541
i= 2 ,x=1.500 --- F(x)=1.947734
i= 3 ,x=1.750 --- F(x)=1.770795
辛普森積分法  S4=2.020651

---------------------------------
i= 1 , x=1.25 ---- f(x)=2.225541
i= 2 , x=1.50 ---- f(x)=1.947734
i= 3 , x=1.75 ---- f(x)=1.770795
梯形積分法  T4=2.0318929

---------------------------------
i= 1 ,x=1.100 --- F(x)=0.909091
i= 2 ,x=1.200 --- F(x)=0.833333
i= 3 ,x=1.300 --- F(x)=0.769231
i= 4 ,x=1.400 --- F(x)=0.714286
i= 5 ,x=1.500 --- F(x)=0.666667
i= 6 ,x=1.600 --- F(x)=0.625000
i= 7 ,x=1.700 --- F(x)=0.588235
i= 8 ,x=1.800 --- F(x)=0.555556
i= 9 ,x=1.900 --- F(x)=0.526316
辛普森積分法  S10=0.693150

---------------------------------
i= 1 , x=1.10 ---- f(x)=0.909091
i= 2 , x=1.20 ---- f(x)=0.833333
i= 3 , x=1.30 ---- f(x)=0.769231
i= 4 , x=1.40 ---- f(x)=0.714286
i= 5 , x=1.50 ---- f(x)=0.666667
i= 6 , x=1.60 ---- f(x)=0.625000
i= 7 , x=1.70 ---- f(x)=0.588235
i= 8 , x=1.80 ---- f(x)=0.555556
i= 9 , x=1.90 ---- f(x)=0.526316
梯形積分法  T10=0.6937714