(* f is the integrand of the integral inside the curly bracket in equation 3.16 of the paper. lambda=Pi, 1/kb=4, and Delta u= dU. I show the variable z with x.*)

f[dU_,x_]:=(x^(4*Pi))*Exp[Pi*dU/2]*Gamma[1-2*Pi,-I*x]*Gamma[1-2*Pi,I*x*Exp[dU/4]]/x

(* IntfK is the integrand of the integral inside the curly bracket in equation 3.16. Lam is the upper cutoff of the integral.*)

IntfK[dU_,Lam_]:=NIntegrate[Re[f[dU,x]],{x,0,Lam}]

(*Here I want approximate the integrand for large x (x is same as z). Hence, I use the series for the incomplete gamma functions for large argument. In Apxp and Appxp, I clean up Mathematica results.*)

Series[Gamma[1-2*Pi,-I*x],{x,Infinity,1}]

Apxp[x_]:=Exp[I*x]*(x^(-2 π)) *((((-I))^(-2 π))-((2 I (((-I))^(-2 π)) π)/x))

Series[Gamma[1-2*Pi,I*x*Exp[dU/4]],{x,Infinity,1}]

Appxp[x_,dU_]:=Exp[-I*x*Exp[dU/4]]*(x^(-2 π)) *((((I (E^(dU/4))))^(-2 π))+((2 I (E^(-dU/4)) (((I (E^(dU/4))))^(-2 π)) π)/x))

(*Then I use these results to find the contribution of the large frequencies to the two-point function (first formulation)*)

Integrate[(x^(4*Pi))*Exp[Pi*dU/2]*(Exp[I*x]*(x^(-2 π)) *((((-I))^(-2 π))-((2 I (((-I))^(-2 π)) π)/x)))*(Exp[-I*x*Exp[dU/4]]*(x^(-2 π)) *((((I (E^(dU/4))))^(-2 π))+((2 I (E^(-dU/4)) (((I (E^(dU/4))))^(-2 π)) π)/x)))/x,x]

(*This is the contribution of the large frequencies to the two-point function (first formulation) including all constant coefficients.*)

G[x_,dU_]:=(1/(x^2))(E^((1/4) dU (-1+2 π)+I ((π^2)+x-(E^(dU/4)) x))) (((I (E^(dU/4))))^(-2 π)) (2 I π ((-1+(E^(dU/4))) x+π (I+(-1+(E^(dU/4))) x))-(E^(I (-1+(E^(dU/4))) x)) (2 π (1+π)+2 (E^(dU/2)) π (1+π)-(E^(dU/4)) (((1+2 π))^2)) (x^2) ExpIntegralEi[-I (-1+(E^(dU/4))) x])

(*Here, I tese what happens to the integral as x goes to the infinity.*)

Limit[Re[G[x,2]],x->Infinity]

Limit[Re[G[x,3]],x->Infinity]

N[Re[G[100,10]-G[50,10]]]

(* These are the final data poinst which have been optained by adding the contribution of the large x (G) and small and intermediate values of x (IntfK).*)

Table[{j/10,(IntfK[j/10,800]-Re[G[800,j/10]])*(1/(2*Pi))},{j,1,50}]

T2={{(1/10),0.20852341357887658`},{(1/5),0.10926081209233297`},{(3/10),0.05856281374426369`},{(2/5),0.028086970896852755`},{(1/2),0.008590198434961688`},{(3/5),-0.004192446310963626`},{(7/10),-0.012594590989644995`},{(4/5),-0.01802785253781675`},{(9/10),-0.021403544438290067`},{1,-0.02333559636935738`},{(11/10),-0.024250048110692778`},{(6/5),-0.024448760177537922`},{(13/10),-0.024148923721832503`},{(7/5),-0.02350803653845114`},{(3/2),-0.022640735362173786`},{(8/5),-0.02163113143186024`},{(17/10),-0.020540451474973548`},{(9/5),-0.01941302841869782`},{(19/10),-0.01828079922708292`},{2,-0.017166356766776695`},{(21/10),-0.01608538450079139`},{(11/5),-0.015048296893953803`},{(23/10),-0.014061545083891952`},{(12/5),-0.013128884299342563`},{(5/2),-0.012251898216667893`},{(13/5),-0.011430592528667927`},{(27/10),-0.010663936554089438`},{(14/5),-0.0099501829027663`},{(29/10),-0.009286959562700255`},{3,-0.008671736593342744`},{(31/10),-0.008101663578803042`},{(16/5),-0.007573934786368284`},{(33/10),-0.0070857062281182385`},{(17/5),-0.006634157747654936`},{(7/2),-0.0062166631782742884`},{(18/5),-0.005830598875776132`},{(37/10),-0.005473586549584412`},{(19/5),-0.005143338236996579`},{(39/10),-0.0048377306466678574`},{4,-0.004554737583497906`},{(41/10),-0.0042925849390256406`},{(21/5),-0.004049545397471967`},{(43/10),-0.0038240547555075445`},{(22/5),-0.0036147063299670846`},{(9/2),-0.003420130788165659`},{(23/5),-0.003239158778249621`},{(47/10),-0.0030706785016410675`},{(24/5),-0.0029136692810896787`},{(49/10),-0.00276720123160325`},{5,-0.002630431313796129`}};

ListLinePlot[T2,PlotStyle->{Thickness[0.003],Black}]

(* Now I compute G_K using the second formulation. Equation 3.18 of the paper*)

