a[n_]:=DiagonalMatrix[Sqrt[Range[n]],1]
ad[n_]:=DiagonalMatrix[Sqrt[Range[n]],-1]
h[n_,x_,y_]:=ArrayFlatten[{{ad[n].a[n]+x^2 IdentityMatrix[n+1],x a[n]+y ad[n]},{x ad[n]+y a[n],ad[n].a[n]+y^2 IdentityMatrix[n+1]}}]
gr[0]=ListPlot[Transpose[Table[{x,#}&/@Sort[Eigenvalues[N[h[50,x,0]],-16]],{x,-1.5,1.5,0.1}],{2,1,3}],Joined->True,PlotRange->All];
e[n_,x_,s_]:=1/2 (x^2+1)+n+s/2 Sqrt[(x^2+1)^2+4 x^2 n]
gr[-1]=Plot[e[#,x,-1]&/@Range[0,6],{x,-1.5,1.5},PlotStyle->Dashed,Evaluated->True];
gr[1]=Plot[e[#,x,1]&/@Range[0,6],{x,-1.5,1.5},PlotStyle->Dotted,Evaluated->True];
Show[{gr[1],gr[0],gr[-1]},PlotRange->{All,{-1,9}}]