hijing/src/gauss1.f

0001
0002 C
0003 C
0004 C*********GAUSSIAN ONE-DIMENSIONAL INTEGRATION PROGRAM*************
0005 C
0006         FUNCTION GAUSS1(F,A,B,EPS)
0007         EXTERNAL F
0008         DIMENSION W(12),X(12)
0009         DATA CONST/1.0E-12/
0010         DATA W/0.1012285,.2223810,.3137067,.3623838,.0271525,
0011      &         .0622535,0.0951585,.1246290,.1495960,.1691565,
0012      &         .1826034,.1894506/
0013         DATA X/0.9602899,.7966665,.5255324,.1834346,.9894009,
0014      &         .9445750,0.8656312,.7554044,.6178762,.4580168,
0015      &         .2816036,.0950125/
0016         DELTA=CONST*ABS(A-B)
0017         GAUSS1=0.0
0018         AA=A
0019 5       Y=B-AA
0020         IF(ABS(Y).LE.DELTA) RETURN
0021 2       BB=AA+Y
0022         C1=0.5*(AA+BB)
0023         C2=C1-AA
0024         S8=0.0
0025         S16=0.0
0026         DO 1 I=1,4
0027         U=X(I)*C2
0028 1       S8=S8+W(I)*(F(C1+U)+F(C1-U))
0029         DO 3 I=5,12
0030         U=X(I)*C2
0031 3       S16=S16+W(I)*(F(C1+U)+F(C1-U))
0032         S8=S8*C2
0033         S16=S16*C2
0034         IF(ABS(S16-S8).GT.EPS*(1.+ABS(S16))) GOTO 4
0035         GAUSS1=GAUSS1+S16
0036         AA=BB
0037         GOTO 5
0038 4       Y=0.5*Y
0039         IF(ABS(Y).GT.DELTA) GOTO 2
0040         WRITE(6,7)
0041         GAUSS1=0.0
0042         RETURN
0043 7       FORMAT(1X,'GAUSS1....TOO HIGH ACURACY REQUIRED')
0044         END