Critical behavior in a chiral molecular model.

Understanding the condensed-phase behavior of chiral molecules is important in biology as well as in a range of technological applications, such as the manufacture of pharmaceuticals. Here, we use molecular dynamics simulations to study a chiral four-site molecular model that exhibits a second-order symmetry-breaking phase transition from a supercritical racemic liquid into subcritical D-rich and L-rich liquids. We determine the infinite-size critical temperature using the fourth-order Binder cumulant, and we show that the finite-size scaling behavior of the order parameter is compatible with the 3D Ising universality class. We also study the spontaneous D-rich to L-rich transition at a slightly subcritical temperature of T = 0.985Tc, and our findings indicate that the free energy barrier for this transformation increases with system size as N2/3, where N is the number of molecules, consistent with a surface-dominated phenomenon. The critical behavior observed herein suggests a mechanism for chirality selection in which a liquid of chiral molecules spontaneously forms a phase enriched in one of the two enantiomers as the temperature is lowered below the critical point. Furthermore, the increasing free energy barrier with system size indicates that fluctuations between the L-rich and D-rich phases are suppressed as the size of the system increases, trapping it in one of the two enantiomerically enriched phases. Such a process could provide the basis for an alternative explanation for the origin of biological homochirality. We also conjecture the possibility of observing nucleation at subcritical temperatures under the action of a suitable chiral external field.

Realizability of iso-g2 processes via effective pair interactions.

An outstanding problem in statistical mechanics is the determination of whether prescribed functional forms of the pair correlation function g2(r) [or equivalently, structure factor S(k)] at some number density ρ can be achieved by many-body systems in d-dimensional Euclidean space. The Zhang-Torquato conjecture states that any realizable set of pair statistics, whether from a nonequilibrium or equilibrium system, can be achieved by equilibrium systems involving up to two-body interactions. To further test this conjecture, we study the realizability problem of the nonequilibrium iso-g2 process, i.e., the determination of density-dependent effective potentials that yield equilibrium states in which g2 remains invariant for a positive range of densities. Using a precise inverse algorithm that determines effective potentials that match hypothesized functional forms of g2(r) for all r and S(k) for all k, we show that the unit-step function g2, which is the zero-density limit of the hard-sphere potential, is remarkably realizable up to the packing fraction ϕ = 0.49 for d = 1. For d = 2 and 3, it is realizable up to the maximum "terminal" packing fraction ϕc = 1/2d, at which the systems are hyperuniform, implying that the explicitly known necessary conditions for realizability are sufficient up through ϕc. For ϕ near but below ϕc, the large-r behaviors of the effective potentials are given exactly by the functional forms exp[ - κ(ϕ)r] for d = 1, r-1/2 exp[ - κ(ϕ)r] for d = 2, and r-1 exp[ - κ(ϕ)r] (Yukawa form) for d = 3, where κ-1(ϕ) is a screening length, and for ϕ = ϕc, the potentials at large r are given by the pure Coulomb forms in the respective dimensions as predicted by Torquato and Stillinger [Phys. Rev. E 68, 041113 (2003)]. We also find that the effective potential for the pair statistics of the 3D "ghost" random sequential addition at the maximum packing fraction ϕc = 1/8 is much shorter ranged than that for the 3D unit-step function g2 at ϕc; thus, it does not constrain the realizability of the unit-step function g2. Our inverse methodology yields effective potentials for realizable targets, and, as expected, it does not reach convergence for a target that is known to be non-realizable, despite the fact that it satisfies all known explicit necessary conditions. Our findings demonstrate that exploring the iso-g2 process via our inverse methodology is an effective and robust means to tackle the realizability problem and is expected to facilitate the design of novel nanoparticle systems with density-dependent effective potentials, including exotic hyperuniform states of matter.

Fluid-fluid phase transitions in a chiral molecular model.

Molecular chirality is a fundamental phenomenon, underlying both life as we know it and industrial pharmaceutical syntheses. Understanding the symmetry breaking phase transitions exhibited by many chiral molecular substances provides basic insights for topics ranging from the origin of life to the rational design of drug manufacturing processes. In this work, we have performed molecular dynamics simulations to investigate the fluid-fluid phase transitions of a flexible three-dimensional four-site chiral molecular model developed by Latinwo et al. [J. Chem. Phys. 145, 154503 (2016)] and Petsev et al. [J. Chem. Phys. 155, 084105 (2021)]. By introducing a bias favoring local homochiral vs heterochiral interactions, the system exhibits a phase transition from a single achiral phase to a single chiral phase that undergoes infrequent interconversion between the two thermodynamically identical chiral states: the L-rich and D-rich phases. According to the phase rule, this reactive binary system has two independent degrees of freedom and exhibits a density-dependent critical locus. Below the liquid-liquid critical locus, there exists a first-order vapor-liquid coexistence region with a single independent degree of freedom. Our results provide basic thermodynamic and kinetic insights for understanding many-body chiral symmetry breaking phenomena.

