Microphysics of liquid water in sub-10 nm ultrafine aerosol particles

. Ultrafine aerosol particles with sizes smaller than 50 nm have been shown in recent studies to serve as a large source of cloud condensation nuclei (CCN) that can promote additional cloud droplet formation under supersaturation conditions. Knowledge of the microphysics of liquid water in these droplets remains limited, particularly in the sub-10 nm particle size range, due to experimental and theoretical challenges associated with the complexity of aerosol components and the small length scales of interest (e.g., difficulty of precisely sampling the 5 liquid-air interface, questionable validity of mean-field theoretical representations). Here, we carried out molecular dynamics (MD) simulations of aerosol particles with diameters between 1 and 10 nm and characterized atomistic-level structure and water dynamics in well-mixed and phase-separated system with different particle sizes, NaCl salinities, and organic surface loadings as a function of distance from the time-averaged Gibbs dividing interface or instantaneous water-air interface. We define a sphericity factor ( ϕ ) that can shed light on the phase-mixing state of 10 nanodroplets, and we reveal an unexpected dependence of mixing state on droplet size. Our results also evidence an ion concentration enhancement in ultrafine aerosols, which should modulate salt nucleation kinetics in ultrafine sea salt aerosols, and provide detailed characterization of the influence of droplet size on surface tension and on water self-diffusivity near the interface. Analysis of water evaporation free energy and water activity demonstrates the validity of the Kelvin equation and Köhler theory at droplet sizes larger than 4 nm under moderate salinities 15 and organic loadings and the need for further extension to account for ion concentration enhancement in sub-10 nm aerosols, droplet-size-dependent phase separation effects, and a sharp decrease in the cohesiveness of liquid water in sub-4 nm droplets. Finally, we show that an idealized fractional surface coating factor ( f s ) can be used to categorize and reconcile water accommodation coefficients ( α ∗ ) observed in MD simulations and experimental results in the presence of organic coatings, and we resolve the droplet-size dependence of α ∗


