Molecular Dynamics Simulations of Substrate Hydrophilicity and Confinement Effects in Capped Nafion Films.

Nafion nanocomposites for energy-related applications are being used extensively because of the attractive properties such as enhanced water retention, low unwanted crossover of electrolytes, and high proton conductivity. We present the results of the molecular dynamics modeling of Nafion films confined between two walls (substrates) of different polymer-wall interaction strengths and of different separation distances to model Nafion nanocomposites. Our goal is to provide insights into the effects of varying hydrophilicity and volume fraction of fillers/nanoparticles on the internal structure and water transport inside the Nafion membrane. The sulfur-sulfur radial distribution function first peak distance and the sulfur-oxygen (water) coordination number in the first hydration shell were negligibly affected by the wall (substrate) hydrophilicity or the film thickness. The Nafion side chains were found to bend toward the substrates with high hydrophilicity which is in qualitative agreement with existing experiments. The amount of bending was observed to reduce with increasing film thickness. However, the side-chain length did not show any noticeable variation with wall (substrate) hydrophilicity or film thickness. The water clusters became smaller and more isolated clusters emerged for highly hydrophilic substrates. In addition, the water cluster sizes showed a decreasing trend with decreasing film thickness in the case of hydrophilic substrates, which has also been observed in experiments of supported Nafion films. The in-plane water diffusion was enhanced considerably for hydrophilic substrates, and this mechanism has also been proposed previously in experiments. The in-plane water diffusion was also found to be a strong function of the substrate selectivity toward the hydrophilic phase. Our simulations can help provide more insights to experimentalists for choosing or modifying nanoparticles for Nafion nanocomposites.

Introduction

Hydrogen, as a fuel for fuel cells,[1] is being increasingly looked as an important alternative source of energy in the transportation sector. There is a lot of research going on in ways to efficiently produce hydrogen by using energy from renewable sources such as solar and wind.[2] In addition, solar and wind energy will become a larger part of the energy mix in household and industrial usage. Batteries are needed to store energy from renewable sources and use it whenever required. In this respect, the flow batteries are being proposed as one of the solutions for large-scale energy storage.[3]

Nafion, shown in Figure , is a widely used polymer electrolyte membrane (PEM) material in PEM fuel cells (PEMFCs) and flow batteries.[4,5] The polymer membrane allows the diffusion of protons and also prevents the crossover of electrolytes in flow batteries and crossover of methanol in direct methanol fuel cells (DMFCs).[6] Nafion nanocomposites are being used to improve proton conductivity[7] and reduce unwanted crossover.[8] The PEM is also found in catalyst layers in PEMFCs along with platinum nanoparticles.[9] It is evident that there will be important interfacial interactions of the PEM with the nanoparticles and/or catalyst layers, which will affect the device performance.[10] A better and more comprehensive theoretical understanding of the effects of such interfacial interactions on phase separation and proton transport will be beneficial for designing improved nanocomposites.

Figure 1

Nafion chain (n = 7, m = 10) for equivalent weight (EW) of 1100; n represents the length of a monomer and m represents the degree of polymerization; the red oval highlights the side-chain protogenic group.

A variety of experiments have been conducted on Nafion-supported films and Nafion nanocomposites to understand and improve the performance of such materials. It has been suggested that Nafion micelles are cylindrical structures containing water and ions and lined by the sulfonic acid groups.[11] Micellar orientation of supported Nafion films was parallel to the substrate in the case of a hydrophilic substrate and it was oriented away from the more hydrophobic substrate.[12] As a consequence, it was proposed that the surface-treated nanopatterned substrates could be used to enhance water transport and ionic conductivity in a desired direction within a Nafion membrane because water and ion mobility takes place mostly along the micelles. Amide thin films on different substrates, such as silica and MgO(110), have shown differences in proton conductivity up to an order of magnitude because of differences in an interfacial structure.[13]

The confinement effect is an important factor, which influences the membrane structure and water transport. Nafion-supported films showed significant reduction in water diffusion at film thickness below 60 nm because of confinement effects.[14] Nafion films supported on silica showed increased phase separation with increasing film thickness.[15] Grazing incidence small-angle X-ray scattering (GISAXS) experiments[16] of supported Nafion films on silica showed increasing d-spacing in both in-plane and perpendicular directions with increasing film thickness with a faster rate of increase in the latter.

Nafion nanocomposites containing highly hydrophilic nanoparticles, such as modified silica, have shown higher proton conductivity as compared to that of bulk Nafion.[7] Crossover of methanol in DMFCs has been reduced by using Nafion nanocomposites.[17] Vanadium ion crossover has also been reduced by doping Nafion with nanoparticles.[8] A previous study[18] has hypothesized the reduction of methanol crossover and increased proton conduction in a Nafion-modified carbon nanotube (CNT) nanocomposite because of the formation of long oriented pathways along the modified CNTs, which were selective to water. The side-chain orientation of Nafion chains has been shown to be affected by the hydrophilicity of the substrate.[19] All of the above-mentioned experiments show that the effect of the film thickness/confinement and substrate hydrophilicity/selectivity to water on the internal hydrated nanostructure of Nafion is significant, which in turn motivates our molecular dynamics (MD) simulation-based study.

Nafion nanocomposites show a large variation in nanoparticle sizes (5–75 nm).[20,21] As a first approximation, a flat substrate of lateral dimensions in the range of 4–6 nm should be suitable for atomistic modeling of the Nafion–substrate interactions. The average interparticle distance in a Nafion titania nanocomposite was shown to be 9 nm.[20] Also, previous supported film experiments have stressed on the importance of a more thorough understanding of Nafion films in thickness less than 10 nm because this is the range of interparticle distance found commonly in catalyst layers.[15] Therefore, in the present study, the film thicknesses in the range of 6–11 nm have been chosen to perform the classical MD simulations.

