Surface viscosity of liquid interfaces from Green-Kubo relations.

The precise determination of surface transport coefficients at liquid interfaces is critical to an array of processes, ranging from atmospheric chemistry to catalysis. Building on our prior results that highlighted the emergence of a greatly reduced surface viscosity in simple liquids via the dispersion relation of surface excitations [Malgaretti et al., J. Chem. Phys. 158, 114705 (2023)], this work introduces a different approach to directly measure surface viscosity. We use modified Green-Kubo relations suitable for inhomogeneous systems to accurately quantify viscosity contributions from fluid slabs of variable thickness through extensive molecular dynamics simulations. This approach distinguishes the viscosity effects of the surface layer vs the bulk, offering an independent measure of surface viscosity and providing a more detailed understanding of interfacial dynamics and its transport coefficients.

Percolation transition and bimodal density distribution in hydrogen fluoride.

Hydrogen-bond networks in associating fluids can be extremely robust and characterize the topological properties of the liquid phase, as in the case of water, over its whole domain of stability and beyond. Here, we report on molecular dynamics simulations of hydrogen fluoride (HF), one of the strongest hydrogen-bonding molecules. HF has more limited connectivity than water but can still create long, dynamic chains, setting it apart from most other small molecular liquids. Our simulation results provide robust evidence of a second-order percolation transition of HF's hydrogen bond network occurring below the critical point. This behavior is remarkable as it underlines the presence of two different cohesive mechanisms in liquid HF, one at low temperatures characterized by a spanning network of long, entangled hydrogen-bonded polymers, as opposed to short oligomers bound by the dispersion interaction above the percolation threshold. This second-order phase transition underlines the presence of marked structural heterogeneity in the fluid, which we found in the form of two liquid populations with distinct local densities.

Local structure of liquid/vapour interfaces approaching the critical point.

Investigating the structure of fluid interfaces at high temperatures is a particularly delicate task that requires effective ways of discriminating liquid from vapour and identifying the location of the liquid phase boundary, thereby allowing to distinguish intrinsic from capillary fluctuations. Several numerical approaches require introducing a coarse-graining length scale, often heuristically chosen to be the molecular size, to determine the location of the liquid phase boundary. Here, we propose an alternative rationale for choosing this coarse-graining length scale; we require the average position of the local liquid phase dividing surface to match its flat, macroscopic counterpart. We show that this approach provides additional insight into the structure of the liquid/vapour interface, suggesting the presence of another length scale beyond the bulk correlation one that plays an important role in determining the interface structure.

Renormalized charge and dielectric effects in colloidal interactions: a numerical solution of the nonlinear Poisson-Boltzmann equation for unknown boundary conditions.

 The Derjaguin-Landau-Verwey-Overbeek (DLVO) theory, introduced more than 70 years ago, is a hallmark of colloidal particle modeling. For highly charged particles in the dilute regime, it is often supplemented by Alexander's prescription (Alexander et al. in J Chem Phys 80:5776, 1984) for using a renormalized charge. Here, we solve the problem of the interaction between two charged colloids at finite ionic strength, including dielectric mismatch effects, using an efficient numerical scheme to solve the nonlinear Poisson-Boltzmann (NPB) equation with unknown boundary conditions. Our results perfectly match the analytical predictions for the renormalized charge by Trizac and coworkers (Aubouy et al. in J Phys A 36:5835, 2003). Moreover, they allow us to reinterpret previous molecular dynamics (MD) simulation results by Kreer et al. (Phys Rev E 74:021401, 2006), rendering them now in agreement with the expected behavior. We furthermore find that the influence of polarization becomes important only when the Debye layers overlap significantly.

High-density liquid (HDL) adsorption at the supercooled water/vapor interface and its possible relation to the second surface tension inflection point.

We investigate the properties of water along the liquid/vapor coexistence line in the supercooled regime down to the no-man's land. Extensive molecular dynamics simulations of the TIP4P/2005 liquid/vapor interface in the range 198-348 K allow us to locate the second surface tension inflection point with a high accuracy at 283 ± 5 K, close to the temperature of maximum density. This temperature also coincides with the appearance of a density anomaly at the interface known as the apophysis. We relate the emergence of the apophysis to the observation of high-density liquid (HDL) water adsorption in the proximity of the liquid/vapor interface.

