Numerical Estimates of the Topological Effects in the Elasticity of

Apr 19, 2019 - Andrei A. Gusev*. Institute of Polymers, Department of Materials, ETH Zürich, CH-8093 Zürich , Switzerland. Macromolecules , Article ...
0 downloads 0 Views 1MB Size
Download PDF
Article Cite This: Macromolecules XXXX, XXX, XXX−XXX

pubs.acs.org/Macromolecules

Numerical Estimates of the Topological Effects in the Elasticity of Gaussian Polymer Networks and Their Exact Theoretical Description Andrei A. Gusev*

Macromolecules Downloaded from pubs.acs.org by UNIV OF LOUISIANA AT LAFAYETTE on 04/19/19. For personal use only.

Institute of Polymers, Department of Materials, ETH Zürich, CH-8093 Zürich, Switzerland ABSTRACT: A general elastic-thread homogenization procedure is introduced and then used to obtain direct numerical estimates of the shear modulus of Gaussian polymer networks on the basis of their representative 3D periodic microstructures. A force relaxation algorithm is employed to find equilibrium positions of the network junctions. It is shown that from a statistical mechanics perspective, this numerical algorithm can simply be regarded as a convenient method for finding the mean (averaged out over thermal fluctuations) positions of the junctions that maximize the entropy of the polymer network and therefore define the observable, equilibrium thermodynamic state of the network. As an illustration, the procedure is used to obtain estimates of the shear modulus of various model 3D periodic end-linked and vulcanized microstructures incorporating up to a million strands. The same computer microstructures are used to extract several specific topological factors operative in the considered theoretical and model descriptions. It is shown that the numerical stiffness estimates always agree exactly with the predictions of the affine network theory. It is found that while the phantom network model is consistently more accurate than the affine network model, neither of these two classical models are really suitable for general purpose quantitative predictions. It is argued that for such reliable predictions, representative 3D network microstructures are generally needed in order to estimate the single key topological factor Γ advocated by the exact affine network theory.

■

constitutive relations25,26 that are required to study and design the large strain deformation behavior of advanced devices and structures from soft elastomeric materials. However, it still remains an open question if and to what extent these two classical models are valid even within their original theoretical framework of phantom Gaussian strands. In this work, representative three-dimensional periodic computer networks with phantom Gaussian strands are built and used to obtain direct numerical estimates of the networks’ shear modulus. For this, an elastic-thread homogenization procedure is used, in which each molecular network strand is replaced by a corresponding Gaussian elastic thread with the same stiffness coefficient.27 These same computer models are used to extract the topological factors operative in the considered classical solutions, and the resulting predictions are compared with the direct numerical estimates in order to assess the validity and the accuracy of these classical models and to demonstrate that it is the affine network theory (ANT) of Gaussian polymer networks that is exact.

INTRODUCTION Understanding the relationships between the topology and the mechanical properties of polymer networks and gels remains one of the biggest challenges in polymer science.1−9 The shear modulus of these systems can vary broadly, being of the order from 10−4 to 1 MPa. Such soft elastomeric materials and their composites are widely employed today not only in their traditional use areas such as automotive tires and various vibration damping and noise cutting applications but also in many other areas such as shape memory polymers,10 drug delivery systems,11,12 selective membranes and gas storage,13 water purification,14,15 tissue engineering scaffolds,16 stretchable electronics,17 soft robotics,18 and additive manufacturing.19 The entropic origin of the elasticity of polymer networks was first understood and described using a theoretical framework of fluctuating phantom Gaussian strands. This approach ignores the presence of the entanglements that can obstruct thermal motion of the neighboring network strands but it has allowed us to formulate and solve classical affine and phantom network models (ANM and PNM, respectively). Their resulting analytical solutions are now commonly employed to analyze and understand the results of various experimental and theoretical studies on the elasticity of polymers networks and gels; these models are used as a basis to develop more realistic polymer network models that account for the presence of entanglements20−24 and also to develop molecular-guided © XXXX American Chemical Society

Received: February 6, 2019 Revised: April 5, 2019

A

DOI: 10.1021/acs.macromol.9b00262 Macromolecules XXXX, XXX, XXX−XXX

Article

Macromolecules

■

THEORETICAL DESCRIPTIONS For an ideal Gaussian network strand of N Kuhn monomers of length b, one formally regards the entropy S(N,R) as a function of the exact value of the mean (averaged out over thermal fluctuations) strand end-to-end vector R and obtains that the entropy change upon deformation is given by20,21,23,24,27−29

The PNM is a widely used non-affine model that assumes that because of thermal fluctuations, the effective stretching of the network strands should differ from the macroscopic stretch applied to the whole network.22−24,27,33,34 The PNM shear modulus can be obtained by using the recurrence relation diagrams24,34 for the effective network strands and it is often written as

S(N , R′) − S(N , R)

G PNM = (νeff − μeff )kT

