The bond rupture force for sulfur chains calculated from quantum chemistry simulations and its relevance to the tensile strength of vulcanized rubber

David E. Hanson * and John L. Barber
Theoretical Division, Los Alamos National Laboratory, P.O. Box 1663, Los Alamos, New Mexico 87545, USA. E-mail:

Received 2nd October 2017 , Accepted 20th November 2017

First published on 20th November 2017

From quantum chemistry simulations using density functional theory, we obtain the total electronic energy of an eight-atom sulfur chain as its end-to-end distance is extended until S–S bond rupture occurs. We find that a sulfur chain can be extended by about 40% beyond its nominally straight conformation, where it experiences rupture at an end-to-end tension of about 1.5 nN. Using this rupture force as the chain failure limit in an explicit polymer network simulation model (EPnet), we predict the tensile failure stress for sulfur crosslinked (vulcanized) natural rubber. Quantitative agreement with published experimental data for the failure stress is obtained in these simulations if we assume that only about 30% of the sulfur chains produce viable network crosslinks. Surprisingly, we also find that the failure stress of a rubber network does not scale linearly with the chain failure force limit.


Although natural rubber was introduced to Europe ca. 1500, its poor durability rendered it unsuitable for economically viable applications and, for several centuries, it was primarily regarded as a scientific curiosity. It wasn’t until 1838, when Charles Goodyear1 found that by adding a few percent sulfur and applying heat it's properties improved dramatically, making it attractive for numerous commercial applications. While sulfur has more allotropes than any other element,2 the most common form is an 8-atom ring and it is likely that this is what Goodyear used in his vulcanization process. Upon heating to about 160 °C for 30 minutes, the sulfur rings open due to thermally induced bond rupture, producing a short chain with radicals at each end. Although the chemistry of rubber crosslinking is complex, likely involving many intermediate reactions that are still not completely understood, the consensus is that the terminal radicals can replace hydrogen atoms on rubber molecules, allowing some of the sulfur chains to form covalent crosslinks between neighboring rubber molecules. The material is transformed from a liquid state into a continuous molecular network.

In time, other compounds, such as dicumyl peroxide, were discovered that could also produce crosslinks, and these were found to produce significantly stronger rubbers. Tensile strengths of peroxide crosslinked rubbers have been measured up to 35 MPa3 but experiments on sulfur crosslinked systems have yielded much lower values. Treloar4 reported about 6.3 MPa using only sulfur as the crosslinking agent and Kluppel et al.5 reported about 5 MPa using sulfur plus several additives. Why the sulfur-crosslinked rubber was so much weaker remained something of a mystery; some researchers suggested that the sulfur somehow can ‘inhibit crystallizations by producing structural modifications’ of the material.6 Here we will present the results of studies showing how the bond rupture strength of the sulfur chains affects the tensile strength of rubber networks.

Quantum chemistry bond rupture simulations

