Example evaluation for package of TMD packages

SetOptions[$Output[[1]],FormatType->StandardForm]

AlphaStrong

SetDirectory["/home/vla18041/WorkingFiles/TMD/Code/FinalCode"];
<<AlphaStrong`

?as

Plot[{as[1,μ],as[2,μ],as[3,μ]},{μ,1.,10},PlotStyle->{Black,Red,Blue}]

?asIt

Plot[{asIt[1,μ],asIt[2,μ],asIt[3,μ]},{μ,1.,10},PlotStyle->{Black,Red,Blue}]

Plot[{asIt[2,μ],as[2,μ]},{μ,1.,10},PlotStyle->{Black,Red,Blue}]

TMDevolution

SetDirectory["/home/vla18041/WorkingFiles/TMD/Code/FinalCode"];
<<TMDevolution`

Anomalous Dimensions

?Γcusp

Γcusp["q",3,Nf]

{γV["g",1,3],γV["g",1,4],γV["g",1,5]}//N//InputForm

γV["g",1,Nf]

?γV

?ddTMD

?dd

ddTMD["q",3,2,Nf]
dd[3,2,Nf]

?ΓLogevolutor

?AsLog

Evolution integrals

?Γevolutor

It is not defined for μ0<Λ. And, thus, gives mistakes.

Plot[{Γevolutor["g",1,10,μ0],Γevolutor["g",2,1,μ0],Γevolutor["g",3,10,μ0]},{μ0,0.3,4}]//Timing

Γevolutor["g",3,10,0.02]

First time evaluation takes some time for tabulating intrinsic variables

Plot[{Γevolutor["g",1,10,μ0],Γevolutor["g",2,1,μ0],Γevolutor["g",3,10,μ0]},{μ0,0.5,4}]//Timing

But next evaluations takes much-much less time.

Plot[{Γevolutor["g",1,10,μ0],Γevolutor["g",2,1,μ0],Γevolutor["g",3,10,μ0]},{μ0,0.5,4}]//Timing

?DTMD

Plot[{DTMD["g",1,10,bt,(C0/bt)],DTMD["g",2,10,bt,(C0/bt)],DTMD["g",3,10,bt,(C0/bt)]},{bt,0,10},PlotStyle->{Black,Red,Blue}]//Timing

Plot[{DTMD["g",1,10,bt,(C0/(bstar[bt,1]))],DTMD["g",2,10,bt,(C0/(bstar[bt,1]))],DTMD["g",3,10,bt,(C0/(bstar[bt,1]))]},{bt,0,10},PlotStyle->{Black,Red,Blue}]//Timing

?ΓLogevolutor

Plot[{Γevolutor["g",1,10,(C0/(bstar[bt,1]))],Γevolutor["g",2,10,(C0/(bstar[bt,1]))],Γevolutor["g",3,10,(C0/(bstar[bt,1]))]},{bt,0,10},PlotStyle->{Black,Red,Blue}]//Timing

?γevolutor

Plot[{γevolutor["q",1,10,(C0/(bstar[bt,1]))],γevolutor["q",2,10,(C0/(bstar[bt,1]))],γevolutor["q",3,10,(C0/(bstar[bt,1]))]},{bt,0,10},PlotStyle->{Black,Red,Blue}]//Timing

?DTMD

?Dboundary

Evolution Kernel R

?RTMD

RTMD["q",3,5,10,9.3,0.9]

Timing[Plot[RTMD["q",1,bt,(10.^2),10.,(C0/(bstar[bt,1.]))],{bt,0,5}]]

Plot[{RTMD["q",1,bt,(10.^2),10.,(C0/(bstar[bt,1.]))],RTMD["q",2,bt,(10.^2),10.,(C0/(bstar[bt,1.]))],RTMD["q",3,bt,(10.^2),10.,(C0/(bstar[bt,1.]))]},{bt,0,5},PlotStyle->{Black,Red,Blue}]

Plot[{RTMD["q",2,bt,(10.^2),10.,(C0/(bstar[bt,0.5]))],RTMD["q",2,bt,(10.^2),10.,(C0/(bstar[bt,1.]))],RTMD["q",2,bt,(10.^2),10.,(C0/(bstar[bt,1.5]))]},{bt,0,5},PlotStyle->{Black,Red,Blue}]

?RTMD2

Plot[{RTMD2["q",1,bt,(10.^2),10.,1.2,(C0/(bstar[bt,1.])),(C0/(bstar[bt,1.]))],RTMD2["q",1,bt,(10.^2),10.,1.2,(C0/(bstar[bt,1.])),(C0/(bstar[bt,1.]))]},{bt,0,4},PlotStyle->{Black,Red},Filling->{1->{2}}]

TMD at small - bt

DO NOT Forget to load MSTW!

SetDirectory["/media/USB DISK/MSTW"];
<<mstwpdf.m

Timing[ReadPDFGrid["/media/USB DISK/MSTW/Grids/mstw2008nnlo",0]]

