1
C Last modification on March 10th 2001 by M.S.
2
C ==================================================================
3
C ================= PROGRAM HDECAY: COMMENTS =======================
4
C ==================================================================
11
C This program calculates the total decay widths and the branching
12
C ratios of the C Standard Model Higgs boson (HSM) as well as those
13
C of the neutral (HL= the light CP-even, HH= the heavy CP-even, HA=
14
C the pseudoscalar) and the charged (HC) Higgs bosons of the Minimal
15
C Supersymmetric extension of the Standard Model (MSSM). It includes:
17
C - All the decay channels which are kinematically allowed and which
18
C have branching ratios larger than 10**(-4).
20
C - All QCD corrections to the fermionic and gluonic decay modes.
21
C Most of these corrections are mapped into running masses in a
22
C consistent way with some freedom for including high order terms.
24
C - Below--threshold three--body decays with off--shell top quarks
25
C or ONE off-shell gauge boson, as well as some decays with one
26
C off-shell Higgs boson in the MSSM.
28
C - Double off-shell decays: HSM,HL,HH --> W*W*,Z*Z* -->4 fermions,
29
C which could be important for Higgs masses close to MW or MZ.
31
C - In the MSSM, the radiative corrections with full squark mixing and
32
C uses the RG improved values of Higgs masses and couplings with the
33
C main NLO corrections implemented (based on M.Carena, M. Quiros and
34
C C.E.M. Wagner, Nucl. Phys. B461 (1996) 407, hep-ph/9508343).
36
C - In the MSSM, all the decays into CHARGINOS, NEUTRALINOS, SLEPTONS
37
C and SQUARKS (with mixing in the stop and sbottom sectors).
39
C - Chargino, slepton and squark loops in the 2 photon decays and squark
40
C loops in the gluonic decays (including QCD corrections).
42
C ===================================================================
43
C This program has been written by A.Djouadi, J.Kalinowski and M.Spira.
44
C For details on how to use the program see: Comp. Phys. Commun. 108
45
C (1998) 56, hep-ph/9704448. For any question, comment, suggestion or
46
C complaint, please contact us at:
47
C djouadi@lpm.univ-montp2.fr
49
C Michael.Spira@cern.ch
52
C ================ IT USES AS INPUT PARAMETERS:
54
C IHIGGS: =0: CALCULATE BRANCHING RATIOS OF SM HIGGS BOSON
55
C =1: CALCULATE BRANCHING RATIOS OF MSSM h BOSON
56
C =2: CALCULATE BRANCHING RATIOS OF MSSM H BOSON
57
C =3: CALCULATE BRANCHING RATIOS OF MSSM A BOSON
58
C =4: CALCULATE BRANCHING RATIOS OF MSSM H+ BOSON
59
C =5: CALCULATE BRANCHING RATIOS OF ALL MSSM HIGGS BOSONS
61
C TGBET: TAN(BETA) FOR MSSM
62
C MABEG: START VALUE OF M_A FOR MSSM AND M_H FOR SM
63
C MAEND: END VALUE OF M_A FOR MSSM AND M_H FOR SM
64
C NMA: NUMBER OF ITERATIONS FOR M_A
65
C ALS(MZ): VALUE FOR ALPHA_S(M_Z)
66
C MSBAR(1): MSBAR MASS OF STRANGE QUARK AT SCALE Q=1 GEV
68
C MB: BOTTOM POLE MASS
72
C ALPH: INVERSE QED COUPLING
78
C VUS: CKM PARAMETER V_US
79
C VCB: CKM PARAMETER V_CB
80
C VUB/VCB: RATIO V_UB/V_CB
81
C 1ST AND 2ND GENERATION:
82
C MSL1: SUSY BREAKING MASS PARAMETERS OF LEFT HANDED SLEPTONS
83
C MER1: SUSY BREAKING MASS PARAMETERS OF RIGHT HANDED SLEPTONS
84
C MQL1: SUSY BREAKING MASS PARAMETERS OF LEFT HANDED SUPS
85
C MUR1: SUSY BREAKING MASS PARAMETERS OF RIGHT HANDED SUPS
86
C MDR1: SUSY BREAKING MASS PARAMETERS OF RIGHT HANDED SDOWNS
88
C MSL: SUSY BREAKING MASS PARAMETERS OF LEFT HANDED STAUS
89
C MER: SUSY BREAKING MASS PARAMETERS OF RIGHT HANDED STAUS
90
C MSQ: SUSY BREAKING MASS PARAMETERS OF LEFT HANDED STOPS
91
C MUR: SUSY BREAKING MASS PARAMETERS OF RIGHT HANDED STOPS
92
C MDR: SUSY BREAKING MASS PARAMETERS OF RIGHT HANDED SBOTTOMS
93
C AL: STAU TRILINEAR SOFT BREAKING TERMS
94
C AU: STOP TRILINEAR SOFT BREAKING TERMS.
95
C AD: SBOTTOM TRILINEAR SOFT BREAKING TERMS.
96
C MU: SUSY HIGGS MASS PARAMETER
97
C M2: GAUGINO MASS PARAMETER.
99
C NNLO (M): =0: USE O(ALPHA_S) FORMULA FOR POLE MASS --> MSBAR MASS
100
C =1: USE O(ALPHA_S**2) FORMULA FOR POLE MASS --> MSBAR MASS
102
C ON-SHELL: =0: INCLUDE OFF_SHELL DECAYS H,A --> T*T*, A --> Z*H,
103
C H --> W*H+,Z*A, H+ --> W*A, W*H, T*B
104
C =1: EXCLUDE THE OFF-SHELL DECAYS ABOVE
106
C ON-SH-WZ: =0: INCLUDE DOUBLE OFF-SHELL PAIR DECAYS PHI --> W*W*,Z*Z*
107
C =1: INCLUDE ONLY SINGLE OFF-SHELL DECAYS PHI --> W*W,Z*Z
109
C IPOLE: =0 COMPUTES RUNNING HIGGS MASSES (FASTER)
110
C =1 COMPUTES POLE HIGGS MASSES
112
C OFF-SUSY: =0: INCLUDE DECAYS (AND LOOPS) INTO SUPERSYMMETRIC PARTICLES
113
C =1: EXCLUDE DECAYS (AND LOOPS) INTO SUPERSYMMETRIC PARTICLES
115
C INIDEC: =0: PRINT OUT SUMS OF CHARGINO/NEUTRALINO/SFERMION DECAYS
116
C =1: PRINT OUT INDIVIDUAL CHARGINO/NEUTRALINO/SFERMION DECAYS
118
C NF-GG: NUMBER OF LIGHT FLAVORS INCLUDED IN THE GLUONIC DECAYS
119
C PHI --> GG* --> GQQ (3,4 OR 5)
121
C =======================================================================
122
C ============== BEGINNING OF THE MAIN PROGRAM ==========================
123
C =======================================================================
127
IMPLICIT DOUBLE PRECISION (A-H,O-Z)
128
PARAMETER(K=6,NI=87,NSA=85,NSB=86,NLA=88,NLB=89,NHA=90,NHB=91,
129
. NHC=92,NAA=93,NAB=94,NCA=95,NCB=96,NRA=97,NRB=98,
130
. NSUSYL=81,NSUSYA=82,NSUSYH=83,NSUSYC=84,NPAR=80,
131
. NSUSYLA=79,NSUSYLB=78,NSUSYLC=77,NSUSYLD=76,NSUSYLE=75,
132
. NSUSYLF=59,NSUSYHF=58,
133
. NSUSYHA=74,NSUSYHB=73,NSUSYHC=72,NSUSYHD=71,NSUSYHE=70,
134
. NSUSYAA=69,NSUSYAB=68,NSUSYAC=67,NSUSYAD=66,NSUSYAE=65,
135
. NSUSYCA=64,NSUSYCB=63,NSUSYCC=62,NSUSYCD=61,NSUSYCE=60)
136
DIMENSION GMN(4),XMN(4),GMC(2),GMST(2),GMSB(2),GMSL(2),
137
. GMSU(2),GMSD(2),GMSE(2),GMSN(2)
138
DIMENSION HLBRSC(2,2),HLBRSN(4,4),HHBRSC(2,2),HHBRSN(4,4),
139
. HABRSC(2,2),HABRSN(4,4),HCBRSU(2,4),
140
. HHBRST(2,2),HHBRSB(2,2),HCBRSTB(2,2)
141
DIMENSION AC1(2,2),AC2(2,2),AC3(2,2),
142
. AN1(4,4),AN2(4,4),AN3(4,4),
143
. ACNL(2,4),ACNR(2,4)
144
DIMENSION GLTT(2,2),GLBB(2,2),GHTT(2,2),GHBB(2,2),GCTB(2,2),
145
. GLEE(2,2),GHEE(2,2),GCEN(2,2)
146
COMMON/MASSES/AMS,AMC,AMB,AMT
148
COMMON/PARAM/GF,ALPH,AMTAU,AMMUON,AMZ,AMW
149
COMMON/CKMPAR/VUS,VCB,VUB
150
COMMON/HMASS/AMSM,AMA,AML,AMH,AMCH,AMAR
151
COMMON/BREAK/AMEL,AMER,AMSQ,AMUR,AMDR,AL,AU,AD,AMU,AM2
152
COMMON/SFER1ST/AMQL1,AMUR1,AMDR1,AMEL1,AMER1
153
COMMON/WZWDTH/GAMC0,GAMT0,GAMT1,GAMW,GAMZ
154
COMMON/COUP/GAT,GAB,GLT,GLB,GHT,GHB,GZAH,GZAL,
155
. GHHH,GLLL,GHLL,GLHH,GHAA,GLAA,GLVV,GHVV,
157
COMMON/ALS/XLAMBDA,AMC0,AMB0,AMT0,N0
158
COMMON/FLAG/IHIGGS,NNLO,IPOLE
159
COMMON/ONSHELL/IONSH,IONWZ,IOFSUSY
161
COMMON/WIDTHSM/SMBRB,SMBRL,SMBRM,SMBRS,SMBRC,SMBRT,SMBRG,SMBRGA,
162
. SMBRZGA,SMBRW,SMBRZ,SMWDTH
163
COMMON/WIDTHA/ABRB,ABRL,ABRM,ABRS,ABRC,ABRT,ABRG,ABRGA,ABRZGA,
165
COMMON/WIDTHHL/HLBRB,HLBRL,HLBRM,HLBRS,HLBRC,HLBRT,HLBRG,HLBRGA,
166
. HLBRZGA,HLBRW,HLBRZ,HLBRA,HLBRAZ,HLBRHW,HLWDTH
167
COMMON/WIDTHHH/HHBRB,HHBRL,HHBRM,HHBRS,HHBRC,HHBRT,HHBRG,HHBRGA,
168
. HHBRZGA,HHBRW,HHBRZ,HHBRH,HHBRA,HHBRAZ,HHBRHW,
170
COMMON/WIDTHHC/HCBRB,HCBRL,HCBRM,HCBRBU,HCBRS,HCBRC,HCBRT,HCBRW,
172
COMMON/WISUSY/HLBRSC,HLBRSN,HHBRSC,HHBRSN,HABRSC,HABRSN,HCBRSU,
173
. HLBRCHT,HHBRCHT,HABRCHT,HLBRNET,HHBRNET,HABRNET,
174
. HCBRCNT,HLBRSL,HHBRSL,HCBRSL,HABRSL,HABRST,HABRSB,
175
. HHBRSQ,HHBRST,HHBRSB,HHBRSQT,HCBRSQ,HCBRSTB,
176
. HCBRSQT,HLBRSQ,HLBRSQT
177
COMMON/WISFER/BHLSLNL,BHLSLEL,BHLSLER,BHLSQUL,BHLSQUR,BHLSQDL,
178
. BHLSQDR,BHLST(2,2),BHLSB(2,2),BHLSTAU(2,2),
179
. BHHSLNL,BHHSLEL,BHHSLER,BHHSQUL,BHHSQUR,BHHSQDL,
180
. BHHSQDR,BHHST(2,2),BHHSB(2,2),BHHSTAU(2,2),
181
. BHASTAU,BHASB,BHAST,
182
. BHCSL00,BHCSL11,BHCSL21,BHCSQ,BHCSTB(2,2)
183
COMMON/SMASS/GMN,XMN,GMC,GMST,GMSB,GMSL,GMSU,GMSD,GMSE,GMSN
185
double precision hdals, hdmhbeg, hdmhend, tgbet
186
common /hdparms/ hdals, hdmhbeg, hdmhend, tgbet
194
c OPEN(NI,FILE='hdecay.in')
195
c OPEN(NPAR,FILE='br.input')
238
c READ(NI,101)IOFSUSY
239
c READ(NI,101)INDIDEC
270
IF(NFGG.GT.5.OR.NFGG.LT.3)THEN
271
C WRITE(6,*)'NF-GG NOT VALID. TAKING THE DEFAULT NF-GG = 3....'
275
100 FORMAT(10X,G30.20)
283
XLAMBDA=XITLA(NLOOP,ALSMZ,ACC)
287
C--INITIALIZE COEFFICIENTS FOR POLYLOGARITHMS
294
AMAR = AMABEG + (AMAEND-AMABEG)/(NMA-1D0)*(II-1D0)
301
***********************************************
310
C =====================================================================
311
C =========== BEGINNING OF THE SUBROUTINE FOR THE DECAYS ==============
312
C !!!!!!!!!!!!!! Any change below this line is at your own risk!!!!!!!!
313
C =====================================================================
315
SUBROUTINE HDEC(TGBET)
316
IMPLICIT DOUBLE PRECISION (A-H,O-Z)
317
DOUBLE PRECISION LAMB
318
DIMENSION XX(4),YY(4)
319
DIMENSION AMCHAR(2),AMNEUT(4),XMNEUT(4),
320
. AC1(2,2),AC2(2,2),AC3(2,2),
321
. AN1(4,4),AN2(4,4),AN3(4,4),
322
. ACNL(2,4),ACNR(2,4),
323
. AMST(2),AMSB(2),AMSL(2),
324
. AMSU(2),AMSD(2),AMSE(2),AMSN(2),
325
. GLTT(2,2),GLBB(2,2),GLEE(2,2),
326
. GHTT(2,2),GHBB(2,2),GHEE(2,2),
327
. GCTB(2,2),GCEN(2,2)
328
DIMENSION GMST(2),GMSB(2),GMSL(2),GMSU(2),GMSD(2),GMSE(2),
330
DIMENSION HLBRSC(2,2),HLBRSN(4,4),HHBRSC(2,2),
331
. HHBRSN(4,4),HABRSC(2,2),HABRSN(4,4),HCBRSU(2,4),
332
. HHBRST(2,2),HHBRSB(2,2),HCBRSTB(2,2)
333
DIMENSION WHLCH(2,2),WHLNE(4,4),WHHCH(2,2),WHHNE(4,4),
334
. WHACH(2,2),WHANE(4,4),WHCCN(2,4),
335
. WHHST(2,2),WHHSB(2,2),WHHSTAU(2,2),WHCSTB(2,2),
336
. WHLST(2,2),WHLSB(2,2),WHLSTAU(2,2)
337
COMPLEX*16 CF,CG,CI1,CI2,CA,CB,CTT,CTB,CTC,CTW,CLT,CLB,CLW,
338
. CAT,CAB,CAC,CAW,CAH,CTH,CLH,CX1,CX2,CAX1,CAX2,CTL,CAL,
339
. CSL,CSQ,CSB1,CSB2,CST1,CST2,CSL1,CSL2,
340
. CXL,CXQ,CXB1,CXB2,CXT1,CXT2,CXL1,CXL2
341
COMPLEX*16 CSEL,CSER,CSUL,CSUR,CSDL,CSDR,
342
. CXEL,CXER,CXUL,CXUR,CXDL,CXDR
343
COMMON/HMASS/AMSM,AMA,AML,AMH,AMCH,AMAR
345
COMMON/MASSES/AMS,AMC,AMB,AMT
346
COMMON/ALS/XLAMBDA,AMC0,AMB0,AMT0,N0
347
COMMON/PARAM/GF,ALPH,AMTAU,AMMUON,AMZ,AMW
348
COMMON/CKMPAR/VUS,VCB,VUB
349
COMMON/BREAK/AMEL,AMER,AMSQ,AMUR,AMDR,AL,AU,AD,AMU,AM2
350
COMMON/WZWDTH/GAMC0,GAMT0,GAMT1,GAMW,GAMZ
351
COMMON/ONSHELL/IONSH,IONWZ,IOFSUSY
353
COMMON/FLAG/IHIGGS,NNLO,IPOLE
354
COMMON/WIDTHSM/SMBRB,SMBRL,SMBRM,SMBRS,SMBRC,SMBRT,SMBRG,SMBRGA,
355
. SMBRZGA,SMBRW,SMBRZ,SMWDTH
356
COMMON/WIDTHA/ABRB,ABRL,ABRM,ABRS,ABRC,ABRT,ABRG,ABRGA,ABRZGA,
358
COMMON/WIDTHHL/HLBRB,HLBRL,HLBRM,HLBRS,HLBRC,HLBRT,HLBRG,HLBRGA,
359
. HLBRZGA,HLBRW,HLBRZ,HLBRA,HLBRAZ,HLBRHW,HLWDTH
360
COMMON/WIDTHHH/HHBRB,HHBRL,HHBRM,HHBRS,HHBRC,HHBRT,HHBRG,HHBRGA,
361
. HHBRZGA,HHBRW,HHBRZ,HHBRH,HHBRA,HHBRAZ,HHBRHW,
363
COMMON/WIDTHHC/HCBRB,HCBRL,HCBRM,HCBRBU,HCBRS,HCBRC,HCBRT,HCBRW,
365
COMMON/WISUSY/HLBRSC,HLBRSN,HHBRSC,HHBRSN,HABRSC,HABRSN,HCBRSU,
366
. HLBRCHT,HHBRCHT,HABRCHT,HLBRNET,HHBRNET,HABRNET,
367
. HCBRCNT,HLBRSL,HHBRSL,HCBRSL,HABRSL,HABRST,HABRSB,
368
. HHBRSQ,HHBRST,HHBRSB,HHBRSQT,HCBRSQ,HCBRSTB,
369
. HCBRSQT,HLBRSQ,HLBRSQT
370
COMMON/WISFER/BHLSLNL,BHLSLEL,BHLSLER,BHLSQUL,BHLSQUR,BHLSQDL,
371
. BHLSQDR,BHLST(2,2),BHLSB(2,2),BHLSTAU(2,2),
372
. BHHSLNL,BHHSLEL,BHHSLER,BHHSQUL,BHHSQUR,BHHSQDL,
373
. BHHSQDR,BHHST(2,2),BHHSB(2,2),BHHSTAU(2,2),
374
. BHASTAU,BHASB,BHAST,
375
. BHCSL00,BHCSL11,BHCSL21,BHCSQ,BHCSTB(2,2)
376
COMMON/SMASS/AMNEUT,XMNEUT,AMCHAR,AMST,AMSB,AMSL,
377
. AMSU,AMSD,AMSE,AMSN
378
COMMON/COUP/GAT,GAB,GLT,GLB,GHT,GHB,GZAH,GZAL,
379
. GHHH,GLLL,GHLL,GLHH,GHAA,GLAA,GLVV,GHVV,
381
HVV(X,Y)= GF/(4.D0*PI*DSQRT(2.D0))*X**3/2.D0*BETA(Y)
382
. *(1.D0-4.D0*Y+12.D0*Y**2)
383
AFF(X,Y)= GF/(4*PI*DSQRT(2.D0))*X**3*Y*(BETA(Y))
384
HFF(X,Y)= GF/(4*PI*DSQRT(2.D0))*X**3*Y*(BETA(Y))**3
385
CFF(Z,TB,X,Y)= GF/(4*PI*DSQRT(2.D0))*Z**3*LAMB(X,Y)
386
. *((1.D0-X-Y)*(X*TB**2+Y/TB**2)-4.D0*X*Y)
387
HV(V)=3.D0*(1.D0-8.D0*V+20.D0*V**2)/DSQRT((4.D0*V-1.D0))
388
. *DACOS((3.D0*V-1.D0)/2.D0/DSQRT(V**3))
389
. -(1.D0-V)*(47.D0/2.D0*V-13.D0/2.D0+1.D0/V)
390
. -3.D0/2.D0*(1.D0-6.D0*V+4.D0*V**2)*DLOG(V)
391
HVH(X,Y)=0.25D0*( (1-X)*(-2+4*X-2*X**2+9*Y+9*X*Y-6*Y**2)
392
. /(3*Y)-2*(1-X-X**2+X**3-3*Y-2*X*Y-3*X**2*Y+3*Y**2
393
. +3*X*Y**2-Y**3)*(-PI/2- DATAN((1-2*X+X**2-Y-X*Y)/
394
. ((1-X)*DSQRT(-1.D0+2*X+2*Y-(X-Y)**2))))/DSQRT(-1.D0
395
. +2*X-(X-Y)**2+2*Y)-(1+X**2-2*Y-2*X*Y+Y**2)*DLOG(X))
396
QCD0(X) = (1+X**2)*(4*SP((1-X)/(1+X)) + 2*SP((X-1)/(X+1))
397
. - 3*DLOG((1+X)/(1-X))*DLOG(2/(1+X))
398
. - 2*DLOG((1+X)/(1-X))*DLOG(X))
399
. - 3*X*DLOG(4/(1-X**2)) - 4*X*DLOG(X)
400
HQCDM(X)=QCD0(X)/X+(3+34*X**2-13*X**4)/16/X**3*DLOG((1+X)/(1-X))
401
. + 3.D0/8/X**2*(7*X**2-1)
402
AQCDM(X)=QCD0(X)/X + (19+2*X**2+3*X**4)/16/X*DLOG((1+X)/(1-X))
404
HQCD(X)=(4.D0/3*HQCDM(BETA(X))
405
. +2*(4.D0/3-DLOG(X))*(1-10*X)/(1-4*X))*ASH/PI
406
. + (29.14671D0 + RATCOUP*(1.570D0 - 2*DLOG(HIGTOP)/3
407
. + DLOG(X)**2/9))*(ASH/PI)**2
408
. + (164.14D0 - 25.77D0*5 + 0.259D0*5**2)*(ASH/PI)**3
409
AQCD(X)=(4.D0/3*AQCDM(BETA(X))
410
. +2*(4.D0/3-DLOG(X))*(1-6*X)/(1-4*X))*ASH/PI
411
. + (29.14671D0 + RATCOUP*(23/6.D0 - DLOG(HIGTOP)
412
. + DLOG(X)**2/6))*(ASH/PI)**2
413
. + (164.14D0 - 25.77D0*5 + 0.259D0*5**2)*(ASH/PI)**3
416
TQCDH(X)=1.D0+4.D0/3*HQCDM(BETA(X))*ASH/PI
417
TQCDA(X)=1.D0+4.D0/3*AQCDM(BETA(X))*ASH/PI
418
QCDC(X,Y)=1.D0+4/3.D0*ASH/PI*(9/4.D0 + (3-2*X+2*Y)/4*DLOG(X/Y)
419
. +((1.5D0-X-Y)*LAMB(X,Y)**2+5*X*Y)/2/LAMB(X,Y)
420
. /(1-X-Y)*DLOG(XI(X,Y)*XI(Y,X))
422
. + ASH/PI*(2*(4/3.D0-DLOG(X))
423
. - (X*2*(4/3.D0-DLOG(X)) + Y*2*(4/3.D0-DLOG(Y)))/(1-X-Y)
424
. - (X*2*(4/3.D0-DLOG(X))*(1-X+Y)
425
. +Y*2*(4/3.D0-DLOG(Y))*(1+X-Y))/LAMB(X,Y)**2)
426
QCDCI(X,Y)=1.D0+4/3.D0*ASH/PI*(3 + (Y-X)/2*DLOG(X/Y)
427
. +(2*(1-X-Y)+LAMB(X,Y)**2)/2/LAMB(X,Y)
428
. *DLOG(XI(X,Y)*XI(Y,X))
430
. + ASH/PI*(2*(4/3.D0-DLOG(X)) + 2*(4/3.D0-DLOG(Y))
431
. - (X*2*(4/3.D0-DLOG(X))*(1-X+Y)
432
. +Y*2*(4/3.D0-DLOG(Y))*(1+X-Y))/LAMB(X,Y)**2)
433
QCDCM(X,Y)=1.D0+4/3.D0*ASH/PI*(9/4.D0 + (3-2*X+2*Y)/4*DLOG(X/Y)
434
. +((1.5D0-X-Y)*LAMB(X,Y)**2+5*X*Y)/2/LAMB(X,Y)
435
. /(1-X-Y)*DLOG(4*X*Y/(1-X-Y+LAMB(X,Y))**2)
437
QCDCMI(X,Y)=1.D0+4/3.D0*ASH/PI*(3 + (Y-X)/2*DLOG(X/Y)
438
. +(2*(1-X-Y)*LAMB(X,Y)**2)/2/LAMB(X,Y)
439
. *DLOG(4*X*Y/(1-X-Y+LAMB(X,Y))**2)
441
CQCD(Z,TB,X,Y)= GF/(4*PI*DSQRT(2.D0))*Z**3*LAMB(X,Y)
442
. *((1.D0-X-Y)*(X*TB**2*QCDC(X,Y)
443
. +Y/TB**2*QCDC(Y,X))
444
. -4.D0*X*Y*QCDCI(X,Y))
445
CQCDM(Z,TB,X,Y)= GF/(4*PI*DSQRT(2.D0))*Z**3*LAMB(X,Y)
446
. *((1.D0-X-Y)*(X*TB**2*QCDCM(X,Y)
447
. +Y/TB**2*QCDCM(Y,X))
448
. -4.D0*X*Y*QCDCMI(X,Y))
449
ELW(AMH,AMF,QF,ACF)=ALPH/PI*3.D0/2*QF**2
450
. *(3.D0/2-DLOG(AMH**2/AMF**2))
451
. +GF/8/DSQRT(2.D0)/PI**2*(ACF*AMT**2
452
. +AMW**2*(3*DLOG(CS)/SS-5)+AMZ**2*(0.5D0
453
. -3*(1-4*SS*DABS(QF))**2))
454
CF(CA) = -CDLOG(-(1+CDSQRT(1-CA))/(1-CDSQRT(1-CA)))**2/4
455
CG(CA) = CDSQRT(1-CA)/2*CDLOG(-(1+CDSQRT(1-CA))/(1-CDSQRT(1-CA)))
456
CI1(CA,CB) = CA*CB/2/(CA-CB)
457
. + CA**2*CB**2/2/(CA-CB)**2*(CF(CA)-CF(CB))
458
. + CA**2*CB/(CA-CB)**2*(CG(CA)-CG(CB))
459
CI2(CA,CB) = -CA*CB/2/(CA-CB)*(CF(CA)-CF(CB))
460
HGGQCD(ASG,NF)=1.D0+ASG/PI*(95.D0/4.D0-NF*7.D0/6.D0)
461
SGGQCD(ASG)=ASG/PI*17.D0/6.D0
462
AGGQCD(ASG,NF)=1.D0+ASG/PI*(97.D0/4.D0-NF*7.D0/6.D0)
463
HFFSELF(AMH)=1.D0+GF*AMH**2/16.D0/PI**2/DSQRT(2.D0)*2.117203D0
464
. -(GF*AMH**2/16.D0/PI**2/DSQRT(2.D0))**2*32.6567D0
465
HVVSELF(AMH)=1.D0+GF*AMH**2/16.D0/PI**2/DSQRT(2.D0)*2.800952D0
466
. +(GF*AMH**2/16.D0/PI**2/DSQRT(2.D0))**2*62.0308D0
472
C--DECOUPLING THE TOP QUARK FROM ALPHAS
475
C--TOP QUARK DECAY WIDTH
476
GAMT0 = GF*AMT**3/8/DSQRT(2D0)/PI*(1-AMW**2/AMT**2)**2
477
. *(1+2*AMW**2/AMT**2)
478
IF(IHIGGS.NE.0.AND.AMT.GT.AMCH+AMB)THEN
479
GAMT1 = GF*AMT**3/8/DSQRT(2D0)/PI*(1-AMCH**2/AMT**2)**2
480
. *((AMB/AMT)**2*TGBET**2 + 1/TGBET**2)
488
C =========================================================
490
C =========================================================
493
C ============= RUNNING MASSES
499
HIGTOP = AMH**2/AMT**2
509
C =============== PARTIAL WIDTHS
515
CTT = 4*AMT**2/AMH**2*DCMPLX(1D0,-EPS)
516
CTB = 4*AMB**2/AMH**2*DCMPLX(1D0,-EPS)
517
CAT = 2*CTT*(1+(1-CTT)*CF(CTT))
518
CAB = 2*CTB*(1+(1-CTB)*CF(CTB))
519
FQCD=HGGQCD(ASG,NFEXT)
520
XFAC = CDABS(CAT+CAB)**2*FQCD
521
HGG=HVV(AMH,0.D0)*(ASG/PI)**2*XFAC/8
523
C H ---> G G* ---> G CC TO BE ADDED TO H ---> CC
526
FQCD=HGGQCD(ASG,NFEXT)
527
XFAC = CDABS(CAT+CAB)**2*FQCD
528
DCC=HVV(AMH,0.D0)*(ASG/PI)**2*XFAC/8 - HGG
529
C H ---> G G* ---> G BB TO BE ADDED TO H ---> BB
532
FQCD=HGGQCD(ASG,NFEXT)
533
XFAC = CDABS(CAT+CAB)**2*FQCD
534
DBB=HVV(AMH,0.D0)*(ASG/PI)**2*XFAC/8 - HGG - DCC
537
HGG = HGG + DBB + DCC
540
ELSEIF(NFGG.EQ.4)THEN
546
IF(AMH.LE.2*AMMUON) THEN
549
HMM=HFF(AMH,(AMMUON/AMH)**2)
550
. *(1+ELW(AMH,AMMUON,-1.D0,7.D0))
554
IF(AMH.LE.2*AMTAU) THEN
557
HLL=HFF(AMH,(AMTAU/AMH)**2)
558
. *(1+ELW(AMH,AMTAU,-1.D0,7.D0))
562
IF(AMH.LE.2*AMS) THEN
565
HS2=3.D0*HFF(AMH,(RMS/AMH)**2)
566
. *QCDH(RMS**2/AMH**2)
567
. *(1+ELW(AMH,RMS,-1.D0/3.D0,7.D0))
569
IF(HS2.LT.0.D0) HS2 = 0
570
HS1=3.D0*HFF(AMH,(AMS/AMH)**2)
571
. *TQCDH(AMS**2/AMH**2)
574
HSS = QQINT(RAT,HS1,HS2)
577
IF(AMH.LE.2*AMC) THEN
580
HC2=3.D0*HFF(AMH,(RMC/AMH)**2)
581
. *QCDH(RMC**2/AMH**2)
582
. *(1+ELW(AMH,RMC,2.D0/3.D0,7.D0))
585
IF(HC2.LT.0.D0) HC2 = 0
586
HC1=3.D0*HFF(AMH,(AMC/AMH)**2)
587
. *TQCDH(AMC**2/AMH**2)
590
HCC = QQINT(RAT,HC1,HC2)
593
IF(AMH.LE.2*AMB) THEN
596
HB2=3.D0*HFF(AMH,(RMB/AMH)**2)
597
. *QCDH(RMB**2/AMH**2)
598
. *(1+ELW(AMH,RMB,-1.D0/3.D0,1.D0))
601
IF(HB2.LT.0.D0) HB2 = 0
602
HB1=3.D0*HFF(AMH,(AMB/AMH)**2)
603
. *TQCDH(AMB**2/AMH**2)
606
HBB = QQINT(RAT,HB1,HB2)
615
IF (AMH.LE.AMT+AMW+AMB) THEN
617
ELSEIF (AMH.LE.XM1) THEN
618
FACTT=6.D0*GF**2*AMH**3*AMT**2/2.D0/128.D0/PI**3
619
CALL HTOTTS(AMH,AMT,AMB,AMW,HTTS)
621
ELSEIF (AMH.LE.XM2) THEN
626
FACTT=6.D0*GF**2*XX(1)**3*AMT**2/2.D0/128.D0/PI**3
627
CALL HTOTTS(XX(1),AMT,AMB,AMW,HTTS)
629
FACTT=6.D0*GF**2*XX(2)**3*AMT**2/2.D0/128.D0/PI**3
630
CALL HTOTTS(XX(2),AMT,AMB,AMW,HTTS)
633
XY2=3.D0*HFF(XX(3),(XMT/XX(3))**2)
634
. *QCDH(XMT**2/XX(3)**2)
636
IF(XY2.LT.0.D0) XY2 = 0
637
XY1=3.D0*HFF(XX(3),(AMT/XX(3))**2)
638
. *TQCDH(AMT**2/XX(3)**2)
641
YY(3) = QQINT(RAT,XY1,XY2)
643
XY2=3.D0*HFF(XX(4),(XMT/XX(4))**2)
644
. *QCDH(XMT**2/XX(4)**2)
646
IF(XY2.LT.0.D0) XY2 = 0
647
XY1=3.D0*HFF(XX(4),(AMT/XX(4))**2)
648
. *TQCDH(AMT**2/XX(4)**2)
651
YY(4) = QQINT(RAT,XY1,XY2)
652
HTT = FINT(AMH,XX,YY)
654
HT2=3.D0*HFF(AMH,(RMT/AMH)**2)
655
. *QCDH(RMT**2/AMH**2)
657
IF(HT2.LT.0.D0) HT2 = 0
658
HT1=3.D0*HFF(AMH,(AMT/AMH)**2)
659
. *TQCDH(AMT**2/AMH**2)
662
HTT = QQINT(RAT,HT1,HT2)
665
IF (AMH.LE.2.D0*AMT) THEN
668
HT2=3.D0*HFF(AMH,(RMT/AMH)**2)
669
. *QCDH(RMT**2/AMH**2)
671
IF(HT2.LT.0.D0) HT2 = 0
672
HT1=3.D0*HFF(AMH,(AMT/AMH)**2)
673
. *TQCDH(AMT**2/AMH**2)
676
HTT = QQINT(RAT,HT1,HT2)
681
CTT = 4*AMT**2/AMH**2*DCMPLX(1D0,-EPS)
682
CTB = 4*AMB**2/AMH**2*DCMPLX(1D0,-EPS)
683
CTC = 4*AMC**2/AMH**2*DCMPLX(1D0,-EPS)
684
CTL = 4*AMTAU**2/AMH**2*DCMPLX(1D0,-EPS)
685
CTW = 4*AMW**2/AMH**2*DCMPLX(1D0,-EPS)
686
CAW = -(2+3*CTW+3*CTW*(2-CTW)*CF(CTW))
687
CAT = 4/3D0 * 2*CTT*(1+(1-CTT)*CF(CTT))
688
CAB = 1/3D0 * 2*CTB*(1+(1-CTB)*CF(CTB))
689
CAC = 4/3D0 * 2*CTC*(1+(1-CTC)*CF(CTC))
690
CAL = 2*CTL*(1+(1-CTL)*CF(CTL))
691
XFAC = CDABS(CAT+CAB+CAC+CAL+CAW)**2
692
HGA=HVV(AMH,0.D0)*(ALPH/PI)**2/16.D0*XFAC
699
FT = -3*2D0/3*(1-4*2D0/3*SS)/DSQRT(SS*CS)
700
FB = 3*1D0/3*(-1+4*1D0/3*SS)/DSQRT(SS*CS)
701
CTT = 4*AMT**2/AMH**2*DCMPLX(1D0,-EPS)
702
CTB = 4*AMB**2/AMH**2*DCMPLX(1D0,-EPS)
703
CTW = 4*AMW**2/AMH**2*DCMPLX(1D0,-EPS)
704
CLT = 4*AMT**2/AMZ**2*DCMPLX(1D0,-EPS)
705
CLB = 4*AMB**2/AMZ**2*DCMPLX(1D0,-EPS)
706
CLW = 4*AMW**2/AMZ**2*DCMPLX(1D0,-EPS)
707
CAT = FT*(CI1(CTT,CLT) - CI2(CTT,CLT))
708
CAB = FB*(CI1(CTB,CLB) - CI2(CTB,CLB))
709
CAW = -1/DSQRT(TS)*(4*(3-TS)*CI2(CTW,CLW)
710
. + ((1+2/CTW)*TS - (5+2/CTW))*CI1(CTW,CLW))
711
XFAC = CDABS(CAT+CAB+CAW)**2
712
ACOUP = DSQRT(2D0)*GF*AMZ**2*SS*CS/PI**2
713
HZGA = GF/(4.D0*PI*DSQRT(2.D0))*AMH**3*(ALPH/PI)*ACOUP/16.D0
714
. *XFAC*(1-AMZ**2/AMH**2)**3
723
CALL HTOVV(AMH,AMW,GAMW,HTWW)
724
HWW = 3D0/2D0*GF*AMW**4/DSQRT(2D0)/PI/AMH**3*HTWW
725
ELSEIF (AMH.LE.XM2) THEN
730
CALL HTOVV(XX(1),AMW,GAMW,HTWW)
731
YY(1)=3D0/2D0*GF*AMW**4/DSQRT(2D0)/PI/XX(1)**3*HTWW
732
CALL HTOVV(XX(2),AMW,GAMW,HTWW)
733
YY(2)=3D0/2D0*GF*AMW**4/DSQRT(2D0)/PI/XX(2)**3*HTWW
734
YY(3)=HVV(XX(3),AMW**2/XX(3)**2)
736
YY(4)=HVV(XX(4),AMW**2/XX(4)**2)
738
HWW = FINT(AMH,XX,YY)
740
HWW=HVV(AMH,AMW**2/AMH**2)
750
ELSE IF (AMH.LE.XM1) THEN
751
CWW=3.D0*GF**2*AMW**4/16.D0/PI**3
752
HWW=HV(AMW**2/AMH**2)*CWW*AMH
753
ELSE IF (AMH.LT.XM2) THEN
754
CWW=3.D0*GF**2*AMW**4/16.D0/PI**3
759
YY(1)=HV(AMW**2/XX(1)**2)*CWW*XX(1)
760
YY(2)=HV(AMW**2/XX(2)**2)*CWW*XX(2)
761
YY(3)=HVV(XX(3),AMW**2/XX(3)**2)
763
YY(4)=HVV(XX(4),AMW**2/XX(4)**2)
765
HWW = FINT(AMH,XX,YY)
767
HWW=HVV(AMH,AMW**2/AMH**2)
778
CALL HTOVV(AMH,AMZ,GAMZ,HTZZ)
779
HZZ = 3D0/4D0*GF*AMZ**4/DSQRT(2D0)/PI/AMH**3*HTZZ
780
ELSEIF (AMH.LE.XM2) THEN
785
CALL HTOVV(XX(1),AMZ,GAMZ,HTZZ)
786
YY(1)=3D0/4D0*GF*AMZ**4/DSQRT(2D0)/PI/XX(1)**3*HTZZ
787
CALL HTOVV(XX(2),AMZ,GAMZ,HTZZ)
788
YY(2)=3D0/4D0*GF*AMZ**4/DSQRT(2D0)/PI/XX(2)**3*HTZZ
789
YY(3)=HVV(XX(3),AMZ**2/XX(3)**2)/2
791
YY(4)=HVV(XX(4),AMZ**2/XX(4)**2)/2
793
HZZ = FINT(AMH,XX,YY)
795
HZZ=HVV(AMH,AMZ**2/AMH**2)/2.D0
805
ELSE IF (AMH.LE.XM1) THEN
806
CZZ=3.D0*GF**2*AMZ**4/192.D0/PI**3*(7-40/3.D0*SS+160/9.D0*SS**2)
807
HZZ=HV(AMZ**2/AMH**2)*CZZ*AMH
808
ELSE IF (AMH.LT.XM2) THEN
809
CZZ=3.D0*GF**2*AMZ**4/192.D0/PI**3*(7-40/3.D0*SS+160/9.D0*SS**2)
814
YY(1)=HV(AMZ**2/XX(1)**2)*CZZ*XX(1)
815
YY(2)=HV(AMZ**2/XX(2)**2)*CZZ*XX(2)
816
YY(3)=HVV(XX(3),AMZ**2/XX(3)**2)/2
818
YY(4)=HVV(XX(4),AMZ**2/XX(4)**2)/2
820
HZZ = FINT(AMH,XX,YY)
822
HZZ=HVV(AMH,AMZ**2/AMH**2)/2.D0
827
C ========== TOTAL WIDTH AND BRANCHING RATIOS
829
WTOT=HLL+HMM+HSS+HCC+HBB+HTT+HGG+HGA+HZGA+HWW+HZZ
850
DOUBLE PRECISION FUNCTION BIJ(X,Y)
851
IMPLICIT DOUBLE PRECISION (A-H,O-Z)
852
DOUBLE PRECISION LAMB
853
BIJ = (1-X-Y)/LAMB(X,Y)*(
854
. 4*SP(XI(X,Y)*XI(Y,X))
855
. - 2*SP(-XI(X,Y)) - 2*SP(-XI(Y,X))
856
. + 2*DLOG(XI(X,Y)*XI(Y,X))*DLOG(1-XI(X,Y)*XI(Y,X))
857
. - DLOG(XI(X,Y))*DLOG(1+XI(X,Y))
858
. - DLOG(XI(Y,X))*DLOG(1+XI(Y,X))
860
. -4*(DLOG(1-XI(X,Y)*XI(Y,X))
861
. +XI(X,Y)*XI(Y,X)/(1-XI(X,Y)*XI(Y,X))*DLOG(XI(X,Y)*XI(Y,X)))
862
. + (LAMB(X,Y)+X-Y)/LAMB(X,Y)*(DLOG(1+XI(X,Y))
863
. - XI(X,Y)/(1+XI(X,Y))*DLOG(XI(X,Y)))
864
. + (LAMB(X,Y)-X+Y)/LAMB(X,Y)*(DLOG(1+XI(Y,X))
865
. - XI(Y,X)/(1+XI(Y,X))*DLOG(XI(Y,X)))
869
DOUBLE PRECISION FUNCTION BETA(X)
870
IMPLICIT DOUBLE PRECISION (A-H,O-Z)
871
BETA=DSQRT(1.D0-4.D0*X)
875
DOUBLE PRECISION FUNCTION LAMB(X,Y)
876
IMPLICIT DOUBLE PRECISION (A-H,O-Z)
877
LAMB=DSQRT((1.D0-X-Y)**2-4.D0*X*Y)
881
DOUBLE PRECISION FUNCTION XI(X,Y)
882
IMPLICIT DOUBLE PRECISION (A-H,O-Z)
883
DOUBLE PRECISION LAMB
884
XI = 2*X/(1-X-Y+LAMB(X,Y))
889
C ===================== THE FUNCTION F0 ===============
890
COMPLEX*16 FUNCTION F0(M1,M2,QSQ)
891
IMPLICIT REAL*8 (A-H,M,O-Z)
892
COMPLEX*16 CD,CR,CQ2,IEPS,CBET,CXX
896
IEPS = DCMPLX(1.D0,1.D-12)
899
CR = CDSQRT((1+CD)**2 - 4*M1SQ/CQ2)
904
F0 = -2.D0 + CR*CDLOG(-(1+CR)/(1-CR))
906
CBET = CDSQRT(1-4*M1*M2/(CQ2 - (M1-M2)**2))
907
CXX = (CBET-1)/(CBET+1)
908
F0 = -1 + ((QSQ+M2SQ-M1SQ)/2/QSQ - M2SQ/(M2SQ-M1SQ))
910
. - (QSQ-(M1-M2)**2)/QSQ*CBET*CDLOG(CXX)
916
C ************************************************************
917
C SUBROUTINE FOR HSM ---> V*V* ---> 4F
918
C ************************************************************
919
SUBROUTINE HTOVV(AMH,AMV,GAMV,HTVV)
920
IMPLICIT DOUBLE PRECISION (A-H,O-Z)
921
COMMON/VVOFF/AMH1,AMV1,GAMV1
933
CALL QGAUS1(FTOVV1,DD,UU,RES)
940
DOUBLE PRECISION FUNCTION FTOVV1(XX)
941
IMPLICIT DOUBLE PRECISION (A-H,O-Z)
951
CALL QGAUS2(FTOVV2,DD,UU,RES)
958
DOUBLE PRECISION FUNCTION FTOVV2(XX)
959
IMPLICIT DOUBLE PRECISION (A-H,O-Z)
968
DOUBLE PRECISION FUNCTION FTOVV(XX)
969
IMPLICIT DOUBLE PRECISION (A-H,O-Z)
970
DOUBLE PRECISION LAMB
972
COMMON/VVOFF/AMH,AMW,GAMW
973
LAMB(X,Y)=DSQRT((1.D0-X-Y)**2-4.D0*X*Y)
978
Y1 = DATAN((YY-AMW**2)/AMW/GAMW)
979
Y2 = -DATAN((AMW**2)/AMW/GAMW)
981
T1 = TAN(Y1*XX(1)+Y2*(1.D0-XX(1)))
982
SP = AMW**2 + AMW*GAMW*T1
983
YY = (AMH-DSQRT(SP))**2
984
Y1 = DATAN((YY-AMW**2)/AMW/GAMW)
985
Y2 = -DATAN((AMW**2)/AMW/GAMW)
987
T2 = TAN(Y1*XX(2)+Y2*(1.D0-XX(2)))
988
SM = AMW**2 + AMW*GAMW*T2
990
GAM = AM2*LAMB(SP/AM2,SM/AM2)*(1+LAMB(SP/AM2,SM/AM2)**2*AMH**4
994
FTOVV = PRO1*PRO2*GAM*DJAC/PI**2
997
SM = (AMH-DSQRT(SP))**2*XX(2)
998
DJAC = AMH**2*(AMH-DSQRT(SP))**2/PI**2
1000
GAM = AM2*LAMB(SP/AM2,SM/AM2)*(1+LAMB(SP/AM2,SM/AM2)**2*AMH**4
1002
PRO1 = SP*GAMW/AMW/((SP-AMW**2)**2+AMW**2*GAMW**2)
1003
PRO2 = SM*GAMW/AMW/((SM-AMW**2)**2+AMW**2*GAMW**2)
1004
FTOVV = PRO1*PRO2*GAM*DJAC
1009
C ************************************************************
1010
C SUBROUTINE FOR HSM ---> TT* ---> TBW
1011
C ************************************************************
1012
SUBROUTINE HTOTTS(AMH,AMT,AMB,AMW,HTTS)
1013
IMPLICIT REAL*8(A-Z)
1017
COMMON/IKSY0/X1,X2,M1,M2,M3,ECM,S
1018
COMMON/TOP0/AMH0,AMT0,AMB0,AMW0
1027
C FIRST INTEGRATE OVER X2, i.e. (1+3) SYSTEM
1028
C CHECK WHETHER ENOUGH PHASE SPACE
1030
IF(MASTOT.GE.AMH) GOTO 12
1041
CALL QGAUS1(FUNSTT1,D,U,SS)
1051
DOUBLE PRECISION FUNCTION FUNSTT1(XL)
1052
IMPLICIT REAL*8(A-Z)
1054
COMMON/IKSY0/X1,X2,M1,M2,M3,ECM,S
1060
E2S=(S-S13-M2**2)/TEM
1061
E3S=(S13+M3**2-M1**2)/TEM
1062
C SECOND INTEGRAL OVER X1, i.e. (2+3) SYSTEM
1063
U1=(E2S+E3S)**2-(DSQRT(E2S**2-M2**2)-DSQRT(E3S**2-M3**2))**2
1064
D1=(E2S+E3S)**2-(DSQRT(E2S**2-M2**2)+DSQRT(E3S**2-M3**2))**2
1071
CALL QGAUS2(FUNSTT2,D,U,SS)
1079
DOUBLE PRECISION FUNCTION FUNSTT2(XK)
1080
IMPLICIT REAL*8(A-Z)
1081
COMMON/IKSY0/X1,X2,M1,M2,M3,ECM,S
1088
SUBROUTINE ELEMSTT(RES)
1089
IMPLICIT REAL*8(A-Z)
1090
COMMON/IKSY0/X1,X2,M1,M2,M3,ECM,S
1091
COMMON/TOP0/AMH,AMT,AMB,AMW
1092
COMMON/WZWDTH/GAMC0,GAMT0,GAMT1,GAMW0,GAMZ0
1093
GAMT=GAMT0**2*AMT**2/AMH**4
1094
GAMW=GAMW0**2*AMW**2/AMH**4
1102
W11=1.D0/((1.D0-X2)**2+GAMT)
1103
W33=1.D0/(W3**2+GAMW**2)
1106
R11=4*T*W-16.*T*W*Y1-4.*T*Y2*Y1+8.*T*Y1+32.*T*W**2-20
1107
. .*T*Y1**2+8.*W*Y2*Y1+4.*W*Y1**2-4.*Y2*Y1**2-16.*T**2*W-
1108
. 32.*T**2*Y1+4.*T**2-16.*T**3-8.*W**2+4.*Y1**2-4.*Y1**3
1109
R33=-4.*T*W+4.*T*W*Y2-2.*T*W*Y2*Y1+4.*T*W*Y1+T*W*Y2**2-
1110
. 3.*T*W*Y1**2+2.*T*Y2*Y1-3.*T*Y2*Y1**2+4.*T*W**2-4.*T*W**3
1111
. +T*Y2**2-3.*T*Y2**2*Y1-T*Y2**3+T*Y1**2-T*Y1**3+4.*T**2
1112
. *W-4.*T**2*W*Y2-4.*T**2*W*Y1-2.*T**2*Y2*Y1-4.*T**2*W**2-
1113
. T**2*Y2**2-T**2*Y1**2+4.*W**2*Y2*Y1-8.*W**3*Y2-8.*W**3*Y1
1115
R13=8.*W-24.*T*W+16.*T*W*Y1 -4.*T*Y2+16.*T*Y2*Y1-4.*T*
1116
. Y1+16.*T*W**2+4.*T*Y2**2+12.*T*Y1**2-8.*W*Y2-12.*W*Y2*Y1
1117
. -8.*W*Y1+4.*W*Y1**2-4.*Y2*Y1+8.*Y2*Y1**2+16.*T**2*W+8.
1118
. *T**2*Y2+8.*T**2*Y1+16.*W**2*Y2+24.*W**2*Y1+4.*Y2**2*Y1-
1119
. 32.*W**3-4.*Y1**2+4.*Y1**3
1120
RES=R11*W11+4.D0*R33*W33/T-2.D0*R13*W13
1124
C **************************************************
1125
C SUBROUTINE FOR A -> TT* -> TBW
1126
C **************************************************
1128
SUBROUTINE ATOTT(AMA,AMT,AMB,AMW,AMCH,ATT0)
1129
IMPLICIT REAL*8(A-Z)
1133
COMMON/IKSY1/X1,X2,M1,M2,M3,ECM,S
1134
COMMON/TOP1/AMA1,AMT1,AMB1,AMW1,AMCH1
1144
C FIRST INTEGRATE OVER X2, i.e. (1+3) SYSTEM
1145
C CHECK WHETHER ENOUGH PHASE SPACE
1147
IF(MASTOT.GE.AMA) GOTO 12
1158
CALL QGAUS1(FUNATT1,D,U,SS)
1168
DOUBLE PRECISION FUNCTION FUNATT1(XL)
1169
IMPLICIT REAL*8(A-Z)
1171
COMMON/IKSY1/X1,X2,M1,M2,M3,ECM,S
1177
E2S=(S-S13-M2**2)/TEM
1178
E3S=(S13+M3**2-M1**2)/TEM
1179
C SECOND INTEGRAL OVER X1, i.e. (2+3) SYSTEM
1180
U1=(E2S+E3S)**2-(DSQRT(E2S**2-M2**2)-DSQRT(E3S**2-M3**2))**2
1181
D1=(E2S+E3S)**2-(DSQRT(E2S**2-M2**2)+DSQRT(E3S**2-M3**2))**2
1188
CALL QGAUS2(FUNATT2,D,U,SS)
1196
DOUBLE PRECISION FUNCTION FUNATT2(XK)
1197
IMPLICIT REAL*8(A-Z)
1198
COMMON/IKSY1/X1,X2,M1,M2,M3,ECM,S
1205
SUBROUTINE ELEMATT(RES)
1206
IMPLICIT REAL*8(A-Z)
1207
COMMON/IKSY1/X1,X2,M1,M2,M3,ECM,S
1208
COMMON/TOP1/AMA,AMT,AMB,AMW,AMCH
1209
COMMON/WZWDTH/GAMC0,GAMT0,GAMT1,GAMW,GAMZ
1210
GAMT=GAMT1**2*AMT**2/AMA**4
1211
GAMC=GAMC0**2*AMCH**2/AMA**4
1220
W22=1.D0/ ((1.D0-X0+W-CH)**2+GAMC)
1221
W11=1.D0/((1.D0-X2)**2+GAMT)
1223
R11=4.D0*T*W-4.D0*T*Y1*Y2+8.D0*T*Y2-4.D0*T*Y2**2+8.D0*W*Y1*Y2+4.D0
1224
. *W*Y2**2-4.D0*Y1*Y2**2+4.D0*T**2-8.D0*W**2+4.D0*Y2**2-4.D0*Y2**3
1225
R22=-16.D0*W+16.D0*T*W-8.D0*T*Y1*Y2-4.D0*T*Y1**2-4.D0*T*Y2**2+16.
1226
.D0*W*Y1+8.D0*W*Y1*Y2+16.D0*W*Y2+4.D0*W*Y1**2+4.D0*W*Y2**2+8.D0*Y1*
1227
. Y2-12.D0*Y1*Y2**2-12.D0*Y1**2*Y2-16.D0*W**2+4.D0*Y1**2-4.D0*Y1**3
1228
. +4.D0*Y2**2-4.D0*Y2**3
1229
R12=16.D0*W-16.D0*T*W-8.D0*T*Y1+16.D0*T*Y1*Y2-8.D0*T*Y2+8.D0*T*Y1
1230
. **2+8.D0*T*Y2**2-16.D0*W*Y1-8.D0*W*Y1*Y2-16.D0*W*Y2-8.D0*W*Y2**2-
1231
. 8.D0*Y1*Y2+16.D0*Y1*Y2**2+8.D0*Y1**2*Y2+16.D0*W**2-8.D0*Y2**2
1233
RES=R11*W11+R22*W22+R12*W12
1237
C ************************************************************
1238
C SUBROUTINE FOR H ---> TT* ---> TBW
1239
C ************************************************************
1240
SUBROUTINE HTOTT(AMH,AMT,AMB,AMW,AMCH,TB,GHT,GAT,GHVV,HTT0)
1241
IMPLICIT REAL*8(A-Z)
1245
COMMON/IKSY2/X1,X2,M1,M2,M3,ECM,S
1246
COMMON/TOP2/AMH2,AMT2,AMB2,AMW2,AMCH2,TB2,GHT2,GAT2,GHVV2
1260
C FIRST INTEGRATE OVER X2, i.e. (1+3) SYSTEM
1261
C CHECK WHETHER ENOUGH PHASE SPACE
1263
IF(MASTOT.GE.AMH) GOTO 12
1274
CALL QGAUS1(FUNHTT1,D,U,SS)
1284
DOUBLE PRECISION FUNCTION FUNHTT1(XL)
1285
IMPLICIT REAL*8(A-Z)
1287
COMMON/IKSY2/X1,X2,M1,M2,M3,ECM,S
1293
E2S=(S-S13-M2**2)/TEM
1294
E3S=(S13+M3**2-M1**2)/TEM
1295
C SECOND INTEGRAL OVER X1, i.e. (2+3) SYSTEM
1296
U1=(E2S+E3S)**2-(DSQRT(E2S**2-M2**2)-DSQRT(E3S**2-M3**2))**2
1297
D1=(E2S+E3S)**2-(DSQRT(E2S**2-M2**2)+DSQRT(E3S**2-M3**2))**2
1304
CALL QGAUS2(FUNHTT2,D,U,SS)
1312
DOUBLE PRECISION FUNCTION FUNHTT2(XK)
1313
IMPLICIT REAL*8(A-Z)
1314
COMMON/IKSY2/X1,X2,M1,M2,M3,ECM,S
1321
SUBROUTINE ELEMHTT(RES)
1322
IMPLICIT REAL*8(A-Z)
1323
COMMON/IKSY2/X1,X2,M1,M2,M3,ECM,S
1324
COMMON/TOP2/AMH,AMT,AMB,AMW,AMCH,TB,GHT,GAT,GHVV
1325
COMMON/WZWDTH/GAMC0,GAMT0,GAMT1,GAMW0,GAMZ0
1326
GAMT=GAMT1**2*AMT**2/AMH**4
1327
GAMC=GAMC0**2*AMCH**2/AMH**4
1328
GAMW=GAMW0**2*AMW**2/AMH**4
1338
W22=1.D0/ ((1.D0-X0+W-CH)**2+GAMC)
1339
W11=1.D0/((1.D0-X2)**2+GAMT)
1340
W33=1.D0/(W3**2+GAMW**2)
1345
R11=4*T*W-16.*T*W*Y1-4.*T*Y2*Y1+8.*T*Y1+32.*T*W**2-20
1346
. .*T*Y1**2+8.*W*Y2*Y1+4.*W*Y1**2-4.*Y2*Y1**2-16.*T**2*W-
1347
. 32.*T**2*Y1+4.*T**2-16.*T**3-8.*W**2+4.*Y1**2-4.*Y1**3
1348
R22=-16.*W+16.*T*W-8.*T*Y2*Y1-4.*T*Y2**2-4.*T*Y1**2+16
1349
. .*W*Y2 + 8.*W*Y2*Y1 + 16.*W*Y1 + 4.*W*Y2**2 + 4.*W*Y1**2+8.*Y2*
1350
. Y1-12.*Y2*Y1**2-12.*Y2**2*Y1-16.*W**2+4.*Y2**2-4.*Y2**3
1351
. +4.*Y1**2-4.*Y1**3
1352
R33=-4.*T*W+4.*T*W*Y2-2.*T*W*Y2*Y1+4.*T*W*Y1+T*W*Y2**2-
1353
. 3.*T*W*Y1**2+2.*T*Y2*Y1-3.*T*Y2*Y1**2+4.*T*W**2-4.*T*W**3
1354
. +T*Y2**2-3.*T*Y2**2*Y1-T*Y2**3+T*Y1**2-T*Y1**3+4.*T**2
1355
. *W-4.*T**2*W*Y2-4.*T**2*W*Y1-2.*T**2*Y2*Y1-4.*T**2*W**2-
1356
. T**2*Y2**2-T**2*Y1**2+4.*W**2*Y2*Y1-8.*W**3*Y2-8.*W**3*Y1
1358
R12=-16.*W+48.*T*W-16.*T*W*Y2+16.*T*W*Y1+8.*T*Y2-32.*T
1359
. *Y2*Y1+8.*T*Y1-8.*T*Y2**2 - 24.*T*Y1**2+16.*W*Y2+8.*W*Y2*
1360
. Y1+16.*W*Y1+8.*W*Y1**2+8.*Y2*Y1-16.*Y2*Y1**2-16.*T**2*Y2
1361
. -16.*T**2*Y1-8.*Y2**2*Y1-16.*W**2+8.*Y1**2-8.*Y1**3
1362
R13=8.*W-24.*T*W+16.*T*W*Y1 -4.*T*Y2+16.*T*Y2*Y1-4.*T*
1363
. Y1+16.*T*W**2+4.*T*Y2**2+12.*T*Y1**2-8.*W*Y2-12.*W*Y2*Y1
1364
. -8.*W*Y1+4.*W*Y1**2-4.*Y2*Y1+8.*Y2*Y1**2+16.*T**2*W+8.
1365
. *T**2*Y2+8.*T**2*Y1+16.*W**2*Y2+24.*W**2*Y1+4.*Y2**2*Y1-
1366
. 32.*W**3-4.*Y1**2+4.*Y1**3
1367
R23=16.*W-16.*T*W+8.*T*W*Y2+8.*T*W*Y1+8.*T*Y2*Y1+4.*T*
1368
. Y2**2+4.*T*Y1**2-16.*W*Y2-16.*W*Y1-4.*W*Y2**2+4.*W*Y1**2
1369
. -8.*Y2*Y1+12.*Y2*Y1**2+8.*W**2*Y2-8.*W**2*Y1+12.*Y2**2*
1370
. Y1-4.*Y2**2+4.*Y2**3-4.*Y1**2+4.*Y1**3
1371
GLVV=DSQRT(1.D0-GHVV**2)
1372
RES=GHT**2*R11*W11+GLVV**2*GAT**2*R22*W22+
1373
. 4.D0*GHVV**2*R33*W33/T+2.D0*GHT*GLVV*GAT*R12*W12+
1374
. 2.D0*GHT*GHVV*R13*W13+2.D0*GHVV*GLVV*GAT*R23*W23
1378
C ************************************************************
1379
C SUBROUTINE FOR H+ ---> BT* ---> BBW
1380
C ************************************************************
1381
SUBROUTINE CTOTT(AMCH,AMT,AMB,AMW,CTT0)
1382
IMPLICIT REAL*8(A-Z)
1386
COMMON/IKSY3/X1,X2,M1,M2,M3,ECM,S
1387
COMMON/TOP3/AMH3,AMT3,AMB3,AMW3
1396
C FIRST INTEGRATE OVER X2, i.e. (1+3) SYSTEM
1397
C CHECK WHETHER ENOUGH PHASE SPACE
1399
IF(MASTOT.GE.AMCH) GOTO 12
1410
CALL QGAUS1(FUNCTT1,D,U,SS)
1420
DOUBLE PRECISION FUNCTION FUNCTT1(XL)
1421
IMPLICIT REAL*8(A-Z)
1423
COMMON/IKSY3/X1,X2,M1,M2,M3,ECM,S
1429
E2S=(S-S13-M2**2)/TEM
1430
E3S=(S13+M3**2-M1**2)/TEM
1431
C SECOND INTEGRAL OVER X1, i.e. (2+3) SYSTEM
1432
U1=(E2S+E3S)**2-(DSQRT(E2S**2-M2**2)-DSQRT(E3S**2-M3**2))**2
1433
D1=(E2S+E3S)**2-(DSQRT(E2S**2-M2**2)+DSQRT(E3S**2-M3**2))**2
1440
CALL QGAUS2(FUNCTT2,D,U,SS)
1448
DOUBLE PRECISION FUNCTION FUNCTT2(XK)
1449
IMPLICIT REAL*8(A-Z)
1450
COMMON/IKSY3/X1,X2,M1,M2,M3,ECM,S
1457
SUBROUTINE ELEMCTT(RES)
1458
IMPLICIT REAL*8(A-Z)
1459
COMMON/IKSY3/X1,X2,M1,M2,M3,ECM,S
1460
COMMON/TOP3/AMCH,AMT,AMB,AMW
1461
COMMON/WZWDTH/GAMC0,GAMT0,GAMT1,GAMW,GAMZ
1462
GAMT=GAMT1**2*AMT**2/AMCH**4
1466
RES=((1.D0-X1-W)*(1.D0-X2-W)+W*(X1+X2-1.D0+W))/
1467
. ((1.D0-X2+B-T)**2+GAMT)
1471
C ***************** INTEGRATION ROUTINE ***********************
1472
C Returns SS as integral of FUNC from A to B, by 10-point Gauss-
1473
C Legendre integration
1474
SUBROUTINE QGAUS1(FUNC,A,B,SS)
1475
IMPLICIT REAL*8(A-Z)
1479
DATA X/.1488743389D0,.4333953941D0,.6794095682D0
1480
. ,.8650633666D0,.9739065285D0/
1481
DATA W/.2955242247D0,.2692667193D0,.2190863625D0
1482
. ,.1494513491D0,.0666713443D0/
1488
SS=SS+W(J)*(FUNC(XM+DX)+FUNC(XM-DX))
1494
C Returns SS as integral of FUNC from A to B, by 10-point Gauss-
1495
C Legendre integration
1496
SUBROUTINE QGAUS2(FUNC,A,B,SS)
1497
IMPLICIT REAL*8(A-Z)
1501
DATA X/.1488743389D0,.4333953941D0,.6794095682D0
1502
. ,.8650633666D0,.9739065285D0/
1503
DATA W/.2955242247D0,.2692667193D0,.2190863625D0
1504
. ,.1494513491D0,.0666713443D0/
1510
SS=SS+W(J)*(FUNC(XM+DX)+FUNC(XM-DX))
1517
C ******************************************************************
1519
C DOUBLE PRECISION FUNCTION RUNP(Q,NF)
1520
C IMPLICIT DOUBLE PRECISION (A-H,O-Z)
1521
C COMMON/RUN/XMSB,XMHAT,XKFAC
1523
C RUNP = RUNM(Q/2.D0,NF)*XKFAC
1527
DOUBLE PRECISION FUNCTION RUNM(Q,NF)
1528
IMPLICIT DOUBLE PRECISION (A-H,O-Z)
1530
PARAMETER (ZETA3 = 1.202056903159594D0)
1531
DIMENSION AM(NN),YMSB(NN)
1532
COMMON/ALS/XLAMBDA,AMCA,AMBA,AMTA,N0A
1533
COMMON/MASSES/AMS,AMC,AMB,AMT
1535
COMMON/RUN/XMSB,XMHAT,XKFAC
1536
COMMON/FLAG/IHIGGS,NNLO,IPOLE
1538
B0(NF)=(33.D0-2.D0*NF)/12D0
1539
B1(NF) = (102D0-38D0/3D0*NF)/16D0
1540
B2(NF) = (2857D0/2D0-5033D0/18D0*NF+325D0/54D0*NF**2)/64D0
1542
G1(NF) = (202D0/3D0-20D0/9D0*NF)/16D0
1543
G2(NF) = (1249D0-(2216D0/27D0+160D0/3D0*ZETA3)*NF
1544
. - 140D0/81D0*NF**2)/64D0
1545
C1(NF) = G1(NF)/B0(NF) - B1(NF)*G0(NF)/B0(NF)**2
1546
C2(NF) = ((G1(NF)/B0(NF) - B1(NF)*G0(NF)/B0(NF)**2)**2
1547
. + G2(NF)/B0(NF) + B1(NF)**2*G0(NF)/B0(NF)**3
1548
. - B1(NF)*G1(NF)/B0(NF)**2 - B2(NF)*G0(NF)/B0(NF)**2)/2D0
1549
TRAN(X,XK)=1D0+4D0/3D0*ALPHAS_HD(X,2)/PI+XK*(ALPHAS_HD(X,2)/PI)**2
1550
CQ(X,NF)=(2D0*B0(NF)*X)**(G0(NF)/B0(NF))
1551
. *(1D0+C1(NF)*X+C2(NF)*X**2)
1557
C--------------------------------------------
1560
IF(ISTRANGE.EQ.0)THEN
1561
C--STRANGE POLE MASS FROM MSBAR-MASS AT 1 GEV
1564
123 AMS = (AMSU+AMSD)/2
1566
XMSB = AMS/CQ(ALPHAS_HD(AMS,2)/PI,3)
1567
. *CQ(ALPHAS_HD(1.D0,2)/PI,3)/TRAN(AMS,0D0)
1568
DD = (XMSB-AMSB)/AMSB
1569
IF(DABS(DD).GE.ACC)THEN
1584
C--------------------------------------------
1591
XK = XK - 1.04D0*(1.D0-AM(I)/AM(NF))
1594
XMSB = AM(NF)/TRAN(AM(NF),0D0)
1595
XMHAT = XMSB/CQ(ALPHAS_HD(AM(NF),2)/PI,NF)
1602
YMSB(4) = YMSB(3)*CQ(ALPHAS_HD(AM(4),2)/PI,3)/
1603
. CQ(ALPHAS_HD(1.D0,2)/PI,3)
1604
YMSB(5) = YMSB(4)*CQ(ALPHAS_HD(AM(5),2)/PI,4)/
1605
. CQ(ALPHAS_HD(AM(4),2)/PI,4)
1606
YMSB(6) = YMSB(5)*CQ(ALPHAS_HD(AM(6),2)/PI,5)/
1607
. CQ(ALPHAS_HD(AM(5),2)/PI,5)
1610
YMSB(5) = YMSB(4)*CQ(ALPHAS_HD(AM(5),2)/PI,4)/
1611
. CQ(ALPHAS_HD(AM(4),2)/PI,4)
1612
YMSB(6) = YMSB(5)*CQ(ALPHAS_HD(AM(6),2)/PI,5)/
1613
. CQ(ALPHAS_HD(AM(5),2)/PI,5)
1616
YMSB(4) = YMSB(5)*CQ(ALPHAS_HD(AM(4),2)/PI,4)/
1617
. CQ(ALPHAS_HD(AM(5),2)/PI,4)
1618
YMSB(6) = YMSB(5)*CQ(ALPHAS_HD(AM(6),2)/PI,5)/
1619
. CQ(ALPHAS_HD(AM(5),2)/PI,5)
1622
YMSB(5) = YMSB(6)*CQ(ALPHAS_HD(AM(5),2)/PI,5)/
1623
. CQ(ALPHAS_HD(AM(6),2)/PI,5)
1624
YMSB(4) = YMSB(5)*CQ(ALPHAS_HD(AM(4),2)/PI,4)/
1625
. CQ(ALPHAS_HD(AM(5),2)/PI,4)
1630
ELSEIF(Q.LE.AMB)THEN
1633
ELSEIF(Q.LE.AMT)THEN
1640
IF(NNLO.EQ.1.AND.NF.GT.3)THEN
1641
XKFAC = TRAN(AM(NF),0D0)/TRAN(AM(NF),XK)
1645
RUNM = YMSB(N0)*CQ(ALPHAS_HD(Q,2)/PI,N0)/
1646
. CQ(ALPHAS_HD(Q0,2)/PI,N0)
1651
DOUBLE PRECISION FUNCTION ALPHAS_HD(Q,N)
1652
IMPLICIT DOUBLE PRECISION (A-H,O-Z)
1654
COMMON/ALSLAM/XLB1(6),XLB2(6)
1655
COMMON/ALS/XLAMBDA,AMC,AMB,AMT,N0
1656
B0(NF)=33.D0-2.D0*NF
1657
B1(NF)=6.D0*(153.D0-19.D0*NF)/B0(NF)**2
1658
ALS1(NF,X)=12.D0*PI/(B0(NF)*DLOG(X**2/XLB(NF)**2))
1659
ALS2(NF,X)=12.D0*PI/(B0(NF)*DLOG(X**2/XLB(NF)**2))
1660
. *(1.D0-B1(NF)*DLOG(DLOG(X**2/XLB(NF)**2))
1661
. /DLOG(X**2/XLB(NF)**2))
1674
ELSEIF(Q.LE.AMB)THEN
1676
ELSEIF(Q.LE.AMT)THEN
1682
ALPHAS_HD=ALS1(NF,Q)
1684
ALPHAS_HD=ALS2(NF,Q)
1689
SUBROUTINE ALSINI(ACC)
1690
IMPLICIT DOUBLE PRECISION (A-H,O-Z)
1692
COMMON/ALSLAM/XLB1(6),XLB2(6)
1693
COMMON/ALS/XLAMBDA,AMC,AMB,AMT,N0
1701
XLB(4)=XLB(3)*(XLB(3)/AMC)**(2.D0/25.D0)
1702
XLB(5)=XLB(4)*(XLB(4)/AMB)**(2.D0/23.D0)
1703
XLB(6)=XLB(5)*(XLB(5)/AMT)**(2.D0/21.D0)
1706
XLB(5)=XLB(4)*(XLB(4)/AMB)**(2.D0/23.D0)
1707
XLB(3)=XLB(4)*(XLB(4)/AMC)**(-2.D0/27.D0)
1708
XLB(6)=XLB(5)*(XLB(5)/AMT)**(2.D0/21.D0)
1711
XLB(4)=XLB(5)*(XLB(5)/AMB)**(-2.D0/25.D0)
1712
XLB(3)=XLB(4)*(XLB(4)/AMC)**(-2.D0/27.D0)
1713
XLB(6)=XLB(5)*(XLB(5)/AMT)**(2.D0/21.D0)
1716
XLB(5)=XLB(6)*(XLB(6)/AMT)**(-2.D0/23.D0)
1717
XLB(4)=XLB(5)*(XLB(5)/AMB)**(-2.D0/25.D0)
1718
XLB(3)=XLB(4)*(XLB(4)/AMC)**(-2.D0/27.D0)
1725
XLB(4)=XLB(3)*(XLB(3)/AMC)**(2.D0/25.D0)
1726
. *(2.D0*DLOG(AMC/XLB(3)))**(-107.D0/1875.D0)
1727
XLB(4)=XITER(AMC,XLB(3),3,XLB(4),4,ACC)
1728
XLB(5)=XLB(4)*(XLB(4)/AMB)**(2.D0/23.D0)
1729
. *(2.D0*DLOG(AMB/XLB(4)))**(-963.D0/13225.D0)
1730
XLB(5)=XITER(AMB,XLB(4),4,XLB(5),5,ACC)
1731
XLB(6)=XLB(5)*(XLB(5)/AMT)**(2.D0/21.D0)
1732
. *(2.D0*DLOG(AMT/XLB(5)))**(-321.D0/3381.D0)
1733
XLB(6)=XITER(AMT,XLB(5),5,XLB(6),6,ACC)
1736
XLB(5)=XLB(4)*(XLB(4)/AMB)**(2.D0/23.D0)
1737
. *(2.D0*DLOG(AMB/XLB(4)))**(-963.D0/13225.D0)
1738
XLB(5)=XITER(AMB,XLB(4),4,XLB(5),5,ACC)
1739
XLB(3)=XLB(4)*(XLB(4)/AMC)**(-2.D0/27.D0)
1740
. *(2.D0*DLOG(AMC/XLB(4)))**(107.D0/2025.D0)
1741
XLB(3)=XITER(AMC,XLB(4),4,XLB(3),3,ACC)
1742
XLB(6)=XLB(5)*(XLB(5)/AMT)**(2.D0/21.D0)
1743
. *(2.D0*DLOG(AMT/XLB(5)))**(-321.D0/3381.D0)
1744
XLB(6)=XITER(AMT,XLB(5),5,XLB(6),6,ACC)
1747
XLB(4)=XLB(5)*(XLB(5)/AMB)**(-2.D0/25.D0)
1748
. *(2.D0*DLOG(AMB/XLB(5)))**(963.D0/14375.D0)
1749
XLB(4)=XITER(AMB,XLB(5),5,XLB(4),4,ACC)
1750
XLB(3)=XLB(4)*(XLB(4)/AMC)**(-2.D0/27.D0)
1751
. *(2.D0*DLOG(AMC/XLB(4)))**(107.D0/2025.D0)
1752
XLB(3)=XITER(AMC,XLB(4),4,XLB(3),3,ACC)
1753
XLB(6)=XLB(5)*(XLB(5)/AMT)**(2.D0/21.D0)
1754
. *(2.D0*DLOG(AMT/XLB(5)))**(-321.D0/3381.D0)
1755
XLB(6)=XITER(AMT,XLB(5),5,XLB(6),6,ACC)
1758
XLB(5)=XLB(6)*(XLB(6)/AMT)**(-2.D0/23.D0)
1759
. *(2.D0*DLOG(AMT/XLB(6)))**(321.D0/3703.D0)
1760
XLB(5)=XITER(AMT,XLB(6),6,XLB(5),5,ACC)
1761
XLB(4)=XLB(5)*(XLB(5)/AMB)**(-2.D0/25.D0)
1762
. *(2.D0*DLOG(AMB/XLB(5)))**(963.D0/14375.D0)
1763
XLB(4)=XITER(AMB,XLB(5),5,XLB(4),4,ACC)
1764
XLB(3)=XLB(4)*(XLB(4)/AMC)**(-2.D0/27.D0)
1765
. *(2.D0*DLOG(AMC/XLB(4)))**(107.D0/2025.D0)
1766
XLB(3)=XITER(AMC,XLB(4),4,XLB(3),3,ACC)
1774
DOUBLE PRECISION FUNCTION XITER(Q,XLB1,NF1,XLB,NF2,ACC)
1775
IMPLICIT DOUBLE PRECISION (A-H,O-Z)
1776
B0(NF)=33.D0-2.D0*NF
1777
B1(NF)=6.D0*(153.D0-19.D0*NF)/B0(NF)**2
1778
ALS2(NF,X,XLB)=12.D0*PI/(B0(NF)*DLOG(X**2/XLB**2))
1779
. *(1.D0-B1(NF)*DLOG(DLOG(X**2/XLB**2))
1780
. /DLOG(X**2/XLB**2))
1781
AA(NF)=12D0*PI/B0(NF)
1782
BB(NF)=B1(NF)/AA(NF)
1783
XIT(A,B,X)=A/2.D0*(1D0+DSQRT(1D0-4D0*B*DLOG(X)))
1788
X=DLOG(Q**2/XLB2**2)
1789
ALP=ALS2(NF1,Q,XLB1)
1793
XLB2=Q*DEXP(-XX/2.D0)
1797
IF(DY.GE.ACC) GOTO 1
1802
DOUBLE PRECISION FUNCTION FINT(Z,XX,YY)
1803
C--ONE-DIMENSIONAL CUBIC INTERPOLATION
1805
C--XX = ARRAY OF 4 DISCRETE X-VALUES AROUND Z
1806
C--YY = ARRAY OF 4 DISCRETE FUNCTION-VALUES AROUND Z
1807
IMPLICIT DOUBLE PRECISION (A-H,O-Z)
1808
DIMENSION XX(4),YY(4)
1818
A0=(X-X1)*(X-X2)*(X-X3)/(X0-X1)/(X0-X2)/(X0-X3)
1819
A1=(X-X0)*(X-X2)*(X-X3)/(X1-X0)/(X1-X2)/(X1-X3)
1820
A2=(X-X0)*(X-X1)*(X-X3)/(X2-X0)/(X2-X1)/(X2-X3)
1821
A3=(X-X0)*(X-X1)*(X-X2)/(X3-X0)/(X3-X1)/(X3-X2)
1822
FINT=DEXP(A0*Y0+A1*Y1+A2*Y2+A3*Y3)
1826
DOUBLE PRECISION FUNCTION SP(X)
1827
C--REAL DILOGARITHM (SPENCE-FUNCTION)
1828
IMPLICIT DOUBLE PRECISION (A-H,O-Z)
1835
COMPLEX*16 FUNCTION LI2(X)
1836
C--COMPLEX DILOGARITHM (SPENCE-FUNCTION)
1837
IMPLICIT DOUBLE PRECISION (A-H,O-Z)
1839
COMMON/CONST/ZETA2,ZETA3
1850
IF(R2.EQ.1.D0.AND.XI.EQ.0.D0)THEN
1854
LI2=-DCMPLX(ZETA2/2.D0)
1857
ELSEIF(R2.GT.1.D0.AND.RR.GT.0.5D0)THEN
1859
LI2=CLI2(Y)+ZETA2-CDLOG(X)*CDLOG(1.D0-X)+0.5D0*CDLOG(X)**2
1861
ELSEIF(R2.GT.1.D0.AND.RR.LE.0.5D0)THEN
1863
LI2=-CLI2(Y)-ZETA2-0.5D0*CDLOG(-X)**2
1865
ELSEIF(R2.LE.1.D0.AND.XR.GT.0.5D0)THEN
1867
LI2=-CLI2(Y)+ZETA2-CDLOG(X)*CDLOG(1.D0-X)
1869
ELSEIF(R2.LE.1.D0.AND.XR.LE.0.5D0)THEN
1876
COMPLEX*16 FUNCTION CLI2(X)
1877
C--TAYLOR-EXPANSION FOR COMPLEX DILOGARITHM (SPENCE-FUNCTION)
1878
IMPLICIT DOUBLE PRECISION (A-H,O-Z)
1880
COMMON/BERNOULLI/B2(18),B12(18),B3(18)
1892
DOUBLE PRECISION FUNCTION FACULT(N)
1893
C--DOUBLE PRECISION VERSION OF FACULTY
1894
IMPLICIT DOUBLE PRECISION (A-H,O-Z)
1898
FACULT=FACULT*DFLOAT(I)
1903
SUBROUTINE BERNINI(N)
1904
C--INITIALIZATION OF COEFFICIENTS FOR POLYLOGARITHMS
1905
IMPLICIT DOUBLE PRECISION (A-H,O-Z)
1906
DIMENSION B(18),PB(19)
1907
COMMON/BERNOULLI/B2(18),B12(18),B3(18)
1908
COMMON/CONST/ZETA2,ZETA3
1925
B(12)=-691.D0/2730.D0
1929
B(16)=-3617.D0/510.D0
1931
B(18)=43867.D0/798.D0
1933
ZETA3=1.202056903159594D0
1936
B2(I)=B(I)/FACULT(I+1)
1937
B12(I)=DFLOAT(I+1)/FACULT(I+2)*B(I)/2.D0
1944
B3(I)=B3(I)+PB(J+1)*PB(I-J+1)/FACULT(I-J)/FACULT(J+1)
1951
DOUBLE PRECISION FUNCTION QQINT(RAT,H1,H2)
1952
IMPLICIT DOUBLE PRECISION (A-H,O-Z)
1954
QQINT = RAT**N * H1 + (1-RAT**N) * H2
1958
DOUBLE PRECISION FUNCTION XITLA(NO,ALP,ACC)
1959
C--ITERATION ROUTINE TO DETERMINE IMPROVED LAMBDA'S
1960
IMPLICIT DOUBLE PRECISION (A-H,O-Z)
1961
COMMON/PARAM/GF,ALPH,AMTAU,AMMUON,AMZ,AMW
1962
B0(NF)=33.D0-2.D0*NF
1963
B1(NF)=6.D0*(153.D0-19.D0*NF)/B0(NF)**2
1964
ALS2(NF,X,XLB)=12.D0*PI/(B0(NF)*DLOG(X**2/XLB**2))
1965
. *(1.D0-B1(NF)*DLOG(DLOG(X**2/XLB**2))
1966
. /DLOG(X**2/XLB**2))
1967
AA(NF)=12D0*PI/B0(NF)
1968
BB(NF)=B1(NF)/AA(NF)
1969
XIT(A,B,X)=A/2.D0*(1D0+DSQRT(1D0-4D0*B*DLOG(X)))
1973
XLB=Q*DEXP(-AA(NF)/ALP/2.D0)
1981
XLB=Q*DEXP(-XX/2.D0)
1985
IF(DY.GE.ACC) GOTO 1