Multi-scale characterization of lignin modified bitumen using experimental and molecular dynamics simulation methods

(cid:1) The molecular model of lignin modiﬁed bitumen was established. (cid:1) The inﬂuence of lignin on the thermodynamics, rheological and adhesion properties of bitumen were investigated. (cid:1) The correlation between microscale and macroscale performance were built. (cid:1) The dynamic diffusion and self-healing properties of lignin modiﬁed bitumen were predicted. Lignin, as a major waste from biofuel and paper industries, could be utilized as a modiﬁer to enhance the relevant performance of bitumen. However, the effects of lignin on the thermodynamics properties and molecular structure of bitumen are rarely studied. Meanwhile, the potential modiﬁcation mechanism of lignin modiﬁed bitumen is still unclear. Molecular dynamics (MD) simulation and laboratory experimental methods are combined to explore the inﬂuence of lignin on the thermodynamics characteristics, rhe- ological properties as well as the molecular structure of bitumen. The lignin modiﬁed bitumen with different dosages of lignin (10, 20 and 30 wt%) were prepared. DSR results from a macroscale view reveal that lignin could signiﬁcantly improve the modulus, elastic recovery and rutting resistance of bitumen, but it adversely affects the fatigue performance. Meanwhile, the MD simulation results from a microscale perspective show that lignin could increase the density, cohesive energy density, shear viscosity, modulus and adhesive strength of bitumen. However, the free volume, diffusion coefﬁcient and self-healing ability of lignin modiﬁed bitumen are weakened with the increase of lignin dosage. The MD simulations results are consistent with the experimental data. Furthermore, the correlations between the microscale and macroscale properties of lignin modiﬁed bitumen indicate that the physical and rheological properties of bitumen both depend on the molecular structure dramatically. The ﬁndings of this research can provide insights for an in-depth understanding of the effect of lignin on bitumen.

Multi-scale characterization of lignin modified bitumen using experimental and molecular dynamics simulation methods Shisong  The molecular model of lignin modified bitumen was established. The influence of lignin on the thermodynamics, rheological and adhesion properties of bitumen were investigated. The correlation between microscale and macroscale performance were built. The dynamic diffusion and selfhealing properties of lignin modified bitumen were predicted. g r a p h i c a l a b s t r a c t a r t i c l e i n f o

Introduction
Bitumen is the most widely used construction material for flexible pavements. However, due to the limited amount of crude oil, bitumen resources are exhaustible [1,2]. Recently, many researchers devoted their efforts to minimize the application of bitumen by implementing techniques to re-use the reclaimed asphalt binder or exploring other alternative materials. Alternative binders, such as bio-based bitumen, could reduce the requirement of petroleumbased bitumen, and asphalt roads with superior performance characteristics can be constructed [3][4][5].
The application of biomass materials in road construction has many advantages in terms of promoting the sustainability in paving industry [6,7]. It has been proved that, lignin as one of the main by-products during the biomass working process can improve the engineering properties of bitumen and can be used as substituent of petroleum-based bituminous binder. Xu et al. found that the addition of lignin was beneficial to enhance the rutting and aging resistance of bitumen [8]. Moreover, others studies pointed out that lignin (with 2, 4, 6 and 8% of the total binder by weight) was an effective modifier that could be used partly to replace the bitumen [9,10]. However, current studies in the field of lignin modified bitumen (LMB) mainly focus on the evaluation of rheological properties at the macroscale. The influence of lignin on the thermodynamics characteristics of bitumen remains unclear, which is difficult to be explored through using macro-scale experimental methods. Moreover, the dosage of lignin into bitumen was limited in previous studies, and the alternative feature of lignin to unrenewable bituminous materials partly has not been estimated yet.
Molecular dynamics (MD) simulation method has been successfully performed to estimate the effects of block copolymers on bitumen as well as the interfacial interaction between bitumen and mineral aggregates [11][12][13]. The most commonly adopted molecular model of bitumen is a twelve-components model, which was originally developed to provide predictions of the chemical and physical characteristics, such as density, solubility parameter and surface energy [14,15]. Through the MD simulations, researchers concentrated on many important issues, such as the bitumen modification, oxidative aging, self-healing and adhesion from a microscale perspective [16][17][18]. Regarding the bitumen modification, Ding et al. reported that the agglomeration of styrene-butadie nestyrene (SBS) copolymers was mainly affected by the molecular structure of asphaltenes in bitumen [19]. Meanwhile, Khabaz et al. also utilized the MD simulations to evaluate the physical characteristics of styrene-butadiene rubber (SBR) modified bitumen, indicating that SBR copolymer intensified the aggregation of asphaltenes and led to reduction of the mobility for all bituminous molecules [20]. At the same time, Guo et al. successfully predicted the physical properties and microstructure of rubber bitumen by using MD simulations [21]. Therefore, MD simulation method has great functions on exploring the thermodynamics of bitumen, as well as fundamentally explaining the potential mechanism of modification, oxidation and rejuvenation of bitumen from atomistic perspective [22][23][24].
The existing studies with regard to the lignin effects on bitumen mostly focused on the material rheological and compositional changes at macro-scale by conducting laboratory experimental programs [4][5][6]. However, the modification mechanism of lignin in bitumen has barely been investigated at micro-scale. In this study, a molecular characterization of LMB is performed through MD simulations [25,26]. The molecular models of LMB are firstly established, and the microscale properties are predicted by MD simulations. At the same time, the influence of lignin on the rheological properties of bitumen at the macroscale perspective is evaluated through laboratory tests. Furthermore, some correlations between thermodynamics and microstructure parameters obtained from MD simulations with the rheological properties at macro-scale are explored to further identify and explain the lignin modification in bitumen [27][28][29][30]. Overall, the main objective of this study is to evaluate: (1) the effects of lignin on the thermodynamic characteristics, rheological properties and molecular structure of bitumen; (2) the lignin modification mechanism in bitumen; (3) the effects of lignin on the self-healing and adhesion of bitumen; and (4) the correlations between microscale and macroscale properties of LMB obtained from the MD simulations and laboratory experiments, respectively.