2 2 2 2 2 2 3 (λx − 1)R x + (λy − 1)R y + (λz − 1)R x =− k 2 Nb2

where μeff is the number density of the elastically effective crosslinks, defined as the ones that connect at least two elastically effective strands. It is then the main focus points of this study to introduce an elastic-thread force-relaxation homogenization procedure and to use it to obtain direct numerical estimates of the shear modulus of a number of 3D periodic Gaussian end-linked and vulcanized network microstructures. The same microstructures are employed to estimate the topological factor Γ and the number density of the elastically effective strands and crosslinks, νeff and μeff, respectively, and these estimates are used together to assess and compare the accuracy and physical validity of the above theoretical predictions 3−5.

(1)

where k is the Boltzmann constant, R = (Rx, Ry, Rz) and R′ = (λxRx, λyRy, λzRz) the mean end-to-end vectors before and after the deformation, respectively, and λx, λy, and λz the stretch factors. Note that in eq 1, a complete thermal averaging over strand conformations is also assumed both before and after the deformation. Consider a network of such Gaussian strands and assume that the network is affine (i.e., that every network strand experiences the same relative deformation). Then, the change in the entropy of the network upon deformation is the sum of the entropy changes of the individual strands, see eq 1. Supposing that the network is isotropic in its undeformed state and that there is no energy change upon deformation, one finds that the change in the specific Helmholtz free energy ΔF as a function of the stretch factors is determined by a single key topological factor Γ ΔF = −T ΔS = Γ=

R2 Nb2

■

END-LINKED NETWORKS To assess and compare the validity of the above theoretical predictions, we first consider a family of model Gaussian networks composed of equal length strands chemically endlinked by tri- or tetrafunctional crosslink agents. For all networks, a Kuhn length of b = 1 nm, a density of ρ = 1 g/cm3, and a molar mass of the Kuhn monomer of M0 = 200 g/mol are assumed, corresponding to a mean monomer number density of n = ρNA/M0 ≈ 3 nm−3, where NA is the Avogadro number. (For reference,24 for a widely used 1,4-polyisoprene rubber b = 0.84 nm, ρ = 0.83 g/cm3, and M0 = 120 g/mol, while for poly(dimethyl siloxane) rubber b = 1.3 nm, ρ = 0.90 g/cm3, and M0 = 380 g/mol.) Periodic cubic computer models with M = 105 strands and K = 2M/ϕ crosslinks are studied, where ϕ is the crosslink functionality defined as the number of chemically active sites per crosslink. The size of the models is determined as L = (MN/n)1/3, where N is the number of Kuhn monomers per network strand. The strands with up to N = 104 monomers are considered, leading to computer models with sizes of up to L = 1 μm. First, the crosslinks are randomly dispersed using a simple Monte Carlo procedure, see Figure 1a.

1 ΓνkT (λx 2 + λy 2 + λz 2 − 3) with 2 (2)

where ν is the number of strands per unit volume, T the temperature, ΔS the change in the specific entropy of the network upon deformation, and the angular brackets denote the ensemble average over all (not necessarily of equal length) network strands. Essentially, the same topological factor was discussed already by James and Guth in their 1947 paper30 in which it was defined in terms of the “vector-mean fractional extensions” of the strands in the undeformed network. But for some reason, this topological factor has not been considered in the contemporary research works on the elasticity of polymer networks and gels. Following standard text-book derivation,20,21,23,24 but without making the usual assumption that the strands are in ideal states with their R2 = Nb2, one finds that according to this above ANT, the shear modulus of an arbitrary isotropic Gaussian network is GANT = ΓνkT

(3)

The simplest way to approximate this topological factor Γ is to assume that in the undeformed state, the strands are in their ideal state where the relation Γ = ⟨R2/Nb2⟩ = 1 holds. This is the classical affine network model (ANM) originally proposed by Kuhn29,31,32 GANM = νeff kT

(5)

Figure 1. Generation of a periodic Gaussian network model with endlinked chains. For visual clarity, a relatively small model with a unit cell consisting of M = 5000 strands is considered, while all numerical calculations are performed using larger models with M = 105 strands. (a) Initial system with 2500 randomly situated tetrafunctional crosslinks (ϕ = 4). (b) Cross-linked network system with a reaction conversion rate of ξ = 99.5%. The network strands are shown by straight cylindrical threads.

(4)

where νeff is the number of elastically effective (active) strands per unit volume, which are defined as the ones that can store elastic energy upon network deformation. B

DOI: 10.1021/acs.macromol.9b00262 Macromolecules XXXX, XXX, XXX−XXX

Article

Macromolecules

