Search

In this paper, we report a new algorithm for creating an adaptive basis set in the Multiconfigurational Ehrenfest (MCE) method, which is termed Full Cloning (FC), and test it together with the existing Multiple Cloning (MC) using the spin-boson model at zero-temperature as a benchmark. The zero-temperature spin-boson regime is a common hurdle in the development of methods that seek to model quantum dynamics. Two versions of MCE exist. We demonstrate that MC is vital for the convergence of MCE version 2 (MCEv2). The first version (MCEv1) converges much better than MCEv2, but FC improves its convergence in a few cases where it is hard to converge it with the help of a reasonably small size of the basis set. [ABSTRACT FROM AUTHOR]

We obtain the phase diagram of the hard core lattice gas with third nearest neighbor exclusion on the triangular lattice using Monte Carlo simulations that are based on a rejection-free flat histogram algorithm. In a recent paper [Darjani et al., J. Chem. Phys. 151, 104702 (2019)], it was claimed that the lattice gas with third nearest neighbor exclusion undergoes two phase transitions with increasing density with the phase at intermediate densities exhibiting hexatic order with continuously varying exponents. Although a hexatic phase is expected when the exclusion range is large, it has not been seen earlier in hard core lattice gases with short range exclusion. In this paper, by numerically determining the entropies for all densities, we show that there is only a single phase transition in the system between a low-density fluid phase and a high density ordered sublattice phase and that a hexatic phase is absent. The transition is shown to be first order in nature, and the critical parameters are determined accurately. [ABSTRACT FROM AUTHOR]

Simulating the unitary dynamics of a quantum system is a fundamental problem of quantum mechanics, in which quantum computers are believed to have significant advantage over their classical counterparts. One prominent such instance is the simulation of electronic dynamics, which plays an essential role in chemical reactions, non-equilibrium dynamics, and material design. These systems are time-dependent, which requires that the corresponding simulation algorithm can be successfully concatenated with itself over different time intervals to reproduce the overall coherent quantum dynamics of the system. In this paper, we quantify such simulation algorithms by the property of being fully-coherent: the algorithm succeeds with arbitrarily high success probability 1 − δ while only requiring a single copy of the initial state. We subsequently develop fully-coherent simulation algorithms based on quantum signal processing (QSP), including a novel algorithm that circumvents the use of amplitude amplification while also achieving a query complexity additive in time t, ln(1/δ), and ln(1/ϵ) for error tolerance ϵ: Θ ‖ H ‖ | t | + ln (1 / ϵ) + ln (1 / δ) . Furthermore, we numerically analyze these algorithms by applying them to the simulation of the spin dynamics of the Heisenberg model and the correlated electronic dynamics of an H2 molecule. Since any electronic Hamiltonian can be mapped to a spin Hamiltonian, our algorithm can efficiently simulate time-dependent ab initio electronic dynamics in the circuit model of quantum computation. Accordingly, it is also our hope that the present work serves as a bridge between QSP-based quantum algorithms and chemical dynamics, stimulating a cross-fertilization between these exciting fields. [ABSTRACT FROM AUTHOR]

Within the context of fewest-switch surface hopping (FSSH) dynamics, one often wishes to remove the angular component of the derivative coupling between states J and K . In a previous set of papers, Shu et al. [J. Phys. Chem. Lett. 11, 1135–1140 (2020)] posited one approach for such a removal based on direct projection, while we isolated a second approach by constructing and differentiating a rotationally invariant basis. Unfortunately, neither approach was able to demonstrate a one-electron operator O ̂ whose matrix element J O ̂ K was the angular component of the derivative coupling. Here, we show that a one-electron operator can, in fact, be constructed efficiently in a semi-local fashion. The present results yield physical insight into designing new surface hopping algorithms and are of immediate use for FSSH calculations. [ABSTRACT FROM AUTHOR]

In a recent paper, we proposed the adaptive shift method for correcting undersampling bias of the initiator-full configuration interaction (FCI) quantum Monte Carlo. The method allows faster convergence with the number of walkers to the FCI limit than the normal initiator method, particularly for large systems. However, in its application to some systems, mostly strongly correlated molecules, the method is prone to overshooting the FCI energy at intermediate walker numbers, with convergence to the FCI limit from below. In this paper, we present a solution to the overshooting problem in such systems, as well as further accelerating convergence to the FCI energy. This is achieved by offsetting the reference energy to a value typically below the Hartree–Fock energy but above the exact energy. This offsetting procedure does not change the exactness property of the algorithm, namely, convergence to the exact FCI solution in the large-walker limit, but at its optimal value, it greatly accelerates convergence. There is no overhead cost associated with this offsetting procedure and is therefore a pure and substantial computational gain. We illustrate the behavior of this offset adaptive shift method by applying it to the N2 molecule, the ozone molecule at three different geometries (an equilibrium open minimum, a hypothetical ring minimum, and a transition state) in three basis sets (cc-pVXZ, X = D, T, Q), and the chromium dimer in the cc-pVDZ basis set, correlating 28 electrons in 76 orbitals. We show that in most cases, the offset adaptive shift method converges much faster than both the normal initiator method and the original adaptive shift method. [ABSTRACT FROM AUTHOR]

In this paper we present the Markovian Multiagent Monte-Carlo Second Order Self-Consistent Field Algorithm (M3-SOSCF). This algorithm provides a highly reliable methodology for converging SCF calculations in single-reference methods using a modified differential evolution approach. Additionally, M3 is embarrassingly parallel and modular in regards to Newton–Raphson subroutines. We show that M3 is able to surpass contemporary SOSCFs in reliability, which is illustrated by a benchmark employing poor initial guesses and a second benchmark with SCF calculations which face difficulties using standard SCF algorithms. Furthermore, we analyse inherent properties of M3 and show that in addition to its robustness and efficiency, it is more user-friendly than current SOSCFs. [ABSTRACT FROM AUTHOR]

