Molecular Dynamics Simulations of Water Sorption in a

Sep 20, 2013 - Kevin B. Daly , Athanassios Z. Panagiotopoulos , Pablo G. Debenedetti , and Jay B. Benziger. The Journal of Physical Chemistry B 2014 1...
2 downloads 0 Views 4MB Size
Subscriber access provided by UNIV OF PITTSBURGH

Article

Molecular Dynamics Simulations of Water Sorption in a Perfluorosulfonic acid Membrane Kevin Brendan Daly, Jay B Benziger, Pablo G. Debenedetti, and Athanassios Z. Panagiotopoulos J. Phys. Chem. B, Just Accepted Manuscript • DOI: 10.1021/jp405440r • Publication Date (Web): 20 Sep 2013 Downloaded from http://pubs.acs.org on September 22, 2013

Just Accepted “Just Accepted” manuscripts have been peer-reviewed and accepted for publication. They are posted online prior to technical editing, formatting for publication and author proofing. The American Chemical Society provides “Just Accepted” as a free service to the research community to expedite the dissemination of scientific material as soon as possible after acceptance. “Just Accepted” manuscripts appear in full in PDF format accompanied by an HTML abstract. “Just Accepted” manuscripts have been fully peer reviewed, but should not be considered the official version of record. They are accessible to all readers and citable by the Digital Object Identifier (DOI®). “Just Accepted” is an optional service offered to authors. Therefore, the “Just Accepted” Web site may not include all articles that will be published in the journal. After a manuscript is technically edited and formatted, it will be removed from the “Just Accepted” Web site and published as an ASAP article. Note that technical editing may introduce minor changes to the manuscript text and/or graphics which could affect content, and all legal disclaimers and ethical guidelines that apply to the journal pertain. ACS cannot be held responsible for errors or consequences arising from the use of information contained in these “Just Accepted” manuscripts.

The Journal of Physical Chemistry B is published by the American Chemical Society. 1155 Sixteenth Street N.W., Washington, DC 20036 Published by American Chemical Society. Copyright © American Chemical Society. However, no copyright claim is made to original U.S. Government works, or works produced by employees of any Commonwealth realm Crown government in the course of their duties.

Page 1 of 39

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

Molecular dynamics simulations of water sorption in a perfluorosulfonic acid membrane Kevin B. Daly, Jay B. Benziger, Pablo G. Debenedetti, and Athanassios Z. Panagiotopoulos∗ Department of Chemical and Biological Engineering, Princeton University, Princeton, New Jersey 08544 E-mail: [email protected]

∗ To

whom correspondence should be addressed

1 ACS Paragon Plus Environment

The Journal of Physical Chemistry

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Abstract Atomistic molecular dynamics simulations are reported over a wide range of water contents and temperatures to obtain a better understanding of the structural and transport aspects of water sorption in Nafion, a perfluorosulfonic acid membrane, under equilibrium conditions. For the short Nafion chains studied, good agreement is found between the water sorption isotherms from simulations and experiments at intermediate hydration (2 . λ . 7, where λ is the number of water molecules per sulfonate group), suggesting that, in that range, the isotherm is insensitive to effects of polymer chain relaxation. If polymer chain relaxation were important for water sorption at these conditions, then the water uptake of experimental membranes, which contain very long chains, might be far from equilibrium, making it difficult to obtain agreement with equilibrated, short-chain simulations. At λ . 7, strong water-sulfonate interactions, rather than chain relaxation, may control water sorption, despite the fact that chain relaxation time increases dramatically with decreasing hydration. Evidence for strong water-sulfonate interactions is found in the observation that sulfonate groups share water molecules in their first coordination shells at λ . 7. Strong water-sulfonate interactions are also observed to influence transport properties like water diffusivity, and are as important for understanding these transport properties as larger-scale phenomena like morphology and percolation transitions. Finally, at low humidity (λ ≈ 1 − 2), rod-like hydrophilic clusters are observed, as well as a mechanism of water diffusion that differs qualitatively from that of water at high hydration (λ & 7) and in the bulk, pure-component phase.

Keywords: Nafion, ionomers, phase equilibrium, diffusion

2 ACS Paragon Plus Environment

Page 2 of 39

Page 3 of 39

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

1

Introduction

Polymer electrolyte membranes are an important component of hydrogen fuel cells, and materials for these membranes are selected to maximize proton conductivity as well as mechanical and chemical stability. Currently, the most widely used material is Nafion, a (poly)perfluorosulfonic acid that was first developed in the late 1960’s by duPont. 1 Nafion consists of comb-like, amphiphilic polymers that are synthesized through a copolymerization of tetrafluoroethylene and perfluorosulfonic vinyl ether monomers, with chemical structure shown in Figure 1. The hydrophobic perfluorocarbon backbones and perfluoroether side chains impart mechanical and chemical stability, while the hydrophilic sulfonic acid groups aggregate into nanoscale clusters through which proton conduction takes place. When Nafion sorbs water, the size of the clusters increases, as does proton conductivity.

Figure 1: Chemical structure of Nafion monomers used in the present work. Blue ellipses denote hydrophobic portions, while the yellow ellipse denotes the hydrophilic portion. In general, the spacing between side-chains may vary among chains in experimental samples. Additionally, the number of perfluoroether moieties may vary in the side chain for other perfluorosulfonic acid membranes, such as Hyflon. 2 Image adapted from Mauritz and Moore 3 .

Given the importance of hydrophilic clusters for proton conductivity, understanding the microscopic structure of Nafion may be a viable route to better membrane materials. However, this understanding so far remains elusive. 3 The clusters have a characteristic length scale that is known to vary with water content between 2.7 and 4.8 nm based on SAXS and SANS data, 4 but their precise morphology is poorly understood. Within these clusters, the coordination of water and sulfonic acid groups is unknown, as well as its relationship to proton transport. The molecular weight of the polymer chains in commercially available Nafion films such as Nafion 117 has proven difficult to measure with standard techniques; the few estimates that exist put it in the range of 10 5 to

3 ACS Paragon Plus Environment

The Journal of Physical Chemistry

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

10 6 Da (> 100 SO3 H groups per chain). 3 Finally, the distribution of side-chains along the backbone is unknown, 5 and is typically assumed to be uniform. 3 These gaps in understanding of Nafion persist despite four decades of experiments, leading researchers more recently to consider modeling approaches. As a result, numerous Nafion-related simulations and analytical models, with resolutions ranging from atomic to mesoscale, have been published over the past decade. 5–20 On the scale probed by to classical, atomistic simulations (O(1 − 10 nm)), pioneering work was done in the early 2000’s by Vishnyakov and Neimark. 12–15 They developed one of the first atomistic force-fields for the side-chain, 12 elements of which were borrowed in later force-fields, such as the force-field of Cui et al. 10 that forms the basis of the present work. Using this sidechain force-field, Vishnyakov and Neimark performed simulations of an isolated side-chain in water, and observed an anisotropic spatial distribution of water around the side-chain, with higher concentrations near the sulfonate group. 12 They also simulated a bulk membrane, and observed microphase segregation, as well as temporary water bridges between hydrophilic clusters. 14 Later in 2004, Jang et al. reported one of the first systematic studies of chain topology, in particular the distribution of side chains along the backbone. 5 For the short chains studied, they found that the extreme case of evenly spaced side-chains gives a modest decrease in water diffusivity (25%) over the opposite extreme of tightly clustered side-chains. More recently, Knox et al. performed very large scale (O(10 nm)) simulations of high molecular weight Nafion with a variety of imposed starting morphologies that have been suggested in the experimental literature. 6 They found that the “ionomer peak” present in SAXS and SANS data at Bragg distances of 2.7-4.8 nm is rather insensitive to morphology. One problem that has not been thoroughly explored with simulations is the non-equilibrated nature of Nafion. Like many polymers, Nafion has relaxation time scales that frequently exceed those of experiments, depending on conditions such as temperature and hydration level. This problem has been particularly apparent in measurements of the water sorption isotherms, in which the reported water content at saturation has been found to vary by nearly a factor of two. 21 Nafion also appears to show behavior consistent with Schroeder’s Paradox, in which water uptake differs

4 ACS Paragon Plus Environment

Page 4 of 39

Page 5 of 39

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