Different varieties of nanoparticles such as silica, zirconia, and modified CNTs[7,18,22] have been used in Nafion nanocomposites. In addition, Nafion can also exist in catalyst layers between the carbon support and platinum nanoparticles.[9,23] All of these nanoparticles and supports have varying levels of hydrophilicity. The polymer material present between nanoparticles in nanocomposites has been modeled previously using capped films.[24] Therefore, the model of Nafion capped between substrates of varying hydrophilicity will be effective to provide insights into the importance of the interfacial interactions in Nafion nanocomposites.

Mashio et al.[25] simulated Nafion-supported films on a graphite sheet and a graphite sheet modified with carboxyl and carboxylate ions using classical MD. The number of water molecules, hydronium ions, and sulfonic acid groups was observed to increase with the presence of the ionic groups in the graphite sheet. The water clusters reduced in size for the functionalized graphite sheet than for the bare graphite sheet. Zhang and Ding[26] simulated Nafion-supported films on a platinum substrate using classical MD. The film thickness variation showed significant diversity in the water cluster morphology. The water diffusion constants varied nonmonotonically with the thickness of the Nafion film. Water diffusion in the thickest film (7.3 nm) was faster than that in a bulk Nafion. Borges et al.[27] performed classical MD simulations of the Nafion (fixed thickness)-supported films on walls of varying hydrophilicity. The films showed changes in phase separation patterns because of water flooding the highly hydrophilic walls (substrates). Water diffusion in the films was found to be greater than that in bulk Nafion at the same hydration level without any noticeable trend with varying wall (substrate) hydrophilicity. Borges et al.[28] also found that varying the hydration levels in such supported Nafion films showed distinct changes in the micellar structure within the films. Dissipative particle dynamics of Nafion films have shown preferential flooding of water at the quartz substrate,[29] which was also in agreement with experiments.[30]

Unlike a supported film, the capped Nafion films have interfaces with two substrates. This can induce additional confinement effects. Also, the side-chain orientation and the sulfonic acid (protogenic) group preferential accumulation in the presence of a substrate would be different than that in the presence of a free interface. The protogenic group locations and side-chain orientations will invariably have an effect on the water clustering within the Nafion-capped film. All these reasons make it necessary to study capped Nafion films. Classical MD simulations allow capturing the effects of deviations of film mass density from bulk mass density in capped films. In addition, atomistic representation of the Nafion molecule allows us to study the important structural properties such as sulfur–sulfur radial distribution functions (RDFs), the side-chain orientations, and the side-chain lengths. The Nafion films were capped by walls of tunable hydrophilicity in our simulations. Such a tunable hydrophilicity allows us to study the effects over a wide range of hydrophilicities as opposed to substrates with fixed chemistry. In our simulations, the side-chain orientations were found to vary with the substrate hydrophilicity and film thickness (confinement effect), while the side-chain lengths did not show any such trends. Water cluster sizes for highly hydrophilic substrates indicated that it was a function of film thickness (confinement effect). In-plane water diffusion for our capped Nafion film simulations was considerably enhanced for hydrophilic substrates, which is different from what was observed for supported Nafion film simulations by Borges et al.[27] This enhancement of in-plane water diffusion occurred despite reduced water cluster sizes for highly hydrophilic substrates. In addition, the in-plane water diffusion was found to be a strong function of the selectivity of the substrate to the hydrophilic phase.

Materials and Methods

Simulation Details

The structure of the Nafion monomer is shown in Figure . The value of n represents the number of repeat units in a monomer. The value of m is the degree of polymerization. n = 7 for this study, which corresponds to an EW of 1100. EW is defined as the weight of the polymer divided by the number of protogenic groups (sulfonic acid groups). The EW of 1100 is a very commonly used variety of Nafion and, hence, has been chosen for this study.[4,31]

The polymer consistent force field (pcff)[32] was used for simulating the polymer matrix, water molecules, and hydronium ions using LAMMPS[33] software. Partial charges for all the atoms were assigned using condensed-phase optimized molecular potentials for atomistic simulation studies (COMPASS) force field.[34]COMPASS charges have been used along with the pcff force field in previous simulations.[35−37] The pcff force field has been used to accurately model polyelectrolytes such as Nafion, SPEEK, sulfonated copolyimides, and other polymers previously.[38−42] Water molecules[38,42−44] and hydronium ions[38,41,42] have also been modeled previously using the pcff force field. Details about the force field validation are provided in the Supporting Information (section III).

Nafion was simulated at one hydration level (λ = 15) and at two different temperatures of T = 300 K and T = 353 K. The hydration level (λ) is defined as the number of water molecules per side chain of Nafion. λ = 15 was chosen because this is a moderate hydration level considering the fact that hydration levels in Nafion can go as high as λ = 30.[45] The sulfonic acid group of Nafion is fully dissociated at λ = 15.[46] Therefore, hydronium ions were introduced into the simulation box to account for this dissociation.

Integrated Lennard-Jones potential,[47]has been used to simulate structureless walls at the top and bottom of the simulation box.[27] The cutoff distance rc is chosen as 15 Å.[27] Two different sets of ϵ values are used in the simulations. ϵphob represents the interaction energy between the wall and the hydrophobic part of the system, which includes all the polymer atoms except the atoms in the sulfonic acid group, water molecules, and hydronium ions. ϵphyl represents the interaction energy between the wall and the hydrophilic part of the system, which includes all the atoms in the sulfonic acid group, water molecules, and hydronium ions.