SetDirectory["/home/vla18041/WorkingFiles/TMD/Code/FinalCode"];
<<TMDmatching`

Looking at TMD

CxF[0,1,"direct",(10^(-3)),1.,10.]
CxF[1,1,"direct",(10^(-3)),1.,10.]
CxF[2,1,"direct",(10^(-3)),1.,10.]

CxF[0,1,"withMIX",(10^(-3)),1.,10.]
CxF[1,1,"withMIX",(10^(-3)),1.,10.]
CxF[2,1,"withMIX",(10^(-3)),1.,10.]

CxF[0,0,"direct",(10^(-3)),1.,10.]
CxF[1,0,"direct",(10^(-3)),1.,10.]
CxF[2,0,"direct",(10^(-3)),1.,10.]

CxF[0,0,"withMIX",(10^(-3)),1.,10.]
CxF[1,0,"withMIX",(10^(-3)),1.,10.]
CxF[2,0,"withMIX",(10^(-3)),1.,10.]

NIntegrate[(-2 (-1+z) z+Lμ (-1+2 z-2 z^2))xf[0,x/z,10,0]/.{x->(10.^(-3)),Lμ->%51},{z,(10^(-3)),1},PrecisionGoal->5]

(-2 (-1+z) z+Lμ (-1+2 z-2 z^2))-(-Lμ+2(1+Lμ)z)//Simplify

Plot[(1/z)xf[0,(0.001/z),10,1],{z,0.001,1},PlotRange->All]

?CxF

CxF[2,1,"withMIX",0.3,0.2,10.]
CxF[2,0,"withMIX",0.3,0.2,10.]

Since evaluation of convolutions is time-consuming, we suggest to evaluate point-by-point.

<<TMDevolution`

Timing[A1=Table[{x,x CxF[2,1,"withMIX",x,1,(C0/(bstar[1,1]))]},{x,(10^(-6)),0.3,0.001}];]
Timing[A2=Table[{x,x CxF[2,2,"withMIX",x,1,(C0/(bstar[1,1]))]},{x,(10^(-6)),0.3,0.001}];]
Timing[A3=Table[{x,x CxF[2,3,"withMIX",x,1,(C0/(bstar[1,1]))]},{x,(10^(-6)),0.3,0.001}];]

ListPlot[{A1,A2,A3},PlotStyle->{Black,Red,Blue},Joined->True]

The effect of mixing is very small, by evaluating it takes tons of time. If “withMIX” is not set, the mixing is not evaluated.

Timing[A1=Table[{x,x CxF[2,1,"direct",x,1,(C0/(bstar[1,1]))]},{x,(10^(-6)),0.3,0.001}];]
Timing[A2=Table[{x,x CxF[2,2,"direct",x,1,(C0/(bstar[1,1]))]},{x,(10^(-6)),0.3,0.001}];]
Timing[A3=Table[{x,x CxF[2,3,"direct",x,1,(C0/(bstar[1,1]))]},{x,(10^(-6)),0.3,0.001}];]

ListPlot[{A1,A2,A3},PlotStyle->{Black,Red,Blue},Joined->True]

Show[%15,%19]

Timing[A1=Table[{bt,CxF[1,1,0.01,bt,(C0/(bstar[bt,1]))]},{bt,(10^(-4)),4,0.1}];]
Timing[A2=Table[{bt,CxF[1,2,0.01,bt,(C0/(bstar[bt,1]))]},{bt,(10^(-4)),4,0.1}];]
Timing[A3=Table[{bt,CxF[1,3,0.01,bt,(C0/(bstar[bt,1]))]},{bt,(10^(-4)),4,0.1}];]

ListPlot[{A1,A2,A3},PlotStyle->{Black,Red,Blue},Joined->True]

Timing[A1=Table[{bt,CxF[0,0,0.01,bt,(C0/(bstar[bt,1]))]},{bt,(10^(-4)),4,0.1}];]
Timing[A2=Table[{bt,CxF[1,0,0.01,bt,(C0/(bstar[bt,1]))]},{bt,(10^(-4)),4,0.1}];]
Timing[A3=Table[{bt,CxF[2,0,0.01,bt,(C0/(bstar[bt,1]))]},{bt,(10^(-4)),4,0.1}];]

ListPlot[{A1,A2,A3},PlotStyle->{Black,Red,Blue},Joined->True]

Timing[A1=Table[{bt,CxF[0,2,0.01,bt,(C0/(bstar[bt,1]))]},{bt,(10^(-4)),4,0.1}];]
Timing[A2=Table[{bt,CxF[1,2,0.01,bt,(C0/(bstar[bt,1]))]},{bt,(10^(-4)),4,0.1}];]
Timing[A3=Table[{bt,CxF[2,2,0.01,bt,(C0/(bstar[bt,1]))]},{bt,(10^(-4)),4,0.1}];]

ListPlot[{A1,A2,A3},PlotStyle->{Black,Red,Blue},Joined->True]

A1=Table[{x,bt,x CxF[1,2,x,bt,(C0/(bstar[bt,1]))]},{x,(10^(-3)),1,0.025},{bt,(10^(-4)),4,0.3}];//Timing
A2=Table[{x,bt,x CxF[2,2,x,bt,(C0/(bstar[bt,1]))]},{x,(10^(-3)),1,0.025},{bt,(10^(-4)),4,0.3}];//Timing

ListPlot3D[{Flatten[A1,1]},AxesLabel->{"x","bt"},PlotRange->{0,2},InterpolationOrder->3]
ListPlot3D[{Flatten[A2,1]},AxesLabel->{"x","bt"},PlotRange->{0,2},InterpolationOrder->3]

ListPointPlot3D[{Flatten[A1,1],Flatten[A2,1]},AxesLabel->{"x","bt"}]

Looking at TMD renormalon correction

