Surface Electrostatic Potential and Water Orientation in the presence

Article pubs.acs.org/Langmuir

Surface Electrostatic Potential and Water Orientation in the presence of Sodium Octanoate Dilute Monolayers Studied by Means of Molecular Dynamics Simulations Kalil Bernardino§ and André F. de Moura*,§ §

Departamento de Química, Universidade Federal de São Carlos, Rodovia Washington Luiz km 235, CP 676, CEP 13565-905, São Carlos, SP Brasil S Supporting Information *

ABSTRACT: A series of atomistic molecular dynamics simulations were performed in the present investigation to assess the spontaneous formation of surfactant monolayers of sodium octanoate at the water−vacuum interface. The surfactant surface coverage increased until a saturation threshold was achieved, after which any further surfactant addition led to the formation of micellar aggregates within the solution. The saturated films were not densely packed, as might be expected for short-chained surfactants, and all films regardless of the surface coverage presented surfactant molecules with the same ordering pattern, namely, with the ionic heads toward the aqueous solution and the tails lying nearly parallel to the interface. The major contributions to the electrostatic surface potential came from the charged heads and the counterion distribution, which nearly canceled out each other. The balance between the oppositely charged ions rendered the electrostatic contributions from water meaningful, amounting to ca. 10% of the contributions arising from the ionic species. And even the aliphatic tails, whose atoms bear relatively small partial atomic charges as compared to the polar molecules and molecular fragments, contributed with ca. 20% of the total electrostatic surface potential of the systems under investigation. Although the aliphatic tails were not so orderly arranged as in a compact film, the C−H bonds assumed a preferential orientation, leading to an increased contribution to the electrostatic properties of the interface. The most prominent feature arising from the partitioning of the electrostatic potential into individual contributions was the long-range ordering of the water molecules. This ordering of the water molecules produced a repulsive dipole−dipole interaction between the two interfaces, which increased with the surface coverage. Only for a water layer wider than 10 nm was true bulk behavior observed, and the repulsive dipole−dipole interaction faded away.

■

INTRODUCTION The main characteristic of aqueous surfactant molecules is their strong tendency to become adsorbed at the different interfaces formed between liquid water and its surroundings. Interfaces as diverse as those formed between water and either air, hydrophobic liquids, or hydrophobic solids always become partially covered by surfactant monolayers, which are highly directional, with the ionic or strongly polar surfactant heads interacting directly with the water molecules while the hydrophobic tails tend to decrease their exposed area to the aqueous environment due to the hydrophobic effect, resulting in an overall reduction of the surface free energy.1,2 In the case of ionic surfactants, the formation of a film generates an electrostatic potential difference that also organizes the counterions and the solvent molecules near the interface, an effect that is described by the electrical double layer models, where the interface is usually taken as a homogeneous charged plan, the counterions as point charges, and the solvent as a continuum dielectric in the analytical derivations.1,2 Although of great importance to understanding the general aspects of the © XXXX American Chemical Society

behavior of charged interfaces as well as the stabilization of colloids in suspension by the so-called double-layer repulsion, the analytical solutions for the electrical double layer ignore the atomistic nature of both the interface and the solvent, the nonhomogeneity of the surfactant distribution at the interface, and the specific interactions between the surfactant and the counterions.1,3 Computer simulations have been used for the last few decades to study molecular-level aspects of both ionic and nonionic surfactants at the water−vapor interface4−8 as well at interfaces between water and hydrophobic surfaces such as graphene9 and trichloroethylene.10 Despite the useful structural properties provided by these and other simulations, to the best of our knowledge, the recent work of Nguyen et al.7 is the only study providing a detailed picture of the electrostatic potential across the water−vapor interface in the presence of a dilute Received: August 5, 2015 Revised: September 15, 2015

A

DOI: 10.1021/acs.langmuir.5b02904 Langmuir XXXX, XXX, XXX−XXX

Article

Langmuir

