Size Effects of Finite Element Model for Three-Dimensional Microstructural Modeling of Asphalt Mixture

Asphalt mixture is a particulate composite material consisting of aggregate, mastic, and air voids.1e computed tomography (CT) image-based finite element approach is used as an effective method to simulate micromechanical response of asphalt mixture. For finite element analysis, the accuracy of the finite results is determined by the size of the finite element. In this paper, a voxel-based three-dimensional (3D) digital reconstruction model of asphalt mixture with the CT images after being processed was proposed. In this 3D model, the aggregate phase was considered as elastic materials while the asphalt mastic phase was considered as linear viscoelastic material. Four micromechanical digital models were generated, whose voxel sizes were 0.5mm, 0.67mm, 1.0mm, and 2.0mm, respectively. 1e four digital models were used to conduct uniaxial creep test for predicting creep stiffness modulus to investigate the effect of voxel size. Simulation results showed that the voxel sizes had a significant effect on creep stiffness modulus. For the creep simulation test, the most appropriate voxel size whose creep stiffness modulus changes within 2.5% is 1.0mm with regard to time steps, computational time, aggregate, and mastic shape representations.


Introduction
1.1.Background.Asphalt mixture is the most widely used road construction material in the paving industry.e traditional trial and error approach in industrial practice has focused on empirical experiments that develop correlations between the macro phenomena and the material characteristics.However, this traditional approach does not provide information that can be used to explain why certain mixtures perform better than others.And also this approach is cost intensive, high-resource consuming, and time consuming.To overcome these shortcomings, some microstructural methods and numerical techniques are developed to gain insight into the performance of asphalt mixture.
e composition of asphalt mixture can be considered to be included in three phases: the mastic phase (asphalt cement with a filler that is smaller than 2.36 mm in size), the aggregate phase, and the air voids.e mechanical behavior of asphalt mixture is complex due to the heterogeneity of the asphalt composite material.Several research studies have investigated the internal properties and performance of asphalt mixture at microstructural scale using the nondestructive X-ray CT technique [1][2][3][4][5][6][7][8][9].ere have also been a number of recent attempts to use X-ray CT images to reconstruct a three-dimensional (3D) finite element specimen to predict the mechanical properties of asphalt mixture [10][11][12][13][14].
e majority of these microstructure numerical samples for asphalt mixture face challenges with regard to accurate geometry of different phases.e geometry accuracy and volume composition of aggregate and mastic phases are directly affected by the resolutions of CT image for asphalt mixture.But most of the previous image-based finite element approaches are limited to low precision numerical models in order to reduce the time steps and computational time [15][16][17].However, it is not sufficient enough to consider the quantitative effects of finite element dimensions on asphalt mixture simulation results.

Objective.
e objective of this study is to investigate size effects of finite element model for three-dimensional microstructural modeling of asphalt mixture.Four imagebased microstructure models with different finite element dimensions for asphalt mixture were established in order to quantify the size effects of finite element model for threedimensional microstructural modeling.Firstly, the nondestructive industrial X-ray CT technique was used to capture the internal microstructure of the asphalt mixture specimen.e grayscale thresholds for dividing aggregate, asphalt mastic, and air voids were chosen based on the Otsu method.A morphological multiscale algorithm was applied to segment the coarse aggregate adhesion images.en, four voxel-based 3D digital reconstruction models of the asphalt mixture specimen is constructed based on the X-ray CT images.After that, the parameters of the asphalt mastic phase are obtained using the dynamic shear rheometer (DSR) test.Finally, the image-based numerical specimens are applied to conduct the creep simulations in order to investigate size effects of finite elements on creep stiffness.