?CxFrenormalon

Timing[R1=Table[{bt,CxFrenormalon[1,0.01,bt,(C0/(bstar[bt,1]))]},{bt,(10^(-4)),4,0.1}];]

Timing[A1=Table[{bt,CxF[1,1,0.01,bt,(C0/(bstar[bt,1]))]},{bt,(10^(-4)),4,0.1}];]

ListPlot[{A1,R1},PlotStyle->{Black,Red,Blue},Joined->True]

Timing[A1=Table[{x,-1.2CxFrenormalon[1,x,1,(C0/(bstar[1,1]))]},{x,(10^(-6)),1,0.001}];]

ListPlot[A1,PlotStyle->{Black,Red,Blue},Joined->True]

?CxFexp

Timing[A1=Table[{x,x CxF[1,2,"withMIX",x,1,(C0/(bstar[1,1]))]},{x,(10^(-6)),1,0.01}];]
Timing[A2=Table[{x,x CxF[2,2,"withMIX",x,1,(C0/(bstar[1,1]))]},{x,(10^(-6)),1,0.01}];]

Timing[B1=Table[{x,x CxFexp[1,2,"withMIX",x,1,(C0/(bstar[1,1])),0.15815688133371278`]},{x,(10^(-6)),1,0.01}];]
Timing[B2=Table[{x,x CxFexp[2,2,"withMIX",x,1,(C0/(bstar[1,1])),0.15815688133371278`]},{x,(10^(-6)),1,0.01}];]

ListPlot[{A1,B1},PlotStyle->{Red,Blue},Joined->True]
ListPlot[{A2,B2},PlotStyle->{Red,Blue},Joined->True]

Timing[A1=Table[{bt,CxFexp[1,2,0.1,bt,(C0/(bstar[bt,1])),0.158]},{bt,(10^(-3)),6,0.1}];]
Timing[A2=Table[{bt,CxFexp[2,2,0.1,bt,(C0/(bstar[bt,1])),0.158]},{bt,(10^(-3)),6,0.1}];]

ListPlot[{A1,A2},PlotStyle->{Red,Blue},Joined->True]

Timing[A1=Table[{bt,CxF[2,2,0.1,bt,(C0/(bstar[bt,1]))]},{bt,(10^(-3)),6,0.1}];]
Timing[A2=Table[{bt,CxFexp[2,2,0.1,bt,(C0/(bstar[bt,1])),0.0]},{bt,(10^(-3)),6,0.1}];]

Total[A1-A2]

Timing[A1=Table[{bt,CxFexp[2,2,0.1,bt,(C0/(bstar[bt,1])),0.0]},{bt,(10^(-3)),6,0.1}];]
Timing[A2=Table[{bt,CxFexp[2,2,0.1,bt,(C0/(bstar[bt,1])),0.05]},{bt,(10^(-3)),6,0.1}];]
Timing[A3=Table[{bt,CxFexp[2,2,0.1,bt,(C0/(bstar[bt,1])),0.1]},{bt,(10^(-3)),6,0.1}];]
Timing[A4=Table[{bt,CxFexp[2,2,0.1,bt,(C0/(bstar[bt,1])),0.15]},{bt,(10^(-3)),6,0.1}];]
Timing[A5=Table[{bt,CxFexp[2,2,0.1,bt,(C0/(bstar[bt,1])),0.2]},{bt,(10^(-3)),6,0.1}];]

ListPlot[{A1,A2,A3,A4,A5},Joined->True]

?CxFrenormalonExp

Timing[B1=Table[{bt,CxFrenormalonExp[2,0.1,bt,(C0/(bstar[bt,1])),-0.1,0.0]},{bt,(10^(-3)),6,0.1}];]
Timing[B2=Table[{bt,CxFrenormalonExp[2,0.1,bt,(C0/(bstar[bt,1])),-0.1,0.01]},{bt,(10^(-3)),6,0.1}];]
Timing[B3=Table[{bt,CxFrenormalonExp[2,0.1,bt,(C0/(bstar[bt,1])),-0.1,0.02]},{bt,(10^(-3)),6,0.1}];]
Timing[B4=Table[{bt,CxFrenormalonExp[2,0.1,bt,(C0/(bstar[bt,1])),-0.1,0.03]},{bt,(10^(-3)),6,0.1}];]
Timing[B5=Table[{bt,CxFrenormalonExp[2,0.1,bt,(C0/(bstar[bt,1])),-0.1,0.04]},{bt,(10^(-3)),6,0.1}];]

ListPlot[{B1,B2,B3,B4,B5}]

Timing[A1=Table[{bt,0.1(CxFexp[2,2,0.1,bt,(C0/(bstar[bt,1])),0.0]+CxFrenormalonExp[2,0.1,bt,(C0/(bstar[bt,1])),-0.1,0.0])},{bt,(10^(-3)),6,0.1}];]
Timing[A2=Table[{bt,0.1(CxFexp[2,2,0.1,bt,(C0/(bstar[bt,1])),0.03]+CxFrenormalonExp[2,0.1,bt,(C0/(bstar[bt,1])),-0.1,0.03])},{bt,(10^(-3)),6,0.1}];]
Timing[A3=Table[{bt,0.1(CxFexp[2,2,0.1,bt,(C0/(bstar[bt,1])),0.06]+CxFrenormalonExp[2,0.1,bt,(C0/(bstar[bt,1])),-0.1,0.06])},{bt,(10^(-3)),6,0.1}];]
Timing[A4=Table[{bt,0.1(CxFexp[2,2,0.1,bt,(C0/(bstar[bt,1])),0.09]+CxFrenormalonExp[2,0.1,bt,(C0/(bstar[bt,1])),-0.1,0.09])},{bt,(10^(-3)),6,0.1}];]
Timing[A5=Table[{bt,0.1(CxFexp[2,2,0.1,bt,(C0/(bstar[bt,1])),0.12]+CxFrenormalonExp[2,0.1,bt,(C0/(bstar[bt,1])),-0.1,0.12])},{bt,(10^(-3)),6,0.1}];]