Attempts to predict molecular structure often try to minimize some potential function over a set of structures. Much effort has gone into creating potential functions and into creating algorithms for minimizing these potential functions. This paper develops a formalism that addresses a complementary question: What are the accuracy requirements for a potential function that predicts molecular structure? The formalism is applied to a simple model of a protein structure potential. The results of this calculation suggest that high accuracy predictions (∼1 Å rms deviation in α-carbon positions) of protein structures require monomer–monomer interaction energies accurate to within 5% to 15%. The paper closes with a discussion of the implications of these results for practical structure prediction. [ABSTRACT FROM AUTHOR]

We report the development of a new Laplace MP2 (second-order Møller–Plesset) implementation using a range separated Coulomb potential, partitioned into short- and long-range parts. The implementation heavily relies on the use of sparse matrix algebra, density fitting techniques for the short-range Coulomb interactions, while a Fourier transformation in spherical coordinates is used for the long-range part of the potential. Localized molecular orbitals are employed for the occupied space, whereas orbital specific virtual orbitals associated with localized molecular orbitals are obtained from the exchange matrix associated with specific localized occupied orbitals. The range separated potential is crucial to achieve efficient treatment of the direct term in the MP2, while extensive screening is employed to reduce the expense of the exchange contribution in MP2. The focus of this paper is on controllable accuracy and linear scaling of the data entering the algorithm. [ABSTRACT FROM AUTHOR]

In a recent paper, Lucco Castello et al. (arXiv:2107.03537) provided an accurate parameterization of classical one-component plasma bridge functions that was embedded in a novel dielectric scheme for strongly coupled electron liquids. Here, this approach is rigorously formulated, its set of equations is formally derived, and its numerical algorithm is scrutinized. A systematic comparison with available and new path integral Monte Carlo simulations reveals a rather unprecedented agreement especially in terms of the interaction energy and the long wavelength limit of the static local field correction. [ABSTRACT FROM AUTHOR]

A new method of the stochastic simulation algorithm (SSA), named the Hashing-Leaping method (HLM), for exact simulations of a class of Markov jump processes, is presented in this paper. The HLM has a conditional constant computational cost per event, which is independent of the number of exponential clocks in the Markov process. The main idea of the HLM is to repeatedly implement a hash-table-like bucket sort algorithm for all times of occurrence covered by a time step with length t. This paper serves as an introduction to this new SSA method. We introduce the method, demonstrate its implementation, analyze its properties, and compare its performance with three other commonly used SSA methods in four examples. Our performance tests and CPU operation statistics show certain advantages of the HLM for large scale problems. [ABSTRACT FROM AUTHOR]

This paper reports the release of PathSum, a new software suite of state-of-the-art path integral methods for studying the dynamics of single or extended systems coupled to harmonic environments. The package includes two modules, suitable for system-bath problems and extended systems comprising many coupled system-bath units, and is offered in C++ and Fortran implementations. The system-bath module offers the recently developed small matrix path integral (SMatPI) and the well-established iterative quasi-adiabatic propagator path integral (i-QuAPI) method for iteration of the reduced density matrix of the system. In the SMatPI module, the dynamics within the entanglement interval can be computed using QuAPI, the blip sum, time evolving matrix product operators, or the quantum-classical path integral method. These methods have distinct convergence characteristics and their combination allows a user to access a variety of regimes. The extended system module provides the user with two algorithms of the modular path integral method, applicable to quantum spin chains or excitonic molecular aggregates. An overview of the methods and code structure is provided, along with guidance on method selection and representative examples., (© 2023 Author(s). Published under an exclusive license by AIP Publishing.)

The standard fewest-switches surface hopping (FSSH) approach fails to model nonadiabatic dynamics when the electronic Hamiltonian is complex-valued and there are multiple nuclear dimensions; FSSH does not include geometric magnetic effects and does not have access to a gauge independent direction for momentum rescaling. In this paper, for the case of a Hamiltonian with two electronic states, we propose an extension of Tully's FSSH algorithm, which includes geometric magnetic forces and, through diabatization, establishes a well-defined rescaling direction. When combined with a decoherence correction, our new algorithm shows satisfying results for a model set of two-dimensional single avoided crossings. [ABSTRACT FROM AUTHOR]

The fruitful interplay of high-resolution spectroscopy and quantum chemistry has a long history, especially in the field of small, semi-rigid molecules. However, in recent years, the targets of spectroscopic studies are shifting toward flexible molecules, characterized by a large number of closely spaced energy minima, all contributing to the overall spectrum. Here, artificial intelligence comes into play since it is at the basis of powerful unsupervised techniques for the exploration of soft degrees of freedom. Integration of such algorithms with a two-stage QM/QM′ (Quantum Mechanical) exploration/refinement strategy driven by a user-friendly graphical interface is the topic of the present paper. We will address in particular: (i) the performances of different semi-empirical methods for the exploration step and (ii) the comparison between stochastic and meta-heuristic algorithms in achieving a cheap yet complete exploration of the conformational space for medium sized chromophores. As test cases, we choose three amino acids of increasing complexity, whose full conformer enumeration has been reached only very recently. Next, we show that systems in condensed phases can be treated at the same level and with the same efficiency when employing a polarizable continuum description of the solvent. Finally, the challenging issue represented by the vibrational circular dichroism spectra of some rhodium complexes with flexible ligands has been addressed, showing that our fully unsupervised approach leads to remarkable agreement with the experiment. [ABSTRACT FROM AUTHOR]

Molecular dynamics simulations require barostats to be performed at a constant pressure. The usual recipe is to employ the Berendsen barostat first, which displays a first-order volume relaxation efficient in equilibration but results in incorrect volume fluctuations, followed by a second-order or a Monte Carlo barostat for production runs. In this paper, we introduce stochastic cell rescaling, a first-order barostat that samples the correct volume fluctuations by including a suitable noise term. The algorithm is shown to report volume fluctuations compatible with the isobaric ensemble and its anisotropic variant is tested on a membrane simulation. Stochastic cell rescaling can be straightforwardly implemented in the existing codes and can be used effectively in both equilibration and production phases. [ABSTRACT FROM AUTHOR]

