Workshop on Modeling Rare Events in Complex Physical Systems - IMS

Workshop on Modeling Rare Events in Complex Physical Systems
(5 - 8 Nov 2013)

~ Abstracts ~

Finite temperature analysis of stochastic networks
Maria Cameron, University of Maryland, USA

The network of potential minima and transition states for the Lennard-Jones-38 (LJ38) cluster created by Wales's group exemplifies a large stochastic network with detailed balance and temperature-dependent pairwise transition rates. The transition process between the two deepest minima of LJ38 is analyzed at the temperature range from zero up to the melting point (T=0.18 e/kT). For T< 0.06 e/kT, the transition process is well-described by a single pathway obtained by means of the zero-temperature asymptotic analysis. However, as the temperature increases, this description gets rapidly inadequate because the process diversifies dramatically. Using the Transition Path Theory as a starting point, we propose an approach to effectively describe and analyze such a network at a finite temperature. It allows us to define and visualize transition tubes in the discrete case and explain the temperature-dependence of the transition rate.

Multiscale-based ideas for path sampling
Hiroshi Fujisaki, Nippon Medical School, Japan

We propose a novel path sampling method based on the Onsager-Machlup action by generalizing the multiscale enhanced sampling technique suggested by Moritsugu and coworkers (J. Chem. Phys. 133, 224105 (2010)). The basic idea of this method is that the system we want to study (for example, some molecular system described by molecular mechanics) is coupled to a coarse-grained (CG) system, which can move more quickly and computed more efficiently than the original system. We simulate this combined system (original + CG system) using Langevin dynamics where different heat baths are coupled to the two systems. When the coupling is strong enough, the original system is guided by the CG system, and able to sample the configuration and path space more efficiency. We need to correct the bias caused by the coupling, however, by employing the Hamiltonian replica exchange where we prepare many path replica with different coupling strengths. As a result, an unbiased path ensemble for the original system can be found in the weakest coupling path ensemble.

We apply this method to a model polymer with Asakura-Oosawa interaction, and compare the results with the conventional transition path sampling method. We also discuss some future applications of this method for large molecular systems.

Fast algorithms for transition state calculations
Weiguo Gao, Fudan University, China

Transition states are fundamental to understanding the reaction dynamics qualitatively. To date various methods of first principle location of transition states have been developed. In the absence of the knowledge of the final structure, the minimal-mode following method climbs up to a transition state without calculating the Hessian matrix. In this talk, we introduce a locally optimal search direction finding algorithm and an iterative minimizing method for the translation. Numerical experiments demonstrate the efficiency of our proposed algorithms. This is joint work with Jing Leng, Zhi-Pan Liu, Cheng Shang and Xiang Zhou.

Sampling saddle points on the free energy surface
Daniel Aitor Holdbrook, A*Star, Bioinformatics Institute

Efficient exploration of the free energy surface of a complex system at finite temperatures, though often desirable, remains a difficult task. Such explorations involve searching for metastable states and bottlenecks for transitions between these metastable states in the space of coarse grained variables. I will show how a set of dynamical equations obtained from a simple modification of Langevin dynamics can be used to determine transition states and escape paths from basins of attraction. The mean forces on the coarse grained variables and the Hessian information are obtained on-the-fly by evolving the underlying microscopic system. The utility of the method is illustrated by studying the saddle points associated with isomerization transition in alanine dipeptide.

Multi-step multi-target umbrella sampling and its application in trans-membrane permeation of ions
Dan Hu, Shanghai Jiao Tong University, China

Multi-target umbrella sampling is an important technique in free energy calculation of phase transitions. In this sampling technique, it is crucial to find a few proper reaction coordinates in obtaining a continuous transition trajectory and an accurate free energy profile. However, when the system if complex, it is usually very hard to find the reaction coordinates. In this talk, we discuss a multi-step multi-target umbrella sampling framework (MMUS), in which local reaction coordinates are designed step by step to perform multi-target umbrella samplings. The difficulty in obtaining proper local reaction coordinates is reduced by making use of the information of hysteresis in samplings. In order to obtain the free energy profile by the sampling results with different reaction coordinates, we also designed a new method based on weighted least squares to take the place of Wham.

