/python/README.rst

.. -*- rst -*- -*- restructuredtext -*-

.. This file should be written using the restructure text
.. conventions.  It will be displayed on the bitbucket source page and
.. serves as the documentation of the directory.

.. default-role:: math

=============================
 Python Implementation Notes
=============================

This file documents some implementation details.  Additional notes, test-cases,
and usage can be found in the documented source-code.

Bases
=====
Several different bases are provided in ``bases.py``.  We generally advocate
using the ``Fourier`` basis as this is faster (kinetic energy uses the FFTW_).
Originally the  ``Dirichlet`` basis was used to provide upper bounds, but
similar bounds can now be obtained with the ``Fourier`` basis using twisted
boundary conditions by applying ``theta_bloch=np.pi``.  This augments all
plane-waves with an overall phase `e^{i \theta_B L x} = e^{i k_B x}`.

The choice ``theta_bloch=np.pi`` thus enforces anti-symmetric boundary
conditions.  In combination with even parity of the ground state this
effectively enforces Dirichlet boundary conditions and one finds that (once UV
convergence is obtained) the same upper bound is achieved.


Triton and Alpha Particles
==========================
The directory ``shell_model`` contains Python code for using DVR bases for
a brute force solution to the few-body problem for three- and four-particle
states with simple attractive s-wave interaction.  We apply a DVR basis to the
relative Jacobi coordinates to find the three-body ("triton") and four-body
("alpha particle") bound states.  This code requires a few more components to
function efficiently – in particular, numexpr_ and anfft_ which wraps the
FFTW_.  (The code should still run, albeit more slowly, if these are not
installed.

Using this code, one can find the lowest energy states using a simple Lanczos_
iteration with lattices of up to `18\times 18\times 18` for the "triton" and
`8\times 8 \times 8` for the "alpha particle".  These correspond to
wavefunctions with `18^6\approx 10^8` and `8^9\approx 10^8` components (1GB per
wavefunction) respectively, requiring about 12GB of RAM, and can be solved in
less than an hour on a modern laptop.  (A more careful implementation caching
matrices on disk and using an in-place FFT could extend these to about `10^9`
components.)  This is a fairly straightforward implementation: we have not
attempted to use any special features – such a symmetries of the ground state –
to reduce the size of the basis.

The python code is organized into the following files:

``shell_model.py``:
   The main code defining the units, problem, Hamiltonian etc.
``dvr.py``:
   Various 1D DVR basis sets including Dirichlet and Sinc bases.
``bases.py``:
   Provides full basis sets for the few-body problem constructed from the DVR
   and Fourier bases.  Also defines the ``RelCoord`` class which specifies the
   appropriate transformations for converting relative Jacobi coordinates to
   absolute coordinates for the three- and four-body systems.
``solvers.py``:
   Code for the Lanczos_ solver that finds the lowest few states of the
   Hamiltonian.  Note: this is a simple implementation that does not implement
   re-orthogonalization schemes (to save memory) but works sufficiently well for
   the problems considered here.
``utils.py``:
   Tries to implement numexpr_ and anfft_, warning the user about potential
   performance problems if these are not provided.
``figures.py``:
   Some code for automating the generation of energies used in the paper.


.. _DVR_Demo.py: \
   https://bitbucket.org/mforbes/paper_dvrvsho/src/tip/DVR_Demo.py
.. _plots.py: \
   https://bitbucket.org/mforbes/paper_dvrvsho/src/tip/plots.py
.. _View the DVR_Demo Notebook: \
   http://nbviewer.ipython.org/urls/bitbucket.org/mforbes/paper_dvrvsho/raw/tip/DVR_Demo.ipynb
.. |EPD| replace:: Enthough Python Distribution (EPD_)
.. _EPD: http://www.enthought.com/products/epd.php
.. _IPython: http://ipython.org/
.. _Ipython notebook: \
   http://ipython.org/ipython-doc/dev/interactive/htmlnotebook.html
.. _SciPy: http://www.scipy.org/
.. _NumPy: http://numpy.scipy.org/
.. _matplotlib: http://matplotlib.org/
.. _SymPy: http://sympy.org/
.. _IPython Notebook Viewer: http://nbviewer.ipython.org/
.. _pymmf: https://bitbucket.org/mforbes/pymmf
.. _numexpr: http://code.google.com/p/numexpr/
.. _anfft: https://code.google.com/p/anfft/
.. _FFTW: http://www.fftw.org
.. _Lanczos: http://en.wikipedia.org/wiki/Lanczos_algorithm