(* Full Screen: Command-Option-F *)
(* Evaluate Cells: Shift-Enter *)
(* Slow evaluations have been pre-computed *)

(* TOWARDS A SPECTRAL PROOF OF RIEMANN'S HYPOTHESIS *)

(* Zeta[s] = Sum[n^{-s},{n,1,Infinity}] | Re[s]>1 *)
(* Gamma[s] = Integrate[x^{s-1}Exp[-x],{x,0,Infinity}] | Re[s]>0 *)

Riemannxi[t_] := -(t^2+1/4) Gamma[I t/2+1/4]Zeta[I t+1/2]/(2 Pi^(I t/2+1/4))
TraditionalForm[Riemannxi[t]]
Plot[Re[Riemannxi[t]],{t,0,100},AxesOrigin->{0,0},PlotRange->All,AxesLabel->{t,ξ}]

ScalingR[t_] := 2^{5/4}Cosh[Pi t/4]/(Pi^{1/4}(t^2+4)^{7/8})
TraditionalForm[ScalingR[t]]
Plot[ScalingR[t]Re[Riemannxi[t]],{t,0,100},AxesLabel->{t,Scaled ξ},PlotRange->All]

Plot[Re[Riemannxi[2Sqrt[e]]],{e,-10,50},AxesLabel->{e,ξ[2Sqrt[e]]}]

Density[e_]:=Log[Sqrt[e]/Pi]/Sqrt[e]
TraditionalForm[Density[e]]

heldint = HoldForm[Integrate[Density[e]/Pi/Sqrt[v-e],{e,0,v},Assumptions->{v>0}]];
int=ReleaseHold[heldint];
TraditionalForm[Row[{heldint," == ",int}]]

Plot[4Pi^2 Exp[2x],{x,0,1},AxesLabel->{x,V},AxesOrigin->{0,0}]

Plot[{ScalingR[t]Re[Riemannxi[t]],ScalingR[t]500BesselK[I t/2,2Pi]},{t,0,100},PlotRange->All,AxesLabel->{t,K[I t/2,2Pi]}]

TraditionalForm[Row[{"Riemannxi[-I tau] ~ ",(tau/(2Pi E))^{tau/2}tau^{7/4}(Pi/2)^{1/4}}]]

Plot[{4Pi^2Exp[2x],4Pi^2Exp[2x]-9Pi Exp[x]},{x,0,1},AxesOrigin->{0,0},AxesLabel->{x,V}]

Plot[{ScalingR[t]Re[Riemannxi[t]],ScalingR[t]2 Pi^{-1/4}WhittakerW[9/4,I t/2,4Pi]},{t,0,100}]

ParametricPlot3D[{Sin[2 Pi x]/y,Cos[2 Pi x]/y,Log[y+Sqrt[y^2-1]]-Sqrt[1-y^{-2}]},{y, 1, 20},{x,0,1},PlotRange->All,PlotPoints->100]

(* Magnetic Laplacian Δ_B = -y^2 (d_x^2 + d_y^2) + 2iBy d_x + B^2 *)
(* ϕ(y)Exp[-2πimx] eigenfunction for eigenvalue (t^2+1)/4 + B^2 *)(* iff (-u^2 d_u^2 + (u^2-1-t^2)/4 \pm Bu)ϕ = 0 with u = 4π|m|y *)
(* Solution decaying as y ->Infinity is WhittakerW[\pm B,it/2,4π|m|y] *)
(* So general solution evaluated at x=0,y=1, is: Sum[a[m]WhittakerW[9/4,it/2,4πm]+a[-m]WhittakerW[-9/4,it/2,4πm],{m,1,Infinity}] *)
(* How to choose the coefficients a? *)

PolyaPhi[v_,NN_] := Sum[(4Pi^2n^4 Exp[9v/2]-6Pi n^2 Exp[5v/2])Exp[-Pi n^2 Exp[2v]],{n,1,NN}] (* even! *)
TraditionalForm[PolyaPhi[v,Infinity]]
Plot[PolyaPhi[v,8],{v,0,1.5}]

RSasymptotics[t_]:=-2Sum[n^{-1/2}Cos[Arg[Gamma[1/4+I t/2]]-t Log[Pi]/2-t Log[n]],{n,1,Sqrt[t/(2Pi)]}] 
Plot[{ScalingR[t]Re[Riemannxi[t]],RSasymptotics[t]},{t,2Pi,100},AxesOrigin->{0,0},AxesLabel->{t,RSasymptotics},PlotPoints->300]

Plot[{ScalingR[t]Re[Riemannxi[t]],ScalingR[t]2Pi^{-1/4}Sum[n^{-3/2}Re[WhittakerW[9/4,I t/2,4Pi n^2]],{n,3}]},{t,0,100}]

p = ParametricPlot3D[{Sin[2 Pi x]/y,Cos[2 Pi x]/y,Log[y+Sqrt[y^2-1]]-Sqrt[1-y^{-2}]},{y, Pi Sqrt[2.5], 20},{x,0,1},PlotRange->All,PlotPoints->100,BoxRatios->{1,1,1}];
c1 = ParametricPlot3D[{0,1/y,Log[y+Sqrt[y^2-1]]-Sqrt[1-y^{-2}]},{y,2 Pi, 20},PlotStyle->Thick];
c2 = ParametricPlot3D[{0,-1/y,Log[y+Sqrt[y^2-1]]-Sqrt[1-y^{-2}]},{y,Pi Sqrt[3], 20},PlotStyle->Thick];
c3 = ParametricPlot3D[{Sin[2 Pi Cos[t]]/(2 Pi Sin[t]),Cos[2 Pi Cos[t]]/(2 Pi Sin[t]),Log[2 Pi Sin[t]+Sqrt[(2 Pi Sin[t])^2 -1]]-Sqrt[1-(2 Pi Sin[t])^{-2}]},{t,Pi/3,2 Pi/3},PlotStyle->Thick];
Show[{p,c1,c2,c3}]