Toward the long-term aging influence and novel reaction kinetics models of bitumen

ABSTRACT

This study aimed to explore the long-term aging influence on chemo-rheological properties and develop novel consecutive models for the long-term aging reaction kinetics of bitumen. The results revealed that the aging index was significantly dependent on the types of selected parameters. The Zero-order model was suitable to describe the long-term aging reaction kinetics of bitumen based on the oxygen-containing functional groups with the reaction rate constants in 0.7–3.3*10⁻⁴ (mol  L⁻¹·h⁻¹). In the SARA-based consecutive reaction model, the most optimum kinetics model for aromatic fraction was the Third-order reaction model and the corresponding reaction kinetics constant (k₁) was 0.02 (mol·L⁻¹)⁻²(h)⁻¹. The Zero-order model could well fit the generation kinetics of asphaltene molecules with the reaction rate constant k₂ of 3.85*10⁻⁴ mol·(L·h)⁻¹. Further, the transformation reaction from the resin to asphaltene molecules was the control step of the whole consecutive reaction model. In this study, when one-unit resin fraction was generated, the consumption amount of aromatic fraction was about 2.82 units. Meanwhile, when one-unit resin fraction was consumed, only 0.58-unit asphaltene could be generated. The developed reaction kinetics models could be beneficial to predict the functional groups distribution and SARA fractions in aged bitumen with different aging degrees.

KEYWORDS:

• Long-term aging
• consecutive kinetics model
• SARA fractions
• quantitative conversion relationship
• chemo-rheological characterization

1. Introduction

Asphalt pavement has been widely applied due to its excellent characteristics of flatness, comfort and flexibility (Yan et al. 2019, Qiu et al. 2020, Ma et al. 2021). However, the service life of asphalt roads is limited because of the external loading and environmental effects as well as the deterioration of material properties (Liu et al. 2019, Oldham et al. 2020). The bitumen plays an important role of binder to connect the aggregates and prevent the asphalt mixture from deformation and cracking. The aging phenomenon of bitumen significantly weakens the low-temperature and fatigue cracking resistance as well as the adhesion properties of asphalt mixture (Liu et al. 2014, Omairey et al. 2019, He et al. 2021). Pavement distresses would occur because of bitumen aging, such as the cracking, moisture damage and raveling. (Petersen and Glaser 2011, Wang et al. 2015a). Thus, it is important to understand the aging influence and mechanism of bitumen, especially for the long-term aging behaviours.

Lots of previous studies focused on the influence of aging on the chemical, physical and rheological properties of bitumen (Das et al. 2015, Yang et al. 2018). The oxidation aging would accelerate the formation of carbonyl and sulfoxide functional groups, and remarkably improve the molecular weight and surface roughness of bitumen. The aging of bitumen significantly increased the softening point, complex modulus, viscosity, rutting parameter and stiffness of bitumen, while the ductility, penetration, m-value, phase angle decreased (Qu et al. 2018). Moreover, the aging effect on the mechanical and rheological properties of bitumen is related to internal change of chemical compositions and microstructure. Lu et al. (2021) analyzed the asphaltenes and maltenes in bitumen after long-term aging, and found that the aging process promoted that some low molecular weight compounds conversed to the part of asphaltene. Moreover, maltene molecules became insoluble in n-heptane after reaction with oxygen. Jing et al. (2021) found a linear correlation between the combined aging index and the crossover map distance, and proposed a chemo-mechanics framework of bitumen aging behaviours. Xu and Wang (2017) conducted the molecular dynamics simulation technology to explore the oxidative aging effects at nanoscale and found that the aging process increased the activation energy barrier for self-healing of bitumen, but distinctly weakened the work of adhesion and the nano-agglomeration behaviour of asphaltene molecules.

The potential oxidation reaction pathway of bitumen molecules has also been studied with the reaction molecular dynamics simulation method. The reaction molecular dynamics simulation and quantum chemistry theory were utilized to detect the oxidative aging of saturate, resin, asphaltene molecules as well as the average molecule of bitumen (Pan et al. 2012, Pan and Tarefder 2016). The results revealed that the aromatic and asphaltene molecules were partially oxidized with the sulfoxidation and ketonization occurring, while the saturate molecule was quite resistant to oxidation. Moreover, the saturate alkane chain would be separated into shorter chains, and the formation of sulfoxide was earlier than the carbonyl term. In addition, Yang et al. (2021) simulated the chemical reaction pathway between 12 types of saturate, aromatic, resin and asphaltene molecules (SARA) in bitumen and molecular oxygen, and the same finding could be observed that the ketones and sulfoxides were the main products of bitumen oxidation aging, while the ketones were more easily generated on the cyclic benzyl carbons than the branched benzyl carbons. Further, it was reported that the molecular or atomic oxygen was attached to the α-C atom to form the hydroperoxides, followed by the formation of sulfoxide and carbonyl functional groups. It indicated that the hydroperioxides (C–O–O·) was the critical intermediate product to determine the occurrence of the oxidative aging of bitumen molecules. Regarding the polymerization mechanism, it was shown that it would occur when the aromatics and free radicals was formed by the ‘dehydrogenation’ oxidation reactions to form asphaltenes (Hu et al. 2020).

Regarding the kinetics model of bitumen aging, there were two periods with the declining fast-rate and constant-rate. Jin et al. (2011) recommended that the combined first-order, zero-order oxidation model was capable of modeling bitumen oxidation kinetics. Liu and Glover (2015) utilized the combined model to study the oxidation kinetics of warm mix asphalt and found that the warm mix additives showed on significant adverse effect on oxidation kinetics of bitumen. Meanwhile, the diffusion model of oxygen was also coupled with the reaction model of bitumen. Cui et al. (2018) considered the whole process that oxygen molecules penetrated into the bitumen film and then reacted with the bitumen molecules and established a diffusion-reaction balance oxidation model of bitumen binder. Additionally, Zhang et al. developed the new short-term and long-term aging models for bitumen binder to predict the viscosity of aged binder accurately (Zhang et al. 2019a, Zhang et al. 2019b). It should be mentioned that the exist oxidation reaction model is based on the functional groups variations from the Fourier transform infrared (FTIR) spectroscopy results, such as the carbonyl and sulfoxide index (Zhao et al. 2003, Wang et al. 2015b). Actually, the oxidation reactions between the oxygen and bitumen molecules include not only the generation of oxygen-containing functional groups, but also the conversion between the chemical compositions (saturate S, aromatic A, resin R and asphaltene A_s) (Fallah et al. 2019). The variation of SARA fractions during the aging process plays a vital role on determining the macroscale rheological and mechanical properties of aged bitumen. However, as far as the authors are aware that the oxidation reaction kinetics model of bitumen from the viewpoint of SARA fractions’ conversion has not studied yet, which obstructs the further understanding of long-term aging mechanism of bitumen.