Materials and preparation
In this study, a 70/100 bituminous binder was used. Table 1 shows the basic properties of bitumen and lignin. The density of lignin is about 1.38 g/cm 3 , which is larger than bitumen. Moreover, the apparent morphology and scanning electron microscopy (SEM) photographs of lignin are displayed in Fig. 1.
Before the preparation of LMB binders, lignin was pre-heated for 2 h in a draft oven at 100℃ to eliminate the moisture influence. After that, lignin (with the dosage of 10, 20 and 30 wt%) was added into the binder, which was also pre-heated at 140℃. The new lignin-bitumen systems were sheared for 2 h by using a highspeed shear mixing device with 3000 rpm shear speed at 160℃ to ensure the sufficient dispersion of lignin in bitumen.

Testing methods
Frequency and temperature sweep tests were conducted by using DSR device to evaluate the impact of lignin on the rheological behaviours of bitumen [29]. During the frequency sweep test, loading frequency varies in the region of 0.1 to 100 rad/s at 25 ℃ [4]. Moreover, both fresh bitumen and LMB binders are subjected to temperature sweep test, in which the constant frequency is 10 rad/s and temperature changes from À5 to 60 ℃. Besides, multiple stress creep and recovery (MSCR) test is performed at 52 ℃ with the stress of 0.1 kPa, according to ASTM D7405. Meanwhile, all bitumen samples are subjected to the Linear Amplitude Sweep (LAS) test at 25 ℃ to assess the influence of lignin on the fatigue life of bitumen, according to AASHTO TP 101. In addition, the low-temperature properties of LMB binders are evaluated with G-R parameter and relaxation test. The G-R parameter is calculated as ð G Ã COS 2 d sind Þ on basis of the complex modulus G* and phase angle d values obtained from frequency sweep test with a frequency of 0.005 rad/s at 15 ℃. Moreover, the relaxation tests of bitumen samples are conducted under the strain-controlled mode at 0℃. The strain of tested sample firstly rises from 0 to 1% in one second, followed by keeping the constant shear strain of 1% with the relaxation time of 100 s. The stress value during relaxation test is recorded to assess the relaxation behavior of bitumen.

Basic theory of MD simulation
The MD simulation is an effective tool to explore the physical properties and the microstructure of bitumen from a nanoscale perspective by solving Newton's law as stated in Eq.1 and 2.
where E represents total energy, U shows the total potential energy, p is the momentum, v refers to the velocity of atomic and m is the mass of one molecule [18]. The software of Material Studio was applied to carry out the simulation in this study. During the MD simulations, the COMPASS Ⅱ force field is utilized to promote the molecular interaction and calculate the potential energies of lignin modified bitumen models, including the non-bonded energy and valence energy as stated in Eq. 3.
where U val represents the covalently bonded interaction composed of bond stretching, angle bending, dihedral torsion and improper terms, and U nb refers to the non-bonded interaction including van der Waals and Coulomb electrostatic interactions [18].

