Volume 145, Issue 7, 21 August 2016

We use a combination of coarsegrained molecular dynamics simulations and theoretical modeling to examine threejunctions in mixed lipid bilayer membranes. These junctions are localized defect lines in which three bilayers merge in such a way that each bilayer shares one monolayer with one of the other two bilayers. The resulting local morphology is nonlamellar, resembling the threefold symmetric defect lines in inverse hexagonal phases, but it regularly occurs during membrane fission and fusion events. We realize a system of junctions by setting up a honeycomb lattice, which in its primitive cell contains two hexagons and four threeline junctions, permitting us to study their stability as well as their line tension. We specifically consider the effects of lipid composition and intrinsic curvature in binary mixtures, which contain a fraction of negatively curved lipids in a curvatureneutral background phase. Threejunction stability results from a competition between the junction and an open edge, which arises if one of the three bilayers detaches from the other two. We show that the stable phase is the one with the lower defect line tension. The strong and opposite monolayer curvatures present in junctions and edges enhance the mole fraction of negatively curved lipids in junctions and deplete it in edges. This lipid sorting affects the two line tensions and in turn the relative stability of the two phases. It also leads to a subtle entropic barrier for the transition between junction and edge that is absent in uniform membranes.
 COMMUNICATIONS


Communication: Fitting potential energy surfaces with fundamental invariant neural network
View Description Hide DescriptionA more flexible neural network (NN) method using the fundamental invariants (FIs) as the input vector is proposed in the construction of potential energy surfaces for molecular systems involving identical atoms. Mathematically, FIs finitely generate the permutation invariant polynomial (PIP) ring. In combination with NN, fundamental invariant neural network (FINN) can approximate any function to arbitrary accuracy. Because FINN minimizes the size of input permutation invariant polynomials, it can efficiently reduce the evaluation time of potential energy, in particular for polyatomic systems. In this work, we provide the FIs for all possible molecular systems up to five atoms. Potential energy surfaces for OH3 and CH4 were constructed with FINN, with the accuracy confirmed by fulldimensional quantum dynamic scattering and bound state calculations.
 Top

 ARTICLES

 Theoretical Methods and Algorithms

Finite state projection based bounds to compare chemical master equation models using singlecell data
View Description Hide DescriptionEmerging techniques now allow for precise quantification of distributions of biological molecules in single cells. These rapidly advancing experimental methods have created a need for more rigorous and efficient modeling tools. Here, we derive new bounds on the likelihood that observations of singlecell, singlemolecule responses come from a discrete stochastic model, posed in the form of the chemical master equation. These strict upper and lower bounds are based on a finite state projection approach, and they converge monotonically to the exact likelihood value. These bounds allow one to discriminate rigorously between models and with a minimum level of computational effort. In practice, these bounds can be incorporated into stochastic model identification and parameter inference routines, which improve the accuracy and efficiency of endeavors to analyze and predict singlecell behavior. We demonstrate the applicability of our approach using simulated data for three example models as well as for experimental measurements of a timevarying stochastic transcriptional response in yeast.

Bootstrap embedding: An internally consistent fragmentbased method
View Description Hide DescriptionStrong correlation poses a difficult problem for electronic structure theory, with computational cost scaling quickly with system size. Fragment embedding is an attractive approach to this problem. By dividing a large complicated system into smaller manageable fragments “embedded” in an approximate description of the rest of the system, we can hope to ameliorate the steep cost of correlated calculations. While appealing, these methods often converge slowly with fragment size because of small errors at the boundary between fragment and bath. We describe a new electronic embedding method, dubbed “Bootstrap Embedding,” a selfconsistent wavefunctioninwavefunction embedding theory that uses overlapping fragments to improve the description of fragment edges. We apply this method to the one dimensional Hubbard model and a translationally asymmetric variant, and find that it performs very well for energies and populations. We find Bootstrap Embedding converges rapidly with embedded fragment size, overcoming the surfaceareatovolumeratio error typical of many embedding methods. We anticipate that this method may lead to a lowscaling, high accuracy treatment of electron correlation in large molecular systems.

Explicitly correlated coupledcluster theory with Brueckner orbitals
View Description Hide DescriptionBrueckner orbitals are the optimal orbitals for use in F12 explicitly correlated coupledcluster (CC) treatments. A novel approach, Brueckner coupledcluster doubles with perturbative triples BCCD(T)(F12*) is presented that avoids the expensive reevaluation of F12 integrals throughout the orbital optimisation and includes a newly derived basis set correction to the Brueckner reference energy. The generalisation of F12 theory to arbitrary nonHartree–Fock references and to Fock operators that include scalar relativistic effects is also presented. The performance of the new Brueckner F12 method is assessed for a test set of 50 open and closedshell reactions and for the ionisation potentials and electron affinities (EAs) of the firstrow transition metal atoms. Benchmark basis set limit coupledcluster singles, doubles and perturbative triples (CCSD(T)) and BCCD(T) values are reported for all energies in the test sets. BCCD(T)(F12*) performs systematically better than CCSD(T)(F12*) for electron affinities where orbital relaxation effects are significant.

