Molecular dynamics finite-temperature

Force field calculations often truncate the non bonded potential energy of a molecular system at some finite distance. Truncation (nonbonded cutoff) saves computing resources. Also, periodic boxes and boundary conditions require it. However, this approximation is too crude for some calculations. For example, a molecular dynamic simulation with an abruptly truncated potential produces anomalous and nonphysical behavior. One symptom is that the solute (for example, a protein) cools and the solvent (water) heats rapidly. The temperatures of system components then slowly converge until the system appears to be in equilibrium, but it is not. 

Because the cohesive energy of the fullerene Cyo is —7.29 eV/atom and that of the graphite sheet is —7.44 eV/atom, the toroidal forms (except torus C192) are energetically stable (see Fig. 5). Finite temperature molecular-dynamics simulations show that all tori (except torus Cm2) are thermodynamically stable. 

The approach to the evaluation of vibrational spectra described above is based on classical simulations for which quantum corrections are possible. The incorporation of quantum effects directly in simulations of large molecular systems is one of the most challenging areas in theoretical chemistry today. The development of quantum simulation methods is particularly important in the area of molecular spectroscopy for which quantum effects can be important and where the goal is to use simulations to help understand the structural and dynamical origins of changes in spectral lineshapes with environmental variables such as the temperature. The direct evaluation of quantum time- correlation functions for anharmonic systems is extremely difficult. Our initial approach to the evaluation of finite temperature anharmonic effects on vibrational lineshapes is derived from the fact that the moments of the vibrational lineshape spectrum can be expressed as functions of expectation values of positional and momentum operators. These expectation values can be evaluated using extremely efficient quantum Monte-Carlo techniques. The main points are summarized below. 

The motion of particles of the film and substrate were calculated by standard molecular dynamics techniques. In the simulations discussed here, our purpose is to calculate equilibrium or metastable configurations of the system at zero Kelvin. For this purpose, we have applied random and dissipative forces to the particles. Finite random forces provide the thermal motion which allows the system to explore different configurations, and the dissipation serves to stabilize the system at a fixed temperature. The potential energy minima are populated by reducing the random forces to zero, thus permitting the dissipation to absorb the kinetic energy. 

The basic principles are described in many textbooks [24, 26]. They are thus only sketchily presented here. In a conventional classical molecular dynamics calculation, a system of particles is placed within a cell of fixed volume, most frequently cubic in size. A set of velocities is also assigned, usually drawn from a Maxwell-Boltzmann distribution appropriate to the temperature of interest and selected in a way so as to make the net linear momentum zero. The subsequent trajectories of the particles are then calculated using the Newton equations of motion. Employing the finite difference method, this set of differential equations is transformed into a set of algebraic equations, which are solved by computer. The particles are assumed to interact through some prescribed force law. The dispersion, dipole-dipole, and polarization forces are typically included whenever possible, they are taken from the literature. 

In some situations we have performed finite temperature molecular dynamics simulations [50, 51] using the aforementioned model systems. On a simplistic level, molecular dynamics can be viewed as the simulation of the finite temperature motion of a system at the atomic level. This contrasts with the conventional static quantum mechanical simulations which map out the potential energy surface at the zero temperature limit. Although static calculations are extremely important in quantifying the potential energy surface of a reaction, its application can be tedious. We have used ah initio molecular dynamics simulations at elevated temperatures (between 300 K and 800 K) to more efficiently explore the potential energy surface. 

Some of the major areas of activity in this field have been the application of the method to more complex materials, molecular dynamics, [28] and the treatment of excited states. [29] We will deal with some of the new materials in the next section. Two major goals of the molecular dynamics calculations are to determine crystal structures from first principles and to include finite temperature effects. By combining molecular dynamics techniques and ah initio pseudopotentials within the local density approximation, it becomes possible to consider complex, large, and disordered solids. 

Figure 1. Density profiles for a droplet confined to a finite volume as predicted by YBG theoiy (solid line) are compared with the results of molecular dynamics simulations (7) ( ). The reduced temperature, kT/e=0.71, is near the triple point. The total number of atoms in each of the tems is indicated.

The need to reliably describe liquid systems for practical purposes as condensed matter with high mobility at a given finite temperature initiated attempts, therefore, to make use of statistical mechanical procedures in combination with molecular models taking into account structure and reactivity of all species present in a liquid and a solution, respectively. The two approaches to such a description, namely Monte Carlo (MC) simulations and molecular dynamics (MD), are still the basis for all common theoretical methods to deal with liquid systems. While MC simulations can provide mainly structural and thermodynamical data, MD simulations give also access to time-dependent processes, such as reaction dynamics and vibrational spectra, thus supplying — connected with a higher computational effort — much more insight into the properties of liquids and solutions. 

The effects of the solvent and finite temperature (entropy) on the Wittig reaction have been studied by using DFT in combination with molecular dynamics and a continuum solvation model.62 The free energy profile has been found to have a significant entropic barrier to the addition step of the reaction where only a small barrier was present in the potential energy curve. 

As molecular packing calculations involve just simple lattice energy minimizations another set of tests have focused on the finite temperature effects. For this purpose, Sorescu et al. [112] have performed isothermal-isobaric Monte Carlo and molecular dynamics simulations in the temperature range 4.2-325 K, at ambient pressure. It was found that the calculated crystal structures at 300 K were in outstanding agreement with experiment within 2% for lattice dimensions and almost no rotational and translational disorder of the molecules in the unit cell. Moreover, the space group symmetry was maintained throughout the simulations. Finally, the calculated expansion coefficients were determined to be in reasonable accord with experiment. 

The difference in the preferred binding mode observed for the Pd- and Ni-based catalysts can be the crucial factor determining activity/inactivity of these two systems in polar copolymerization. However, the question arises about the stability of the alternative binding modes at finite temperature. If the minima were separated by relatively low barriers and fast interconversion between the two isomer complexes could occur, then this difference would be of minor importance. In order to check the stability of the two modes and get the insight into the mechanism of possible interconversions, a series of molecular dynamics simulations was performed. 

In order to advance the usefulness of theoretical calculations in broader applications, the computations should be carried out for a system that describes experimental conditions as accurately as possible. To do this, consider the following topics for a description of a system beyond using a single static structure at 0 K conformational averaging of static gas phase structures (this has been partially addressed already in Sect. 3), solvent effects on static structures, zero-point and finite temperature vibrational averaging, and molecular dynamics (MD) or Monte-Carlo (MC) sampling without and with solvation. 

As described briefly in Chapter 3, a promising new method of electronic structure calculation utilizing combined molecular-dynamics and density-functional theory has recently been developed by Car and Pari-nello (1985). This approach has recently been applied to cristobalite, yielding equilibrium lattice constants within 1% of experiment (Allan and Teter, 1987), as shown in Table 7.2. New oxygen nonlocal pseudopotentials were also an important part of this study. Such a method is a substantial advance upon density-functional pseudopotential band theory, since it can be efficiently applied both to amorphous systems and to systems at finite temperature. 

---