(* Analytical solution to y[z]. *) 
(z/2)(∫_0^1)(Exp[-z^2 mu/2]/((1-mu^2)^(1/4)))\[DifferentialD]mu

(* No good for values of z>7, so we seek numerical approach *)
y[z_]:=(√(π/2) z Gamma[3/4] (BesselI[1/4,(z^2)/2]-StruveL[1/4,(z^2)/2]))/(2 ((((z^2)))^(1/4)))
Plot[y[z],{z,-8.5,8.5}]

(* Numerical representation, which we compare with 1/z *)
ClearAll[yN];
yN[z_?NumericQ]:=NIntegrate[(z/2)(Exp[-z^2 mu/2]/((1-mu^2)^(1/4))),{mu,0,1}]
Plot[{yN[z],1/z},{z,-10,10},Frame->True,PlotStyle->{Red,Blue},FrameLabel->{{"y(z)",""},{"z",""}},PlotLegends->Placed[{"y(z)","1/z"},{Left,Top}]]

(* Plot y[z] and compare it with (∂_z)(∂_z)y(z).  The latter is obtained by the definition of the differential equation (∂_z)(∂_z)y(z) = z^2 y[z] - z  *) 

ypp[z_]:=z^2 y[z] - z 


Plot[{y[z],-ypp[z]},{z,-5,5},Frame->True,FrameLabel->{{"y(z), and (∂_z)(∂_z)y(z)",""},{"z",""}},PlotLegends->Placed[{"y(z)","(∂_z)(∂_z)y(z)"},{Left,Top}]]

(* Now calculate Delta'(X) in the layer, where X would be taken to infinity.  Look for Delta'(x) to saturate to correct value, 2π Gamma[3/4]/Gamma[1/4, for large X *)


ClearAll[Delp];
Delp[X_?NumericQ]:=NIntegrate[1-x yN[x],{x,-X,X}]
Plot[Delp[X],{X,1,10},FrameLabel->{{"Δ'(X)",""},{"X",""}},PlotLegends->Placed[{"(∫_(-X)^X)[1-z y(z)]\[DifferentialD]z"},{Right,Bottom}],PlotRange->{{0,10},{1.7,2.2}},Frame->True]

(* Compare Delp[X], with X=10 with 2π Gamma[3/4]/Gamma[1/4].  Within 0.06 percent: *)

(Delp[10]-2π Gamma[3/4]/Gamma[1/4])/Delp[10]