between Nafion equilibrated with 100% relative humidity (RH) vapor and with liquid water, a violation of thermodynamic equilibrium considerations. 21 Onishi and coworkers showed that this apparent paradox could be traced back to inconsistent thermal histories of the experimental samples, confirming that the samples were not truly at equilibrium. 21 Since water uptake is strongly correlated with proton conductivity, understanding the non-equilibrated nature of Nafion is of practical interest for fuel cell operation. An open question is the extent to which this non-equilibrated behavior manifests itself at humidities lower than 100%. Also of interest is the question of whether the rate of desorption could be lowered so that the membrane retains more water at low humidity, thereby improving performance. Operating at low humidity is of practical interest, because it enables high temperature operation without the need to pressurize the feed to the fuel cell. High temperatures in turn prevent poisoning of the anode catalyst. 22 A starting point for answering these questions is finding stable equilibrium states of Nafion/water systems with respect to water sorption, a task that is beyond the current capabilities of experiments because the chains used in actual experiments are too long. Synthesizing shorter chains with low polydispersity is difficult because the standard synthesis involves free radical polymerization. 23 In addition, the highly electronegative fluorine atoms preclude anionic polymerization, which would yield more monodisperse chains. The experimental literature provides information on additional competing effects that determine how fast Nafion approaches equilibrium. If chain relaxation were the only effect, then one would expect the rate of sorption of water to decrease with decreasing water content. After all, mechanical measurements show that the membrane becomes dramatically stiffer as water is removed, with the exception of anomalous behavior at high temperatures and ultra-low humidity (λ = Nw /NSO3 H < 1, where Nw is the number of waters and NSO3 H is the number of sulfonic acid groups). 24 However, experiments show that the rate of sorption actually increases with decreasing water content. 25 For example, when a dry membrane with a specific thermal history was exposed to a 99% RH environment, it took approximately 1 min to go from a water content of λ = 0 to λ = 4, 1 hour to go from λ = 4 to λ = 8, and 1,000 hours to go from λ = 8 to λ = 12. 25 These dra-

5 ACS Paragon Plus Environment

The Journal of Physical Chemistry

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 6 of 39

matically different rates suggest that in the limit of low water content, the activity gradient between the membrane and its surroundings overpowers the resistance of the chains to swelling, leading to a high overall flux of water. Experiments have also shown that the rate of sorption can be limited by interfacial transport resistance, particularly at high water content. 26 Finally, the self-diffusivity of water within the bulk membrane decreases by at least two orders of magnitude with decreasing water content. 27 All of these experimental complications make it difficult to determine how far an experimental membrane is from its true equilibrium state with respect to water sorption. The goal of the present computational work is to better understand equilibrium states, with respect to water sorption, of Nafion-water systems and to find the point at which non-equilibrated effects become significant. Unlike experiments, simulations allow chain length and structure to be controlled precisely, and the chains can be made short enough to reach equilibrium even on simulation time scales. We can also perform simulations with periodic boundaries to eliminate interfacial resistance to mass transfer. Because we are restricted to simulating short chains, we necessarily focus on properties that are insensitive to chain length. Our approach is to calculate these properties for our well-equilibrated simulations and compare them to experimental measurements performed on possibly non-equilibrated Nafion samples. By finding the point at which the two data sets deviate significantly, we can identify the approximate hydration threshold above which equilibration effects are important. To differentiate between error in the model and genuine non-equilibrated effects, we attempt to reproduce a variety of experimental properties. If we can demonstrate good agreement between simulations and experiments at low to moderate water content, where equilibration effects are not very significant in experiments, 25 then we can reasonably expect the model to perform well under conditions where non-equilibrium effects do matter. Although non-equilibrium effects appear in experimental measurements of both the proton conductivity and water sorption isotherm, the latter property is more accessible to classical simulations, and will therefore be the focus of the current work. We demonstrate that the sorption isotherm is indeed not sensitive to chain length by computing the water activity for a fixed water

6 ACS Paragon Plus Environment

Page 7 of 39

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

content and a range of polymer chain lengths. Details of this calculation are given in the Supporting information. To the best of our knowledge, only two prior simulation papers exist in which at least part of the sorption isotherm of Nafion was calculated. 28 The focus of the first work was predicting water content at saturation. 28 To that end, the authors calculated the water activity over a range of imposed hydration levels, and were able to find water activities both less than and greater than one. Activities greater than one imply metastability, and therefore a phase transition. An example of such a transition would be a polymer-solvent system with χ  1/2, where χ is the Flory-Huggins interaction parameter. 29 However, there is currently no experimental evidence for such a transition in Nafion membranes. Alternatively, the simulations may not have reached equilibrium, since the equilibration time following annealing was only 5 ns. The second simulation paper to calculate the sorption isotherm modeled Nafion as a network of rigid slit pores of two different sizes, and imposed a range of water activities using grand canonical Monte Carlo simulations. 30 Although this model gives a simple, conceptual understanding of water sorption and transport in Nafion, it significantly underpredicts water activity at low to moderate hydration levels. We believe this may be a result of the rigidity of the pores in the model, which contain large voids at low hydration. This model also depends on a careful choice of pore size distribution and sulfonate grafting density. Therefore, the predictions of this model warrant a follow-up simulation study that includes the flexible polymer backbone and does not impose a pore morphology. In this paper, classical molecular dynamics was used to simulate the equilibrium states of Nafion samples of low molecular weight across a wide range of water contents and temperatures. The sorption isotherm was calculated, as well as transport and structural properties that are experimentally measurable. The microscopic detail of the simulations was used to gain insights on the mechanism of water transport, as well as the morphology of the hydrophilic clusters. The paper is organized as follows: the setup, execution, and analysis of the simulations are described in Section 2; the results are discussed in Section 3, and the main conclusions are summarized in Section 4.

7 ACS Paragon Plus Environment

The Journal of Physical Chemistry

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

2

Methods

2.1

Force field

Nafion chains were modeled as a combination of united-atom sites for CFx moieties and atomistic sites for the remaining atoms. The sulfonate groups were assumed to be fully deprotonated so that each side chain has a net charge of -1 that is balanced by a free hydronium ion in the system. For sites in the same chain separated by three covalent bonds or less, interactions were represented through a combination of flexible bond, angle, and dihedral potentials, the parameters and functional forms of which were obtained from Cui et al. 10 . All other intra- and intermolecular interactions were represented with a combination of Lennard Jones and Coulombic terms, the latter acting between partial charges on each site. The parameters characterizing Coulombic and dispersive ineractions were also taken from Cui et al. 10 . In addition, a cut-off distance of 1 nm was used for Lennard-Jones interactions, 10 and standard long-range corrections were added to both the pressure and potential energy. 31 Electrostatic interactions were computed using the Smooth Particle Mesh Ewald (SPME) method 32 with fourth-order interpolation, a short-ranged cut-off of 1.2 nm, a κ of 2.6 nm−1 , and a grid-spacing of no more than 0.1 nm. Water was represented with the rigid, four-site TIP4P/2005 model, 33 which reproduces several thermodynamic, structural, and transport properties of liquid and crystalline water with good accuracy. 34 For hydronium, a four-site model was also used, albeit with flexible bonds and angles. The parameters for this model come from Cui et al. 10 , the same source as for the parameters of the Nafion chain. Since the degree of microphase separation depends on interactions between water molecules and CFx groups, we felt that it was necessary to test the force-field parameters for a simple system that was well characterized by experiments. Therefore, the solubility of water was calculated in perfluorohexane using a combination of Molecular Dynamics and Widom test particle insertions, 35 and compared to the experimental data of Freire et al. 36 at several temperatures. In addition to the

8 ACS Paragon Plus Environment

Page 8 of 39

Page 9 of 39

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

base model of Cui et al. 10 , a variety of water and polymer force-fields, both atomistic and unitedatom, were tested. Neither the water nor the perfluorohexane models have been parameterized to reproduce properties of perfluorocarbon-water mixtures. For interactions between unlike particles, both commonly used Lorentz-Berthelot and Kong 37 rules were tested. All model combinations were found to systematically underpredict the mole fraction of water in liquid perfluorohexane by factors between 2 and 8.5. We deemed this level of agreement unsatisfactory, and optimized the unlikeinteraction parameters of the TIP4P/2005 force field paired with the united-atom perfluoroalkane force field used in our simulations of Nafion. This optimization consisted of adjusting the LennardJones σ and ε for interactions between water oxygens and CFx groups. After the optimization, the solubility of water in perfluorohexane could be predicted to within 10% of experiments over a range of temperatures. The details of these force-field parameter determination simulations have been included in the Supporting Information.

2.2

Molecular dynamics simulations

The majority of our simulations consisted of 64 Nafion chains with a molecular weight of 3450 Da. Each chain had three evenly-spaced side chains terminating with sulfonic acid moieties. Although the spacing was uniform, the placement of each set of three side-chains along the backbone was asymmetric: the first side-chain was located at the second carbon from the nearest terminus and the last was located at the fifteenth carbon from the nearest terminus. The case of symmetric placement was also considered, and the pair correlation functions between hydrophilic species in equilibrated membranes were found to be similar to those of the asymmetric case. The additional results for symmetric placement of chains are reported in the Supporting Information. The simulation boxes had periodic boundaries and lengths ranging from 5 to 9 nm after equilibration, depending on the hydration level. Systems for which calculations are reported in Section 3 can be assumed to meet all of these specifications unless otherwise noted. For each simulation, an artificial initial configuration was constructed in which chains were 9 ACS Paragon Plus Environment

