# Molecular Dynamics Simulations of Polymer Grouting Material: Its Mechanical Behavior under Uniaxial Tension, Cyclic Tensile Loading, and Stress Relaxation

Authors: Mingrui Du; Haoyue Zhang; Yan Zhao; Zongliang Huang; Chengzhi Zheng; Rui Wu; Xinlei Gao
Journal: Geofluids
 (2022)

---

## Abstract

At present, polyurethane (PU) has been extensively used as a grouting material in civil engineering. The mechanical properties of PU are the key to achieving the desirable grouting effect. This study presents the research results of the mechanical behavior of PU matrix under tensile, successive cyclic tensile, and stress relaxation at the nanoscale, using the coarse-grained molecular dynamics simulation method. The influences of the number of molecule chains and strain rate on the tensile mechanical properties are discussed, and the tensile deformation mechanism of PU matrix is revealed. The tensile strength of PU matrix is independent of loading path, and after yielding, the strain of PU matrix contains the elastic strain, plastic strain, and viscous strain. In the stress relaxation process, the evolution of the axial stress is mainly caused by the varied van der Waals interactions. The stress relaxation behavior of PU matrix can be described by the viscoelastic model consisting of one elastic element in parallel with one Maxwell element.


---

## Body

## 1. Introduction

Polyurethane (PU) [1], the product of polymerization reaction between isocyanate and polyols, is a versatile class of copolymer used in coatings [2], armors [3], and composites [4] for its properties of high durability, high toughness, high capacity to absorb and dissipate energy, and low density. At present, by virtue of the characteristics of the strong permeability resistance [5], fast-setting [6], high expansion force [7], and eco-friendliness [8], the foamed polyurethane has been increasingly adopted as grouting material in civil engineering. As shown in Figure 1, the representative applications of the polymer grouting material include filling soil cavities to lift the sedimentary underground pipelines [9, 10] or foundations [11, 12] and constructing barriers to prevent water seepage [13, 14] or the leakage of pollutants [14, 15]. In these applications, the two raw materials, isocyanate and polyols, are simultaneously injected into soil at the constant mass ratio as alternatives to cement-based grouts [16, 17], and after polymerization reaction, the PU matrices will be retained in the soil mass to play the reinforcing role [9, 10, 16]. For example, it has been reported that when the void surrounding the underground pipe was filled through PU grouting, under the same vertical load, the strain and bending moment in the pipe was declined by about 76.9% and 67.9% [10]. The mechanical properties of PU matrices are the key to achieving the desirable grouting effect, and the polymer grouting materials with insufficient strength may cause the refailure of the reinforced engineering structures [9–11, 18].Figure 1
The representative applications of the foamed polyurethane (PU) as grouting materials in geotechnical engineering.To date, numerous experimental studies have been carried out to investigate the tensile, compressive, flexural, and shear strength of polymer grouting materials [19–23]. For instance, the tensile strength and compressive strength of the plain PU composites are about 44.4 MPa and 15.5 MPa, and the appropriate addition of PVA fibers can improve the mechanical properties of PU composites by about 4.5-22.3% and 9.0-12.2% [19]; in Hussain et al.’s study [21], it was found that the experimentally tested flexural strength of PU composites could be about 11.0 MPa; in the study of Zheng et al. [22], it was found that the flexural strength and compressive strength of PU grouting materials with density of 0.1-0.5 g/cm3 were in the range of about 1.2-18.1 MPa and 1.0-7.5 MPa, respectively. It has been reported that, according to the micro morphology observations [6], PU grouting materials are porous materials containing numerous hollow closed foams, making them perform with obvious discontinuity characteristics. The foams are packed together, and their walls are made of amorphous PU [6], and under the external loading, it is the deformation of the foam walls, i.e., the PU matrix, that cause the failure of polymer grouting materials. By incorporating the beam deformation theory and the membrane stretching theory, the Gibson-Ashby equation that reflects the mechanical behavior of PU grouting materials has been proposed [24]. The Gibson-Ashby equation shows that knowing the mechanical properties of the PU matrix from the molecular scale is the foundation of the theoretical research on the mechanical properties of polymer grouting materials [24], which, however, is rarely seen at present.In the past decades, the molecular dynamics (MD) simulation has been proven to be a reliable and promising tool to explore the nanoscale mechanical behavior of various materials. For example, studies on the mechanical behavior of amorphous polymer materials [25, 26], the microscopic deformation characteristics of inorganic minerals or metals [27, 28], the tensile strength of a single carbon nanotube [29, 30] or one layer of graphene [31, 32], and the nanoscale shear enhancing mechanism of carbon nanotubes in cement [33], all using the MD simulation method, have been reported. A large amount of significant efforts to develop the united atom [34, 35] and coarse-grained [36–38] have been carried out. When the united atom or coarse-grained MD simulation methods are used, groups of atoms are lumped into single interaction sites with the simplified interaction potentials, and thus, the total number of degrees of freedom of the simulated models can be reduced. The united atom and coarse-grained MD methods require less computation, due to which, the time and length scales of the models can be increased to simulate the bulk behavior. With the united atom or coarse-grained MD simulation methods, the nanoscale mechanical properties and deformation mechanism of amorphous polyethylene under the tensile-compressive load [25] or pure tensile load [35] have been understood. In the study of Park et al. [38], the thermomechanical behavior of the shape-memory polyurethane copolymers was also studied using the coarse-grained MD simulation method. Studies on the nanoscale modulus [39], cohesion strength [40], and self-assembly behavior [41] of PU material with coarse-grained MD simulation method have also been reported. The coarse-grained MD simulation has provided a new approach to understand the mechanical behavior of PU matrix from the molecular perspective.In this study, using the coarse-grained MD simulation method, the mechanical behavior of the amorphous PU matrix under tensile was firstly studied, and then, both the successive cyclic tensile simulation and stress relaxation simulation were carried out. This study is aimed at addressing the mechanical properties of the amorphous PU matrices at the nanoscale under different types of loadings that are commonly seen in geotechnical engineering [9, 11, 18, 42]. The reason why the tensile simulation was firstly carried out is that under external load, the failure of the PU foams is mainly caused by the tension deformation [6, 24, 25]. This study is organized as follows: the construction of the coarse-grained model of amorphous PU and the simulation details are included in Section 2. In Section 3, both the influences of the number of the molecular chains and the strain rate on the tensile strength of PU matrix, as well as the tensile deformation mechanism of PU matrix, are discussed. The mechanical behavior of PU matrix under cyclic tensile loading and the stress relaxation characteristics of PU matrix under different tensile strains are included in Section 3. The main conclusions are summarized in Section 4.


## 2. Simulation Models and Methods

### 2.1. Coarse-Grained Method and Force Field

