Monte Carlo simulation

As we have seen in a previous chapter, the motion of a system of $N$ particles with fixed $N,E,V$ can be simulated integrating the equations of motion performing a MD simulation and time averaging the quantities of interest. How can we do an ensemble average? One way would enumerating all the microstates, but this approach is not practical since the number of configurations is generally too large. In the spirit of Monte Carlo, we wish to develop a practical method to obtain representative sample of the total number of microstates. An obvious procedure is to fix $N$ and $V$, change the positions and velocities of individual particles at random, and retain the configuration of it has the desired total energy. However, this procedure is very inefficient, since in general most configuration will not have the desired total energy and must be discarded.

Imagine a macroscopic system which is divided into two ``subsystems'': the original system of interest, referred to as the system, and a subsystem which has one element. For historical reasons this extra degree of freedom is called a ``demon''. The demon travels about the system transferring energy as it attempts to change the dynamical variables of the system. If the demon has sufficient energy in its sack, it gives energy to any element of the system who requires energy to make the desired change. Conversely if the desired change lowers the energy of the system, the excess energy is given to the demon. The only limitation on the demon is that it cannot have negative energy. The algorithm is summarized in the following:

1. choose a particle at random and make a trial change of its coordinates;

2. compute the change in the energy of the system due to the change in coordinates;

3. if the trial change decreases the energy of the system, the system gives its energy to the demon, and the new configuration is accepted;

4. if the trial change increases the energy of the system, the new configuration is accepted if the demon has sufficient energy to give to the system. Otherwise, the new configuration is rejected and the particle retains its old coordinates;

5. if the trial change does not change the energy of the system, the new coordinate is accepted.

The above procedure is repeated until a representative sample of states is obtained. After a period of time for equilibration, the demon and the system will reach a compromise and agree on an average energy for each. The total energy remains constant, and since the demon is only one degree of freedom in comparison to the many degrees of freedom of the system, we expect that the energy fluctuations of the system will be small.

How do we know that his Monte Carlo simulation of the microcanonical ensemble will yield results equivalent to the time-average results of MD simulations? The assumption that the two averages yield equivalent results is called the ``ergodic hypothesis'' (or more accurately the quasi-ergodic hypothesis) Although these two averages cannot be shown to be identical in general, they have been found to yield equivalent results in all cases of practical interest.