1 subroutine CTEQ5evolve(x,Q,pdf)
2 implicit real*8(a-h,o-z)
3 include 'parmsetup.inc'
4 character*16 name(nmxset)
5 integer nmem(nmxset),ndef(nmxset),mmem
6 common/NAME/name,nmem,ndef,mmem
10 PARAMETER (MXX = 105, MXQ = 25, MXF = 6)
11 PARAMETER (MXPQX = (MXF *2 +2) * MXQ * MXX)
12 PARAMETER (M= 2, M1 = M + 1)
14 > / CtqPar1 / Al, XV(0:MXX), QL(0:MXQ), UPD(MXPQX)
15 > / CtqPar2 / Nx, Nt, NfMx
16 > / XQrange / Qini, Qmax, Xmin
17 > / QCDtable / Alambda, Nfl, Iorder
18 > / Masstbl / Amass(6)
19 data pi / 3.141592653589793d0 /
22 U = X * CtLhCtq5Pdf(1,X,Q)
23 D = X * CtLhCtq5Pdf(2,X,Q)
24 USEA = X * CtLhCtq5Pdf(-1,X,Q)
25 DSEA = X * CtLhCtq5Pdf(-2,X,Q)
26 STR = X * CtLhCtq5Pdf(3,X,Q)
27 CHM = X * CtLhCtq5Pdf(4,X,Q)
28 BOT = X * CtLhCtq5Pdf(5,X,Q)
29 GLU = X * CtLhCtq5Pdf(0,X,Q)
54 read(1,*)nmem(nset),ndef(nset)
57 Read (1, *) Dr, Fl, Al, (Amass(I),I=1,6)
63 Read (1, *) NX, NT, NfMx
66 Read (1, *) QINI, QMAX, (QL(I), I =0, NT)
69 Read (1, *) XMIN, (XV(I), I =0, NX)
72 QL(Iq) = Log (QL(Iq) /Al)
75 C Since quark = anti-quark for nfl>2 at this stage,
76 C we Read out only the non-redundent data points
77 C No of flavors = NfMx (sea) + 1 (gluon) + 2 (valence)
79 Nblk = (NX+1) * (NT+1)
80 Npts = Nblk * (NfMx+3)
83 Read (1, *, IOSTAT=IRET) (UPD(I), I=1,Npts)
87 entry CTEQ5alfa(alfas,Qalfa)
89 alfas = pi*CtLhALPI(Qalfa)
93 entry CTEQ5init(Eorder,Q2fit)
97 entry CTEQ5pdf(mem) !error corrected (jcp)
100 call setnmem(iset,mem)
105 FUNCTION CtLhPartonX5 (IPRTN, X, Q)
107 C Given the parton distribution function in the array UPD in
108 C COMMON / CtqPar1 / , this routine fetches u(fl, x, q) at any value of
109 C x and q using Mth-order polynomial interpolation for x and Ln(Q/Lambda).
111 IMPLICIT DOUBLE PRECISION (A-H, O-Z)
113 PARAMETER (MXX = 105, MXQ = 25, MXF = 6)
114 PARAMETER (MXPQX = (MXF *2 +2) * MXQ * MXX)
115 PARAMETER (M= 2, M1 = M + 1)
118 > / CtqPar1 / Al, XV(0:MXX), QL(0:MXQ), UPD(MXPQX)
119 > / CtqPar2 / Nx, Nt, NfMx
120 > / XQrange / Qini, Qmax, Xmin
122 Dimension Fq(M1), Df(M1)
124 c note, this routine doesn't have the features of not recalculating
125 c the x or q point values if they have not changed since the last call.
126 c that makes it slower than cteq6, but since cteq5 is already obsolete,
127 c speed is not so important for it.
131 C Find lower end of interval containing X
134 11 If (JU-JL .GT. 1) Then
136 If (X .GT. XV(JM)) Then
145 c crude treatment if outside the defined interval.
148 Elseif (Jx .GT. Nx-M) Then
151 C Find the interval where Q lies
154 12 If (JU-JL .GT. 1) Then
156 If (QG .GT. QL(JM)) Then
167 c If (Q .lt. Qini) Print '(A, 2(1pE12.4))',
168 c > ' WARNING: Q << Qini, extrapolation used; Q,Qini=',Q,Qini
169 Elseif (Jq .GT. Nt-M) Then
171 c If (Q .gt. Qmax) Print '(A, 2(1pE12.4))',
172 c > ' WARNING: Q > Qmax, extrapolation used; Q, Qmax =', Q, Qmax
175 If (Iprtn .GE. 3) Then
180 C Find the off-set in the linear array Upd
182 J0 = (JFL * (NT+1) + Jq) * (NX+1) + Jx
184 C Now interpolate in x for M1 Q's
186 J1 = J0 + (Nx+1)*(Iq-1) + 1
187 Call CtLhPolint3 (XV(Jx), Upd(J1), M1, X, Fq(Iq), Df(Iq))
189 C Finish off by interpolating in Q
190 Call CtLhPolint3 (QL(Jq), Fq(1), M1, QG, Ftmp, Ddf)
196 C ****************************
199 Function CtLhCtq5Pdf (Iparton, X, Q)
200 Implicit Double Precision (A-H,O-Z)
203 > / CtqPar2 / Nx, Nt, NfMx
204 > / QCDtable / Alambda, Nfl, Iorder
209 If (X .lt. 0D0 .or. X .gt. 1D0) Then
210 Print *, 'X out of range in CtLhCtq5Pdf: ', X
213 If (Q .lt. Alambda) Then
214 Print *, 'Q out of range in CtLhCtq5Pdf: ', Q
217 If ((Iparton .lt. -NfMx .or. Iparton .gt. NfMx)) Then
219 C put a warning for calling extra flavor.
221 Print *, 'Warning: Iparton out of range in CtLhCtq5Pdf: '
228 CtLhCtq5Pdf = CtLhPartonX5 (Iparton, X, Q)
229 if(CtLhCtq5Pdf.lt.0.D0) CtLhCtq5Pdf = 0.D0
233 C ********************