Topological Structure of Networks Formed from ... - ACS Publications

Oct 30, 2017 - ABSTRACT: Gels formed by coupling two different four-arm star polymers lead to polymer networks with high strength and low spatial hete...
0 downloads 3 Views 1MB Size
Article Cite This: Macromolecules XXXX, XXX, XXX−XXX

Topological Structure of Networks Formed from Symmetric FourArm Precursors Tzyy-Shyang Lin,† Rui Wang,† Jeremiah A. Johnson,‡ and Bradley D. Olsen*,† †

Department of Chemical Engineering and ‡Department of Chemistry, Massachusetts Institute of Technology, 77 Massachusetts Avenue, Cambridge, Massachusetts 02139, United States S Supporting Information *

ABSTRACT: Gels formed by coupling two different four-arm star polymers lead to polymer networks with high strength and low spatial heterogeneity. However, like all real polymer networks, these materials contain topological defects which affect their properties. In this study, kinetic graph theory and Monte Carlo simulation are used to investigate the structure and cyclic defects formed via A−B type end-linking of symmetric tetra-arm star polymer precursors. While loops constituting of odd number of junctions are forbidden by precursor chemistry, the amount and the correlation of secondary loops are found to increase with decreasing precursor concentration. This suppresses gelation, and the delay of gel point is quantitatively predicted by the topological simulations. Furthermore, comparison with network formed with asymmetric bifunctional−tetrafunctional precursors revealed that the behavior of loops consisting of 2n junctions in the symmetric system is analogous to the behavior of loops consisting of n junctions in the asymmetrical system, suggesting analogies between chemically distinct networks. and Temple,27 Stepto and co-workers,30,31 and Lang et al.33 to predict the formation of loops. Besides theoretical approaches, Monte Carlo simulations have also been extensively used to study loop kinetics. Eichinger and co-workers developed a dynamic Monte Carlo algorithm that captures the kinetics of loop formation for end-linked elastomers.34,35 Stepto and coworkers also developed a Monte Carlo algorithm that takes into account the effect of enhanced internal concentrations between reactive groups on the same molecule.32 In addition, Rankin et al. also studied the formation of gel networks with extensive cyclic defects with Monte Carlo methods.36 More recently, Schwenke et al. applied the bond-fluctuating model (BFM) to simulate the formation of homopolymer star networks and copolymer star networks.37 While these studies successfully capture many aspects of loop formation in end-linking gels, many of them are not directly comparable to experimental findings or rely on fitting parameters to quantitatively match experimental results. The recent invention of network disassembly spectroscopy (NDS) by Johnson and co-workers has enabled the quantitative characterization of loop defects in polymer networks without making any assumptions regarding network structure or dynamics,38−40 providing a powerful basis for testing theories of network topology. Inspired by the early work of Stepto,30,31

I. INTRODUCTION Polymer networks have been used in a wide range of applications, including superabsorbers, high impact rubbers, shape memory polymers, selective membranes, drug delivery devices, and tissue engineering scaffolds.1−15 A substantial number of these systems are end-linked polymer networks or polycondensation networks;16,17 the well-defined chemistry of these systems also makes them ideal models for understanding network physics. Despite extensive studies on these systems, much of our fundamental knowledge is still based on the ideal tree model, which neglects the effects of topological defects (loops) in the networks.18−22 Loops weaken the material by making part of the polymer network elastically ineffective in the sense that looped chains make a reduced contribution to the modulus of the polymer network.23,24 Loops also decrease the connectivity of the network, delaying the gel point (critical conversion).25−27 Because of the importance of these properties, methods to accurately predict the real topological structure of polymer networks and to quantify the effect of this structure on network properties are critical to the design of new materials. Theories and simulation methods have been developed to quantify the formation of cyclic defects. Tonelli and Helfand developed a chain statistics-based theory that predicts the formation of elastically ineffective loops during gelation.25,26 On the other hand, rate theory based models, where loop statistics are tracked by specifying the reactions and interconversions between different species, have also been adopted by Gordon © XXXX American Chemical Society

Received: August 24, 2017 Revised: October 30, 2017

A

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

Article

Macromolecules

concentrations for all experimental data to be compared in this work are beyond the dilute regime. While the chain statistics used in the simulation can be readily replaced by other chain statistics, this approximation has proven quantitatively accurate in models of previous gels.23,38−41 The kinetic graph theory used in this study is a modified rate theory based on the work of Stepto and co-workers,30−32 Zhou et al.,38 and Wang et al.41 Kinetic graph theory tracks the kinetics of the formation of infinitely large network through a set of differential equations describing interconversion between a series of finite subgraphs. Herein, the subgraphs are limited to a maximal size of two junctions, and junctions are assumed to be uncorrelated beyond the maximal size. Two junctions are required within the finite subgraph to describe all possible local connection structures that a tetra-arm precursor can have. Each functional group on the junction can be either dangling (unreacted), connected to the other junction within the junction pair, or connected outside of the finite subgraph. A major difference with the original KGT detailed in previous report38,41 is that for this more complex system the exact type of bridging connection between the junction and the underlying network needs to be explicitly taken into account. In the previous study, it is assumed that only one type of bridging connection to the network is possible, whereas in this work, the bridging connection to the network is categorized into several connection types. For a junction that has x connections to the underlying network, there are Px possible types of connection to the network, where Px is the number of unrestricted partitions of integer x. For example, a junction with three connections to the network has P3 = 3 types of connection, representing the case where each of the three groups on the junction connects to distinct junctions in the network, with two groups connecting to the same junction and the other connecting to another and with all groups connecting to the same junction. Given this difference in the type of connection beyond the two junction closure, there are a total of 15 possible two junction structures at full conversion; these species are illustrated in Figure 1c. A full accounting of all possible subgraphs is provided in the Supporting Information. For functional groups belonging to different molecules or different subgraphs, the rate of formation of the bridging connection (intermolecular reaction) is

