Chapt6Exercise4c.zip_mechanics_quantum_quantum mechanics

%Chapt6Exercise4c.m %simple harmonic oscillator wave function psi_n(xi)and %modulus wavefunction squared |psi_n(xi)|^2 %xi=(m*w/hbar)*x %using relation %psi_n=((2/n)^0.5)*((xi*psi_n-1)-((((n-1)/2)^0.5)*psi_n-2)) clear; clf; npoints=500; %(number of data points in plot) - 1 nlim=100; %arbitrary limit to value of n that can be plotted n=input('Input quantum index n = ');%Input from keyboard value of quantum index n %of wave function to be plotted (n=0,1,2,3, ... ) if n < 0; error('minimum value of n must be greater or equal to 0'); end; if n > nlim; error('maximum value of n must be less than or equal to 100'); end; ximax=sqrt((2*n)+1); %classical turning point xiplot=(3/((n+1)^(1/3)))+(1.2*ximax); %plot range of x-axis deltaxi=2*xiplot/npoints; %increment in xi Ao=(1/pi)^(1/4); %Normalization amplitude for n=0 An=1.1*Ao/((n+1)^0.1); %fix vertical scale xi(1)=-xiplot; %first value of xi for j=2:1:npoints+1 xi(j)=xi(j-1)+deltaxi; psi1(j)=Ao*exp((-xi(j)^2)/2); %known n=0 ground-state wave function psi2(j)=(sqrt(2))*xi(j)*psi1(j); %known n=1 first excited-state wave function psi(j)=psi2(j); end if n < 1; psi=psi1; end; if n>=2 for ni=2:1:n nin=double(ni); for j=2:1:npoints+1 psi(j)=(sqrt(2/nin))*((xi(j)*psi2(j))-((sqrt((nin-1)/2))*psi1(j))); psi1(j)=psi2(j); %update new value of psi_(n-2) psi2(j)=psi(j); %update new value of psi_(n-1) end end end figure(1); subplot(1,2,1),plot(xi,psi); axis([-xiplot,xiplot,-An,An]); xlabel('Position, \xi (m)'),ylabel('Wave function, \psi(\xi)'); ttl = sprintf('Chapt6Exercise4c, xi =(m*omega/hbar)*x, n=%3.0f',n); title (ttl); subplot(1,2,2),plot(xi,abs(psi.^2)); axis([-xiplot,xiplot,0,An^2]); xlabel('Position, \xi (m)'); ylabel('Modulus wave function squared, |\psi(\xi)| ^2'); ttl2=sprintf('Classical turning-point=+/-%5.3f',ximax); title (ttl2);