We have successfully applied the new framework in the study of trans-membrane permeation process of ions. Proper reaction coordinates are designed and continuous transition trajectories are obtained in a relatively easier way. Our work makes clear the trans-membrane permeation mechanism of ions for the first time, in which a water chain is formed across the membrane to bridge the permeation of ions. With the presence of the water chain, the free energy barrier is significantly reduced for thick membranes.

Investigation conformational changes of biological macromolecules using kinetic network models
Xuhui Huang, The Hong Kong University of Science and Technology, Hong Kong

Simulating biologically relevant timescales at atomic resolution is a challenging task since typical atomistic simulations are at least two orders of magnitude shorter. Markov State Models (MSMs), a kinetic network model, built from molecular dynamics (MD) simulations provide one means of overcoming this gap without sacrificing atomic resolution by extracting long time dynamics from short MD simulations. In this talk, I will demonstrate the power of MSMs by applying it to simulate the complex conformational changes, that occurs at tens of microsecond timescales for a large RNA transcription complex (close to half million atoms). In the second part of my talk, I will introduce a new efficient dynamic clustering algorithm for the automatic construction of MSMs for multi-body systems. We have successfully applied this new algorithm to model the protein-ligand recognition and hydrophobic collapse processes that occur at a mixture of different timescales.

The application of the string method to the dislocation dynamics
Congming Jin, Zhejiang Sci-tech University, China

We apply the string method to the dislocation dynamics. First, we present a numerical method to compute the transition rates of thermally activated events in dislocation dynamics on the atomistic scale, based on the formula from the transition state theory. The method is applied to the migration of kinks in 30^o partial dislocations in silicon. To our knowledge, this is the first time that the contribution of entropy to the transition rate of such events has been calculated using reliable atomistic models. Then we systematically study the cross-slip process of the screw dislocation in aluminum using the string method, focusing on the dependence of the cross-slip mechanism on the dislocation length and the relative position of the dislocation at the primary and the cross-slip plane. We find that when the dislocation segment is short, the energetically most favorable cross-slip process follows the Fleicher model uniformly. On the other hand, for a longer dislocation segment, we have identified two cross-slip mechanisms with very similar energy barriers. In both mechanisms, the cross-slip of the entire segment is achieved by the motion of two constriction points following the Friedel-Escaig model. However, when the constriction starts to form, the cross-slip proceeds via either the symmetric transfer of the dislocation from the primary to the cross-slip plane without further recombination, or via the Fleicher model as in the short segment with the formation of an asymmetric constriction pair. The cross-slip energy barrier for the short and longer dislocation segment is 0.028 eV/b and 0.60eV, respectively.

Constructing the energy landscape of the genetic switching driven by intrinsic noise
Tiejun Li, Peking University, China

Genetic switching driven by noise is a fundamental cellular process in genetic regulatory networks. Quan- titatively characterizing this switching and its fluctuation properties is a key problem in computational biology. With the autoregulatory dimer model as a specific example, we design a general methodology to quantitatively understand the metastability in gene expressions perturbed by the intrinsic noise. Based on the large deviation theory, we develop new analytical techniques to characterize and compute the op- timal transition paths between the on and off states. We also construct the global quasi-potential energy landscape for the dimer model. From the obtained quasi-potential, we can extract quantitative results such as the stationary distributions of mRNA, protein and dimer, the noise strength of expression state, and the mean switching time starting from either stable state. In the final stage, we apply this procedure to a transcriptional cascade model. Our results suggest that the quasi-potential energy landscape and the proposed methodology are general to understand the metastability in other biological systems with intrinsic noise.

Analysis and simulation of multiscale intracellular bio-chemical reacting networks
Richard Di Liu, Michigan State University, USA

Intracellular reacting networks involving gene regulation often exhibits multiscale properties. That includes multiple reacting rates, multiple population magnitudes and multi-stability. Direct Stochastic Simulation Algorithm (SSA) would turn out to be inefficient dealing with such systems. Schemes such as Nested SSA and Tau-leaping method have proved to be effective for certain asymptotic regimes. I will discuss recent progress on applying the Transition Path Theory to study metastable reacting networks.