In this paper, we have analyzed the time series associated with the iterative scheme of a double similarity transformed coupled cluster theory. The coupled iterative scheme to solve the ground state Schrödinger equation is cast as a multivariate time-discrete map, and the solutions show the universal Feigenbaum dynamics. Using recurrence analysis, it is shown that the dynamics of the iterative process is dictated by a small subgroup of cluster operators, mostly those involving chemically active orbitals, whereas all other cluster operators with smaller amplitudes are enslaved. Using synergetics, we will indicate how the master-slave dynamics can suitably be exploited to develop a novel coupled-cluster algorithm in a much reduced dimension. [ABSTRACT FROM AUTHOR]

We present NECI, a state-of-the-art implementation of the Full Configuration Interaction Quantum Monte Carlo (FCIQMC) algorithm, a method based on a stochastic application of the Hamiltonian matrix on a sparse sampling of the wave function. The program utilizes a very powerful parallelization and scales efficiently to more than 24 000 central processing unit cores. In this paper, we describe the core functionalities of NECI and its recent developments. This includes the capabilities to calculate ground and excited state energies, properties via the one- and two-body reduced density matrices, as well as spectral and Green's functions for ab initio and model systems. A number of enhancements of the bare FCIQMC algorithm are available within NECI, allowing us to use a partially deterministic formulation of the algorithm, working in a spin-adapted basis or supporting transcorrelated Hamiltonians. NECI supports the FCIDUMP file format for integrals, supplying a convenient interface to numerous quantum chemistry programs, and it is licensed under GPL-3.0. [ABSTRACT FROM AUTHOR]

We have developed an algorithm to detect holes in multi-dimensional real-valued surfaces—such as the potential energy surfaces (PESs) that describe the nuclear motion of molecules in the context of the Born–Oppenheimer approximation. For our purposes, a PES "hole" is defined as an unphysical saddle point, beyond which the potential energy drops (typically) without limit to negative infinity. PES holes are numerical artifacts that can arise when fitting PES functional forms to discrete ab initio data—even when the data is of high quality, and/or for comparatively few degrees of freedom (DOF). Often undetected, PES holes can have devastating effects on subsequent dynamical calculations, especially if they occur at low energies. In this paper, we present a highly efficient algorithm designed to systematically identify hole configurations and energies. The method is applied to a variety of molecular PESs ranging up to 30 DOF. A number of evidently previously undetected PES holes are reported here—surprisingly, even for PESs that have been available for decades. The code itself (Crystal) is presented together with a user manual. These tools may be of great benefit for PES developers, who can use the information they provide to fix holes, once identified. More generally, the methodology can be applied in any context involving multi-dimensional surfaces. [ABSTRACT FROM AUTHOR]

We introduce a new, topologically-based method for coarse-graining polymer chains. Based on the wavelet transform, a multiresolution data analysis technique, this method assigns a cluster of particles to a coarse-grained bead located at the center of mass of the cluster, thereby reducing the complexity of the problem by dividing the simulation into several stages, each with a fraction of the number of beads as the overall chain. At each stage, we compute the distributions of coarse-grained internal coordinates as well as potential functions required for subsequent simulation stages. In this paper, we present the basic algorithm, and apply it to freely jointed chains; the companion paper describes its applications to self-avoiding chains. [ABSTRACT FROM AUTHOR]

In Paper I of this series the formulations of the replica-exchange multicanonical algorithm and the multicanonical replica-exchange method for Monte Carlo versions have been presented. The effectiveness of these algorithms were then tested with the system of a penta peptide, Met-enkephalin, in the gas phase. In this article the detailed comparisons of performances of these algorithms together with the regular replica-exchange method are made, taking a more complex system of a 17-residue helical peptide. It is shown that these two new algorithms are more efficient than the regular replica-exchange method. [ABSTRACT FROM AUTHOR]

A good understanding and determination of colloidal interactions is paramount to comprehend and model the thermodynamic and structural properties of colloidal suspensions. In concentrated aqueous suspensions of colloids with a titratable surface charge, this determination is, however, complicated by the density dependence of the effective pair potential due to both the many-body interactions and the charge regulation of the colloids. In addition, colloids generally present a size distribution which results in a virtually infinite combination of colloid pairs. In this paper, we develop two methods and describe the corresponding algorithms to solve this problem for arbitrary size distributions. An implementation in Nim is also provided. The methods, inspired by the seminal work of Torres et al., [J. Chem. Phys. 128, 154906 (2008)] are based on a generalization of the cell and renormalized jellium models to polydisperse suspensions of spherical colloids with a charge regulating boundary condition. The latter is described by the one-pK-Stern model. The predictions of the models are confronted to the equations of state of various commercially available silica dispersions. The renormalized Yukawa parameters (effective charges and screening lengths) are also calculated. The importance of size and charge polydispersity as well as the validity of these two models is discussed in light of the results. [ABSTRACT FROM AUTHOR]

This paper presents in detail our fast semistochastic heat-bath configuration interaction (SHCI) method for solving the many-body Schrödinger equation. We identify and eliminate computational bottlenecks in both the variational and perturbative steps of the SHCI algorithm. We also describe the parallelization and the key data structures in our implementation, such as the distributed hash table. The improved SHCI algorithm enables us to include in our variational wavefunction two orders of magnitude more determinants than has been reported previously with other selected configuration interaction methods. We use our algorithm to calculate an accurate benchmark energy for the chromium dimer with the X2C relativistic Hamiltonian in the cc-pVDZ-DK basis, correlating 28 electrons in 76 spatial orbitals. Our largest calculation uses two billion Slater determinants in the variational space and semistochastically includes perturbative contributions from at least trillions of additional determinants with better than 10−5 Ha statistical uncertainty. [ABSTRACT FROM AUTHOR]

