780.20: 1094 Session 13
Handouts: Excerpts on autocorrelation functions
and variational Monte Carlo, printouts
of codes
In this session, we'll continue our survey of Monte Carlo computational
methods with a look at autocorrelation
and variational Monte Carlo applied to simple problems.
The discussion is easily generalized to the more complex problems
for which these approaches are well suited.
Session 12 (cont.)
Continue with "The TwoD Ising Model", but move on to
Session 13 "Autocorrelation" after an hour. If you have time later
(possibly next week),
go back and finish it.
Autocorrelation
When evaluating an average over Monte Carlo configurations,
we want to skip the first n_{1} steps and then use data taken
every n0 Monte Carlo steps. How do we determine how large
to take n0 and n1?
We'll use a simple integration problem to explore these issues,
by integrating x^{2}e^{x2}
from 0 to 10 and dividing by the integral of e^{x2}
(so this is like < x^{2} > with
a probability distribution proportional to e^{x2}).
 Calculate the (normalized) integral in Mathematica for
reference.
 Consider first random sampling as a review/recap of our discussions.
 Create an autocorrelation_test project
and then run the autocorrelation_test.cpp
code a number of times (at least 10)
and see how the random sampling estimates
vary. Is the calculation of errors in the code (e.g.
error_random) consistent
with the discussion in the background notes from Session 11e?

Are the estimated errors consistent with the actual deviations
of the estimates from the exact result?
(E.g., if the error is one standard deviation, about what fraction
of the results should lie between the exact result +/ the error?
What fraction do you observe?)
 Vary the number N of configurations used to make the
estimates to see how the error scales.
(The number of configurations
is just the number of x's sampled in this case.)
With what power of N does it scale?
Is this result consistent with expectations from the notes?
 Next consider the Metropolis algorithm. There are three
parameters for you to adjust (besides the total number of iterations):
"max_step", "initial_skip" and "skip".
 What do each of these affect?
(Why do we skip steps at the beginning? Why do we skip
configurations? How do we know how much to change x?
See the Pang excerpt for answers, but discuss it first.)
 Adjust "max_step" so the acceptance rate is around 50% (or
less). [If stepping a vector X, we test every component to form a
Monte Carlo step and look for choosing a step size so that the
acceptance rate is around 50%.] What "max_step" did you
use?
 Devise a way to decide on the "initial_skip" that includes a
plot. Set "initial_skip" in the code to the value you deduce.
(Note that your result will depend on your choice of "max_step".)
Explain your value.
 To figure out a good value for "skip", we'll calculate the
"autocorrelation function". This is defined in the Pang handout
as equation (9.21). For us, "A" is "x^{2}". Your task
is to generate the analog of Fig. 9.1 for the current problem.
 What is the value of
the autocorrelation function C(l) at l=0?
 What is the condition that makes C(l) tend to zero (based
on the defining equation)?
 You have everything you need in the code already for (9.21)
except for the average in (9.22). Modify the code to
calculate that and then print out and plot C(l). Attach the plot.
[Hint: An easy way (although not the most efficient way!)
is to introduce an array "result_save[i]" in the Monte Carlo
step loop and save all of the "result" evaluations.
Then after this loop you can step through values of l (say from
0 to 100) and calculate (9.22) for each l (you should only use
iterations  l configurations to do this calculation; why?).
Then you can calculate (9.21) and output it.]
 Now that the Metropolis results are reliable, repeat step 2.
above, but now for the Metropolis sampling.
Variational Monte Carlo
We'll do a simple example of variational Monte Carlo to illustrate
the basic idea. Generalizing to more complex systems is straightforward
(but takes a lot longer to run!).
 Take a look at the "variational_SHO.cpp" code and how it
implements variational Monte Carlo for a onedimensional harmonic
oscillator using the VariationalMC class.
Make a variational_SHO project and add variational_SHO.cpp, random_seed.cpp,
VaritionalMC.h, and VariationalMC.cpp.
Run it. You'll be asked to supply
a range and step size for the variational parameter "a". This will
require some experimentation to make sure the minimum with respect to
"a" is in the interval you select.
 The variational_SHO.dat file is suitable for plotting in gnuplot.
The plotfile "variational_SHO.plt" is provided to illustrate
how to plot it with error bars.
Try it out (and attach a plot).
Does the graph make sense? How might you modify the code to
find the minimum automatically (rather than
graphically)?
 Does this code implement the features explored
in the "Autocorrelation" section? If not, how would you improve
this code?
 Try modifying the trial function to another reasonable form
with one variational parameter "a" (use your
imagination!). You'll have to modify the "Psi" and "E_local"
functions. Find the minimum graphically and compare to the
exact solution.
780.20: 1094 Session 13.
Last modified: 12:56 pm, May 05, 2008.
furnstahl.1@osu.edu