- #1

fluidistic

Gold Member

- 3,794

- 172

I forgot to say, my program may be a little hard to understand. When you execute it, first it will ask the bounds of the integral, just put 0,1. Then chose 100 for n and 200 for m and we're done.

Program Simpson

implicit none

Real(8) :: x,s,x_i,x_f,r,q

Integer :: k,n,m

Write(*,*)'Enter x_i y x_f'

Read(*,*)x_i,x_f

Write(*,*)'Enter n'

Read(*,*)n

Write(*,*)'Enter m'

Read(*,*)m

Call Simp(x_i,x_f,s,n)

Write(*,*)'When n=100, s is worth',s

Call Simp2(x_i,x_f,r,m)

Write(*,*)'When n=200, s is worth',r

Q=(s-(1-exp(-1.)))/(r-(1-exp(-1.)))

Write(*,*)'The quotient of precision is',Q

Contains

Subroutine Simp(x_i,x_f,s,n)

implicit none

Real(8), intent(in) :: x_i,x_f

Real(8), Intent(out) :: s

Real(8) :: dx

Integer :: n

dx=(x_f-x_i)/n

s=(dx/3.)/(f(x_i)+f(x_f))

Do k=2,n-2,2

x=x_i+k*dx

s=s+2*f(x)

end do

do k=1,n-1,2

x=x_i+k*dx

s=s+4*f(x)

end do

s=(dx/3.)*s

end subroutine

Subroutine Simp2(x_i,x_f,r,m)

implicit none

Real(8), Intent(in) :: x_i,x_f

Real(8), Intent(out) :: r

Real(8) :: dx

Integer :: m

dx=(x_f-x_i)/m

r=(dx/3.)/(f(x_i)+f(x_f))

Do k=2,m-2,2

x=x_i+k*dx

r=r+2*f(x)

end do

do k=1,m-1,2

x=x_i+k*dx

r=r+4*f(x)

end do

r=(dx/3.)*r

end subroutine

Real(8) Function f(x)

Real(8), Intent(in) :: x

f=exp(-x)

end function

end program