Global Optimization: Quantum Thermal Annealing with Path Integral Monte Carlo

Quantum polyamorphism in compressed distinguishable helium-4

We demonstrate that two amorphous solid states can exist in ⁴He consisting of distinguishable Boltzmann atoms under compressed conditions. The isothermal compression of normal or supercritical fluid ⁴He was conducted at 3-25 K using the isobaric-isothermal path integral centroid molecular dynamics simulation. The compression of fluid first produced the low-dispersion amorphous (LDA) state possessing modest extension of atomic necklaces. Further isothermal compression up to the order of 10 kbar to 1 Mbar or an isobaric cooling of LDA induced the transition to the high-dispersion amorphous (HDA) state. The HDA was characterized by long quantum wavelengths of atoms extended over several Angstroms and the promotion of atomic residual diffusion. They were related to the quantum tunneling of atoms bestriding the potential saddle points in this glass. The change in pressure or temperature induced the LDA-HDA transition reversibly with hysteresis, while it resembled the coil-globule transition of classical polymers. The HDA had lower kinetic and higher Gibbs free energies than the LDA at close temperature. The HDA was absent at T ≥ 13 K, while the LDA-HDA transition pressure significantly decreased with lowering temperature. The LDA and HDA correspond to the trapped and tunneling regimes proposed by Markland et al. [J. Chem. Phys. 136, 074511 (2012)], respectively. The same reentrant behavior as they found was observed for the expansion factor of the quantum wavelength as well as for atomic diffusivity.

Electrostatics-Assisted Building-Up Procedure for Capturing Energy Minima of Metal Clusters: Test Case of Agn Clusters

Global geometry optimization of metal clusters is an important problem in nanophysics. The starting geometries of the clusters generated with empirical or other model potentials are generally optimized further by density functional theory (DFT)-based energy minimization. For this purpose, several algorithms such as simulated annealing, genetic algorithms, basin hopping, etc. are used. Our building-up procedure generates putative lower-energy structures of metal (M) clusters, M_n+1, M_n+2, etc., by anchoring one or more metal atoms in the vicinity of the minima of the molecular electrostatic potential (MESP) of M_n. Here, we report an application of this method to Ag_n clusters, for 5 ≤ n ≤ 20, followed up by DFT-based geometry optimization, generating several lower-energy structures than those reported in the literature. New low-energy isomers are obtained by applying the same procedure to the test case of mixed-metal clusters, Ni_nAg_m, for n + m = 4 and 5. In conclusion, our MESP-based building-up procedure offers a new general methodology for generating lower-energy geometries of metal clusters.

Development of a generalized hybrid Monte Carlo algorithm to generate the multicanonical ensemble with applications to molecular systems

In the present paper, a generalized hybrid Monte Carlo method to generate the multicanonical ensemble has been developed, which is a generalization of the multicanonical hybrid Monte Carlo (HMC) method by Hansmann and co-workers [Chem. Phys. Lett. 259, 321 (1996)]. The generalized hybrid Monte Carlo (GHMC) method is an equations-of-motion guided Monte Carlo combined with partial momentum refreshment. We successfully applied our multicanonical GHMC to dense Lennard-Jones fluids and a coarse grained protein model. It is found that good computational efficiency can be gained in the case of the acceptance ratio around 60% for the models examined. While a large number of molecular dynamics (MD) steps in a single GHMC cycle is needed to yield good computational efficiency at a large mixing ratio of momenta with thermal noise vectors, corresponding to the original multicanonical HMC method, a small number of MD steps are enough to achieve good efficiency at a small mixing ratio. This property is useful to develop a composite algorithm combining the present GHMC method with other Monte Carlo moves.

Robust quantum optimizer with full connectivity

Quantum phenomena have the potential to speed up the solution of hard optimization problems. For example, quantum annealing, based on the quantum tunneling effect, has recently been shown to scale exponentially better with system size than classical simulated annealing. However, current realizations of quantum annealers with superconducting qubits face two major challenges. First, the connectivity between the qubits is limited, excluding many optimization problems from a direct implementation. Second, decoherence degrades the success probability of the optimization. We address both of these shortcomings and propose an architecture in which the qubits are robustly encoded in continuous variable degrees of freedom. By leveraging the phenomenon of flux quantization, all-to-all connectivity with sufficient tunability to implement many relevant optimization problems is obtained without overhead. Furthermore, we demonstrate the robustness of this architecture by simulating the optimal solution of a small instance of the nondeterministic polynomial-time hard (NP-hard) and fully connected number partitioning problem in the presence of dissipation.

Steered quantum dynamics for energy minimization