The main objective of this study was to investigate the influence of long-term aging of bitumen on chemo-rheological properties and develop the reaction kinetics models. It should be noted that the reaction kinetics models play an important role on predicting the functional group distribution and SARA fractions of aged binders with different long-term aging durations. To this end, Figure 1 illustrates the research protocol. Firstly, the aged binders with different long-term aging degrees were prepared. Afterward, the aging influence on the variation of chemical compositions, functional groups and rheological properties of bitumen were investigated. Moreover, the long-term aging sensitivity of different chemical and rheological parameters were analyzed and compared. Importantly, the functional groups-based and SARA fractions-based long-term aging kinetics models of bitumen were established and discussed. Lastly, the quantitative conversion relationship between the aromatic, resin and asphaltene fractions in bitumen during long-term aging was speculated.

Figure 1. The research protocol of this study.

2. Materials and methods

2.1. Materials

The virgin bitumen used in the study is the Total PEN 70/100, which is commonly applied in the Netherlands. Table 1 presents the physical properties and elements concentration of the virgin bitumen.

Table 1. The physical and chemical properties of PEN 70/100.

Properties	Unit	PEN 70/100	Test standard
Penetration at 25°C	1/10 mm	91	ASTM D5
Softening point	°C	48	ASTM D36
135°C Dynamic viscosity	Pa·s	0.8	AASHTO T316
Density at 25°C	g/cm^3	1.017	EN 15326
Density at 60°C	g/cm^3	0.996
Carbon C	wt%	84.06	ASTM D 7343
Hydrogen H	wt%	10.91
Oxygen O	wt%	0.62
Sulphur S	wt%	3.52
Nitrogen N	wt%	0.90

2.2. Aging procedures and sample preparation

To simulate the short-term and long-term aging of bitumen, the thin film oven test (TFOT) (ASTM D1754) and pressure aging vessel (PAV) (AASHTO R28) were employed to manufacture the aged binders with different aging degrees. The virgin bitumen was subjected to the TFOT test at 163°C for 5 h to obtain the short-term aged bitumen. Afterward, TFOT-aged samples were further aged in the PAV device. The temperature and pressure for PAV test was 100°C and 2.1 MPa, and the aging time of PAV test was selected as 20, 40 and 80 h to prepare different long-term aged bitumen. For convenience, the virgin, short-term aged and long-term aged for 20, 40 and 80 h were simplified as VB, SAB, LAB20, LAB40 and LAB80, respectively.

2.3. Characterization methods

2.3.1. Thin-layer chromatography with flame ionization detection (TLC-FID)

The bitumen is the residue material derived from the distillation procedure of crude oil, and it is an extremely complex material and composed by thousands of different molecules (Xu and Wang 2017, Fallah et al. 2019, Liu et al. 2021). Based on the solubility and polarity difference in solvents, bitumen is further divided into different species, including the saturate, aromatic, resin and asphaltene constituents (SARA fractions). In this study, the SARA fractions distributions of virgin and different aged bitumen were detected with the TLC-FID technology (IP 469 2001). Firstly, a solution involving 0.1 g bitumen and 10 ml toluene was prepared, which then went through a silica chrome rod with three types of solvents: n-heptane, toluene/n-heptane (80:20) and dichloromethane/methanol (95:5). In detail, the saturate constituents were separated chromatographically in n-heptane, while the aromatic and resin fractions could be further detached according to the solubility difference in toluene/n-heptane and dichloromethane/methanol solvents. Lastly, the asphaltene components were the residual in the sampling spot. After the separation procedure, the chromatogram columns of the chrome rod were monitored through the FID device to determine the mass fractions of SARA components in virgin and aged bitumen. All tests were conducted at least three repetitions for each sample.

2.3.2. Attenuated total reflectance-Fourier transform infrared (ATR-FTIR) spectroscopy

It was reported that the functional groups distribution in bitumen was strongly associated with its physical and rheological properties (Liu et al. 2019, He et al. 2021). In this study, the influence of short-and-long terms aging on the functional groups distribution of bitumen was detected with the ATR-FTIR spectrometer from PerkinElmer together with a single-point ATR fixture (Waltham, MA, USA). The scanning range of wavenumber was 600–4000 cm⁻¹, and the number of scan cycles was 12 times with a fixed instrument resolution of 4 cm⁻¹.

2.3.3. Rheological tests

The rotational viscometer (RV) was employed to test the dynamic viscosity of virgin and aged bitumen at the temperature region of 80–170°C with an increment of 10°C (AASHTO T316-13). The spinning speed of rotor was 20 rad/s, and the holding time for each temperature was longer than 30 min to ensure the data stability.

The frequency sweep tests were conducted with a dynamic shear rheometer (DSR) to measure the complex modulus G* and phase angle δ of virgin and aged bitumen. The frequency sweep region increased from 0.1 to 100 rad/s, and the temperature varied from 0 to 70°C with an increment of 10°C (AASHTO T315).

3. Results and discussion

3.1. SARA fractions distribution

The performance loss of bitumen during the aging process is significantly due to the variation of chemical compositions. Figure 2 shows the SARA fractions distributions of virgin and different aged bitumen. The saturate dosage in all binders is the lowest, followed by the asphaltene and resin components, while the aromatic content is the largest except when the long-term aging time exceeds 40 h. The saturate dosage in virgin and aged bitumen keeps constant, indicating that aging has no obvious influence on the saturate concentration in bitumen. However, the aromatic, resin and asphaltene fractions in bitumen show a significant dependent on the aging degree, especially during the long-term aging procedure. With the aging degree deepens, the aromatic fraction dosage decreases, while the resin and asphaltene fractions increase dramatically. It implies that the oxidative aging would accelerate the inter-substance transformation from aromatic to resin and asphaltene fractions. However, the conversion relationship and aging kinetics between the SARA fractions in bitumen during the long-term aging process are still unclear, which will be further studied and discussed in following sections.

Figure 2. Aging effects on the SARA fractions and colloidal index CI of bitumen.

