d3r9.f90

Visual Fortran常用数值算法集

PROGRAM D3R9
!Driver for routine QUAD3D
COMMON XMAX
PARAMETER (PI=3.1415926,NVAL=10)
WRITE(*,*) 'Integral of r^2 over a spherical volume'
WRITE(*,'(/4X,A,T14,A,T24,A)') 'Radius','QUAD3D',&
        'Actual'
DO I=1,NVAL
    XMAX=0.1*I
    XMIN=-XMAX
    CALL QUAD3D(XMIN,XMAX,S)
    WRITE(*,'(1X,F8.2,2F10.4)')&
                       XMAX,S,4.0*PI*(XMAX**5)/5.0
END DO
END PROGRAM
FUNCTION FUNC(X,Y,Z)
  FUNC=X**2+Y**2+Z**2
END FUNCTION FUNC
FUNCTION Z1(X,Y)
COMMON XMAX
  Z1=-SQRT(ABS(XMAX**2-X**2-Y**2))
END FUNCTION Z1
FUNCTION Z2(X,Y)
COMMON XMAX
  Z2=SQRT(ABS(XMAX**2-X**2-Y**2))
END FUNCTION Z2
FUNCTION Y1(X)
COMMON XMAX
  Y1=-SQRT(ABS(XMAX**2-X**2))
END FUNCTION Y1
FUNCTION Y2(X)
COMMON XMAX
  Y2=SQRT(ABS(XMAX**2-X**2))
END FUNCTION Y2
SUBROUTINE QGAUSX(FUNC,A,B,SS)
DIMENSION X(5),W(5)
DATA X/.1488743389,.4333953941,.6794095682,&
       .8650633666,.9739065285/
DATA W/.2955242247,.2692667193,.2190863625,&
       .1494513491,.0666713443/
XM=0.5*(B+A)
XR=0.5*(B-A)
SS=0.
DO J=1,5
    DX=XR*X(J)
    SS=SS+W(J)*(FUNC(XM+DX)+FUNC(XM-DX))
END DO
SS=XR*SS
END SUBROUTINE QGAUSX
SUBROUTINE QGAUSY(FUNC,A,B,SS)
DIMENSION X(5),W(5)
DATA X/.1488743389,.4333953941,.6794095682,&
       .8650633666,.9739065285/
DATA W/.2955242247,.2692667193,.2190863625,&
       .1494513491,.0666713443/
XM=0.5*(B+A)
XR=0.5*(B-A)
SS=0.
DO J=1,5
    DX=XR*X(J)
    SS=SS+W(J)*(FUNC(XM+DX)+FUNC(XM-DX))
END DO
SS=XR*SS
RETURN
END SUBROUTINE QGAUSY
SUBROUTINE QGAUSZ(FUNC,A,B,SS)
DIMENSION X(5),W(5)
DATA X/.1488743389,.4333953941,.6794095682,&
       .8650633666,.9739065285/
DATA W/.2955242247,.2692667193,.2190863625,&
       .1494513491,.0666713443/
XM=0.5*(B+A)
XR=0.5*(B-A)
SS=0.
DO J=1,5
    DX=XR*X(J)
    SS=SS+W(J)*(FUNC(XM+DX)+FUNC(XM-DX))
END DO
SS=XR*SS
END SUBROUTINE QGAUSZ