Electron correlation within the relativistic nopair approximation
View Description Hide DescriptionThis paper addresses the definition of correlation energy within 4component relativistic atomic and molecular calculations. In the nonrelativistic domain the correlation energy is defined as the difference between the exact eigenvalue of the electronic Hamiltonian and the HartreeFock energy. In practice, what is reported is the basis set correlation energy, where the “exact” value is provided by a full Configuration Interaction (CI) calculation with some specified oneparticle basis. The extension of this definition to the relativistic domain is not straightforward since the corresponding electronic Hamiltonian, the DiracCoulomb Hamiltonian, has no bound solutions. Presentday relativistic calculations are carried out within the nopair approximation, where the DiracCoulomb Hamiltonian is embedded by projectors eliminating the troublesome negativeenergy solutions. HartreeFock calculations are carried out with the implicit use of such projectors and only positiveenergy orbitals are retained at the correlated level, meaning that the HartreeFock projectors are frozen at the correlated level. We argue that the projection operators should be optimized also at the correlated level and that this is possible by full Multiconfigurational SelfConsistent Field (MCSCF) calculations, that is, MCSCF calculations using a nopair full CI expansion, but including orbital relaxation from the negativeenergy orbitals. We show by variational perturbation theory that the MCSCF correlation energy is a pure MP2like correlation expression, whereas the corresponding CI correlation energy contains an additional relaxation term. We explore numerically our theoretical analysis by carrying out variational and perturbative calculations on the twoelectron rare gas atoms with specially tailored basis sets. In particular, we show that the correlation energy obtained by the suggested MCSCF procedure is smaller than the nopair full CI correlation energy, in accordance with the underlying minmax principle and our theoretical analysis. We also show that the relativistic correlation energy, obtained from nopair full MCSCF calculations, scales at worst as X ^{−2} with respect to the cardinal number X of our correlationconsistent basis sets optimized for the twoelectron atoms. This is better than the X ^{−1} scaling suggested by previous studies, but worse than the X ^{−3} scaling observed in the nonrelativistic domain. The wellknown 1/Z expansion in nonrelativistic atomic theory follows from coordinate scaling. We point out that coordinate scaling for consistency should be accompanied by velocity scaling. In the nonrelativistic domain this comes about automatically, whereas in the relativistic domain an explicit scaling of the speed of light is required. This in turn explains why the relativistic correlation energy to the lowest order is not independent of nuclear charge, in contrast to nonrelativistic theory.

On the description of Brownian particles in confinement on a nonCartesian coordinates basis
View Description Hide DescriptionWe developed a theoretical framework to study the diffusion of Brownian pointlike particles in bounded geometries in two and three dimensions. We use the FrenetSerret moving frame as the coordinate system. For narrow tubes and channels, we use an effective onedimensional description reducing the diffusion equation to a FickJacobslike equation. From this last equation, we can calculate the effective diffusion coefficient applying Neumann boundary conditions. On one hand, for channels with a straight axis our theoretical approximation for the effective coefficient does coincide with the reported in the literature [D. Reguera and J. M. Rubí, Phys. Rev. E 64, 061106 (2001) and P. Kalinay and J. K. Percus, ibid. 74, 041203 (2006)]. On the other hand, for tubes with a straight axis and circular crosssection our analytical expression does not coincide with the reported by Rubí and Reguera and by Kalinay and Percus, although it is practically identical.

Acceleration of saddlepoint searches with machine learning
View Description Hide DescriptionIn atomistic simulations, the location of the saddle point on the potentialenergy surface (PES) gives important information on transitions between local minima, for example, via transitionstate theory. However, the search for saddle points often involves hundreds or thousands of ab initio force calls, which are typically all done at full accuracy. This results in the vast majority of the computational effort being spent calculating the electronic structure of states not important to the researcher, and very little time performing the calculation of the saddle point state itself. In this work, we describe how machine learning (ML) can reduce the number of intermediate ab initio calculations needed to locate saddle points. Since machinelearning models can learn from, and thus mimic, atomistic simulations, the saddlepoint search can be conducted rapidly in the machinelearning representation. The saddlepoint prediction can then be verified by an ab initio calculation; if it is incorrect, this strategically has identified regions of the PES where the machinelearning representation has insufficient training data. When these training data are used to improve the machinelearning model, the estimates greatly improve. This approach can be systematized, and in two simple example problems we demonstrate a dramatic reduction in the number of ab initio force calls. We expect that this approach and future refinements will greatly accelerate searches for saddle points, as well as other searches on the potential energy surface, as machinelearning methods see greater adoption by the atomistics community.

