FUNCTION: dmc - computes ground state energy and wave function for
particle in a box of length 2.
CALLING SEQUENCE:
dmc(tau,n_equil_blks, n_meas_blks, n_sweeps,n_wlks)
PARAMETERS:
tau - (imaginary) time step size (Del_tau in class)
n_equil_blks - number of blocks used to equilibrate the random walk
n_meas_blks - number of blocks used to measure energy
Must be greater than 1; see aver_error in output.
n_sweeps - number of steps -- sweeps thru walkers -- in each block
n_wlks - number of walkers used to represent wave function
in each step
SYNOPSIS:
- Scheme of code is as follows:
1. Randomly initialize the walkers and set up histogram to store them.
2. Equilibrate -- let random walk converge to ground state.
3. Measure the ground-state energy and wavefunction.
4. Output energy, error estimate and wavefunction in form of histogram
in a file histdata.
5. There is a program plt which will plot histogram of walkers and
compare it to the exact wavefunction.
- Running time grows as (n_equil_blks + n_meas_blks) * n_sweeps * n_wlks;
Second example below is 12 times as long as first example.
EXAMPLES:
dmc(0.05, 1, 2, 5, 10);
dmc(0.05, 4, 5, 10, 20);
FUNCTION: gaussian - generates a Gaussian random number
with zero mean and unit variance.
CALLING SEQUENCE:
gaussian()
PARAMETERS:
None
SYNOPSIS:
- Based on the Box-Mueller algorithm discussed in class (Mon, 1/27).
- The scheme is as follows:
1. Generate two random numbers x and y between [-1..1].
2. Exclude those (x,y) such r = x^2 + y^2 >= 1.
3. From the "inverse factor" fac = sqrt(-2 ln(r))/r, we can produce
two gaussian distributed numbers: (x*fac, y*fac).
Note that the code avoids direct computations of 'cos' and 'sin',
which are expensive.
- This scheme is efficient; A VERY GOOD EXERCISE would be to
compare it with the Metropolis algorithm.
FUNCTION: initialize - initializes walkers and histogram.
CALLING SEQUENCE:
initialize(tau, n_equil_blks, n_meas_blks, n_sweeps, n_wlks)
PARAMETERS:
Same as used by main program dmc.
SYNOPSIS:
- Dmc calls this function automatically.
- Scheme of code is as follows:
1. Generates initial walker distribution (e.g., uniform inside box).
2. Initailizes value of histogram (wave function) in each bin
to zero.
3. Creates functions - midpt - for computing midpoint position
of each histogram bin - and bin - for computing which bin a
given position x is in.
4. Also evaluates tauinv = 1/tau and rrtau = sqrt(tau) that are
used elsewhere.
FUNCTION: output - puts successively numbered runs into file "histdata"
The principal data list is called dataN where N is run number.
The list contains histogram data for sampled wave function
as well as exact wave function for plotting.
CALLING SEQUENCE:
output()
PARAMETERS:
None actually, but global varables are passed and created.
the most important of which are:
histarray - where output puts the histrogram data prior to output.
histdata - name of file where output stores results for plotting
dataN - data structue with wave function histogram suitable for
plotting. The number N is assigned automatically and
is printed on screen for your information.
SYNOPSIS:
- Prints on screen at conclusion of run: energy, error,
run number (N), execution time
- Rest of data is put in file histdata. Some data is printed as a
comment such as the details of run parameters. These run parameters
are also printed in the first two lines of the file dataN but are not
presently used by programs.
- Most of data structure dataN contains the histogram: bins, computed
wavefunction and exact wavefunction for use by plotting routine.
FUNCTION: step - does one step of random walk (sweep thru all walkers)
CALLING SEQUENCE:
step()
PARAMETERS:
None, but there are global variables shared with other routines:
rr() - generates uniform random numbers between 0 and 1.
gaussian() - generates gaussian random numbers with mean of zero
and variance of unity.
wlks - array storing the walkers.
tauinv - 1/tau, used in computing image correction of the
free Green's function, due to boundary condition
rttau - sqrt(tau), used in generating free Green's function.
SYNOPSIS:
- Used inside the main program [dmc] to generate new set of walkers
for each sweep. But idea of code is reusable.
- Scheme of code is as follows: [ going right to left in Eq. (3) ]
(a). Pick at random a walker whose value is a position -1=<dnow=<1
(b). Generate potential new walker from diffusion with free Green's
function g_0, namely, dnxt := dnow + rttau*gaussian();
(c). Now deal with K(dnxt,dnow):
1). Test if dnxt is consistent with being inside length-two box.
Free Green's function g_0 extends to infinity.
2). P in eq. (4) is exp(-2(1-abs(dnow))(1-abs(dnxt))/tau).
If this is less than a random number, we accept dnxt.
Otherwise return to step (a).
BUG CHECK.
- If program does not work, set printlevel to at least 5. And then type
step(); If output not informative, increase printlevel to 10.