This method uses a standard Monte Carlo method for a system coupled to a suitably chosen finite bath (H ¨uller, 1993). We consider an N-particle system
01:16:24
with energy E coupled to a finite reservoir which is an ideal gas with M degrees of freedom and kinetic energy k. One then studies the micro-canonical ensemble of the total, coupled system with fixed total energy G. An analysis of detailed balance shows that the ratio of the transition probabilities between two states is then
Wb →a
Wa→b = (G − Ea)N−2N (G − Eb)N−2N ≈ e−ζ (Eb−Ea) (4.93) where ζ = (N − 2)2Nkb and where kb = (G − Eb)N is the mean kinetic energy per particle in the bath. The only difference in the Monte Carlo method is that the effective inverse temperature ζ is adjusted dynamically during the course of the simulation. Data are then obtained by computing the mean value of the energy on the spin system E and the mean value of the temperature from kb. This method becomes accurate in the limit of large system size.
Plots of E vs. T then trace out the complete ‘van der Waals loop’ at a first order phase transition.
4.5.3 Q2R
The Q2R cellular automaton has been proposed as an alternative, microcanon-ical method for studying the Ising model. In a cellular automaton model the state of each spin in the system at each time step is determined completely by consideration of its nearest neighbors at the previous time step. The Q2R rule states that a spin is flipped if, and only if, half of its nearest neighbors are up and half down. Thus, the local (and global) energy change is zero. A starting spin configuration of a given energy must first be chosen and then the Q2R rule applied to all spins; this method is thus also deterministic after the ini-tial state is chosen. Thermodynamic properties are generally well reproduced, although the susceptibility below Tcis too low. (Other cellular automata will be discussed in Chapter 8.)
Problem 4.11 Simulate an L = 10, q = 10 Potts model square lattice using a microcanonical method and estimate the transition temperature. How does your answer compare with that obtained in Problem 4.8?
4 . 6 G E N E R A L R E M A R K S , C H O I C E O F E N S E M B L E
We have already indicated how models may be studied in different ensembles by different methods. There are sometimes advantages in using one ensem-ble over the other. In some cases there may be computational advantages to choosing a particular ensemble, in other situations there may be a symmetry which can be exploited in one ensemble as opposed to the other. One of the simplest cases is the study of a phase diagram of a system with a tricritical point.
Here there are both first order and second order transitions. As shown in Fig.
4.25 the phase boundaries look quite different when shown in the canonical and grand-canonical ensembles. Thus, for low ‘density’ (or magnetization in
01:16:24
Fig. 4.25 Phase diagram for an Ising antiferromagnet with nearest and next-nearest neighbor couplings with a tricritical point: (top) canonical ensemble, the shaded area is a region of two-phase coexistence; (bottom) grand canonical ensemble.
magnetic language) two phase transitions are encountered as the temperature is increased whereas if the ‘field’ is kept fixed as the temperature is swept only a single transition is found. Of course, to trace out the energy–field relation in the region where it is double valued, it is preferable to use a microcanonical ensemble (as was described in the previous section) or even other ensem-bles, e.g. a Gaussian ensemble (Challa and Hetherington, 1988). The use of a microcanonical ensemble for the study of protein folding will be given in a later chapter.
A situation in which it is preferable to use an ensemble which differs from those available to the experimentalist is the case of fluid or solid binary (A, B) mixtures. In the laboratory, for a given volume V and temperature T, the parti-cle numbers NAand NBwill be fixed (i.e. the relevant ensemble is the canonical ensemble). In a simulation it is often preferable to work in the ‘semi-grand canonical ensemble’ in which only the total number of particles N = NA+ NB
is fixed and an additional intensive variable, the chemical potential difference μ enters the Boltzmann factor in the transition probability. Due to the difference in chemical potential μ, ‘identity switches’, A → B or B → A may occur as attempted Monte Carlo moves. This is not ‘alchemy’ (like medieval chemists trying to transform lead into gold) but a valid method in statistical mechanics that is preferred because of its faster equilibration (in particular for solid alloys where ordered superstructures like in β-brass or copper–gold alloys occur).
An example for the use of this ensemble will be given in Chapter 6 of this book.
01:16:24
Fig. 4.26 Dynamic Monte Carlo algorithms for SAWs on a simple cubic lattice: (a) generalized Verdier–Stockmayer algorithm; (b) slithering snake algorithm; (c) pivot algorithm.
4 . 7 S TAT I C S A N D DY N A M I C S O F P O LY M E R M O D E L S O N L AT T I C E S
4.7.1 Background
Real polymers are quite complex and their simulation is a daunting task (Binder, 1995). There are a number of physically realistic approximations which can be made, however, and these enable us to construct far simpler models which (hopefully) have fundamentally the same behavior. First we recognize that the bond lengths of polymers tend to be rather fixed as do bond angles. Thus, as a more computationally friendly model we may construct a ‘polymer’ which is made up of bonds which connect nearest neighbor sites (monomers) on a lattice and which obey an excluded volume constraint. The sites and bonds on the lattice do not represent individual atoms and molecular bonds but are rather the building blocks for a coarse-grained model. Even within this simplified view of the physical situation simulations can become quite complicated since the chains may wind up in very entangled states in which further movement is almost impossible.