model system with pure water was simulated as well to provide a comparison framework for the electrostatic potential arising from the monolayers. The interfacial area was kept constant throughout the simulations (constant (N, V, T) ensemble). In practice, this prescription creates two vacuum−solution interfaces, similar to a planar soap bubble, which may be regarded as a good approximation to the typical air−solution interfaces in experiments because the number of collisions of gas molecules from the atmosphere on a nanometer-scale surface during tens to hundreds of nanoseconds is rather small and can be safely disregarded. Since the z axis is large enough to avoid any interaction between the two interfaces through the intervening vacuum, there are no artifacts arising from the application of periodic boundary conditions in the direction perpendicular to the interface. Bulk surfactant molecules spontaneously migrated to form dilute films at both water−vacuum interfaces, which in turn led to the buildup of a different counterion distribution in the z direction for each system. (For instance, Figure 1 presents the

surfactant monolayer. Their investigation assessed the effects of two isomeric alcohols adsorbed at the air−water interface on the electrostatic surface potential, both in the presence and in the absence of an inert electrolyte. On the other hand, there are a significant number of simulations studying the electrostatic surface potential and the counterion distribution of lipid bilayers in aqueous solution.11−13 For instance, Yi et al.14 achieved good agreement between their computer modeling results and those predicted by the standard Gouy−Chapman theory for the ion distribution across the interfacial region. Despite the similarities between the lipid aggregation in bilayers and the surfactant adsorption at air−water interfaces with regard to the distribution of either counterions or inert electrolytes, the case of a dilute interfacial film is expected to have a more complex electrostatic behavior. The exposed water molecules at the clean water−air interface should change their orientation differently from those just below surfactant-covered regions. As a consequence, the electrostatic surface potential is expected to depend strongly on the degree of surface coverage, which may span a wide range of values between the clean surface and the maximum coverage, as opposed to the rather narrow range of surface area per molecule in any typical lipid bilayer. The present investigation is aimed at filling this gap in our understanding of these complex interfacial systems by means of a series of molecular dynamics simulations of the spontaneous adsorption and self-organization of sodium octanoate at the water−vacuum interface. Several concentrations and system sizes have been studied, providing a comprehensive molecular-level picture of this interfacial model system. This manuscript is organized as follows: we will survey the protocol used for the model system setup and the molecular dynamics simulations conditions in the Methodology section. In the Results and Discussion section we will in turn present and discuss the charge density and the electrostatic potential across the water−vacuum interface as the surfactant concentration changes, splitting the latter into individual contributions arising from each chemical species comprising the system, with especial emphasis on the analysis of water orientation across the systems. Finally, we will describe a water-mediated repulsion energy between charged interfaces that arises from the reorientation of water molecules in response to the electrostatic potential at the interfaces and is amenable to an approximate analytical model.

Figure 1. Initial (top) and final (bottom) structures of the system with 40 octanoate ions (surface concentration of 0.833 molecules nm−2). The water molecules are represented by the translucent red surface.

■

initial and final structures for the N = 40 system; the interested reader is referred to the accompanying graphic animation in the Supporting Information for a typical structural evolution.) The structural relaxation was quantitatively determined using both the convergence of the energy components as well as the octanoate and sodium density profiles across the interfaces. (See Supporting Information for more details, especially Figure S1.) The time scale for the full structural relaxation increased with the number of surfactant units, varying from 5 ns for the system with 10 ions to 30 ns for the system with 50 ions. The total sampling time used for analyses (discarding the time needed for the relaxation) for each system was 9 ns for pure water, 45 ns for N = 10, 40 ns for N = 20, 35 ns for N = 30 and 40, and 40 ns for N = 50 (the latter system reached a total time of 70 ns). It is important to acknowledge the fact that sodium octanoate has a high cmc value3 and thus an appreciable solubility in water. As a consequence, it is impossible to form an

METHODOLOGY Model Systems. There were two reasons to choose sodium octanoate: first, short-chained surfactants present fast structural and energetic relaxation in bulk simulations,15−18 allowing a proper sampling of the equilibrium states to be achieved within affordable computer times, and second, sodium octanoate has a large cmc value,19 rendering realistic model systems much smaller than those describing surfactants with longer hydrophobic tails. Model systems comprised either 10, 20, 30, 40, or 50 octanoate anions along with the same number of sodium counterions and 3300 water molecules. All ions and molecules were randomly placed in the cubic boxes with a density of ca. 1 g cm−3, and then the edges of the boxes in the z direction were increased to 10 nm in order to form slabs with two liquid− vacuum interfaces and periodic boundary conditions in the xy plane, with an interfacial area of 24.01 nm2 on each side. A B


Article

Langmuir

Figure 2. Graphical representation of the octanoate tails at the water−vacuum interface. The red surface stands for exposed water whereas the bluish surfaces stand for the hydrophobic tails pointing outward the aqueous phase. The figures on the top left of each structure stand for the surface concentration of octanoate ions in each system (nm−2).

octanoate compact film even under lateral compression, leading to a considerable extension of exposed solvent, even for the more concentrated systems. It is not possible to form a compact monolayer by increasing the bulk surfactant concentration, as demonstrated by the extra simulations we performed with the number of octanoate ions up to 100, which simultaneously formed both interfacial films and micellar aggregates within the solution for all systems with N > 50. Typically, the surfactant concentration at the interfaces did not change upon addition of extra octanoate molecules, which selfassembled into micelles. These concentrated systems presented a complex, hybrid structure with both interfacial films and micellar aggregates, rendering them topologically difficult to analyze. Besides, the more concentrated systems were not thermodynamically equilibrated after 100 ns of simulation. If we consider that systems with N ⩽ 50 were fully relaxed at time scales shorter than 100 ns and that bulk sodium octanoate micelles would have relaxed as well,15−18 then it is clear that these hybrid systems have not attained equilibrium due to the film−micelle interaction. This complex behavior is beyond our present scope, so the study of these concentrated systems (N > 50) will be postponed until a forthcoming investigation. (Some qualitative findings may be found in the Supporting Information.) The appreciable solubility of octanoate anions in water allowed some degree of surfactant exchange between the interfaces and the bulk solution even for the most dilute systems. Nonetheless, after equilibrium was attained almost all octanoate ions were located at the interfaces, with nearly the same number of surfactant molecules at each water−vapor interface, as expected for a system in thermodynamic equilibrium. For instance, the system closer to interface saturation (N = 50) reached an average concentration of octanoate at the midplane of the simulation box of ca. 0.03 mol L−1 (about 10-fold lower than the cmc),19 in contrast to the much larger concentration at the interfaces of ca. 4 mol L−1 (about 10-fold higher than the cmc). Considering that all of the surfactant became adsorbed at the two interfaces and that surface area amounted to 24.01 nm2 for all systems, the surface