Timing[B1=Table[{bt,0.9(CxFexp[2,2,0.9,bt,(C0/(bstar[bt,1])),0.0]+CxFrenormalonExp[2,0.9,bt,(C0/(bstar[bt,1])),-0.1,0.0])},{bt,(10^(-3)),6,0.1}];]
Timing[B2=Table[{bt,0.9(CxFexp[2,2,0.9,bt,(C0/(bstar[bt,1])),0.03]+CxFrenormalonExp[2,0.9,bt,(C0/(bstar[bt,1])),-0.1,0.03])},{bt,(10^(-3)),6,0.1}];]
Timing[B3=Table[{bt,0.9(CxFexp[2,2,0.9,bt,(C0/(bstar[bt,1])),0.06]+CxFrenormalonExp[2,0.9,bt,(C0/(bstar[bt,1])),-0.1,0.06])},{bt,(10^(-3)),6,0.1}];]
Timing[B4=Table[{bt,0.9(CxFexp[2,2,0.9,bt,(C0/(bstar[bt,1])),0.09]+CxFrenormalonExp[2,0.9,bt,(C0/(bstar[bt,1])),-0.1,0.09])},{bt,(10^(-3)),6,0.1}];]
Timing[B5=Table[{bt,0.9(CxFexp[2,2,0.9,bt,(C0/(bstar[bt,1])),0.12]+CxFrenormalonExp[2,0.9,bt,(C0/(bstar[bt,1])),-0.1,0.12])},{bt,(10^(-3)),6,0.1}];]

Timing[C1=Table[{bt,0.9(CxFexp[2,2,0.01,bt,(C0/(bstar[bt,1])),0.0]+CxFrenormalonExp[2,0.01,bt,(C0/(bstar[bt,1])),-0.1,0.0])},{bt,(10^(-3)),6,0.1}];]
Timing[C2=Table[{bt,0.9(CxFexp[2,2,0.01,bt,(C0/(bstar[bt,1])),0.03]+CxFrenormalonExp[2,0.01,bt,(C0/(bstar[bt,1])),-0.1,0.03])},{bt,(10^(-3)),6,0.1}];]
Timing[C3=Table[{bt,0.9(CxFexp[2,2,0.01,bt,(C0/(bstar[bt,1])),0.06]+CxFrenormalonExp[2,0.01,bt,(C0/(bstar[bt,1])),-0.1,0.06])},{bt,(10^(-3)),6,0.1}];]
Timing[C4=Table[{bt,0.9(CxFexp[2,2,0.01,bt,(C0/(bstar[bt,1])),0.09]+CxFrenormalonExp[2,0.01,bt,(C0/(bstar[bt,1])),-0.1,0.09])},{bt,(10^(-3)),6,0.1}];]
Timing[C5=Table[{bt,0.9(CxFexp[2,2,0.01,bt,(C0/(bstar[bt,1])),0.12]+CxFrenormalonExp[2,0.01,bt,(C0/(bstar[bt,1])),-0.1,0.12])},{bt,(10^(-3)),6,0.1}];]

ListPlot[{B1,B2,B3,B4,B5},Joined->True]
ListPlot[{A1,A2,A3,A4,A5},Joined->True]
ListPlot[{C1,C2,C3,C4,C5},Joined->True]

Timing[B1=Table[{x,x(CxF[1,2,x,1,(C0/(bstar[1,1]))]+1.2CxFrenormalon[2,x,1,(C0/(bstar[1,1]))])},{x,(10^(-3)),1,0.01}];]
Timing[B2=Table[{x,x(CxF[1,2,x,1,(C0/(bstar[1,1]))]+0CxFrenormalon[2,x,1,(C0/(bstar[1,1]))])},{x,(10^(-3)),1,0.01}];]
Timing[B3=Table[{x,x(CxF[1,2,x,1,(C0/(bstar[1,1]))]-1.2CxFrenormalon[2,x,1,(C0/(bstar[1,1]))])},{x,(10^(-3)),1,0.01}];]

ListPlot[{B1,B2,B3},Joined->True,PlotStyle->{{Thickness[0.002],Red},{Thick,Black},{Thickness[0.002],Red}}]