ϵphob has been fixed at 0.25 kcal/mol and five different values of ϵphyl = 0.25, 0.50, 1.20, 1.50, and 2.00 kcal/mol have been used to simulate the effects of varying hydrophilicity of nanoparticles.[27] The paper shows results for these set of values unless mentioned otherwise. Additional simulations have been performed in which ϵphyl has been fixed at 2.00 kcal/mol and five different values of ϵphob = 0.25, 0.50, 1.20, 1.50, and 2.00 kcal/mol have been used to understand the effect of contrast between the ϵphyl and ϵphob on the water transport within the Nafion-capped film.

A nanocomposite has fillers/nanoparticles dispersed inside the matrix (polymer). The matrix material present between any two nanoparticles is the representative volume element (RVE) being modeled in our simulations. This RVE was modeled by confining 17 Nafion chains along with water molecules and hydronium ions between structureless walls of tunable hydrophilicity,[27,28] as shown in Figure a. The walls represent the nanoparticle surfaces of variable hydrophilicity. This representation has been used to model nanocomposites previously.[24,48] Henceforth, ϵphyl = 0.25 and 0.50 kcal/mol walls will be referred to as low hydrophilicity (LH) walls and ϵphyl = 1.20, 1.50, and 2.00 kcal/mol walls will be referred to as high hydrophilicity (HH) walls in what follows. In both these cases, ϵphob has been fixed at 0.25 kcal/mol unless mentioned otherwise.

Figure 2

(a) Hydrated Nafion film between two structureless walls. Z-axis is the direction perpendicular to the walls, and X and Y axes are parallel to the walls. Z-direction has fixed boundaries, and the film is periodic in X and Y directions. Blue color represents water molecules and hydronium ions, and orange color is used for Nafion molecules. (b) Side-chain vector (vector connecting the first carbon in the backbone to the sulfur in the sulfonic acid group) orientation, that is, the angle between the side-chain vector and the Z-axis. The simulation box is divided into three equal layers as shown in (a), and the side-chain orientation was computed in these three layers.

Three different film thickness values of 6.3, 8.7, and 11.5 nm were simulated for each of the wall hydrophilicity values. The film thickness was varied in the Z-direction (Figure ). The thickness variation represented effectively the variation of the filler fraction in a nanocomposite, that is, higher film thickness corresponds to lower filler fraction and vice versa. The simulations were run for a total of 8 ns, and the last 3 ns of the production runs was used for analysis. The density, with a variation less than 0.05%, had stabilized after 2.5 ns from the start of the simulation and the energy was also stable. The average water cluster size showed variations less than 1% during the production run. A detailed description of the model construction and simulation protocol have been presented in the Supporting Information (sections I and II). Note that all the results shown in this paper are for T = 353 K. Qualitatively similar results were obtained for T = 300 K. The duration of our simulations and implemented system sizes are consistent with previous simulation studies.[4,49−51] Each film simulation consumed around 120–170 CPU hours on 32 cores of the Lisa computing cluster in SURFsara (Amsterdam).

Calculation Methods

From the production runs, structural and dynamic characteristics such as RDFs, side-chain orientations, cluster distribution of water molecules and/or hydronium ions, and diffusion coefficients of water molecules were calculated.

The number density is defined as the number of atoms of a particular type divided by the total number of atoms of the same type in a layer of thickness 0.2 Å at a particular distance from either of the walls in the Z-direction. The total mass density is defined as the total mass present in a 0.2 Å thick layer, at a particular distance (in Z-direction) from either of the walls, by the volume of the layer.

The RDF g(r) is proportional to the probability of finding an atom B at a distance r from the reference atom A inside a shell of thickness dr.[31] The coordination number (CN) is the average number of atoms of a particular type found at a certain distance from a particular central atom of a certain type. The sulfur–sulfur RDF has been analyzed to check for any significant changes in the distance between the side-chain protogenic groups. The sulfur–sulfur and sulfur–oxygen (water) CNs have also been analyzed.

The side-chain orientation is defined as the angle between the side-chain vector and the Z-axis, as shown in Figure b. The side-chain vector is defined as the vector from the carbon connecting the side chain to the backbone toward the sulfur in the sulfonic acid group. The simulation box was divided into three equal layers from top to bottom. The angle between the side-chain vector and the +Z-axis was computed for the top two layers, and the angle between the side-chain vector and the −Z-axis was computed for the bottom layer. Analysis was done using a custom MATLAB script.

The cluster distribution of water molecules was computed for the different hydration levels (λ) using the OVITO software.[52] A cluster is defined as a group of atoms in which each atom is within a particular predefined cutoff distance of at least another atom within that group. The oxygen atom in the water molecule was used for computing cluster sizes, that is, cluster of 10 oxygen atoms is assumed to represent the cluster of 10 water molecules. The cluster distribution plots number of clusters, averaged over a time interval, versus the cluster size. Python scripts were used to access OVITO API, and MATLAB scripts were used for further postprocessing.

Water channel sizes were computed using the Zeo++ software[53] which uses Voronoi tessellation for its internal calculations. All the atoms associated with the water molecules were removed, and the remaining atom positions and types were provided as input to this software for channel size computation. The water channel sizes for the HH walls were estimated from the water number density profiles, which has been explained in Water Cluster Distribution section later.

The translational diffusion coefficients for water molecules were computed by analyzing their mean square displacement (MSD) using the Einstein relation in the diffusive regime.[54] These diffusion coefficients were computed as an average for the entire Nafion film between the walls. The simulation box was divided into five equal layers from top to bottom in the Z-direction, and the layer-resolved diffusion coefficients were also computed in these layers using a custom MATLAB script.

Results and Discussion

Snapshots and Density Profiles