Zhou et al. developed Monte Carlo simulations that quantitatively capture the formation of primary loops.38 Later, Wang et al. developed kinetic graph theory41 and a revised topological Monte Carlo simulation that does not require any fitting parameter.41,42 Their results show quantitative agreement with the results of NDS, and they show that a universal cyclic topology, depending only on one dimensionless parameter, exists for systems prepared via end-linking bifunctional (A2) polymer precursors and multifunctional (Bf) junctions. Using this information on defect densities, Zhong et al. developed a “real elastic network theory” (RENT) that corrects the classical phantom network theory for the presence of defects, quantitatively predicting the network modulus.23 In contrast to the A2 + Bf systems, the effect of primary loops can be totally eliminated by the clever design of polymer precursors. Gels formed through A−B type end-linking reactions between multifunctional precursors Ag + Bf with f, g > 2 can eliminate all odd order loops within the network structure. Sakai and co-workers have developed tetra-PEG gels formed from symmetric tetra-arm poly(ethylene glycol) (PEG) precursors43,44 as such model hydrogels for studying the fundamental physics of polymer networks. SANS studies on tetra-PEG gels showed that they can achieve a higher level of spatial homogeneity as compared to A2 + Bf gels.45−49 In addition, tetra-PEG gels exhibit higher breaking strength under compression than other hydrogels, with breaking strength comparable to native articular cartilage.44,50 Sakai and coworkers hypothesize that this is due to reduced network heterogeneity, and the structure of tetra-PEG polymer network was shown to be near ideal in some conditions.44,49,50 However, while it was initially proposed that the structure of tetra-PEG gels is a homogeneous diamond-like structure,44,49,50 later studies on local structure of tetra-PEG gels show that microscopic inhomogeneity exists, and the actual structure of the gel is far from the ideal diamond lattice.51,52 Beyond tetraPEG gels, Ossipov and Hilborn found that poly(vinyl alcohol) (PVA) gels formed from polyfunctionalized precursors have larger gel fractions and higher elastic moduli compared to the gel formed with bifunctional cross-linker.53 This work develops kinetic graph theory and Monte Carlo theory for multiarm star polymer precursors and applies these methods to characterize network topology and its impact on gelation and linear mechanics. While the majority of the analysis is focused on tetra-arm precursors, all the simulation and analytical methods presented in this study can be readily generalized to study precursors with arbitrary functionality. The topological structure, modulus, and gel point of the end-linked polymer networks are characterized to reveal the detailed structure of the network formed, and the results are compared to experimental data. Finally, the star polymer results are compared to asymmetric network systems, revealing important similarities.

II. SIMULATION METHOD The structure of networks formed via end-linking of stoichiometric amounts of symmetric tetra-arm precursors is studied with kinetic graph theory (KGT) and kinetic Monte Carlo (MC) methods. All reaction kinetics are assumed to be kinetically controlled without diffusion limitations, and the equal reactivity assumption is invoked. This corresponds to the kinetics found for tetra-PEG gelation by Nishi et al.54 Finally, the arms of the tetra-functional star polymers are modeled as monodisperse flexible Gaussian chains, since the initial polymer

Figure 1. Schematic illustration of KGT calculation. (a) Starting precursors, (b) intermediate products (see Supporting Information for all 72 subgraphs considered), and (c) terminal products that exist at full conversion. Solid lines indicate the cutoff of a connection to the underlying network. Multiple arms connected to the same line represents the connection of multiple arms to the same junction in the network. B

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

Article

Macromolecules R ij ,bridge = kABNA, iNB, jCiCj

Unless otherwise specified, MC simulation results are obtained from systems consisting of 10 000 B4 precursors and a stoichiometric amount of A precursors. In this study, at least 20 independent simulation runs are carried out for every condition, and the results shown are the ensemble average of the different runs. All MC results of average number of loops per junction have errors smaller than the size of symbols.

(1)

where kAB is the second-order rate constant, NA,i and NB,j are the number of functional groups A and B on subgraph i and j and Ci, Cj are the concentration of species (subgraph) i and j. This is equivalent to the more widely used form of the rate law which has the concentration of A group on subgraph type i and the concentration of B group on subgraph type j. The meansquare end-to-end distance between different functional groups on the same junction is given by ⟨R2⟩ = (2Ma/m)b2, where Ma is the molar mass of each arm of the star polymer, m is the molar mass of a Kuhn monomer, and b is the Kuhn length. The reaction rate of functional groups on the same molecule depends on the probability of two unreacted functional groups encountering each other. For two functional groups separated by n − 1 junctions, the probability of end-to-end interaction is Pn = (3/2πn⟨R2⟩)3/2; hence, the rate of two functional groups on the same subgraph forming an intermolecular loop of order n is 3/2 n n ⎞⎛ ⎛ NA, 3 ⎞ iNB, i R i ,loop = kAB⎜ ⎟⎜ ⎟ Ci ⎝ NAv ⎠⎝ 2πn⟨R2⟩ ⎠

