Theoretical Studies of Proton Transfer in Water and Model Polymer

Results from small-angle X-ray measurement of Nafion10 have verified the ... Density functional theory (DFT) studies of proton transfer from solution ...
0 downloads 0 Views 208KB Size
Download PDF
Ind. Eng. Chem. Res. 2001, 40, 4789-4800

4789

Theoretical Studies of Proton Transfer in Water and Model Polymer Electrolyte Systems Tao Li, Aaron Wlaschin, and Perla B. Balbuena* Department of Chemical Engineering, University of South Carolina, Columbia, South Carolina 29208

The interactions of H3O+ with water and model Nafion structures are studied using ab initio, density functional theory, and molecular dynamics simulations. In the gas phase, H3O+ is solvated by a well-defined first solvation shell with three water molecules, while an additional water molecule locates in the second shell. This same structure is preferred in the liquid phase. In the absence of an electric field, proton transfer between H3O+ and a water molecule proceeds without barrier if the hydrogen bond length O-H‚‚‚O is smaller than 2.7 Å. When an electric field is applied in a direction opposite to that of the system dipole, the activation barrier for proton transport is significantly reduced. The Nafion side chain, ending in a sulfonic group, is folded in the gas phase, whereas in aqueous media, it becomes stretched, and proton transfer takes place from the acid group to the water molecules. Results from MD simulations indicate that the contact ion pair formed between the sulfonic acid anion -SO3- and the hydronium ion is very stable. Introduction Nafion, a perfluorosulfonic acid (PFSA) membrane developed by DuPont, has been reported as chemically stable, structurally durable, and proton-conductive.1 One of its well-known applications is in low-temperature-operation polymer membrane electrolyte fuel cell systems. A typical PFSA has a hydrophobic backbone that consists of poly(tetrafluoroethylene) (PTFE) and hydrophilic perfluorovinyl ether pendant side chains terminated by sulfonate ionic groups.2 The ionic groups comprise up to 15% of the polymer.3 The molecular structure of H-Nafion is4

where n can vary between 5 and 13.5 and 100 < m < 1000.1 The polymer equivalent weight, defined as the weight of dry polymer in grams containing 1 mol of unit charge, is usually in the range of 1000-2000. The -SO3H groups tend to dissociate in polar solvents, and they protonate the solvent. This protonated solvent species serves as the major charge carrier in the membrane. Experimental data suggest that the negatively charged -SO3- groups, along with water molecules and positively charged counterions such as H+, Na+, or K+ depending on the pretreatment, tend to aggregate and form hydrophilic clusters.5,6 These clusters are dispersed in the continuous hydrophobic fluorocarbon matrix. The cluster network model has been proposed by Hsu and Gierke to describe the membrane microstructure.7-9 According to this model, the membrane is formed by clusters and channels, with the ionic groups located at the inner surface of the 3-5-nmdiameter clusters. As the degree of hydration increases, * Corresponding author: Perla B. Balbuena. Phone: (803) 777-8022. Fax: (803) 777- 8265. E-mail: [email protected].

the clusters swell and increase in size. Results from small-angle X-ray measurement of Nafion10 have verified the existence of clusters of about 4-nm diameter, with channels of 1-nm diameter and length that limit the transport properties of water and ions. The cluster model also implies a rather curled polymer side chain in the membrane. Other researchers have proposed that the sulfonic acid domains are parallel to each other, that is, that a lamellar geometry exists in the membrane, and that the side chain is much more stretched.11 Yeager et al.12 suggested that the hydrated membrane is formed by two continuous interpenetrating phases with irregularly shaped clusters. Spectroscopic studies suggest that the size of the clusters should be smaller than that proposed by Gierkes’ model.13 The actual microstructure has not been definitively resolved, and significant debate still exists about it.14 In fuel cells that utilize hydrogen as a source of fuel, solvated proton clusters serve as the major charge carriers. The proton transfer inside the membrane is usually described by the vehicle mechanism15-17 and by the hopping mechanism,15,18,19 also known as the Grotthuss model.20 As opposed to proton hopping, the vehicle mechanism proposes that protons bounded with a fraction of water molecules travel through the electrolyte via mutual diffusion. To some extent, proton solvation/desolvation can be treated as a special case of the proton-transfer reaction. Many experiments have been performed to determine the proton solvation energy,21 and studies agree that this energy release is mostly due to the formation of the hydronium ion.22 Three coordination water molecules were found by ab initio calculations (HF/6-31G* and MP2/6-31G*),23 mass spectroscopy,24 and ultrafast solvation dynamics experiments,25,26 whereas some earlier diffraction and X-ray experiments reported the presence of a fourth water molecule in the first ionic shell of the hydronium ion.27 The hopping mechanism suggests that the proton propagates along the O-H‚‚‚O hydrogen bond of the H3O+(H2O) complex. Results from a mixed molecular dynamics/Monte Carlo (MD/MC) algorithm used to

10.1021/ie010467y CCC: $20.00 © 2001 American Chemical Society Published on Web 10/09/2001

4790

Ind. Eng. Chem. Res., Vol. 40, No. 22, 2001

describe the proton transfer by an instantaneous hopping mechanism19 indicate that the proton migration rate is mainly determined by the lifetime of the solvation complexes. Density functional theory (DFT) studies of proton transfer from solution to a metal cluster18 reveal that the presence of an external electric field could effectively accelerate the proton-transfer rate. The proton transport properties are strongly coupled with the water content and the microstructure of the material. The experimentally determined proton diffusion coefficient in a Nafion-117 membrane is 1.4 × 10-9 m2/s,2 whereas in solution at infinite dilution, it is 9.31 × 10-9 m2/s.28 The presence of the hydrophilic pendant chain -OCF2CF(CF3)OCF2CF2SO3- imposes a dramatic influence on the solvent structure and the associated transport properties; thus, structural information regarding the pendant chain is important. Because of the importance of proton transfer in fuel cell design, studies of proton solvation and transport processes at the microscopic level are essential. MD simulations29 have been used to study the motion of water molecules and protons in a Nafion membrane,30 modeled as a large number of identical cylindrical pores. The calculated electro-osmotic drag coefficients were overestimated because of the simplified physical model used. Only when the wall charge density was comparatively low was the proton diffusion coefficient obtained found to be of the same order of magnitude as the experimental value. Applications of a statistical mechanical model31 and a phenomenological model32 were also reported recently. In this work, we first examine the minimum-energy structure and solvation energies of the H+(H2O)n (n ) 1-5) complexes using DFT. Next, we investigate the proton-transfer mechanism by analyzing the potential energy surface and its corresponding activation barrier, as well as the effect of an external electric field on proton transfer between water molecules, as a function of the oxygen-oxygen distance. The proton-transfer reaction from the molecular unit CF2dCFOCF2CF(CF3)OCF2CF2SO3H to the water molecule is investigated. The fully optimized geometry of the ionized molecular unit CF2dCFOCF2CF(CF3)OCF2CF2SO3- is used in MD simulations, where we examine the membrane fragment conformation and proton-transport process in bulk water and in a model polymer/water system. 2. Computational Details 2.1. Ab Initio Calculations. Ab initio and DFT calculations were carried out using Gaussian 9433and Gaussian 98.34 The basis sets used in the calculations included both 6-31+G (d, p) and 6-311++G (d, P), which have polarization and diffusion functions on the oxygen atoms. Full geometry optimizations of the solvation clusters were carried out using the B3P8635 and B3PW9136,37 functionals. The calculations of potential energy curves for proton transfer were carried out only at the B3P86/6-31+G (d, p) level. A finite external electric field (EF) was added to the calculation. The field was generated by a dipole, and the sign of the field was chosen so that a positive sign was in the direction parallel to the original dipole direction of the complex. The EF intensities were in the range of 0.001-0.02 au, where 1 au corresponds to 5.14 × 10 9 V/cm. Interactions of the H-terminated side chain [represented by CF2dCFOCF2CF(CF3)OCF2CF2SO3H] with water were studied at both the HF/6-31+G** and

