LHAPDF/lhapdf-5.9.1/src/wrapgjr.f

   1       subroutine GJRevolve(xin,qin,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
   7       double precision gridx(nmxgridx),gridq(nmxgridq)
   8       integer ngridx,ngridq,jx,jq
   9       CHARACTER*80 LINE
  10       dimension pdf(-6:6)
  11       integer ng(2),init,set,i,j,k,l,nset,iset
  12       double precision fgrid(118,99,-5:3,0:26),grid(217)
  13       !double precision fgrid(118,99,-5:3,0:0),grid(217)
  14 !      common/fgridc/fgrid
  15       double precision upv,dnv,usea,dsea,str,chm,bot,glu
  16       double precision arg(2)
  17       double precision lha_dfint
  18       double precision lha_gjr08
  19       data ng /118,99/
  20
  21       data grid &
  22      & /1d-9,1.25d-9,1.6d-9,2d-9,2.5d-9,3.16d-9,4d-9,5d-9,6.3d-9,8d-9, &
  23      &  1d-8,1.25d-8,1.6d-8,2d-8,2.5d-8,3.16d-8,4d-8,5d-8,6.3d-8,8d-8, &
  24      &  1d-7,1.25d-7,1.6d-7,2d-7,2.5d-7,3.16d-7,4d-7,5d-7,6.3d-7,8d-7, &
  25      &  1d-6,1.25d-6,1.6d-6,2d-6,2.5d-6,3.16d-6,4d-6,5d-6,6.3d-6,8d-6, &
  26      &  1d-5,1.25d-5,1.6d-5,2d-5,2.5d-5,3.16d-5,4d-5,5d-5,6.3d-5,8d-5, &
  27      &  1d-4,1.25d-4,1.6d-4,2d-4,2.5d-4,3.16d-4,4d-4,5d-4,6.3d-4,8d-4, &
  28      &  1d-3,1.25d-3,1.6d-3,2d-3,2.5d-3,3.16d-3,4d-3,5d-3,6.3d-3,8d-3, &
  29      &  1d-2,1.25d-2,1.6d-2,2d-2,2.5d-2,3.16d-2,4d-2,5d-2,6.3d-2,8d-2, &
  30      &  0.10d0,0.125d0,0.15d0,0.175d0,0.20d0,0.225d0,0.25d0,0.275d0, &
  31      &  0.30d0,0.325d0,0.35d0,0.375d0,0.40d0,0.425d0,0.45d0,0.475d0, &
  32      &  0.50d0,0.525d0,0.55d0,0.575d0,0.60d0,0.625d0,0.65d0,0.675d0, &
  33      &  0.70d0,0.725d0,0.75d0,0.775d0,0.80d0,0.825d0,0.85d0,0.875d0, &
  34      &  0.9d0,0.920d0,0.94d0,0.960d0,0.98d0,1d0, &
  35      &  0.3d0,0.31d0,0.35d0,0.375d0,0.4d0,0.45d0,0.5d0,0.51d0,0.525d0, &
  36      &  0.55d0,0.575d0,0.6d0,0.65d0,0.7d0,0.75d0,0.8d0,0.85d0,0.9d0, &
  37      &  1d0,1.25d0,1.6d0,2d0,2.5d0,3.16d0,4d0,5d0,6.3d0,8d0, &
  38      &  1d1,1.25d1,1.6d1,2d1,2.5d1,3.16d1,4d1,5d1,6.3d1,8d1, &
  39      &  1d2,1.25d2,1.6d2,2d2,2.5d2,3.16d2,4d2,5d2,6.3d2,8d2, &
  40      &  1d3,1.25d3,1.6d3,2d3,2.5d3,3.16d3,4d3,5d3,6.3d3,8d3, &
  41      &  1d4,1.25d4,1.6d4,2d4,2.5d4,3.16d4,4d4,5d4,6.3d4,8d4, &
  42      &  1d5,1.25d5,1.6d5,2d5,2.5d5,3.16d5,4d5,5d5,6.3d5,8d5, &
  43      &  1d6,1.25d6,1.6d6,2d6,2.5d6,3.16d6,4d6,5d6,6.3d6,8d6, &
  44      &  1d7,1.25d7,1.6d7,2d7,2.5d7,3.16d7,4d7,5d7,6.3d7,8d7,1d8/
  45
  46       save
  47       x=xin
  48       q2=qin*qin
  49       call getnset(iset)
  50       call getnmem(iset,imem)
  51        upv =  LHA_GJR08(x,Q2,grid,fgrid,ng,1,imem)
  52       dnv =  LHA_GJR08(x,Q2,grid,fgrid,ng,2,imem)
  53       usea = LHA_GJR08(x,Q2,grid,fgrid,ng,-1,imem)
  54       dsea = LHA_GJR08(x,Q2,grid,fgrid,ng,-2,imem)
  55       str =  LHA_GJR08(x,Q2,grid,fgrid,ng,-3,imem)
  56       glu =  LHA_GJR08(x,Q2,grid,fgrid,ng,0,imem)
  57       pdf(-6) = 0.0d0
  58        pdf(6) = 0.0d0
  59       pdf(-5) = 0.0d0
  60        pdf(5) = 0.0d0
  61       pdf(-4) = 0.0d0
  62        pdf(4) = 0.0d0
  63       pdf(-3) = str
  64        pdf(3) = str
  65       pdf(-2) = usea
  66        pdf(2) = upv+usea
  67       pdf(-1) = dsea
  68        pdf(1) = dnv+dsea
  69        pdf(0) = glu
  70       if(name(iset)(1:7).eq.'GJR08VF'.or. &
  71      &   name(iset)(1:7).eq.'GJR08LO') then
  72         chm =  LHA_GJR08(x,Q2,grid,fgrid,ng,-4,imem)
  73        bot =  LHA_GJR08(x,Q2,grid,fgrid,ng,-5,imem)
  74        pdf(-5) = bot
  75         pdf(5) = bot
  76        pdf(-4) = chm
  77         pdf(4) = chm
  78       endif
  79       return
  80 !cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
  81 !
  82       entry GJRgetgrid(nset,ngridx,ngridq,gridx,gridq)
  83
  84       ngridx=118
  85       do jx=1,118
  86           gridx(jx)=grid(jx)
  87       enddo
  88
  89       ngridq=99
  90       do jq=1,99
  91           gridq(jq)=grid(118+jq)
  92       enddo
  93
  94       return
  95
  96 !cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
  97       entry GJRread(nset)
  98 !
  99       call getnmem(nset,imem)
 100       read(1,*)nmem(nset),ndef(nset)
 101       lstart = -3
 102       if(name(nset)(1:7).eq.'GJR08VF') lstart=-5
 103       do i=0,nmem(nset)
 104       !do ii=0,nmem(nset)
 105 !       i = 0
 106        do j=1,118
 107         do k=1,99
 108          if(name(nset)(1:7).eq.'GJR08VF'.or. &
 109      &      name(nset)(1:7).eq.'GJR08LO') then
 110            read(1,*) fgrid(j,k,-5,i),fgrid(j,k,-4,i), &
 111      &             fgrid(j,k,-3,i),fgrid(j,k,-2,i),fgrid(j,k,-1,i), &
 112      &             fgrid(j,k,0,i),fgrid(j,k,1,i),fgrid(j,k,2,i), &
 113      &             fgrid(j,k,3,i)
 114          else
 115           read(1,*) fgrid(j,k,-3,i),fgrid(j,k,-2,i),fgrid(j,k,-1,i), &
 116      &             fgrid(j,k,0,i),fgrid(j,k,1,i),fgrid(j,k,2,i), &
 117      &              fgrid(j,k,3,i)
 118          endif
 119         if(name(nset)(1:7).ne.'GJR08LO') then
 120           do  l=-lstart,3
 121             if (grid(118+k) < 0.5d0) then
 122               fgrid(j,k,l,i)=0d0
 123               fgrid(j,k,l,i)=0d0/fgrid(j,k,l,i)
 124             endif
 125           enddo
 126          endif
 127         enddo
 128        enddo
 129       enddo
 130
 131       return
 132 !
 133
 134 !cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
 135       entry GJRalfa(alfas,qalfa)
 136       call getnset(iset)
 137       call getnmem(iset,imem)
 138       arg(1) = 1d-9
 139       arg(2) = qalfa*qalfa
 140 !      imem = 0
 141       alfas = lha_dfint(2,arg,ng,grid,fgrid(1,1,3,imem))
 142       return
 143 !
 144 !cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
 145       entry GJRinit(Eorder,Q2fit)
 146       return
 147 !
 148 !cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
 149       entry GJRpdf(mem)
 150       imem = mem
 151       call getnset(iset)
 152        call setnmem(iset,imem)
 153       return
 154 !
 155  1000 format(5e13.5)
 156       end
 157 !
 158 !cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
 159
 160       double precision function LHA_GJR08(x,Q2,grid,fgrid,ng,n,set)
 161       implicit none
 162       integer ng(2),n,set
 163       double precision grid(217),arg(2),x,Q2
 164       double precision lha_dfint
 165       double precision fgrid(118,99,-5:3,0:26)
 166 !      common/fgridc/fgrid
 167       arg(1) = x
 168       arg(2) = Q2
 169        LHA_GJR08 = lha_dfint(2,arg,ng,grid,fgrid(1,1,n,set))
 170       end
 171
 172 !cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
 173
 174
 175
 176 !! CERNLIB E104 modified to be used with GJR08 GRIDS:
 177 !! Name changed from fint to dfint.
 178 !! Name changed from dfint to lha_dfint.
 179 !! Real variables changed to double precision.
 180 !! External references to CERNLIB (error handling) routines removed.
 181           DOUBLE PRECISION FUNCTION LHA_DFINT(NARG,ARG,NENT,ENT,TABLE)
 182           INTEGER   NENT(9), INDEX(32)
 183           DOUBLE PRECISION ARG(9),   ENT(9),   TABLE(9), WEIGHT(32)
 184           LHA_DFINT  =  0d0
 185           IF(NARG .LT. 1  .OR.  NARG .GT. 5)  GOTO 300
 186           LMAX      =  0
 187           ISTEP     =  1
 188           KNOTS     =  1
 189           INDEX(1)  =  1
 190           WEIGHT(1) =  1d0
 191           DO 100    N  =  1, NARG
 192              X     =  ARG(N)
 193              NDIM  =  NENT(N)
 194              LOCA  =  LMAX
 195              LMIN  =  LMAX + 1
 196              LMAX  =  LMAX + NDIM
 197              IF(NDIM .GT. 2)  GOTO 10
 198              IF(NDIM .EQ. 1)  GOTO 100
 199              H  =  X - ENT(LMIN)
 200              IF(H .EQ. 0.)  GOTO 90
 201              ISHIFT  =  ISTEP
 202              IF(X-ENT(LMIN+1) .EQ. 0d0)  GOTO 21
 203              ISHIFT  =  0
 204              ETA     =  H / (ENT(LMIN+1) - ENT(LMIN))
 205              GOTO 30
 206   10         LOCB  =  LMAX + 1
 207   11         LOCC  =  (LOCA+LOCB) / 2
 208 !             IF(X-ENT(LOCC))  12, 20, 13
 209              IF(X-ENT(LOCC).lt.0)  goto 12
 210              IF(X-ENT(LOCC).eq.0)  goto 20
 211              IF(X-ENT(LOCC).gt.0)  goto 13
 212   12         LOCB  =  LOCC
 213              GOTO 14
 214   13         LOCA  =  LOCC
 215   14         IF(LOCB-LOCA .GT. 1)  GOTO 11
 216              LOCA    =  MIN0( MAX0(LOCA,LMIN), LMAX-1 )
 217              ISHIFT  =  (LOCA - LMIN) * ISTEP
 218              ETA     =  (X - ENT(LOCA)) / (ENT(LOCA+1) - ENT(LOCA))
 219              GOTO 30
 220   20         ISHIFT  =  (LOCC - LMIN) * ISTEP
 221   21         DO 22  K  =  1, KNOTS
 222                 INDEX(K)  =  INDEX(K) + ISHIFT
 223   22            CONTINUE
 224              GOTO 90
 225   30         DO 31  K  =  1, KNOTS
 226                 INDEX(K)         =  INDEX(K) + ISHIFT
 227                 INDEX(K+KNOTS)   =  INDEX(K) + ISTEP
 228                 WEIGHT(K+KNOTS)  =  WEIGHT(K) * ETA
 229                 WEIGHT(K)        =  WEIGHT(K) - WEIGHT(K+KNOTS)
 230   31            CONTINUE
 231              KNOTS  =  2*KNOTS
 232   90         ISTEP  =  ISTEP * NDIM
 233  100         CONTINUE
 234           DO 200    K  =  1, KNOTS
 235              I  =  INDEX(K)
 236              LHA_DFINT  =  LHA_DFINT + WEIGHT(K) * TABLE(I)
 237  200         CONTINUE
 238           RETURN
 239  300      WRITE(*,1000) NARG
 240           STOP
 241 1000      FORMAT( 7X, 24HFUNCTION DFINT... NARG =,I6, &
 242      &              17H NOT WITHIN RANGE)
 243           END
 244