How can I animate quantum tunneling using Python and pseudocode?

[Avatrin]

Hi

Let's say I am supposed to animate the movement of some one-dimensional wavefunction as it hits a finite potential barrier like this:


I was supposed to know how to do this over a year ago, but I am never certain how to write the algorithm. What confuses me the most is how to show the incident wave "split into" the transmitted and reflected waves.

So, let's say I am given some wavefunction Ψ(x,t) that is moving towards a potential barrier or well. How should I write the algorithm? Is there some book you guys recommend that teaches this?

Also, I use Python. So, I would prefer any example code to be written in either Python or pseudocode. I do the plotting and animation with Scitools and Easyviz (I will take any recommendations if somebody recommend something else).

 

[vanhees71]

I did this for a QM 1 lecture once. I just did the brute-force Fourier decomposition of a wave packet as the one shown in your youtube stream. You just solve the energy-eigenvalue problem (time-independent Schrödinger equation) for your potential well and then Fourier transform the initial wave packet in terms of energy eigenfunction. I just did the one-dimensional integral for each time with a Gauss quadrature. It's not that demanding on CPU time. So you don't need more fancy algorithms for this.

If you like to see the animations, you can find them here (however the text around it is in German):

Free particle:
http://theory.gsi.de/~vanhees/faq/quant/node33.html

total reflection on a finite potential well (step function):
http://theory.gsi.de/~vanhees/faq/quant/node35.html

with a small part running through that well:
http://theory.gsi.de/~vanhees/faq/quant/node35.html

some general case (different particle energies):
http://theory.gsi.de/~vanhees/faq/quant/node37.html
http://theory.gsi.de/~vanhees/faq/quant/node38.html
http://theory.gsi.de/~vanhees/faq/quant/node39.html

Rectangular finite-depth potential pot:

energy close to a resonance value (asymptotic final state is the wave packet nearly completely run through the well but with some time delay)
http://theory.gsi.de/~vanhees/faq/quant/node44.html

and at some other energy:
http://theory.gsi.de/~vanhees/faq/quant/node45.html

wave packet as superposition of some bound states:
http://theory.gsi.de/~vanhees/faq/quant/node46.html

 

[Avatrin]

The specific issue I am struggling with is tunneling. Animating that is what I am struggling with. The most confusing aspect is how I get two probability distributions from one. Does that come "automatically" when solving the problem numerically?

Also, almost all of my knowledge on how to solve the one-dimensional Schrodingers equation comes from the second chapter of David J. Griffiths book. So, the way you suggest that I use the Fourier transform is a little foreign to me, and I have never heard of Gaussian quadrature.

 

[HomogenousCow]

The specific issue I am struggling with is tunneling. Animating that is what I am struggling with. The most confusing aspect is how I get two probability distributions from one. Does that come "automatically" when solving the problem numerically?

Also, almost all of my knowledge on how to solve the one-dimensional Schrodingers equation comes from the second chapter of David J. Griffiths book. So, the way you suggest that I use the Fourier transform is a little foreign to me, and I have never heard of Gaussian quadrature.

You're only supposed to get one probability distribution.
The energy eigenvectors will be these piecewise function for the three regions (left, inside, and right of the barrier), they're supposed to be orthogonal so you can just sum a bunch of them and make a gaussian to the right of the barrier.

 

[vanhees71]

Of course, the wave packets shown in my animations are just solutions of the time-dependent Schrödinger equation. I'm plotting $|\psi(t,x)|^2$, i.e., the probability distribution. It's done by numerical integration of the explicit solution in terms of energy eigenfunctions. In the general case you have a discrete spectrum (usually eigenvalues $E_n<0$) and a continuous one for $E>0$. Then you take an incoming wave packet far away from the support of the potential or, in the case where the potential is non-zero also at infinity like in the case of the potential step in a region of constant potential ("asymptotic free state") and let it run towards the region, where the potential is active. This produces the scattering solutions shown in the animations. Another posibility is that a particle is trapped in a potential valley as is exemplified for the potential pot. Given the initial wave packet as $\psi(t=0,x=\psi_0(x)$, the time evolution follows from the expansion in terms of energy eigenfunctions $u_n$ in the discrete and $u_E$ in the continuous spectrum, i.e.,
$$\psi(t,x)=\sum_{n} u_n(x) \psi_n \exp(-\mathrm{i} E_n t/\hbar)+\int_0^{\infty} \mathrm{d} E \psi_E u_E(x) \exp(-\mathrm{i} E t/\hbar),$$
where
$\psi_n=\int_{\mathbb{R}} \mathrm{d} x u_n^*(x) \psi_0(x), \quad \psi_E=\int_{\mathbb{R}} \mathrm{d} x u_E^*(x) \psi_0(x),$
are the matrix elements of the initial wave function wrt. the energy eigenfunctions in the discrete and continuous spectrum of the Hamiltonian.

Gauss quadrature is just a numerical integration technique. You should find tons of such standard routines for integration/quadrature in open-source software libraries.

 

[Avatrin]

I need to resurrect this thread. The approaches above have not really helped me. The kind of exercises I am struggling with seem to be simplifications rather than accurate one-dimensional representations of the solution of Schrodingers equation. Let's construct an example:

I start off with a wave function of the form:
$$\psi (x,0) = Ae^{kx}$$
This function is going to move toward the right towards a, let's say, potential step or a potential hill. I do not know how to approach this.

In my textbook, I learn how to solve the time independent solutions to such problems, but I do not know what to do when time is a variable. One of the problems I have is that the wave is constantly "colliding" with the potential hill since the function is non-zero everywhere.

I would like recommendations for textbooks that teach this (since Griffiths do not).

 

[vanhees71]

What's unclear with my posting from Dec 28? This really works nicely on a usual personal computer. It's no big deal with a simple numerical-integration routine.

 

[Avatrin]

What's unclear with my posting from Dec 28? This really works nicely on a usual personal computer. It's no big deal with a simple numerical-integration routine.

Well, it assumes more knowledge than is required in the course I am going to take (energy eigenfunctions). So, I guess, part of it is that I am looking for an easier approach.

However, the bigger problem is that I cannot wrap my head around the time-dependent solutions of Schrödinger's equation. For instance, in the case of the finite potential well, I will have three solutions; One outside the left side of the well, one inside and one on the right side, outside, the well. Am I just supposed to multiply all of these with $e^{i\frac{E_n}{t}}$? And, run the algorithm for all of the three solutions simultaneously? I think I need to see an example...

 

[HomogenousCow]

Shankar has an excellent section on this.

 

[vanhees71]

You don't have three solutions but energy eigenstates, defined in the various regions corresponding to your potential. They are all properly matched together by the appropriate boundary conditions and finally make up one function. I don't know, what your course is supposed to teach you, but in any case you should have some idea what the energy eigenfunctions are good for, namely to be able to solve the initial-value problem for the time-dependent Schrödinger equation. After all that's what a big deal of QT is about!