Intrinsic Stress Field for Liquid Surfaces.

The microscopic stress field inhomogeneity in the interfacial region adjacent to the liquid surface is the fundamental origin of the liquid surface tension, but because of broadening due to capillary fluctuations, a detailed molecular level understanding of the stress field remains elusive. In this work, we deconvolute the capillary fluctuations to reveal the intrinsic stress field and show that the atomic-level contributions to the surface tension are similar in functional form across a variety of monatomic systems. These contributions are confined to an interfacial region approximately 1.5±0.1 times the particle diameter for all systems studied. In addition, the intrinsic density and stress profiles show a strong spatial correlation that should be useful in the development of a statistical mechanical theory for the prediction of surface stress and surface tension.

Analysis of probability of inserting a hard spherical particle with small diameter in hard-sphere fluid.

The probability of inserting, without overlap, a hard spherical particle of diameter σ in a hard-sphere fluid of diameter σ0 and packing fraction η determines its excess chemical potential at infinite dilution, µex(σ, η). In our previous work [R. L. Davidchack and B. B. Laird, J. Chem. Phys. 157, 074701 (2022)], we used Widom's particle insertion method within molecular dynamics simulations to obtain high precision results for µex(σ, η) with σ/σ0 ≤ 4 and η ≤ 0.5. In the current work, we investigate the behavior of this quantity at small σ. In particular, using the inclusion-exclusion principle, we relate the insertion probability to the hard-sphere fluid distribution functions and thus derive the higher-order terms in the Taylor expansion of µex(σ, η) at σ = 0. We also use direct evaluation of the excluded volume for pairs and triplets of hard spheres to obtain simulation results for µex(σ, η) at σ/σ0 ≤ 0.2247 that are of much higher precision than those obtained earlier with Widom's method. These results allow us to improve the quality of the small-σ correction in the empirical expression for µex(σ, η) presented in our previous work.

Generation of Amorphous Silica Surfaces with Controlled Roughness.

Amorphous silica (a-SiO2) surfaces, when grafted with select metals on the active sites of the functionalized surfaces, can act as useful heterogeneous catalysts. From a molecular modeling perspective, one challenge has been generating a-SiO2 slab models with controllable surface roughness to facilitate the study of the effect of surface morphology on the material properties. Previous computational methods either generate relatively flat surfaces or periodically corrugated surfaces that do not mimic the full range of potential surface roughness of the amorphous silica material. In this work, we present a new method, inspired by the capillary fluctuation theory of interfaces, in which rough silica slabs are generated by cleaving a bulk amorphous sample using a cleaving plane with Fourier components randomly generated from a Gaussian distribution. The width of this Gaussian distribution (and thus the degree of surface roughness) can be tuned by varying the surface roughness parameter α. Using the van Beest, Kramer, and van Santen (BKS) force field, we create a large number of silica slabs using cleaving surfaces of varying roughness (α) and using two different system sizes. These surfaces are then characterized to determine their roughness (mean-squared displacement), density profile, and ring size distribution. This analysis shows a higher concentration of surface defects (under-/overcoordinated atoms and strained rings) as the surface roughness increases. To examine the effect of the roughness on surface reactivity, we re-equilibriate a subset of these slabs using the reactive force field ReaxFF and then expose the slabs to water and observe the formation of surface silanols. We observe that the rougher surfaces exhibit higher silanol concentrations as well as bimodal acidity.

Atomistic characterization of the SiO2 high-density liquid/low-density liquid interface.

