← Back to branch summary

~mg5core1/mg5amcnlo/2.6.4

« back to all changes in this revision

Viewing changes to DECAY/alfas_functions.f

move ./decay to ./mg5decay; resolve unit tests (n.b. __init__ does not check keys of input dictionaries, followed last revision)

Show diffs side-by-side

added added

removed removed

Lines of Context:
DECAY/alfas_functions.f
 
1
C
 
2
C-----------------------------------------------------------------------------
 
3
C
 
4
      double precision function alfa(alfa0,qsq )
 
5
C
 
6
C-----------------------------------------------------------------------------
 
7
C
 
8
C       This function returns the 1-loop value of alpha.
 
9
C
 
10
C       INPUT: 
 
11
C               qsq   = Q^2
 
12
C
 
13
C-----------------------------------------------------------------------------
 
14
C
 
15
      implicit none
 
16
      double precision  qsq,alfa0
 
17
c
 
18
c constants
 
19
c
 
20
      double precision  One, Three, Pi,zmass
 
21
      parameter( One = 1.0d0, Three = 3.0d0 )
 
22
      parameter( Pi = 3.14159265358979323846d0 )
 
23
      parameter( zmass = 91.188d0 )
 
24
cc
 
25
      alfa = alfa0 / ( 1.0d0 - alfa0*dlog( qsq/zmass**2 ) /Three /Pi )
 
26
ccc
 
27
      return
 
28
      end
 
29
 
 
30
C
 
31
C-----------------------------------------------------------------------------
 
32
C
 
33
      double precision function alfaw(alfaw0,qsq,nh )
 
34
C
 
35
C-----------------------------------------------------------------------------
 
36
C
 
37
C       This function returns the 1-loop value of alpha_w.
 
38
C
 
39
C       INPUT: 
 
40
C               qsq = Q^2
 
41
C               nh  = # of Higgs doublets
 
42
C
 
43
C-----------------------------------------------------------------------------
 
44
C
 
45
      implicit none
 
46
      double precision  qsq, alphaw, dum,alfaw0
 
47
      integer  nh, nq
 
48
c
 
49
c         include
 
50
c
 
51
          
 
52
c
 
53
c constants
 
54
c
 
55
      double precision  Two, Four, Pi, Twpi, zmass,tmass
 
56
      parameter( Two = 2.0d0, Four = 4.0d0 )
 
57
      parameter( Pi = 3.14159265358979323846d0 )
 
58
      parameter( Twpi = 3.0d0*Four*Pi )
 
59
      parameter( zmass = 91.188d0,tmass=174d0 )
 
60
cc
 
61
      if ( qsq.ge.tmass**2 ) then
 
62
         nq = 6
 
63
      else
 
64
         nq = 5
 
65
      end if
 
66
      dum = ( 22.0d0 - Four*nq - nh/Two ) / Twpi
 
67
      alfaw = alfaw0 / ( 1.0d0 + dum*alfaw0*dlog( qsq/zmass**2 ) )
 
68
ccc
 
69
      return
 
70
      end
 
71
 
 
72
C-----------------------------------------------------------------------------
 
73
C
 
74
      DOUBLE PRECISION FUNCTION ALPHAS(Q)
 
75
c
 
76
c     Evaluation of strong coupling constant alpha_S
 
77
c     Author: R.K. Ellis
 
78
c
 
79
c     q -- scale at which alpha_s is to be evaluated
 
80
c
 
81
c-- common block alfas.inc
 
82
c     asmz -- value of alpha_s at the mass of the Z-boson
 
83
c     nloop -- the number of loops (1,2, or 3) at which beta 
 
84
c
 
85
c     function is evaluated to determine running.
 
86
c     the values of the cmass and the bmass should be set
 
87
c     in common block qmass.
 
88
C-----------------------------------------------------------------------------
 
89
 
 
90
      IMPLICIT NONE
 
91
c
 
92
      include 'alfas.inc'
 
93
      DOUBLE PRECISION Q,T,AMZ0,AMB,AMC
 
94
      DOUBLE PRECISION AS_OUT
 
95
      INTEGER NLOOP0,NF3,NF4,NF5
 
96
      PARAMETER(NF5=5,NF4=4,NF3=3)
 
97
C
 
98
      REAL*8       CMASS,BMASS
 
99
      COMMON/QMASS/CMASS,BMASS
 
100
      DATA CMASS,BMASS/1.42D0,4.7D0/  ! HEAVY QUARK MASSES FOR THRESHOLDS
 
101
C
 
102
      REAL*8 ZMASS
 
103
      DATA ZMASS/91.188D0/
 
104
C
 
