**Spring 2015**

Be able to apply two central methods, Variational Monte Carlo and Hartree-Fock (and density functional theory) to properties of atoms and molecules.

Understand how to simulate qauntum mechanical systems with many interacting particles. The methods are relevant for atomic, molecular, solid state, materials science, nanotechnology, quantum chemistry and nuclear physics.

Monday 9.15am-12pm: First lecture: Presentation of the course, aims and content

Monday 9.15am-12pm: Introduction to Monte Carlo methods and basic many-body physics

Tuesday 9.15am-12pm: Basic many-body physics and variational Monte Carlo methods

Wednesday 9.15am-12pm: Variational Monte Carlo methods and presentation of project 1

Thursday 2.15pm-4pm: Continued discussion of first project, variational Monte Carlo methods and codes

Computer lab: Thursday 4pm-7pm. First time: Thursday this week at room FV329.

Lectures: Thursday (2.15pm-4pm), remotely via adove connect every second week.

Computerlab: Thursday (2.15pm-7pm), first time January 22, last lab session May 14.

Weekly plans and all other information are on the webpage, see URL: ''

Second intensive week starts April 7 and ends April 10.

First project to be handed in February 27 at noon.

Second project to be handed in March 27 at noon.

Third and final project to be handed in June 1 at noon.

Three compulsory projects. Electronic reports only. You are free to choose your format.

Evaluation and grading: The first two projects count 25% while the last project counts 50% of the final grade. There is no oral or written exam.

The computer lab (room FV329) consists of 16 Linux PCs, but many prefer own laptops. C/C++ is the default programming language, but Fortran2008 and Python are also used. All source codes discussed during the lectures can be found at the webpage of the course. We recommend either C/C++, Fortran2008 or Python as languages.

Parallelization (MPI and OpenMP), high-performance computing topics. Choose between Fortran2008 and/or C++ as programming languages. Python also possible as programming language.

Algorithms for Monte Carlo Simulations (multidimensional integrals), Metropolis-Hastings and importance sampling algorithms. Improved Monte Carlo methods.

Statistical analysis of data from Monte Carlo calculations, blocking method. (exercise 1 of part 1)

Eigenvalue solvers, efficient computations of integrals

Search for minima in multidimensional spaces (conjugate gradient method, steepest descent method, quasi-Newton-Raphson, Broyden-Jacobian).

Iterative methods for solutions of non-linear equations.

Object orientation

Variational Monte Carlo for 'ab initio' studies of quantum mechanical many-body systems.

Simulation of three-dimensional systems like atoms or molecules.

Hartree-Fock method to study atoms and molecules

Most quantum mechanical problems of interest in for example atomic, molecular, nuclear and solid state physics consist of a large number of interacting electrons and ions or nucleons.

The total number of particles $N$ is usually sufficiently large that an exact solution cannot be found.

Typically, the expectation value for a chosen hamiltonian for a system of $N$ particles is

$$
\langle H \rangle =
\frac{\int d\boldsymbol{R}_1d\boldsymbol{R}_2\dots d\boldsymbol{R}_N
\Psi^{\ast}(\boldsymbol{R_1},\boldsymbol{R}_2,\dots,\boldsymbol{R}_N)
H(\boldsymbol{R_1},\boldsymbol{R}_2,\dots,\boldsymbol{R}_N)
\Psi(\boldsymbol{R_1},\boldsymbol{R}_2,\dots,\boldsymbol{R}_N)}
{\int d\boldsymbol{R}_1d\boldsymbol{R}_2\dots d\boldsymbol{R}_N
\Psi^{\ast}(\boldsymbol{R_1},\boldsymbol{R}_2,\dots,\boldsymbol{R}_N)
\Psi(\boldsymbol{R_1},\boldsymbol{R}_2,\dots,\boldsymbol{R}_N)},
$$

an in general intractable problem.

This integral is actually the starting point in a Variational Monte Carlo calculation. **Gaussian quadrature: Forget it**! Given 10 particles and 10 mesh points for each degree of freedom
and an
ideal 1 Tflops machine (all operations take the same time), how long will it take to compute the above integral? The lifetime of the universe is of the order of $10^{17}$ s.