Crystal Prediction via Genetic Algorithms in a Model Chiral System.

Chiral crystals and their constituent molecules play a prominent role in theories about the origin of biological homochirality and in drug discovery, design, and stability. Although the prediction and identification of stable chiral crystal structures is crucial for numerous technologies, including separation processes and polymorph selection and control, predictive ability is often complicated by a combination of many-body interactions and molecular complexity and handedness. In this work, we address these challenges by applying genetic algorithms to predict the ground-state crystal lattices formed by a chiral tetramer molecular model, which we have previously shown to exhibit complex fluid-phase behavior. Using this approach, we explore the relative stability and structures of the model's conglomerate and racemic crystals, and present a structural phase diagram for the stable Bravais crystal types in the zero-temperature limit.

Interconversion-controlled liquid-liquid phase separation in a molecular chiral model.

Liquid-liquid phase separation of fluids exhibiting interconversion between alternative states has been proposed as an underlying mechanism for fluid polyamorphism and may be of relevance to the protein function and intracellular organization. However, molecular-level insight into the interplay between competing forces that can drive or restrict phase separation in interconverting fluids remains elusive. Here, we utilize an off-lattice model of enantiomers with tunable chiral interconversion and interaction properties to elucidate the physics underlying the stabilization and tunability of phase separation in fluids with interconverting states. We show that introducing an imbalance in the intermolecular forces between two enantiomers results in nonequilibrium, arrested phase separation into microdomains. We also find that in the equilibrium case, when all interaction forces are conservative, the growth of the phase domain is restricted only by the system size. In this case, we observe phase amplification, in which one of the two alternative phases grows at the expense of the other. These findings provide novel insights on how the interplay between dynamics and thermodynamics defines the equilibrium and steady-state morphologies of phase transitions in fluids with interconverting molecular or supramolecular states.

Effect of configuration-dependent multi-body forces on interconversion kinetics of a chiral tetramer model.

We describe a reformulation of the four-site molecular model for chiral phenomena introduced by Latinwo et al. ["Molecular model for chirality phenomena," J. Chem. Phys. 145, 154503 (2016)]. The reformulation includes an additional eight-body force that arises from an explicit configuration-dependent term in the potential energy function, resulting in a coarse-grained energy-conserving force field for molecular dynamics simulations of chirality phenomena. In this model, the coarse-grained interaction energy between two tetramers depends on their respective chiralities and is controlled by a parameter λ, where λ < 0 favors local configurations involving tetramers of opposite chirality and λ > 0 gives energetic preference to configurations involving tetramers of the same chirality. We compute the autocorrelation function for a quantitative chirality metric and demonstrate that the multi-body force modifies the interconversion kinetics such that λ ≠ 0 increases the effective barrier for enantiomer inversion. Our simulations reveal that for λ > 0 and temperatures below a sharply defined threshold value, this effect is dramatic, giving rise to spontaneous chiral symmetry breaking and locking molecules into their chiral identity.

Effects of Trehalose on Lipid Membranes under Rapid Cooling using All-Atom and Coarse-Grained Molecular Simulations.

We investigate the effect of the cryopreservative α-α-trehalose on a model 1,2-dimyristoyl-sn-glycero-3-phosphocholine (DMPC) lipid membrane undergoing cooling from 350 to 250 K using all-atom (AA) and coarse-grained (CG) molecular dynamics simulation. In the AA simulations, we find that the addition of trehalose alters the Lα (liquid crystalline) to Pß (ripple) phase transition, suppressing the major domain of the Pß phase and increasing the degree of leaflet interdigitation (the minor domain) which yields a thinner membrane with a higher area per lipid. Calculation of dihedral angle distributions for the lipid tails shows a greater fraction of gauche angles in the Pß phase as trehalose concentration is increased, indicating that trehalose increases lipid disorder in the membrane. In contrast, the CG simulations transition directly from the Lα to the Lß (gel) phase upon cooling without exhibiting the Pß phase (likely due to increased lipid mobility in the CG system). Even so, the CG simulations show that the addition of trehalose clearly suppresses the Lα to Lß phase transition, demonstrating that trehalose increases lipid disorder at low temperatures for the CG system, similar to the AA. Analysis using a two-state binding model provides net affinity coefficients between trehalose and the membrane as well as trehalose partition coefficients between the membrane interface and the bulk solution for both the AA and CG systems.

