J. Phys. Chem. B 2007, 111, 6357-6363
6357
Proton Hopping in Phosphoric Acid Solvated Nafion Membrane: A Molecular Simulation Study Liuming Yan,*,† Suhua Zhu,†,‡ Xiaobo Ji,†,‡ and Wencong Lu† Department of Chemistry, College of Sciences, and School of Material Science and Engineering, Shanghai UniVersity, 99 Shangda Road, Shanghai 200444, China ReceiVed: February 5, 2007; In Final Form: April 10, 2007
Molecular simulation studies of the microstructure and of the proton transport properties of phosphoric acid solvated Nafion membrane are carried out. The ab initio calculations show that the phosphoric acid is a good solvent to promote the proton ionization of the sulfonic acid group, and only two phosphoric acid molecules are necessary for the dissociation of one sulfonic acid group. A mechanism of proton hopping between phosphoric acid and protonated phosphoric acid cation in the hydrophilic subphase is also elucidated by ab initio calculations. The molecular dynamics simulations, conducted at a phosphoric acid concentration of 25.4% (wt) which is slightly lower than that of phosphoric acid swollen Nafion, show that the phosphoric acid exists in subphases and that it cannot develop into a continuous subphase. Thus, proton-hopping pathways are interrupted, and the conductivity is expected to be lower than that for pure phosphoric acid. The molecular dynamics simulations, conducted at a phosphoric acid concentration of 45.1% (wt) which corresponds to an unstable state, show that the hydrophobic poly(tetrafluoroethylene) backbones trend to gather together forming hydrophobic clusters and that the phosphoric acid forms a continuous subphase with the sulfonic acid groups located at the hydrophobic/hydrophilic interface. Thus, proton-hopping pathways can develop uninterruptedly like the pure phosphoric acid, and high conductivity is expected. The molecular dynamics study also shows that the hydrogen-bonding characteristics of phosphoric acid and sulfonate anion are similar regardless of the factor that the former can move freely while the latter is attached to Nafion backbone.
1. Introduction
TABLE 1: Nonbonded Force Parameters for the Phosphoric Acid and Its Protonated Cation
The high-temperature proton exchange membrane (PEM) for hydrogen fuel cell has attracted numerous research interests because of the fast oxygen electroreduction rate at cathode as well as of the high CO tolerance to the fed gas when the fuel cell operates at temperatures well above 100 °C.1-3 The most common form of membrane for the proton exchange is the Nafion membrane. Its proton conductivity performance is closely related to the degree of hydration of the membrane. In general, the operating temperature of Nafion membrane should be lower than the boiling point of water, because above this temperature, most of the water in the Nafion membrane will evaporate, and this lowers the proton conductivity. The relative low-operating temperature would lead to several significant drawbacks, including the lower oxygen electroreduction rate, the inactivation of precious metal catalysts, the “crossover” of significant amounts of water or methanol through permeation, and the electro-osmotic drag. To keep high conductivity at elevated temperatures, much effort has been focused on the following two main strategies. One is the increasing of the proton conductivity at low relative humidity and high temperature, which is achieved by the addition of a high boiling point liquid, such as the heterocyclic-structured (1-butyl, 3-methyl)imidazolium triflate4 or other ionic liquids.5-7 The other one is to form an organic-inorganic composite membrane with higher mechanical, thermal, and chemical stabilities by the addition * To whom correspondence should be addressed. Tel: 8621-66132405. E-mail:
[email protected]. † Department of Chemistry. ‡ School of Material Science and Engineering.
q/e P OM OA H
σ(Å)
(kcal‚mol-1)
H3PO4
P(OH)4+
4.595 3.111 2.955
0.0936 0.1491 0.2029
1.093 -0.601 -0.605 0.441
1.184 -0.567 0.521
of an inorganic ingredient, thus improving the water retention at elevated temperatures.8-12 As one of the proton conduction media with high efficiency, aqueous solution of phosphoric acid has been used as the sole proton conduction media in the phosphoric acid fuel cell operating at temperature between 150 and 200 °C for a long time.13 Phosphoric acid is an amphoteric molecular liquid with boiling point of 213 °C and high dielectric constant. The first attempt to dope Nafion with phosphoric acid was made by Savinell et al. and achieved a conductivity of 0.05 S/cm at 150 °C.14 Since phosphoric acid is a weak acid compared to sulfonic acid, it acts as a Brønsted base and solvates the proton from the strong sulfonic acid in the same way as water does.15,16 Phosphoric acid has also been doped to poly(benzimidazole) (PBI) and PBI composite membrane to improve the proton conductivity.17,18 The most magnificent advantage of the introduction of phosphoric acid is that it can conduct proton even in an anhydrous form18,19 because of its proton-solvating ability and self-ionization behavior.20 Meanwhile, immobilized phosphoric acids, such as zirconium phosphate or phosphonic acid, arouse tremendous research interest in the field of composite proton conduction membrane.5,11,12,17,18,21 Motivated
10.1021/jp071005m CCC: $37.00 © 2007 American Chemical Society Published on Web 05/23/2007
6358 J. Phys. Chem. B, Vol. 111, No. 23, 2007
Yan et al.
SCHEME 1: Molecular Model of the Phosphoric Acid Solvated Nafion Membrane; 1, the Oligomer Consisting of a Poly(tetrafluoroethylene) Backbone and Four Pendant Side Chains with Sulfonic Acid End Groups; 2, 3, 4, Atomic Site Name of the Side Chain, the Phosphoric Acid, and the Protonated Phosphoric Acid Cation
by the fascinating self-ionization behavior, phosphoric acid is expected to be a good candidate for the enhancement of proton conductivity of Nafion membrane on a high temperature and a low water content.22 Herein, we demonstrate a theoretical study on the Nafion membrane solvated by phosphoric acid. The proton affinity energy of a single phosphoric acid molecule as well as phosphoric acid clusters was calculated. Subsequently, a mechanism of proton hopping between a phosphoric acid molecule and a protonated phosphoric acid cation was proposed by ab initio calculations. Moreover, we also study the structural and transport properties of the phosphoric acid solvated Nafion membrane by the use of molecular dynamics simulations. These theoretical works would benefit the further development of the high-temperature proton conduction membrane with improved performance.
Figure 1. Optimized conformations for clusters of triflic acid with (a) one, (b) two, and (c) three phosphoric acid molecules. The first dissociation of the acidic proton occurs when two phosphoric acid molecules are added, and the protonated phosphoric acid cation has two hydrogen bonds with the sulfonate anion. In the cluster with three phosphoric acid molecules, the ionized proton is shared by two of the phosphoric acid molecules at distances of 1.113 and 1.276 Å, respectively. The third phosphoric acid molecule is close to the triflate anion at a distance of hydrogen bonding, thus, the hydrogen-bonding energy can compensate the energy needed in the ionization of the triflic acid. Color code: O, red; H, white; P, tan; F, green; C, cyan; and S, yellow.
2. The Calculation Methods and the Molecular Models 2.1. The ab Initio Calculation Methods. All the ab initio calculations are based on density functional theory at level of B3PW91/6-31G(d,p)23-28 as implemented in GAUSSIAN 03 program.29 All the structures are first optimized and then verified to be local minimum by the calculation of second derivatives. The calculated energies are calibrated by zero-point energies from the frequencies of vibration. Partial charges are calculated with the ChelpG method30 which reproduces the electrostatic potential of the isolated molecule or ion. 2.2. The Molecular Models and Force Field Parameters. The Nafion membrane is modeled by oligomers which consists of a poly(tetrafluoroethylene) backbone with a length of 64 carbons and 4 pendant side chains terminated by sulfonic acid groups (Scheme 1).31 All the sulfonic acid groups are supposed to ionize leaving a negatively charged SO3- group when solvated by phosphoric acid. Protonated phosphoric acid cations P(OH)4 + are added to ensure the neutrality of the membrane. In our molecular dynamics simulation, all-atom models are applied to all the ingredients composing the phosphoric acid solvated Nafion membrane. The total energy is a sum over all intramolecular and intermolecular interactions. The intramolecular interactions, or bonding interactions, consist of bond stretching ub ) ∑i1/2Kb(b - b0)2, bond angle bending
Figure 2. Proton-hopping mechanism: (a) phosphoric acid molecule P1 and (b) its protonated cation P2. (c) When P1 approaches P2 from the left, the proton is shared between P1 and P2 with equal distance of 1.190 Å, which is about 0.2 Å farther than the usual O-H bond distance in a phosphoric acid molecule but much nearer than that of the hydrogen bond distance (about 1.6 Å). (d) When a third phosphoric acid molecule P3 approaches P1 from the left, the shared proton between P1 and P2 is attracted by P3, and the distance between the proton and P2 increases to 1.413 Å. As a result, the proton hops from P2 to P1. All the systems are optimized to local minimum, and the energies noted below each picture are relative to isolated phosphoric acid. Color code: O, red; H, white; and P, tan.
uθ ) ∑i1/2Kθ(θ - θ0)2, and dihedral torsion uφ ) ∑i1/2Kφ(1 + cos 3φ). The intermolecular interactions, or nonbonding interactions, are the sum over all pairs of van der Waals interactions and electrostatic interactions. The partial charges for phosphoric acid and its protonated cation are derived from density functional theory calculation. The partial charges evaluated in this way are different for each oxygen atom or each hydrogen atom because of the lack of symmetry and the intramolecular hydrogen-oxygen interaction.
Proton Hopping in Phosphoric Acid Nafion Membrane
J. Phys. Chem. B, Vol. 111, No. 23, 2007 6359
TABLE 2: Bonded Force Field Parameters for Phosphoric Acid and Its Protonated Cation P(OH)4+
H3PO4 Bond Stretching (Units: b0 Å; Kb kcal ‚mol-1‚Å-2) b0 P-OM P-OA OA-H
1.471 1.610 0.966
Kb
b0
Kb
900 600 750
1.559 0.972
600 750
θ0
Kθ
109.5 119.0
110 95
Bond Angle Bending (Units: θ0 Degree; Kθ kcal‚mol-1‚rad-2) θ0 Kθ OM-P-ΟΑ OM-P-ΟΑ P-ΟΑ-Η OM-P-OA-H or OA-P-OA-H
115.8 102.5 111.4
95 110 95
Dihedral Torsion (Unit: uφ kcal‚mol-1) uφ ) 0.75(1 + cos 2φ) + 0.25(1 + cos 3φ)
TABLE 3: Model Molecules Used for the Evaluation of Force Field Parameters 5 6 7 8 9 10
CF3CF3 CF3CF2CF3 CF3CF2CF2CF3 CF3OCF3 CF3CF2OCF3 CF3CF2CF2OCF3
11 12 13 14 15 16
CF3OCF2CF2OCF3 CF3SO3CF3CF2SO3CF3OCF2CF2SO3CF3CF2CF2CF2CF2CF2CF3 CF3OCF2CF(CF3)OCF2CF2SO3-
TABLE 4: van der Waals Parameters and Partial Charges of the Nafion Oligomera C (side chain) F (side chain) S OS O2 C (backbone, terminal) F (backbone, terminal) C (backbone, inner) F (backbone, inner) a
σ (Å)
(kcal‚mol-1)
q/e
3.473 3.093 3.550 3.070 3.150 3.473 3.093 3.473 3.093
0.0951 0.0725 0.2499 0.1700 0.2000 0.0951 0.0725 0.0951 0.0725
0.298 -0.121 0.923 -0.305 -0.561 0.390 -0.130 0.180 -0.090
Vishnyakov and Neimark.
In our simulations, the calculated partial charges are averaged for each group of atoms (Table 1). These values are slightly smaller compared to that used in the literature.32 The van der Waals parameters are adopted from Spieser and Leeflang,32 except for that of hydrogen atom which is set to zero for both attraction and repulsion interactions instead. We simulate the pristine phosphoric acid using these parameters, and we get a liquid density of 1.844 g/cm3 at 298 K and 1 atm which is only slightly smaller than the experimental value of 1.868 g/cm3 at the same temperature and pressure.33 If the Spieser parameters32 are used, the simulated density is 1.807 g/cm3. The equilibrium bond lengths b0 and bond angles θ0 are derived from density functional theory calculation and are averaged for each set of parameters (Table 2). All the bonded force constants, including bond stretching, bond angle bending, and dihedral torsion, are adopted from the literature.32 United atom models are used by most of the simulations of Nafion in literature. The force field parameters are usually derived from small molecules of perfluoroethers,34,35 hydrofluorocarbons,36 or even generic parameters such as the Cornell,37 the DREIDING,38 Amber 99,39 and so forth. Although a united atom force field model predicts accurately the thermodynamic and structural properties and the computational cost is only about one-ninth compared to the computational cost that an all-atom model is used.35 However, a united atom model might not be expected to be quite so successful for predicting transport properties.35 Therefore, an all-atom model is used for the Nafion oligomer.
TABLE 5: Bonded Force Parameters for the Nafion Oligomer Bond Stretching (Units: b0 Å; Kb kcal‚mol-1‚Å-2) C-C C-F C-OS
b0
Kb
1.558 1.342 1.406
624 976 851
C-S S-O2
b0
Kb
1.878 1.475
361 1302
Bond Angle Bending (Units: θ0 degree; Kθ kcal‚mol-1‚rad-2) C-C-C F-C-C F-C-F C-C-OS F-C-OS
θ0
Kθ
116.7 109.3 107.7 109.1 110.6
233 172 206 211 285
C-OS-C S-C-C S-C-F C-S-O2 O2-S-O2
θ0
Kθ
121.7 117.1 109.7 102.3 115.6
170 150 194 177 319
Dihedral Torsion (Unit: Kφ kcal‚mol-1) Kφ C-C-C-Ca OS-C-C-C OS-C-C-OS C-OS-C-Cb S-C-C-OS O2-S-C-C
0.4809 0.3649 0.5476 0.0000
Kφ O2-S-C-F F-C-C-F F-C-C-C OS-C-C-F C-OS-C-F S-C-C-F
0.3579 0.3245 0.2378 0.4535 0.2594 0.3185
a uφ ) 0.458 + 3.192 cos φ - 0.649 cos2 φ - 7.462 cos3 φ + 2.134 cos4 φ + 6.292 cos5 φ. b uφ ) 1.962 + 4.967 cos φ - 1.481 cos2 φ - 9.141 cos3 φ + 1.992 cos4 φ + 6.756 cos5 φ.
We derive the force field parameters from a set of model molecules (Table 3). The equilibrium bond lengths and bond angles for the oligomer backbone are derived from molecule 15, and those for the pendant are from 16. Then, each bond length or bond angle is perturbed from its equilibrium value (not the averaged value, the geometrically optimized value for each molecule listed in Table 3) one by one while keeping the rest of the geometry fixed. The energy differences between the perturbed and equilibrium geometries are used to evaluate the bond stretching or bond angle bending constant. The dihedral torsion parameters are evaluated by partial optimization of the molecules 5-14. During the partial optimization, the target dihedral is frozen at a predetermined value and the rest of the freedoms are allowed to relax since the torsional motions have much lower frequencies compared to bond stretching and bond angle bending, thus, there is enough time for their relaxation to take place. The torsional potential profile as a function of dihedral is calculated as the energy differences between these partially optimized geometries and the fully optimized geometry and is fitted to the empirical formula of torsional potential. The partial charges of the side chains and the backbone are derived from molecules 15 and 16, respectively. The van der Waals parameters are adopted from Vishnyakov and Neimark.40,41
6360 J. Phys. Chem. B, Vol. 111, No. 23, 2007
Yan et al.
Figure 3. Snapshots of the Nafion oligomers solvated by the phosphoric acid: one of the eight oligomers is shown in a ball-packing model and the other seven oligomers are shown in a stick-and-ball model; the phosphoric acid molecules and the protonated cations are shown in licorice model. (a) At phosphoric acid concentration of 45.1%, the phosphoric acid molecules gather together and form continuous subphase, in which proton-hopping pathways develop and proton transport is enhanced. (b) At concentration of 25.4%, no continuous phosphoric acid subphase forms and only discontinuous clusters exist, and thus, proton-hopping pathways are interrupted. At both cases, the hydrophobic poly(tetrafluoroethylene) backbones fold together and form hydrophobic clusters, while the hydrophilic sulfonate anions stretch out to the interface between the hydrophobic cluster and the hydrophilic subphase or clusters. Color code: C, cyan; F, green; O, red; S, yellow; P, tan; and H, white.
2.3. Simulation Details. The molecular dynamics simulations are carried out by use of the DL_POLY program.42 LennardJones interactions are cut off at 8.0 Å. The Ewald summation algorithm is applied to the calculation of electrostatic potential energy in account for long-range interactions. The temperature and pressure are maintained by the Berendsen thermostat.43 The equations of motion are solved by the Verlet scheme with time step of 1 fs.44 In the simulation, eight Nafion oligomer molecules are placed in a rectangular simulation cell as well as 32 protonated phosphoric acid cations to compensate the charges of the oligomers. In addition, 277 or 96 phosphoric acid molecules are inserted into the simulation cell, corresponding to a phosphoric acid concentration of 45.1% and 25.4%. That is to say, our simulated systems have a mole ratio of phosphoric acid to sulfonic acid of about 9.6 to 1 and of about 4 to 1. The simulated system with high phosphoric acid concentration does represent an unstable system, since the phosphoric acid equilibrated Nafion 117 contains about 5 mol of phosphoric acid per mole of sulfonic acid in the polymer.14 However, our purpose to exploit this unstable system is to explain (see Conclusions) why the conductivity for phosphoric acid swollen Nafion is lower than that of pure phosphoric acid.14 The density of the initial configuration is about 0.5 g/cm3. During the simulation, the simulation cell is gradually squeezed; the equilibrium densities are about 1.820 and 1.887 g/cm3 for these systems, respectively. Since the dry density of Nafion is about 2.0 g/cm3,45 and that of phosphoric acid is 1.868 g/cm3, the theoretical densities of these two systems are 1.923 and 1.956 g/cm3. The simulated densities are in good agreement with theoretical densities. Two hundred thousand time steps are run in the first period to ensure that the system reaches equilibration. After equilibration, 40 000 more time steps are carried out, and a trajectory file is recorded every time step
with 40 000 configurations in total. This trajectory file is used to generate the structural information as well as transport properties. 3. Results and Discussion 3.1. Energetic Consideration. An energetic consideration is carried out by ab initio calculation. The proton affinity energy of a single phosphoric acid is 8.87 eV, which is about 1.4 eV larger than the proton affinity energy of water (7.49 eV) calculated at the same level of theory.31 This proton affinity energy cannot compensate the total energy needed for the proton ionization of the sulfonic acid; for example, the proton ionization energy of a triflic acid is 13.79 eV and similar values are reported for other types of sulfonic acid.31,46 However, the protonation energy of phosphoric acid cluster composed of two phosphoric acid molecules can compensate the proton ionization energy (Figure 1) compared to the factor that three water molecules are needed to compensate this proton ionization energy. For example, the calculated proton affinity energies are 10.43 and 11.43 eV for clusters with two and three phosphoric acid molecules. Comparing these with that of water clusters, the proton affinity energies are 9.23 and 12.03 eV for clusters of two and three water molecules.31 By the consideration that the charged species are only separated by a finite distance, the Coulomb attraction between a protonated phosphoric acid cluster and a sulfonate anion will compensate another 1 or 2 eV. In addition, there is more hydrogen bonding between a phosphoric acid molecule and a sulfonate anion than there is between a water molecule and a sulfonate anion. Therefore, proton ionization can spontaneously occur with the existence of two phosphoric acid molecules for each sulfonic acid. From this consideration, it can be deduced that phosphoric acid and water share some common properties of proton solvating; both of them can promote the proton ionization of the sulfonic acid side chains of Nafion. Figure 1 shows the optimized conformations of
Proton Hopping in Phosphoric Acid Nafion Membrane clusters composed of one triflic acid molecule and one to three phosphoric acid molecules. It is clearly shown that the triflic acid ionized in case that two phosphoric acid molecules exist in its vicinity. 3.2. Proton-Hopping Mechanism. It can be expected that both diffusion and hopping contribute to the proton transport in the phosphoric acid solvated Nafion membrane. Since a nondissociated harmonic bond model is used for the O-H bond in phosphoric acid and its protonated cation, the hopping mechanism, similar to a Grotthuss mechanism in water,47 is not allowed in our molecular dynamics simulation. However, the possibility of a hopping mechanism can be verified by ab initio calculations. From Figure 2c, it is found that a proton is shared by two phosphoric acid molecules, P1 and P2, at equal distance of 1.190 Å in a cluster composed of two acid molecules. This distance is much longer than the usual O-H bond distance in phosphoric acid molecules; thus, the O-H bonds are weakened because of the existence of a neighboring acid molecule. When a third phosphoric acid molecule, P3, approaches P1 from the left (Figure 2d), the shared proton between P1 and P2 is attracted by P3 and the H-O distance in P2 increases to 1.413 Å. Therefore, the proton hops from P2 to P1. If this process continues, a proton-hopping pathway develops and the proton transports from right to left without any transport of phosphoric acid molecules in the system. 3.3. Microstructure of the Phosphoric Acid Solvated Nafion Membrane. The characteristics that the backbone is extremely hydrophobic and that the side chains are extremely hydrophilic shade great influences on the microstructure of the Nafion membrane. At phosphoric acid concentration of 45.1%, the phosphoric acid molecules trend to gather together and form continuous hydrophilic subphase. This is a decisive structural factor that governs the proton transport characteristics since uninterrupted proton-hopping pathways could develop in the hydrophilic subphase. At phosphoric acid concentration of 25.4%, the phosphoric acid molecules also trend to gather together forming hydrophilic clusters; however, no continuous subphase forms; therefore, proton-hopping pathways are interrupted. On the other hand, the oligomer backbones trend to fold and form hydrophobic clusters with the hydrophilic side chains stretching out on the interface between the hydrophobic cluster and hydrophilic subphase or cluster (Figure 3). In addition, a connectivity analysis is carried out. If two phosphoric acid molecules are at a distance of less than 5.0 Å (P-P distance), a connection between these two molecules is drawn. If two phosphoric acid molecules can be connected by a series connection, they are supposed to be in the same connective set. For the system with phosphoric acid concentration of 45.1%, the connectivity analysis reveals that almost all the phosphoric acid molecules including the protonated cations belong to the same connective set, and only a few (140 °C) proton exchange membrane fuel cells. J. Power Sources 2005, 142 (1-2), 223-237. (22) Li, Q. F.; He, R. H.; Jensen, J. O.; Bjerrum, N. J. Approaches and recent development of polymer electrolyte membranes for fuel cells operating above 100 °C. Chem. Mater. 2003, 15 (26), 4896-4915. (23) Perdew, J. P.; Wang, Y. Accurate and simple analytic representation of the electron-gas correlation energy. Phys. ReV. B 1992, 45 (23), 1324413249. (24) Perdew, J. P. Unified Theory of Exchange and Correlation beyond the Local Density Approximation. In Electronic Structure of Solids; Ziesche, P., Eschrig, H., Eds.; Akademie Verlag: Berlin, 1991; pp 11-20. (25) Perdew, J. P.; Chevary, J. A.; Vosko, S. H.; Jackson, K. A.; Pederson, M. R.; Singh, D. J.; Fiolhais, C. Atoms, molecules, solids, and surfaces: Applications of the generalized gradient approximation for exchange and correlation. Phys. ReV. B 1992, 46 (11), 6671-6687. (26) Becke, A. D. Density-functional thermochemistry II. The effect of the Perdew-Wang generalized-gradient correlation correction. J. Chem. Phys. 1992, 97 (12), 9173-9177. (27) Hehre, W. J.; Ditchfield, R.; Pople, J. A. Self-consistent molecular orbital methods. XII. Further extensions of Gaussian-type basis sets for use in molecular orbital studies of organic molecules. J. Chem. Phys. 1972, 56 (5), 2257-2261. (28) Gordon, M. S. The isomers of silacyclopropane. Chem. Phys. Lett. 1980, 76 (1), 163-168. (29) Frisch, M. J.; Trucks, G. W.; Schlegel, H. B.; Scuseria, G. E.; Robb, M. A.; Cheeseman, J. R.; Montgomery, J. A.; Vreven, T., Jr.; Kudin, K. N.; Burant, J. C.; Millam, J. M.; Iyengar, S. S.; Tomasi, J.; Barone, V.; Mennucci, B.; Cossi, M.; Scalmani, G.; Rega, N.; Petersson, G. A.; Nakatsuji, H.; Hada, M.; Ehara, M.; Toyota, K.; Fukuda, R.; Hasegawa, J.; Ishida, M.; Nakajima, T.; Honda, Y.; Kitao, O.; Nakai, H.; Klene, M.; Li, X.; Knox, J. E.; Hratchian, H. P.; Cross, J. B.; Adamo, C.; Jaramillo, J.; Gomperts, R.; Stratmann, R. E.; Yazyev, O.; Austin, A. J.; Cammi, R.; Pomelli, C.; Ochterski, J. W.; Ayala, P. Y.; Morokuma, K.; Voth, G. A.; Salvador, P.; Dannenberg, J. J.; Zakrzewski, V. G.; Dapprich, S.; Daniels, A. D.; Strain, M. C.; Farkas, O.; Malick, D. K.; Rabuck, A. D.; Raghavachari, K.; Foresman, J. B.; Ortiz, J. V.; Cui, Q.; Baboul, A. G.; Clifford, S.; Cioslowski, J.; Stefanov, B. B.; Liu, G.; Liashenko, A.; Piskorz, P.; Komaromi, I.; Martin, R. L.; Fox, D. J.; Keith, T.; Al-Laham, M. A.; Peng, C. Y.; Nanayakkara, A.; Challacombe, M.; Gill, P. M. W.; Johnson, B.; Chen, W.; Wong, M. W.; Gonzalez, C.; Pople, J. A. Gaussian-2003, Revision C.2; Gaussian, Inc.: Pittsburgh, PA, 2003. (30) Breneman, C. M.; Wiberg, K. B. Determining atom-centered monopoles from molecular electrostatic potentials. The need for high sampling density in formamide conformational analysis. J. Comp. Chem. 1990, 11 (3), 361-373. (31) Yan, L.; Balbuena, P. B.; Seminario, J. M. Perfluorobutane sulfonic acid hydration: A model for a proton conductive polymer fragment. J. Phys. Chem. A 2006, 110 (13), 4574-4581. (32) Spieser, S. A. H.; Leeflang, B. R. A force field for phosphoric acid: Comparison of simulated with experimental data in the solid and liquid state. J. Phys. Chem. A 2000, 104 (31), 7333-7338. (33) Dean, J. A. Lange’s handbook of chemistry; McGraw-Hill: New York, 1999; Vol. 3.
Proton Hopping in Phosphoric Acid Nafion Membrane (34) Li, H.-C.; McCabe, C.; Cui, S. T.; Cummings, P. T.; Cochran, H. D. Development of a Force Field for Molecular Simulation of the Phase Equilibria of Perfluoroethers. Mol. Phys. 2003, 101 (14), 2157-2169. (35) Li, H.-C.; McCabe, C.; Cui, S. T.; Cummings, P. T.; Cochran, H. D. Development of a force field for molecular simulation of the phase equilibria of perfluoromethylpropyl ether. Mol. Phys. 2002, 100 (2), 265272. (36) Fermeglia, M.; Ferrone, M.; Pricl, S., Development of an all-atoms force field from ab initio calculations for alternative refrigerants. Fluid Phase Equilib. 2003, 210, 105-116. (37) Cornell, W. D.; Cieplak, P.; Bayly, C. I.; Gould, I. R.; Merz, K. M.; Ferguson, D. M.; Spellmeyer, D. C.; Fox, T.; Caldwell, J. W.; Kollman, P. A. A second generation force field for the simulation of proteins, nucleic acids, and organic molecules. J. Am. Chem. Soc. 1995, 117 (19), 51795197. (38) Mayo, S. L.; Olafson, B. D.; Goddard, W. A., III. DREIDING: A generic force field for molecular simulations. J. Phys. Chem. 1990, 94 (26), 8897-8909. (39) The amber force field file. http://amber.scripps.edu/ (accessed Dec 2006). (40) Vishnyakov, A.; Neimark, A. V. Molecular dynamics simulation of microstructure and molecular mobilities in swollen Nafion membranes. J. Phys. Chem. B 2001, 105 (39), 9586-9594.
J. Phys. Chem. B, Vol. 111, No. 23, 2007 6363 (41) Vishnyakov, A.; Neimark, A. V. Molecular dynamics simulation of Nafion oligomer solvation in equimolar methanol-water mixture. J. Phys. Chem. B 2001, 105 (32), 7830-7834. (42) Smith, W.; Leslie, M.; Forester, T. R. Computer code DL_POLY_2.14; CCLRC, Daresbury Laboratory: Daresbury, U.K., 2003. (43) Berendsen, H. J. C.; Postma, J. P. M.; van Gunsteren, W. F.; DiNola, A.; Haak, J. R. Molecular dynamics with coupling to an external bath. J. Chem. Phys. 1984, 81 (8), 3684-3690. (44) Verlet, L. Computer “experiments” on classical fluids. I. Thermodynamical properties of Lennard-Jones molecules. Phys. ReV. 1967, 159 (1), 98-103. (45) Nandan, D.; Mohan, H.; Iyer, R. M. Methanol and water uptake, densities, equivalental volumes and thicknesses of several uni- and divalent ionic perfluorosulphonate exchange membranes (Nafion-117) and their methanol-water fractionation behaviour at 298 K. J. Membr. Sci. 1992, 71 (1-2), 69-80. (46) Paddison, S. J. Proton conduction mechanisms at low degrees of hydration in sulfonic acid-based polymer electrolyte membranes. Annu. ReV. Mater. Res. 2003, 33 (1), 289-319. (47) Agmon, N. The Grotthuss mechanism. Chem. Phys. Lett. 1995, 244 (5-6), 456-462. (48) Tsuchida, E. Ab initio molecular-dynamics simulation of concentrated phosphoric acid. J. Phys. Soc. Jpn. 2006, 75 (5), 054801-054805.
