Comparing the efficiency of Metropolis Monte Carlo and molecular-dynamics methods for configuration space sampling

The investigations repor ted at this sympos ium are in the process of publication elsewhere [1, 2]; since the la t te r papers are r a the r detailed, we will content ourselves here with a ra ther brief survey of the work, with only a few addit ional results. The t e rm (~ Monte Carlo ,) has come into use to designate numerical methods in which specifically stochastic elements are introduced, in contras t to the whole body of classical numerical techniques which consist of numerical evaluations of complete ly de te rminan t algebraic expressions. The par t icular Monte Carlo me thod used here was devised by METROPOHS et al. [3]; its essential fea ture is t ha t i t produces a Markov chain [4] in which the individual 5 Ia rkov s tates are points in the usual configuration space of s tat is t ical mechanics, for a sys tem of N molecules confined a t a t e m p e r a t u r e T in ~ volume Y. The t ransi t ion probabil i t ies character izing the Markov chain are de termined in such a way tha t the va lue of any funct ion of configuration state, averaged over all s tates developed in the Markov cha.in, approaches the pe t i t canonical ensemble average value of the same function. I t should be noticed t ha t the present me thod does no t consist in an evaluat ion of the configurational phase in tegral b y means of a random, selection of points in configuration space. A me thod of the l a t t e r sort, which m a y be compared with a sequence of independent throws of a die, is indeed

A novel method was developed to enhance canonical sampling. A system is divided into virtually introduced sub-states, called "virtual states," which does not exist in reality. The configuration sampling is achieved by a standard canonical sampling method, the Metropolis Monte Carlo method, and confined in a virtual state for a while. In contrast, inter-virtual state motions are controlled by transition probabilities, which can be set arbitrarily. A simple recursive equation was introduced to determine the inter-virtual state transition probabilities, by which the sampling is enhanced considerably. We named this method "virtual-system coupled canonical Monte Carlo (VcMC) sampling." A simple method was proposed to reconstruct a canonical distribution function at a certain temperature from the resultant VcMC sampling data. Two systems, a one-dimensional double-well potential and a three-dimensional ligand-receptor binding/unbinding model, were examined. VcMC produced an accurate canonical distribution much more quickly than a conventional canonical Monte Carlo simulation does.

The problem of the global Lagrangian description of a dynamical system which has a global Hamiltonian description is considered; the motion of a charged particle in the field of a static magnetic monopole is a prototype of such dynamical problems. The symplectic structure on the (cotangent bundle of the) manifold is associated with a closed but not necessarily exact two-form. In conventional language this means that the tensor field inverse to the symplectic metric is divergence free but not necessarily expressible globally as the curl of a vector field; and hence the usual passage to Hamilton's principle and the Lagrangian can be carried out in patches but not necessarily globally. Generalizing techniques developed elsewhere, we show that the fundamental expression for the action integral is like the surface integral of the closed but not necessarily exact antisymmetric tensor field (two-form) rather than like the line integral of a vector field (one-form). This naturally leads to the path-space formalism in which the action functional depends not only on a point in configuration space but also a path from a chosen fixed point terminating in the relevant point in configuration space. We then show how this is a universal unfolding and how the topological obstructions in the configuration-space description are circumvented in the path-space description. The path-space description contains redundancies and they can be reduced; locally the new manifold is obtained by associating an angle of rotation [a U(1) fiber] with each configuration-space point, but globally there would be a nontrivial fiber-bundle structure. These questions naturally carry over to the quantum theory of such systems and it is shown how quantization conditions arise. The structure of the quantum-theory formalism for such systems is analyzed. The considerations are applied to the monopole problem where the enlarged space in which the Lagrangian description is possible is explicitly identified as the manifold of the group SU(2).

