midinf.f90

来自「FORTRANvisualfortran常用数值算法集及源码」· F90 代码 · 共 28 行

SUBROUTINE midinf(funk,aa,bb,s,n)
INTEGER n
REAL aa,bb,s,funk
EXTERNAL funk
INTEGER it,j
REAL a,b,ddel,del,sum,tnm,func,x
func(x)=funk(1./x)/x**2
b=1./aa
a=1./bb
if (n==1) then
  s=(b-a)*func(0.5*(a+b))
else
  it=3**(n-2)
  tnm=it
  del=(b-a)/(3.*tnm)
  ddel=del+del
  x=a+0.5*del
  sum=0.
  do j=1,it
    sum=sum+func(x)
    x=x+ddel
    sum=sum+func(x)
    x=x+del
  end do
  s=(s+(b-a)*sum/tnm)/3.
endif
END SUBROUTINE midinf