concentration of the surfactant was computed to be 0, 0.208, 0.417, 0.625, 0.833, and 1.024 molecules nm−2, which are typical values for dilute monolayers. We will use these figures to refer to each system instead of using the total number of octanoate ions comprising the systems since the former definition is an intensive quantity and is numerically equal to the charge density per unity area (in e nm−2) for each system, allowing one to compare the results arising from our modeling with those from future investigations of other interfacial systems with similar surface charges densities. Two larger systems were studied to assess the effects of the size on the interfacial electrostatic properties, one comprising only pure water and the other with 40 octanoate anions and 40 sodium counterions. The number of water molecules was increased to 11 600 while the interfacial area was kept fixed, increasing the separation between the two interfaces from ca. 5 nm to ca. 17 nm. The production runs reached 9 ns for pure water and 30 ns for the N = 40 system. These extra simulations are meant to provide a comparison framework for the changes in the electrostatic potential, counterion distribution, and water orientation as two charged interfaces approach each other in solution. Simulations Conditions and Parameters. All of the simulations were performed using the GROMACS suite20−22 (version 4.5.5) with the leapfrog algorithm and an integration step of 1 fs. The total volume of the system was kept constant, and the temperature was controlled with the weak coupling scheme of Berendsen23 (T = 300 K, τT = 0.1 ps). The potential energy was described using the OPLS-AA parameters24 for the octanoate anion, the Åqvist parameters25 for sodium, and the SPC model26 for water. Nonbonded interactions were truncated at 1.5 nm, with the particle-mesh Ewald correction (PME)27 for the long-range Coulomb interactions and a shift function for the Lennard-Jones interactions. The standard 3D PME has been used instead of its corrected version for slab systems since these two schemes yielded average correction energies differing typically by less than 1 kJ mol−1. These small differences also indicated that a vacuum region of 10 nm was enough to prevent any appreciable interaction between periodic C


Article

Langmuir images in the z direction. VMD (version 1.8.7)28 was used to visualize the trajectories and render the graphical representations of the systems.

■

RESULTS AND DISCUSSIONS Electrostatic Surface Potential and Charge Distribution. In a classical molecular dynamics simulation, each atom in the system has a well-defined partial charge, enabling the direct calculation of the local charge density from the distribution of each atomic type comprising the system. Each simulation box was divided in slices along the z direction, and the average charge density σ(z) was calculated for each slice along each trajectory, leading straightforwardly to the one-dimensional Poisson equation relating the electrostatic potential V(z) across the interface to the charge density σ(z) in eq 1,1,7 where ε0 is the permittivity of the vacuum. It should be borne in mind that the relative permittivity need not be included in the Poisson equation since both the solvent and the ionic screening naturally arise from their charge distributions, as discussed at length in the remaining two sections. d2V (z) σ (z ) =− 2 ε0 dz

Figure 3. Local charge density σ(z) (top panel) and electrostatic potential V(z) (bottom panel) along the direction perpendicular to the interface. The different colors stand for the octanoate surface concentration in each system (nm−2). The insets expand the regions of interest for each plot.

(1)

The Poisson equation can yield meaningful information only if the system is well equilibrated. This was not the case for the more concentrated systems (N > 50), which were not fully relaxed after 100 ns of trajectory, mostly due to the complex interaction between the micelles and the surfactant films and the slow diffusion of the clusters. With regard to the current investigation, only the more diluted systems will be considered (N≤50), which were equilibrated within the time scale of the simulations. (The interested reader is referred to the Supporting Information for a qualitative description of the systems in which micelles were observed along with the films.) Even for the systems which attained equilibrium within the time scale under consideration, the instantaneous distribution of octanoate anions at the interface was not uniform in the xy plane (Figure 2). Since the Poisson equation assumes the presence of a uniform charge density, the time average over tens of nanoseconds has to be considered, rendering the interfaces homogeneous (results not shown). The charge density and the electrostatic potential across the interface (Figure 3) have been symmetrized by taking the average charge density between each pair of slices equidistant from the center of mass of the system. (Symmetrizing σ(z) leads to symmetrical V(z) as a consequence of eq 1.) This procedure is justifiable by considering that both water−vacuum interfaces must be equivalent for any system under thermodynamic equilibrium, meaning that the small differences in charge density which were observed between the two interfaces are due to poor sampling and would eventually vanish if much longer simulations were affordable. The results for the larger systems may be found in the Supporting Information. The electrostatic potential drop of ca. 570 mV for pure water between the vacuum region and the minimum in the interfacial region is comparable to values reported by other investigators who also employed nonpolarizable water models to study pure water interfaces. The magnitude of the potential drop does not vary appreciably for different models of liquid water, except in the cases of either polarizable or quantum chemical potentials.29,30 The minimum electrostatic potential at the

interface became more negative with the increase in octanoate concentration, reaching −594 mV for the two moreconcentrated systems, with octanoate surface concentrations of 0.833 and 1.024 molecules nm−2, consistent with the saturation of the octanoate surface concentration. (The interested reader is referred to the discussion in the Supporting Information.) The electrostatic potential difference between the minimum at the interface and the bulk solution, taken as z = 0 nm (inset of Figure 3) changed linearly with the octanoate surface concentration (Figure 4), as expected from the electric double layer models for the plane capacitors, which predict potential differences varying linearly with the surface charge density.1,7 However, even in the absence of ions, there is a nonzero electrostatic potential difference between the interface and the