Stochastic surface walking method and beyond
Zhipan Liu, Fudan University, China

We propose an unbiased general-purpose potential energy surface (PES) searching method for both the structure and the pathway prediction of complex system from cluster to crystal, namely, stochastic surface walking method (SSW) [1,2]. The method is based on the idea of bias-potential-driven dynamics and Metropolis Monte Carlo, which is a natural extension of automated transition state location approach developed in the group via bias-potential driven constrained-Broyden dimer method [3,4,5]. A central feature of the SSW method is able to perturb smoothly a structural configuration toward a new configuration and simultaneously has the ability to surmount the high barrier in the path. We apply the method for locating the global minimum (GM) of model potential clusters (Lennard-Jones, Short-ranged Morse), carbon clusters up to 100 atoms from random structures. In addition to GM searching, the method can identify the pathways for chemical reactions with large dimensionality by connecting SSW trajectories with a newly-developed double-ended method within the same theoretical framework to locate the transition state (TS). We compare the efficiency of this method with nudged elastic band method, and show that our method is more efficient in TS location for a wide range of reactions in the gas phase and on the surfaces.

[1]Zhang, X.-J.; Shang, C. and Liu, Z-P.*, "From Atoms to Fullerene: Stochastic Surface Walking Solution for Automated Structure Prediction of Complex Material", J. Chem. Theory Comput, 2013, 9, 3252
[2]Shang, C. and Liu, Z.-P.*;"Stochastic Surface Walking Method for Structure Prediction and Pathway Searching", J. Chem. Theory Comput, 2013, 9, 183
[3]Shang, C. and Liu, Z.-P.*;"Constrained Broyden Dimer Method with Bias Potential for Exploring Potential Energy Surface of Multistep Reaction Process", J. Chem. Theory Comput, 2012, 8, 2215
[4]Shang, C.; Liu, Z.-P.*, "Constrained-Broyden-Minimization Combined with the Dimer Method for Locating Transition State of Complex Reactions" ,J. Chem. Theory Comput, 2010, 6, 1136
[5]Wang, H.-F.; Liu, Z.-P.* "Comprehensive Mechanism and Structure-Sensitivity of Ethanol Oxidation on Platinum: New Transition-State Searching Method for Resolving Complex Reaction Network", J. Am. Chem. Soc. 2008, 130,10996

Theory of transition path process
Jianfeng Lu, Duke University, USA

Understanding rare events like transitions of chemical system from reactant to product states is a challenging problem due to the time scale separation. In this talk, we will discuss some recent progress in mathematical theory of transition paths. In particular, we identify and characterize the stochastic process corresponds to transition paths. The study of transition path process helps to understand the transition mechanism and provides a framework to analyze numerical approaches for rare event sampling and simulation.

Simulating rare events in biomolecules via temperature accelerated molecular dynamics
Luca Maragliano, Italian Institute of Technology, Italy

I will present a method designed to explore the free energy landscape of a system associated with a large set of collective variables. An extended system is introduced where the collective variables and the physical ones are evolved concurrently in effective adiabatic separation. This allows to explore directly the free energy space of the collective variables and to accelerate sampling via an artificially high temperature acting on them. I will discuss the technique and its implementation, and illustrate its potential via applications to ligand diffusion and conformational changes in proteins and RNA molecules.

Droplet motion driven by liquid-vapor transition at three-phase contact line
Tiezheng Qian, The Hong Kong University of Science and Technology, Hong Kong

Recently, the dynamic van der Waals theory (DvdWT) has been presented for the study of hydrodynamics in one-component fluids with liquid-vapor transition in inhomogeneous temperature fields [Onuki A 2005 Phys. Rev. Lett. 94 054501]. We first derive the hydrodynamic boundary conditions at the fluid-solid interface for the DvdWT using conservation laws and the positive definiteness of entropy production together with the Onsager reciprocal relation. We then apply the DvdWT to the study of droplet motion driven by thermal gradients on solid substrates. The effect of thermal singularity at the liquid-vapor-solid three phase contact line is investigated. The droplet motion predicted by the continuum hydrodynamic model is also observed and semi-quantitatively verified by performing molecular dynamics simulations for confined one-component two-phase fluids. The effect of hydrodynamic fluctuations on the droplet motion may be further explored using the so-called Model H in the Hohenberg-Halperin notation.