Kinetic Frustration Effects on Dense Two-Dimensional Packings of Convex Particles and Their Structural Characteristics.

The study of hard-particle packings is of fundamental importance in physics, chemistry, cell biology, and discrete geometry. Much of the previous work on hard-particle packings concerns their densest possible arrangements. By contrast, we examine kinetic effects inevitably present in both numerical and experimental packing protocols. Specifically, we determine how changing the compression/shear rate of a two-dimensional packing of noncircular particles causes it to deviate from its densest possible configuration, which is always periodic. The adaptive shrinking cell (ASC) optimization scheme maximizes the packing fraction of a hard-particle packing by first applying random translations and rotations to the particles and then isotropically compressing and shearing the simulation box repeatedly until a possibly jammed state is reached. We use a stochastic implementation of the ASC optimization scheme to mimic different effective time scales by varying the number of particle moves between compressions/shears. We generate dense, effectively jammed, monodisperse, two-dimensional packings of obtuse scalene triangle, rhombus, curved triangle, lens, and "ice cream cone" (a semicircle grafted onto an isosceles triangle) shaped particles, with a wide range of packing fractions and degrees of order. To quantify these kinetic effects, we introduce the kinetic frustration index K, which measures the deviation of a packing from its maximum possible packing fraction. To investigate how kinetics affect short- and long-range ordering in these packings, we compute their spectral densities χ̃V(k) and characterize their contact networks. We find that kinetic effects are most significant when the particles have greater asphericity, less curvature, and less rotational symmetry. This work may be relevant to the design of laboratory packing protocols.

Thermodynamics of DNA Hybridization from Atomistic Simulations.

Studying DNA hybridization equilibrium at atomistic length scales, either via molecular dynamics (MD) or through commonly used advanced sampling approaches, is notoriously difficult. In this work, we describe an order-parameter-based advanced sampling technique to calculate the free energy of hybridization, and estimate the melting temperature of DNA oligomers at atomistic resolution. The free energy landscapes are reported as a function of a native-topology-based order parameter for the Drew-Dickerson dodecamer and for a range of DNA decamer sequences of different GC content. Our estimated melting temperatures match the experimental numbers within ±15 °C. As a test of the numerical reliability of the procedures employed, it was verified that the predicted free energy surfaces and melting temperatures of the d- and l-enantiomers of the Drew-Dickerson dodecamer were indistinguishable within numerical accuracy.

Genetic Algorithm Approach for the Optimization of Protein Antifreeze Activity Using Molecular Simulations.

Antifreeze proteins (AFPs) are of much interest for their ability to inhibit ice growth at low concentrations. In this work, we present a genetic algorithm for the in silico design of AFP mutants with improved antifreeze activity, measured as the predicted thermal hysteresis at a fixed concentration, ΔTC. Central to the algorithm is our recently developed neural network method for predicting ΔTC from molecular simulations [Kozuch et al., PNAS, 115, 13252 (2018)]. Applying the algorithm to three structurally diverse AFPs, wfAFP, rQAE, and RiAFP, we find that significantly improved mutants are discovered for rQAE and RiAFP. Testing of the optimized mutants shows an increase in ΔTC of 0.572 ± 0.11 K (262 ± 50.6%) and 1.33 ± 0.14 K (39.9 ± 4.19%) over the native structures for rQAE and RiAFP, respectively. Structural analysis of the optimized mutants reveals that the algorithm is able to exploit two pathways for enhancing the predicted antifreeze activity of the mutants: (1) increasing the local order of surface waters by encouraging the formation of internal water channels in the protein and (2) increasing the total ice-binding area by improving the planar structure of the ice-binding surface. Additionally, analysis of all mutants explored by the algorithm reveals that a subset of residues, mainly nonpolar, are particularly helpful in improving antifreeze activity at the ice-binding surface.

Sensitivity of pair statistics on pair potentials in many-body systems.

We study the sensitivity and practicality of Henderson's theorem in classical statistical mechanics, which states that the pair potential v(r) that gives rise to a given pair correlation function g2(r) [or equivalently, the structure factor S(k)] in a classical many-body system at number density ρ and temperature T is unique up to an additive constant. While widely invoked in inverse-problem studies, the utility of the theorem has not been quantitatively scrutinized to any large degree. We show that Henderson's theorem has practical shortcomings for disordered and ordered phases for certain densities and temperatures. Using proposed sensitivity metrics, we identify illustrative cases in which distinctly different potential functions give very similar pair correlation functions and/or structure factors up to their corresponding correlation lengths. Our results reveal that due to a limited range and precision of pair information in either direct or reciprocal space, there is effective ambiguity of solutions to inverse problems that utilize pair information only, and more caution must be exercised when one claims the uniqueness of any resulting effective pair potential found in practice. We have also identified systems that possess virtually identical pair statistics but have distinctly different higher-order correlations. Such differences should be reflected in their individually distinct dynamics (e.g., glassy behaviors). Finally, we prove a more general version of Henderson's theorem that extends the uniqueness statement to include potentials that involve two- and higher-body interactions.