For the PU copolymers, different coarse-grained methods and force fields have been developed [38, 39, 43]. In Park et al.’s study [38], the hard segment of the PU molecular was divided into two different beads, one of which includes two halves of aromatic rings, and the remaining component of the benzenes and the urethane linkage were represented with another bead. As for the coarse-grained PU model proposed by Uddin et al. [39, 43], it has the physical and mechanical properties that are a good match with the experimental results; therefore, it was adopted in this study. More details are included as follows.As shown in Figure2(a), a single PU molecule chain is composed of the hard and soft segments, where the formers come from diphenylmethane diisocyanate (MDI) and butanediol [44]. The hard segments were treated as separate beads to obtain interactions between them. The segments from MDI consist of two phenyl rings connected by a -CH2 group, and the internal interactions are difficult to capture; therefore, the hard segment was divided into three types of beads [39, 43]. The –(CH2)4–O group, –(CONH)–C6H6–CH2 group, and the –(C6H6)–NHCOO group are represented by type A bead, type B bead, and type C bead, respectively. The molecular weights of A, B, and C beads are 72.1 g/mol, 135.2 g/mol, and 135.1 g/mol, respectively. The hard segments and soft segments in the PU copolymers often have a degree of polymerization in the range of 1–10 and 15–30, respectively [44]. Based on the commonly seen polymerization degrees of the hard segments and soft segments, the coarse-grained PU model (as shown in Figure 2(b)) was built. Figure 2(b) shows one-fifth of the coarse-grained PU molecule chain, and an entire coarse-grained PU chain contains 280 beads in total, composed of 230 A beads, 25 B beads, and 25 C beads. The polymerization degrees of the hard and soft segments in the coarse-grained model are 2 and 14, respectively, and the molecular weight is 23340 g/mol. The full-atomic polyurethane model mapped from the coarse-grained model is consistent with the typical PU molecule formula presented in literature [44], and thus, it is representative. With this coarse-grained method, for a typical PU matrix with the polymerization of 100, the total number of particles is reduced to 28000 from 389000, making it possible to compute relatively large systems at a physically meaningful simulation scale.Figure 2
(a) Construction of the coarse-grained PU model with three different kinds of beads; (b) 1/5 of the coarse-grained PU molecule chain.
(a)(b)The potential energy of the coarse-grained PU model contains three terms, the bond stretching, the bond angle bending, and the van der Waals interactions, whereas the term associated with the dihedral angle torsion is negligible [37]; therefore, the potential energy can be expressed as follows:
(1)Eptotal=Ebond+Eangle+Enonbond.The function forms ofEbond, Eangle, and Enonbond in Equation (1) are as follows:
(2)Ebondr=12Kbondr−r02,(3)Eangleθ=12Kangleθ−θ02,(4)Enonbondr=4εσr12−σr6,r≤rc,where Kbond and Kangle are the stiffness constants for the bond length and bond angle potential energy, respectively; r0 and θ0 are the equilibrium bond length and equilibrium bond angle, respectively; r and θ are the distance between two atoms and the angle between two bonds, respectively; ε is the energy well depth, σ is the distance between two nonbonded atoms at which the potential energy is zero; and rc is the cutoff distance which is defined as 2.5σ [44].In coarse-grained MD simulations [45, 46], Equation (2) and Equation (3) are the most commonly used functions to describe the variations of potential energy associated with bond length and bond angle, and the van der Waals interactions are often expressed by the Lennard-Jones 12-6 function (Equation (4)). Based on the trajectories of the full-atomic polyurethane equilibrium model developed using the polymer consistent force field, the probability distribution functions among the beads can be derived with the inverse Boltzmann (IBM) method [39, 43]. The center of mass of each group of atoms is defined as the center of mass of the corresponding bead. Then, the potentials on the bond stretching and the angle bending are fitted with Equation (2) and Equation (3), and the pair potentials are fitted with Equation (4), through which the coarse-grained force field parameters can be obtained. The coarse-grained force field parameters are listed in Tables 1–3.Table 1
Parameters of bond stretching potential energy.
Bond typer0 (ANG)Kbond (kcal/mol/ANG2)A-A5.211.53A-B6.111.59B-C7.115.36C-A6.261.34Table 2
Parameters of bond angle bending potential energy.
Bond angle typeθ (degree)Kangle (kcal/mol/rad2)A-A-A177.80.83A-A-B148.91.46A-B-C171.60.87B-C-A167.80.66C-A-B173.01.58C-A-A159.42.35Table 3
Parameters of nonbonded interactions.
Bead iBead jσ (ANG)ϵ (kcal/mol)AA5.620.31BB5.840.29CC5.660.32AB6.060.21BC6.030.69CA4.461.24


### 2.2. Constructing the Coarse-Grained Model of Amorphous PU Matrix

Figure3 shows the construction procedure of the coarse-grained model. First, an initially large size model was built by replicating the coarse-grained molecule chain, and then, the initial unstable model was relaxed at 600 K using a Langevin thermostat and the NVE ensemble for 10 ns. After that, the model was relaxed for another 10 ns at 600 K and then cooled to 298 K at a rate of 25 K per 50 ps, followed by a steady-state relaxation at 298 K for 10 ns, both with the NPT ensemble and the external pressure of one atmosphere (0.1 MPa). The simulation box was deformed to obtain a model with the desired density, after which, the model was equilibrated for 10 ns while keeping the density unchanged, both using the NVT ensemble at 298 K. The Berendsen thermostat was used to reset the temperature in both the NPT and NVT simulations, and the Berendsen barostat was used to reset the pressure of the system in NPT simulation.Figure 3
Construction of the coarse-grained model of amorphous PU matrix. The different molecular chains were painted with different colors, and the unit of side length is Å.The finally obtained coarse-grained PU models have a density of about 0.86 g/cm3, which is comparable with the density obtained from the full-atomic simulation (about 0.90 g/cm3) [39, 43] and the experimentally obtained density (0.97-1.24 g/cm3) [47]. In addition, the nanoscaled mechanical properties of the coarse-grained PU model are also comparable with the experimental ones [39, 43]. The coarse-grained models are large collections of PU chains that are randomly jumbled together. The cubic simulation cells with the side lengths of 76.9, 102.3, and 162.5 Å were established, and each was composed of 10, 25, and 100 PU chains, respectively (Figure 3). With these three different models, the influence of the number of chains on the tensile strength can be studied.


### 2.3. Simulation Details

For conducting the tensile simulation, a constant engineering strain rate was applied to the coarse-grained models with the NVT ensemble at the temperature of 298 K. The periodic boundary condition was adopted, and the Berendsen thermostat was used to reset the temperature of the system. According to the symmetric pressure tensor and the geometric size of the models, the tensile stress (σ) and engineering strain (ε) were computed, and the potential energy contributions associated with different terms were tracked as a function of ε. When presenting our results, true stress (σT) and true strain (εT) were used as alternatives of σ and ε, respectively, and they were calculated as follows:
(5)εT=ln1+ε,σT=σ1+ε.The strain rate of5×1010/s was applied to the three models with different side lengths until ε reached 1.0, and the tensile strain was along the x-axis. For the configuration containing 100 PU molecule chains, the uniaxial tension simulation at two other strain rates (2×1011/s and 1×1010/s) was also conducted to study the influences of strain rate. The tensile rate adopted here is in the range suggested by previous studies [25, 35].The successive cyclic tension loading simulation was conducted on the 100-chain coarse-grained model. Both the loading and unloading rates of the successive cyclic tensile simulation were5×1010/s, and a total of 7 cycles were simulated. For each cycle, the model was unloaded to zero stress at first and then reloaded to a larger strain than the previous cycle. εT at each unloading point was 0.04, 0.09, 0.18, 0.24, 0.32, 0.45, and 0.60, respectively.Finally, the stress relaxation simulation was conducted on the 100-chain coarse-grained model. The model was stretched to different true strains (0.05, 0.1, 0.2, 0.4, 0.6, and 0.8) at the rate of5×1010/s and then relaxed for a period of time while keeping the strain constant. The evolutions of the axial stress and potential energy of the model with the relaxation time were recorded. Both the model establishment and the tensile simulations were carried out by the LAMMPS package [48], and the open visualization tool OVITO was adopted for molecular visualization [49].