We introduce a quantum optimal control algorithm for energy minimization that combines the diffeomorphic modulation under observable response preserving homotopy (D-MORPH) gradient and the Broyden Fletcher Goldfarb Shanno (BFGS) iterative scheme for nonlinear optimization. An extended set of controls defining the time-dependent mass, dipole moment, and external perturbational field are optimized to find an effective Hamiltonian that steers the dynamics of the system into the global minimum without getting trapped into local minima. The algorithm is illustrated as applied to energy minimization on rugged surfaces and golf potentials comparable to those previously explored for testing quantum annealing methodologies.

Replica exchange statistical temperature molecular dynamics algorithm

The replica exchange statistical temperature molecular dynamics (RESTMD) algorithm is presented, designed to alleviate an extensive increase of the number of replicas required as system size increases in the conventional temperature replica exchange method (tREM), and to obtain improved sampling in individual replicas. RESTMD optimally integrates multiple STMD (Phys. Rev. Lett. 2006, 97, 050601) runs with replica exchanges, giving rise to a flat energy sampling in each replica with a self-adjusting weight determination. The expanded flat energy dynamic sampling range allows the use of significantly fewer STMD replicas while maintaining the desired acceptance probability for replica exchanges. The computational advantages of RESTMD over conventional REM and single-replica STMD are explicitly demonstrated with an application to a coarse-grained protein model. The effect of two different kinetic temperature control schemes on the sampling efficiency is explored for diverse simulation conditions.

Ab initio random structure searching

It is essential to know the arrangement of the atoms in a material in order to compute and understand its properties. Searching for stable structures of materials using first-principles electronic structure methods, such as density-functional-theory (DFT), is a rapidly growing field. Here we describe our simple, elegant and powerful approach to searching for structures with DFT, which we call ab initio random structure searching (AIRSS). Applications to discovering the structures of solids, point defects, surfaces, and clusters are reviewed. New results for iron clusters on graphene, silicon clusters, polymeric nitrogen, hydrogen-rich lithium hydrides, and boron are presented.

Quantum-thermal annealing with a cluster-flip algorithm

A quantum-thermal annealing method using a cluster-flip algorithm is studied in the two-dimensional spin-glass model. The temperature (T) and the transverse field (Gamma) are decreased simultaneously with the same rate along a linear path on the T-Gamma plane. We found that the additional pulse of the transverse field to the frozen local spins produces a good approximate solution with a low computational cost.

Statistical temperature molecular dynamics: Application to coarse-grained β-barrel-forming protein models

Recently the authors proposed a novel sampling algorithm, "statistical temperature molecular dynamics" (STMD) [J. Kim et al., Phys. Rev. Lett. 97, 050601 (2006)], which combines ingredients of multicanonical molecular dynamics and Wang-Landau sampling. Exploiting the relation between the statistical temperature and the density of states, STMD generates a flat energy distribution and efficient sampling with a dynamic update of the statistical temperature, transforming an initial constant estimate to the true statistical temperature T(U), with U being the potential energy. Here, the performance of STMD is examined in the Lennard-Jones fluid with diverse simulation conditions, and in the coarse-grained, off-lattice BLN 46-mer and 69-mer protein models, exhibiting rugged potential energy landscapes with a high degree of frustration. STMD simulations combined with inherent structure (IS) analysis allow an accurate determination of protein thermodynamics down to very low temperatures, overcoming quasiergodicity, and illuminate the transitions occurring in folding in terms of the energy landscape. It is found that a thermodynamic signature of folding is significantly suppressed by accurate sampling, due to an incoherent contribution from low-lying non-native IS in multifunneled landscapes. It is also shown that preferred accessibility to such IS during the collapse transition is intimately related to misfolding or poor foldability.

Quantum annealing of an Ising spin-glass by Green's function Monte Carlo

We present an implementation of quantum annealing (QA) via lattice Green's function Monte Carlo (GFMC), focusing on its application to the Ising spin glass in transverse field. In particular, we study whether or not such a method is more effective than the path-integral Monte Carlo- (PIMC) based QA, as well as classical simulated annealing (CA), previously tested on the same optimization problem. We identify the issue of importance sampling, i.e., the necessity of possessing reasonably good (variational) trial wave functions, as the key point of the algorithm. We performed GFMC-QA runs using such a Boltzmann-type trial wave function, finding results for the residual energies that are qualitatively similar to those of CA (but at a much larger computational cost), and definitely worse than PIMC-QA. We conclude that, at present, without a serious effort in constructing reliable importance sampling variational wave functions for a quantum glass, GFMC-QA is not a true competitor of PIMC-QA.

Inherent Structure Analysis of Protein Folding