Three-Dimensional Microstructural Reconstruction
2.1.CT Image Processing.At present, there is no image analysis software especially for material segmentation of CT slice images of asphalt mixture.e usual method involves segmentation of the internal materials according to the gray distribution characteristics of asphalt mixture.e most widely used approach is the Otsu method [18][19][20].
Because the difference of material density is not very obvious in the areas where the particles are in close contact with each other, cohesion of aggregate particles occurs when using the Otsu method to segment the aggregate and asphalt mastic phases.So, the adhesive aggregates are identified as individual particles by computer in the subsequent mechanical simulation, resulting in a calculation result which is inconsistent with the actual situation.erefore, before constructing the numerical model, adhesive aggregates in the images must be separated.At present, coarse aggregate adhesive separation is mainly based on manual processing, which is convenient for single image processing [21][22][23].However, it is not conducive for rapid and efficient image processing if the number of images is too large.erefore, an automatic separation algorithm for adhesive aggregates is necessary.
is section presents a morphological multiscale separation algorithm for adhesive aggregates based on different structural element radius.e algorithm first classifies the binary CT image in MATLAB before eroding in order to prevent smaller aggregates from being filtered out in the subsequent segmentation process.en, structural element with a minimum size is applied to perform ultimate erosion on all the labeled binary images.In the ultimate erosion process, if two or more divided objects appeared in the binary image, then they are classified as possible adhesive aggregates for the next analysis step.If two or more split objects are not observed, then the object could be considered as independent aggregate.Structural elements with different sizes are then selected as the ultimate erosion structural elements for the morphological algorithm of the adhesion particle binary map.According to the adhesive separation map based on morphological algorithm under different structural elements, the separation number with the largest number of division lines is taken as the actual separation of the binary image, and the separation image corresponding to the smallest structural element is taken as the actual separation image.Finally, all the separated marked images are superimposed to form a complete glue separation image.e algorithm flow chart is shown in Figure 1, which was implemented in the MATLAB.

Numerical Model Reconstructing Method.
e threedimensional microstructural reconstruction of asphalt mixture based on voxel is a reconstruction of complete asphalt mixture specimen from a series of asphalt mixture CT images. is reconstruction is a reverse solution process.For an M × N binary image, the basic constituent unit is a pixel.e expansion of a pixel in three dimensions is a voxel, as shown in Figure 2. A voxel is the most basic unit in the reconstruction of three-dimensional numerical samples.
After CT scanning of the asphalt mixture specimen, a series of continuous CT images can be obtained.ese continuous pseudocolor CT images were converted into a matrix of gray images, and the adjacent images were superimposed on each other, thereby expanding into threedimensional space to form a three-dimensional gray unit space network, as shown in Figure 3.
For finite element mesh models, the element, node number, and node coordinates must be specified.For a series of consecutive CT images, the unit and node number of the gray pixel matrix of the first CT image were defined as shown in Figure 4(a), and the unit and node number of the gray pixel matrix of no.k CT image were defined as shown in Figures 4(b) .Figure 5 shows a diagram of an initial threedimensional numerical model obtained in accordance with the above algorithm.
Figure 5 is a numerical model containing the background and mixture.e background, void, asphalt mastic, and aggregate of the CT image were segmented according to the relevant method in Sections 2.1 and 2.2.e segmentation process is shown in Figure 6.

Reconstruction of Numerical Model of AC-16 Asphalt
Mixture.In this paper, an AC-16 asphalt mixture test specimen was used for analysis.According to the above numerical modeling method, four microstructure numerical models with voxel lengths of 0.5 mm, 0.67 mm, 1 mm, and 2 mm, respectively, were reconstructed for AC-16 asphalt mixture.
e numerical models of four different voxel lengths for AC-16 asphalt mixture are shown in Figure 7.
e eight-node 3D solid reduction integration elements 2 Advances in Materials Science and Engineering (C3D8R) were adopted in constructing the mesh for the four types of numerical models.Numbers of units and nodes are shown in Table 1.

