Daan Frenkel

Configurational bias Monte Carlo: a new sampling scheme for flexible chains

We propose a novel approach that allows efficient numerical simulation of systems consisting of flexible chain molecules. The method is especially suitable for the numerical simulation of dense chain systems and monolayers. A new type of Monte Carlo move is introduced that makes it possible to carry out large scale conformational changes of the chain molecule in a single trial move. Our scheme is based on the selfavoiding random walk algorithm of Rosenbluth and Rosenbluth. As an illustration, we compare the results of a calculation of mean-square end to end lengths for single chains on a two-dimensional square lattice with corresponding data gained from other simulations.

New Monte Carlo method to compute the free energy of arbitrary solids. Application to the fcc and hcp phases of hard spheres

We present a new method to compute the absolute free energy of arbitrary solid phases by Monte Carlo simulation. The method is based on the construction of a reversible path from the solid phase under consideration to an Einstein crystal with the same crystallographic structure. As an application of the method we have recomputed the free energy of the fcc hard‐sphere solid at melting. Our results agree well with the single occupancy cell results of Hoover and Ree. The major source of error is the nature of the extrapolation procedure to the thermodynamic limit. We have also computed the free energy difference between hcp and fcc hard‐sphere solids at densities close to melting. We find that this free energy difference is not significantly different from zero: −0.001

Molecular Dynamics Simulations

Molecular Dynamics (MD) simulations are in many respects very similar to real experiments. In MD, first, sample is prepared, a model system consisting of N particles is selected, and then Newtons equations of motion are solved for the system until the properties of the system no longer change with time. To measure an observable quantity in a MD simulation, one must first of all be able to express this observable as a function of the positions and momenta of the particles in the system. The best introduction to MD simulations is to consider a simple program. To start the simulation, one should assign initial positions and velocities to all particles in the system. The particle positions should be chosen compatible with the structure that one is aiming to simulate. A good MD program requires a good algorithm to integrate Newtons equations of motion. Accuracy for large time steps is more important because the longer the time step that one can use, the fewer evaluations of the forces are needed per unit of simulation time. For most MD applications, Verlet-like algorithms are perfectly adequate. However, sometimes it is convenient to employ a higher-order algorithm.

Prediction of absolute crystal-nucleation rate in hard-sphere colloids

Crystal nucleation is a much-studied phenomenon, yet the rate at which it occurs remains difficult to predict. Small crystal nuclei form spontaneously in supersaturated solutions, but unless their size exceeds a critical value—the so-called critical nucleus—they will re-dissolve rather than grow. It is this rate-limiting step that has proved difficult to probe experimentally. The crystal nucleation rate depends on Pcrit, the (very small) probability that a critical nucleus forms spontaneously, and on a kinetic factor (κ) that measures the rate at which critical nuclei subsequently grow. Given the absence of a priori knowledge of either quantity, classical nucleation theory is commonly used to analyse crystal nucleation experiments, with the unconstrained parameters adjusted to fit the observations. This approach yields no ‘first principles’ prediction of absolute nucleation rates. Here we approach the problem from a different angle, simulating the nucleation process in a suspension of hard colloidal spheres, to obtain quantitative numerical predictions of the crystal nucleation rate. We find large discrepancies between the computed nucleation rates and those deduced from experiments: the best experimental estimates of Pcrit seem to be too large by several orders of magnitude.

Tracing the phase boundaries of hard spherocylinders

We have mapped out the complete phase diagram of hard spherocylinders as a function of the shape anisotropy L/D. Special computational techniques were required to locate phase transitions in the limit L/D→∞ and in the close-packing limit for L/D→0. The phase boundaries of five different phases were established: the isotropic fluid, the liquid crystalline smectic A and nematic phases, the orientationally ordered solids—in AAA and ABC stacking—and the plastic or rotator solid. The rotator phase is unstable for L/D⩾0.35 and the AAA crystal becomes unstable for lengths smaller than L/D≈7. The triple points isotropic-smectic-A-solid and isotropic-nematic-smectic-A are estimated to occur at L/D=3.1 and L/D=3.7, respectively. For the low L/D region, a modified version of the Gibbs–Duhem integration method was used to calculate the isotropic-solid coexistence curves. This method was also applied to the I-N transition for L/D>10. For large L/D the simulation results approach the predictions of the Onsager theory. ...

