Intro to Mathematica

Press 'Shift' + 'Enter' to calculate cell

2+2

<<PhysicalConstants`

(Subscript[m, p])=ProtonMass

ℏ=PlanckConstantReduced

Time Independent Potential

V[x_]:=0

DSolve[Energy*ψ[x]==((-(ℏ^2))/(2*(Subscript[m, p])))ψ''[x]+V[x]*ψ[x],ψ[x],x]

ψ[0]

Apply Boundary Conditions (BC's)

DSolve[{Energy*ψ[x]==((-(ℏ^2))/(2*(Subscript[m, p])))ψ''[x],ψ[0]==0},ψ[x],x]

ψ[x_]:=C[1]* Sin[(((Sqrt[2]) (Sqrt[Energy]) x (Sqrt[(Subscript[m, p])]))/ℏ)]

Solve[ ψ[x]==0,x]

What about the other bc?

Reduce[{ Sin[(((Sqrt[2]) (Sqrt[Energy])*(Sqrt[(Subscript[m, p])])*L )/ℏ)]==0,Energy>0,L>0},Energy,Reals]

(Kilogram==0&&L>0&&Energy>0)||(C[1]∈Integers&&((C[1]≤-1.`&&((Joule<0&&Kilogram>0&&L>0&&Second>0&&(Energy==((1.312454319956745`*^-40 (Joule^2) (Second^2) ((C[1])^2))/(Kilogram (L^2)))||Energy==(1/(Kilogram (L^2)))3.324485629362923`*^-42 (9.869604401089358` (Joule^2) (Second^2)+39.47841760435743` (Joule^2) (Second^2) C[1]+39.47841760435743` (Joule^2) (Second^2) ((C[1])^2))))||(Joule>0&&Kilogram>0&&L>0&&Second<0&&(Energy==((1.312454319956745`*^-40 (Joule^2) (Second^2) ((C[1])^2))/(Kilogram (L^2)))||Energy==(1/(Kilogram (L^2)))3.324485629362923`*^-42 (9.869604401089358` (Joule^2) (Second^2)+39.47841760435743` (Joule^2) (Second^2) C[1]+39.47841760435743` (Joule^2) (Second^2) ((C[1])^2))))))||(C[1]≥1.`&&((Joule<0&&Kilogram>0&&L>0&&Second<0&&Energy==((1.312454319956745`*^-40 (Joule^2) (Second^2) ((C[1])^2))/(Kilogram (L^2))))||(Joule>0&&Kilogram>0&&L>0&&Second>0&&Energy==((1.312454319956745`*^-40 (Joule^2) (Second^2) ((C[1])^2))/(Kilogram (L^2))))))||(C[1]≥0&&((Joule<0&&Kilogram>0&&L>0&&Second<0&&Energy==(1/(Kilogram (L^2)))3.324485629362923`*^-42 (9.869604401089358` (Joule^2) (Second^2)+39.47841760435743` (Joule^2) (Second^2) C[1]+39.47841760435743` (Joule^2) (Second^2) ((C[1])^2)))||(Joule>0&&Kilogram>0&&L>0&&Second>0&&Energy==(1/(Kilogram (L^2)))3.324485629362923`*^-42 (9.869604401089358` (Joule^2) (Second^2)+39.47841760435743` (Joule^2) (Second^2) C[1]+39.47841760435743` (Joule^2) (Second^2) ((C[1])^2)))))))

Energy==(1/(Kilogram (L^2)))3.324485629362923`*^-42 *9.869604401089358` (Joule^2) (Second^2) for n=1
Energy==((1.312454319956745`*^-40 (Joule^2) (Second^2) ((C[1])^2))/(Kilogram (L^2)))  is a factor of 4 larger than the ground state for n larger than 1.

Normalize

(∫_0^L)(((Abs[ψ[x]]))^2)\[DifferentialD]x

Note that  Sin[((1.0969018417211581`*^21 (Sqrt[Energy]) (Sqrt[Kilogram]) L)/(Joule Second))=0   Mathematica does not know that.

NSolve[((C[1])^2) (0.5` L)==1,C[1]]

Let's just pick   C[1]->((1.4142135623730951`)/(Sqrt[L]))   since it doesn't matter if C[1] is positive or negative.

Now, we can define our wave function:

L=1

ψ[n_,x_]:=((1.4142135623730951`)/(Sqrt[L]))* Sin[n*Pi*x/L]

Plot[Evaluate@
  Append[Table[UnitStep[x]*(ψ[n,x]+n^2)*UnitStep[-x+L],{n,1,7}],x^2/2],{x,-2L,2L},Filling->Table[n->n^2,{n,1,7}],PlotRange->{0,58}]

Plot[Evaluate@
  Append[Table[UnitStep[x]*(ψ[n,x]+n^2)*UnitStep[L-x],{n,1,7}],x^2/2],{x,-2L,2L},Filling->Table[n->n^2,{n,1,7}],PlotRange->{0,58}]

Parabolic Well

Define a function

f[n_,x_]:= Abs[((1/Pi)^(1/4)HermiteH[n,x])/(E^(x^2/2)Sqrt[2^n n!])]^2

Plot[Evaluate@
  Append[Table[f[n,x]+n+1/2,{n,0,7}],x^2/2],{x,-4,4},Filling->Table[n->n-1/2,{n,1,8}]]

ProtonMass