The Journal of Physical Chemistry

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 10 of 39

extended and spaced evenly apart in an 18 nm box. A slab of water molecules from an equilibrated configuration of liquid water was also inserted into the box. The system was then allowed to spontaneously find its stable equilibrium state by evolving via molecular dynamics in the NPT ensemble using Gromacs 4.5.3. 38 The temperature was fixed with the Nosé-Hoover thermostat 39 and the pressure with the Parinello-Rahman barostat, 40 both with τ = 10 ps. The box length shrunk from 18 nm to its final value of 5-9 nm during the course of the equilibration phase. The densities at various hydration levels are included in the Supporting Information. The time step was set at 0.5 fs and the equilibration phase was continued until the total energy of the system appeared to plateau. Equilibration took anywhere from 10 to 120 ns to reach, with longer times required for lower temperatures and hydration levels. Note that some previous studies using united atom models were equilibrated for only 1-2 ns. 10,41 After equilibration was complete, a production run was performed for 30-200 ns, again with longer times required for lower temperatures and hydration levels.

2.3

Thermodynamic and structural properties

The chemical potential of water was calculated in the production MD period of the simulations using a GPU implemention of Bennett’s method, 42 as described in Ref. 43. Configurations were taken at 1 ps intervals, and in each configuration 2,560,000 insertions and the maximum number of deletions (equal to the number of water molecules in the system) for the given hydration level were performed. A number of insertions and deletions in this range have been shown to optimally balance sampling and computational efficiency for aqueous systems of this size. 43 The structure factor, S(k), is a measure of spatial correlations between particles that is proportional to the intensity of a scattered wave in small-angle scattering experiments. The calculation of S(k) can be performed in two ways that vary in efficiency and accuracy. 44 The first way, which is refered to as the “direct method”, involves explicitly the particle coordinates ri , : * S(k) =

1 N

"

N

#2 +

∑ sin(k · ri)

i=1

* +

1 N

"

N

#2 +

∑ cos(k · ri)

i=1

10 ACS Paragon Plus Environment

(1)

Page 11 of 39

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

where k = |k| = |2π/L · (nx , ny , nz )| and L is the box length. nx , ny , and nz are integers ranging from 1 to ∞. The orientation of k is defined with respect to an arbitrary reference coordinate system. The advantage of this method is that it can access the low-k region down to k = 2π/L. The disadvantage of this method is that it becomes computationally intractable at high k, since the number of corresponding k’s is proprotional to k2 . Therefore, in the high-k region, the Fourier transform method becomes advantageous. This method involves the pair-correlation function:

S(k) = 1 + 4πρ

Z ∞ sin(kr)

r

0

k

[g(r) − 1]dr

(2)

where ρ is number density and g(r) is the pair correlation function. The disadvantage of this method is that it is limited to k ≥ 4π/L, and is very sensitive to error in g(r) at low k. In Section 3 both methods are used to estimate the structure factor and demonstrate that they are consistent where their k-ranges overlap. + The surface area of clusters of hydrophilic species (SO− 3 , H2 O, and H3 O ) is estimated by

computing the Connolly surface obtained by rolling a probe molecule over their van der Waals spheres. 45 All other species were deleted from the simulation trajectory. The radius of the probe molecule was defined to be that of a backbone CF2 group. The Connolly surface is computed by first discretizing the van der Waals spheres of hydrophilic particles and then performing discrete image processing operations, as described in detail in the Supporting Information. Using this procedure, not only could the surface area of the hydrophilic phase be estimated, but also the volume contained within it, allowing us to study the scaling of surface area to volume ratio, as discussed in Section 3.

11 ACS Paragon Plus Environment

The Journal of Physical Chemistry

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

2.4

Page 12 of 39

Transport properties

The diffusivity of water D was computed by measuring the mean-squared displacement hr2 (t)i and applying the Einstein relationship: 1 dhr2 (t)i t→∞ 6 dt

D = lim

(3)

To restrict our calculations to the diffusive regime, linear regression of lnhr2 (t)i vs. lnt was performed over successive logarithmic decades to estimate the scaling constant n, where hr2 (t)i ∼ t n . D was then estimated from the slope of hr2 (t)i vs. t only if |1 − n| < 0.05. Although diffusivity is the most relevant measure of water or hydronium translational mobility, this quantity is difficult to calculate at low hydration because molecules at these conditions reach the diffusive regime only at very long times. Therefore, a translational relaxation time was computed. The translational relaxation time is defined to be the average time needed for a molecule to move a root mean-squared distance equal to its size, considered to be 3 Å for both water and hydronium molecules. To measure rotational mobility, the autocorrelation function of a specific unit vector u was computed in the reference frame of the molecule:

C(t) = hu(t) · u(0)i

(4)

For water, u was chosen to be parallel to the dipole. For side chains, u was chosen to be parallel to the difference in the position vectors of the S atom and the CF group that anchors the chain to the backbone. The autocorrelation function C(t) is a meaningful measure of rotational mobility since it is also the first Legendre polynomial of cos θ , where θ is the angle between the unit vector at time t and at time zero. This Legendre polynomial arises as the first term of a series representation of the probability of angular displacements P(θ ,t) that also contains the rotational diffusion coefficient. 46

12 ACS Paragon Plus Environment

Page 13 of 39

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

Moreover, the characteristic time of C(t) can be measured experimentally by infrared absorption, Raman scattering, and NMR. In simulations, this characteristic time τ was estimated by fitting C(t) to a stretched exponential function: 47 i h C(t) = exp −(t/τ)β

(5)

where β is a parameter obtained by fitting the data and 0 ≤ β ≤ 1. The fitting procedure consisted of taking the logarithm of Eq. 5 twice and performing a linear regression. The relaxation time of a monomer was measured by computing the autocorrelation function of the normal modes of the backbone X p : 48

C(t) =

hX p (t) · X p (0)i hX p (0) · X p (0)i

(6)

The set of {X p } was obtained from discrete Fourier transforms of monomer coordinates from all the chains. 48 A value of the mode index p was selected to have a corresponding wavelength N/p equal to the length of a monomer, where N is the degree of polymerization. We defined a monomer as 16 CFx beads. We then estimated a characteristic time τ by fitting C(t) to a stretched exponential function using the procedure discussed in Section 2.4.

3 3.1

Results and Discussion Water sorption isotherm

Figure 2 shows isotherms at 353 K from simulations and experiments. 27 Between λ ≈ 2 and λ ≈ 7, the isotherms agree quantitatively, deviating at most by ca. 10%. By comparison, the isotherm of Hwang et al. 30 deviated by up to 600% in the same range. The agreement of our results with experiments is notable given that the simulated systems are well-equilibrated with respect to chain relaxation, while the experimental membranes are decidedly non-equilibrated with respect to chain

13 ACS Paragon Plus Environment

The Journal of Physical Chemistry

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

relaxation. Slowly relaxing chains resist swelling of the hydrophilic clusters, lowering the rate of sorption, and by extension, the apparent, non-equilibrated water uptake. Experimentally, this effect is very pronounced at saturation (aw = 1). 21 However, the experimental rate of sorption also accelerates by orders of magnitude with decreasing hydration, 25 indicating that chain resistance to swelling may become overpowered by some other effect. At a low enough hydration level, the water uptake should be able to reach equilibrium values on typical experimental time scales. This hydration level appears to be in the neighborhood of λ ≈ 7, based on the agreement seen between isotherms from experiments and equilibrium simulations.

Figure 2: Water activity, aw = fw / fw,pure (T, p), as a function of water content, λ = Nw /NSO3 .  , experimental isotherm at T = 353 K (from 27); , simulation isotherm at T = 298 K; N, simulation isotherm at T = 353 K; , simulation isotherm at T = 353 K with unoptimized CFx /H2 O LJ mixing parameters. −, isotherm fitted to simulation data at T = 353 K with DGWC theory (see text for details). 49 −− −−, isotherm fitted to experimental data at T = 353 K with DGWC theory. Dashed line marks the water content at which SO− 3 groups no longer share water molecules in their first coordination shells (see Figure 3). Insert shows the data in greater detail for 0 ≤ λ ≤ 10.