This review provides a comprehensive description of the self-guided Langevin dynamics (SGLD) and the self-guided molecular dynamics (SGMD) methods and their applications. Example systems are included to provide guidance on optimal application of these methods in simulation studies. SGMD/SGLD has enhanced ability to overcome energy barriers and accelerate rare events to affordable time scales. It has been demonstrated that with moderate parameters, SGLD can routinely cross energy barriers of 20 kT at a rate that molecular dynamics (MD) or Langevin dynamics (LD) crosses 10 kT barriers. The core of these methods is the use of local averages of forces and momenta in a direct manner that can preserve the canonical ensemble. The use of such local averages results in methods where low frequency motion "borrows" energy from high frequency degrees of freedom when a barrier is approached and then returns that excess energy after a barrier is crossed. This self-guiding effect also results in an accelerated diffusion to enhance conformational sampling efficiency. The resulting ensemble with SGLD deviates in a small way from the canonical ensemble, and that deviation can be corrected with either an on-the-fly or a post processing reweighting procedure that provides an excellent canonical ensemble for systems with a limited number of accelerated degrees of freedom. Since reweighting procedures are generally not size extensive, a newer method, SGLDfp, uses local averages of both momenta and forces to preserve the ensemble without reweighting. The SGLDfp approach is size extensive and can be used to accelerate low frequency motion in large systems, or in systems with explicit solvent where solvent diffusion is also to be enhanced. Since these methods are direct and straightforward, they can be used in conjunction with many other sampling methods or free energy methods by simply replacing the integration of degrees of freedom that are normally sampled by MD or LD.

The paper aims at a comparison of techniques based on the kinetic Monte Carlo (kMC) and the conventional Metropolis Monte Carlo (MC) methods as applied to the hard-sphere (HS) fluid and solid. In the case of the kMC, an alternative representation of the chemical potential is explored [E. A. Ustinov and D. D. Do, J. Colloid Interface Sci. 366, 216 (2012)], which does not require any external procedure like the Widom test particle insertion method. A direct evaluation of the chemical potential of the fluid and solid without thermodynamic integration is achieved by molecular simulation in an elongated box with an external potential imposed on the system in order to reduce the particle density in the vicinity of the box ends. The existence of rarefied zones allows one to determine the chemical potential of the crystalline phase and substantially increases its accuracy for the disordered dense phase in the central zone of the simulation box. This method is applicable to both the Metropolis MC and the kMC, but in the latter case, the chemical potential is determined with higher accuracy at the same conditions and the number of MC steps. Thermodynamic functions of the disordered fluid and crystalline face-centered cubic (FCC) phase for the hard-sphere system have been evaluated with the kinetic MC and the standard MC coupled with the Widom procedure over a wide range of density. The melting transition parameters have been determined by the point of intersection of the pressure-chemical potential curves for the disordered HS fluid and FCC crystal using the Gibbs-Duhem equation as a constraint. A detailed thermodynamic analysis of the hard-sphere fluid has provided a rigorous verification of the approach, which can be extended to more complex systems.

Robotic manipulators designed for home assistance and new surgical procedures often have significant uncertainty in their actuation due to compliance requirements, cost constraints, and size limits. We introduce a new integrated motion planning and control algorithm for robotic manipulators that makes safety a priority by explicitly considering the probability of unwanted collisions. We first present a fast method for estimating the probability of collision of a motion plan for a robotic manipulator under the assumptions of Gaussian motion and sensing uncertainty. Our approach quickly computes distances to obstacles in the workspace and appropriately transforms this information into the configuration space using a Newton method to estimate the most relevant collision points in configuration space. We then present a sampling-based motion planner based on executing multiple independent rapidly exploring random trees that returns a plan that, under reasonable assumptions, asymptotically converges to a plan that minimizes the estimated collision probability. We demonstrate the speed and safety of our plans in simulation for (1) a 3-D manipulator with 6 DOF, and (2) a concentric tube robot, a tentacle-like robot designed for surgical applications.

