Effect of Neutralization on the Structure and Dynamics of Model

Jan 3, 2018 - These stickers are similar to backbone monomers but have higher cohesive interaction strength with each other and with the ionic groups,...
0 downloads 8 Views 4MB Size
Article Cite This: Macromolecules XXXX, XXX, XXX−XXX

Effect of Neutralization on the Structure and Dynamics of Model Ionomer Melts Janani Sampath and Lisa M. Hall* William G. Lowrie Department of Chemical and Biomolecular Engineering, The Ohio State University, 151 W. Woodruff Ave., Columbus, Ohio 43210, United States S Supporting Information *

ABSTRACT: Using coarse-grained molecular dynamics simulations, we investigate how ionomer behavior depends on ion and acid content. Our coarse-grained model builds on a prior model that was successful in capturing the low wavevector ionomer peak observed in X-ray scattering. As in prior work, we consider a linear backbone of neutral monomer beads with charged pendant beads and unbound counterions; we now also include “stickers” pendant to the backbone to represent unneutralized acid groups. These stickers are similar to backbone monomers but have higher cohesive interaction strength with each other and with the ionic groups, chosen to approximately match structural features of prior atomistic simulations. This allows us to simulate partially neutralized materials over longer time and length scales than is possible atomistically. We find that stickers can aid in ionic aggregate rearrangements such that samples containing sticker groups can show a lower viscosity than samples with the same ion content but without stickers; this agrees with experimental results of analogous materials.



deformation on aggregate morphology.18 All pendant moieties of the coarse-grained model ionomers were charged (representing 100% neutralization with sodium); however, typical experimental systems contain both acid groups (−COOH) and anions (−COO−) pendant to the polymer backbone.19−22 The level of neutralization is an important determinant of ionomer properties, affecting both dynamics and aggregate morphology. In fact, the atomistic MD work of Bolinteneneau et al.8,9 showed that both the charge−charge interactions (between carboxylate anions and sodium counterions) and the hydrogen bonding between carboxyl groups were important factors setting the overall aggregate structure. This work showed a wide variety of aggregate morphologies, vastly different from the predominantly spherical aggregates assumed by some scattering models.23,24 Recent work by Winey et al. showed that the simulated scattering is in excellent agreement with experimental scattering.10 Regarding the impact of ionic groups on dynamic properties, Register and co-workers have extensively studied the rheological behavior of ethylene−methacrylate ionomers at different ion concentrations and neutralization levels.20,22,25−28 These ionomers are similar to PEAA but include an additional methyl group bonded to the same carbon as the carboxyl group, as in the commercially important ionomer Surlyn. The Register group’s work showed that polyethylene-co-methacrylic acid (PEMAA) materials (with carboxyl groups but no ionic

INTRODUCTION Ionomers have been the subject of ongoing study for decades due to their mechanical strength and potentially tunable properties, with applications ranging from water treatment to energy storage.1−4 Their desirable mechanical properties are a result of the presence of a small fraction of charged groups (lower than 15%) covalently bound to a nonpolar backbone, which, along with their counterions, form nanoscale ionic aggregates. Because of the importance of these aggregates in determining the overall material properties, detailed knowledge of the aggregate microstructure and its effect on dynamic properties can facilitate future design and utilization of ionomeric materials. Experimental and simulation work focused on polyethyleneco-acrylic acid (PEAA) polymers and their neutralized ionomer counterparts has provided significant insight into the aggregate structure and resulting bulk properties of associating polymers.5−10 Of particular interest here is the recent experimental work by Winey and others in characterizing PEAA and partially neutralized PEAA with precisely spaced AA groups.5,11−13 There has been significant coarse-grained and atomistic molecular dynamics (MD) modeling of analogous materials by Frischknecht and co-workers.6,8−10,14A coarsegrained, fully neutralized ionomer model developed by Hall et al. was able to replicate the major trends in the low wavevector scattering peak seen in small-angle X-ray scattering (SAXS) experiments as a function of spacing,15−17 which is related to the aggregate size and interaggregate ordering. Using this fully neutralized coarse-grained model, we previously reported rheological properties and the effect of uniaxial tensile © XXXX American Chemical Society

Received: September 24, 2017 Revised: November 28, 2017

A

DOI: 10.1021/acs.macromol.7b02073 Macromolecules XXXX, XXX, XXX−XXX

Article

Macromolecules

report rheological properties and compare to prior experimental work of Register and co-workers.20,22

groups) have a higher viscosity than their esterified derivatives (in which the acid groups have been converted to esters which do not strongly associate), as expected due to the hydrogen bonding between the acid groups. However, the viscosity of partially neutralized materials (with both acid and ionic groups) is lower than that of their esterified derivatives (in which the same ionic groups are present but the acid groups are esterified). In other words, acid groups increased viscosity in uncharged materials but decreased viscosity when ionic groups were also present. They described this as a “plasticization” of the ionic associations by acid groups, which accelerates the overall chain relaxation. In the current work, we update the prior coarse-grained model to account for partial neutralization, allowing efficient study of the long-time dynamics of these materials including the effects of the carboxyl groups. Specifically, we include “sticker” groups that represent carboxyl moieties; these groups have additional cohesive interactions with other sticker groups and with ions (relative to standard monomer−monomer interactions). While the hydrogen-bonding interaction between carboxyl groups is directional, the sticker−sticker interactions are radially symmetric for simplicity. Note that some prior work on associating polymer melts has used “sticker” groups to account for the association of ionic groups by using stronger interactions but without the long-ranged Coulombic interaction.29,30 While such models can provide a general understanding of association effects, they do not reproduce local structures of alternating positive and negative charges. In contrast, the current model includes long-ranged Coulombic interactions between the ionic groups and also relatively short-ranged sticker−sticker and sticker−ion interactions. The prior coarse-grained work focused on ionomers with precisely spaced ions along the chain backbone; the structural regularity of these precise ionomers allows for interesting properties, and characterization of these in comparison to pseudorandom materials has improved our understanding of structure−morphology relationships in ionomeric materials. However, there is significant technological interest in the more common and commercially available partially neutralized random ionomers such as Surlyn, which is used in a range of applications such as in packaging and as the outer covering of many golf balls. In the current work, we thus consider random ionomers (which were also considered in parts of the prior work, but without carboxyl groups15,16) that are partially neutralized with sodium. To better map to PEMAA (as in Surlyn) rather than PEAA (as in prior work), the uncharged bead attached to the charged pendant bead should include an additional methyl group. This methyl group is small relative to our coarse-grained bead size, and we suggest this methyl group should simply be included within the CG bead to which the pendant group is bonded (which also contains another backbone carbon, for the same total of three carbons as in other uncharged backbone beads). Though commercial materials are typically branched and contain a range of molecular weights, for simplicity and ease of comparison with prior work, we consider linear, monodisperse, relatively low molecular weight systems. Below, after discussing the methodology, we report the effect of stickers on the local aggregate structure as well as the overall aggregate order for different chain architectures. We study the dynamics by calculating the mean-squared displacement (MSD) and the polymer end-toend vector and cluster autocorrelation functions. Finally, we



