************************************************************************ * * Tabulation.f: * * Example program used for the LH benchmark. * ************************************************************************ program Tabulation * implicit none * integer ilha double precision Q0,Q double precision Q02,Q2 double precision AlphaQCD,AlphaQED double precision xPDF,xgamma double precision eps double precision xlha(11) parameter(eps=1d-10) data xlha / 1d-7, 1d-6, 1d-5, 1d-4, 1d-3, 1d-2, 1 1d-1, 3d-1, 5d-1, 7d-1, 9d-1 / * * Some examples ... * c call SetFastEvolution(.true.) c call LockGrids(.true.) c call EnableEvolutionOperator(.true.) c call SetFFNS(4) c call SetTheory("QavDS") c call SetTheory("QUniD") c call SetPerturbativeOrder(0) c call SetPDFEvolution("exactalpha") c call SetPDFSet("NNPDF23_nlo_as_0119_qed") c call SetPDFSet("MRST2004qed") call SetPDFSet("external") c call SetNumberOfGrids(1) c call SetGridParameters(1,30,3,1d-5) c call SetGridParameters(2,30,3,2d-1) c call SetGridParameters(3,30,3,8d-1) c call SetPDFSet("NNPDF30_nnlo_as_0118") c call SetAlphaQCDRef(0.118d0,91.2d0) c call SetAlphaEvolution("expanded") c call SetPDFEvolution("expandalpha") c call SetPoleMasses(1.275d0,4.18d0,173.03d0) c call SetMaxFlavourPDFs(5) c call SetMaxFlavourAlpha(5) * * Initializes integrals on the grids * call InitializeAPFEL * * Evolve PDFs on the grids * write(6,*) "Enter initial and final scale in GeV^2" read(5,*) Q02,Q2 * Q0 = dsqrt(Q02) - eps Q = dsqrt(Q2) call EvolveAPFEL(Q0,Q) * * Tabulate PDFs for the LHA x values * write(6,*) "alpha_QCD(mu2F) =",AlphaQCD(Q) write(6,*) "alpha_QED(mu2F) =",AlphaQED(Q) write(6,*) " " * write(6,*) "Standard evolution:" write(6,'(a5,2a12,a14,a10,2a12)') "x", 1 "u-ubar","d-dbar","2(ubr+dbr)","c+cbar","gluon","photon" do ilha=3,11 write(6,'(es7.1,6es12.4)') 1 xlha(ilha), 2 xPDF(2,xlha(ilha)) - xPDF(-2,xlha(ilha)), 3 xPDF(1,xlha(ilha)) - xPDF(-1,xlha(ilha)), 4 2d0 * ( xPDF(-1,xlha(ilha)) + xPDF(-2,xlha(ilha)) ), 5 xPDF(4,xlha(ilha)) + xPDF(-4,xlha(ilha)), 6 xPDF(0,xlha(ilha)), 7 xgamma(xlha(ilha)) enddo write(*,*) " " * end * ************************************************************************ * * External PDF set that the user can feed to APFEL to use as initial * scale set. * The external set must have the name "ExternalPDFSetAPFEL" and the * following structure: * * call ExternalSetAPFEL(x,xf) * * where "x" is a double corresponding to the Bjorken's variable and * xf(-6:7) is an array of double that corresponds to the set of * distributions in the physical basis (tbar, bbar, ..., g, ..., b, t) * and where the component xf(7) is the photon distribution. * This subroutine is called only if the user sets: * * call SetPDFSet("external") * ************************************************************************ subroutine ExternalSetAPFEL(x,Q,xf) * implicit none ** * Input Variables * double precision x, Q ** * Internal Variables * integer ipdf double precision xuv,xdv,xg,xdbar,xubar,xs,xsbar double precision N_uv,auv,buv,N_dv,adv,bdv,N_g,ag double precision bg,N_db,adb,bdb,fs ** * Output Variables * double precision xf(-6:7) * * Parameters of the User defined PDFs * N_uv = 5.107200d0 auv = 0.8d0 buv = 3d0 N_dv = 3.064320d0 adv = 0.8d0 bdv = 4d0 N_g = 1.7d0 ag = -0.1d0 bg = 5d0 N_db = 0.1939875d0 adb = -0.1d0 bdb = 6d0 fs = 0.2d0 * * User defined PDFs * xuv = N_uv * x**auv * ( 1d0 - x )**buv xdv = N_dv * x**adv * ( 1d0 - x )**bdv xg = N_g * x**ag * ( 1d0 - x )**bg xdbar = N_db * x**adb * ( 1d0 - x )**bdb xubar = xdbar * ( 1d0 - x ) xs = fs * ( xdbar + xubar ) xsbar = xs * * Initialize PDFs to zero * do ipdf=-6,7 xf(ipdf) = 0d0 enddo * if(x.gt.1d0) return * xf(3) = xs xf(2) = xuv + xubar xf(1) = xdv + xdbar xf(0) = xg xf(-1) = xdbar xf(-2) = xubar xf(-3) = xsbar * return end