Figure 4. Electrostatic potential difference between the minimum at the interface (Figure 3) and the midpoint of the system (z = 0 nm) as a function of the octanoate surface concentration. The solid straight line stands for the linear regression of the data (R2 = 0.97). D


Article

Langmuir bulk (Figure 4) due to the water organization at the interface, which produced a nonzero charge distribution (Figure 3). It is advisible, if not mandatory, to assess the finite-size effects on the electrostatic properties of the systems. The larger systems had the same interfacial area and the same surfactant surface concentrations but had a nearly 3-fold-thicker water layer separating the two vacuum−water interfaces, as described at length in the Supporting Information. On the one hand, only minor changes were observed for the electrostatic potential difference between the vacuum region and the interface, in both the presence and absence of the surfactant. On the other hand, the electrostatic potential difference between the interface and the bulk, which was farther now, changed appreciably from −0.11 V for the smaller system to −0.16 V for the larger system, both with an octanoate surface concentration of 0.833 molecules nm−2. This difference may be mainly ascribed to the solvent orientation, as discussed in the Conclusions, and this is an example of the uncertainty with respect to the proper definition of where in the solution the bulk behavior is actually observed for systems with charged interfaces. Surfactant and Counterion Contributions to the Electrostatic Potential. Although partial atomic charges are ill-defined quantities both in experimental measurements and in quantum chemical calculations, classical force fields usually rely on parametrized atomic partial charges which are fitted to reproduce a few relevant properties of the systems of interest. These charges are the origin of the electrostatic properties which have been discussed in the previous section, and further insight may be obtained if the charge densities and the electrostatic potential (Figures 3 and 4) are split into contributions arising from the relevant molecular fragments comprising the systems, namely, the octanoate charged head, defined by the group CH2CO2−, the octanoate tail C6H13, and the sodium counterions (the water contribution will be discussed in the next section). The first CH2 group of octanoate is considered to be part of the charged head since its atoms bear different partial charges as compared to all of the other CH2 groups, and they contribute to the 1− charge of the head along with the carboxylate group.4 Despite the fact that octanoate is fairly water-soluble, essentially all of the surfactant remained at the interfaces, with only a few rare events of octanoate ions crossing the bulk solution after the systems reached equilibrium. As a consequence, narrow charge distributions arising from the octanoate heads were observed at both interfaces (Figure 5), which are similar to the uniformly charged planes from the electrical double layer models.17 It is well known from basic electrostatics that the space between two equally charged parallel plates must have a constant electrostatic potential, so the large electrostatic potential differences arising from the charged heads remained constant inside the aqueous solution. The large electrostatic potential differences produced by the charged head between the bulk solution and the vacuum regions are due to the nonzero charge arising from these groups, but these contributions were counteracted by the electrostatic potential difference arising from the oppositely charged counterions (Figure 6), resulting in an overall null potential in the vacuum regions outside the bulk solutions (Figure 3). With regard to the bulk region, however, the continuous distribution of sodium ions, as described by diffuse layer models, cannot fully neutralize the contributions from the octanoate anions, yielding

Figure 5. Contribution from octanoate charged heads (CH2CO2−) to the charge density (top) and electrostatic potential (bottom). The different colors stand for different octanoate surface concentrations (in molecules nm−2).

Figure 6. Contribution from sodium counterions to the charge density (top) and electrostatic potential (bottom). The different colors stand for different octanoate surface concentrations (in molecules nm−2).

nonzero electrostatic potential profiles across the interfacial regions. Despite being electrically neutral, the surfactant aliphatic tails (C6H13) assumed a preferential orientation at the interfaces in order to keep the charged heads pointing toward the solution (Figure 1), as expected for amphiphilic molecules. The ordering of the tails produced a charge distribution with a local excess of either positive charges arising from the hydrogen atoms or negative charges arising from the carbon atoms (Figure 7), which were an order of magnitude smaller than those arising from the charged species (Figures 5 and 6). The excess hydrogen atoms exposed to the vacuum regions produced an electrostatic potential drop that is proportional to the surfactant concentration in the concentration range under investigation, without any appreciable change in the tail orientation with the surfactant concentration (results not shown). The tail E


Article

Langmuir

Figure 7. Contribution from octanoate aliphatic tails (C6H13) to the charge density (top) and electrostatic potential (bottom). The different colors stand for different octanoate surface concentrations (in molecules nm−2).

Figure 8. Contribution of water to the charge density (top) and electrostatic potential (bottom). The different colors stand for different octanoate surface concentrations (in molecules nm−2).

contribution corresponded to ca. 20% of the total surface potential drop in the more concentrated system (yellow curves in Figures 3 and 7), demonstrating that, although small, the electrostatic contributions arising from the tails cannot be disregarded, consistent with the empirical models describing the electrostatic potential of monolayers as either two-layer31 or three-layer capacitors,32 both of which explicitly consider the dipole contributions from the organic monolayer−vacuum interface. On the other hand, united-atom and coarse-grained models cannot account for this electrostatic contribution from the surfactant tails, since the bonds comprising the CH2 and CH3 groups are not treated explicitly and consequently their local dipole moment contributions are also missing. For instance, Nguyen et al.7 observed an unexpected behavior for the total electrostatic potential variation with the concentration of long-chain alcohol surfactants, which they assumed to be probably an artifact of the simplified description of the aliphatic groups as united atoms. It is also important to stress that the chemical modification of the hydrophobic tails may be used to fine tune the electrostatic potential of surfactant monolayers, as already observed experimentally,33 and computer simulations might play an important role in the assessment of these effects. Contributions of Water Orientation to the Electrostatic Potential. The water contribution for the charge density and the electrostatic potential was more complex and size-sensitive, so the rest of the paper is dedicated to the definition and the calculation of a water-mediated repulsion between charged interfaces. The charge density profile for the pure water (Figure 8, brown curve) presented a positive peak, indicating excess hydrogen atoms pointing toward the vacuum regions, followed by a negative peak corresponding to the excess oxygen pointing to the liquid phase without any further electric structure inside the liquid, resulting in a sudden decrease in the electrostatic potential at the interface, followed by a constant negative potential in the bulk. In the presence of the surfactant monolayers, on the other hand, more complex patterns were observed, with two closely lying peaks of positive charge at each

