Balancing the sustainable component of ethylene-vinyl acetate for achieved better compatibility improvement of wax-based warm mix additives in bitumen

.

The ethylene-vinyl acetate (EVA) polymers are always doped into waxy bitumen to inhibit network of wax crystals in bitumen. However, the compatibility improvement behaviors between wax-based warm mix (WWM) additives and bitumen by EVA are not clear, and the sustainable components of EVA for corresponding WWM additives to achieve better compatibility improvement are also not determined. This paper investigated compatibility improvement behaviors between commonly used WWM additives and bitumen after the addition of EVA to obtain sustainable components of EVA through experimental method of activation energy of viscous flow (AEVF) and density function theory-molecular dynamic (DFT-MD) calculations. The results show that the repulsions between the end of main-chain with highest electronegativity in WWM additives and polar molecules of EVA can alleviate the aggregation behaviors of WWM additives and EVA displays the best and worst compatibility improvements for additives with shortest and longest carbon chains, respectively. The dispersed asphaltenes combined with EVA can form the composite wax inhibitors (WIs) systems to increase diffusion

Introduction
In January 2020, the European Commission had passed the European Climate Law, which proposed to achieve the goals of carbon neutrality in 2050 [1][2][3][4]. The warm mixing technology was firstly developed by Shell and Kolo-veidekke companies in 1995, which had attracted a lot of attention of many researchers and engineers now. It also had the advantages of green environmental protection, energy saving, emission reduction and easy construction [5][6][7][8]. The application of warm mixing technology in the modified asphalt mixtures, reclaimed asphalt mixtures and other asphalt mixtures can reduce the viscosity and the mixing and compaction temperatures during construction, and these bring good environmental, economic and production benefits for the road construction [9][10][11][12]. Therefore, the warm mixing technology is an inevitable trend to achieve goal of carbon neutrality, which is also one of the important ways to achieve the development of green pavement in the worldwide.
However, the popularization and application of wax-based warm mix (WWM) additives in warm mixing technology are limited mainly due to the inferior low-temperature cracking resistance of warm mix asphalts (WMAs) [13][14][15][16]. According to the double edge notched tension (DENT) and extend bending beam rheometer (Ex-BBR) tests, the intermediate-temperature ductile and low-temperature cracking resistances of bitumen become worse after doping with the WWM additives [13,15], and these additives also increase the thermal stresses in bitumen [16]. Therefore, the addition of WWM additives that damages the low-temperature performance of bitumen becomes a concern of researchers. In order to solve this problem, some studies have determined the optimal average carbon numbers of WWM additives from the components perspective of additives through high-temperature gas chromatography (HTGC), multiple stress creep recovery (MSCR), DENT and Ex-BBR experiments, which achieves the better low-temperature cracking resistance of WMAs [13], and other studies introduce new kinds of WWM additives (polyethylene (PE)wax [17][18][19][20] and Siligate [21]) or combine the use of rubber crumb modifier (CRM) and WWM additives [22][23][24][25] to improve the cracking resistance of bitumen at low temperatures compared with the existed WWM additives. Interestingly, some researchers propose that the WWM additives have similar structures of the waxes in bitumen, which are mixtures of n-alkanes, and have concluded that poor compatibility between WWM additives and bitumen is the core factor for the inferior low-temperature cracking resistance of WMAs according to the wax precipitation theory, which leads to the phase separation and uneven state of WWM additives in WMAs [13,15,26,27]. Therefore, it is necessary to take measures to improve the compatibility between additives and bitumen from the essence of phase separation for WWM additives.
At present, there have been several investigations on the addition of wax inhibitors (WIs) into waxy bitumen to increase the compatibility between the wax and other components in bitumen [28][29][30][31]. The commonly used WIs including ethylene-vinyl acetate (EVA) copolymer [32,33] and Poly (methyl methacrylate) (PMMA) [34,35] are added to waxy bitumen, and the wax precipitation temperature (WPT) is characterized by the activation energy of viscous flow (AEVF) experimental method derived from dynamic shear rheometer (DSR) tests [29][30][31]. The addition of WIs can decrease the WPT of waxy bitumen, which indicates that the WIs can increase the compatibility of wax in bitumen [29][30][31]. Based on this, some researchers prepare nano-SiO 2 hybrid PMMA from the aspect of nanocomposite materials, which is added into the waxy bitumen [29]. The experimental results show that the nanocomposite materials can further reduce WPT to improve compatibility of wax compared with PMMA. Other researchers also prepare the ionic liquid (IL) modified with graphene oxide from the new perspective of dispersing asphaltenes, and the compatibility of wax in bitumen can also be further improved compared with the pure polymer [31]. Therefore, the addition of WIs can increase the compatibility of waxes in bitumen.
According to the previous publications, few investigations have focused on the methods of compatibility improvement for WWM additives in bitumen, and the researches on the compatibility improvement behaviors at the atom level are also blank. In this paper, the compatibility improvement of WWM additives induced by EVA as WIs is verified by the reductions of WPT from AEVF experimental method, and the charge properties of WWM additives and EVA monomer are characterized by density function theory (DFT) calculation. Based on these, the molecular models of the WWM additives-asphaltenes-EVA-solvent system are established, and the thermodynamic parameters about the interaction between EVA, asphaltenes and WWM additives are investigated to obtain the compatibility improvement behaviors by molecular dynamic (MD) simulation. The investigation can provide the new chemical method to improve the compatibility of WWM additives with bitumen and give insight into the compatibility improvement behaviors induced by EVA at the molecular level, which can provide the inspirations on how to choose the sustainable components of EVA to achieve high-efficiency compatibility improvement for the corresponding WMAs with different average carbon numbers.