A model sulfur (S) crosslink in vulcanized rubber was constructed by bonding the ends of an eight-atom sulfur chain to short alkane chains, composed of three carbon (C) atoms, terminated with hydrogen atoms. The sulfur chain was attached to the alkane moities as shown in Fig. 1a. All calculations were performed with the Gaussian09 suite of electronic structure codes7 using the hybrid B3LYP functional, in conjunction with a standard 6-31g(d) basis set. The spin unrestricted self-consistent-field (SCF) calculation of the ground state at each point along the extension curve is straightforward. All energies are expressed with respect to the optimized initial conformation energy, which has been set to zero. The simulation methodology was essentially the same as that used in a previous study8 of mechanically induced bond rupture in polyisoprene. The chain was stretched to bond rupture in the simulation by moving apart the carbon atoms (labeled C1 and C2 in Fig. 1a), in steps of 0.025 nm, and optimizing the geometry, while holding the end-to-end distance between carbon atoms C1 and C2 fixed. The electronic structure optimization necessarily imposes a temperature of 0 K. While this may seem inappropriate if we are making comparisons to experiments carried out at room temperature, we shall see below that the total electronic energy stored in the 7 S–S bonds of the chain at the point of bond rupture is much greater than the thermal energy kBT, where T is the temperature and kB is Boltzmann's constant. Consequently, we are justified in ignoring thermal contributions. The initial end-to-end distance (C1–C2 in Fig. 1a) was 0.6158 nm and bond rupture, shown in Fig. 1b, occurred at the center of the chain between atoms S1 and S2 in Fig. 1b, at an end-to-end distance of 1.5908 nm. The change in electronic energy remains near zero until an end-to-end distance of about 1.25 nm, when bond distortions begin to occur. We find that an S–S bond can stretch by nearly a factor of two before rupturing. The initial distance between the center S atom pair was 0.2135 nm and increased to 0.4045 nm just before bond rupture. The total electronic energy as a function of end-to-end distance is shown in Fig. 2. The abrupt decrease in energy, occurring at an end-to-end distance of about 1 nm, is due to an increase in the bond angles for S atoms 2–3–4 and 3–4–5, from about 103 deg to the nominal value of 109 deg, possibly an artifact of the optimization algorithm in the electronic structure code. At finite temperature, bond rotations would likely preclude such behavior. The tensile force, obtained by numerically differentiating the energy curve in Fig. 2 with respect to the chain end-to-end distance, is shown in Fig. 3. We believe that the abrupt dips in the force are attributable to an artifact of the optimization algorithm of the electronic structure code. Small variations in the electronic energy become more pronounced when the derivative is taken. The maximum force of about 1.5 nN occurs slightly before an end-to-end distance of 1.52 nm. Since this force is much smaller than the value calculated for bond rupture in a polyisoprene chain (6.8 nN),8 we expect that it can account for the low failure stress observed in vulcanized rubber, about 6.3 MPa as shown in Fig. 4.
image file: c7cp06730e-f1.tif
Fig. 1 (a) Model sulfur crosslink: optimized eight-atom sulfur chain terminated by short sections of 3 carbon alkane chains. The sulfur chain is extended by increasing the distance between carbon atoms C1 and C2 in increments of 0.025 nm. (b) Sulfur chain at point of bond rupture at center S1–S2 bond. The missing bonds between other S atom pairs are an artifact of the molecule rendering program.

image file: c7cp06730e-f2.tif
Fig. 2 Quantum chemistry total optimized electronic energy (relative to initial state) as eight-atom sulfur chain between C atoms C1 and C2 (Fig. 1) is extended to bond rupture. S–S bond distortions begin at an end-to-end distance of ∼1.2 nm.

image file: c7cp06730e-f3.tif
Fig. 3 Quantum chemistry calculated tensile force as sulfur chain extended to bond rupture. Force is calculated as the numerical derivative of the curve in Fig. 3.

image file: c7cp06730e-f4.tif
Fig. 4 EPnet simulations of stress vs. tensile strain for dicumyl peroxide crosslinking (dashed line, open red circles), sulfur chain crosslinking (solid line, open green triangles, grey error bars) and experimental data of Treloar4 (blue filled circles) for natural rubber crosslinked with 8 pph sulfur (crosslink density: 5 × 1019 cm−3).

Network morphology

Not surprisingly, the crosslinking agent used affects the morphology of a natural rubber network. The crosslinking mechanism for both sulfur and dicumyl peroxide are thought to be similar: thermal decomposition of each agent produces two radicals that abstract hydrogen atoms from neighboring polyisoprene chains. In one case, the sulfur chain becomes the crosslink node, in the other the two chains become covalently bonded via a C–C bond and the network crosslink is considered to be the covalent bond. Schematically, we may represent the C–C crosslink node as the letter ‘X’, two carbon atoms connected by a single covalent bond with two network chains emanating from each (for a functionality of four). Previous quantum chemistry studies8 have shown that a polyisoprene chain can be stretched about 40% beyond its initial contour length before rupture, a condition that we shall refer to as ‘taut’. Simulations showed that the weakest bond along a polyisoprene chain is between the two carbon atoms, each of which has two hydrogen atoms bonded to it. The rupture force was calculated to be about 6.8 nN.

In the case of sulfur crosslinking, we can visualize the node as the letter ‘H’, with the horizontal line representing a sulfur chain. (See Fig. 5. for schematic illustrations.) Since the rupture force for a sulfur bond is much lower than for an isoprene chain (1.5 nN vs. 6.8 nN) in rubber networks, it is the crosslink itself that is the weakest point for chain rupture. At rupture, the crosslink node essentially disappears and the local network reverts back to its un-crosslinked state, i.e., replaced by two longer chains (Fig. 5c). It is possible that the network path containing a taut chain does not include the sulfur chain at either of its terminal nodes (Fig. 5b). At each of its crosslink nodes, the taut network path could lie along the original chains before the crosslink occurred, i.e., not along the cross of the ‘H’. However, this is unlikely since a chain has three possible network paths at each end and, of the nine possible paths through both end crosslink nodes, only one avoids both sulfur chains. For eight of the nine possible paths, at least one sulfur chain will lie along the taut network path and it will be the weakest link for chain rupture.