interface (Figure 8). The intensity of the second positive peak depended on the surfactant concentration, and these peaks were located close to the regions of maximum probability of finding the negatively charged octanoate heads (Figure 5), suggesting the presence of water layers polarized by the surfactant monolayers, with the water hydrogen atoms oriented toward the surfactant heads. With regard to the first positive peak of charge density (Figure 8), water molecules also become polarized by the existence of the liquid−vacuum interface with water, and the intensity of these peaks was not sensitive to the surfactant concentration at the interfaces. The charge density and electrostatic potential profiles for the systems with surfactant did not form plateaus in the middle of the system (Figure 8), meaning that the thickness of the water slab between the two interfaces was not enough to ensure that the orientational correlation between solvent molecules and the interfaces was lost in a true bulk region inside the liquid phase. The minimum in the electrostatic potential arising from water at z = 0 nm for all systems with surfactant implies an inversion of the algebraic sign of the electric field produced by water in the middle of the slab. Since water is a dipolar molecule, the inversion of the direction of the electric field must arise from the inversion of the average orientation of the water molecules with respect to the z axis at z = 0 nm. The average angle θ between the dipole moment vector μ of the water molecules and the z axis provided a quantitative estimate of this reorientation of the water molecules (Figure 9a). We also computed the orientation distributions in two regions with different average orientations, taking the slices with the z coordinates ranging from −2.3 to −2.0 nm and from −1.7 to −1.4 nm, corresponding approximately to the first and second positive peaks of the charge density profiles for the systems with octanoate, respectively (Figure 9b,c). The thickness of each slab was chosen to be close to the diameter of a water molecule (ca. 0.3 nm). For z < 0, an average angle smaller than 90° stands for a preferential orientation with the hydrogen atoms pointing to the closest interface, whereas for z > 0 the same interpretation holds for an average larger than 90° (an average of 90°, for any F


Article

Langmuir

preferential orientation of water molecules with the hydrogen atoms pointing to the negatively charged interfaces. The water slab between the two octanoate-rich interfaces was not thick enough to allow the full relaxation of the water dipole moment in the z direction, so at any given position z inside the solution, the water molecules were oriented in response to the closest interface. As a consequence, a sudden inversion of the water dipole moment orientation near the middle of the system was observed, which must produce a dipole−dipole repulsion between the water molecules on different sides of the system. This water-mediated repulsion between two interfaces must have a significant impact on the repulsion between charged surfaces in water and other polar solvents, in addition to the well-known double-layer interaction. Guldbrand et al.34 reported a series of Monte Carlo simulations of charged planar interfaces interacting by means of the intervening ionic strength. Although the interfaces had the same charge, the interaction eventually became attractive due to instantaneous charge fluctuations in the two halves of the system. Each half of the model system was further divided into three parallel layers, and the forces between layers were computed, with the largest contributions arising from the interaction between the layers closer to the plane in the center of the system. These conclusions refer only to the interactions between the ionic clouds, without any treatment of the dipolar solvent, which was treated as a continuum dielectric medium. The treatment that follows is complementary to their approach in the sense that we shall provide a description of the solventmediated dipolar interaction between the planar charged interfaces. At any given time, a slice of the water slab with width dz at position z has an instantaneous resultant dipole moment μ(z) arising from the vector sum of the dipole moment of all water molecules comprising the slice. The dipole vector components which are parallel to the interface (xy plane components, μx(z) and μy(z)) are not relevant to the present discussion and should average to zero if proper sampling is carried out, whereas the component along the z direction (μz(z)) should become null only if we move far enough into the bulk solution. The dipole moment component μz(z) of any slice may be conveniently expressed as the algebraic sum of the z component of the dipole vectors of water molecules in the slice, which amounts to μwcos(θ), where μw is the singlemolecule water dipole moment, hereby assumed to be constant and equal to the dipole moment of the SPC water model in the equilibrium geometry (2.27 D),6 and θ is the angle between the water dipole moment and the z axis, as defined in Figure 9a. The dipole moment of each slice is now defined by eq 2, where the brackets stand for the instantaneous average over the molecules comprising the slice, not the time average, and Nw(z) is the number of water molecules in the slice, as computed during the simulations.

Figure 9. (a) Mean value of the angle θ between the water dipole moment and the normal to the interface as a function of the oxygen atom position along the z direction. (b) Orientation distribution for the slices with a z coordinate ranging from −2.3 to −2.0 nm and (c) from −1.7 to −1.4 nm. The different colors stand for different octanoate surface concentrations (molecules nm−2), and the dotted black curves stand for the ideal distribution (no preferential orientation).

position along the z direction, stands for the ideal distribution with no preferential alignment of the dipole moment vectors). With regard to pure water (Figure 9a brown curve), there is a preferential orientation only at the water−vacuum interfaces, with molecules pointing preferably with their hydrogen atoms toward the vacuum regions, as already indicated by the charge density profiles (Figure 8). Besides producing an average larger than 90° in the region with z coordinates from −2.3 to −2.0 nm, the distribution was narrower than the ideal distribution. This narrowing is consistent with enthalpic contributions acting on the preferential orientation of the water molecules besides the entropic forces which are responsible for the ideal distribution. On the other hand, the bulk region away from the interfaces had no preferential orientation of the dipole vectors with respect to the z axis, and the angle distribution exactly matched the ideal distribution, indicating that only entropic forces are governing the orientation of the water molecule dipole vectors. The octanoate concentration did not seem to affect the orientation distribution of the water molecules appreciably; these were directly exposed to the vacuum interface, with only a minimal shifting of the maximum probability in favor of the hydrogen atoms pointing to the vapor phase (Figure 9b), which is coherent with the small changes in the charge density and electrostatic potential at the interface (Figure 8). On the other hand, water molecules deeper in the solution presented distorted orientation distributions, and the distortions from the ideal distribution increased with the surfactant concentration at the interface (Figure 9c), consistent with the