III. RESULTS AND DISCUSSION Local Topological Structure of A4 + B4 End-Linking Network. The average number of loops per junction (including both A and B junctions) ni of secondary and higher-order loops (fourth- and sixth-order loops) obtained from KGT and MC simulation is shown in Figure 2. As for

(2)

where NnA,i and NnB,i are the number of functional groups A and B in subgraph i which can form nth-order loops when reacted, and NAv is Avogadro’s number. In this study, the maximal size is set to two junctions, which is the minimal size to model topological defects in a tetra-PEG network. Because network chemistry only permits even ordered loops in tetra-PEG gels, these are the only defects captured in the two junction approximation. At this level of approximation, there are a total of 1520 intermolecular reactions considered in the A4 + B4 KGT. Of these, 21 reactions are intramolecular (looping) reactions. To have a more accurate description of local structures and higher order loops, kinetic graph theory with larger maximum size would be needed; e.g., to account for fourth-order loops, the maximum graph size must be no less than four. The Monte Carlo simulations used in this study are based on the algorithm developed by Stepto 30−32 and adapted accordingly by Zhou et al.38 One A group and B group are selected sequentially for each reaction with probability

PAB =

1 V

+

(

⎡1 ∑ij ⎢ V + ⎣

3 2πn⟨R2⟩

(

Figure 2. Results of KGT and MC simulation. (a) KGT (lines) and MC (symbols) results of secondary loop fraction at different concentrations, (b) MC results for fourth- and sixth-order loop fraction, (c) KGT results for fraction of different local structures at different concentrations, and (d) comparison between the fraction local structures found by KGT (lines) and experimental results obtained by Lange et al.52

3/2

)

3 2πn⟨R2⟩

asymmetric systems,41,42 the topology of end-linked star polymer networks is uniquely determined by the dimensionless concentration cb3(M/m)3/2, where c is the initial concentration of functional groups and M = 2Ma. The dimensionless concentration represents the effective ratio of the molecular volume to the volume of solution per molecule. The number of secondary loops per junction monotonically decreases with increasing dimensionless concentration, and 1/n2 is found to be a linear function of dimensionless concentration for cb3(M/ m)3/2 > 0.5. This suggests that the concentration dependence of secondary loop fraction taking the form of the Langmuir isotherm

3/2 ⎤

)

⎥ ⎦

(3)

where V is system volume and the sum in the denominator is over all remaining unreacted functional group A, B pairs. The 1/V term represents the probability of a pair of reactive groups finding each other, while the latter term represents the effect of a locally enhanced concentration for intramolecular reactions, detailed by Stepto.30−32 The algorithm adopted in this study has been experimentally validated by comparison to experimental measurements of loop densities for asymmetric systems.23,38−42 In addition, reactions between any pair of functional groups are allowed to take place, not restricting to neighboring reactive groups.34,35 Furthermore, the purely topological nature of the algorithm greatly speeds up the simulation, allowing simulations on larger systems with improved accuracy.

n2 =

n2max 1 + kcb3(M /m)3/2

(4)

nmax 2

where = 1.5 is the maximal value of n2 in this network, corresponding to the formation of isolated spindle doublets (the species on the last row and last column of Figure 1c) at infinite dilution that has three secondary loops per two C

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

Article

Macromolecules junctions. Furthermore, it is observed that the number of higher order loops (fourth and sixth) is significantly less than the number of secondary loops over all dimensionless concentrations, while the numbers of fourth-order and sixthorder loops show peaks at intermediate value of the dimensionless concentration (and intermediate values of n2). With increasing dimensionless concentration, the structure of individual junctions transitions from dominantly associating with a single partner into a minimum sized cluster (spindles at dilute concentration) to associating with four distinct partner junctions (tree at high concentration), as shown in Figure 2c. The fraction of tree structures increased monotonically with increasing dimensionless concentration, while the fraction of spindle structures decreased monotonically. As this transition occurs, three different regimes of concentration are observed. In region 1 (cb3(M/m)3/2 < 5.50 × 10−5), where f tree/fspindle < 10−4, most of the products are spindles. In the infinite dilution limit, it is expected that all tetra-PEG molecules will form spindles, which can satisfy all covalent bonds with just two isolated molecules. Gel does not form in this regime; the network products are sol. In contrast, in region 3 (cb3(M/m)3/2 > 4.47), where f tree/fspindle > 104, tree species dominated and loops are rare. The system gets closer to the ideal tree-like structure limit. In the intermediate regime, loops are loosely correlated and separated by tree structures. In this regime, three other structures form in significant amounts, and their quantities peak in intermediate values of dimensionless concentration. With decreasing dimensionless concentration, the amount of triconnected structure peaks first, followed by the adjacent double loops and then contacted double loops structures. This shows that with decreasing dimensionless concentration not only the does the number of secondary loops increase but also the loops become correlated and secondary loops tend to form in vicinity of other loops. The predictions of the KGT and MC simulation match closely to experimental results by Lange et al. using 1H MQ NMR data to estimate the fraction of secondary loops,52 as shown in Figure 2d. This agreement illustrates the predictive power of these simulation methods in quantifying network topology. Notably, there is no fitting parameter involved in this comparison. The parameters for PEG used in relating experimental and simulation results are b = 0.76 nm and m = 119 g/mol, as measured by Kienberger et al.55 This same set of parameters is used for all experimental comparisons throughout this study. Comparison of data for the symmetric precursors (A4 + B4) and asymmetric bifunctional-tetrafunctional precursors (A2 + B4) studied by Wang et al.41,42 shows that there is a strong analogy between the 2n-th order loops in symmetric system and the n-th order loop in the asymmetric system. As observed, the behavior of secondary and higher order loops bears close resemblance to the behavior of primary and higher order loops in an asymmetric system, with the lowest order loop taking a Langmuir isotherm form and the number of higher order loops peaking at an intermediate number of lowest order loops. If the symmetric system with each arm of the four-arm star having N/ 2 Kuhn monomers is compared to a system with A2 strands of 2N monomers, the networks show a strong topological similarity. The number of primary loops normalized by the maximal number of primary loops in the asymmetric system is nearly identical to the normalized number of secondary loops in the symmetric system, as shown in Figure 3, since the relative reaction rate of looping to branching is identical in both