a classical elastic-thread network model are exactly the mean positions of the fluctuating crosslinks in the comparable molecular network of ideal phantom Gaussian chains. To show this, they have first written the number of configurations (microstates) Ω available for a Gaussian polymer network with its junctions fixed at some given positions as ÅÄÅ ÑÉ ÅÅ 3 (R α − R β)2 ÑÑÑÑ Å Ω(R1, ..., R K ) = const × expÅÅÅ− ∑ ÑÑ ÅÅ 2 α > β Nαβb2 ÑÑÑÑ ÅÇ (7) Ö

The network strands are added sequentially, by selecting one of the remaining active crosslink sites at random and, with a given probability ξ, connecting the strand to this site and deactivating it, or leaving the site active otherwise. A set of all remaining active sites within a cut off distance of Nb (contour length of the strands) is formed and one of them is selected at random assuming a Gaussian end-to-end distance distribution with a mean-square value of Nb2. The second strand end is connected to the selected site with the same probability ξ and the site is deactivated, or it remains active otherwise. The procedure is continued until all M network strands have been considered, see Figure 1b. The resulting reaction conversion ratio is then simply given by ξ and the networks have on average (1 − ξ)2 free strands, 2ξ(1 − ξ) dangling strands, and ξ2 connected strands. However among those connected strands, there are still some elastically non-effective strands, for example, those that are involved in the various higher-order dangling (pending) and free (soluble) structures and non-loadbearing topological loops.1−9 As seen from the above, the end-linked networks studied in this work are generated by simply using phantom strands to connect randomly distributed crosslinks solely depending on their availability and distance. As compared to real end-linked polymer networks, no additional packing constraints, composition fluctuations, and mixing related problems are considered. These are among the reasons why such relatively high conversion ratios of up to 99.5% can be reached (though by allowing for some highly forced strands at the later stages of cross-linking when it becomes increasingly difficult to find a nearby active site), whereas in practice conversion ratios of up to only 95% are typically reached.35 As a further simplification, no network relaxation is allowed during cross-linking, while in practice both these processes occur simultaneously. Nevertheless, despite all these simplifications, such generated computer microstructures are fully suitable to be used to demonstrate below that as compared to direct numerical shear modulus estimates; the ANT formula 3 always gives exact predictions, whereas both classical ANM and PNM solutions 4 and 5 only approximate ones.

where vectors Rα and Rβ are the fixed positions of junctions α and β, respectively, Nαβ the number of Kuhn monomers in the strand joining these two junctions and the sum is over all strands of the network. Equation 7 merely states that the number of configurations available for a Gaussian network with its junctions fixed at some given positions is the product of the numbers of the configurations available for its individual network strands with their ends fixed at the same corresponding positions. Note that in eq 7 a complete thermal relaxation (averaging out over all possible conformations) of the individual strands is assumed for each given set of the positions of the junctions. This assumption is commonly used in theoretical analysis of the elasticity of Gaussian networks.3,24,36 It has then been proven by James and Guth that remarkably,27 for Gaussian polymer networks it is the mean (over thermal fluctuations) positions of the junctions that make Ω a maximum Ω max = Ω(R̅ 1, ..., R̅ K ) ÄÅ ÉÑ ÅÅ (R̅ α − R̅ β)2 ÑÑÑÑ ÅÅ 3 ÑÑ = const × expÅÅÅ− ∑ 2 ÑÑ ÅÅ 2 α > β N b ÑÑÖ αβ ÅÇ

where the overbar stands for the mean positions of the junctions. The equivalent thread network has the same connectivity as the original molecular chain network, with the αβ-strand of the original network replaced by a classical elastic thread with the force-elongation relation [cf. eq 6]

■

HOMOGENIZATION PROCEDURE A force relaxation method is used to find mean equilibrium positions of the crosslinks, by iteratively moving all of them simultaneously along the forces acting on them R α ← R α + ε fα,

fα =

fαβ =

∑ fαβ,

καβ =

3kT Nαβb2

3kT R αβ = καβ R αβ Nαβb2

(9)

The potential energy of the thread network is then given by

β

fαβ = καβ R αβ ,

(8)

V (R1, ..., R K ) =

(6)

where Rα and fα are the position vector of crosslink α and the force acting on it, respectively, Rαβ = Rβ − Rα and καβ the endto-end vector and the stiffness coefficient of strand αβ, fαβ the force exerted by crosslink β on crosslink α, and the summation is performed over all crosslinks β that are connected to crosslink α. It is assumed that the strands may have a variable number of monomers Nαβ. The relaxation parameter ε is selected to be large enough but still within the stability region of this simple explicit relaxation scheme. As a stopping criterion, a 108 times reduction of the first norm of the global force vector (f1, ..., fK) is used. It was established by James and Guth in their seminal 1943 paper27 that the equilibrium positions of the crosslinks in such

1 ∑ καβ(R α − R β)2 2 α>β

(10)

