[Fortran] 纯文本查看 复制代码
external f
double precision fc,a,b,s
a=0.0
b=1.0
eps=0.000001
call simp(a,b,fc,eps,s)
write(*, '(1x, "s=", d15.6)')s
end
function fc(x)
double precision fc,x
fc=log(1.0d0+x)/(1.0d0+x*x)
return
end
subroutine simp(a,b,fc,eps,t)
double precision a,b,fc,t,h,t1,s1,p,x,t2,s2
n=1
h=b-a
t1=h*(fc(a)+fc(b))/2.0d0
s1=t1
10 p=0.0d0
do 20 k=0,n-1
x=a+(k+0.5d0)*h
p=p+fc(x)
20 continue
t2=(t1+h*p)/2.0d0
s2=(4.0d0*t2-t1)/3.0d0
if (abs(s2-s1) .ge. eps)then
t1=t2
n=n+n
h=h/2.0d0
s1=s2
goto 10
endif
t=s2
return
end