]>
Commit | Line | Data |
---|---|---|
0795afa3 | 1 | #include "isajet/pilot.h" |
2 | REAL FUNCTION SSXINT(XL,F,XR) | |
3 | C----------------------------------------------------------------------- | |
4 | C Integrate F over (XL,XR) using adaptive Gaussian quadrature. | |
5 | C TOLABS = 1e-35: absolute error for convergence. | |
6 | C TOLREL = 5e-5: relative error for convergence. | |
7 | C TOLBIN = 1e-3: relative bin size limit. Set contribution to | |
8 | C zero if bin falls below this. | |
9 | C----------------------------------------------------------------------- | |
10 | #if defined(CERNLIB_IMPNONE) | |
11 | IMPLICIT NONE | |
12 | #endif | |
13 | #include "isajet/sslun.inc" | |
14 | EXTERNAL F | |
15 | INTEGER NMAX | |
16 | REAL TOLABS,TOLREL,TOLBIN,XMIN,XLIMS(100) | |
17 | REAL R(93),W(93) | |
18 | INTEGER PTR(4),NORD(4) | |
19 | INTEGER ICOUNT,IBAD | |
20 | REAL XL,XR,F | |
21 | REAL AA,BB,TVAL,VAL,TOL | |
22 | INTEGER NLIMS,I,J | |
23 | C | |
24 | DATA PTR,NORD/4,10,22,46, 6,12,24,48/ | |
25 | DATA R/.2386191860,.6612093865,.9324695142, | |
26 | 1 .1252334085,.3678314990,.5873179543,.7699026742,.9041172563, | |
27 | 1 .9815606342,.0640568929,.1911188675,.3150426797,.4337935076, | |
28 | 1 .5454214714,.6480936519,.7401241916,.8200019860,.8864155270, | |
29 | 1 .9382745520,.9747285560,.9951872200,.0323801710,.0970046992, | |
30 | 1 .1612223561,.2247637903,.2873624873,.3487558863,.4086864820, | |
31 | 1 .4669029048,.5231609747,.5772247261,.6288673968,.6778723796, | |
32 | 1 .7240341309,.7671590325,.8070662040,.8435882616,.8765720203, | |
33 | 1 .9058791367,.9313866907,.9529877032,.9705915925,.9841245837, | |
34 | 1 .9935301723,.9987710073,.0162767488,.0488129851,.0812974955, | |
35 | 1 .1136958501,.1459737146,.1780968824,.2100313105,.2417431561, | |
36 | 1 .2731988126,.3043649444,.3352085229,.3656968614,.3957976498, | |
37 | 1 .4254789884,.4547094222,.4834579739,.5116941772,.5393881083, | |
38 | 1 .5665104186,.5930323648,.6189258401,.6441634037,.6687183100, | |
39 | 1 .6925645366,.7156768123,.7380306437,.7596023411,.7803690438, | |
40 | 1 .8003087441,.8194003107,.8376235112,.8549590334,.8713885059, | |
41 | 1 .8868945174,.9014606353,.9150714231,.9277124567,.9393703398, | |
42 | 1 .9500327178,.9596882914,.9683268285,.9759391746,.9825172636, | |
43 | 1 .9880541263,.9925439003,.9959818430,.9983643759,.9996895039/ | |
44 | DATA W/.4679139346,.3607615730,.1713244924, | |
45 | 1 .2491470458,.2334925365,.2031674267,.1600783285,.1069393260, | |
46 | 1 .0471753364,.1279381953,.1258374563,.1216704729,.1155056681, | |
47 | 1 .1074442701,.0976186521,.0861901615,.0733464814,.0592985849, | |
48 | 1 .0442774388,.0285313886,.0123412298,.0647376968,.0644661644, | |
49 | 1 .0639242386,.0631141923,.0620394232,.0607044392,.0591148397, | |
50 | 1 .0572772921,.0551995037,.0528901894,.0503590356,.0476166585, | |
51 | 1 .0446745609,.0415450829,.0382413511,.0347772226,.0311672278, | |
52 | 1 .0274265097,.0235707608,.0196161605,.0155793157,.0114772346, | |
53 | 1 .0073275539,.0031533461,.0325506145,.0325161187,.0324471637, | |
54 | 1 .0323438226,.0322062048,.0320344562,.0318287589,.0315893308, | |
55 | 1 .0313164256,.0310103326,.0306713761,.0302999154,.0298963441, | |
56 | 1 .0294610900,.0289946142,.0284974111,.0279700076,.0274129627, | |
57 | 1 .0268268667,.0262123407,.0255700360,.0249006332,.0242048418, | |
58 | 1 .0234833991,.0227370697,.0219666444,.0211729399,.0203567972, | |
59 | 1 .0195190811,.0186606796,.0177825023,.0168854799,.0159705629, | |
60 | 1 .0150387210,.0140909418,.0131282296,.0121516047,.0111621020, | |
61 | 1 .0101607705,.0091486712,.0081268769,.0070964708,.0060585455, | |
62 | 1 .0050142027,.0039645543,.0029107318,.0018539608,.0007967921/ | |
63 | C | |
64 | DATA TOLABS,TOLREL,TOLBIN,NMAX/1.E-35,5.E-5,1E-3,100/ | |
65 | C | |
66 | SSXINT=0 | |
67 | NLIMS=2 | |
68 | XLIMS(1)=XL | |
69 | XLIMS(2)=XR | |
70 | ICOUNT = 0 | |
71 | IBAD=0 | |
72 | XMIN=TOLBIN*ABS(XR-XL) | |
73 | C | |
74 | 10 AA=(XLIMS(NLIMS)-XLIMS(NLIMS-1))/2 | |
75 | BB=(XLIMS(NLIMS)+XLIMS(NLIMS-1))/2 | |
76 | TVAL=0 | |
77 | DO 15 I=1,3 | |
78 | 15 TVAL=TVAL+W(I)*(F(BB+AA*R(I))+F(BB-AA*R(I))) | |
79 | TVAL=TVAL*AA | |
80 | DO 25 J=1,4 | |
81 | VAL=0 | |
82 | DO 20 I=PTR(J),PTR(J)-1+NORD(J) | |
83 | ICOUNT = ICOUNT + 1 | |
84 | IF (ICOUNT .GT. 1E5) THEN | |
85 | WRITE(LOUT,*) 'ERROR IN SSXINT: SET SSXINT = 0' | |
86 | SSXINT=0 | |
87 | RETURN | |
88 | ENDIF | |
89 | 20 VAL=VAL+W(I)*(F(BB+AA*R(I))+F(BB-AA*R(I))) | |
90 | VAL=VAL*AA | |
91 | TOL=MAX(TOLABS,TOLREL*ABS(VAL)) | |
92 | IF (ABS(TVAL-VAL).LT.TOL) THEN | |
93 | SSXINT=SSXINT+VAL | |
94 | NLIMS=NLIMS-2 | |
95 | IF (NLIMS.NE.0) GO TO 10 | |
96 | GO TO 999 | |
97 | ELSEIF(ABS(AA).LT.XMIN.AND.J.EQ.4) THEN | |
98 | C Bin is too small -- set VAL = 0. -- FEP | |
99 | IBAD=IBAD+1 | |
100 | NLIMS=NLIMS-2 | |
101 | IF (NLIMS.NE.0) GO TO 10 | |
102 | GO TO 999 | |
103 | ENDIF | |
104 | 25 TVAL=VAL | |
105 | IF (NMAX.EQ.2) THEN | |
106 | SSXINT=VAL | |
107 | GO TO 999 | |
108 | END IF | |
109 | IF (NLIMS.GT.(NMAX-2)) THEN | |
110 | WRITE(LOUT,50) SSXINT,NMAX,BB-AA,BB+AA | |
111 | 50 FORMAT (' ERROR IN SSXINT, SSXINT,NMAX,XL,XR=',G15.7,I5,2G15.7) | |
112 | RETURN | |
113 | END IF | |
114 | XLIMS(NLIMS+1)=BB | |
115 | XLIMS(NLIMS+2)=BB+AA | |
116 | XLIMS(NLIMS)=BB | |
117 | NLIMS=NLIMS+2 | |
118 | GO TO 10 | |
119 | C | |
120 | 999 IF(IBAD.GT.0) THEN | |
121 | WRITE(LOUT,*) 'WARNING IN SSXINT: BAD CONVERGENCE FOR ',IBAD, | |
122 | $ ' INTERVALS.' | |
123 | ENDIF | |
124 | RETURN | |
125 | END |