Surface viscosity in simple liquids.

The response of Newtonian liquids to small perturbations is usually considered to be fully described by homogeneous transport coefficients like shear and dilatational viscosity. However, the presence of strong density gradients at the liquid/vapor boundary of fluids hints at the possible existence of an inhomogeneous viscosity. Here, we show that a surface viscosity emerges from the collective dynamics of interfacial layers in molecular simulations of simple liquids. We estimate the surface viscosity to be 8-16 times smaller than that of the bulk fluid at the thermodynamic point considered. This result can have important implications for reactions at liquid surfaces in atmospheric chemistry and catalysis.

Capillary interactions between soft capsules protruding through thin fluid films.

When a suspension dries, the suspending fluid evaporates, leaving behind a dry film composed of the suspended particles. During the final stages of drying, the height of the fluid film on the substrate drops below the particle size, inducing local interface deformations that lead to strong capillary interactions among the particles. Although capillary interactions between rigid particles are well studied, much is still to be understood about the behaviour of soft particles and the role of their softness during the final stages of film drying. Here, we use our recently-introduced numerical method that couples a fluid described using the lattice Boltzmann approach to a finite element description of deformable objects to investigate the drying process of a film with suspended soft particles. Our measured menisci deformations and lateral capillary forces, which agree well with previous theoretical and experimental works in case of rigid particles, show that the deformations become smaller with increasing particle softness, resulting in weaker lateral interaction forces. At large interparticle distances, the force approaches that of rigid particles. Finally, we investigate the time dependent formation of particle clusters at the late stages of the film drying.

Ab initio structure and thermodynamics of the RPBE-D3 water/vapor interface by neural-network molecular dynamics.

Aided by a neural network representation of the density functional theory potential energy landscape of water in the Revised Perdew-Burke-Ernzerhof approximation corrected for dispersion, we calculate several structural and thermodynamic properties of its liquid/vapor interface. The neural network speed allows us to bridge the size and time scale gaps required to sample the properties of water along its liquid/vapor coexistence line with unprecedented precision.

PDADMAC/PSS Oligoelectrolyte Multilayers: Internal Structure and Hydration Properties at Early Growth Stages from Atomistic Simulations.

We analyze the internal structure and hydration properties of poly(diallyl dimethyl ammonium chloride)/poly(styrene sulfonate sodium salt) oligoelectrolyte multilayers at early stages of their layer-by-layer growth process. Our study is based on large-scale molecular dynamics simulations with atomistic resolution that we presented recently [Sánchez et al., Soft Matter 2019, 15, 9437], in which we produced the first four deposition cycles of a multilayer obtained by alternate exposure of a flat silica substrate to aqueous electrolyte solutions of such polymers at 0.1M of NaCl. In contrast to any previous work, here we perform a local structural analysis that allows us to determine the dependence of the multilayer properties on the distance to the substrate. We prove that the large accumulation of water and ions next to the substrate observed in previous overall measurements actually decreases the degree of intrinsic charge compensation, but this remains as the main mechanism within the interface region. We show that the range of influence of the substrate reaches approximately 3 nm, whereas the structure of the outer region is rather independent from the position. This detailed characterization is essential for the development of accurate mesoscale models able to reach length and time scales of technological interest.

Atomistic simulation of PDADMAC/PSS oligoelectrolyte multilayers: overall comparison of tri- and tetra-layer systems.

By employing large-scale molecular dynamics simulations of atomistically resolved oligoelectrolytes in aqueous solutions, we study in detail the first four layer-by-layer deposition cycles of an oligoelectrolyte multilayer made of poly(diallyl dimethyl ammonium chloride)/poly(styrene sulfonate sodium salt) (PDADMAC/PSS). The multilayers are grown on a silica substrate in 0.1 M NaCl electrolyte solutions and the swollen structures are then subsequently exposed to varying added salt concentration. We investigated the microscopic properties of the films, analyzing in detail the differences between three- and four-layer systems. Our simulations provide insights into the early stages of growth of a multilayer, which are particularly challenging for experimental observations. We found rather strong complexation of the oligoelectrolytes, with fuzzy layering of the film structure. The main charge compensation mechanism is for all cases intrinsic, whereas extrinsic compensation is relatively enhanced for the layer of the last deposition cycle. In addition, we quantified other fundamental observables of these systems, such as the film thickness, water uptake, and overcharge fractions for each deposition layer.