The λ ≈ 7 transition region may have a microscopic basis, as suggested by Figure 3. In the inset, the sulfur-oxygen (H2 O and H3 O+ ; data for S/O(H3 O+ ) can be found in the Supporting Information) pair-correlation functions exhibit a clear peak over a range of hydration levels at a distance r < 4.65 Å from the sulfur atom, indicating the presence of a well-defined first coordination shell. The average number of water molecules, or occupancy (Nw,1 ), in this shell can be found by integrating the area under the corresponding peak. Figure 3 shows how at λ . 7, Nw,1 is actually 14 ACS Paragon Plus Environment

Page 14 of 39

Page 15 of 39

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

higher than the total number of waters per sulfonate in the system, λ . In other words, the computed curves in this region lie above the line Nw,1 = λ in Figure 3. Therefore, in order the achieve the observed occupancy, the sulfonate groups must share water molecules in their first coordination shell. These shared water molecules interact more strongly with the membrane than water molecules near only one sulfonate group, perhaps strongly enough to overpower the resistance of chains to cluster swelling. The significance of λ ≈ 7 has been identified by both experiments and other simulations. Pivovar and Pivovar observed through quasielastic neutron scattering that the fraction of non-diffusing hydrogen atoms decreases with increasing λ , and levels off at λ ≈ 7. 50 Devanathan et al. performed atomistic molecular dynamics simulations and reported that the average number of sulfonates coordinated to a hydronium ion decreases sharply from roughly 2.5 at λ ≈ 1 to 1 at λ ≈ 7, and then slowly approaches 0.5 (λ ≈ 20). 16 A coordination number greater than one implies that sulfonate groups are sharing hydroniums, much as they share waters in Figure 3. Notably, the simulations results of Devanathan et al. were obtained using a different force field, which includes a fully atomistic representation of the backbone and a different water model. This consistency between their force field and the force field used in the present work implies that the precise representation of the Nafion backbone does not strongly influence the hydration of the sulfonate groups. Consequently, model simplifications like the united-atom backbone used in this work are not only justified, but also advantageous since they enable longer runs and larger systems.

15 ACS Paragon Plus Environment

The Journal of Physical Chemistry

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Figure 3: Number of water molecules in the first solvation shell of SO− 3 , Nw,1 , where the shell is defined by the first peak in g(r) ending at r = 4.65 Å. , simulations at T = 298 K; , simulations at T = 353 K. Solid line represents Nw,1 = λ . Dashed line marks approximate intersection between the data and the solid line; this is the water content at which SO− 3 groups no longer share water molecules in their first coordination shells. Inset: pair correlation function for pairs of sulfur and oxygen (H2 O and H3 O+ ) at T = 353 K and λ = 1.3, 2.7, 8.8, and 16, where the height of the first peak decreases with λ . Solid line: definition of the first hydration shell.

Beyond λ ≈ 7, the occupancy increases slowly, indicating that the first coordination shell is near full capacity. At this point, water added to the system would begin to resemble bulk water (aw = 1), no longer strongly coordinated with sulfonate groups. This behavior explains why the simulation isotherm slowly approaches aw = 1 after λ ≈ 7, as shown in Figure 2. On the other hand, the experimental isotherm reaches aw = 1 at only λ ≈ 15, far away from the dilute, waterrich limit. At this hydration level, water sorption in experiments is known to be extremely slow to the extent that membranes continue to sorb water after one month of equilibration. 24 Therefore, the equilibrium water uptake may be difficult to access experimentally in this high hydration region (λ & 7). Experimental rates of water sorption become very fast at λ . 2, suggesting that chain relaxation is even less important for water sorption than it is at 2 . λ . 7. We would expect even better agreement with simulations, but in fact, simulation activities were found to be systematically lower than those of experiments in this region, as shown in Figure 2. We attribute this deviation to a known limitation of our force field, namely the assumption of complete dissociation of SO–3 and 16 ACS Paragon Plus Environment

Page 16 of 39

Page 17 of 39

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

H3 O+ . This assumption is not accurate for λ < 3 according to ab initio simulations of triflic acid and water. 51 Deprotonated sulfonate groups have a net charge and interact more strongly with polar water molecules, leading to artificially high heats of sorption, as we will show. The success of our model in the region of 2 . λ . 7 can be partly attributed to the optimized CFx /H2 O LJ mixing parameters discussed in Section 2.1. Figure 2 shows an isotherm computed from simulations without the optimized parameters. Not suprisingly, the activities are substantially higher than the experimental activities; in Section 2.1 we explained how standard mixing rules greatly underpredict the solubility of water in pure perfluoroalkanes. Despite the lack of quantitative agreement with experiments, the unoptimized isotherm still levels off around λ ≈ 7, suggesting that the isotherm is controlled primarily by water/side-chain interactions. Figure 2 shows an isotherm at 298 K from simulations with optimized parameters, in addition to the isotherm at 353 K already discussed. Both isotherms have a sigmoidal shape characteristic of Type II isotherms, 52 in qualitative agreement with the experimental isotherm at T = 353 K. Type II isotherms can be fit well by BET theory, in which the absorbate can form multiple layers on the absorbent, with the first layer having a heat of sorption higher than the heats of subsequent layers, which in turn are set equal to the heat of vaporization of the bulk absorbate qvap . The heat of sorption, qst , of simulated Nafion is also larger than the heat of vaporization of water, especially at low water content, since aw increases with T , and 

∂ (ln aw ) R ∂ (1/T )

 = −(qst − qvap )

(7)

Nw ,Nn ,p

Here the heat of sorption is defined to be the difference between water partial molar enthalpies in Nafion and in the vapor phase. aw is the activity of water, or fw / fw,pure , where fw is the fugacity of water at the given temperature and the pure water reference state is at p = 1 bar. Note that at this low pressure, to an excellent approximation, the fugacity is equal to the vapor pressure, fw,pure ≈ pw,sat , so our activities can be compared to experimental activities where the reference state is saturation. Nw is the number of moles of water, Nn is the number of moles of Nafion, and p

17 ACS Paragon Plus Environment

The Journal of Physical Chemistry

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

is the pressure. qvap is the heat of vaporization of pure water. qst was calculated at various values of λ = Nw /NSO3 by plotting ln 1/ f versus 1/T for data at T = 298, 314, 332, and 353 K, as shown in the Supporting material. A line was then fitted to data at each composition, and from the slope of the line, qst was obtained. qst is shown in Figure 4, and is a monotically decreasing function of λ that approaches qvap as λ → ∞. If we interpret each sulfonate group to be a sorption site in the context of BET theory, 53 then it is clear that this theory is too simple to fit both the isotherm data and the heat of sorption data. In the BET theory, only the first monolayer has a higher heat of sorption than the subsequent layers, whereas in our data, qst is still appreciably greater than qvap past the first monolayer (λ = 1).

Figure 4: Isosteric heat of sorption, qst , normalized by the heat of vaporization, qvap , of bulk water. , heats obtained from ln 1/ f vs. 1/T plots of simulation isotherms at T = 298, 314, 332, and 353 K. , heats predicted by the fit of the simulation isotherm at T = 353 K to DGWC theory. 49 4 , heats predicted by the fit of the experimental isotherm at T = 353 K to DGWC theory. Solid lines are meant to guide the eye. Dashed line marks the water content at which SO− 3 groups no longer share water molecules in their first coordination shells (see Figure 3). The BET theory was recently extended by Dutcher and coworkers, who formulated the DGWC theory, 49 that can accommodate an arbitrary number of layers (n) with qst > qvap . Introducing more layers adds more adjustable parameters, but the advantage of these BET-like theories is that the adjustable parameters can be readily compared to the heats of sorption obtained empirically from plots of ln 1/ f versus 1/T . In principle, other theories based on reaction equilibria could be used to fit the isotherm, 54 however the adjustable parameters (equilibrium constants) do not 18 ACS Paragon Plus Environment

Page 18 of 39

Page 19 of 39

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

lend themselves to easy comparison with heats of sorption. In the DGWC theory, an additional parameter, KA , controls the extent to which the heat of sorption in the outer layers differs from the bulk heat of vaporization. Setting KA = 1 and n = 2 recovers BET theory. KA was kept equal to 1 when modeling Nafion with DGWC theory, since Figure 4 clearly shows that qst has practically converged to qvap at high hydrations. n was increased to a value of 5 until a satisfactory fit was obtained. The DGWC isotherm was then fitted to the simulation isotherm at 353 K and λ ≤ 7, since this is the region where water/membrane interactions are strongest, as demonstrated in Figure 3. Details of the fitting procedure are provided in the Supporting Information. Given the close agreement between the model predictions for qst and the corresponding simulation values, it is clear that DGWC theory is an accurate physical description of the system for λ . 7. The same fitting procedure was performed for the experimental data, except n was set equal to 4 to obtain a good fit. The predicted values of qst were found to be qualitatively similar to the simulation values, but are significantly lower for λ < 2. This deviation confirms our hypothesis about the deviation in Figure 2; namely, the assumption of complete deprotonation leads to artificially strong water/sulfonate interactions at very low hydration.