B3PW91/6-311++G** levels of theory. The ionic side chain, represented by CF2dCFOCF2CF(CF3)OCF2CF2SO3-, was fully optimized at the HF/6-31+G** level. 2.2. Molecular Dynamics Simulations. MD simulations were carried out using Cerius2.38 Except for the bulk water simulations, all of the other systems under consideration contain 500 water molecules, x hydronium ions, and x SO3--terminated Nafion oligomers, where x varies between 0 and 2. The Nafion monomer and oligomer were constructed using the Cerius 2 polymer builder, and an initial minimization was done using a combination of the steepest-descent, Newton-Raphson, quasi-Newton, and truncated-Newton methods.38 We have chosen n ) 7, which gives an equivalent weight (EW) of 1165 g of polymer per equivalent of ionic groups in the H+ form. This structure corresponds to Nafion1200, one of the most extensively studied membranes. Unless otherwise specified, the simulations were carried out in the canonical ensemble, i.e., at constant number of molecules N, volume V, and temperature T. The force field used has the functional form of the Dreiding force field39 that has been reported to be effective for polymer materials. Long-range contributions of the Coulomb interaction are calculated with the Ewald summation. The nonbonding potential includes also the short-range Lennard-Jones terms. Intramolecular terms such as bond stretching, bending, and torsion are included. The equilibrium bond distance for the Dreiding force field is assumed to follow a simple additivity of bond radii, where the bond radii are based on structural data of standard reference molecules. The equilibrium bond angle is assumed to depend only on the central atom.39 Periodic boundary conditions are applied. The simulations are performed at 298.15 and 353.15 K. The temperature is maintained using a T-damping thermostat.40 The Dreiding force field39 has been also applied to water, using the SPC/E charges.41 To test the validity of the Dreiding force field in reproducing the thermodynamic and structural properties of bulk water at ambient conditions, molecular dynamics simulations of 250 water molecules in the NPT ensemble were carried out. The temperature is controlled at 298.15 K, and the pressure is 1 bar. The total length of these simulations was 100 ps, and the last 50 ps of data were collected for statistical analysis. Another set of MD simulations was done to investigate the dynamics of H3O+. To obtain an accurate result for the proton dynamics, a small time step, 1 fs, was used in the simulations. The first 50 ps were usually discarded, and then an additional 200 or 250 ps were run for the production stage. The hydronium ion diffusion coefficient is calculated on the basis of the mean square displacement (MSD) as29

〈|r b(t) - b r (0)|2〉 tf∞ 6t

DMSD ) lim

where r is the position vector and t is the time. As this equation indicates, the calculated diffusion coefficients become more accurate with long simulation times. For short simulation times, the MSD shows a nonlinearity as a function of time, and this leads to an overestimation of the diffusion coefficient. If the time evolution of the MSD is plotted on a logarithmic scale, the slope for Einstein diffusion should be 1.42 Therefore, for each data analysis, we plot the change of the MSD with time on a logarithmic scale, and determine the region where the

Ind. Eng. Chem. Res., Vol. 40, No. 22, 2001 4791

Figure 1. Calculated gas-phase (B3PW91/6-311++G**) proton solvation complexes H+(H2O)n and corresponding atomic Mulliken charges. (a) n ) 1, (b) n ) 2, (c) n ) 3, (d) n ) 4, (e) n ) 5.

slope is 1. The diffusion coefficient is calculated using that range of data. 3. Results and Discussions 3.1. Proton Solvation in Water. Figure 1 displays DFT-optimized structures and representative Mulliken atomic charges of the H+(H2O)n clusters with n ) 1-5; the corresponding equilibrium bond lengths are reported in Table 1. The rO‚‚‚O parameter is defined as the distance between the H3O+ oxygen and the oxygen of the hydrogen-bonded water molecule, rO-H is the O-H bond length in H3O+ or H2O, and rH‚‚‚O is the H‚‚‚O distance in the OH‚‚‚O hydrogen bond between H3O+ and a water molecule. In general, the optimized structures are in good agreement with other calculated and experimental results (Table 1). The two different levels of theory (functional/basis set) tested in this study, B3P86/6-31+G** and B3PW91/6-311++G**, give comparable geometries, although the O-H bond lengths are slightly longer at the B3P86 level, whereas the O‚‚‚O distances are a little bit shorter than those from B3PW91. For both methods, the calculated OH bond length in H3O+ is 0.017 Å longer than the OH bond length in H2O. The H+(H2O)2 structure (Figure 1b) is highly symmetric, with the proton located exactly between the two water molecules. This is because the proton is forming overlaps between the two lone-pair electron orbitals belonging to two distinct water oxygen atoms separated

by 2.39 Å from each other. The linear H bond formed between H3O+ and H2O lengthens the O-H bond in H3O+. Ojama¨e and co-workers43,44 reported two different conformations for H5O2+ based on ab initio MP2 calculations: a structure of C2 symmetry with a symmetric hydrogen bond, (H2O‚‚‚H‚‚‚OH2)+ and a structure of Cs symmetry with an asymmetric hydrogen bond, H3O+‚‚‚H2O. The first conformation is viewed as the global minimum. Our minimum-energy structure reported in Figure 1b is in agreement with this global minimum: all bond distances are in close agreement with those from the MP2 results. Our calculations also indicate that, at the B3PW91 level, the H3O+‚‚‚H2O complex of Cs symmetry is not stable. To locate the minimum-energy structure for the H+(H2O)3 cluster, the two water molecules are initially located in symmetric positions with respect to H3O+. After geometry optimization, this symmetry with respect to the hydronium ion remains (Figure 1c). Each O-H bond in H3O+ that is part of a H bond with a water molecule is stretched to 1.05 Å, and this length is 0.14 Å smaller than the corresponding bond length in the H+(H2O)2 cluster (Table 1, B3PW91/6-311++G**). The optimized geometry of the H+(H2O)4 complex is shown in Figure 1d. This structure is close to the global minimum reported by Oja¨mae et al.43 Our rH‚‚‚O distance (1.51 Å) is slightly shorter than the value reported by Oja¨mae (1.56 Å). In Figure 1e, the minimum-energy structure of the H+(H2O)5 complex is shown; however,

4792

Ind. Eng. Chem. Res., Vol. 40, No. 22, 2001

Table 1. Key Structural Parameters for the Proton Hydrated Clustersa system

rO-H (Å)

rO‚‚‚O (Å)

rH‚‚‚O (Å)

B3P86/6-31+G** calculations H2O

0.963

H+(H2O)

0.980

H+(H

0.969 (O1-H3)

2O)2

2.388

1.192 (O1-H7) 1.192 (H7-O4)

H+(H2O)3 0.970 (O1-H2) 2.458 (O1-O8) 1.051 (H7-O8) 0.967 (O8-H10)

1.409 (H7-O1)

