6810: 1094 Activities 14

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 that are easily generalized.

Autocorrelation

When evaluating an average over Monte Carlo configurations, we want to skip the first n1 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²e^-x2 from 0 to 10 and dividing by the integral of e^-x2 (so this is like < x² > with a probability distribution proportional to e^-x2).

1. Calculate the (normalized) integral in Mathematica for reference. Answer:
2. Consider first random sampling as a review/recap of our discussions.
   • 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 12e? Explain.
     
     
     
     
   • 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?
     
     
     
3. 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 Session 14 notes for answers.)
     
     
     
     
   • 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" (use a plot). Set "initial_skip" in the code to this value. (Note that your result will depend on your choice of "max_step".) Explain your value.
   
   
   
   
4. To figure out a good value for "skip", we'll calculate the "autocorrelation function". This is defined in the Session notes as equation (14.2). For us, "A" is "x²". Your task is to generate a figure like the one in the notes.
   • 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 (14.2) except for the rightmost average. 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 (14.2) for each l (you should only use [iterations - l] configurations to do this calculation; why?). Then you can calculate (14.2) and output it.]
5. Now that the Metropolis results are reliable, repeat step 2. above, but now for the Metropolis sampling. What do you find?
   
   

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!).

1. Take a look at the "variational_SHO.cpp" code and see how it implements variational Monte Carlo for a one-dimensional harmonic oscillator using the VariationalMC class. Compile, link, and run it (there is a makefile). What is the trial wave function? Why is this a reasonable one to use?
   
   
   
   
   
2. 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.
3. The variational_SHO.dat file can be plotted in gnuplot. The plotfile "variational_SHO.plt" is provided to plot it with error bars. Try it out a few times each with 100, 1000, 10000 steps (attach a plot with a sample result). Does the graph make sense for a variational calculation? What about the error bars? How might you modify the code to find the minimum automatically (rather than graphically)?
   
   
   
   
   
4. Does this code implement all the features explored in the "Autocorrelation" section? If not, how would you improve this code?
   
   
   
   
   
5. Bonus: The trial function is too good a guess. Try modifying it 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.
   
   
   
   

Checking for Bad Input Values

You may have noticed that the nonlinear pendulum codes with the menus went crazy if you accidently typed a letter or a symbol rather than an integer. Here we see how to avoid this.

1. Use make_input_check to create input_check and run it. First try some integers and then -1 to see how it is supposed to work. Then restart it and put in a letter or symbol (or whatever you want). To stop it, use Ctrl-c.
2. A fix is given in input_check_fixed (which has its own makefile). Try it out. Look at the printout or on the screen to see what was done. What was done to fix it?