He Ma, Marco Govoni, and Giulia Galli

*J. Open Source Softw.* 5(47), 2160 (2020) DOI:10.21105/joss.02160

In this work we describe the code PyZFS for the calculation of the ZFS tensor D of molecules and solids, based on wavefunctions obtained from density functional theory (DFT) calculations.

Hythem Sidky, Wei Chen, and Andrew L. Ferguson

*Mol. Phys.* 118, 1-22 (2020). DOI: 10.1080/00268976.2020.1737742

Classical molecular dynamics simulates the time evolution of molecular systems through the phase space spanned by the positions and velocities of the constituent atoms. Molecular-level thermodynamic, kinetic, and structural data extracted from the resulting trajectories provide valuable information for the understanding, engineering, and design of biological and molecular materials. The cost of simulating many-body atomic systems makes simulations of large molecules prohibitively expensive, and the high-dimensionality of the resulting trajectories presents a challenge for analysis. Driven by advances in algorithms, hardware, and data availability, there has been a flare of interest in recent years in the applications of machine learning – especially deep learning – to molecular simulation. These techniques have demonstrated great power and flexibility in both extracting mechanistic understanding of the important nonlinear collective variables governing the dynamics of a molecular system, and in furnishing good low-dimensional system representations with which to perform enhanced sampling or develop long-timescale dynamical models. It is the purpose of this article to introduce the key machine learning approaches, describe how they are married with statistical mechanical theory into domain-specific tools, and detail applications of these approaches in understanding and accelerating biomolecular simulation.

Emre Sevgen, Ashley Z. Guo, Hythem Sidky, Jonathan K. Whitmer, and Juan J. de Pablo

*J. Chem. Theory Comput.* 16, 1448-1455 (2020). DOI: 10.1021/acs.jctc.9b00883

An adaptive, machine learning-based sampling method is presented for simulation of systems having rugged, multidimensional free energy landscapes. The method’s main strength resides in its ability to learn both from the frequency of visits to distinct states and the generalized force estimates that arise in a system as it evolves in phase space. This is accomplished by introducing a self-integrating artificial neural network, which generates an estimate of the free energy directly from its derivatives. The usefulness of the proposed combined approach is examined in the context of two concrete examples, namely, an alanine dipeptide molecule in water and a polymer diffusing through a narrow pore. This new method is found to be robust, faster, and more accurate than approaches that rely only on frequency-based or generalized force-based estimations. After combining the proposed approach with overfill protection and support for sparse data storage and training, the method is shown to be more effective than comparable, previously available techniques and capable of scaling efficiently to larger numbers of collective variables.

Raul Laasner, Xiaochen Du, Aditya Tanikanti, Connor Clayton, Marco Govoni, Giulia Galli, Matti Ropo, and Volker Blum

*J. Open Source Softw.* 5(45), 1945 (2020). DOI: 10.21105/joss.01945

MatD3 is an open-source, dedicated database and web application framework designed tostore, curate and disseminate experimental and theoretical materials data generated by indi-vidual research groups or research consortia. A research group can set up its own instance ofMatD3and publish scientific results or simply use an existing online MatD3instance. Dissem-inating research data in this form enables broader access, reproducibility, and repurposing ofscientific products. MatD3is a general purpose database that does not focus on any specificlevel of theory or experimental method. Instead, the focus is on storing and making accessiblethe data and making it straightforward to curate them.

Nicholas E. Jackson, Alec S. Bowen, and Juan J. de Pablo

*Macromolecules* 53, 482-490 (2020). DOI: 10.1021/acs.macromol.9b02020

Conjugated polymers represent a high-potential material class due to the tunability of their optoelectronic properties via straightforward processing protocols. However, to in silico tailor these optoelectronic properties, one must compute the conformationally dependent electronic structure over mesoscopic length and time scales. This task represents a challenging multiscale computational problem in which quantum chemistry and atomistic or coarse-grained molecular dynamics must be integrated. Recently, we introduced a methodology, called electronic coarse graining (ECG), that utilizes artificial neural networks to compute ab initio quality electronic structures using only the reduced degrees of freedom of coarse-grained models. Here, we adapt ECG to variable-molecular weight conjugated oligomers by casting the problem in the framework of sequence detection. A machine-learning (ML) architecture employing one-dimensional convolutional neural networks and long short-term memory networks is utilized to compute ground-state orbital energies, charge density distributions, and optical spectra solely from the coarse-grained model’s configurational degrees of freedom. Robust molecular weight transferability for ECG is established via a Δ-ML approach that leverages model electronic Hamiltonians for ground and excited-state property determination. The accuracy and transferability of the ECG methodology opens the door for scalable optoelectronic property prediction in conjugated polymers directly from coarse-grained degrees of freedom.

Han Yang, Marco Govoni, and Giulia Galli

*J. Chem. Phys.* 151, 224102 (2019). DOI: 10.1063/1.5126214

Recently it was shown that the calculation of quasiparticle energies using the G0W0 approximation can be performed without computing explicitly any virtual electronic states, by expanding the Green function and screened Coulomb interaction in terms of the eigenstates of the static dielectric matrix. Avoiding the evaluation of virtual electronic states leads to improved efficiency and ease: of convergence of G0W0 calculations. Here we propose a further improvement of the efficiency of these calculations, based on an approximation of density-density response functions of molecules and solids. The approximation relies on the calculation of a subset of eigenvectors of the dielectric matrix using the kinetic operator instead of the full Hamiltonian, and it does not lead to any substantial loss of accuracy for the quasiparticle energies. The computational savings introduced by this approximation depend on the system, and they become more substantial as the number of electrons increases.

Anne C. Leonhard and Jonathan K. Whitmer

*J. Chem. Theory Comput.* 15, 5761-5768 (2019). DOI: 10.1021/acs.jctc.9b00348

Calculation of binding free energies, which characterize the stability and specificity of molecular interactions (e.g., drug formulations), has been established as an important application of modern computational chemistry. Though several methodologies currently exist within the literature, there is always a need for improvements in speed and accuracy that allow for efficient and accurate chemical screening. In this paper we explore the application of an unrestrained advanced sampling method, adaptive biasing force (ABF), to the calculation of standard binding free energies for a set of molecular host–guest systems for which the binding free energy is well-known: the cucurbit[7]uril-hosted systems utilized in the HYDROPHOBE challenge [Assaf et al. J. Phys. Chem. B 2017, 121, 11144−11162]. We demonstrate that the use of the ABF advanced sampling method yields values systematically closer to experiment, with a lower uncertainty in the calculation, than thermodynamic integration and restrained sampling techniques commonly used in the literature. We hypothesize that this is due to the algorithm’s ability to drive multiple binding and dissociation events through applied bias forces. These results show the promise of unrestrained advanced sampling methods along well-defined reaction coordinates to the calculation of binding energies in host–guest systems.

Michael D. LaCount and Francois Gygi

*J. Chem. Phys.* 151, 164101 (2019), *Editor's Pick*. DOI: 10.1063/1.5124957

We present an ensemble of 16 independent first-principles molecular dynamics simulations of water performed using the Strongly Constrained and Appropriately Normed (SCAN) meta-generalized gradient approximation exchange-correlation functional. These simulations were used to compute the structural and electronic properties of liquid water, as well as polarizabilities, Raman and infrared spectra. Overall, we find that the SCAN functional used at a simulation temperature of 330 K provides an accurate description of the structural and electronic properties of water while incurring a moderate computational cost. The availability of an ensemble of independent simulations provides a quantitative estimate of the uncertainty in computed structural and electronic properties. Results are also compared with a similar dataset generated using the Perdew, Burke, and Ernzerhof exchange-correlation functional at a temperature of 400 K. All simulation data and trajectories are available at http://quantum-simulation.org.

Jiale Shi, Hythem Sidky, and Jonathan K. Whitmer

*Soft Matter* 15, 8219-8226 (2019). DOI: 10.1039/C9SM01395D

Bent-shaped liquid crystals have attracted significant attention recently due to their novel mesostructure and the intriguing behavior of their elastic constants, which are strongly anisotropic and have an unusual temperature dependence. Though theories explain the onset of the twist-bend nematic phase (NTB) through spontaneous symmetry breaking concomitant with transition to a negative bend (K3) elastic constant, this has not been observed as yet in experiments. There, the small bend elastic constant has a strongly non-monotonic temperature dependence, which first increases after crossing the isotropic (I)–nematic (N) transition, then dips near the nematic (N)-twist-bend (NTB) transition before it increases again as the transition is crossed. The molecular mechanisms responsible for this exotic behavior are unclear. Here, we utilize density of states algorithms in Monte Carlo simulation applied to a variant of the Lebwohl–Lasher model which includes bent-shaped-like interactions to analyze the mechanism behind elastic response in this novel mesostructure and understand the temperature dependence of its Frank–Oseen elastic constants.

Abhijit Hazarika, Igor Fedin, Liang Hong, Jinglong Guo, Vishwas Srivastava, Wooje Cho, Igor Coropceanu, Joshua Portner, Benjamin T. Diroll, John P. Philbin, Eran Rabani, Robert Klie, and Dmitri V. Talapin

*J. Amer. Chem. Soc.* 141, 13487-13496 (2019). DOI: 10.1021/jacs.9b04866

In contrast to molecular systems, which are defined with atomic precision, nanomaterials generally show some heterogeneity in size, shape, and composition. The sample inhomogeneity translates into a distribution of energy levels, band gaps, work functions, and other characteristics, which detrimentally affect practically every property of functional nanomaterials. We discuss a novel synthetic strategy, colloidal atomic layer deposition (c-ALD) with stationary reactant phases, which largely circumvents the limitations of traditional colloidal syntheses of nano-heterostructures with atomic precision. This approach allows for significant reduction of inhomogeneity in nanomaterials in complex nanostructures without compromising their structural perfection and enables the synthesis of epitaxial nano-heterostructures of unprecedented complexity. The improved synthetic control ultimately enables bandgap and strain engineering in colloidal nanomaterials with close to atomic accuracy. To demonstrate the power of the new c-ALD method, we synthesize a library of complex II–VI semiconductor nanoplatelet heterostructures. By combining spectroscopic and computational studies, we elucidate the subtle interplay between quantum confinement and strain effects on the optical properties of semiconductor nanostructures.

Trung Dac Nguyen, Honghao Li, Debarshee Bagchi, Francisco J. Solis, and Monica Olvera de la Cruz

*Comp. Phys. Commun.* 241, 80-91 (2019). DOI: 10.1016/j.cpc.2019.03.006

We implement and extend three methods for incorporating surface polarization effects into coarse-grained Molecular Dynamics (MD) simulations for LAMMPS: the boundary element method using the generalized minimum residual (GMRES) solver for solving the Poisson equation, the induced charge computation method using the successive over relaxations scheme, and the direct optimization of an energy functional of induced charge density. Our implementations are validated against the analytical results for several cases and against published results on electrolyte and polyelectrolyte adsorption on charged interfaces. We have also employed the mentioned methods to examine the effects of surface polarization due to dielectric mismatch on the conformational behavior and relaxation times of the Rouse modes of a highly charged polyelectrolyte with explicit multivalent counterions in a spherical droplet surrounded by a low dielectric medium. We have found that dielectric mismatch coupled with spatial confinement results in noticeable changes in the radius of gyration of the collapsed conformation as well as the work required to compress and stretch the chain. Dielectric confinement is therefore believed to play an important role in the packaging of polyelectrolytes, which was generally overlooked in previous studies. Finally, the parallel performance of the implemented methods for practical simulations is characterized for a model electrolyte system. For polarizable interfaces composed of more than 1000 mesh points, the induced charge computation method is found to be the most efficient.

Andrey Sokolov, Ali Mozaffari, Rui Zhang, Juan J. de Pablo, and Alexey Snezhko

*Phys. Rev. X* 9, 031014 (2019). DOI: 10.1103/PhysRevX.9.031014

Living liquid crystals, a realization of active nematics where a lyotropic liquid crystal is combined with active bacteria, exhibit a plethora of out-of-equilibrium phenomena that range from active turbulence and dynamic spatiotemporal patterns to the creation and annihilation of motile topological defects. Experiments and hydrodynamic simulations are used here to report on the emergence of bend stripes, which arise as spontaneous undulations of the director field in circularly aligned lyotropic liquid crystals doped with bacteria. The interplay between bacterial-induced hydrodynamic flows and elastic forces in the material induces remarkable deformation patterns consisting of branched, radially elongated bands of a high curvature of the director field. The average number of such branches increases with the distance from the center of the circular alignment, leading to the formation of a radial tree of bands that is reminiscent of a snowflake structure. Hydrodynamic simulations, which are in agreement with the experiments, are used to explain the origin of such structures and to provide additional insights into regimes that are beyond the limit of experimental measurements. In particular, it is found that when activity is switched off in the early stages of pattern formation, a pronounced decay of bend-distortion energy ensues, with little change of the splay energy, serving to confirm that the bend stripes are an outcome of activity-driven bend-instability phenomena. Taken together, experiments and simulations demonstrate a system in which strain and geometry can be combined to dynamically manipulate pattern formation in active matter, paving the way to a deeper understanding and finer control of active colloidal systems.

Huihuo Zheng, Marco Govoni, and Giulia Galli

*Phys. Rev. Mat.* 3, 073803 (2019). DOI: 10.1103/PhysRevMaterials.3.073803