Computing rare event properties by discretization and optimal control
Marco Sarich, Freie Universität Berlin, Germany

Markov processes are widely used to model physical, chemical, or biological systems, and these processes often exhibit metastability. This leads to the presence of rare events. That is, the computational effort for direct sampling is typically too high in order to achieve good estimates for intersting properties of these events. We will discuss computational approaches to time-reversible problems that are based on discretizaion and optimal control. Further, we will discuss possible extensions to the non-reversible case using a directed network example.

Modeling of reactive events
Eric Vanden-Eijnden, New York University, USA

Dynamics in nature often proceed in the form of reactive events. The system under study spends very long periods of time at various metastable states; only very rarely it hops from one metastable state to another. Understanding the dynamics of such systems requires us to study the ensemble of transition paths between the different metastable states. Transition path theory (TPT) is a general mathematical framework developed for this purpose. It is also the foundation for developing modern numerical algorithms such as the string method for finding the transition pathways or milestoning to calculate the reaction rate, and it can also be used in the context of Markov State Models (MSMs). In these lectures, we will review the basic ingredients of the transition path theory and discuss connections with transition state theory (TST) as well as approaches to metastability based on spectral and potential theory. We also discuss how the string method arises in order to find approximate solutions in the framework of the transition path theory, the connections between milestoning and TPT, and the way the theory help building MSMs.

Exploring energy landscapes: from molecules to nanodevices
David Wales, University of Cambridge, UK

Coarse-graining the potential energy surface into the basins of attraction of local minima provides a computational framework for investigating structure, dynamics and thermodynamics in molecular science. Steps between local minima form the basis for global optimisation via basin-hopping and for calculating thermodynamic properties using the superposition approach and basin-sampling. To treat global dynamics we must include transition states of the potential energy surface, which link local minima via steepest-descent paths. We may then apply the discrete path sampling method, which provides access to rate constants for rare events. In large systems the paths between minima with unrelated structures may involve hundreds of stationary points of the potential energy surface. New algorithms have been developed for both geometry optimisation and finding connections between distant local minima, which allow us to treat such systems. A graph transformation approach enables rate constants and committor probabilities to be extracted from kinetic transition networks containing over a million states. Applications will be presented for a range of different examples, including atomic and molecular clusters, biomolecules, condensed matter, and coarse-grained models of mesoscopic structures.

Selected Publications:
D.J. Wales, Curr. Op. Struct. Biol., 20, 3-10 (2010)
D.J. Wales, J. Chem. Phys., 130, 204111 (2009)
B. Strodel and D.J. Wales, Chem. Phys. Lett., 466, 105-115 (2008)
D.J. Wales and T.V. Bogdan, J. Phys. Chem. B, 110, 20765-20776 (2006)
D.J. Wales, Int. Rev. Phys. Chem., 25, 237-282 (2006)
D.J. Wales, "Energy Landscapes", Cambridge University Press, Cambridge, 2003

Hybrid parallel minimum action method and its applications
Xiaoliang Wan, Louisiana State University, USA

In this work, we report a hybrid (MPI/OpenMP) parallelization strategy for the minimum action method. The crucial part of the minimum action method is to minimize the Freidlin-Wentzell action functional. Due to the fact that the corresponding Euler-Lagrange equation is, in general, highly nonlinear and of high order, we solve the optimization problem directly instead of discretizing the Euler-Lagrange equation to provide a general but equivalent numerical framework. To enhance the efficiency of the minimum action method for general dynamical systems we consider parallel computing. In particular, we present a hybrid parallelization strategy based on MPI and OpenMP. We apply the hybrid parallelization strategy to construct a minimum action method for small perturbations of Navier-Stokes equations. We will discuss our numerical experiments to study the nonlinear instability of two dimensional Poiseuille flows from the large deviation point of view.