H+(H2O)4 0.965 (O1-H2) 0.965 (O4-H5) 1.020 (O8-H9) [1.017]27

1.508 (H9-O4)

2.527 (O4-O8) [2.52]27

H+(H2O)5 0.964 (O1-H2) 2.566 (O1-O8) 1.559 (O4-H9) 1.009 (H7-O8) 2.446 (O8-O11) 1.556 (H7-O1) 1.063 (O8-H10) 2.626 (O11-O14) 1.384 (H11-O10) H2O

B3PW91/6-311++G** calculations 0.960

H+(H2O)

0.978

H+(H2O)2 0.967 (O1-H3)

2.389

1.192 (O1-H7) 1.196 (H7-O4)

H+(H2O)3 0.967 (O1-H2) 2.478 (O1-O8) 1.042 (H7-O8) 0.965 (O8-H10)

1.436 (H7-O1)

H+(H2O)4 0.962 (O1-H2) 0.963 (O4-H5) 1.013 (O8-H9) [1.017]27

1.538 (H9-O4)

2.547 (O4-O8) [2.52]27

H+(H2O)5 0.962 (O1-H2) 2.584 (O1-O8) 1.585 (O4-H9) 1.003 (H7-O8) 2.460 (O8-O11) 1.582 (H7-O1) 1.053 (O8-H10) 2.662 (O11-O14) 1.408 (H11-O10) a Atom numbers in parentheses correspond to Figure 1. Experimental results are in brackets.

it cannot be conclusively identified as the global energy minimum. The optimized structures of H+(H2O)4 and H+(H2O)5 suggest that H3O+ is solvated by three water molecules in the first shell, with a fourth water molecule located in the second shell. This observation agrees with an ab initio study23 at a Moller-Plesset level (MP2/ 6-31G*)45 that uses a continuum model46 to represent the hydronium hydration, whereas a Hartree-Fock (HF) study47 indicated the possibility of the fourth water being in the first solvation shell. Recent ab initio MD calculations suggest that the most favorable structures for hydrated protons in water would be H+(H2O)2 and H+(H2O)3.48 With the increase of the solvation cluster size, the hydronium ion rO-H values decrease, whereas the rO‚‚‚O distance increases. The rH‚‚‚O distances also increase, except when a water molecule is located in the second solvation shell. Therefore, the H-bond strength becomes weaker when the cluster size increases. This change of bond strength can be related to the amount of charge transferred between the proton and the solvation water molecules. In Figure 1b, each solvation water molecule bears a charge of 0.21 au, whereas the charges are reduced to 0.07 and 0.05 au for Figure 1c and d, respectively. For Figure 1e, barely zero charge transfer occurs between the fourth solvation water and the H(H2O)3+ ion. Thus, the attractive interactions of the proton with surrounding water molecules diminish as more water molecules are added. The atomic charges in H3O+ also change significantly as the number of solvation water molecules increases. The charge on the oxygen atom becomes strongly negative, changing from

-0.20 au in H3O+ to -0.69 au in H(H2O)5+, whereas the charge on the H atom interacting with water decreases from 0.57 au in H(H2O)2+ to 0.49 au in H(H2O)5+. The charge on the noninteracting hydronium H atom (such as H10 in Figure 1c) also decreases from 0.40 au in H3O+ to 0.32 au in H(H2O)5+. Table 2 lists the proton solvation energies ∆E, proton solvation enthalpies ∆H, and proton solvation free energies ∆G, defined as the corresponding energy differences between the product complexes and the reactants for the reactions

H+ + nH2O S H+(H2O)n All of the energies and energy differences are slightly higher at the B3P86 theory level. The solvation reactions are all spontaneous, as indicated by their ∆G values. Literature-reported values for proton solvation enthalpies are between -254 and -260.5 kcal/ mol,21,22,49,50 which correspond to the calculated values for n ) 4-5 solvation water molecules (Table 2). It can also been seen that the initial formation of the hydronium ion contributes most to the total solvation energy of H+(H2O)n. The other part of the solvation energy is due to interactions of the hydronium ions with the solvation water molecules. From our calculations, the enthalpy for hydronium ion formation is -165.07 kcal/ mol, and Meot-Ner et al.22,50 reported a measured value of -166.3 kcal/mol. The energy difference between H+(H2O)n and H+(H2O)n-1 becomes smaller with increasing number of hydration molecules. However, even the second hydration shell produces a significant increase (∼13 kcal/mol) of the solvation energy. Energies of formation for complexes such as H+(H2O)n are calculated by subtracting the energies of the isolated monomers from the energy of the complex. In the calculation of the complex energy, the bonded monomers are described with the basis set of the whole complex, while the energies of the isolated monomers are calculated using the monomer basis sets. The energies of formation will thus be overestimated. This type of error is called the basis set superposition error (BSSE), which is known to be relatively important in systems containing hydrogen bonds.51 The standard Boys and Bernardi counterpoise procedure,52 where the energies of the isolated monomer are calculated with the basis set of the complex, has been followed. The BSSE-corrected free energies of solvation (B3PW91/6-311++G**) are listed in Table 2. Note that, without the BSSE correction, the calculated solvation free energies are comparatively higher than the experimental results,22 with errors within 5% of the experimental values. 3.2. Proton Transfer between Two Water Molecules: The Hopping Mechanism. According to the proton hopping mechanism,15 the proton propagates along the O-H‚‚‚O hydrogen bond direction of the H3O+(H2O)n complexes. This process is shown in Figure 2a, and the formal reaction is

H3O+ + H2O T H2O + H3O+ The coordinates used to construct the potential energy curves for the above reaction are the O-O distance R1 and the proton-oxygen distance R2 (Figure 2b). For each curve, R1 is kept constant while R2 is increased in increments of 0.01 Å. The calculated potential energy curves for several values of R1 are displayed in Figure

Ind. Eng. Chem. Res., Vol. 40, No. 22, 2001 4793 Table 2. Summary of Thermodynamic Properties for the Proton Hydration Clustersa E + ZPE (hartree)

system

E (hartree)

ZPE (hartree)

H2O H+(H2O)

-76.610 55 -76.887 38

0.021 49 0.034 59

H+(H2O)2 H+(H2O)3 H+(H2O)4

-153.559 73 -230.210 98 -306.854 66

0.057 95 0.083 43 0.108 90

B3P86/6-31+G** calculations -76.589 07 -76.852 79 -165.49 -165.43 [-166.3]22 -153.501 78 -203.09 -204.25 -230.127 55 -226.12 -227.72 -306.745 76 -244.41 -246.35

H+(H2O)5

-383.489 75

0.133 26

-383.356 49

H2O H+(H2O)

-76.428 22 -76.704 45

0.021 48 0.034 53

H+(H2O)2 H+(H2O)3 H+(H2O)4

-153.191 21 -229.658 03 -306.117 83

0.057 60 0.083 15 0.108 24

-153.133 61 -229.574 89 -306.009 60

-200.85 -222.50 -240.03

-201.97 -223.95 -241.68

H+(H2O)5

-382.568 48

0.132 74

-382.435 74

-252.19

-254.35

∆E (kcal/mol)

-258.01

∆H (kcal/mol)

-260.39

B3PW91/6-311++G** calculations -76.406 771 -76.669 92 -165.13 -165.07

∆G (kcal/mol)

BSSE (kcal/mol)

∆G + BSSE (kcal/mol)

0.66

-165.41

