\documentclass[11pt]{article}
\usepackage[T1]{fontenc}
\usepackage{graphicx}
\hoffset=-1.3cm
\voffset=-3.01cm
\textheight=26.9cm
\textwidth=14.5cm
\topsep=0.3cm
\newcommand{\simgeq}{\; \raisebox{-0.4ex}{\tiny$\stackrel
{{\textstyle>}}{\sim}$}\;}
\newcommand{\simleq}{\; \raisebox{-0.4ex}{\tiny$\stackrel
{{\textstyle<}}{\sim}$}\;}
\newcommand{\beq}{\begin{equation}}
\newcommand{\beqar}{\begin{eqnarray}}
\newcommand{\eeq}[1]{\label{#1} \end{equation}}
\newcommand{\eeqar}[1]{\label{#1} \end{eqnarray}}
\pagestyle{empty}

\begin{document}
\noindent
{\large {\bf Beräkningsverktyg, lösning till övn. 3.4c i `Kvantvärldens fenomen'}}\\
\vspace{-3cm}

\parbox[t]{0.45\textwidth}
{\noindent
In regions 1 and 3:
\[
\frac{-\hbar^2}{2m}\phi'' + 0\phi = E\phi \;\;\;\;\Longrightarrow\;\;\;\;
\]
\[
\phi'' = - \frac{2mE}{\hbar^2} \phi = -k^2 \phi\;\;\;\; 
(k^2 = \frac{2mE}{\hbar^2})
\]
}
%\vspace{20mm}
\hspace{0.05\textwidth}
%\parbox[b]{70mm}
\parbox{0.35\textwidth}
{\vspace{35mm}
\centering
{\includegraphics[clip=true,width=0.32\textwidth]{well.eps}}
}

%In regions 1 and 3:
%\[
%\frac{-\hbar^2}{2m}\phi'' + 0\phi = E\phi \;\;\;\;\Longrightarrow\;\;\;\;
%\phi'' = - \frac{2mE}{\hbar^2} \phi = -k^2 \phi\;\;\;\; 
%(k^2 = \frac{2mE}{\hbar^2})
%\]
\noindent
Solutions of the form (with no returned wave in region 3):
\begin{eqnarray*}
\phi_1 & = & A e^{ikx} + B e^{-ikx} \;\;\;\;\Longrightarrow\;\;\;\;
\phi_1' = A ik e^{ikx} + B(-ik) e^{-ikx} \\
\phi_3 & = & F e^{ikx}  \;\;\;\;\Longrightarrow\;\;\;\;
\phi_1' = F (ik) e^{ikx} 
\end{eqnarray*}
In region 2:
\[
\frac{-\hbar^2}{2m}\phi'' - V_0\phi = E\phi \;\;\;\;\Longrightarrow\;\;\;\;
\phi'' = - \frac{2m(E+V_0}{\hbar^2} \phi = -K^2 \phi\;\;\;\; 
(k^2 = \frac{2m(E+V_0)}{\hbar^2})
\]
with the same type of solution as in regions 1 and 3:
\[
\phi_2  =  C e^{iKx} + D e^{-iKx} \;\;\;\;\Longrightarrow\;\;\;\;
\phi_2' = C (iK) e^{iKx} + D(-iK) e^{-iKx} \\
\]
Continuity of the wave-function and its first derivative at $x=0$ and $x=a$:
\begin{eqnarray*}
A+B & = & C+D \\
kA -kB & = & kC - kD \\
C e^{iKa} + D e^{-iKa} & = & F e^{ika} \\
C (iK) e^{iKa} + D(-iK) e^{-iKa} & = & F (ik) e^{ika}
\end{eqnarray*}
With the amplitudes in units of $A$
\[
\left\{\begin{array}{lllll}
+\;(B/A)   & -\; (C/A)             & -\;(D/A)              &                  & =-1 \\
-\;k(B/A) & -\; K(C/A)          & +\;K(D/A)           &                  & = -k \\
        & +\; e^{iKa}(C/A)    & +\; e^{-iKa}(D/A)   & -\; e^{ika}(F/A)   & = 0 \\
        & +\; iK e^{iKa}(C/A) & -\;iK e^{-iKa}(C/A) & -\;ik e^{ika}(E/A) & = 0
\end{array}
\right.
\]
corresponding to the system 
\[
\mathbf{A}\mathbf{x} = \mathbf{y}
\]
with
\[
\mathbf{A} = \left[\begin{array}{cccc}
1   &  -1             & -1            &  0       \\
-\;k & -\; K          & K           &        0     \\
   0     &  e^{iKa}    &  e^{-iKa}   & -e^{ika}   \\
   0     & +\; K e^{iKa}) & -\;K e^{-iKa} & -\;k e^{ika}
\end{array}
\right],
\;\;\;
\mathbf{x} = \left[\begin{array}{c}
B/A \\ C/A \\ D/A \\ F/A
\end{array}
\right],
\;\;\;
\mathbf{y} = \left[\begin{array}{c}
-1 \\ k \\ 0 \\ 0
\end{array}
\right]
\]
The probability current of the incoming wave in region 1 is given as
$(\hbar k/m)|A|^2$ and for the outcoming wave in region 3 $(\hbar
k/m)|F|^2$. This gives the transmission coefficient:
\[
T = T(E) = \frac{(\hbar k/m)|F|^2}{(\hbar k/m)|A|^2} = \left|\frac{F}{A}\right|^2
\]
Constants:
\[
k = \sqrt{\frac{2mE}{\hbar^2}} = \sqrt{\frac{2mc^2E}{(\hbar c)^2}}
\]
With $\hbar c = 197.3271$ eV nm and $m_e c^2 = 510.999 10^3$ eV and with
$E$ in eV, $k$ will be in units of [1/nm]. The only difference for $K$
compared with $k$ is that $E$ is replaced by $(E+V_0)$, where $V_0 = 5$ eV.  

It is possible to find an analytic solution (see p. 201, G. Ohl\'en,
Kvantvärldens fenomen:
\[
T= \left(1 + \frac{V_0}{4E\left(E+V_0\right)} \sin^2{(ka)} \right)^{-1}
\]

\noindent
Matlab-fil för lösning av uppgiften:
\begin{verbatim}
% transmission.m  inlamningsuppgift 2, deluppgift 1, ht 2005

clear all
close all

L=50; % antal punkter i energivektorn.
E_vector = linspace(eps,10,L);
T_vector=zeros(L,1);

a = 1; % [nm]
mC2=510999; % [eV]
hbarC=197.3271; % [eV nm]
V0=5; % [eV]

for j = 1:L
    E=E_vector(j);
    K = sqrt( (2*mC2*(E+V0))/(hbarC^2) );  
    k = sqrt( (2*mC2*(E))/(hbarC^2) );
    
    y = [-1 -k 0 0]'; % ekvationssystemets högerled
    
    A = [            % matris for koefficienterna till de fyra obekanta "B, C, D, E".
        1 -1 -1 0;
        -k -K K 0;
        0 exp(i*K*a) exp(-i*K*a) -exp(i*k*a);
        0 K*exp(i*K*a) -K*exp(-i*K*a) -k*exp(i*k*a)
    ];
   % x = A\y; % losningen till ekvation systemet. 
   % x = inv(A)*y;
   x = A\y
    
    T_vector(j) = abs(x(4))^2; 
   % 4.e positionen i x vektorn är amplituden för den transmitterade vagen.
end

plot (E_vector, T_vector,'-')
hold on

% Jämförelse med analytisk beräkning
Ein_vector = linspace(eps,10,200);
k_vector=sqrt( 2*mC2*(Ein_vector+V0)/hbarC^2 );
analytic=1./(1 + V0^2./(4*Ein_vector.*(Ein_vector+V0)).* sin(k_vector*a).^2 ); 
% sidan 201 i 'Kvantvärldens fenomen'.
plot(Ein_vector,analytic,'r--')
axis([0 10 0 1.4])
legend('Numeriskt beraknad','Analytisk formel, sid. 201','Location','NW')
xlabel('Energi [e.V.]','fontsize',15)
ylabel('Transmitans','fontsize',15)
\end{verbatim}
\end{document}
