Effect of cross-linking degree on thermal-mechanical properties of vulcanized SBR

Tian Yan 1,2,3

Kejian Wang 1✉ Emailwangkj@mail.buct.edu.cn

Xiuying Zhao 2

1 Mechanical and Electronical colleges Beijing University of Chemical Technology 100029 Beijing China

2 Key Laboratory of Beijing City on Preparation and Processing of Novel Polymer Materials Beijing University of Chemical Technology 100029 Beijing China

3 Jiangsu Hengtong Intelligent Equipment Co., LTD 215214 Suzhou China

Tian Yan^{1, 2, 3}, Kejian Wang^1*, Xiuying Zhao²

1. Mechanical and Electronical colleges, Beijing University of Chemical Technology, 100029, Beijing, China

2. Key Laboratory of Beijing City on Preparation and Processing of Novel Polymer Materials, Beijing University of Chemical Technology, 100029, Beijing, China

3. Jiangsu Hengtong Intelligent Equipment Co., LTD, Suzhou, 215214, China

Corresponding authors:

Kejian Wang: wangkj@mail.buct.edu.cn

Abstract

$\dot {\gamma }$

) enhances the susceptibility of molecular chains to orientation, thereby decreasing η by 83.4–93.2%. Moreover, greater

$\dot {\gamma }$

resulted in heightened intermolecular forces, constraining the movement of the molecular chains and consequently increasing η_b by 88.6-100.9%. Conversely, higher Dc intensifies intermolecular constraints, elevating η and η_b by 65.5–23.3% and 6.6–11.1%, respectively, within the same shear regime. Thermally, a 48.4% increase in κ is observed at Dc = 8.0 compared to Dc = 1.0, attributed to amplified low-frequency phonon density in tightly cross-linking networks. Meanwhile, restricted chain mobility at elevated Dc suppresses heat absorption, decreasing C_p and C_v by 8.9% and 8.6%, respectively, across 300-330K. Critically, these insights provide actionable guidelines for tailoring vulcanization parameters to balance processability and service performance in SBR-based products. The elucidated micro-mechanisms that spanning phonon transport modulation and chain dynamics restriction offer fundamental principles for designing cross-linked elastomers with targeted multifunctional properties.

Keywords:

vulcanized SBR

thermal properties

mechanical properties

MD simulations

micro mechanisms

1. Introduction

Rubber has found extensive applications in diverse fields like construction, medicine, and aerospace ^[1–3], owing to its favorable cost-effectiveness and outstanding damping characteristics. The growing interest in styrene-butadiene rubber (SBR) is attributed to its affordability and ease of processing. It is widely recognized that the thermal and mechanical properties of rubber significantly influence its overall performance. Specifically, in applications where rubber serves as a thermal conductor, superior thermal conductivity is typically sought after, whereas for applications in shock absorbers, greater emphasis is placed on its mechanical properties. Studying the effects and mechanisms of various factors on the thermal and mechanical properties of rubber is crucial for the development of high-performance rubber matrices. Research on the thermal and mechanical properties of rubber matrices primarily focuses on fillers ^[4–6] and blend modifications ^[7–8]. It is essential to note that uncured mixed rubber cannot fully meet practical application requirements and thus requires vulcanization. The characteristics of vulcanization significantly influence the thermal and mechanical properties of rubber, directly impacting its usability. Therefore, investigating the impact of vulcanization characteristics on the thermal and mechanical properties of rubber, and elucidating the underlying micro mechanisms, is vital for the targeted advancement of high-performance rubber materials.

Molecular dynamics (MD) simulations offer a superior method compared to traditional experiments for elucidating the microdynamic mechanisms underlying macroscopic phenomena ^[9–11]. Zhang et al ^[12] utilized MD simulation techniques to analyze the deformation behavior of Styrene-Butadiene Rubber (SBR) with varying degrees of cross-linking during uniaxial tensile processes, examining the correlation between molecular chain structure and tensile properties. Their findings revealed that cross-linking SBR initiates plastic deformation when strain surpasses 220%, with no overall fracture observed within a strain range of 400%. In the realm of thermal property exploration, Bhowmik et al ^[13] calculated the specific heat capacity (Cp) of typical polymers, including polytetrafluoroethylene rubber, at room temperature using MD simulations. Their study highlighted that Cp values derived from MD simulations exceeded experimental results, attributable to the lack of precision in the selected force field. In a study by Cai et al ^[6], a blend of MD simulations and experimental data was employed to evaluate the impact of graphene with diverse functional groups on the thermal conductivity and dielectric properties of carboxylated Acrylonitrile Butadiene Rubber (NBR). Results suggested that enhancements in the dielectric and thermodynamic properties of composite materials are predominantly influenced by interfacial behaviors. The scrutiny of vulcanization properties on the thermal and mechanical characteristics of vulcanized Styrene-Butadiene Rubber (SBR) through MD simulations presents a novel investigative angle.

This study investigates the impact of changes in the degree of cross-linking (Dc) on the thermal properties and viscoelastic behavior of SBR using MD simulation and experimental techniques. A molecular dynamics model of vulcanized SBR was developed, enabling exploration of viscoelastic responses at different shear rates and assessment of how Dc influences thermal conductivity and specific heat capacity of the vulcanized SBR system. Furthermore, the study reveals the micro-mechanisms through which Dc affects the thermal properties and viscoelastic behavior of vulcanized SBR. Experimental validation was performed to confirm the accuracy of the molecular dynamic simulation. The objective of this chapter is to introduce a novel approach for investigating the influence and micro-mechanisms of Dc on the thermal and mechanical properties of vulcanized SBR.

2. Simulation details and experiment methods

2.1 Simulation part

2.1.1 Establishment of MD model and selection of force field

The SBR matrix utilized in this study was a random copolymer comprising styrene and butadiene monomers. The butadiene monomer is comprised of 1,2-butadiene, cis-1,4-butadiene and trans-1,4-butadiene. Details regarding the monomer composition in the individual SBR chain can be found in Table 1. Each SBR polymer chain consists of 60 repeating units and the initial SBR structure was generated using the Amorphous Cell module within Materials Studio 7.0 software.