10.23 13.29 13.31

-183.96 -195.56 -203.91

11.94

-211.62

-166.43 [-158.9]22 -196.4022 -212.12 -222.76 [-205.8]22 -228.94 [-211.4]22

-166.07 [-158.9]22 -194.19 -208.85 -217.23 [-205.8]22 -223.55 [-211.4]22

a Enthalpies and free energies are at 298.15 K. Note that absolute energies (E) and zero-point energies (ZPE) are in hartrees, whereas relative energies (∆E), enthalpies (∆H), and Gibbs free energies (∆G), as well as basis set superposition error (BSSE) corrections, are in kilocalories per mole.

Figure 2. (a) Schematic illustration of Grotthuss proton hopping mechanism. (b) Coordinates used for potential energy surface scan.

3. Only when R1 is small enough (R1 ) 2.4 Å) does the the potential energy curve show a single minimum; all other curves have symmetric double-well shapes. This observation is in agreement with MP2 and B3LYP results.43,44 The activation energy (Table 3), defined as the energy difference between the maximum and minimum points of the curve, decreases significantly as R1 decreases. At fixed R1, the left minimum (Figure 3) corresponds to the equilibrium of a left water molecule with a right H3O+, whereas the symmetric minimum located at the right corresponds to a left H3O+ in equilibrium with a right water molecule. A transition state for proton transfer corresponds to the maximum point in Figure 3, i.e., a proton is located equidistant between the two water molecules. From Figure 3, we observe that, with the increase of the R1 distance, the activation barrier for proton transfer increases sharply, and at a certain small R1 value, the transfer is barrierfree. As pointed out earlier,19 because the zero-point energy E0 of a O-H vibration is 5.3 kcal/mol,16 as long as the O-O distance is smaller than 2.7 Å, the barrier will be under E0, and the proton will be able to move

Figure 3. Potential energy surface for proton transfer between two water molecules as a function of the OH distance (R2) and parametric in the O-O distance (R1), as defined in Figure 2b. The energy E is the energy difference between the energy of the H+(H2O+)n complex and the energies of the isolated monomers (H3O+ and H2O).

freely between the two solvation water molecules. In these calculations, we do not consider the dynamics of the process, which can also influence the proton transfer. Other studies suggest that the rate-determining step in proton diffusion in water is the reorientation of the H2O molecules in the first hydration shell of H3O+ to allow proton motion through the water molecules.16 The proton transfer between the H3O+‚‚‚H2O complex then would take place by quantum mechanical tunneling.16 To investigate the effect of an electric field (EF) on proton transfer, calculations were repeated at O-H‚‚‚O distances R1 of 2.7 and 2.8 Å (Figure 4a and b, respectively). The corresponding activation barriers are listed in Table 3. From Figure 4a, it can be seen that even a low-intensity EF destroys the symmetry of the PES curves reported in Figure 3. In our model, an electrical field on the x axis is created by a dipole that

4794

Ind. Eng. Chem. Res., Vol. 40, No. 22, 2001

Table 3. Activation Barriers Estimated from a Potential Energy Surface Scan for Proton Transfer between Two Water Molecules R1 (O-O distance)a (Å)

electric field (V/cm)

barrier (kcal/mol)

2.40 2.55 2.60 2.65 2.70 2.80 2.90 3.10 3.50 2.70 2.70 2.80 2.80 2.80 2.80

0 0 0 0 0 0 0 0 0 -5.14 × 106 -25.70 × 106 -5.14 × 107 -10.28 × 107 +5.14 × 107 +10.28 × 107

0 0.25 0.94 1.88 3.21 6.72 10.86 20.75 42.44 2.52 0.55 1.02 0 15.17 -

a The O-O distance is one of the scan coordinates, shown in Figure 2b.

When a negative EF is applied (Figure 4a), both PES minima are shifted to larger R2 distances, and the maximum is shifted to smaller R2 distances. As schematically shown in Figure 4a, a negative field favors the proton transfer, with the left minimum becoming less favorable as the field EF increases because of the proximity of H3O+ to the positive pole. Note that the corresponding energy increases as the field intensity increases. Therefore, after proton transfer, the system stabilizes at a much lower energy minimum, and consequently, the activation barrier for proton transfer decreases significantly as the field intensity increases, thus favoring the hopping transfer mechanism (Table 3). For example, at R1 ) 2.8 Å (Figure 4b), when the EF intensity is -10.28 × 107 V/cm, the barrier for proton transfer disappears, and the proton can move freely from H3O+ to the other solvation water. When the field direction is reversed (EF > 0), the energy of the system increases, and the barrier increases sharply (Table 3), with a very stable minimum corresponding to the H2O‚ ‚‚H3O+ complex and a large activation barrier for proton transfer, which increases with the field intensity. Therefore, an appropriate external electric field lowers the activation barrier for the proton-transfer reaction and increases its reaction rate. The reorientation of the H3O+(H2O) complex as a result of the applied field was not considered in the calculations, and inclusion of this effect will possibly change the absolute value of the barriers. It should also be pointed out that, although our calculations of the EF effect have been performed only at the O-H‚‚‚O distances R1 of 2.7 and 2.8 Å, a monotonic increase of the proton-transfer barrier with the R1 distance can be expected. 3.3. Structural Properties of the Nafion-Hydronium-Water System. The presence of water causes dissociation of the -SO3H group of H-Nafion and water protonation, i.e., protons are transferred from H-Nafion to water, thus forming the -SO3-‚‚‚H3O+ ion pair

-SO3H + H2O f -SO3-‚‚‚H3O+

Figure 4. (a) Potential energy surface for proton transfer between two water molecules in the presence of an external electric field (EF) at R1 ) 2.70 Å. The energies are normalized with respect to the minima at R1 ) 2.70 Å and EF ) 0. (b) Same at R1 ) 2.80 Å, with energies normalized with respect to the minima at R1 ) 2.80 Å and EF ) 0.

is aligned with this x axis. We define a positive field as one that goes from positive to negative along the positive x-axis direction. Conversely, when the field is negative, it goes from positive to negative along the negative x-axis direction, i.e., EF always points in the direction opposite to the dipole that originates it. In our system, the positive x axis has the direction O1O4 (Figure 2b).

The solvation of the ion pair to form -SO3-‚‚‚H3O+(H2O)n complexes, the natural next step in aqueous media, is reported to be irreversible.53,54 The ion-pair stability was investigated by optimizing (HF/6-31+G**) the neutral molecular complex structure -OCF2CF(CF3)OCF2CF2SO3H‚‚‚H2O attached to a CF2-CF2 group representing a small portion of the backbone, as shown in Figure 5a. In the initial configuration for geometry optimization, the water molecule is located close to the -SO3H group, as previous ab initio studies reported water to be preferentially located at the terminal part of the side chain.55 It was found that the complex, formed mainly by H-bond interactions, is more stable than the -SO3-‚‚‚H3O+ ion pair; therefore, at this level of theory, there is no evidence of proton transfer from the side chain to the water molecule. The O-H bond (O21-H27) of the -SO3H group is stretched by about 0.02 Å because of the H-bond interaction with the water molecule. The (-SO3)H‚‚‚O(H2) distance is 1.743 Å. The binding energy of the neutral complex is -10.28 kcal/ mol. The enthalpy of complex formation is -10.59 kcal/ mol, and the free energy change is -2.0 kcal/mol, which indicates that the formation of the molecular complex is thermodynamically favorable. The entropic term is negative and large, T∆S ) -8.59 kcal/mol, indicating the formation of a highly ordered molecular complex.