Timing[A1=Table[{bt,(CxFexp[2,2,0.1,bt,(C0/(bstar[bt,1])),0.0]+CxFrenormalonExp[2,0.1,bt,(C0/(bstar[bt,1])),-0.1,0.0])},{bt,(10^(-3)),10,0.1}];]
Timing[A2=Table[{bt,(CxFexp[2,2,0.1,bt,(C0/(bstar[bt,1])),0.03]+CxFrenormalonExp[2,0.1,bt,(C0/(bstar[bt,1])),-0.1,0.03])},{bt,(10^(-3)),10,0.1}];]
Timing[A3=Table[{bt,(CxFexp[2,2,0.1,bt,(C0/(bstar[bt,1])),0.06]+CxFrenormalonExp[2,0.1,bt,(C0/(bstar[bt,1])),-0.1,0.06])},{bt,(10^(-3)),10,0.1}];]
Timing[A4=Table[{bt,(CxFexp[2,2,0.1,bt,(C0/(bstar[bt,1])),0.09]+CxFrenormalonExp[2,0.1,bt,(C0/(bstar[bt,1])),-0.1,0.09])},{bt,(10^(-3)),10,0.1}];]
Timing[A5=Table[{bt,(CxFexp[2,2,0.1,bt,(C0/(bstar[bt,1])),0.12]+CxFrenormalonExp[2,0.1,bt,(C0/(bstar[bt,1])),-0.1,0.12])},{bt,(10^(-3)),10,0.1}];]

ListPlot[{A1,A2,A3,A4,A5},Joined->True]

Timing[A1=Table[{bt,(CxFexp[2,2,0.1,bt,(C0/(bstar[bt,1])),0.0]+CxFrenormalonExp[2,0.1,bt,(C0/(bstar[bt,1])),0.1,0.0])},{bt,(10^(-3)),10,0.1}];]
Timing[A2=Table[{bt,(CxFexp[2,2,0.1,bt,(C0/(bstar[bt,1])),0.03]+CxFrenormalonExp[2,0.1,bt,(C0/(bstar[bt,1])),0.1,0.03])},{bt,(10^(-3)),10,0.1}];]
Timing[A3=Table[{bt,(CxFexp[2,2,0.1,bt,(C0/(bstar[bt,1])),0.06]+CxFrenormalonExp[2,0.1,bt,(C0/(bstar[bt,1])),0.1,0.06])},{bt,(10^(-3)),10,0.1}];]
Timing[A4=Table[{bt,(CxFexp[2,2,0.1,bt,(C0/(bstar[bt,1])),0.09]+CxFrenormalonExp[2,0.1,bt,(C0/(bstar[bt,1])),0.1,0.09])},{bt,(10^(-3)),10,0.1}];]
Timing[A5=Table[{bt,(CxFexp[2,2,0.1,bt,(C0/(bstar[bt,1])),0.12]+CxFrenormalonExp[2,0.1,bt,(C0/(bstar[bt,1])),0.1,0.12])},{bt,(10^(-3)),10,0.1}];]

ListPlot[{A1,A2,A3,A4,A5},Joined->True]

Timing[R0=Table[{bt,(CxFexp[0,2,0.1,bt,(C0/(bstar[bt,1])),0.09]+CxFrenormalonExp[2,0.1,bt,(C0/(bstar[bt,1])),-0.1,0.09])},{bt,(10^(-3)),10,0.1}];]
Timing[R1=Table[{bt,(CxFexp[1,2,0.1,bt,(C0/(bstar[bt,1])),0.09]+CxFrenormalonExp[2,0.1,bt,(C0/(bstar[bt,1])),-0.1,0.09])},{bt,(10^(-3)),10,0.1}];]
Timing[R2=Table[{bt,(CxFexp[2,2,0.1,bt,(C0/(bstar[bt,1])),0.09]+CxFrenormalonExp[2,0.1,bt,(C0/(bstar[bt,1])),-0.1,0.09])},{bt,(10^(-3)),10,0.1}];]

ListPlot[{R0,R1,R2},Joined->True,PlotStyle->{Black,Red,Blue}]

Timing[R0=Table[{bt,(CxFexp[0,2,0.1,bt,2,0.09]+CxFrenormalonExp[2,0.1,bt,2,-0.1,0.09])},{bt,(10^(-3)),10,0.1}];]

Timing[R1=Table[{bt,(CxFexp[1,2,0.1,bt,2,0.09]+CxFrenormalonExp[2,0.1,bt,2,-0.1,0.09])},{bt,(10^(-3)),10,0.1}];]
Timing[R2=Table[{bt,(CxFexp[2,2,0.1,bt,2,0.09]+CxFrenormalonExp[2,0.1,bt,2,-0.1,0.09])},{bt,(10^(-3)),10,0.1}];]

ListPlot[{R0,R1,R2},Joined->True,PlotRange->{-3,20},PlotStyle->{Black,Red,Blue}]

Timing[R0=Table[{bt,(CxFexp[0,2,0.1,bt,(C0/(bstar[bt,1])),0.1]+CxFrenormalonExp[2,0.1,bt,(C0/(bstar[bt,1])),-0.02,0.1])},{bt,(10^(-3)),10,0.1}];]
Timing[R1=Table[{bt,(CxFexp[1,2,0.1,bt,(C0/(bstar[bt,1])),0.1]+CxFrenormalonExp[2,0.1,bt,(C0/(bstar[bt,1])),-0.02,0.1])},{bt,(10^(-3)),10,0.1}];]
Timing[R2=Table[{bt,(CxFexp[2,2,0.1,bt,(C0/(bstar[bt,1])),0.1]+CxFrenormalonExp[2,0.1,bt,(C0/(bstar[bt,1])),-0.02,0.1])},{bt,(10^(-3)),10,0.1}];]

ListPlot[{R0,R1,R2},Joined->True,PlotRange->{-3,10},PlotStyle->{Black,Red,Blue}]

?CxF