Material Parameters
Material properties are required for numerical modeling before simulation.e coarse aggregate can be considered to be linear elastic material due to its small deformation.In this paper, modulus and Poisson's ratio of the coarse aggregate were assumed to be 25 GPa and 0.25, respectively.Asphalt mastic was set as a linear viscoelastic material, and the generalized Maxwell model is often used as its constitutive model [24][25][26].e material properties of asphalt mastic in this numerical model adopted the 8-element generalized Maxwell model, and the parameters were obtained by       Advances in Materials Science and Engineering dynamic shear rheological test.Initially, the dynamic shear rheometer (DSR) test was performed at five temperatures (15, 40, 50, 60, and 70 °C), and the shear strain level was controlled within 0.15%.en, the Williams-Landel-Ferry (WLF) equation was applied to shift the measured shear modulus at test temperature to the reference temperature (25 °C in this study).e curve of the real part of the complex modulus can be described by the sigmoidal mathematical model: where G ′ is the real part of the complex modulus; δ, α, β, and c are the parameters of the equation; and t r is the loading time at the reference temperature.e real part of the complex modulus data was fitted by equation ( 1), as shown in Figure 8. e generalized Maxwell model parameters can be obtained (Table 2) according to the Prony series representation of the real part of the complex modulus as shown in equation ( 2) using the collocation method.
where m is the number of elements of the generalized Maxwell model; G e is the equilibrium modulus; G i is the relaxation strength of no.i; and ρ i is the shear relaxation time.
Given that the viscoelastic material was set to the generalized Maxwell model, which was the line viscoelastic constitutive, the load magnitude did not affect the analysis results.
erefore, the uniaxial creep simulation test was carried out with a dead load of 80 kPa.e displacementbased boundary conditions were applied to the asphalt mixture digital sample.e test specimen was completely fixed at one end to prevent the rigid body from being displaced, and the load was applied to the other end of the test specimen.For simplicity, the rigid body was used to simulate the loading plates at top and bottom of the specimen.A rigid body is a collection of nodes and elements whose motion is governed by the motion of a reference node.In this study, the dead load was imposed to the reference node, and the rigid body distributed the concentrated load on the top surface of the asphalt mixture digital model.e viscoelastic analysis step was used in this paper.e total analysis step time was 10 and the incremental step time was 0.1.

Computational Time Analysis.
As can be seen in Table 1, the length of voxel became smaller; the number of elements of the model and the number of nodes increased.Taking a model of about 1.7 million elements as an example, the minimum memory required for computer access is approximately 244 GB.So, it is difficult for a PC or a workstation computer to conduct the simulation.erefore, using supercomputer for massive parallel simulation is necessary.Simulation and analysis of the four numerical models in this study were based on the supercomputer Milkyway-2 in the National Supercomputer Center in Guangzhou.In this paper, massive parallel simulation experiments were conducted on the numerical models with a total of 5 nodes having cumulative capacity of 120 processors and 600 GB of memory.e computational time of each numerical model is shown in Table 3. e computational time consumed by different voxel length models varied greatly, and the maximum calculation time was 225.92 min.us, the computational time can be greatly reduced by increasing the number of applied nodes.

Creep Stiffness Modulus Analysis.
Creep stiffness modulus was used as the result analysis index in this paper.e   9. When the voxel lengths were 0.5, 0.67, and 1 mm, the creep stiffness modulus curves were very close.
Since Figure 9 cannot clearly demonstrate differences among the four numerical models' creep stiffness modulus, stiffness ratios were used to quantify their differences.e ratio of the creep stiffness modulus for each numerical model shown in equation ( 3) is defined as a fraction of the reference creep stiffness modulus (creep stiffness modulus with a voxel length of 0.5 mm). is is taken as a variation in the evaluation of each creep stiffness modulus curve, as shown in Figure 10.