3.2

Cluster Morphology

The only direct experimental measurement of cluster morphology in Nafion is provided by smallangle scattering spectra, which are related to pair correlation functions between particles in the system. Between small-angle neutron scattering (SANS) and small-angle x-ray scattering (SAXS) spectra, the SANS spectrum appears to be more directly comparable to simulations, since it can be transformed into a structure factor without accounting for the atomic scattering function, which for SAXS depends on the wave vector, k = 2π/d, where d is the Bragg distance. 55 Moreover, the neutron scattering cross-section of hydrogen is roughly an order of magnitude higher than that of any other element in Nafion, so the SANS spectra can be readily interpreted as the structure factor of the hydrogen-hydrogen pair correlation function. The experimental SANS spectrum contains a broad peak at k∗ ≈ 1.3 − 2.3 nm−1 (d = 4.8 − 2.7 19 ACS Paragon Plus Environment

The Journal of Physical Chemistry

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

nm), with k∗ inversely proportional to hydration level. 4 The structure factor from simulations also contains a broad peak within this range, as shown in Figure 5. The exact location of this peak at a water content of λ = 7.1 is d = 2.9 nm (k = 2.15 nm–1 ), which is slightly lower than the experimental value at the same water content, 3.8 nm (k = 1.7 nm–1 ). Note that d, the Bragg distance, is often interpreted to be the distance between two equivalent scattering centers, e.g. spherical water clusters. 3 The small deviation between experiment and simulations does not appear to be a consequence of finite system size or chain length, as demonstrated in Figure 5. However, it could be attributed to the lack of polymer crystallinity in the simulations. 3 Nevertheless, the agreement between simulations and experiments is notable considering that the predicted morphology is spontaneous, and not imposed at the beginning of the simulation.

Figure 5: Structure factor of hydrogen atoms in simulations of Nafion at T = 353 K and λ = 7.1. In all cases, error bars are smaller than the size of the symbols. , S(k) for system with L = 11.65 nm and number of sulfonate monomers per chain, N, equal to 3 calculated by the direct method (Eq. 1); accompanying solid line is meant to guide the eye. , S(k) for the same system calculated by the Fourier transform method (Eq. 2); N, S(k) for system with L = 5.83 nm and N = 3 calculated by the direct method. , S(k) for system with L = 5.83 nm and N = 24 calculated by the direct method. Dash-dotted lines: positions of the “ionomer peak” 3 of a dry membrane (larger k) and an intact membrane with aw = 1 (smaller k) from SANS experiments. 4 Dotted line: position of the ionomer peak of a moderately hydrated membrane (λ = 7.1) from SANS experiments. 4

The simulations predict another feature of the experimental SANS spectra, the Porod region starting at k > 1.5 nm−1 . In this region, the intensity is proportional to k−4 at length scales less than the spacing between mesoscopic domains but greater than the thickness of the domain interface. 56 20 ACS Paragon Plus Environment

Page 20 of 39

Page 21 of 39

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

Such behavior is evident in the S(k) from simulations when plotted as S(k)k4 vs. k, as done in Figure 6. The Porod region appears as a horizontal segment of the rescaled S(k). The existence of this region in the simulations suggests that there is a sharp interface between hydrophobic and hydrophilic regions.

Figure 6: Rescaled structure factor of hydrogen atoms in simulations of Nafion at T = 353 K and λ = 7.1. , S(k) for system with L = 11.65 nm calculated by the direct method. The horizontal line marks the Porod region where S(k) ∼ k−4 .

Since the structure factor reveals little about the shape of the clusters, a snapshot of Nafion at λ = 1.3 is shown in Figure 7. This snapshot shows narrow, rod-like clusters that span the length of the box, suggesting that the clusters percolate even at nearly dry conditions. The shape of these clusters is consistent with analysis of experimental proton conductivities based on percolation theory. 54 This analysis suggests that in order for Nafion to be conductive at low hydration, the clusters must likely have a high aspect ratio consistent with rods or lamellae, but not spheres. Upon closer inspection, the clusters observed in Figure 7 appear to be composed of interdigitated sulfonate and hydronium molecules, as shown in Figure 8. This arrangement of sulfonate and hydronium moieties is consistent with X-ray diffraction measurements of triflic acid monohydrate. 57 Triflic acid has a chemical structure very similar to that of the side chains in Nafion. Although the consistency with experimental data is encouraging, these observed configurations at low hydration must be further verified with a more detailed model that explicitly incorporates protonation/deprotonation

21 ACS Paragon Plus Environment

The Journal of Physical Chemistry

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

and more accurately predicts the sorption isotherm at low hydration, as discussed in Section 3.1.

Figure 7: Snapshot of Nafion at λ = 1.3 and T = 353.15 K. Only hydrophilic species (SO− 3, + H2 O, H3 O ) are shown. A segment of one cluster is wrapped in a Connolly surface (purple) to emphasize its cylindrical shape.

Figure 8: Top: snapshot of Nafion at λ = 1.3 and T = 353.15 K. Only hydrophilic species (SO− 3, − + + H2 O, H3 O ) are shown. A segment of a cluster consisting of H3 O and SO3 species is colored to emphasize its molecular structure. Bottom: idealized cartoon of the observed structure with H3 O+ and SO− 3. For fuel cell applications, it is particularly important to know how the clusters change with 22 ACS Paragon Plus Environment

Page 22 of 39

Page 23 of 39

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

hydration, since hydration is strongly correlated with proton conductivity. Changes in cluster morphology can be quantified by computing the surface area to volume ratio, S/V . One simple case we might expect for rod-like clusters would be swelling in only the radial direction and no appearance of new clusters with hydration. This behavior corresponds to two-dimensional scaling, S/V ∼ λ −1/2 , as explained in the Supporting Information. In fact, S/V ∼ λ −1/3 is observed, indicating that the swelling behavior is complex and three-dimensional.

Figure 9: Water content dependence of the surface area to volume ratio, normalized by its value at an arbitrary reference condition (S0 /V0 ) taken here to be λ0 = 1.3. , simulations at T = 298 K; N, simulations at T = 353 K. Dash-dotted line: expected scaling for 3D swelling (S/V ∼ λ −1/3 ); solid line: expected scaling for 2D swelling (S/V ∼ λ −1/2 ); dotted line: expected scaling for 1D swelling (S/V ∼ λ −1 ). Dashed vertical line marks the water content at which the first solvation shell of SO− 3 groups is approximately full.

3.3

Water/proton dynamics

In Figure 10, various measures of characteristic times are plotted, in particular water translational relaxation time, inverse diffusivity, and hydronium translational relaxation time. These characteristic times can also be found in a table in the Supporting Information, along with their values at very high hydration (λ ≈ 100). Figure 10 shows that all of these characteristic times behave in qualitatively the same way with respect to water content (λ ). Starting from dry conditions, characteristic times drop by several orders of magnitude. They are within one order of magnitude of their bulk values when there is enough water for sulfonate groups to no longer share water molecules in 23 ACS Paragon Plus Environment

The Journal of Physical Chemistry

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

their first coordination shell. We identified this threshold as λ ≈ 7 in Section 3.1. It is remarkable that characteristic times obtained from small-displacement transport properties, such as the translational relaxation time, exhibit the same general trend as those obtained from large-displacement properties like diffusivity, although the exact hydration threshold for divergence varies from property to property. This diverging behavior also appears in other smalldisplacement properties defined in different studies, such as the mean-residence times of H2 O and H3 O computed by Devanathan et al. 17 In general, diffusivities in porous media are sensitive to both the Å-scale molecular environment and nm-scale morphology. 58 On the other hand, translational relaxation times are not sensitive nm-scale morphology. Therefore, it appears from Figure 10 that changes in the Å-scale molecular environment with hydration significantly contribute to water transport in Nafion.

Figure 10: Comparison of various characteristic times related to water and proton transport as function as water content (λ ). , water translational relaxation time normalized by the corresponding value in the pure water (τ/τb ) from simulations at T = 353 K; , hydronium translational relaxation time normalized by the corresponding value in the dilute limit (isolated chains in water). (τ/τb ) from simulations at T = 353 K; N, inverse diffusivity of water normalized by the corresponding value in pure water (Db /D) from simulations at T = 353 K; ♦ , inverse diffusivity of water normalized by the corresponding value for pure water (Db /D) from experiments at T = 343 K; 27 O , inverse conductivity for triflic acid normalized by the inverse conductivity of pure water (σ /σb ) from experiments at T = 298 K; 59,60 D, inverse conductivity for Nafion normalized by the inverse conductivity of pure water (σ /σb ) from experiments at T = 353 K. 60,61 Dashed line marks the water content at which SO− 3 groups no longer share water molecules in their first coordination shells (see Figure 3). Solid lines are meant to guide the eye.