?CxFrenormalonExp

(CxFexp[1,2,0.8,1,10.,0.1]+CxFrenormalonExp[2,0.8,1.,1.,0.3,0.1])

This is model for u quark from the paper

Needs["PlotLegends`"]

MODELu[n_,x_,bt_,μ_,gb_,gq_]:=CxFexp[n,2,x,bt,μ,gb]+CxFrenormalonExp[2,x,bt,μ,0.25,gq,gb]

FIGURE 3 LEFT

Timing[R0=Join[Table[{bt,MODELu[0,0.1,bt,(C0/(bstar[bt,1.5])),0.2,0.01]},{bt,(10^(-3)),0.2,0.01}],Table[{bt,MODELu[0,0.1,bt,(C0/(bstar[bt,1.5])),0.2,0.01]},{bt,0.2,8,0.1}]];]
Timing[R1=Join[Table[{bt,MODELu[1,0.1,bt,(C0/(bstar[bt,1.5])),0.2,0.01]},{bt,(10^(-3)),0.2,0.01}],Table[{bt,MODELu[1,0.1,bt,(C0/(bstar[bt,1.5])),0.2,0.01]},{bt,0.2,8,0.1}]];]
Timing[R2=Join[Table[{bt,MODELu[2,0.1,bt,(C0/(bstar[bt,1.5])),0.2,0.01]},{bt,(10^(-3)),0.2,0.01}],Table[{bt,MODELu[2,0.1,bt,(C0/(bstar[bt,1.5])),0.2,0.01]},{bt,0.2,8,0.1}]];]

PLOT1=ListPlot[{R0,R1,R2},Joined->True,PlotStyle->{
{Thick,Dashing[{0.02,0.02,0.001,0.02}],Black},
{Thick,Dashing[{0.02,0.02}],Red},
{Thick,Dashing[{1,0}],Blue}
}]

Show[Graphics[{
Text[Style["b((GeV^(-1)))",FontSize->14],{8.5,0.4}],
{Thick,Dashing[{0.02,0.02,0.001,0.02}],Black,Line[{{5,6},{6,6}}]},
{Thick,Dashing[{0.02,0.02}],Red,Line[{{5,5.2},{6,5.2}}]},
{Thick,Dashing[{1,0}],Blue,Line[{{5,4.4},{6,4.4}}]},
Text[Style["LO",FontSize->14],{6.2,6},{Left,Center}],
Text[Style["NLO",FontSize->14],{6.2,5.2},{Left,Center}],
Text[Style["NNLO",FontSize->14],{6.2,4.4},{Left,Center}]
}]
,PLOT1,
Axes->True,AspectRatio->1/1.6,
AxesLabel->{"",Style["(Subscript[F, u ← N])",FontSize->14]}]

FIGURE 3 RIGHT

Timing[W1=Join[Table[{bt,MODELu[2,0.1,bt,(C0/(bstar[bt,0.5])),0.2,0.01]},{bt,(10^(-3)),0.2,0.01}],Table[{bt,MODELu[2,0.1,bt,(C0/(bstar[bt,0.5])),0.2,0.01]},{bt,0.2,8,0.1}]];]
Timing[W2=Join[Table[{bt,MODELu[2,0.1,bt,(C0/(bstar[bt,1.0])),0.2,0.01]},{bt,(10^(-3)),0.2,0.01}],Table[{bt,MODELu[2,0.1,bt,(C0/(bstar[bt,1.0])),0.2,0.01]},{bt,0.2,8,0.1}]];]
Timing[W3=Join[Table[{bt,MODELu[2,0.1,bt,(C0/(bstar[bt,1.5])),0.2,0.01]},{bt,(10^(-3)),0.2,0.01}],Table[{bt,MODELu[2,0.1,bt,(C0/(bstar[bt,1.5])),0.2,0.01]},{bt,0.2,8,0.1}]];]
Timing[W4=Join[Table[{bt,MODELu[2,0.1,bt,(C0/(bstar[bt,2.])),0.2,0.01]},{bt,(10^(-3)),0.2,0.01}],Table[{bt,MODELu[2,0.1,bt,(C0/(bstar[bt,2.])),0.2,0.01]},{bt,0.2,8,0.1}]];]

PLOT2=ListPlot[{W1,W2,W3,W4},Joined->True,PlotStyle->{
{Thick,Dashing[{0.02,0.02,0.001,0.02}],Black},
{Thick,Dashing[{0.02,0.02}],Red},
{Thick,Dashing[{0.03,0.01}],Darker@Green},
{Thick,Dashing[{1,0}],Blue}
}]

Show[Graphics[{
Text[Style["b((GeV^(-1)))",FontSize->14],{8.5,0.4}],
{Thick,Dashing[{0.02,0.02,0.001,0.02}],Black,Line[{{5,6},{6,6}}]},
{Thick,Dashing[{0.02,0.02}],Red,Line[{{5,5.2},{6,5.2}}]},
{Thick,Dashing[{0.03,0.01}],Darker@Green,Line[{{5,4.4},{6,4.4}}]},
{Thick,Dashing[{1,0}],Blue,Line[{{5,3.6},{6,3.6}}]},
Text[Style["(Subscript[b, max])=0.5",FontSize->14],{6.2,6},{Left,Center}],
Text[Style["(Subscript[b, max])=1.0",FontSize->14],{6.2,5.2},{Left,Center}],
Text[Style["(Subscript[b, max])=1.5",FontSize->14],{6.2,4.4},{Left,Center}],
Text[Style["(Subscript[b, max])=2.0",FontSize->14],{6.2,3.6},{Left,Center}]
}]
,PLOT2,
Axes->True,AspectRatio->1/1.6,
AxesLabel->{"",Style["(Subscript[F, u ← N])",FontSize->14]}]