We derive a dielectric-dependent hybrid functional which accurately describes the electronic properties of heterogeneous interfaces and surfaces, as well as those of three- and two-dimensional bulk solids. The functional, which does not contain any adjustable parameter, is a generalization of self-consistent hybrid functionals introduced for homogeneous solids, where the screened Coulomb interaction is defined using a spatially varying, local dielectric function. The latter is determined self-consistently using density functional calculations in finite electric fields. We present results for the band gaps and dielectric constants of 3D and 2D bulk materials, and band offsets for interfaces, showing an accuracy comparable to that of GW calculations.

Marcello Puligheddu, Yi Xia, Maria K. Y. Chan, and Giulia Galli

*Phys. Rev. Mat.* 3, 085401 (2019). DOI: 10.1103/PhysRevMaterials.3.085401.

The predictive modeling of lattice thermal conductivity is of fundamental importance for the understanding and design of materials for a wide range of applications. Two major approaches, namely molecular dynamics (MD) simulations and calculations solving approximately the Boltzmann transport equation (BTE), have been developed to compute the lattice thermal conductivity. We present a detailed direct comparison of these two approaches, using as prototypical cases MgO and PbTe. The comparison, carried out using empirical potentials, takes into account the effects of fourth order phonon scattering, temperature-dependent phonon frequencies (phonon renormalization), and investigates the effects of quantum vs. classical statistics. We clarify that equipartition, as opposed to Maxwell Boltzmann, govern the statistics of phonons in MD simulations. We find that lattice thermal conductivity values from MD and BTE show an apparent, satisfactory agreement; however such an agreement is the result of error cancellations. We also show that the primary effect of statistics on thermal conductivity is via the scattering rate dependence on phonon populations.

Joshua Lequieu, Andrés Córdoba, Joshua Moller, and Juan J. de Pablo

*J. Chem. Phys.* 150, 215102 (2019). DOI: 10.1063/1.5092976

A central question in epigenetics is how histone modifications influence the 3D structure of eukaryotic genomes and, ultimately, how this 3D structure is manifested in gene expression. The wide range of length scales that influence the 3D genome structure presents important challenges; epigenetic modifications to histones occur on scales of angstroms, yet the resulting effects of these modifications on genome structure can span micrometers. There is a scarcity of computational tools capable of providing a mechanistic picture of how molecular information from individual histones is propagated up to large regions of the genome. In this work, a new molecular model of chromatin is presented that provides such a picture. This new model, referred to as 1CPN, is structured around a rigorous multiscale approach, whereby free energies from an established and extensively validated model of the nucleosome are mapped onto a reduced coarse-grained topology. As such, 1CPN incorporates detailed physics from the nucleosome, such as histone modifications and DNA sequence, while maintaining the computational efficiency that is required to permit kilobase-scale simulations of genomic DNA. The 1CPN model reproduces the free energies and dynamics of both single nucleosomes and short chromatin fibers, and it is shown to be compatible with recently developed models of the linker histone. It is applied here to examine the effects of the linker DNA on the free energies of chromatin assembly and to demonstrate that these free energies are strongly dependent on the linker DNA length, pitch, and even DNA sequence. The 1CPN model is implemented in the LAMMPS simulation package and is distributed freely for public use.

Ngoc Linh Nguyen, He Ma, Marco Govoni, Francois Gygi, and Giulia Galli

*Phys. Rev. Lett.* 122, 237402 (2019). DOI: 10.1103/PhysRevLett.122.237402

We present a method to compute optical spectra and exciton binding energies of molecules and solids based on the solution of the Bethe-Salpeter equation (BSE) and the calculation of the screened Coulomb interaction in finite field. The method does not require the explicit evaluation of dielectric matrices nor of virtual electronic states, and can be easily applied without resorting to the random phase approximation. In addition it utilizes localized orbitals obtained from Bloch states using bisection techniques, thus greatly reducing the complexity of the calculation and enabling the efficient use of hybrid functionals to obtain single particle wavefunctions. We report exciton binding energies of several molecules and absorption spectra of condensed systems of unprecedented size, including water and ice samples with hundreds of atoms.

Jarol Molina, Juan J. de Pablo, and Juan P. Hernández-Ortiz

*Phys. Chem. Chem. Phys.* 21, 9362 (2019). DOI: 10.1039/c9cp00598f

The design of polymeric membranes for proton or ionic exchange highly depends on the fundamentalunderstanding of the physical and molecular mechanisms that control the formation of the conductionchannels. There is an inherent relation between the dynamical structure of the polymeric membraneand the electrostatic forces that drive membrane segregation and proton transport. Here, we used amulti-scale computational approach to analyze the morphology of sulfonated poly(ether ether ketone)membranes at the mesoscale. A self-consistent description of the electrostatic phenomenon wasadopted, where discrete polymer chains and a continuum proton field were embedded in a continuumfluid. Brownian dynamics was used for the evolution of the suspended polymer molecules, while aconvection–diffusion transport equation, including the Nernst–Planck diffusion mechanism, accountedfor the dynamics of the proton concentration field. We varied the polymer concentration, the degree ofsulfonation and the level of confinement to find relationships between membrane structure and protonconduction. Our results indicate that the reduced mobility of polymer chains, at concentrations aboveoverlap, and a moderate degree of sulfonation –i.e., 30% – are essential elements for membrane segregationand proton domain connectivity. These conditions also ensure that the membrane structure is not affectedby size or by potential gradients. Importantly, our analysis shows that membrane conductivity and current arelinearly dependent on polymer concentration and quadratically dependent on the degree of sulfonation. Wefound that the optimal polymeric membrane design requires a polymer concentration above overlap and adegree of sulfonation around 50%. These conditions promote a dynamical membrane morphology with aconstant density of proton channels. Our results and measurements agree with previous experimental works,thereby validating our model and observations.

Krishnendu Ghosh, He Ma, Vikram Gavini, and Giulia Galli

*Phys. Rev. Mat.* 3, 043801 (2019). DOI: 10.1103/PhysRevMaterials.3.043801

The interaction between electronic and nuclear spins in the presence of external magnetic fields can be described by a spin Hamiltonian, with parameters obtained from first principles, electronic structure calculations. We describe an approach to compute these parameters, applicable to both molecules and solids, which is based on Density Functional Theory (DFT) and real-space, all-electron calculations using finite elements (FE). We report results for hyperfine tensors, zero field splitting tensors (spin-spin component) and nuclear quadrupole tensors of a series of molecules and of the nitrogen-vacancy center in diamond. We compare our results with those of calculations using Gaussian orbitals and plane-wave basis sets, and we discuss their numerical accuracy. We show that calculations based on FE can be systematically converged with respect to the basis set, thus allowing one to establish reference values for the spin Hamiltonian parameters, at a given level of DFT.

Nicholas E. Jackson, Michael A. Webb, and Juan J. de Pablo

*Curr. Opin. Chem. Eng.* 23, 106-114, (2019). DOI: 10.1016/j.coche.2019.03.005

The multiscale design of soft materials requires an ensemble of computational techniques spanning quantum-chemistry to molecular dynamics to continuum modeling. The recent emergence of machine-learning (ML) and modern optimization algorithms has accelerated material property prediction, as well as stimulated the development of hybrid ML/molecular modeling methodologies capable of providing physical insights unobtainable from purely physics-based modeling and intuition. Such hybrid techniques also have important ramifications for the ML-enhanced interpretation of results from simulations and experiments alike. Leveraging ML techniques for the design of chemical or morphological structures based on a target property or functionality represents an exciting goal for the general area of soft materials, including polymers, liquid crystals, colloids, or biomolecules, to name a few representative classes of systems. Here, we provide a perspective on recent work using ML techniques of relevance for the multiscale design of soft materials and outline potential future directions of interest to the soft materials community.

Ye Zhou, Bohdan Senyuk, Rui Zhang, Ivan I. Smalyukh, and Juan J. de Pablo

*Nat. Commun.* 10, 1000, (2019). DOI: 10.1038/s41467-019-08645-9

The defect structure associated with a colloid in a nematic liquid crystal is dictated bymolecular orientation at the colloid surface. Perpendicular or parallel orientations to thesurface lead to dipole-like or quadrupole-like defect structures. However, the so-called elastichexadecapole discovered recently, has been assumed to result from a conic anchoringcondition. In order to understand it at a fundamental level, a model for this anchoring isintroduced here in the context of a Landau-de Gennes free energy functional. We investigatethe evolution of defect configurations, as well as colloidal interactions, by tuning the preferredtilt angle (θe). The model predicts an elastic dipole whose stability decreases asθeincreases,along with a dipole-hexadecapole transformation, which are confirmed by our experimentalobservations. Taken together, our results suggest that previously unanticipated avenues mayexist for design of self-assembled structures via control of tilt angle.

Igor Coropceanu, Michael A. Boles, and Dmitri V. Talapin

*J. Amer. Chem. Soc.* 141, 5728 (2019). DOI: 10.1021/jacs.8b12539

The self-assembly of two sizes of spherical nanocrystals has revealed a surprisingly diverse library of structures. To date, at least 15 distinct binary nanocrystal superlattice (BNSL) structures have been identified. The stability of these binary phases cannot be fully explained using the traditional conceptual framework treating the assembly process as entropy-driven crystallization of rigid spherical particles. Such deviation from hard sphere behavior may be explained by the soft and deformable layer of ligands that envelops the nanocrystals, which contributes significantly to the overall size and shape of assembling particles. In this work, we describe a set of experiments designed to elucidate the role of the ligand corona in shaping the thermodynamics and kinetics of BNSL assembly. Using hydrocarbon-capped Au and PbS nanocrystals as a model binary system, we systematically tuned the core radius (R) and ligand chain length (L) of particles and subsequently assembled them into binary superlattices. The resulting database of binary structures enabled a detailed analysis of the role of effective nanocrystal size ratio, as well as softness expressed as L/R, in directing the assembly of binary structures. This catalog of superlattices allowed us to not only study the frequency of different phases but to also systematically measure the geometric parameters of the BNSLs. This analysis allowed us to evaluate new theoretical models treating the cocrystallization of deformable spheres and to formulate new hypotheses about the factors affecting the nucleation and growth of the binary superlattices. Among other insights, our results suggest that the relative abundance of the binary phases observed may be explained not only by considerations of thermodynamic stability, but also by a postulated preordering of the binary fluid into local structures with icosahedral or polytetrahedral symmetry prior to nucleation.

Nicholas E. Jackson, Alec S. Bowen, Lucas W. Antony, Michael A. Webb, Venkatram Vishwanath, and Juan J. de Pablo

*Sci. Adv.* 5, eaav1190 (2019). DOI: 10.1126/sciadv.aav1190

Computational studies aimed at understanding conformationally dependent electronic structure in soft materials require a combination of classical and quantum-mechanical simulations, for which the sampling of conformational space can be particularly demanding. Coarse-grained (CG) models provide a means of accessing relevant time scales, but CG configurations must be back-mapped into atomistic representations to perform quantum-chemical calculations, which is computationally intensive and inconsistent with the spatial resolution of the CG models. A machine learning approach, denoted as artificial neural network electronic coarse graining (ANN-ECG), is presented here in which the conformationally dependent electronic structure of a molecule is mapped directly to CG pseudo-atom configurations. By averaging over decimated degrees of freedom, ANN-ECG accelerates simulations by eliminating backmapping and repeated quantum-chemical calculations. The approach is accurate, consistent with the CG spatial resolution, and can be used to identify computationally optimal CG resolutions.

Yamil J. Colon, Ashley Z. Guo, Lucas W. Antony, Kyle Q. Hoffmann, and Juan J. de Pablo

*J. Chem. Phys* 150, 104502 (2019). DOI: 10.1063/1.5063588

Metal-organic frameworks (MOFs) represent an important class of materials. Careful selection of building blocks allows for tailoring of the properties of the resulting framework. The self-assembly process, however, is not understood, and without detailed knowledge of the underlying molecular mechanism, it is difficult to anticipate whether a particular design can be realized, or whether the material adopts a metastable, kinetically arrested state. We present a detailed examination of early-stage self-assembly pathways of the MOF-5. Enhanced sampling techniques are used to model a self-assembly in an explicit solvent (dimethylformamide, DMF). We identify several free energy barriers encountered during the assembly of the final MOF, which arise from structural rearrangements preceding MOF formation and from disrupted MOF-solvent interactions as formation proceeds. In all cases considered here, MOFs exhibit favorable entropic gains during the assembly. More generally, the strategy presented provides a step toward the experimental design characterizing the formation of ordered frameworks and possible sources of polymorphism.

Marco Govoni, Milson Munakami, Aditya Tanikanti, Jonathan H. Skone, Hakizumwami B. Runesha, Federico Giberti, Juan J. de Pablo, and Giulia Galli

*Sci. Data* 6, 190002 (2019). DOI: 10.1038/sdata.2019.2

We propose a strategy and present a simple tool to facilitate scientific data reproducibility by making available, in a distributed manner, all data and procedures presented in scientific papers, together with metadata to render them searchable and discoverable. In particular, we describe a graphical user interface (GUI), Qresp, to curate papers (i.e. generate metadata) and to explore curated papers and automatically access the data presented in scientific publications.

He Ma, Marco Govoni, Francois Gygi, and Giulia Galli

*J. Chem. Theory Comput.* 15, 154-164 (2019). DOI: 10.1021/acs.jctc.8b00864

We describe a finite-field approach to compute density response functions, which allows for efficient G0W0 and G0W0Γ0 calculations beyond the random phase approximation. The method is easily applicable to density functional calculations performed with hybrid functionals. We present results for the electronic properties of molecules and solids and we discuss a general scheme to overcome slow convergence of quasiparticle energies obtained from G0W0Γ0 calculations, as a function of the basis set used to represent the dielectric matrix.

