Assessing long-range contributions to the charge asymmetry of ion adsorption at the air–water interface†

Anions generally associate more favorably with the air–water interface than cations. In addition to solute size and polarizability, the intrinsic structure of the unperturbed interface has been discussed as an important contributor to this bias. Here we assess quantitatively the role that intrinsic charge asymmetry of water's surface plays in ion adsorption, using computer simulations to compare model solutes of various size and charge. In doing so, we also evaluate the degree to which linear response theory for solvent polarization is a reasonable approach for comparing the thermodynamics of bulk and interfacial ion solvation. Consistent with previous works on bulk ion solvation, we find that the average electrostatic potential at the center of a neutral, sub-nanometer solute at the air–water interface depends sensitively on its radius, and that this potential changes quite nonlinearly as the solute's charge is introduced. The nonlinear response closely resembles that of the bulk. As a result, the net nonlinearity of ion adsorption is weaker than in bulk, but still substantial, comparable to the apparent magnitude of macroscopically nonlocal contributions from the undisturbed interface. For the simple-point-charge model of water we study, these results argue distinctly against rationalizing ion adsorption in terms of surface potentials inherent to molecular structure of the liquid's boundary.


I. ELECTROSTATIC CONTRIBUTIONS FROM NEAR AND FAR
The challenge of identifying and interpreting a potential drop across the liquid-vapor interface can be viewed as an issue of partitioning molecules between distinct regions of space.
Consider a macroscopic droplet of liquid bounded by an interface S with the vapor phase (as illustrated in Fig. S1). The origin of our coordinate system lies deep within the bulk liquid phase. We will aim to calculate the average electric potential φ at the origin, distinguishing contributions of molecules that are far from the probe (including those at the phase boundary) from those that lie nearer the origin. Specifically, we will divide the two populations at an imaginary surface B that is also deep within the bulk liquid. We will take B to be distant enough from the origin that liquid structure on this surface is bulk in character, even if the microscopic vicinity of the origin is complicated by a solute's excluded volume.

A. Partitioning schemes
The vast majority of molecules in the droplet are unambiguously located either outside B where N is the total number of molecules in the droplet and α indexes charged sites within each molecule. Here, p(r, Ω) = δ(r − r (0) )δ(Ω − Ω (0) ) is the joint probability distribution of a molecule's position (i.e., r (0) ) and intramolecular configuration Ω (0) (specified relative to the reference position r (0) , as indicated by the superscript). 2 By Δr α = r α − r (0) we denote the displacement of charge q α from the reference point r (0) . This intramolecular displacement is entirely determined by Ω (0) .
For the P-scheme, each charge α contributes to φ P far if r α lies outside B. The corresponding far-field potential is where p α (r, Ω) is the joint probability distribution for site position r α and intramolecular configuration of a solvent molecule. In Eq. S3 we have made use of the connection between the distributions p and p α .

B. Multipole expansion
Since the entire "far" region is macroscopically distant from the origin, small-Δr α expansions of |r + Δr α | −1 and p(r − Δr α , Ω) are well justified. These yield where we have omitted leading terms proportional to α q α , which vanish by molecular charge neutrality. When carried through subsequent calculations, terms beyond quadrupole order in these expansions would vanish due either to symmetry or to the macroscopic scale of the droplet.
Defining dipole and quadrupole densities as and and Integrating by parts, and noting that m(r) and ∇ : Q(r) vanish both on B and at infinity, Using the divergence theorem, where R is a point on B andn(R) is the corresponding local inward-pointing normal vector.
Since B lies within the bulk liquid, where the average quadrupole density Q liq is isotropic, Q(r) = I (TrQ liq /3) everywhere on this surface. As a result, These two measures of the far-field potential are thus different. Moreover, the quadrupole trace that determines this difference depends on the choice of r (0) . This ambiguity is a well-known feature of the so-called Bethe potential −(4π/3)TrQ liq . 3-8