Figure 3. Comparison between normalized number of loops per junction n of symmetric and asymmetric systems calculated from MC simulations. The arms of precursors of symmetric system are consisted of N/2 Kuhn monomers, while the A2 strands consist of 2N monomers.

systems. This is despite the fact that the order of loops formed is not equal for both systems. The behavior of higher order loops is also similar for the two systems, as shown in the inset of Figure 3. This result shows that the secondary loops in the A4 + B4 system are essentially analogous to the primary loops in the A2 + B4 system, and the higher order 2n-th order loops in the symmetric system are analogous to the n-th order loops in the asymmetric system. This finding shows that the symmetric network is topologically analogous to the asymmetric networks studied previously. In the limit of large concentration, loops are rare and the network structure is essentially a local perturbation of the ideal tree structure. In the limiting case of infinite dilution, the network becomes a loop-free tree. These results also favor a random network structure over a diamond-like lattice structure for tetra-PEG gels, as the diamond lattice contains large amounts of sixth-order loops, inconsistent with the predictions of the theory and simulation. On the basis of these results, it is possible to imagine that different network chemistries may fall into sets that are topologically similar and that these topologically similar groups may in the future provide a basis for the classification and characterization of polymer networks. Gel Point of the A4 + B4 Network. Using the Monte Carlo simulations for gelation, the gel point can be predicted. It is well-known that loops in real networks lead to the delay of gelation. Various experimental studies have already demonstrated such effects in asymmetric systems,27−29 and rigorous theoretical studies on the effects of loops on gel point delay have also been conducted.22,56 Wang et al.56 show that by adopting the procedure proposed by Somvarsky and Dusek and Rankin et al.36,57 to estimate the real gel point of infinitely large system from finite size simulations, the MC simulation algorithm used in this study can accurately predict the gel point of A2 + B4 systems.56 Adopting the same procedure as Wang et al.,56 the gel point of each MC simulation run is extracted by locating the peak of the reduced weight-average degree of polymerization (DPrw, the weight-average degree of polymerization excluding the largest molecule in the simulation system). Several representative reduced weight-average degree of polymerization as a function of conversion p curves are plotted in Figure 4a. The position of D

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

Article

Macromolecules

Figure 4. Effect of loops on gel point. (a) Schematic illustration of gel point delay. (b) Extrapolation method for estimating the real gel point. The error bar shown is the standard error of the mean. (c) Dependence of gel point on concentration. The error bars of the prediction are the error estimated from χ2 fitting from (b) with error for each system size estimated from standard error of the mean of 20 simulation runs. Experimental data are extracted from Sakai et al.51

Elasticity of the Network. Using the loop densities from KGT and Monte Carlo simulation as inputs, the elasticity of the gel network can be predicted. Classical network theories (the affine network theory and the phantom network theory) are commonly used to link the network topology to the gel elasticity. While there are other models, such as the constraint junction models or slip-link models,59−61 that refine the prediction of network elasticity, they generally require information that is not captured in the simplistic simulation used in this study. Therefore, the following discussion will be focused on the prediction of the affine and phantom network theories. If the gel is assumed to follow affine network model, then all strands in loops with order 2 and larger are totally elastically effective. Since primary loops are completely eliminated by chemistry, the modulus of an affine network can be calculated as Gaffine = Gideal. On the other hand, if the polymer network follows the phantom network model, RENT23 may be used to predict the linear elastic modulus of the polymer network. RENT theory predicts that the elasticity of a phantom network is given by the equation