image file: c7cp06730e-f5.tif
Fig. 5 Schematic diagrams illustrating network morphology. Here solid lines represent polyisoprene chains, and dashed lines indicate sulfur crosslinks. (a) A portion of a taut network path which passes through a sulfur crosslink. (b) A portion of a taut network path which does not pass through a sulfur crosslink. (c) Changes in the local morphology of a network due to the break of a sulfur crosslink. Here chain segments 1–4, initially connected to the same crosslink, combine upon crosslink rupture to form two longer chains (1 + 3 and 2 + 4).

To quantitatively relate the lower rupture force to the macroscopic failure stress of a network requires numerical simulations of a rubber network as it is subjected to high tensile strain.

Network simulations

Numerical network simulations of tensile strain were carried out using the EPnet code,9 with the same methodology used previously10 for poly-isoprene (natural rubber) networks. The code constructs a representative volume element of an explicit, 3-dimensional network that respects both the Random Walk distribution of chain end-to-end distances and the exponentially decreasing probability for longer chain lengths. It uses the Molecular Kink9 network chain models to compute chain forces in response to imposed strains. There are three chain extension force models, depending on the degree of strain in individual chains. Two are based on local entropy changes in the chain as it becomes straight, and the third treats chain extension between taut and bond rupture. Used here, the term ‘taut’ refers to a chain whose end-to-end distance is equal to or greater than its initial contour length, where the contour length is defined as the sum of the end-to-end distances of its constituent isoprene units. The force model for taut chain extensions is obtained from quantum chemistry simulations8 of isoprene chains stretched to bond rupture and it has no free parameters. It is this force model that is most relevant for network strains leading to material failure. Numerical simulations show that, at high tensile strain, network paths emerge, comprised of sequential taut chains that span the entire simulation cell and actually support nearly all of the tensile stress.

Model parameter values and the network construction algorithms used for this study are the same as those found10 to agree with tensile and compressive strain experiments on natural rubber crosslinked with dicumyl peroxide.11 These parameter values were found to lie within the range of uncertainty of estimates obtained by independent analyses of Molecular Dynamics12 and Statistical Mechanics studies of polyisoprene.13 The tensile stress calculated during a simulation is determined only by the network morphology and the chain extension force models.9 The only difference in the input parameters to the EPnet code for this study was that the force threshold for chain failure (bond rupture) was set to 1.5 nN rather than the value for dicumyl peroxide crosslinked poly-isoprene of 6.8 nN.8

The simulation cell was a 55.3 nm cube with periodic boundary conditions, containing approximately 10[thin space (1/6-em)]500 randomly placed network nodes (crosslinks) and 21[thin space (1/6-em)]000 chains. Network chains were constructed between nodes using an algorithm that respects both a Random Walk probability distribution for end-to-end distance and a chain length probability distribution based on the assumption that all isoprene units have an equal probability to participate as a crosslink node. The cell was large enough to ensure that, on average, there were at least 10 sequential network chains along any axis.9 As macroscopic strain (the change in cell length divided by the initial length) is imposed, network node coordinates follow an affine trajectory until the net vector force on a node exceeds a threshold of 10 pN. Nodes are maintained at force equilibrium (net force below 10 pN) by an iterative relaxation procedure at the end of each strain step. Since the length of a sulfur crosslink chain is much less than that of a network chain, it is assumed to have zero length in the network simulations.