Molecular components of bitumen and lignin
Bituminous binder is composed of complicated components with different molecular weight and structure, such as the alkanes, cycloalkanes and aromatic hydrocarbons [31]. In accordance of the polarity difference, bitumen could be separated into four main fractions dissolved in different solvents, which are called as the saturate, aromatic, resin and asphaltene [12,32,33]. Li and Green-field proposed the twelve-components molecular model for different kinds of virgin bitumen and found that some predicted properties of binder from molecular simulation, such as density and solubility parameter, were consistent with experimental results [34]. In this study, the twelve-components molecular model is used to establish the molecular model of bitumen and these twelve-components molecules are presented in Fig. 2.
The infrared (IR) spectra of lignin is illustrated in Fig. 3(a), in which some peaks could be remarkably detected at the wavenumber of 2938 cm À1 , 1696 cm À1 , 1595 cm À1 , 1217 cm À1 as well as 1029 cm À1 , which represents the existence of several typical functional groups of -OH, C = O, benzene ring, C-O and C-H bonds, respectively. Based on the IR spectra and previous literature [35], the representative molecular structure of lignin is drawn in Fig. 3 (b). The molecule formula of lignin is C 47 H 49 O 18 , and its molecular weight is 901 g/mol. It could be seen that the molecular structure of lignin contains the benzene ring, phenolic, hydroxide, aldehyde and methoxy groups, which have been proved that are resistant against oxidation [36].

Molecular model of bitumen and LMB binders
The molecular model of bitumen is built in accordance with the detailed components shown in Table 2. The initial density is set as 0.1 g/cm 3 to prevent the overlapping of atoms. After establishing the initial molecular model of bitumen, geometry optimization is performed to minimize the total energy. Following that, molecular dynamics simulation at 298.15 K is conducted to promote the system reaching an equilibrium state, which is composed of the canonical ensemble (NVT) for 200 ps as well as the isothermalisobaric (NPT) for 50 ps with the time step of 1 fs. The Nose thermostat and Andersen barostat are employed to control the temperature and pressure during MD simulations, respectively. Importantly, the electrostatic and Van der Waals interaction are set as Ewald and Atom-based with a cut-off distance of 12.5 Å. The stable molecular model of virgin bitumen after MD simulation is illustrated in Fig. 4.
To investigate the effects of lignin on the thermodynamics characteristics and molecular structure of bitumen, the molecular models of LMB with various lignin contents are also established. Based on the lignin dosage of 10, 20 and 30 wt%, the amount of lignin molecule in the molecular models of LMB is four, nine and fifteen, respectively. After being built, the initial molecular models of LMB binders are subjected to the NVT and NPT dynamics equilibrium simulations.
To check the equilibrium state of the aforementioned molecular models, the energy variation as a function of the simulation time during the NPT simulation is recorded and presented in Fig. 5, including the potential energy, kinetic energy, non-bond energy and total energy. All types of energy show stable during the NPT simulations at 298 K when the simulation time exceeds to 15 ps. Take the 10 wt% lignin modified bitumen as an example, the difference between its molecular structure before and after NPT equilibrium process can be seen in Fig. 6. The dynamics equilibrium operation could make the whole structure denser. The equilibrium molecular models of LMB with lignin dosage of 10, 20 and 30 wt% are obtained, as illustrated in Fig. 7, of which the thermodynamics properties can be further predicted.

Effect of lignin on viscoelasticity
The master curves of complex modulus and phase angle of LMB binders are displayed in Fig. 8. The complex modulus values of bitumen increase remarkably with the frequency increasing, showing a linear correlation in double logarithmic coordinates. Meanwhile, the addition of lignin enhances the complex modulus and improves the deformation resistance of bitumen dramatically. It should be noted that adding lignin could decline the phase angle of binder significantly. According to these results, it seems that lignin operates as modifier in bitumen improving the overall material elasticity, mainly at high temperatures. Additionally, with the lignin dosage increases, the viscous characteristic of LMB tends to more obvious. Hence, lignin is a sustainable modifier to strengthen the high-temperature performance and reduce bitumen consumption.