the peak is a function of the system size Ns, the total number of precursors initially placed in the simulation. Somvarsky and Dusek57 suggested that the real gel point for an infinitely large system should be extrapolated by fitting the observed gel point −1/2 , where pNc is the gel point extracted with pNc = p∞ c + aNs from finite system with size Ns (the system size is defined to be the number of B4 precursors initially placed in the simulation system), p∞ c is the gel point for infinitely large system, and a is a fitting parameter. As pointed out by Wang et al.,56 the Ns−1/2 dependence of gel point on system size originates from the percolation of the mean-field Bethe lattice. According to percolation theory, the largest cluster grows as N* ∼ |p − pc|−1/σ. However, for a finite size system, the system is observed to percolate when Ns ∼ N*.58 Equating the two relations, the observed percolation point, or the gel point, should be related to the real gel point for an infinitely large system in the form pNc −σ = p∞ c + bNs . For the ideal system, the Bethe lattice has an critical exponent of σ = 1/2, and the expression reduces into the original form used by Somvarsky and Dusek.57 For a network with loops, the critical exponent is not guaranteed to have the same value as that of the Bethe lattice. However, to avoid the error in nonlinear fitting, the size dependence of gel −1/2 to extrapolate and point is still fitted to pNc = p∞ c + aNs extract the gel point. Sensitivity analysis shows that for 0.3 < σ < 0.7 the extracted value of p∞ c is insensitive on the value of σ used (see Supporting Information). The dependence of the observed gel point on system size and the corresponding fit are shown in Figure 4b. Gel points obtained from the extrapolation are plotted as a function of dimensionless concentration in Figure 4c. The lower concentration regime corresponds to the cases with larger loop fractions, and the gel point is observed to be extensively delayed (higher critical conversion for gelation). The predicted gel point delay agrees well with the experimental measurements of Sakai et al.51 This demonstrates that the simulation methods used in this study are not only predictive in the local topologies of the gel network but also a powerful tool in predicting the global properties and chemistry of end-linking gel networks. This ability to predict the gel point without the need of fitting parameters is important to industrial and medical applications, where means of precise prediction of scorch time is essential.

⎞ 1⎛ G RENT f − 2⎛ 8(f − 1) 3 ⎞ ⎜1 − 2 n2, f ⎟ = ⎜1 − n2⎟ ≈ ⎝ Gideal f ⎝ 2 4 ⎠ f (f − 2) ⎠ (5)

where f = 4 is the functionality of the junctions. Since all odd order loops are vanishing and loops with order 4 and larger have a negligible effect on elasticity, only the correction term for secondary loops is retained. Although the topological connection of the tetra-arm system is analogous to the traditional asymmetric system, the elimination of primary loops and third-order loops has a dramatic effect on the elastic modulus of the system. By eliminating the odd-order loops through clever precursor design in tetra-PEG gels, the modulus of the network is dramatically increased. Figure 5 shows the predicted modulus of tetra-PEG network. The predicted elasticity is compared to the modulus of tetraPEG gel measured by Akagi et al.62 Nishi et al.54 have shown that over the range of precursor concentration used by Akagi et al. the diffusion of the reactants and the adjustment of chain configurations are much faster than the reaction rate during the end-linking formation of tetra-PEG gels. This justifies the direct comparison between the experimental data with the results of a kinetically controlled algorithm. Furthermore, while in many E

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

Article

Macromolecules

simulation. Direct comparison with NMR data validates the predictive power of the simulation methods to model network topology without any adjustable parameters. In addition, the delay of gel point is captured well with the pure topological simulation. It was shown that in the limit of very large precursor concentration (cb3(M/m)3/2 ≫ 1) the network is nearly loop free and behaves asymptotically tree-like; at moderately high precursor concentrations the network is a perturbation of the ideally branched network, and defects are generally isolated and uncorrelated. At lower concentrations, abundant loops exist and become correlated, strongly delaying the gel point. Finally, at very low concentration, the prevalence of strongly correlated cyclic defects prevents the formation of gel altogether. Careful comparison of the statistics of topological defects and the kinetics of the percolation process reveals that the behavior of lowest order loops and the relation between higher order loops and the lowest order loops are analogous in the tetra-arm symmetric A4 + B4 system and the asymmetric A2 + B4 system with corresponding polymer strand length. Finally, due to defect correlation, affine/phantom crossover, and the effects of unreacted dangling chains, RENT is shown not applicable for predicting the linear elastic modulus of tetra-PEG gels.

Figure 5. Comparison with elasticity calculated by loop-corrected affine network model and modified RENT and elasticity of A4 + B4 tetra-PEG gel with Mn = 5 kDa (◇), 10 kDa (○), and 20 kDa (△). Data are extracted from Akagi et al.62

model gels the effect of trapped entanglements are important,63−65 Akagi et al. showed that trapped entanglements are essentially absent in their tetra-PEG model network with tearing tests.62 This enables the direct comparison between the experimentally measured modulus and the prediction of RENT. In spite of the fact that the precursors used by Akagi et al. have large end-group functionality (X0 > 0.94), there will still be considerable amount of dangling ends present in the real experimental systems due to nonideal end-group functionality X0 < 1 and conversion p < 1. To determine the effect of unreacted functional groups, RENT is extended to account for the effects of dangling chains (see Supporting Information for derivation). For a system of conversion p, the modulus is predicted by the modified RENT to be G RENT 1 − p⎞ 1⎛ 3 ≈ ⎜1 − n2 − ⎟ Gideal 2⎝ 4 p ⎠



