ECE_456_Reports/PS2/problem_3_part_b_v.m
2021-03-05 15:25:31 -07:00

50 lines
1.3 KiB
Matlab

% A simple MATLAB code to plot the momentum
% wave function for a particle in an escape-proof "box"
clear all
hbar = 1.054e-34;
L = 101e-10;
% Set up k' and p' axes
k_prime = linspace(-2*pi/L,2*pi/L,100);
dk_prime = k_prime(2) - k_prime(1);
p_prime = hbar * k_prime;
dp_prime = p_prime(2) - p_prime(1);
% k1 and p1
k1 = (pi/L); p1 = hbar * k1;
% Fourier transform A(k')
A1 = 0.5*sqrt(L/pi)*( sinc((k_prime+k1)*L/(2*pi)) );
A2 = 0.5*sqrt(L/pi)*( sinc((k_prime-k1)*L/(2*pi)) );
A = A1 + A2;
A_sum = sum(abs(A).^2)*dk_prime;
% Momentum wave function Phi(p')
Phi1 = 0.5*sqrt(L/(pi*hbar))*( sinc((p_prime+p1)*L/(2*pi*hbar)) );
Phi2 = 0.5*sqrt(L/(pi*hbar))*( sinc((p_prime-p1)*L/(2*pi*hbar)) );
Phi = Phi1 + Phi2;
Phi_sum = sum(abs(Phi).^2)*dp_prime;
% Normalized p' axis for plotting purposes; the points of classical
% momenta occur when the variable pp_N is plus or minus unity
pp_N = p_prime / p1;
figure(1); clf;
h = plot(pp_N,Phi1,'kx',pp_N,Phi2,'ko',pp_N,Phi,'k--');
grid on; set(h,'linewidth',[2.0]); set(gca,'Fontsize',[18]);
xlabel('NORMALIZED MOMENTUM [p^{\prime}/p_1]');
ylabel('WAVE FUNCTION [sqrt(s / kg m)]');
legend('Phi1','Phi2','Phi');
figure(2); clf; h = plot(pp_N,abs(Phi).^2,'k-');
grid on; set(h,'linewidth',[2.0]); set(gca,'Fontsize',[18]);
xlabel('NORMALIZED MOMENTUM [p^{\prime}/p_1]');
ylabel('PROBABILITY DENSITY [s / kg m]');