fortran 寫的 simpson 積分程式範例

[sjgau02]

採用 Compaq Visual Fortran 95
free format 格式
為了教學使用
程式不是很精簡，不是很有效率

有問題的朋友，歡迎跟我用 e-mail 討論
e-mail: sjga...@gmail.com

program VF0734
implicit none

real x1,x2, y1,y2, ans1, dx
real u0,u1,u2, v0,v1,v2, sum
integer i,j,k, no

! 你可以改變 x1, x2 的值 來試試看 程式的正確性
x1= 0.5
x2= 2.9

! sin() 的積分是 -cos()
! 你可以改變成 不同的函數 來測試看看
y1= -cos(x1)
y2= -cos(x2)
ans1= y2 - y1
write(*,*)'x1, x2= ', x1, x2
write(*,*)'ans1= ', ans1
pause
! -------------------------------------------

! for no=2, no<=1000, no*= 2
! no= 32 的時候 就已經 達到所需要的接確度
! 可見，這個積分程式的效率非常高

no= 2
do while (no .LE. 1000)
dx= (x2 - x1)/no

sum= 0.0
do i=0, (no - 1), 2
u0= x1 + i*dx
v0= sin(u0)

u1= u0 + dx
v1= sin(u1)

u2= u0 + 2*dx
v2= sin(u2)

sum= sum + v0 + 4*v1 + v2
! write(*,*)'i, u0, v0= ', i, u0, v0
end do

sum= (dx/3)*sum
write(*,*)'no, ans1, sum= ', no, ans1, sum
pause

no= no*2
end do
end program VF0734