ASSOCIATED CONTENT

S Supporting Information *

The Supporting Information is available free of charge on the ACS Publications website at DOI: 10.1021/acs.macromol.7b01829. Details on KGT subgraphs, MC simulation algorithm, sensitivity analysis of p∞ c (σ), system size effects on ni, and sources of errors in the calculation of γ (PDF)



(6)

AUTHOR INFORMATION

Corresponding Author

*E-mail: [email protected] (B.D.O).

The predicted moduli at p = 0.8 and p = 0.9 are plotted in Figure 5. Compared to the experimentally measured value, RENT does not quantitatively capture the measured elastic modulus. While the qualitative trend for the effect of loops is correct, the quantitative value of the modulus is not predicted. The disagreement is likely to result from the failure of several assumptions used in the derivation of RENT. In RENT, it is assumed that the loops are not correlated, and all nodes are reacted. However, as shown in Figure 2c, in the crossover region from cb3(M/m)3/2 = 0.1 to 10, there are plenty of “contacted double loops” and a few “adjacent double loop” structures which contain secondary loops that are directly correlated with each other. These structures are not modeled in the current manifestation of RENT. In addition, the modified RENT predicts the modulus of a network with isolated defects, and the cooperative effect of unreacted dangling chains and loops has not been considered. Finally, it is believed that the gel makes a transition from phantom-like behavior to affine-like behavior as a function of increasing concentration,62 which is not captured by RENT theory. Clearly, these issues motivate the development of more advanced versions of RENT to incorporate these network effects.

ORCID

Jeremiah A. Johnson: 0000-0001-9157-6491 Bradley D. Olsen: 0000-0002-7272-7140 Notes

The authors declare no competing financial interest.

■ ■

ACKNOWLEDGMENTS This work was supported by the National Science Foundation (Award CHE-1629358). REFERENCES

(1) Flory, P. J. Principles of Polymer Chemistry; Cornell University Press: Ithaca, NY, 1953. (2) Colby, R. H.; Rubinstein, M. Polymer Physics; Oxford University Press: Oxford, 2003. (3) de Gennes, P. G. Scaling Concepts in Polymer Physics, 1st ed.; Cornell University Press: Ithaca, NY, 1979. (4) Lee, K. Y.; Mooney, D. J. Hydrogels for tissue engineering. Chem. Rev. 2001, 101, 1869−1880. (5) Zohuriaan-Mehr, M. J.; Kabiri, K. Superabsorbent polymer materials: a review. Iranian Polym. J. 2008, 17, 451. (6) Laftah, W. A.; Hashim, S.; Ibrahim, A. N. Polymer hydrogels: A review. Polym.-Plast. Technol. Eng. 2011, 50, 1475−1486. (7) Chen, J.; Park, H.; Park, K. Synthesis of superporous hydrogels: hydrogels with fast swelling and superabsorbent properties. J. Biomed. Mater. Res. 1999, 44, 53−62. (8) Ratna, D.; Karger-Kocsis, J. Recent advances in shape memory polymers and composites: a review. J. Mater. Sci. 2008, 43, 254−269.

IV. CONCLUSION This study shows that the structure of the polymer networks formed via end-linking of tetra-arm star polymer precursors can be modeled through kinetic graph theory and Monte Carlo F

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

Article