## 2.1. Coarse-Grained Method and Force Field

For the PU copolymers, different coarse-grained methods and force fields have been developed [38, 39, 43]. In Park et al.’s study [38], the hard segment of the PU molecular was divided into two different beads, one of which includes two halves of aromatic rings, and the remaining component of the benzenes and the urethane linkage were represented with another bead. As for the coarse-grained PU model proposed by Uddin et al. [39, 43], it has the physical and mechanical properties that are a good match with the experimental results; therefore, it was adopted in this study. More details are included as follows.As shown in Figure2(a), a single PU molecule chain is composed of the hard and soft segments, where the formers come from diphenylmethane diisocyanate (MDI) and butanediol [44]. The hard segments were treated as separate beads to obtain interactions between them. The segments from MDI consist of two phenyl rings connected by a -CH2 group, and the internal interactions are difficult to capture; therefore, the hard segment was divided into three types of beads [39, 43]. The –(CH2)4–O group, –(CONH)–C6H6–CH2 group, and the –(C6H6)–NHCOO group are represented by type A bead, type B bead, and type C bead, respectively. The molecular weights of A, B, and C beads are 72.1 g/mol, 135.2 g/mol, and 135.1 g/mol, respectively. The hard segments and soft segments in the PU copolymers often have a degree of polymerization in the range of 1–10 and 15–30, respectively [44]. Based on the commonly seen polymerization degrees of the hard segments and soft segments, the coarse-grained PU model (as shown in Figure 2(b)) was built. Figure 2(b) shows one-fifth of the coarse-grained PU molecule chain, and an entire coarse-grained PU chain contains 280 beads in total, composed of 230 A beads, 25 B beads, and 25 C beads. The polymerization degrees of the hard and soft segments in the coarse-grained model are 2 and 14, respectively, and the molecular weight is 23340 g/mol. The full-atomic polyurethane model mapped from the coarse-grained model is consistent with the typical PU molecule formula presented in literature [44], and thus, it is representative. With this coarse-grained method, for a typical PU matrix with the polymerization of 100, the total number of particles is reduced to 28000 from 389000, making it possible to compute relatively large systems at a physically meaningful simulation scale.Figure 2
(a) Construction of the coarse-grained PU model with three different kinds of beads; (b) 1/5 of the coarse-grained PU molecule chain.
(a)(b)The potential energy of the coarse-grained PU model contains three terms, the bond stretching, the bond angle bending, and the van der Waals interactions, whereas the term associated with the dihedral angle torsion is negligible [37]; therefore, the potential energy can be expressed as follows:
(1)Eptotal=Ebond+Eangle+Enonbond.The function forms ofEbond, Eangle, and Enonbond in Equation (1) are as follows:
(2)Ebondr=12Kbondr−r02,(3)Eangleθ=12Kangleθ−θ02,(4)Enonbondr=4εσr12−σr6,r≤rc,where Kbond and Kangle are the stiffness constants for the bond length and bond angle potential energy, respectively; r0 and θ0 are the equilibrium bond length and equilibrium bond angle, respectively; r and θ are the distance between two atoms and the angle between two bonds, respectively; ε is the energy well depth, σ is the distance between two nonbonded atoms at which the potential energy is zero; and rc is the cutoff distance which is defined as 2.5σ [44].In coarse-grained MD simulations [45, 46], Equation (2) and Equation (3) are the most commonly used functions to describe the variations of potential energy associated with bond length and bond angle, and the van der Waals interactions are often expressed by the Lennard-Jones 12-6 function (Equation (4)). Based on the trajectories of the full-atomic polyurethane equilibrium model developed using the polymer consistent force field, the probability distribution functions among the beads can be derived with the inverse Boltzmann (IBM) method [39, 43]. The center of mass of each group of atoms is defined as the center of mass of the corresponding bead. Then, the potentials on the bond stretching and the angle bending are fitted with Equation (2) and Equation (3), and the pair potentials are fitted with Equation (4), through which the coarse-grained force field parameters can be obtained. The coarse-grained force field parameters are listed in Tables 1–3.Table 1
Parameters of bond stretching potential energy.
Bond typer0 (ANG)Kbond (kcal/mol/ANG2)A-A5.211.53A-B6.111.59B-C7.115.36C-A6.261.34Table 2
Parameters of bond angle bending potential energy.
Bond angle typeθ (degree)Kangle (kcal/mol/rad2)A-A-A177.80.83A-A-B148.91.46A-B-C171.60.87B-C-A167.80.66C-A-B173.01.58C-A-A159.42.35Table 3
Parameters of nonbonded interactions.
Bead iBead jσ (ANG)ϵ (kcal/mol)AA5.620.31BB5.840.29CC5.660.32AB6.060.21BC6.030.69CA4.461.24


## 2.2. Constructing the Coarse-Grained Model of Amorphous PU Matrix

Figure3 shows the construction procedure of the coarse-grained model. First, an initially large size model was built by replicating the coarse-grained molecule chain, and then, the initial unstable model was relaxed at 600 K using a Langevin thermostat and the NVE ensemble for 10 ns. After that, the model was relaxed for another 10 ns at 600 K and then cooled to 298 K at a rate of 25 K per 50 ps, followed by a steady-state relaxation at 298 K for 10 ns, both with the NPT ensemble and the external pressure of one atmosphere (0.1 MPa). The simulation box was deformed to obtain a model with the desired density, after which, the model was equilibrated for 10 ns while keeping the density unchanged, both using the NVT ensemble at 298 K. The Berendsen thermostat was used to reset the temperature in both the NPT and NVT simulations, and the Berendsen barostat was used to reset the pressure of the system in NPT simulation.Figure 3
Construction of the coarse-grained model of amorphous PU matrix. The different molecular chains were painted with different colors, and the unit of side length is Å.The finally obtained coarse-grained PU models have a density of about 0.86 g/cm3, which is comparable with the density obtained from the full-atomic simulation (about 0.90 g/cm3) [39, 43] and the experimentally obtained density (0.97-1.24 g/cm3) [47]. In addition, the nanoscaled mechanical properties of the coarse-grained PU model are also comparable with the experimental ones [39, 43]. The coarse-grained models are large collections of PU chains that are randomly jumbled together. The cubic simulation cells with the side lengths of 76.9, 102.3, and 162.5 Å were established, and each was composed of 10, 25, and 100 PU chains, respectively (Figure 3). With these three different models, the influence of the number of chains on the tensile strength can be studied.


