From dab9e9f1df53e7dacad92f88ae9d0ab0e7c0cebb Mon Sep 17 00:00:00 2001 From: scemama Date: Thu, 4 Feb 2021 12:16:50 +0000 Subject: [PATCH] deploy: 8a8ae73f0360dbbe50604685f4c832b649603ee8 --- index.html | 511 ++++++++++++++++++++++++++++++++++++++++------------- 1 file changed, 385 insertions(+), 126 deletions(-) diff --git a/index.html b/index.html index 040badc..196cd95 100644 --- a/index.html +++ b/index.html @@ -3,7 +3,7 @@ "http://www.w3.org/TR/xhtml1/DTD/xhtml1-strict.dtd"> - + Quantum Monte Carlo @@ -329,92 +329,116 @@ for the JavaScript code in this tag.

Table of Contents

-
-

1 Introduction

+
+

1 Introduction

This website contains the QMC tutorial of the 2021 LTTC winter school @@ -454,8 +478,8 @@ coordinates, etc).

-
-

1.1 Energy and local energy

+
+

1.1 Energy and local energy

For a given system with Hamiltonian \(\hat{H}\) and wave function \(\Psi\), we define the local energy as @@ -538,8 +562,8 @@ energy computed over these configurations:

-
-

2 Numerical evaluation of the energy of the hydrogen atom

+
+

2 Numerical evaluation of the energy of the hydrogen atom

In this section, we consider the hydrogen atom with the following @@ -568,8 +592,8 @@ To do that, we will compute the local energy and check whether it is constant.

-
-

2.1 Local energy

+
+

2.1 Local energy

You will now program all quantities needed to compute the local energy of the H atom for the given wave function. @@ -596,8 +620,8 @@ to catch the error.

-
-

2.1.1 Exercise 1

+
+

2.1.1 Exercise 1

@@ -641,10 +665,51 @@ and returns the potential.

+ +
+
2.1.1.1 Solution   solution2
+
+

+Python +

+
+
#!/usr/bin/env python3
+import numpy as np
+
+def potential(r):
+    distance = np.sqrt(np.dot(r,r))
+    assert (distance > 0)
+    return -1. / distance
+
-
-

2.1.2 Exercise 2

+

+Fortran +

+
+
double precision function potential(r)
+  implicit none
+  double precision, intent(in) :: r(3)
+
+  double precision             :: distance
+
+  distance = dsqrt( r(1)*r(1) + r(2)*r(2) + r(3)*r(3) )
+
+  if (distance > 0.d0) then
+     potential = -1.d0 / distance
+  else
+     stop 'potential at r=0.d0 diverges'
+  end if
+
+end function potential
+
+
+
+
+
+ +
+

2.1.2 Exercise 2

@@ -678,10 +743,37 @@ input arguments, and returns a scalar.

+ +
+
2.1.2.1 Solution   solution2
+
+

+Python +

+
+
def psi(a, r):
+    return np.exp(-a*np.sqrt(np.dot(r,r)))
+
-
-

2.1.3 Exercise 3

+

+Fortran +

+
+
double precision function psi(a, r)
+  implicit none
+  double precision, intent(in) :: a, r(3)
+
+  psi = dexp(-a * dsqrt( r(1)*r(1) + r(2)*r(2) + r(3)*r(3) ))
+end function psi
+
+
+
+
+
+ +
+

2.1.3 Exercise 3

@@ -761,10 +853,51 @@ Therefore, the local kinetic energy is

+ +
+
2.1.3.1 Solution   solution2
+
+

+Python +

+
+
def kinetic(a,r):
+    distance = np.sqrt(np.dot(r,r))
+    assert (distance > 0.)
+
+    return a * (1./distance - 0.5 * a)
+
-
-

2.1.4 Exercise 4

+

+Fortran +

+
+
double precision function kinetic(a,r)
+  implicit none
+  double precision, intent(in) :: a, r(3)
+
+  double precision             :: distance
+
+  distance = dsqrt( r(1)*r(1) + r(2)*r(2) + r(3)*r(3) ) 
+
+  if (distance > 0.d0) then
+
+     kinetic =  a * (1.d0 / distance - 0.5d0 * a)
+
+  else
+     stop 'kinetic energy diverges at r=0'
+  end if
+
+end function kinetic
+
+
+
+
+
+ +
+

2.1.4 Exercise 4

@@ -822,23 +955,73 @@ are calling is yours.

+ +
+
2.1.4.1 Solution   solution2
+
+

+Python +

+
+
def e_loc(a,r):
+    return kinetic(a,r) + potential(r)
+
-
-

2.1.5 Exercise 5

+

+Fortran +

+
+
double precision function e_loc(a,r)
+  implicit none
+  double precision, intent(in) :: a, r(3)
+
+  double precision, external :: kinetic
+  double precision, external :: potential
+
+  e_loc = kinetic(a,r) + potential(r)
+
+end function e_loc
+
+
+
+
+
+ +
+

2.1.5 Exercise 5

Find the theoretical value of \(a\) for which \(\Psi\) is an eigenfunction of \(\hat{H}\).

+
+
+ +
+
2.1.5.1 Solution   solution2
+
+\begin{eqnarray*} +E &=& \frac{\hat{H} \Psi}{\Psi} = - \frac{1}{2} \frac{\Delta \Psi}{\Psi} - +\frac{1}{|\mathbf{r}|} \\ + &=& -\frac{1}{2}\left(a^2 - \frac{2a}{\mathbf{|r|}} \right) - +\frac{1}{|\mathbf{r}|} \\ + &=& +-\frac{1}{2} a^2 + \frac{a-1}{\mathbf{|r|}} +\end{eqnarray*} + +