Table 1

Monomer components in a single SBR chain
Monomer	Styrene	1, 2 Butadiene	cis-1,4-Butadiene	trans-1,4-Butadiene
Content (wt. %)	23.5	14	12.5	50

The pcff force field was employed to characterize atomic interactions, Eq. 1 represents the formulation of the pcff force field.

$\begin{gathered} {E_{pot}}=\sum\limits_{b} {[{k_2}{{(b - {b_0})}^2}+{k_3}{{(b - {b_0})}^2}+{k_4}{{(b - {b_0})}^2}]} \hfill \\ +\sum\limits_{\theta } {[{H_2}{{(\theta - {\theta _0})}^2}+{H_3}{{(\theta - {\theta _0})}^3}+{H_4}{{(\theta - {\theta _0})}^4}]} \hfill \\ +\sum\limits_{\varphi } {[{V_1}(1 - \cos \varphi )+{V_2}(1 - \cos 2\varphi )+{V_3}(1 - \cos 3\varphi )]} +\sum\limits_{\chi } {{K_\chi }{\chi ^2}} \hfill \\ +\sum\limits_{{b,{b^,}}} {{F_{bb}}(b - {b_0})({b^,} - {b^,}_{0})+\sum\limits_{{b,\theta }} {{F_{b\theta }}(b - {b_0})(\theta - {\theta _0})} } +\sum\limits_{{{\theta ^,}{\theta ^,}}} {{F_{\theta {\theta ^,}}}(\theta - {\theta ^,})} (\theta - {\theta ^,}) \hfill \\ +\sum\limits_{{b,\varphi }} {(b - {b_0})(V\cos \varphi +{V_2}\cos 2\varphi +{V_3}\cos 3\varphi )} \hfill \\ +\sum\limits_{{\theta ,\varphi }} {(\theta - {\theta _0})(V\cos \varphi +{V_2}\cos 2\varphi +{V_3}\cos 3\varphi )} \hfill \\ +\sum\limits_{{b,\theta ,\varphi }} {{F_{b\theta \varphi }}(b - {b_0})(\theta - {\theta _0})\cos \varphi } \hfill \\ +\sum\limits_{{i,j}} {\frac{{{q_i}{q_j}}}{{\varepsilon {r_{ij}}}}} \hfill \\ \end{gathered}$

The van der Waals interaction among atoms was described by the LJ 9 − 6 potential as illustrated in Eq. 2.