24 ACS Paragon Plus Environment

Page 24 of 39

Page 25 of 39

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

This hypothesis is further strengthened by the qualitative agreement of the diffusivity computed from simulations and that from experiments where their respective λ ranges overlap. At low water contents, the system could not be simulated long enough to reach the diffusive regime, and at high water contents, reliable experimental data is not available. Nevertheless, the agreement at moderate water contents proves that our simulations capture the essential features of water transport in Nafion. The agreement also reinforces our conclusion from Section 3.1 that slow chain relaxation is not important for water properties at moderate to low water contents (λ . 7). Having established the accuracy of water diffusion in our simulations, we demonstrate that the mechanism of diffusion undergoes qualitative changes as a function of hydration. Figure 11 shows the ratio of rotational to translational relaxation time, defined in Section 2. Clearly, the motion of water molecules resembles motion in bulk water at λ & 7, while at low hydration levels, rotation becomes progressively slower relative to translation. Remarkably, this behavior persists across a 55 K temperature interval. This behavior is qualitatively similar to that of supercooled liquids, when the conventional Debye model is used to quantify rotational mobility. 46 This unusual behavior can be interpreted as follows for Nafion. In Figure 12, water molecules inhabit the interface of the hydrophilic clusters, with their dipoles oriented perpendicular to the interface. As they diffuse along the interface, their orientation changes little relative to their position.

25 ACS Paragon Plus Environment

The Journal of Physical Chemistry

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Figure 11: Ratio of rotational to translation relaxation time for water normalized by its bulk value as a function of water content (λ ). , simulations at T = 298 K; N, simulations at T = 314 K; , simulations at T = 332 K; , simulations at T = 353 K. Dashed line marks the water content at which SO− 3 groups no longer share water molecules in their first coordination shells (see Figure 3). Solid black line is meant to guide the eye.

Figure 12: Top: snapshots of Nafion at λ = 1.3 and T = 353.15 K. Only hydrophilic species + [SO− 3 (orange), H2 O (purple), H3 O (cyan)] are shown. A green dot is used to track a single water molecule. Bottom: cartoons emphasizing how water molecules sit at the interface of the hydrophilic and hydrophobic phases at low water content. Left: snapshot at arbitrary time. Right: snapshot taken after approximately 100 ps have elapsed. Note that the simulations do not incorporate Grotthuss-type diffusion of protons. 62

The orientation of water molecules with respect to the interface of hydrophilic clusters can also be observed by calculating spatial distribution functions. Figure 13 shows the average spatial 26 ACS Paragon Plus Environment

Page 26 of 39

Page 27 of 39

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

distribution of different species around a water molecule at λ = 1.3 and λ = 16. At the lower hydration level, hydroniums are located preferentially near the oxygen atom, while hydrophobic CFx groups are located preferentially around the hydrogens, confirming that water molecules are oriented with their dipoles perpendicular to the interface. At the higher hydration level, hydroniums are distributed uniformly in the vicinity of all hydrogen bonding sites, while hydrophobic species tend to avoid hydrogen bonding sites, a picture that we would expect in bulk solution. Not surprisingly, at these conditions, the mechanism of water diffusion is practically bulk-like, according to Figure 11. The spatial distribution functions suggest that the arrangement of water molecules in Figure 12 is indeed representative of the ensemble.

Figure 13: Spatial distribution of either hydrophobic species (CFx , perfluoroether oxygens) or water/hydronium around a water molecule at T = 353.15 K. All visible bins have an occupancy (number density) above some arbitrary threshold. Top Left: water/hydronium distribution in Nafion at λ = 1.3. Threshold is set so that the total occupancy is 1 in the highlighted volume. Top Right: hydrophobic distribution in Nafion at λ = 1.3. Threshold is set so that the total occupancy is 5 in the highlighted volume. Bottom Left: water/hydronium distribution in Nafion at λ = 16. Threshold is set so that the total occupancy is 1 in the highlighted volume. Bottom Right: hydrophobic distribution in Nafion at λ = 16. Threshold is set so that the total occupancy is 5 in the highlighted volume. Although Figure 11 to Figure 13 suggest that the mechanism of water transport resembles that of the bulk at λ & 7, the magnitude of transport properties are still at least a factor of two greater than that of the bulk, even at λ ≈ 16, as shown in Figure 10. This effect manifests itself in the translational relaxation time, indicating that it is related to the Å-scale molecular environment of 27 ACS Paragon Plus Environment

The Journal of Physical Chemistry

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

water, rather than nm-scale morphology. In the Supporting Information, we show that all of these properties asymptotically approach the their bulk values as λ → ∞. We speculate that this unexpectedly high translational relaxation time at high water content is related to the sluggish motion of the polymer chains. It is known that this type of motion creates the Å-scale density fluctuations necessary for diffusion in homogeneous polymer/solvent systems. 63 We examine this motion in detail by first considering just the relaxation of the side chains. Figure 14 shows how relaxation times for water and the side chain at a variety of water contents and temperatures of 298-332 K collapse onto a single curve at low hydration values when plotted against each other, proving that they are strongly correlated. Intriguingly, this correlation breaks down at high hydration (λ & 7; see Figure 14), and at high enough temperatures (data not shown). Therefore, we speculate that at high hydration, the backbone rather than the side-chain slows down water dynamics. In the Supporting Information, we include a plot of relaxation times for backbone monomer segments as a function of water translational relaxation times. Although the plot suggests that these two properties are correlated at high hydration, there is too much noise in the backbone monomer relaxation times to reach a definitive conclusion, so further simulations are needed.

28 ACS Paragon Plus Environment

Page 28 of 39

Page 29 of 39

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

Figure 14: Correlation between the translational relaxation time of water, τw , and the rotational relaxation time of the side-chain, τS , at various temperatures and water contents (λ ). Quantities are normalized by their values at the same temperature at infinite dilution (λ → ∞). , T = 298 K; N, T = 314 K; , T = 332 K; Dotted line marks the value of τw /τw,b at which T = 298 and SO− 3 groups no longer share water molecules in their first coordination shells (see Figure 3). Solid lines are meant to guide the eye.

Notably, the relaxation times of both the side chains and the backbone are also far more sensitive to hydration than the relaxation time of water; the side chain and backbone data span six and four decades, respectively, as shown in Figure 14 and the Supporting Information. On the other hand, the water data span only one decade in both plots. Figure 15 illustrates this sensitivity in terms of snapshots of polymer configurations at an arbitrary time and 20 ns later. At low hydration, the chains barely move, whereas at high hydration, the chains at 20 ns are uncorrelated with those in the original configuration.

29 ACS Paragon Plus Environment

The Journal of Physical Chemistry

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Figure 15: Snapshot of Nafion backbones at arbitrary time (orange) and after 20 ns have elapsed (cyan) for a system at 353 K and different water contents. Left: λ = 1.3; right: λ = 16. The simulation box lengths are 5.53 nm and 6.28 nm, respectively. Box boundaries were periodic, but particle coordinates have been unwrapped for visualization. Interestingly, the hydronium relaxation time decreases dramatically with increasing hydration, especially at λ < 7, as shown in Figure 10, similar to the behavior of the water translational relaxation time. We speculate that like water dynamics, proton dynamics at low hydration may be strongly influenced by the Å-scale molecular environment, in addition to nm-scale morphology. It is worth noting that the experimental proton conductivity shown in Figure 10 behaves the same way, as does the experimental conductivity of concentrated solutions of triflic acid, a compound very similar to the sulfonic acid groups in Nafion that lacks nanoscale morphology and dependence on thermal history. Direct calculations of conductivity at low hydration that incorporate Grotthuss-type diffusion of protons 62 might further clarify the relative importance of Å-scale molecular environment versus nm-scale morphology in proton transport.

4

Conclusions

Molecular dynamics simulations were performed for a system of short, easily equilibrated Nafion chains across a range of temperatures and water contents, and a variety of experimentally accessible quantities were computed, in particular the sorption isotherms, diffusivities, and the hydrogenhydrogen structure factor. In the case of the sorption isotherm, good agreement was obtained with experiments at intermediate hydration (2 . λ . 7). At high hydration, the simulation isotherm 30 ACS Paragon Plus Environment

Page 30 of 39

Page 31 of 39

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