Figure shows the snapshots at the end of the production runs for five different values of wall hydrophilicity. For low values of wall hydrophilicity (ϵphyl = 0.25 and 0.50 kcal/mol), there is negligible accumulation of water molecules near the walls. However, for the HH walls (ϵphyl = 1.20, 1.50, and 2.00 kcal/mol), there is a considerable accumulation of water near the walls.

Figure 3

Snapshots for (a) ϵphyl = 0.25 kcal/mol, (b) ϵphyl = 0.50 kcal/mol, (c) ϵphyl = 1.20 kcal/mol, (d) ϵphyl = 1.50 kcal/mol, and (e) ϵphyl = 2.00 kcal/mol where blue color shows the water molecules and hydronium ions, and orange color shows the Nafion atoms.

Figure shows the water number density profiles for different wall hydrophilicity values for the 6.3 nm film. The water number density profile for the ϵphyl = 0.25 kcal/mol wall was very uniform throughout the thickness of the film. The water number density shows small peaks near the walls for a slightly higher hydrophilicity wall (ϵphyl = 0.50 kcal/mol). Previous simulations of Nafion supported on a primary hydrophobic graphite have shown emergence of small peaks in the water density profiles when the graphite was modified by hydrophilic carboxylate ions.[25] The water number density near the walls is much higher for HH walls than for the LH walls because of a considerable accumulation of water close to the walls, as shown in Figure . Also, the water number density near the center of the film is higher for the LH walls than for the HH walls.

Figure 4

Water (oxygen) number density profiles for the 6.3 nm film at different wall hydrophilicity (ϵphyl) values.

Figure shows the water, carbon, and sulfur number density profiles for the lowest and highest hydrophilicity wall for the 6.3 nm film. The carbon number density shows small peaks near the walls in the LH case. These peaks disappear for the HH wall and a smooth profile appears, which reaches a stable value at a distance further away from the walls as compared to the LH case. This suggests that carbon atoms are moving away from the walls in the highest hydrophilicity case because of water accumulation near the walls. The sulfur number density is similar to the water number density profile, that is, number density near the walls is considerably higher for the HH wall than for the LH wall. This indicates that both sulfur and water show preferential accumulation at the HH wall. Similar trends in the carbon, sulfur, and water density profiles are also observed for other film thicknesses.

Figure 5

Water (oxygen), carbon, and sulfur number density profiles for the 6.3 nm film for (a) LH film, ϵphyl = 0.25 kcal/mol and (b) HH film, ϵphyl = 2.00 kcal/mol.

Figure shows the total mass density profiles, normalized by the bulk density, for different film thickness values for the lowest and highest hydrophilicity walls. The X-axis in Figure is the relative distance (t/T), defined as the distance (t) from a wall divided by the film thickness (T). The bulk domain is defined as the space where the normalized mass density is equal to 1. For both the LH and HH cases, a broadening of the bulk domains can be observed with increasing film thickness. However, there is an important difference between these two cases. The LH wall shows an almost uniform density profile throughout the film thickness, whereas the HH wall shows high density values near the walls. These high density values are due to the preferential accumulation of the hydrophilic components, such as water and sulfur, near the walls.

Figure 6

Total film mass density, normalized by the bulk density = 1.79 g/cm3, profiles for different film thickness values for (a) LH film, ϵphyl = 0.25 kcal/mol and (b) HH film, ϵphyl = 2.00 kcal/mol.

RDFs and CNs

The distance between protogenic sulfonic acid groups is an important characteristic to probe the internal structure of the membrane. Previous simulation studies have shown that the sulfur–sulfur distance (interprotogenic group distance) less than 6.5 Å increased water binding to sulfonic acid groups and also affected the ease of the proton dissociation.[55] Hence, the sulfur–sulfur RDF at small atomic separations (<8 Å) and sulfur–sulfur CN have been analyzed in this study to check for any significant changes. The water structure around the sulfur atom is also important for proton dissociation.[46] The sulfur–oxygen (water) CN has also been analyzed to check for any significant changes in the first hydration shell (∼4.7 Å, Figure S1).

Figure a shows the S–S RDF values for different levels of the wall hydrophilicity for the 6.3 nm film. The RDF plots for ϵphyl = 0.25, 0.50, and 1.20 kcal/mol have their first maximum at almost the same distance as the bulk Nafion first peak distance of 4.3 Å. RDF plots for ϵphyl = 1.50 and 2.00 kcal/mol have the first peak at a slightly higher distance of 4.5 Å. The Nafion chains tend to move away from the HH walls because of the preferential accumulation of water and are packed into a more confined space toward the center of the film. This could increase the repulsion between the negatively charged sulfonic acid groups, which can explain the slightly higher distance of the first peak for ϵphyl = 1.50 and 2.00 kcal/mol walls. This effect can also be seen for the ϵphyl = 2.00 kcal/mol wall for higher film thicknesses (Figure S2). In conclusion, the negligible increase in the position of the first peak implies that the hydrophilicity of the substrate has no considerable effect on the S–S distance.

Figure 7

(a) Sulfur–sulfur RDFs and (b) S–S CNs for the 6.3 nm film and the different wall hydrophilicity values (ϵphyl). RDF and CN for bulk Nafion have also been shown.

However, the values of the S–S RDFs for the HH walls (ϵphyl = 1.20, 1.50, and 2.00 kcal/mol) were visibly higher than those for the LH walls (ϵphyl = 0.25 and 0.50 kcal/mol) up to a distance of 8 Å (Figure a). A similar trend is also observed for thicker films (Figure S2). The accumulation of the Nafion chains near the center of the film would increase the probability of finding a sulfur atom within the close proximity of another sulfur atom. This would explain the rise in the RDF values for the HH walls.