${V_{lj}}(r)=\left\{ \begin{gathered} 4\varepsilon {\left( {\frac{\sigma }{r}} \right)^9} - {\left( {\frac{\sigma }{r}} \right)^6},(r \leqslant {r_c}) \hfill \\ 0,(r>{r_c}) \hfill \\ \end{gathered} \right.$

where

$\varepsilon$

is the bond energy at the equilibrium position(kcal/mol),

$\sigma$

is the collision diameter(Å), r_c is cutoff radius(Å), this work takes 12.5 Å, r is the distance between two atoms (Å).

The vulcanized SBR models with different Dc were established according to Eq. 3. The Dc in our work were 0.0, 1.0, 2.0, 5.0, 6.0and 8.0, respectively.

$D{\text{c}}=Ns/Mc$

where Ns is the number of cross-linking bonds, Mc is the number of SBR single chains.

The rules of molecular chain cross-linking are as follows: if the distance between two carbon (C) atoms on the main chain of SBR is in the range of 1.6–3.45 Å, then a sulfur (S) atom will be added at the geometric center of these two C atoms to form a C-S-C bond ^[14]. Considering the phenomenon of self-cross-linking in the actual vulcanization process of SBR, every two C atoms on a single chain can also form a C-S-C bond ^[15]. Additionally, once a C atom has formed a cross-linking bond, it cannot form new cross-linking bonds with other S atoms. Cross-linking is implemented using Perl scripts in Materials Studio 7.0 software. The four monomers composing SBR, SBR single chain, the molecular formula of cross-linking and the cross-linking results are shown in Fig. 1.

In this study, unless indicated otherwise, the temperature and pressure of the system in all MD simulations were regulated by the Anderson thermostat and the Berendsen barostat, respectively. Dynamics integration was carried out utilizing the Velocity-Verlet algorithm. The MD simulations were executed using the LAMMPS software ^[16] on a high-performance supercomputer with 120 cores, and the outcomes are visualized using The Open Visualization Tool (OVITO) software ^[17].

Fig. 1

Modeling details of MD models including four monomers that make up SBR chain (a), single chain containing 60 monomers (b), cross-linking model (c), the molecular formula of the cross-linking (d) and a local magnification (e) of red box in (c). In (c) and (e), the magenta beads represent carbon atoms of main chain, the gray beads represent carbon atoms of benzene ring, the white beads represent hydrogen atoms and the yellow beads represent sulfur atoms

We established two models of different sizes to simulate different properties of vulcanized SBR with varying Dc. The small model comprised 50 SBR single chains, whereas the large model comprised 120 single chains. Annealing relaxation was conducted on both cross-linking models using specific methods. For the small model, a 5 ns NPT relaxation at 300 K and 1 atm was followed by increasing the system temperature to 600 K under the NVT ensemble, then cooling it back to 300 K for annealing, with each annealing cycle lasting 4 ns and repeated 4 times. Subsequently, a 5 ns final NPT relaxation was carried out at 300 K and 1 atm. In the case of the large model, a 10 ns NPT relaxation at 300 K and 1 atm was followed by a similar annealing process but with an 8 ns duration for each cycle, repeated 4 times. Finally, a 10 ns final NPT relaxation was performed at 300 K and 1 atm for the large model. Notably, both the large and small models used a simulation time step of 1 fs for relaxation and annealing. During the last 5 ns of the final NPT relaxation, mean square displacement (MSD) and radial distribution function (RDF) were recorded for the small model. Further details on the model construction and validation of the approach's accuracy can be found in our previous work ^[18]. For further details on the MD simulation techniques, such as NVT ensemble, NPT ensemble and velocity-Verlet algorithm, see Appendix Ⅰ.

Upon completion of the final NPT relaxation of the model, the thermodynamic properties were computed, focusing primarily on thermal conductivity and specific heat, alongside mechanical properties containing shear viscosity and bulk viscosity.

2.1.2 Simulation of shear viscosity

The Non-equilibrium Molecular Dynamics (NEMD) method was employed, various Dc SBR systems were simulated at different strain rates (

$\dot {\gamma }$

) using the SLLOD motion equations to model the cyclic shearing process for shear viscosity (

$\eta$

) calculations. The cyclic shearing was applied in the XY plane of small model with the strain rates of 5 ps^− 1, 10 ps^− 1 and 20 ps^− 1, respectively, and the time step was 0.5 fs. Each SBR system with different Dc underwent 10 complete shearing cycles at varying rates. Subsequently, the shear viscosity was computed for the diverse Dc SBR systems based on Eq. 4.

$\eta {\text{=}}\frac{{{p_{XY}}}}{{\nabla {v_{XY}}}}$

where

${p_{XY}}$

is the pressure tensor along XY direction (kcal/mol) and

$\nabla {v_{XY}}$

is the gradient of the velocity of the particle in X direction along Y direction (s^− 1),

$\nabla {v_{XY}}{\text{=}}d{v_X}/dY$

${v_X}$

is the velocity of the particle in X direction (m/s) and

is the length of simulation box along Y direction.

2.1.3 Simulation of bulk viscosity

The simulation method for bulk viscosity (

${\eta _{\text{b}}}$

) aligns with that for

$\eta$

, however, their calculation methods differ. During computation, the non-diagonal elements of the pressure tensor in Eq. 4 must be substituted with diagonal elements ^[19]. A detailed explanation of the simulation method for

${\eta _{\text{b}}}$

is omitted, with the calculation approach elucidated in Eq. 5.

${\eta _{\text{b}}}=\frac{{{p_{XX}}+{p_{YY}}}}{{2\nabla {v_{XY}}}}$

where

${p_{XX}}$

and

${p_{YY}}$

are the pressure tensors along the X direction and Y direction, respectively.

2.1.4 Simulation of thermal conductivity

Reverse Non-Equilibrium Molecular Dynamics (rNEMD) technique developed by Müller-Plathe ^[20], also referred to as the M-P method, was employed to calculate the thermal conductivity (

$\kappa$

) of various Dc SBR systems. Renowned for its superior computational accuracy and practicality, this approach finds wide application in the investigation of material thermal conductive properties. Figure 2 illustrates a schematic diagram derived from simulations using M-P method. The method establishes a consistent temperature gradient by exchanging atomic momentum for computations of

$\kappa$

. Momentum exchange instances the transfer of greater energy towards the proximity of the heat source, reciprocally. Over time, the system attains a stable temperature gradient.

Fig. 2

Diagram of simulating thermal conductivity by M-P method

The specific simulation method for

$\kappa$

is outlined as follows: The large model is partitioned into 100 layers along the direction of heat flow. The central two layers are designated as the hot region, while the outermost layers (comprising four layers at both ends) are designated as the cold region. The simulation utilizes a time step of 1 fs, encompassing a total simulation duration of 5 ns. Momentum exchanges occur at intervals of 0.2 ps in NVT ensemble. Subsequent to the simulation, the

$\kappa$

of various Dc SBR systems was computed in accordance with the Fourier law of heat conduction. It is crucial to highlight that the model is subject to periodic boundary conditions in the direction of heat flow, resulting in bidirectional heat flow. This necessitates dividing the effective heat flow in one direction by 2 to arrive at the ultimate value, as delineated in Eq. 6.

$\kappa =\frac{{QdZ}}{{2tSdT}}$

where Q is the heat flow generated during the exchange of momentum(W), t is total simulation time(fs), S is the cross-sectional area through which the heat flow flows(Å²),

$dT$

and

$dZ$

temperature difference (K) and distance (A) between hot zone and cold zone, respectively.

2.1.4 Simulation of specific heat capacity

The specific heat capacity can be categorized into two types: specific heat capacity at constant pressure (C_p) and specific heat capacity at constant volume (C_v). These parameters represent the material's ability to absorb heat under consistent pressure and consistent volume conditions. As reported in reference ^[13], the C_p and C_v values of polymers exhibit relatively stable behavior within the temperature range of 300–330 K. To determine C_p, the following approach is employed: Various Dc SBR systems were brought to equilibrium in an NPT ensemble within the temperature range of 300–330 K, using a temperature increment of 5 K. The simulation employs a time step of 1 fs and a pressure of 1 atm, with each temperature increment having an equilibration duration of 1 ns. The averaged enthalpy of system is determined in final 400 ps of the equilibration and a graphical representation illustrating the changes in enthalpy values with temperature between 300 and 330 K is generated based on the simulation outcomes. The C_p value can be derived using Eq. 7.

${C_{\text{p}}}=\frac{{\Delta H}}{{\Delta T \times m}}$

where

$\frac{{\Delta H}}{{\Delta T}}$

is the enthalpy change rate of the system with temperature(J/K), that is the slope of the enthalpy-temperature curve, m is the quality of the system.

The C_v of a substance refers to its ability to absorb heat under constant volume, while enthalpy H is defined as the sum of a system's internal energy U and the product of pressure P and volume V, expressed as H = U + PV, where U, P, and V signify the system's internal energy, pressure, and volume, respectively. Consequently, employing enthalpy for calculation of C_v can lead to inaccurate outcomes. To mitigate the impact of pressure and volume on the calculations, the NVT ensemble is utilized to stabilize the system and uphold a constant volume. The parameters such as temperature range, temperature step size, time step and equilibration time required for C_v mirror those needed for C_p calculations. In the concluding 400 ps of the simulation, the average internal energy (as opposed to enthalpy) of system is computed for each temperature step. The C_v value can be computed using Eq. 8.

${C_{\text{v}}}=\frac{{\Delta U}}{{\Delta T \times m}}$

where

$\frac{{\Delta U}}{{\Delta T}}$

is the internal energy change rate with temperature in the system (J/K), that is the slope of the internal energy-temperature curve.

2.2 Experiment part

2.2.1 Preparation of mixed SBR

Styrene-butadiene rubber (SBR) serves as the matrix material, obtained through blending and vulcanization in an open mixer. Table 2 presents the raw materials and composition of the SBR compound, including the sulfur content calculated using the MD model for various Dc SBR systems. Following calculations, it was determined that the sulfur contents for various Dc SBR systems were 1.04 phr, 1.99 phr, 5.03 phr, 6.29 phr and 8.28 phr.

Table 2

Raw materials and contents of mixed SBR
Raw material	Mark	Manufacturer	Content (phr)
Solution-polymerized SBR	1502	China Petrochemical Qilu Branch	100
Insoluble sulfur	IS60	Avalible in the market	variable
Stearic acid(SA)	SA1801	Avalible in the market	1.0
Zinc oxide(ZnO)	--	Hebei Sanshi (Group) Co., LTD	4.0
Acclerator DZ	--	Avalible in the market	1.0
Acclerator TT	--	Avalible in the market	0.1

Firstly, the pure SBR was placed on an open mill with a roll diameter of 152.4 mm and a small roll spacing for plasticization at room temperature. After achieving uniform stretching, ZnO, SA, accelerator TT, accelerator DZ and insoluble sulfur are added sequentially. To ensure the even dispersion of each additive, the rubber underwent repeated cutting and plasticization after the addition of each additive. Subsequently, following the addition of all additives, the mixed rubber underwent a process of alternating between kneading and pressing into a triangular shape three times to further guarantee the even dispersion of the additives. Finally, the roll spacing is increased to sheet the rubber.

2.2.2 Vulcanization characteristic

Fig. 3

Vulcanization characteristic curves of SBR with different sulfur content

After allowing the mixed rubber to stand for 24 hours, its vulcanization characteristic was tested using rotorless vulcanizer manufactured by GOTHCH Company and following the standard GB/T 16584 − 1996. About 5g of the mixed rubber was then placed into the mold when the mold's temperature reaches 150°C and the testing commences within 5 seconds of mold closure. The test was terminated once the vulcanization curve stabilizes and maintains equilibrium for a minimum of 10 minutes. The T₉₀, which represents the time taken for the rubber to achieve 90% vulcanization, can be derived from the vulcanization curve.

The vulcanization characteristic curves of mixed SBR with various sulfur contents are depicted in Fig. 3. It can be observed from Fig. 3 that as the sulfur content increases, the torque of SBR also rises. At a sulfur content of 8.28 phr, the torque alteration is approximately 10N·m, nearly five times higher than that of mixed SBR with 1.04 phr sulfur content. This occurrence is primarily attributed to the escalation in sulfur content, which results in the formation of a more compact cross-linking network, consequently enhancing the strength of vulcanized SBR.

2.2.3 Cross-linking degree test

After allowing the vulcanized SBR to stand for 24 h, the cross-linking density analyzer with the model XLDS-15 is used to test the cross-linking density by the method of nuclear magnetic resonance (NMR). The SBR samples are tested under specific conditions of 90℃ for the temperature, 0.35 T for the magnetic field intensity and 15 MHz for the resonance frequency. The cross-linking density is calculated by measuring the relaxation time (T₂) of the SBR samples with varying sulfur contents. The same formulae are tested five times and the average values are adopted as the final result.

2.2.4 Thermal conductive test

The mixed SBR with varying sulfur contents were vulcanized using a flat vulcanizer to fabricate samples for thermal conductivity assessment. These samples are cylindrical, with diameters and thicknesses measuring 38 mm and 10 mm, respectively. Figure 4 illustrates the samples utilized in the thermal conductivity analysis. The thermal conductivity of the vulcanized SBR with varying sulfur contents was measured using a DTC-300 thermal conductivity meter following the standard of ASTM E1530.

Fig. 4

Samples for thermal conductivity testing

2.2.5 Heat capacity test

The C_p of vulcanized SBR with varying sulfur contents at 20°C was determined using the indirect test method. This method involves utilizing a standard substance with consistent physical and chemical properties within the test temperature range and a known C_p value. In this research, sapphire was selected as the standard substance. The DSC curves of the empty crucible, sapphire and the test sample were recorded at a heating rate of 10 K/min. From the test outcomes, the disparity between the enthalpy change rates of the test sample, sapphire and the enthalpy change rate of the empty crucible was computed and designated as y and y', respectively. The equation for the test sample is C_p= C'×((m'×y)/(m×y')), where C', m and m' represent the C_p values at 20°C for sapphire, the mass of test sample and sapphire, respectively. However, current research in the field of C_v predominantly depends on MD simulation, theoretical derivation or indirect conversion using C_p and other parameters to determine the substance's C_v, with minimal literature on direct experimental measurement of C_v due to the difficulty in maintaining constant volume conditions in laboratory setups during temperature variations. Consequently, this study does not include experimental measurement of C_v.

3. Results and discussion

3.1 Cross-linking degree analysis

The cross-linking degree test results are presented in Fig. 5. It is evident from Fig. 5 that the cross-linking degree calculated using various Dc MD models exhibits minimal error when compared to the test results of vulcanized SBR with corresponding sulfur content. This suggests that the model can effectively replace experimental procedures.

Fig. 5

The cross-linking degree calculated by different Dc models and the cross-linking degree of vulcanized SBR with corresponding sulfur content are measured

3.2 Analysis of equilibrium process

The time evolution of MSD and RDF for the C atoms on the main chain of SBR in various Dc SBR small models during 5 ns of final NPT relaxation process is illustrated in Fig. 6 and Fig. 7.

Fig. 6

MSD of SBR with different Dc during NPT relaxation stage

Fig. 7

RDF of SBR with different Dc after 5 ns of NPT relaxation (a) and (b) is local amplification at red box in (a)

By analyzing the movement of molecular chains and the arrangement of atoms within these chains, the influence on micro structure of vulcanized SBR is examined through Fig. 6 and Fig. 7. The MSD decreases as Dc increases during the final NPT relaxation process. Following 5 ns of NPT equilibrium, the MSD for the SBR system is approximately six times higher when Dc = 1.0 compared to Dc = 8.0. This disparity is primarily due to the increased presence of cross-linking in higher Dc systems, which restricts the mobility of the molecular chains. Additionally, Fig. 7 illustrates that the density of carbon atoms along the SBR main chain rises with higher Dc. Specifically, the packing efficiency of C atoms on the main chain within the SBR system with Dc = 8.0 is approximately 6% greater than that of the system with Dc = 1.0. This phenomenon is attributed to the abundance of cross-linking in the SBR system with higher Dc, resulting in more compact arrangement of the SBR molecular chains.

3.3 Analysis of thermal properties

3.3.1 Thermal conductivity property

Take the SBR system with Dc = 8.0 serves as an example, Fig. 8 shows the extracted temperature gradient between the cold and hot layers along the direction of heat flow within the simulation, as well as diverse Dc SBR systems derived from simulations and experiments. Furthermore, Fig. 8 illustrates the temperature discrepancy (ΔT) alongside the corresponding distance (Δl). Subsequent to ample momentum interchange between the cold and hot regions, a conspicuous temperature gradient emerges (Fig. 8(a)) and κ escalating with the augmented Dc of the SBR system (Fig. 8(b)). In contrast with the SBR system at Dc = 1.0, the κ of SBR system at Dc = 8.0 exhibits a 48.76% improvement. This phenomenon primarily stems from the denser molecular chain arrangement in SBR system with higher Dc, facilitated by the presence of cross-linking bond (Fig. 7). The condensed molecular chains engender additional pathways for thermal conduction. Conversely, in the SBR system with lower Dc results in larger intermolecular gaps, impeding heat transfer and causing a decline in κ. This assertion finds validation in the research conducted by Kikugawa ^[22] examining the impact of cross-linking on the thermal conductivity of amorphous polymers.

Furthermore, Fig. 8(b) illustrates that as the Dc rises from 1.0 to 8.0, the discrepancies between the simulation and experiment results of the SBR system are 19.54%, 15.94%, 13.05%, 16.13% and 10.33% correspondingly. This robustly validates the precision of forecasting κ of the SBR systems with different Dc using MD simulation.

Fig. 8

Temperature difference along direction of heat flow when calculating the κ (taking the SBR system with Dc = 8.0 as an example) (a) and the κ of SBR with different Dc obtained from MD simulations and experiments(b)

Understanding the phonon vibrations within vulcanized SBR is crucial for deciphering the micro-mechanisms that influence its thermal conductivity. Phonon Density of State (PDOS) obtained through the Fourier transformation of the velocity autocorrelation function during atomic momentum exchange serves as a vital tool in analyzing heat conduction properties. The calculation method is detailed in Eq. 9.

$I(\omega )=\operatorname{Re} \left[ {\frac{1}{{2\pi }}\int_{0}^{\infty } {\left\langle {\vec {v}\left( t \right)\cdot \vec {v}\left( 0 \right)} \right\rangle } {e^{i\omega t}}dt} \right]$

where

$I(\omega )$

is Phonon density corresponding to frequency, Re represents a real number,

$\omega$

is the corresponding frequency (Hz), t is the momentum exchange time (fs),

$\vec {v}\left( 0 \right)$

and

$\vec {v}\left( t \right)$

is the velocity vector of the atom at the initial and final time (Å/fs). According to the trajectories of each atom in the momentum exchange process, Eq. 9 can be rewritten as Eq. 10.

$I(\omega ) \propto \mathop {\lim }\limits_{{N \to \infty }} \sum\limits_{{n=0}}^{N} {{C_{\text{v}}}n\Delta t\cos (\omega n\Delta t)} \Delta t$

where N is the total steps of simulation, n is the current number of simulations, C_v is the specific heat capacity at constant volume (J/K·g) and Δt is time difference (fs).

The PDOS in various Dc SBR systems during the momentum exchange process is illustrated in Fig. 9. To facilitate analysis, the peaks in PDOS of each Dc SBR system in Fig. 9 are indicated by circles. As Dc increases, the respective PDOS peaks gradually shift towards lower frequencies. Specifically, for Dc = 1.0, the peak in PDOS is observed around 37.5 THz. In contrast, for Dc values of 2.0 and 5.0, the peak corresponds to a frequency of approximately 7.5 THz, additionally, the SBR system with Dc = 5.0 exhibits a peak near 42 THz. Notably, as Dc values reach 6.0 and 8.0, the maximum PDOS value emerges at around 0.8 THz.

The primary function of the PDOS is to highlight that higher phonon density within a specific frequency range results in more noticeable thermal diffusion and enhanced thermal transport effects, especially in the low-frequency spectrum (less than 14 THz) ^[23]. Analysis of Fig. 9 reveals that SBR systems with elevated Dc exhibit increased phonon density at lower frequencies, thereby amplifying thermal diffusion and enhancing overall thermal conductivity. This enhancement in thermal conductivity as Dc rises can be attributed to the presence of cross-linking in SBR, which not only promote a more tightly-packed molecular chain configuration but also elevate the PDOS within the low-frequency range.

Fig. 9

PODS of SBR with different Dc during the momentum exchange process

3.3.2 Heat absorption property

Illustrating with the calculation of C_p for the Dc = 8.0 SBR system, this study extracted the relationship between enthalpy and temperature throughout the heating process and presented the C_p and C_v values for various Dc SBR systems in Fig. 10. The enthalpy of the vulcanized SBR system exhibits a linear increase with temperature. The curve fit yielded an R² value of 0.998 (Fig. 10(a)), affirming the reliability of using enthalpy's temperature change rate for C_p determination. Furthermore, both C_p and C_v exhibit a gradual decline as Dc increases (Fig. 10(b)), consistent with findings by Weston ^[24] regarding the influence of cross-linking density on specific heat capacity. In comparison to the SBR system with Dc = 1.0, the system with Dc = 8.0 saw reductions of 8.92% and 8.58% in C_p and C_v, respectively. This phenomenon stems from the cross-linking converting SBR chains into a stiff network structure, restricting chain mobility as depicted in Fig. 6. For SBR systems with lower Dc, chains can rotate and vibrate freely, storing energy as heat and yielding higher specific heat capacities. Conversely, higher Dc constrains rotational and vibrational energies, lowering the specific heat capacity.

Fig. 10

Relationship between temperature and enthalpy during the heating process when calculating C_p (Taking the SBR with Dc = 8.0 as an example) (a) and the simulated results of C_p and C_v of SBR with different Dc (b)

Another notable observation in Fig. 10(b) is that C_v consistently remains lower than C_p for vulcanized SBR. This difference arises from the distinct thermodynamic processes governing heat absorption under constant volume and constant pressure conditions. Under constant volume, the absorbed heat primarily increases the internal energy of the system by enhancing molecular thermal motion (e.g., vibrational, rotational, and translational modes), without performing work on the system’s volume. In contrast, under constant pressure, the absorbed heat not only increases the internal energy but also performs work by expanding the volume of the material. This additional work term, accounted for by enthalpy changes, results in a higher heat capacity for C_p compared to C_v. Therefore, C_p calculations must consider both internal energy changes and volume work, whereas C_v calculations focus solely on internal energy changes. This distinction underscores the fundamental thermodynamic principles underlying the difference between C_p and C_v in vulcanized SBR systems.

Figure 11 shows the C_p obtained through experiments and MD simulations, as well as the ratio of the simulation results to the experimental results, in order to determine their error. When Dc is 1.0, 2.0, 5.0, 6.0 and 8.0, the simulated C_p are 1.78, 1.8, 1.84, 1.88, and 1.85 times the experimental results, respectively. According to the comparative standards reported in Reference ^[13], the error between the experiment and simulation results is small. This also demonstrates the accuracy of using MD simulation methods to predict the C_p of different Dc SBR systems.

Fig. 11

C_p obtained from experiments and MD simulations, as well as the ratio between them

3.4 Analysis of rheological properties

3.4.1 Shear viscosity

The η and η_b of SBR systems with different Dc under varied

$\dot {\gamma }$

loads are calculated by Eq. 4 and Eq. 5, as depicted in Fig. 12. As the

$\dot {\gamma }$

increases, a shear thinning effect is observed in SBR systems with varying Dc. Specifically,

$\eta$

decreases with escalating

$\dot {\gamma }$

(Fig. 12(a)), reflecting a typical rheological response where higher Dc values in the SBR system correspond to increased

$\eta$

. Conversely,

${\eta _{\text{b}}}$

in diverse Dc SBR systems increases with rising

$\dot {\gamma }$

(Fig. 12(b)) and elevated SBR system Dc results in higher

${\eta _{\text{b}}}$

. Subsequent analysis will delve into the underlying microscopic mechanisms of these observations.

The increase in η and η_b with Dc can be elucidated through Fig. 6. As Dc increases, the motion of SBR chains becomes restricted, impeding molecular chain movement during shear and enhancing system rigidity, consequently leading to the escalation of η and η_b. The phenomenon of shear thinning is attributed to the orientation of polymer chains during shearing process. The orientation of polymer chain can be determined using Eq. 11 ^[25], with the results illustrated in Fig. 13. Additionally, Fig. 13 presents the average Rg value for SBR systems with different Dc under varying levels of

$\dot {\gamma }$

$\left\langle {\tan 2\theta } \right\rangle =\frac{{2\left\langle {{G_{xy}}} \right\rangle }}{{\left\langle {{G_{xx}}} \right\rangle - \left\langle {{G_{yy}}} \right\rangle }}$

where θ is the orientation angle of the molecular chain (°), <G_xx>, <G_yy> and < G_xy> are the squared components of the rotation tensor of the molecular chain in different directions, respectively. Additionally, it should be noted that < G_xx>, <G_yy> and < G_xy> are obtained by measuring the C-C bonds on the SBR main chain during the simulated shearing process.

Fig. 12

η (a) and η_b (b) of SBR with different Dc under different

$\dot {\gamma }$

The radius of gyration (Rg) quantifies the spatial compactness of polymer chains, reflecting their conformational flexibility, i.e., the higher the Rg value, the poorer the flexibility of the polymer system. In cross-linking SBR systems, Rg increases with rising cross-linking degree (Dc) due to increased chain confinement. Rg can be calculated from molecular dynamics (MD) trajectories using Eq. 12.

$R{\text{g}}=\sqrt {\frac{1}{N}{{\sum\limits_{{i=1}}^{N} {\left( {{{\mathbf{r}}_i} - {{\mathbf{r}}_{{\text{cm}}}}} \right)} }^2}}$

where

${{\mathbf{r}}_i}$

is the position of the i-th monomer,

${{\mathbf{r}}_{{\text{cm}}}}$

is the chain’s center of mass, and N is the number of monomers per chain.

Examination of Fig. 13 (a) reveals a clear inverse relationship between the orientation degree of SBR molecular chains and the increasing

$\dot {\gamma }$

, aligning closely with the trend observed in η. A lower orientation degree signifies a more orderly arrangement of molecular chains, resulting in improved flexibility. This orderly arrangement is a key factor contributing to the shear thinning behavior in the SBR system as

$\dot {\gamma }$

increases. Moreover, Fig. 13 (b) illustrates that as

$\dot {\gamma }$

rises, the average value of R_g decreases, indicating enhanced molecular chain flexibility under higher

$\dot {\gamma }$

conditions. Furthermore, Fig. 13 indicate that SBR systems with smaller Dc exhibit a higher tendency to orient during shearing process. This propensity is primarily due to the fact that increasing Dc in the SBR system corresponds to a higher number of cross-linking, impeding molecular chain movement and hindering orientation. Apart from molecular chain alignment, the breakage of cross-linking bond during shearing is another factor contributing to shear thinning, as it promotes the flowability of molecular chains. Additionally, the acceleration of cross-linking rupture with increasing

$\dot {\gamma }$

intensifies the flowability of molecular chains, consequently resulting in a reduction of η.

Fig. 13

Variation of < tan 2θ > and average value of Rg with

$\dot {\gamma }$

of SBR with different Dc

3.4.2 Bulk viscosity

The correlation between the increase in η_b and the corresponding increase in

$\dot {\gamma }$

can be elucidated by considering two key factors during the shearing procedure: the intermolecular interaction forces and the dynamics of molecular chains. For instance, in the SBR system with Dc = 1.0, Fig. 14 illustrates the variations in intermolecular interaction energy (E_pair) and mean MSD under various shear load conditions of

$\dot {\gamma }$

When

$\dot {\gamma }$

is low, the molecular chains of SBR have sufficient time for reorganization, resulting in a lower η_b in the vulcanized SBR system. With an increase in

$\dot {\gamma }$

, the speed of molecular chain movement accelerates, impeding reorganization and leading to a rise in η_b. The vigorous movement of the molecular chains inevitably enhances the intermolecular interaction forces, as evident in Fig. 14(a): E_pair significantly increases with higher

$\dot {\gamma }$

. Strengthening intermolecular interaction forces impact the fluidity of molecular chains within the system. In particular, stronger intermolecular interaction forces correspond to diminished fluidity of the molecular chains. As illustrated in Fig. 14(b), increasing

$\dot {\gamma }$

causes a gradual decrease in MDS, indicating that higher

$\dot {\gamma }$

diminishes the fluidity of the molecular chains, ultimately resulting in an increase in η_b as

$\dot {\gamma }$

rises.

Fig. 14

Variation of E_pair (a) and MSD (b) with shearing period under different

$\dot {\gamma }$

(taking the SBR with Dc = 1.0 as an example)

In conclusion, as

$\dot {\gamma }$

increases, the SBR molecular chains undergo orientation, leading to a more orderly arrangement and enhanced flow characteristics, thereby displaying rheological properties. Furthermore, the rupture of cross-linking can further enhance the flow properties of the molecular chains; a larger

$\dot {\gamma }$

accelerates the rate of cross-linking breakage, thereby facilitating shear thinning. Additionally, the rise in component η_b with the increase in

$\dot {\gamma }$

is primarily attributed to the larger

$\dot {\gamma }$

strengthening intermolecular interactions and decreasing the flow characteristics of the molecular chains.

4. Conclusions

This work combines MD simulation with experiments to investigate the thermal and mechanical properties of different Dc vulcanized SBR. Compared to the SBR system with Dc = 1.0, the system with Dc = 8.0 showed an increase of 48.76% in some properties, while C_p and C_v decreased by 8.92% and 8.58% respectively. The micro-mechanism behind this is that the low frequency PDOS in the SBR system increases with the increase of Dc, and the system's heat absorption capacity decreases with the increase of Dc. Additionally, different Dc vulcanized SBR systems exhibited distinct rheological behaviors under shearing loads with different

$\dot {\gamma }$

. The fundamental reasons for the decrease and increase in η and η_b with increasing

$\dot {\gamma }$

are the molecular chain alignment and enhanced intermolecular forces during the shearing process. Ultimately, the accuracy of the MD simulation results was confirmed by experimental results.

Funding

The above work has not been given any kind of financial aid.

Data Availability

Data will be provided on the request.

Author Contribution

Tian Yan: simulation, experiments and writing.Kejian Wang: Supervisor, framework of the manuscript, central idea, most of financial support. Xiuying Zhao: financial support for experiments.

Declarations

of interest conflict

The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.

Electronic Supplementary Material

Below is the link to the electronic supplementary material

Supplementary Material 1

Additional Files

Additional file 1

References

X. Li, Y. Dong, Z. R. Li, et al. Experimental study on the temperature dependence of hyperelastic behavior of tire rubbers under moderate finite deformation[J]. Rubber Chem. Technol., 2011, 84: 215–228. DOI: https://doi.org/10.5254/1.3577534.

S.-J. Moon., J.-S. Kwak, J. -H. Chung, et al. A study on the hybrid mount against vibration and shock for naval ships[J]. Shock Vib., 2010, 17: 269–283. DOI: https://doi.org/10.3233/SAV-2010-0511.

D. Petkovic, M. Issa, N. D. Pavlovic, et al. Electrical properties estimation of conductive silicone rubber for tactile sensing structure[J]. J. Mech. Sci. Technol., 2011, 25: 1159–1165. DOI: https://doi.org/10.1108/02602281311299653.

W. Hong, J. P. Lin, X. H. Tian, et al. Viscoelasticity of nanosheet-filled polymer composites: Three regimes in the enhancement of moduli[J]. J. Phys. Chem. B, 2020, 124: 6437–6447. DOI: https://doi.org/10.1021/acs.jpcb.0c04235.

Y. Y. Gao, F. Y. Hu, Y. P. Wu, et al. Understanding the structural evolution under the oscillatory shear field to determine the viscoelastic behavior of nanorod filled polymer nanocomposites[J]. Comp. Mater. Sci., 2018, 142: 192–199. DOI: https://doi.org/10.1016/j.commatsci.2017.09.051.

F. Cai, Y. L. Luo, W. Yang, et al. Study on the thermal and dielectric properties of covalently modified GO/XNBR composites[J] Mater. Design, 2021, 198: 109335. DOI: https://doi.org/10.1016/j.matdes.2020.109335.

L. W. Mu, J. He, Y. F. Li, et al. Molecular origin of efficient phonon transfer in modulated polymer blends: Effect of hydrogen bonding on polymer coil size and assembled microstructure[J]. J. Phys. Chem. C, 2017, 121(26): 14204–14212. DOI: https://doi.org/10.1021/acs.jpcc.7b03726.

J. S. Seo, H. T. Jeon, T. H. Han. Rheological investigation of relaxation behavior of polycarbonate/acrylonitrile-butadiene-styrene blends[J]. Polymers, 2020, 12(9): 1916. DOI:

Y. Bai, Y. Z. Wang, X. Li, et al. Interaction mechanism and binding mode of phycocyanin to lysozyme: Molecular docking and molecular dynamics simulation[J]. Food Chem., 2024, 438: 138001. DOI: https://doi.org/10.1016/j.foodchem.2023.138001.

10.

N. R. Sun, Y. J. Li, N. X. Qiu, et al. Adsorption behaviors for clathrate hydrates of CO2 with mixed gases[J]. Fuel, 2024, 358: 130265. DOI: https://doi.org/10.1016/j.fuel.2023.130265.

11.

Z. H. Zhang, M. F. Long, N. Zheng, et al. Microstructural, physicochemical properties, and interaction mechanism of hydrogel nanoparticles modified by high catalytic activity transglutaminase crosslinking[J]. Food Hydrocolloid., 2024, 147: 109384. DOI: https://doi.org/10.1016/j.foodhyd.2023.109384.

12.

J. Zhang, W. Wang, Y. B. Wang, et al. Effect of cross-linked structures on mechanical properties of styrene-butadiene rubber via molecular dynamics simulation[J]. Macromol. Theor. Simul., 2022, 31: 2100054. DOI: https://doi.org/10.1002/mats.202100054.

13.

R. Bhowmik, S. Sihn, V. Varshney, et al. Calculation of specific heat of polymers using molecular dynamics simulations[J]. Polymer, 2019, 167: 176–181. DOI: https://doi.org/10.1016/j.polymer.2019.02.013.

14.

T. Zhang, H. B. Huang, W. Li, et al. Vulcanization modeling and mechanism for improved tribological performance of styrenebutadiene rubber at the atomic scale[J]. Tribol. Lett., 2020, 68(3): 83. DOI: https://doi.org/10.1007/s11249-020-01321-w.

15.

J. Kruzelak, I. Hudec. Vulcanization systems for rubber compounds based on IIR and halogenated IIR: An overview[J]. Rubber Chem. Technol, 2018, 91(1):167–183. DOI: https://doi.org/10.5254/rct-18-82609.

16.

S. Plimpton. Fast parallel algorithms for short-range molecular dynamics[J]. J. Comput. Phys., 1995, 117(1): 1–19.

17.

A. Stukowski. Visualization and analysis of atomistic simulation data with OVITO – the Open Visualization Tool[J]. Model. Simul. Mater. Sci., 2010, 18(1): 015012. DOI: https://doi.org/10.1088/0965-0393/18/1/015012.

18.

T. Yan, K. J. Wang, X. Y. Zhao, et al. Effect of cross-linking density on non-linear viscoelasticity of vulcanized SBR: A MD simulation and experimental study[J]. Int. J. Mol. Sci., 2023, 24: 9970. DOI: https://doi.org/10.3390/ijms24129970.

19.

D. I. Dimitrov, A. Milchev, K. Binder. Local viscosity in the vicinity of a wall coated by polymer brush from Green-Kubo relations[J]. Macromol. Theor. Simul., 2008, 17(6): 313–318. DOI: https://doi.org/10.1002/mats.200800038.

20.

F. Müller-Plathe. A simple nonequilibrium molecular dynamics method for calculating the thermal conductivity[J]. J Chem. Phys., 1997, 106(14): 6082–6085. DOI: http://jcp.aip.org/resource/1/JCPSA6/v106/i14?ver=pdfcov.

21.

J. Nichele, A. B. de Oliveira, L. S. D. Alves, et al. Accurate calculation of near-critical heat capacities Cp and Cv of argon using molecular dynamics[J]. J. Mol. Liq., 2017, 237: 65–70. DOI: https://doi.org/10.1016/j.molliq.2017.03.120.

22.

G. Kikugawa, T. G. Desai, P. Keblinski, et al. Effect of crosslink formation on heat conduction in amorphous polymers[J]. J. Appl. Phys., 2013, 114(3): 034302. DOI: https://doi.org/10.1063/1.4813505.

23.

A. Vasilev, T. Lorenz, C. Breitkopf. Prediction of thermal conductivities of rubbers by MD simulations—new insights[J]. Polymers, 2022, 14: 2046. DOI: https://doi.org/10.3390/polym14102046.

24.

D. A. Weston. Quantitative evaluation of rubber/silica particle interphase through measurement of heat capacity[D]. South Carolina, Clemson University, 2004.

25.

J. S. Yang, D. H. Huang. Rheological properties of ring and linear polymers under start-up shear by molecular dynamics simulations[J]. Acta Phys. Sin., 2019, 68(13): 138301. DOI: https://doi.org/10.7498/aps.68.20190403.

Yes

Abstract

Styrene-butadiene rubber (SBR) has attracted widespread attention due to its excellent thermal-mechanical properties, which critically govern its processing and service performance. While existing research predominantly focuses on blending and filling modifications, the influence of vulcanization characteristics on these properties and their underlying micro-mechanisms remain underexplored. This work systematically investigates how cross-linking degree (Dc) modulates the thermal-mechanical properties of vulcanized SBR, including shear viscosity η, bulk viscosity ηb, specific heat capacity(Cp, Cv), as well as the thermal conductivity(κ). Key findings reveal that increasing shear rate ( ) enhances the susceptibility of molecular chains to orientation, thereby decreasing η by 83.4-93.2%. Moreover, greater resulted in heightened intermolecular forces, constraining the movement of the molecular chains and consequently increasing ηb by 88.6-100.9%. Conversely, higher Dc intensifies intermolecular constraints, elevating η and ηb by 65.5–23.3% and 6.6-11.1%, respectively, within the same shear regime. Thermally, a 48.4% increase in κ is observed at Dc =8.0 compared to Dc = 1.0, attributed to amplified low-frequency phonon density in tightly cross-linking networks. Meanwhile, restricted chain mobility at elevated Dc suppresses heat absorption, decreasing Cp and Cv by 8.9% and 8.6%, respectively, across 300-330K. Critically, these insights provide actionable guidelines for tailoring vulcanization parameters to balance processability and service performance in SBR-based products. The elucidated micro-mechanisms that spanning phonon transport modulation and chain dynamics restriction offer fundamental principles for designing cross-linked elastomers with targeted multifunctional properties.