## 2.3. Simulation Details

For conducting the tensile simulation, a constant engineering strain rate was applied to the coarse-grained models with the NVT ensemble at the temperature of 298 K. The periodic boundary condition was adopted, and the Berendsen thermostat was used to reset the temperature of the system. According to the symmetric pressure tensor and the geometric size of the models, the tensile stress (σ) and engineering strain (ε) were computed, and the potential energy contributions associated with different terms were tracked as a function of ε. When presenting our results, true stress (σT) and true strain (εT) were used as alternatives of σ and ε, respectively, and they were calculated as follows:
(5)εT=ln1+ε,σT=σ1+ε.The strain rate of5×1010/s was applied to the three models with different side lengths until ε reached 1.0, and the tensile strain was along the x-axis. For the configuration containing 100 PU molecule chains, the uniaxial tension simulation at two other strain rates (2×1011/s and 1×1010/s) was also conducted to study the influences of strain rate. The tensile rate adopted here is in the range suggested by previous studies [25, 35].The successive cyclic tension loading simulation was conducted on the 100-chain coarse-grained model. Both the loading and unloading rates of the successive cyclic tensile simulation were5×1010/s, and a total of 7 cycles were simulated. For each cycle, the model was unloaded to zero stress at first and then reloaded to a larger strain than the previous cycle. εT at each unloading point was 0.04, 0.09, 0.18, 0.24, 0.32, 0.45, and 0.60, respectively.Finally, the stress relaxation simulation was conducted on the 100-chain coarse-grained model. The model was stretched to different true strains (0.05, 0.1, 0.2, 0.4, 0.6, and 0.8) at the rate of5×1010/s and then relaxed for a period of time while keeping the strain constant. The evolutions of the axial stress and potential energy of the model with the relaxation time were recorded. Both the model establishment and the tensile simulations were carried out by the LAMMPS package [48], and the open visualization tool OVITO was adopted for molecular visualization [49].


## 3. Results and Discussion

### 3.1. Effects of the Number of PU Chains on the Tensile Mechanical Properties

Figures4(a)–4(c) show the variation laws of the stress of the coarse-grained PU models with different numbers of molecular chains under uniaxial tension, from which one can see that when εT was lower than 0.03, the stress increased almost linearly with the increasing strain, corresponding to the elastic deformation behavior. When εT was over 0.03, the variation rate of the tensile stress started to decline, implying that plastic deformation had occurred to the model, and at εT of about 0.20 (yield strain), the tensile stress achieved a plateau. It can be found that the number of molecular chains has little influence on the yield strain or yield strength (both are about 0.20 and 30 MPa). The stress-strain curve becomes smoother with a higher number of PU chains, suggesting that the simulations moved closer to the bulk behavior. Therefore, in all the subsequent simulations, the coarse-grained PU model with 100 molecular chains was used.Figure 4
(a–c) Stress-strain curves of coarse-grained PU models with different numbers of chains under tensile at a strain rate of5×1010/s and (d) the evolution of molecular configuration containing 10 PU chains at different strains.
(a)(b)(c)(d)The snapshots of the deformation process of the configurations containing 10 PU molecule chains are presented in Figure4(d). It can be seen that the deformation of the model mainly consists of two parts: the slippage between the PU chains and the stretching of PU chains themselves. Before the yielding, the stretching deformation of each independent PU molecule chain was less significant than the slippage between them, whereas it became more and more obvious with the increasing strain after yielding. The slippage deformation occurred mainly along the direction of tensile strain.


### 3.2. Effects of Strain Rate on the Tensile Mechanical Properties

Figure5 shows the variations of the stress of the coarse-grained PU matrix containing 100 molecule chains under uniaxial tension at the strain rate of 2×1011/s and 1×1010/s. The curves follow the same laws as discussed in Section 3.1 except for the different yield strengths and different yield strains. According to Figures 4(c) and 5, one can see that the yield strength of the model declines with the declined strain rate. At the strain rates of 2×1011/s, 5×1010/s, and 1×1010/s, the yield strength was about 70, 30, and 16 MPa, respectively. On the contrary, with the decreasing strain rate, the yield strain increased gradually. For example, the yield strain was about 0.20 and 0.22 at the strain rate of 5×1010/s and 1×1010/s, respectively, which were about double that at the strain rate of 2×1011/s. Moreover, it can also be seen that when a higher strain rate was adopted, the stress-strain curves became smoother.Figure 5
Stress-strain curves of the coarse-grained PU matrix containing 100 chains under tensile tension at the strain rate of (a)2×1011/s and (b) 1×1010/s.
(a)(b)Figure6 illustrates the evolution of potential energy (PEtotal) of the coarse-grained PU matrix at different tensile strain rates, corresponding to the stress-strain behavior presented in Figures 4(c) and 5. The evolutions of individual components, i.e., the bond potential energy (PEbond), bond angle potential energy (PEangle), and van der Waals interactions (PEnonbond), are also included in Figure 6. E0 represents the potential energy at initial time.Figure 6
Variations of the potential energy of the coarse-grained PU matrix versus true strain for a strain rate of (a)2×1011/s, (b) 5×1010/s, and (c) 1×1010/s.
(a)(b)(c)From Figure6(a), it can be seen that when the axial strain increased at a rate of 2×1011/s, both PEtotal and PEnonbond increased gradually, and the increasing trends are curves with a downward opening, and that the PEtotal curve and PEnonbond curves almost coincided until εT reached to about 0.10 (the yield strain). PEtotal was about 6000 kcal/mol at the yield strain, and after yielding, PEtotal became lower than PEnonbond; meanwhile, the PEangle term started to decline in an almost linear way. PEangle was about -6000 kcal/mol at εT of 0.70. During the tensile procedure, PEbond increased with the increasing strain, but its variation was lower than those of the other two terms (with the maximum increment of lower than 1500 kcal/mol).The changedPEnonbond, PEangle, and PEbond are primarily attributed to the slippage between different molecules, bond angle bending behavior, and bond stretching behavior, respectively. Figure 6(a) shows that, at the strain rate of 2×1011/s, under uniaxial tension, only the slippage between different PU molecules occurred before yielding, but after yielding, in addition to the slippage between different PU molecules that still remained, on each PU molecule, the bond angles started to bend, and the bonds were stretched to synthetically accommodate deformation. The difference between PEnonbond and PEtotal after yielding was dissipated by the bond angle bending behavior and bond stretching behavior, but primarily via the former. Figures 6(b) and 6(c) show that whether at the strain rate of 5×1010/s or 1×1010/s, there was little change in the overall evolution laws of PEtotal and PEnonbond, but PEangle started to decline from the very beginning of the simulation, and that PEbond was always equal to zero. This indicates that at a lower strain rate, as soon as the tensile load was applied, the bond angles displayed bending behavior, but the bond lengths were little influenced.Figure6 also shows that the change of total potential energy at yielding is independent of the strain rate. At the three loading rates, both PEtotal at yielding are about 6000 kcal/mol, but after yielding, the PEtotal of the model sustaining a lower strain rate has a lower value, which is due to that PEnonbond has significantly declined. With the decreasing strain rate, the increased yield strain (Figures 4(c) and 5) leads to the increased volume strain of the model. Considering that the yield strength of the model is positively dependent on the potential energy but negatively influenced by the volume [50], when the strain rate declines, PEtotal is almost unchanged, but the increased volume strain of the coarse-grained PU models at yielding is the main reason for the lower yield strength.Figures4(d) and 6 show that the slippage between different molecules mainly dominates the deformation of the PU matrix. The molecules in the original model are tightly compacted together, and the relative slippage between them would occur once the strain was applied (Figure 4(d)), and this would increase the average spacing between different molecules and thus change the van der Waals interactions (Figure 6). The increased molecular spacing also provides some space for the restoration of the deformed bond angles; therefore, the potential energy associated with the bond angles would decline gradually with the increasing strain. When the model is deformed to a larger size, the average distance between different molecules would be increased, which would result in the less and less pronounced van der Waals interactions (Figure 6). At a lower strain rate, the amount of slippage per unit time between different molecules is smaller, and the time for the variation of bond angles is longer. The bond angles have more time to revert to their equilibrium, and thus, the PEangle begins to decline as soon as the strain appears. However, most of the bond angles are still not in their equilibrium state. This can be illustrated by Figure 7, from which one can see that even at εT of 0.60, for these 6 types of different bond angles, over 90% of them were smaller than the equilibrium bond angles. Therefore, the PEangle term continuously and steadily declines.Figure 7
Probability fractions of values of different bond angles of the 100-chain coarse-grained PU atεT of 0.60. Over 90% of the bond angles were smaller than the equilibrium bond angles (θ0).To find out which types of bond angles reverted the most, the evaluations ofPEangle of each type of bond angle of the model at a tension rate of 2×1011/s were calculated, and they are drawn in Figure 8.Figure 8
Evolution of total potential energy of different types of bond angles for a coarse-grained PU matrix containing 100 molecular chains for a strain rate of2×1011/s.From Figure8, one can see that the potential energy of bond angles A-A-B, A-B-C, B-C-A, C-A-B, and C-A-A varied little under the tension, whereas that of bond angles A-A-A and B-C-A declined gradually with the increasing strain when the strain was over 0.10. The variations of PEangle of bond angle A-A-A and bond angle B-C-A correspond well with the result shown in Figure 6(a). Moreover, even the B-C-A bond angles have the minimum Kangle than the others, the total change in PEangle of A-A-A bond angles is more significant than that of the B-C-A bond angles. For example, at εT of 0.40, the total potential energy of the A-A-A bond angles decreased by about 2400 kcal/mol, which is over 6 times that of B-C-A bond angles (about 380 kcal/mol). This is because the total number of the A-A-A bond angles is the maximum (190 per PU chain); therefore, the total variation of potential energy of the A-A-A bond angles is the most significant. This implies that the change in PEangle depends on both the bond angle stiffness and the total number of bond angles [35]. The A-A-A bond angles are mapped to the bond angles in the soft segment in the all-atom model (Figure 2(a)), which implies that the deformation of the polyurethane molecule mainly occurs in the soft segment.