Vikramjit S. Rathee, Hythem Sidky, Benjamin J. Sikora, and Jonathan K. Whitmer

*Polymers* 11, 183 (2019). DOI: 10.3390/polym11010183

The titration behavior of weak polyelectrolytes is of high importance, due to their uses in new technologies including nanofiltration and drug delivery applications. A comprehensive picture of polyelectrolyte titration under relevant conditions is currently lacking, due to the complexity of systems involved in the process. One must contend with the inherent structural and solvation properties of the polymer, the presence of counterions, and local chemical equilibria enforced by background salt concentration and solution acidity. Moreover, for these cases, the systems of interest have locally high concentrations of monomers, induced by polymer connectivity or confinement, and thus deviate from ideal titration behavior. This work furthers knowledge in this limit utilizing hybrid Monte Carlo–Molecular Dynamics simulations to investigate the influence of salt concentration, pKa, pH, and counterion valence in determining the coil-to-globule transition of poorly solvated weak polyelectrolytes. We characterize this transition at a range of experimentally relevant salt concentrations and explicitly examine the role multivalent salts play in determining polyelectrolyte ionization behavior and conformations. These simulations serve as an essential starting point in understanding the complexation between weak polyelectrolytes and ion rejection of self-assembled copolymer membranes

Michael Webb, Jean-Yves Delannoy, and Juan J. de Pablo

*J. Chem. Theory Comput.* 15, 1199-1208 (2019). DOI: 10.1021/acs.jctc.8b00920

A novel methodology is introduced here to generate coarse-grained (CG) representations of molecular models for simulations. The proposed strategy relies on basic graph-theoretic principles, and is referred to as graph-based coarse-graining (GBCG). It treats a given system as a molecular graph, and derives a corresponding CG representation by using edge contractions to combine nodes in the graph, which correspond to atoms in the molecule, into CG sites. A key element of this methodology is that the nodes are combined according to well defined protocols that rank-order nodes based on the underlying chemical connectivity. By iteratively performing these operations, successively coarser representations of the original atomic system can be produced to yield a systematic set of CG mappings with hierarchical resolution in an automated fashion. These capabilities are demonstrated in the context of several test systems, including toluene, pentadecane, a polysaccharide dimer, and a rhodopsin protein. In these examples, GBCG yields multiple, intuitive structures that naturally preserve the chemical topology of the system. Importantly, these representations are rendered from algorithmic implementation rather than an arbitrary ansatz, which, until now, has been the conventional approach for defining CG mapping schemes. Overall, the results presented here indicate that GBCG is efficient, robust, and unambiguous in its application, making it a valuable tool for future CG modeling.

T.J. Smart, F. Wu, M. Govoni, and Y. Ping

*Phys. Rev. Mat.* 2, 124002 (2018). DOI: 10.1103/PhysRevMaterials.2.124002

Defects in 2D materials are becoming prominent candidates for quantum emitters and scalable optoelectronic applications. However, several physical properties that characterize their behavior, such as charged defect ionization energies, are difficult to simulate with conventional first-principles methods, mainly because of the weak and anisotropic dielectric screening caused by the reduced dimensionality. We establish fundamental principles for accurate and efficient calculations of charged defect ionization energies and electronic structure in ultrathin 2D materials. We propose to use the vacuum level as the reference for defect charge transition levels (CTLs) because it gives robust results insensitive to the level of theory, unlike commonly used band edge positions. Furthermore, we determine the fraction of Fock exchange in hybrid functionals for accurate band gaps and band edge positions of 2D materials by enforcing the generalized Koopmans’ condition of localized defect states. We found the obtained fractions of Fock exchange vary significantly from 0.2 for bulk h-BN to 0.4 for monolayer h-BN, whose band gaps are also in good agreement with experimental results and calculated GW results. The combination of these methods allows for reliable and efficient prediction of defect ionization energies (difference between CTLs and band edge positions). We motivate and generalize these findings with several examples including different defects in monolayer to few-layer hexagonal boron nitride (h-BN), monolayer MoS2 and graphane. Finally, we show that increasing the number of layers ofh-BN systematically lowers defect ionization energies, mainly through CTLs shifting towards vacuum, with conduction band minima kept almost unchanged.

Carl S. Adorf, Vyas Ramasubramani, Joshua A. Anderson, and Sharon C. Glotzer

*Comput. Sci. Eng.* 21, 66-79 (2018). DOI: 10.1109/mcse.2018.2882355

A critical challenge in scientific computing is balancing developing high quality software with the need for immediate scientific progress. We present a flexible approach that emphasizes writing specialized code that is refactored only when future immediate scientific goals demand it. Our lazy refactoringtechnique, which calls for code with clearly defined interfaces and sharply delimited scopes to maximize reuse and integrability, helps reduce total development time and accelerates the production of scientific results. We offer guidelines for how to implement such code, as well as criteria to aid in the evaluation of existing tools. To demonstrate their application, we showcase the development progression of tools for particle simulations originating from the Glotzer Group at the University of Michigan. We emphasize the evolution of these tools into a loosely integrated software stack of highly reusable software that can be maintained to ensure the long-term stability of established research workflows.

James M Hodges, Shiqiang Hao, Jann A Grovogui, Xiaomi Zhang, Trevor Paul Bailey, Xiang Li, Zhehong Gan, Yan-Yan Hu, Ctirad Uher, Vinayak P. Dravid, Christopher Wolverton, and Mercouri G. Kanatzidis

*J. Am. Chem. Soc.* 140, 18115-18123 (2018). DOI: 10.1021/jacs.8b11050

Thermoelectric generators can convert heat directly into usable electric power but suffer from low efficiencies and high costs, which have hindered wide-scale applications. Accordingly, an important goal in the field of thermoelectricity is to develop new high performance materials that are composed of more earth-abundant elements. The best systems for midtemperature power generation rely on heavily doped PbTe, but the Te in these materials is scarce in the Earth’s crust. PbSe is emerging as a less expensive alternative to PbTe, although it displays inferior performance due to a considerably smaller power factor S2σ, where S is the Seebeck coefficient and σ is electrical conductivity. Here, we present a new p-type PbSe system, Pb0.98Na0.02Se–x%HgSe, which yields a very high power factor of ∼20 μW·cm–1·K–2 at 963 K when x = 2, a 15% improvement over the best performing PbSe–x%MSe materials. The enhancement is attributed to a combination of high carrier mobility and the early onset of band convergence in the Hg-alloyed samples (∼550 K), which results in a significant increase in the Seebeck coefficient. Interestingly, we find that the Hg2+ cations sit at an off-centered position within the PbSe lattice, and we dub the displaced Hg atoms “discordant”. DFT calculations indicate that this feature plays a role in lowering thermal conductivity, and we believe that this insight may inspire new design criteria for engineering high performance thermoelectric materials. The high power factor combined with a decrease in thermal conductivity gives a high figure of merit ZT of 1.7 at 970 K, the highest value reported for p-type PbSe to date.

Vikramjit S. Rathee, Hythem Sidky, Benjamin J. Sikora, and Jonathan K. Whitmer

*J. Am. Chem. Soc.* 140, 15319–15328 (2018). DOI: 10.1021/jacs.8b08649

Polyelectrolytes may be classified into two primary categories (strong and weak) depending on how their charge state responds to the local environment. Both of these find use in many applications, including drug delivery, gene therapy, layer-by-layer films, and fabrication of ion filtration membranes. The mechanism of polyelectrolyte complexation is, however, still not completely understood, though experimental investigations suggest that entropy gain due to release of counterions is the key driving force for strong polyelectrolyte complexation. Here we perform a comprehensive thermodynamic investigation through coarse-grained molecular simulations permitting us to calculate the free energy of complex formation. Importantly, our expanded-ensemble methods permit the explicit separation of energetic and entropic contributions to the free energy. Our investigations indicate that entropic contributions indeed dominate the free energy of complex formation for strong polyelectrolytes, but are less important than energetic contributions when weak electrostatic coupling or weak polyelectrolytes are present. Our results provide a new view of the free energy of polyelectrolyte complex formation driven by polymer association, which should also arise in systems with large charge spacings or bulky counterions, both of which act to weaken ion–polymer binding.

Ryan McAvoy, Marco Govoni, and Giulia Galli

*J. Chem. Theory Comput.* 14, 6269–6275 (2018). DOI: 10.1021/acs.jctc.8b00728

We report first principles calculations of electronic gaps, lifetimes and photoelectron spectra of a series of molecules, performed by efficiently combining the computation of electron-electron and electron-phonon self-energies. The dielectric matrix is represented in terms of dielectric eigenpotentials, utilized for both the calculation of $G_0W_0$ quasi-particle energies and the diagonalization of the dynamical matrix; virtual electronic states are never explicitly computed and all self-energies are evaluated over the full frequency spectrum. Our formulation enables electronic structure calculations at the many-body perturbation theory level, inclusive of electron-phonon coupling, for systems with hundreds of electrons.

Daniel M. Kroupa, Gregory F. Pach, Márton Vörös, Federico Giberti, Boris D Chernomordik, Ryan W. Crisp, Arthur J Nozik, Justin C. Johnson, Rohan Singh, Victor I. Klimov, Giulia Galli, and Matthew C. Beard

*ACS Nano* 12, 10084–10094 (2018). DOI: 10.1021/acsnano.8b04850

Generating multiple excitons by a single high-energy photon is a promising third generation solar energy conversion strategy. We demonstrate that multiple exciton generation (MEG) in PbS|CdS Janus-like hetero-nanostructures is enhanced over that of single-component and core/shell nanocrystal architectures, with an onset close to two times the PbS band gap. We attribute the enhanced MEG to the asymmetric nature of the hetero-nanostructure that results in an increase in the effective Coulomb interaction that drives MEG and a reduction of the competing hot exciton cooling rate. Slowed cooling occurs through effective trapping of hot- holes by a manifold of valence band interfacial states having both PbS and CdS character, as evidenced by photoluminescence studies and by ab initio calculations. Using transient photocurrent spectroscopy, we find that the MEG characteristics of the individual nanostructures are maintained in conductive arrays and demonstrate that these quasi-spherical PbS|CdS nanocrystals can be incorporated as the main absorber layer in functional solid-state solar cell architectures. Finally, based upon our analysis, we provide design rules for the next generation of engineered nanocrystals to further improve the MEG characteristics.

Matteo Gerosa, Francois Gygi, Marco Govoni, and Giulia Galli

*Nat. Mater.* 17, 1122-1127 (2018). DOI: 10.1038/s41563-018-0192-4

Computational screening of materials for solar to fuel conversion technologies have mostly focused on bulk properties, thus neglecting the structure and chemistry of surfaces and interfaces with water. We report a finite temperature study of WO3, a promising anode for photo-electrochemical cells, carried out using first principles molecular dynamics simulations coupled with many body perturbation theory. We identified three major factors determining the chemical reactivity of the material in- terfaced with water: the presence of surface defects, the dynamics of excess charge at the surface, and finite temperature fluctuations of the surface electronic orbitals. These general descriptors are essential for the understanding and prediction of opti- mal oxide photo-absorbers for water oxidation.

James M Hodges, Yi Xia, Christos D. Malliakas, Grant C. B. Alexander, Maria K. Y. Chan, and Mercouri G. Kanatzidis

*Chem. Mater.* 30, 7245-7254 (2018). DOI: 10.1021/acs.chemmater.8b03306

Metal chalcogenides underpin a wide variety of energy-related applications and are ideal systems for probing lattice dynamics and fundamental transport phenomena. Here, we describe the synthesis and transport properties of CsAg5TeS2 and its solid solution CsAg5Te3-xSx (x = 1-2), new semiconductors with tunable bandgaps ranging from 0.17 to 0.30 eV. CsAg5TeS2 has a fully ordered two-dimensional structure that includes a group of Ag atoms in a heteroleptic, tetrahedral coordination geometry (AgTe2S2). Single crystal X-ray diffraction indicates that the compounds crystallize in the tetragonal space group P4/mmm, while pair-distribution function (PDF) analysis reveals off-centering at the heteroleptic Ag sites, signifying a lower symmetry I4/mcm space group. The underlying disorder acts as a phonon-blocking mechanism that helps facilitate ultralow lattice thermal conductivity below 0.40 W⋅m-1⋅K-1 at ~300 K and highlights the importance of local disorder in thermal transport. Density functional theory (DFT) provides additional insight into the electronic and thermal properties of the materials, which are good candidates for p-type thermoelectrics.

Vyas Ramasubramani, Carl S. Adorf, Paul M. Dodd, Bradley D. Dice, Sharon C. Glotzer

*Proc. 17th Python Sci. Conf. * 152–159 (2018). DOI: 10.25080/Majora-4af1f417-016

Computational research requires versatile data and workflow man- agement tools that can easily adapt to the highly dynamic requirements of scientific investigations. Many existing tools require strict adherence to a par- ticular usage pattern, so researchers often use less robust ad hoc solutions that they find easier to adopt. The resulting data fragmentation and methodological incompatibilities significantly impede research. Our talk showcases signac, an open-source Python framework that offers highly modular and scalable solutions for this problem. Named for the Pointillist painter Paul Signac, the framework’s powerful workflow management tools enable users to construct and automate workflows that transition seamlessly from laptops to HPC clusters. Crucially, the underlying data model is completely independent of the workflow. The flexible, serverless, and schema-free signac database can be introduced into other workflows with essentially no overhead and no recourse to the signac workflow model. Additionally, the data model’s simplicity makes it easy to parse the underlying data without using signac at all. This modularity and simplicity eliminates significant barriers for consistent data management across projects, facilitating improved provenance management and data sharing with minimal overhead.