The particle-particle random phase approximation (pp-RPA) has been used to investigate excitation problems in our recent paper [Y. Yang, H. van Aggelen, and W. Yang, J. Chem. Phys. 139, 224105 (2013)]. It has been shown to be capable of describing double, Rydberg, and charge transfer excitations, which are challenging for conventional time-dependent density functional theory (TDDFT). However, its performance on larger molecules is unknown as a result of its expensive O(N6) scaling. In this article, we derive and implement a Davidson iterative algorithm for the pp-RPA to calculate the lowest few excitations for large systems. The formal scaling is reduced to O(N4), which is comparable with the commonly used configuration interaction singles (CIS) and TDDFT methods. With this iterative algorithm, we carried out benchmark tests on molecules that are significantly larger than the molecules in our previous paper with a reasonably large basis set. Despite some self-consistent field convergence problems with ground state calculations of (N - 2)-electron systems, we are able to accurately capture lowest few excitations for systems with converged calculations. Compared to CIS and TDDFT, there is no systematic bias for the pp-RPA with the mean signed error close to zero. The mean absolute error of pp-RPA with B3LYP or PBE references is similar to that of TDDFT, which suggests that the pp-RPA is a comparable method to TDDFT for large molecules. Moreover, excitations with relatively large non-HOMO excitation contributions are also well described in terms of excitation energies, as long as there is also a relatively large HOMO excitation contribution. These findings, in conjunction with the capability of pp-RPA for describing challenging excitations shown earlier, further demonstrate the potential of pp-RPA as a reliable and general method to describe excitations, and to be a good alternative to TDDFT methods. [ABSTRACT FROM AUTHOR]

In the present paper, thermal conductivity coefficient of high-temperature molecular and atomic gases with excited electronic states is studied using both the kinetic theory algorithm developed by authors earlier and the well known simple expression for the thermal conductivity coefficient proposed by Eucken and generalized by Hirschfelder. The influence of large collision diameters of excited states on the thermal conductivity is discussed. The limit of validity of the Eucken correction is evaluated on the basis of the kinetic theory calculations; an improved model suitable for air species under high-temperature conditions is proposed. [ABSTRACT FROM AUTHOR]

In recent years, the machine learning techniques have shown great potent1ial in various problems from a multitude of disciplines, including materials design and drug discovery. The high computational speed on the one hand and the accuracy comparable to that of density functional theory on another hand make machine learning algorithms efficient for high-throughput screening through chemical and configurational space. However, the machine learning algorithms available in the literature require large training datasets to reach the chemical accuracy and also show large errors for the so-called outliers—the out-of-sample molecules, not well-represented in the training set. In the present paper, we propose a new machine learning algorithm for predicting molecular properties that addresses these two issues: it is based on a local model of interatomic interactions providing high accuracy when trained on relatively small training sets and an active learning algorithm of optimally choosing the training set that significantly reduces the errors for the outliers. We compare our model to the other state-of-the-art algorithms from the literature on the widely used benchmark tests. [ABSTRACT FROM AUTHOR]

The reaction-diffusion master equation (RDME) is a model that allows for efficient on-lattice simulation of spatially resolved stochastic chemical kinetics. Compared to off-lattice hard-sphere simulations with Brownian dynamics or Green's function reaction dynamics, the RDME can be orders of magnitude faster if the lattice spacing can be chosen coarse enough. However, strongly diffusion-controlled reactions mandate a very fine mesh resolution for acceptable accuracy. It is common that reactions in the same model differ in their degree of diffusion control and therefore require different degrees of mesh resolution. This renders mesoscopic simulation inefficient for systems with multiscale properties. Mesoscopic-microscopic hybrid methods address this problem by resolving the most challenging reactions with a microscale, off-lattice simulation. However, all methods to date require manual partitioning of a system, effectively limiting their usefulness as "black-box" simulation codes. In this paper, we propose a hybrid simulation algorithm with automatic system partitioning based on indirect a priori error estimates.We demonstrate the accuracy and efficiency of the method on models of diffusion-controlled networks in 3D. [ABSTRACT FROM AUTHOR]

Molecular Dynamics-Green's Function Reaction Dynamics (MD-GFRD) is a multiscale simulation method for particle dynamics or particle-based reaction-diffusion dynamics that is suited for systems involving low particle densities. Particles in a low-density region are just diffusing and not interacting. In this case, one can avoid the costly integration of microscopic equations of motion, such as molecular dynamics (MD), and instead turn to an event-based scheme in which the times to the next particle interaction and the new particle positions at that time can be sampled. At high (local) concentrations, however, e.g., when particles are interacting in a nontrivial way, particle positions must still be updated with small time steps of the microscopic dynamical equations. The efficiency of a multiscale simulation that uses these two schemes largely depends on the coupling between them and the decisions when to switch between the two scales. Here we present an efficient scheme for multi-scale MD-GFRDsimulations. It has been shownthatMD-GFRDschemes are more efficient than brute-force molecular dynamics simulations up to a molar concentration of 10² µM. In this paper, we showthat the choice of the propagation domains has a relevant impact on the computational performance. Domains are constructed using a local optimization of their sizes and a minimal domain size is proposed. The algorithm is shown to be more efficient than brute-force Brownian dynamics simulations up to a molar concentration of 10³ µM and is up to an order of magnitude more efficient compared with previous MD-GFRD schemes. [ABSTRACT FROM AUTHOR]