The equilibrium silica liquid-liquid interface between the high-density liquid (HDL) phase and the low-density liquid (LDL) phase is examined using molecular-dynamics simulation. The structure, thermodynamics, and dynamics within the interfacial region are characterized in detail and compared with previous studies on the liquid-liquid phase transition (LLPT) in bulk silica, as well as traditional crystal-melt interfaces. We find that the silica HDL-LDL interface exhibits a spatial fragile-to-strong transition across the interface. Calculations of dynamics properties reveal three types of dynamical heterogeneity hybridizing within the silica HDL-LDL interface. We also observe that as the interface is traversed from HDL to LDL, the Si/O coordination number ratio jumps to an unexpectedly large value, defining a thin region of the interface where HDL and LDL exhibit significant mixing. In addition, the LLPT phase coexistence is interpreted in the framework of the traditional thermodynamics of alloys and phase equilibria.

Local collective dynamics at equilibrium BCC crystal-melt interfaces.

We present a classical molecular-dynamics study of the collective dynamical properties of the coexisting liquid phase at equilibrium body-centered cubic (BCC) Fe crystal-melt interfaces. For the three interfacial orientations (100), (110), and (111), the collective dynamics are characterized through the calculation of the intermediate scattering functions, dynamical structure factors, and density relaxation times in a sequential local region of interest. An anisotropic speedup of the collective dynamics in all three BCC crystal-melt interfacial orientations is observed. This trend differs significantly from the previously observed slowing down of the local collective dynamics at the liquid-vapor interface [del Rio and González, Acta Mater. 198, 281 (2020)]. Examining the interfacial density relaxation times, we revisit the validity of the recently developed time-dependent Ginzburg-Landau theory for the solidification crystal-melt interface kinetic coefficients, resulting in excellent agreement with both the magnitude and the kinetic anisotropy of the crystal-melt interface kinetic coefficients measured from the non-equilibrium molecular-dynamics simulations.

Inside and out: Surface thermodynamics from positive to negative curvature.

To explore the curvature dependence of solid-fluid interfacial thermodynamics, we calculate, using Grand Canonical Monte Carlo simulation, the surface free energy for a 2d hard-disk fluid confined in a circular hard container of radius R as a function of the bulk packing fraction η and wall curvature C̄=-1/R. (The curvature is negative because the surface is concave.) Combining this with our previous data [Martin et al., J. Phys. Chem. B 124, 7938-7947 (2020)] for the positive curvature case (a hard-disk fluid at a circular wall, C̄=+1/R), we obtain a complete picture of surface thermodynamics in this system over the full range of positive and negative wall curvatures. Our results show that γ is linear in C̄ with a slope that is the same for both positive and negative wall curvatures, with deviations seen only at high negative curvatures (strong confinement) and high density. This observation indicates that the surface thermodynamics of this system is consistent with the predictions of so-called morphometric thermodynamics at both positive and negative curvatures. In addition, we show that classical density functional theory and a generalized scaled particle theory can be constructed that give excellent agreement with the simulation data over most of the range of curvatures and densities. For extremely high curvatures, where only one or two disks can occupy the container at maximum packing, it is possible to calculate γ exactly. In this limit, the simulations and density functional theory calculations are in remarkable agreement with the exact results.

Chemical potential and surface free energy of a hard spherical particle in hard-sphere fluid over the full range of particle diameters.

The excess chemical potential µex(σ, η) of a test hard spherical particle of diameter σ in a fluid of hard spheres of diameter σ0 and packing fraction η can be computed with high precision using Widom's particle insertion method [B. Widom, J. Chem. Phys. 39, 2808 (1963)] for σ between 0 and just larger than 1 and/or small η. Heyes and Santos [J. Chem. Phys. 145, 214504 (2016)] analytically showed that the only polynomial representation of µex consistent with the limits of σ at zero and infinity has a cubic form. On the other hand, through the solvation free energy relationship between µex and the surface free energy γ of hard-sphere fluids at a hard spherical wall, we can obtain precise measurements of µex for large σ, extending up to infinity (flat wall) [R. L. Davidchack and B. B. Laird, J. Chem. Phys. 149, 174706 (2018)]. Within this approach, the cubic polynomial representation is consistent with the assumptions of morphometric thermodynamics. In this work, we present the measurements of µex that combine the two methods to obtain high-precision results for the full range of σ values from zero to infinity, which show statistically significant deviations from the cubic polynomial form. We propose an empirical functional form for the µex dependence on σ and η, which better fits the measurement data while remaining consistent with the analytical limiting behavior at zero and infinite σ.

