7754
J. Phys. Chem. B 2008, 112, 7754–7761
Orientational Dynamics of Water in the Nafion Polymer Electrolyte Membrane and Its Relationship to Proton Transport Matt K. Petersen,† Alison J. Hatt,‡ and Gregory A. Voth* Department of Chemistry and Center for Biophysical Modeling and Simulation, UniVersity of Utah, Salt Lake City, Utah 84112-0850 ReceiVed: January 9, 2008; ReVised Manuscript ReceiVed: March 11, 2008
Orientational anisotropies are calculated from molecular dynamics simulations of bulk water and the Na+ and H+ forms of hydrated Nafion and then compared with corresponding experimental values. The extended jump model of Laage and Hynes is applied to water reorientations for each system, and the anisotropies are explored as a product of hydrogen bond restricted “wobble-in-a-cone” reorientations and that due to the discrete jumps of hydrogen bond reorganization. Additionally, the timescales of hydrogen bond switching and proton transport are presented for bulk water and the H+ form of hydrated Nafion. The short time scale of proton hopping is found to be independent of Nafion water loading, suggesting the short time dynamics of proton hopping are relatively insensitive to the level of hydration. Furthermore, the long time decay for the forward rate of hydrogen bond switching is shown to be identical to the long time decay in the forward rate of proton hopping, for bulk water and all water loadings of Nafion investigated, suggesting a unified process. Introduction 1,2
experimental3
and interest in The significant theoretical DuPont’s Nafion has been driven by the successful application of Nafion (and related perfluorosulfonic acid membranes: Asahi Chemical’s Aciplex and Flemion and Solvay-Solexis’ Hyflon, formerly Dow Chemical’s Dow Membrane), as the electrolyte in polymer electrolyte membrane fuel cells (PEMFC). PEMFC, in particular, have weight and operational temperature characteristics suited to applications ranging from ubiquitous portable electronics to a potential replacement for the internal combustion engine in automobiles.4–6 PEMFC are capable of using renewable fuels such as methanol and ethanol, with few harmful emissions relative to the current technologies. However, for PEMFCs to become a competitive energy delivery technology, the efficiency-limiting properties (and, conversely, the desirable properties) need to be well characterized with a view toward development of novel membrane replacements. Nafion (Figure 1) is composed of a fluorocarbon backbone with periodically spaced sulfonic acid terminated perfluoroalkyl ether pendants. When hydrated, the polymer is characterized by the formation of phase-separated domains. The extreme polarity mismatch of the ionomer drives the organization of hydrophobic regions containing the polymer backbone and hydrophilic regions that contain the sulfonate and counterions solvated by constituent water. The high proton mobility exploited in PEMFC is facilitated by an extended network formed by the hydrophilic phase of the hydrated polymer. Within these hydrophilic regions the cation is free to diffuse and more interestingly, in the case of H+, diffuse through structural rearrangements in the extended hydrogen bond network of water consistent with the Grotthuss mechanism.7,8 We have previously investigated the mechanism of proton * To whom correspondence should be addressed. E-mail:
[email protected]. † Center for Integrated Nanotechnology, Sandia National Laboratories, Albuquerque, New Mexico 87185-1415. ‡ Materials Department, University of California, Santa Barbara, California 93106-5050.
transport in hydrated Nafion with respect to ion pair solvation structures found in these hydrophilic regions and their effect on proton transport dynamics.9,10 In this work we focus on the relationship between proton shuttling and reorientation of the hydrogen bond network. We present results obtained through molecular dynamics simulations for both the Na+ and acidic form of Nafion as well as bulk H2O for comparison. The related dynamics of proton transport and hydrogen bond rearrangements are investigated through an extension of the Ivanov jump mechanism11–14 and compared to recently reported and similarly analyzed experimental reorientational anisotropies.15,16 This article is organized as follows: the first section contains a short description of the methods used to generate the data for each individual simulation set. This methods section is followed by a discussion section covering the water reorientational dynamics in the simulated Na+-Nafion systems with a comparison to recent experimental results.15 The discussion then covers water reorientation in H+-Nafion and relates this reorganization to proton hopping rates. Concluding remarks follow the discussion section. Methods The results presented below were gathered from five distinct sets of molecular dynamics (MD) simulations: (1) The first set consisting of Nafion at hydration levels such that λ ) [SO3-]/[H2O] ) 1, 5, 3, 7.5 wherein Na+ was the counterion. These simulation afford direct comparison to the experimental treatment of Moilanen et al.15 (2) The second set, previously presented,9 consists of the acidic form of Nafion at hydration λ ) 15 and 7. The preponderance of the sulfonic acid protons were treated with a classical hydronium, or “single state” approximation, while a single proton of 1 of the 40 ion pairs was treated with the multistate empirical valence bond method (MS-EVB).17,18 (3) The third set, also previously presented elsewhere (with the exception of the lowest water loading),10 parallels the second set of simulations. However, although the water loading is the
10.1021/jp800221x CCC: $40.75 2008 American Chemical Society Published on Web 06/06/2008
Orientational Dynamics of Water
Figure 1. Chemical structure for the acid form of Nafion. The number of monomers in a given chain is on the order of several hundreds, while the length of the fluorinated backbone of each unit can vary from 5 to 15. Nafion of an equivalent weight of 1100 (used in this study) has on average m ) 7.
same, all of the excess proton ion pairs were treated with the self-consistent multi state empirical valence bond (SCI-MSEVB) method.19 (4) The fourth set consists of a single excess proton in bulk water treated with the MS-EVB method. (5) And the last, likewise, parallels the fourth set with a classical treatment of a single excess proton in bulk water. It is worth noting that all of the trajectories analyzed, with the exception of the Na+-Nafion and lowest water loading SCIMS-EVB H+-Nafion systems, were part of previous studies of proton transport.9,10,20 The advent of the Moilanen et al. publication15 of the extended jump model as applied to the orientational correlation function clearly warranted revisiting each system. MS-EVB and SCI-MS-EVB Simulations. Proton diffusion in water is distinct from that of other monovalent cations in that an excess proton may either diffuse through typical “vehicular” diffusion or by shuttling along the extended hydrogen bond network through successive covalent and hydrogen bond formation/cleavage consistent with the Grotthuss mechanism.7,8,18,20,21 Typical empirical potentials used in classical MD simulations do not model this fluctuating bonding topology, and hence omit both the Grotthuss shuttling process and the delocalization of the excess proton charge defect. One computationally tractable approach, capable of modeling this explicit proton physics, is the MS-EVB and the related multiproton SCI-MS-EVB methods. In short, the MS-EVB method (well described elsewhere, in refs 18 and 22–24 and references cited therein) assumes that the ground-state potential energy surface can be described through a linear combination of limiting bonding topologies. At each time-step in the dynamical trajectory, a restricted set of possible bond permutations is identified and the associated potential energy for each configuration is calculated. The Hamiltonian for the system is constructed such that the diagonal elements are the potential energy expressions for these limiting states. Additionally, off-diagonal elements are included, which allow for coupling between the limiting configurations and further provide for electrostatic interactions between the transferring moiety and the surrounding solvent. A straightforward extension of the MS-EVB method to more than just a few excess protons results in an intractably large number of basis states. An approximation capable of reducing the number of states is the SCI-MS-EVB method.10,19 The SCIMS-EVB method achieves a tractable set of bonding topologies by deconstructing the multiproton Hamiltonian into a set of single proton “EVB centers”. Each EVB center is comprised of a single excess proton and the associated solvating water molecules. In any given time-step, an initial solution to each single proton EVB matrix is calculated according to that from the previous step. From this initial guess, subsequent solutions are calculated within the effective field of all other EVB centers wherein the interaction parameters are scaled according to the EVB coefficients calculated from the previous iteration. The
J. Phys. Chem. B, Vol. 112, No. 26, 2008 7755 process is repeated until the calculated potential energy converges within a selected tolerance. Na+-Nafion. Four Nafion chains of 1100 EW, each containing 10 equally spaced monomers, were simulated with the potential first employed by Jang et al.25 The water portion of the potential was the flexible TIP3P model, nominally modified for use in the MS-EVB2 force-field,17 and the Na+ ions parameters were taken from the OPLS-AA force-field. The potentials were mixed with the standard Lorentz-Berthelot mixing rules, and the Ewald summation method was used for all electrostatic interactions. The polymer, ions, and water were randomly placed at approximately 50% of the experimental density. A constant temperature equilibration was carried out at 400 K for 500 ps, followed by 1 ns of constant temperature/pressure equilibration at 300 K and 1 atm. The equilibrated densities were all within 5% of experimental densities26 and were found to agree with previously published computer simulations.27 The ending configurations of these equilibration trajectories were then used to generate 500 ps constant volume trajectories over which data was collected, wherein the temperature was held at 300 K with a Nose´-Hoover thermostat.28 The same procedure was used for each water loading, that is λ ) 1, 5, 3, and 7.5. H+-Nafion (MS-EVB). The details of the H+-Nafion MSEVB simulations have been presented in detail elsewhere.9 In brief, the simulations were constructed as in the previous set with the exception of using a single state classical hydronium cation as the counterion in the construction and equilibration process. As mentioned earlier, a classical hydronium cation model does not include Grotthuss shuttling and proton charge delocalization. Like the Na+-Nafion simulations, the densities were in reasonable agreement with previous simulations25 and experiment.29 Following the initial low density/high temperature run, and subsequent 300 K/1 atm annealing, 10 distinct trajectories were created by randomly replacing one of the classical hydronium cations with one to be treated with the full MS-EVB2 potential. Each configuration was used to generate a 200 ps constant NVT trajectory over which data was collected, the temperature being held constant with a Nose´-Hoover thermostat. H+-Nafion (SCI-MS-EVB). This set of simulations, again detailed elsewhere,10 was constructed like the first two with a single state hydronium cation used as the counterion throughout the construction and pre-equilibration. The resulting densities were identical to those of the H+-Nafion MS-EVB simulations. Following the low density/high temperature run, and subsequent 300 K/1 atm annealing, 5 seed configurations were selected at evenly spaced intervals along a 300 K constant temperature/ volume trajectory. These configurations were then equilibrated with the full SCI-MS-EVB potential for 200 ps at 300 K and constant volume. Lastly, the ending configurations of the previous trajectories were used to start five, 500 ps microcanonical (constant NVE) trajectories over which data was collected for analysis. H+-H2O (MS-EVB). A single state classical hydronium cation was added to a previously equilibrated configuration of 256 water molecules. The system was further equilibrated at constant temperature/volume with a Nose´-Hoover thermostat for 500 ps. Twenty unique configurations were collected at equal spaced intervals from an additional 500 ps of a constant temperature/volume trajectory. From these 20 configurations, 500 ps microcanonical trajectories were generated over which data was collected. The excess proton in bulk water has been extensively studied with the MS-EVB2 potential.17,20,30,31 Proton
7756 J. Phys. Chem. B, Vol. 112, No. 26, 2008
Petersen et al.
Figure 2. Reorientational anisotropies calculated from the simulated Na+-Nafion systems for a series of water loadings, λ. The short time, fast decay is shown in the inset.
and water diffusion, pair distribution functions, density, etc. were all found to be in agreement with these previous studies. H3O+-H2O. This set of simulations was generated with the same procedure as that for the H+-H2O MS-EVB simulations with the exception that single-state classical hydronium approximation for the hydronium ion was employed. In this approximation the same parameters used in the MS-EVB2 and SCI-MS-EVB force field for the hydronium ion are retained; however, no alternate bonding topologies are considered in the potential, thus eliminating the possibility of excess proton charge delocalization or Grotthuss shuttling.
Figure 3. Comparison of the simulated (broken line) and experimental (solid line) Na+-Nafion anisotropy for a series of water loadings, λ. The experimental anisotropies were constructed from the population relaxation times reported by Moilanen et al.15
include this regime in our fitting and comparison with experiment to avoid possible confusion. Figure 3 depicts the orientational anisotropy for the simulated Na+-Nafion systems, as well as experimental curves recreated from the published biexponential parameters of Moilanen et al.15 It has been shown35 that water confined on the nanometer length scale exhibits biexponential anisotropic decay such that
r(t) ) a1 exp(-t/τ1) + a2 exp(-t/τ2)
(3)
al.35
Discussion Water Reorientation in Bulk H2O and Na+-Nafion. The orientational anisotropy
r(t) )
S|(t) - S⊥(t) ) 0.4C2(t) S|(t) + 2S⊥(t)
(1)
is measured from the normalized parallel and perpendicular transient absorption in time-resolved polarization selective pump-probe experiments.32,33 It is readily related to the autocorrelation function of the transition dipole through the second order Legendre polynomial such that
C2(t) ) 〈P2[µˆ (0) · µˆ (t)]〉
(2)
The H2O orientational anisotropy of bulk water, and likewise in hydrated Nafion, is characterized by a fast (