The chemical master equation (CME) is frequently used in systems biology to quantify the effects of stochastic fluctuations that arise due to biomolecular species with low copy numbers. The CME is a system of ordinary differential equations that describes the evolution of probability density for each population vector in the state-space of the stochastic reaction dynamics. For many examples of interest, this state-space is infinite, making it difficult to obtain exact solutions of the CME. To deal with this problem, the Finite State Projection (FSP) algorithm was developed by Munsky and Khammash [J. Chem. Phys. 124(4), 044104 (2006)], to provide approximate solutions to the CME by truncating the state-space. The FSP works well for finite time-periods but it cannot be used for estimating the stationary solutions of CMEs, which are often of interest in systems biology. The aim of this paper is to develop a version of FSP which we refer to as the stationary FSP (sFSP) that allows one to obtain accurate approximations of the stationary solutions of a CME by solving a finite linear-algebraic system that yields the stationary distribution of a continuous-time Markov chain over the truncated state-space.We derive bounds for the approximation error incurred by sFSP and we establish that under certain stability conditions, these errors can be made arbitrarily small by appropriately expanding the truncated state-space. We provide several examples to illustrate our sFSP method and demonstrate its efficiency in estimating the stationary distributions. In particular, we show that using a quantized tensor-train implementation of our sFSP method, problems admitting more than 100 x 106 states can be efficiently solved. [ABSTRACT FROM AUTHOR]

In the quest towards coarse-grained potentials and new water models, we present an extension of the force matching technique to parameterize an all-atom force field for rigid water. The methodology presented here allows to improve the matching procedure by first optimizing the weighting exponents present in the objective function. A new gauge for unambiguously evaluating the quality of the fit has been introduced; it is based on the root mean square difference of the distributions of target properties between reference data and fitted potentials. Four rigid water models have been parameterized; the matching procedure has been used to assess the role of the ghost atom in TIP4P-like models and of electrostatic damping. In the former case, burying the negative charge inside the molecule allows to fit better the torques. In the latter, since short-range interactions are damped, a better fit of the forces is obtained. Overall, the best performing model is the one with a ghost atom and with electrostatic damping. The approach shown in this paper is of general validity and could be applied to any matching algorithm and to any level of coarse graining, also for non-rigid molecules. [ABSTRACT FROM AUTHOR]

Two iterative procedures have been proposed recently to calculate the chemical potentials corresponding to prescribed concentrations from grand canonical Monte Carlo (GCMC) simulations. Both are based on repeated GCMC simulations with updated excess chemical potentials until the desired concentrations are established. In this paper, we propose combining our robust and fast converging iteration algorithm [Malasics, Gillespie, and Boda, J. Chem. Phys. 128, 124102 (2008)] with the suggestion of Lamperski [Mol. Simul. 33, 1193 (2007)] to average the chemical potentials in the iterations (instead of just using the chemical potentials obtained in the last iteration). We apply the unified method for various electrolyte solutions and show that our algorithm is more efficient if we use the averaging procedure. We discuss the convergence problems arising from violation of charge neutrality when inserting/deleting individual ions instead of neutral groups of ions (salts). We suggest a correction term to the iteration procedure that makes the algorithm efficient to determine the chemical potentials of individual ions too. [ABSTRACT FROM AUTHOR]

The weighted stochastic simulation algorithm (wSSA) recently introduced by Kuwahara and Mura [J. Chem. Phys. 129, 165101 (2008)] is an innovative variation on the stochastic simulation algorithm (SSA). It enables one to estimate, with much less computational effort than was previously thought possible using a Monte Carlo simulation procedure, the probability that a specified event will occur in a chemically reacting system within a specified time when that probability is very small. This paper presents some procedural extensions to the wSSA that enhance its effectiveness in practical applications. The paper also attempts to clarify some theoretical issues connected with the wSSA, including its connection to first passage time theory and its relation to the SSA. [ABSTRACT FROM AUTHOR]

While closed-loop control of quantum dynamics in the laboratory is proving to be broadly successful, the control mechanisms induced by the fields are often left obscure. Hamiltonian encoding (HE) was originally introduced as a method for understanding mechanisms in quantum dynamics in the context of computational simulations, based on access to the system wavefunction. As a step towards laboratory implementation of HE, this paper addresses the issues raised by the use of observables rather than the wavefunction in HE. The goal of laboratory based HE is to obtain an understanding of control mechanism through a sequence of systematic control experiments, whose collective information can identify the underlying control mechanism defined as the set of significant amplitudes connecting the initial and final states. Mechanism is determined by means of observing the dynamics of special sequences of system Hamiltonians encoded through the control field. The proposed algorithm can handle complex systems, operates with no recourse to dynamical simulations, and functions with limited understanding of the system Hamiltonian. As with the closed-loop control experiments, the HE control mechanism identification algorithm performs a new experiment each time the dynamical outcome from an encoded Hamiltonian is called for. This paper presents the basic HE algorithm in the context of physical systems described by a finite dimensional Hilbert space. The method is simulated with simple models, and the extension to more complex systems is discussed. [ABSTRACT FROM AUTHOR]

Local coupled-cluster theory provides an algorithm for measuring electronic correlation quickly, using only the spatial locality of localized electronic orbitals. Previously, we showed [J. Subotnik et al., J. Chem. Phys. 125, 074116 (2006)] that one may construct a local coupled-cluster singles-doubles theory which (i) yields smooth potential energy surfaces and (ii) achieves near linear scaling. That theory selected which orbitals to correlate based only on the distances between the centers of different, localized orbitals, and the approximate potential energy surfaces were characterized as smooth using only visual identification. This paper now extends our previous algorithm in three important ways. First, locality is now based on both the distances between the centers of orbitals as well as the spatial extent of the orbitals. We find that, by accounting for the spatial extent of a delocalized orbital, one can account for electronic correlation in systems with some electronic delocalization using fast correlation methods designed around orbital locality. Second, we now enforce locality on not just the amplitudes (which measure the exact electron-electron correlation), but also on the two-electron integrals themselves (which measure the bare electron-electron interaction). Our conclusion is that we can bump integrals as well as amplitudes, thereby gaining a tremendous increase in speed and paradoxically increasing the accuracy of our LCCSD approach. Third and finally, we now make a rigorous definition of chemical smoothness as requiring that potential energy surfaces not support artificial maxima, minima, or inflection points. By looking at first and second derivatives from finite difference techniques, we demonstrate complete chemical smoothness of our potential energy surfaces (bumping both amplitudes and integrals). These results are significant both from a theoretical and from a computationally practical point of view. [ABSTRACT FROM AUTHOR]