### 3.3. Mechanical Behavior under Successive Cyclic Tension Loading

Figure9 shows the stress-strain curve of the coarse-grained PU matrix under successive tensile loading-unloading cycles. It can be seen that, generally, the contour of the stress-strain curve of the coarse-grained PU sustaining the successive loading-unloading cycles follows that of the tensile stress-strain curve. The strength after each reloading was almost equal to the tensile strength at the same strain. After having passed through the unloading point, the reloading curves precisely follow the original stress-strain curve, which is known as “stress memory” behavior. The “stress memory” phenomenon has also been observed in other polymeric materials like high-density polyethylene [25]. This observation is consistent with the experimental results [51], and it indicates that the strength of polyurethane under tension is independent of the loading path.Figure 9
The stress-strain curve of the coarse-grained PU matrix under successive cyclic tensile loading (the colored one). The stress-strain curve under tensile was also drawn as the black one for comparison. The strain rate is5×1010/s.Figure9 also shows that, before yielding, the stress curves of unloading and subsequent reloading change linearly, and they are almost coincident, but after yielding, they display the nonlinear variation law. More specifically, the stress curves of reloading increase sharply and linearly within a small strain range, and then, they turn to increase at a much lower rate up to the value of the unloading point. Taking the fourth loop as an example, the slope of the reloading stress curve was about 1850 MPa per strain when εT increased from 0.22 to 0.24, but it declined to about 90 MPa prestrain in the εT range of 0.24 to 0.26. The unloading curves vary similarly but in a reverse direction. Due to the nonlinear variation characteristic, after yielding, hysteresis loops are seen on the stress curves of all the unloading-reloading cycles. The hysteresis loops indicate that the strain of the unloading points contains the viscous strain [25]. The viscous strain (εv) of polymeric materials, calculated by subtracting the elastic strain (εe) and plastic strain (εp) from the strain of the unloading point (as shown in Figure 9), originates from the chain dynamics that are determined by the friction between polymer atoms, driven by entropic effects. When the polymer is subjected to relatively small or moderate strain, the mobility of tightly bonded molecules is constrained by the van der Waals interactions that cause friction between polymer atoms. When the deformation increased, the molecules are gradually separated in the space, and the accumulated plastic deformation causes some damage which weakens the van der Waals interactions between the molecules, thus resulting in viscous strain of the polymers. From Figure 9, it can also be read that when the strain of the unloading point increases the viscous strain in PU matrix increases gradually.


### 3.4. Mechanical Behavior under the Stress Relaxation