Molecular Simulations of Phase Equilibria and Transport Properties in a Model CO2-Expanded Lithium Perchlorate Electrolyte.

Carbon-dioxide (CO2)-expanded liquids, in which a significant mole fraction of CO2 is dissolved into an organic solvent, have been of significant interest, especially as catalytic support media. Because the CO2 mole fraction and density can be controlled over a significant range by changing the CO2 partial pressure, the transport properties of these solvents are highly tunable. Recently, these liquids have garnered interest as potential electrolyte solutions for catalytic electrochemistry; however, little is currently known about the influence of the electrolyte on CO2 expansion. In the present work, we use molecular-dynamics simulations to study diffusion and viscosity in a model lithium perchlorate electrolyte in CO2-expanded acetonitrile and demonstrate that these properties are highly dependent on the concentration of the electrolyte. Our present results indicate that the electrolyte slows down diffusion of both CO2 and MeCN, and that the slowed diffusion in the former is driven by changes in the activation entropy, whereas slowed diffusion in the latter is driven by changes in the activation energy.

Surface Free Energy of a Hard-Disk Fluid at Curved Hard Walls: Theory and Simulation.

In this work, we examine the surface thermodynamics of a hard-disk fluid at curved hard walls using Monte Carlo (MC) simulation and a generalized scaled particle theory (gSPT). The curved walls are modeled as hard disks of varying radii, R. The surface free energy, γ, and excess surface volume, vex, for this system are calculated as functions of both the fluid packing fraction and the wall radius. The simulation results are used to test, for this system, the assumptions of morphometric thermodynamics (MT), which predicts that both γ and vex are linear functions of the surface curvature, 1/R, for a two-dimensional system. In addition, we compare the simulation results to the gSPT developed in this work, as well as with virial expansions derived from the known virial coefficients of the binary hard-sphere fluid. At low to intermediate packing fractions, the non-MT terms (terms of higher order than 1/R in a expansion of γ and vex) of γ are zero within the simulation error; however, at the highest densities, deviations from MT become significant, similar to what was seen in our earlier simulation work on the three-dimensional hard-sphere/hard-wall system. In addition, the new gSPT gives improved results for both γ and vex over standard scaled particle theory (SPT) but underestimates the deviations from MT at high density.

Pressure and Temperature Tuning of Gas-Expanded Liquid Structure and Dynamics.

Carbon dioxide-expanded liquids (CXLs) represent an important class of reaction media that provide tunability of mass transport, solvation, and solubility. Their properties have been demonstrated to provide advantages over traditional organic solvents. However, the molecular-level effects of the CO2 expansion on the structure and dynamics of the liquid that lead to this result have not been fully explored. To address this question, we have used molecular simulations to examine the behavior of two CXLs relevant to the hydroformylation of 1-octene, which has been demonstrated to benefit from the use of gas-expanded reaction media. Specifically, the phase equilibrium properties of CO2-expanded 1-octene and nonanal are calculated as functions of temperature and pressure using Gibbs ensemble Monte Carlo simulations to determine the pressure-composition phase diagrams and volume expansion. In addition, molecular dynamics (MD) simulations were conducted to compute the liquid structure, diffusion coefficients, and shear viscosities. The simulated phase diagrams are in excellent agreement with previous experimental data when available, validating the models used. The MD simulations reveal a direct, linear relationship between the liquid viscosity and the volume expansion, which has not been previously reported. In contrast, deviations from such a relationship are observed for the diffusion coefficient at large volume expansion, indicating that a single Stokes-Einstein relation cannot describe the behavior at all pressures.