f[d_]:=2*Pi*I*((1/(4 (π^2) (-1+4 (π^2))))I (E^(-d/4)) (Hypergeometric2F1[2,1-2 π,2-2 π,(E^(-d/4))]+2 π Hypergeometric2F1[2,1-2 π,2-2 π,(E^(-d/4))]-Hypergeometric2F1[2,1+2 π,2+2 π,(E^(-d/4))]+2 π Hypergeometric2F1[2,1+2 π,2+2 π,(E^(-d/4))]))

g[d_]:=2*Pi*I*Exp[2*I*(Pi)^2]*Exp[-Pi*d/2]/(2*Sinh[2*Pi*I*Pi])

p1=Plot[(1/(4*Pi))*(f[d]+Re[g[d]]),{d,1/10,10},PlotRange->{{1/10,10},{-0.2,0.25}},PlotStyle->Darker[Green]]

p1=Plot[(1/(4*Pi))(f[d]+Re[g[d]]),{d,1/10,10},Frame->True,FrameLabel->{"Δu",""},ImageSize->300,LabelStyle->{FontSize->12,FontFamily->"Courier New",Black},PlotRange->{{0,10},{-0.03,0.025}},PlotStyle->{{Black,Thickness[0.003]}},FrameTicks->Automatic,Axes->False]

p2=ListPlot[T2,Frame->True,FrameLabel->{"(Subscript[Δu, s])",""},ImageSize->300,LabelStyle->{FontSize->12,FontFamily->"Courier New",Black},PlotRange->{{0,10},{-0.2,0.2}},PlotStyle->{{Black,Thickness[0.001]}},PlotMarkers->{Automatic,1.5},FrameTicks->Automatic,Axes->False]

Show[p1,p2]

(* So they have agreement. Now I make a nicer plot*)

(*The following table is for the two-point function in the first formulation. This is only the integral I numerically computed before. *)

Table[{j/2,IntfK[j/2,800]-Re[G[800,j/2]]},{j,1,20}]

T9={{19/50,0.20805220319369944`},{21/50,0.1475764314762373`},{23/50,0.09679214294651874`},{(1/2),0.05397380859230835`},{27/50,0.017781673520523917`},{(3/5),-0.026341999421578867`},
{34/50,-0.07029454605834501`},
{(4/5),-0.11327205911233235`},{1,-0.1466218762422196`},{(3/2),-0.14225593577135162`},{2,-0.1078594006148142`},{(5/2),-0.07698094686002747`},{3,-0.054486127951022686`},{(7/2),-0.03906044674141736`},{4,-0.028618260262692694`},{(9/2),-0.021489315516835007`},{5,-0.016527487382388933`},{(11/2),-0.012988962874729384`},{6,-0.010399761166469386`},{(13/2),-0.008457589245118899`},{7,-0.0069673585810162255`},{(15/2),-0.005800857455409442`},{(81/10),-0.004709143092991647`},{(43/5),-0.003989033026100046`},{(46/5),-0.0032931523266987554`},{(48/5),-0.0029091945981811317`},{10,-0.0025768198299440086`}};

p9=ListPlot[T9,Frame->True,FrameLabel->{"Δu",""},ImageSize->300,LabelStyle->{FontSize->12,FontFamily->"Courier New",Black},PlotRange->{{0,10},{-0.2,0.25}},PlotStyle->{{Black,Thickness[0.001]}},PlotMarkers->{Automatic,3},FrameTicks->Automatic,Axes->False]

(*But I need to multiply the vertical axis by (1/2pi)*)

T10={{19/50,(1/(2*Pi))*0.20805220319369944`},{21/50,0.1475764314762373`*(1/(2*Pi))},{23/50,0.09679214294651874`*(1/(2*Pi))},{(1/2),0.05397380859230835`*(1/(2*Pi))},{27/50,0.017781673520523917`*(1/(2*Pi))},{(3/5),-0.026341999421578867`*(1/(2*Pi))},
{34/50,-0.07029454605834501`*(1/(2*Pi))},
{(4/5),-0.11327205911233235`*(1/(2*Pi))},{1,-0.1466218762422196`*(1/(2*Pi))},{(3/2),-0.14225593577135162`*(1/(2*Pi))},{2,-0.1078594006148142`*(1/(2*Pi))},{(5/2),-0.07698094686002747`*(1/(2*Pi))},{3,-0.054486127951022686`*(1/(2*Pi))},{(7/2),-0.03906044674141736`*(1/(2*Pi))},{4,-0.028618260262692694`*(1/(2*Pi))},{(9/2),-0.021489315516835007`*(1/(2*Pi))},{5,-0.016527487382388933`*(1/(2*Pi))},{(11/2),-0.012988962874729384`*(1/(2*Pi))},{6,-0.010399761166469386`*(1/(2*Pi))},{(13/2),-0.008457589245118899`*(1/(2*Pi))},{7,-0.0069673585810162255`*(1/(2*Pi))},{(15/2),-0.005800857455409442`*(1/(2*Pi))},{(81/10),-0.004709143092991647`*(1/(2*Pi))},{(43/5),-0.003989033026100046`*(1/(2*Pi))},{(46/5),-0.0032931523266987554`*(1/(2*Pi))},{(48/5),-0.0029091945981811317`*(1/(2*Pi))},{10,-0.0025768198299440086`*(1/(2*Pi))}};

p10=ListPlot[T10,Frame->True,FrameLabel->{"Δu",""},ImageSize->300,LabelStyle->{FontSize->12,FontFamily->"Courier New",Black},PlotRange->{{0,10},{-0.2,0.25}},PlotStyle->{{Black,Thickness[0.001]}},PlotMarkers->{Automatic,3},FrameTicks->Automatic,Axes->False]

Show[p1,p10]