This commit was generated by cvs2svn to compensate for changes in r2,
[u/mrichter/AliRoot.git] / MINICERN / mathlib / gen / d / radapt.F
1 *
2 * $Id$
3 *
4 * $Log$
5 * Revision 1.1.1.1  1996/04/01 15:02:13  mclareni
6 * Mathlib gen
7 *
8 *
9 #include "gen/pilot.h"
10       SUBROUTINE RADAPT(F,A,B,NSEG,RELTOL,ABSTOL,RES,ERR)
11  
12 C     RES = Estimated Integral of F from A to B,
13 C     ERR = Estimated absolute error on RES.
14 C     NSEG  specifies how the adaptation is to be done:
15 C        =0   means use previous binning,
16 C        =1   means fully automatic, adapt until tolerance attained.
17 C        =n>1 means first split interval into n equal segments,
18 C             then adapt as necessary to attain tolerance.
19 C     The specified tolerances are:
20 C            relative: RELTOL ;  absolute: ABSTOL.
21 C        It stops when one OR the other is satisfied, or number of
22 C        segments exceeds NDIM.  Either TOLA or TOLR (but not both!)
23 C        can be set to zero, in which case only the other is used.
24  
25       DOUBLE PRECISION TVALS,TERSS
26       EXTERNAL F
27  
28       PARAMETER (NDIM=100)
29       PARAMETER (R1 = 1, HF = R1/2)
30  
31       DIMENSION XLO(NDIM),XHI(NDIM),TVAL(NDIM),TERS(NDIM)
32       SAVE XLO,XHI,TVAL,TERS,NTER
33       DATA NTER /0/
34  
35       IF(NSEG .LE. 0)  THEN
36        IF(NTER .EQ. 0) THEN
37         NSEGD=1
38         GO TO 2
39        ENDIF
40        TVALS=0
41        TERSS=0
42        DO 1 I = 1,NTER
43        CALL RGS56P(F,XLO(I),XHI(I),TVAL(I),TE)
44        TERS(I)=TE**2
45        TVALS=TVALS+TVAL(I)
46        TERSS=TERSS+TERS(I)
47     1  CONTINUE
48        ROOT= SQRT(2*TERSS)
49        GO TO 9
50       ENDIF
51       NSEGD=MIN(NSEG,NDIM)
52     2 XHIB=A
53       BIN=(B-A)/NSEGD
54       DO 3 I = 1,NSEGD
55       XLO(I)=XHIB
56       XLOB=XLO(I)
57       XHI(I)=XHIB+BIN
58       IF(I .EQ. NSEGD) XHI(I)=B
59       XHIB=XHI(I)
60       CALL RGS56P(F,XLOB,XHIB,TVAL(I),TE)
61       TERS(I)=TE**2
62     3 CONTINUE
63       NTER=NSEGD
64       DO 4 ITER = 1,NDIM
65       TVALS=TVAL(1)
66       TERSS=TERS(1)
67       DO 5 I = 2,NTER
68       TVALS=TVALS+TVAL(I)
69       TERSS=TERSS+TERS(I)
70     5 CONTINUE
71       ROOT= SQRT(2*TERSS)
72       IF(ROOT .LE. ABSTOL .OR. ROOT .LE. RELTOL*ABS(TVALS)) GO TO 9
73       IF(NTER .EQ. NDIM) GO TO 9
74       BIGE=TERS(1)
75       IBIG=1
76       DO 6 I = 2,NTER
77       IF(TERS(I) .GT. BIGE) THEN
78        BIGE=TERS(I)
79        IBIG=I
80       ENDIF
81     6 CONTINUE
82       NTER=NTER+1
83       XHI(NTER)=XHI(IBIG)
84       XNEW=HF*(XLO(IBIG)+XHI(IBIG))
85       XHI(IBIG)=XNEW
86       XLO(NTER)=XNEW
87       CALL RGS56P(F,XLO(IBIG),XHI(IBIG),TVAL(IBIG),TE)
88       TERS(IBIG)=TE**2
89       CALL RGS56P(F,XLO(NTER),XHI(NTER),TVAL(NTER),TE)
90       TERS(NTER)=TE**2
91     4 CONTINUE
92     9 RES=TVALS
93       ERR=ROOT
94       RETURN
95       END