Resonance energy transfer: The unified theory via vector spherical harmonics
View Description Hide DescriptionIn this work, we derive the wellestablished expression for the quantum amplitude associated with the resonance energy transfer (RET) process between a pair of molecules that are beyond wavefunction overlap. The novelty of this work is that the field of the mediating photon is described in terms of a spherical wave rather than a plane wave. The angular components of the field are constructed in terms of vector spherical harmonics while Hankel functions are used to define the radial component. This approach alleviates the problem of having to select physically correct solution from nonphysical solutions, which seems to be inherent in plane wave derivations. The spherical coordinate system allows one to easily decompose the photon’s fields into longitudinal and transverse components and offers a natural way to analyse near, intermediate, and farzone RET within the context of the relative orientation of the transition dipole moments for the two molecules.

Real space electrostatics for multipoles. III. Dielectric properties
View Description Hide DescriptionIn Papers I and II, we developed new shifted potential, gradient shifted force, and Taylor shifted force realspace methods for multipole interactions in condensed phase simulations. Here, we discuss the dielectric properties of fluids that emerge from simulations using these methods. Most electrostatic methods (including the Ewald sum) require correction to the conducting boundary fluctuation formula for the static dielectric constants, and we discuss the derivation of these corrections for the new real space methods. For quadrupolar fluids, the analogous material property is the quadrupolar susceptibility. As in the dipolar case, the fluctuation formula for the quadrupolar susceptibility has corrections that depend on the electrostatic method being utilized. One of the most important effects measured by both the static dielectric and quadrupolar susceptibility is the ability to screen charges embedded in the fluid. We use potentials of mean force between solvated ions to discuss how geometric factors can lead to distancedependent screening in both quadrupolar and dipolar fluids.

Systemsize corrections for selfdiffusion coefficients calculated from molecular dynamics simulations: The case of CO2, nalkanes, and poly(ethylene glycol) dimethyl ethers
View Description Hide DescriptionMolecular dynamics simulations were carried out to study the selfdiffusion coefficients of CO2, methane, propane, nhexane, nhexadecane, and various poly(ethylene glycol) dimethyl ethers (glymes in short, CH3O–(CH2CH2O)n–CH3 with n = 1, 2, 3, and 4, labeled as G1, G2, G3, and G4, respectively) at different conditions. Various system sizes were examined. The widely used Yeh and Hummer [J. Phys. Chem. B 108, 15873 (2004)] correction for the prediction of diffusion coefficient at the thermodynamic limit was applied and shown to be accurate in all cases compared to extrapolated values at infinite system size. The magnitude of correction, in all cases examined, is significant, with the smallest systems examined giving for some cases a selfdiffusion coefficient approximately 15% lower than the infinite systemsize extrapolated value. The results suggest that finite size corrections to computed selfdiffusivities must be used in order to obtain accurate results.

Relativistic equationofmotion coupledcluster method using openshell reference wavefunction: Application to ionization potential
View Description Hide DescriptionThe openshell reference relativistic equationofmotion coupledcluster method within its fourcomponent description is successfully implemented with the consideration of single and double excitation approximations using the DiracCoulomb Hamiltonian. At the first attempt, the implemented method is employed to calculate ionization potential value of heavy atomic (Ag, Cs, Au, Fr, and Lr) and molecular (HgH and PbF) systems, where the effect of relativity does really matter to obtain highly accurate results. Not only the relativistic effect but also the effect of electron correlation is crucial in these heavy atomic and molecular systems. To justify the fact, we have taken two further approximations in the fourcomponent relativistic equationofmotion framework to quantify how the effect of electron correlation plays a role in the calculated values at different levels of theory. All these calculated results are compared with the available experimental data as well as with other theoretically calculated values to judge the extent of accuracy obtained in our calculations.