Ind. Eng. Chem. Res., Vol. 40, No. 22, 2001 4795

Figure 6. HF/6-31+G**-optimized structure of Nafion pendant anion chain. Table 4. Geometric Parameters for Nafion Monomer Anion Depicted in Figure 6 bond length (Å)

Figure 5. (a) HF/6-31+G**-optimized structure of H-Nafion/ water complex. (b) B3PW91/6-311++G**-optimized structure of H-Nafion/water complex indicating proton transfer from the acid group to the water molecule.

It is known that, although the HF geometric parameters (bond lengths, angles, dihedral angles) are in excellent agreement with experimental results for many systems,56 the HF method cannot capture some states where electron correlation effects are important. For this reason, we re-optimized the same system using DFT. Figure 5b illustrates the results of a geometry optimization (B3PW91/6-311++G**) of the same system in Figure 5a, and in its initial configuration, the water molecule is located close to the -SO3H group. The O-H bond O21-H27 of the -SO3H group is now stretched to 1.009 Å because of the H-bond interaction, and although the H27-O28 distance is still large (1.616 Å), the proton transfer from the -SO3H group to H2O can be clearly seen. The charge of H27 increases to 0.51 au, compared to a value of 0.31 au when no water is present. This result illustrates the effects of electron correlation and large basis sets on ion-pair association. To understand the solvent effect on the chain conformation, we analyzed the anion structure. The fully optimized (HF/6-31+G**) geometry of the ionized Nafion side chain CF2CF-OCF2CF(CF3)OCF2CF2SO3- is depicted in Figure 6. Corresponding geometric parameters are summarized in Table 4. Most of the parameters are in agreement with the reported calculated data.14 The distance of the -SO3- group to the backbone (CF2-CF2) measured between S18 and C22 (Figure 6) is 4.33 Å, and overall, the side chain is quite folded, with an estimated diameter of around 8 Å, in agreement with previous B3LYP/6-31G** calculations.55 On the basis of this gas-phase folded structure, it has been argued that

S-O C-S C-C C-O C-F

1.44 1.85 1.54 1.36 1.33

angle (°) OSO OSC SCF CCS COC FCO CCF CCC FCF

115.72 102.65 110.30 116.81 125.23 109.57 110.18 113.87 106.71 108.64

Table 5. Thermodynamic Properties of Bulk Water Determined From NPT Simulations at 298 K and 1 bar Using the Dreiding Force Field density Evdw Ecoulomb EH bond ∆hvap D (g/cm3) (kcal/mol) (kcal/mol) (kcal/mol) (kcal/mol) (10-9 m2/s) 0.987

2.616

-12.242

-4.777

14.403

3.350

Gierke’s model,57 a 30-50 Å cluster consisting of 70 sulfonated side chains and 1000 water molecules in a fully hydrated Nafion membrane, was impossible.55 In the next section, we report the results of MD simulations to elucidate this point. 3.4. Dynamic Properties of the Nafion-Hydronium-Water System. To verify the accuracy of the Dreiding force field in representing liquid water, simulations are first performed for ambient water in the NPT ensemble. The properties obtained are summarized in Table 5. The average density of the system is 0.987 g/cm3, within 1% of the experimental value.58 The heat of vaporization59,60 is calculated as

∆hvap ) - Epot(inter) where Epot(inter) includes all contributions from the nonbonded intermolecular interactions to the potential energy. The calculated value of ∆hvap is 14.403 kcal/mol, which is 3.893 kcal/mol higher than the experimental value. As seen from Table 5, the inclusion of H-bond terms in the potential energy accounts for this overestimation. The calculated diffusion coefficient is 3.35

4796

Ind. Eng. Chem. Res., Vol. 40, No. 22, 2001

Table 6. Comparison of the Calculated O-O, O-H, and H-H Radial Distribution Functions of Water at 298 K and 1 bar Using the Dreiding Force Field with Those of the SPC/E Model and Experimental Values first peak position

first minimum position

second peak position

second minimum position

model

rOO (Å)

rOH (Å)

rHH (Å)

rOO (Å)

rOH (Å)

rHH (Å)

rOH (Å)

rHH (Å)

rOH (Å)

rHH (Å)

Dreiding SPC/E 64 expt64

2.89 2.90 2.80

1.95 1.95 1.85

2.67 2.50 2.40

3.41 3.40 3.45

2.39 2.40 2.45

3.13 3.10 3.00

3.31 3.30 3.30

3.92 3.85 3.85

4.30 4.40 4.40

4.70 4.70 4.72

Table 7. Mulliken Atomic Charges (HF/6-31+G**) Used in MD Simulations

a

molecule

atom

H2O

O H

-0.84 0.42

q (e)

H3O+

O H

-0.49 0.49

Nafiona

O19 O20 O21 O1 O11 S18 F8 F15 F17 C2 C5 C6 C12 C13

-0.93 -0.81 -0.83 -0.53 -0.19 2.08 -0.40 -0.39 -0.45 0.94 -0.11 1.70 0.83 0.58

Nafion atom numbers correspond to Figure 6.

× 10-9 m2/s, overestimated in comparison with the 2.3 × 10-9 m2/s experimental result61 at the same temperature. We have also computed the radial distribution function (RDF), which provides structural information. The first peak and the first minimum positions (Table 6) of the distribution functions gOO(r), gOH(r), and gHH(r) agree well with the experimental RDFs and SPC/E water model results.62-64 The most important structural discrepancy between the computed and experimental RDF is that the peak height of the computed gOO(r) is overestimated with respect to the experimental value.64 Overall, we consider the performance of the Dreiding model fair for our purposes of testing the relative mobility of the hydronium in water with and without the presence of Nafion functional groups. Thus, to further investigate the proton-transport process and conformational changes of the hydrophilic Nafion side chain in the bulk solvent, MD simulations for both the water/H3O+ and the water/SO3--terminated Nafion monomer/H3O+ systems were performed. H3O+ is considered rather than H+ assuming that H+ motion in water takes place either as H3O+ or as H3O+(H2O)n. Atomic charges for water, H3O+, and the Nafion side chain extracted from our ab initio results in Figure 6 are included in Table 7. The calculated diffusion coefficients are collected in Table 8, along with the MD simulation conditions, dimensions of the simulation box, and simulation times. Because of the limited size of the simulation box, the diffusion coefficients reported hereafter are local, i.e., they reflect the influence of the neighboring functional groups. The calculated H3O+ diffusion coefficient in water (2 × 10-9 m2/s) is slightly lower than that for the self-diffusion of pure water, 3.35 × 10-9 m2/s (Table 5). Note that we report the diffusion coefficient of H3O+, instead of that of H+. For that reason, the value is comparable to those of a larger ion such as Rb+,65

Figure 7. Snapshot of simulation 1 (Table 8) illustrating a H3O+ solvation configuration in the liquid phase, where the background water molecules have been removed for clarity.