Definitions Symmetry operations leave a set of objects in indistinguishable configurations which are said to be equivalent . A set of symmetry operators always contains at least one element, the identity operator E . When operating with E the final configuration is not only indistinguishable from the initial one, it is identical to it. A proper rotation , or simply rotation , is effected by the operator R ( φ n ), which means “carry out a rotation of configuration space with respect to fixed axes through an angle φ about an axis along some unit vector n .” The range of φ is − π φ ≤ π . Configuration space is the three-dimensional (3-D) space ℛ 3 of real vectors in which physical objects such as atoms, molecules, and crystals may be represented. Points in configuration space are described with respect to a system of three space-fixed right-handed orthonormal axes x , y , z , which are collinear with OX, OY, OZ (Figure 2.1(a)). (A right-handed system of axes means that a right-handed screw advancing from the origin along OX would rotate OY into OZ, or advancing along OY would rotate OZ into OX, or along OZ would rotate OX into OY.) The convention in which the axes x , y , z remain fixed, while the whole of configuration space is rotated with respect to fixed axes, is called the active representation .

A theory is developed for diffusion in supercooled liquids. The viewpoint taken is that the dynamics of system, represented by a point in configuration space, consists of long sojourns in the vicinity of the local minima of the total potential, interrupted by relatively infrequent barrier crossings to adjacent minima. An equation due to Zwanzig, plus a new relation for the hopping rate among the local minima, are combined in a theory for the self-diffusion constant D. All the quantities entering the two equations are shown to be derivable from the configuration-averaged densities of vibrational states, 〈ρ(ω)〉, introduced by us previously; the contribution of unstable modes plays a crucial role. Solution of the coupled equations reveals the existence of a critical temperature below which D vanishes. The theory is in excellent agreement with molecular dynamics values for D, over a wide temperature range in which D displays activated behavior, in Lennard-Jones argon at a reduced density of unity.

Uniform acceptance force biased Monte Carlo (UFMC) simulations have previously been shown to be a powerful tool to simulate atomic scale processes, enabling one to follow the dynamical path during the simulation. In this contribution, we present a simple proof to demonstrate that this uniform acceptance still complies with the condition of detailed balance, on the condition that the characteristic parameter λ = 1/2 and that the maximum allowed step size is chosen to be sufficiently small. Furthermore, the relation to Metropolis Monte Carlo (MMC) is also established, and it is shown that UFMC reduces to MMC by choosing the characteristic parameter λ = 0 [Rao, M. et al. Mol. Phys.1979, 37, 1773]. Finally, a simple example compares the UFMC and MMC methods.

1,3,5-tri-amino-2,4,6-tri-nitrobenzene (TATB) is a typical and wide- studied explosive molecular and its single crystal is a typical triclinic lattice with a, b, c respectively 9.010, 9.028, 6.812 angstrom and α, β, γ 108.59, 91.82, 119.97. This explosive crystal is a stiff and anisotropic material. Packing of polymer on its crystal surface can obviously improve its mechanics properties. As a powerful tool, molecular dynamics (MD) simulation can be used to calculate the mechanics properties mainly elastic properties. But when applied the software package of Materials Studio (MS) to carry out MD, only the isotropic elastic properties which were averaged in x, y and z directions of materials can be obtained. To calculate the elastic properties in one direction of anisotropic materials, we developed a method which is similar to the experimental determination of elastic properties in one direction. Firstly, a P1 periodic super cell of TATB with 34.06×36.04 ×28.838 angstroms which ab plane was designated as (0 1 0) planar of the crystal was constructed. After several fixing and relaxing steps, the cell was pre-equilibrated 500 ps and performed 100 ps MD at 298 K within NVT ensemble. The averaged isotropic tensile modulus and Poisson' ratio were calculated by MS analysis module of elastic properties (static) and they were 1796±90 GPa and 0.231 respectively. Secondly, to obtain the elastic properties in each direction, NPT ensemble was chosen and different stresses in six different directions were added through many tries to keep the cell parameters fluctuating around those values in NVT ensemble while performing constant pressure MD. The MS averaged isotropic tensile modulus and Poisson' ratio of well equilibrated NPT's system were then 2361 ±30 GPa and 0.273 respectively. Finally, while carrying out subsequent constant pressure MD, different magnitude of compressive or tensile stress was applied to the cell in one direction. The elastic properties were then obtained via the strain-stress profile. The calculated anisotropic tensile modulus was separately 1423, 2558 and 1955 GPa in x, y and z direction and the averaged value was 1979 GPa while averaged Poisson's Ratio was 0.217. These two averaged values were well agreed with the MS calculated ones within about 20% deviation. The result revealed that our method can be applied to calculate elastic properties of anisotropic materials.