μz (z) = μw, z (z) Nw(z) = μw cos θ Nw(z)

(2)

The interaction between the resultant dipoles of two slices at positions z1 and z2 can now be computed straightforwardly by the well-known expression for the interaction between two parallel or antiparallel dipoles (eq 3). Just like in the case of the Poisson eq (eq 1), this expression assumes that a homogeneous distribution of dipoles exists within any slice at any moment, which cannot be granted. On the one hand, this limitation renders eq 3 less accurate than the actual calculation between G


Article

Langmuir

limited width of the water slab (Figure 9a). The larger systems had exactly the same interfacial area (24.01 nm2) and surfactant concentrations (0 and 0.833 molecules nm−2), but the solventmediated repulsion between the interfaces was negligible as compared to that of the smaller systems (Table 1). The main structural difference of the larger systems was the presence of a relatively large region with randomly oriented water molecules between the two interfacial regions (results shown in the Supporting Information), which had not been observed for the smaller systems (Figure 9a).

all pairs of dipoles in the two slices, but on the other hand the interaction between these idealized average dipoles is physically sound and provides a better framework in which to discuss the solvent-mediated interaction between interfaces. In practice, the interaction energy contributions were calculated with all of the details available from our atomistic force field, including the coulomb energy between partial atomic charges and dispersion interactions, but discussions followed this orientation-based reasoning. E1,2 = −

■

μz (z1) μz (z 2) 2πε0 |z 2 − z1|3

CONCLUSIONS The role played by surface science analytical models, e.g., the DLVO theory,1 in our understanding of the interactions between charged interfaces cannot be regarded as purely historical since we still rely on them to foresee and to explain the stability of a wide range of colloidal systems. Nonetheless, the ever-increasing computing power has had (and will probably continue to have) a tremendous impact on the more detailed description of interfacial phenomena. For instance, earlier simulations had already included more realistic descriptions of the ionic species comprising the electric double layer, but only in the last couple of decades has it become feasible to study model systems with explicit solvent instead of a dielectric continuum. We have moved one step forward in the direction of the realistic modeling of charged interfaces by letting the model systems evolve spontaneously from random distributions of surfactant molecules into orderly interfacial films. The spontaneous adsorption of the anionic surfactant at water− vacuum interfaces and the counterion relaxation into a diffuse layer could be observed on time scales suitable for atomistic computer simulations, and the results provided a rich description of the electrostatic properties of dilute surfactant monolayers on the molecular level. The surfactant monolayers produced two interfacial regions with distinct orientation of the water molecules, one oriented mainly by the presence of the interface itself and the other oriented in response to the surfactant electrostatic charge. The fact that water orientation changes in response to the surface charge cannot be dealt with by models without an explicit description of the solvent molecules on the atomic level. Even coarse-grained models would fail, especially if multiple, close-lying solvation layers with varying degrees of orientation are formed. This finding from our atomistic simulations enriches the general discussions on colloidal stability, adding a new repulsive, solvent-mediated interaction between equally charged interfaces, in addition to other wellknown interactions, e.g., the double-layer repulsion1,2 and depletion forces.35,36 These dipole−dipole interactions should be considered to be short-ranged as compared to plain Coulomb interactions between ionic species, becoming effective only for surface separations in the range of a few nanometers. It is expected that for oppositely charged surfaces this dipole−dipole force should become attractive instead, due to the inversion of the relative orientation of the water molecules with respect to the interface. Besides providing a detailed description of the water orientation at charged interfaces, the simulations also demonstrated that CH2 and CH3 groups comprising the aliphatic hydrophobic tails may be responsible for at least 20% of the total electrostatic potential at the interfaces. This is an obvious drawback of coarse-grained models in general, which

(3)

Equation 3 can be applied only to separations between slices larger than the charge separation in the water dipole moment, so they cannot be used to estimate the interaction between neighboring slices. However, the repulsive dipole−dipole interaction between water molecules on opposite sides of the system may be conveniently calculated as the sum of interactions between each possible pair of slices involving one slice from the first half of the system and one from the second half (as defined by the inversion of the water dipole orientation in Figure 9a at position zinv, which may fluctuate from one system to another). Erepulsion = −