An analysis in terms of the inherent structures (IS, local minima) of the multidimensional potential energy landscape is applied to proteins. Detailed calculations are performed for the 46 bead BLN model, which folds into a four-stranded beta-barrel. Enhanced sampling has allowed determination of 239 199 IS states, believed to encompass nearly all the compact, low-energy states, and of well-averaged thermodynamic quantities at low temperature. The density of states shows distinct lobes for compact and extended states, and entropic barriers for the collapse and local ordering transitions. A two-dimensional scatterplot or density of states clearly shows the multifunnel structure of the energy landscape. The anharmonic vibrational free energy is found to play a crucial role in protein folding. The problem of determining the folding transition in a multifunnel system is discussed, and novel indicators of folding are introduced. A particularly clear picture is obtained through the occupation probabilities, pi, of individual low-lying IS, which become finite below the collapse temperature; it is suggested that poor foldability corresponds to a large "misfolding interval" where the excited state pi>0 exceeds that of the native state p0.

Statistical-temperature Monte Carlo and molecular dynamics algorithms

A simulation method is presented that achieves a flat energy distribution by updating the statistical temperature instead of the density of states in Wang-Landau sampling. A novel molecular dynamics algorithm (STMD) applicable to complex systems and a Monte Carlo algorithm are developed from this point of view. Accelerated convergence for large energy bins, essential for large systems, is demonstrated in tests on the Ising model, the Lennard-Jones fluid, and bead models of proteins. STMD shows a superior ability to find local minima in proteins and new global minima are found for the 55 bead AB model in two and three dimensions. Calculations of the occupation probabilities of individual protein inherent structures provide new insights into folding and misfolding.

Local-structural diversity and protein folding: Application to all-β off-lattice protein models

Global measures of structural diversity within a distribution of biopolymers, such as the radius of gyration and percent native contacts, have proven useful in the analysis of simulation data for protein folding. In this paper we describe a statistical-based methodology to quantify the local structural variability of a distribution of biopolymers, applied to 46- and 69-"residue" off-lattice, three-color model proteins. Each folds into beta-barrel structures. First we perform a principal component analysis of all interbead distance variables for a large number of independent, converged Boltzmann-distributed samples of conformations collected at each of a wide range of temperatures. Next, the principal component vectors are subjected to orthogonal (varimax) rotation. The results are displayed on so-called "squared-loading" plots. These provide a quantitative measure of the contribution to the sample variance of the position of each residue relative to the others. Dominant structural elements, those having the largest structural diversity within the sampled distribution, are responsible for peaks and shoulders observed in the specific heat versus temperature curves, generated using the weighted histogram analysis method. The loading plots indicate that the local-structural diversity of these systems changes gradually with temperature through the folding transition but radically changes near the collapse transition temperature. The analysis of the structural overlap order statistic suggests that the 46-mer thermodynamic folding transition involves the native state and at least three other nearly native intermediates. In the case of the 46-mer protein model, data are generated at sufficiently low temperatures that squared-loading plots, coupled with cluster analysis, provide a local and energetic description of its glassy state.

Optimization by quantum annealing: lessons from hard satisfiability problems

The path integral Monte Carlo simulated quantum annealing algorithm is applied to the optimization of a large hard instance of the random satisfiability problem (N = 10,000). The dynamical behavior of the quantum and the classical annealing are compared, showing important qualitative differences in the way of exploring the complex energy landscape of the combinatorial optimization problem. At variance with the results obtained for the Ising spin glass and for the traveling salesman problem, in the present case the linear-schedule quantum annealing performance is definitely worse than classical annealing. Nevertheless, a quantum cooling protocol based on field-cycling and able to outperform standard classical simulated annealing over short time scales is introduced.

Energy landscape distortions and the mechanical unfolding of proteins

Molecular simulations and an energy landscape analysis are used to examine the stretching of a model protein. A mapping of the energy landscape shows that stretching the protein causes energy minima and energy barriers to flatten out and disappear, and new energy minima to be created. The implications of these landscape distortions depend on the timescale regime under which the protein is stretched. When the timescale for thermally activated processes is longer than the timescale of stretching, the disappearances of energy barriers provide the mechanism for protein unfolding. When the timescale for thermally activated processes is shorter than the timescale of stretching, the landscape distortions influence the stretching process by changing the number and types of energy minima in which the system can exist.

Quantum annealing of the traveling-salesman problem

We propose a path-integral Monte Carlo quantum annealing scheme for the symmetric traveling-salesman problem, based on a highly constrained Ising-like representation, and we compare its performance against standard thermal simulated annealing. The Monte Carlo moves implemented are standard, and consist in restructuring a tour by exchanging two links (two-opt moves). The quantum annealing scheme, even with a drastically simple form of kinetic energy, appears definitely superior to the classical one, when tested on a 1002-city instance of the standard TSPLIB.

Long time molecular dynamics for enhanced conformational sampling in biomolecular systems

