LHAPDF/lhapdf5.2.2/wrapabfkwpi.f

   1       subroutine ABFKWPevolve(xin,qin,pdf)
   2       include 'parmsetup.inc'
   3       PARAMETER(NX=50)
   4       PARAMETER(NQ=19)
   5       real*8 xin,qin,pdf(-6:6),xval(45),qcdl4,qcdl5
   6       real*8 upv,dnv,usea,dsea,str,chm,bot,top,glu
   7       real*8 calcpi(8,20,25,3),calcpio(8,20,25),parpi(40,3)
   8       common /ABFKWP/ CALCPI,CALCPIO,PARPI,lastmem
   9       character*16 name(nmxset)
  10       integer nmem(nmxset),ndef(nmxset),mmem
  11       common/NAME/name,nmem,ndef,mmem
  12       integer nset
  13       save
  14
  15       iimem = imem
  16       if(iimem.eq.0) iimem = 1
  17       if(iimem.le.3) then
  18        call ABFKWxx(iimem,xin,qin,upv,dnv,usea,dsea, str,chm,glu)
  19       endif
  20
  21
  22       pdf(-6)= 0.0d0
  23       pdf(6)= 0.0d0
  24       pdf(-5)= 0.0d0
  25       pdf(5 )= 0.0d0
  26       pdf(-4)= chm
  27       pdf(4 )= chm
  28       pdf(-3)= str
  29       pdf(3 )= str
  30       pdf(-2)= usea
  31       pdf(2 )= upv+usea
  32       pdf(-1)= dsea
  33       pdf(1 )= dnv+dsea
  34       pdf(0 )= glu
  35
  36       return
  37 ccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
  38       entry ABFKWPread(nset)
  39       read(1,*)nmem(nset),ndef(nset)
  40 c      print *,nmem,ndef
  41       lastmem = -999
  42       do j=1,3
  43         read(1,*)(parpi(k,j),k=1,4)
  44         read(1,*)(parpi(k,j),k=5,8)
  45         read(1,*)(parpi(k,j),k=9,12)
  46         read(1,*)(parpi(k,j),k=13,16)
  47         read(1,*)(parpi(k,j),k=17,20)
  48         read(1,*)(parpi(k,j),k=21,24)
  49         read(1,*)(parpi(k,j),k=25,28)
  50         read(1,*)(parpi(k,j),k=29,32)
  51         read(1,*)(parpi(k,j),k=33,36)
  52         read(1,*)(parpi(k,j),k=37,40)
  53         do l=1,25
  54           do k=1,20
  55             read(1,*)(CALCPI(m,k,l,j),m=1,4)
  56             read(1,*)(CALCPI(m,k,l,j),m=5,8)
  57           enddo
  58         enddo
  59       enddo
  60       return
  61 c
  62 ccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
  63       entry ABFKWPalfa(alfas,qalfa)
  64         call getnset(iset)
  65         call GetOrderAsM(iset,iord)
  66         call Getlam4M(iset,imem,qcdl4)
  67         call Getlam5M(iset,imem,qcdl5)
  68         call aspdflib(alfas,Qalfa,iord,qcdl5)
  69       return
  70 c
  71 ccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
  72       entry ABFKWPinit(Eorder,Q2fit)
  73       return
  74 c
  75 ccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
  76       entry ABFKWPpdf(mem)
  77       imem = mem
  78       return
  79 c
  80  1000 format(5e13.5)
  81       end
  82 ccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
  83 *
  84 * $Id$
  85 *
  86 * $Log$
  87 * Revision 1.2  2005/10/07 15:15:05  whalley
  88 * Changes to most files for V5 - multiset initializations
  89 *
  90 * Revision 1.1.1.1  2005/05/06 14:54:43  whalley
  91 * Initial CVS import of the LHAPDF code and data sets
  92 *
  93 * Revision 1.1.1.2  1996/10/30 08:27:26  cernlib
  94 * Version 7.04
  95 *
  96 * Revision 1.1.1.1  1996/04/12 15:28:53  plothow
  97 * Version 7.01
  98 *
  99 *
 100       SUBROUTINE ABFKWxx(imem,DX,DQ,DUPV,DDNV,DUSEA,DDSEA,DSTR,DCHM,DGL)
 101       double precision
 102      +       PARPI(40,3),CALCPI(8,20,25,3),CALCPIO(8,20,25),ZEROD,
 103      +       DX,DQ,DUPV,DDNV,DUSEA,DDSEA,DSTR,DCHM,DGL
 104       REAL    X, Q, UPV, DNV, USEA, DSEA, STR, CHM, GL
 105
 106       common /ABFKWP/CALCPI,CALCPIO,PARPI,lastmem
 107
 108 c      COMMON/W5051Ixx/CALCPIO
 109       REAL   XPDF(7)
 110       DATA ZEROD/0.D0/
 111 C----------------------------------------------------------------------
 112        DATA ISTART/0/
 113        SAVE ISTART,OWLAM2,Q02PI
 114 C
 115       if(imem.ne.lastmem) then
 116          istart = 0
 117          lastmem = imem
 118       endif
 119       IF (ISTART.EQ.0) THEN
 120         ISTART=1
 121         DO 11 K=1,25
 122         DO 11 I=1,20
 123         DO 11 M=1,8
 124    11   CALCPIO(M,I,K) = CALCPI(M,I,K,imem)
 125            OWLAM=PARPI(1,imem)
 126            OWLAM2=OWLAM**2
 127            Q02PI=PARPI(39,imem)
 128            Q2MAX=PARPI(40,imem)
 129          ENDIF
 130 C
 131 C the conventions are : q(1)=x*u, q(2)=x*d, q(3)=x*str, q(4)=x*usea,
 132 C                       q(5)=x*dsea, q(6)=x*charm, q(7)=x*gluon
 133 C
 134       X = DX
 135       Q = DQ
 136       Q2 = Q*Q
 137       IDQ2=2
 138       SB=0.
 139       IF(Q2-Q02PI) 1,1,2
 140     2 IF(IDQ2-1) 1,1,3
 141     3 SB= LOG( LOG( MAX(Q02PI,Q2)/OWLAM2)/ LOG(Q02PI/OWLAM2))
 142     1 CALL AURPIx(1,0,X,SB,XPDF(1))
 143       CALL AURPIx(2,0,X,SB,XPDF(2))
 144       CALL AURPIx(3,0,X,SB,XPDF(3))
 145       CALL AURPIx(4,0,X,SB,XPDF(4))
 146       CALL AURPIx(5,0,X,SB,XPDF(5))
 147       CALL AURPIx(8,0,X,SB,XPDF(6))
 148       CALL AURPIx(7,0,X,SB,XPDF(7))
 149 C
 150       DUPV=XPDF(1) - XPDF(4)
 151       DDNV=XPDF(2) - XPDF(5)
 152       DUSEA=XPDF(4)
 153       DDSEA=XPDF(5)
 154       DSTR=XPDF(3)
 155       DCHM=XPDF(6)
 156       DGL =XPDF(7)
 157 C
 158       RETURN
 159       END
 160 c==============================================================
 161 *
 162 * $Id$
 163 *
 164 * $Log$
 165 * Revision 1.2  2005/10/07 15:15:05  whalley
 166 * Changes to most files for V5 - multiset initializations
 167 *
 168 * Revision 1.1.1.1  2005/05/06 14:54:43  whalley
 169 * Initial CVS import of the LHAPDF code and data sets
 170 *
 171 * Revision 1.1.1.2  1996/10/30 08:27:36  cernlib
 172 * Version 7.04
 173 *
 174 * Revision 1.1.1.1  1996/04/12 15:29:03  plothow
 175 * Version 7.01
 176 *
 177 *
 178 C
 179       SUBROUTINE AURPIx(I,NDRV,X,S,ANS)
 180       double precision
 181      +       CALCPI(8,20,25,3),CALCPIO(8,20,25),parpi(40,3)
 182       common /ABFKWP/CALCPI,CALCPIO,parpi,lastmem
 183 c      COMMON/W5051I4/CALCPIO
 184       REAL   F1(25),F2(25)
 185       DATA DELTA/.10/
 186       ANS=0.
 187       IF(X.GT.0.9985) RETURN
 188       IF(I.EQ.3.AND.X.GT.0.95) RETURN
 189       IF(I.EQ.8.AND.X.GT.0.95) RETURN
 190       IS=S/DELTA+1
 191       IS1=IS+1
 192       DO 1 L=1,25
 193       KL=L+NDRV*25
 194       F1(L)=CALCPIO(I,IS,KL)
 195       F2(L)=CALCPIO(I,IS1,KL)
 196     1 CONTINUE
 197       A1=AUGETFV(X,F1)
 198       A2=AUGETFV(X,F2)
 199       S1=(IS-1)*DELTA
 200       S2=S1+DELTA
 201       ANS=A1*(S-S2)/(S1-S2)+A2*(S-S1)/(S2-S1)
 202       RETURN
 203       END
 204 c===============================================================
 205 *
 206 * $Id$
 207 *
 208 * $Log$
 209 * Revision 1.2  2005/10/07 15:15:05  whalley
 210 * Changes to most files for V5 - multiset initializations
 211 *
 212 * Revision 1.1.1.1  2005/05/06 14:54:43  whalley
 213 * Initial CVS import of the LHAPDF code and data sets
 214 *
 215 * Revision 1.1.1.2  1996/10/30 08:27:34  cernlib
 216 * Version 7.04
 217 *
 218 * Revision 1.1.1.1  1996/04/12 15:29:02  plothow
 219 * Version 7.01
 220 *
 221 *
 222 C
 223       FUNCTION AUGETFV(X,FVL)
 224 C  LOGARITHMIC INTERPOLATOR - WATCH OUT FOR NEGATIVE
 225 C  FUNCTIONS AND/OR X VALUES OUTSIDE THE RANGE 0 TO 1.
 226 C  NOTE: DIMENSION OF FVL IS OVERWRITTEN BY VALUE USED
 227 C  IN MAIN ROUTINE.
 228       DIMENSION FVL(25),XGRID(25)
 229       DATA NX,XGRID/25,.001,.002,.004,.008,.016,.032,.064,.1,.15,
 230      *.2,.25,.3,.35,.4,.45,.5,.55,.6,.65,.7,.75,.8,.85,.9,.95/
 231       AUGETFV=0.
 232       DO 1 I=1,NX
 233       IF(X.LT.XGRID(I)) GO TO 2
 234     1 CONTINUE
 235     2 I=I-1
 236       IF(I.EQ.0) THEN
 237          I=I+1
 238       ELSE IF(I.GT.23) THEN
 239          I=23
 240       ENDIF
 241       J=I+1
 242       K=J+1
 243       AXI= LOG(XGRID(I))
 244       BXI= LOG(1.-XGRID(I))
 245       AXJ= LOG(XGRID(J))
 246       BXJ= LOG(1.-XGRID(J))
 247       AXK= LOG(XGRID(K))
 248       BXK= LOG(1.-XGRID(K))
 249       FI= LOG(ABS(FVL(I)) +1.E-15)
 250       FJ= LOG(ABS(FVL(J)) +1.E-16)
 251       FK= LOG(ABS(FVL(K)) +1.E-17)
 252       DET=AXI*(BXJ-BXK)+AXJ*(BXK-BXI)+AXK*(BXI-BXJ)
 253       ALOGA=(FI*(AXJ*BXK-AXK*BXJ)+FJ*(AXK*BXI-AXI*BXK)+FK*(AXI*BXJ-AXJ*
 254      $ BXI))/DET
 255       ALPHA=(FI*(BXJ-BXK)+FJ*(BXK-BXI)+FK*(BXI-BXJ))/DET
 256       BETA=(FI*(AXK-AXJ)+FJ*(AXI-AXK)+FK*(AXJ-AXI))/DET
 257       IF(ABS(ALPHA).GT.99..OR.ABS(BETA).GT.99..OR.ABS(ALOGA).GT.99.)
 258      1RETURN
 259 C      IF(ALPHA.GT.50..OR.BETA.GT.50.) THEN
 260 C         WRITE(6,2001) X,FVL
 261 C 2001    FORMAT(8E12.4)
 262 C         WRITE(6,2001) ALPHA,BETA,ALOGA,DET
 263 C      ENDIF
 264       AUGETFV=EXP(ALOGA)*X**ALPHA*(1.-X)**BETA
 265       RETURN
 266       END
 267 c============================================================