105
      SAVE AMZ0,NLOOP0,AMB,AMC
 
106
      DATA AMZ0,NLOOP0/0D0,0/
 
107
      IF (Q .LE. 0D0) THEN 
 
108
         WRITE(6,*) 'q .le. 0 in alphas'
 
109
         WRITE(6,*) 'q= ',Q
 
110
         STOP
 
111
      ENDIF
 
112
      IF (asmz .LE. 0D0) THEN 
 
113
         WRITE(6,*) 'asmz .le. 0 in alphas',asmz
 
114
c         WRITE(6,*) 'continue with asmz=0.1185'
 
115
         STOP
 
116
         asmz=0.1185D0
 
117
      ENDIF
 
118
      IF (CMASS .LE. 0.3D0) THEN 
 
119
         WRITE(6,*) 'cmass .le. 0.3GeV in alphas',CMASS
 
120
c         WRITE(6,*) 'continue with cmass=1.5GeV?'
 
121
         STOP
 
122
         CMASS=1.42D0
 
123
      ENDIF
 
124
      IF (BMASS .LE. 0D0) THEN 
 
125
         WRITE(6,*) 'bmass .le. 0 in alphas',BMASS
 
126
         WRITE(6,*) 'COMMON/QMASS/CMASS,BMASS'
 
127
         STOP
 
128
         BMASS=4.7D0
 
129
      ENDIF
 
130
c--- establish value of coupling at b- and c-mass and save
 
131
      IF ((asmz .NE. AMZ0) .OR. (NLOOP .NE. NLOOP0)) THEN
 
132
         AMZ0=asmz
 
133
         NLOOP0=NLOOP
 
134
         T=2D0*DLOG(BMASS/ZMASS)
 
135
         CALL NEWTON1(T,asmz,AMB,NLOOP,NF5)
 
136
         T=2D0*DLOG(CMASS/BMASS)
 
137
         CALL NEWTON1(T,AMB,AMC,NLOOP,NF4)
 
138
      ENDIF
 
139
 
 
140
c--- evaluate strong coupling at scale q
 
141
      IF (Q  .LT. BMASS) THEN
 
142
           IF (Q  .LT. CMASS) THEN
 
143
             T=2D0*DLOG(Q/CMASS)
 
144
             CALL NEWTON1(T,AMC,AS_OUT,NLOOP,NF3)
 
145
           ELSE
 
146
             T=2D0*DLOG(Q/BMASS)
 
147
             CALL NEWTON1(T,AMB,AS_OUT,NLOOP,NF4)
 
148
           ENDIF
 
149
      ELSE
 
150
      T=2D0*DLOG(Q/ZMASS)
 
151
      CALL NEWTON1(T,asmz,AS_OUT,NLOOP,NF5)
 
152
      ENDIF
 
153
      ALPHAS=AS_OUT
 
154
      RETURN
 
155
      END
 
156
 
 
157
 
 
158
      SUBROUTINE NEWTON1(T,A_IN,A_OUT,NLOOP,NF)
 
159
C     Author: R.K. Ellis
 
160
 
 
161
c---  calculate a_out using nloop beta-function evolution 
 
162
c---  with nf flavours, given starting value as-in
 
163
c---  given as_in and logarithmic separation between 
 
164
c---  input scale and output scale t.
 
165
c---  Evolution is performed using Newton's method,
 
166
c---  with a precision given by tol.
 
167
 
 
168
      IMPLICIT NONE
 
169
      INTEGER NLOOP,NF
 
170
      REAL*8 T,A_IN,A_OUT,AS,TOL,F2,F3,F,FP,DELTA
 
171
      REAL*8 B0(3:5),C1(3:5),C2(3:5),DEL(3:5)
 
172
      PARAMETER(TOL=5.D-4)
 
173
C---     B0=(11.-2.*NF/3.)/4./PI
 
174
      DATA B0/0.716197243913527D0,0.66314559621623D0,0.61009394851893D0/
 
175
C---     C1=(102.D0-38.D0/3.D0*NF)/4.D0/PI/(11.D0-2.D0/3.D0*NF)
 
176
      DATA C1/.565884242104515D0,0.49019722472304D0,0.40134724779695D0/
 
177
C---     C2=(2857.D0/2.D0-5033*NF/18.D0+325*NF**2/54)
 
178
C---     /16.D0/PI**2/(11.D0-2.D0/3.D0*NF)
 
179
      DATA C2/0.453013579178645D0,0.30879037953664D0,0.14942733137107D0/
 
180
C---     DEL=SQRT(4*C2-C1**2)
 
181
      DATA DEL/1.22140465909230D0,0.99743079911360D0,0.66077962451190D0/
 
