Google
 

Trailing-Edge - PDP-10 Archives - decuslib10-02 - 43,50145/dqg16.ssp
There are 2 other files named dqg16.ssp in the archive. Click here to see a list.
C                                                                       DG16  10
C     ..................................................................DG16  20
C                                                                       DG16  30
C        SUBROUTINE DQG16                                               DG16  40
C                                                                       DG16  50
C        PURPOSE                                                        DG16  60
C           TO COMPUTE INTEGRAL(FCT(X), SUMMED OVER X FROM XL TO XU)    DG16  70
C                                                                       DG16  80
C        USAGE                                                          DG16  90
C           CALL DQG16 (XL,XU,FCT,Y)                                    DG16 100
C           PARAMETER FCT REQUIRES AN EXTERNAL STATEMENT                DG16 110
C                                                                       DG16 120
C        DESCRIPTION OF PARAMETERS                                      DG16 130
C           XL     - DOUBLE PRECISION LOWER BOUND OF THE INTERVAL.      DG16 140
C           XU     - DOUBLE PRECISION UPPER BOUND OF THE INTERVAL.      DG16 150
C           FCT    - THE NAME OF AN EXTERNAL DOUBLE PRECISION FUNCTION  DG16 160
C                    SUBPROGRAM USED.                                   DG16 170
C           Y      - THE RESULTING DOUBLE PRECISION INTEGRAL VALUE.     DG16 180
C                                                                       DG16 190
C        REMARKS                                                        DG16 200
C           NONE                                                        DG16 210
C                                                                       DG16 220
C        SUBROUTINES AND FUNCTION SUBPROGRAMS REQUIRED                  DG16 230
C           THE EXTERNAL DOUBLE PRECISION FUNCTION SUBPROGRAM FCT(X)    DG16 240
C           MUST BE FURNISHED BY THE USER.                              DG16 250
C                                                                       DG16 260
C        METHOD                                                         DG16 270
C           EVALUATION IS DONE BY MEANS OF 16-POINT GAUSS QUADRATURE    DG16 280
C           FORMULA, WHICH INTEGRATES POLYNOMIALS UP TO DEGREE 31       DG16 290
C           EXACTLY. FOR REFERENCE, SEE                                 DG16 300
C           V.I.KRYLOV, APPROXIMATE CALCULATION OF INTEGRALS,           DG16 310
C           MACMILLAN, NEW YORK/LONDON, 1962, PP.100-111 AND 337-340.   DG16 320
C                                                                       DG16 330
C     ..................................................................DG16 340
C                                                                       DG16 350
      SUBROUTINE DQG16(XL,XU,FCT,Y)                                     DG16 360
C                                                                       DG16 370
C                                                                       DG16 380
      DOUBLE PRECISION XL,XU,Y,A,B,C,FCT                                DG16 390
C                                                                       DG16 400
      A=.5D0*(XU+XL)                                                    DG16 410
      B=XU-XL                                                           DG16 420
      C=.49470046749582497D0*B                                          DG16 430
      Y=.13576229705877047D-1*(FCT(A+C)+FCT(A-C))                       DG16 440
      C=.47228751153661629D0*B                                          DG16 450
      Y=Y+.31126761969323946D-1*(FCT(A+C)+FCT(A-C))                     DG16 460
      C=.43281560119391587D0*B                                          DG16 470
      Y=Y+.47579255841246392D-1*(FCT(A+C)+FCT(A-C))                     DG16 480
      C=.37770220417750152D0*B                                          DG16 490
      Y=Y+.62314485627766936D-1*(FCT(A+C)+FCT(A-C))                     DG16 500
      C=.30893812220132187D0*B                                          DG16 510
      Y=Y+.7479799440828837D-1*(FCT(A+C)+FCT(A-C))                      DG16 520
      C=.22900838882861369D0*B                                          DG16 530
      Y=Y+.8457825969750127D-1*(FCT(A+C)+FCT(A-C))                      DG16 540
      C=.14080177538962946D0*B                                          DG16 550
      Y=Y+.9130170752246179D-1*(FCT(A+C)+FCT(A-C))                      DG16 560
      C=.47506254918818720D-1*B                                         DG16 570
      Y=B*(Y+.9472530522753425D-1*(FCT(A+C)+FCT(A-C)))                  DG16 580
      RETURN                                                            DG16 590
      END                                                               DG16 600