Figure b shows the S–S CNs for the 6.3 nm film for different wall hydrophilicity levels. The HH walls show a slightly higher CN as compared to the LH walls for all distances. This pattern is seen for higher film thicknesses as well (Figure S3). These slightly higher CNs for HH walls combined with no noticeable change in S–S RDF for close range distances (<8 Å) imply that there should not be any detrimental effect on S–S close range ordering in the HH wall films.

Figure shows the sulfur–oxygen (water) (S–Ow) CNs for the 6.3 nm film for different wall hydrophilicity levels. The S–Ow CNs are slightly lower for HH walls than for the LH walls, and this pattern is seen for higher film thicknesses as well (Figure S4). This is due to the fact that a considerable amount of water accumulates near the HH walls, which reduces the average number of water molecules near the sulfur atoms. However, the CNs for both the LH and HH walls are quite similar up to 4.7 Å (first coordination shell of S–Ow, Figure S1), which implies that the high wall hydrophilicity does not have a detrimental effect on close range water solvation structure around the sulfur atoms in the sulfonic acid group.

Figure 8

S–Ow CNs for the 6.3 nm film and the different wall hydrophilicity values (ϵphyl). CN for bulk Nafion has also been shown.

Side-Chain Orientation

The catalyst layer in fuel cells can have Nafion present in between platinum nanoparticles.[9] Also, Nafion can have an interface with carbon in the catalyst layer and in the electrodes in a fuel cell.[9] The side-chain orientation with respect to the nanoparticles will have an impact on the compatibility of such nanocomposites.[19] Therefore, the effects of varying substrate hydrophilicity and filler fraction on the side-chain orientation have been investigated.

Figure shows the layer-resolved orientation of side chains for the 6.3 nm film for different values of wall hydrophilicity in three equal film layers. As can be seen, the angle between the side chain and the Z-axis does not show noticeable variation across the layers for the LH walls. However, the angle value reduces considerably in the top and bottom layers for the HH walls. A similar pattern was also observed for higher film thickness values (Figure S5).

Figure 9

Side-chain orientation with respect to the Z-axis in three layers for the 6.3 nm film for different wall hydrophilicity (ϵphyl) values.

The simulations for each of the different hydrophilicity values started from the same initial configuration. Therefore, the decrease in the angle in the top and bottom layers for the HH walls is purely due to the effect of the walls. This same effect was also seen for a different set of initial configuration values, which further confirms our hypothesis. It can be concluded that the HH walls tend to bend the side chains toward them, which agrees qualitatively with experimental observations.[19]

Figure shows the side-chain orientation for the lowest and highest hydrophilicity wall. The angles for the lowest hydrophilicity wall in all three layers show a negligible variation with varying film thicknesses. A similar trend is also observed for the other LH wall (ϵphyl = 0.50 kcal/mol) (Figure S5a). In contrast, the value of the angles in the top and bottom layer increases progressively with increasing film thickness for the highest hydrophilicity wall. A similar trend is observed for another HH wall (ϵphyl = 1.50 kcal/mol) as well (Figure S5c). For the remaining HH wall (ϵphyl = 1.20 kcal/mol), the amount of side chain bending toward the walls was subdued for film thicknesses of 8.7 and 11.5 nm (Figure S5b). This implies that the effect of side chains bending toward the walls for the HH walls reduces with the increasing film thickness, at least in the film thickness range investigated. In conclusion, varying the nanoparticle hydrophilicity and/or nanocomposite filler fraction can be used to alter the side-chain orientation with respect to the nanoparticle surface.

Figure 10

Side-chain orientation with respect to Z-axis in three layers for different film thickness values for (a) ϵphyl = 0.25 kcal/mol LH wall and (b) ϵphyl = 2.00 kcal/mol HH wall.

Side-Chain Length

The differences in the side-chain lengths have been shown to affect the diffusion of water within the PEMs.[50] Therefore, the side-chain lengths have been analyzed for different wall hydrophilicities and the film thicknesses, and the same has been shown in Figure . The bulk Nafion side-chain length was on average about 7.3 Å. The capped Nafion film side-chain lengths did not show any considerable deviation from the bulk value. There was no noticeable trend for the side-chain lengths with wall hydrophilicity and/or film thickness. Therefore, it can be concluded that changes in water diffusion amount (if any) should not be due to the changes in side-chain lengths in capped Nafion films.

Figure 11

Nafion side-chain lengths for different wall hydrophilicity (ϵphyl) values and different film thicknesses.

Water Cluster Distribution

Water clusters present in the hydrated Nafion nanostructure form percolated channels at sufficiently high hydration levels.[56] These percolated channels allow the transport of protons across the membrane, which allows the fuel cells and flow batteries to function. However, these percolated water channels can also allow unwanted crossover of methanol and vanadium ions. Nafion nanocomposites have been shown to reduce the crossover of methanol[18] and vanadium ions.[8] Therefore, it is important to understand the effect of nanoparticle hydrophilicity and filler fraction on the water cluster distribution.

All the water cluster analysis shown here are for a cutoff distance of 3.7 Å (see the Calculation Methods section earlier) averaged over 3 ns of simulated physical time. This cutoff distance was chosen because it is close to the first coordination shell (Figure S6) of water, and hence, this distance will encompass a majority of the water molecules. A single large cluster is observed for this cutoff distance of 3.7 Å for bulk Nafion (Figure ). No such large cluster is observed for a cutoff distance of 3 Å for bulk Nafion because this distance is close to the first peak of oxygen (water)–oxygen (water) RDF (Figure S6) and, hence, encompasses few water molecules.

Figure 12

Water cluster distribution for the 6.3 nm film at different wall hydrophilicity values and also for bulk Nafion. The cluster distribution shown is for the cluster sizes from 100 to 2380.