182
      F2(AS)=1D0/AS+C1(NF)*LOG((C1(NF)*AS)/(1D0+C1(NF)*AS))
 
183
      F3(AS)=1D0/AS+0.5D0*C1(NF)
 
184
     & *LOG((C2(NF)*AS**2)/(1D0+C1(NF)*AS+C2(NF)*AS**2))
 
185
     & -(C1(NF)**2-2D0*C2(NF))/DEL(NF)
 
186
     & *ATAN((2D0*C2(NF)*AS+C1(NF))/DEL(NF))
 
187
 
 
188
           
 
189
      A_OUT=A_IN/(1D0+A_IN*B0(NF)*T)
 
190
      IF (NLOOP .EQ. 1) RETURN
 
191
      A_OUT=A_IN/(1D0+B0(NF)*A_IN*T+C1(NF)*A_IN*LOG(1D0+A_IN*B0(NF)*T))
 
192
      IF (A_OUT .LT. 0D0) AS=0.3D0
 
193
 30   AS=A_OUT
 
194
 
 
195
      IF (NLOOP .EQ. 2) THEN
 
196
      F=B0(NF)*T+F2(A_IN)-F2(AS)
 
197
      FP=1D0/(AS**2*(1D0+C1(NF)*AS))
 
198
      ENDIF
 
199
      IF (NLOOP .EQ. 3) THEN
 
200
      F=B0(NF)*T+F3(A_IN)-F3(AS)
 
201
      FP=1D0/(AS**2*(1D0+C1(NF)*AS+C2(NF)*AS**2))
 
202
      ENDIF
 
203
      A_OUT=AS-F/FP
 
204
      DELTA=ABS(F/FP/AS)
 
205
      IF (DELTA .GT. TOL) GO TO 30
 
206
      RETURN
 
207
      END
 
208
 
 
209
 
 
210
C-----------------------------------------------------------------------------
 
211
C
 
212
      double precision function mfrun(mf,scale,asmz,nloop)
 
213
C
 
214
C-----------------------------------------------------------------------------
 
215
C
 
216
C       This function returns the 2-loop value of a MSbar fermion mass
 
217
C       at a given scale.
 
218
C
 
219
C       INPUT: mf    = MSbar mass of fermion at MSbar fermion mass scale 
 
220
C              scale = scale at which the running mass is evaluated
 
221
C              asmz  = AS(MZ) : this is passed to alphas(scale,asmz,nloop)
 
222
C              nloop = # of loops in the evolution
 
223
C       
 
224
C
 
225
C
 
226
C       EXTERNAL:      double precision alphas(scale,asmz,nloop)
 
227
C                      
 
228
C-----------------------------------------------------------------------------
 
229
C
 
230
      implicit none
 
231
C
 
232
C     ARGUMENTS
 
233
C
 
234
      double precision  mf,scale,asmz
 
235
      integer           nloop
 
236
C
 
237
C     LOCAL
 
238
C
 
239
      double precision  beta0, beta1,gamma0,gamma1
 
240
      double precision  A1,as,asmf,l2
 
241
      integer  nf
 
242
C
 
243
C     EXTERNAL
 
244
C
 
245
      double precision  alphas
 
246
      external          alphas
 
247
c
 
248
c     CONSTANTS
 
249
c
 
250
      double precision  One, Two, Three, Pi
 
251
      parameter( One = 1.0d0, Two = 2.0d0, Three = 3.0d0 )
 
252
      parameter( Pi = 3.14159265358979323846d0) 
 
253
      double precision tmass
 
254
      parameter(tmass=174d0)
 
255
cc
 
256
C
 
257
C
 
258
      if ( mf.gt.tmass ) then
 
259
         nf = 6
 
260
      else
 
261
         nf = 5
 
262
      end if
 
263
 
 
264
      beta0 = ( 11.0d0 - Two/Three *nf )/4d0
 
265
      beta1 = ( 102d0  - 38d0/Three*nf )/16d0
 
266
      gamma0= 1d0
 
267
      gamma1= ( 202d0/3d0  - 20d0/9d0*nf )/16d0
 
268
      A1    = -beta1*gamma0/beta0**2+gamma1/beta0
 
269
      as    = alphas(scale)
 
270
      asmf  = alphas(mf)
 
271
      l2    = (1+ A1*as/Pi)/(1+ A1*asmf/Pi)
 
272
      
 
273
      
 
274
      mfrun = mf * (as/asmf)**(gamma0/beta0)
 
275
 
 
276
      if(nloop.eq.2) mfrun =mfrun*l2
 
277
ccc
 
278
      return
 
279
      end
 
280