For high precision in source reconstruction of magnetoencephalography (MEG) or electroencephalography data, high accuracy of the coregistration of sources and sensors is mandatory. Usually, the source space is derived from magnetic resonance imaging (MRI). In most cases, however, no quality assessment is reported for sensor-to-MRI coregistrations. If any, typically root mean squares (RMS) of point residuals are provided. It has been shown, however, that RMS of residuals do not correlate with coregistration errors. We suggest using target registration error (TRE) as criterion for the quality of sensor-to-MRI coregistrations. TRE measures the effect of uncertainty in coregistrations at all points of interest. In total, 5544 data sets with sensor-to-head and 128 head-to-MRI coregistrations, from a single MEG laboratory, were analyzed. An adaptive Metropolis algorithm was used to estimate the optimal coregistration and to sample the coregistration parameters (rotation and translation). We found an average TRE between 1.3 and 2.3 mm at the head surface. Further, we observed a mean absolute difference in coregistration parameters between the Metropolis and iterative closest point algorithm of and m. A paired sample t-test indicated a significant improvement in goal function minimization by using the Metropolis algorithm. The sampled parameters allowed computation of TRE on the entire grid of the MRI volume. Hence, we recommend the Metropolis algorithm for head-to-MRI coregistrations.

It is a commonly observed phenomenon that spherical particles with inertia in an incompressible fluid do not behave as ideal tracers. Due to the inertia of the particle, the planar dynamics are described in a four-dimensional phase space and thus can differ considerably from the ideal tracer dynamics. Using finite-time Lyapunov exponents, we compute the sensitivity of the final position of a particle with respect to its initial velocity, relative to the fluid, and thus partition the relative velocity subspace at each point in configuration space. The computations are done at every point in the relative velocity subspace, thus giving a sensitivity field. The Stokes number, being a measure of the independence of the particle from the underlying fluid flow, acts as a parameter in determining the variation in these partitions. We demonstrate how this partition framework can be used to segregate particles by Stokes number in a fluid. The fluid model used for demonstration is a two-dimensional cellular flow.

By combining a well‐known theorem of Longuet‐Higgins on the phase changes of electronic wave functions along loops encircling a conical intersection of potential surfaces, and a recently proven “vertical symmetry theorem” on the relative distributions of point symmetry groups and potential surface critical points in configuration space, a new result is obtained, suitable for the prediction of the presence of a critical point within a domain of configurations.

Given an eigenstate of a time-independent Hamiltonian, an equation is derived which generates the corresponding energy by means of a pointwise application of only the kinetic operator. Satisfaction of the equation at all points in configuration space, for arbitrary values of a scale factor, constitutes a necessary eigenstate condition. Explicit knowledge of the potential operator, other than knowledge of its homogeneity, is not required. In addition, integration of the pointwise equation yields a generalization of the familiar virial theorem. 2 references.

It is shown that it is possible to define ‘true local entropy’ in velocity space, in an approximate version of the two-fluid formulation of quantum theory introduced by the present author in earlier papers. Using this definition, it is then shown that it is possible to define finite forms for total entropy at all points in configuration space. This important step is achieved by the introduction of a ‘responding’ velocity space. The use of a basis system which responds to occupation number density, makes possible a clear separation of the statistics and the dynamics of the underlying quantum process, and also makes possible the unambiguous use of certain divergent and oscillatory integrals.