In this paper, we present a path integral hybrid Monte Carlo (PIHMC) method for rotating molecules in quantum fluids. This is an extension of our PIHMC for correlated Bose fluids [S. Miura and J. Tanaka, J. Chem. Phys. 120, 2160 (2004)] to handle the molecular rotation quantum mechanically. A novel technique referred to be an effective potential of quantum rotation is introduced to incorporate the rotational degree of freedom in the path integral molecular dynamics or hybrid Monte Carlo algorithm. For a permutation move to satisfy Bose statistics, we devise a multilevel Metropolis method combined with a configurational-bias technique for efficiently sampling the permutation and the associated atomic coordinates. Then, we have applied the PIHMC to a helium-4 cluster doped with a carbonyl sulfide molecule. The effects of the quantum rotation on the solvation structure and energetics were examined. Translational and rotational fluctuations of the dopant in the superfluid cluster were also analyzed. [ABSTRACT FROM AUTHOR]

In this paper, we introduce algorithms based on Fourier synthesis for designing phase compensating rf pulse sequences for high-resolution nuclear magnetic resonance (NMR) spectroscopy in an inhomogeneous B0 field. We show that using radio frequency pulses and time varying linear gradients in three dimensions, it is possible to impart the transverse magnetization of spins phase, which is a desired function of the spatial (x,y,z) location. Such a sequence can be used to precompensate the phase that will be acquired by spins at different spatial locations due to inhomogeneous magnetic fields. With this precompensation, the chemical shift information of the spins can be reliably extracted and high resolution NMR spectrum can be obtained. [ABSTRACT FROM AUTHOR]

The explicit tau-leaping procedure attempts to speed up the stochastic simulation of a chemically reacting system by approximating the number of firings of each reaction channel during a chosen time increment τ as a Poisson random variable. Since the Poisson random variable can have arbitrarily large sample values, there is always the possibility that this procedure will cause one or more reaction channels to fire so many times during τ that the population of some reactant species will be driven negative. Two recent papers have shown how that unacceptable occurrence can be avoided by replacing the Poisson random variables with binomial random variables, whose values are naturally bounded. This paper describes a modified Poisson tau-leaping procedure that also avoids negative populations, but is easier to implement than the binomial procedure. The new Poisson procedure also introduces a second control parameter, whose value essentially dials the procedure from the original Poisson tau-leaping at one extreme to the exact stochastic simulation algorithm at the other; therefore, the modified Poisson procedure will generally be more accurate than the original Poisson procedure. [ABSTRACT FROM AUTHOR]

Genetic algorithm global optimization of Ag-Pd, Ag-Au, and Pd-Pt clusters is performed. The 34- and 38-atom clusters are optimized for all compositions. The atom-atom interactions are modeled by a semiempirical potential. All three systems are characterized by a small size mismatch and a weak tendency of the larger atoms to segregate at the surface of the smaller ones. As a result, the global minimum structures exhibit a larger mixing than in Ag-Cu and Ag-Ni clusters. Polyicosahedral structures present generally favorable energetic configurations, even though they are less favorable than in the case of the size-mismatched systems. A comparison between all the systems studied here and in the previous paper (on size-mismatched systems) is presented. [ABSTRACT FROM AUTHOR]

This paper proposes methods for calculating the derivative couplings between adiabatic states in density-functional theory (DFT) and compares them with each other and with multiconfigurational self-consistent field calculations. They are shown to be accurate and, as expected, the costs of their calculation scale more favorably with system size than post-Hartree-Fock calculations. The proposed methods are based on single-particle excitations and the associated Slater transition-state densities to overcome the problem of the unavailability of multielectron states in DFT which precludes a straightforward calculation of the matrix elements of the nuclear gradient operator. An iterative scheme employing linear-response theory was found to offer the best trade-off between accuracy and efficiency. The algorithms presented here have been implemented for doublet-doublet excitations within a plane-wave-basis and pseudopotential framework but are easily generalizable to other excitations and basis sets. Owing to their fundamental importance in cases where the Born-Oppenheimer separation of motions is not valid, these derivative couplings can facilitate, for example, the treatment of nonadiabatic charge transfers, of electron-phonon couplings, and of radiationless electronic transitions in DFT.© 2005 American Institute of Physics. [ABSTRACT FROM AUTHOR]

Most of the numerical simulation in quantum (bilinear) control have used one of the monotonically convergent algorithms of Krotov (introduced by Tannor et al.) or of Zhu and Rabitz. However, until now no explicit relationship has been revealed between the two algorithms in order to understand their common properties. Within this framework, we propose in this paper a unified formulation that comprises both algorithms and that extends to a new class of monotonically convergent algorithms. Numerical results show that the newly derived algorithms behave as well as (and sometimes better than) the well-known algorithms cited above. © 2003 American Institute of Physics. [ABSTRACT FROM AUTHOR]