Macromolecules (9) McKeown, N. B.; Budd, P. M. Polymers of intrinsic microporosity (PIMs): organic materials for membrane separations, heterogeneous catalysis and hydrogen storage. Chem. Soc. Rev. 2006, 35, 675−683. (10) Hoare, T. R.; Kohane, D. S. Hydrogels in drug delivery: progress and challenges. Polymer 2008, 49, 1993−2007. (11) Dragan, E. S. Design and applications of interpenetrating polymer network hydrogels. A review. Chem. Eng. J. 2014, 243, 572− 590. (12) Ashley, G. W.; Henise, J.; Reid, R.; Santi, D. V. Hydrogel drug delivery system with predictable and tunable drug release and degradation rates. Proc. Natl. Acad. Sci. U. S. A. 2013, 110, 2318−2323. (13) Stuart, M. A. C.; Huck, W. T.; Genzer, J.; Müller, M.; Ober, C.; Stamm, M.; Sukhorukov, G. B.; Szleifer, I.; Tsukruk, V. V.; Urban, M.; Winnik, F.; Zauscher, S.; Luzinov, I.; Minko, S. Emerging applications of stimuli-responsive polymer materials. Nat. Mater. 2010, 9, 101−113. (14) Sun, J. Y.; Zhao, X.; Illeperuma, W. R.; Chaudhuri, O.; Oh, K. H.; Mooney, D. J.; Vlassak, J. J.; Suo, Z. Highly stretchable and tough hydrogels. Nature 2012, 489, 133−136. (15) Shin, H.; Jo, S.; Mikos, A. G. Biomimetic materials for tissue engineering. Biomaterials 2003, 24, 4353−4364. (16) Hild, G. Model networks based onendlinking’processes: synthesis, structure and properties. Prog. Polym. Sci. 1998, 23, 1019− 1149. (17) Mark, J. E. The use of model polymer networks to elucidate molecular aspects of rubberlike elasticity. Adv. Polym. Sci. 1982, 44, 1− 26. (18) Stockmayer, W. H. Theory of molecular size distribution and gel formation in branched-chain polymers. J. Chem. Phys. 1943, 11, 45− 55. (19) Stockmayer, W. H. Theory of molecular size distribution and gel formation in branched polymers II. General cross linking. J. Chem. Phys. 1944, 12, 125−131. (20) James, H. M.; Guth, E. Statistical thermodynamics of rubber elasticity. J. Chem. Phys. 1953, 21, 1039−1049. (21) Mark, J. E.; Sullivan, J. L. Model networks of end-linked polydimethylsiloxane chains. I. Comparisons between experimental and theoretical values of the elastic modulus and the equilibrium degree of swelling. J. Chem. Phys. 1977, 66, 1006−1011. (22) Suematsu, K. Theory of gel point in real polymer solutions. Eur. Phys. J. B 1998, 6, 93−100. (23) Zhong, M.; Wang, R.; Kawamoto, K.; Olsen, B. D.; Johnson, J. A. Quantifying the impact of molecular defects on polymer network elasticity. Science 2016, 353, 1264−1268. (24) Gordon, M.; Scantlebury, G. R. Statistical kinetics of polyesterification of adipic acid with pentaerythritol or trimethylol ethane. J. Chem. Soc. B 1967, 1−13. (25) Tonelli, A. E.; Helfand, E. Elastically ineffective crosslinks in rubbers. Rubber Chem. Technol. 1974, 47, 1116−1126. (26) Helfand, E.; Tonelli, A. E. Elastically ineffective polymer chains in rubbers. Macromolecules 1974, 7, 832−834. (27) Gordon, M.; Temple, W. B. Ring-chain competition kinetics in linear polymers. Makromol. Chem. 1972, 152, 277−289. (28) Ahmad, Z.; Stepto, R. F. T. Approximate theories of gelation. Colloid Polym. Sci. 1980, 258, 663−674. (29) Stepto, R. F. T. Dependence of the gel point on molecular structure and reaction conditions. Faraday Discuss. Chem. Soc. 1974, 57, 69−79. (30) Stanford, J. L.; Stepto, R. F. T. Rate theory of irreversible linear random polymrisation. Part 1.Basic theory. J. Chem. Soc., Faraday Trans. 1 1975, 71, 1292−1307. (31) Stanford, J. L.; Stepto, R. F. T.; Waywell, D. R. Rate theory of irreversible linear random polymerisation. Part 2.Application to intramolecular reaction in AA+ BB type polymerisations. J. Chem. Soc., Faraday Trans. 1 1975, 71, 1308−1326. (32) Rolfes, H.; Stepto, R. F. Network formation and properties: Rate theory description of effects of ring formation on elastic shear modulus of RA2+ RB3 networks. Makromol. Chem., Macromol. Symp. 1990, 40, 61−79.