As an example from the nuclear many-body problem, we have Schroedinger's equation as a differential equation

$$
\hat{H}\Psi(\boldsymbol{r}_1,..,\boldsymbol{r}_A,\alpha_1,..,\alpha_A)=E\Psi(\boldsymbol{r}_1,..,\boldsymbol{r}_A,\alpha_1,..,\alpha_A)
$$

where

$$
\boldsymbol{r}_1,..,\boldsymbol{r}_A,
$$

are the coordinates and

$$
\alpha_1,..,\alpha_A,
$$

$$
2^A\times \left(\begin{array}{c} A\\ Z\end{array}\right)
$$

coupled second-order differential equations in $3A$ dimensions.

For a nucleus like ${}^{10}\mbox{Be}$ this number is **215040**.
This is a truely challenging many-body problem.

Methods like partial differential equations can at most be used for 2-3 particles.

Monte-Carlo methods

Renormalization group (RG) methods, in particular density matrix RG

Large-scale diagonalization (Iterative methods, Lanczo's method, dimensionalities $10^{10}$ states)

Coupled cluster theory, favoured method in quantum chemistry, molecular and atomic physics. Applications to ab initio calculations in nuclear physics as well for large nuclei.

Perturbative many-body methods

Green's function methods

Density functional theory/Mean-field theory and Hartree-Fock theory

The physics of the system hints at which many-body methods to use.

**Pros and Cons of Monte Carlo.**

Is physically intuitive.

Allows one to study systems with many degrees of freedom. Diffusion Monte Carlo (DMC) and Green's function Monte Carlo (GFMC) yield in principle the exact solution to Schroedinger's equation.

Variational Monte Carlo (VMC) is easy to implement but needs a reliable trial wave function, can be difficult to obtain. This is where we will use Hartree-Fock theory to construct an optimal basis.

DMC/GFMC for fermions (spin with half-integer values, electrons, baryons, neutrinos, quarks) has a sign problem. Nature prefers an anti-symmetric wave function. PDF in this case given distribution of random walkers ($p\ge 0$).

The solution has a statistical error, which can be large.

There is a limit for how large systems one can study, DMC needs a huge number of random walkers in order to achieve stable results.

Obtain only the lowest-lying states with a given symmetry. Can get excited states.

**Where and why do we use Monte Carlo Methods in Quantum Physics.**

Quantum systems with many particles at finite temperature: Path Integral Monte Carlo with applications to dense matter and quantum liquids (phase transitions from normal fluid to superfluid). Strong correlations.

Bose-Einstein condensation of dilute gases, method transition from non-linear PDE to Diffusion Monte Carlo as density increases.

Light atoms, molecules, solids and nuclei.

Lattice Quantum-Chromo Dynamics. Impossible to solve without MC calculations.

Simulations of systems in solid state physics, from semiconductors to spin systems. Many electrons active and possibly strong correlations.

Given a hamiltonian $H$ and a trial wave function $\Psi_T$, the variational principle states that the expectation value of $\langle H \rangle$, defined through

$$
E[H]= \langle H \rangle =
\frac{\int d\boldsymbol{R}\Psi^{\ast}_T(\boldsymbol{R})H(\boldsymbol{R})\Psi_T(\boldsymbol{R})}
{\int d\boldsymbol{R}\Psi^{\ast}_T(\boldsymbol{R})\Psi_T(\boldsymbol{R})},
$$

is an upper bound to the ground state energy $E_0$ of the hamiltonian $H$, that is

$$
E_0 \le \langle H \rangle .
$$

In general, the integrals involved in the calculation of various expectation values are multi-dimensional ones. Traditional integration methods such as the Gauss-Legendre will not be adequate for say the computation of the energy of a many-body system.

The trial wave function can be expanded in the eigenstates of the hamiltonian since they form a complete set, viz.,

$$
\Psi_T(\boldsymbol{R})=\sum_i a_i\Psi_i(\boldsymbol{R}),
$$

and assuming the set of eigenfunctions to be normalized one obtains

$$
\frac{\sum_{nm}a^*_ma_n \int d\boldsymbol{R}\Psi^{\ast}_m(\boldsymbol{R})H(\boldsymbol{R})\Psi_n(\boldsymbol{R})}
{\sum_{nm}a^*_ma_n \int d\boldsymbol{R}\Psi^{\ast}_m(\boldsymbol{R})\Psi_n(\boldsymbol{R})} =\frac{\sum_{n}a^2_n E_n}
{\sum_{n}a^2_n} \ge E_0,
$$

where we used that $H(\boldsymbol{R})\Psi_n(\boldsymbol{R})=E_n\Psi_n(\boldsymbol{R})$. In general, the integrals involved in the calculation of various expectation values are multi-dimensional ones. The variational principle yields the lowest state of a given symmetry.

In most cases, a wave function has only small values in large parts of configuration space, and a straightforward procedure which uses homogenously distributed random points in configuration space will most likely lead to poor results. This may suggest that some kind of importance sampling combined with e.g., the Metropolis algorithm may be a more efficient way of obtaining the ground state energy. The hope is then that those regions of configurations space where the wave function assumes appreciable values are sampled more efficiently.

The tedious part in a VMC calculation is the search for the variational minimum. A good knowledge of the system is required in order to carry out reasonable VMC calculations. This is not always the case, and often VMC calculations serve rather as the starting point for so-called diffusion Monte Carlo calculations (DMC). DMC is a way of solving exactly the many-body Schroedinger equation by means of a stochastic procedure. A good guess on the binding energy and its wave function is however necessary. A carefully performed VMC calculation can aid in this context.

- Construct first a trial wave function $\psi_T(\boldsymbol{R},\boldsymbol{\alpha})$, for a many-body system consisting of $N$ particles located at positions

$\boldsymbol{R}=(\boldsymbol{R}_1,\dots ,\boldsymbol{R}_N)$. The trial wave function depends on $\alpha$ variational parameters $\boldsymbol{\alpha}=(\alpha_1,\dots ,\alpha_M)$.

- Then we evaluate the expectation value of the hamiltonian $H$

$$
E[H]=\langle H \rangle =
\frac{\int d\boldsymbol{R}\Psi^{\ast}_{T}(\boldsymbol{R},\boldsymbol{\alpha})H(\boldsymbol{R})\Psi_{T}(\boldsymbol{R},\boldsymbol{\alpha})}
{\int d\boldsymbol{R}\Psi^{\ast}_{T}(\boldsymbol{R},\boldsymbol{\alpha})\Psi_{T}(\boldsymbol{R},\boldsymbol{\alpha})}.
$$

$$
P(\boldsymbol{R})= \frac{\left|\psi_T(\boldsymbol{R})\right|^2}{\int \left|\psi_T(\boldsymbol{R})\right|^2d\boldsymbol{R}}.
$$

$$
E[H(\boldsymbol{\alpha})] =
\frac{\int d\boldsymbol{R}\Psi^{\ast}_T(\boldsymbol{R},\boldsymbol{\alpha})H(\boldsymbol{R})\Psi_T(\boldsymbol{R},\boldsymbol{\alpha})}
{\int d\boldsymbol{R}\Psi^{\ast}_T(\boldsymbol{R},\boldsymbol{\alpha})\Psi_T(\boldsymbol{R},\boldsymbol{\alpha})}.
$$

$$
E_L(\boldsymbol{R},\boldsymbol{\alpha})=\frac{1}{\psi_T(\boldsymbol{R},\boldsymbol{\alpha})}H\psi_T(\boldsymbol{R},\boldsymbol{\alpha}),
\label{eq:locale1} \tag{1}
$$

called the local energy, which, together with our trial PDF yields

$$
E[H(\boldsymbol{\alpha})]= = \int P(\boldsymbol{R})E_L(\boldsymbol{R}) d\boldsymbol{R}\approx \frac{1}{N}\sum_{i=1}^NP(\boldsymbol{R_i},\boldsymbol{\alpha})E_L(\boldsymbol{R_i},\boldsymbol{\alpha})
\label{eq:vmc1} \tag{2}
$$

with $N$ being the number of Monte Carlo samples.