Tensile strain is imposed on the simulation cell in increments of 0.1 and the engineering stress (tensile force divided by the un-strained cross sectional area of the simulation cell) is computed at each step, up to a maximum strain of 10. This is repeated for ten statistically independent representative volume elements of the network to obtain an average stress vs. strain curve. The average stress vs. tensile strain from these simulations is shown in Fig. 4 with the experimental data of Treloar4 for comparison. Since the cross linking efficiency of sulfur in natural rubber is not known, we chose a crosslink density for the simulations that produced acceptable agreement with experiment over both high and low strain regions. Using our value of 1.5 nN for the chain rupture force limit and a crosslink density of 5 × 1019 cm−3 produced acceptable agreement with the experimental data. Engineering stress vs. strain calculated from network simulations are shown in Fig. 4 (open green triangles, solid line). Error bars, shown in grey, are the standard deviation of the 10 independent simulations. The large statistical variation of the stress for strains greater than 9 is due to the small number of taut network chains near failure. While the agreement between theory and experiment for the stress vs. strain is not perfect, we believe that it is not unreasonable for an initial study. We have only a single experiment data set (without error bars) and the disagreement could be due to deficiencies in the theory or experimental artifacts. The failure stress (the average over the stress plateau region) from the simulations was 7.3 MPa, is about 15% greater than the experimental value of 6.3 MPa. Treloar commented that his specimen would have had a higher tensile limit if a ‘dumbell’ shaped sample had been used. For comparison, we also include the results from simulations with the chain rupture force limit set to the value for dicumyl peroxide crosslinked natural rubber of 6.8 MPa8 (open red circles and solid line).

The network crosslink density that we used for our simulations is much less than the stoichiometric value to which a concentration of 8 pph (parts per hundred, by weight) of sulfur corresponds. The molecular weight of an eight-atom sulfur chain (256 Da) is nearly the same as that of dicumyl peroxide (270 Da) so, if every sulfur chain produced a network node, a sulfur concentration of 8 pph would produce a very high crosslink density (∼16 × 1019 cm−3) and this would make the calculated stress vs. strain curve very inconsistent with the experiment. Evidently, the crosslinking reaction of sulfur is significantly less efficient than unity, with only about 30% of the sulfur chains able to produce viable network crosslinks. The 30% factor for the crosslinking efficiency for sulfur crosslinked rubber does affect the failure stress to some extent because it partially determines the crosslink density of the network. The morphology of the network is determined entirely by the density of crosslink nodes (the crosslinker concentration times the crosslinking efficiency) and this morphology affects both the maximum stress and the shape of the stress vs. strain curve, especially where the stress upturn occurs. The failure stress is determined by two parameters: (1) the chain force rupture limit and (2) the density of network chains. For the network simulations to be consistent with the stress vs. strain experiment, both parameters must be set independently. So far as we know, the 30% crosslinking efficiency only pertains to sulfur. For dicumyl peroxide, the cross linking efficiency is thought to be unity and this value is consistent with our previous studies on polyisoprene at low and intermediate strains.

One can conjecture several possible reasons why the crosslinking efficiency might be so low: both ends of a sulfur chain could react with the same polyisoprene chain to produce an 8-atom loop instead of a crosslink; two 8-atom sulfur chains might bond to produce a 16-atom chain or an abstracted H atom might terminate the end of the sulfur chain, removing the S radical.

The EPnet code uses a simple algorithm to detect the failure of network chains due to bond rupture: a chain is considered to have ruptured when either the chain tension or strain (ratio of its end-to-end distance to the un-strained contour length) exceeds its bond rupture limit. Network chains that are flagged as ruptured are ignored in further simulation strain steps. This algorithm was designed for natural rubber networks crosslinked with dicumyl peroxide, i.e., having a strong crosslink consisting of a C–C bond. In these networks, bond rupture can occur at any isoprene unit along a chain if the chain tension exceeds 6.8 nN. However, for sulfur crosslinked networks, bond rupture occurs on the connecting sulfur chain and this gives rise to an inconsistency in the way EPnet treats chain failure. Because the EPnet code assumes that only individual chains can rupture, a rupture event results in the loss of only one chain at each of the end nodes, leaving three remaining chains. But for sulfur crosslinked networks, a rupture event removes the crosslink entirely and only two network chains remain. We would expect that this inconsistency would cause EPnet to over estimate the network stress near network failure. We attempted to estimate the effect of this inconsistency by modifying the EPnet code to delete not only the ruptured chain but also an additional chain having the highest tension at one of its crosslink nodes, leaving only two chains at one crosslink node instead of three. However, we observed no significant change in the failure stress, so we believe that the inconsistency is not significant.