Phase stability of the ice XVII-based CO2 chiral hydrate from molecular dynamics simulations.

We computed the phase diagram of CO2 hydrates at high pressure (HP), from 0.3 to 20 kbar, by means of molecular dynamics simulations. The two CO2 hydrates known to occur in this pressure range are the cubic structure I (sI) clathrate and the HP hydrate, whose water framework is the recently discovered ice XVII. We investigated the stability of both hydrates upon heating (melting) as well as the phase changes upon compression. The CO2-filled ice XVII is found to be more stable than the sI clathrate and than the mixture of ice VI and dry ice at pressure values ranging from 6 to 18 kbar and in a wide temperature range, although a phenomenological correction suggests that the stability should more realistically range from 6.5 to 13.5 kbar. Our simulation results support the current hypothesis that the HP hydrate is stable at temperatures above the melting curve of ice VI.

Pytim: A python package for the interfacial analysis of molecular simulations.

Pytim is a versatile python framework for the analysis of interfacial properties in molecular simulations. The code implements several algorithms for the identification of instantaneous interfaces of arbitrary shape, and analysis tools written specifically for the study of interfacial properties, such as intrinsic profiles. The code is written in the python language, and makes use of the numpy and scipy packages to deliver high computational performances. Pytim relies on the MDAnalysis library to analyze the trajectory file formats of popular simulation packages such as gromacs, charmm, namd, lammps or Amber, and can be used to steer OpenMM simulations. Pytim can write information about surfaces and surface atomic layers to vtk, cube, and pdb files for easy visualization. The classes of Pytim can be easily customized and extended to include new interfacial algorithms or analysis tools. The code is available as open source and is free of charge. © 2018 The Authors. Journal of Computational Chemistry published by Wiley Periodicals, Inc.

Electrokinetic droplet transport from electroosmosis to electrophoresis.

Droplet transport in microfluidic channels by electrically induced flows often entails the simultaneous presence of electroosmosis and electrophoresis. Here we make use of coupled lattice-Boltzmann/molecular dynamics simulations to compute the mobility of a droplet in a microchannel under the effect of an external electric field. By varying the droplet solvation free energy of the counterions released at the channel walls, we observe the continuous transition between the electroosmotic and electrophoretic regime. We show that it is possible to describe the mobility of a droplet in a unified, consistent way, by combining the theoretical description of the electroosmotic flow with, in this case, the Hückel limit of electrophoresis, modified in order to take into account the Hadamard-Rybczynski droplet drag.

The impact of tensorial temperature on equilibrium thermodynamics.

Thermodynamic temperature is a scalar. However, the connection with the kinetic energy tensor in statistical mechanics leaves open the possibility to define a tensorial temperature. This concept has sometimes been used to simulate isothermal conditions in out-of-equilibrium systems. Here, we show, by studying a sessile water droplet, that a tensorial temperature leads to the wrong thermodynamics, or, in other words, the equilibrium isothermal ensemble generated using a tensorial temperature is not the canonical one, with interfacial free energies that can differ up to 40% from the correct ones.

ForConX: A forcefield conversion tool based on XML.

The force field conversion from one MD program to another one is exhausting and error-prone. Although single conversion tools from one MD program to another exist not every combination and both directions of conversion are available for the favorite MD programs Amber, Charmm, Dl-Poly, Gromacs, and Lammps. We present here a general tool for the force field conversion on the basis of an XML document. The force field is converted to and from this XML structure facilitating the implementation of new MD programs for the conversion. Furthermore, the XML structure is human readable and can be manipulated before continuing the conversion. We report, as testcases, the conversions of topologies for acetonitrile, dimethylformamide, and 1-ethyl-3-methylimidazolium trifluoromethanesulfonate comprising also Urey-Bradley and Ryckaert-Bellemans potentials. © 2017 Wiley Periodicals, Inc.