Figure10(a) shows the variations of the stress of the coarse-grained PU matrix with the running time at different strains of relaxation, from which one can see that when the stress relaxation begins, the axial stress declines sharply within a small time range, but not to zero. When the strain of relaxation is equal to 0.05 or 0.10, the axial stress declines linearly to a certain value, but when the strain of relaxation is equal to or higher than 0.20, the axial stress declines in a nonlinear way. Besides, when the strain of relaxation is lower than 0.20, after relaxation, the final values of the axial stress increase with the increasing strain of relaxation. For example, when the stress relaxation happened at εT of 0.20, after relaxation, the axial stress of the model was about 10 MPa, about 33% higher than that relaxed at εT of 0.10. However, when the strain of relaxation was equal to or higher than 0.40, after stress relaxation, all the axial stresses were almost equal to 15 MPa. The stress variation curves of the PU model during the stress relaxation process are similar to the experimental results of the macroscopic polyurethane [52, 53], indicating the viscoelastic behavior of the material.Figure 10
Evolutions of (a) stress with the relaxation time at different strains of relaxation and (b) potential energy with the relaxation time at the strain of relaxation of 0.40.
(a)(b)Figure10(b) presents the evolutions of the potential energy of the model during tension and subsequent stress relaxation (the strain of stress relaxation was 0.40). From Figure 10(b), it can be seen that during tension, the variations of the potential energy associated with the van der Waals interactions, bond lengths, and bond angles are similar to those shown in Figure 6. In the stress relaxation process, the PEangle term was maintained at about -2600 kcal/mol, but the PEnonbond declined gradually from about 14000 kcal/mol to 12000 kcal/mol, which lasted for around 8 ps. The PEnonbond was maintained at 12000 kcal/mol in the following time. Besides, in the stress relaxation process, the variation of PEtotal was equal to that of PEnonbond. The variation of PEnonbond and PEtotal corresponds well with that of the axial stress presented in Figure 10(b), indicating that in the stress relaxation, the variation of axial stress of the PU model was mainly dominated by the varied van der Waals interactions.According to Figure10(a), in the stress relaxation, Et, the relaxation modulus of the coarse-grained PU model, defined as the stress corresponding to unit strain, can be calculated. Et is a parameter of time dependence. Considering that the general Maxwell model illustrated in Figure 11(a) can be used to describe the stress relaxation property of viscoelastic solid materials [53, 54], Et can be expressed as the following function:
(6)Et=E∞+∑i=1NEi⋅exp−tλi,where E∞ represents the average of Et when the relaxation time (t) reaches infinity, Ei and λi represent the real-time relaxation modulus and relaxation time constant of element i in the general Maxwell model, respectively, and N is the total number of the Maxwell elements.Figure 11
(a) Schematic diagram of the general Maxwell model and (b) the comparison of the nondimension relaxation modulus between simulation results (the scatter points) and theoretical values (the dashed lines) atεR. Here, the starting time of stress relaxation was zeroed out.
(a)(b)In Equation (6), E∞ can be calculated using the following function:
(7)E∞=E0−∑i=1NEi,where E0 represents the relaxation modulus at the initial time.Therefore, Equation (6) can be rewritten as the following:
(8)Et=E0−∑i=1NEi⋅1−exp−tλi.To analyze the relaxation modulus,Et was transformed to the dimensionless form by dividing Equation (8) by E0, and the nondimension relaxation modulus can be expressed as the following:
(9)et=1−∑i=1Nei⋅1−exp−tλi.The variations ofet of the model in the stress relaxation at different εT are plotted in Figure 11(b) as the scatter points. The time at the beginning of stress relaxation in Figure 10(a) has been shifted to zero. The variations of et are fitted with Equation (9), and it was found that the best fit can be obtained for the model at different εT when N is equal to 1. Therefore, the viscoelastic model for PU matrix during the stress relaxation process is described by the following:
(10)et=1−e⋅1−exp−tλ.The parameters of the viscoelastic model and the goodness of fit (R2) are listed in Table 4, and the theoretical values of et calculated by Equation (10) have been drawn in Figure 11(b). It can be seen that even though a difference exists, the theoretical values of et are close to the simulation results. The theoretical et declines gradually with the increasing relaxation time but not to zero. With the increasing strain of relaxation, the decline of theoretical et becomes more nonlinear. The theoretical stress value can be obtained by multiplying the theoretical et with εT. Since εT is constant, it can be deduced that the theoretical stress will vary according to the same law as the theoretical et. Therefore, the viscoelastic model consisting of one elastic element in parallel with one Maxwell element can be used to describe the mechanical behavior of an amorphous polyurethane material under stress relaxation.Table 4
Parameters of the viscoelastic model for coarse-grained PU matrix at the different strains of stress relaxation.
εTeλR20.050.8120.2380.7510.100.7090.2170.7700.200.5890.3030.7830.400.5490.4460.7290.600.5210.9540.7700.800.5350.9160.735In practical engineering, polymer grouting material refers to the product of the reaction between isocyanate and polyols with a constant mass ratio [5, 9]. Therefore, in this study, only the influences of different loading types were investigated, and the influences of different raw materials or mass ratios on the mechanical properties of PU matrix will be studied in the further study.


## 3.1. Effects of the Number of PU Chains on the Tensile Mechanical Properties

Figures4(a)–4(c) show the variation laws of the stress of the coarse-grained PU models with different numbers of molecular chains under uniaxial tension, from which one can see that when εT was lower than 0.03, the stress increased almost linearly with the increasing strain, corresponding to the elastic deformation behavior. When εT was over 0.03, the variation rate of the tensile stress started to decline, implying that plastic deformation had occurred to the model, and at εT of about 0.20 (yield strain), the tensile stress achieved a plateau. It can be found that the number of molecular chains has little influence on the yield strain or yield strength (both are about 0.20 and 30 MPa). The stress-strain curve becomes smoother with a higher number of PU chains, suggesting that the simulations moved closer to the bulk behavior. Therefore, in all the subsequent simulations, the coarse-grained PU model with 100 molecular chains was used.Figure 4
(a–c) Stress-strain curves of coarse-grained PU models with different numbers of chains under tensile at a strain rate of5×1010/s and (d) the evolution of molecular configuration containing 10 PU chains at different strains.
(a)(b)(c)(d)The snapshots of the deformation process of the configurations containing 10 PU molecule chains are presented in Figure4(d). It can be seen that the deformation of the model mainly consists of two parts: the slippage between the PU chains and the stretching of PU chains themselves. Before the yielding, the stretching deformation of each independent PU molecule chain was less significant than the slippage between them, whereas it became more and more obvious with the increasing strain after yielding. The slippage deformation occurred mainly along the direction of tensile strain.


## 3.2. Effects of Strain Rate on the Tensile Mechanical Properties