Research objectives and methods
Through the combinations of compatibility experiment and DFT-MD simulation method, this paper investigates the compatibility improvement behaviors between WWM additives and bitumen after the addition of EVA to achieve the sustainable use of EVA, and the diagram of research objectives and methods are shown in Fig. 1. There are four objectives in this paper as follows.
(1) To determine the compatibility of WWM additives with virgin bitumen by the AEVF method from dynamic shear rheometer (DSR) tests. (2) To characterize the atom charges distributions of WWM additives and EVA monomer by DFT calculation (3) To establish the molecular models of WWM additivesasphaltenes-EVA -solvent system in MD simulation. (4) To investigate the effect of temperatures and asphaltenes contents on the crystallization behaviors of WWM additives at molecular level. (5) To characterize the compatibility improvement behaviors of WWM additives by EVA at molecular level through thermodynamic parameters including diffusion coefficient, binding energy and cohesive energy density.

Virgin bitumen
According to the previous literature, the virgin bitumen with different wax contents had different compatibilities of wax with other components of bitumen [36,37]. Therefore, two represent virgin bitumen with significant differences in wax contents were chosen to investigate the compatibility behaviors of WWM additives with virgin bitumen in this paper. One of virgin bitumen with less wax content from Sichuan, China had been widely used in Western China, which was labeled as SC, and other virgin bitumen labeled as GS with more wax content was purchased from GS Caltex Corporation, South Korea. The basic physical characteristics of virgin bitumen are shown in Table 1.

