\documentclass[12pt]{article} \textwidth 16cm \textheight 23cm \topmargin 0cm \evensidemargin 0cm \oddsidemargin 0cm \parindent=0pt \parskip=10pt \begin{document} \begin{center} \Large \bf Computational Methods and Data Analysis 2003\\ \large \end{center} \vspace{0.4cm} \noindent {\bf Ex.~4: Bound states in a finite square well}\\ This exercise requires the implementation of algorithms for finding the zero's of nonlinear equations and for numerical integration. Consider a square potential well of finite height V$_0$ and width $2a$, as in figure ~\ref{fig:put}. \begin{figure}[hbt] \centering \setlength{\unitlength}{0.00050000in}% % \begin{picture}(6066,5493)(2968,-6700) \thicklines \put(6001,-6361){\line( 0, 1){4800}} \put(3001,-2761){\line( 1, 0){1200}} \put(4201,-2761){\line( 0,-1){3600}} \put(4201,-6361){\line( 1, 0){3600}} \put(7801,-6361){\line( 0, 1){3600}} \put(7801,-2761){\line( 1, 0){1200}} \put(6001,-1561){\line( 0, 1){ 75}} \put(6001,-1486){\vector( 0, 1){ 75}} \put(6001,-6661){\small 0} \put(3901,-6661){\large $-a$} \put(7801,-6661){\large $a$} \put(6226,-1411){\large $V(x)$} \put(6301,-2761){\large $V_{0}$} \end{picture} \caption[1]{Finite potential well} \label{fig:put} \end{figure} The Schr\"odinger equation: \[ \frac{d^2\Psi}{dx^2}-\frac{2m}{\hbar^2}\left[V(x)-E\right]\Psi=0 \] with $V(x)$ given by: \begin{eqnarray*} V(x)&=&V_0 \qquad x<-a \mbox{~~~~~~~~~~~~~~~~~~~~~region $I$~~~}\\ V(x)&=&0 \qquad -aa \mbox{~~~~~~~~~~~~~~~~~~~~~~~~region $III$~} \end{eqnarray*} has general solutions for the bound states at energy less than $V_0$ in the three regions given by: \begin{eqnarray*} \Psi_{I}(x)&=&A~\exp^{\beta x}+B~\exp^{-\beta x} \qquad x<-a\\ \Psi_{II}(x)&=&C~\sin\alpha x+D~\cos\alpha x \qquad -aa\\[1ex] \alpha&=&\sqrt{\frac{2mE}{\hbar^2}}\\ \beta&=&\sqrt{2m(V_0-E)/\hbar^2} \end{eqnarray*} The boundary conditions require continuity of the wavefunction and of its derivative at the boundary, whence \\[1ex] Even states: $B=E=0, C=0, D\neq0, A=F \qquad \alpha \tan \alpha a = \beta$\\ Odd states: $B=E=0 , C\neq0, D=0, A=-F \qquad \alpha \cot \alpha a = -\beta$\\ Therefore eigenvalues of the Schr\"odinger equation are found by solving numerically: \begin{center} \begin{tabular}{ccc} $\alpha tg(\alpha a)-\beta$=0 &for even states\\ $\alpha ctg(\alpha a)+\beta$=0 &for odd states \end{tabular} \end{center} A more convenient form is: \begin{center} \begin{tabular}{ccc} $\alpha sin(\alpha a)-\beta cos(\alpha a)$=0 &for even states\\ $\alpha cos(\alpha a)+\beta sin(\alpha a)$=0 &for odd states \end{tabular} \end{center} Use $2a$=10~\AA, $m=1~m_e$. It is useful to use $\hbar^2=7.6199682 m_e$~eV~\AA$^2$ \begin{enumerate} \item Bracket the range in which solutions exist for even and odd states. This can be done by noting that $\beta/\alpha$ is always positive so that, for instance, solutions for the even state will occur in ranges of $\alpha a$ such that the tangent is positive. Remember that, in one dimension, a square well has always at least one bound state at $E\neq 0$. \item Find all bound states with precision 10$^{-4}$~eV by use of the bisection and Newton-Raphson methods for several values of $V_0$ between $V_0$=0.2~eV and $V_0$=5~eV. \item Make a plot of the eigenvalues as a function of $V_0$, comparing also to the limiting case of an infinite well: \begin{equation} E_n= {\hbar^2\over 2m_e}({n\pi\over 2a})^2. \end{equation} Give some of the data also in a Table in order to check that the results are correct within the desired precision. Do not forget to label the axis of the figures and to make clear table and figure captions, giving the used units and for which values of the parameters (e.g which value of $a$) the calculation is done. The idea is that, even without the text of the exercise, by reading the report one can understand and reproduce the calculation. Comment briefly the results. %If you use gnuplot, you could plot the eigenvalues with linespoints so that you% can trace back for which values of $V_0$ the calculations has been performed. \item Facultative: Make a plot of the eigenvalues as a function of $a$ for $V_0$=3 eV. \item Calculate (and plot) the normalized wavefunction corresponding to the first two bound states for two values of $V_0$ and comment the results. To perform the integral needed to normalize the wavefunction you need to specify the range of integration. This obviously requires an approximation, the desired range is in principle from minus to plus infinity. However, once you know the energy of a given state you can estimate the decay length of the wavefunctions in the forbidden barrier regions in order to define a sensible range for integration. Describe your reasoning in doing this. You can either perform the integrals by use of the built-in matlab functions or by constructing your own code for trapezoidal and Simpson integration. In defining the grid of points on which the wavefunction is calculated to perform the integration choose a spacing which fits into the well. \item Check numerically that, as expected from quantum mechanics, the eigenfunctions corresponding to different eigenvalues are orthogonal. \item Facultative: Integrate the wavefunction of one of the states by means of Monte Carlo integration by writing \begin{equation} \int_a^b{f(x) dx}= (b-a) {1\over N} \Sigma_{i=1}^N f(x_i) \end{equation} where $x_i$ are $N$ random numbers uniformely distributed between $a$ and $b$. Study the dependence on $N$. %\item Could you find out which is the value of $V_0$ which gives the first 3 %bound states at energy $E_1$=0.26230 $E_2$=1.03633 $E_3$= 2.27006 eV? \end{enumerate} \end{document}