Figure shows the water cluster distribution at different wall hydrophilicity values for the 6.3 nm film. The cluster distributions for the LH walls (ϵphyl = 0.25 and 0.50 kcal/mol) are very close to the bulk cluster distribution. The largest clusters (cluster sizes close to 2400) for the HH walls (ϵphyl = 1.20, 1.50, and 2.00 kcal/mol) are lesser in number as compared to those for the LH walls. Also, there is an emergence of clusters in the size range of 900–1500 for the HH walls. This shows that the cluster sizes decrease considerably for the HH walls for a fixed film thickness. Similar effects are also seen at higher film thicknesses (Figure S7).

The insets in Figure show a continuous percolating cluster for the LH wall films spanning the whole box in all three dimensions and isolated clusters for the HH wall films near the center of the box (film). The HH wall films also form two percolating roughly cuboidal water channels along the walls, which are reflected in the two peaks in the cluster distribution at around 800 and 1400. The thickness of these channels in the Z-direction for the 6.3 nm film is around 9–10 Å as estimated from the water number density plots (Figure ), that is, difference between the distance where the number density reaches its minimum after the 3rd peak from any wall (box edge) and the distance where the number density first acquires a nonzero value. The only water channel in the lowest hydrophilicity wall 6.3 nm film had a maximum channel diameter of 13.6 Å and a minimum diameter of 6.5 Å. The corresponding quantities for the single channel in bulk Nafion were 11.5 and 5.3 Å, respectively. It is clear that water, which has an average van der Waals diameter of 2.8 Å, can diffuse through both the LH and HH wall films. However, the HH wall film provides uniformly wide and straight water channels along the walls, whereas the LH wall film water channel has bottlenecks (minimum diameter) and is more tortuous (extends through the box in all three dimensions). As a result, water is observed to diffuse noticeably faster in the in-plane direction for the HH wall films than for the LH wall films, which has been discussed later in the Water Transport section.

Figure shows the water cluster count normalized by the bulk water cluster count for different wall hydrophilicity and film thickness values. All the normalized cluster counts are larger than 1, which implies a larger number of water clusters for all the wall hydrophilicity values and the film thickness, as compared to bulk Nafion water cluster count. This indicates a more dispersed water cluster network in the Nafion films as compared to the bulk Nafion.

Figure 13

Cluster count, normalized by the bulk cluster count, vs different wall hydrophilicity (ϵphyl) values for different film thickness values.

The water cluster count is higher for the HH walls than for the LH walls for all three different film thicknesses (Figure ). This effect is universal and is weakly dependent on the film thickness. The higher cluster count indicates a more dispersed water cluster network for the HH wall films than for the LH wall films, which can also be seen in the inset for the HH wall films in Figure . The existing experiments have shown that unwanted crossover reduces because of the highly hydrophilic nanoparticles such as silica, clay, and so forth[17] added to Nafion. In fact, the existence of long-range oriented pathways along the modified CNTs was the proposed mechanism for the observed enhanced proton transport and reduced methanol crossover in a Nafion-modified CNT nanocomposite.[18] Our simulations also show the preferential accumulation of water along the HH walls and a concomitant increase in the water cluster count because of the emergence of a more dispersed water phase and isolated water clusters. It is likely that less polar molecules such as methanol will move away from the highly hydrophilic nanoparticles similar to carbon moving away from the HH walls as seen in Figure b. This will increase the chances of such molecules being trapped in the isolated clusters, which are found at larger distances from the HH walls as seen in the inset for HH wall films in Figure .

Figure shows the average water cluster size for different wall hydrophilicity and film thickness values. The average cluster sizes are almost constant with increasing film thickness values for the LH walls. However, there is a distinct pattern for the HH walls, which shows that the average water cluster sizes show an increasing trend with increasing film thickness. The water channels along the HH walls become more connected through the center of the film (box) with increasing film thickness. This is evidenced by the increasing average number of clusters in the 2200–2330 size range and a concomitant decrease in the 800–1400 size range with increasing film thickness (Figures and S7). This behavior for the HH walls indicates higher phase separation for increasing film thickness. Previous transmission electron microscopy images[15] and GISAXS experiments[16] also show similar trends versus film thickness for Nafion films supported on hydrophilic silica substrates.

Figure 14

Average cluster size vs different film thickness values for different wall hydrophilicity (ϵphyl) values. Average cluster size for bulk Nafion is also shown. The dashed–dotted line shows the trend for LH walls, and the dashed line shows the trend for HH walls. These lines are not numerical fits.

Water Transport

Water diffusion through the Nafion nanostructure plays an important role in the working of fuel cells and flow batteries. In these devices, the protons attach themselves to water molecules and are transported from anodic to cathodic side or vice versa via the so-called vehicular transport mechanism. There is also an alternative method of proton transport in which the proton hops across hydrogen bonds in the water phase. Nanoparticles are added to Nafion to enhance water retention and proton conductivity.[7,57] Therefore, it is important to study the effect on water transport because of the nanoparticle hydrophilicity and filler fraction.

Water diffusion constants (D, D, and D) have been computed in the X, Y, and Z directions using the Einstein relation for diffusive motion. Diffusion constants have been calculated from the time period where water transport is in a diffusive regime (Figure S8). Water diffusion in the films’ XY-plane is studied using the in-plane diffusivity (D)and was compared to the analogous (two-third of total water diffusion coefficient) values for Nafion bulk (Dbulk)

Henceforth, the diffusion in the XY-plane will be referred to as in-plane water transport. The total water diffusion coefficient (1.5 × Dbulk) in bulk Nafion for λ = 15 at T = 353 K was found to be 1.93 × 10–5 cm2/s from our simulations.a