The Handedness of DNA Assembly around Carbon Nanotubes Is Determined by the Chirality of DNA.

Single-stranded DNA chains enable postsynthesis sorting of single-walled carbon nanotubes (CNTs) according to their diameter and helicity by wrapping helically around CNT surfaces. Both DNA chains and CNTs in these CNT-DNA conjugates are intrinsically chiral. Using a single-stranded DNA chain in both of its chiral realizations, we systematically study cross-chiral interactions between DNA and CNTs by varying the helicity of CNTs within a relatively narrow range of diameters. We find that regardless of the helicity or handedness of the carbon nanotube, the chirality of DNA dictates the handedness of its predominant helical wrap around carbon nanotubes.

Computational investigation of retro-isomer equilibrium structures: Intrinsically disordered, foldable, and cyclic peptides.

We use all-atom modeling and advanced-sampling molecular dynamics simulations to investigate quantitatively the effect of peptide bond directionality on the equilibrium structures of four linear (two foldable, two disordered) and two cyclic peptides. We find that the retro forms of cyclic and foldable linear peptides adopt distinctively different conformations compared to their parents. While the retro form of a linear intrinsically disordered peptide with transient secondary structure fails to reproduce a secondary structure content similar to that of its parent, the retro form of a shorter disordered linear peptide shows only minor differences compared to its parent.

Low temperature protein refolding suggested by molecular simulation.

The function of critical biological materials, such as proteins, is intrinsically tied to their structure, and this structure is in turn heavily dependent on the properties of the solvent, most commonly water or dilute aqueous solutions. As water is known to exhibit anomalous properties, especially at supercooled temperatures, it is natural to ask how these properties might impact the thermodynamics of protein folding. To investigate this question, we use molecular simulation to explore the behavior of a model miniprotein, Trp-cage, as low as 70 K below the freezing point of the solvent at ambient pressure. Surprisingly, we find that while the expected cold denaturation of the protein is observed at moderate supercooling, further cooling to more than 55 K below the freezing point leads to cold refolding of the protein. Structural and hydrogen bonding analysis suggests that this refolding is driven by the desolvation of the protein's hydrophobic core, likely related to the pronounced decrease in density at this temperature. Beyond their intrinsic fundamental interest, these results have implications for cryomicroscopy and cryopreservation, where biological materials are often transiently subjected to these extreme conditions.

Structural degeneracy in pair distance distributions.

A traditional basic descriptor of many-particle systems has been the distribution of interparticle pair distances. In the case of structureless particles at thermal equilibrium, with just additive pair interactions, this suffices to determine pressure and mean energy. However, it is usually the case that a given set of pair distances can emerge from a multiplicity of distinguishable many-particle configurations. This paper focuses on the ways in which such a configuration detail can be overlooked. After providing some elementary small-system examples in which full pair distance specification still permits distinct configurational pattern ambiguity, subsequent analysis concentrates on large-system classical canonical ensembles. In that context, configurational degeneracy is analyzed in two-dimensional systems for the shape distribution of triangles, whose chirality occurrence can be controlled by suitable three-particle interactions. For many-particle systems in three dimensions, the possibility is explored that a set of three-particle "pair-invisible" interactions can exist which modify the three-particle distribution function, but which have no effect on the pair distribution function, and thus remain undetected by conventional diffraction experiments. For illustration, a specific mathematical example is presented, applicable to the case where two-particle interactions vanish.

Effect of heterochiral inversions on the structure of a ß-hairpin peptide.

We study computationally a family of ß-hairpin peptides with systematically introduced chiral inversions, in explicit water, and we investigate the extent to which the backbone structure is able to fold in the presence of heterochiral perturbations. In contrast to the recently investigated case of a helical peptide, we do not find a monotonic change in secondary structure content as a function of the number of L- to D-inversions. The effects of L- to D-inversions are instead found to be highly position-specific. Additionally, in contrast to the helical peptide, some inversions increase the stability of the folded peptide: in such cases, we compute an increase in ß-sheet content in the aqueous solution equilibrium ensemble. However, the tertiary structures of the stable (folded) configurations for peptides for which inversions cause an increase in ß-sheet content show differences from one another, as well as from the native fold of the nonchirally perturbed ß-hairpin. Our results suggest that although some chiral perturbations can increase folding stability, chirally perturbed proteins may still underperform functionally, given the relationship between structure and function.

