Monte Carlo method in statistical physics
Encyclopedia
Monte Carlo in statistical physics refers to the application of the Monte Carlo method
to problems in statistical physics
, or statistical mechanics
.
, PS for simplicity, the mean value of A using the Boltzmann distribution:
.
where
is the energy of the system for a given state defined by
- a vector with all the degrees of freedom (for instance, for a mechanical system, 
),
and
is the partition function
.
One possible approach to solve this multivariable integral is to exactly enumerate all possible configurations
of the system, and calculate averages at will. This is actually done in
exactly solvable systems, and in simulations of simple systems with few
particles. In realistic systems, on the other hand, even an exact enumeration
can be difficult to implement.
For those systems, the Monte Carlo integration (and not to be confused with Monte Carlo method
,
which is used to simulate molecular chains) is generally employed. The main motivation for its use is the fact that, with the Monte Carlo integration, the error goes as
, independently of the dimension of the integral.
Another important concept related to the Monte Carlo integration is the importance sampling
, a technique that improves the computational time of the simulation.
On the following sections, the general implementation of the Monte Carlo integration for solving this kind of problems is discussed.
is
where
are uniformly obtained from all the phase space (PS) and N is the number of sampling points (or function evaluations).
From all the phase space, some zones of it are generally more important to the mean of the variable A than others. In particular, those that have the value of
sufficiently high when compared to the rest of the energy spectra are the most relevant for the integral. Using this fact, the natural question to ask is: is it possible to choose, with more frequency, the states that are known to be more relevant to the integral? The answer is yes, using the Importance sampling
technique.
Lets assume
is a distribution that chooses the states that are known to be more relevant to the integral.
The mean value of
can be rewritten as
,
where
are the sampled values taking into account the importance probability 
. This integral can be estimated by
where
are now randomly generated using the 
distribution. Since most of the times it is not easy to find a way of generating states with a given distribution, the Metropolis algorithm must be used.
, to choose for the importance sampling is the Boltzmann distribution or micro-canonic distribution. Let
be the distribution to use. Substituting on the previous sum,
.
So, the procedure to obtain a mean value of a given variable, using metropolis algorithm, with the micro-canonical distribution, is to use the Metropolis algorithm to generate states given by the distribution
and perform means over 
.
One important issue must be considered when using the metropolis algorithm with the micro-canonical distribution: when performing a given measure, i.e realization of
, one must ensure that that realization is not correlated with the previous state of the system (otherwise the states are not being "randomly" generated). On systems with relevant energy gaps, this is the major drawback of the use of the micro-canonical distribution because the time needed to the system de-correlate from the previous state can tend to infinity.
The multi-canonic approach uses a different choice for importance sampling:
where
is the density of states
of the system. The major advantage of this choice is that the energy histogram is flat, i.e. the generated states are equally distributed on energy. This means that, when using the Metropolis algorithm, the simulation doesn't see the "rough energy landscape", because every energy is treated equally.
The major drawback of this choice is the fact that, on most systems,
is unknown. To overcome this, the Landau algorithm
is normally used to obtain the DOS during the simulation. Note that after the DOS is known, the mean values of every variable can be calculated for every temperature, since the generation of states does not depend on
.
. Lets consider a 2 dimensional spin network, with L spins (lattice sites) on each side. There are naturally
spins, and so, the phase space is discrete and is characterized by N spins, 
where 
is the spin of each lattice site. the system's energy is given by 
, where 
are the group of first neighborhood spins of i and J is the interaction matrix (for a ferromagnetic ising model, J is the identity matrix). The problem is stated.
On this example, the objective is to obtain
and 
(for instance, to obtain the magnetic susceptibility
of the system) since it's straightforward to generalize to other observables. According to the definition,
.
be the system's Boltzmann temperature and initialize the system with an initial state (which can be anyone since the final result should not depend on it).
With micro-canonic choice, the metropolis method must be employed. Because there is no right way of choosing which state is to be picked, one can particularize and choose to try to flip one spin at the time. This choice is usually called single spin flip. The following steps are to be made to perform a single measurement.
step 1: generate a state that follows the
distribution:
step 1.1: Perform TT times the following iteration:
step 1.1.1: pick a lattice site at random (with probability 1/N), which will be called i, with spin
.
step 1.1.2: pick a random number
.
step 1.1.3: calculate the energy change of trying to flip the spin i:
and its magnetization change:
step 1.1.4: if
, flip the spin (
), otherwise, don't.
step 1.1.5: update the several macroscopic variables in case the spin flipped:
, 
after TT times, the system is considered to be not correlated from its previous state, which means that, at this moment, the probability of the system to be on a given state follows the Boltzmann distribution, which is the objective proposed by this method.
step 2 -> perform the measurement:
step 2.1: save, on an histogram, the values of M and M^2.
As a final note, one should note that TT is not easy to estimate because it is not easy to say when the system is de-correlated from the previous state. To surpass this point, one generally do not use a fixed TT, but TT as a tunneling time. One tunneling time is defined as the number of steps 1. the system needs to make to go from the minimum of its energy to the maximum of its energy and return back.
A major drawback of this method with the single spin flip choice in systems like Ising model is that the tunneling time scales as a power law as
where z is greater than 0.5, phenomenon known as critical slowing down.
, lack a dynamical description and are only defined by an energy prescription; for these the Monte Carlo approach is the only one feasible.
for optimization, in which a fictitious temperature is introduced and then gradually lowered.
Monte Carlo method
Monte Carlo methods are a class of computational algorithms that rely on repeated random sampling to compute their results. Monte Carlo methods are often used in computer simulations of physical and mathematical systems...
to problems in statistical physics
Statistical physics
Statistical physics is the branch of physics that uses methods of probability theory and statistics, and particularly the mathematical tools for dealing with large populations and approximations, in solving physical problems. It can describe a wide variety of fields with an inherently stochastic...
, or statistical mechanics
Statistical mechanics
Statistical mechanics or statistical thermodynamicsThe terms statistical mechanics and statistical thermodynamics are used interchangeably...
.
Overview
The general motivation to use the Monte Carlo method in statistical physics is to evaluate a multivariable integral. The typical problem begins with a system of which the Hamiltonian is known, it is at a given temperature and it follows the Boltzmann statistics. To obtain the mean value of some macroscopic variables, say A, the general approach is to compute, over all the phase spacePhase space
In mathematics and physics, a phase space, introduced by Willard Gibbs in 1901, is a space in which all possible states of a system are represented, with each possible state of the system corresponding to one unique point in the phase space...
, PS for simplicity, the mean value of A using the Boltzmann distribution:
.where
is the energy of the system for a given state defined by - a vector with all the degrees of freedom (for instance, for a mechanical system, ), andis the partition function
Partition function (statistical mechanics)
Partition functions describe the statistical properties of a system in thermodynamic equilibrium. It is a function of temperature and other parameters, such as the volume enclosing a gas...
.
One possible approach to solve this multivariable integral is to exactly enumerate all possible configurations
of the system, and calculate averages at will. This is actually done in
exactly solvable systems, and in simulations of simple systems with few
particles. In realistic systems, on the other hand, even an exact enumeration
can be difficult to implement.
For those systems, the Monte Carlo integration (and not to be confused with Monte Carlo method
Monte Carlo method
Monte Carlo methods are a class of computational algorithms that rely on repeated random sampling to compute their results. Monte Carlo methods are often used in computer simulations of physical and mathematical systems...
,
which is used to simulate molecular chains) is generally employed. The main motivation for its use is the fact that, with the Monte Carlo integration, the error goes as
, independently of the dimension of the integral.Another important concept related to the Monte Carlo integration is the importance sampling
Importance sampling
In statistics, importance sampling is a general technique for estimating properties of a particular distribution, while only having samples generated from a different distribution rather than the distribution of interest. It is related to Umbrella sampling in computational physics...
, a technique that improves the computational time of the simulation.
On the following sections, the general implementation of the Monte Carlo integration for solving this kind of problems is discussed.
Importance sampling
The estimative, under Monte Carlo integration, of an integral defined asis
where
are uniformly obtained from all the phase space (PS) and N is the number of sampling points (or function evaluations).From all the phase space, some zones of it are generally more important to the mean of the variable A than others. In particular, those that have the value of
sufficiently high when compared to the rest of the energy spectra are the most relevant for the integral. Using this fact, the natural question to ask is: is it possible to choose, with more frequency, the states that are known to be more relevant to the integral? The answer is yes, using the Importance samplingImportance sampling
In statistics, importance sampling is a general technique for estimating properties of a particular distribution, while only having samples generated from a different distribution rather than the distribution of interest. It is related to Umbrella sampling in computational physics...
technique.
Lets assume
is a distribution that chooses the states that are known to be more relevant to the integral.The mean value of
can be rewritten as,where
are the sampled values taking into account the importance probability . This integral can be estimated bywhere
are now randomly generated using the distribution. Since most of the times it is not easy to find a way of generating states with a given distribution, the Metropolis algorithm must be used.Canonical
Because it is known that the most likelihood states are those that maximize the Boltzmann distribution, a good distribution,, to choose for the importance sampling is the Boltzmann distribution or micro-canonic distribution. Letbe the distribution to use. Substituting on the previous sum,
.So, the procedure to obtain a mean value of a given variable, using metropolis algorithm, with the micro-canonical distribution, is to use the Metropolis algorithm to generate states given by the distribution
and perform means over .One important issue must be considered when using the metropolis algorithm with the micro-canonical distribution: when performing a given measure, i.e realization of
, one must ensure that that realization is not correlated with the previous state of the system (otherwise the states are not being "randomly" generated). On systems with relevant energy gaps, this is the major drawback of the use of the micro-canonical distribution because the time needed to the system de-correlate from the previous state can tend to infinity.Multi-canonical
As stated before, micro-canonical approach has a major drawback, which becomes relevant in most of the systems that use Monte Carlo Integration. For those systems with "rough energy landscapes", the multi-canonical approach can be used.The multi-canonic approach uses a different choice for importance sampling:
where
is the density of statesDensity of states
In solid-state and condensed matter physics, the density of states of a system describes the number of states per interval of energy at each energy level that are available to be occupied by electrons. Unlike isolated systems, like atoms or molecules in gas phase, the density distributions are not...
of the system. The major advantage of this choice is that the energy histogram is flat, i.e. the generated states are equally distributed on energy. This means that, when using the Metropolis algorithm, the simulation doesn't see the "rough energy landscape", because every energy is treated equally.
The major drawback of this choice is the fact that, on most systems,
is unknown. To overcome this, the Landau algorithmWang and Landau algorithm
The Wang and Landau algorithm proposed by Fugao Wang and David P. Landau is an extension of Metropolis Monte Carlo sampling. It is designed to calculate the density of states of a computer-simulated system, such as an Ising model of spin glasses, or model atoms in a molecular force field...
is normally used to obtain the DOS during the simulation. Note that after the DOS is known, the mean values of every variable can be calculated for every temperature, since the generation of states does not depend on
.Implementation
On this section, the implementation will focus on the Ising modelIsing model
The Ising model is a mathematical model of ferromagnetism in statistical mechanics. The model consists of discrete variables called spins that can be in one of two states . The spins are arranged in a graph , and each spin interacts with its nearest neighbors...
. Lets consider a 2 dimensional spin network, with L spins (lattice sites) on each side. There are naturally
spins, and so, the phase space is discrete and is characterized by N spins, where is the spin of each lattice site. the system's energy is given by , where are the group of first neighborhood spins of i and J is the interaction matrix (for a ferromagnetic ising model, J is the identity matrix). The problem is stated.On this example, the objective is to obtain
and (for instance, to obtain the magnetic susceptibilityMagnetic susceptibility
In electromagnetism, the magnetic susceptibility \chi_m is a dimensionless proportionality constant that indicates the degree of magnetization of a material in response to an applied magnetic field...
of the system) since it's straightforward to generalize to other observables. According to the definition,
.Canonical
First, the system must be initialized: let be the system's Boltzmann temperature and initialize the system with an initial state (which can be anyone since the final result should not depend on it).With micro-canonic choice, the metropolis method must be employed. Because there is no right way of choosing which state is to be picked, one can particularize and choose to try to flip one spin at the time. This choice is usually called single spin flip. The following steps are to be made to perform a single measurement.
step 1: generate a state that follows the
distribution:step 1.1: Perform TT times the following iteration:
step 1.1.1: pick a lattice site at random (with probability 1/N), which will be called i, with spin
.step 1.1.2: pick a random number
.step 1.1.3: calculate the energy change of trying to flip the spin i:
and its magnetization change:
step 1.1.4: if
, flip the spin ( ), otherwise, don't.step 1.1.5: update the several macroscopic variables in case the spin flipped:
, after TT times, the system is considered to be not correlated from its previous state, which means that, at this moment, the probability of the system to be on a given state follows the Boltzmann distribution, which is the objective proposed by this method.
step 2 -> perform the measurement:
step 2.1: save, on an histogram, the values of M and M^2.
As a final note, one should note that TT is not easy to estimate because it is not easy to say when the system is de-correlated from the previous state. To surpass this point, one generally do not use a fixed TT, but TT as a tunneling time. One tunneling time is defined as the number of steps 1. the system needs to make to go from the minimum of its energy to the maximum of its energy and return back.
A major drawback of this method with the single spin flip choice in systems like Ising model is that the tunneling time scales as a power law as
where z is greater than 0.5, phenomenon known as critical slowing down.Applicability
The method thus neglects dynamics, which can be a major drawback, or a great advantage. Indeed, the method can only be applied to static quantities, but the freedom to choose moves makes the method very flexible. An additional advantage is that some systems, such as the Ising modelIsing model
The Ising model is a mathematical model of ferromagnetism in statistical mechanics. The model consists of discrete variables called spins that can be in one of two states . The spins are arranged in a graph , and each spin interacts with its nearest neighbors...
, lack a dynamical description and are only defined by an energy prescription; for these the Monte Carlo approach is the only one feasible.
Generalisations
The great success of this method in statistical mechanics has led to various generalizations such as the method of simulated annealingSimulated annealing
Simulated annealing is a generic probabilistic metaheuristic for the global optimization problem of locating a good approximation to the global optimum of a given function in a large search space. It is often used when the search space is discrete...
for optimization, in which a fictitious temperature is introduced and then gradually lowered.
See also
- Monte Carlo integration
- Metropolis algorithm
- Importance samplingImportance samplingIn statistics, importance sampling is a general technique for estimating properties of a particular distribution, while only having samples generated from a different distribution rather than the distribution of interest. It is related to Umbrella sampling in computational physics...
- Quantum Monte CarloQuantum Monte CarloQuantum Monte Carlo is a large class of computer algorithms that simulate quantum systems with the idea of solving the quantum many-body problem. They use, in one way or another, the Monte Carlo method to handle the many-dimensional integrals that arise...
- Monte Carlo molecular modelingMonte Carlo molecular modelingMonte Carlo molecular modeling is the application of Monte Carlo methods to molecular problems. These problems can also be modeled by the molecular dynamics method. The difference is that this approach relies on statistical mechanics rather than molecular dynamics. Instead of trying to reproduce...