It can be easily seen that the equilibrium (minimum potential energy) positions of the junctions of this classical thread network are just the mean positions R̅ α of the junctions of the original molecular chain network, which becomes obvious at once when, using eqs 9 and 10, one rewrites eq 7 as É ÄÅ ÅÅ V (R1, ..., R K ) ÑÑÑ ÑÑ Å Å Ω(R1, ..., R K ) = const × expÅ− ÑÑ ÅÅ kT ÑÑÖ (11) ÅÇ Therefore, the condition that V shall be a minimum is just the condition that Ω shall be a maximum. The equilibrium positions of the junctions can be determined from C

DOI: 10.1021/acs.macromol.9b00262 Macromolecules XXXX, XXX, XXX−XXX

Article

Macromolecules fα = − =

minutes are typically needed to perform a single force relaxation run.

∂V (R1, ..., R K ) = − ∑ καβ(R α − R β) ∂R α α>β

■

RESULTS FOR END-LINKED NETWORKS The topological factor Γ of eq 2 can be written as

∑ καβ R αβ = 0 (12)

α>β

M

and hence, the iterative force relaxation method 6 can conveniently be used to find the mean positions of the junctions that make Ω a maximum for the original polymer network. This classical reduction is exact for any molecular Gaussian network.27,37 One can also take an alternative perspective and formally regard the network entropy as a function of the exact positions of the junctions28,38

Γ=

(R α − R β)2 3 k∑ 2 α>β Nαβb2

(13)

where a complete thermal relaxation of the individual strands (averaging out over all possible conformations) is assumed for each given set of the positions of the junctions. Since logarithm is a monotonous function of its argument, the condition that Ω shall be a maximum is also the condition that S shall be a maximum, and hence [cf. eq 6] ∂S(R1, ..., R K ) 3k = −∑ (R α − R β ) 2 ∂R α N αβb α>β =

1 ∑ καβ R αβ = 1 fα = 0 T α>β T

(14)

Thus, from a statistical mechanics perspective, one can view the force relaxation algorithm 6 simply as a convenient numerical method for finding the mean positions of the junctions that maximize the entropy of the polymer network and therefore define the observable, equilibrium thermodynamic state of the network [cf. eq 8]. After initial force relaxation, six pure shear deformations are applied to every relaxed structure, two in each of the three orthogonal coordinate planes. These volume-preserving deformations are imposed by changing two of the three box dimensions while keeping the remaining one unchanged. For example, in the xy-plane the two deformations are prescribed by putting λx = λ, λy = λ−1, λz = 1 and λx = λ−1, λy = λ, λz = 1, respectively. After six force relaxation runs, the resulting deformation free energies are evaluated and their average value ⟨ΔF⟩ is used to obtain an estimate of the structure’s shear modulus as [cf. eq 2] Gnum = 2⟨ΔF ⟩/(λ 2 + λ−2 − 2)

(16)

where subscript τ labels the strands of the network, Rτ is the mean (over thermal fluctuations) end-to-end vector of strand τ, and Nτ the number of Kuhn monomers in it. Since the mean end-to-end vector Rτ of the dangling and free strands is zero, these strands do not contribute to the topological factor Γ, and hence, they do not contribute the network shear modulus GANT, see eq 3. For the simplest topological loops consisting of a single network strand with both ends attached to the same junction, the mean end-to-end vector Rτ is obviously zero and such trivial loops also do not contribute to the network shear modulus. Equation 16 states that the contribution of any network strand (including those involved in higher-order loops with multiple strands and junctions) is completely determined by the square of its mean end-to-end vector Rτ, which suggests a practical approach to quantify the “elastic effectiveness” of network strands. For example,3 given a detailed molecular model of a network, one can use molecular dynamics or Monte Carlo sampling to evaluate the mean end-to-end vectors of the individual strands. Those with zero values of Rτ can then be classified as elastically ineffective, whereas the remaining strands as elastically effective with their quantitatively defined contributions to the network shear modulus as given by eq 16. This allows for rational topological design of ideal Gaussian polymer networks and can also be used to approach approximate design of non-Gaussian networks and help establish link to experimental NMR studies of real polymer networks.39 In this work, a coarse-grained elastic-thread representation is employed to study ideal Gaussian networks. After force relaxation, the resulting positions of the junctions of such networks exactly correspond to the mean (over thermal fluctuations) positions of the junctions in the comparable molecular networks,27 and these positions are used to calculate the mean strand end-to-end vectors Rτ and then the corresponding strand contributions to the topological factor Γ as defined by eq 16. Figure 2 shows the distribution of these strand contributions for one of the studied end-linked networks. One can see that the distribution is relatively broad and that practically all network strands are in non-ideal states with their R2/Nb2 < 1. The topological factor Γ of this specific force-relaxed equilibrium network has a value of Γ ≈ 0.12, which is significantly below unity, whereas initially before force relaxation it had a value of Γ ≈ 0.97. This clearly demonstrates the significance of force relaxation for accurate shear modulus predictions. Figure 2 shows that for this specific studied network, about a quarter of all network strands have zero end-to-end vectors. (This fraction can actually vary strongly depending on ϕ and ξ.) These strands do not store elastic energy upon network deformation, and hence, they are elastically non-effective. The remaining strands are elastically effective with a broad distribution of the elastic effectiveness that is quantified by the corresponding strand contributions to Γ as defined by eq