Effect of lignin on rutting and fatigue performance
The effects of lignin dosage on the rutting and fatigue resistance of fresh or aged bitumen are evaluated, which are presented in Fig. 9. With the increase of temperature, both rutting and fatigue parameters of all bitumen samples decline significantly, which is associated with the weakening of molecular interaction and internal friction, as well as the improvement of molecular mobility. Thus, the high temperature would weaken the rutting resistance of bitumen, whereas it would be beneficial to form systems with improved fatigue cracking resistance. The correlation equations between testing temperature and rutting as well as fatigue parameters are listed in Fig. 9. The temperature sensitivity of bitumen could be characterized based on the slope value. The performance of binder with a larger slope value is more sensitive to temperature variation. It is interesting that the addition of lignin would slightly increase the temperature sensitivity of bitumen at first. However, when lignin content is 30 wt%, the temperature sensitivity of lignin modified bitumen is better than that of virgin bitumen, which means that high lignin dosage would enhance the thermal stability of bitumen dramatically.
In addition, both rutting and fatigue parameters of bituminous binders increase as the lignin dosage increases, indicating that lignin can significantly improve the rutting resistance of bitumen, while it would lead to binders prone to fatigue damage. To quantitatively evaluate the lignin influence on the rutting and fatigue performance of bitumen, the rutting and fatigue failure temperatures of virgin and LMB binders are calculated when the rutting and fatigue parameter is equal to 1.0 and 5000 kPa, respectively. Fig. 10(a) indicates that adding lignin could increase both rutting and fatigue failure temperatures of binder. Compared with virgin binder, when the lignin dosage is 10, 20 and 30 wt%, the rutting failure temperature of bitumen enhances by 1.7, 5.0 and 12.5%, respectively. At the same time, the fatigue failure temperature rises by 1.9, 3.2 and 5.6°C. Therefore, lignin dosage should be controlled to balance the rutting and fatigue performance of LMB.
Moreover, the LAS test is also conducted to calculate the fatigue life (N f ) of bitumen binders according to the viscoelastic continuum damage. The fatigue life at two stress levels of 2.5% and 5.0% is displayed in Fig. 10b. The fatigue life of bitumen is weakened with the applied stress increasing from 2.5 to 5.0%. As expected, adding lignin reduces the fatigue life of bitumen regard- less of loading stress. As a result, the lignin content would adversely deteriorate the fatigue resistance and service life of bitumen, which is consistent with the aforementioned fatigue parameter results.

Effect of lignin on creep and elastic recovery behaviours
The multiple stress creep and recovery (MSCR) test is performed to further investigate the effects of lignin on the elastic behaviour

Effect of lignin on low temperature performance
The effects of lignin on the low temperature properties of bitumen are evaluated with the G-R parameter and relaxation behaviour. Fig. 12(a) shows the G-R parameter value variation of different lignin modified bitumen before and after long-term aging, which can be used to predict the crack potential of lignin modified bitumen. The G-R value for representing the onset of cracking and surface cracking is 180 and 450 kPa. As expected, the G-R value increases as the aging degree deepens dramatically, indicating that the crack of bitumen binder occurs easily when it suffers from aging. Lignin has significant influence on the G-R parameter of bitumen, which enhances with the increase of lignin dosage dramatically. That means that the addition of lignin adversely increases the cracking propensity of bitumen at low temperatures. When the lignin dosage increases to 30 wt%, the G-R parameter of lignin modified binder after 20 h long-term aging is higher than the surface cracking standard of 450 kPa. And the G-R values of all LMB binders after 40 h long-term aging are larger than surface cracking point. Hence, the lignin content should be controlled to guarantee the low temperature cracking resistance.
In addition, the relaxation test is employed to further investigate the low-temperature relaxation ability of LMB. The detailed relaxation curves of virgin and LMB are displayed in Fig. 12(b). With the relaxation time prolongs, the residue stress of bitumen reduces dramatically. It is clear that the lignin has great influence on the relaxation behaviour of binder. Fig. 12(c) illustrates the maximum stress, residue stress at 100 s and the residue stress ratio of bitumen binders. With the increase of lignin dosage, the maximum stress and residue stress at 100 s both significantly enhance.
Meanwhile, the residue stress ratio increases remarkably of LMB as the lignin content increases, indicating that the addition of lignin could reduce the relaxation ability of bitumen and increase the cracking propensity. Further, Fig. 12(d) displays the relaxation time of bitumen binders when the stress reduction is 25% and 50%, respectively. Similarly, the relaxation time of LMB binder extends as the lignin dosage increasing. According to the results from G-R parameter and relaxation test, it is clear that excessive lignin damages the low temperature properties of bitumen. Hence, the lignin dosage should be limited to fulfil the requirements of low temperature performance for LMB binders.

Effect of lignin on density
The density values of virgin and LMB binders can be acquired from the MD simulation (see Fig. 13(a)). It can be found that the density increases dramatically, which is related to the decline of volume during the NPT simulations. After about 20 ps, the density value keeps constant, indicating that the molecular structures of binders reach to equilibrium state. The density values of all stable molecular models are shown in Fig. 13(b), and the density of virgin bitumen is approximately 0.98 g/cm 3 , which is consistent with the previous studies [37 -39]. As illustrated in Table 1, the density of lignin is 1.38 g/cm 3 , which is larger than that of virgin bitumen. The density of LMB rises with the increase of lignin dosage. When the lignin content varies from 10, 20 to 30 wt%, the density rises slightly from 0.98 to 1.00, 1.04 and 1.05 g/cm 3 , respectively. The additional lignin may occupy the free volume and absorb the light components in the virgin bitumen, which hinders the molecular mobility and enhances the density of binder. The free volume results will be presented and discussed in a later section.