Figure5 shows the variations of the stress of the coarse-grained PU matrix containing 100 molecule chains under uniaxial tension at the strain rate of 2×1011/s and 1×1010/s. The curves follow the same laws as discussed in Section 3.1 except for the different yield strengths and different yield strains. According to Figures 4(c) and 5, one can see that the yield strength of the model declines with the declined strain rate. At the strain rates of 2×1011/s, 5×1010/s, and 1×1010/s, the yield strength was about 70, 30, and 16 MPa, respectively. On the contrary, with the decreasing strain rate, the yield strain increased gradually. For example, the yield strain was about 0.20 and 0.22 at the strain rate of 5×1010/s and 1×1010/s, respectively, which were about double that at the strain rate of 2×1011/s. Moreover, it can also be seen that when a higher strain rate was adopted, the stress-strain curves became smoother.Figure 5
Stress-strain curves of the coarse-grained PU matrix containing 100 chains under tensile tension at the strain rate of (a)2×1011/s and (b) 1×1010/s.
(a)(b)Figure6 illustrates the evolution of potential energy (PEtotal) of the coarse-grained PU matrix at different tensile strain rates, corresponding to the stress-strain behavior presented in Figures 4(c) and 5. The evolutions of individual components, i.e., the bond potential energy (PEbond), bond angle potential energy (PEangle), and van der Waals interactions (PEnonbond), are also included in Figure 6. E0 represents the potential energy at initial time.Figure 6
Variations of the potential energy of the coarse-grained PU matrix versus true strain for a strain rate of (a)2×1011/s, (b) 5×1010/s, and (c) 1×1010/s.
(a)(b)(c)From Figure6(a), it can be seen that when the axial strain increased at a rate of 2×1011/s, both PEtotal and PEnonbond increased gradually, and the increasing trends are curves with a downward opening, and that the PEtotal curve and PEnonbond curves almost coincided until εT reached to about 0.10 (the yield strain). PEtotal was about 6000 kcal/mol at the yield strain, and after yielding, PEtotal became lower than PEnonbond; meanwhile, the PEangle term started to decline in an almost linear way. PEangle was about -6000 kcal/mol at εT of 0.70. During the tensile procedure, PEbond increased with the increasing strain, but its variation was lower than those of the other two terms (with the maximum increment of lower than 1500 kcal/mol).The changedPEnonbond, PEangle, and PEbond are primarily attributed to the slippage between different molecules, bond angle bending behavior, and bond stretching behavior, respectively. Figure 6(a) shows that, at the strain rate of 2×1011/s, under uniaxial tension, only the slippage between different PU molecules occurred before yielding, but after yielding, in addition to the slippage between different PU molecules that still remained, on each PU molecule, the bond angles started to bend, and the bonds were stretched to synthetically accommodate deformation. The difference between PEnonbond and PEtotal after yielding was dissipated by the bond angle bending behavior and bond stretching behavior, but primarily via the former. Figures 6(b) and 6(c) show that whether at the strain rate of 5×1010/s or 1×1010/s, there was little change in the overall evolution laws of PEtotal and PEnonbond, but PEangle started to decline from the very beginning of the simulation, and that PEbond was always equal to zero. This indicates that at a lower strain rate, as soon as the tensile load was applied, the bond angles displayed bending behavior, but the bond lengths were little influenced.Figure6 also shows that the change of total potential energy at yielding is independent of the strain rate. At the three loading rates, both PEtotal at yielding are about 6000 kcal/mol, but after yielding, the PEtotal of the model sustaining a lower strain rate has a lower value, which is due to that PEnonbond has significantly declined. With the decreasing strain rate, the increased yield strain (Figures 4(c) and 5) leads to the increased volume strain of the model. Considering that the yield strength of the model is positively dependent on the potential energy but negatively influenced by the volume [50], when the strain rate declines, PEtotal is almost unchanged, but the increased volume strain of the coarse-grained PU models at yielding is the main reason for the lower yield strength.Figures4(d) and 6 show that the slippage between different molecules mainly dominates the deformation of the PU matrix. The molecules in the original model are tightly compacted together, and the relative slippage between them would occur once the strain was applied (Figure 4(d)), and this would increase the average spacing between different molecules and thus change the van der Waals interactions (Figure 6). The increased molecular spacing also provides some space for the restoration of the deformed bond angles; therefore, the potential energy associated with the bond angles would decline gradually with the increasing strain. When the model is deformed to a larger size, the average distance between different molecules would be increased, which would result in the less and less pronounced van der Waals interactions (Figure 6). At a lower strain rate, the amount of slippage per unit time between different molecules is smaller, and the time for the variation of bond angles is longer. The bond angles have more time to revert to their equilibrium, and thus, the PEangle begins to decline as soon as the strain appears. However, most of the bond angles are still not in their equilibrium state. This can be illustrated by Figure 7, from which one can see that even at εT of 0.60, for these 6 types of different bond angles, over 90% of them were smaller than the equilibrium bond angles. Therefore, the PEangle term continuously and steadily declines.Figure 7
Probability fractions of values of different bond angles of the 100-chain coarse-grained PU atεT of 0.60. Over 90% of the bond angles were smaller than the equilibrium bond angles (θ0).To find out which types of bond angles reverted the most, the evaluations ofPEangle of each type of bond angle of the model at a tension rate of 2×1011/s were calculated, and they are drawn in Figure 8.Figure 8
Evolution of total potential energy of different types of bond angles for a coarse-grained PU matrix containing 100 molecular chains for a strain rate of2×1011/s.From Figure8, one can see that the potential energy of bond angles A-A-B, A-B-C, B-C-A, C-A-B, and C-A-A varied little under the tension, whereas that of bond angles A-A-A and B-C-A declined gradually with the increasing strain when the strain was over 0.10. The variations of PEangle of bond angle A-A-A and bond angle B-C-A correspond well with the result shown in Figure 6(a). Moreover, even the B-C-A bond angles have the minimum Kangle than the others, the total change in PEangle of A-A-A bond angles is more significant than that of the B-C-A bond angles. For example, at εT of 0.40, the total potential energy of the A-A-A bond angles decreased by about 2400 kcal/mol, which is over 6 times that of B-C-A bond angles (about 380 kcal/mol). This is because the total number of the A-A-A bond angles is the maximum (190 per PU chain); therefore, the total variation of potential energy of the A-A-A bond angles is the most significant. This implies that the change in PEangle depends on both the bond angle stiffness and the total number of bond angles [35]. The A-A-A bond angles are mapped to the bond angles in the soft segment in the all-atom model (Figure 2(a)), which implies that the deformation of the polyurethane molecule mainly occurs in the soft segment.


## 3.3. Mechanical Behavior under Successive Cyclic Tension Loading

Figure9 shows the stress-strain curve of the coarse-grained PU matrix under successive tensile loading-unloading cycles. It can be seen that, generally, the contour of the stress-strain curve of the coarse-grained PU sustaining the successive loading-unloading cycles follows that of the tensile stress-strain curve. The strength after each reloading was almost equal to the tensile strength at the same strain. After having passed through the unloading point, the reloading curves precisely follow the original stress-strain curve, which is known as “stress memory” behavior. The “stress memory” phenomenon has also been observed in other polymeric materials like high-density polyethylene [25]. This observation is consistent with the experimental results [51], and it indicates that the strength of polyurethane under tension is independent of the loading path.Figure 9
The stress-strain curve of the coarse-grained PU matrix under successive cyclic tensile loading (the colored one). The stress-strain curve under tensile was also drawn as the black one for comparison. The strain rate is5×1010/s.Figure9 also shows that, before yielding, the stress curves of unloading and subsequent reloading change linearly, and they are almost coincident, but after yielding, they display the nonlinear variation law. More specifically, the stress curves of reloading increase sharply and linearly within a small strain range, and then, they turn to increase at a much lower rate up to the value of the unloading point. Taking the fourth loop as an example, the slope of the reloading stress curve was about 1850 MPa per strain when εT increased from 0.22 to 0.24, but it declined to about 90 MPa prestrain in the εT range of 0.24 to 0.26. The unloading curves vary similarly but in a reverse direction. Due to the nonlinear variation characteristic, after yielding, hysteresis loops are seen on the stress curves of all the unloading-reloading cycles. The hysteresis loops indicate that the strain of the unloading points contains the viscous strain [25]. The viscous strain (εv) of polymeric materials, calculated by subtracting the elastic strain (εe) and plastic strain (εp) from the strain of the unloading point (as shown in Figure 9), originates from the chain dynamics that are determined by the friction between polymer atoms, driven by entropic effects. When the polymer is subjected to relatively small or moderate strain, the mobility of tightly bonded molecules is constrained by the van der Waals interactions that cause friction between polymer atoms. When the deformation increased, the molecules are gradually separated in the space, and the accumulated plastic deformation causes some damage which weakens the van der Waals interactions between the molecules, thus resulting in viscous strain of the polymers. From Figure 9, it can also be read that when the strain of the unloading point increases the viscous strain in PU matrix increases gradually.