whereas the experimental value15 reported for H+ diffusion in water is 9.31 × 10-9. Figure 7, a snapshot of simulation 1 (Table 8), illustrates the dynamic behavior of the first shell of water molecules surrounding the hydronium ion. Compared to the results shown in Figure 1e, the distances from the H atoms in H3O+ to the water oxygen atoms have been stretched by more than 0.3 Å. On the other hand, the three-water structure of Figure 1d and e is distorted in Figure 7, indicating an exchange of the water molecules located at the left of the figure, probably because of a thermal hydrogenbond breaking in the second solvation shell, as found recently by ab initio MD simulations.48 According to our model, at 298.15 K, H3O+ diffuses 4-4.5 times faster in pure water than in the presence in a Nafion/water environment (Table 8). Although the diffusion coefficients in simulations 1-3 (Table 8) were calculated for hydronium ions not attached to the anion group, the electrostatic interaction between the membrane and the hydronium ion, and the possible structure change of the solvation waters because of the presence of the ionic chain, might account for the smaller diffusion coefficient in the Nafion membrane. The proton diffusion values in simulations 1 and 3 differ very slightly because of the negligible difference in hydration conditions for the Nafion chains between the two cases. Simulation 2 was conducted at a higher temperature, 353.15 K, an upper limit for fuel cell operation. At this temperature, the model predicts a 4-fold increase in the H3O+ diffusion coefficient with respect to that at room temperature (Table 8). Experimentally, Kreuer4 has reported a 2-fold increase in the proton diffusion coefficients in highly hydrated membranes (Table 8). Thus, our simulation results agree qualitatively with the proton-transport behavior observed experimentally.

Ind. Eng. Chem. Res., Vol. 40, No. 22, 2001 4797 Table 8. Calculated and Experimental H3O+ Diffusion Coefficients in Bulk Water and Nafion/Water Systems system

simulation no.

unit cell (Å)

temperature (K)

total simulation time (ps)

diffusion coefficient (10-9 m2/s)

exptl value (10-9 m2/s)

500 water molecules

0

24.64

298.15

100

2.00

Nafion (1 monomer)/ 500 water molecules

1 2

24.77 24.77

298.15 353.15

250 300

0.45 1.78

1.40 at 298 K2 0.65 at 298 K4

Nafion (2 monomers)/ 500 water molecules

3

25.74

298.15

250

0.50

1.2 at 353 K4

Changes of the side-chain conformation in the hydrated state of Nafion membrane with respect to the gas-phase structure were investigated by analyzing the distance between the first backbone carbon atom (C22 in Figure 5) and the sulfur atom (S18) in the pendant chain. The variation of this distance during the simulation is displayed in Figure 8. Two measurements were evaluated under the conditions of simulations 1 and 3 indicated in Table 8. The calculated gas-phase C22S18 distance (Figure 6) is 4.33 Å, whereas in simulation 1, the average distance is 7.50 Å. In simulation 2, the distance increases to 8.15 Å. Thus, the side chain stretches dramatically because of the interactions of the -SO3- group with solvation water molecules and H3O+. It can be concluded, in agreement with previous reports, that the unfolding of side chain is favored in the aqueous environment.55 For the two-unit Nafion oligomer, the distance between the two neighboring -SO3- groups during the simulation is around 13 Å. Consequently, the side chains can be considered independent of each other. Aldebert et al.66 reported an estimated experimental distance of 18 Å between two -SO3- groups in a Nafion membrane with EW ) 1100. Thus, it can be inferred by comparison with our simulation results that the polymeric backbone between the -SO3- groups is not very rigid. Figure 9a and b illustrates the radial distribution functions gSF(r) calculated at 298.15 and 353.15 K, respectively. These functions provide the position of the F atoms surrounding the sulfonic group, therefore indicating the structural changes of the side chain as a function of temperature. Comparing Figure 9a and b, we observe that the positions of the first two peaks, corresponding to the first two sets of F atoms in Figure 6, remain unchanged from 298 to 353.15 K. The

next two peaks show only slight changes, but the last peak is shifted to higher separations, because of the high-temperature-induced chain elongation illustrated in Figure 8. At the beginning of simulation 3 (two-unit oligomer in water at 298.15 K), one H3O+ was deliberately located close to one of the -SO3- groups, to study the coordination of H3O+ with the sulfonic acid anion. Figure 10 illustrates that H3O+ remained coordinated to the -SO3- group during the 250-ps simulation. On average, H3O+ is 3.75 Å away from the -SO3- group, and this distance shows only very slight fluctuations. Figure 11 illustrates the structure of the solvated hydronium ion in contact with the sulfonic group, where all of the background water molecules and part of the polymer

Figure 8. Length of the side chain measured as the distance between C22 and S18 in Figure 6 during MD simulation; simulation states specified in Table 8.

Figure 9. (a) Radial distribution function gSF(r) illustrating the F atoms surrounding a central S atom in the Nafion chain. (a) At 298.15 K, (b) at 353.15 K.

4798

Ind. Eng. Chem. Res., Vol. 40, No. 22, 2001

Figure 10. Time evolution of the distance between the -SO3group and H3O+ from simulation 3 (Table 8).

and the energy required for additional water molecules attached by linear hydrogen bonds reduced as n increases. In the gas phase, when the O-O distance is less than 2.7 Å, the proton-transfer barrier is lower than the zeropoint energy for vibration of the OH bond, and the transfer can take place readily. In the presence of an external electric field, the symmetry of the double-wellshaped PES curve is destroyed. When the applied field is opposite to the dipole direction of the complex, the activation barrier diminishes dramatically. Calculated H3O+ diffusion coefficients in bulk water are found to be 4-4.5 times faster than those in model Nafion/water systems. At 353.15 K, the diffusion coefficient in model Nafion/water system increases 4 times with respect to that at ambient temperature conditions. The Nafion sulfonic acid-terminated side chain is folded in the gas phase, whereas in aqueous environment, MD results reveal that much more stretched side chains are favored and that the stretching increases when the temperature increases from 298.15 to 353.15 K. DFT calculations indicate that proton transfer from the acid group is favored and that an ion complex -SO3-‚‚‚H3O+ is formed, whereas results from MD simulations show that this contact ion pair formed is very stable and that the proton-transfer reaction from the side chain to water is irreversible. Acknowledgment This work is supported by the National Science Foundation Career Award Grant CTS-9876065 and by the Army Research Office Grant DAAD19-00-1-0087. A.W. was an undergraduate participant, supported under Grant NSF/REU EEC-9732345. The use of computational facilities at the National Center for Supercomputing Applications (NCSA) at the University of Illinois, Urbana-Champaign, and at NERSC is gratefully acknowledged. Literature Cited

Figure 11. Snapshot from simulation 3 illustrating the structure and solvation of the ion pair SO3-H3O+.

structure have been removed for clarity. Although the relevant distances between the sulfonic group and the hydronium ion are elongated with respect to those in Figure 5b, the ion pair remains in close contact. It can be concluded that the contact ion pair formed between the cation and the acid anion is very stable and that the protonation of the solvent by H-Nafion is preferred and irreversible, in agreement with experimental studies.54 4. Conclusions Stable gas-phase proton hydration complexes H+(H2O)n with n ) 1-5 are formed by exothermic reactions. The major contribution to the enthalpy and free energy of formation comes from the formation of H3O+,

