Harmonic Oscillator

eqn=-ℏ^2/(2μ)ψ''[x]+k x^2/2 ψ[x]==Ev ψ[x]

Ev:=(-ℏ^2/(2μ)ψ''[x]+k x^2/2 ψ[x])/ψ[x]

α=Sqrt[k μ]/ℏ

ψ[x_]=Exp[-α x^2 / 2]

Ev

Simplify[%,{k > 0, μ > 0, ℏ > 0}]

PowerExpand[%]

ψ[x_]=x Exp[-α x^2 / 2]

Simplify[Ev,{k > 0, μ > 0, μ > 0}]

PowerExpand[%]

Cnorm[v_]:=Sqrt[(1/(2^v v!))](α/Pi)^(1/4)

Cnorm[0]

Herm[0,y_]=1

Herm[1,y_]:= 2y

Herm[v_,y_]:=2 y Herm[v-1,y]-2(v-1)Herm[v-2,y]

Herm[3,y]

Herm[2,y]

ψv[v_,y_]:=Cnorm[v] Herm[v,y] Exp[-y^2/2]

ψv[0,y]

ψv[1,y]

Integrate[ψv[0,y]^2 / Sqrt[α],{y,-Infinity,Infinity}]

Integrate[ψv[1,y]^2 / Sqrt[α],{y,-Infinity,Infinity}]

Integrate[ψv[2,y]^2 / Sqrt[α],{y,-Infinity,Infinity}]

Integrate[ψv[3,y]^2 / Sqrt[α],{y,-Infinity,Infinity}]

Plot[{Herm[0,y] Exp[-y^2/2],Herm[1,y] Exp[-y^2/2],Herm[2,y] Exp[-y^2/2]},{y,-3,3}]

Evib[v_]:=(-ℏ^2/(2μ)D[ψv[v,Sqrt[α]x],{x,2}]+k x^2/2 ψv[v,Sqrt[α]x])/ψv[v,Sqrt[α]x]/(ℏ Sqrt[k/μ])

condition= {k > 0, μ > 0, ℏ > 0}

Do[Print["v = ",n,"  E = ",Simplify[Evib[n],condition]," hν"],{n,0,3}]

avex=Integrate[y ψv[0,y]^2 / Sqrt[α],{y,-Infinity,Infinity}] / Sqrt[α]

avex2=Integrate[y^2 ψv[0,y]^2 / Sqrt[α],{y,-Infinity,Infinity}] / α

△x=Sqrt[avex2-avex^2]

avep=-ℏ Integrate[ ψv[0,y] D[ψv[0,y],y] ,{y,-Infinity,Infinity}]

avep2=-ℏ^2 Integrate[ ψv[0,y] D[ψv[0,y],{y,2}]  Sqrt[α],{y,-Infinity,Infinity}]

△p=Sqrt[avep2-avep^2]

△x △p

Simplify[△x △p,condition]

eqny=1/2 k (y/Sqrt[α])^2==(v + 1/2)ℏ Sqrt[k/μ]

Solve[eqny,y]

Simplify[%,condition]

ytp[v_]:= Sqrt[1+ 2 v]

Pforb[v_]:=2 Integrate[ψv[v,y]^2 ,{y,ytp[v],Infinity}] /Sqrt[α]

Pforb[0]

N[%]

Pforb[1]

N[%]

N[Pforb[2]]

Pforb[1]

N[Pforb[3]]

Pforb[3]

Pforb[4]

N[%]

Do[Pf[n]=N[Pforb[n]],{n,0,6}]

DiscretePlot[Pf[n],{n,0,6},PlotRange->{{0,6},{0,0.2}}]