## 3.4. Mechanical Behavior under the Stress Relaxation

Figure10(a) shows the variations of the stress of the coarse-grained PU matrix with the running time at different strains of relaxation, from which one can see that when the stress relaxation begins, the axial stress declines sharply within a small time range, but not to zero. When the strain of relaxation is equal to 0.05 or 0.10, the axial stress declines linearly to a certain value, but when the strain of relaxation is equal to or higher than 0.20, the axial stress declines in a nonlinear way. Besides, when the strain of relaxation is lower than 0.20, after relaxation, the final values of the axial stress increase with the increasing strain of relaxation. For example, when the stress relaxation happened at εT of 0.20, after relaxation, the axial stress of the model was about 10 MPa, about 33% higher than that relaxed at εT of 0.10. However, when the strain of relaxation was equal to or higher than 0.40, after stress relaxation, all the axial stresses were almost equal to 15 MPa. The stress variation curves of the PU model during the stress relaxation process are similar to the experimental results of the macroscopic polyurethane [52, 53], indicating the viscoelastic behavior of the material.Figure 10
Evolutions of (a) stress with the relaxation time at different strains of relaxation and (b) potential energy with the relaxation time at the strain of relaxation of 0.40.
(a)(b)Figure10(b) presents the evolutions of the potential energy of the model during tension and subsequent stress relaxation (the strain of stress relaxation was 0.40). From Figure 10(b), it can be seen that during tension, the variations of the potential energy associated with the van der Waals interactions, bond lengths, and bond angles are similar to those shown in Figure 6. In the stress relaxation process, the PEangle term was maintained at about -2600 kcal/mol, but the PEnonbond declined gradually from about 14000 kcal/mol to 12000 kcal/mol, which lasted for around 8 ps. The PEnonbond was maintained at 12000 kcal/mol in the following time. Besides, in the stress relaxation process, the variation of PEtotal was equal to that of PEnonbond. The variation of PEnonbond and PEtotal corresponds well with that of the axial stress presented in Figure 10(b), indicating that in the stress relaxation, the variation of axial stress of the PU model was mainly dominated by the varied van der Waals interactions.According to Figure10(a), in the stress relaxation, Et, the relaxation modulus of the coarse-grained PU model, defined as the stress corresponding to unit strain, can be calculated. Et is a parameter of time dependence. Considering that the general Maxwell model illustrated in Figure 11(a) can be used to describe the stress relaxation property of viscoelastic solid materials [53, 54], Et can be expressed as the following function:
(6)Et=E∞+∑i=1NEi⋅exp−tλi,where E∞ represents the average of Et when the relaxation time (t) reaches infinity, Ei and λi represent the real-time relaxation modulus and relaxation time constant of element i in the general Maxwell model, respectively, and N is the total number of the Maxwell elements.Figure 11
(a) Schematic diagram of the general Maxwell model and (b) the comparison of the nondimension relaxation modulus between simulation results (the scatter points) and theoretical values (the dashed lines) atεR. Here, the starting time of stress relaxation was zeroed out.
(a)(b)In Equation (6), E∞ can be calculated using the following function:
(7)E∞=E0−∑i=1NEi,where E0 represents the relaxation modulus at the initial time.Therefore, Equation (6) can be rewritten as the following:
(8)Et=E0−∑i=1NEi⋅1−exp−tλi.To analyze the relaxation modulus,Et was transformed to the dimensionless form by dividing Equation (8) by E0, and the nondimension relaxation modulus can be expressed as the following:
(9)et=1−∑i=1Nei⋅1−exp−tλi.The variations ofet of the model in the stress relaxation at different εT are plotted in Figure 11(b) as the scatter points. The time at the beginning of stress relaxation in Figure 10(a) has been shifted to zero. The variations of et are fitted with Equation (9), and it was found that the best fit can be obtained for the model at different εT when N is equal to 1. Therefore, the viscoelastic model for PU matrix during the stress relaxation process is described by the following:
(10)et=1−e⋅1−exp−tλ.The parameters of the viscoelastic model and the goodness of fit (R2) are listed in Table 4, and the theoretical values of et calculated by Equation (10) have been drawn in Figure 11(b). It can be seen that even though a difference exists, the theoretical values of et are close to the simulation results. The theoretical et declines gradually with the increasing relaxation time but not to zero. With the increasing strain of relaxation, the decline of theoretical et becomes more nonlinear. The theoretical stress value can be obtained by multiplying the theoretical et with εT. Since εT is constant, it can be deduced that the theoretical stress will vary according to the same law as the theoretical et. Therefore, the viscoelastic model consisting of one elastic element in parallel with one Maxwell element can be used to describe the mechanical behavior of an amorphous polyurethane material under stress relaxation.Table 4
Parameters of the viscoelastic model for coarse-grained PU matrix at the different strains of stress relaxation.
εTeλR20.050.8120.2380.7510.100.7090.2170.7700.200.5890.3030.7830.400.5490.4460.7290.600.5210.9540.7700.800.5350.9160.735In practical engineering, polymer grouting material refers to the product of the reaction between isocyanate and polyols with a constant mass ratio [5, 9]. Therefore, in this study, only the influences of different loading types were investigated, and the influences of different raw materials or mass ratios on the mechanical properties of PU matrix will be studied in the further study.


## 4. Conclusions

The mechanical behavior of amorphous PU matrix under uniaxial tension, successive cyclic tension, and stress relaxation at nanoscale was investigated with coarse-grained MD simulation method. The conclusions are drawn as follows:(1)
The mechanical properties of PU matrix under tensile are little affected by the number of molecules but change significantly when the strain rate varies. With the increasing strain rate, the yield strength increases, and the yield strains decrease. At yield, the total change in the potential energy is independent of the strain rate, and it is about 6000 kcal/mol. When the strain rate increases, the declined volume strain of the PU matrix is the main reason for the increased yield strength(2)
Under tensile, the potential energy associated with the bond stretching varies less significantly than that associated with the bond angle bending or van der Waals interactions, indicating that the deformation of the model is mainly caused by the slippage between different molecule chains and the bending of the bond angles. The bending of bond angles mainly occurs in the soft segments of the molecule chains. After yielding, most of the bond angles are still not in equilibrium(3)
The contour of the stress-strain curve of the model under cyclic tension is consistent with that under uniaxial tension, indicating that the model performs with the stress memory behavior, and thus, its tensile strength is independent of loading path. After yielding, the stress curves of the unloading-reloading cycles contain the hysteresis loops, which imply that the strain of the unloading points contains elastic strain, plastic strain, and viscous strain(4)
During the stress relaxation process, the evolution of the axial stress is primarily caused by the varied van der Waals interactions. The mechanical behavior of the amorphous PU matrix under stress relaxation at the nanoscale can be described by the viscoelastic model consisting of one elastic element in parallel with one Maxwell element



---

*Source: 4923385-2022-09-21.xml*