S(R1, ..., R K ) = k ln Ω(R1, ..., R K ) = k ln(const) −

R 2 1 ∑ τ2 M τ = 1 Nτb

(15)

A stretch factor of λ = 1.1 is used in the calculations but it was checked that the resulting estimates of Gnum are independent of λ. Since a complete thermal relaxation (averaging out) is assumed for the network strands and junctions [see eqs 7 and 8], the so-obtained shear modulus corresponds to the equilibrium shear modulus of the network at long times (or low frequencies). Numerical simulations are performed with in-house C code compiled and run on the ETH Zürich mainframe Euler cluster (https://scicomp.ethz.ch/wiki/Euler) using Hewlett-Packard 1215 m710x compute nodes equipped with quad-core Intel Xeon E3-1285Lv5 processors. Several (single processor) CPU D

DOI: 10.1021/acs.macromol.9b00262 Macromolecules XXXX, XXX, XXX−XXX

Article

Macromolecules fα =

∑ fαβ = 0 (17)

β

Suppose an affine transformation is applied to the system (both to the periodic cell and the crosslink coordinates, see Figure 1) and let Λ be the related position-independent transformation matrix. The crosslink positions are transformed as R′α = ΛRα, where the prime refers to the positions after the deformation. It then follows immediately that for an arbitrary equilibrium Gaussian polymer network, any affine transformation does not perturb the total forces acting on the crosslinks, that is, they all remain zero. Indeed, using eqs 6 and 17, one finds Figure 2. Distribution of the elastic effectiveness of the strands of a force-relaxed network with ν = 0.5 nm−3, ϕ = 3, and ξ = 90%. The histogram gives the probability density function and the solid squares the cumulative distribution function. For this network, Γ = ⟨R2/Nb2⟩ ≈ 0.120. The open circle shows the combined fraction of the dangling and free strands (i.e., 1 − ξ2 = 0.19 as discussed in the paragraph after Figure 1), while the open square shows the total fraction of the strands with a zero mean end-to-end vector. This total fraction (≈0.25) includes all elastically non-effective strands including those involved in the various network dangling and free structures and loops. It is this total fraction that is used in this work to define the fraction of elastically effective strands (νeff/ν ≈ 0.75 for this considered network).

f ′α =

∑ καβ R′αβ = ∑ καβ ΛR αβ = Λ∑ καβ R αβ = Λfα β

=0

β

β

(18)

which proves that the mean positions of the fluctuating crosslink junctions of equilibrium Gaussian networks deform affinely. To my knowledge, this important result was already first established by James and Guth in their seminal 1943 paper.27 Figure 3b shows that it is the topological factor Γ that is key to the accurate prediction of the shear modulus. As followed from eqs 2 and 16, this factor does not receive any contribution from the elastically non-effective strands that after force relaxation all have a zero end-to-end distance. And the ANM correctly discards these strands by considering only the elastically effective strands, see eq 4. However, the ANM formula 4 still greatly overestimates the network shear modulus, see Figure 4a. The underlying reason is that the principal ANM assumption that all elastically effective strands are in ideal states with their R2/Nb2 = 1 is simply unjustified. In reality, there exists a broad distribution and for the vast majority of the strands of the studied networks, their actual R2/ Nb2 contributions to Γ are significantly smaller than 1, as illustrated by Figure 2. Figure 4b shows that the PNM predictions are more accurate than the ANM predictions. Interestingly, for the highly cross-linked networks with a conversion rate of ξ = 99.5%, the PNM predictions 5 can become near exact.

16. These elastically effective strands do store elastic energy upon deformation and it is their number density νeff that determines the ANM shear modulus by eq 4. It is also these same elastically effective strands that are used to define the density of elastically effective crosslinks, μeff, required for the PNM predictions, see eq 5. Figure 3a shows that for the shear modulus of the many different Gaussian networks studied, the agreement between the ANT predictions and the direct numerical estimates is nothing but exact. This exact agreement evidently requires physical explanation, which is provided by the remarkable fact that the Gaussian polymer networks always deform affinely. Indeed, after force relaxation, the total forces acting on each of the crosslinks are zero [cf. eq 6]

Figure 3. Model polymer networks with end-linked strands. (a) After force relaxation, topological factor Γ is estimated based on the equilibrium positions of the crosslinks and eq 3 is used to obtain the network shear modulus GANT. Each symbol corresponds to a specific network with its own combination of ϕ, ξ, and ν = n/N. The resulting shear moduli cover a broad range representative of all commonly used polymer networks and gels. The dash line shows exact agreement, GANT = Gnum. (b) Topological factors for the networks with ϕ = 3 and ξ = 90%. After force relaxation, the elastically effective strands have a non-zero end-to-end distance, while the ends of elastically non-effective strands come together to the same positions. One can see that both νeff/ν and μeff/ν are noticeably smaller than their respective limiting theoretical values of 1 and 2/ϕ that are representative of the defect free networks. E

DOI: 10.1021/acs.macromol.9b00262 Macromolecules XXXX, XXX, XXX−XXX

Article

Macromolecules

Figure 4. Shear modulus predictions for the end-linked polymer networks. The dash lines show exact agreement. (a) ANM. (b) PNM.

Figure 5. Vulcanized networks. (a) Solid triangles are the ANM and PNM predictions. Solid circles are the ANT predictions while the dash line corresponds to exact agreement. (b) Topological factors. The horizontal dash line drawn at 0.5 is a guide for eyes to emphasize that in these studied networks, practically all their crosslinks are elastically effective.

repeatedly attempting to place each of them at a randomly selected grid site, until a site is found that has at least two remaining nonterminal uncross-linked monomers. Two of such monomers are selected at random and permanently connected, thus creating a junction linking together four network strands, and the placement of the next crosslink is performed until all the crosslinks are placed. The grid is then discarded and the above force relaxation method is used to find equilibrium offgrid positions of the junctions and to determine the shear modulus of the network using eq 15. Note that there are at most only two dangling and no free strands in such vulcanized networks so their shear modulus is solely determined by some other topological structures such as, for example, different order topological loops. Figure 5a shows that the ANT predictions are again exact while the ANM predictions are significantly larger (by ca. 250%) than the numerical estimates. Figure 5b indicates that the reason for this poor performance is the fact that the fraction of the elastically effective strands does not constitute any accurate approximation to the key topological factor Γ. As for the PNM predictions, similar to the polymer networks with end-linked strands studied above (see Figure 4), they are much more accurate than the ANM ones, being only 10% larger than the numerical estimates. Figure 5b shows that for these studied vulcanized networks, the fraction of the elastically effective strands is rather close to νeff/ν = 1, while the fraction of the elastically effective crosslinks is practically indistinguishable from μeff/ν = 0.5. As a consequence, for such networks without dangling and free strands some first quick PNM estimates (accurate to about 30%) can already be

Classical PNM formula 5 is derived assuming a defect-free tree-like topological structure,22,24 and our results show that for the highly cross-linked networks studied this assumption is indeed appropriate. However, for imperfect networks, the PNM predictions 5 are less accurate (see Figure 4b) and they overestimates the network shear modulus by up to an order of magnitude. It has been shown in recent works3,6,36 that various topological defects such as dangling strands and loops, which commonly occur in imperfect networks, can be incorporated in the underlying PNM tree-like structure and so the effect of these defects on the shear modulus of the network can be theoretically quantified by using the real elastic network theory.6,36 It would be interesting to see if and how one could use our detailed 3D periodic computer microstructures to assess the physical validity and the accuracy of these recent theoretical results. This will be a topic of our separate future work.

■

VULCANIZED NETWORKS To further assess the accuracy of the considered theoretical approaches, model polymer networks obtained by intrachain cross-linking (vulcanization) are studied. A single phantom chain with N = 106 Kuhn monomers is randomly grown on a periodic simple cubic grid with a spacing of b = 1 nm. The size of the grid is selected as L = 3 N /n , where, as before, a mean monomer number density of n = 3 nm−3 is assumed. At each growth step, one of the six nearest neighbor sites is selected using periodic boundary conditions. After the growth procedure is completed, a required number of tetrafunctional crosslinks (ϕ = 4) are sequentially placed on the grid, by F

DOI: 10.1021/acs.macromol.9b00262 Macromolecules XXXX, XXX, XXX−XXX

Article

Macromolecules

Brownian dynamics, or elastic-thread studies), there seems to be no obvious reason for not adopting the exact ANT. These two classical models are widely used today to estimate the shear modulus of real, non-Gaussian polymer networks and gels and to relate its measured values to the topology of the networks, they are frequently used as a basis for more advanced theoretical descriptions of the elasticity of entangled polymer networks, and they are also often employed as a starting point to develop molecular-guided constitutive relations for the simulation and design of the large strain deformation responses of advanced engineering devices and structures from such elastomeric polymer materials. The ANT advocates that the shear modulus of an arbitrary fluctuating Gaussian polymer network is completely determined by a single key topological factor Γ and it is suggestive to also consider this exact theory as a possible, firm theoretical basis to approach the above problems. Alternatively, the elastic-thread numerical method employed in this work can in its own right be adopted and directly used to study these problems.