Yi Xia

*App. Phys. Lett.* 113, 073901 (2018). DOI: 10.1063/1.5040887

We perform a first-principles study of lattice thermal transport in PbTe by explicitly considering anharmonicity up to 4th order. To determine the temperature-dependent lattice constant of PbTe beyond quasiharmonic approximation, we introduce a simple yet effective scheme to account for anharmonic phonon renormalization at finite temperature. Moreover, we explicitly compute mode-resolved phonon lifetimes by including both three- and four-phonon scatterings. We find that (1) anharmonic phonon renormalization leads to strong vibrational frequency shifts which improve the agreement between simulated and experimental lattice constants; (2) these frequency shifts lead to a significant increase in lattice thermal conductivity (κl) because of reduced phonon scattering phase space; and (3) four-phonon scatterings are responsible for severe reduction in κl on top of three-phonon scatterings, making κl consistent with experiments. Our results suggest that the predicted κl and its temperature dependence without considering thermal expansion, anharmonic phonon renormalization and four-phonon scatterings could accidentally agree with experiments due to error cancellation. Our study not only deepens the understanding of lattice thermal transport in PbTe but also exemplifies a widely applicable approach to investigate lattice dynamics and thermal transport properties from first-principles calculations including high-order anharmonicity.

Malte Lange and Timothy Berkelbach

*J. Chem. Theory Comput.* 14, 4224-4236 (2018). DOI: 10.1021/acs.jctc.8b00455

We discuss the analytic and diagrammatic structure of ionization potential (IP) and electron affinity (EA) equation-of-motion coupled-cluster (EOM-CC) theory, in order to put it on equal footing with the prevalent GW approximation. The comparison is most straightforward for the time-ordered one-particle Green’s function, and we show that the Green’s function calculated by EOM-CC with single and double excitations (EOM-CCSD) includes fewer ring diagrams at higher order than does the GW approximation, due to the former’s unbalanced treatment of time-ordering. However, the EOM-CCSD Green’s function contains a large number of vertex corrections, including ladder diagrams, mixed ring-ladder diagrams, and exchange diagrams. By including triple excitations, the EOM-CCSDT Green’s function includes all diagrams contained in the GW approximation, along with many high-order vertex corrections. In the same language, we discuss a number of common approximations to the EOM-CCSD equations, many of which can be classified as elimination of diagrams. Finally, we present numerical results by calculating the principal charged excitations energies of the molecules contained in the so-called GW100 test set [ J. Chem. Theory Comput. 2015, 11, 5665−5687]. We argue that (in molecules) exchange is as important as screening, advocating for a Hartree–Fock reference and second-order exchange in the self-energy.

Ding Pan, Marco Govoni, and Giulia Galli

*J. Chem. Phys.* 149, 051101 (2018). DOI: 10.1063/1.5044636

The dielectric properties of molecules and nanostructures are usually modified in a complex manner, when assembled into a condensed phase. We propose a first-principles method to compute polarizabilities of sub-entities of solids and liquids, which accounts for multipolar interactions at all orders and is applicable to semiconductors and insulators. The method only requires the evaluation of induced fields in the condensed phase, with no need of multiple calculations for each constituent. As an example, we present results for the molecular polarizabilities of water in a wide pressure and temperature range. We found that at ambient conditions, the dipole-induced-dipole approximation is sufficiently accurate and the Clausius-Mossotti relation may be used, e.g., to obtain molecular polarizabilities from experimental refractive indexes. However with increasing pressure, this approximation becomes unreliable and in the case of ice X the Clausius-Mossotti relation is not valid.

Nicholas E. Jackson, Michael A. Webb, and Juan J. de Pablo

*J. Chem. Phys.* 149, 072326 (2018). DOI: 10.1063/1.5030531

A configurational sampling algorithm based on nested layerings of Markov chains (Layered Nested Markov Chain Monte Carlo or L-NMCMC) is presented for simulations of systems characterized by rugged free energy landscapes. The layerings are generated using a set of auxiliary potential energy surfaces. The implementation of the method is demonstrated in the context of a rugged, two-dimensional potential energy surface. The versatility of the algorithm is next demonstrated on a simple, many-body system, namely, a canonical Lennard-Jones fluid in the liquid state. In that example, different layering schemes and auxiliary potentials are used, including variable cutoff distances and excluded-volume tempering. In addition to calculating a variety of properties of the system, it is also shown that L-NMCMC, when combined with a free-energy perturbation formalism, provides a straightforward means to construct approximate free-energy surfaces at no additional computational cost using the sampling distributions of each auxiliary Markov chain. The proposed L-NMCMC scheme is general in that it could be complementary to any number of methods that rely on sampling from a target distribution or methods that exploit a hierarchy of time scales and/or length scales through decomposition of the potential energy.

Xikai Jiang, Jiyuan Li, Victor Lee, Heinrich M. Jaeger, Olle G. Heinonen, and Juan J. de Pablo

*J. Chem. Phys.* 148, 234302 (2018). DOI: 10.1063/1.5027435

We report a computational strategy to obtain the charges of individual dielectric particles from experimental observation of their interactions as a function of time. This strategy uses evolutionary optimization to minimize the difference between trajectories extracted from experiment and simulated trajectories based on many-particle force fields. The force fields include both Coulombic interactions and dielectric polarization effects that arise due to particle-particle charge mismatch and particle-environment dielectric contrast. The strategy was applied to systems of free falling charged granular particles in vacuum, where electrostatic interactions are the only driving forces that influence the particles' motion. We show that when the particles' initial positions and velocities are known, the optimizer requires only an initial and final particle configuration of a short trajectory in order to accurately infer the particles' charges; when the initial velocities are unknown and only the initial positions are given, the optimizer can learn from multiple frames along the trajectory to determine the particles' initial velocities and charges. While the results presented here offer a proof-of-concept demonstration of the proposed ideas, the proposed strategy could be extended to more complex systems of electrostatically charged granular matter.

Emilio Scalise, Vishwas Srivastava, Eric M. Janke, Dmitri Talapin, Giulia Galli, and Stefan Wipperman,

*Nat. Nanotech.* 13, 841-848 (2018). DOI: 10.1038/s41565-018-0189-9

Semiconducting nanomaterials synthesized using wet chemical techniques play an important role in emerging optoelectronic and photonic technologies. Controlling the surface chemistry of the nano building blocks and their interfaces with ligands is one of the outstanding challenges for the rational design of these systems. We present an integrated theoretical and experimental approach to characterize, at the atomistic level, buried interfaces in solids of InAs nanoparticles capped with Sn2S64– ligands. These prototypical nanocomposites are known for their promising transport properties and unusual negative photoconductivity. We found that inorganic ligands dissociate on InAs to form a surface passivation layer. A nanocomposite with unique electronic and transport properties is formed, that exhibits type II heterojunctions favourable for exciton dissociation. We identified how the matrix density, sulfur content and specific defects may be designed to attain desirable electronic and transport properties, and we explain the origin of the measured negative photoconductivity of the nanocrystalline solids.

Viktor Rozsa, Ding Pan, Federico Giberti, and Giulia Galli

*PNAS* 115, 6952-6957 (2018). DOI: 10.1073/pnas.1800123115

The phase diagram of water at extreme conditions plays a critical role in Earth and planetary science, yet remains poorly understood. Here we report a first-principles investigation of the liquid at high temperature, between 11 GPa and 20 GPa—a region where numerous controversial results have been reported over the past three decades. Our results are consistent with the recent estimates of the water melting line below 1,000 K and show that on the 1,000-K isotherm the liquid is rapidly dissociating and recombining through a bimolecular mechanism. We found that short-lived ionic species act as charge carriers, giving rise to an ionic conductivity that at 11 GPa and 20 GPa is six and seven orders of magnitude larger, respectively, than at ambient conditions. Conductivity calculations were performed entirely from first principles, with no a priori assumptions on the nature of charge carriers. Despite frequent dissociative events, we observed that hydrogen bonding persists at high pressure, up to at least 20 GPa. Our computed Raman spectra, which are in excellent agreement with experiment, show no distinctive signatures of the hydronium and hydroxide ions present in our simulations. Instead, we found that infrared spectra are sensitive probes of molecular dissociation, exhibiting a broad band below the OH stretching mode ascribable to vibrations of complex ions.

Trung Dac Nguyen, Baofu Qiao, and Monica Olvera de la Cruz

*PNAS* 115, 6578-6583 (2018). DOI: 10.1073/pnas.1806207115

Membraneless organelles are aggregates of disordered proteins that form spontaneously to promote specific cellular functions in vivo. The possibility of synthesizing membraneless organelles out of cells will therefore enable fabrication of protein-based materials with functions inherent to biological matter. Since random copolymers contain various compositions and sequences of solvophobic and solvophilic groups, they are expected to function in nonbiological media similarly to a set of disordered proteins in membraneless organelles. Interestingly, the internal environment of these organelles has been noted to behave more like an organic solvent than like water. Therefore, an adsorbed layer of random copolymers that mimics the function of disordered proteins could, in principle, protect and enhance the proteins’ enzymatic activity even in organic solvents, which are ideal when the products and/or the reactants have limited solubility in aqueous media. Here, we demonstrate via multiscale simulations that random copolymers efficiently incorporate proteins into different solvents with the potential to optimize their enzymatic activity. We investigate the key factors that govern the ability of random copolymers to encapsulate proteins, including the adsorption energy, copolymer average composition, and solvent selectivity. The adsorbed polymer chains have remarkably similar sequences, indicating that the proteins are able to select certain sequences that best reduce their exposure to the solvent. We also find that the protein surface coverage decreases when the fluctuation in the average distance between the protein adsorption sites increases. The results herein set the stage for computational design of random copolymers for stabilizing and delivering proteins across multiple media.

Alex Gaiduk, Jeffrey Gustafson, Francois Gygi, and Giulia Galli

*J. Phys. Chem. Lett.* 9, 3068-3073 (2018). DOI: 10.1021/acs.jpclett.8b01017

We carried out first-principles simulations of liquid water under ambient conditions using a dielectric-dependent hybrid functional, where the fraction of exact exchange is set equal to the inverse of the high-frequency dielectric constant of the liquid. We found excellent agreement with experiment for the oxygen–oxygen partial correlation function at the experimental equilibrium density and 311 ± 3 K. Other structural and dynamical properties, such as the diffusion coefficient, molecular dipole moments, and vibrational spectra, are also in good agreement with experiment. Our results, together with previous findings on electronic properties of the liquid with the same functional, show that the dielectric-dependent hybrid functional accurately describes both the structural and electronic properties of liquid water.

Yi Xia, James M. Hodges, Mercouri G. Kanatzidis, and Maria K. Y. Chan

*Appl. Phys. Lett. * 112, 181906 (2018). DOI: 10.1063/1.5002587

PbTe, one of the most promising thermoelectric materials, has recently demonstrated a thermoelectric figure of merit (ZT) of above 2.0 when alloyed with group II elements. The improvements are due mainly to significant reduction of lattice thermal conductivity (κl), which was in turn attributed to nanoparticle precipitates. However, a fundamental understanding of various phonon scattering mechanisms within the bulk alloy is still lacking. In this work, we apply the newly-developed density-functional-theory-based compressive sensing lattice dynamics approach to model lattice heat transport in PbTe, MTe, and Pb0.94M0.06Te (M = Mg, Ca, Sr, and Ba) and compare our results with experimental measurements, with focus on the strain effect and mass disorder scattering. We find that (1) CaTe, SrTe, and BaTe in the rock-salt structure exhibit much higher κl than PbTe, while MgTe in the same structure shows anomalously low κl; (2) lattice heat transport of PbTe is extremely sensitive to static strain induced by alloying atoms in solid solution form; (3) mass disorder scattering plays a major role in reducing κl for Mg/Ca/Sr-alloyed PbTe through strongly suppressing the lifetimes of intermediate- and high-frequency phonons, while for Ba-alloyed PbTe, precipitated nanoparticles are also important.

Emre Sevgen, Federico Giberti, Hythem Sidky, Jonathan Whitmer, Giulia Galli, Francois Gygi, and Juan de Pablo

*J. Chem. Theory Comput.* 14, 2881–2888 (2018). DOI: 10.1021/acs.jctc.8b00192

We present a seamless coupling of a suite of codes designed to perform advanced sampling simulations, with a first-principles molecular dynamics (MD) engine. As an illustrative example, we discuss results for the free energy and potential surfaces of the alanine dipeptide obtained using both local and hybrid density functionals (DFT), and we compare them with those of a widely used classical force field, Amber99sb. In our calculations, the efficiency of first-principles MD using hybrid functionals is augmented by hierarchical sampling, where hybrid free energy calculations are initiated using estimates obtained with local functionals. We find that the free energy surfaces obtained from classical and first-principles calculations differ. Compared to DFT results, the classical force field overestimates the internal energy contribution of high free energy states, and it underestimates the entropic contribution along the entire free energy profile. Using the string method, we illustrate how these differences lead to different transition pathways connecting the metastable minima of the alanine dipeptide. In larger peptides, those differences would lead to qualitatively different results for the equilibrium structure and conformation of these molecules.

Ashley Z. Guo, Emre Sevgen, Hythem Sidky, Jonathan K. Whitmer, Jeffrey A. Hubbell, and Juan J. de Pablo

