Solution Structure Effects on the Properties of Electric Double Layers

Oct 24, 2018 - Center for Micro-Engineered Materials and Department of Chemical and Biological Engineering, University of New Mexico , Albuquerque , N...
0 downloads 0 Views 4MB Size
Invited Feature Article pubs.acs.org/Langmuir

Cite This: Langmuir 2018, 34, 13808−13820

Solution Structure Effects on the Properties of Electric Double Layers with Surface Charge Regulation Assessed by Density Functional Theory Frank B. van Swol and Dimiter N. Petsev*

Downloaded via TULANE UNIV on December 16, 2018 at 17:06:38 (UTC). See https://pubs.acs.org/sharingguidelines for options on how to legitimately share published articles.

Center for Micro-Engineered Materials and Department of Chemical and Biological Engineering, University of New Mexico, Albuquerque, New Mexico 87131, United States ABSTRACT: The structure of electrolyte solutions in electric double layers is primarily determined by the solvent, despite the fact that it is usually neutral and not subject to Coulombic interactions. The number of solvent molecules in a typical electrolyte solution may be significantly greater that the number of ions. Hence, the charged species compete for space with a much larger number of neutral molecules, which has a strong effect on the density distributions near charged surfaces. Even for very dilute electrolyte solutions, the density profiles resemble liquidlike structure, which is entirely due to the presence of the dense solvent. Our work demonstrates that the solvent structural effect strongly couples to the surface chemistry, which governs the charge and potential. We argue that a comprehensive statistical−mechanical approach, such as classical density functional theory that explicitly includes all solution species, in combination with a surface charge regulation condition at the interface, provides an excellent approach for describing charged interfaces. It allows for revealing important physical features and includes non-Coulombic contributions such as ionic and surface solvation.

1. INTRODUCTION

∇2 Ψ = −

Interfaces between electrolyte solutions and other substances occur in many systems of fundamental and practical importance. Examples include surface and colloid science, electrochemistry and corrosion, soft materials and biomaterials. These interfaces are subject to complex physical and chemical interactions with fluids, such as chemical reactions between surface terminal groups and certain solution components, physical adsorption, and desorption.1,2 All of these interactions may lead to the accumulation of charge at the interface. As a result, the mobile charges in the vicinity redistribute in accordance with the electrostatic field generated by the surface charge. Hence, the effect propagates in the neighboring phases and is characterized by spatial charge and electrostatic potential distributions. The whole system, consisting of a charged interface and locally distributed mobile charges and a potential, is often referred to as the electric double layer (EDL). The name has historical roots, reflecting the initial notion that the charged interface with the electrolyte solution can be represented by a simple capacitor in which the surface charges and the solution countercharges are placed on two well-defined planes in space.3 This is certainly not an accurate representation of the physical situation because all mobile charges (e.g, ions in the solution) are subject to Brownian motion at finite temperatures. An alternative model has been suggested,4−6 on the basis of Maxwell’s theory of electrodynamics, relating the potential Ψ(r) to the charge density distribution ρe(r) © 2018 American Chemical Society

ρe εε0

(1)

where ε0 is the dielectric constant of vacuum and ε is the relative dielectric permittivity of the solvent. Assuming an equilibrium consisting of the Boltzmann distribution of all ionic species ρe =