Introduction
Water is a predominant component of atmospheric aerosol, with mass fraction typically larger than 70% at moderate to high relative humidities (RH > 40%), and plays crucial roles in the formation and growth of aerosol particles (Pierce et al., 2012;Karlsson et al., 2022). For example, aerosol liquid water facilitates the partitioning of gas-phase bulk liquid water.
Molecular dynamics (MD) simulations have been extensively used to gain insight into the properties of liquid water near interfaces at the atomistic level. In the case of liquid water-air interfaces, most of these studies have focused on flat interfaces relevant to large aerosol particles (Brown et al., 2005;Thomas et al., 2007;D'Auria and Tobias, 2009;Tang et al., 2020;Luo et al., 2020). Comparatively few have examined the properties of liquid water 65 in sub-10 nm aerosol particles (Chowdhary and Ladanyi, 2009;Li et al., 2011;Ma et al., 2011;Lovrić et al., 2016;Sun et al., 2013;Vardanega and Picaud, 2014;Radola et al., 2017;Karadima et al., 2019;Li et al., 2013;Radola et al., 2015Radola et al., , 2021. These simulations have shown that interfacial water adopts distinct properties from those of bulk liquid water and that the surface tension of water droplets decreases with increasing interfacial curvature at the nanometer scale (Bourg and Steefel, 2012;Li et al., 2010). In addition, interfacial water molecules have been shown 70 to have highly distinct dynamics compared to bulk liquid water, including slower hydrogen bond dynamics and faster diffusion dynamics (Ni et al., 2013;Liu et al., 2004;Wick and Dang, 2005), but the apparent discrepancy between these observations and their potential implications have not been resolved. Finally, a complicating factor revealed by these simulations and by experimental observations is that the shape of nanodroplets often exhibits significant fluctuations and is rarely a perfect sphere, which complicates the determination of interface thickness and water 75 uptake properties (Li et al., 2016;Garrett et al., 2006).
Here, we use MD simulations to characterize nano-droplet structure, water diffusion, and water evaporation energetics in well-mixed and phase separated system with different sizes, NaCl salinities, and organic surface loadings.
Our approach builds on previous MD simulation studies of nano-aerosol droplets but includes two distinct features that have the potential to help resolve the challenges identified in previous studies. First, we use two methods 80 of interface characterization -the time-averaged Gibbs dividing interface (Lau et al., 2015) and the instantaneous interface (Willard and Chandler, 2010) -to compare different perspectives on the properties of the interface. Second, we characterize a broad range of structural, dynamic, and energetic properties to help reveal links between these properties. We re-examine the validity of the Kelvin equation and Köhler theory and, more broadly, of treating nano-droplet water as a bulk aqueous phase. We find that the instantaneous interface characterization enables a 85 novel characterization of the droplet phase state and helps resolve apparent contradictions between experimental observation and MD simulation results on interfacial water dynamics. We also find that the Kelvin equation and Köhler theory accurately describe water microphysics under moderate salinity and organic loadings in nanometerscale droplets. However, these theories underestimate water activity at droplet diameters below 4 nm and cannot account for the droplet size effect on ion concentration enhancement in NaCl-bearing droplets and on phase separation 90 in droplets bearing surface-active organic compounds. Overall, our results indicate that the properties of nano-aerosol water deviate strongly from those of bulk water in particles smaller than 4 nm. One implication of this finding is that molecular-scale insight is required to accurately represent the initial stages of droplet nucleation and growth. 95 Molecular dynamics simulations of droplets containing 100, 500, 1000, 3000, or 5000 water molecules (N w ) with NaCl molality of 0, 1, 2, or 4 m were used to examine the effect of droplet size and salinity on water microphysics in nano-aerosol particles. In addition, droplets with N w = 100, 500, 1000, 3000, or 5000 with low or high organic surface loadings (0.7 or 3.3 pimelic acid molecules nm −2 , corresponding to surface loadings of 0.3 or 1.5 monolayers) were used to examine the impact of organic compounds. Sodium chloride is used as a representative inorganic 100 salt in aerosol particles for two reasons: first, it accounts for a significant fraction of sea-spray and mineral dust aerosols, which constitute the largest mass of particulate matter emitted to the atmosphere (Lewis and Schwartz, 2004;Wise et al., 2007) and, second, relative to other inorganic ions present in aerosol particles, such as NH + 4 and SO 2-4 , alkali halide salts such as NaCl have comparatively well-tested and accurate interatomic interaction potential models, the main input controlling the accuracy of MD simulation predictions (Smith and Dang, 1994;Joung and 105 Cheatham III, 2008III, , 2009. Pimelic acid (PML, C 7 H 12 O 4 ) was used as a representative organic component, because dicarboxylic acids contribute a large fraction of the total identifiable resolved organic mass in fine aerosols and have been identified and quantified in both rural and urban air (Hyder et al., 2012;Kawamura and Yasui, 2005;Tran et al., 2000). Furthermore, the O/C ratio of PML (0.57) lies near the midpoint of the range commonly observed for aerosol organic materials (0.2 to 1.0) (Zhang et al., 2007;Hallquist et al., 2009;Song et al., 2018) and near the value below 110 which liquid-liquid phase separation is commonly observed in aerosol particles (∼ 0.7 to 0.8 for organic-salt-water aerosols, ∼ 0.6 for organic-water aerosols) (Ciobanu et al., 2009;Song et al., 2012;You et al., 2014;Renbaum-Wolff et al., 2016;Rastak et al., 2017;Song et al., 2018). These observations suggest that PML should have the potential to mimic key properties of organic substances in nano-aerosol droplets. The sizes of droplets studied in this work range from 2 to 8 nm (Appendix Table A1).

System preparation and MD simulations
Molecular dynamics simulations were performed with the Gromacs program (Abraham et al., 2015) using cubic simulation cells subject to periodic boundary conditions in all three directions. The edge length of the simulation cell is set to be at least 6 nm larger than the droplet diameter to avoid interactions between the droplet and its periodic images (Appendix Figure A1). Interatomic interactions were represented using the SPC/E model for water 120 molecules (Berendsen et al., 1987), the Joung-Cheatham model for Na + and Clions (Joung and Cheatham III, 2008), and the OPLS-AA model for pimelic acid molecules (Jorgensen et al., 1996). In all simulations, bonds involving H atoms were restrained using the LINCS algorithm (Hess, 2008). Simulations were carried out in the NVT ensemble with temperature maintained at 298.15 K using the Nosé-Hoover thermostat (Evans and Holian, 1985) and with the overall translational and rotational momentum of the simulated system set to zero. Coulomb and van der Waals 125 interactions were resolved using a real space cutoff of 1.2 nm with a particle mesh Ewald sum treatment of long-range comprised of 5 ns for equilibration and 100 ns for data analysis. We note that the Ewald sum treatment of long-range Coulomb interactions used in our study is commonly used in MD simulations, including in simulations of aerosol particles, but care should be taken that it introduces artifacts that remain incompletely examined, particularly in 130 the dynamics of charged and dipolar species (Frenkel and Smit, 2001;Hub et al., 2014).

Droplet structure
To calculate the radius of the water droplets in most systems (except at high organic loading) the time-averaged water number density profile ρ(r) was obtained for each droplet as a function of distance from the droplet center of mass. Results were fitted to the relation (Thompson et al., 1984), where ρ l is the liquid number density, ρ v is the vapor number density, and d i is the interfacial width. The density profile was used to locate the equimolar (Gibbs) dividing surface, which corresponds to a vanishing adsorption in a one-component system. The radius R e of the equimolar dividing surface was calculated using numerical quadrature, For droplets with a high organic loading, the radius was calculated instead based on the molecular volumes of water and pimelic acid: where N w is the number of water molecules, N org is the number of organic molecules, ρ w is the bulk water density (34.02 ± 0.19 molecules nm −3 ) and ρ org is the bulk organic density (4.373±0.003 molecules nm −3 ) calculated in the 145 conditions of our simulations.

Diffusion coefficient
The self-diffusion coefficient of water molecules was calculated by analyzing water trajectories at 1 ps interval during 25 ns simulations for all droplets. To analyze the dependence of diffusion on position relative to the interface, we binned water molecules into different 3 Å thick layers based on their distance relative to the interface. Any molecule 150 initially present in a specified layer that passed through the layer boundaries no longer contributed to the calculation.
The self-diffusion of water molecules was investigated in directions parallel and normal to the interface.
The parallel diffusion coefficient (D ∥ ) was calculated using the Einstein relation (Liu et al., 2004) where l 2 ∥,i is the mean-square displacement (MSD) of water molecules in directions parallel to the interface, averaged 155 over all molecules of interest remaining in i th layer during the entire time interval τ . The infinite-time limit in Eq. 4 was approximated by evaluating the slope of MSD vs. τ for τ = 2 to 6 ps based on the linear relationship between MSD and τ on this time scale; use of longer time-scales yielded less precise results as fewer molecules remained consistently located in the i th layer (Appendix Figure A6).
The diffusion coefficient normal to the interface (D ⊥ ) was evaluated using the anisotropic Smoluchowski equation 160 in the r direction: where p = p(r, t|r 0 , t 0 ) is the conditional probability distribution function, D ⊥ = D rr is the diffusion tensor perpen- ∂r is the mean force acting on the molecule, and F (r) = −k B T ln ρ(r) ρ l is the potential of mean force (PMF) calculated from the water 165 density profile. D ⊥ values in different water layers were determined by matching the water survival probabilities calculated from MD simulation and from the numerical solution to the Smoluchowski equation as described by Wick and Dang (2005).

Surface tension
Surface tension was calculated by analyzing the normal component of the Irving-Kirkwood pressure tensor (Irving 170 and Kirkwood, 1950) at 1 ps intervals during 25 ns simulations for all droplets. According to the formula given by Thompson et al. (1984), the normal pressure tensor P N (r) can be written as where P K (r) = k B T ρ(r) is the kinetic term, P U (r) = S −1 k f k is the configuration term, S is the surface area of the spherical surface of radius r, and f k is the normal component of the force between a pair of molecules acting 175 across the surface S. The surface tension is evaluated as where W is the work of formation as described by Zakharov et al. (1997): In Eq. 8, P β is the vapor pressure of the droplet, and R β is the radius of a sphere in the vapor region, which can be 180 approximated as the equimolar radius R e . We note that the surface tension should be rigorously calculated at the radius R s where the surface tension applies (McGraw and Laaksonen, 1997). Our use of R e in Eq. 7 might introduce systematic error of σ calculation, i.e. less than 8% for droplet with N w ≥ 500 (see Appendix A for more details).
This potential source of systematic error is commensurate with the statistical error of our predicted σ values (as shown by our reported error bars of ∼ ±8%).

Free energy of evaporation
The potential of mean force (PMF) associated with the evaporation of water molecules as a function of the distance r from the droplet center, F (r), was calculated using two methods. The first method consisted in calculating the PMF based on the water density profile during 100 ns equilibrium simulations by The second method is the well-known umbrella sampling technique (Torrie and Valleau, 1977) and the F (r) is calculated by successively tethering a single water molecule at different distances from the droplet center of mass using a harmonic potential with a force constant of 1000 kJ mol −1 nm −1 and 100 ns of simulation time at each tethering distance. The evaporation free energy of water molecules was calculated as the free energy difference of water molecules in the vapor phase and the liquid phase as

Accommodation coefficient
To calculate the accommodation coefficient for a single droplet, the phase-state of every water molecule in the system was recorded during a 100 ns trajectory with a timestep of 1 fs. The phase state of individual water molecules was divided into three categories (Varilly and Chandler, 2013): (1) liquid phase (state=-1): any water molecule that is 200 hydrogen-bonded to at least one water molecule and is at most 4 Å away from the nearest water molecules (where the "position" of a water molecule refers to the position of its center of mass unless otherwise stated); (2) vapor phase (state=1): any water molecule that has no hydrogen-bonded water neighbors and is more than 8 Å away from its nearest water neighbor; (3) intermediate phase (state=0): any water molecule that is neither in liquid nor vapor phase. An accommodation event is identified if the state of a water molecule changes from 1 at timestep t 1 to -1 at 205 timestep t 2 , without adopting a non-zero state value during the time interval (t 1 , t 2 ). A reflection event is identified if the state of a single water molecule changes from 1 at t 1 to 1 at t 2 , and stays in state 0 during the time interval (t 1 , t 2 ). The accommodation coefficient α was calculated as where N acc is the number of accommodation events and N ref is the number of reflection events. The reported 210 normalized accommodation coefficient α * was calculated as α * = α/α 0,w , where α 0,w (the value for the flat surface of pure water) was approximated as the α value of the largest simulated pure water droplet with N w = 5000. Insets (a)-(e) in Figure 1 show representative snapshots of simulated droplets that illustrate the variety of predicted 215 morphologies and phase mixing states. Additional droplet snapshots are presented in Appendix Figure A3. For NaCl droplets, Na + and Clions are repelled from the water-air interface and distributed inside the droplets in all simulations. These ions are crystallized in the center of the droplet at 4 m concentration, forming a phase state illustrated in Figure 1(e), whereas they are well-mixed with water molecules at 1 m and 2 m concentrations with phase state illustrated in Figure 1(a). The NaCl solubility implied by our results (between 2 and 4 m) is consistent 220 with that observed in previous MD studies that examined bulk NaCl solutions using the same interatomic potential (Mester and Panagiotopoulos, 2015a, b). For droplets with PML, the arrangement of organic molecules varies with droplet size and organic surface loading. At low organic loading, PML molecules are predominantly adsorbed at the water-air interface with minimal clustering, as illustrated in Figure 1(b). At high organic loading, when droplet size is the droplet surface, forming the 'engulfed' morphology shown in Figure 1(c). In smaller droplets with a high organic loading, however, the organic molecules cluster in the center of the droplet and displace water molecules towards the surface as shown in Figure 1(d). This transition suggests the existence of a threshold, at a droplet diameter of 3 to 4 nm, in the competition between water hydrogen bonding, which promotes water clustering, and hydrophobic attraction between PML molecules, which promotes organic clustering. Our results are consistent with previous 230 findings that the hydrogen bond network formed by bulk liquid water involves cooperativity over length scales of 1 to 2 nm and becomes significantly weaker upon attenuation of long-range water-water orientational correlations (Alper and Levy, 1989;Ohmine and Tanaka, 1993).

Morphology and phase state characterization
To characterize the various phase states identified above, we propose a quantitative variable, referred to hereafter as the sphericity factor (ϕ), defined as the ratio of the surface area of a perfect spherical water droplet, S w,sp , to the 235 instantaneous interface surface area formed by water molecules with the same amount of water, S w,inst : The ϕ-value is in the range of (0, 1] and reflects the tendency towards water clustering. In particular, ϕ = 1 indicates that water forms a perfect spherical droplet with no capillary waves or disruption due to phase separation.  Figure 1) distinguishes between well-mixed (above the threshold) and phase separated droplets (below the threshold). For well-mixed systems, at fixed ion concentration or organic surface excess, ϕ decreases slightly with increasing droplet size, indicating that 245 larger droplets have a more fluctuating water-air interface. This observation is consistent with our results showing that the water-air interface width observed in time-averaged density profiles increases with droplet size (Appendix Table A1) and agrees with the expectation that the capillary wave contribution to interfacial fluctuation increases with interfacial area (Ismail et al., 2006;Lau et al., 2015). We also find that ϕ increases with ion concentration but decreases with organics at the interface, which indicates that NaCl reduces, but PML enhances, interfacial 250 fluctuations. This observation is consistent with the variations in interfacial width with salinity and organic surface loading (Appendix Table A1) and with the expectation that NaCl and PML should respectively increase or decrease surface tension because of their negative (NaCl) or positive (organic) adsorption at the water-air interface, as confirmed below.
Phase-separated systems (i.e., droplets with either 4 m NaCl or a high organic loading) also show a small decrease 255 in ϕ with increasing droplet size in most cases. A significant exception is presented by the smallest droplet with a high organic loading, which shows a sharp decrease in ϕ with decreasing droplet size between D = 4.1 nm and 2.7 nm. This decrease in ϕ coincides with a change in droplet phase-mixing state from a water-core structure at D ≥ 4.1 nm to an organic-core structure at D = 2.7 nm and shows that the phase-mixing state of nanodroplets can abruptly vary with particle size. Organic particles with sizes below 5 nm and different organic-water ratios, which lie in the 260 above mentioned phase-mixing transition region and may shed light on the early stages of hygroscopic growth of secondary organic aerosols , will be examined in more detail in an upcoming study.

Interfacial structure of nanodroplets
To gain additional insight into the phase state behaviors outlined above, the radial density profiles of water, ions, and PML molecules in our nano-droplets were calculated relative to either the time-averaged interface (Gibbs dividing 265 surface) (Lau et al., 2015) or the instantaneous interface (Willard and Chandler, 2010). As shown in previous studies of flat water-air interfaces (Willard and Chandler, 2010) and in Appendix Figure A4, the atomic density profile obtained under the instantaneous interface scheme provides more detailed information about molecular structure at the interface (e.g., water and solute density layering). Figure 2(a) shows ion density relative to the instantaneous interface in 1 m NaCl droplets with different numbers of water molecules. Our results show that both Na + and Cl -

270
ions are negatively adsorbed at the water-air interface, but Na + is more strongly excluded from the interface than Cl -, such that the interface is negatively charged, in agreement with previous studies of flat water-air interfaces (Jungwirth and Tobias, 2000;Underwood and Greenwell, 2018) (see Appendix Figure A4 for more density profiles). Although all droplets highlighted in Figure 2(a) were initialized as 1 m uniformly mixed NaCl solutions, ion concentrations are strongly enhanced (up to >2 times the initial concentration) in the core of the droplets at equilibrium due to the 275 negative adsorption of ions at the water-air interface. This ion concentration enhancement increases with decreasing droplet size, particularly in droplets smaller than 4 nm (N w < 1000). We note that the interface propensity of halide anions varies with anion size and polarizability and that the simple interatomic potential models used here may underestimate the affinity of halide anions for the interface (Levin et al., 2009;Jungwirth and Tobias, 2006;Caleman et al., 2011). However, the negative adsorption of alkali cations predicted here is well-established (Underwood and 280 Greenwell, 2018; Levin et al., 2009;Jungwirth and Tobias, 2006;Caleman et al., 2011), and the overall negative adsorption of NaCl salt at the water-air interface is correctly predicted as shown by previous MD simulation results on the salinity-dependence of water-air interfacial tension Tobias, 2001, 2006;D'Auria and Tobias, 2009;Horinek et al., 2009;Sun et al., 2012).
To illustrate the expected effect of droplet size on ion concentration, we define a salt concentration enhancement 285 factor (ϵ) as the ratio between salt concentration in the core of the droplet and average salt concentration in the entire droplet. Predicted ϵ-values for different alkali-halide salts are shown in Figure 2(b) based on the simplifying assumption that the salt is excluded from the vicinity of the water surface, in a region with a thickness equal to 1.2 times the first hydration radius of alkali cations (see Appendix B for sensitivity analysis), quantified as the first minimum of the ion-water radial distribution function predicted using the SPC/E water model and compatible ions 290 models (the first hydration shell radii of alkali cations Li + , Na + , and K + are 0.279, 0.318, and 0.359 nm, respectively) (Joung and Cheatham III, 2009). Our simplistic model, shown as solid curves in Figure 2(b) for different salts, is triangles are previous MD simulation results for NaCl (Bahadur and Russell, 2008); and the orange dashed line is the NaCl DRH prediction based on a thermodynamic model (Russell and Ming, 2002). (d) Atomic density profiles of water molecules (ρw) and PML molecules (ρorg) relative to the Gibbs dividing surface under the time-averaged interface scheme at low organic loading.
consistent with our MD simulation predictions (blue circles, shown for our systems with 1 m NaCl), indicating that the length scale associated with salt exclusion at the water-air interface is invariant with surface curvature. Horizontal lines in Figure 2(b) show the ϵ-values necessary to achieve crystallization of different salts based on their experimental 295 solubility (Joung and Cheatham III, 2009) in droplets where the average concentration is 1 m. Overall, both our simplistic model and our detailed atomistic simulation predictions indicate that the concentration enhancement of simple alkali-halide salts is well above 1.25 for sub-10 nm droplets, and may be more than 10 for droplets with a diameter below 1.5 nm. In other words, in droplets with common sea-salt ionic species, crystallization should occur more readily than expected based on the average salt concentration, particularly at small droplet sizes. Our results 300 suggest that the tendency for salt deliquescence should be diminished and salt efflorescence should be enhanced in ultrafine nanoparticles. The ion concentration enhancement discussed here is consistent with experimental findings that both the deliquescence and efflorescence relative humidities (DRH, ERH) of NaCl and KCl particles increase as droplet size decreases from 60 nm to less 10 nm, with a sharper decrease at smaller droplet sizes (Biskos et al., 2006a, b;Giamarelou et al., 2018;Hämeri et al., 2001). Shaded areas in Figure 2(c) show the KCl and NaCl DRH as 305 a function of nanoparticle diameter predicted from the ion concentration enhancement model discussed above (see Appendix C for calculation details). Our predictions are consistent with both the size dependence of NaCl DRH reported in previous studies (Biskos et al., 2006a;Bahadur and Russell, 2008;Russell and Ming, 2002) and with previous observations that KCl has a larger DRH than NaCl (Giamarelou et al., 2018). Overall, the concentration enhancement discussed above may help explain experimental observations of a size-dependence of aerosol shape, in 310 which nano-aerosol particles generated from well-controlled sea spray chambers are more cubic at smaller sizes but more spherical at larger sizes (Zieger et al., 2018). Our observations indicate that the ϵ factor becomes particularly significant for sub-10 nm particles and should be taken into consideration in predictions of salt microphyics in ultra-fine particles. sharp peak near the water-air interface, confirming the high affinity of the organic compound for the interface. The location of the PML density peak shifts into the liquid water region with decreasing droplet size, particularly at sizes below 3 nm (N w < 500). This observation suggests that the hydrophilic character of the water phase decreases with decreasing number of water molecules in the droplet, as already noted above in the case of the smallest droplet 320 at high organic loading. The PML density curves at high organic loading (Appendix Figure A5(d)) reveal a more complex and variable pattern of organic-water mixing phase states, in agreement with the transition between waterand organic-core morphologies identified above.

Water diffusion kinetics in nanodroplets
Water diffusion in nanodroplets is a fundamental phenomenon that influences interfacial mass fluxes and chemical 325 transformation rates (Chan et al., 2014;Slade and Knopf, 2014;Julin et al., 2013;Evoy et al., 2020). Details of this phenomenon are complicated by the inherently anisotropic and non-uniform nature of water diffusion near the interface. In this study, the water self-diffusion coefficient is calculated in discrete 3 Å thick water layers in directions parallel (D ∥ ) or perpendicular (D ⊥ ) to the interface. As explained in the Methods section, we calculated D ∥ using the Einstein relation and D ⊥ using the anisotropic Smoluchowski equation (Liu et al., 2004). Diffusion coefficients 330 were calculated relative to the center of mass of the droplet ('com' scheme) or to the instantaneous interface ('inst' scheme). As expected, both schemes yield identical results for D ∥ , while for D ⊥ the 'com' scheme yields larger values due to the presence of capillary waves at the interface (Appendix D2.1). The following discussion focuses primarily on our predicted D ∥ values. Our D ⊥ values are somewhat less precise and are discussed in the supplementary information (Appendix D2.2).

335
Figure 3(a) shows D ∥ in different simulated systems relative to the distance from the Gibbs dividing surface. Our results show that D ∥ values are bulk-liquid-like at distances > 9 Å (i.e. > 3 water layers) from the interface for all well-mixed systems. The D ∥ values at the interface are 2 to 3 times larger than in the core of the droplet, with potential implications for mass transfer and reaction kinetics at the interface. In addition, D ∥ in the droplet core decreases from 0.25 Å 2 ps −1 to 0.15 Å 2 ps −1 as ion concentration increases from 0 m to 2 m, while D ∥ at the interface decreases with increasing organic surface loading. The D ∥ difference between interfacial layer and bulk liquid region, and the D ∥ change with salinity and organic addition observed in our results are consistent with the trend reported in previous studies (Kim et al., 2012;Paul and Chandra, 2005;Tsimpanogiannis et al., 2019) as shown in Figure   3(a). In all well-mixed NaCl droplets (1, 2 m), we consistently observe a droplet size dependence of water diffusion near the interface as illustrated in Figure 3(b). More precisely, D ∥ at the interface decreases with droplet size, while D ∥ in the few monolayers below the interface increases with droplet size. This observation is consistent with our results showing that for droplets containing water and NaCl, smaller droplets have lower ion density at the interface and greater ion density in the droplet core. The minor offset between our D ∥ values in the bulk-liquid-like region and those reported by Kim et al. (2012) at the same salinity (calculated in bulk liquid NaCl solutions) as shown in i.e., water diffusivity normal to the interface is smaller than in the droplet core, in agreement with studies showing slower reorientational dynamics of water near hydrophobic surfaces (Laage et al., 2009;Fayer et al., 2010;Fayer, 2012).

Surface tension of nanodroplets
Surface tension is an important factor that governs the hygroscopic growth and the activation of aerosol particles and is incorporated in climate models through Köhler theory. However, the surface tension of nano-aerosol particles, especially at the sub-10 nm scale, is not well known and remains arduous to characterize experimentally (Noziere et al., 2014;Bzdek et al., 2016Bzdek et al., , 2020bLee and Tivanski, 2021). Appendix Figure A12 shows surface tensions (σ) 360 calculated for our different droplets based on the Irving-Kirkwood pressure tensor (Irving and Kirkwood, 1950;Thompson et al., 1984). A representative example of the instantaneous normal pressure tensor is shown in Appendix Figure A11. All calculated σ values are also tabulated in Appendix Table A1. We note that σ values reported in different MD simulation studies differ significantly even for pure liquid water system (by up to ≈ 20 mN m −1 ), depending on the interatomic potential models and surface tension determination methods employed (Bahadur et al.,365 2007; Bahadur and Russell, 2008;Li et al., 2010;Lau et al., 2015). For well-mixed NaCl droplets with 0, 1, and 2 m NaCl concentrations, our results show that surface tension increases with droplet size and ion concentration and approaches to 55, 56, and 58 mN m −1 respectively when N w ≥ 3000 (i.e., D > 5 nm). The predicted increase in surface tension with NaCl salinity is consistent with experimental values reported for flat water-air interfaces (1.3 to 1.8 mN m −1 per 1 m of salinity) (Hård and Johansson, 1977;Washburn et al., 1926;Abramzon and Gaukhberg, 370 1993;Tuckermann, 2007).  (Ismail et al., 2006). The solid light blue line in Figure 4 shows the fit of our σ/σ 0 results for 0 m droplets to the 375 second-order Helfrich expansion, i.e. σ/σ 0 = 1 − 2δ 0 /r + c/r 2 , where δ 0 is the Tolman length, r is the droplet radius, and c is the rigidity coefficient (Wilhelmsen et al., 2015). The best fit is achieved using δ 0 = 0.06 nm (95% confidence interval [0.032, 0.094] nm) and c = −0.025 nm 2 (95% confidence interval [-0.096, 0.045] nm 2 ). The positive value of δ 0 in this study, which indicates that σ increases with droplet size, agrees with the trends observed in other studies as shown in Figure 4. In addition, although the absolute σ values reporting in different studies differ depending on the 380 choice of simulation methods, σ/σ 0 data from previous studies (Bahadur and Russell, 2008;Li et al., 2010;Lau et al., 2015) agree well with our simulation results and collapse near the light blue fitting curve. This observation shows that, for pure water droplets, σ decreases sharply as droplet size decreases below 4 nm, and can even be less than half of σ 0 at sizes below 1 nm. For droplets with an organic surface loading, surface tension shows the same trend with particle size as for pure water droplet, but also decreases with organic surface loading for all systems where the 385 droplet has a water core (N w ≥ 500), in agreement with experiment observations that water surface tension decreases with organic adsorption at the water-air interface (Tuckermann, 2007;Kiss et al., 2005). An exception is observed in Figure 4 in the case of the smallest droplets with low and high organic loadings (left-most yellow and orange circles in Figure 4), which both contain the same amount of water (N w = 100), yet the droplet with high organic loading has a higher surface tension. The surface tension of the organic-core configuration (i.e., the smallest droplet with 390 a high organic loading) is close to that observed for pure droplets of PML (on the order of 30 mN m −1 at D ≈ 3 nm), suggesting that the unexpected increase in surface tension with organic loading is associated with the observed morphology reversal from a water core to an organic core. In short, the common assumption in aerosol microphysics studies that surface tension invariably decreases with increasing organic loading in aqueous aerosols appears invalid, at least in some cases, for droplets smaller than 4 nm. This finding may have important implications, in particular,

395
in understanding the initial nucleation and growth of secondary organic aerosols (Enghoff and Svensmark, 2008;Kürten et al., 2018;Guo et al., 2020). The evaporation free energy δF , which is the free energy difference associated with transferring a single water molecules from the droplet phase to the vapor phase, is an important factor that governs water evaporation kinetics 400 through the relation (Seinfeld and Pandis, 2008):

Validity of the Kelvin equation
where j is the evaporation rate, c is a constant factor, R is the ideal gas constant, and T is absolute temperature.
From classical thermodynamic Kelvin theory, the water evaporation free energy for a droplet with radius r, δF (r), relates to the value for a bulk liquid system with a flat interface and with the same aqueous chemistry, δF b , through 405 the relation where P v (r) is the saturated vapor pressure of the droplet, P v,b is the saturated vapor pressure of the For each simulated system, we calculated δF using two different approaches. The first method consisted in determining the free energy of water evaporation using the umbrella sampling method (Torrie and Valleau, 1977), a biased MD simulation method in which a single water molecule is progressively steered from the particle to the 415 vapor phase. The second method consisted in determining the equilibrium water density in the liquid and vapor phases, ρ l and ρ v , from which the evaporation free energy was obtained as δF = −RT ln ρ v /ρ l . The second method was applied only for systems where a well-defined aqueous phase (and corresponding ρ l value) could be identified.
From Figure 5(a), we observe that the umbrella sampling and density profile methods yield consistent δF values for all well-mixed systems (0 m, 1 m, 2 m and low organic loading). We note that there are limited δF literature 420 data even for pure water system. Even for MD simulations of the flat water-air interface carried out with the same water model as in the present study (SPC/E), δF values reported in previous studies vary significantly, as shown in Figure 5(a) (Garrett et al., 2006;Li et al., 2018;Varilly and Chandler, 2013). The comprehensive δF dataset reported in our study, with longer simulation times, multiple analysis approaches, and a range of droplet sizes and water chemistry conditions, provide a useful guide for future examinations of the evaporation free energy of water 425 molecules from aerosol particles.
Calculated δF values in our study increase with droplet size (i.e., water evaporates more readily from smaller droplets) as expected from Eq. 14. Although the impact of ion concentration on evaporation energy cannot be precisely resolved, a linear relationship between δF and 2σ/D is clearly observed in the case of well-mixed droplets, in agreement with the Kelvin equation. The theoretical slope for pure water droplets according to the Kelvin equation 430 prediction (k = −36 × 10 −6 m 3 mol −1 , calculated using the molar volume of liquid water predicted by the SPC/E water model) is indicated by the black line in Figure 5(a). This theoretical slope lies within the 95% confidence interval of the slope obtained by linear regression of the data for well-mixed NaCl droplets with sizes larger than 3.5 nm (N w ≥ 1000) (blue solid line in Figure 5(a) with slope k f = −41 ± 13 × 10 −6 m 3 mol −1 ). However, if the evaporation data obtained for droplets smaller than 3 nm (N w ≤ 500) are included in the regression, the best fitting 435 line is deflected downward and k falls outside of k f 's 95% confidence interval. Our observation that the Kelvin equation is valid at droplet sizes larger than 4 nm under moderate salinities, has also been revealed by experimental studies (Winkler et al., 2012).
For droplets with a low organic loading, the trend shown in Figure 5(a) is generally consistent with that observed for well-mixed water and NaCl droplets. However, the smallest low organic loading droplet represents a clear outlier 440 in Figure 5(a), as also observed in Figure 4, suggesting a possible link between its low surface tension and its low δF value. For the 4 m NaCl droplet system, the trend in δF (dark purple circles) is similar to that observed for the well-mixed NaCl systems, suggesting that although the ions in center of the droplet form a solid, the water surface is similar to that of the well-mixed ion-water systems. In contrast, for the phase separated high organic loading system (orange circles), predicted δF values are consistently much smaller than those of the well-mixed systems.

445
This indicates that at high organic loading, the organic coating strongly enhances water evaporation and invalidates Eq. 14. For the smallest droplet with a high organic loading, δF is only 24 kJ mol −1 , deviating significantly even from the values obtained for larger droplets with high organic surface loading, an observation that may be related to the reversal from water-core to organic-core morphology identified above.
To conclude, the Kelvin equation accurately predicts water evaporation energetics for salty NaCl droplets with 450 sizes larger than 3.5 nm, even at salinities up to 4 m where salt precipitation is observed. However, it overestimates the stability of droplets smaller than 3 nm and of droplets with a high organic surface loading, in agreement with the expectation that water becomes a somewhat distinct (less hydrophilic) fluid when the collective nature of the water hydrogen bond network is disrupted.

455
Water activity a w , which is the ratio of the partial pressure of water in solution and the saturated water vapor pressure of bulk pure water, is an important factor that governs the hygroscopic growth and mass concentration of aerosol particles at a given relative humidity (Seinfeld and Pandis, 2008). According to Köhler theory (Seinfeld and Pandis, 2008), the water activity of an aerosol droplet can be expressed as where x w is the effective mole fraction of water molecules in the solution, which accounts for the effect of solutes as represented using Raoult's law, and P v,0 is the saturated vapor pressure of pure water system with a flat interface.
To evaluate Köhler theory at the nanometer scale, we calculated a w in all simulated droplet systems and compared our results with Eq. 15 as shown in Figure 5(b). For each system, we calculated the water activity as where ρ v,0 is the density of water in the vapor phase in a system with pure water and a flat water-air interface.

465
The theoretical prediction was determined by fitting a w values obtained for droplets with sizes larger than 3.5 nm using Eq. 15, with σ and V m obtained from simulations of systems with bulk-liquid-like water and a flat water-air interface (with the same salinity or organic surface loading) and x w fitted to minimize the root mean square deviation from simulation results. From Figure 5(b), we see that a w increases with decreasing droplet size. For pure water droplets, within the precision of our results, the simulation data match the theoretical prediction well for droplets 470 larger than 3.5 nm in diameter. For droplets smaller than 3 nm, water activity is underestimated by Köhler theory, indicating that small nano-droplets are less stable than expected (i.e., they exhibit higher saturated vapor pressure than expected, as noted above). For 2 m NaCl droplets, within the precision of our results, the simulation results are consistent with Eq. 15. However, we note that under the ideal solution assumption that underlies Köhler theory, 2 m NaCl droplet should have x w = 0.93, whereas our result are consistent with Eq. 15 with an effective water mole 475 fraction of x w = 0.86 ± 0.02, which indicates that the impact of NaCl on the activity of water is roughly twice as large as predicted by Raoult's law. Since the interatomic potential models used in this study accurately represent experimental water activity in bulk liquid water up to 3 m NaCl salinity (Mester and Panagiotopoulos, 2015a), this observation likely relates to the presence of the water-air interface. The magnitude of the observed deviation is consistent with the ion concentration enhancement effect discussed above. For droplets containing organic matter, as 480 the mixing state of water and organic matter is often poorly defined, the water effective mole fraction is often fitted to match measured values of a w (Seinfeld and Pandis, 2008). From our fitting results, the effective x w value in our low organic loading system equals 0.996±0.006, indicating that at low surface coverage PML has essentially no impact on water activity. In contrast, in our high organic loading systems, the effective x w value equals 0.890±0.011, indicating that PML can lower water activity by 10% in engulfed structure droplets. Overall, our results are consistent with

485
Köhler theory for droplets with sizes larger than 3.5 nm. However, Köhler theory overestimates the stability of water in droplets smaller than 3 nm (a phenomenon that may be related to the lower surface tension and water cohesion in small droplets identified above) and underestimates the impact of ions on water activity (a phenomenon that may be related to the ion concentration enhancement noted above).

Accommodation coefficient of water and condensation kinetics 490
The mass accommodation coefficient of water molecules on nano-aerosol particles, defined as the probability that a water vapor molecule that impinges on the particle surface becomes incorporated in the particle, is a fundamental factor that governs the hygroscopic growth kinetics of atmospheric aerosols (Diveky et al., 2021;Seinfeld and Pandis, 2008). Its precise value in atmospheric conditions remains highly uncertain (within a range of 0.01 to 1) despite decades of study because of its complex dependence on particle composition, phase-mixing state, and surface curva-495 ture (Diveky et al., 2021;Shiraiwa and Pöschl, 2021;Barclay and Lukes, 2019;Davis, 2006;Winkler et al., 2004).
In our study, the water accommodation coefficient (α) was calculated by tracking the trajectory and phase-state of water molecules during a 100 ns trajectory for each simulated system. We note that even based on atomistic trajectories determined using MD simulations, there are several possible ways of defining whether an impinging water molecule is absorbed or reflected and, as a result, different reported α-values (ranging from 0.9 to 1.0) even  (Ma et al., 2011;Chakraborty and Zachariah, 2008;Li et al., 2019;Takahama and Russell, 2011) and experimental results (Diveky et al., 2019;Duffey et al., 2013;Diveky et al., 2021) obtained at ambient temperature categorized by their fs values.
for the flat surface of pure liquid water (Julin et al., 2013;Barclay and Lukes, 2019;Skarbalius et al., 2021;Tsuruta and Nagayama, 2004;Ohashi et al., 2020). To highlight the impacts of curvature, salinity, and organic coatings on α and compare our results with those of other studies, we report normalized values as α * = α/α 0,w , where α 0,w is the value reported for the flat surface of pure liquid water in the same study. For studies that did not report α 0,w , we used α 0,w = 1 for simplicity. We convert organic concentrations into a surface coating factor, f s , defined as the 505 ratio of the surface area occupied by an idealized monolayer of organic molecules (with the thickness of a single carbon monolayer, 0.35 nm) to the surface area of an idealized spherical water droplet. For consistency, we restrict our comparison to studies that considered semi-volatile or non-volatile organic molecules with chain-like structures and > 5 carbon atoms. Figure 6 shows α * as a function of 1/D based on our results (circles) and previous studies of aerosols containing 510 liquid water and organic matter (stars) (Ma et al., 2011;Chakraborty and Zachariah, 2008;Li et al., 2019;Takahama and Russell, 2011;Diveky et al., 2019;Duffey et al., 2013;Diveky et al., 2021). The α * -values of the droplets containing pure or salty water (f s = 0) show no apparent dependence on salinity but an increase with droplet size.
This size-dependence of α * is consistent with previous observations that the mass accommodation coefficient of small water droplets depends on their surface curvature, particularly at droplet diameters below 3 nm (Barclay and Lukes, . Although the origin of this size dependence is unclear, we note that previous studies have demonstrated that the angle of incidence of impinging molecules can significantly impact α * , with more tangential collisions resulting in smaller α * -values (Nagayama and Tsuruta, 2003;Garrett et al., 2006). For an idealized geometry (with point-sized water molecules impacting a perfectly spherical droplet), the distribution of incidence angles of impinging water vapor molecules should be invariant with droplet size. The size-dependence of α * , suggests that deviations from 520 ideality, such as particle surface roughness or the finite size of incident water molecule, may impact the outcome of particle-water collision events.
For our low organic loading systems (f s = 0.30), shown as yellow circles, calculated α * -values show similar sizedependence as for NaCl droplets and converge to 0.70 in droplets larger than 3 nm. This observation is consistent with the expectation that for droplets with a sub-monolayer organic coating, the accommodation coefficient can be 525 roughly approximated as the fraction of the droplet surface that is not coated by organic matter, i.e. α * ≃ 1 − f s , as expected based on the simplistic assumption that α * ≃ 1 on the surface of liquid water and α * ≃ 0 on the organiccoated surface. This relationship between f s and α * is also consistent with previous MD simulation results obtained with water droplets containing azelaic acid, a C9-dicarboxylic acid (f s = 0.35, α * = 0.60) (Ma et al., 2011). For our high organic loading systems (f s = 1.50), shown as orange circles, predicted α * values are significantly smaller than 530 for our low organic loading systems, as expected, with α * ≃ 0.23 at the largest droplet sizes. Again, our α * results for high organic loading systems are consistent with values obtained in previous MD simulation and experimental studies for systems with f s > 1, which range from 0 to 0.3 as shown in Figure 6. Interestingly, α * now decreases with increasing droplet size, displaying a size dependence opposite to that observed for the other systems discussed above.
This unexpected behavior appears to be caused by the distinct size-dependence of the phase-mixing state observed 535 for organic-rich droplets, with smaller droplets displaying an organic-core, water-coated morphology, whereas larger droplets have an 'engulfed' structure that exposes less water at its surface.

Conclusions
To the best of our knowledge, this work presents an unprecedentedly broad MD simulation examination of the properties of liquid water in sub-10 nm aerosol particles. Specifically, we characterize droplet shape, water and solute 540 distribution, water diffusivity parallel and normal to the interface, surface tension, water evaporation free energy, water activity, and the water accommodation coefficient in aerosol particles with six different aqueous chemistries for a range of particle sizes. We build upon previous MD simulation studies that have carried out more piecemeal characterizations of water structure, dynamics, and energetics in aerosol particles (Karadima et al., 2019;Chowdhary and Ladanyi, 2009;Li et al., 2011Li et al., , 2010. We propose a sphericity factor (ϕ) that enables the quantitive character-545 ization of the mixing structure of water-, organic-, and salt-rich phases and helps reveal the links between aerosol particle morphology and the distinct energetic and kinetic properties of water molecules at nanometer scale.
A key finding of this study is that the size dependence of water microphysics in ultrafine aerosol particles is relatively mild and can be represented using continuum scale models in aerosol particles larger than 4 nm, whereas abrupt changes that require more detailed understanding of molecular level interactions are observed in sub-4 nm 550 diameter particles. For particles consisting of a water-rich phase, smaller droplets display less water cohesion, stronger ion concentration enhancement in the droplet core, faster interfacial water diffusion, smaller surface tension, smaller water accommodation coefficient, and higher water activity. For organic-rich aerosols, a micelle structure (with an organic core and water coating) can form at small particle sizes, which results in distinct water microphysics compared to water-core aerosols at similar sizes.

555
Finally, we reiterate that our simulations use a highly simplified proxy for aerosol organic matter as a single compound (PML). This compound was selected for its similarity to the compounds most abundantly observed in organic aerosols in terms of molecular weight (Thornton et al., 2020), number of C atoms (Zheng et al., 2021), O:C ratio (Zheng et al., 2021;Song et al., 2018), and functional groups (Zhao et al., 2018;Qi et al., 2022). Although the use of such an idealized proxy facilitates our effort to evaluate the sensitivity of aerosol properties to organic 560 loading, future studies should examine whether our predictions can be generalized to other organic compounds (or mixtures of organic compounds) abundantly found in organic aerosols, such as alcohols or other dicarboxylic acids.
Data availability. All data relevant to this research are available in the maintext and SI Appendix.

Appendix A: Discussion of using equimolar radius Re into surface tension calculation
As noted in the main text, the surface tension should be rigorously calculated at the radius Rs, where the surface water droplets, bubbles, water cavitation in minerals and water with hydrophobic solute) (Min and Berkowitz, 2019;Wilhelmsen et al., 2015;Kim and Jhe, 2020;Burian et al., 2017;Menzl et al., 2016;Azouzi et al., 2013;Joswiak et al., 2013;Sedlmeier and Netz, 2012). Among these studies, Kim and Jhe directly calculated (Re − Rs) from MD simulations of pure water droplet with radii ranging from 0.5 to 1.5 nm. Their results show that for nanodroplets with radii larger than 0.9 nm, (Re − Rs) is in the range of [-0.05, 0.05] nm. Our results in Section 3.4 similarly show 575 that δ is very small (∼0.06 nm). The sensitivity of Eq. 7 to the choice of R value should be given by the expression: For (R s − R e ) < 0.06 nm, which is the case in our study, the maximum error introduced by the use of R e in σ calculation is less than 8% for droplets with N w ≥500. This systematic error is commensurate with the statistical error of our predicted σ values (as shown by our reported error bars of ∼ ±8%).
(2) In MD simulation practice, 580 as noted in our manuscript, the droplets simulated here do not have a perfectly spherical shape. As a result, it is impossible to accurately determine the surface tension surface Rs in a traditional manner. The surface tension calculation scheme used in our study (Eq. 6-Eq. 8), which uses Re for approximation, is widely used in previous MD simulation studies for nanodroplet surface tension calculation (Li et al., 2010;Zhao et al., 2019;Li et al., 2013;Zakharov et al., 1997).

Appendix B: Sensitivity analysis of ion concentration enhancement to ion interfacial exclusion width
The first hydration shell radii (r hyd ) of alkali cations Li + , Na + , and K + are 0.279, 0.318, and 0.359 nm, respectively, as reported by Joung and Cheatham III (2009). The interfacial width parameter (exclusion width d exc ) was determined following two steps: (1) from the ion density profile relative to the Gibbs dividing interface ( Figure A4e), we can see 590 that the distance where Na+ density reaches the well-mixed density 1m, ranges from 0.3 to 0.5 nm for all simulated droplets, which provides estimates of the upper and lower bounds of d exc ; (2) we tested different values of the ion exclusion width and calculated the corresponding ion concentration enhancement factor and, then, compared the calculated results with our MD simulation results. Our results show that setting d exc to 1.2r hyd , yielded the best match to our MD simulation results. The sensitivity of ϵ values of NaCl droplets to d exc are illustrated in Figure A2.

Appendix C: Deliquesence relative humidity as a function of nanoparticle size
To evaluate the DRH based on our simulation results, we assume that for a nanoparticle at the DRH, ion concentration in the droplet core equals the solubility of the salt C sol . We denote the wet particle diameter (including water molecules) as D wet and the cation hydration radius as r hyd . If we follow the same assumptions as in Section 3.2 (i.e., that the salt is excluded from the vicinity of the water surface in a region with a thickness equal to 1.2 r hyd ), the 600 ion concentration corresponding to a well-mixed system, C, can be calculated from The dry diameter of this deliquescent particle D can be approximately evaluated as falling within the range [D salt , And the expression of a w (C) is taken from the experimental study of bulk-liquid solutions (El Guendouzi and Dinane, 2000).
We tested the impact of different thermostat settings -NVT, NVE, NVT with fixed system angular momentum, and Langevin thermostat (LD) with friction constant γ = 0.5 ps −1 (this value results in a friction that is lower than the internal friction of water, while still providing efficient thermostatting) -on the value of the water diffusion coefficient D ∥ using the Einstein relation for a water droplet with N w = 1000. Resulting plots of water mean square displacement (MSD) versus time are shown in Figure A6. Different symbols in Figure A6 represents the MSD values 615 of water molecules in different 3 Å layers categorized based on distance from the time-averaged Gibbs dividing interface. As shown in Figure A6, the main text are based on simulations carried out using the NVT thermostat with fixed system angular momentum.
D2 Sensitivity to the frame of reference ('com' and 'inst' scheme)

D2.1 Diffusion coefficient parallel to the interface (D ∥ )
The Einstein relation was used to calculate D ∥ under both the 'com' and 'inst' schemes, i.e., relative to the average or instantaneous interface. Both schemes are expected to provide mutually consistent results for D ∥ , irrespective of 635 distance from the interface, because capillary waves exist only at the interface and are mostly transverse waves that should have negligible impact on water displacement parallel to the interface. A comparison of D ∥ values calculated under both schemes is presented in Figure A8(a) for droplets with N w = 3000. As expected, the two schemes provide mutually consistent values. D ∥ calculated for different simulated systems under the 'com' scheme are shown in Figure   A8(b).

D2.2 Diffusion coefficient normal to the interface (D ⊥ )
As discussed in the main text, D ⊥ values were determined as a function of distance from the interface by comparing of the interfacial water dynamics, whereas the 'inst' scheme excludes the interface fluctuation and only quantifies water displacement relative to the instantaneous interface. Since the magnitude of capillary waves increases with droplet size (as observed based on both ϕ values and interfacial width as noted above), the difference between the D ⊥ values at the interface calculated with the two schemes should increase with droplet size.
As shown in Figure A9, we obtained mutually consistent D ⊥ values near 0.25 Å 2 ps −1 under both 'com' and 'inst'  (Liu et al., 2004;Wick and Dang, 2005). This discrepancy might be due to different water models, different interfacial layer definitions, or the differences in the curve fitting procedure required to determine D ⊥ . As noted in the main text, our observations of relatively low D ⊥ values near the interface are consistent with reports of slower water reorientational kinetics near hydrophobic surfaces. The discrepancies between different reported values of D ⊥ for interfacial water should be carefully examined in future studies.  ∥ Figure A7. Diffusion coefficient D ∥ of water molecules (i.e., in the direction parallel to the interface) as a function of distance from the interface in simulations of a pure water droplet with Nw = 1000 carried out using different thermostats (calculated from the slope, from 2 to 6 ps, of the MSD data shown in Figure A4).