************************************************************************ * * DISObservables.f: * * Test program for the DIS observables. * ************************************************************************ program DISObservables * implicit none ** * Variables * integer ilha double precision Q0,Q double precision Q02,Q2 double precision AlphaQCD,AlphaQED double precision F2light,F2charm,F2bottom,F2total double precision FLlight,FLcharm,FLbottom,FLtotal double precision F3light,F3charm,F3bottom,F3total double precision StructureFunctionxQ double precision eps double precision xlha(11) character*2 proc 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 / * * Settings * proc = "NC" c call SetMassScheme("ZM-VFNS") call SetMassScheme("FONLL-C") call SetProcessDIS(proc) c call SetQLimits(1.4d0,250d0) c call EnableDynamicalScaleVariations(.true.) c call SetRenQRatio(0.5d0) c call SetFacQRatio(0.5d0) c call SetEWCouplings(1d0,1d0,1d0,1d0) c call SetPropagatorCorrection(0.1d0) c call SetPolarizationDIS(0d0) c call SetProjectileDIS("electron") c call SetTargetDIS("proton") c call EnableTargetMassCorrections(.false.) c call EnableDampingFONLL(.true.) c call SetFastEvolution(.true.) c call LockGrids(.true.) c call EnableEvolutionOperator(.true.) c call SetFFNS(3) c call SetTheory("QUniD") c call EnableNLOQEDCorrections(.true.) c call EnableLeptonEvolution(.true.) c call SetTauMass(1d10) c call SetPerturbativeOrder(0) c call SetPDFEvolution("exactalpha") c call SetPDFSet("NNPDF23_nlo_as_0119_qed") c call SetPDFSet("MRST2004qed") 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(dsqrt(2d0),4.5d0,100000d0) c call SetMaxFlavourPDFs(5) c call SetMaxFlavourAlpha(5) c call SetQGridParameters(100,2) * * Initializes integrals on the grids * call InitializeAPFEL_DIS * * 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 ComputeStructureFunctionsAPFEL(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,'(a5,4a12)') "x","F2light","F2charm","F2bottom","F2total" do ilha=3,11 write(6,'(es7.1,9es12.4)') 1 xlha(ilha), 2 F2light(xlha(ilha)), 3 F2charm(xlha(ilha)), 4 F2bottom(xlha(ilha)), 5 F2total(xlha(ilha)) enddo write(*,*) " " * write(6,'(a5,4a12)') "x","FLlight","FLcharm","FLbottom","FLtotal" do ilha=3,11 write(6,'(es7.1,9es12.4)') 1 xlha(ilha), 2 FLlight(xlha(ilha)), 3 FLcharm(xlha(ilha)), 4 FLbottom(xlha(ilha)), 5 FLtotal(xlha(ilha)) enddo write(*,*) " " * write(6,'(a5,4a12)') "x","F3light","F3charm","F3bottom","F3total" do ilha=3,11 write(6,'(es7.1,9es12.4)') 1 xlha(ilha), 2 F3light(xlha(ilha)), 3 F3charm(xlha(ilha)), 4 F3bottom(xlha(ilha)), 5 F3total(xlha(ilha)) enddo write(*,*) " " * * Cache Structure functions * call CacheStructureFunctionsAPFEL(Q0) * write(6,'(a5,4a12)') "x","F2light","F2charm","F2bottom","F2total" do ilha=3,11 write(6,'(es7.1,9es12.4)') 1 xlha(ilha), 2 StructureFunctionxQ(proc,"F2","light",xlha(ilha),Q), 3 StructureFunctionxQ(proc,"F2","charm",xlha(ilha),Q), 4 StructureFunctionxQ(proc,"F2","bottom",xlha(ilha),Q), 5 StructureFunctionxQ(proc,"F2","total",xlha(ilha),Q) enddo write(*,*) " " * write(6,'(a5,4a12)') "x","FLlight","FLcharm","FLbottom","FLtotal" do ilha=3,11 write(6,'(es7.1,9es12.4)') 1 xlha(ilha), 2 StructureFunctionxQ(proc,"FL","light",xlha(ilha),Q), 3 StructureFunctionxQ(proc,"FL","charm",xlha(ilha),Q), 4 StructureFunctionxQ(proc,"FL","bottom",xlha(ilha),Q), 5 StructureFunctionxQ(proc,"FL","total",xlha(ilha),Q) enddo write(*,*) " " * write(6,'(a5,4a12)') "x","F3light","F3charm","F3bottom","F3total" do ilha=3,11 write(6,'(es7.1,9es12.4)') 1 xlha(ilha), 2 StructureFunctionxQ(proc,"F3","light",xlha(ilha),Q), 3 StructureFunctionxQ(proc,"F3","charm",xlha(ilha),Q), 4 StructureFunctionxQ(proc,"F3","bottom",xlha(ilha),Q), 5 StructureFunctionxQ(proc,"F3","total",xlha(ilha),Q) enddo write(*,*) " " * end