i −qi Ψ yz zz zz k T k B {

∑ ρi0 qi expjjjjj i

(2)

and inserting it into eq 1 leads to the very popular Poisson− Boltzmann (PB) model of the EDL. qi is the ion charge, and kBT is the thermal energy. The linearized version of the PB equation was used by Debye and Huckel to develop their celebrated theory of electrolyte solutions.7 The PB model is based on a second-order differential equation (eq 1), and as such, it requires boundary conditions. The simplest examples are the Dirichlet (fixed potential) and Neumann (fixed charge) conditions or a combination of the two. This approach has the great advantage of being mathematically easy to manipulate, and despite the nonlinearity of the PB model (eqs 1 and 2), it allows us to obtain analytical solutions for simple geometries.4−6 However, a close inspection of the charging mechanisms of the interface with the Received: July 19, 2018 Revised: October 23, 2018 Published: October 24, 2018 13808

DOI: 10.1021/acs.langmuir.8b02453 Langmuir 2018, 34, 13808−13820

Langmuir

Invited Feature Article

2. DENSITY FUNCTIONAL THEORY OF ELECTRIC DOUBLE LAYERS WITH SURFACE CHARGE REGULATION 2.1. Grand Thermodynamic Potential. EDLs involve multicomponent solutions containing charged (ionic) and neutral (solvent) species. These species are subjected to the effects of the field exerted by the interface with the substrate. The field has electrostatic and van der Waals components, which lead to a variety of interactions with the fluid. In addition, the solution components (M) in the fluid interact with each other. All of these interactions are easily incorporated into a cDFT model. The electrolyte solution is described in terms of a grand thermodynamic potential functional, which for single flat EDL reads27−31,34,35

electrolyte solution indicates that neither the Dirichlet nor the Neumann boundary conditions are justified. The surface charge is due to a thermodynamic equilibrium between the surface reactive groups and the solution components. A boundary condition that correctly accounts for the physical situation was first proposed by Ninham and Parsegian8 and is often referred to as charge regulation. The surface charge regulation reflects the mechanism of charge formation and therefore can distinguish between different realistic cases.9−16 An excellent overview of the PB analysis of EDLs using charge regulation boundary conditions was recently published by Trefalt et al.2 The PB method suffers from a major deficiency. It assumes that the electrolyte solution is a structureless continuum and that all effects stemming from the actual size and various interaction between the ions and the solvent molecules are lost. This is particularly inaccurate at small (nanometer) length scales. Hence, very important fundamental and applied problems related to the properties of electrolyte solutions in confined spaces very close to the charged interface cannot be adequately analyzed. The PB approach has been modified to account for the finite ion size and steric interactions near the charged interface,17−19 but this approach remains still quite limited. A more accurate alternative is offered by modern statistical mechanics. However, its methods can be applied at different levels of detail. A popular model (sometimes referred to as “primitive”) takes into account the ion sizes and interactions but continues to consider the solvent as structureless. For low to moderate electrolyte concentrations, the primitive model was shown to be equivalent to the modified PB approach mentioned above.20 A better model (known as “civilized”) should include the explicit structural contribution due to the solvent molecules and all interactions associated with them.21−31 Note that the civilized model would never approach the modified PB limit, to reduce the ion concentration, because of the structural effect of the solvent. The solvent density is much greater than the densities of all ionic species in the electrolyte solution. Hence, it dominates the overall liquid structure. Recently we examined the interplay between the local solution structure near the EDL charged interface and the charge regulation chemical equilibrium that governs the surface charge using classical density functional theory (cDFT).28−31 This approach is well suited to account for a wide range of molecular and ionic interactions in the fluid as well as with the adjacent substrate. In addition, cDFT has been extensively tested against computer simulations, and it was demonstrated that the resultant density distributions from both approaches are in excellent agreement.32 Our cDFT analysis indicates that the solution structure is dominated by the presence of the neutral solvent, and it strongly couples to the surface charge regulation chemistry. Hence, it has a significant impact on the EDL electrostatic properties. This coupling effect is the main focus of the present study. A rigorous statistical mechanical approach allows us to better interpret established quantities in electrochemistry and surface and colloid science like the Helmholtz planes3 and the Stern layer33 and offers a pathway for the quantitative modeling of charged interfaces that involve electrolyte solutions.

ex ex Ω[{ρi (z)}] = F id[{ρi (z)}] + FHS [{ρi (z)}] + Flong [{ρi (z)}] M

+ 2π ∑

∫ 9 d9 ∫ dz ρi (9, z)[Viext(9, z) − μi ]

i=1

(3)

The radial coordinate 9 runs parallel to the EDL interface with the substrate. The distance coordinate z is normal to the interface. The first term on the right-hand side of eq 3 corresponds to the ideal contribution to the free energy, M

F id[{ρi (z)}] = 2πkBT ∑

∫ 9 d9

i=1

∫ dz ρi (9, z){ln[λi3ρi (9, z)] − 1}

(4)

where λ = h2/(2πmikBT ) is the thermal de Broglie wavelength, h is Planck’s constant, mi is the mass of species i, and ρi(z) is the local density of component i along the z coordinate normal to the wall. The excess free energy consists of hard-sphere and longrange parts. The hard-sphere contribution is based on the derivation of Rosenfeld36,37 and reads ex FHS [{ρi (z)}] = 2πkBT

∫ 9 d9 ∫ dz ΦHS{nα(9, z)}

(5)

ΦHS{nα(9, z)} is the hard-sphere reduced free energy, and nα(9, z′) = nα(r) (with r being the position vector) is the weighted local density. ΦHS{nα(r)} = −n0 log(1 − n3) + +

n1n2 − n1·n2 1 − n3

n2 3 − 3n2 n2 ·n2 24(1 − n3)2

(6)

Remarkably, the functional form of the reduced free energy ΦHS{nα(9, z)} is independent of the number of components M. The weighted densities nα(r) exhibit such a dependence and are defined by M

nα(r) =

∑ ∫ d3r′ρi (r′) ωαi (r − r′) i=1

(7)

and the weighting functions ωiα(r) are 13809

DOI: 10.1021/acs.langmuir.8b02453 Langmuir 2018, 34, 13808−13820

ÄÅ É ÅÅ 2 i d y9 i d y3ÑÑÑ ÅÅ jj s zz Ñ s j z ΦLJ(z) = ϵs − w ÅÅ jj zz − jjj zzz ÑÑÑ, z > ds/2 ÅÅ 15 k z { ÑÑ z k { ÑÖ ÇÅ

Langmuir

Invited Feature Article

ω3i(r) = Θ(R i − r ), ω2i (r) = δ(R i − r ), ω1i (r) =

ωi2(r) i ωi (r) , ω0(r) = 2 2 , 4πR i 4πR i

is used in the analysis. The electrostatic interaction of the ions with the interface is q σz Φel (z) = i , z > di /2 2εε0 (15)

ωi2(r)

r Θ(R i − r ), ω1i (r) = 4πR i r

ωi2(r) =

(8)

θ(Ri-r) and δ(Ri-r) are the step, and Dirac delta functions respectively. The long-range contribution to the free energy functional is ex Flong [{ρi (z)}] =

π 2

M

Using a solution model based on a combination of Coulombic and LJ interactions is clearly an approximation when polar or polarizable solvents (i.e., water) are involved. Nevertheless, it captures important features of electrolyte solutions, such as attractive and repulsive interactions between all solution components.28 This approximation works reasonably well for low to moderate ionic concentrations.22,29,41−43 In such cases, the average distance between the charged species is sufficiently large so that any ion-induced local perturbation of the solvent molecular structure relaxes to the bulk state before reaching another ion.44 For example, in our case, the overall ionic concentration corresponds to 0.01 M. Assuming that the solvent molarity is 55.5 M, the average distance between two ions is ∼30 solvent molecular diameters. Even close to the charged surface, where the ion densities are greater (Figure 1), the distance separating two ions is ∼10 molecular diameters. As Carnie and Chan showed,41 the local polarization perturbation decays to the bulk value over distances of the order of 3 to 4 molecular diameters. The overall free energy may also be reduced due to local polarization response of the solvent near the charged interface. We can mimic that effect by varying the solvent−wall interaction. While empirical and therefore not ideal, this approach has the advantage of being considerably simpler and computationally efficient. The Coulombic interactions are included in the term Fex long[{ρi(z)}] (eq 9) and are therefore treated within a meanfield approximation. cDFT allows us to include higher-order correlations that are beyond the electrostatic mean-field approximation, but they lead to noticeable contributions at much higher ionic concentrations.45 2.2. Charge Regulation and Charge Formation at the Electric Double Layer Interface. The electrostatic properties (charge and potential) at the interface of the EDL are determined by a charge regulation boundary condition.8,9 This condition derives the surface charge from the specific chemical equilibrium between the surface reactive groups and certain ionic species in the solutions: the potential-determining ions (PDIs). This representation is much more realistic than fixing the surface charge or potential as external parameters. We use the model of Chan et al.,10 which considers the surface reaction equilibria

M

∑ ∑ ∫ 9 d9 ∫ dz × i=1 j=1

∫ dz′ρi (9, z) ρj (9, z′) ΦLR (9, |z − z′|)

(9)

where ΦLR (9, |z − z′|) is the contribution of the long-range interactions of the reduced free energy. Equation 9 indicates that the long-range interactions are accounted for in a meanfield limit.38 The last term in eq 3 is the Lagrangian constraint, which accounts for the external fields Viext and fixed chemical potentials μi for all species. The densities of all components ρi(z) are found by minimizing the functional39 δ Ω[{ρi (z)}] δρi (z)

=0 (10)

The charge density distribution in the EDL is derived by summing over the individual charged species, taking into account the charge numbers with their sign M

ρe (z) =

∑ qiρi (z) i=1

(11)

This approach includes the contributions of the interactions between the solutions species, which are present in the term 36,37 Fex accounting for the excluded volume effects, HS[{ρi(z)}], ex and Flong[{ρi(z)}], which captures the long-range interactions. In our model, these interactions are of the Lennard-Jones (LJ)40 ÅÄÅ ÑÉ ÅÅij d yz12 ij d yz6ÑÑÑ ij ÅÅjj ij zz Ñ ΦLJ(rij) = 4ϵijÅÅÅjj zz − jjjj zzzz ÑÑÑÑ, rij > dij j z j z r Ñ ÅÅÅk rij { k ij { ÑÑÑÖ (12) ÅÇ and Coulombic Φel (rij) =

qiqje

AH 2+ + BH F AH + BH 2+, pK + = −log10 K+

2

4πεε0rij

(14)

, rij > dij

AH + BH (13)

F A− + BH 2+, pK − = − log10 K −

(16)

in conjunction with the bulk reaction

types27−31 where dij = (di + dj)/2, di is the diameter of component i, and rij(R, z) is the distance between species i and j. All non-Coulombic interactions of a molecule (or ion) of type i with the interface are assumed to be of the “hard wall” type (i.e., only the excluded volume effects are taken into account). The only exception is our analysis of the solvent− wall interactions presented in Section 3.3 below, where the LJ (9-3) potential40

2BH F B− + BH 2+

(17)

where K+ and K− are the equilibrium constants and AH is a surface chemical group that can attach or release a proton depending on the local densities of PDIs species BH2+. This model is relatively simple, but it still allows us to examine subtle effects such as pH variation around the surface isoelectric point. 13810

DOI: 10.1021/acs.langmuir.8b02453 Langmuir 2018, 34, 13808−13820

Langmuir

Invited Feature Article

chemical equilibrium at the EDL interface (e.g., eq 19) is equivalent to adding an appropriate surface functional term to eq 3.46 The grand thermodynamic potential (eq 10) is minimized using the free-ware Tramonto code.47 The parameters governing the surface chemistry are Γd2/e = 0.66, pK− = 6, and pK+ = −2, which are common to oxide surfaces (e.g., silica).28,29,48 This choice of parameters assigns a negative charge to the surface for all values of the solution pH greater than pI = 2. The diameters of all solution species d are the same (except in cases where indicated otherwise). In order to ensure a correct conversion between scaled densities and molarity, we select d = 0.288 nm.28,29 Thus, 1 L of solvent contains 55.5 M. The Debye screening parameter is defined as7

The densities of all ions in the immediate vicinity of the charged interfaces are given by28,29 d /2

s

ρi =

∫0 i dz ρi (z) d /2

∫0 i dz

(18) 10

The surface charge is σ

= eΓ

s s ρ AH + − ρ − A 2

s s s ρAH + ρ AH + + ρ − A 2

= eΓ

δ sinh[e(ΨN − Ψs)/kBT ] 1 + δ cosh[e(ΨN − Ψs)/kBT ]

(19)

ij 1 κ = jjjj j εε0kBT k

where δ = 2 K −/K+ , Γ is the surface concentration of ionizable groups, and Ψs is the surface potential. The latter is directly related to the overall potential distribution in the EDL. Hence, even ions that are chemically inert and are not involved in the chemical reactions (eq 16) affect the surface charge. This is due to simple charge screening or, more subtly, the effect of the chemical potential and therefore the local concentration of potential-determining ions. Quantity ΨN = [ln(10)kBT/e](pI − pHb) is the Nernst potential, where pI = (pK− + pK+)/2 is the surface isoelectric point and pHb is the pH value in the bulk solution, far from the charged interface. The surface charge depends on the surface potential and the Nernst potential and on equilibrium constants K+ and K−. The Nernst potential depends on the pH in the solution bulk and on the isoelectric point, which in turn depends only on K+ and K−. Equation 19 corresponds to the charge at the physical interface between the substrate and the electrolyte solution, also known as the inner Helmholtz plane (IHP).3 The charge of the Stern layer is located at the outer Helmholtz plane (OHP).3,33 The actual physical situation could be more complicated, as we discuss in Section 3.3. Still, for clarity we refer to the adjacent surface layer as the Stern layer when discussing some of the results below. Finding the density distributions of all ionic species in the solution allows for determining the overall charge density distribution more accurately than using the approximate equation (eq 2). The charge density can then be introduced into eq 1 to find the electrostatic potential distribution in the solution and at the surface. The surface charge is also equal (and opposite in sign) to the net fluid charge, or M

σ = −∑ i=1

∫0

2z zz

1/2

(21)

The Debye parameter has the dimension of inverse characteristic length, which is relevant to the asymptotic charge and potential decay far from the interface. In the present study, (κd)−1 = 10.56, and we use that number to determine the size of the system, which we set to 100d. It is also instructive to compare this long-ranged electrostatic length to the characteristic length scale associated with the structural decay, which was found to be ∼1.3d.28,29 The boundary condition for the surface charge is given by eq 19. The densities, far away from the interface, are equal to their fixed bulk values. More details of the computational procedure are provided elsewhere.28,29

3. SOLUTION STRUCTURE AND SURFACE CHARGE REGULATION The solvent in an electrolyte solution is involved in various interactions with the solutes (ions) and with the surface where the charge regulation reactions occur. These interactions are not of the simple Coulombic type, but they are still very important. The solvent occupies volume, which is unavailable to the ions. It may also be involved in attractive interactions with itself, the ions, and the surface. Therefore, it is not surprising that the solution structure is significantly affected by these non-Coulombic solvent effects. The structure then determines the local densities of all species. Among those species are the PDIs, which take part in the surface chemical equilibria.16 Hence, despite the neutral nature of the solvent molecules, they participate in a range of complex interactions, and because of their great number, they may have a dominant effect on the properties of the EDLs. Below we discuss a few typical examples that are related to ionic and surface solvation. 3.1. Ion Density Distributions in the Electric Double Layer. The ion density distributions in an EDL near the charged surface are shown in Figure 1. The plots are normalized by the corresponding bulk densities ρbi . The surface charge is governed by the regulation condition (eq 19). We explore the combination of two pH values and two different ionic strengths. The pH variation affects the surface charge (i.e., a greater difference between the solution pH and the surface pI leads to a larger resultant surface charge). The variation of the bulk density of ions, or the ionic strength, primarily has an effect on the range of the electrostatic potential in the solution phase. However, since the ionic



qiρi (z) dz

y ∑ ρi qi zzz i { 0

(20)

Since the solvent is neutral, it is formally treated as a component with zero charge (qs = 0). We consider a flat EDL, which means that all charge distributions are averaged over planes parallel to the interface. This is evident from the integration over the radial coordinate 9 in eqs 4, 5, and 9. The charge regulation model, described above, also results in a surface charge which is uniformly distributed over the interface. In that respect, the surface charge treatment is in agreement with the cDFT representation of the adjacent fluid. Equations 1, 19, 11, and 20 are then iteratively solved to obtain the surface charges and potentials as well as the corresponding distributions in the fluid phase. Defining a 13811

DOI: 10.1021/acs.langmuir.8b02453 Langmuir 2018, 34, 13808−13820

Langmuir

Invited Feature Article

Figure 1. Ion density profiles ρi(z) near a charge-regulating interface. ρbi is the bulk density of component i, far from the surface. The solid curves are for solvent diameter equal to that of the ions, d. The dashed curves are for solvent diameter equal to 0.01d. All blue curves correspond to the positive non-PDIs, while the red curves correspond to the negative ions. (a) pH = 4, I = 0.01 M; (b) pH = 4, I = 0.001 M; (c) pH = 5, I = 0.01 M; and (d) pH = 5, I = 0.001 M. The remaining parameters are given in the text.

involving the solvent molecules are not ignored. At the same time, the ion densities are too low to produce any noticeable structural peaks near the charged surface. Another observation is that the curves with and without the solvent contributions do not perfectly match, even at distances where the structure is no longer present. This is due to the solvent-induced structural effect on the surface charge regulation (eq 19). The surface charges and potentials are slightly more negative when the solution does not exhibit any structuring, as shown in Table 1.

strength modifies the local potential distribution, it leads to slight variations in the surface charge (eq 19). The selected surface parameters and pH values result in a negative surface charge. Hence, the positive ions (blue curves) are attracted to the interface, and their density locally increases. The negative ions (red curves) are repelled, and hence their density near the charged interface is lower. The solid curves in Figure 1 show the density distributions obtained from the cDFT, explicitly taking into account the presence of the solvent molecules. The diameter of all species is the same. All interactions between the solution components include excluded volume and LJ attraction with energy parameter ϵij = kBT (eq 12). In addition, all ions experience Coulombic interactions in accordance with eq 13. It is important to emphasize that the ionic density is rather low (gaslike), but the density profiles exhibit well-defined (liquidlike) structural peaks. This structuring is due to the explicit presence of the solvent, which dominates the underlying liquid structure. The ions can be placed only in a void space that has not been occupied by a solvent molecule. This means that they effectively replicate the dense liquid structure. The electrostatic field exerted by the charged interface contributes to the slopes of the curves at larger distances. These slopes are overall negative for the positive ions and positive for the negative ions. The dashed curves illustrate the electrolyte behavior in the case where the neutral solvent is practically structureless (i.e.,the primitive model limit). The detailed interactions

Table 1. Surface Potential and Surface Charge Variations with the Solution pH and Electrolyte Concentration structureless solvent pH 4 4 5 5

I(M) −2

10 10−3 10−2 10−3

explicit solvent

eΨs/kBT

2

σd /e

eΨs/kBT

σd 2/e

−0.883 −1.440 −2.222 −2.872

−0.00264 −0.00147 −0.00766 −0.00369

−0.721 −1.239 −1.980 −2.634

−0.00215 −0.00122 −0.00666 −0.00322

This finding demonstrates that ignoring the solvent-induced solution structure leads to overestimated results for the magnitude of the surface charges and potentials. Using traditional boundary conditions such as constant charge or potential at the surface would also miss this effect. Hence, it is the coupling between solution structure and surface chemistry that reveals this detail. 13812

DOI: 10.1021/acs.langmuir.8b02453 Langmuir 2018, 34, 13808−13820

Langmuir

Invited Feature Article

Figure 1a shows the ionic density profiles for pH = 4 and I = 0.01 M. The surface charge and potentials are low (Table 1), and the Debye screening length is equal to κ−1 = 10.56d. For an ionic diameter of 0.288 nm, the dimensional screening corresponds to 3.04 nm. Reducing the ionic strength to I = 0.001 M increases the screening length to κ−1 = 33.38d, or 9.61 nm. This is the reason for the more gradual decrease in the local density with distance (Figure 1b). In addition, the background electrolyte concentration reduction also leads to an increase in the magnitude of the surface potential and a decrease in the magnitude of the surface charge (Table 1). Changing the solution pH to 5 (while keeping I = 0.01 M) increases the magnitudes of both the surface potential and charge (Figure 1c and Table 1). As a result, the local density of positive ions near the surface increases while that of the negative ions decreases. Finally, Figure 1d depicts the case where pH = 4 and I = 0.001 M. The surface potential magnitude increases further, while the surface charge magnitude drops in comparison with the higher electrolyte concentration case shown in Figure 1c. The local density of positive ions near the surface is even higher, and that for the negative ions is lower. The propagation of the charge and potential into the fluid is greater because of the larger screening length. An important conclusion stemming from the above analysis is that the explicit account for the neutral solvent effects leads to a proper description of the solution structure in the EDL. This structure is important for the correct computation of the surface potential and charge in the presence of surface charge regulation (eq 19). The surface potential and charge then determine the properties of the whole EDL, including the far regions where structure has already disappeared. The solution structure varies with the distance from the surface as ∼exp(−βz) with characteristic length β−1 = 1.3d.28,29 Taking into account the structuring near the surface allows for the determination of the composition of the Stern layer30,31 and the location of the electrokinetic shear plane.29 3.2. Effects Due to Ionic Solvation. An important advantage of solvent-explicit molecular models of EDLs is that they can naturally incorporate the variety of non-Coulombic interactions that involve the neutral solvent molecules and the ionic species in the solution. These effects can be extremely important for the properties of electrolyte solutions and EDLs.49−51 The solvation interactions can be relatively simply taken into account in the framework of an LJ solvent. Clearly this is an approximation for the case of polar solvents such as water, but it is a reasonable one, provided that the ionic concentration is low to moderate.28,29 We have recently studied the effect of weak and moderate solvent−solvent and solvent−ion interactions and found that they have a significant effect of properties of charged interfaces.30 Below we offer new results that extend the analysis to the strong ionic solvation limit. 3.2.1. Solvation of the Potential-Determining Ions. Figure 2 shows the variation of local PDI density with distance and the strength of their solvation energy (i.e., the LJ parameter, ϵs−PDI). As the PDIs become more solvated, their density near the interface with the substrate decreases and the peaks move farther into the solution. The reduction of PDI concentration shifts the chemical equilibrium (eq 17) toward a more negative surface. (See also eq 19 and Figure 3a,b.) This leads to a strong attraction of positive charge in the fluid in the immediate vicinity of the negatively charged interface, as shown by the

Figure 2. PDI density profiles for varying solvation LJ energy ϵs−PDI/ kBT.

peaks in the fluid charge density distribution in Figure 3c. The charge distribution in the EDL is naturally accompanied by an increase in the potential distribution (Figure 3d). The surface charge and potential exhibit a particularly abrupt change at around ϵs−PDI/kBT ≈ 3.4. This change correlates with a similarly abrupt shift in the PDI density distributions shown in Figure 2. For the selected parameters listed above, σe/d2 = 0.664 corresponds to a fully ionized surface (Figure 3a). The surface potential magnitudes that correspond to a fully ionized surface are nonphysically high, which implies that such a scenario is unlikely. Still, it clearly shows how the surface responds to the extreme solvation of the PDIs through the surface charge regulation mechanism. The electrostatic potential distributions in the EDL (Figure 3d) exhibit a linear change between z/d = 0 and 0.5. The plane z/d = 0.5 corresponds to the location of the center of charges in the first layer of fluid next to the charged interface (Figure 3c). This plane defines the location of the Stern layer.33 More discussion of the Stern layer in the framework of cDFT applied to charge-regulating EDLs is provided below. (See also refs 29−31.) 3.2.2. Solvation of the Positive Ions. The presence of background electrolyte (consisting of ions that are not directly involved in chemical and potential-determining interactions with the interface) may still have a strong effect on the properties of EDLs. These ions may still be involved in physical (Coulombic or LJ) interactions with the surface as well as with all species in the fluid including the PDIs. The positive ions carry the same charge as the PDIs; therefore, they are attracted by the negative surface. Hence, they compete for space in the Stern layer, which affects the surface charge by shifting the equilibrium (eq 17). The density profiles for the positive ions near the charged interface and for different solvation energies ϵs−Pos are shown in Figure 4. The positive ion density in the Stern layer decreases with the increase in the solvation energy. For ϵs−Pos ≳ 2kBT, the positive ions are practically absent from the layer in the immediate vicinity of the charged surface. Instead, they form layers of ions farther into the fluid as indicated by the structural peaks. This observation is in qualitative agreement with X-ray measurements of the effect of hydration on the local ionic densities near muscovite−water interfaces.52 The solvation of the positive ions has an effect on the charge and potential at the surface and in the fluid (Figure 5). For ϵs−Pos = kBT, the interaction is the same as for the remaining fluid components, including the PDIs. The positive ions and the PDIs tend to occupy the Stern layer, and because of their identical charge, they experience mutual repulsion. Hence, the positive ions create a barrier for the PDIs, which is both steric and electrostatic, and reduce their concentration in the layer 13813

DOI: 10.1021/acs.langmuir.8b02453 Langmuir 2018, 34, 13808−13820

Langmuir

Invited Feature Article

Figure 3. Effect of the PDI solvation on the (a) surface charge, (b) surface potential, (c) total charge distribution in the fluid near the charged surface (the inset provides an enlarged view), and (d) potential distributions in the fluid (the potentials become more negative with ϵs−PDI; see panel b). The linear potential regions between z/d = 0 and 0.5 correspond to the distance between the charged interface and the centers of fluid molecules in the first adjacent layer. See also the first peak in the charge distribution in panel c. Adapted from J. Chem. Phys. 2017, 147, 214704, DOI: 10.1063/1.5005060.

3.2.3. Solvation of the Negative Ions. Figure 6 shows the density profile for the negative ions in the EDL for varying solvation energy ϵs−Neg. Similar to the cases for the PDIs and the positive non-PDIs above, the local negative ion density shifts away from the surface into the fluid as the solvation energy increases. The overall effect of the negative ions’ solvation on the surface charge and potential, as well as on their distributions in the fluid, is relatively weak (Figure 7). The surface charge (Figure 7a) becomes slightly more negative as the negative ions become more solvated. At the same time, the electrostatic surface potential (Figure 7b) in the fluid becomes (again slightly) less negative. Consequently, positive ions move in close to the surface, which results in a net positive fluid charge in the Stern layer and beyond (Figure 7c). The potential distribution in the EDL (Figure 7d) follows the trend determined by the surface values and the screening due to the dissolved ions. The linear portions of the curves between z/d = 0 and 0.5 are again due to the potential change inside the Stern layer. 3.3. Effects Due to Solvent−Surface Interactions. The interface between a charged wall and the fluid in an EDL can be considered to be a system component similar to those in the fluid. It interacts with all solution species, including the solvent. Here we describe the effect of long-ranged solvent interaction with the wall in the framework of the LJ model (eq 15). Below we examine the effect of the solvent−wall interactions of the LJ type, while the ions experience only Coulombic and excluded volume forces. An analysis of the effects of the LJ ion−wall attractive interactions and their competition with the solvent on the properties of EDLs is presented elsewhere.31 Figure 8 shows the density profiles for the solvent near the charged interface. The different curves are for varying attractive

Figure 4. Positive ion density profiles for varying LJ solvation energy ϵs−Pos/kBT.

next to the surface. (See also refs 29 and 30.) This leads to a highly negative charge at the surface. As the positive ions become more solvated, they are removed from the layer that is in immediate contact with the surface. Hence, the PDIs have better access to the reactive groups and reduce the negative charge by means of the charge regulation reaction (eq 17), as shown in Figure 5a. However, as the ions are removed from the region near the charged wall, the potential becomes more negative (Figure 5b). The charge distribution in the fluid is shown in Figure 5c. As the solvation of the positive ions increases (with ϵs−Pos), the Stern layer becomes more populated with negative ions. Farther away from the surface, however, the fluid charges exhibit positive peaks. Hence, the solvation effects lead to strong structuring of both the liquid and the charge near the surface. The electrostatic potential distribution in the EDL again exhibits a linear portion between z/d = 0 and 0.5 and monotonically increases beyond the Stern layer (Figure 5d). 13814

DOI: 10.1021/acs.langmuir.8b02453 Langmuir 2018, 34, 13808−13820

Langmuir

Invited Feature Article

Figure 5. Effect of positive ions’ solvation on the (a) surface charge, (b) surface potential, (c) total charge distribution in fluid near the charged surface, and (d) potential distributions in the fluid (the potentials become more negative with ϵs−Pos; see panel b). The linear potential regions between z/d = 0 and 0.5 correspond to the distance between the charged interface and the centers of fluid molecules in the first adjacent layer. (See also the first peak in the charge distribution in panel c.) Adapted from J. Chem. Phys. 2017, 147, 214704, DOI: 10.1063/1.5005060.

(mostly positive) ions form structured layers, which are again due to the attractive and long-ranged Coulombic force (Figure 9c, eq 15). The potential distribution curves in the EDL and for different values of ϵs−w are plotted in Figure 9d. The potentials are negative and exhibit a well-defined linear region between z/ d = 0 and 0.5. A close inspection of the potential curves reveals that the changes between the regions corresponding to the position of the main charge peaks in Figure 9c are also very close to linear. These linear portions of the potential distributions correspond to the changes between the layer ions near the surface. 3.4. Molecular Interpretation of the Stern Layer and Helmholtz Planes. It is useful to relate the results obtained by the cDFT approach to established concepts from traditional electrochemistry and colloid and interface science. According to Stern,33 the first layer that is adjacent to the charged surface has special significance and merits a separate analysis. This physical representation is very similar to the first Helmholtz model of the EDL.3 Stern, however, considered the region outside of the adjacent layer to be described by the Gouy− Chapman model.4−6 Hence, the Stern layer is modeled as an adsorption, while the remaining EDL is analyzed by means of the Poisson−Boltzmann equation (i.e., the combination of eqs 1 and 2). Grahame53 refined the concept by defining a plane where the physically adsorbed or covalently bound ions are located, often referred to as the inner Helmholtz plane. Fully solvated ions form the next, outer Helmholtz plane. Clearly, these ideas are oversimplified. While correctly assuming that the ion dimensions are important, they are not taken into account in the diffuse part of the EDL. In addition, factors such as surface roughness or ion-size variations do not allow the placement of the ions on well-defined planes. Finally, the effect

Figure 6. Negative ion density profiles for varying LJ solvation energy ϵs−Neg/kBT.

energies with the wall. The first peak in the solvent density exhibits a noticeable increase with LJ energy parameter ϵs−w. The second, third, and remaining distant peaks also become more pronounced, which indicates that the fluid becomes more structured (i.e., layered) near the wall. This has an effect on the local density distribution of all ionic species, which in turn affects the charges and potentials both at the surface and in the fluid as shown in Figure 9. The surface charge (Figure 9a) increases with the local density of the solvent molecules, which displace and effectively dilute the concentration of PDIs. Thus, the surface equilibrium (eq 16) is shifted toward greater dissociation of the reactive groups. The surface potential (Figure 9b) also becomes more negative with LJ parameter ϵs−w. The charge distribution in the fluid is shown in Figure 9c. Interestingly, it indicates some positive charge density in the layer adjacent to the wall but a lot more charge in the next few layers. Hence, while the ions are mostly displaced near the surface, some remain there, attracted by the Coulombic force exerted by the oppositely charged interface. As the solvent density drops with the distance from the wall (Figure 8), the 13815

DOI: 10.1021/acs.langmuir.8b02453 Langmuir 2018, 34, 13808−13820

Langmuir

Invited Feature Article

Figure 7. Effect of the negative ion solvation on the (a) surface charge, (b) surface potential, (c) total charge distribution in the fluid near the charged surface, and (d) potential distributions in the fluid (the potentials become less negative with ϵs−Neg; see panel b). The linear potential regions between z/d = 0 and 0.5 correspond to the distance between the charged interface and the centers of fluid molecules in the first adjacent layer. (See also the first peak in the charge distribution in panel c.) Adapted from J. Chem. Phys. 2017, 147, 214704, DOI: 10.1063/1.5005060.

support such a distinction. An oxide surface in contact with aqueous solution will regulate its charge by binding or releasing protons. However, these are usually carried around as hydronium ions, which are much larger. Therefore, the plane the hydronium ions reach as they touch the surface is not the same as the reaction plane, where the chemical interaction involving the much smaller hydrogen ion occurs. In addition, the layered structure can persist for more than a single layer of ions and molecules. Hence, the realistic physical picture can be more complex that the simple Stern, or even Graham, model and is entirely governed by the particular type and strength of the interactions between all solutions species and with the charged surface. (See also refs 30 and 31.) In light of this discussion, it is useful to consider the relationship between the charge and potential distributions in the fluid. While the fluid charge distributions exhibit multiple peaks near the surface, the respective potential distributions are mostly monotonic. This is not surprising since eq 1 states that the charge distribution is related to the curvature of the potential function. Hence, as long as the charge distribution does not change sign, the curvature of the electrostatic potential distribution functions will not change either. Obviously, a constant curvature requirement prevents the potential distribution from exhibiting the oscillatory behavior that is observed for the charge or the density profiles. There are examples, as in Figure 5c above, where the fluid charge density changes sign and yet the corresponding potential curves, shown in Figure 5d, remain monotonic. (More examples can be found in refs 30 and 31.) A close inspection of the potential curves, however, shows that there are variations in the curvatures, but their magnitude is too small to be visible on this scale. This is also why the potential curves do not exhibit minima and maxima like the fluid charge distributions. In fact, twice differentiating the potential curves

Figure 8. Solvent density profiles for varying LJ solvent-charged wall interactions ϵs−w/kBT.

of the solvent is underestimated because it plays a dominant role in the local structuring near the charged interface. A statistical mechanical model with molecular resolution, like the cDFT discussed and applied in this article, does not require an ad hoc introduction of particular structure such as layers and planes. The molecular organization of the solution in the EDL comes out naturally as a result of the minimization of the thermodynamic potential (eq 10). Still, a layered structure of all of the species (charged ions and neutral solvent) emerges, undoubtedly enhanced by the fact that we also assumed a molecularly smooth interface between the substrate and solution54 as well as because all solution species are assumed to have the same size. Our analysis suggests that indeed there are certain important planes where species tend to locate. The first one is the plane where a surface charge regulation chemical reaction occurs. The model we used above distinguishes this plane from the one where the physical adsorption of ions or solvent molecules may occur. This is different from Grahame’s representation, where both chemical bonding and physical bonding occur on the same plane. There are common real-world examples which 13816

DOI: 10.1021/acs.langmuir.8b02453 Langmuir 2018, 34, 13808−13820

Langmuir

Invited Feature Article

Figure 9. Effect of the solvent wall−surface interactions on the (a) surface charge, (b) surface potential, (c) total charge distribution in the fluid near the charged surface, and (d) potential distribution in the fluid (the potentials become more negative with ϵs−w; see panel b). The linear potential regions between z/d = 0 and 0.5 correspond to the distance between the charged interface and the centers of fluid molecules in the first adjacent layer. (See also the first peak in the charge distribution in panel c.) Adapted from J. Chem. Phys. 2018, 148, 044702, DOI: 10.1063/ 1.5012090.

agreement with the results for the fluid charge. (See also eq 1.) Note that eq 1 is as accurate as charge density distribution ρe on its right-hand side. Combining eq 1 with eq 2 leads to the rather approximate Poisson−Boltzmann model. Having a more accurate distribution that accounts for size and molecular and ionic interactions leads to a much more realistic relationship between the fluid charge and potential. Such a framework allows one to analyze the entire EDL as a single entity instead of rather artificially breaking it into layers that are governed by different models, as suggested by Stern33 and Grahame.53

in Figure 5d recovers the charge density curves in Figure 5c. Another way to demonstrate the curvature variation in the potential distribution is by subtracting the linear contribution, which dominates for small z/d, and examining the nonlinear residual, ΔΨNL, defined by ÄÅ ÉÑ ÅÅÅ ÑÑ i y ∂Ψ j z zz z ÑÑÑÑ ΔΨNL(z) = Ψ(z) − ÅÅÅΨs + jj k ∂z { z = 0 ÑÑÑÖ ÅÅÅÇ (22) The quantity ΨNL is plotted in Figure 10. The plotted curve for ϵs−Pos = 2kBT, which is shown in Figure 5c, corresponds to

4. CONCLUSIONS Electric double layers are complex and diverse systems which present many interesting fundamental physical and chemical problems. At the same time, they appear as important elements in many established and emerging practical applications. Therefore, it is not surprising that electric double layers have been extensively studied for more that a century. Still, quite often the theoretical models that are used are rather oversimplified. In a recent Feature Article, Trefal at al.2 presented a strong case in favor of using a surface charge regulation condition at the substrate−electrolyte interface instead of setting a constant potential or charge as an external parameter. In this article, we argue that an account of the precise solution structure is extremely important. Furthermore, we demonstrate (see also refs 28−31) that this structure is largely determined by the neutral solvent, despite the fact that it does not participate in any Coulombic interactions. Typically, the number of solvent molecules greatly exceeds that of the ions, which is the reason for its strong effect. The charged ionic species are forced to assume locations in a structural matrix created by the dense solvent. That is why

Figure 10. Nonlinear residual of the potential distribution (eq 22) in the EDL for ϵs−Pos = 2kBT.

a case where the fluid charge density exhibits a negative minimum between z/d = 0.5 and 1.35. For z/d > 1.35, the charge density is always positive. The plot in Figure 10 shows that the potential curvature is positive between z/d = 0.5 and 1.35 and negative at larger distances. Hence, the curvature variations in the potential distribution function are in perfect 13817

DOI: 10.1021/acs.langmuir.8b02453 Langmuir 2018, 34, 13808−13820

Langmuir

Invited Feature Article

sometimes useful, are usually oversimplified, and the real physical picture is more complex. The charge layering has been known for a few decades, and attempts were made to incorporate it into a theoretical model. Many of these models, however, ignore the effect of the solvent contribution to the solution structure. Some of these models are based on modifying the Poisson−Boltzmann model17−19,55,56 or employed integral equations approaches while assuming that the solvent is a structureless continuum.57−62 Other statistical− mechanical analyses did include the solvent molecular contributions but did not account for the surface charge regulation.21−26,63 The message of the present study is that the solvent structure and the surface chemistry are coupled, and taking one into account while ignoring the other may result in missing crucially important properties of electric double layers.

their density profiles near the surface resemble a liquidlike structure while their actual concentration is gaslike. The solution structure near the charged interface affects the surface chemistry, and hence there is strong coupling between the solvent-induced effects and the charge regulation. (For more details, see refs 28 and 29.) The surface charge then determines the electrostatic force experienced by the ions in the electric double layer. The intermolecular interactions inside an electrolyte solution are complex. Traditionally, in an effort to improve our understanding of solution phenomena we break these interactions up into a variety of contributions (e.g., we speak of ionic, polar, dipolar, induced dipolar, and hydrogen bonding interactions). It is understood, but not always stressed, that the division of these contributions is arbitrary. That is, experimentally there are no ways to consider physical quantities, such as the heat of dissolution of a salt, and assign the ionic and polarization contributions to that quantity. The situation is different for modeling approaches, such as MD and cDFT, because there we can specify a Hamiltonian and allow ourselves to perform thought experiments. Traditionally, we attempt to keep the design of a Hamiltonian as simple as possible in order to facilitate the understanding and physical intuition that is to be generated by the study of model systems. The celebrated van der Waals model for fluids and fluid mixtures is a good example in this context. It consists of just two kinds of contributions: excluded volume (hard sphere interactions) and mean field attractions. Despite its simplicity, the van der Waals model has been extremely effective in the description of phase diagrams and interfacial phenomena. However, the study of critical phenomena demonstrated that the values of the critical exponents were not quite right. Our approach to electrolyte solutions proceeds in a similar vein. Clearly, we recognize that the inclusion of the solvent introduces a plethora of interactions, as indicated above. However, at this stage we have sought to keep the computational challenges manageable by simplifying the interactions. It is in this spirit that we consider representing the solvent by a simple Lennard-Jones interaction, which is therefore to be considered an effective representation of water. This, we believe, captures the most important aspects: excluded volume and attractive interactions. The ion−solvent interactions are generally on the order of the thermal energy. Poorly solvated ions have a parameter that is less than kBT, while the parameter for strongly solvated ions is larger than kBT. This is system-specific, and that is why we covered a range of values to examine how different cases would affect properties such as the charge and potential distributions. A model that explicitly includes the solvent at the molecular level allows one to examine the effect of solvation nonCoulombic interactions. These are very common and are responsible for the variations in the electric double layer properties containing ions with the same ionic charge but being different otherwise.49−51 In this article, we demonstrated the power of cDFT to model the effects of the solvation interactions in combination with the surface chemical equilibrium on the charge and potential distributions at the surface and in the fluid. More examples are available elsewhere.30,31 A detailed molecular model of the electric double layer offers a better physical interpretation of established concepts such as the Stern layer and the Helmholtz planes. An elaborate molecular theory demonstrates that these concepts, while



AUTHOR INFORMATION

Corresponding Author

*E-mail: [email protected]. ORCID

Dimiter N. Petsev: 0000-0003-3257-4853 Notes

The authors declare no competing financial interest. Biographies

Frank B. van Swol is a research professor in the Department of Chemical and Biological Engineering at the University of New Mexico (Albuquerque, NM, United States). He received his Ph.D. in physical chemistry from the University of Amsterdam in 1988. He was a Ramsay Memorial Fellow at the Physical Chemistry Laboratory at Oxford University (Oxford, U.K.). He then moved to the United States as a visiting faculty member in the Chemical Engineering Department of Cornell University and subsequently joined the Chemical Engineering Department at the University of Illinois at UrbanaChampaign. He joined Sandia National Laboratories in 1994 and retired from there in 2015. His research focus concentrates on interfacial phenomena including colloids and electrolytes. He is also interested in additive manufacturing processes such as selective laser melting. 13818

DOI: 10.1021/acs.langmuir.8b02453 Langmuir 2018, 34, 13808−13820

Langmuir

Invited Feature Article

Dissimilar Amphoteric Surfaces. J. Chem. Soc., Faraday Trans. 1 1976, 72, 2844−2865. (11) Chan, D. Y. C.; Healy, T. W.; Supasiti, T.; Usui, S. Electrical double layer interactions between dissimilar oxide surfaces with charge regulation and Stern-Grahame layers. J. Colloid Interface Sci. 2006, 296, 150−158. (12) Behrens, S. H.; Grier, D. G. The charge of glass and silica surfaces. J. Chem. Phys. 2001, 115, 6716−6721. (13) Borkovec, M. Origin of 1-pK and 2-pK Models for Ionizable Water-Solid Interfaces. Langmuir 1997, 13, 2608−2613. (14) Behrens, S. H.; Borkovec, M. Electric double layer interaction of ionizable surfaces: Charge regulation for arbitrary potentials. J. Chem. Phys. 1999, 111, 382−385. (15) Behrens, S. H.; Borkovec, M. Electrostatic Interaction of Colloidal Surfaces with Variable Charge. J. Phys. Chem. B 1999, 103, 2918−2928. (16) Popa, I.; Sinha, P.; Finessi, M.; Maroni, P.; Papastavrou, G.; Borkovec, M. Importance of Charge Regulation in Attractive DoubleLayer Forces between Dissimilar Surfaces. Phys. Rev. Lett. 2010, 104, 228301. (17) Borukhov, I.; Andelman, D.; Orland, H. Steric effects in electrolytes: A modified Poisson-boltzmann equation. Phys. Rev. Lett. 1997, 79, 435−438. (18) Ben-Yaakov, D.; Andelman, D.; Podgornik, R. Beyond standard Poisson-Boltzmann theory: Ion-specific interactions in aqueous solutions. J. Phys.: Condens. Matter 2009, 21, 424106. (19) Ben-Yaakov, D.; Andelman, D.; Podgornik, R.; Harries, D. Ionspecific hydration effects: Extending the Poisson-Boltzmann theory. Curr. Opin. Colloid Interface Sci. 2011, 16, 542−550. (20) Attard, P. Electrolytes and the Electric Double Layer. Adv. Chem. Phys. 2007, 92, 1−159. (21) Adelman, S. A.; Deutch, J. M. Exact solution of the mean spherical model for strong electrolytes in polar solvents. J. Chem. Phys. 1974, 60, 3935−3949. (22) Chan, D. Y. C.; Mitchell, D. J.; Ninham, B. W.; Pailthorpe, B. A. On the theory of dipolar fluids and ion-dipole mixtures. J. Chem. Phys. 1978, 69, 691−696. (23) Blum, L.; Henderson, D. Mixtures of hard ions and dipoles against a charged wall: The Ornstein-Zernike equation, some exact results, and the mean spherical approximation. J. Chem. Phys. 1981, 74, 1902−1910. (24) Dong, W.; Rosinberg, M. L.; Perera, A.; Patey, G. N. A theoretical study of the solid-electrolyte solution interface. I. Structure of a hard sphere ion-dipole mixture near an uncharged hard wall. J. Chem. Phys. 1988, 89, 4994−5009. (25) Groot, R. D. Density-functional theory for inhomogeneous fluids. Phys. Rev. A: At., Mol., Opt. Phys. 1988, 37, 3456−3464. (26) Perera, L.; Berkowitz, M. L. Dynamics of ion solvation in a Stockmayer fluid. J. Chem. Phys. 1992, 96, 3092−3101. (27) Frink, L. J. D.; van Swol, F. Solvation Forces and Colloidal Stability - A Combined Monte-Carlo and Density-Functional Theory Approach. J. Chem. Phys. 1994, 100, 9106−9116. (28) Fleharty, M. E.; van Swol, F.; Petsev, D. N. Solvent Role in the Formation of Electric Double Layers with Surface Charge Regulation: A Bystander or a Key Participant? Phys. Rev. Lett. 2016, 116, 048301. (29) Vangara, R.; Brown, D. C. R.; van Swol, F.; Petsev, D. N. Electrolyte solution structure and its effect on the properties of electric double layers with surface charge regulation. J. Colloid Interface Sci. 2017, 488, 180−189. (30) Vangara, R.; van Swol, F.; Petsev, D. N. Solvation effects on the potential and charge distributions in electric double layers. J. Chem. Phys. 2017, 147, 214704. (31) Vangara, R.; van Swol, F.; Petsev, D. N. Solvophilic and solvophobic surfaces and non-Coulombic surface interactions in charge regulating electric double layers. J. Chem. Phys. 2018, 148, 044702. (32) Lee, J. W.; Nilson, R. H.; Templeton, J. A.; Griffiths, S. K.; Kung, A.; Wong, B. M. Comparison of Molecular Dynamics with Classical Density Functional and Poisson?Boltzmann Theories of the

Dimiter N. Petsev is a professor in the Department of Chemical and Biological Engineering at the University of New Mexico. He received his Ph.D. in physical chemistry at the University of Sofia (Sofia, Bulgaria) in 1996. After graduating from the University of Sofia, he spent a year and a half as postdoctoral fellow at Purdue University, followed by research positions at the Center for Microgravity and Materials Research at the University of Alabama in Huntsville and the Chemical Engineering Department at the University of Houston. His research interests include surface and interfacial phenomena, structure and dynamics of charged colloidal suspensions, soft matter, nanomaterials, electrokinetic phenomena, and capillary flows.



ACKNOWLEDGMENTS This research was supported by the Air Force Academy under cooperative agreement FA7000-14-2-0017 through funding from the Corrosion Policy and Oversight Office and the United States Department of Energy, Office of Basic Energy Sciences, Division of Materials Sciences and Engineering. The authors are indebted to Professor David G. Whitten for his encouragement and comments on an earlier version of the manuscript. We are also thankful to the UNM Center for Advanced Research Computing for providing support and access to computational resources.



REFERENCES

(1) Dukhin, S. S. In Surface and Colloid Science; Matijevic, E., Ed.; Wiley Interscience: New York, 1974; Chapter 1, pp 1−49. (2) Trefalt, G.; Behrens, S. H.; Borkovec, M. Charge Regulation in the Electrical Double Layer: Ion Adsorption and Surface Interactions. Langmuir 2016, 32, 380−400. (3) Helmholtz, H. Ueber einige Gesetze der Vertheilung elektrischer Ströme in körperlichen Leitern mit Anwendung auf die thierischelektrischen Versuche. Ann. Phys. 1853, 165, 211−233. (4) Gouy, G. Sur la constitution de la charge electrique ‘a la surface dun electrolyte. J. Phys. Theor. Appl. 1910, 9, 457−468. (5) Gouy, G. Sur la fonction lectrocapillaire. Ann. Phys. 1917, 9, 129−184. (6) Chapman, D. L. A contribution to the theory of electrocapillarity. Philos. Mag. 1913, 25, 475−481. (7) Debye, P.; Huckel, E. Zur theorie der elektrolyte. I. Gefrierpunktserniedrigung und verwandte erscheinungen. Phys. Z. 1923, 24, 185−206. (8) Ninham, B. W.; Parsegian, V. A. Electrostatic Potential between Surface Bearing Ionizable Groups in Ionnic Equilibrium with Physiologic Saline Solution. J. Theor. Biol. 1971, 31, 405−428. (9) Ettelaie, R.; Buscall, R. Electrical double layer interactions for spherical charge regulating colloidal particles. Adv. Colloid Interface Sci. 1995, 61, 131−160. (10) Chan, D. Y. C.; Healy, T. W.; White, L. R. Electrical Double Layer Interactions under Regulation by Surface Ionization Equilibria13819

DOI: 10.1021/acs.langmuir.8b02453 Langmuir 2018, 34, 13808−13820

Langmuir

Invited Feature Article

Electric Double Layer in Nanochannels. J. Chem. Theory Comput. 2012, 8, 2012−2022. (33) Stern, O. Zur theorie der elektrolytischen doppelschicht. Z. Electrochem 1924, 30, 508−516. (34) Evans, R. In Fundamentals of Inhomogeneous Fluids; Henderson, D., Ed.; Marcel Dekker: New York, 1992; Chapter 3, pp 85−175. (35) Wu, J. Density Functional Theory for Chemical Engineering: From Capillarity to Soft Materials. AIChE J. 2006, 52, 1169−1191. (36) Rosenfeld, Y. Free-energy model for the inhomogeneous hardsphere fluid mixture and density functional theory of freezing. Phys. Rev. Lett. 1989, 63, 980−983. (37) Roth, R. Fundamental measure theory for hard-sphere mixtures: a review. J. Phys.: Condens. Matter 2010, 22, 063102. (38) Archer, A. J.; Chacko, B.; Evans, R. Spinodal decomposition in a Lennard-Jones fluid. arXiv:1706.08744v1, 2017. (39) Morse, P.; Feshbach, H. Methods of Theoretical Physics; Feshbach Publishing, 1953; Vol. 1. (40) Hill, T. L. An Introduction to Statistical Thermodynamics; Dover: New York, 1986. (41) Carnie, S.; Chan, D. Y. C. The structure of electrolytes at charged surfaces: Ion-dipole mixtures. J. Chem. Phys. 1980, 73, 2949− 2957. (42) Ballenegger, V.; Hansen, J.-P. Dielectric permittivity profiles of confined polar fluids. J. Chem. Phys. 2005, 122, 114711. (43) Gillespie, D.; Nonner, W.; Eisenberg, R. S. Density functional theory of charged, hard-sphere fluids. Phys. Rev. E: Stat. Phys., Plasmas, Fluids, Relat. Interdiscip. Top. 2003, 68, 031503. (44) Barthel, J. M. G.; Krienke, H. Physical Chemistry of Electrolyte Solutions: Modern Aspects; Springer: Darmstadt, 1998. (45) Valisko, M.; Boda, D.; Gillespie, D. Selective adsorption of ions with different diameter and valence at highly charged interfaces. J. Phys. Chem. C 2007, 111, 15575−15585. (46) Bozic, A. L.; Podgornik, R. Anomalous Multipole Expansion: Charge Regulation of Patchy Inhomogeneously Charged Spherical Particles. J. Chem. Phys. 2018, 149, 163307. (47) Sandia National Laboratories Tramonto Software: https:// software.Sandia.Gov/tramonto/. (48) van Hal, R. E. G.; Eijkel, J. C. T.; Bergveld, P. A general model to describe the electrostatic potential at electrolyte oxide interfaces. Adv. Colloid Interface Sci. 1996, 69, 31−62. (49) Tadros, T. F.; Lyklema, J. Adsorption of potential-determining ions at the silica-aqueous electrolyte interface and the role of some cations. J. Electroanal. Chem. Interfacial Electrochem. 1968, 17, 267− 275. (50) Kunz, W.; Henle, J.; Ninham, B. W. The present state of affairs with Hofmeister effects. Curr. Opin. Colloid Interface Sci. 2004, 9, 1− 18. (51) Kunz, W.; Henle, J.; Ninham, B. W. ’Zur Lehre von der Wirkung der Salze’ (about the science of the effect of salts): Franz Hofmeister’s historical papers. Curr. Opin. Colloid Interface Sci. 2004, 9, 19−37. (52) Lee, S. S.; Fenter, P.; Park, C.; Sturchio, N. C.; Nagy, K. L. Hydrated Cation Speciation at the Muscovite (001)-Water Interface. Langmuir 2010, 26, 16647−16651. (53) Grahame, D. C. The Electrical Double Layer and the Theory of Electrocapillarity. Chem. Rev. 1947, 41, 441−501. (54) Frink, L. J. D.; van Swol, F. Solvation forces between rough surfaces. J. Chem. Phys. 1998, 108, 5588−5598. (55) Biesheuvel, P. M.; v.d. Veen, M.; Norde, W. A modified Poisson-Boltzmann model including charge regulation for the adsorption of ionizable polyelectrolytes to charged interfaces, applied to lysozyme adsorption on silica. J. Phys. Chem. B 2005, 109, 4172− 4180. (56) López-García, J. J.; Horno, J. Poisson-Boltzmann Description of the Electrical Double Layer Including Ion Size Effects. Langmuir 2011, 27, 13970−13974. (57) Plischke, M. Pair correlation functions and density profiles in the primitive model of the electric double layer. J. Chem. Phys. 1988, 88, 2712−2718.

(58) Yu, Y. X.; Wu, J. Z.; Gao, G. H. Density-functional theory of spherical electric double layers and zeta potentials of colloidal particles in restricted-primitive-model electrolyte solutions. J. Chem. Phys. 2004, 120, 7223−7233. (59) Tang, Z.; Striven, L. E.; Davis, H. T. Interactions between primitive electrical double layers. J. Chem. Phys. 1992, 97, 9258−9266. (60) Pizio, O.; Sokolowski, S. On the effects of ion-wall chemical association on the electric double layer: A density functional approach for the restricted primitive model at a charged wall. J. Chem. Phys. 2006, 125, 024512. (61) Jiang, J.; Cao, D.; Henderson, D.; Wu, J. Revisiting density functionals for the primitive model of electric double layers. J. Chem. Phys. 2014, 140, 044714. (62) Blum, L.; Henderson, D. In Fundamentals of Inhomogeneous Fluids; Henderson, D., Ed.; : New York, 1992; Chapter 6, pp 239− 276. (63) Grimson, M. J.; Rickayzen, G. Forces between Surfaces in Electrolyte Solutions. Chem. Phys. Lett. 1982, 86, 71−75.

13820

DOI: 10.1021/acs.langmuir.8b02453 Langmuir 2018, 34, 13808−13820