From Figure 2(a), the generation rate of resin and asphaltene fractions are similar, which are both lower than the consumption rate of aromatic component. It should be mentioned that the generation and consumption reactions for resin fractions occur simultaneously, hence its aging mechanism and kinetics are complex. In addition, when the long-term aging duration is 80 h, the resin fraction is larger than aromatic and asphaltene fractions. It indicates that the generation rate of resin molecules is faster than its consumption rate, which is associated with the stoichiometric numbers of aromatic and resin molecules in the oxidation reaction formula. Overall, during the long-term aging procedure, the ratio of heavy-weight components (resin and asphaltene) continues to increase, while the light-weight aromatic fraction decreases distinctly. Table 2 lists the corresponding changes of high-polarity (resin and asphaltene) and light-molecular weight (saturate and aromatic) fractions. After the short-term and long-term aging for 20, 40 and 80 h, the light-weight fraction decreases by 3.0%, 16.5%, 25.3% and 36.4%, while the heavy-polarity component increases by 3.9%, 21.8%, 33.4% and 48.0%, respectively. When the long-term aging time is 80 h, the polarity ratio increases from 0.757 to 1.762. In summary, the long-term aging significantly alters the relative ratio of light/heavy weight and low/high polarity components, which breaks the original colloidal balance and promotes the agglomeration of high-polarity fractions for reducing the polar molecule surface energy and adapting the new solvent environment with limited saturate and aromatic fractions.

Table 2. The light and heavy-weight fractions in bitumen.

Binder	S + A (wt%)	R + As (wt%)	Polarity ratio (R + As)/(S + A)
VB	56.9	43.1	0.757
SAB	55.2	44.8	0.812
LAB20	47.5	52.5	1.105
LAB40	42.5	57.5	1.353
LAB80	36.2	63.8	1.762

The parameter of colloidal index is normally utilized to evaluate the variation of colloidal structure of bitumen, which is calculated as follows: (1) $\mathrm{C}\mathrm{o}\mathrm{l}\mathrm{l}\mathrm{o}\mathrm{i}\mathrm{d}\mathrm{a}\mathrm{l}\mathrm{i}\mathrm{n}\mathrm{d}\mathrm{e}\mathrm{x}\mathrm{C}\mathrm{I}={\displaystyle \frac{R+A}{A_S+S}}$(1) where R, A, A_s and S represent the mass fraction of resin, aromatic, asphaltene and saturate, respectively. Figure 2(b) demonstrates the CI values of virgin and aged binders, and the bitumen with higher CI value would show the more stable colloidal structure. It can be found that with the increase of aging degree, the colloidal index of bitumen remarkably decreases, which drops down to about half when the long-term aging time is 80 h. Overall, the colloidal structure of bitumen tends to be more unstable after aging. From the viewpoint of colloid chemistry, the unstable colloidal structure is mainly due to the imbalance of micelles (asphaltene and resin) and dispersion medium (aromatic and saturate). The increase of micelle number and decrease of dispersion medium both accelerate the micelle aggregation and settlement.

3.2. Functional groups analysis

The aging of bitumen not only changes the SARA fractions distribution, but also introduces the oxygen-containing functional groups into aromatic, resin and asphaltene molecules, such as the carbonyl and sulfoxide groups (Pan et al. 2012, Qu et al. 2018, Fallah et al. 2019). The depolymerization and condensation polymerization reactions of bitumen molecules would occur simultaneously, influencing the saturation and aromaticity of bitumen molecules significantly. The variation of functional groups distribution in bitumen is closely related to its rheological and viscoelastic behaviours, and it is essential to detect the chemical structure changes of bitumen molecules during the aging process to further establish more exact molecular model of aged bitumen.

Figure 3 depicts the normalized FTIR results of virgin, short-term aged and long-term aged binders. The virgin bitumen has several characteristic peaks at 2952, 2862, 1600, 1460, 1375, 864, 814 and 743 cm⁻¹, and the corresponding functional groups of these peaks are listed in Table 3. During the aging process, the peak absorbance at 1700, 1600 and 1030cm⁻¹ increases dramatically, indicating that the amount of C = O, C = C and S = O functional groups in bitumen rise gradually. To quantitatively assess the aging influence on the chemical groups in bitumen molecules, the Carbonyl Index (CI), Sulfoxide Index (SI), Aromatic Index (ARI) and Aliphatic Index (AII) are proposed here: (2) $\mathrm{C}\mathrm{a}\mathrm{r}\mathrm{b}\mathrm{o}\mathrm{n}\mathrm{y}\mathrm{l}\mathrm{i}\mathrm{n}\mathrm{d}\mathrm{e}\mathrm{x}\mathrm{C}\mathrm{I}={\displaystyle \frac{A_{1700}}{\sum A}}$(2)

(3) $\mathrm{S}\mathrm{u}\mathrm{l}\mathrm{f}\mathrm{o}\mathrm{x}\mathrm{i}\mathrm{d}\mathrm{e}\mathrm{i}\mathrm{n}\mathrm{d}\mathrm{e}\mathrm{x}\mathrm{S}\mathrm{I}={\displaystyle \frac{A_{1030}}{\sum A}}$(3) (4) $\mathrm{A}\mathrm{r}\mathrm{o}\mathrm{m}\mathrm{a}\mathrm{t}\mathrm{i}\mathrm{c}\mathrm{i}\mathrm{n}\mathrm{d}\mathrm{e}\mathrm{x}\mathrm{A}\mathrm{R}\mathrm{I}={\displaystyle \frac{A_{1600}}{\sum A}}$(4) (5) $\mathrm{A}\mathrm{l}\mathrm{i}\mathrm{p}\mathrm{h}\mathrm{a}\mathrm{t}\mathrm{i}\mathrm{c}\mathrm{i}\mathrm{n}\mathrm{d}\mathrm{e}\mathrm{x}\mathrm{A}\mathrm{I}\mathrm{I}={\displaystyle \frac{A_{1460}+A_{1376}}{\sum A}}$(5) (6) $\sum \mathrm{A}={\mathrm{A}}_{\left(2952,2862\right)}+{\mathrm{A}}_{1700}+{\mathrm{A}}_{1600}+{\mathrm{A}}_{1460}+{\mathrm{A}}_{1375}+{\mathrm{A}}_{1030}+{\mathrm{A}}_{864}+{\mathrm{A}}_{814}+{\mathrm{A}}_{743}+{\mathrm{A}}_{724}$(6) where $\sum A$ is the sum of peak areas, and A_1700, A_1030, A₁₆₀₀, A₁₄₆₀ and A₁₃₇₅ is the peak area at the wavenumber position of 1700, 1030, 1060, 1460 and 1376 cm⁻¹, respectively.