Effect of lignin on cohesive energy density and solubility parameter
The cohesion energy density (CED) is an effective indication to characterize the attractive interaction between molecules and assess the cohesive ability of substances. Bitumen with higher CED value would possess superior cohesive strength. The CED values of bitumen binders at 298 K are illustrated in Fig. 14a. The Van der Waals attraction between bitumen molecules is approximately 10-100 times higher than the electrostatic interaction. The addition of lignin remarkably strengthens the CED value and enhances the cohesive strength of bitumen. As shown in Fig. 3, there are many oxygen atoms in the lignin molecule, which in turn is associated with the high-polarity of lignin. Thus, the lignin molecules would have stronger interaction with the bituminous components. Compared to virgin bitumen, the CED value of LMB binder with the lignin dosage of 10, 20 and 30 wt% is enhanced by 36.2, 45.4 and 47.4%, respectively. The effect of lignin on improving the cohesive performance of bitumen is limited when the lignin dosage exceeds to 20 wt%.
Additionally, the solubility parameters l of binders are also obtained from simulation results, which can be calculated from CED values, according to Eq.4.
where l is the total solubility parameter, while l vdw and l elc refers to the solubility parameter contributed from van der Waals and electrostatic interaction, respectively. The solubility parameter of LMB was presented in Fig. 14b at 298 K. It is worth noting that the van der Waals interaction has a greater effect on the solubility parameter of bitumen than electrostatic interaction. The total solubility parameter of virgin bitumen is about 15.4 (J/cm 3 ) 0.5 , which is close to the experimental result. The addition of lignin increases the solubility parameter of bitumen, and the maximal increase is 21.4%, when the lignin dosage is 20 wt%. As the lignin dosage increases, the solubility parameter difference between virgin bitumen and LMB binder becomes larger, indicating that the LMB with high lignin content would occur phase separation easily. Interestingly, when the lignin content exceeds to 10 wt%, the influence of lignin on electrostatic interaction tends to be more obvious. For instance, the van der Waals solubility parameter starts to decline when the lignin dosage increases from 20 to 30 wt%, while the electrostatic solubility parameter increases continuously, which may be attributed with the high-polarity characteristic of lignin.

Effect of lignin on free volume
The free volume and occupied volume of LMB can be determined based on the molecular structure obtained from MD simu-  lations, which is shown in Fig. 15. The free volume is an important parameter to connect the physical and rheological properties of material. The lower the free volume of bitumen is, the smaller the molecular mobility is, which would result in a phase transition from viscoelastic state to glass state. Obviously, bitumen with glass state would have lower free volume, and the molecules move difficultly, which is beneficial to enhance the modulus but has an adverse effect on the thermal cracking resistance of bitumen. The occupied and free volume values of LMB are displayed in Fig. 16a. In regards to virgin bitumen, the occupy volume is two times larger than free volume. The former is contributed to the stiffness, while the latter is related to the flexibility of bitumen. The addition of lignin remarkably increases the occupied volume and reduces the free volume of bitumen, which is associated with the filling contribution of lignin. With the lignin dosage increases, free volume of LMB binder continues to decrease and the compressive possibility reduces, which leads to the increase of density. Furthermore, the relationship between free volume and lignin concentration is shown in Fig. 16b. There is a linear relation between lignin dosage and free volume fraction, which is helpful to predict the free volume of modified binders with various lignin proportions.

Effect of lignin on inter-diffusion coefficient
The rheological properties of bitumen are highly related to the self-diffusion ability of molecules [31,33]. Previous studies also reported that there was significant interaction between SARA fraction in form of colloidal structure of bituminous material [40][41][42]. To evaluate the influence of lignin on the diffusion coefficient of bitumen molecules, the mean space displacement (MSD) analysis of SARA fractions is performed, including saturate, aromatic, resin as well as asphaltene molecules. During the MD simulations, the MSD values of molecules are obtained by measuring the mean square distance between the initial and final position at the simulation time t, which represents the movement displacement of bitumen molecules. The MSD of bitumen molecules overtime in an equilibrium model can be calculated following Eq.5.

