C$TEST QGGD
C TO RUN AS A MAIN PROGRAM REMOVE NEXT LINE
      SUBROUTINE QGGD
C***********************************************************************
C
C  TEST OF THE PORT PROGRAMS FOR GAUSS QUADRATURE (DOUBLE-PRECISION)
C
C***********************************************************************
      COMMON/CSTAK/D(500)
      DOUBLE PRECISION D
C
      CALL ISTKIN(500,4)
C
      CALL TDG1
      CALL TDGEX
      CALL TDGEX2
      CALL TDGEXA
      CALL TDGLOG
      CALL TDGXA
      CALL TDGXAB
C
      STOP
      END
      SUBROUTINE TDGEX2
      DOUBLE PRECISION X(5),W(5),FEX2,SUM,TRUE,PI,ERR
      DOUBLE PRECISION DATAN2, DSQRT, DEXP
C
      CALL DGQEX2(5,X,W)
      IOUT=I1MACH(2)
      WRITE(IOUT,1)
      DO 10 J=1,5
   10    WRITE(IOUT,2) J, X(J),W(J)
      SUM=0.D0
      DO 20 J=1,5
   20    SUM=SUM+W(J)*FEX2(X(J))
      PI=2.D0*DATAN2(1.D0,0.D0)
      TRUE=DSQRT(PI)*DEXP(-0.25D0)
      ERR=TRUE-SUM
      WRITE(IOUT,3) TRUE,SUM,ERR
      RETURN
    1 FORMAT(///15H TEST OF DGQEX2//30H0ABSCISSAS AND WEIGHTS FOR N=5)
    2 FORMAT(I4,0P2D18.9)
    3 FORMAT(15H0SAMPLE PROBLEM/6H TRUE=,1PD18.9/
     X   6H CALC=,1PD18.9/6H ERR =,1PD11.2)
      END
      DOUBLE PRECISION FUNCTION FEX2(X)
      DOUBLE PRECISION X
      DOUBLE PRECISION DCOS
      FEX2=DCOS(X)
      RETURN
      END
      SUBROUTINE TDGEXA
      DOUBLE PRECISION X(5),W(5),FEXA,SUM,TRUE,PI,ERR
      DOUBLE PRECISION DATAN2, DSQRT
C
      CALL DGQEXA(5,-0.5D0,X,W)
      IOUT=I1MACH(2)
      WRITE(IOUT,1)
      DO 10 J=1,5
   10    WRITE(IOUT,2) J, X(J),W(J)
      SUM=0.D0
      DO 20 J=1,5
   20    SUM=SUM+W(J)*FEXA(X(J))
      PI=2.D0*DATAN2(1.D0,0.D0)
      TRUE=0.5D0*DSQRT(PI)*(1.D0-1.D0/DSQRT(3.D0))
      ERR=TRUE-SUM
      WRITE(IOUT,3) TRUE,SUM,ERR
      RETURN
    1 FORMAT(///15H TEST OF DGQEXA//30H0ABSCISSAS AND WEIGHTS FOR N=5)
    2 FORMAT(I4,0P2D18.9)
    3 FORMAT(15H0SAMPLE PROBLEM/6H TRUE=,1PD18.9/
     X   6H CALC=,1PD18.9/6H ERR =,1PD11.2)
      END
      DOUBLE PRECISION FUNCTION FEXA(X)
      DOUBLE PRECISION X
      DOUBLE PRECISION DEXP
      FEXA=0.5D0*(1.D0-DEXP(-2.D0*X))
      RETURN
      END
      SUBROUTINE TDGLOG
      DOUBLE PRECISION X(5),W(5),FLOG,SUM,TRUE,PI2,ERR
      DOUBLE PRECISION DATAN2
C
      CALL DGQLOG(5,X,W)
      IOUT=I1MACH(2)
      WRITE(IOUT,1)
      DO 10 J=1,5
   10    WRITE(IOUT,2) J, X(J),W(J)
      SUM=0.D0
      DO 20 J=1,5
   20    SUM=SUM+W(J)*FLOG(X(J))
      PI2=DATAN2(1.D0,0.D0)
      TRUE=-(PI2**2/3.D0)
      ERR=TRUE-SUM
      WRITE(IOUT,3) TRUE,SUM,ERR
      RETURN
    1 FORMAT(///15H TEST OF DGQLOG//30H0ABSCISSAS AND WEIGHTS FOR N=5)
    2 FORMAT(I4,0P2D18.9)
    3 FORMAT(15H0SAMPLE PROBLEM/6H TRUE=,1PD18.9/
     X   6H CALC=,1PD18.9/6H ERR =,1PD11.2)
      END
      DOUBLE PRECISION FUNCTION FLOG(X)
      DOUBLE PRECISION X
      FLOG=-1.D0/(1.D0+X)
      RETURN
      END
      SUBROUTINE TDGXA
      DOUBLE PRECISION X(5),W(5),FXA,SUM,TRUE,B(1),PI2,ERR
      DOUBLE PRECISION DATAN2
C
      CALL DGQXA(5,-0.5D0,X,W)
      IOUT=I1MACH(2)
      WRITE(IOUT,1)
      DO 10 J=1,5
   10    WRITE(IOUT,2) J, X(J),W(J)
      SUM=0.D0
      DO 20 J=1,5
   20    SUM=SUM+W(J)*FXA(X(J))
      CALL DBESRJ(1.D0,1,B)
      PI2=DATAN2(1.D0,0.D0)
      TRUE=PI2*B(1)
      ERR=TRUE-SUM
      WRITE(IOUT,3) TRUE,SUM,ERR
      RETURN
    1 FORMAT(///14H TEST OF DGQXA//30H0ABSCISSAS AND WEIGHTS FOR N=5)
    2 FORMAT(I4,0P2D18.9)
    3 FORMAT(15H0SAMPLE PROBLEM/6H TRUE=,1PD18.9/
     X   6H CALC=,1PD18.9/6H ERR =,1PD11.2)
      END
      DOUBLE PRECISION FUNCTION FXA(X)
      DOUBLE PRECISION X
      DOUBLE PRECISION DCOS, DSQRT
      FXA=DCOS(1.D0-X)/DSQRT(2.D0-X)
      RETURN
      END
      SUBROUTINE TDGXAB
      DOUBLE PRECISION X(5),W(5),FXAB,SUM,TRUE,PI,ERR
      DOUBLE PRECISION DATAN2, DSQRT
C
      CALL DGQXAB(5,-0.5D0,0.5D0,X,W)
      IOUT=I1MACH(2)
      WRITE(IOUT,1)
      DO 10 J=1,5
   10    WRITE(IOUT,2) J, X(J),W(J)
      SUM=0.D0
      DO 20 J=1,5
   20    SUM=SUM+W(J)*FXAB(X(J))
      PI=2.D0*DATAN2(1.D0,0.D0)
      TRUE=PI*(1.D0-1.D0/DSQRT(3.D0))
      ERR=TRUE-SUM
      WRITE(IOUT,3) TRUE,SUM,ERR
      RETURN
    1 FORMAT(///15H TEST OF DGQXAB//30H0ABSCISSAS AND WEIGHTS FOR N=5)
    2 FORMAT(I4,0P2D18.9)
    3 FORMAT(15H0SAMPLE PROBLEM/6H TRUE=,1PD18.9/
     X   6H CALC=,1PD18.9/6H ERR =,1PD11.2)
      END
      DOUBLE PRECISION FUNCTION FXAB(X)
      DOUBLE PRECISION X
      FXAB=1.D0/(2.D0+X)
      RETURN
      END
      DOUBLE PRECISION FUNCTION DGAMMA(X)
      DOUBLE PRECISION X
      DOUBLE PRECISION Y, DSQRT, DATAN2
C J. L. BLUE, 16 DEC 77
C DUMMY VERSION - ONLY GOOD FOR INTEGERS OR HALF-INTEGERS
      Y = X
      DGAMMA = 1.D0
   1  IF (Y .LE. 1.D0) GOTO  2
         Y = Y-1.D0
         DGAMMA = Y*DGAMMA
         GOTO  1
   2  IF (Y .NE. 1.D0) GOTO 3
         GOTO  6
   3     IF (Y .NE. 0.5D0) GOTO 4
            DGAMMA = DGAMMA*DSQRT(2.D0*DATAN2(1.D0, 0.D0))
            GOTO  5
C/6S
   4        CALL SETERR(16HDGAMMA - NOT YET, 16, 1, 2)
C/7S
C  4        CALL SETERR('DGAMMA - NOT YET', 16, 1, 2)
C/
   5  CONTINUE
   6  RETURN
      END
      SUBROUTINE TDGEX
      DOUBLE PRECISION X(5),W(5),FEX,SUM,TRUE,PI,ERR
      DOUBLE PRECISION DATAN2
C
      CALL DGQEX(5,X,W)
      IOUT=I1MACH(2)
      WRITE(IOUT,1)
      DO 10 J=1,5
   10    WRITE(IOUT,2) J, X(J),W(J)
      SUM=0.D0
      DO 20 J=1,5
   20    SUM=SUM+W(J)*FEX(X(J))
      PI=2.D0*DATAN2(1.D0,0.D0)
      TRUE=PI**2/6.D0
      ERR=TRUE-SUM
      WRITE(IOUT,3) TRUE,SUM,ERR
      RETURN
    1 FORMAT(///14H TEST OF DGQEX//30H0ABSCISSAS AND WEIGHTS FOR N=5)
    2 FORMAT(I4,0P2D18.9)
    3 FORMAT(15H0SAMPLE PROBLEM/6H TRUE=,1PD18.9/
     X   6H CALC=,1PD18.9/6H ERR =,1PD11.2)
      END
      DOUBLE PRECISION FUNCTION FEX(X)
      DOUBLE PRECISION X
      DOUBLE PRECISION DEXP
      FEX=X/(1.D0-DEXP(-X))
      RETURN
      END
      SUBROUTINE TDG1
      DOUBLE PRECISION X(5),W(5),F1,SUM,TRUE,ERR
      DOUBLE PRECISION DLOG
C
      CALL DGQ1(5,X,W)
      IOUT=I1MACH(2)
      WRITE(IOUT,1)
      DO 10 J=1,5
   10    WRITE(IOUT,2) J, X(J),W(J)
      SUM=0.D0
      DO 20 J=1,5
   20    SUM=SUM+W(J)*F1(X(J))
      TRUE=DLOG(3.D0)
      ERR=TRUE-SUM
      WRITE(IOUT,3) TRUE,SUM,ERR
      RETURN
    1 FORMAT(///13H TEST OF DGQ1//30H0ABSCISSAS AND WEIGHTS FOR N=5)
    2 FORMAT(I4,0P2D18.9)
    3 FORMAT(15H0SAMPLE PROBLEM/6H TRUE=,1PD18.9/
     X   6H CALC=,1PD18.9/6H ERR =,1PD11.2)
      END
      DOUBLE PRECISION FUNCTION F1(X)
      DOUBLE PRECISION X
      F1=1.D0/(2.D0+X)
      RETURN
       END
C OUTPUT FROM TEST....
C
C TEST OF DGQ1
C ABSCISSAS AND WEIGHTS FOR N=5
C   1  -0.906179852D 00   0.236926886D 00
C   2  -0.538469312D 00   0.478628672D 00
C   3   0.                0.568888889D 00
C   4   0.538469312D 00   0.478628672D 00
C   5   0.906179852D 00   0.236926886D 00
C SAMPLE PROBLEM
C TRUE=   1.098612288D 00
C CALC=   1.098609242D 00
C ERR =   3.05D-06
C TEST OF DGQEX
C ABSCISSAS AND WEIGHTS FOR N=5
C   1   0.263560320D 00   0.521755609D 00
C   2   0.141340307D 01   0.398666812D 00
C   3   0.359642579D 01   0.759424499D-01
C   4   0.708580998D 01   0.361175868D-02
C   5   0.126408009D 02   0.233699724D-04
C SAMPLE PROBLEM
C TRUE=   1.644934067D 00
C CALC=   1.644924428D 00
C ERR =   9.64D-06
C TEST OF DGQEX2
C ABSCISSAS AND WEIGHTS FOR N=5
C   1  -0.202018288D 01   0.199532420D-01
C   2  -0.958572467D 00   0.393619324D 00
C   3   0.                0.945308723D 00
C   4   0.958572467D 00   0.393619324D 00
C   5   0.202018288D 01   0.199532420D-01
C SAMPLE PROBLEM
C TRUE=   1.380388448D 00
C CALC=   1.380390076D 00
C ERR =  -1.63D-06
C TEST OF DGQEXA
C ABSCISSAS AND WEIGHTS FOR N=5
C   1   0.117581320D 00   0.122172528D 01
C   2   0.107456201D 01   0.480277222D 00
C   3   0.308593744D 01   0.677487891D-01
C   4   0.641472972D 01   0.268729148D-02
C   5   0.118071896D 02   0.152808657D-04
C SAMPLE PROBLEM
C TRUE=   3.745635718D-01
C CALC=   3.753077036D-01
C ERR =  -7.44D-04
C TEST OF DGQLOG
C ABSCISSAS AND WEIGHTS FOR N=5
C   1   0.291344723D-01   0.297893472D 00
C   2   0.173977214D 00   0.349776227D 00
C   3   0.411702521D 00   0.234488291D 00
C   4   0.677314176D 00   0.989304601D-01
C   5   0.894771366D 00   0.189115521D-01
C SAMPLE PROBLEM
C TRUE=  -8.224670336D-01
C CALC=  -8.224670176D-01
C ERR =  -1.59D-08
C TEST OF DGQXA
C ABSCISSAS AND WEIGHTS FOR N=5
C   1   0.221635689D-01   0.591048448D 00
C   2   0.187831568D 00   0.538533440D 00
C   3   0.461597363D 00   0.438172723D 00
C   4   0.748334630D 00   0.298902697D 00
C   5   0.948493926D 00   0.133342689D 00
C SAMPLE PROBLEM
C TRUE=   1.201969715D 00
C CALC=   1.201969712D 00
C ERR =   3.93D-09
C TEST OF DGQXAB
C ABSCISSAS AND WEIGHTS FOR N=5
C   1  -0.841253536D 00   0.906757696D-01
C   2  -0.415415014D 00   0.333914163D 00
C   3   0.142314838D 00   0.652488710D 00
C   4   0.654860736D 00   0.945254246D 00
C   5   0.959492979D 00   0.111925977D 01
C SAMPLE PROBLEM
C TRUE=   1.327793289D 00
C CALC=   1.327791436D 00
C ERR =   1.85D-06
