6. Quantum Mechanics

Legal Notice, Disclaimer, and Book Ordering Information:

Null

Setup

\(Clear["Global`*"]; 
Off[General::spell]; 
Off[General::spell1]; 
    \)

Needs["Calculus`VectorAnalysis`"]

Needs["Graphics`PlotField3D`"]

<< "D:\Mma\Inits\ch6init.m"

Clear["Global`*"];

6.2 Problems

6.2.1  One Dimensional Schrödinger Equation

Problem 2:  Bound particle in an finite potential well

Suppose a particle is bound by the finite potential well V[x] = -V0 if -a<x<a and =0 otherwise.

Part (a)  Use the condition that the wave function must vanish at x->∞ or x->-∞ and write the general form of the wave function in the regions to the left, right, and center of the well.

ψW[x_] = waveTrig[x,ki,A,B]

kRule = {k -> Sqrt[2m(En-V)]/ℏ}//PowerExpand

kiRule = {ki -> k}/.kRule//.{V->-V0, En->-Wn}

waveExp[x,k,c1,c2]

ψL[x_] = waveExp[x,I k0,0,C]  (* x < -a  *)

ψR[x_] = waveExp[x,I k0,D,0]  (* x > a  *)

k0Rule = {k0 -> -I k}/.kRule//.{V->0, En->-Wn}//PowerExpand

Part (b)  Require that the wave function and its derivative be continuous across the boundaries of the well and show that there are two classes of solutions, odd and even.  Find the relation between the wave function constants and derive an equation for the energy eigenvalues.

eq2 = {(ψL[ x] - ψW[ x] == 0)/.{x -> -a}
	  ,(ψW[ x] - ψR[ x] == 0)/.{x -> +a}
	  ,(ψL'[x] - ψW'[x] == 0)/.{x -> -a}
	  ,(ψW'[x] - ψR'[x] == 0)/.{x -> +a}}

eq2a = {(ψL'[x]/ψL[x] == ψW'[x]/ψW[x])/.{x -> -a}
	  ,(ψW'[x]/ψW[x] == ψR'[x]/ψR[x])/.{x -> +a}}

sSol=Solve[eq2a,{k0,A}]//Flatten

symEn=Thread[(sSol[[1]]/.{Rule->Equal})/ki,Equal]

aSol=Solve[eq2a,{k0,B}]//Flatten

asymEn=Thread[(aSol[[1]]/.{Rule->Equal})/ki,Equal]

eq3s = eq2//.sSol