MSD t ð Þ ¼< ½r
where MSD is the mean space distance of molecules, r i (t) and r i (0) represent the position vector of particle i at time 0 (initial state) and t (equilibrium state), respectively. Moreover, the angular brackets are the average of such particles in the system over the whole simulation process [18]. These equilibrium molecular models of virgin and LMB binders are subjected to additional MD simulation using microcanonical ensemble (NVE) for 180 ps at 298.15 K. The MSD values of SARA fractions for binders are recorded and shown in Fig. 17. It is clear that the MSD values of SARA molecules increase with the simulation time prolongs. Meanwhile, saturate fraction with the lowest molecular weight has the highest MSD value, while the MSD value of asphaltene is the smallest, which is related to its largest molecular weight. In other words, the high molecular weight of molecules would weaken the mobility. Regarding the effect of lignin, the MSD values of all SARA fractions decline with the increase of lignin dosage, which indicates that lignin would deteriorate the molecular mobility of bitumen. On the one hand, the addition of lignin reduces the free volume of bitumen dramatically. On the other hand, there are many polar functional groups in the lignin molecule, which remarkably strengthen the lignin-bitumen interaction.
Based on the MSD results, the diffusion coefficient parameter can be obtained by Eq. (6) [18].  where N shows the total number of particles and a represents the slope value of the MSD-time fitting curve. It can be found that the slope value of MAD-time curve is highly related to the diffusion coefficient of molecules. The variation of MSD curve of SARA fractions as a function of simulation time is shown in Fig. 18. The initial nonlinear and final noisy parts of the MSD curves are removed and the related simulation time is utilized in the region from 20 to 100 ps. There is a significant linear relationship between MSD value of bitumen molecules and simulation time. These MSD curves fitting equations of SARA fractions are displayed in Table 3. Further, corresponding diffusion coefficient values are calculated and listed in Fig. 19.   As shown in Fig. 19, saturate molecules have the highest diffusion coefficient, which is approximately 2-3 times higher than that of other molecules in bitumen. Moreover, resin and asphaltene molecules exhibit the lowest diffusion ability, which is associated with their high molecular weight. The SARA fractions of LMB binders have lower diffusion coefficient values than virgin bitumen, demonstrating that adding lignin reduces the diffusion coefficient of bitumen molecules, especially for the saturate. Thus, the mixing temperature of lignin modified bitumen needs to rise to improve the diffusion coefficient of bitumen molecules and ensure the homogenous status of asphalt mix. The diffusion coefficient of bitumen molecules in LMB enhances as the lignin dosage decreasing. To the resin and asphaltene molecules, the influence of lignin concentration on the diffusion coefficient is limited.

Effect of lignin on mechanical properties
In this study, a mechanical simulation of LMB is subjected with a small strain of 0.1% by using NVT-equilibrium program [30]. Elastic modulus (E), bulk modulus (K) as well as shear modulus (G) of binders are calculated. The effects of lignin on the moduli of bitumen at 298 K are illustrated in Fig. 20(a). It is clear that the order of modulus magnitude is E > K > G for all bitumen samples. The addition of lignin has a positive effect on enhancing the moduli of bitumen, which enhances gradually with the increase of lignin dosage.
Noteworthy, when lignin content reaches to 30 wt%, the improvement effect of lignin on the modulus of bitumen is more obvious, which agrees with the aforementioned experimental results. If only considering the modulus of bitumen, the higher lignin content is recommended.
On the other hand, the lignin influence on the shear viscosity of bitumen at 433.15 K is explored and shown in Fig. 20(b), which is obtained from the equilibrium molecular model of bitumen after 300 ps NPT ensemble simulation. It is revealed that the shear viscosity of virgin bitumen is about 5.0cP. As expected, the addition of lignin could increase the viscosity of bitumen dramatically. With the increase of lignin dosage, the shear viscosity of LMB enhances dramatically, which results from the increase of lignin-bitumen interaction. The polar functional groups in the lignin molecule as well as the decreasing free volume reduce the molecular mobility of LMB. It is interesting to see that the viscosity of lignin modified bitumen would not increase when the lignin dosage is higher than 20 wt%, which can be found in the aforementioned simulation results in terms of the lignin influence on density, cohesive energy density as well as solubility parameter of bitumen.