Although molecular dynamics methods are commonly used to drive biomolecular simulations, the technique provides insufficient sampling to impact studies of the 200-300 residue proteins of greatest interest. One severe limitation of molecular dynamics is that the integrators are restricted by resonance phenomena to small time steps (Delta t<8 fs) much slower then the time scales of important structural and solvent rearrangements. Here, a novel set of equations of motion and a reversible, resonance-free, integrator are designed which permit step sizes on the order of 100 fs to be used.

Ligand-protein docking using a quantum stochastic tunneling optimization method

A novel hybrid optimization method called quantum stochastic tunneling has been recently introduced. Here, we report its implementation within a new docking program called EasyDock and a validation with the CCDC/Astex data set of ligand-protein complexes using the PLP score to represent the ligand-protein potential energy surface and ScreenScore to score the ligand-protein binding energies. When taking the top energy-ranked ligand binding mode pose, we were able to predict the correct crystallographic ligand binding mode in up to 75% of the cases. By using this novel optimization method run times for typical docking simulations are significantly shortened.

Application of a chaperone-based refolding method to two- and three-dimensional off-lattice protein models

A model of protein-chaperone interaction as a two-phase (unfolding/refolding) iterative annealing mechanism able to promote structural segregation of hydrophobic and hydrophilic monomers and thereby facilitate access to nativelike states has recently been applied successfully to two 22-mers of the Honeycutt and Thirumalai BLN (hydrophobic, hydrophilic, neutral) heteropolymer model. This technique is here applied to a much wider data set: 94 8-mers of the off-lattice protein model originally presented in two dimensions by Stillinger and Head-Gordon, and later extended into three dimensions by Irbäck and Potthast; the model chaperone is shown to be equally successful, and by progressive elaboration of the chaperone model as in the earlier BLN model work, to be utilizing very similar underlying mechanisms. It is demonstrated that on average, contacts with the model chaperone give rise to a consistent movement in structure space in the direction of more nativelike structures; this method of global minimization does not therefore rely fundamentally on random search. Insofar as the responses to the chaperone of the two- and three-dimensional forms of the substrate model do differ, this can be interpreted as reflecting the different handling of hydrophilic monomers in the models-in particular, whether there is active repulsion between these and monomers of hydrophobic character. The chaperone-induced refolding method is also tested on a set of 220 9-mer chains of each version of the substrate model, where it is seen that the two-dimensional model, with its more clearly distinguished roles for the hydrophobic and hydrophilic monomers, shows a more favorable scaling behavior.

Theory of quantum annealing of an Ising spin glass

Probing the lowest energy configuration of a complex system by quantum annealing was recently found to be more effective than its classical, thermal counterpart. By comparing classical and quantum Monte Carlo annealing protocols on the two-dimensional random Ising model (a prototype spin glass), we confirm the superiority of quantum annealing relative to classical annealing. We also propose a theory of quantum annealing based on a cascade of Landau-Zener tunneling events. For both classical and quantum annealing, the residual energy after annealing is inversely proportional to a power of the logarithm of the annealing time, but the quantum case has a larger power that makes it faster.

Using novel variable transformations to enhance conformational sampling in molecular dynamics

One of the computational "grand challenges" is to develop methodology capable of sampling conformational equilibria in systems with rough energy landscapes. Here, a significant step forward is made by combining molecular dynamics with a novel variable transformation designed to enhance sampling by reducing barriers without introducing bias and, thus, to preserve, perfectly, equilibrium properties.

Global minimization of an off-lattice potential energy function using a chaperone-based refolding method

A global energy minimization method based on what is known about the mechanisms of the GroEL/GroES chaperonin system is applied to two 22-mers of an off-lattice protein model whose native states are beta-hairpins and which have structural similarity to short peptides known to interact strongly with the GroEL substrate binding domain. These model substrates have been used by other workers to test the effectiveness of a number of global minimization techniques, and are regarded as providing a significant challenge. The minimization method developed here is progressively elaborated from an initial simple form that targets exposed hydrophobic regions for unfolding to include a refolding phase that encourages the later recompactification of partly unfolded substrate; this refolding phase is seen to be crucial in the successful application of the method. The optimal handling of hydrophilic monomers within the model is also systematically explored, and it is seen that the best interpretation of their role is one that allows the chaperonin model to operate in "proofreading" mode whereby misfolded substrates are recognized by their surface exposure of a large proportion of hydrophobic monomers. The final version of the model allows native-like structures to be found quickly, on average for the two 22-mer substrates after 6 or 7 chaperone contacts. These results compare very favorably with those that have been obtained elsewhere using generic global minimization methods such as those based on thermal annealing. The paper concludes with a discussion of the place of the technique within the general category of hypersurface deformation methods for global minimization, and with suggestions as to how the chaperone-based method developed here could be elaborated so as to be effective on longer substrate chains that give rise to more complex tertiary structures in their native states.