symSol=Join[Solve[eq3s,{B,C}]//Flatten,{A->0}]

eq3a = eq2//.aSol

asymSol=Join[Solve[eq3a,{A,C}]//Flatten,{B->0}]

Part (c)    Express the energy eigenvalues Wn in a form similar to those of the infinite well bymeasuring the energy levels from the bottom of the well and parametrizing the values with n, where Wn = -(n^2)(ℏ^2)(π^2)/(8(a^2) m).  Solve for the allowed n values by plotting two sides of the eigenvalue equation and observing where the curves intersect.  Verify that the values of n approach the integer values of the infinite well as V->∞.

kRule2 = Join[kiRule,k0Rule]

nRule={Wn->V0 - (n ℏ π/a)^2/(8 m)}

eq6=symEn/.kRule2/.nRule//PowerExpand

values={a->1,m->1,ℏ->1,V0->{100,200,500,∞}}

text1={
	 Text[FontForm["V0=∞"  ,{"Courier-Bold",10}],{1.8,1.65}]
	,Text[FontForm["V0=500",{"Courier-Bold",10}],{5.85,1.4}]
	,Text[FontForm["V0=200",{"Courier-Bold",10}],{7.75,1.05}]
	,Text[FontForm["V0=100",{"Courier-Bold",10}],{7.75,0.75}]};

plot1=Plot[{ArcTan[eq6[[1]]],ArcTan[eq6[[2]]]}
	/.values//Flatten//Evaluate
	,{n,0.01,9}
	,Frame->True
	,FrameTicks->{{{1,"n=1"},{3,"n=3"},{5,"n=5"}
				  ,{7,"n=7"},{9,"n=9"}},{0,0.5,1,1.5}}
	
	,PlotRange->{0,2}
	,Epilog->Graphics[text1][[1]]];

eq7=asymEn/.kRule2/.nRule//PowerExpand

text2={
	 Text[FontForm["V0=∞"  ,{"Courier-Bold",10}],{1.8,1.65}]
	,Text[FontForm["V0=500",{"Courier-Bold",10}],{4.8,1.43}]
	,Text[FontForm["V0=200",{"Courier-Bold",10}],{6.75,1.13}]
	,Text[FontForm["V0=100",{"Courier-Bold",10}],{6.75,0.87}]};

plot2=Plot[{ArcTan[eq7[[1]]],ArcTan[eq7[[2]]]}
	/.values//Flatten//Evaluate
	,{n,0.01,9}
	,Frame->True
	,FrameTicks->{{{2,"n=2"},{4,"n=4"},{6,"n=6"}
				  ,{8,"n=8"},{10,"n=10"}},{0,0.5,1,1.5}}
	
	,PlotRange->{0,2}
	,Epilog->Graphics[text2][[1]]];

Part (d)  Solve for the n values in a second way by finding the roots of the eigenvalue equation using the command FindRoot.

nValues[eq_,potential_,guess_]:=
	(n/.FindRoot[eq/.{m->1,a->1,ℏ->1,V0->potential}
	//Evaluate
	,{n,guess}][[1]])

symGuess={.9,2.9,4.9,6.9,8.9};

symValues = {nValues[eq6,100,#]
			,nValues[eq6,200,#]
			,nValues[eq6,500,#]}& /@ symGuess;

asymGuess={1.9,3.9,5.9,7.9.9.9};

asymValues = {nValues[eq7,100,#]
			,nValues[eq7,200,#]
			,nValues[eq7,500,#]}& /@ asymGuess;

(Partition[Sort[{symValues,asymValues}
	//Flatten//N[#,4]&],3]
	//TableForm[#,TableSpacing->{1,3}
				 ,TableHeadings->{
				 	{"n=1","n=2","n=3","n=4","n=5"
				 		  ,"n=6","n=7","n=8","n=9"},
				 	{"V0=100","V0=200","V0=500"}}]&)

Part (e)  Plot the various wave functions for a given value of the potential.

symRules={kRule2,nRule,symSol
		,{a->1,ℏ->1,m->1,V0->100,D->100}}//Flatten;

asymRules={kRule2,nRule,asymSol
		 ,{a->1,ℏ->1,m->1,V0->100,D->100}}//Flatten;

Clear[Ψs,Ψa]

(*Symmetric*)
Ψs[x_/;   x<-1, n0_]:= (ψL[x]//.symRules//.{n->n0});
Ψs[x_/;-1<=x<1, n0_]:= (ψW[x]//.symRules//.{n->n0});
Ψs[x_/;   x>=1, n0_]:= (ψR[x]//.symRules//.{n->n0});

(*Asymmetric*)
Ψa[x_/;   x<-1, n0_]:= (ψL[x]//.asymRules//.{n->n0});
Ψa[x_/;-1<=x<1, n0_]:= (ψW[x]//.asymRules//.{n->n0});
Ψa[x_/;   x>=1, n0_]:= (ψR[x]//.asymRules//.{n->n0});

symEnergy = nValues[eq6,100,#]& /@ symGuess

plotSym = Plot[Ψs[x,#]//.symRules//Evaluate
			,{x,-2,2}
			,PlotLabel->N[#,3]
			,Axes->None
			,DisplayFunction->Identity]& /@ symEnergy;

Show[GraphicsArray[{plotSym[[{1,2,3}]],plotSym[[{4,5}]]}]
			,DisplayFunction->$DisplayFunction];

asymEnergy = nValues[eq7,100,#]& /@ asymGuess

plotAsym = Plot[Ψa[x,#]//.asymRules//Evaluate
			,{x,-2,2}
			,PlotLabel->N[#,3]
			,Axes->None
			,DisplayFunction->Identity]& /@ asymEnergy;

Show[GraphicsArray[{plotAsym[[{1,2}]],plotAsym[[{3,4}]]}]
			,DisplayFunction->$DisplayFunction];