1 2πε0

∑∑

μz (z1) μz (z 2) |z 2 − z1|3

(4)

The solvent-mediated repulsion was calculated along the simulations using a slice thickness of dz = 0.2 nm. The average values are reported along with the total water−water interaction, which includes both Coulomb and dispersion interactions (Table 1). Both quantities are reported per unit area since their magnitude should scale with the interfacial area of the system. Table 1. Total Water−Water Interaction Energy (Ewater−water) and the Dipole Repulsion Energy between the Two Halves of the System (Erepulsion) octanoate surface concentration (nm−2)

Ewater−water (kJ mol−1 nm−2)

Erepulsion (kJ mol−1 nm−2)

0 0.208 0.417 0.625 0.833 1.042 0 (larger system) 0.833 (larger system)

−6580 −6413 −6259 −6127 −5993 −5876 −23 470 −22 870

0.02 3.2 11.6 24.2 31.2 49.0 0.00 0.04

The cohesive energy between water molecules becomes less negative with the surfactant surface concentration (Table 1), a fact that cannot be ascribed to the formation of the interfaces themselves since the reference system is pure water with the two water−vacuum interfaces already there. This energy decrease must be ascribed instead to distortions caused by the ionic species on the water structure, including but not limited to the water dipole polarization in the z direction. The repulsive dipole−dipole interaction between the two halves of the system, on the other hand, increased with the surfactant surface concentration, but this interaction was meaningful only for the smaller systems, whose water orientation could not reach the random pattern expected for bulk solution due to the H


Article

Langmuir

(10) Shi, W.-X.; Guo, H.-X. Structure, interfacial properties, and dynamics of the sodium alkyl sulfate type surfactant monolayer at the water/trichloroethylene interface: a molecular dynamics simulation study. J. Phys. Chem. B 2010, 114, 6365−6376. (11) Gurtovenko, A. A.; Vattulainen, I. Calculation of the electrostatic potential of lipid bilayers from molecular dynamics simulations: methodological issues. J. Chem. Phys. 2009, 130, 215107. (12) Zarzycki, P. Interfacial water screens the protein induced transmembrane voltage. J. Phys. Chem. B 2015, 119, 1474−1482. (13) Högberg, C. J.; Lyubartsev, A. P. Effect of local anesthetic lidocaine on electrostatic properties of a lipid bilayer. Biophys. J. 2008, 94, 525−531. (14) Yi, M.; Nymeyer, H.; Zhou, H.-X. Test of Gouy-Chapman theory for a charged lipid membrane against explicit-solvent molecular dynamics simulations. Phys. Rev. Lett. 2008, 101, 038103. (15) de Moura, A. F.; Freitas, L. C. G. Molecular dynamics simulation of the sodium octanoate micelle in aqueous solution: comparison of force field parameters and molecular topology effects on the micellar structure. Braz. J. Phys. 2004, 34, 64−72. (16) de Moura, A. F.; Freitas, L. C. G. Molecular dynamics simulation of the sodium octanoate micelle in aqueous solution. Chem. Phys. Lett. 2005, 411, 474−478. (17) de Moura, A. F.; Bernardino, K.; de Oliveira, O. V.; Freitas, L. C. G. Solvation of sodium octanoate micelles in concentrated urea solution studied by means of molecular dynamics simulations. J. Phys. Chem. B 2011, 115, 14582−14590. (18) Bernardino, K.; de Moura, A. F. Aggregation thermodynamics of sodium octanoate micelles studied by means of molecular dynamics simulations. J. Phys. Chem. B 2013, 117, 7324−7334. (19) Rosenholm, J. B. The structure and properties of medium-chain surfactant solutions: a case study of sodium octanoate. Adv. Colloid Interface Sci. 1992, 41, 197−239. (20) Berendsen, H. J. C.; van der Spoel, D.; van Drunen, R. GROMACS: a message-passing parallel molecular dynamics implementation. Comput. Phys. Commun. 1995, 91, 43−56. (21) Lindahl, E.; Hess, B.; van der Spoel, D. GROMACS 3.0: a package for molecular simulation and trajectory analysis. J. Mol. Mod. 2001, 7, 306−317. (22) van der Spoel, D.; Lindahl, E.; Hess, B.; Groenhof, G.; Mark, A. E.; Berendsen, H. J. C. GROMACS: fast, flexible, and free. J. Comput. Chem. 2005, 26, 1701−1718. (23) Berendsen, H. J. C.; Postma, J. P. M.; DiNola, A.; Haak, J. R. Molecular dynamics with coupling to an external bath. J. Chem. Phys. 1984, 81, 3684−3690. (24) Jorgensen, W. L.; Maxwell, D. S.; Tirado-Rives, J. Development and testing of the OPLS all-atom force field on conformational energetics and properties of organic liquids. J. Am. Chem. Soc. 1996, 118, 11225−11236. (25) Åqvist, J. Ion-water interaction potentials derived from free energy perturbation simulations. J. Phys. Chem. 1990, 94, 8021−8024. (26) Berendsen, H. J. C.; Postma, J. P. M.; van Gusteren, W. F.; Hermans, J. In Intermolecular Forces; Pullman, B., Ed.; Reidel: Dordrecht, 1981. (27) Essmann, U.; Perera, L.; Berkowitz, M. L.; Darden, T.; Lee, H.; Pedersen, L. G. A Smooth particle mesh Ewald method. J. Chem. Phys. 1995, 103, 8577−8593. (28) Humphrey, W.; Dalke, A.; Schulten, K. VMD Visual Molecular Dynamics. J. Mol. Graphics 1996, 14 (1), 33−38. (29) Kathmann, S. M.; Kuo, I-F. W.; Mundy, C. J. Electronic effects on the surface potential at vapor-liquid interface of water. J. Am. Chem. Soc. 2008, 130, 16556−16561. (30) Sokhan, V. P.; Tildesley, D. J. Molecular physics lecture: the free surface of water: molecular orientation, surface potential and nonlinear susceptibility. Mol. Phys. 1997, 92, 625−640. (31) Vogel, V.; Möbius, D. Local surface potentials and electric dipole moments of lipid monolayers: contributions of the water/lipid and the lipid/air interfaces. J. Colloid Interface Sci. 1988, 126, 408−420.