Figure 3. The normalized FTIR spectra of virgin and aged bitumen.

Table 3. The characteristic peaks of different functional groups.

Wavenumber (cm^-1)	Functional groups
2952 and 2862	Asymmetric and symmetric stretching vibrations of C–H on aliphatic hydrogen
1700	Stretching vibration of C = O functional group
1600	Stretching vibration of C = C bond in aromatic ring
1460	Bending vibration of C–H bond methylene
1375	Bending vibration of C–H bond in methyl
1030	Stretching vibration of S = O functional group
864, 814 and 743	Bending vibration of C–H bond on the substituents in aromatic rings

Figure 4 depicts the CI and SI parameters of bitumen binders, which both have an obvious increasing trend as the aging duration prolonging. It implies that the amount of C = O and S = O functional groups in bitumen molecules enlarges. In addition, the generation stage of sulfoxide group is earlier than the carbonyl group, which is due to the lower formation energy of S = O bond. The influence of short-term aging on sulfoxide index is larger than carbonyl index, while the carbonyl groups is more dependent on the long-term aging duration than sulfoxide groups because of the limited amount of sulfur atoms in bitumen molecules. When the long-term aging time is shorter than 80 h, the carbonyl index of aged bitumen is lower than sulfoxide index. In following Section 4, the reaction rate of carbonyl and sulfoxide index will be further analyzed to explore the functional groups-based long-term aging kinetics of bitumen.

Figure 4. The carbonyl index CI and sulfoxide index SI of virgin and aged binders.

The aromaticity and aliphatic indices of virgin and aged binders are shown in Figure 5. As the aging degree depeening, the aromaticity index of bitumen significantly increases, while the aliphatic index decreases correspondingly. In summary, during the aging process, the aromaticity and saturation of bitumen molecules increase and decrease, respectively. The conclusion is in good agreement with the SARA fractions distribution results. The increased asphaltene and resin molecules show the higher aromaticity and lower aliphatic degree than the reduced aromatic fraction. The generation of polycyclic aromatic hydrocarbons and consumption of aliphatic hydrocarbons both result in the improvement of molecular interaction and bitumen stiffness.

Figure 5. The aromaticity and aliphatic indices of virgin and aged binders.

3.3. Viscosity and workability

The influence of short-and-long terms aging on the viscosity and workability of bitumen is estimated using the rotational viscometer. Figure 6(a) displays the dynamic viscosity values of virgin and aged binders at different testing temperatures. As expected, the viscosity of all binders decreases as the testing temperature rising. The reason is that at high temperature, the molecular energy and Brownian motion both intensify, while the molecular interaction and intermolecular friction weaken. Meanwhile, with the aging time prolongs, the viscosity of bitumen significantly increases, which leads to the worse workability and insufficient adhesion to aggregates. From the perspective of molecular level, the aging process promotes the increase of polar functional groups and molecules, which is beneficial to decrease the molecular distance and enhance the intermolecular force. Moreover, the reduction of light-weight fractions would shrink the free volume in bitumen system, which increases the difficulty for the movement of bitumen molecules.

Figure 6. The influence of aging on the viscosity of bitumen at various temperatures.

There is a linear relationship between the logarithm value of viscosity (Lnη) and reciprocal of temperature (1/T). With the 1/T value increases, the Lnη value raises linearly. The correlation equations of virgin and aged binders are recorded in Figure 6(b). As the aging degree deepens, the slope value of bitumen in correlation equation increases gradually, which means the viscosity of aged bitumen with higher aging degree is more sensitive to temperature variation.

The Arrhenius formula is utilized to describe the viscosity-temperature correlation: (7) $\mathrm{L}\mathrm{n}\left(\eta \right)={\displaystyle \frac{E_{\mathrm{v}\mathrm{i}\mathrm{s}}}{RT}+\mathrm{L}\mathrm{n}A}$(7) where η is the rotational viscosity of bitumen, Pa·s; T refers to the testing temperature, K; E_vis and A represents the flow activation energy and pre-exponential factor, respectively; and R is the gas constant, 8.314 J/(mol·K). By definition, the flow activation energy E_vis is associated with the required energy for bitumen molecules, and the sample with larger E_vis value requires more energy for flow.

Table 4 lists the E_vis and A values of virgin and aged binders. After the bitumen is subjected to short-term and long-term aging for 20, 40 and 80 h, the flow activation energy increases by 13.9%, 24.5%, 29.9% and 57.7%, respectively. Therefore, the flow hindrance of bitumen markedly improves during the aging process, which deteriorates the workability and the high temperature is required to ensure the satisfied adhesion of bitumen on aggregates (Lin et al. 2021, Ren et al. 2021).

Table 4. The flow energy and pre-exponential factor of virgin and aged binders.

	VB	SAB	LAB20	LAB40	LAB80
Flow activation energy Evis (J)	64489	73430	80300	83742	101681
Pre-exponential factor A	3.14E-9	6.15E-10	5.98E-11	8.65E-11	1.36E-12

3.4. Rheological properties

The rheological parameters, complex modulus G* and phase angle δ, of virgin and aged binders are evaluated. The corresponding master curves are established with a generalized logistic formula and Williams-Landel-Ferry (WLF) function, shown in Equations (8) and (9): (8) $\mathrm{l}\mathrm{o}\mathrm{g}P=v+{\displaystyle \frac{\alpha }{{\left[1+\lambda e^{\left(\text{β}+\gamma \left(\mathrm{l}\mathrm{o}\mathrm{g}\omega \right)\right)}\right]}^{1/\lambda }}}$(8) (9) $\mathrm{l}\mathrm{o}\mathrm{g}{\alpha }_T\left(T\right)={\displaystyle \frac{-C_1\left(T-T_f\right)}{C_2+\left(T-T_f\right)}}$(9) where P refers the complex modulus or phase angle of bitumen; ω shows the loading frequency; v is the lower asymptote; α represents the difference between the values of upper and lower asymptote. Moreover, ${\alpha }_T\left(T\right)$ denotes the shifting factor; T and T_f are the test temperature and reference temperature, respectively; and β, γ, λ, C₁ and C₂ are all constants.