obtained by simply using these limiting theoretical values of 1 and 2/ϕ = 0.5, respectively, with no need for any further detailed topological analysis. However, for accurate quantitative predictions, one should still use representative forcerelaxed 3D network models to estimate the key topological factor Γ advocated by the exact ANT.

■

CONCLUSIONS To conclude, a force relaxation homogenization procedure has been introduced and then used to obtain direct numerical estimates of the shear modulus of various model, both endlinked and vulcanized networks of ideal Gaussian chains. For this, an exact classical analogue to the fluctuating network of molecular Gaussian strands is used in which each strand is replaced by a corresponding elastic thread with the same elastic stiffness coefficient. It is found that the resulting numerical stiffness estimates always agree exactly with the predictions of the ANT but not with the predictions of classical ANM and PNM. It has been seen, in line with recent literature studies,3,6,7 that while it is relatively easy and straightforward to identify and discard the contributions of the dangling and free network strands to the overall network shear modulus, it is quantitatively nontrivial to determine the contributions of the remaining strands. It has been shown in this work that for any force-relaxed Gaussian polymer network, it is the ANT key topological factor Γ that quantitatively captures the elastic effectiveness of the network strands, which can help study and understand the importance of the various specific topological features of the network (such as, e.g., different order loops and various dangling structures). In this work, force-relaxed 3D periodic microstructures are used to estimate the mean (averaged out over thermal fluctuations) end-to-end vectors of the strands that are required to calculate the ANT topological factor Γ. As stated by eq 16, this factor is defined by a sum of the contributions of the individual strands, which may generally depend on the local surroundings of the strands. Here, a multibody optimization problem is solved numerically (using the force relaxation method) to estimate the elastic effectiveness of the individual strands as directly determined by their contributions to factor Γ. It is then an interesting problem to identify those specific topological features of the surroundings that control the elastic effectiveness of the strands and to also relate the obtained shear modulus estimates to the underlying chemical parameters such as the average molecular weight and polydispersity of the strands, junction functionality, and conversion of the cross-linking reaction. One can as well use our force-relaxed 3D periodic microstructures to systematically assess the physical validity and the accuracy of the branchingtheory shear modulus prediction40 and also recent theoretical solutions for the elastic effectiveness of the strands near isolated topological defects,3,6,36 which were obtained by including isolated dangling strands and loops into the underlying tree-like PNM topological structure.24,34 All this remains to be elaborated in a future work. For ideal Gaussian networks, classical ANM and PNM models are nothing but approximations to the exact ANT. Their use is justified in the situations where a quick estimate of the networks’ shear modulus is needed but there is not enough topological information to estimate the required key factor Γ. However, in the situations where such topological information can be available (e.g., in molecular dynamics, Monte Carlo,

■

AUTHOR INFORMATION

Corresponding Author

*E-mail: [email protected]. ORCID

Andrei A. Gusev: 0000-0002-7075-8111 Notes

The author declares no competing financial interest.

■

ACKNOWLEDGMENTS The author is thankful to I. Tsimouri for her help with the visualization of the computer models.

■

REFERENCES

(1) Wang, J.; Lin, T.-S.; Gu, Y.; Wang, R.; Olsen, B. D.; Johnson, J. A. Counting Secondary Loops Is Required for Accurate Prediction of End-Linked Polymer Network Elasticity. ACS Macro Lett. 2018, 7, 244−249. (2) Gu, Y.; Schauenburg, D.; Bode, J. W.; Johnson, J. A. Leaving Groups as Traceless Topological Modifiers for the Synthesis of Topologically Isomeric Polymer Networks. J. Am. Chem. Soc. 2018, 140, 14033−14037. (3) Lang, M. Elasticity of Phantom Model Networks with Cyclic Defects. ACS Macro Lett. 2018, 7, 536−539. (4) Lin, T.-S.; Wang, R.; Johnson, J. A.; Olsen, B. D. Topological Structure of Networks Formed from Symmetric Four-Arm Precursors. Macromolecules 2018, 51, 1224−1231. (5) 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. (6) 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. (7) Wang, R.; Alexander-Katz, A.; Johnson, J. A.; Olsen, B. D. Universal Cyclic Topology in Polymer Networks. Phys. Rev. Lett. 2016, 116, 188302. (8) 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. (9) 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. (10) Hu, J.; Zhu, Y.; Huang, H.; Lu, J. Recent advances in shapememory polymers: Structure, mechanism, functionality, modeling and applications. Prog. Polym. Sci. 2012, 37, 1720−1763.

G

DOI: 10.1021/acs.macromol.9b00262 Macromolecules XXXX, XXX, XXX−XXX

Article

Macromolecules

(38) Landau, L. D.; Lifshitz, E. M. Statistical Physics, 3rd ed.; Butterworth-Heinemann Publications: China, 2009; Part 1, Volume 5 of Course of Theoretical Physics, Chapter XII: Fluctuations, pp 333− 396. (39) Lang, M. Monomer Fluctuations and the Distribution of Residual Bond Orientations in Polymer Networks. Macromolecules 2013, 46, 9782−9797. (40) Miller, D. R.; Macosko, C. W. A New Derivation of Post Gel Properties of Network Polymers. Macromolecules 1976, 9, 206−211.

(11) Hoare, T. R.; Kohane, D. S. Hydrogels in drug delivery: Progress and challenges. Polymer 2008, 49, 1993−2007. (12) Dragan, E. S. Design and applications of interpenetrating polymer network hydrogels. A review. Chem. Eng. J. 2014, 243, 572− 590. (13) 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. (14) Shannon, M. A.; Bohn, P. W.; Elimelech, M.; Georgiadis, J. G.; Mariñas, B. J.; Mayes, A. M. Science and technology for water purification in the coming decades. Nature 2008, 452, 301−310. (15) Geise, G. M.; Lee, H.-S.; Miller, D. J.; Freeman, B. D.; McGrath, J. E.; Paul, D. R. Water Purification by Membranes: The Role of Polymer Science. J. Polym. Sci., Part B: Polym. Phys. 2010, 48, 1685−1718. (16) Sun, J.-Y.; Zhao, X.; Illeperuma, W. R. K.; Chaudhuri, O.; Oh, K. H.; Mooney, D. J.; Vlassak, J. J.; Suo, Z. Highly stretchable and tough hydrogels. Nature 2012, 489, 133−136. (17) Rogers, J. A.; Someya, T.; Huang, Y. Materials and Mechanics for Stretchable Electronics. Science 2010, 327, 1603−1607. (18) Rus, D.; Tolley, M. T. Design, fabrication and control of soft robots. Nature 2015, 521, 467−475. (19) Jungst, T.; Smolan, W.; Schacht, K.; Scheibel, T.; Groll, J. Strategies and Molecular Design Criteria for 3D Printable Hydrogels. Chem. Rev. 2016, 116, 1496−1539. (20) Doi, M. Introduction to Polymer Physics; Oxford University Press Inc.: New York, 2001. (21) Treloar, L. R. G. The Physics of Rubber Elasticity; Oxford University Press Inc.: New York, 2005. (22) Rubinstein, M.; Panyukov, S. Elasticity of polymer networks. Macromolecules 2002, 35, 6670−6686. (23) Mark, J. E.; Erman, B. Rubberlike Elasticity: A Molecular Primer; Cambridge University Press: New York, 2007. (24) Rubinstein, M.; Colby, R. H. Polymer Physics; Oxford University Press Inc.: New York, 2017. (25) Davidson, J. D.; Goulbourne, N. C. A nonaffine network model for elastomers undergoing finite deformations. J. Mech. Phys. Solids 2013, 61, 1784−1797. (26) Xiang, Y.; Zhong, D.; Wang, P.; Mao, G.; Yu, H.; Qu, S. A general constitutive model of soft elastomers. J. Mech. Phys. Solids 2018, 117, 110−122. (27) James, H. M.; Guth, E. Theory of the elastic properties of rubber. J. Chem. Phys. 1943, 11, 455−481. (28) Einstein, A. Theorie der Opaleszenz von homogenen Flüssigkeiten und Flüssigkeitsgemischen in der Nähe des kritischen Zustandes. Ann. Phys. 1910, 338, 1275−1298. (29) Kuhn, W. Gestalt und Eigenschaften fadenförmiger Moleküle in Lösungen (und im elastisch festen Zustande). Angew. Chem. 1936, 49, 858−862. (30) James, H. M.; Guth, E. Theory of the increase in rigidity of rubber during cure. J. Chem. Phys. 1947, 15, 669−683. (31) Wall, F. T. Statistical Thermodynamics of Rubber. II. J. Chem. Phys. 1942, 10, 485−488. (32) Treloar, L. R. G. The elasticity of a network of long-chain molecules. I. Trans. Faraday Soc. 1943, 39, 36−41. (33) Flory, P. J.; Gordon, M.; McCrum, N. G. Statistical Thermodynamics of Random Networks [and Discussion]. Proc. R. Soc. London, Ser. A 1976, 351, 351−380. (34) Rubinstein, M.; Panyukov, S. Nonaffine deformation and elasticity of polymer networks. Macromolecules 1997, 30, 8036−8044. (35) Gottlieb, M.; Macosko, C. W.; Benjamin, G. S.; Meyers, K. O.; Merrill, E. W. Equilibrium modulus of model poly(dimethylsiloxane) networks. Macromolecules 1981, 14, 1039−1046. (36) Lin, T.-S.; Wang, R.; Johnson, J. A.; Olsen, B. D. Revisiting the Elasticity Theory for Real Gaussian Phantom Networks. Macromolecules 2019, 52, 1685−1694. (37) James, H. M. Statistical properties of networks of flexible chains. J. Chem. Phys. 1947, 15, 651−668. H

DOI: 10.1021/acs.macromol.9b00262 Macromolecules XXXX, XXX, XXX−XXX