C. Dipole surface potential
To simplify the result for φ M far , note that Q(r) is isotropic everywhere outside B, except in the microscopic vicinity of S. In the bulk regions of the far domain, we then have Q(r) : ∇∇r −1 ∝ δ(r) = 0. The final term in Eq. S10 therefore has nonzero contributions only from a thin shell whose volume is proportional to L 2 , where L is the macroscopic scale of the droplet. Since ∇∇r −1 ∼ L −3 in this shell, the quadrupolar contribution to φ M far has a negligible magnitude, L −1 . As a result, This integral similarly has nonzero contributions only from a microscopically thin shell of broken symmetry, centered on the phase boundary S. Since the macroscopic surface is very smooth on this scale, and because the average dipole density points normal to the locally planar interface, the far-field potential may be written where R is the point on S nearest to r, the coordinate z = (r−R)·n(R) is the perpendicular displacement from the liquid-vapor interface,n(R) is the outward-pointing normal of S, and m ⊥ (z)n(R) is the average dipole field at r. Neglecting contributions of O(z/L), we may replace r −1 by R −1 , and easily evaluate the surface integral, yielding where the integral is performed in the direction from liquid (z liq < 0) to vapor (z vap > 0).

D. Near-field potential
In evaluating φ far , we have made no assumptions about the liquid's structure near the probe. If the origin lies inside a solute's excluded volume, then the near-field potential is complicated by the microscopically heterogeneous arrangement of solvent molecules in its vicinity. If, however, the probe is simply a point within the isotropic bulk liquid, then φ near can be easily determined.
For a probe that resides in uniform bulk liquid, m(r) = 0 and Q(r) = Q liq everywhere inside B. In the P-scheme we can conclude immediately from the analogue of Eq. S11 that φ near = 0. In the M-scheme we have In either case the total potential sums to

A. Outline
Here we present details of the piecewise linear response (PLR) model discussed in the main article. The PLR model is based on the observation that solvent response to charging a solute is linear for both anions and cations, but differs between the two cases. [13][14][15] In such a model, the average electrostatic potential due to the solvent at the center of a charged cavity can be written as where q c is the value of the 'crossover charge' between the two linear regimes, (δφ solv ) 2 + is the variance of φ solv for q ≥ q c , and (δφ solv ) 2 − is the variance of φ solv for q < q c . (As written, it is implicitly assumed that q c ≤ 0, as suggested by simulations.) Let us define and ψ is, In general, φ neut , q c and J will depend upon solute size, and whether or not the solute is located in bulk or at the interface.  in bulk, but some small deviations are seen. These deviations are more pronounced when the solute is at the interface. For R = 0.240 nm, the above PLR model breaks down at large negative q, but it remains reasonable for smaller values of the absolute charge. By fitting straight lines to the anion and cation response, we can obtain values for q c , (δφ solv ) 2 + and (δφ solv ) 2 − . The results from using these in Eq. S26 to compute Δ ads ψ (PLR) are presented in Fig. 4b in the main article. Results for R = 0.75 nm and R = 1 nm are not shown because, while anion and cation response do still differ, the degree of nonlinearity is much less on an absolute scale than for the smaller solutes. This makes it challenging to reliably obtain q c .  IV. CONSTRUCTING P 0 (φ solv )

Figures
In order to compute F chg (q) from Eq. 7, we require P 0 (φ solv ), the probability distribution of φ solv . For the range of q of interest, i.e. −1 ≤ q/e ≤ 1, sampling P 0 directly (in the absence of solute charge) would yield grossly insufficient data in the extreme wings of the distribution. Instead, we obtain P 0 by histogram reweighting using MBAR. 17 As an illustration, Fig. S10 (a) shows probability distributions P q (φ solv ) of φ solv at the center of the solute (R = 0.240 nm) with different values of q. Using data from simulations across the full range of q, we then construct P 0 (φ solv ), as shown in Fig. S10 (b). Gaussian distribution with mean and variance obtained from the simulation at q/e = 0. Note that finite size corrections have not been applied to these plots.