Numerical calculation of the rate of crystal nucleation in a Lennard‐Jones system at moderate undercooling

We report a computer‐simulation study of the rate of homogeneous crystal nucleation and the structure of crystal nuclei in a Lennard‐Jones system at moderate undercooling. The height of the nucleation barrier has been determined using umbrella sampling, whereas the barrier crossing rate is calculated using molecular dynamics simulation. The simulations clearly show that the barrier crossing is a diffusive process. Nevertheless, the kinetic prefactor in the nucleation rate is found to be some two orders of magnitude larger than predicted by classical nucleation theory. The height of the barrier is in good agreement with the theoretical prediction. Although the Lennard‐Jones system has a stable face‐centered cubic (fcc) phase below the melting line, the precritical nuclei are found to be mainly body‐centered cubic (bcc) ordered. As they grow to their critical size, they become more fcc ordered in the core. However, the critical and postcritical nuclei retain a high degree of bcc ordering in the interface. F...

Accelerating Monte Carlo Sampling

This chapter discusses various advanced Monte Carlo (MC) techniques. The method of parallel tempering provides good sampling of systems that have a free energy landscape with many local minima. It resembles the technique of simulating annealing and is related to several other schemes such as the extended-ensemble method, simulated tempering, and J-walking. The idea of parallel tempering is to include MC trial moves that attempt to “swap” systems that belong to different thermodynamic states. The basic idea behind the hybrid MC scheme is that one can use Molecular Dynamics (MD) to generate MC trial moves. For every trial move, the particle velocities are chosen at random from a Maxwell distribution. It is often advantageous to construct a trial move that consists of a sequence of MD steps. Yet, one cannot make the time step for a single hybrid MC move too long because then the acceptance would become very small. It is also interesting to use hybrid MC on models that have an expensive (many-body) potential energy function that may, to a first approximation, be modeled using a cheap (pair) potential. One of the differences of simulations of models with continuous interactions compared to those of models with hard-core potentials is the way in which MC moves are optimized.

Computer simulation study of gas-liquid nucleation in a Lennard-Jones system

We report a computer-simulation study of homogeneous gas–liquid nucleation in a Lennard-Jones system. Using umbrella sampling, we compute the free energy of a cluster as a function of its size. A thermodynamic integration scheme is employed to determine the height of the nucleation barrier as a function of supersaturation. Our simulations illustrate that the mechanical and the thermodynamical surfaces of tension and surface tension differ significantly. In particular, we show that the mechanical definition of the surface tension cannot be used to compute this barrier height. We find that the relations recently proposed by McGraw and Laaksonen [J. Chem. Phys. 106, 5284 (1997)] for the height of the barrier and for the size of the critical nucleus are obeyed.

Novel scheme to study structural and thermal properties of continuously deformable molecules.

The authors present a method for calculating the chemical potential of arbitrary chain molecules in a computer simulation. The method is based on a generalization of Siepmanns method for calculating the chemical potential of chain molecules with a finite number of conformations. Next, the authors show that it is also possible to extend the configurational-bias Monte Carlo scheme developed recently by Siepmann and Frenkel (1992) to continuously deformable molecules. The utility of their technique for computing the chemical potential of chain molecules is demonstrated by computing the chemical potential of a fully flexible chain consisting of 10-20 segments in a moderately dense atomic fluid. Under these conditions the conventional particle-insertion schemes fail completely. In addition, they show that their novel configurational-bias Monte Carlo scheme compares favourably with conventional Monte Carlo procedures for chain molecules.

Fluid-fluid coexistence in colloidal systems with short-ranged strongly directional attraction

We present a systematic numerical study of the phase behavior of square-well fluids with a “patchy” short-ranged attraction. In particular, we study the effect of the size and number of attractive patches on the fluid–fluid coexistence. The model that we use is a generalization of the hard sphere square well model. The systems that we study have a stronger tendency to form gels than the isotropic square-well system. For this reason, we had to combine Gibbs ensemble simulations of the fluid–fluid coexistence with a parallel tempering scheme. For moderate directionality, changes of the critical density and the width of coexistence curves are small. For strong directionality, however, we find clear deviations from the extended law of corresponding states: in contrast to isotropic attractions, the critical point is not characterized by a universal value of the reduced second virial coefficient. Furthermore, as the directionality increases, multiparticle bonding affects the critical temperature. We discuss imp...