(33) Lang, M.; Schwenke, K.; Sommer, J. U. Short cyclic structures in polymer model networks: A test of mean field approximation by Monte Carlo simulations. Macromolecules 2012, 45, 4886−4895. (34) Leung, Y.-K.; Eichinger, B. E. Computer simulation of endlinked elastomers. I. Trifunctional networks cured in the bulk. J. Chem. Phys. 1984, 80, 3877−3884. (35) Leung, Y.-K.; Eichinger, B. E. Computer simulation of endlinked elastomers. II. Bulk cured tetrafunctional networks. J. Chem. Phys. 1984, 80, 3885−3891. (36) Rankin, S. E.; Kasehagen, L. J.; McCormick, A. V.; Macosko, C. W. Dynamic Monte Carlo simulation of gelation with extensive cyclization. Macromolecules 2000, 33, 7639−7648. (37) Schwenke, K.; Lang, M.; Sommer, J. U. On the structure of star−polymer networks. Macromolecules 2011, 44, 9464−9472. (38) Zhou, H.; Woo, J.; Cok, A. M.; Wang, M.; Olsen, B. D.; Johnson, J. A. Counting primary loops in polymer gels. Proc. Natl. Acad. Sci. U. S. A. 2012, 109, 19119−19124. (39) Zhou, H.; Schön, E. M.; Wang, M.; Glassman, M. J.; Liu, J.; Zhong, M.; Diaz, D. D.; Olsen, B. D.; Johnson, J. A. Crossover experiments applied to network formation reactions: Improved strategies for counting elastically inactive molecular defects in PEG gels and hyperbranched polymers. J. Am. Chem. Soc. 2014, 136, 9464− 9470. (40) Kawamoto, K.; Zhong, M.; Wang, R.; Olsen, B. D.; Johnson, J. A. Loops versus branch functionality in model click hydrogels. Macromolecules 2015, 48, 8980−8988. (41) Wang, R.; Alexander-Katz, A.; Johnson, J. A.; Olsen, B. D. Universal cyclic topology in polymer networks. Phys. Rev. Lett. 2016, 116, 188302. (42) Wang, R.; Johnson, J. A.; Olsen, B. D. Odd−even effect of junction functionality on the topology and elasticity of polymer networks. Macromolecules 2017, 50, 2556−2564. (43) Sakai, T. Experimental verification of homogeneity in polymer gels. Polym. J. 2014, 46, 517−523. (44) Sakai, T.; Matsunaga, T.; Yamamoto, Y.; Ito, C.; Yoshida, R.; Suzuki, S.; Sasaki, N.; Shibayama, M.; Chung, U. I. Design and fabrication of a high-strength hydrogel with ideally homogeneous network structure from tetrahedron-like macromonomers. Macromolecules 2008, 41, 5379−5384. (45) Matsunaga, T.; Sakai, T.; Akagi, Y.; Chung, U. I.; Shibayama, M. Structure characterization of tetra-PEG gel by small-angle neutron scattering. Macromolecules 2009, 42, 1344−1351. (46) Matsunaga, T.; Sakai, T.; Akagi, Y.; Chung, U. I.; Shibayama, M. SANS and SLS studies on tetra-arm PEG gels in as-prepared and swollen states. Macromolecules 2009, 42, 6245−6252. (47) Saffer, E. M.; Lackey, M. A.; Griffin, D. M.; Kishore, S.; Tew, G. N.; Bhatia, S. R. SANS study of highly resilient poly (ethylene glycol) hydrogels. Soft Matter 2014, 10, 1905−1916. (48) Matsunaga, T.; Asai, H.; Akagi, Y.; Sakai, T.; Chung, U. I.; Shibayama, M. SANS studies on Tetra-PEG Gel under uniaxial deformation. Macromolecules 2011, 44, 1203−1210. (49) Asai, H.; Fujii, K.; Ueki, T.; Sakai, T.; Chung, U. I.; Watanabe, M.; Han, Y. S.; Kim, T. H.; Shibayama, M. Structural analysis of high performance ion-gel comprising tetra-peg network. Macromolecules 2012, 45, 3902−3909. (50) Sugimura, A.; Asai, M.; Matsunaga, T.; Akagi, Y.; Sakai, T.; Noguchi, H.; Shibayama, M. Mechanical properties of a polymer network of Tetra-PEG gel. Polym. J. 2013, 45, 300−306. (51) Sakai, T.; Katashima, T.; Matsushita, T.; Chung, U. I. Sol-gel transition behavior near critical concentration and connectivity. Polym. J. 2016, 48, 629−634. (52) Lange, F.; Schwenke, K.; Kurakazu, M.; Akagi, Y.; Chung, U. I.; Lang, M.; Sommer, J. U.; Sakai, T.; Saalwächter, K. Connectivity and structural defects in model hydrogels: A combined proton NMR and Monte Carlo simulation study. Macromolecules 2011, 44, 9666−9674. (53) Ossipov, D. A.; Hilborn, J. Poly (vinyl alcohol)-based hydrogels formed by “click chemistry. Macromolecules 2006, 39, 1709−1718. G

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

Article

Macromolecules (54) Nishi, K.; Fujii, K.; Katsumoto, Y.; Sakai, T.; Shibayama, M. Kinetic aspect on gelation mechanism of tetra-PEG hydrogel. Macromolecules 2014, 47, 3274−3281. (55) Kienberger, F.; Pastushenko, V. P.; Kada, G.; Gruber, H. J.; Riener, C.; Schindler, H.; Hinterdorfer, P. Static and dynamical properties of single poly (ethylene glycol) molecules investigated by force spectroscopy. Single Mol. 2000, 1, 123−128. (56) Wang, R.; Lin, T. S.; Johnson, J. A.; Olsen, B. D. Kinetic Monte Carlo Simulation for Quantification of the Gel Point of Polymer Networks. ACS Macro Lett. 2017, 6, 1414−1419. (57) Somvarsky, J.; Dusek, K. Kinetic Monte-Carlo Simulation of Network Formation. 2. Effect of System Size. Polym. Bull. 1994, 33, 377−384. (58) Christensen, K. Percolation Theory; Imperial College: London, 2002; Vol. 1. (59) Flory, P. J. Theory of elasticity of polymer networks. The effect of local constraints on junctions. J. Chem. Phys. 1977, 66, 5720−5729. (60) Edwards, S. F.; Vilgis, T. A. The tube model theory of rubber elasticity. Rep. Prog. Phys. 1988, 51, 243−297. (61) Vilgis, T. A.; Erman, B. Comparison of the constrained junction and the slip-link models of rubber elasticity. Macromolecules 1993, 26, 6657−6659. (62) Akagi, Y.; Gong, J. P.; Chung, U. I.; Sakai, T. Transition between phantom and affine network model observed in polymer gels with controlled network structure. Macromolecules 2013, 46, 1035−1040. (63) Campise, F.; Agudelo, D. C.; Acosta, R. H.; Villar, M. A.; Vallés, E. M.; Monti, G. A.; Vega, D. A. Contribution of Entanglements to Polymer Network Elasticity. Macromolecules 2017, 50, 2964−2972. (64) Lang, M.; Sommer, J. U. Analysis of entanglement length and segmental order parameter in polymer networks. Phys. Rev. Lett. 2010, 104, 177801. (65) Schlögl, S.; Trutschel, M. L.; Chassé, W.; Riess, G.; Saalwächter, K. Entanglement effects in elastomers: Macroscopic vs microscopic properties. Macromolecules 2014, 47, 2759−2773.

H

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