(1) Heitner-Wirguin, C. Recent Advances in Perfluorinated Ionomer Membranes: Structure, Properties and Applications. J. Membr. Sci. 1996, 120, 1. (2) Verbrugge, M. W.; Hill, R. F. Analysis of Promising Perfluorosulfonic Acid Membranes for Fuel-cell Electrolytes. J. Electrochem. Soc. 1990, 137, 3770. (3) Schlick, S. Ionomers: Characterization, Theory and Applications; CRC Press: Boca Raton, FL, 1996. (4) Kreuer, K. D. On the Development of Proton Conducting Materials for Technological Applications. Solid State Ionics 1997, 97, 1. (5) James, P. J.; Elliott, J. A.; McMaster, T. J.; Newton, J. M.; Elliott, A. M. S.; Hanna, S.; Miles, M. J. Hydration of Nafion Studied by AFM and X-ray Scattering. J. Mater. Sci. 2000, 35, 5111. (6) Porat, Z.; Fryer, J. R.; Huxham, M.; Rubinstein, I. Electron Microscopy Investigation of the Microstructure of Nafion Films. J. Phys. Chem. 1995, 99, 4667. (7) Gierke, T. D.; Hsu, W. J. The Cluster-Network Model of Ion Clustering in Perfluorosulfonated Membranes; ACS Symposium Series No. 180; American Chemical Society: Washington, D.C., 1982. (8) Gierke, T. D.; Mumn, G. E.; Wilson, F. C. Morphology of Perfluorosulfonated Membrane Products. J. Polym. Sci. 1981, 19, 1687. (9) Hsu, W. J.; Barkley, J. R.; Meakin, P. Ion Percolation and Insulator-to-Conductor Transition in Nafion Perfluorosulfonic Acid Membranes. Macromolecules 1980, 13, 198.

Ind. Eng. Chem. Res., Vol. 40, No. 22, 2001 4799 (10) Yeo, S. C.; Eisenberg, A. Physical Properties and Supermolecular Structure of Perfluorinated Ion-Containing (Nafion) Polymers. J. Appl. Polym. Sci. 1977, 21, 875. (11) Rebrov, A. V.; Ozerin, A. N.; Svergun, D. I.; Bobrva, L. P.; Bakeev, N. F. Study of Aggregation of Macromolecules of Perfluorosulfinated Ionomer in Solution by the Small-Angle X-rayScattering Methodology. Polym. Sci. U.S.S.R. 1990, 32, 1593. (12) Yeager, H. L.; Steck, A. Cation and Water Diffusion in Nafion Ion Exchange Membranes: Influence of Polymer Structure. J. Electrochem. Soc. 1981, 128, 1880. (13) Falk, M. An Infrared Study of Water in Perfluorosulfonate (Nafion) Membranes. Can. J. Chem. 1980, 58, 1495. (14) Vishnyakov, A.; Neimark, A. Molecular Simulation Study of Nafion Solvation in Water and Methanol. J. Phys. Chem. B 2000, 104, 4471. (15) Atkins, P. W. Physical Chemistry, 5th ed.; Oxford University Press: Oxford, U.K., 1994. (16) Conway, B. E.; Bockris, J. O. M.; Linton, H. Proton Conductance and the Existence of the H3O‚ ion. J. Phys. Chem. 1956, 24, 834. (17) Eikerling, M.; Kornyshev, A. A.; Stimming, U. Electrophysical Properties of Polymer Electrolyte Membranes: A Random Network Model. J. Phys. Chem. B 1997, 101, 10807. (18) Ohwaki, T.; Yamashita, K. Electric Field Effects on Proton Transfer Reactions at Metal Electrodes: A DFT Study on H+(H2O)2/Pt(111) and Ag(111). Electrochemistry 1999, 67, 1214. (19) Schmidt, R. G.; Brickmann, J. Molecular Dynamics Simulations Study of a Hydronium Ion in Liquid Water with Implementation of the Proton Transfer by Means of a Hopping Mechanism. Solid Ionics 1995, 77, 3. (20) Grotthuss, C. J. T. d. Sur la decomposition de l’eau et des corps qu’elle tient en dissolution a’l’aide de l’electricite galvanique. Ann. Chim. 1806, LVIII, 54. (21) Lim, C.; Bashford, D.; Karplus, M. Absolute pKa Calculations with Continuum Dilelectric Methods. J. Phys. Chem. 1991, 95, 5610. (22) Meot-Ner, M.; Speller, C. V. Filling of Solvent Shells About Ions. 1. Thermochemical Criteria and the Effects of Isometric Clusters. J. Phys. Chem. 1986, 90, 6616. (23) Tunon, I.; Silla, E.; Bertran, J. Proton Solvation in Liquid Water. An Ab Initio Study Using Continuum Model. J. Phys. Chem. 1993, 97, 5547. (24) Lau, Y. K.; Ikuta, S.; Kebarle, P. Thermodynamics and Kinetics of the Gas-phase Reactions: H3O+(H2O)n-1 + H2O ) H3O+(H2O)n. J. Am. Chem. Soc. 1982, 104, 1462. (25) Lee, J.; Robinson, G. W.; Webb, S. P.; Philips, L. A.; Clark, J. H. Hydration Dynamics of Protons from Photon Initiated Acids. J. Am. Chem. Soc. 1986, 108, 6538. (26) Robinson, G. W. Proton Charge Transfer Involving the Water Solvent. J. Phys. Chem. 1991, 95, 10386. (27) Triolo, R.; Narten, A. H. Diffraction Pattern and Structure of Aqueous Hydrochloric Acid Solutions at 20 °C. J. Chem. Phys. 1975, 63, 3624. (28) Newman, J. S. Electrochemical Systems; Prentice Hall: Englewood Cliffs, NJ, 1973. (29) Allen, M. P.; Tildesley, D. J. Computer Simulation of Liquids; Oxford University Press: Oxford, U.K., 1990. (30) Din, X.-D.; Michaelides, E. E. Transport Processes of Water and Protons Through Micropores. AIChE J. 1998, 44, 35. (31) Paddison, S. J.; Paul, R.; Zawodzinski, T. A. A Statistical Mechanical Model of Proton and Water Transport in a Proton Exchange Membrane. J. Electrochem. Soc. 2000, 147, 617. (32) Eikerling, M.; Kornyshev, A. A.; Kuznetsov, A. M.; Ulstrup, J.; Walbran, S. Mechanisms of Proton Conductance in Polymer Electrolyte Membranes. J. Phys. Chem. B 2001, 105, 3646. (33) Frisch, M. J.; Trucks, G. W.; Schlegel, H. B.; Gill, P. M. W.; Johnson, B. G.; Robb, M. A.; Cheeseman, J. R.; Keith, T.; Petersson, G. A.; Montgomery, J. A.; Raghavachari, K.; Al-Laham, M. A.; Zakrzewski, V. G.; Ortiz, J. V.; Foresman, J. B.; Peng, C. Y.; Ayala, P. Y.; Chen, W.; Wong, M. W.; Andres, J. L.; Replogle, E. S.; Gomperts, R.; Martin, R. L.; Fox, D. J.; Binkley, J. S.; Defrees, D. J.; Baker, J.; Stewart, J. P.; Head-Gordon, M.; Gonzalez, C.; Pople, J. A. Gaussian 94, revision E.1; Gaussian Inc.: Pittsburgh, PA, 1997. (34) Frisch, M. J.; Trucks, G. W.; Schlegel, H. B.; Scuseria, G. E.; Robb, M. A.; Cheeseman, J. R.; Zakrzewski, V. G.; Montgomery, J. A.; Stratmann, R. E.; Burant, J. C.; Dapprich, S.; Millam, J. M.; Daniels, A. D.; Kudin, K. N.; Strain, O. F. M. C.; Tomasi, J.; Barone, B.; Cossi, M.; Cammi, R.; Mennucci, B.; Pomelli, C.;