The fourfold way is a general algorithm for generating diabatic electronic wave functions that span the same space as a small set of variationally optimized adiabatic electronic wave functions and for using the resulting diabatic wave functions to generate diabatic potential energy surfaces and their couplings. In this paper we extend the fourfold way so it is applicable to more complex polyatomic systems and in particular to the calculation of global potential energy surfaces for such systems. The extension involves partitioning the active space into three blocks, introducing restricted orbital rotation within two of the blocks, introducing a specific resolution of the subspace containing molecular orbitals that are doubly occupied in all dominant configuration state functions, and introducing specific orientations of the coordinate systems for reference molecular orbitals and resolution molecular orbitals. The major strength of the improved method presented in this paper is that it allows the diabatic molecular orbitals to exhibit a gradual change of chemical character with smooth deformation along the reaction coordinate for a change of chemical arrangement while preserving the orbital character required for a physical ordering of the orbitals. This feature is required for the convenient construction of global potential energy surfaces for non-BornOppenheimer rearrangements. The resulting extended algorithm is illustrated by calculating diabatic potential energy surfaces and couplings for the two lowest singlet potential energy surfaces of HNCO. [ABSTRACT FROM AUTHOR]

Molecular dynamics simulations are useful in obtaining thermodynamic and kinetic properties of bio-molecules, but they are limited by the time scale barrier. That is, we may not obtain properties' efficiently because we need to run microseconds or longer simulations using femtosecond time steps. To overcome this time scale barrier, we can use the weighted ensemble (WE) method, a powerful enhanced sampling method that efficiently samples thermodynamic and kinetic properties. However, the WE method requires an appropriate partitioning of phase space into discrete macrostates, which can be problematic when we have a high-dimensional collective space or when little is known a priori about the molecular system. Hence, we developed a new WE-based method, called the "Concurrent Adaptive Sampling (CAS) algorithm," to tackle these issues. The CAS algorithm is not constrained to use only one or two collective variables, unlike most reaction coordinate-dependent methods. Instead, it can use a large number of collective variables and adaptive macrostates to enhance the sampling in the high-dimensional space. This is especially useful for systems in which we do not know what the right reaction coordinates are, in which case we can use many collective variables to sample conformations and pathways. In addition, a clustering technique based on the committor function is used to accelerate sampling the slowest process in the molecular system. In this paper, we introduce the new method and show results from two-dimensional models and bio-molecules, specifically penta-alanine and a triazine trimer. [ABSTRACT FROM AUTHOR]

Massively parallel algorithms are presented in this paper to reduce the computational burden associated with quantum transport simulations from first-principles. The power of modern hybrid computer architectures is harvested in order to determine the open boundary conditions that connect the simulation domain with its environment and to solve the resulting Schrödinger equation. While the former operation takes the form of an eigenvalue problem that is solved by a contour integration technique on the available central processing units (CPUs), the latter can be cast into a linear system of equations that is simultaneously processed by SplitSolve, a two-step algorithm, on general-purpose graphics processing units (GPUs). A significant decrease of the computational time by up to two orders of magnitude is obtained as compared to standard solution methods. [ABSTRACT FROM AUTHOR]

Expanding the set of stable, accurate, and scalable methods for simulating molecular quantum dynamics is important for accelerating the computational exploration of molecular processes. In this paper, we adapt the signed particles Monte Carlo algorithm for solving the transient Wigner equation to scenarios of chemical interest. This approach was used in the past to study electronic processes in semi-conductors, but to the best of our knowledge, it had never been applied to molecular modeling. We present the algorithm and demonstrate its excellent performance on harmonic and double well potentials for electronic and nuclear systems. We explore the stability of the algorithm, discuss the choice of hyper-parameters, and cautiously speculate that it may be used in quantum molecular dynamics simulations. [ABSTRACT FROM AUTHOR]

Nuclear long-lived spin states represent spin density operator configurations that are exceptionally well protected against spin relaxation phenomena. Their long-lived character is exploited in a variety of Nuclear Magnetic Resonance (NMR) techniques. Despite the growing importance of long-lived spin states in modern NMR, strategies for their identification have changed little over the last decade. The standard approach heavily relies on a chain of group theoretical arguments. In this paper, we present a more streamlined method for the calculation of such configurations. Instead of focusing on the symmetry properties of the relaxation superoperator, we focus on its corresponding relaxation algebra. This enables us to analyze long-lived spin states with Lie algebraic methods rather than group theoretical arguments. We show that the centralizer of the relaxation algebra forms a basis for the set of long-lived spin states. The characterization of the centralizer, on the other hand, does not rely on any special symmetry arguments, and its calculation is straightforward. We outline a basic algorithm and illustrate advantages by considering long-lived spin states for some spin-1/2 pairs and rapidly rotating methyl groups. [ABSTRACT FROM AUTHOR]

In this paper we optimize several algorithms for the computation of reaction rates based on information calculated along minimum energy reaction paths and we evaluate the efficiencies of the optimized algorithms. The investigations are based on the calculation of chemical reaction rate constants using variational transition state theory and multidimensional semiclassical transmission coefficients including reaction path curvature. Several methods are evaluated and compared by a systematic set of applications to test cases involving the hydrogen-atom transfer reactions CH3+H2→CH4+H and OH+H2→H2O+H. For each method we present general recommendations for all algorithmic choices other than gradient step size so that future calculations may be carried out reasonably efficiently by varying only one parameter. In the process of these optimizations we have found that the accuracy of the Euler stabilization method can be significantly increased by choosing the auxiliary parameters differently than in previous work; the optimized algorithm is called ES1*. Our final recommendations for future work are (i) when the Hessian/gradient computational cost ratio is low (<=3): the Page–McIver algorithm with the Hessian recalculated at every step, with a cubic starting step, and with curvature calculated from the derivative of the gradient, and (ii) when the Hessian/gradient computational cost ratio is moderate or large: the ES1* algorithm with a Hessian step size three times larger than the gradient step size, with a quadratic starting step, and with curvature calculated from the derivative of the gradient. [ABSTRACT FROM AUTHOR]