+\(a=1\) cancels the \(1/|r|\) term, and makes the energy constant and +equal to -0.5 atomic units. +

-
-

2.2 Plot of the local energy along the \(x\) axis

+
+

2.2 Plot of the local energy along the \(x\) axis

The program you will write in this section will be written in @@ -869,8 +1052,8 @@ In Fortran, you will need to compile all the source files together:

-
-

2.2.1 Exercise

+
+

2.2.1 Exercise

@@ -963,11 +1146,87 @@ plot './data' index 0 using 1:2 with lines title 'a=0.1', \

+ +
+
2.2.1.1 Solution   solution2
+
+

+Python +

+
+
#!/usr/bin/env python3
+
+import numpy as np
+import matplotlib.pyplot as plt
+
+from hydrogen import e_loc
+
+x=np.linspace(-5,5)
+plt.figure(figsize=(10,5))
+
+for a in [0.1, 0.2, 0.5, 1., 1.5, 2.]:
+  y=np.array([ e_loc(a, np.array([t,0.,0.]) ) for t in x])
+  plt.plot(x,y,label=f"a={a}")
+
+plt.tight_layout()
+plt.legend()
+plt.savefig("plot_py.png")
+
+
+ + +
+

plot_py.png +

+
+ +

+Fortran +

+
+
program plot
+  implicit none
+  double precision, external :: e_loc
+
+  double precision :: x(50), energy, dx, r(3), a(6)
+  integer :: i, j
+
+  a = (/ 0.1d0, 0.2d0, 0.5d0, 1.d0, 1.5d0, 2.d0 /)
+
+  dx = 10.d0/(size(x)-1)
+  do i=1,size(x)
+     x(i) = -5.d0 + (i-1)*dx
+  end do
+
+  r(:) = 0.d0
+
+  do j=1,size(a)
+     print *, '# a=', a(j)
+     do i=1,size(x)
+        r(1) = x(i)
+        energy = e_loc( a(j), r )
+        print *, x(i), energy
+     end do
+     print *, ''
+     print *, ''
+  end do
+
+end program plot
+
+
+ + +
+

plot.png +

+
+
+
-
-

2.3 Numerical estimation of the energy

+
+

2.3 Numerical estimation of the energy

If the space is discretized in small volume elements \(\mathbf{r}_i\) @@ -997,8 +1256,8 @@ The energy is biased because:

-
-

2.3.1 Exercise

+
+

2.3.1 Exercise

@@ -1071,8 +1330,8 @@ To compile the Fortran and run it:

-
-

2.4 Variance of the local energy

+
+

2.4 Variance of the local energy

The variance of the local energy is a functional of \(\Psi\) @@ -1099,8 +1358,8 @@ energy can be used as a measure of the quality of a wave function.

-
-

2.4.1 Exercise (optional)

+
+

2.4.1 Exercise (optional)

@@ -1111,8 +1370,8 @@ Prove that :

-
-

2.4.2 Exercise

+
+

2.4.2 Exercise

@@ -1191,8 +1450,8 @@ To compile and run:

-
-

3 Variational Monte Carlo

+
+

3 Variational Monte Carlo

Numerical integration with deterministic methods is very efficient @@ -1208,8 +1467,8 @@ interval.

-
-

3.1 Computation of the statistical error

+
+

3.1 Computation of the statistical error

To compute the statistical error, you need to perform \(M\) @@ -1249,8 +1508,8 @@ And the confidence interval is given by

-
-

3.1.1 Exercise

+
+

3.1.1 Exercise

@@ -1292,8 +1551,8 @@ input array.

-
-

3.2 Uniform sampling in the box

+
+

3.2 Uniform sampling in the box

We will now perform our first Monte Carlo calculation to compute the @@ -1354,8 +1613,8 @@ compute the statistical error.

-
-

3.2.1 Exercise

+
+

3.2.1 Exercise

@@ -1459,8 +1718,8 @@ well as the index of the current step.

-
-

3.3 Metropolis sampling with \(\Psi^2\)

+
+

3.3 Metropolis sampling with \(\Psi^2\)

We will now use the square of the wave function to sample random @@ -1579,8 +1838,8 @@ All samples should be kept, from both accepted and rejected moves.

-
-

3.3.1 Optimal step size

+
+

3.3.1 Optimal step size

If the box is infinitely small, the ratio will be very close @@ -1615,8 +1874,8 @@ the same variable later on to store a time step.

-
-

3.3.2 Exercise

+
+

3.3.2 Exercise

@@ -1727,8 +1986,8 @@ Can you observe a reduction in the statistical error?

-
-

3.4 Generalized Metropolis algorithm

+
+

3.4 Generalized Metropolis algorithm

One can use more efficient numerical schemes to move the electrons by choosing a smarter expression for the transition probability. @@ -1849,8 +2108,8 @@ The algorithm of the previous exercise is only slighlty modified as:

-
-

3.4.1 Gaussian random number generator

+
+

3.4.1 Gaussian random number generator

To obtain Gaussian-distributed random numbers, you can apply the @@ -1914,8 +2173,8 @@ In Python, you can use the -

3.4.2 Exercise 1

+
+

3.4.2 Exercise 1

@@ -1958,8 +2217,8 @@ Write a function to compute the drift vector \(\frac{\nabla \Psi(\mathbf{r})}{\P

-
-

3.4.3 Exercise 2

+
+

3.4.3 Exercise 2

@@ -2058,8 +2317,8 @@ Modify the previous program to introduce the drift-diffusion scheme.

-
-

4 Diffusion Monte Carlo

+
+

4 Diffusion Monte Carlo

As we have seen, Variational Monte Carlo is a powerful method to @@ -2076,8 +2335,8 @@ finding a near-exact numerical solution to the Schrödinger equation.

-