*J. Chem. Phys.* 148, 134108 (2018). DOI: 10.1063/1.5020733

A machine learning assisted method is presented for molecular simulation of systems with rugged free energy landscapes. The method is general and can be combined with other advanced sampling techniques. In the particular implementation proposed here, it is illustrated in the context of an adaptive biasing force approach where, rather than relying on discrete force estimates, one can resort to a self-regularizing artificial neural network to generate continuous, estimated generalized forces. By doing so, the proposed approach addresses several shortcomings common to adaptive biasing force and other algorithms. Specifically, the neural network enables (1) smooth estimates of generalized forces in sparsely sampled regions, (2) force estimates in previously unexplored regions, and (3) continuous force estimates with which to bias the simulation, as opposed to biases generated at specific points of a discrete grid. The usefulness of the method is illustrated with three different examples, chosen to highlight the wide range of applicability of the underlying concepts. In all three cases, the new method is found to enhance considerably the underlying traditional adaptive biasing force approach. The method is also found to provide improvements over previous implementations of neural network assisted algorithms.

Bozhi Tian, Shuai Xu, John A Rogers, Stefano Cestellos-Blanco, Peidong Yang, João L Carvalho-de-Souza, Francisco Bezanilla, Jia Liu, Zhenan Bao, Martin Hjort, Yuhong Cao, Nicholas Melosh, Guglielmo Lanzani, Fabio Benfenati, Giulia Galli, Francois Gygi, Rylan Kautz, Alon A Gorodetsky, Samuel S Kim, Timothy K Lu, Polina Anikeeva, Michal Cifra, Ondrej Krivosudský, Daniel Havelka, and Yuanwen Jiang

*Phys. Biol.* 15, 031002 (2018). DOI: 10.1088/1478-3975/aa9f34

This roadmap outlines the role semiconductor-based materials play in understanding the complex biophysical dynamics at multiple length scales, as well as the design and implementation of next-generation electronic, optoelectronic, and mechanical devices for biointerfaces. The roadmap emphasizes the advantages of semiconductor building blocks in interfacing, monitoring, and manipulating the activity of biological components, and discusses the possibility of using active semiconductor–cell interfaces for discovering new signaling processes in the biological world.

William Dawson, Francois Gygi

*J. Chem. Phys.* 148, 124501 (2018). DOI: 10.1063/1.5018116

First-principles molecular dynamics (FPMD) simulations based on density functional theory are becoming increasingly popular for the description of liquids. In view of the high computational cost of these simulations, the choice of an appropriate equilibration protocol is critical. We assess two methods of estimation of equilibration times using a large dataset of first-principles molecular dynamics simulations of water. The Gelman-Rubin potential scale reduction factor [A. Gelman and D. B. Rubin, Stat. Sci. 7, 457 (1992)] and the marginal standard error rule heuristic proposed by White [Simulation 69, 323 (1997)] are evaluated on a set of 32 independent 64-molecule simulations of 58 ps each, amounting to a combined cumulative time of 1.85 ns. The availability of multiple independent simulations also allows for an estimation of the variance of averaged quantities, both within MD runs and between runs. We analyze atomic trajectories, focusing on correlations of the Kohn-Sham energy, pair correlation functions, number of hydrogen bonds, and diffusion coefficient. The observed variability across samples provides a measure of the uncertainty associated with these quantities, thus facilitating meaningful comparisons of different approximations used in the simulations. We find that the computed diffusion coefficient and average number of hydrogen bonds are affected by a significant uncertainty in spite of the large size of the dataset used. A comparison with classical simulations using the TIP4P/2005 model confirms that the variability of the diffusivity is also observed after long equilibration times. Complete atomic trajectories and simulation output files are available online for further analysis.

Hythem Sidky, Jonathan Whitmer

*J. Chem. Phys.* 148, 104111 (2018). DOI: 10.1063/1.5018708

Existing adaptive bias techniques, which seek to estimate free energies and physical properties from molecular simulations, are limited by their reliance on fixed kernels or basis sets which hinder their ability to efficiently conform to varied free energy landscapes. Further, user-specified parameters are in general non-intuitive yet significantly affect the convergence rate and accuracy of the free energy estimate. Here we propose a novel method, wherein artificial neural networks (ANNs) are used to develop an adaptive biasing potential which learns free energy landscapes. We demonstrate that this method is capable of rapidly adapting to complex free energy landscapes and is not prone to boundary or oscillation problems. The method is made robust to hyperparameters and overfitting through Bayesian regularization which penalizes network weights and auto-regulates the number of effective parameters in the network. ANN sampling represents a promising innovative approach which can resolve complex free energy landscapes in less time than conventional approaches while requiring minimal user input.

Carl S. Adorf, Paul M. Dodd, Vyas Ramasubramani, Sharon C. Glotzer

*Comp. Mater. Sci.* 146, 220-229 (2018). DOI: 10.1016/j.commatsci.2018.01.035com

Researchers in the fields of materials science, chemistry, and computational physics are regularly posed with the challenge of managing large and heterogeneous data spaces. The amount of data increases in lockstep with computational efficiency multiplied by the amount of available computational resources, which shifts the bottleneck in the scientific process from data acquisition to data processing and analysis. We present a framework designed to aid in the integration of various specialized data formats, tools and workflows. The signac framework provides all basic components required to create a well-defined and thus collectively accessible and searchable data space, simplifying data access and modification through a homogeneous data interface that is largely agnostic to the data source, i.e., computation or experiment. The framework’s data model is designed to not require absolute commitment to the presented implementation, simplifying adaption into existing data sets and workflows. This approach not only increases the efficiency with which scientific results can be produced, but also significantly lowers barriers for collaborations requiring shared data access.

Vikramjit S. Rathee, Aristotle J. Zervoudakis, Hythem Sidky, Benjamin J. Sikora and Jonathan K. Whitmer

*J. Chem. Phys.* 148, 114901 (2018). DOI: 10.1063/1.5017941

Weak polyelectrolytes are relevant for a wide range of fields; in particular, they are useful in the creation of “smart” materials for chemical separations and drug delivery. The charges on weak polyelectrolytes are dynamic and sensitive to the surroundings, causing polymer chains to adopt environmentally-responsive conformations. Currently, there exists no comprehensive picture of this behavior, particularly where polymer–polymer interactions have the potential to affect charging properties significantly. In this study, we elucidate the novel interplay between weak polyelectrolyte charging and complexation behavior through coupled molecular dynamics and Monte Carlo simulations. Specifically, we investigate a model of two equal-length and oppositely charging polymer chains in implicit salt solution represented through Debye–Hückel interactions. The charging tendency of each chain, along with the salt concentration, are varied to determine the existence and extent of cooperativity in charging and complexation. Strong cooperation in the charging of these chains is observed at large Debye lengths, corresponding to low salt concentrations, while at lower Debye lengths (higher salt concentrations), the chains behave in apparent isolation. When the electrostatic coupling is long-ranged, we find that a highly charged chain strongly promotes the charging of its partner chain, even if the environment is unfavorable for an isolated version of that partner chain. Evidence of this phenomenon is supported by a drop in the potential energy of the system, which does not occur at the lower Debye lengths where both potential energies and charge fractions converge for all partner chain charging tendencies. The discovery of this cooperation will be helpful in developing “smart” drug delivery mechanisms by allowing for better predictions for the dissociation point of delivery complexes.

Hythem Sidky, Juan J. de Pablo, and Jonathan K. Whitmer

*Phys. Rev. Lett.* 120, 107801 (2018). DOI: 10.1103/PhysRevLett.120.107801

Experiments on confined droplets of the nematic liquid crystal 5CB have questioned long- established bounds imposed on the elastic free energy of nematic systems. This elasticity, which derives from molecular alignment within nematic systems, is quantified through a set of moduli which can be difficult to measure experimentally and, in some cases, can only be probed indirectly. This is particularly true of the surface-like saddle-splay elastic term, for which the available ex- perimental data indicate values on the cusp of stability, often with large uncertainties. Here, we demonstrate that all nematic elastic moduli, including the saddle-splay elastic constant k24, may be calculated directly from atomistic molecular simulation. Importantly, results obtained through in silico measurements of the 5CB elastic properties demonstrate unambiguously that saddle-splay elasticity alone is unable to describe the observed confined morphologies.

Weiwei Chu, Jian Qin, and Juan J. de Pablo

*Macromolecules* 51, 1986–1991 (2018). DOI: 10.1021/acs.macromol.7b02508

We examine the distribution of ionic species in charged, microphase-separated diblock copolymers by relying on coarse-grained simulations of the underlying materials. The model adopted here has been particularly useful in understanding the behavior of neutral block polymer systems. In this work, it is extended to describe charged molecules. A simulation methodology is proposed in which dielectric inhomogeneities within the phase-separated systems are taken into account by solving on-the-fly Poisson’s equation for spatially varying dielectric permittivity and calculating the corresponding electric fields. The Green’s function appropriate for a periodic and ionic system is explicitly calculated and used in simulations to arrive at phase diagrams as a function of salt concentration and copolymer composition. The systems considered here are representative of lithium salt-doped diblock copolymers, which could be of potential use for solid-electrolyte batteries.

Marco Govoni and Giulia Galli

*J. Chem. Theory Comput.* 14, 1895–1909 (2018). DOI: 10.1021/acs.jctc.7b00952

The reproducibility of calculations carried out within many-body perturbation the- ory at the G0W0 level is assessed for 100 closed shell molecules and compared to that of density functional theory. We consider vertical ionization potentials (VIP) and electron affinities (VEA) obtained with five different codes: BerkeleyGW, FHI-aims, TURBOMOLE, VASP, and WEST. We review the approximations and parameters that control the accuracy of G0W0 results in each code and we discuss in detail the effect of extrapolation techniques for the parameters entering the WEST code. Differences between the VIP and VEA computed with the various codes are within ∼ 60 meV and ∼ 120 meV respectively, which is up to four times larger than in the case of the best results obtained with DFT codes. Vertical ionization potentials are validated against experiment and CCSD(T) quantum chemistry results showing a mean absolute relative error of ∼ 4% for data obtained with WEST. Our analysis of the differences between localized orbitals and plane-wave implementations points out at molecules containing Cu, I, Ga, and Xe as major sources of discrepancies, which call for a re-evaluation of the pseudopotentials used for these systems in G0W0 calculations.

Alex Gaiduk, T. Anh Pham, Marco Govoni, Francesco Paesani and Giulia Galli

*Natur. Commun.* 9, 237 (2018). DOI: 10.1038/s41467-017-02673-z

Understanding redox and photochemical reactions in aqueous environments requires a precise knowledge of the ionization potential and electron affinity of liquid water. The former has been measured, but not the latter. We predict the electron affinity of liquid water and of its surface from first principles, coupling path-integral molecular dynamics with ab initio potentials, and many-body perturbation theory. Our results for the surface (0.8 eV) agree well with recent pump-probe spectroscopy measurements on amorphous ice. Those for the bulk (0.1–0.3 eV) differ from several estimates adopted in the literature, which we critically revisit. We show that the ionization potential of the bulk and surface are almost identical; instead their electron affinities differ substantially, with the conduction band edge of the surface much deeper in energy than that of the bulk. We also discuss the significant impact of nuclear quantum effects on the fundamental gap and band edges of the liquid.

Michael Quevillon and Jonathan K. Whitmer

*Materials* 11, 64 (2018). DOI: 10.3390/ma11010064

Ionic liquid crystals occupy an intriguing middle ground between room-temperature ionic liquids and mesostructured liquid crystals. Here, we examine a non-polarizable, fully atomistic model of the 1-alkyl-3-methylimidazolium nitrate family using molecular dynamics in the constant pressure–constant temperature ensemble. These materials exhibit a distinct “smectic” liquid phase, characterized by layers formed by the molecules, which separate the ionic and aliphatic moieties. In particular, we discuss the implications this layering may have for electrolyte applications.

Vikramjit S. Rathee, Benjamin J. Sikora, Hythem Sidky and Jonathan K. Whitmer

*Materials Research Express* 5, 014010 (2018). DOI: 10.1088/2053-1591/aaa049

The coil–globule transition in weak (annealed) polyelectrolytes involves a subtle balance of pH, charge strength, and solvation forces. In this work, we utilize a coarse-grained hybrid grand-canonical Monte Carlo and molecular dynamics approach to explore the swelling behavior of weak linear and star polyelectrolytes under different ionic screening conditions and pH. Importantly, we are able to quantify topology-dependent effects in charging which arise at the core of star polymers. Our results are suggestive of suppression of charging in star weak polyelectrolytes in comparison to linear weak polyelectrolytes. Furthermore, we characterize the coil–globule transition in linear and star weak polyelectrolyte through expanded ensemble density-of-states simulations which suggest a change from a first order to second order phase transition moving from linear to star polyelectrolytes. Lastly, we characterize the inhomogeneous charging across the weak star polyelectrolyte through observed shifts in ${\rm{\Delta }}{{\rm{pK}}}_{{\rm{o}}},$ and compare with experimental work. We discuss these results in relation to surfaces functionalized by weak polyelectrolyte brushes and weak polyelectrolyte-based drug delivery applications

Beatriz Seoane, Daniel R. Reid, Juan J. de Pablo, and Francesco Zamponi

*Phys. Rev. Materials* 2, 015602 (2018). DOI: 10.1103/PhysRevMaterials.2.015602