This paper presents a new, quantum mechanical, time dependent approach to the photodissociation of triatomic molecules in Jacobi coordinates. The algorithm is based on a nested interaction representation, designed to make the representation of the time evolving wave packet as compact as possible. The new equations of motion are solved numerically using a synthesis of grid techniques: the fast Fourier transform (FFT) method is used in Cartesian-like coordinates, and the discrete variable representation (DVR) method in the angular or bending coordinate. A variant on the short iterative Lanczos (SIL) procedure is used for the temporal propagation of the wave packet. Rotational state distributions obtained from this new algorithm are presented for the single surface photodissociation of ClCN and for the two surface photodissociation of ICN. The ClCN results are in good agreement with the semiclassical results of Barts and Halpern [J. Phys. Chem. 93, 7346 (1989)] and in excellent agreement with the time independent quantum results of Schinke [J. Chem. Phys. 92, 2397 (1990)]. Rotational state distributions for the two electronic surface photodissociation of ICN are in good agreement with the time independent quantum results of Guo and Schatz [J. Chem. Phys. 92, 1634 (1990)] and illustrate the flexibility of the method for dealing with nonadiabatic processes. The numerical efficiency of the method is comparable with standard time independent techniques, but has the attractive feature of yielding final state distributions at all energies from a single wave packet propagation. [ABSTRACT FROM AUTHOR]

Variational basis-set formulations of the quantum mechanical reactive scattering problem lead to large, dense sets of equations. In previous work, we showed that the generalized minimal residual (GMRes) algorithm is sometimes competitive in terms of computer time with direct methods for these dense matrices, even when complex-valued boundary conditions are used, leading to non-Hermitian matrices. This paper presents a preconditioning scheme to accelerate convergence and improve performance. We block the potential energy coupling into a series of distortion blocks, and we employ the outgoing wave variational principle with nonorthogonal basis functions, including both dynamically adapted Green’s functions for the distortion blocks and also square integrable functions. The coefficient matrix of the resulting linear system couples the blocks. We have found that preconditioners formed from diagonal blocks of the coefficient matrix corresponding to the distortion blocks and vibrational blocks are effective at accelerating the iterative method in every test case, by factors of 2.9–20, with an average speedup of a factor of 6.5. The storage requirements and computational efficiency of the new scheme compare favorably to those for preconditioners based on banded matrices of variable bandwidth. The new preconditioners yield converged transition probabilities in less computer time than a direct solver even in cases which do not converge in a reasonable amount of time without preconditioning, and the average speedup compared to the direct solution is a factor of 7.6. [ABSTRACT FROM AUTHOR]

In the present study, variational path integral molecular dynamics and associated hybrid Monte Carlo (HMC) methods have been developed on the basis of a fourth order approximation of a density operator. To reveal various parameter dependence of physical quantities, we analytically solve one dimensional harmonic oscillators by the variational path integral; as a byproduct, we obtain the analytical expression of the discretized density matrix using the fourth order approximation for the oscillators. Then, we apply our methods to realistic systems like a water molecule and a para-hydrogen cluster. In the HMC, we adopt two level description to avoid the time consuming Hessian evaluation. For the systems examined in this paper, the HMC method is found to be about three times more efficient than the molecular dynamics method if appropriate HMC parameters are adopted; the advantage of the HMC method is suggested to be more evident for systems described by many body interaction. [ABSTRACT FROM AUTHOR]

The frequently used reduction technique is based on the chemical master equation for stochastic chemical kinetics with two-time scales, which yields the modified stochastic simulation algorithm (SSA). For the chemical reaction processes involving a large number of molecular species and reactions, the collection of slow reactions may still include a large number of molecular species and reactions. Consequently, the SSA is still computationally expensive. Because the chemical Langevin equations (CLEs) can effectively work for a large number of molecular species and reactions, this paper develops a reduction method based on the CLE by the stochastic averaging principle developed in the work of Khasminskii and Yin [SIAM J. Appl. Math. 56, 1766-1793 (1996); ibid. 56, 1794-1819 (1996)] to average out the fast-reacting variables. This reduction method leads to a limit averaging system, which is an approximation of the slow reactions. Because in the stochastic chemical kinetics, the CLE is seen as the approximation of the SSA, the limit averaging system can be treated as the approximation of the slow reactions. As an application, we examine the reduction of computation complexity for the gene regulatory networks with two-time scales driven by intrinsic noise. For linear and nonlinear protein production functions, the simulations show that the sample average (expectation) of the limit averaging system is close to that of the slow-reaction process based on the SSA. It demonstrates that the limit averaging system is an efficient approximation of the slow-reaction process in the sense of the weak convergence. [ABSTRACT FROM AUTHOR]

In this paper, we present a novel hybrid method to simulate discrete stochastic reaction-diffusion models arising in biochemical signaling pathways. We study moderately stiff systems, for which we can partition each reaction or diffusion channel into either a slow or fast subset, based on its propensity. Numerical approaches missing this distinction are often limited with respect to computational run time or approximation quality. We design an approximate scheme that remedies these pitfalls by using a new blending strategy of the well-established inhomogeneous stochastic simulation algorithm and the tau-leaping simulation method. The advantages of our hybrid simulation algorithm are demonstrated on three benchmarking systems, with special focus on approximation accuracy and efficiency. [ABSTRACT FROM AUTHOR]

In this paper, we describe a new and fully coherent stochastic surface hopping method for simulating mixed quantum-classical systems. We illustrate the approach on the simple but unforgiving problem of quantum evolution of a two-state quantum system in the limit of unperturbed pure state dynamics and for dissipative evolution in the presence of both stationary and nonstationary random environments. We formulate our approach in the Liouville representation and describe the density matrix elements by ensembles of trajectories. Population dynamics are represented by stochastic surface hops for trajectories representing diagonal density matrix elements. These are combined with an unconventional coherent stochastic hopping algorithm for trajectories representing off-diagonal quantum coherences. The latter generalizes the binary (0,1) "probability" of a trajectory to be associated with a given state to allow integers that can be negative or greater than unity in magnitude. Unlike existing surface hopping methods, the dynamics of the ensembles are fully entangled, correctly capturing the coherent and nonlocal structure of quantum mechanics. [ABSTRACT FROM AUTHOR]