continues to sorb water up to the dilute, water-rich limit, with a corresponding activity that slowly approaches one. On the other hand, the experimental isotherm appears to saturate at a water content of λ ≈ 15, far from the dilute limit. We speculate that these differences between simulations and experiments at high activity are related to experimental difficulties reaching equilibrium near saturation. 21,25 In experiments at high activity, the slowly relaxing, high molecular weight chains resist swelling of the hydrophilic clusters, lowering the apparent water uptake. Intriguingly, although the relaxation time of the chains increases with decreasing hydration in both simulation and experiments, 24 agreement between simulation and experiments improves. We propose that the chain resistance to swelling may be overpowered by local, Å-scale water-sulfonate interactions, which become stronger with decreasing hydration. To further support this hypothesis, calculations were performed to obtain several Å-scale structural properties, as well as small-displacement (Å-scale) transport properties. In particular, the average number of waters in the first coordination shell of sulfonate groups was computed, and was found to sharply increase with increasing hydration up to intermediate hydration (λ ≈ 7), in a manner similar to the water activity. Moreover, below λ ≈ 7, nearby sulfonate groups share water molecules in their first coordination shells, a phenomenon that may lead to strong water-sulfonate interactions. At λ ≈ 7, there is enough water in the system for sulfonate groups to no longer share waters. The following small-displacement (Å-scale) transport properties were also computed: the translational and rotational relaxation time of water, the translational relaxation time of hydronium, and the rotational relaxation time of side chains. All these quantities undergo the same dramatic changes at low to intermediate water content as the water sorption isotherm and coordination shell occupancy do. Therefore, they may also reflect the strong water-sulfonate interactions at low to intermediat hydration. Their dependence on water content resembles that of experimental diffusivities, suggesting that Å-scale water-sulfonate interactions are as important for water transport as larger scale phenomena like morphology or percolation. Finally, rod-like water clusters are observed at very low humidity that consist of interdigitated

31 ACS Paragon Plus Environment

The Journal of Physical Chemistry

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

sulfonate and hydronium molecules. Water molecules inhabit the interface of these clusters, and diffuse through a mechanism in which rotations are much slower than translations relative to the mechanism in bulk water. It remains to be seen whether this behavior at low humidity can be reproduced with a more detailed model that explicitly incorporates protonation/deprotonation and more accurately predicts the sorption isotherm under these conditions. As humidity is increased, the mechanism for water diffusion becomes more bulk-like as humidity is increased, while the characteristic length scale of the clusters scales as the 1/3 power of water content, indicating complex, three-dimensional behavior.

Acknowledgements P. G. D. gratefully acknowledges the support of the National Science Foundation (Grant No. CHE1213343). A. Z. P. would like to acknowledge support for this work from the Department of Energy, Office of Basic Energy Sciences, under grant DE-SC0002128.

References (1) Costamagna, P.; Srinivasan, S. Quantum jumps in the PEMFC science and technology from the 1960s to the year 2000: Part II. Engineering, technology development and application aspects. Journal of Power Sources 2001, 102, 253–269. (2) Arcella, V.; Troglia, C.; Ghielmi, A. Hyflon Ion Membranes for Fuel Cells. Industrial & Engineering Chemistry Research 2005, 44, 7646–7651. (3) Mauritz, K. A.; Moore, R. B. State of Understanding of Nafion. Chemical Reviews 2004, 104, 4535–4586. (4) Gebel, G. Structural evolution of water swollen perfluorosulfonated ionomers from dry membrane to solution. Polymer 2000, 41, 5829–5838.

32 ACS Paragon Plus Environment

Page 32 of 39

Page 33 of 39

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

(5) Jang, S. S.; Molinero, V.; Çagin, T.; Goddard III, W. A. Effect of monomeric sequence on nanostructure and water dynamics in Nafion 117. Solid State Ionics 2004, 175, 805–808. (6) Knox, C. K.; Voth, G. A. Probing Selected Morphological Models of Hydrated Nafion Using Large-Scale Molecular Dynamics Simulations. The Journal of Physical Chemistry B 2010, 114, 3205–3218. (7) Choe, Y.-K.; Tsuchida, E.; Ikeshoji, T.; Yamakawa, S.; Hyodo, S.-a. Nature of proton dynamics in a polymer electrolyte membrane, nafion: a first-principles molecular dynamics study. Physical Chemistry Chemical Physics 2009, 11, 3892. (8) Jorn, R.; Voth, G. A. Mesoscale Simulation of Proton Transport in Proton Exchange Membranes. J. Phys. Chem. C 2012, 10476–10489. (9) Paddison, S. J.; Paul, R.; Kreuer, K.-D. Theoretically computed proton diffusion coefficients in hydrated PEEKK membranes. Physical Chemistry Chemical Physics 2002, 4, 1151–1157. (10) Cui, S.; Liu, J.; Selvan, M. E.; Keffer, D. J.; Edwards, B. J.; Steele, W. V. A Molecular Dynamics Study of a Nafion Polyelectrolyte Membrane and the Aqueous Phase Structure for Proton Transport. The Journal of Physical Chemistry B 2007, 111, 2208–2218. (11) Devanathan, R.; Venkatnathan, A.; Rousseau, R.; Dupuis, M.; Frigato, T.; Gu, W.; Helms, V. Atomistic Simulation of Water Percolation and Proton Hopping in Nafion Fuel Cell Membrane. The Journal of Physical Chemistry B 2010, 114, 13681–13690. (12) Vishnyakov, A.; Neimark, A. V. Molecular Simulation Study of Nafion Membrane Solvation in Water and Methanol. The Journal of Physical Chemistry B 2000, 104, 4471–4478. (13) Vishnyakov, A.; Neimark, A. V. Molecular Dynamics Simulation of Nafion Oligomer Solvation in Equimolar Methanol-Water Mixture. The Journal of Physical Chemistry B 2001, 105, 7830–7834.

33 ACS Paragon Plus Environment

The Journal of Physical Chemistry

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

(14) Vishnyakov, A.; Neimark, A. V. Molecular Dynamics Simulation of Microstructure and Molecular Mobilities in Swollen Nafion Membranes. The Journal of Physical Chemistry B 2001, 105, 9586–9594. (15) Rivin, D.; Meermeier, G.; Schneider, N. S.; Vishnyakov, A.; Neimark, A. V. Simultaneous Transport of Water and Organic Molecules through Polyelectrolyte Membranes. The Journal of Physical Chemistry B 2004, 108, 8900–8909. (16) Devanathan, R.; Venkatnathan, A.; Dupuis, M. Atomistic Simulation of Nafion Membrane: I. Effect of Hydration on Membrane Nanostructure. The Journal of Physical Chemistry B 2007, 111, 8069–8079. (17) Devanathan, R.; Venkatnathan, A.; Dupuis, M. Atomistic Simulation of Nafion Membrane. 2. Dynamics of Water Molecules and Hydronium Ions. J. Phys. Chem. B 2007, 111, 13006– 13013. (18) Calvo-Muñoz, E. M.; Selvan, M. E.; Xiong, R.; Ojha, M.; Keffer, D. J.; Nicholson, D. M.; Egami, T. Applications of a general random-walk theory for confined diffusion. Physical Review E 2011, 83, 011120. (19) Petersen, M. K.; Hatt, A. J.; Voth, G. A. Orientational Dynamics of Water in the Nafion Polymer Electrolyte Membrane and Its Relationship to Proton Transport. The Journal of Physical Chemistry B 2008, 112, 7754–7761. (20) Liu, J.; Suraweera, N.; Keffer, D. J.; Cui, S.; Paddison, S. J. On the Relationship between Polymer Electrolyte Structure and Hydrated Morphology of Perfluorosulfonic Acid Membranes. The Journal of Physical Chemistry C 2010, 114, 11279–11292. (21) Onishi, L. M.; Prausnitz, J. M.; Newman, J. Water-Nafion Equilibria. Absence of Schroeder’s Paradox. The Journal of Physical Chemistry B 2007, 111, 10166–10173.

34 ACS Paragon Plus Environment

Page 34 of 39

Page 35 of 39

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