lack dipole vectors from the individual C−H bonds. Taking standard bond distances and partial atomic charges from OPLSAA or another similar atomistic force field, each C−H bond contributes with a dipole moment that is 15-fold smaller than that of a typical O−H bond from water molecules, but the ordering of the bonds with respect to the interfacial plane matters, and C−H bonds were arranged in an orderly manner, which increased their importance with respect to the electrostatic properties of the surfactant films.

■

ASSOCIATED CONTENT

* Supporting Information S

The Supporting Information is available free of charge on the ACS Publications website at DOI: 10.1021/acs.langmuir.5b02904. Two videos with the initial relaxation (first 25 ns) for the N = 40 and 80 systems (AVI) (AVI) Detailed description of the simulation protocols and the model systems, along with a discussion of the structural relaxation and the contributions comprising the total electrostatic potential (PDF)

■

AUTHOR INFORMATION

Corresponding Author

*E-mail: [email protected]. Phone: +55-16-3351-8090. Fax: +55-16-3351-8350. Author Contributions

The manuscript was written through the contributions of all authors. All authors have given approval to the final version of the manuscript. Notes

The authors declare no competing financial interest.

■

ACKNOWLEDGMENTS We are indebted to FAPESP (grant number 2012/15147-4) and CAPES for financial support. K.B. thanks CNPq and FAPESP (2011/14562-5) for scholarships. A.F.M. thanks MEC/PET for a fellowship.

■

REFERENCES

(1) Israelachvili, J. N. Intermolecular and Surface Forces, 2nd ed.; Academic Press Limited: London, 1998. (2) Hiemenz, P. C.; Rajagopalan, R. Principles of Colloid and Surface Chemistry, 3rd ed.; Taylor & Francis Group: New York, 1997. (3) Israelachvili, J. Self-assembly in two dimensions: surface micelles and domain formation in monolayers. Langmuir 1994, 10, 3774−3781. (4) Abranko-Ridge, N.; Darvas, M.; Horvai, G.; Jedlovszky, P. Immersion depth of surfactants at the free water surface: a computer simulation and ITIM analysis study. J. Phys. Chem. B 2013, 117, 8733− 8746. (5) Shelley, M. Y.; Sprik, M.; Shelley, J. C. Pattern formation in a selfassembled soap monolayer on the surface of water: a computer simulation study. Langmuir 2000, 16, 626−630. (6) Shi, L.; Tummala, N. R.; Striolo, A. C12E6 and SDS surfactants simulated at vacuum-water interface. Langmuir 2010, 26, 5462−5474. (7) Nguyen, C. V.; Phan, C. M.; Ang, H. M.; Nakahara, H.; Shibata, O.; Moroi, Y. Molecular dynamics investigation on adsorption layers of alcohols at the air/brine interface. Langmuir 2015, 31, 50−56. (8) Gang, H. Z.; Liu, J. F.; Mu, B. Z. Monolayer of aerosol-OT surfactants adsorbed at the air/water interface: an atomistic computer simulation study. J. Phys. Chem. B 2011, 115, 12770−12777. (9) Domínguez, H. Self-aggregation of the SDS surfactant at a solidliquid interface. J. Phys. Chem. B 2007, 111, 4054−4059. I


Article

Langmuir (32) Demchak, R. J.; Fort, T., Jr. Surface dipole moments of closepacked un-ionized monolayers at the air-water interface. J. Colloid Interface Sci. 1974, 46, 191−202. (33) Petrov, J. G.; Polymeropoulos, E. E.; Möhwald, H. Threecapacitor model for surface potential of insoluble monolayers. J. Phys. Chem. 1996, 100, 9860−9869. (34) Guldbrand, L.; Jönsson, B.; Wennerström, H.; Linse, P. Electrical double layer forces. a Monte Carlo study. J. Chem. Phys. 1984, 80, 2221−2228. (35) Dalmaschio, C. J.; Firmiano, E. G.; Pinheiro, A. N.; Sobrinho, D. G.; de Moura, A. F.; Leite, E. R. Nanocrystals self-assembled in superlattices directed by the solvent-organic capping interaction. Nanoscale 2013, 5, 5602−5610. (36) de Moura, A. F.; Bernardino, K.; Dalmaschio, C. J.; Leite, E. R.; Kotov, N. A. Thermodynamic insights into the self-assembly of capped nanoparticles using molecular dynamic simulations. Phys. Chem. Chem. Phys. 2015, 17, 3820−3831.

J


Surface Electrostatic Potential and Water Orientation in the presence

Recommend Documents