Figure a shows the in-plane water diffusion, normalized by the bulk water diffusion constant for different wall hydrophilicity values and for the different film thicknesses. The in-plane water diffusion is noticeably higher for the HH walls (ϵphyl = 1.20, 1.50, and 2.00) than for the LH walls (ϵphyl = 0.25, 0.50), for all three different film thicknesses. The existing experiments have shown micellar orientation in Nafion-supported films along hydrophilic substrates and away from hydrophobic substrates.[12] It was further proposed that treated nanopatterned substrates can be used to enhance the directional transport of water within the Nafion membrane because the water transport is mostly along the micelles.[12] A similar enhancement of water transport for the HH walls (substrates) is also observed in our simulations.

Figure 15

Film-averaged (a) in-plane water diffusion constants (D) normalized by the corresponding two-dimensional (2D) water diffusion (Dbulk) constant at λ = 15 for bulk Nafion. (b) Water diffusion anisotropy ratio values vs wall hydrophilicity (ϵphyl) for different film thicknesses.

It is important to keep in mind that the water cluster sizes showed a significant size reduction for the HH walls. The bulk classical MD simulations of PEMs such as Nafion,[49] SPEEK,[51] and PFIA[58] have shown the water diffusion to increase with increasing water cluster sizes. However, the capped Nafion films show a decrease in water diffusion despite larger water cluster sizes for the LH walls. This is due to the formation of water channels parallel to the HH walls with a uniform width as compared to the long tortuous water channel with bottlenecks in the LH wall film.

Previous simulations done for a supported Nafion film did not show any noticeable distinction between water diffusion constants for less and more hydrophilic substrates.[27] However, the capped Nafion films, simulated in this paper, show a clear difference in water diffusion rates between the LH and HH walls. The possible reason for this behavior has been explained later by analyzing the layer-resolved in-plane diffusion. There was no monotonicity observed in the film-averaged diffusion coefficients with respect to the film thickness in the thickness range investigated. This nonmonotonicity observation agrees well with the previous experimental conductivity[15] and simulated water diffusion rates[26] in the film thickness range investigated.

Anisotropy in water diffusion is defined as the in-plane water diffusion constant (D) divided by twice the Z-direction diffusion constant (D). Figure b shows the anisotropy in water diffusion for different hydrophilicity walls and different film thicknesses. Anisotropy for the HH walls is higher than that for the LH walls. This effect is due to the high in-plane diffusion near the walls in the HH walls. Anisotropy in water diffusion, for the HH walls, decreases considerably on increasing the film thickness from 6.3 to 8.7 nm. This can be explained possibly by the strong confinement in the Z-direction in the thinnest film.

In-plane water transport was resolved in five equal layers in the Z-direction. Figure shows the layer-resolved in-plane water diffusion constant, normalized by the bulk 2D water diffusion constant, for the 6.3 nm film. The water diffusion constants are slightly smaller than bulk values throughout the thickness of the film for the LH walls. For the HH walls, the diffusion constant near the center of the film is close to that for the bulk, but the diffusion increases considerably on moving closer to the walls. Similar trends are also observed for higher film thicknesses (Figure S9). It is the presence of such highly mobile water layers near both the walls in a capped Nafion film that can explain the noticeably high film-averaged in-plane water diffusion constant for the HH walls. Previously simulated supported Nafion films were shown to have considerably less in-plane water diffusion near the free interface as compared to that near the highly hydrophilic substrates.[27] In contrast, we observe occurrence of highly mobile layers at both the walls (substrates) for the HH cases across all the film thicknesses. This fact can explain the considerably higher film-averaged in-plane diffusion for the capped Nafion films confined by HH walls than by LH walls unlike the previously simulated[27] supported Nafion films.

Figure 16

Layer-resolved in-plane water diffusion constants (D) normalized by the 2D water diffusion constant (Dbulk) at λ = 15 for bulk Nafion. Results are shown for the 6.3 nm film for varying wall hydrophilicity (ϵphyl).

Figure shows the layer-resolved in-plane water transport for the lowest and highest hydrophilicity wall for different film thicknesses. For the lowest hydrophilicity wall, the in-plane water diffusion constants start to deviate more from the bulk diffusion constant with increasing film thickness. Water is confined towards the center of the film in the LH wall films. Increasing film thickness allows water more space to move inside the film away from the walls. This could be the reason for the slightly higher deviation from bulk diffusion values for the higher film thicknesses. For the highest hydrophilicity wall, water diffusion is not affected noticeably by the film thickness. Water is mostly concentrated near the walls for these cases, and so, increasing film thickness plays a negligible role in in-plane water diffusion in the thickness range investigated.

Figure 17

Layer-resolved in-plane water diffusion constants (D) normalized by the 2D water diffusion constant (Dbulk) at λ = 15 for bulk Nafion. Results are shown for varying film thicknesses for (a) ϵphyl = 0.25 kcal/mol wall and (b) ϵphyl = 2.00 kcal/mol wall.

Enhanced film averaged in-plane transport of water has been observed for the HH walls. It is important to ascertain whether this high in-plane transport is due to just the  high hydrophilicity of the walls or due to the contrast between the ϵphyl and ϵphob. Simulations were run by fixing ϵphyl = 2.00 kcal/mol and varying ϵphob = 0.25, 0.50, 1.20, 1.50, 2.00 kcal/mol to understand the effect of this contrast on the in-plane water transport. Figure shows the in-plane water diffusion constants, normalized by the corresponding bulk values, for different ϵphob values for the 6.3 nm film. In-plane water diffusion constants decreased with increasing ϵphob values. This implies that the water diffusion is a function of the contrast between the ϵphyl and ϵphob values. Higher contrast results in more in-plane water diffusion. The insets also show a preferential accumulation of water at low ϵphob values or high contrast between the ϵphyl and ϵphob values. Nanoparticles can be made more selective toward water by modifying their surface.