Necessary conditions for the emergence of homochirality via autocatalytic selfreplication
View Description Hide DescriptionWe analyze a recent proposal for spontaneous mirror symmetry breaking based on the coupling of firstorder enantioselective autocatalysis and direct production of the enantiomers that invokes a critical role for intrinsic reaction noise. For isolated systems, the racemic state is the unique stable outcome for both stochastic and deterministic dynamics when the system is in compliance with the constraints dictated by the thermodynamics of chemical reaction processes. In open systems, the racemic outcome also results for both stochastic and deterministic dynamics when driving the autocatalysis unidirectionally by external reagents. Nonracemic states can result in the latter only if the reverse reactions are strictly zero: these are kinetically controlled outcomes for small populations and volumes, and can be simulated by stochastic dynamics. However, the stability of the thermodynamic limit proves that the racemic outcome is the unique stable state for strictly irreversible externally driven autocatalysis. These findings contradict the suggestion that the inhibition requirement of the Frank autocatalytic model for the emergence of homochirality may be relaxed in a noiseinduced mechanism.

An efficient algorithm for finding the minimum energy path for cation migration in ionic materials
View Description Hide DescriptionThe Nudged Elastic Band (NEB) is an established method for finding minimumenergy paths and energy barriers of ion migration in materials, but has been hampered in its general application by its significant computational expense when coupled with density functional theory (DFT) calculations. Typically, an NEB calculation is initialized from a linear interpolation of successive intermediate structures (also known as images) between known initial and final states. However, the linear interpolation introduces two problems: (1) slow convergence of the calculation, particularly in cases where the final path exhibits notable curvature; (2) divergence of the NEB calculations if any intermediate image comes too close to a nondiffusing species, causing instabilities in the ensuing calculation. In this work, we propose a new scheme to accelerate NEB calculations through an improved path initialization and associated energy estimation workflow. We demonstrate that for cation migration in an ionic framework, initializing the diffusion path as the minimum energy path through a static potential built upon the DFT charge density reproduces the true NEB path within a 0.2 Å deviation and yields up to a 25% improvement in typical NEB runtimes. Furthermore, we find that the locally relaxed energy barrier derived from this initialization yields a good approximation of the NEB barrier, with errors within 20 meV of the true NEB value, while reducing computational expense by up to a factor of 5. Finally, and of critical importance for the automation of migration path calculations in highthroughput studies, we find that the new approach significantly enhances the stability of the calculation by avoiding unphysical image initialization. Our algorithm promises to enable efficient calculations of diffusion pathways, resolving a longstanding obstacle to the computational screening of intercalation compounds for Liion and multivalent batteries.

Generalized average local ionization energy and its representations in terms of Dyson and energy orbitals
View Description Hide DescriptionRyabinkin and Staroverov [J. Chem. Phys. 141, 084107 (2014)] extended the concept of average local ionization energy (ALIE) to correlated wavefunctions by defining the generalized ALIE as , where λj are the eigenvalues of the generalized Fock operator and fj(r) are the corresponding eigenfunctions (energy orbitals). Here we show that one can equivalently express the generalized ALIE as , where Ik are singleelectron removal energies and dk(r) are the corresponding Dyson orbitals. The two expressions for emphasize different physical interpretations of this quantity; their equivalence enables one to calculate the ALIE at any level of ab initio theory without generating the computationally expensive Dyson orbitals.

Variational path integral molecular dynamics and hybrid Monte Carlo algorithms using a fourth order propagator with applications to molecular systems
View Description Hide DescriptionIn 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 parahydrogen 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.

A multiscale transport model for LennardJones binary mixtures based on interfacial friction
View Description Hide DescriptionWe propose a onedimensional isothermal hydrodynamic transport model for nonreacting binary mixtures in slit shaped nanochannels. The coupled species momentum equations contain viscous dissipation and interspecies friction term of MaxwellStefan form. Species partial viscosity variations in the confinement are modeled using the van der Waals one fluid approximation and the local average density method. Species specific macroscopic friction coefficient based Robin boundary conditions are provided to capture the species wall slip effects. The value of this friction coefficient is computed using a species specific generalized Langevin formulation. Gravity driven flow of methanehydrogen and methaneargon mixtures confined between graphene slit shaped nanochannels are considered as examples. The proposed model yields good quantitative agreement with the velocity profiles obtained from the nonequilibrium molecular dynamics simulations. The mixtures considered are observed to behave as single species pseudo fluid, with the interfacial friction displaying linear dependence on molar composition of the mixture. The results also indicate that the different species have different slip lengths, which remain unchanged with the channel width.