FIGURE 4 LEFT

Timing[Y1=Join[Table[{bt,MODELu[2,0.1,bt,(C0/(bstar[bt,1.5])),0.2,0.1]},{bt,(10^(-3)),0.2,0.01}],Table[{bt,MODELu[2,0.1,bt,(C0/(bstar[bt,1.5])),0.2,0.1]},{bt,0.2,8,0.1}]];]
Timing[Y2=Join[Table[{bt,MODELu[2,0.1,bt,(C0/(bstar[bt,1.5])),0.2,0.001]},{bt,(10^(-3)),0.2,0.01}],Table[{bt,MODELu[2,0.1,bt,(C0/(bstar[bt,1.5])),0.2,0.001]},{bt,0.2,8,0.1}]];]
Timing[Y3=Join[Table[{bt,MODELu[2,0.1,bt,(C0/(bstar[bt,1.5])),0.2,-0.1]},{bt,(10^(-3)),0.2,0.01}],Table[{bt,MODELu[2,0.1,bt,(C0/(bstar[bt,1.5])),0.2,-0.1]},{bt,0.2,8,0.1}]];]

PLOT3=ListPlot[{Y1,Y2,Y3},Joined->True,PlotStyle->{
{Thick,Dashing[{0.02,0.02,0.001,0.02}],Black},
{Thick,Dashing[{1,0}],Blue},
{Thick,Dashing[{0.02,0.02}],Red}
}]

Show[Graphics[{
Text[Style["b((GeV^(-1)))",FontSize->14],{8.5,0.4}],
{Thick,Dashing[{0.02,0.02,0.001,0.02}],Black,Line[{{5,6},{6,6}}]},
{Thick,Dashing[{1,0}],Blue,Line[{{5,5.2},{6,5.2}}]},
{Thick,Dashing[{0.02,0.02}],Red,Line[{{5,4.4},{6,4.4}}]},
Text[Style["(Subscript[g, q])=+0.1",FontSize->14],{6.2,6},{Left,Center}],
Text[Style["(Subscript[g, q])=0",FontSize->14],{6.2,5.2},{Left,Center}],
Text[Style["(Subscript[g, q])=-0.1",FontSize->14],{6.2,4.4},{Left,Center}]
}]
,PLOT3,
Axes->True,AspectRatio->1/1.6,
AxesLabel->{"",Style["(Subscript[F, u ← N])",FontSize->14]}]

FIGURE 4 RIGHT

Timing[U1=Join[Table[{bt,MODELu[2,0.1,bt,(C0/(bstar[bt,1.5])),0.3,0.01]},{bt,(10^(-3)),0.2,0.01}],Table[{bt,MODELu[2,0.1,bt,(C0/(bstar[bt,1.5])),0.3,0.01]},{bt,0.2,8,0.1}]];]
Timing[U2=Join[Table[{bt,MODELu[2,0.1,bt,(C0/(bstar[bt,1.5])),0.2,0.01]},{bt,(10^(-3)),0.2,0.01}],Table[{bt,MODELu[2,0.1,bt,(C0/(bstar[bt,1.5])),0.2,0.01]},{bt,0.2,8,0.1}]];]
Timing[U3=Join[Table[{bt,MODELu[2,0.1,bt,(C0/(bstar[bt,1.5])),0.1,0.01]},{bt,(10^(-3)),0.2,0.01}],Table[{bt,MODELu[2,0.1,bt,(C0/(bstar[bt,1.5])),0.1,0.01]},{bt,0.2,8,0.1}]];]
Timing[U4=Join[Table[{bt,MODELu[2,0.1,bt,(C0/(bstar[bt,1.5])),0.01,0.01]},{bt,(10^(-3)),0.2,0.01}],Table[{bt,MODELu[2,0.1,bt,(C0/(bstar[bt,1.5])),0.01,0.01]},{bt,0.2,8,0.1}]];]

PLOT4=ListPlot[{U1,U2,U3,U4},Joined->True,PlotStyle->{
{Thick,Dashing[{0.02,0.02,0.001,0.02}],Black},
{Thick,Dashing[{0.02,0.02}],Red},
{Thick,Dashing[{1,0}],Blue},
{Thick,Dashing[{0.03,0.01}],Darker@Green}
}]

Show[Graphics[{
Text[Style["b((GeV^(-1)))",FontSize->14],{8.5,0.4}],
{Thick,Dashing[{0.03,0.01}],Darker@Green,Line[{{5+3,6},{6+3,6}}]},
{Thick,Dashing[{1,0}],Blue,Line[{{5+3,5.2},{6+3,5.2}}]},
{Thick,Dashing[{0.02,0.02}],Red,Line[{{5+3,4.4},{6+3,4.4}}]},
{Thick,Dashing[{0.02,0.02,0.001,0.02}],Black,Line[{{5+3,3.6},{6+3,3.6}}]},
Text[Style["(Subscript[g, b])=0.01",FontSize->14],{6.2+3,6},{Left,Center}],
Text[Style["(Subscript[g, b])=0.1",FontSize->14],{6.2+3,5.2},{Left,Center}],
Text[Style["(Subscript[g, b])=0.2",FontSize->14],{6.2+3,4.4},{Left,Center}],
Text[Style["(Subscript[g, b])=0.3",FontSize->14],{6.2+3,3.6},{Left,Center}]
}]
,PLOT4,
Axes->True,AspectRatio->1/1.6,
AxesLabel->{"",Style["(Subscript[F, u ← N])",FontSize->14]}]