We investigate the low-temperature properties of two-dimensional Lennard-Jones glass films, prepared in silico both by liquid cooling and by physical vapor deposition. We identify deep in the solid phase a crossover temperature T∗, at which slow dynamics and enhanced heterogeneity emerge. Around T∗, localized defects become visible, leading to vibrational anomalies as compared to standard solids. We find that, on average, T∗ decreases in samples with lower inherent structure energy, suggesting that such anomalies will be suppressed in ultrastable glass films, prepared both by very slow liquid cooling and vapor deposition.

Hythem Sidky, Yamil J. Colon, Julian Helfferich, Benjamin J. Sikora, Cody Bezik, Weiwei Chu, Federico Giberti, Ashley Z. Guo, Xikai Jiang, Joshua Lequieu, Jiyuan Li, Joshua Moller, Michael J. Quevillon, Mohammad Rahimi, Hadi Ramezani-Dakhel, Vikramjit S. Rathee, Daniel R. Reid, Emre Sevgen, Vikram Thapar, Michael A. Webb, Jonathan K. Whitmer, and Juan J. de Pablo

*J. Chem. Phys.* 148, 044104 (2018). DOI: 10.1063/1.5008853

Molecular simulation has emerged as an essential tool for modern-day research, but obtaining proper results and making reliable conclusions from simulations requires adequate sampling of the system under consideration. To this end, a variety of methods exist in the literature that can enhance sampling considerably, and increasingly sophisticated, effective algorithms continue to be developed at a rapid pace. Implementation of these techniques, however, can be challenging for experts and non-experts alike. There is a clear need for software that provides rapid, reliable, and easy access to a wide range of advanced sampling methods and that facilitates implementation of new techniques as they emerge. Here we present SSAGES, a publicly available Software Suite for Advanced General Ensemble Simulations designed to interface with multiple widely used molecular dynamics simulations packages. SSAGES allows facile application of a variety of enhanced sampling techniques—including adaptive biasing force, string methods, and forward flux sampling—that extract meaningful free energy and transition path data from all-atom and coarse-grained simulations. A noteworthy feature of SSAGES is a user-friendly framework that facilitates further development and implementation of new methods and collective variables. In this work, the use of SSAGES is illustrated in the context of simple representative applications involving distinct methods and different collective variables that are available in the current release of the suite. The code may be found at: https://github.com/MICCoM/SSAGES-public.

Hadi Ramezani-Dakhel, Monirosadat Sadati, Rui Zhang, Mohammad Rahimi, Khia Kurtenbach, Benoît Roux, and Juan J. de Pablo

*ACS Cent. Sci.* 3, 1345–1349 (2017). DOI: 10.1021/acscentsci.7b00495

It is well understood that the adsorption of solutes at the interface between a bulk liquid crystal phase and an aqueous phase can lead to orientational or anchoring transitions. A different principle is introduced here, whereby a transient reorientation of a thermotropic liquid crystal is triggered by a spontaneous flux of water across the interface. A critical water flux can be generated by the addition of an electrolyte to the bulk aqueous phase, leading to a change in the solvent activity; water is then transported through the liquid crystal phase and across the interface. The magnitude of the spontaneous water flux can be controlled by the concentration and type of solutes, as well as the rate of salt addition. These results present new, previously unappreciated fundamental principles that could potentially be used for the design of materials involving transient gating mechanisms, including biological sensors, drug delivery systems, separation media, and molecular machines.

Hosung Seo, He Ma, Marco Govoni, and Giulia Galli

*Phys. Rev. Materials* 1, 075002 (2017). DOI: 0.1103/PhysRevMaterials.1.075002

The development of novel quantum bits is key to extending the scope of solid-state quantum-information science and technology. Using first-principles calculations, we propose that large metal ion–vacancy pairs are promising qubit candidates in two binary crystals: 4H-SiC and w-AlN. In particular, we found that the formation of neutral Hf- and Zr-vacancy pairs is energetically favorable in both solids; these defects have spin-triplet ground states, with electronic structures similar to those of the diamond nitrogen-vacancy center and the SiC divacancy. Interestingly, they exhibit different spin-strain coupling characteristics, and the nature of heavy metal ions may allow for easy defect implantation in desired lattice locations and ensure stability against defect diffusion. To support future experimental identification of the proposed defects, we report predictions of their optical zero-phonon line, zero-field splitting, and hyperfine parameters. The defect design concept identified here may be generalized to other binary semiconductors to facilitate the exploration of new solid-state qubit

Kristy L. Kounovsky-Shafer, Juan P. Hernandez-Ortiz, Konstantinos Potamousis, Gene Tsvid, Michael Place, Prabu Ravindran, Kyubong Jo, Shiguo Zhou, Theo Odijk, Juan J. de Pablo and David C. Schwartz

*PNAS* 114, 13400-13405 (2017). DOI: 10.1073/pnas.1711069114

Very large DNA molecules enable comprehensive analysis of complex genomes, such as human, cancer, and plants because they span across sequence repeats and complex somatic events. When physically manipulated, or analyzed as single molecules, long polyelectrolytes are problematic because of mechanical considerations that include shear-mediated breakage, dealing with the massive size of these coils, or the length of stretched DNAs using common experimental techniques and fluidic devices. Accordingly, we harness analyte “issues” as exploitable advantages by our invention and characterization of the “molecular gate,” which controls and synchronizes formation of stretched DNA molecules as DNA dumbbells within nanoslit geometries. Molecular gate geometries comprise micro- and nanoscale features designed to synergize very low ionic strength conditions in ways we show effectively create an “electrostatic bottle.” This effect greatly enhances molecular confinement within large slit geometries and supports facile, synchronized electrokinetic loading of nanoslits, even without dumbbell formation. Device geometries were considered at the molecular and continuum scales through computer simulations, which also guided our efforts to optimize design and functionalities. In addition, we show that the molecular gate may govern DNA separations because DNA molecules can be electrokinetically triggered, by varying applied voltage, to enter slits in a size-dependent manner. Lastly, mapping the Mesoplasma florum genome, via synchronized dumbbell formation, validates our nascent approach as a viable starting point for advanced development that will build an integrated system capable of large-scale genome analysis.

Marcello Puligheddu, Francois Gygi, and Giulia Galli

*Phys. Rev. Materials* 1, 060802(R) (2017). DOI: 10.1103/PhysRevMaterials.1.060802

Advances in understanding heat transport in solids were recently reported by both experiment and theory. However an efficient and predictive quantum simulation framework to investigate thermal properties of solids, with the same complexity as classical simulations, has not yet been developed. Here we present a method to compute the thermal conductivity of solids by performing ab initio molecular dynamics at close to equilibrium conditions, which only requires calculations of first-principles trajectories and atomic forces, thus avoiding direct computation of heat currents and energy densities. In addition the method requires much shorter sequential simulation times than ordinary molecular dynamics techniques, making it applicable within density functional theory. We discuss results for a representative oxide, MgO, at different temperatures and for ordered and nanostructured morphologies, showing the performance of the method in different conditions.

Yaohua Li, Martin Girard, Meng Shen, Jaime Andres Millan, and Monica Olvera de la Cruz

*PNAS* 114, 11838 (2017). DOI: 10.1073/pnas.1713168114

Controlling interactions between proteins and nanoparticles in electrolyte solutions is crucial for advancing biological sciences and biotechnology. The assembly of charged nanoparticles (NPs) and proteins in aqueous solutions can be directed by modifying the salt concentration. High concentrations of monovalent salt can induce the solubilization or crystallization of NPs and proteins. By using a multiscale coarse-grained molecular dynamics approach, we show that, due to ionic correlations in the electrolyte, NPs pairs at high monovalent salt concentrations interact via remarkably strong long-range attractions or repulsions, which can be split into three regimes depending on the surface charge densities of the NPs. NPs with zero-to-low surface charge densities interact via a long-range attraction that is stronger and has a similar range to the depletion attraction induced by polymers with radius of gyrations comparable to the NP diameter. On the other hand, moderately charged NPs with smooth surfaces as well as DNA-functionalized NPs with no possibility of hybridization between them interact via a strong repulsion of range and strength larger than the repulsion predicted by models that neglect ionic correlations, including the Derjaguin–Landau–Vervey–Overbeek (DLVO) model. Interactions between strongly charged NPs (>2 e/nm2), both types smooth and DNA-functionalized NPs, show an attractive potential well at intermediate-to-high salt concentrations, which demonstrates that electrolytes can induce aggregation of strongly charged NPs. Our work provides an improved understanding of the role of ionic correlations in NP assembly and design rules to utilize the salting-out process to crystallize NPs.

Joshua Lequieu, David C. Schwartz and Juan J. de Pablo

*PNAS* 114, E9197 (2017). DOI: 10.1073/pnas.1705685114

Nucleosomes represent the basic building block of chromatin and provide an important mechanism by which cellular processes are controlled. The locations of nucleosomes across the genome are not random but instead depend on both the underlying DNA sequence and the dynamic action of other proteins within the nucleus. These processes are central to cellular function, and the molecular details of the interplay between DNA sequence and nucleosome dynamics remain poorly understood. In this work, we investigate this interplay in detail by relying on a molecular model, which permits development of a comprehensive picture of the underlying free energy surfaces and the corresponding dynamics of nucleosome repositioning. The mechanism of nucleosome repositioning is shown to be strongly linked to DNA sequence and directly related to the binding energy of a given DNA sequence to the histone core. It is also demonstrated that chromatin remodelers can override DNA-sequence preferences by exerting torque, and the histone H4 tail is then identified as a key component by which DNA-sequence, histone modifications, and chromatin remodelers could in fact be coupled.

Meng Shen, Honghao Li, and Monica Olvera de la Cruz

*Phys. Rev. Lett.* 119, 138002 (2017). DOI: 10.1103/PhysRevLett.119.138002

Surface polarization in ion-containing heterogeneous dielectric media such as cell media and emulsions is determined by and determines the positions of the ions. We compute the surface polarization self-consistently as the ions move and analyze their effects on the interactions between electro-neutral, ion-containing droplets using coarse-grained molecular dynamics simulations based on the true energy functional. For water droplets immersed in oil, the interdroplet interaction is attractive, and the surface polarization makes the major contribution. By contrast, for oil droplets in water, the ion-surface induced charge interaction is repulsive and counteracts the attraction between the ions, leading to a small attractive interaction between the droplets. This research improves our understanding of self-assembly in mixed phases such as metal extraction for recovering rare earth elements and nuclear waste as well as water purification.

Peter Scherpelz and Giulia Galli

*Phys. Rev. Materials* 1, 021602(R) (2017). DOI: 10.1103/PhysRevMaterials.1.021602

Surface defects created and probed with scanning tunneling microscopes are a promising platform for 1 atomic-scale electronics and quantuminformation technology applications. Using first-principles calculations we demonstrate how to engineer dangling bond (DB) defects on hydrogenated Si(100) surfaces, which give rise to isolated impurity states that can be used in atomic-scale devices. In particular, we show that sample thickness and biaxial strain can serve as control parameters to design the electronic properties of DB defects. While in thick Si samples the neutral DB state is resonant with bulk valence bands, ultrathin samples (1–2 nm) lead to an isolated impurity state in the gap; similar behavior is seen for DB pairs and DB wires. Strain further isolates the DB from the valence band, with the response to strain heavily dependent on sample thickness. These findings suggest new methods for tuning the properties of defects on surfaces for electronic and quantum information applications. Finally, we present a consistent and unifying interpretation of many results presented in the literature for DB defects on hydrogenated silicon surfaces, rationalizing apparent discrepancies between different experiments and simulations.

Xujun Zhao, Jiyuan Li, Xikai Jiang, Dmitry Karpeev, Olle Heinonen, Barry Smith, Juan P. Hernandez-Ortiz, and Juan J. de Pablo

*J. Chem. Phys.* 146, 244114 (2017). DOI: 10.1063/1.4989545

An efficient parallel Stokes’s solver has been developed for complete description of hydrodynamic interactions between Brownian particles in bulk and confined geomteries. A Langevin description of the particle dynamics is adopted, where the long-range interactions are included using a Green’s function formalism. A scalable parallel computational approach is presented, where the general geometry Stokeslet is calculated following a matrix-free highly-efficient iterative finite-element Stokes’ solver for the accurate treatment of long-range hydrodynamic interactions in arbitrary confined geometries. A combination of mid-point time integration of the Brownian stochastic differential equation, the parallel Stokes’ solver, and a Chebyshev polynomial approximation for the fluctuation-dissipation theorem lead to an O(N) parallel algorithm. We illustrate the new algorithm in the context of the dynamics of confined polymer solutions under equilibrium and non-equilibrium conditions. The method is then extended to treat suspended finite size particles of arbitrary shape in any geometry using an Immersed Boundary approach.

Hythem Sidky and Jonathan K. Whitmer

*J. Phys. Chem. B.* 121, 6691-6698 (2017). DOI: 10.1021/acs.jpcb.7b03113

Chromonic liquid crystals exhibit a unique self-assembly process which is of both theoretical and practical interest. A characteristic feature of chromonics is the occurrence of molecular association through stacking at extremely low concentrations. Experimental evidence has suggested that this process is approximately isodesmic across a broad concentration range. Another important aspect of these phases is a separation of energy scales for elastic deformation that leads to novel bulk morphologies and defect arrangements. To date, only a handful of computational studies have managed to reproduce crucial aspects of chromonic phases, with ionic chromonics treated only by expensive fully atomistic simulations. Here, we present a simple model based on the competition between long-range repulsions and short-range anisotropic attractions capable of capturing all features of the chromonic phase. Molecular simulations of coarse-grained mesogens are used to map out the phase behavior and explore how structural and energetic anisotropies influence their ordering and response. This work presents the first computational investigation into the nematic elasticity of these phases and demonstrates key correlations between elastic response and stack growth.