Combined molecular dynamics and neural network method for predicting protein antifreeze activity.

Antifreeze proteins (AFPs) are a diverse class of proteins that depress the kinetically observable freezing point of water. AFPs have been of scientific interest for decades, but the lack of an accurate model for predicting AFP activity has hindered the logical design of novel antifreeze systems. To address this, we perform molecular dynamics simulation for a collection of well-studied AFPs. By analyzing both the dynamic behavior of water near the protein surface and the geometric structure of the protein, we introduce a method that automatically detects the ice binding face of AFPs. From these data, we construct a simple neural network that is capable of quantitatively predicting experimentally observed thermal hysteresis from a trio of relevant physical variables. The model's accuracy is tested against data for 17 known AFPs and 5 non-AFP controls.

Computational Investigation of the Effect of Backbone Chiral Inversions on Polypeptide Structure.

Studying a set of helix-folding polyalanine peptides with systematically inserted chiral inversions in explicit water, we investigate quantitatively the effect of chiral perturbations on the structural ensembles of the peptides, thereby assessing the extent to which the backbone structure is able to fold in the presence of systematic heterochiral perturbations. Starting from the homochiral l-Ala20 peptide, we invert the backbone chiralities of Ala residues one by one along a specific perturbation pathway, until reaching the homochiral d-Ala20 peptide. Analysis of the helical contents of the simulated structural ensembles of the peptides shows that even a single inversion in the middle of the peptide completely breaks the helical structure in its vicinity and drastically reduces the helical content of the peptide. Further inversions in the middle of the peptide monotonically decrease the original helical content, that is, the right-handed helical content for l-Ala, and increase the helical content of the opposite chirality. Further analysis of the peptide ensembles using several size- and shape-related order parameters also indicate the drastic global changes in the peptide structure due to the local effects caused by the chiral inversions, such as formation of a reverse turn. However, the degree of the structural changes introduced by opposite chirality substitutions depends on the position of the inversion.

Rational design of stealthy hyperuniform two-phase media with tunable order.

Disordered stealthy hyperuniform materials are exotic amorphous states of matter that have attracted recent attention because of their novel structural characteristics (hidden order at large length scales) and physical properties, including desirable photonic and transport properties. It is therefore useful to devise algorithms that enable one to design a wide class of such amorphous configurations at will. In this paper, we present several algorithms enabling the systematic identification and generation of discrete (digitized) stealthy hyperuniform patterns with a tunable degree of order, paving the way towards the rational design of disordered materials endowed with novel thermodynamic and physical properties. To quantify the degree of order or disorder of the stealthy systems, we utilize the discrete version of the τ order metric, which accounts for the underlying spatial correlations that exist across all relevant length scales in a given digitized two-phase (or, equivalently, a two-spin state) system of interest. Our results impinge on a myriad of fields, ranging from physics, materials science and engineering, visual perception, and information theory to modern data science.

Cavitation transition in the energy landscape: Distinct tensile yielding behavior in strongly and weakly attractive systems.

Particles with cohesive interactions display a tensile instability in the energy landscape at the Sastry density ρS. The signature of this tensile limit is a minimum in the landscape equation of state, the pressure-density relationship of inherent structures sampled along a liquid isotherm. Our previous work [Y. E. Altabet, F. H. Stillinger, and P. G. Debenedetti, J. Chem. Phys. 145, 211905 (2016)] revisited the phenomenology of Sastry behavior and found that the evolution of the landscape equation of state with system size for particles with interactions typical of molecular liquids indicates the presence of an athermal first-order phase transition between homogeneous and fractured inherent structures, the latter containing several large voids. Here, we study how this tensile limit manifests itself for different interparticle cohesive strengths and identify two distinct regimes. Particles with sufficiently strong cohesion display an athermal first-order phase transition, consistent with our prior characterization. Weak cohesion also displays a tensile instability. However, the landscape equation of state for this regime is independent of system size, suggesting the absence of a first-order phase transition. An analysis of the voids suggests that yielding in the energy landscape of weakly cohesive systems is associated with the emergence of a highly interconnected network of small voids. While strongly cohesive systems transition from exclusively homogeneous to exclusively fractured configurations at ρS in the thermodynamic limit, this interconnected network develops gradually, starting at ρS, even at infinite system size.