Spotting the difference in molecular dynamics simulations of biomolecules
View Description Hide DescriptionComparing two trajectories from molecular simulations conducted under different conditions is not a trivial task. In this study, we apply a method called Linear Discriminant Analysis with ITERative procedure (LDAITER) to compare two molecular simulation results by finding the appropriate projection vectors. Because LDAITER attempts to determine a projection such that the projections of the two trajectories do not overlap, the comparison does not suffer from a strong anisotropy, which is an issue in protein dynamics. LDAITER is applied to two test cases: the T4 lysozyme protein simulation with or without a point mutation and the allosteric protein PDZ2 domain of hPTP1E with or without a ligand. The projection determined by the method agrees with the experimental data and previous simulations. The proposed procedure, which complements existing methods, is a versatile analytical method that is specialized to find the “difference” between two trajectories.

Membrane curvature generated by asymmetric depletion layers of ions, small molecules, and nanoparticles
View Description Hide DescriptionBiomimetic and biological membranes consist of molecular bilayers with two leaflets that are typically exposed to different aqueous solutions. We consider solutions of “particles” that experience effectively repulsive interactions with these membranes and form depletion layers in front of the membrane leaflets. The particles considered here are watersoluble, have a size between a few angstrom and a few nanometers as well as a rigid, more or less globular shape, and do neither adsorb onto the membranes nor permeate these membranes. Examples are provided by ions, small sugar molecules, globular proteins, or inorganic nanoparticles with a hydrophilic surface. We first study depletion layers in a hardcore system based on ideal particle solutions as well as hardwall interactions between these particles and the membrane. For this system, we obtain exact expressions for the coverages and tensions of the two leaflets as well as for the spontaneous curvature of the bilayer membrane. All of these quantities depend linearly on the particle concentrations. The exact results for the hardcore system also show that the spontaneous curvature can be directly deduced from the planar membrane geometry. Our results for the hardcore system apply both to ions and solutes that are small compared to the membrane thickness and to nanoparticles with a size that is comparable to the membrane thickness, provided the particle solutions are sufficiently dilute. We then corroborate the different relationships found for the hardcore system by extensive simulations of a softcore particle system using dissipative particle dynamics. The simulations confirm the linear relationships obtained for the hardcore system. Both our analytical and our simulation results show that the spontaneous curvature induced by a single particle species can be quite large. When one leaflet of the membrane is exposed, e.g., to a 100 mM solution of glucose, a lipid bilayer can acquire a spontaneous curvature of ±1/(270 nm). Our theoretical results can be scrutinized by systematic experimental studies using a large variety of different types of particles.
 Advanced Experimental Techniques

A CMOS millimeterwave transceiver embedded in a semiconfocal FabryPerot cavity for molecular spectroscopy
View Description Hide DescriptionThe extension of radio frequency complementary metal oxide semiconductor (CMOS) circuitry into millimeter wavelengths promises the extension of spectroscopic techniques in compact, power efficient systems. We are now beginning to use CMOS millimeter devices for lowmass, lowpower instrumentation capable of remote or in situ detection of gas composition during space missions. We have chosen to develop a FlygareBalle type spectrometer, with a semiconfocal FabryPerot cavity to amplify the pump power of a mmwavelength CMOS transmitter that is directly coupled to the planar mirror of the cavity. We have built a pulsed transceiver system at 92105 GHz inside a 3 cm base length cavity and demonstrated quality factor up to 4680, allowing for modes with 20 MHz bandwidth, with a sufficient cavity amplification factor for mW class transmitters. This work describes the initial gas measurements and outlines the challenges and next steps.
 Atoms, Molecules, and Clusters

Unifying the rotational and permutation symmetry of nuclear spin states: SchurWeyl duality in molecular physics
View Description Hide DescriptionIn modern physics and chemistry concerned with manybody systems, one of the mainstays is identicalparticlepermutation symmetry. In particular, both the intramolecular dynamics of a single molecule and the intermolecular dynamics associated, for example, with reactive molecular collisions are strongly affected by selection rules originating in nuclearpermutation symmetry operations being applied to the total internal wavefunctions, including nuclear spin, of the molecules involved. We propose here a general tool to determine coherently the permutation symmetry and the rotational symmetry (associated with the group of arbitrary rotations of the entire molecule in space) of molecular wavefunctions, in particular the nuclearspin functions. Thus far, these two symmetries were believed to be mutually independent and it has even been argued that under certain circumstances, it is impossible to establish a onetoone correspondence between them. However, using the SchurWeyl duality theorem we show that the two types of symmetry are inherently coupled. In addition, we use the ingenious representationtheory technique of Young tableaus to represent the molecular nuclearspin degrees of freedom in terms of welldefined mathematical objects. This simplifies the symmetry classification of the nuclear wavefunction even for large molecules. Also, the application to reactive collisions is very straightforward and provides a much simplified approach to obtaining selection rules.