Adamo, C.; Clifford, S.; Ochterski, J.; Petersson, G. A.; Ayala, P. Y.; Cui, Q.; Morokuma, K.; Malick, D. K.; Rabuck, A. D.; Raghavachari, K.; Foresman, J. B.; Ciolovski, J.; Ortiz, J. V.; Stefanov, V. V.; Liu, G.; Liashensko, A.; Piskorz, P.; Komaroni, I.; Gomperts, R.; Martin, R. L.; Fox, D. J.; Keith, T.; Al-Laham, M. A.; Peng, C. Y.; Nanayakkara, A.; Gonzalez, C.; Challacombe, M.; Gill, P. M. W.; Johnson, B.; Chen, W.; Wong, M. W.; Andres, J. L.; HeadGordon, M.; Replogle, E. S.; Pople, J. A. Gaussian 98, revision 0.1; Gaussian Inc.: Pittsburgh, PA, 1998. (35) Perdew, J. P.; Wang, Y. Accurate and Simple Density Functional for the Electronic Exchange Energy Generalized Gradient Approximation. Phys. Rev. B 1986, 33, 8800. (36) Becke, A. D. Density functional thermochemsitry II. The effect of the Perdew-Wang generalized gradient correlation correction. J. Chem. Phys. 1992, 97, 9173. (37) Becke, A. D. Density-functional thermochemistry. III. The role of exact exchange. J. Chem. Phys. 1993, 98, 5648. (38) Cerius2 User’s Manual; Molecular Simulations Inc.: San Diego, CA, 1997. (39) Mayo, S. L.; Olafson, B. D.; Goddard, W. A. Dreiding: A generic force field for molecular simulations. J. Phys. Chem. 1990, 94, 8897. (40) Berendsen, H. J. C.; Postma, J. P. M.; Gunsteren, W. F. v.; DiNola, A.; Haak, J. R. Molecular Dynamics with Coupling to an External Bath. J. Chem. Phys. 1984, 81, 3684. (41) Berendsen, H. J. C.; Postma, J. P. M.; Gunsteren, W. F. v.; Hermans, J. Intermolecular Forces; Reidel: Dordrecht, The Netherlands, 1981. (42) Ennari, J.; Elomaa, M.; Sundholm, F. Modeling a Polyelectrolyte System in Water to Estimate the Ion Conductivity. Polymer 1999, 40, 5035. (43) Oja¨mae, L.; Shavitt, I.; Singer, S. J. Potential Energy Surfaces and Vibrational Spectra of H5O2+ and Larger Hydrated Proton Complexes. Int. J. Quantum Chem: Quantum Chem. Symp. 1995, 29, 657. (44) Oja¨mae, L.; Shavitt, I.; Singer, S. J. Potential Models for Simulations of the Solvated Proton in Water. J. Chem. Phys. 1998, 109, 5547. (45) Moller, C.; Plesset, M. S. Note on an Approximation Treatment for Many-Electron Systems. Phys. Rev. 1934, 46, 618. (46) Tomasi, J.; Persico, M. Molecular Interactions in Solution: An Overview of Methods Based on Continuous Distributions of the Solvent. Chem. Rev. 1994, 94, 2027. (47) Kalstro¨m, G. Proton Transport in Water Modeled by a Quantum Chemical Dielectric Cavity Model. J. Phys. Chem. 1988, 92, 1318. (48) Marx, D.; Tuckerman, M. E.; Hutter, J.; Parrinello, M. The Nature of the Hydrated Excess Proton in Water. Nature 1999, 397, 601. (49) Noyes, R. M. Thermodynamics of Ion Hydration as a Measure of Effective Dielectric Properties of Water. J. Am. Chem. Soc. 1962, 84, 513. (50) Meot-Ner, M. Comparative Stabilities of Cationic and Anionic Hydrogen-Bonded Networks. Mixed Clusters of WaterMethanol. J. Am. Chem. Soc. 1986, 108, 6189. (51) Jensen, F. Introduction to Computational Chemistry; Wiley: New York, 1999. (52) Boys, S. F.; Bernardi, F. The Calculation of Small Molecular Interaction by the Difference of Separate Total Energies. Some Procedures with Reduced Error. Mol. Phys. 1970, 19, 553. (53) Laporta, M.; Pegoraro, M.; Zanderighi, L. Perfluorosulfonated Membrane (Nafion): FT-IR Study of the State of Water with Increasing Humidity. Phys. Chem. Chem. Phys. 1999, 1, 4619. (54) Buzzoni, R.; Bordiga, S.; Ricchiardi, G.; Spoto, G.; Zecchina, A. Interaction of H2O, CH3OH, (CH3)2O, CH3CN, and Pyridine with the Superacid Perofluorosulfonic Membrane Nafion: An IR and Raman Study. J. Phys. Chem. 1995, 99, 11937. (55) Paddison, S. J.; Zawodzinski, T. A. Molecular Modeling of the Pendant Chain in Nafion. Solid State Ionics 1998, 113-115, 333. (56) Hehre, W. J.; Radom, L.; Schleyer, P. v. R.; Pople, J. A. Ab Initio Molecular Orbital Theory; John Wiley & Sons: New York, 1986. (57) Hsu, W. Y.; Gierke, T. D. Ion Transport and Clustering in Nafion Perfluorinated Membranes. J. Membr. Sci. 1983, 13, 307. (58) Handbook of Chemistry and Physics, 77th. ed.; Lide, D. R., Ed.; CRC Press: Boca Raton, FL, 1997.

4800

Ind. Eng. Chem. Res., Vol. 40, No. 22, 2001

(59) Fuller, T. F.; Newman, J. Experimental Determination of the Transport Number of Water in Nafion 117 Membrane. J. Electrochem. Soc. 1992, 139, 1332. (60) Mu¨ller-Plathe, F. An All-Atom Force Field for Liquid EthanolsProperties of Ethanol-Water Mixtures. Mol. Simul. 1996, 19, 133. (61) Colonomos, P.; Wolynes, P. G. Molecular Theory of Solvated Ion Dynamics. II. Fluid Structure and Ionic Mobilities. J. Chem. Phys. 1979, 71, 2644. (62) Koneshan, S.; Rasaiah, J. C.; Dang, L. X. Computer Simulation Studies of Aqueous Solutions at Ambient and Supercritical Conditions Using Effective Pair Potential and Polarizable Potential Models for Water. J. Chem. Phys. 2001, 114, 7544. (63) Dang, L. X. Importance of Polarization Effects in Modeling the Hydrogen Bond in Water Using Classical Molecular Dynamics Techniques. J. Phys. Chem. B 1998, 102, 620.

(64) Chialvo, A. A.; Cummings, P. T. Microstructure of Ambient and Supercritical Water. Direct Comparison between Simulation and Neutron Scattering Experiments. J. Phys. Chem. 1996, 100, 1309. (65) Oelkers, E. H.; Helgeson, H. C. Calculation of the Transport Properties of Aqueous Species at Pressures to 5 kB and Temperatures to 1000 °C. J. Solution Chem. 1989, 18, 601. (66) Aldebert, P.; Dreyfus, B.; Pineri, M. Small-Angle Neutron Scattering of Perfluorosulfonated Ionomers in Solution. Macromolecules 1986, 19, 9, 2651.

Received for review May 29, 2001 Revised manuscript received August 24, 2001 Accepted August 28, 2001 IE010467Y