Figure 7 shows the master curves of G* and δ for virgin and aged binders with different aging degrees. As expected, the G* and δ parameters of bitumen show significant frequency dependence. As the frequency raising, the G* value increases, while the δ reduces dramatically. It is attributed to the enhanced elastic response at high frequency. The long-term aging influence on the viscoelastic behaviours of bitumen is more obvious at low frequency region. With the long-term aging duration lengthens, the G* and δ of aged bitumen increases and decreases, respectively. It indicates that the stiffness and elasticity of bitumen both significantly improve. During the long-term aging, the molecular interaction enhances and the free volume of whole bitumen system decreases, which results in the increase of cohesive energy and shear modulus. The black curves of virgin and aged binders are illustrated in Figure 7(c). The G* of bitumen decreases with the increase of δ value. It can be seen that the G* of bitumen at the same δ level decreases as the long-term aging degree prolonging. Meanwhile, when the G* parameter is fixed, the δ value of binder decreases. Overall, with the aging degree increases, the black curve of bitumen moves to the upper right corner and its scope becomes narrower.

Figure 7. Master curves and black curves of bitumen binders.

The crossover modulus and crossover frequency are two important parameters to assess the macroscale viscoelastic fluid-solid transitional behaviours as well as the microscale molecular weight distribution and polydispersity characteristics of bitumen (Jing et al. 2021). Definitionally, the crossover frequency is the loading frequency value when the storage and loss modulus are same, while the crossover modulus is the complex modulus at the crossover frequency. In this study, the crossover frequency values of virgin and aged binders are determined when the phase angle is 45°, while the corresponding crossover modulus values are identified from black diagrams.

Figure 8 illustrates the crossover modulus and crossover frequency of virgin and aged binders. The crossover modulus and crossover frequency values of bitumen both significantly decrease as the long-term aging degree deepening. In detail, the crossover modulus of bitumen decreases by 37.9, 63.9, 76.0 and 87.0% after short-term and long-term aging for 20, 40 and 80 h, respectively. At the same time, the corresponding crossover frequency decreases from 65.74 to 15.61, 1.26, 0.17 and 0.009 Hz. It implies that the aged binder shows the higher and wider molecular mass distribution, which results in the increase of stiffness, elastic components and relaxation time. Additionally, Figure 8(c) depicts that with the crossover frequency raises, the crossover modulus increases linearly, and the virgin bitumen has the highest crossover frequency and crossover modulus. With the increase of aging degree, the coordinate point (crossover frequency, crossover modulus) of aged binder would move down right gradually.

Figure 8. The crossover modulus and crossover frequency of virgin and aged binders.

3.5. Aging index from different parameters

The aging index is employed to quantitatively describe the aging degree of bitumen, which is strongly dependent on the selected parameter type. In this study, different kinds of parameters, including chemical components (aromatic, resin, asphaltene and colloidal index), functional group (combined index) and rheological properties (viscosity, flow activation energy, complex modulus, phase angle and crossover modulus), are used to calculate the aging index as Equation (10). (10) $\mathrm{A}\mathrm{g}\mathrm{i}\mathrm{n}\mathrm{g}\mathrm{i}\mathrm{n}\mathrm{d}\mathrm{e}\mathrm{x}\mathrm{A}\mathrm{I}={\displaystyle \frac{|P_{\mathrm{a}\mathrm{g}\mathrm{e}\mathrm{d}}-P_{\mathrm{v}\mathrm{i}\mathrm{r}\mathrm{g}\mathrm{i}\mathrm{n}}|}{P_{\mathrm{v}\mathrm{i}\mathrm{r}\mathrm{g}\mathrm{i}\mathrm{n}}}}$(10) where the P_aged and P_virgin are the selected parameter value of aged and virgin bitumen.

Figure 9 displays the chemical components-based aging index of aged binders. As expected, the aged bitumen with longer aging duration has the larger aging index. It is noted that different chemical components show the various aging index values, and all aging indices are lower than 1.0. The asphaltene-based aging index is higher than that of aromatic and resin-based regardless of the aging duration. Moreover, the resin-based aging index is the lowest. Interestingly, the colloidal index-based aging index is between the aromatic and asphaltene-based ones.

Figure 9. The SARA and colloidal index-based aging index of bitumen.

The combined index is the sum of carbonyl and sulfoxide index, which represents the oxidation degree of bitumen molecules. Figure 10(a) shows the combined index of different aged binders, and the long-term oxidation aging remarkably increases the combined index of bitumen. Moreover, the combined index-based aging index is presented in Figure 10(b). The magnitude of combined index-based aging index is larger than chemical components-based aging index. Thus, the functional groups are more sensitive to long-term aging than the SARA fractions.

Figure 10. The various aging indices of bitumen.

In addition, Figure 10(c,d) depict the viscosity and flow activation energy E_vis-based aging index of bitumen. It can be found that the high aging degree enlarges the aging index of bitumen distinctly, especially when the long-term aging time is 80 h. It is interesting to mention that the temperature shows notable influence on the magnitude of viscosity-based aging index. With the increase of testing temperature, the viscosity-based aging index of bitumen drops down dramatically. It is associated with the magnitude difference of viscosity value, which reduces significantly as the testing temperature rising. Compared to the viscosity-based aging index, the magnitude of E_vis-based aging index is much lower, which is even smaller than 1.0. Therefore, the viscosity is more suitable as an aging index parameter than the flow activation energy, but the temperature influence should be considered.

Lastly, the aging indices are calculated using rheological parameters of complex modulus, phase angle and crossover modulus, which are illustrated in Figure 10(e,f). Compared to the G*-based aging index, the δ-based AI is much lower, while the crossover modulus-based AI values are in between. Herein, the complex modulus is more susceptible to long-term aging time than the phase angle and crossover modulus. To compare the long-term aging sensitivity of different parameters, Table 5 summaries the aging indices of 80 h long-term aged bitumen calculated with various chemical and rheological parameters. Three levels of aging index are determined, including the Level I (AI > 1.0), Level II (0.5 < AI < 1.0) and Level III (AI < 0.5). Moreover, the sensitivity ranking is also listed in Table 5. In this study, the parameters in Level I are the rotational viscosity (80°C, 135°C and 160°C), complex modulus and combined index, while the asphaltene dosage, crossover modulus, aromatic dosage, flow activation energy and colloidal index all belong to Level II. The phase angle and resin dosage of bitumen are the least sensitive to long-term aging, while the 80°C viscosity is the most susceptible parameter. The low-temperature viscosity, complex modulus and combined index are recommended as the key indicators for calculating the long-term aging index of bitumen.

