Phys690: Homework Assignment 3

3/2

==> Due on Mon 3/19 in class <==


1-d harmonic oscillator in thermal equilibrium.

Classical statistical mechanics tells us that equilibrium properties of a system are determined by the so-called partition function exp(-V(x)/kT), where V(x) is the total potential energy of configuration x. We will work in units such that kT=1. Consider a particle in a one-dimensional harmonic potential V(x)=x^2/2; here x is simply the one-dimensional coordinate of the particle. The average potential energy of the particle is given by

        INT_{-infty, infty} [ x^2/2 exp(-x^2/2) dx ]
  Ev = ---------------------------------------------- .
          INT_{-infty, infty} [ exp(-x^2/2) dx ]
(INT_{-infty, infty} [....] stands for ``integrating [....] from minus infinity to plus infinity''.)

Write a program to calculate this integral by the Metropolis algorithm. Equilibrate the random walk for several hundred steps before collecting samples. After equilibration, divide the random walk into n_blk blocks of n_smpls each, like in HW2. Compute Ev, (including error bar) and also estimate an average acceptance ratio from the program. Start the random walk in any way you think is reasonable.

  1. Start with a maximum step size of 1. Use n_blk=20 and n_smpls=100. Compare your computed Ev with exact (analytical) result. What is the acceptance ratio? Adjust the maximum step size if necessary.

  2. Use a step size of 0.1. After equilibration, make n_smpls=1, i.e., each block consisting of only one sample. Use 10,000 blocks. What is your error estimate? Given this error bar, does the computed Ev agree with the exact result? Now use n_smpls=100 and n_blk=100. What happens? Explain. Study some other combinations of n_smpls and n_blk while keeping their product fixed at 10,000. Plot the error estimate as a function of n_smpls. Explain the curve.

  3. As we learned in class, the Gaussian distribution can be sampled directly by the Box-Muller method. If we had used this method rather than Metropolis to compute Ev, what would the plot in 2 look like? Given 10,000 samples from the Box-Muller method, how would the accuracy compare with that of the 100*100 combination in 2 (lager, smaller, same)? When in doubt, do the numerical experiment.