FIGURE 5

Xvalues=Table[Exp[-x]-(10^(-6)),{x,0.,7,0.01}];

Timing[X1=Table[{x,x MODELu[2,x,1.5,(C0/(bstar[1.5,1.5])),0.2,0.1]},{x,Xvalues}];]
Timing[X2=Table[{x,x MODELu[2,x,1.5,(C0/(bstar[1.5,1.5])),0.2,0.001]},{x,Xvalues}];]
Timing[X3=Table[{x,x MODELu[2,x,1.5,(C0/(bstar[1.5,1.5])),0.2,-0.1]},{x,Xvalues}];]

PLOT5=ListLogLinearPlot[{X1,X2,X3},Joined->True,PlotStyle->{
{Thick,Dashing[{0.02,0.02,0.001,0.02}],Black},
{Thick,Dashing[{1,0}],Blue},
{Thick,Dashing[{0.02,0.02}],Red}
},
AxesLabel->{Style["x",FontSize->14],Style["(Subscript[xF, u ← N])",FontSize->14]},
Ticks->{{1,0.5,{0.1,"(10^(-1))"},{0.01,"(10^(-2))"},{0.001,"(10^(-3))",0.03}},Automatic},
 PlotLegend->{"(Subscript[g, q])=0.1", "(Subscript[g, q])=0.", "(Subscript[g, q])=-0.1"},
LegendPosition->{.7,0.2},LegendSize->0.5, LegendShadow->None,LegendBorder->None,
AspectRatio->1/1.4
]

Show[Graphics[{
Text[Style["x",FontSize->14],{8.5,0.4}],
{Thick,Dashing[{0.02,0.02,0.001,0.02}],Black,Line[{{5+4,6},{6+4,6}}]},
{Thick,Dashing[{1,0}],Blue,Line[{{5+4,5.2},{6+4,5.2}}]},
{Thick,Dashing[{0.02,0.02}],Red,Line[{{5+4,4.4},{6+4,4.4}}]},
Text[Style["(Subscript[g, q])=+0.1",FontSize->14],{6.2+4,6},{Left,Center}],
Text[Style["(Subscript[g, q])=0",FontSize->14],{6.2+4,5.2},{Left,Center}],
Text[Style["(Subscript[g, q])=-0.1",FontSize->14],{6.2+4,4.4},{Left,Center}]
}],AspectRatio->1/1.6]

Show[%107,%105]

Plot for QCD evol 2016 proceedings

Timing[R0=Join[Table[{bt,MODELu[0,0.1,bt,(C0/(bstar[bt,1.5])),0,0.00000000001]},{bt,(10^(-3)),0.2,0.01}],Table[{bt,MODELu[0,0.1,bt,(C0/(bstar[bt,1.5])),0,0.0000000001]},{bt,0.2,8,0.1}]];]
Timing[R1=Join[Table[{bt,MODELu[1,0.1,bt,(C0/(bstar[bt,1.5])),0,0.00000000001]},{bt,(10^(-3)),0.2,0.01}],Table[{bt,MODELu[1,0.1,bt,(C0/(bstar[bt,1.5])),0,0.00000000001]},{bt,0.2,8,0.1}]];]
Timing[R2=Join[Table[{bt,MODELu[2,0.1,bt,(C0/(bstar[bt,1.5])),0,0.00000000001]},{bt,(10^(-3)),0.2,0.01}],Table[{bt,MODELu[2,0.1,bt,(C0/(bstar[bt,1.5])),0,0.00000000001]},{bt,0.2,8,0.1}]];]

PLOT1=ListPlot[{R0,R1,R2},Joined->True,PlotStyle->{
{Thick,Dashing[{0.02,0.02,0.001,0.02}],Black},
{Thick,Dashing[{0.02,0.02}],Red},
{Thick,Dashing[{1,0}],Blue}
},PlotRange->{0,15}]

With[{a=1,b=2},Show[Graphics[{
Text[Style["b((GeV^(-1)))",FontSize->14],{8.5,0.4}],
{Thick,Dashing[{0.02,0.02,0.001,0.02}],Black,Line[{{5+a,6+b},{6+a,6+b}}]},
{Thick,Dashing[{0.02,0.02}],Red,Line[{{5+a,5.2+b},{6+a,5.2+b}}]},
{Thick,Dashing[{1,0}],Blue,Line[{{5+a,4.4+b},{6+a,4.4+b}}]},
Text[Style["LO",FontSize->14],{6.2+a,6+b},{Left,Center}],
Text[Style["NLO",FontSize->14],{6.2+a,5.2+b},{Left,Center}],
Text[Style["NNLO",FontSize->14],{6.2+a,4.4+b},{Left,Center}]
}]
,PLOT1,
Axes->True,AspectRatio->1/1.6,
AxesLabel->{"",Style["(Subscript[F, u ← N])",FontSize->14]},PlotRange->{0,10}]]