Table 5. The difference and ranking of different parameters for aging index.

	Parameters	AI of LAB80	Ranking
Level I (>1.0)	80°C viscosity	46.49	1
	Complex modulus (1 Hz)	10.45	2
	135°C viscosity	14.03	3
	160°C viscosity	5.50	4
	Combined index	2.68	5
Level II (0.5–1.0)	Asphaltene dosage	0.98	6
	Crossover modulus	0.87	7
	Aromatic dosage	0.64	8
	Flow activation energy	0.58	9
	Colloidal index	0.52	10
Level III (<0.5)	Phase angle	0.48	11
	Resin dosage	0.21	12

4. Aging reaction kinetics models

During the long-term aging process, different chemical reactions of bitumen molecules happen simultaneously, such as the molecular cracking, oxidation and polymerization, etc., which are accompanied with the physical phenomenon like the volatilization and agglomeration (Fallah et al. 2019). Although the aging mechanism of bitumen is complex, it is interesting and meaningful to fundamentally understand and establish the aging reaction kinetics models of bitumen, which is beneficial to predict the chemical components in aged bitumen with different aging degrees. In this section, the aging reaction kinetics models of bitumen are proposed for the first time in accordance with the variation of function groups and SARA fractions of bitumen during the long-term aging process.

4.1. Functional groups-based models

Figure 11 shows the variation of carbonyl and sulfoxide indices as a function of long-term aging time. With the aging time prolongs, the index difference increases linearly. Meanwhile, previous literatures reported that the increase rate of carbonyl and sulfoxide groups kept constant during the long-term aging process of bitumen (Liu et al. 2019, Jing et al. 2021). It indicates that the generation rate of functional groups is independent on their initial concentration. Therefore, the Zero-order model is suitable to describe the long-term aging reaction kinetics of bitumen from the viewpoint of the carbonyl and sulfoxide functional groups. The Zero-order models are listed as Equations (12)–(14). (11) $CI=CI+SI$(11) (12) $\mathrm{\Delta }CI=C{I}_t-C{I}_0=k_{CI}\cdot t$(12) (13) $\mathrm{\Delta }SI=S{I}_t-S{I}_0=k_{SI}\cdot t$(13) (14) $\mathrm{\Delta }CAI=CA{I}_t-CA{I}_0=k_{CAI}\cdot t$(14) where the CI₀, CI_t, SI₀, SI_t and CAI₀, CAI_t are the carbonyl, sulfoxide and combined index at reaction time t = 0 and t = t, respectively, while the k_CI, k_SI and k_CAI refer to the corresponding reaction rate constants.

Figure 11. Overall reaction models of bitumen from FTIR test results.

From Figure 11, it is depicted that the Zero-order model can well fit the correlation curve between the functional group index difference and aging reaction time with the correlation coefficient R² higher than 0.994, listed in Table 6. The functional group-based reaction rate constants are located in the region of 0.7–3.3*10⁻⁴ (mol·L⁻¹·h⁻¹). Besides, the k_CI and k_CAI values have the same magnitude, while the k_SI value is the lowest. Hence, the functional group-based aging reaction model of bitumen is significantly dependent on the type of selected functional groups.

Table 6. The reaction rate constant of functional groups in bitumen.

Aging index	CI	SI	CAI
k (mol·L^-1·h^-1)	2.598E-4	7.089E-5	3.307E-4
R^2	0.9947	0.9999	0.9965

4.2. SARA fractions-based models

The two-steps consecutive reaction model is adopted to quantitatively describe the aging reaction kinetics of bitumen for the first time. During the aging of bitumen, the aromatic fraction would convert to the resin fraction, which further changes into the asphaltene components. There are different complex chemical reaction models with several-steps reactions, including the chain reaction, opposing reaction, parallel reaction as well as the consecutive reaction (Guo and Lua 2001, Tjahjono et al. 2009, Leroy et al. 2013). By definition, during the consecutive reaction model, there are at least two-steps reactions, in which the reactant in the next reaction is the product of the previous step.

In this study, it is assumed that aging mechanism of the component’s transformation between the aromatics, resin and asphaltene fractions of bitumen belongs to the consecutive reaction. It should be mentioned that the air pressure is constant to the 2.1 MPa during the whole PAV test, indicating the variation of oxygen concentration is slight. To simplify the complex aging reaction, the oxygen reactant is omitted here. And the influence of oxygen concentration on the reaction rate is considered into the reaction rate constant k₁ and k₂. The detailed consecutive reaction equation is shown as follows: (15) $\left(\text{β}\right)A\stackrel{k_1}{\to }R\stackrel{k_2}{\to }\left(\gamma \right)A_S$(15) $t=0A_0{R}_0{A}_{S0}$ $t=t{A}_t{R}_t{A}_{St}$where A refers to the aromatic molecule, R represents the resin molecule and $A_S$ is the asphaltene molecule. The parameters of β and γ are the stoichiometric number. It means that during the long-term aging reaction, when one resin molecule is generated, there are β aromatic molecules consumed. Similarly, when one resin molecule is depleted, the number of new asphaltene molecules is γ. Moreover, the k₁ and k₂ is the reaction rate constant of the first and second step reactions, respectively.

Before the aging process (at t = 0), the concentration value of aromatic, resin and asphaltene fraction is A₀, R₀ and A_S0, respectively. When the aging reaction time is at t = t, the corresponding concentration is A_t, R_t and A_st. The reaction order can be different from the sum of stoichiometric number of reactants. To explore the aging reaction kinetics, the reaction order of first and second reactions is assumed as n and m, respectively. Thus, the consumption rate of aromatic fraction and the generation rate of asphaltene components is calculated as Equations (16) and (17). (16) ${\displaystyle \frac{\mathrm{d}A}{\mathrm{d}t}=-k_1\cdot A^n}$(16) (17) ${\displaystyle \frac{\mathrm{d}{A}_S}{\mathrm{d}t}=k_2\cdot R^m}$(17) where A and A_s are the concentration of aromatic and asphaltene fractions, respectively; t refers to the reaction time; k₁ and k₂ are the corresponding reaction rate constants.