�
Variation of stiffness modulus Stiffness modulus of the voxel side length of 0.5 mm . ( Figure 10 shows that a longer length of the voxel results in a higher rate of variation.When the length of the voxel was 0.67 mm, the ratio of the creep stiffness modulus from the first second to the tenth second was 0.55%-1.86%.When the length of the voxel was 1 mm, the ratio of the creep stiffness modulus from the first second to the tenth second was 2.14%-4.15%.When the length of the voxel was 2 mm, the ratio of the creep stiffness modulus from the first second to the tenth second was 6.56%-15.58%.us, when the accuracy of the CT image of the asphalt mixture reached a certain level, that is, when the length of the voxel is less than or equal to 1 mm, the ratio of the creep stiffness modulus within 10 s is less than 5%, that is, the length of the voxel had no significant effect on the creep test results.When the length of the voxel was 2 mm, the creep stiffness modulus curve was different from the other dimensions, with the maximum change of 15.58%.A smaller voxel size provides more accurate angular information of the aggregate.Voxel size has a great influence on the three-dimensional visualization of the asphalt mixture and the simulation results at local mesolevel.However, remarkably small voxel size showed no significant influence on the creep result.To reduce the computational time of simulation, a numerical model reconstructed with a pixel size of 1 mm can be recommended as an analytical model when performing the creep numerical simulation test.

Summary and Conclusions.
e size effect of finite elements in the three-dimensional microstructural numerical model of asphalt mixture was studied in this paper.Firstly, the industrial X-ray CT technique was used to obtain the CT image of the internal composition of the asphalt mixture.
e Otsu method was used to segment the CT image material phases, and a morphological multiscale segmentation algorithm was proposed to separate the adhesive aggregates after material classification.Finally, a three-dimensional numerical model reconstruction method for asphalt mixture based on voxel was proposed from a series of discrete asphalt mixture CT images.Four microstructural numerical models with voxel lengths of 0.5, 0.67, 1, and 2 mm were reconstructed.Uniaxial creep virtual test was implemented on the numerical model with four voxel sizes using the supercomputer Milkyway-2 in the National Supercomputer Center in Guangzhou to analyze the influence of the size effect on the simulation results, and the following conclusions were obtained: (1) A morphological multiscale algorithm for ultimate erosion of coarse aggregate adhesion image of asphalt mixture was proposed based on structural  (4) When the length of the voxel was less than 1 mm, the creep stiffness modulus changed within 2.5%. is finding suggests that when the voxel size is less than 1 mm, the influence of the voxel size on creep stiffness modulus is not significant.Moreover, voxel size mainly affects the computational time and the angular information of the aggregate.A numerical model reconstructed with 1 mm voxel size is recommended as the analytical model.

Recommendations for Future Work
(1) Many researchers have developed the low-resolution digital samples to simulate the behavior of asphalt mixture subjected to specific loading conditions, and this approach has been verified as a relatively reliable method to conduct the virtual test compared with the experimental test results.It is vital to validate the high-resolution simulation results with a series of experimental tests, and it will be the focus in our future papers.(2) e aggregate-mastic interface is a very challenging task in the asphalt mixture virtual simulation.For simplicity, the mastic phase was assumed to be perfectly bonded to the aggregate phase ignoring the internal cohesive behavior of the two phases in this paper.In order to obtain a more accurate digital sample, the interfacial conditions will be considered as part of our ongoing modeling efforts.(3) Asphalt mixture is a highly heterogeneous material consisting of aggregate, mastic, and air voids.e asphalt mastic has very complex temperature-, time-, and rate-dependent three-dimensional viscoelastic property.In this paper, the mastic is considered as a one-dimensional linear viscoelastic material.e underway paper will be focused on a three-dimensional viscoelastic constitutive model developed for asphalt mastic.

Figure 4 :Figure 6 :
Figure 4: (a) Element and node label for the first CT image pixel matrix; (b) element and node label for no.k CT image pixel matrix.

Figure 8 :
Figure 8: Fitting of the sigmoidal function.

Table 2 :
Prony series constants computed for the shear modulus.

Table 1 :
Four types of numerical model element and node information.

Table 3 :
Calculating time of different models.