Surface free energy of a hard-sphere fluid at curved walls: Deviations from morphometric thermodynamics.

We report molecular-dynamics (MD) simulation results for the surface free energy of a hard-sphere fluid at cylindrical and spherical hard walls of different radii. The precision of the results is much higher than that in our previous study [B. B. Laird et al., Phys. Rev. E 86, 060602 (2012)], allowing us to estimate the size of deviations from the predictions of Morphometric Thermodynamics (MT). We compare our results to the analytical expressions for the surface energy as a function of wall radius R and fluid density derived from the White Bear II variant of the density functional theory, as well as to the leading terms of the virial expansion. For the cylindrical wall, we observe deviations from MT proportional to R -2 and R -3, which are consistent with the available virial expressions. For the spherical wall, while the precision is not sufficient to detect statistically significant deviations from MT, the MD results indicate the range of densities for which the truncated virial expansions are applicable.

Thermodynamics of the hard-disk fluid at a planar hard wall: Generalized scaled-particle theory and Monte Carlo simulation.

A generalized scaled-particle theory for the uniform hard-disk mixture is derived in the spirit of the White Bear II free energy of the hard-sphere fluid [H. Hansen-Goos and R. Roth, J. Phys. C: Condens. Matter 18, 8413 (2006)]. The theory provides a very simple result for the interfacial free energy γ of the hard-disk fluid at a planar hard wall (which in d = 2 is a line) in terms of the equation of state. To complement and assess the theory, we perform Monte Carlo simulations from which we obtain γ using Gibbs-Cahn integration. While we find excellent overall agreement between theory and simulation, it also becomes apparent that the set of scaled-particle variables available in d = 2 is too limited, prohibiting a quasi-exact result for γ. Furthermore, this is reflected in the mixture equation of state resulting from our theory, which, similar to a previous attempt by Santos et al. [Mol. Phys. 96, 1 (1999)], displays a small but systematic deviation from simulations.

Properties of the hard-sphere fluid at a planar wall using virial series and molecular-dynamics simulation.

We study the hard-sphere fluid in contact with a planar hard wall. By combining the inhomogeneous virial series with simulation results, we achieve a new benchmark of accuracy for the calculation of surface thermodynamics properties such as surface adsorption Γ and the surface free energy (or surface tension), γ. We briefly introduce the problem of choosing a position for the dividing surface and avoid it by proposing the use of alternative functions to Γ and γ that are independent of the adopted frame of reference. Finally, we present analytic expressions for the dependence of system surface thermodynamic properties on packing fraction, ensuring the high accuracy of the parameterized functions for any frame of reference. The proposed parametric expressions for both, Γ and γ, fit the accurate simulation results within the statistical error.

Electric potential calculation in molecular simulation of electric double layer capacitors.

For the molecular simulation of electric double layer capacitors (EDLCs), a number of methods have been proposed and implemented to determine the one-dimensional electric potential profile between the two electrodes at a fixed potential difference. In this work, we compare several of these methods for a model LiClO4-acetonitrile/graphite EDLC simulated using both the traditional fixed-charged method (FCM), in which a fixed charge is assigned a priori to the electrode atoms, or the recently developed constant potential method (CPM) (2007 J. Chem. Phys. 126 084704), where the electrode charges are allowed to fluctuate to keep the potential fixed. Based on an analysis of the full three-dimensional electric potential field, we suggest a method for determining the averaged one-dimensional electric potential profile that can be applied to both the FCM and CPM simulations. Compared to traditional methods based on numerically solving the one-dimensional Poisson's equation, this method yields better accuracy and no supplemental assumptions.

Orientation dependence of heterogeneous nucleation at the Cu-Pb solid-liquid interface.