Effect of lignin on self-healing behavior
Bitumen has the specific function of self-healing, which effectively relieves the cracking damage and prolongs the service life of asphalt road [31]. Meanwhile, previous studies also revealed that bitumen components or additives have great effects on its thermomechanical and self-healing properties [43][44][45][46]. In this research, the self-healing behaviors of studied binders at a microscale perspective are explored through MD simulation method, and the influence of lignin on the self-healing property of bitumen is evaluated. Firstly, the self-healing molecular models of LMB binders are built by using layer construction tool, and the representative self-healing model is illustrated in Fig. 21(a). In this model, the vacuum space of 10 Å is set in the middle of two bitumen layers, which represents the microcrack width of bituminous material. The geometry optimization and 200 NPT ensemble dynamics simulation are performed to simulate the self-healing process.
The density variation of binders as a function of simulation time can be tracked during the self-healing MD simulation. Take the virgin bitumen as an example, the density curve shown in Fig. 21(b) reveals that the initial density of the cracking model is about 0.5 g/ cm 3 due to the existence of microcrack. With the growth of healing time, the two separated layers close to each other, which results in the reduction of microcrack width and increase of the density of  whole system gradually. After NPT simulation for 50 ps, the density value keeps stabilized at 0.98 g/cm 3 , which is close to that of bulk model. In other words, morphological healing phenomena can be completed in a relative short-time scale and bitumen has abilities in terms of self-healing as well as restoration of physical properties. To quantify the influence of lignin on the self-healing ability of bitumen, the wetting time is calculated, which represents the time consumed when the microcrack between the separated layers of bitumen disappears. Fig. 21(c) presents the wetting time of all binders. With the increase of lignin dosage, the wetting time of LMB increases dramatically, which reveals that the addition of lignin weakens the self-healing speed of bitumen. When the lignin is added into bitumen, the molecular mobility and free space of bitumen model are both declined significantly, which hinders the selfdiffusion ability of molecules in binders. Compared to virgin bitumen, the wetting time of LMB binders with lignin dosage of 10, 20 and 30 wt% increases by 50, 120 and 210%, respectively.

Effect of lignin on adhesion behavior
The adhesion interaction between bitumen and aggregates significantly affects the durability and service life of asphalt pavement. The influence of lignin on the adhesion performance of bitumen is investigated with MD simulations at a microscale perspective. The silica unit (SiO 2 , see Fig. 22(a)) is utilized to create the molecular model of aggregate. Firstly, the silica unit is cleaved in the (0,0,1) direction, which is followed by geometry optimization and establishment of supercell structure, which is presented in Fig. 22(b). In order to evaluate the adhesion interaction between bitumen and mineral aggregate, the aggregate-bitumen layer model is built by putting the bitumen molecular model on the aggregate model. At the same time, a 50 Å vacuum layer is added upper to the bitumen molecular model to prevent the periodic boundary. After that, a 200 ps NVT dynamics simulation is performed to make the interfacial system stable. The molecular structures of the aggregate-bitumen layer model before and after MD simulations are displayed in Fig. 23. The bitumen molecules tend to adhere to the aggregate surface when the layer system reaches to equilibrium. To quantitatively assess the influence of lignin on the adhesion ability of bitumen, the interaction energy between aggregate and bitumen is calculated using Eq.7.
where E inter represents the interaction energy between the bitumen and aggregate, E Total refers to the total potential energy of the whole bitumen-aggregate interface model, E bitumen is the potential energy of the stable sole bitumen layer without aggregate molecules and E aggregate shows the potential energy of individual aggregate layer without the existence of bitumen molecules [18]. Obviously, the bitumen-aggregate layer with higher interaction energy E inter value would possess better adhesion property, and the larger outdoor energy is needed to separate the binder from aggregate. Fig. 24 indicates the interaction energy of bitumen-aggregate interfacial models of LMB with different lignin dosages. It is clear that the addition of lignin has a positive impact on increasing the interaction energy of bitumen-aggregate layer models and improving the adhesion ability of bitumen. To be more specific, when the lignin is added with the content of 10, 20 and 30 wt%, the interaction energy between aggregate and bitumen is enhanced by 30.8, 115 and 196%, respectively. From the perspectives of high temper-  ature properties and adhesion performance, lignin is a suitable modifier and substitution of bitumen. Meanwhile, the moisture, temperature as well as aggregate component also significant influence the adhesion property of LMB, which should be further investigated in the future work.
In addition, to validate the aforementioned MD simulation results and further assess the influence of lignin on the adhesion performance of bitumen, the experimental direct tension test is employed to measure the parameters of maximum tension stress and strain energy of binders by using the universal test machine (UTM). The preparation of adhesion samples and testing program are shown in the Fig. 25(a). A certain amount of flow bitumen is added into the middle gap between two limestone pillars to prepare the bitumen-aggregate interfacial samples. The thickness and diameter value of bitumen film is 0.2 mm and 20 mm, respectively. Before the tension test, the interfacial sample is conditioned in the UTM bin at 25℃ for at least four hours to guarantee the constant temperature. Then, the sample is installed in the fixture and subjected to the direct tension test with the stretching speed of 0.004 mm/s and strain rate of 2%/s, separately. The clear interfacial damage of bitumen-aggregate systems can be found after the tension test, which is also presented in Fig. 25(a).
The stress variation as a function of increasing strain for interfacial sample is recorded to investigate the effects of lignin dosage on the adhesion property of bitumen-aggregate systems, which is displayed in Fig. 25(b). The local fluctuation in the stress-strain curves is related to the instrument accuracy. It can be found that with the increase of strain, the stress value of bitumen-aggregate interfacial sample increases firstly and then decreases when the stress reaches to the maximum value. The smooth stress-strain curves of interfacial samples are illustrated in Fig. 25(c). To quantitatively evaluate the influence of lignin dosage on the adhesion properties of bitumen-aggregate interfacial systems, the maximum value and strain energy are calculated, which is associated with the peak value and integral area value of stress-strain curve, respectively. The influence of lignin on the adhesion parameters can be seen in Fig. 25(d). Compared to the virgin bitumen, the addition of lignin remarkably increases the maximum stress and strain energy value, indicating lignin is beneficial to improve the adhesion property of bitumen-aggregate interfacial system, which agrees to the MD simulation result. Meanwhile, with the increase of lignin dosage, the maximum stress and strain energy both enhances gradually. When the lignin content is 10, 20 and 30 wt %, the maximum stress value increases by 1.5, 3.1 and 8.8%, respectively, while the strain energy improves by 2.9, 33.8 and 55.2%,    separately. To this end, the MD simulation and experimental results both show the positive effect of lignin on improving the adhesion performance of bitumen-aggregate interfacial system. The high polarity characteristic of lignin contributes to the enlargement of molecular interaction of bulk bitumen and bitumenaggregate interfacial systems, which significantly intensifies the cohesion and adhesion properties of LMB binders.