Commonly used WWM additives
The addition of WWM additives can reduce the construction temperatures and energy consumptions, which can weaken the chemical aging of bitumen during construction and effectively enhance the service life of asphalt pavement [5][6][7]. Four typical WWM additives with different average carbon numbers were chosen to add into the virgin bitumen to investigate the compatibility behaviors in this paper, and the detailed information and preparation process of these WWM additives could be found in the previous literature [16,[38][39][40]. In order to avoid the suspicion of commercialization, four commonly used WWM additives were labeled as LA wax, PE wax, FA wax, and FT wax, respectively. Fig. 2 demonstrates the appearance of these WWM additives, and the pertinent properties of additives in this paper are shown in Table 2. All commonly used WWM additives were purchased from Shanghai Aladdin Chemical Co., Ltd. (China), which were of analytical grade and directly applied without further purification.

Wax inhibitors
At present, the commonly used commercial WIs in bitumen industry are mainly EVA copolymers, which is always doped into the pure bitumen to disperse the wax crystals by changing the crystallization behaviors, which can effectively inhibit the formation of the threedimensional network structures of the wax crystals [32,33]. This will also decrease the viscosity and WPT of bitumen and increase the compatibility of waxes with bitumen. The EVA copolymers consist of the polyethylene (PE) segments along the backbone, and vinyl acetate separated the crystalline domains [41,42]. Therefore, the crystallinity is characterized by how often ethyl acetate branched off the backbone. The EVA was purchased from Shanghai Aladdin Chemical Co., Ltd. (China), which was of analytical grade and directly applied without further purification.
In order to investigate the effect of EVA on the compatibility of WWM additives, the EVA copolymers were doped into WMAs. According to the discussion of the optimal contents in the existing literature results, the content of EVA was set as 200 ppm in this paper [29,30].

Sample preparation
Firstly, two kinds of virgin bitumen were heated to 165 • C in an oven to fully flow. Then, the four WWM additives were added to the virgin bitumen at 1.5% and 3% by the weight of bitumen [13,15,16], and a high-speed shearing rate of 1000 revolutions per minute (rpm) was used to stir thoroughly for 30 min at high temperatures, which made the additives to completely dissolve in the virgin bitumen to form the homogeneous WMAs. According to the AASHTO T240 standard, all the WMAs were aged by the rolling thin-film oven (RTFO) at 163 • C for 85 min to reflect the short-term aging of bitumen at the stage of the pavement construction, and then aged for 20 h with the pressure aging vessel (PAV) according to the AASHTO R28 standard, which aimed to  Note: a PG=Performance grade.
simulate the long-term aging at the operation of pavement [43,44]. Finally, the corresponding contents of EVA were also added to the hot WMAs and stirred evenly into the bitumen at high temperatures.

Experimental characterization of compatibility improvement
It was reported that waxes in bitumen consisted of n-alkanes in the C 20 to C 40 range and large iso-alkanes and cycloalkanes soluble in nheptane [45]. The wax precipitation theory could be illustrated as the crystallization process of alkane molecules upon the cooling of bitumen. Initially, most alkane molecules were dissolved in bitumen at higher temperatures, which could lead to the good compatibility of alkanes with bitumen [29][30][31]. In the cold region, the cooling of bitumen would lead to the continuous decreases in the solubility of alkanes in bitumen, and the alkanes could precipitate from the bitumen to form the wax crystals, which resulted in the phase separation and damaged the compatibility of the waxes in bitumen [30,31].
The WWM additives had the similar structures of the wax in bitumen, which were the mixtures of n-alkanes [13,15,26,27]. Therefore, it was essential to use WPT to characterize the precipitation process of WWM additives, and the smaller WPT indicated the better compatibility of additives based on the wax precipitation theory in bitumen [29][30][31]. It was reported that the AEVF method was high-precision and low-cost to characterize the WPT of bitumen, which was conducted by the dynamic shear rheometer (DHR-3) in this paper, and the detailed procedures of theory and experiment could also be referred to previous literature [46][47][48]. On this basis, the diagram of experimental compatibility characterization for WWM additives is demonstrated in Fig. 3, and the compatibility improvement of WWM additives can be characterized as the following Eq. (1).
Where WPT 1 is the WPT of the WMAs before the addition of EVA, and WPT 2 is the WPT after doped with EVA. The bigger △WPT indicates that the better compatibility improvement and the experiments of each sample were repeated three times to obtain the average values of the results as the experimental values.

DFT theory
DFT was first proposed by Hohn P. and Kohn W. in 1964, which was a set of theoretical framework in the field of solid physics and chemistry, and it was also a theoretical method to describe the electronic structures and characteristics of solid materials [49][50][51]. It had gradually become one of the most important methods for calculating electronic structures and material properties in the field of theoretical computing materials science.
The core principle of DFT was to introduce the electron density in the space, replace the ground state wave function with the ground state electron density, and transform the complex problems of electron motion into the electron density for research. The total energy (E) of multielectron systems can be represented as function of the electron density (ρ) [50,51].
Where ρ(r) is the function of electron density, T, U, V xc are function of kinetic, coulomb and exchange-correlation energies, respectively, r is the atom radius, ▽', U'and V xc ' are the Hermitian operators of kinetic, coulomb and exchange-correlation energies, respectively, φ i is the eigenfunction of the operators, and ε i is the eigenvalue. In this paper, the DFT calculation was used to mainly obtain the atom charge distributions of EVA monomer and additives.

Calculation details and methods
Demol 3 is a module for calculating the molecular geometric configuration, energy, electronic structure and other properties based on quantum mechanics [52,53], and the DFT calculation was performed on Demol 3 module of Materials Studio 2020 in this paper. Firstly, the models of WWM additives were established through the n-alkane structures based on average carbon number (Table 2), and the structure of EVA monomer was referred to existed publications [32,33]. Then, the structures of WWM additives and EVA monomer were optimized through the hybrid functions of Becke's three-parameter exchange and Lee, Yang, and Parr's local density approximation (B3LYP), which is demonstrated in Fig. 4. The maximum convergence criteria for energy and displacement were set as 5.0 × 10 − 6 Hartree and 1.0 × 10 − 5 Å, respectively. Finally, the atom charges distributions of the entire molecules were calculated in the energy task of Demol 3 module. In addition, the quality, integration accuracy, self-consistent field were chosen as "fine" grid to achieve the high accuracy.

Molecular structures establishment
In order to characterize the compatibility improvement behaviors of WWM additives in bitumen induced by EVA, the mixing systems of WWM additives-EVA-asphaltenes-solvent were established by the Materials Studio 2020. According to the AEVF experimental method, the WWM additives also consider the n-alkanes with different average carbon numbers as shown in Table 2. The EVA with different carbon numbers of main-chain have different performance, and the carbon numbers were selected as 10, 40, 70, and 100, respectively to investigate the effect of carbon numbers on the compatibility improvement for the WWM additives with different average carbon numbers [32,33], which was labeled as EVA-10, EVA-40, EVA-70 and EVA-100, respectively. The contents of asphaltenes were set as 500 ppm, 1500 ppm, 2500 ppm and 3500 ppm, respectively to distinguish the existence states of asphaltenes in bitumen [54][55][56], and the molecular structure of asphaltenes was referred to the previous publication [57]. In addition, the methanol was selected as the solvent to consider asphaltenes aggregation. The molecular models of asphaltenes, WWM additives, EVA and solvent are illustrated in Fig. 5.
The mixing systems of WWM additives-EVA-asphaltenes-solvent were established in the amorphous cell module of Materials Studio 2020, and Fig. 6 demonstrates the amorphous cell models of various mixing systems according to the different kinds of WWM additives and EVA.

Force field selection
The Condensed-Phase Optimized Molecular Potentials for Atomistic Simulation Studies (COMPASS) is one of the common molecular force fields, which is often used in the calculation of MD simulation [58]. The COMPASS II force field is an important extension of the COMPASS force field based on the molecular mechanics, which can well consider the characteristics of electron cloud distributions, bond lengths and bond angles of atoms in molecules to calculate the interaction energy between atoms in molecules and simulate the mechanical properties of molecules [59]. Therefore, the interactions between the WWM additives, EVA and asphaltenes were characterized by the COMPASSII force field in this paper, and the theoretical method and parameter determination of the COMPASSII force field could be referred to the previous publications [60,61].

Simulation processes
The entire processes of MD calculation were carried out by the  Material Studio 2020, and the detailed steps of MD simulation were as follows.
(1) The components of WWM additives, asphaltenes, EVA and solvent were modeled in the Visualizer module, and then the singlemolecule structure was geometrically optimized in the Forcite module.
(2) Based on the single-molecular model, the COMPASS II force field was selected to establish the three-dimensional models of WWM additives-EVA-asphaltenes-solvent through the Amorphous Cell module. (3) The established preliminary three-dimensional models were also geometrically optimized through the COMPASS II force field and the Smart algorithm, and the maximum numbers of iterations were 20,000 times to minimize the energy of the models. (4) After the geometry optimization, the models were annealed to achieve the global energy re-minimization and more balanced geometry. The temperatures of the models were raised to 1000 K under the isothermal-isobaric ensemble (NPT ensemble, the unchanged atoms number, total pressure and temperature) and then gradually quenched to 298.15 K. The totals of 5 such annealing cycles were performed to eliminate local unreasonable structures in the model. (5) After annealing, the model was dynamically relaxed under the canonical ensemble (NVT ensemble, the unchanged atoms number, total volume and temperature). The 20 ns NVT simulation could make the simulation system obtain the equilibrium, which was consistent with previous publication, and the energy fluctuation satisfied the convergence criterion within 5% [31]. Therefore, the equilibrium phase simulation of 20 ns was carried out under the NVT ensemble to ensure the rationalizations of molecular structures at 298. 15

Thermodynamic parameters of MD models
In order to quantitatively characterize the crystallization behaviors of WWM additives and the compatibility improvement behaviors for WWM additives induced by EVA, the thermodynamic parameters including relative mass density, diffusion coefficient, binding energy and cohesive energy density were selected to perform the analysis.

Relative mass density
The relative mass density was an important parameter to characterize the local mass density of molecules in MD simulations [64]. Therefore, the aggregation states of WWM additives in bitumen could be well characterized by the local mass density.

Diffusion coefficient
The WWM additives, asphaltenes and EVA molecules were always in a state of dynamic equilibrium, that is, these molecules would diffuse into each other. The molecular diffusion could be used to evaluate the intermolecular compatibility, and the greater diffusion degree indicated the better compatibility between the molecules [65]. The mean square displacement (MSD) was used to characterize the molecular diffusivity of the WWM additives, which is demonstrated in Eq. (4) [65].
Where MSD is the mean square displacement of calculated particles, N is the number of calculated particles, t is the statistic time, and r j (t) and r j (0) are the positions of the particles' centroid at the initial moment and time t, respectively. In order to further quantify the molecular diffusivity of WWM additives, a relatively stable interval in the MSD curve was selected for fitting, which introduced the diffusion coefficient (D) of WWM additives. The calculation formula is shown in Eq. (5), and the D was used to quantitatively evaluate the compatibility behaviors of WWM additives from the perspective of molecules diffusions in this paper [65].

Binding energy
The interaction between molecules was also a key factor in the compatibility between molecules. According to the established models of EVA, asphaltenes and WWM additives, the interaction strengths between WWM additives and asphaltenes and EVA were calculated. In this paper, the binding energy E bind was used to characterize the interaction strength between molecules (Eq. (6)) [66]. Where E inter is the interaction energy of each simulated system, E ab is the total energy of the combined simulated systems a and b, E a is the total energy of the simulated systems a, E b is the total energy of the simulated systems b. When E bind was greater than zero, it represented that the force between the two molecules was attractive, and the larger values indicated the stronger intermolecular interaction and better compatibility. When E bind was less than zero, it meant that the force between the two molecules was repulsive, and the larger absolute values Fig. 7. WPT values of (a) SC and (b) GS virgin bitumen with different contents of WWM additives before and after additions of EVA.
represented the weaker intermolecular interaction and worse compatibility [66]. The binding energy could evaluate the intermolecular compatibility from the perspective of molecules binding.

Cohesive energy density
The cohesive energy density (CED) would represent the average energy required to separate molecule to an infinite distance, which was deduced by Eq. (7) [67].
Where E coh , E vdw , E elect and V are the total cohesive energy, van der Waals cohesive energy, electrostatic cohesive energy and volume of the individual simulated system, respectively.
It was reported that the CED of crystals in the amorphous state was smaller than that of the crystals at the oriented order, and the increases of CED were the dominant factor for the transition of crystals from the amorphous state to ordered state [68]. Therefore, the percentage increment of CED could characterize the crystallization behavior changes of WWM additives before and after the additions of EVA from the perspective of energy. Fig. 7 illustrates the WPT values of WMAs with and without the additions of EVA. It can be seen from Fig. 7 that the WPT values of GS bitumen system are higher than that of SC bitumen system, which can attributed to the higher wax contents in GS bitumen and is consistent with the previous literature [29][30][31]. The WWM additives can increase the WPT of bitumen and the higher contents of WWM additives will perform this more obvious effect, which leads to the uneven distribution of WWM additives in bitumen. The WPT values rankings for the four WWM additives modified bitumen is FT wax > FA wax> PE wax> LA wax, and this is also fit with previous publication [15] and can be explained that the WWM additives with higher average carbon numbers will have poor compatibility with bitumen. Interestingly, the WPT values of WMAs all decrease after the addition of EVA due to the repulsion of the EVA and WWM additives, which indicate that the EVA can enhance the compatibility of WWM additives with bitumen.

Experimental results of compatibility improvement for WWM additives
According to Fig. 7, it can been seen that the WPT values of SC pure bitumen increase by 0.8 • C from 18.2 • C to 19.0 • C, while the WPT values of GS pure bitumen change by 4.7 • C from 27.6 • C to 32.3 • C. This indicates that the EVA has the better compatibility improvement of waxes in bitumen with more wax contents, which is consistent with the previous results [29]. Fig. 8 demonstrates the WPT values changes of WMAs after the addition of EVA. The EVA presents the best and worst compatibility improvement for the LA wax and FT wax, respectively. Meanwhile, the compatibility improvement effects of PE wax and FA wax are in the middle, and the latter is worse than the former. This can be attributed to the stronger interaction of EVA on the WWM additives with smaller average carbon numbers. Interestingly, the increased contents of LA wax, PE wax and FA wax can enhance the compatibility improvement of EVA, but the compatibility improvement will be inhibited with the higher contents of FT wax. This can be explained that the supersaturated status of FT wax in bitumen will cause the crystal disorder to inhibit the compatibility improvement of EVA.

Atom charge distributions of WWM additives and EVA monomer
The atom charge distributions for WWM additives and EVA monomer are demonstrated in Fig. 9. According to Fig. 9, the two carbon atoms at the end of main-chain for the four WWM additives have the highest electronegativity values compared with other carbon atoms, and the six hydrogen atoms attached to these carbon atoms also have the highest electronegativity values. This is consistent with the results of previous publication [69]. Interestingly, the electronegativity values of carbon atoms other than terminal carbon atoms are significantly reduced, which can be attributed to the transfers of the covalent electron clouds and the torsion of bond angles. From the aspect of EVA monomer, the two oxygen atoms of acetate have the higher electronegativity values, which can form the polar covalent interactions with the higher electronegative carbon atoms between them. This will lead to the stronger polarity of acetate [69].
Interestingly, the carbon atoms at both ends of the double-bond exhibit two opposite charge properties of the electropositivity and electronegativity values, which can provide the polarity to the doublebond. This can be explained by the fact that the electropositive carbon atoms are influenced by the close electronegative oxygen atoms. It was reported that the polar and non-polar molecules had the strong repulsion reactions [69]. Therefore, the interactions between the end of carbon chain with highest electronegativity in WWM additives and polar molecules of EVA are the largest contributors to the repulsion of the EVA and WWM additives molecules, which can alleviate the aggregation behaviors of WWM additives.

MD model effectiveness validation
It was reported that the density differences between the experiment and simulation were within 10%, which could verify the accuracy of molecular models and force fields [62]. Table 3 illustrates the experimental and simulated density of n-C 28 and n-C 36 systems at different temperatures and pressures. It can be seen from the Table 3 that the density of bitumen decreases with the increases of temperatures, and the variation trends of density for the simulation and experiment are the consistent. In addition, the density obtained by simulation is all within 5% of the experiment, which validates the effectiveness of force field parameters and simulation methods.

Influence of temperatures on crystallization behaviors
The relative mass density curves of WWM additives-asphaltenessolvent mixing systems at different temperatures are shown in Fig. 10. According to Fig. 10, the relative mass density values of the WWM additives molecules are more evenly distributed along the path at 175 • C compared with other temperatures. This indicates that WWM additives are randomly and uniformly distributed at high temperatures, which can explain why the WWM additives are mixed with bitumen at high temperatures in the WMAs experiments [13]. At the same time, the relative mass density values of additives molecules reach the maximum at the center of the cell as the temperatures decrease from 175 • C to 25 • C, which indicates that the additives gradually aggregate towards the center of cell. The reductions of temperatures also increase the relative mass density values, and this can present that the WWM additives will aggregate together at low temperatures. In addition, the FT wax always has the highest relative mass density values, followed by FA wax and PE wax, and LA wax has the lowest relative mass density values, which is directly related to FT waxes with the largest average carbon numbers. Fig. 11 demonstrates the diffusion coefficient values of WWM additives-asphaltenes-solvent mixing systems at different temperatures. It is expected that the diffusion coefficient values of WWM additives decrease as the temperatures change from 175 • C to 25 • C due to the aggregations of WWM additives at low temperatures, and this result is also confirmed by the experimental conclusions that the diffusion coefficient values of pure alkanes reduce with decreasing temperatures [70].The LA wax displays the maximum diffusion coefficient values, which is fit with the conclusion from the relative mass density. Compared to the WWM additives at 25 • C, the diffusion coefficient values of LA wax, PE wax, FA wax and FT wax molecules at 175 • C increase by 3.25, 3.81, 3.92 and 3.52 times, respectively. Therefore, the temperatures have the greater influence on the diffusion coefficient of WWM additives with larger average carbon numbers.
The CED parameters of WWM additives-asphaltenes-solvent mixing systems at different temperatures are demonstrated in Fig. 12. As the temperatures decrease from 25 • C to 175 • C, the CED parameters of the mixing system all increase, which is consistent with the results of previous literature [65,66] and this also indicates that the CED values are required as the driving forces for WWM additives to change from the amorphous state to the ordered state. The CED parameters of the LA wax, PE wax, FA wax and FT wax mixing systems at 25 • C reduce by 26.17%, 19.38%, 17.24% and 15.01% compared with the mixing systems at 175 • C, respectively. Thus, the larger CED values changes will lead to the more severe aggregations of WWM additives with higher average carbon numbers at low temperatures. Fig. 13 shows the relative mass density curves of WWM additivesasphaltenes-solvent mixing systems at different contents of asphaltenes. It can be seen that the relative mass density values of WWM additives show a decreasing trend as the contents of asphaltenes increase from 500 ppm to 1500 ppm and then increase as the contents of asphaltenes continue to increase, which indicates that the asphaltenes contents less than 1500 ppm can inhibit the aggregations of WWM additives, and the asphaltenes contents greater than 1500 ppm will lead to the aggregations of WWM additives. The previous literature showed that the lower asphaltenes contents displayed the dispersed status of asphaltenes, but the asphaltenes would aggregate together as the asphaltenes contents increased to a certain amount [54]. Therefore, the dispersed asphaltenes can act as the WIs to hinder the crystal networks of WWM additives, and the aggregated asphaltenes will promote the  clusters of WWM additives. This can explain the interesting experimental phenomenon that the WPT decreases after the addition of asphaltene dispersants. The optimal asphaltenes contents for best improving the compatibility of WWM additives is 1500 ppm among these aspaltenes contents in this paper. The diffusion coefficient values of WWM additives-asphaltenessolvent mixing systems at different contents of asphaltenes are demonstrated in Fig. 14. It can be seen from Fig. 14 that the diffusion coefficient values of WWM additives increase as the asphaltenes contents increase from 500 ppm to 1500 ppm while the asphaltenes contents changing from 1500 ppm to 3500 ppm can lead to the reduction of diffusion coefficient values, which is consistent with conclusions from the above relative mass density. Compared to the diffusion coefficient values at the asphaltenes content of 500 ppm, the diffusion coefficient values of FT wax, FA wax, PE wax and LA wax increase by 25.88%,    Fig. 15 illustrates the diffusion coefficient values changes of WWM additives-asphaltenes-solvent mixing systems after the addition of EVA with different main-chain carbon numbers. It is expected that the additions of EVA can increase the diffusion coefficient values of WWM additives to improve the compatibility of additives. As for the EVA with different carbon numbers of main-chain, the EVA-10 can increase the diffusion coefficient values of LA wax, PE wax and FA wax most significantly among these EVA, and the diffusion coefficient value change of LA wax is the largest, followed by the PE wax and FA wax. Interestingly, the diffusion coefficient values changes of LA wax, PE wax and FA wax decrease as the carbon numbers of main-chain for EVA increase from 10 to 100. Therefore, the carbon numbers of the main-chain for EVA that are slightly less than the average carbon numbers of the WWM additive will help to improve the compatibility of the WWM additives, but a large gap between the carbon numbers of the main-chain for EVA and the average carbon numbers of the additives will lead to the opposite effect. This conclusion can be validated from FT wax modified bitumen. According to Fig. 15, the EVA-70 has the greatest compatibility improvement of FT wax among these EVA. This can also be attributed to the fact that the main-chain carbon number of EVA-70 is slightly smaller than the average carbon number of FT wax.

Binding energy results
In order to investigate whether EVA and asphaltenes have the synergistic effects, the binding energy values of asphaltenes and WWM additives, asphaltenes and EVA and EVA and WWM additives are demonstrated in Fig. 16. According to Fig. 16 (a), the asphaltenes-WWM additives has the binding energy to inhibit the aggregations of WWM additives, and the binding energy ranking of asphaltenes-WWM additives is consistent with the results of relative mass density and diffusion coefficient. Interestingly, there is also an interaction between EVA and asphaltenes by referring to the binding energy of asphaltenes-EVA, and the binding energy of asphaltenes-EVA increases from 21.97 kJ/mol to 36.24 kJ/mol when the carbon number of main-chain for EVA changes from 10 to 100. It can be seen from Fig. 16 (b) that the EVA-WWM additives also have the binding energy to inhibit the aggregation of WWM additives, which is higher than the binding energy of asphaltenes-WWM additives. This indicates that the EVA and asphaltenes can form the composite WIs system to synergistically inhibit the aggregation of WWM additives. In addition, the binding energy ranking of EVA-WWM additives is consistent with the results of above diffusion coefficient.

CED results
The percentage increment values of CED for WWM additivesasphaltenes-solvent mixing systems before and after the addition of EVA with different main-chain carbon numbers are illustrated in Fig. 17. According to Fig. 17, the EVA can reduce the percentage increment values of CED to disrupt the ordered degree of WWM additives at low temperatures from the aspect of energy, which makes the WWM additives more difficult to transform from the amorphous phase to the ordered crystalline phase according to the previous literature [29]. It can be also found that the EVA-10 can achieve the better compatibility improvement for LA wax, PE wax and FA wax, and the better compatibility improvement of FT wax can be obtained from the EVA-70, which is consistent with results of above diffusion coefficient.

Determining sustainable component of EVA for corresponding WWM additives
Based on the combinations between the results of experiments and DFT-MD simulation, the proposed compatibility improvement mechanism of WWM additives by EVA and asphaltenes is demonstrated in Fig. 18. The dispersed asphaltenes combined with EVA can form the composite WIs system through the binding interactions to repel the WWM additives molecules, which can inhibit the aggregation network of WWM additives and improve the compatibility of WWM additives with bitumen.   According to the effects of EVA with different main-chain carbon numbers on the crystallization behaviors of WWM additives in Section 5.4 and the principle of the minimum main-chain carbon numbers of EVA required by the average carbon numbers of wax molecules [71], the sustainable and effective components of EVA for corresponding WWM additives are determined, which is demonstrated in Fig. 19. The sustainable ranges of acetate percent for EVA for the corresponding LA wax, PE wax, FA wax and FT wax are referred to 25.88%− 31.72%, 17.35%− 21.86%, 14.89%− 18.91% and 5.01%− 6.56%, respectively, which can achieve the better compatibility improvement for the individual WWM additives.

Conclusions and Recommendations
In this paper, the experimental compatibility improvement was characterized by the AEVF method, and the atom charges distributions of WWM additives and EVA monomer and the crystallization behaviors of WWM additives without and with the additions of EVA were characterized through DFT-MD methods. On this basis, the compatibility improvement mechanism of WWM additives induced by EVA was proposed at the atom level to achieve the sustainable use of EVA. The conclusions can be drawn as follows.
(1) The interactions between the end of carbon chain with highest electronegativity in WWM additives and polar molecules of EVA are the largest contributors to the repulsion of the EVA and WWM additives molecules, which can alleviate the aggregation behaviors of WWM additives. (2) The WWM additives with longer carbon chains have the poor compatibility with bitumen, and the addition of EVA can improve the compatibility of WWM additives. The EVA presents the best and worst compatibility improvement for the additives with shortest and largest carbon chain, respectively. (3) The increases of temperatures can decrease the relative mass density and increase the diffusion coefficient of WWM additives. The dispersed asphaltenes can act as the WIs to hinder the crystal network of WWM additives, and the aggregated asphaltenes will promote the cluster of WWM additives. (4) The dispersed asphaltenes combined with EVA can form the composite WIs systems to increase the diffusion coefficient and Fig. 18. Proposed compatibility improvement mechanism of WWM additives by EVA and asphaltenes. Fig. 19. Selection of EVA with the sustainable acetate percent according to corresponding WWM additives.
H. Zhang et al. reduce the percentage increment values of CED to disrupt the ordered degree of WWM additives at low temperatures, which can synergistically improve the compatibility of WWM additives. (5) The sustainable carbon numbers of the main-chain for EVA that are slightly less than the average carbon numbers of the WWM additives will help to improve the compatibility of the WWM additives, but a large gap between the carbon numbers of the main-chain for EVA and the average carbon numbers of the additives will lead to the opposite effect.
This paper gives new insight into the microscopic behaviors of compatibility improvement for WWM additives with bitumen induced by EVA, which can provide the inspiration on how to choose the sustainable components of EVA to achieve high-efficiency compatibility improvement for the corresponding WMAs with different average carbon numbers. Our future work will further investigate the dynamic changes of the interface adhesion behaviors between the aggregates and WMAs and microscopic images before and after the additions of EVA through MD simulation and polarized light microscopy, respectively.

CRediT authorship contribution statement
Haopeng Zhang: Data collection, Data processing, Software, Methodology, Writing of draft manuscript. Shisong Ren: Software, Methodology, Revision of draft manuscript, Funds Support. Yanjun Qiu: Conceptualization, Supervision, Funds Support.

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.

Data Availability
Data will be made available on request.