Phase and interface determination in computer simulations of liquid mixtures with high partial miscibility.

Partially miscible solutions can represent a challenge from the computer simulation standpoint, especially if the mutual solubility of the components is so large that their concentrations do not change much from one phase to another. In this case, identifying which molecules belong to which phase becomes a complicated task. Here, we propose a density-based clustering approach with self-tuning capabilities and apply it to the case of the mixture of an ionic liquid with benzene. The almost linear scaling of the algorithm makes it suitable for the analysis of long Molecular Dynamics or Monte Carlo trajectories.

The role of a small-scale cutoff in determining molecular layers at fluid interfaces.

The existence of molecular layers at liquid/vapour interfaces has been a long debated issue. More than ten years ago it was shown, using computer simulations, that correlations at the liquid/vapour interface resemble those of bulk liquids, even though they can be detected in experiments only in a few cases, where they are so strong that they cannot be concealed by the geometrical smearing of capillary fluctuations. The results of the intrinsic analysis techniques used in computer experiments, however, are still often questioned because of their dependence on a free parameter that usually represents a small-scale cutoff used to determine the interface. In this work I show that there is only one value of the cutoff that can ensure a quantitative explanation of the intrinsic density correlation peaks in terms of successive layer contributions. The value of the cutoff coincides, with a high accuracy, with the molecular diameter.

Kinetic dielectric decrement revisited: phenomenology of finite ion concentrations.

With the help of a recently developed non-equilibrium approach, we investigate the ionic strength dependence of the Hubbard-Onsager dielectric decrement. We compute the depolarization of water molecules caused by the motion of ions in sodium chloride solutions from the dilute regime (0.035 M) up close to the saturation concentration (4.24 M), and find that the kinetic decrement displays a strong non-monotonic behavior, in contrast to the prediction of available models. We introduce a phenomenological modification of the Hubbard-Onsager continuum theory, which takes into account the screening due to the ionic cloud at the mean-field level and, which is able to describe the kinetic decrement at high concentrations including the presence of a pronounced minimum.

The effect of anaesthetics on the properties of a lipid membrane in the biologically relevant phase: a computer simulation study.

Molecular dynamics simulations of the fully hydrated neat dipalmitoylphosphatidylcholine (DPPC) membrane as well as DPPC membranes containing four different general anaesthetic molecules, namely chloroform, halothane, diethyl ether and enflurane, have been simulated at two different pressures, i.e., at 1 bar and 1000 bar, at the temperature of 310 K. At this temperature the model used in this study is known to be in the biologically most relevant liquid crystalline (Lα) phase. To find out which properties of the membrane might possibly be related to the molecular mechanism of anaesthesia, we have been looking for properties that change in the same way in the presence of any general anaesthetic molecule, and change in the opposite way by the increase of pressure. This way, we have ruled out the density distribution of various groups along the membrane normal axis, orientation of the lipid heads and tails, self-association of the anaesthetics, as well as the local order of the lipid tails as possible molecular reasons of anaesthesia. On the other hand, we have found that the molecular surface area, and hence also the molecular volume of the membrane, is increased by the presence of any anaesthetic molecule, and decreased by the pressure, in accordance with the more than half a century old critical volume hypothesis. We have also found that anaesthetic molecules prefer two different positions along the membrane normal axis, namely the middle of the membrane and the outer edge of the hydrocarbon region, close to the polar headgroups. The increase of pressure is found to decrease the former, and increase the latter preference, and hence it might also be related to the pressure reversal of anaesthesia.

Layer-by-layer and intrinsic analysis of molecular and thermodynamic properties across soft interfaces.

Interfaces are ubiquitous objects, whose thermodynamic behavior we only recently started to understand at the microscopic detail. Here, we borrow concepts from the techniques of surface identification and intrinsic analysis, to provide a complementary point of view on the density, stress, energy, and free energy distribution across liquid ("soft") interfaces by analyzing the respective contributions coming from successive layers.