Tuan Anh Pham, Marco Govoni, Robert Seidel, Stephen E. Bradforth, Eric Schwegler, and Giulia Galli

*Sci. Adv.* 3, 1603210 (2017). DOI: 10.1126/sciadv.1603210

Predicting the electronic properties of aqueous liquids has been a long-standing challenge for quantum mechanical methods. Yet it is a crucial step in understanding and predicting the key role played by aqueous solutions and electrolytes in a wide variety of emerging energy and environmental technologies, including battery and photoelectrochemical cell design. Here we propose an efficient and accurate approach to predict the electronic properties of aqueous solutions, based on the combination of first-principles methods and experimental validation using state-of-the-art spectroscopic measurements. We present results of the photoelectron spectra of a broad range of solvated ions, showing that first-principles molecular dynamics simulations and electronic structure calculations using dielectric hybrid functionals provide a quantitative description of the electronic properties of the solvent and solutes, including excitation energies. The proposed computational framework is general and applicable to other liquids, thereby offering great promise in understanding and engineering solutions and liquid electrolytes for a variety of important energy technologies.

Nicholas P. Brawand, Marco Govoni, Márton Vörös, and Giulia Galli

J. Chem. Theory Comput. 13, 3318 (2017), DOI: 10.1021/acs.jctc.7b00368

We analyze the performance of the recently proposed screened exchange constant functional (SX) (Phys. Rev. X **2016**, 10, 642) on the GW100 test set, and we discuss results obtained at different levels of self-consistency. The SX functional is a generalization of dielectric dependent hybrid functionals to finite systems, it is non-empirical and depends on the average screening of the exchange interaction. We compare results for ionization potentials obtained with SX to those of CCSD(T) calculations and experiments and we find excellent agreement, on par with recent state of the art methods based on many body perturbation theory. Applying SX perturbatively to correct PBE eigenvalues yields improved results in most cases, except for ionic molecules, for which wave function self-consistency is instead crucial. Calculations where wave functions and the screened exchange constant (α^{SX}) are determined self-consistently, and those where α^{SX} is fixed to the value determined within PBE, yield results of comparable accuracy. Perturbative G_{0}W_{0} corrections of eigenvalues obtained with self-consistent α^{SX} are small on average, for all molecules in the GW100 test set.

Matthew Goldey, Nicholas Brawand, Márton Vörös, and Giulia Galli

*J. Chem. Theory Comput.*, 13, 2581-2590 (2017). DOI: 10.1021/acs.jctc.7b00088

The in silico design of novel complex materials for energy conversion requires accurate, ab initio simulation of charge transport. In this work, we present an implementation of constrained density functional theory (CDFT) for the calculation of parameters for charge transport in the hopping regime. We verify our implementation against literature results for molecular systems, and we discuss the dependence of results on numerical parameters and the choice of localization potentials. In addition, we compare CDFT results with those of other commonly used methods for simulating charge transport between nanoscale building blocks. We show that some of these methods give unphysical results for thermally disordered configurations, while CDFT proves to be a viable and robust approach.

Lucas W. Antony, Nicholas E. Jackson, Ivan Lyubimov, Venkatram Vishwanath, Mark D. Ediger, and Juan J. de Pablo

*ACS Central Science*, 3, 415-424 (2017). DOI: 10.1021/acscentsci.7b00041

Organic glass films formed by physical vapor deposition exhibit enhanced stability relative to those formed by conventional liquid cooling and aging techniques. Recently, experimental and computational evidence has emerged indicating that the average molecular orientation can be tuned by controlling the substrate temperature at which these “stable glasses” are grown. In this work, we present a comprehensive all-atom simulation study of ethylbenzene, a canonical stable-glass former, using a computational film formation procedure that closely mimics the vapor deposition process. Atomistic studies of experimentally formed vapor-deposited glasses have not been performed before, and this study therefore begins by verifying that the model and method utilized here reproduces key structural features observed experimentally. Having established agreement between several simulated and experimental macroscopic observables, simulations are used to examine the substrate temperature dependence of molecular orientation. The results indicate that ethylbenzene glasses are anisotropic, depending upon substrate temperature, and that this dependence can be understood from the orientation present at the surface of the equilibrium liquid. By treating ethylbenzene as a simple model for molecular semiconducting materials, a quantum-chemical analysis is then used to show that the vapor-deposited glasses exhibit decreased energetic disorder and increased magnitude of the mean-squared transfer integral relative to isotropic, liquid-cooled films, an effect that is attributed to the anisotropic ordering of the molecular film. These results suggest a novel structure–function simulation strategy capable of tuning the electronic properties of organic semiconducting glasses prior to experimental deposition, which could have considerable potential for organic electronic materials design.

Federico Giberti, Márton Vörös, and Giulia Galli

*Nano Lett.* 17, 2547–2553 (2017). DOI: 10.1021/acs.nanolett.7b00283

Heterogeneous nanostructures, such as quantum dots (QDs) embedded in solid matrices or core–shell nanoparticles, are promising platforms for a wide variety of applications, including phosphors with increased quantum yield, photocatalysis, and solar energy conversion. However, characterizing and controlling their interfacial morphology and defects, which greatly influence their electronic properties, have proven difficult in numerous cases. Here we carried out atomistic calculations on chalcogenide nanostructured materials, i.e., PbSe QDs in CdSe matrices and CdSe embedded in PbSe, and we established how interfacial and core structures affect their electronic properties. In particular, we showed that defects present at interfaces of PbSe nanoparticles and CdSe matrices give rise to detrimental intragap states, degrading the performance of photovoltaic devices. Instead, the electronic gaps of the inverted system (CdSe dots in PbSe) are clean, indicating that this material has superior electronic properties for solar applications. In addition, our calculations predicted that the core structure of CdSe and in turn its band gap may be tuned by applying pressure to the PbSe matrix, providing a means to engineering the properties of new functional materials.

Hadi Ramezani-Dakhel, Monirosadat Sadati, Mohammad Rahimi, Abelardo Ramirez-Hernandez, Benoît Roux, and Juan J. de Pablo

*J. Chem. Theory Comput.* 13, 237–244 (2017). DOI: 10.1021/acs.jctc.6b00844

The ordered environment presented by liquid crystals at interfaces enables a range of novel functionalities that is only now beginning to be exploited in applications ranging from light focusing devices to biosensors. One key feature of liquid crystals is that molecular events occurring at an interface propagate over large distances through the bulk. In spite of their importance, our fundamental understanding of liquid crystal–water and liquid crystal–air interfaces remains limited. In this work, we present results from large-scale atomistic molecular dynamics simulations on the organization of the nematic and isotropic phases of the nitrile-containing mesogenic molecule 4-cyano-4′-pentylbiphenyl (5CB) in the vicinity of vacuum and aqueous interfaces. Hybrid boundary conditions are imposed by confining 5CB films between vacuum and an aqueous medium to examine how those two types of interfaces influence the specific structural arrangement and ordering of 5CB. Consistent with experiments, our results indicate that 5CB exhibits homeotropic anchoring at the vacuum interface, and planar alignment at aqueous interfaces. Two-dimensional molecular dynamics potential of mean force calculations and average polarization densities show that the polar nitrile group of 5CB remains hydrated near the aqueous interface, where it modulates the orientation of water molecules. Estimates of the anchoring strength reveal an oscillatory decay and a semilinear decay with distance from the interface in vacuum and water, respectively.

Monirosadat Sadati, Hadi Ramezani-Dakhel, Wei Bu, Emre Sevgen, Zhu Liang, Cem Erol, Mohammad Rahimi, Nader Taheri Qazvini, Binhua Lin, Nicholas L. Abbott, Benoı̂t Roux, Mark L. Schlossman, and Juan J. de Pablo

*J. Am. Chem. Soc.* 139, 3841–3850 (2017). DOI: 10.1021/jacs.7b00167

Numerous applications of liquid crystals rely on control of molecular orientation at an interface. However, little is known about the precise molecular structure of such interfaces. In this work, synchrotron X-ray reflectivity measurements, accompanied by large-scale atomistic molecular dynamics simulations, are used for the first time to reconstruct the air-liquid crystal interface of a nematic material, namely, 4-pentyl-4′-cyanobiphenyl (5CB). The results are compared to those for 4-octyl-4′-cyanobiphenyl (8CB) which, in addition to adopting isotropic and nematic states, can also form a smectic phase. Our findings indicate that the air interface imprints a highly ordered structure into the material; such a local structure then propagates well into the bulk of the liquid crystal, particularly for nematic and smectic phases.

Alex P. Gaiduk and Giulia Galli

*J. Phys. Chem. Lett.* 8, 1496–1502 (2017). DOI: 10.1021/acs.jpclett.7b00239

Determining how the structure of water is modified by the presence of salts is instrumental to understanding the solvation of biomolecules and, in general, the role played by salts in biochemical processes. However, the extent of hydrogen bonding disruption induced by salts remains controversial. We performed extensive first-principles simulations of solutions of a simple salt (NaCl) and found that, while the cation does not significantly change the structure of water beyond the first solvation shell, the anion has a further reaching effect, modifying the hydrogen-bond network even outside its second solvation shell. We found that a distinctive fingerprint of hydrogen bonding modification is the change in polarizability of water molecules. Molecular dipole moments are instead insensitive probes of long-range modifications induced by Na^{+} and Cl^{–} ions. Though noticeable, the long-range effect of Cl^{–} is expected to be too weak to affect solubility of large biomolecules.

Margaret Hudson, Dmitriy Dolzhnikov, Alexander S. Filatov, Eric M. Janke, Jaeyoung Jang, Byeongdu Lee, Chengjun Sun, and Dmitri Talapin

*J. Am. Chem. Soc.* 139, 3368-3377 (2017). DOI: 10.1021/jacs.6b10077

This work investigates the structure and properties of soluble chalcogenidocadmates, a molecular form of cadmium chalcogenides with unprecedented one-dimensional bonding motifs. The single crystal X-ray structure reveals that sodium selenocadmate consists of infinite one-dimensional wires of (Cd2Se3)n2n– charge balanced by Na+ and stabilized by coordinating solvent molecules. Exchanging the sodium cation with tetraethylammonium or didodecyldimethylammonium expands the versatility of selenocadmate by improving its solubility in a variety of polar and nonpolar solvents without changing the anion structure and properties. The introduction of a micelle-forming cationic surfactant allows for the templating of selenocadmate, or the analogous telluride species, to create ordered organic–inorganic hybrid CdSe or CdTe mesostructures. Finally, the interaction of selenocadmate “wires” with Cd2+ ions creates an unprecedented gel-like form of stoichiometric CdSe. We also demonstrate that these low-dimensional CdSe species show characteristic semiconductor behavior, and can be used in photodetectors and field-effect transistors.

Yuping He and Giulia Galli

*Chem. Mater.* 29, 682–689 (2017). DOI: 10.1021/acs.chemmater.6b04300

We carried out calculations based on density functional theory to investigate the electronic, vibrational, and dielectric properties of mixed halide perovskites CH_{3}NH_{3} AI_{3–x}Cl_{x} with A = Pb and Sn. Computed free energies indicated that Cl mixed systems may be formed only for Cl concentrations not exceeding 10^{19} cm^{–3}, and phonon calculations showed that the disorder induced in the host lattice by the presence of a smaller halogen is responsible for mechanical instabilities. However, we found that the presence of chloride may be beneficial to the electronic properties of the perovskites. Chloride anions cause the organic cations to be displaced from the center of the cage; such a displacement induces preferential orientations of the cation dipole, which in turn are responsible for notable changes in the dielectric properties of the material and possibly for the formation of local ferroelectric domains. The latter are instrumental in separating electron hole pairs and hence in contributing to long charge-carrier diffusion lengths, in spite of polarons being more likely formed in mixed perovksites than in CH_{3}NH_{3} AI_{3}.

Daniel Reid, Ivan Lyubimov, M. D. Ediger, and Juan J. de Pablo

*Nat. Commun.* 7, 13062 (2016). DOI: 10.1038/ncomms13062

Glass films prepared by a process of physical vapour deposition have been shown to have thermodynamic and kinetic stability comparable to those of ordinary glasses aged for thousands of years. A central question in the study of vapour-deposited glasses, particularly in light of new knowledge regarding anisotropy in these materials, is whether the ultra-stable glassy films formed by vapour deposition are ever equivalent to those obtained by liquid cooling. Here we present a computational study of vapour deposition for a two-dimensional glass forming liquid using a methodology, which closely mimics experiment. We find that for the model considered here, structures that arise in vapour-deposited materials are statistically identical to those observed in ordinary glasses, provided the two are compared at the same inherent structure energy. We also find that newly deposited hot molecules produce cascades of hot particles that propagate far into the film, possibly influencing the relaxation of the material.

Nicholas P. Brawand, Márton Vörös, Marco Govoni, and Giulia Galli

*Phys. Rev. X* 6, 041002 (2016). DOI: 10.1103/PhysRevX.6.041002

The accurate prediction of electronic and optical properties of molecules and solids is a persistent challenge for methods based on density functional theory. We propose a generalization of dielectric-dependent hybrid functionals to finite systems where the definition of the mixing fraction of exact and semilocal exchange is physically motivated, nonempirical, and system dependent. The proposed functional yields ionization potentials, and fundamental and optical gaps of many, diverse molecular systems in excellent agreement with experiments, including organic and inorganic molecules and semiconducting nanocrystals. We further demonstrate that this hybrid functional gives the correct alignment between energy levels of the exemplary TTF-TCNQ donor-acceptor system.

