Define
(1.1)
The expectation value of the Hamiltonian is
(1.2)
Define
(1.3)
Then change variables in the r integration to
(1.4)
So that
(1.5)
..\ImportanceSampling.doc .htm
Making
(1.6)
Define
(1.7)
So that
(1.8)
Each configuration involves N choices of yi, μi and ji so that
(1.9)
A Monte Carlo estimate of <H> is
(1.10)
The standard deviation in (1.10) for an eigenfunction is zero (VarianceFromZeroMCError.doc .htm), but is nonero for the numerator and denominator separately. (Random integral error analysis.doc .htm) For an eigenfunction, these both become the same. Thus the desired guiding function is such that
(1.11)
Let
(1.12)
For m = 0
(1.13)
For m > 0
(1.14)
This is faster when explicit. Code for finding and testing G(r) and g(r) as given by (1.12) is in Gofr.zip.
The hard part is solving (1.4) for r which can be rewritten as
(1.15)
This will normally be tabulated as a series of values. The tabulated value will give r0. Then Newton’s method (1.16) ..\..\solving\Newton.doc .htm should yield r(y) in just a few iterations.
(1.16)
Equation (1.16) is coded into ConGAtom.for in the subroutine XNEFIN about line 46.
In principle, this code gives
(1.17)
This in principle is then ready to place into the denominator of (1.10) to yield zero Monte Carlo error in the evaluation of the normalization integral for the correct set of C’s and Y. In practice, there is a problem with a 1/rij term that can slip into the wave function and that is always present in the Hamiltonian.
In systems with two singular centers v(r1)
and v(r2), both of which go as 1/r, w can be averaged over the
possible ways in which a point could have been selected. This is rigorous when done correctly because
w is the exact correction for the probability of selecting the point in the way
chosen. Thus if there are two ways of
selecting the same point, the probability of selecting it in either way is
exactly (w1+w2)/2. In the two center
case the point (x,y,z)i can be chosen either with respect to R1
with probability wi1
given by
(1.18)
Or it could be selected with respect to R2 with
(1.19)
The average probability of finding the point (x,y,z)i is
(1.20)
Note that when the Monte-Carlo picks the same f, but with different w’s, The average is <1/w> accounting of course for the relative frequency of selecting the different (x,y,z)’s. In this case where I definitely did not wait for an unlikely, small w, event, but noted analytically that it should be there, the average is of <w>. Note that when the argument of v(r1-R1) is small that there is a 1/(r1-R1)2 in the in the denominator of w1 which becomes the numerator of the sum term and multiplies the singularity to zero and that the same also happens for r1-R2.
In the case of an atom, the last electron is selected with probability pele with respect to one of the other electrons.
This means that with probability pele that any other electron distance could be chosen according to FMonte\eegfun.doc .htm with
(1.21)
This is averaged over all electron pair distances to give
(1.22)
Define
(1.23)
The final wi for the r defined in (1.1) is
(1.24)
The 1/w(rNe) in the second term in (1.24) corrects for the fact that the Ne’th electron in is selected with respect to the (Ne-1)th electron rather than with respect to the nucleus. This can be averaged over the fact that any of the Ne-1 electrons previously selected could be considered to be the second to last selected
(1.24)
In addition the order can be re-done so that any of the Ne electrons can be selected last
(1.24)
Note that the 1/w(rj) prevents the usual i < j simplification of the i,j sum.
The code for an atom is in ConGAtom.for. This uses random from ..\..\definitions\Random\random.for and locate and bli from ..\..\interpolation\locate.for and ..\..\interpolation\bli.for. The code assumes that unit 8 has been opened and in generation mode writes X(3,Nelec),w to this unit. In config mode it reads from this unit. The code is tested in FMonte\The rij tests.doc .htm. (The last line in the last table is the results from ConGenAtom). There are links there to the test code and to a zip of the complete set of Fortran codes used in the test. Internally the code allocates a vector X(3,Nelec) and W(Nelec) in two different routines. The code ConGAt77.for uses a parameter Maxele=100 to permanently over-allocate these. This code is also included in the zip.