The ratio of bond rupture forces for sulfur chains vs. C–C bonds is 0.22 (1.5 nN/6.8 nN) but the ratio of the failure stresses from our simulations for the two cases is only 0.14 (7.3 MPa/52 MPa10). Clearly the simulated failure stress does not scale linearly with the chain rupture force limit. To explore this further, we performed an additional series of EPnet simulations with chain force limits ranging between 1.5 and 6.8 nN. Fig. 6 shows that the network failure stress increases with the rupture force limit faster than linearly. The dashed straight line between the first two points is included to show the reader the expected behavior if the scaling were linear. Fig. 7 shows the calculated average failure strain and fraction of network chains that are taut. As the chain rupture limit increases, the network is able to support a higher strain because of the increase in the number of taut chains and it is this effect that leads to a higher failure stress. Values for the average failure stress and strain for the different chain force limits were obtained visually from plots of stress vs. strain averaged over 5 statistically independent EPnet simulations, analogous to those described above for 1.5 nN, the bond rupture force for sulfur. The accuracy for these estimates is probably comparable to the error bars shown in Fig. 4.

image file: c7cp06730e-f6.tif
Fig. 6 Failure stress calculated from EPnet network simulations vs. chain failure force (solid line, blue triangles) and extrapolated straight line from first two points (dashed line) for comparison. The points labeled ‘S–S’ and ‘C–C’ correspond to the bond rupture forces for sulfur and polyisoprene, respectively.

image file: c7cp06730e-f7.tif
Fig. 7 Simulated network failure strain (solid line, filled triangles) and fraction of network chains that are taut (dashed line, open circles) vs. chain failure force.


From Quantum chemistry simulations, we have calculated the tensile force and strain of an eight-atom sulfur chain as it is stretched to its breaking (rupture) point. The bond rupture force limit that we obtained (1.5 nN) was then used as the network chain failure limit in numerical simulations of explicit rubber networks using the EPnet code. All other input parameter values for these simulations were taken from previous studies10 of dicumyl peroxide crosslinked natural rubber networks that were found to be in agreement with experiments.11 Quantitative agreement with experimental stress vs. strain data for vulcanized rubber4 is achieved if we assume that only about 30% of the sulfur chains produce viable network crosslink nodes, i.e., producing a crosslink density of 5 × 1019 cm−3. The failure stress obtained from our simulations, about 7.3 MPa, is close to the experimental value of 6.3 MPa but much less than the value calculated for dicumyl peroxide crosslinked systems, 52 MPa.10 We believe that the lower failure stress can be completely attributed to the lower tensile strength of the sulfur chain crosslinks. We also find that the failure stress does not scale linearly with the bond rupture force (chain failure force limit). This work provides additional direct confirmation for the use of a numerical explicit polymer network model (EPnet) to study the mechanical response and failure of molecular networks.

Conflicts of interest

There are no conflicts to declare.


This work was partially performed under the auspices of Los Alamos National Laboratory, which is operated by Los Alamos National Security, LLC, for the National Nuclear Security Administration of the U.S. Department of Energy under contract DE-AC52-06NA25396.


  1. C. Goodyear, Gum-elastic and its Varieties: With a Detailed Account of its Applications and Uses, and of the Discovery of Vulcanization, Author, New Haven, CT, 1853 Search PubMed .
  2. R. Steudel and B. Eckert, Top. Curr. Chem., 2003, 230, 1–80 CrossRef CAS .
  3. C. Bell, D. Stinson and A. Thomas, Rubber Chem. Technol., 1982, 55, 66–75 CrossRef .
  4. L. R. G. Treloar, Trans. Faraday Soc., 1944, 40, 0059 RSC .
  5. H. Kluppel, H. Menge, H. Schmidt, H. Schneider and R. H. Schuster, Macromolecules, 2001, 34, 8107–8116 CrossRef .
  6. G. Gee, J. Polym. Sci., 1947, 2, 452–462 CrossRef .
  7. M. J. Frisch, et al., Gaussian 09, Revision A.1, Gaussian, Inc., Wallingford CT, 2009 Search PubMed.
  8. D. E. Hanson and R. L. Martin, J. Chem. Phys., 2009, 130, 064903 CrossRef PubMed .
  9. D. E. Hanson and J. L. Barber, Contemp. Phys., 2015, 56, 319–337 CrossRef CAS .
  10. D. E. Hanson and J. L. Barber, Modell. Simul. Mater. Sci. Eng., 2013, 21, 075013 CrossRef .
  11. P. H. Mott and C. M. Roland, Macromolecules, 1996, 29, 6941 CrossRef CAS .
  12. D. E. Hanson and R. L. Martin, J. Chem. Phys., 2010, 133, 084903 CrossRef PubMed .
  13. D. E. Hanson, J. L. Barber and G. Subramanian, J. Chem. Phys., 2013, 139, 224906 CrossRef PubMed .

This journal is © the Owner Societies 2018