(22) Baschuk, J. J.; Li, X. Carbon monoxide poisoning of proton exchange membrane fuel cells. International Journal of Energy Research 2001, 25, 695–713. (23) Hickner, M. A.; Ghassemi, H.; Kim, Y. S.; Einsla, B. R.; McGrath, J. E. Alternative polymer systems for proton exchange membranes (PEMs). Chemical Reviews-Columbus 2004, 104, 4587–4612. (24) Majsztrik, P. W.; Bocarsly, A. B.; Benziger, J. B. Viscoelastic Response of Nafion. Effects of Temperature and Hydration on Tensile Creep. Macromolecules 2008, 41, 9849–9862. (25) Benziger, J.; Bocarsly, A.; Cheah, M. J.; Majsztrik, P.; Satterfield, B.; Zhao, Q. Mechanical and Transport Properties of Nafion: Effects of Temperature and Water Activity. In Fuel Cells and Hydrogen Storage; Bocarsly, A., Mingos, D. M. P., Eds.; Structure and Bonding 141; Springer Berlin Heidelberg, 2011; pp 85–113. (26) Majsztrik, P.; Bocarsly, A.; Benziger, J. Water Permeation through Nafion Membranes: The Role of Water Activity. The Journal of Physical Chemistry B 2008, 112, 16280–16289. (27) Zhao, Q.; Majsztrik, P.; Benziger, J. Diffusion and Interfacial Transport of Water in Nafion. J. Phys. Chem. B 2011, 115, 2717–2727. (28) Li, X.; Li, F.; Shi, Y.; Chen, Q.; Sun, H. Predicting water uptake in poly(perfluorosulfonic acids) using force field simulation methods. Physical Chemistry Chemical Physics 2010, 12, 14543–14552. (29) Rubinstein, M.; Colby, R. H. Polymer Physics; Oxford University Press, USA, 2003. (30) Hwang, G. S.; Kaviany, M.; Gostick, J. T.; Kientiz, B.; Weber, A. Z.; Kim, M. H. Role of water states on water uptake and proton transport in Nafion using molecular simulations and bimodal network. Polymer 2011, 52, 2584–2593. (31) Frenkel, D.; Smit, B. Understanding molecular simulation: from algorithms to applications; Academic Press, 2002. 35 ACS Paragon Plus Environment

The Journal of Physical Chemistry

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

(32) Essmann, U.; Perera, L.; Berkowitz, M. L.; Darden, T.; Lee, H.; Pedersen, L. G. A smooth particle mesh Ewald method. Journal of Chemical Physics 1995, 103, 8577–8593. (33) Abascal, J. L. F.; Vega, C. A general purpose model for the condensed phases of water: TIP4P/2005. The Journal of Chemical Physics 2005, 123, 234505. (34) Rozmanov, D.; Kusalik, P. G. Transport coefficients of the TIP4P-2005 water model. The Journal of Chemical Physics 2012, 136, 044507–044507–13. (35) Widom, B. Some Topics in the Theory of Fluids. The Journal of Chemical Physics 1963, 39, 2808. (36) Freire, M. G.; Gomes, L.; Santos, L. M. N. B. F.; Marrucho, I. M.; Coutinho, J. A. P. Water Solubility in Linear Fluoroalkanes Used in Blood Substitute Formulations. The Journal of Physical Chemistry B 2006, 110, 22923–22929. (37) Kong, C. L. Combining rules for intermolecular potential parameters. II. Rules for the Lennard-Jones (12-6) potential and the Morse potential. The Journal of Chemical Physics 1973, 59, 2464. (38) Hess, B.; Kutzner, C.; van der Spoel, D.; Lindahl, E. GROMACS 4: Algorithms for Highly Efficient, Load-Balanced, and Scalable Molecular Simulation. Journal of Chemical Theory and Computation 2008, 4, 435–447. (39) Nosé, S. A unified formulation of the constant temperature molecular dynamics methods. The Journal of Chemical Physics 1984, 81, 511. (40) Parrinello, M. Polymorphic transitions in single crystals: A new molecular dynamics method. Journal of Applied Physics 1981, 52, 7182. (41) Urata, S.; Irisawa, J.; Takada, A.; Shinoda, W.; Tsuzuki, S.; Mikami, M. Molecular Dynamics Simulation of Swollen Membrane of Perfluorinated Ionomer. The Journal of Physical Chemistry B 2005, 109, 4269–4278. 36 ACS Paragon Plus Environment

Page 36 of 39

Page 37 of 39

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

(42) Charles H, B. Efficient estimation of free energy differences from Monte Carlo data. Journal of Computational Physics 1976, 22, 245–268. (43) Daly, K. B.; Benziger, J. B.; Debenedetti, P. G.; Panagiotopoulos, A. Z. Massively parallel chemical potential calculation on graphics processing units. Computer Physics Communications 2012, 183, 2054–2062. (44) Sedlmeier, F.; Horinek, D.; Netz, R. R. Spatial Correlations of Density and Structural Fluctuations in Liquid Water: A Comparative Simulation Study. Journal of the American Chemical Society 2011, 133, 1391–1398. (45) Connolly, M. L. Analytical molecular surface calculation. Journal of Applied Crystallography 1983, 16, 548–558. (46) Lombardo, T. G.; Debenedetti, P. G.; Stillinger, F. H. Computational probes of molecular motion in the Lewis-Wahnström model for ortho-terphenyl. The Journal of Chemical Physics 2006, 125, 174507–174507–10. (47) Williams, G.; Watts, D. C. Non-symmetrical dielectric relaxation behaviour arising from a simple empirical decay function. Transactions of the Faraday Society 1970, 66, 80–85. ˘ A molecularâA ˘ Rˇ (48) Kremer, K.; Grest, G. S. Dynamics of entangled linear polymer melts:âAL’ dynamics simulation. The Journal of Chemical Physics 1990, 92, 5057–5086. (49) Dutcher, C. S.; Ge, X.; Wexler, A. S.; Clegg, S. L. Statistical Mechanics of Multilayer Sorption: Extension of the Brunauer-Emmett-Teller (BET) and Guggenheim-Anderson-de Boer (GAB) Adsorption Isotherms. The Journal of Physical Chemistry C 2011, 115, 16474–16487. (50) Pivovar, A. M.; Pivovar, B. S. Dynamic Behavior of Water within a Polymer Electrolyte Fuel Cell Membrane at Low Hydration Levels. The Journal of Physical Chemistry B 2005, 109, 785–793.

37 ACS Paragon Plus Environment

The Journal of Physical Chemistry

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

(51) Paddison, S. J. Proton Conduction in PEMs: Complexity, Cooperativity and Connectivity. DeVice and Materials Modeling in PEM Fuel Cells 2009, 385–412. (52) Sing, K. S. W. Reporting physisorption data for gas/solid systems with special reference to the determination of surface area and porosity (Recommendations 1984). Pure and Applied Chemistry 1985, 57, 603–619. (53) Hill, T. L. An introduction to statistical thermodynamics; Courier Dover Publications, 1987. (54) Wu, X.; Wang, X.; He, G.; Benziger, J. Differences in water sorption and proton conductivity between Nafion and SPEEK. Journal of Polymer Science Part B: Polymer Physics 2011, 49, 1437–1445. (55) Gusev, A. I.; Rempel, A. A.; Magerl, A. J. Disorder and Order in Strongly Nonstoichiometric Compounds: Transition Metal Carbides, Nitrides and Oxides; Springer, 2001. (56) Lifshin, E. X-ray Characterization of Materials; John Wiley & Sons, 2008. (57) Spencer, J. B.; Lundgren, J. O. Hydrogen Bond Studies. LXXIII.* The Crystal Structure of Tritluoromethanesulphonie Acid Monohydrate, H3O+CF3SO-, at 298 and 83 K. Acta Cryst 1973, 29, 1923. (58) Mitra, P. P.; Sen, P. N.; Schwartz, L. M.; Le Doussal, P. Diffusion propagator as a probe of the structure of porous media. Physical Review Letters 1992, 68, 3555–3558. (59) Corkum, R.; Milne, J. The density, electrical conductivity, freezing point, and viscosity of mixtures of trifluoromethanesulfonic acid and water. Canadian Journal of Chemistry 1978, 56, 1832–1835. (60) Light, T. S.; Licht, S. L. Conductivity and resistivity of water from the melting to critical point. Analytical Chemistry 1987, 59, 2327–2330.

38 ACS Paragon Plus Environment

Page 38 of 39

Page 39 of 39

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

(61) Yang, C.; Srinivasan, S.; Bocarsly, A. B.; Tulyani, S.; Benziger, J. B. A comparison of physical properties and fuel cell performance of Nafion and zirconium phosphate/Nafion composite membranes. Journal of Membrane Science 2004, 237, 145–161. (62) Marx, D.; Tuckerman, M. E.; Hutter, J.; Parrinello, M. The nature of the hydrated excess proton in water. Nature 1999, 397, 601–604. (63) Duda, J. L.; Vrentas, J. S.; Ju, S. T.; Liu, H. T. Prediction of diffusion coefficients for polymersolvent systems. AIChE Journal 1982, 28, 279–285.

Figure 16: Table of Contents Image

Supporting Information Available: Details of the force field; details of the methodology for obtaining the Connolly surface and heats of sorption; polymer backbone relaxation times; densities; water/hydronium transport properties; radial distribution functions; additional simulation results for varying chain lengths and side chain placements. This material is available free of charge via the Internet at http://pubs.acs.org.

39 ACS Paragon Plus Environment