METHODS Ionomer Model. The coarse-grained (CG) ionomer model we employ is based on that of refs 15−18, with added sticker groups. Our simulation box contains linear bead−spring backbone beads similar to those employed by Kremer and Grest.31 Along these uncharged backbones, pendant anions are placed randomly. Free counterions are added to make the system charge neutral. The mapping to experimental system of PEMAA is as shown in Figure 1. Three backbone carbons map

Figure 1. Schematic showing (a) a PEMAA-based ionomer partially neutralized with sodium, with the acid groups placed randomly on the backbone, and (b) the corresponding model system; uncharged polymer beads are shown in gray, anionic pendant in red, counterion in blue, and sticker in green.

to one uncharged backbone bead, which is the same size as a pendant bead which maps to a COO− anion and is twice the size of the counterion which maps to Na+. Monomers along a backbone are bonded according to a nonlinear spring (FENE) potential, defined by ÄÅ É l 2Ñ o ÅÅÅ ÑÑÑ o i y r o j z 2 Å Ñ o o o−0.5kR 0 lnÅÅÅÅ1 − jjjj zzzz ÑÑÑÑ, r < R 0 R UFENE(r ) = m Å ÑÑÖ 0 k { o ÅÇ o o o o o r ≥ R0 o n∞ , The maximum bond length R0 is set to 1.5σ to avoid chain crossing,31 and k = 30ε/σ2 is the spring constant, where σ = 1.0 is the Lennard-Jones (LJ) diameter of a polymer bead (used as the unit of length) and ε = 1.0 is the LJ interaction strength of backbone monomers (used as the unit of energy). All beads interact via the truncated LJ potential given by Ä ÉÑ l o ÅÅÅÅi σij y12 i σij y6 ÑÑ o o j z j z o 4εijÅÅÅjj zz − jj zz + S ÑÑÑÑ, r ≤ rc o o ÅÅk r { ÑÑ ULJ, ij(r ) = m kr { ÑÖ o o ÅÇ o o o o r > rc n 0, For bonded beads, the cutoff is at the minimum, 21/6σ, while for nonbonded beads the attractive tail of the potential is included to a cutoff rc = 2.5σ, and the potential is shifted to 0 at the cutoff by a factor S (which depends on the cutoff distance). The polymer beads have a diameter of σij = 1.0σ, and counterions have diameter σij = 0.5σ; all beads have unit mass. The interaction strength εij = 1 for all interactions other than sticker−sticker and sticker−ion interactions discussed below. B

DOI: 10.1021/acs.macromol.7b02073 Macromolecules XXXX, XXX, XXX−XXX

Article

Macromolecules

We found appropriate values for εss, εsa, and εsc by performing short simulations with a range of values and comparing the aggregate structures to those previously obtained through atomistic simulations (containing explicit hydrogen bonds).8,9,37 We used precisely spaced pendant groups with the same average backbone spacer length between pendant groups (mapping three backbone C atoms to one CG bead) and the same levels of neutralization as used in the prior atomistic study; whether a pendant group was a sticker or an ion was chosen randomly according to the neutralization level. Specifically, the average spacer length is 3, 5, or 7 backbone beads per pendant bead (corresponding to 9, 15, and 21 backbone carbons), and we consider four levels of neutralization: 10%, 25%, 43%, and 75%. We adopted a naming convention to describe systems’ architectures, similar to that of prior work: system names start with r or p, where r refers to random and p to precise spacing of pendant groups, are followed by Nbb#- where the # refers to the average spacer length, and end with the percentage of sodium neutralization (e.g., for the system pNbb7-43%Na, pendant groups are placed precisely along the chain every 7 backbone beads and 43% of these pendant groups are charged, with the remainder being sticker groups). The appropriate number of counterions are also added to create a charge neutral system. To quantify the structural similarity of aggregates between the CG and atomistic models, we compared the height of the first peak in the sticker−sticker, sticker−anion, and sticker− counterion radial distribution functions, gss(r*), gsa(r*), and gsc(r*) where r* is the location of the first peak, and the r number of nearest neighbors n(rc,n), defined as 4π∫ 0c,nρg(r)r2 dr, up to a neighbor cutoff radius rc,n. For the atomistic simulations, the sticker and anion locations were defined to be at the center of mass of the atoms of the group (COOH or COO−). We used the first minimum in the CG g(r) as the cutoff radius to calculate n(r) for both the CG and atomistic systems; to map distance from CG to atomistic units, we used 1σ = 0.39 nm as this makes the CG and atomistic first g(r) peaks coincide. This was done for consistency in the comparison and because the atomistic simulations have multiple features just after the first g(r) peak (meaning that the appropriate neighbor cutoff radius to choose based on the atomistic g(r) would be unclear). While we report comparisons of CG and atomistic results for both g(r*) and n(rc,n) below, we consider that approximately reproducing the n(rc,n) values is more important, as the g(r) peak height depends on very local atomistic packing/ordering (detailed features around the first peak in g(r)) that are not expected to be reproduced by our CG model. In initial tests, discussed here to provide context, we considered that all three new interaction parameterssticker− sticker (εss), sticker−anion (εsa), and sticker−counterion (εsc)were the same for simplicity. All parametrization tests were performed as described in the Simulation Details section below except that the simulations were shorter; 50000τ each with structural data saved after 10000τ. On the basis of our mapping to real units and two hydrogen bonds’ interaction strength (approximately 14 kcal/mol)38 expected when two sticker groups come together to form a dimer, we briefly attempted to simulate at the very high LJ interaction parameter of ε = 22.0 (this single test was for the pNbb7-43%Na system only). The ion and sticker motion was slow enough that full equilibration would be very expensive or intractable (the anions’ MSD only reached 11σ2 by 50000τ), so the aggregate

The counterions and charged pendant beads experience a long-range Coulombic interaction of the form Uc(r ) =

q1q2 4πε0εrr

where q1 and q2 represent charges of +1e or −1e depending on ion type, ε0 is the vacuum permittivity, and εr is the background dielectric constant that accounts for the chemical nature of the atomistic system not represented by our neutral coarse-grained beads. We choose the effective dielectric constant εr = 4 and convert between electron charge and our reduced units following prior work in refs 15−18. This background dielectric reduces Coulombic interactions uniformly throughout the system, as the simplest approximate way of reproducing the effect of having the background of many partially charged atoms as in the atomistic system (e.g., each C and H on the polymer backbone would have a partial charge and these would surround ionic groups and partially screen their interactions).15−18,32 As ions are added, the dielectric constant of the resulting material also changes, though the effects of all ion−ion Coulomb interactions are explicitly included. Adding acid groups, which have somewhat more pronounced partial charge separation than the other polymer backbone groups, would further increase the static dielectric constant in a way that is not accounted for in our model. We expect this is a relatively small increase for the acid contents considered here and have chosen to keep a simplified model with a single, constant background dielectric. We map to real length units as in prior work by considering that σ = 0.4 nm, as this corresponds approximately to the size of three carbons along the backbone (one of our CG monomers) or a COO− ion (one of our pendant ions), which is also approximately the same size as a COOH group.33 The counterion diameter 0.5σ then corresponds to the size of an Na+ ion (approximately 0.2 nm16). This mapping means that to give εr = 4, we set the Coulombic interaction strength of 1/0.028 between unit charges separated by 1σ. Sticker Parametrization. Sticker groups representing COOH were modeled simply; stickers are the same as other uncharged monomers except with different, more favorable pairwise LJ interaction strengths with other stickers and with ionic groups. This adds three new parameters to our model versus that of prior work: the sticker−sticker, sticker−anion, and sticker−counterion LJ interaction strengths, or εss, εsa, and εsc. These are meant to account for the additional cohesive chemical interactions of these groups (primarily hydrogen bonding). While hydrogen bonding is strongly orientation dependent, for simplicity, all interactions are radially symmetric in the current model. Therefore, the interaction strengths chosen from the parametrization process described below should not be considered to represent the strength of the hydrogen bonds but simply the values that gave the desired overall properties for our model. There are many other options to model ionomers at finer or coarser length scales; for instance, Maranas and co-workers have developed useful coarse-grained models of PEO-based ionomer materials where the backbone is treated implicitly and only the ions are considered.34−36 We choose to model hydrogen-bonding groups in a simple way at the same level as monomer beads and with a few additional parameters, as this choice allows us to study the general physics of associating polymers while remaining comparable to prior ionomer model of Hall et al.17 C

DOI: 10.1021/acs.macromol.7b02073 Macromolecules XXXX, XXX, XXX−XXX

Article

Macromolecules

Table 1. Pairwise n(rc,n) and g(r*) Values from CG Simulations with the Chosen Interaction Parameters as Well as the Percentage Error between These Values and Those of Atomistic Simulations, Calculated as (Atomistic Value − CG Value)/ Atomistic Value (Percentages Are Rounded Off to the Nearest Whole Number); Values fromNbb = 5 Are Italicized To Emphasize Intermediate System Which Is Most Closely Matched no. of nearest neighbors, n(rc,n) Nbb = 7 Nbb = 5 Nbb = 3

first peak in radial distribution function, g(r*)

% error in n(rc,n)

% error in g(r*)

sticker−sticker

sticker−anion

sticker−cation

S−S

S−A

S−C

sticker−sticker

sticker−anion

sticker−cation

S−S

S−A

S−C

0.91 1.29 1.83

1.12 1.73 2.21

0.73 1.12 1.36

26 −3 −25

33 −1 −20

40 −6 −11

2.54 2.86 2.57

4.83 6.32 4.84

9.31 8.29 5.25

1 −36 −70

63 32 12

55 41 25

reasonable simulation time; further discussion can be found in ref 18. Periodic boundary conditions and a velocity-Verlet algorithm with a time step of 0.005τ are used in all our simulations. Each simulated system contains 800 ionomers, unless otherwise stated. We placed pendant groups fully randomly along the backbone, allowing any backbone bead to carry a pendant group (either charged or a sticker, representing COO− or COOH, chosen randomly according to the neutralization level). This is a closer representation of the commercially available ionomer Surlyn than our prior model (although unlike Surlyn, our chains are short and we do not include branching). All the systems have 35−36 backbone beads per chain, and this is lower than the entanglement limit for typical Kremer−Grest bead−spring chains.31,41 Following prior work, we report pendant group concentration in terms of the ratio of backbone beads per pendant bead, Nbb. To calculate the model’s MAA mole fraction (to consider acid content in the manner it is typically reported for Surlyn), one can use 1/ (((Nbb − 1) × 3/2) + 1) (the ratio of MAA group to the total groups of MAA and ethylene). For this conversion, we considered that the MAA group is represented by the pendant and the backbone bead to which it is connected, and each uncharged bead not connected to the pendant contains three −CH2 groups (3/2 ethylene groups). From this, the three MAA mole fractions considered here are 0.25, 0.143, and 0.1, corresponding to the average number of backbone beads per pendant bead of Nbb = 3, Nbb = 5, and Nbb = 7, respectively. The three systems with varying pendant mole fractions (equal to acid content before neutralization) have been analyzed at the same neutralization of 43%, and the Nbb = 3 system was analyzed at three additional neutralizations of 10, 25, and 75%. The chosen set of systems allows us to explore both changes in pendant content (which we also refer to as acid content) at the same neutralization and changes in neutralization levels for the same acid content. Each system was equilibrated for 2 × 105τ. Over this time, the mean-squared displacement of the polymers was greater than 3Rg2, and the polymer end-to-end vector autocorrelation (normalized to 1 at time 0) decayed to less than 1/e. The average density after the first 10000τ of an NPT run (during which time the density remained relatively constant) was used to set the density for the NVT ensemble runs. For collecting rheological data, the systems were equilibrated for 2 × 105τ in the NVT ensemble, then the thermostat was removed, and a brief further equilibration was run in the NVE ensemble for 20000τ; then data were collected in the NVE ensemble for 105τ. The slowest system (rNbb3-75%Na) took 240 h to run for 2 × 105τ on 48 processors. Relatedly, ref 9 details the challenges of simulating similar materials with a fully atomistic model; their work used a smaller number of polymers in the simulation box and advanced methods to help equilibrate the

structures are likely kinetically trapped on the time scale of our test, especially with regard to their longer-ranged ordering. Additionally, for this interaction parameter set, both the g(r*) and n(rc,n) were very different from those of the atomistic system. We found the simulation was more tractable with the lower values εss = εsa = εsc = 10.0, and the g(r*) and n(rc,n) also were relatively closer to the atomistic results. The large differences in gss(r*), gsa(r*), and gsc(r*) and corresponding differences in the three n(rc,n) values from atomistic simulations suggest that the interaction parameters should be different for the different interacting pairs. We varied all three parameters in increments of 0.5, searching for those parameters that led to a good matching of g(r*) and especially n(rc,n) values versus the atomistic results. After finding the values chosen below, we varied each parameter independently to ensure that this set of parameters gave a local minimum with respect to the sum of the squared error between the coarsegrained and atomistic nss(rc,n), nsa(rc,n), and nsc(rc,n) values for the middle-ranged ion content system pNbb5-43%Na; the error was also acceptable for the Nbb = 3 and 7 systems. Specifically, the chosen parameters had the lowest total weighted error, calculated by summing the squared errors in n(rc,n) and g(r*) values across all systems, after scaling the n(rc,n) error values by 20 to weight them more than g(r*) errors. When this chosen system is equilibrated for a longer period, as discussed below, the n(rc,n) and g(r*) values are constant within 4% of their values during the test simulation. A complete set of all interaction parameters tested and their errors can be found in the Supporting Information. We note that for the Nbb = 7 system there was better correlation with the atomistic structures with higher εsc interaction strengths. However, not only were the dynamics considerably slower for εsc ≥ 5.0 but also the errors in the n(rc,n) for Nbb = 3 and Nbb = 5 systems were higher, leading εsc = 4.5 to be the best overall choice given the criterion mentioned above. Table 1 shows the g(r*) and the n(rc,n) values, and the difference in these values from those of atomistic simulations, for the chosen LJ interaction strengths of εss = 2.0, εsa = 1.5, and εsc = 4.5. The stickers interact with uncharged monomers with ε = 1.0. Simulation Details. All simulations were carried out using LAMMPS molecular dynamics package.39,40 A Nosé−Hoover thermostat and barostat were employed to carry out simulations in the NPT ensemble, which was used to set the box size to simulate the system in NVT (with the same thermostat and no barostat). A constant energy simulation (NVE) was performed to collect stress−stress autocorrelation data by turning off the thermostat after equilibrating in NVT. Specifically, a temperature damping parameter of 1.0τ was used to maintain the temperature at 1.25T*, where the reduced temperature T* = 1 kT/ε, and a pressure damping parameter of 100τ was used to maintain the pressure at 0. This temperature choice allows the system to equilibrate in a D

DOI: 10.1021/acs.macromol.7b02073 Macromolecules XXXX, XXX, XXX−XXX

Article

Macromolecules system, and these simulations were not able to reach the diffusive regime.



RESULTS Structure. We analyzed ionic aggregates (clusters) via a similar cluster analysis as in prior work,15 except we included sticker groups and ions in the clusters. This is similar to the inclusion of oxygen−carboxyl hydrogen pairing in prior atomistic simulations’ analysis,8,9 which is important especially because hydrogen-bonded networks are formed at low levels of neutralization in these systems. To set the cutoff radius for clustering (particles within this distance are defined to be in the same cluster), we used 105% of the radius of the gij peak, r* (setting different distances for different pairwise interactions). Specifically, this was 0.79σ for counterion−sticker and counterion−anion pairs and 1.05σ for sticker−anion and sticker−sticker pairs (anion−anion and counterion−counterion interactions are not used to define clusters, as these are rarely the nearest neighbors within the cluster). We included only those clusters that have more than two counterions. Note that for the atomistic simulations of ref 9 the clustering cutoff was defined to capture the entire first peak in g(r), but these pairwise g(r) consider atom−atom interactions and contain more features than in our coarse-grained model. The cutoffs used here make the clusters appear distinct across the three acid concentrations for the CG model, as was found in the atomistic clustering results. All images were created using VMD.42 Many models of ionomer scattering predict the presence of spherical clusters with liquidlike dispersion;23,24 however, we find the actual cluster morphologies are more complex and depend on both neutralization level and acid content, in line with the conclusions drawn in prior atomistic work. Snapshots of clusters with different acid contents at 43% neutralization are shown in Figure 2. It is apparent that the density of the clusters decreases as the acid content decreases, and the morphology of the clusters also changes significantly. For the rNbb3-43%Na system with the highest mole fraction of acid groups, the clusters are long and branched. It is expected that as the neutralization level increases, these branched clusters will join and percolate through the simulation box. At 43% neutralization, however, the clusters are discrete in the Nbb = 3 system. The rNbb5-43%Na clusters are typically long and stringy, as seen in Figure 2b, whereas the clusters in the rNbb743%Na system are compact (although a few appear more extended/stringlike). At all three mole fractions, ions tend to closely pack and form the central regions of the clusters, while sticker groups decorate the outside of the clusters. The cluster morphologies are visually similar to those of prior atomistic simulations of the similar but precisely spaced material.8,9 To quantify the long-range order of the ionic aggregates, we calculated the ion−ion structure factor given by Sion − ion(k) =

1 N

N

Figure 2. Snapshots of all ionic aggregates (left) and a representative aggregate (right) with anions in red, counterions in blue, and stickers in green for (a) rNbb3-43%Na, (b) rNbb5-43%Na, and (c) rNbb7-43% Na; ions are colored by cluster with the largest 20 clusters being assigned unique colors and smaller clusters assigned different shades, so that the clusters can be seen as distinct.

N

∑ ∑ ⟨e−ik(r − r )⟩ i

j

Figure 3. Ion−ion structure factors for three systems with different average acid content at 43% neutralization, as labeled.

i=1 j=1

where N is the total number of ions in the system, ri,j are the vector positions of the ionic species, and k is the wave vector discretized based on the size of the simulation box. The sticker−sticker structure factor can be found in the Supporting Information. From Figure 3, we observe a sharp peak at low k for all three systems, signifying mesoscale ordering of ionic aggregates, as

established in prior work.16,17 There is a shift of the ionomer peak to lower wavevectors as the acid content decreases, which in agreement with prior experimental, CG, and atomistic results. This corresponds to an increase in the average separation of the aggregates in real space. The increase in E

DOI: 10.1021/acs.macromol.7b02073 Macromolecules XXXX, XXX, XXX−XXX

Article

Macromolecules

directions, but there are many ions that do not belong to this aggregate. From prior work it is known that at 100% neutralization with this ion content a single aggregate will percolate through the box in all three directions with most ions belonging to this largest aggregate.16,18 The structure factor of these Nbb = 3 systems with varying levels of neutralization is shown in Figure 5. (In order to

the peak intensity as the acid content decreases shown in Figure 3 is also in line with experimental10 and prior CG simulation results.16 The reverse trend in peak intensities observed in atomistic simulations was attributed to the small box size of those simulations, which restricted the sampling of longer ranged interaggregate distances that contribute to the ionomer peak.8 Figure 4 shows the aggregate structure of the highest acid content (Nbb = 3) system with different levels of neutralization.

Figure 5. Ion−ion structure factors for the rNbb3 systems at varying levels of neutralization, as labeled.

improve statistics for the lowest ion concentration system, the results shown for the rNbb3-10%Na system are from a simulation of 1600 polymers instead of 800 as in the other systems.) There is only a slight shift in the ionomer peak position to lower wavevectors as neutralization is increased from 10% to 75%, as the ion concentration increases. The peak intensity increases with neutralization, which is expected as neutralization increases the overall number of scatterers contributing to the ion−ion S(k) peak. (For reference, prior work16 considering the same architecture, but with fully repulsive LJ interactions and run at a lower density, reported the peak at 100% neutralization to have a height of 10 at 1.2k[σ−1], which would be a further increase from our 75% neutralized system.) These trends in structure factor are in agreement with prior atomistic simulation results.8 Dynamics. To show the ionic aggregate and polymer dynamics, we calculated the mean-squared displacement (MSD), polymer end-to-end autocorrelation function (ACFee), and the cluster autocorrelation function (ACFcluster, discussed further below). All dynamic data were calculated from frames logarithmically spaced in time (at time steps 20, 21, 23, ..., 220) starting every 220 steps, allowing us to observe the behavior at short and long times while reducing the number frames saved. All data shown are the average over four data series that are started 5200τ apart. The MSD is calculated after subtracting the small center-of-mass displacement of the overall system. The MSD of polymers’ centers of mass (COM) and counterions are shown in Figures 6 and 7 for varying acid contents and neutralization levels, respectively. MSDs of anions and stickers can be found in the Supporting Information. At very short times, both monomers and counterions move ballistically, with a log(MSD) vs log(τ) slope of 2. This ballistic motion is followed by a subdiffusive regime (where the log−

Figure 4. Snapshots of all ionic aggregates (above) and a representative aggregate (below) with anions in red, counterions in blue, and stickers in green for (a) rNbb3-10%Na, (b) rNbb3-25%Na, (c) rNbb3-43%Na, and (d) rNbb3-75%Na; ions are colored by cluster with the largest 20 clusters being assigned unique colors and smaller clusters assigned different shades, so that the clusters can be seen as distinct.

At a neutralization level of 10%, the clusters are small and stringy, often with just a few ions connected to other small groups of nearby ions by a thin string of stickers, as seen in the example cluster shown below the simulation box in Figure 4a. This is in line with results of prior atomistic simulations at low neutralizations, where aggregates propagate through hydrogen bonding.9 With an increase in the neutralization level to 25%, the clusters become more dense and longer but remain stringlike with little branching. Again, what would have been isolated ionic clusters (if only ions were used to define the clusters) are connected into larger aggregates via stickers as seen in Figure 4b. At higher neutralizations, the clusters more often have a continuous core of ions, and their morphology is often branched (Figure 4c). At 75% neutralization, some clusters become very large and branched as seen in the yellow aggregate in Figure 4d; this aggregate spans the box in two F

DOI: 10.1021/acs.macromol.7b02073 Macromolecules XXXX, XXX, XXX−XXX

Article

Macromolecules

Figure 6. MSD of the (a) polymer centers of mass and (b) counterions for three systems with different acid content at 43% neutralization, as labeled; the thin line at long times indicates a log− log slope of 1, and the horizontal dotted line in (a) indicates average Rg2 of the three systems; data are connected by thin lines, and only 3% of data points are marked with a symbol for clarity.

Figure 7. MSD of the (a) polymer centers of mass and (b) counterions for Nbb = 3 systems at various % neutralization, as labeled; the thin line at long times indicates a log−log slope of 1, and the horizontal dotted line in (a) indicates the average Rg2 of the four systems; data are connected by thin lines, and only 3% of data points are marked with a symbol for clarity.

log slope is less than 1). At very long times, all curves are expected to reach the Fickian diffusion regime with a slope of 1. Long time log−log slopes reported below are based on a least-squares regression of data from 1 × 105τ to 2 × 105τ (beyond the final time shown in Figure 6). Initially (at times