Michael A. Boles, Michael Engel, and Dmitri V. Talapin

*Chem. Rev.* 116, 11220–11289 (2016). DOI: 0.1021/acs.chemrev.6b00196

Chemical methods developed over the past two decades enable preparation of colloidal nanocrystals with uniform size and shape. These Brownian objects readily order into superlattices. Recently, the range of accessible inorganic cores and tunable surface chemistries dramatically increased, expanding the set of nanocrystal arrangements experimentally attainable. In this review, we discuss efforts to create next-generation materials via bottom-up organization of nanocrystals with preprogrammed functionality and self-assembly instructions. This process is often driven by both interparticle interactions and the influence of the assembly environment. The introduction provides the reader with a practical overview of nanocrystal synthesis, self-assembly, and superlattice characterization. We then summarize the theory of nanocrystal interactions and examine fundamental principles governing nanocrystal self-assembly from hard and soft particle perspectives borrowed from the comparatively established fields of micrometer colloids and block copolymer assembly. We outline the extensive catalog of superlattices prepared to date using hydrocarbon-capped nanocrystals with spherical, polyhedral, rod, plate, and branched inorganic core shapes, as well as those obtained by mixing combinations thereof. We also provide an overview of structural defects in nanocrystal superlattices. We then explore the unique possibilities offered by leveraging nontraditional surface chemistries and assembly environments to control superlattice structure and produce nonbulk assemblies. We end with a discussion of the unique optical, magnetic, electronic, and catalytic properties of ordered nanocrystal superlattices, and the coming advances required to make use of this new class of solids.

Wooram Park, Sin-Jung Park, Soojeong Cho, Heejun Shin, Young-Seok Jung, Byeongdu Lee, Kun Na, and Dong-Hyun Kim

*J. Am. Chem. Soc.* 138, 10734-10737 (2016). DOI: 10.1021/jacs.6b04875

We developed a thermoswitchable polymeric photosensitizer (T-PPS) by conjugating PS (Pheophorbide-a, PPb-a) to a temperature-responsive polymer backbone of biocompatible hydroxypropyl cellulose. Self-quenched PS molecules linked in close proximity by π–π stacking in T-PPS were easily transited to an active monomeric state by the temperature-induced phase transition of polymer backbones. The temperature-responsive intermolecular interaction changes of PS molecules in T-PPS were demonstrated in synchrotron small-angle X-ray scattering and UV–vis spectrophotometer analysis. The T-PPS allowed switchable activation and synergistically enhanced cancer cell killing effect at the hyperthermia temperature (45 °C). Our developed T-PPS has the considerable potential not only as a new class of photomedicine in clinics but also as a biosensor based on temperature responsiveness.

Joshua Lequieu, Andrés Córdoba, Daniel Hinckley, and Juan J. de Pablo

*ACS Cent. Sci.* 2, 614-620 (2016). DOI: 10.1021/acscentsci.6b00170

The self-assembly of DNA-conjugated nanoparticles represents a promising avenue toward the design of engineered hierarchical materials. By using DNA to encode nanoscale interactions, macroscale crystals can be formed with mechanical properties that can, at least in principle, be tuned. Here we present in silico evidence that the mechanical response of these assemblies can indeed be controlled, and that subtle modifications of the linking DNA sequences can change the Young’s modulus from 97 kPa to 2.1 MPa. We rely on a detailed molecular model to quantify the energetics of DNA–nanoparticle assembly and demonstrate that the mechanical response is governed by entropic, rather than enthalpic, contributions and that the response of the entire network can be estimated from the elastic properties of an individual nanoparticle. The results here provide a first step toward the mechanical characterization of DNA–nanoparticle assemblies, and suggest the possibility of mechanical metamaterials constructed using DNA.

Joshua Lequieu, Andrés Córdoba, David C. Schwartz, and Juan J. de Pablo

*ACS Cent. Sci.* 2, 660-666 (2016). DOI: 10.1021/acscentsci.6b00201

Nucleosomes form the basic unit of compaction within eukaryotic genomes, and their locations represent an important, yet poorly understood, mechanism of genetic regulation. Quantifying the strength of interactions within the nucleosome is a central problem in biophysics and is critical to understanding how nucleosome positions influence gene expression. By comparing to single-molecule experiments, we demonstrate that a coarse-grained molecular model of the nucleosome can reproduce key aspects of nucleosome unwrapping. Using detailed simulations of DNA and histone proteins, we calculate the tension-dependent free energy surface corresponding to the unwrapping process. The model reproduces quantitatively the forces required to unwrap the nucleosome and reveals the role played by electrostatic interactions during this process. We then demonstrate that histone modifications and DNA sequence can have significant effects on the energies of nucleosome formation. Most notably, we show that histone tails contribute asymmetrically to the stability of the outer and inner turn of nucleosomal DNA and that depending on which histone tails are modified, the tension-dependent response is modulated differently.

Xikai Jiang, Jiyuan Li, Xujun Zhao, Jian Qin, Dmitry Karpeev, Juan Hernandez-Ortiz, Juan J. de Pablo,and Olle Heinonen

*J. Chem. Phys.* 145, 064307 (2016). DOI: 10.1063/1.4960436

Large classes of materials systems in physics and engineering are governed by magnetic and electrostatic interactions. Continuum or mesoscale descriptions of such systems can be cast in terms of integral equations, whose direct computational evaluation requires O(N^{2}) operations, where N is the number of unknowns. Such a scaling, which arises from the many-body nature of the relevant Green’s function, has precluded wide-spread adoption of integral methods for solution of large-scale scientific and engineering problems. In this work, a parallel computational approach is presented that relies on using scalable open source libraries and utilizes a kernel-independent Fast Multipole Method (FMM) to evaluate the integrals in O(N) operations, with O(N) memory cost, thereby substantially improving the scalability and efficiency of computational integral methods. We demonstrate the accuracy, efficiency, and scalability of our approach in the context of two examples. In the first, we solve a boundary value problem for a ferroelectric/ferromagnetic volume in free space. In the second, we solve an electrostatic problem involving polarizable dielectric bodies in an unbounded dielectric medium. The results from these test cases show that our proposed parallel approach, which is built on a kernel-independent FMM, can enable highly efficient and accurate simulations and allow for considerable flexibility in a broad range of applications.

Vikramjit Singh Rathee, Siyi Qu, William A. Phillip, and Jonathan K. Whitmer

*Mol. Syst. Des. Eng.* 1, 301-312 (2016). DOI: 10.1039/C6ME00045B

Membrane-based filtration has shown promise in the treatment of municipal and industrial water supplies. Most existing filtration membranes rely on size-selective rejection mechanisms and do not offer a systematic means to tune the chemistry of the membrane and, in turn, the rejection of solutes based on chemical factors. Technological advancements in copolymer self-assembly have permitted the creation of highly porous membranes through arrested microphase separation. The resulting structures exhibit nanoscale pores whose size and chemistry are controlled by the macromolecular architecture of the copolymer precursors. Here, we examine in detail how the structure and chemistry of charge-functionalized copolymers affect their performance as nanofiltration membranes through experiments and coarse-grained molecular simulation of membranes formed from poly(acrylonitrile-r-oligo(ethylene glycol) methyl ether methacrylate-r-glycidyl methacrylate) (P(AN-r-OEGMA-r-GMA)). Charge-selective moieties were introduced into the membranes through reactions with diamine salts of varying length, enabling these membranes to offer preferential rejection of multivalent ionic species. Using a minimal model, we explore the essential thermodynamics features of salt rejection in these membranes, and develop a model capable of linking the performance of the membranes to their molecular character through measurement of the underlying free energy profile. Our key result is an effective membrane porosity which quantitatively describes experiments for a wide variety of situations after calibrating a single parameter. These results demonstrate that at moderate ion concentrations and flow rates, the selection of ions is driven predominantly by equilibrium thermodynamics. Going forward, the ability of the model to capture the rejection of dissolved solutes can be leveraged to design optimized membranes for targeted performance profiles in silico before undertaking time consuming batteries of Edisonian experiments.

Hythem Sidky and Jonathan K. Whitmer

*Liquid Crystals* 43, 2285-2299 (2016). DOI: 10.1080/02678292.2016.1201869

A free energy perturbation method is used to systematically study the elastic properties of four common Gay–Berne nematogenic models; two with a length-to-diameter ratio κ = 3 [(3, 5, 1, 2) and (3, 5, 1, 3)], a model with κ = 4.4 parameterised for p-terphenyl (4.4, 20, 1, 1), and a discogen with κ = 0.345 (0.345, 0.2, 1, 2). In accordance with previous measurements, we find that for κ = 3, models, . We additionally find the latter two models in particular accurately capture the experimentally measured elastic ratios in apolar achiral systems. The (4.4, 20, 1, 1) model reproduces the elastic constant ratios of p-azoxyanisole remarkably well, and maps to within 30% of the absolute. The (0.345, 0.2, 1, 2) model elastic constants exhibit an unusual temperature dependence similar to recent experimental studies. Here we find , in line with theoretical predictions. All models deviate from the mean-field expectation kii ∝ S2. These results represent a crucial first step towards quantitatively accurate coarse-grained liquid crystalline models of self-assembly and response, enabling one to choose a Gay–Berne model based on its measured elastic ratios rather than just its shape and energy anisotropy.

Peter Scherpelz, Marco Govoni, Ikutaro Hamada, and Giulia Galli

*J. Chem. Theory Comput.* 12, 3523-3544 (2016). DOI: 10.1021/acs.jctc.6b00114

We present an implementation of G_{0}W_{0} calculations including spin–orbit coupling (SOC) enabling investigations of large systems, with thousands of electrons, and we discuss results for molecules, solids, and nanocrystals. Using a newly developed set of molecules with heavy elements (called GW-SOC81), we find that, when based upon hybrid density functional calculations, fully relativistic (FR) and scalar-relativistic (SR) G_{0}W_{0} calculations of vertical ionization potentials both yield excellent performance compared to experiment, with errors below 1.9%. We demonstrate that while SR calculations have higher random errors, FR calculations systematically underestimate the VIP by 0.1 to 0.2 eV. We further verify that SOC effects may be well approximated at the FR density functional level and then added to SR G_{0}W_{0} results for a broad class of systems. We also address the use of different root-finding algorithms for the G_{0}W_{0} quasiparticle equation and the significant influence of including d electrons in the valence partition of the pseudopotential for G_{0}W_{0} calculations. Finally, we present statistical analyses of our data, highlighting the importance of separating definitive improvements from those that may occur by chance due to a limited number of samples. We suggest the statistical analyses used here will be useful in the assessment of the accuracy of a large variety of electronic structure methods.

Hythem Sidky and Jonathan K. Whitmer

*Soft Matter* 12, 4489-4498 (2016). DOI: 10.1039/C5SM03107A

Utilizing density-of-states simulations, we perform a full mapping of the phase behavior and elastic responses of binary liquid crystalline mixtures represented by the multicomponent Lebwohl–Lasher model. Our techniques are able to characterize the complete phase diagram, including nematic–nematic phase separation predicted by mean-field theories, but previously not observed in simulations. Mapping this phase diagram permits detailed study of elastic properties across the miscible nematic region. Importantly, we observe for the first time local phase separation and disordering driven by the application of small linear perturbations near the transition temperature and more significantly through nonlinear stresses. These findings are of key importance in systems of blended nematics which contain particulate inclusions, or are otherwise confined.

Alex P. Gaiduk, Marco Govoni, Robert Seidel, Jonathan Skone, Bernd Winter, and Giulia Galli

*J. Am. Chem. Soc.* 138, 6912-6915 (2016). DOI: 10.1021/jacs.6b00225

We present a combined computational and experimental study of the photoelectron spectrum of a simple aqueous solution of NaCl. Measurements were conducted on microjets, and first-principles calculations were performed using hybrid functionals and many-body perturbation theory at the G0W0 level, starting with wave functions computed in ab initio molecular dynamics simulations. We show excellent agreement between theory and experiments for the positions of both the solute and solvent excitation energies on an absolute energy scale and for peak intensities. The best comparison was obtained using wave functions obtained with dielectric-dependent self-consistent and range-separated hybrid functionals. Our computational protocol opens the way to accurate, predictive calculations of the electronic properties of electrolytes, of interest to a variety of energy problems.

Rui Zhang, Tyler Roberts, Igor S. Aranson, and Juan J. de Pablo

*J. Chem. Phys.* 144, 084905 (2016). DOI: 10.1063/1.4940342

Liquid crystals (LCs) display many of the flow characteristics of liquids but exhibit long range orientational order. In the nematic phase, the coupling of structure and flow leads to complex hydrodynamic effects that remain to be fully elucidated. Here, we consider the hydrodynamics of a nematic LC in a hybrid cell, where opposite walls have conflicting anchoring boundary conditions, and we employ a 3D lattice Boltzmann method to simulate the time-dependent flow patterns that can arise. Due to the symmetry breaking of the director field within the hybrid cell, we observe that at low to moderate shear rates, the volumetric flow rate under Couette and Poiseuille flows is different for opposite flow directions. At high shear rates, the director field may undergo a topological transition which leads to symmetric flows. By applying an oscillatory pressure gradient to the channel, a net volumetric flow rate is found to depend on the magnitude and frequency of the oscillation, as well as the anchoring strength. Taken together, our findings suggest several intriguing new applications for LCs in microfluidic devices.