The integration formula of Equation (16) can be rewritten as follows: (18) $A_t^{\left(1-n\right)}-A_0^{\left(1-n\right)}=\left(n-1\right)\cdot k_1\cdot t,\left(n\ne 1\right)$(18) For the Zero-order reaction model, the reaction rate is independent on the initial concentration of reactants. The Equation (16) can be rewritten as follows: (19) ${\displaystyle \frac{\mathrm{d}A}{\mathrm{d}t}=-k_1\cdot A^0=-k_1}$(19)

By further simultaneous integration of both sides, Equation (19) can be expressed as shown: (20) $A_t=A_0-k_1\cdot t$(20) In addition, when the reaction model is First-order, the integral form of Equation (16) is as follows: (21) $A_t=A_0\cdot e^{\left(-k_1\cdot t\right)}$(21)

Similarly, regarding the Second-order reaction model, the relationship between the reactant concentration and time can be described as: (22) ${\displaystyle \frac{1}{A_t}={\displaystyle \frac{1}{A_0}+k_{1}t}}$(22) To investigate the reaction kinetics rate k₁ from the aromatic to resin molecules, the molar concentration of aromatic fraction as a function of long-term aging time is displayed in Figure 15. With the aging time prolongs, the aromatic molarity in bitumen decreases dramatically. Moreover, the decreasing rate seems to slow down gradually, which may be related to the reduced concentration of reactant and increased product dosage during the oxidative reaction process. The common reaction kinetics models, Zero-order, First-order and Second order, are firstly used to describe the relationship between the aromatic molarity and aging time. The correlation curves with different kinetics models are also shown in Figure 12. The reaction rate constant k₁ and correlation coefficient R² values can be obtained and listed in Table 7. It can be seen that the order of R² value is Second-order > First-order > Zero-order, indicating that the Second-order reaction model is more accurate to describe the oxidative aging kinetics of aromatic fractions. However, these reaction models are assumed to fit the reaction curve of aromatic molecules, and the optimum reaction model is still unclear, although the Second-order model seems to be more reasonable.

Figure 12. The molar concentration of aromatic fractions.

Table 7. The parameters of aromatic oxidative reaction in different kinetics models.

Kinetics models	Zero-order	First-order	Second-order	Third-order
k1	0.00313	0.00592	0.0109	0.02
R^2	0.928	0.967	0.992	0.999

The general n-order reaction model (shown in Equation (18)) is then adopted to fit the correlation curve between the aromatic molarity and reaction time, which is also illustrated in Figure 12. The correlation coefficient R² value is 0.999, and the reaction order of aromatic molecules is determined as 3.045. Hence, it can be deduced that the optimum kinetics model for aromatic fraction is the Third-order reaction model. Table 7 displays the reaction rate constant k₁ of aromatic with the different kinetics models during the long-term aging process. With the reaction order increasing from zero to three, the k₁ and R² value both increase. When the Third-order reaction model is applied, the calculated k₁ value is 0.02 (mol·L⁻¹)⁻²(h)⁻¹ and the R² value is 0.999.

In order to estimate the second step reaction kinetics from the resin to asphaltene molecules, the asphaltene molarity as a function of long-term aging time is drawn in Figure 13. The linear correlation between asphaltene molarity and reaction time is significant, and the Zero-order kinetics model is applied here to describe the generation rate k₂ of asphaltene molecules. In Equation (23), t is the aging time, A_st and A_s0 are the asphaltene concentration at t = t and t = 0, and k₂ represents the reaction rate constant regarding the conversion from resin to asphaltene fractions. (23) $A_{st}=A_{s0}+k_2\cdot t$(23) It means that the generation rate of asphaltene molecules during the long-term aging process is constant and independent on the reactant concentration (resin molecules), which is similar to that of aforementioned functional groups from FTIR test. The reaction rate constant k₂ and correlation coefficient R² for the second step reaction are 3.85E-4 mol·(L·h)⁻¹ and 0.9987, respectively. It indicates that the Zero-order model can well fit the aging reaction kinetics of asphaltene molecules.

Figure 13. Reaction models from resin to asphaltene fractions of bitumen.

It is worth mentioning that the second step reaction rate constant k₂ is much lower than the k₁ in the first step in the consecutive reaction model. Hence, the reaction rate from resin to asphaltene is much slower than that from aromatic to resin, which leads to the large consumption of aromatic fraction and the increase of resin concentration. Moreover, the rate control step of the whole consecutive reaction is the second-step, and the reaction rate constant k of the whole consecutive reaction is close to k₂. From the viewpoint of reaction kinetics knowledge, the increase of asphaltene fractions would be beneficial to hinder the aging reaction rate of bitumen. Interestingly, the reaction rate constant k₂ is very close to the value of k_CI and k_CAI obtained from FTIR test. It further validates that the transformation from the resin to asphaltene molecules is the control step of the whole consecutive reaction model. Further, the k_SI value calculated from the change rate of sulfoxide functional groups is the lowest. Therefore, the parameters of k_CI and k_CAI are more appropriate than the k_SI to describe the long-term aging reaction kinetics of bitumen.

4.3. Quantitative conversion relationship between aromatic, resin and asphaltene

The final step of the consecutive reaction model is to determine the stoichiometric number β and γ, which is of significance to understanding the material transformation relationship between the aromatic, resin and asphaltene components during the artificial long-term aging process from the molecular level. From the Equation (15), when the aging time is t, the consumption of aromatic is (A₀−A_t), while the yield of asphaltene is (A_st−A_s0). According to the consecutive reaction model, the concentration difference of resin fraction at time t comes from two parts, including the generation in the first step reaction and the consumption in the second step reaction. Thus, the concentration difference of resin fraction can be introduced as follows: (24) $R_t-R_0={\displaystyle \frac{A_0-A_t}{\text{β}}-{\displaystyle \frac{A_{st}-A_{s0}}{\gamma }}}$(24) where the R₀, A₀, and A_s0 are the initial molarity of resin, aromatic and asphaltene fractions at t = 0, respectively; while the R_t, A_t, and A_st are their corresponding molarity value at reaction time t; β and γ are the stoichiometric number.

We assume that the first step is the Three-order reaction, thus the reaction rate equation of aromatic fraction is as follows: (25) ${\displaystyle \frac{1}{A_t^2}-{\displaystyle \frac{1}{A_0^2}=2\cdot k_1\cdot t}}$(25) Similarly, the second step reaction is assumed as the Zero-order model, and the reaction rate equation of asphaltene fraction is: (26) $A_{st}-A_{s0}={\mathrm{k}}_2\cdot t$(26)