Rare events in dislocation dynamics using generalized Peierls-Nabarro model
Yang Xiang, The Hong Kong University of Science and Technology, Hong Kong

Dislocations are one-dimensional topological defects in crystalline solids. The motion of dislocations is directly responsible for the plastic deformation of these materials. The Peierls-Nabarro model is a hybrid framework that incorporates atomic features into the continuum theory. Such models are able to provide quantitative understanding of dislocation core-related properties at the continuum level, and a reliable alternative to the expensive, empirical atomistic simulations. We present a generalized Peierls-Nabarro model for curved dislocations with applications to dislocation kink migration and other rare events in dislocation dynamics.

A large deviation framework to analyze metastable behavior in climate systems
Xu Yang, University of California, Santa Barbara, USA

We studied the dynamic transition/excursion phenomena in climate systems. We built a framework using large deviation theory, in which different climate regimes are represented by the statistical most likely states and the transition is described by the most likelihood pathways connecting either metastable states or target sets in the small noise limit. Specifically we considered the energy-constrained stochastic dynamics (equilibrium statistical system), the most likely states of whose invariant measure coincide with the selective decay states. We compute the transition pathways using a constrained String method. Nonequilibrium statistical climate systems were also analyzed where the transition pathways were computed by the geometric minimum action method.

Repository based adaptive umbrella sampling
Yingkai Zhang, New York University, USA

A widely employed strategy to determine free energy surfaces is the introduction of an external biasing potential into the Hamiltonian which forces the system to visit the high barrier regions. An ideal choice of the biasing potential would be the negative of the free energy surfaces so that uniform sampling along chosen reaction coordinates can be achieved. Unfortunately, such info is exactly what we try to obtain from simulations and is not known in the first place. To overcome this inherent challenge, we have developed a new adaptive sampling approach, called as "repository based adaptive umbrella sampling" (RBAUS). Its main idea is that a sampling repository is continuously updated based on the latest simulation data, and the accumulated knowledge and sampling history is then employed to determine whether and how to update the biasing umbrella potential for subsequent simulations. In comparison with other adaptive methods, a unique and attractive feature of the RBAUS approach is that the frequency for updating the biasing potential depends on the sampling history and is adaptively simulated on-the-fly, which makes it possible to smoothly bridge non-equilibrium and quasi-equilibrium simulations. Such an adaptive updating is achieved by employing the following general principle: the biasing potential needs to be updated if the subsequent simulations still explores the previously oversampled region more than the previously undervisited; on the other hand, the biasing potential update is not needed if the subsequent simulations is achieving more uniform sampling. Thus we can see that a key element of the the RBAUS approach is the determination of sampling uniformity. Here we have demonstrated sampling entropy (SE) as an excellent indicator for uniform sampling as well as for the convergence of free energy simulations. By introducing sampling entropy and concentration theorem into the biasing potential updating scheme, the adaptivity, robustness and applicability of the RBAUS approach has been significantly improved, including the determination of two dimensional free energy surfaces.

Computation of saddle point and its application on materials and biology
Lei Zhang, Peking University, China

Saddle point search on an energy surface has attracted much attention in various areas. In this talk, I will first introduce efficient numerical methods for finding saddle points. Then I will combine the developed methods with the phase field model to show two applications: one is the nucleation during phase transition in FeCr alloys and the other is neuroblast delamination in Drosophila.

An iterative minimization scheme for saddle search
Xiang Zhou, City University of Hong Kong, Hong Kong

The gentlest ascent dynamics (GAD, E and ZHOU, 2011 Nonlinearity) transforms saddles of energy potential into a stable fixed point. Inspired by GAD, in this talk, I introduce a new formulation of iteratively minimizing a sequence of modified potential to find the saddles of the original function. We show that the iteration converges quadratically. An 175-atom example is illustrated as an application . This is the joint work with Weiguo Gao and Jing Leng.

Best viewed with IE 7 and above