Correlation analysis of micro-and macro-properties of LMB binders
It is clear that the macroscale mechanical properties mostly depend on the thermodynamics characteristics, chemical components as well as the microstructure of bituminous materials. Thus, it is of great importance to find the potential relationships between micro and macro properties, which is beneficial to predict the macro-properties of bitumen on basis of the chemical components and thermodynamics performance. In this study, the correlation analysis between thermodynamics characteristics from MD simulations (e.g. density, free volume fraction, cohesive energy density, solubility parameter and diffusion coefficient) and rheological properties (complex modulus, phase angle, recovery percent, creep compliance and fatigue life) of LMB binders are performed to find the potential relationship between micro-and macro-performance of bitumen samples with the Statistical Product and Service Solutions (SPSS) software. The detailed correlation coefficient values are calculated and shown in Table 4. When the correlation coefficient value between the two parameters is larger than 0.95, it means that there is a good relationship between these properties.
As illustrated in Table 4, the correlation coefficient values between density and recovery percent as well as non-creep compliance are both higher than 0.95, which indicates that there are great correlations between density and MSCR results. Regarding the free volume fraction at micro scale, some macroscale indexes (complex modulus and fatigue life) have great relationships with it. In addition, the correlation curves between density and recovery percent as well as creep compliance are displayed in Fig. 26(a) and (b), respectively. The linear correlation functions between density with logarithm value of recovery and creep compliance are shown in corresponding curves. For all binders, with the increase of density, the recovery percent enhances gradually, while the creep compliance reduces, which is related to the strengthened molecular interaction.
Moreover, it is obvious that there are great linear relationships between cohesive energy density with logarithm values of the selfdiffusion coefficient and complex modulus. The increased cohesive energy density of bitumen samples is beneficial to improve the complex modulus and has an adverse influence on molecular mobility and diffusion capacity. Additionally, Fig. 26(e) and (f) indicate that the free volume fraction of bitumen is closely associated with the parameters of complex modulus and fatigue life. With the free volume fraction increasing, the complex modulus of binder enhances dramatically, while its fatigue life is weakened. From the perspective of polymer science, when the free volume of a system starts to decline, the density and molecular interaction increases, which contributes to the improvement of modulus. However, due to few free spaces for molecules to move and relax, both stress concentration and fatigue damage development under the cyclic loading occur.

Conclusions and recommendations
This paper aims to investigate the effects of lignin on the thermodynamics and rheological properties of bitumen through the combination of MD simulation and experimental methods, and also preliminarily explores the correlation functions between micro and macroscale performance. The main conclusions are as follows: Lignin shows significant effects on improving the complex modulus and rutting factor of bitumen, while adversely affects the fatigue resistance. Meanwhile, lignin would increase the polar components and reinforce the molecular interaction, resulting in the modulus increase of bitumen. Molecular simulation results show that adding lignin would increase the density, cohesive energy density, modulus, viscosity and adhesion properties of bitumen. However, the free volume, diffusion coefficient and self-healing ability are deteriorated. There are great relationships between microscale and macroscale properties of LMB binders. Some macro-performance of LMB could be predicted and determined based on the thermodynamics properties obtained from the MD simulation directly.
In this study, the molecule model of lignin is limited, and the effects of various lignin components on the microstructure and physical properties of bitumen need to be considered. Moreover, the effects of lignin on the performance of asphalt mixtures would be studied in the future.

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