By substituting the Equations (25) and (26), the Equation (24) can be rewritten as: (27) $\mathrm{\Delta }R=R_t-R_0={\displaystyle \frac{A_0-{\displaystyle \frac{A_0}{\sqrt{1+2\cdot k_1\cdot t\cdot A_0^2}}}}{\text{β}}-{\displaystyle \frac{k_2\cdot t}{\gamma }}}$(27) where ΔR refers to the molarity difference of resin fraction, k₁ and k₂ are the reaction rate constant of the first and second-step reactions, respectively.

The Equation (27) introduces the relationship between the molarity difference of resin fraction ΔR and the long-term aging time t, which is used to fit the curve in Figure 14. It should be noted that in this study, the value of A₀, k₁ and k₂ is 0.694 mol·L⁻¹, 0.02(mol·L⁻¹)⁻²(h)⁻¹ and 3.85E-4 mol·(L·h)⁻¹. The predicted value of β and γ is approximately 2.82 and 0.58, respectively, and the R² value is 0.992. To this end, the detailed consecutive reaction of bitumen during the long-term aging process is shown as follows: (28) $\left(2.82\right)A\stackrel{k_1}{\to }R\stackrel{k_2}{\to }\left(0.58\right)A_S$(28) According to the consecutive reaction Equation (28), it can be interpreted that when one resin molecule is generated, the consumption number of aromatic molecules is 2.82. At the same time, when one resin molecule is consumed, only 0.58 asphaltene molecule can be obtained. The detailed illustration of the consecutive reaction model is presented in Figure 14.

Figure 14. Reaction model of resin fraction.

It should be mentioned that the all aromatic, resin and asphaltene fraction are the components groups containing amounts of various molecule types. Meanwhile, the SARA fractions of bitumen here is measured through the thin-layer chromatography, and the aging reaction is composed to the component conversion and molecular agglomeration. For instance, due to the molecular interaction and existence of polar functional groups in aromatic, resin and asphaltene molecules, the molecular agglomeration would occur, which would increase the difficulty of distinguishing the high-weight molecules with the aggregates of light-weight molecules. Thence, the chemical reaction and physical aggregation are both considered in the consecutive reaction model during the long-term aging process of bitumen.

Figure 15 illustrates the long-term aging reaction kinetics models of bitumen based on variation rate of functional groups (C = O or S = O) and chemical compositions (aromatic, resin and asphaltene). From the viewpoint of functional groups, the average reaction rates of carbonyl, sulfoxide and combined indices are calculated, which is the common-used model of bitumen oxidative aging. Although the functional groups-based model can reflect the oxidative rate of bitumen molecules with oxygen to a certain degree, it fails to elaborate the conversion mechanism between chemical compositions of bitumen during the long-term aging process. Hence, one novel two-step consecutive model is established according to the conversion correlation between aromatic, resin and asphaltene fractions. The reaction rate constants in SARA-based kinetics models are calculated, and the corresponding stoichiometric numbers of aromatic, resin and asphaltene are firstly determined. It should be mentioned that the increase of resin and asphaltene dosage are due to the molecular transformation, but also related to the molecular agglomeration. During the aging process, the bitumen molecules would aggregate together due to the increased polarity and intermolecular attraction, which is difficultly separated by solvents in SARA fractions measurement and leads to the increased number of large-scale molecules (resins and asphaltenes) in SARA results. Hence, more works should be done to eliminate the influence of molecular agglomeration on the SARA fractions distribution with the help of functional solvents and effective separation methods.

Figure 15. The illustration graph of long-term aging reaction kinetics models.

5. Conclusions and recommendations

In this study, the functional groups-based model and a novel SARA fractions-based long-term aging two-step consecutive reaction kinetics models were established. Meanwhile, the quantitative conversion relationships between the aromatic, resin and asphaltene molecules were derived for the first time. The main conclusions are listed as follows:

1. During the long-term aging of bitumen, the heavy-weight components (resin and asphaltene) continued to increase. Moreover, both C = O and S = O functional groups in bitumen molecules enlarged, while the aromaticity and saturation of bitumen molecules increased and decreased, respectively.

2. The asphaltene-based aging index was higher than that of aromatic and resin-based regardless of the aging duration. The functional groups were more sensitive to long-term aging than the SARA fractions. Additionally, the viscosity was more suitable as an aging index parameter than the flow activation energy, but the temperature influence should be considered.

3. The Zero-order model was suitable to describe the long-term aging reaction kinetics of bitumen from the viewpoint of the carbonyl and sulfoxide functional groups, and the functional group-based reaction rate constants were located in the region of 0.7–3.3*10⁻⁴ (mol·L⁻¹·h⁻¹). Besides, the k_CI and k_CAI values had the same magnitude, while the k_SI value was the lowest.

4. The aromatic fraction would convert to the resin fraction and further change into the asphaltene components. In the consecutive reaction model, there were at least two-step reactions and the reactant in the next reaction was the product of the previous step. The Second-order reaction model was more accurate to describe the oxidative aging kinetics of aromatic fractions than the Zero-order and First-order models. However, the most optimum kinetics model for aromatic fraction was the Third-order reaction model and the corresponding k₁ value was 0.02 (mol·L^-1)⁻²(h)⁻¹.

5. The Zero-order model could well fit the aging reaction kinetics of asphaltene molecules and the related reaction rate constant k₂ was 3.85E-4 mol·(L·h)⁻¹. The reaction rate from resin to asphaltene was much slower than that from aromatic to resin. Moreover, the transformation from the resin to asphaltene molecules was the control step of the whole consecutive reaction model.

6. When one resin molecule was generated, the consumption number of aromatic molecules was about 2.82. At the same time, when one resin molecule was consumed, only 0.58 asphaltene molecule could be generated.

In this study, the functional groups and SARA fractions-based long-term aging reaction kinetics models of bitumen were established for the first time. However, more types and long-term aging conditions of bitumen should be conducted to validate the reasonability of the proposed kinetics models in future works. Meanwhile, the dynamic correlation model regarding the chemo-rheological properties of bitumen binder during long-term aging process will be further studied. More importantly, the developed reaction kinetics models should be further validated by the aged binders from the thin-film oven aging and field aging procedures.

Disclosure statement

No potential conflict of interest was reported by the author(s).

Additional information

Funding

The first author would like thanks for the financial support from China Scholarship Council [grant number 201906450025].