In this work, we examine the effect of surface structure on the heterogeneous nucleation of Pb crystals from the melt at a Cu substrate using molecular-dynamics (MD) simulation. In a previous work [Palafox-Hernandez et al., Acta Mater. 59, 3137 (2011)] studying the Cu/Pb solid-liquid interface with MD simulation, we observed that the structure of the Cu(111) and Cu(100) interfaces was significantly different at 625 K, just above the Pb melting temperature (618 K for the model). The Cu(100) interface exhibited significant surface alloying in the crystal plane in contact with the melt. In contrast, no surface alloying was seen at the Cu(111) interface; however, a prefreezing layer of crystalline Pb, 2-3 atomic planes thick and slightly compressed relative to bulk Pb crystal, was observed to form at the interface. We observe that at the Cu(111) interface the prefreezing layer is no longer present at 750 K, but surface alloying in the Cu(100) interface persists. In a series of undercooling MD simulations, heterogeneous nucleation of fcc Pb is observed at the Cu(111) interface within the simulation time (5 ns) at 592 K-a 26 K undercooling. Nucleation and growth at Cu(111) proceeded layerwise with a nearly planar critical nucleus. Quantitative analysis yielded heterogeneous nucleation barriers that are more than two orders of magnitude smaller than the predicted homogeneous nucleation barriers from classical nucleation theory. Nucleation was considerably more difficult on the Cu(100) surface-alloyed substrate. An undercooling of approximately 170 K was necessary to observe nucleation at this interface within the simulation time. From qualitative observation, the critical nucleus showed a contact angle with the Cu(100) surface of over 90°, indicating poor wetting of the Cu(100) surface by the nucleating phase, which according to classical heterogeneous nucleation theory provides an explanation of the large undercooling necessary to nucleate on the Cu(100) surface, relative to Cu(111), whose surface is more similar to the nucleating phase due to the presence of the prefreezing layer.

Evaluation of the constant potential method in simulating electric double-layer capacitors.

A major challenge in the molecular simulation of electric double layer capacitors (EDLCs) is the choice of an appropriate model for the electrode. Typically, in such simulations the electrode surface is modeled using a uniform fixed charge on each of the electrode atoms, which ignores the electrode response to local charge fluctuations in the electrolyte solution. In this work, we evaluate and compare this Fixed Charge Method (FCM) with the more realistic Constant Potential Method (CPM), [S. K. Reed et al., J. Chem. Phys. 126, 084704 (2007)], in which the electrode charges fluctuate in order to maintain constant electric potential in each electrode. For this comparison, we utilize a simplified LiClO4-acetonitrile/graphite EDLC. At low potential difference (ΔΨ â©½ 2 V), the two methods yield essentially identical results for ion and solvent density profiles; however, significant differences appear at higher ΔΨ. At ΔΨ â©¾ 4 V, the CPM ion density profiles show significant enhancement (over FCM) of "inner-sphere adsorbed" Li(+) ions very close to the electrode surface. The ability of the CPM electrode to respond to local charge fluctuations in the electrolyte is seen to significantly lower the energy (and barrier) for the approach of Li(+) ions to the electrode surface.

Thermodynamics and intrinsic structure of the Al-Pb liquid-liquid interface: a molecular dynamics simulation study.

We examine the thermodynamics and intrinsic structure of the Al-Pb liquid-liquid interface using molecular dynamics simulation and embedded atom method potentials. The instantaneous interfacial positions, from which the intrinsic structure and the capillary fluctuation spectrum are determined, are calculated using a grid-based method. The interfacial free energy extracted from the capillary fluctuation spectrum is shown to be in excellent agreement with that calculated mechanically by integrating the stress profile. The intrinsic liquid-liquid interfacial density profile shows structural oscillations in the liquid phases in the interfacial region that are shown to be quantitatively similar to the radial distribution functions of the bulk liquid, consistent with theoretical predictions from classical density functional theory and with earlier simulations on liquid-liquid and liquid-vapor interfaces. In addition, we show the mean interfacial density profile for this system is well described as a convolution of the intrinsic density profile and the probability distribution of interfacial position.