Figure 18

Film-averaged in-plane water diffusion constants (D) normalized by 2D water diffusion (Dbulk) constant at λ = 15 for bulk Nafion. Results are shown for the 6.3 nm film. ϵphyl = 2.00 kcal/mol is kept fixed, and ϵphob is varying.

Conclusions

Capped Nafion films were simulated at a moderate hydration level (λ = 15) at T = 353 K to model the interactions present in a RVE consisting of the matrix (hydrated polymer) confined between any two nanoparticles in a Nafion nanocomposite. The Nafion films were capped by walls of different hydrophilicities to study the effect of nanoparticle hydrophilicity on the Nafion nanostructure. The film thickness was varied to study the effect of nanoparticle filler fraction on the Nafion nanostructure and water transport.

The simulated sulfur–sulfur RDFs indicated that there was a negligible effect on the close range sulfur–sulfur distance because of the wall hydrophilicity and the film thickness. Although the RDFs and the CNs suggested that sulfur atoms were more likely to be near each other in a distance less than 8 Å for the higher hydrophilicity walls. The number of water molecules around the sulfonic acid group in the first solvation shell (up to a distance of 4.7 Å) also showed negligible differences with varying wall hydrophilicity. Therefore, nanoparticle hydrophilicity and filler fraction should not have a detrimental effect on the sulfur–sulfur close range distance or the close range hydration structure around the sulfonic acid group.

The Nafion side-chain lengths did not show any noticeable trend with wall hydrophilicity and/or film thickness. However, the side chains were found to bend toward the HH walls. Also, the amount of bending reduced with increased film thickness for the HH walls. Experiments have also shown increased bending of the side chains toward highly hydrophilic substrates.[19] In effect, nanoparticle hydrophilicity and filler fraction could be used to control side-chain orientation with respect to the nanoparticles.

Reduced crossover of methanol has been observed in experiments for Nafion doped with hydrophilic nanoparticles.[17,18] The emergence of isolated water clusters as indicated by the higher cluster count for the highly hydrophilic substrates could explain such experimental observations. On average, the water cluster sizes increased with increasing film thickness for the HH walls, which indicates stronger phase separation with increasing thickness. Qualitatively similar experimental observations have been seen for supported Nafion films on a silica substrate.[15,16]

The water in-plane transport was enhanced considerably by the HH walls, in spite of lower water cluster sizes for the HH wall films. This effect was observed for all the different film thicknesses. Layer-resolved in-plane transport indicated a very highly mobile water layer near both the walls for the HH wall films. The LH wall (ϵphyl = 0.25 kcal/mol) 6.3 nm film had a single water channel with a maximum and minimum diameter of 13.6 and 6.5 Å, respectively. The water channels in the HH wall films were roughly cuboidal blocks along the walls with a width of 9–10 Å. The water channels in the HH wall films had no bottlenecks (minimum diameter) and were visibly less tortuous than the water channel in the LH wall film. These observations explain the enhanced in-plane film-averaged transport for the HH wall films. Previous experiments have also proposed directional enhancement of water transport by altering the hydrophilicity of substrates.[12] The water diffusion anisotropy ratio, ratio of in-plane water diffusion to the water diffusion in the perpendicular direction, was noticeably higher for HH wall films than for the LH wall films. Water diffusion anisotropy ratio appeared to reduce with increasing film thickness.

The in-plane water transport was also examined for different values of the contrast between ϵphyl and ϵphob values. In-plane diffusion of water was enhanced for the larger contrast between ϵphyl and ϵphob. In effect, the enhanced in-plane transport of water was found to be a function of the contrast between ϵphyl and ϵphob values. This fact can be of use in designing nanoparticles by increasing the selectivity to the hydrophilic phase.

To summarize, our simulations showed that high selectivity of walls (nanoparticles surfaces) toward the water phase results in water channels forming along the walls and isolated water clusters emerging at distances further away from the walls. Less polar molecules such as methanol are likely to move away from the selectively hydrophilic surfaces (walls) and get trapped in these isolated clusters. Our simulations show that average water cluster sizes for hydrated Nafion films confined by highly hydrophilic surfaces (HH walls) increase with increasing film thickness. This is due to the increasing connectivity between the water channels, which form along the highly hydrophilic surfaces, through the center of the film. Therefore, it can be suggested that increasing the filler fraction (reducing film thickness) will lead to lesser connectivity of these water channels, which form along the hydrophilic surfaces, at larger distances from these surfaces, which in turn can lead to lower crossover of low polarity molecules such as methanol. Our simulations also showed that the side chains of Nafion can be made to orient toward highly hydrophilic surfaces (HH walls) and that the amount of orientation can be increased by reducing the film thickness (increasing filler percentage of nanoparticles). Our simulations show that water is preferentially accumulated near the highly hydrophilic surfaces (HH walls). In addition, the close range solvation structure near the sulfonic acid group was minimally affected with varying wall hydrophilicity. Therefore, the orientation of side chains toward the water-rich environment near these surfaces (hydrophilic nanoparticle surfaces) can be advantageous in a high-temperature environment. Our simulations also showed directional enhancement of water transport for highly hydrophilic surfaces (HH wall films) because of highly mobile water layers along these surfaces.

In the present study, we did not explicitly model curvature effects of the surface (wall) on the transport of water. Nevertheless, these effects can be important. Experiments have shown that one-dimensional and 2D nanoparticles such as modified CNT and graphene oxide result in higher proton conductivity[18,59,60] because of long-range transport along these nanoparticles and the ordering of these nanoparticles themselves. Larger-scale simulations incorporating explicit fillers into the simulation box will provide more insights into these effects.