Numerical Study on the Asphalt Concrete Structure for Blast and Impact Load Using the Karagozian and Case Concrete Model

The behaviour of an asphalt concrete structure subjected to severe loading, such as blast and impact loadings, is becoming critical for safety and anti-terrorist reasons. With the development of high-speed computational capabilities, it is possible to carry out the numerical simulation of an asphalt concrete structure subjected to blast or impact loading. In the simulation, the constitutive model plays a key role as the model defines the essential physical mechanisms of the material under different stress and loading conditions. In this paper, the key features of the Karagozian and Case concrete model (KCC) adopted in LSDYNA are evaluated and discussed. The formulations of the strength surfaces and the damage factor in the KCC model are verified. Both static and dynamic tests are used to determine the parameters of asphalt concrete in the KCC model. The modified damage factor is proposed to represent the higher failure strain that can improve the simulation of the behaviour of AC material. Furthermore, a series test of the asphalt concrete structure subjected to blast and impact loadings is conducted and simulated by using the KCC model. The simulation results are then compared with those from both field and laboratory tests. The results show that the use of the KCC model to simulate asphalt concrete structures can reproduce similar results as the field and laboratory test.


Introduction
The behaviour of structures or infrastructures under extreme loadings has become a hot topic in the area of civil, mechanical and material engineering.Critical infrastructures, such as runway pavement designed for normal aircraft landing and taking off, are expected to have adequate resistance when subjected to extreme loadings, such as impact or blast loadings (e.g., heavy airplane landing or taking off, air plane crash or terrorist attack).Asphalt concrete (AC) is made of bitumen binder and coarse aggregate.It is usually used as the surface course for both highway and runway flexible pavement [1].The dynamic load in the daily application for AC pavement normally corresponds to a strain rate less than 10 −1 s −1 .Tashman et al. [2] conducted experiments of AC under triaxial compressive loading at the strain rate from 10 −6 s −1 to 10 −3 s −1 .The results showed that the failure stress increased with the increase of the applied strain rate.Seibi et al. [3] studied AC subjected to uniaxial compressive loading with strain rate from 0.064 s −1 to 0.28 s −1 .It was found that the yield stress was significantly dependent on the strain rate.Park et al. [4] carried out tests on AC under uniaxial and triaxial compression with the strain rate changing from 10 −4 s −1 to 0.07 s −1 .The results showed that with the increase of the applied strain rates, the yield stress and failure stress increased, and the strain rate dependency was clearly showing up at the higher strains.It was also found that the viscous behaviour of AC decreased with the increase of the strain rate.However, when the pavement structure is subjected to the impact loading from further heavy traffic loading or aircraft loading, the corresponding strain rate exceeds 10 −1 s −1 , especially for heavy aircraft loading, the related strain rate reaches around 100 s −1 .Based on previous research [5,6], the compressive strength of AC material can be enhanced with the increase of strain rate, and AC material exhibits high plastic behaviour at the high strain rate.However, it is very expensive to conduct field test to investigate the actual behaviour of AC under severe dynamic loading, the numerical simulation is an effective alternative solution.
There are many factors that influence the reliability of numerical simulation.Among these factors, the material model plays a key role because it should reproduce the essential physical mechanisms of the material under various loading conditions.Seibi et al. [3] and Park et al. [4] used the Drucker-Prager yield function to simulate the compressive behaviour of AC under dynamic strain loading (strain rate from 0.0001 s −1 to 0.0701 s −1 ).Tashman et al. [2] developed a microstructure-based viscoplastic continuum model to take into account the effect of temperature in AC material with the strain rate ranging from 10 −6 s −1 to 10 −3 s −1 .Since the late 1990s, microstructural-based discrete element models (DEM) have been used for better understanding of asphalt pavement concrete (e.g., [7][8][9][10]).In DEM models, the assumption of elastic or viscoelastic behaviour was employed for the simulation of the static or creep behaviour of AC under normal traffic loading.However, it was found that when the pavement structure was under blast or high impact load, plastic deformation and severe damage would occur.Thus, a robust material model should be developed to consider the strain rate effect, strain hardening, strain softening and damage of the AC material under the severe dynamic loading.Recently, several concrete-like material models subjected to dynamic loadings have been developed, such as the Riedel-Hiermaier-Thoma (RHT) model [11], the Advanced Fundamental Concrete (AFC) model [12,13], the Karagozian and Case concrete model (KCC) [14] and the Holmquist-Johnson-Cook (HJC) concrete model [15,16].These robust material models are capable of capturing varying concrete-like materials' behaviour under different loading conditions.When subjected to dynamic loading, such as blast loading or high impact loading, concrete-like materials show a highly non-linear response.Besides, due to the general complexity of the constitutive models, the determination of the parameters (i.e., residual strength, failure strain and failure criteria in model) also plays an important role in achieving the actual performance of the concrete-like materials.This requires sufficient understanding of the modelling formulation and the associated considerations.In the current study, the KCC model [14] is used to simulate the AC material.This model is capable of capturing the varied concrete-like material behaviours under different loading conditions.However, it should be noticed that this model cannot consider the temperature effect.In this study, the dynamic behaviour of AC under a high strain rate is investigated, and the temperature effect is not considered.Hence, the KCC can be used; otherwise, the temperature issue should come up.
In this paper, several key features of the KCC model are firstly discussed, and then, the determination of the parameters in the KCC model for AC is provided.An application example on the AC structure under blast and impact loading is also illustrated and validated based on a field blast test and a laboratory drop weight impact test.

Review on KCC Model
When subjected to blast loading or high impact loading, concrete or other concrete-like materials shows a highly non-linear response.They usually exhibit pressure hardening and strain hardening under static loading and strain rate hardening in tension and compression under dynamic loading.A number of material models has been developed to model concrete-like materials recently [11][12][13][14][15][16].Among them, the KCC model is widely used to analyse concrete-like materials' response to blast and impact loading due to its simple implementation.In addition, the KCC model can capture the non-linear behaviour of the material under dynamic loading [14].The key features of the model are discussed briefly in the following section.

Strength Surface in KCC Model
The KCC model decouples stress into the hydrostatic pressure and deviatoric stress as shown in Equation (1): where σ ij is the stress tensor, s ij is the deviatoric stress tensor and σ ii is the hydrostatic pressure tensor.It should be noted that stress is positive in tension, and pressure is positive in compression.
The hydrostatic pressure is related to the volumetric change of material, while the deviatoric stress is related to shear resistance of the material and is usually expressed by the second invariant of the deviatoric stress tensor, J 2 : The KCC model has three independent strength surfaces: maximum strength surface, yield surface and residual strength surface, which are shown in Figure 1.The general formation of strength surfaces can be written as: in which ∆σ is the principal stress difference and p is the hydrostatic pressure.Usually, the above Equation (3) refers to the compressive meridian.The whole failure curve can be obtained through rotation of the compressive meridian around the hydrostatic pressure axis by multiplying r 3 (θ L ), which has the formation: where ψ = r t/r c and r t and r c are the radius of tensile and compressive meridian, respectively.According to Equation ( 5), it can be found that the r 3 (θ L ) depended on ψ and θ L .The parameter ψ in turn relies on the hydrostatic pressure.For the concrete-like material, the value of ψ varies from 1 2 at negative (tensile) pressures to unity at high compressive pressures [14].The value of Lode angle θ L can be obtained from: During the initial increase of hydrostatic pressure p, the deviatoric stresses ∆σ remains in the elastic region until the yield surface is reached.Deviatoric stress can be further developed until the maximum strength surface is reached, and the material will subsequently start to fail.After failure is initiated, the material will gradually lose its load-carrying capacity and reaches its residual strength surface.The formations of these three surfaces are also given in Equations ( 7)- (9).Yield surface ∆σ y = a 0y + p a 1y + a 2y p Maximum strength surface where the eight parameters, namely, a 0 , a 1 , a 2 , a 1 f , a 2 f , a 0y , a 1y and a 2y , for three surfaces can be obtained from the experimental data (e.g., triaxial compression test, biaxial compression test or uniaxial tension/compression test) [14].Residual strength surface where the eight parameters, namely, 0 1 2 1 2 0 1 , , , , , , a a a a a a a and 2 y a , for three surfaces can be obtained from the experimental data (e.g., triaxial compression test, biaxial compression test or uniaxial tension/compression test) [14].

Damage Factor in the KCC Model
After reaching the initial yield surface, but before reaching the maximum strength surface, the current surface can be obtained as a linear interpolation between yield surface ∆σy and maximum strength surface ∆σm: ( ) After reaching the maximum strength surface, the current failure is interpolated between the maximum strength surface ∆σm and the residual strength surface ∆σr, which is similar to the above computation: where η varies from zero to one depending on the accumulated effective plastic strain parameter λ.The value of η normally starts at zero and increases to unity at λ = λm, then decreases back to zero at some larger value of λ. λm is the plastic strain at maximum strength surface.The accumulated effective plastic strain λ can be expressed as follows: where ft is the quasi-static tensile strength, p dε is effective plastic strain increment and rf is the dynamic increase factor (DIF) of the material under dynamic loading.The damage factors b1 and b2 define the softening behaviour due to compression ( 0 P ≥ ) and tension ( 0 P < ), respectively.Parameter b1 can be determined by considering compressive energy Gc (area under the compressive stress-strain curve) obtained from the uniaxial compression test in single element simulation.It is obtained iteratively until the area under the stress-stain curve from single element simulation

Damage Factor in the KCC Model
After reaching the initial yield surface, but before reaching the maximum strength surface, the current surface can be obtained as a linear interpolation between yield surface ∆σ y and maximum strength surface ∆σ m : ∆σ = η(∆σ m − ∆σ y ) + ∆σ y (10) After reaching the maximum strength surface, the current failure is interpolated between the maximum strength surface ∆σ m and the residual strength surface ∆σ r , which is similar to the above computation: where η varies from zero to one depending on the accumulated effective plastic strain parameter λ.The value of η normally starts at zero and increases to unity at λ = λ m , then decreases back to zero at some larger value of λ. λ m is the plastic strain at maximum strength surface.The accumulated effective plastic strain λ can be expressed as follows: where f t is the quasi-static tensile strength, dε p is effective plastic strain increment and r f is the dynamic increase factor (DIF) of the material under dynamic loading.The damage factors b 1 and b 2 define the softening behaviour due to compression (P ≥ 0) and tension (P < 0), respectively.Parameter b 1 can be determined by considering compressive energy G c (area under the compressive stress-strain curve) obtained from the uniaxial compression test in single element simulation.It is obtained iteratively until the area under the stress-stain curve from single element simulation coincides with G c /h, where h is the element size.The damage factor b 2 in Equation ( 13) is related to tensile softening of the material and determined from experimental data.The fracture energy G f can be obtained from the uniaxial tensile test or three-point notched beam test.Then, the single element simulation of uniaxial tensile test is employed to obtain the stress-strain curve from the numerical analysis.Hence, the value of b 2 can be obtained until the area under the tensile stress-stain curve from a single element coincides with G f /w c , where w c is the localization width, and typically, w c is normally taken as one-to six-times the maximum aggregate size [14].
Based on Equations ( 10) and (11), the stress softening factors η and λ are governed by the accumulation of effective plastic strain.However, when the stress path is very close to the negative hydrostatic pressure axis, i.e., isotropic tension, wherein the hydrostatic pressure would decrease from zero to −f t , where no deviatoric stress occurred, no damage accumulation would occur based on these equations.However, in such concrete-like materials, damage cannot be avoided even at this state.Therefore, the above condition has to be modified by including pressure-softening effects near or after tensile failures.In this case, a volumetric damage increment is calculated and added to the total damage factor λ whenever the stress path is close to the triaxial tensile path.

Strain Rate Effect
The material model KCC also includes a radial rate enhancement on the material failure surface.This is because experimental data for concrete-like materials are typically obtained along radial paths from the origin in deviatoric stresses versus hydrostatic pressure via unconfined compressive and tensile tests.It is well known that the strain rate effect is important for concrete-like materials under severe dynamic loading.A typical DIF-strain rate curve for concrete-like materials can be obtained from the servo hydraulic fast loading tests and the split Hopkinson pressure bar (SHPB) test.

Equation of State
In addition to the strength surface model, an equation of state (EOS) is needed to describe the relationship between hydrostatic pressure and volume change of the material subject to dynamic load.EOS is usually determined using a fly impact (i.e., for steel) test or triaxial compressive test (i.e., for concrete or geomaterials).The isotropic compression portion of the KCC material model consists of pairs of hydrostatic pressure P and corresponding volume strain µ.It is implemented as a piece-wise curve in this model.

Determination of Parameters for Asphalt Concrete Material
The KCC model is employed to simulate AC material to capture dynamic response under impact and blast loading.This section describes the determination of the model parameter for the asphalt concrete with the compressive strength of 4.6 MPa, which is used in the following application example.For the asphalt concrete with a different compressive strength, this method can be used as a reference.

Strength Surface
As mentioned in Section 2.1, the KCC material model has three strength surfaces: strength, residual strength and yield surfaces.These three surfaces can be obtained through curve fitting of suitable experimental data.In this study, due to the few data for triaxial compressive tests of asphalt concrete, available data are extracted from Park et al. [4] with the compressive strength f c = 0.311 MPa for AC. Figure 2 presents the determination of the three surfaces by curve fitting for AC with f c = 0.311 MPa.The intersection point of maximum strength surface and residual strength surface is the so-called brittle-to-ductile point.This point can be determined by experimental data under high confining pressure.However, it is difficult to determine this point in the strength surface as no experimental data are available for AC materials.Based on the experimental data for concrete [17], this point was usually taken as p/f c = 3.878.Considering that the size and strength of aggregates used in AC and concrete material are almost the same, the brittle-to-ductile point for AC is taken to be the same as that for concrete.This value might be conservative for AC due to the higher content of coarse aggregate in AC.However, for simulation purpose, this value is acceptable.
aggregates used in AC and concrete material are almost the same, the brittle-to-ductile point for AC is taken to be the same as that for concrete.This value might be conservative for AC due to the higher content of coarse aggregate in AC.However, for simulation purpose, this value is acceptable.If new asphalt concrete with known unconfined compression strength , ' c new f is to be modelled, but its strength surfaces are otherwise unknown, then one way of scaling data from a known material is proposed as follows [14]: where ' , c old f is the unconfined compressive strength for a previously modelled AC.Then, the new material strength surface can be taken as: a a r a a a a r = = = .The new asphalt concrete with unconfined compressive strength fc = 0.8 MPa [2] is used to validate the parameters obtained from the scaling method.Figure 3 shows the maximum strength surface determined by the scaling method.It can be seen that the maximum strength surface fit very well with the experimental data, and thus, it can be concluded that the parameters for AC with different compressive strengths can be obtained by the scaling method.
In this study, the unconfined compressive strength for AC material is 4.6 MPa, and the tensile strength is 0.7 MPa at 35 °C.By scaling the data from the established curves given in Figure 2 [11], the appropriate strength surface of the current materials can be determined; the strength parameters are given in Table 1.If new asphalt concrete with known unconfined compression strength f c,new is to be modelled, but its strength surfaces are otherwise unknown, then one way of scaling data from a known material is proposed as follows [14]: where f c,old is the unconfined compressive strength for a previously modelled AC.Then, the new material strength surface can be taken as: in which a 0n = a 0 r, a 1n = a 1 , a 2n = a 2 /r.The new asphalt concrete with unconfined compressive strength f c = 0.8 MPa [2] is used to validate the parameters obtained from the scaling method.Figure 3 shows the maximum strength surface determined by the scaling method.It can be seen that the maximum strength surface fit very well with the experimental data, and thus, it can be concluded that the parameters for AC with different compressive strengths can be obtained by the scaling method.In this study, the unconfined compressive strength for AC material is 4.6 MPa, and the tensile strength is 0.7 MPa at 35 • C. By scaling the data from the established curves given in Figure 2 [11], the appropriate strength surface of the current materials can be determined; the strength parameters are given in Table 1.

Damage Factor
The strain hardening and softening pairs (η, λ) in Equations ( 10) and ( 11) describe the material behaviour transmitted from the yield surface to the maximum strength surface and from the maximum strength surface to the residual strength surface, respectively.During the transmission, parameter η varies from zero to one, depending on the accumulated effective plastic strain parameter λ.However, it is found that the original damage factor pairs (η, λ) in the KCC model are only suitable for concrete and not for the AC material due to AC having higher plastic failure strain.Thus, the input for accumulated effective plastic strain λ should be modified.Based on the uniaxial compressive test for AC, it was found that at peak stress, the corresponding strain was approximately 0.018, and the final failure strain was about 0.1; while for normal concrete, the corresponding strain at peak stress was around 0.0022.Hence, the λ is modified to give the high failure strain for AC in the current study.Additionally, it is found that when λ is adjusted to 10-times the original λ value, the numerical results seemed to show good agreement with the experimental results from the unconfined compressive test for AC. Figure 4 shows the modified and original series of (η, λ) pairs.It can be seen that the modified damage factor provided smoother descending than the original damage factor and had a higher failure strain that matched the behaviour of AC very well.

Equation of State
There are limited EOS data for AC material.The available EOS data are for AC with compressive strength of fc = 3.8 MPa [5].In this study, the compressive strength for AC is fc = 4.6 MPa.Thus, the pressure-volume pairs can be calculated using the volumetric scaling method [18].In this method, assuming that new data are obtained at the same volumetric strains, thus, the new corresponding pressure (pcnew) can be:

Equation of State
There are limited EOS data for AC material.The available EOS data are for AC with compressive strength of f c = 3.8 MPa [5].In this study, the compressive strength for AC is f c = 4.6 MPa.Thus, the pressure-volume pairs can be calculated using the volumetric scaling method [18].In this method, assuming that new data are obtained at the same volumetric strains, thus, the new corresponding pressure (pc new ) can be: and the new corresponding unloading bulk modulus (kµ new ) is: where r is the scaling factor, which is the ratio of compression strength for new material to the compression strength of the previous material modelled.Hence, the EOS data for f c = 4.6 MPa are calculated based on Equations ( 16) and ( 17), and the EOS inputted in the numerical model is shown in Figure 5.

Equation of State
There are limited EOS data for AC material.The available EOS data are for AC with compressive strength of fc = 3.8 MPa [5].In this study, the compressive strength for AC is fc = 4.6 MPa.Thus, the pressure-volume pairs can be calculated using the volumetric scaling method [18].In this method, assuming that new data are obtained at the same volumetric strains, thus, the new corresponding pressure (pcnew) can be: and the new corresponding unloading bulk modulus (kμnew) is: where r is the scaling factor, which is the ratio of compression strength for new material to the compression strength of the previous material modelled.Hence, the EOS data for fc = 4.6 MPa are calculated based on Equations ( 16) and ( 17), and the EOS inputted in the numerical model is shown in Figure 5.

Softening Parameters b 1 , b 2
The softening parameters (b 1 , b 2 ) shown in Equations ( 12) and ( 13) control the material softening behaviour after peak stress.These parameters are obtained from experiments, as detailed below.

•
Value of b 1 from the uniaxial compressive test: The uniaxial compressive test was conducted for AC according to ASTM 1074.The detailed test results and setup can be further referred to [11].Based on the test results, it is found that the corresponding strain at peak stress (f c = 4.6 MPa) was about 0.018, and the final failure strain was about 0.1, which was higher than that of concrete.The Young's modulus obtained from stain gauges attached at the middle height of the sample was 598 MPa.Based on experimental results, the compressive energy G c was calculated at 15.1 MPa•mm.Hence, for example, the b 1 value for a 10-mm mesh size was calculated as 3.45 using the method stated in Section 2.1.

•
Value of b 2 from the fractural test: The value of b 2 is determined by fracture energy G f , which can be obtained from the uniaxial tensile test or the three-point single-edge notched beam test (SNB).In the current study, the SNB test was carried out to evaluate fracture energy G f for the AC material.The detailed theory about the SNB test can be found in the established literature [19].Therefore, only the test result is presented here.In the SNB test, the compacted AC beam was fabricated with a dimension of 400 × 100 × 100 mm 3 depth.A mechanical notch was sawn with a depth of 20 mm, which gave a ratio of notch to beam depth of 0.2.The simply supported sample with a span length of 340 mm was tested under a 35 • C temperature.From the test, fracture toughness K IC can be obtained according to the formula suggested by Karihaloo and Nallathambi [19].Then, the fracture energy G f is calculated using: in which E is the elastic modulus and ν is Poisson's ratio.
The parameter b 2 is further determined by assigning fracture energy G f in the use of single element simulation of the uniaxial tensile test.The b 2 is then obtained via an iterative procedure until the area under the stress-stain curve from the single element simulation coincides with the value of G f /w c .The parameters obtained from SNB and single element simulation for AC (f c = 4.6 MPa) are summarized in Table 2.
Table 2. Parameters from the single-edge notched beam test (SNB) and single element simulation.

Parameters
Unit Value

Strain Rate Effect
The DIF curve for AC under different strain rates was obtained using servo hydraulic fast loading tests and the split Hopkinson pressure bar (SHPB) test in the current study.The strain rate produced by the servo hydraulic machine was approximately 10 −5 to 1 s −1 , and the higher strain rate loading was obtained through SHPB testing.The detailed setup and procedure for SHPB and the hydraulic test for AC can be referred to Wu [20].The DIF value for compressive (DIF c ) and tensile (DIF t ) strength obtained from the test are given as: However, for the compressive DIF curve, a numerical modelling of the SHPB test adopting this DIF curve found that the initial segment of this curve matched the experimental results very well, while the numerical model results for strain rate larger than 100 s −1 seemed to overestimate the stress.This can be due to the "double counting" of the inertia effect in the numerical modelling when the strain rate exceeded 100 s −1 .Hence, in the current model, the second segment in compressive DIF (Equation ( 19)) is ignored when the strain rate exceeds 100 s −1 .Beyond this, the DIF is assumed to remain a constant value.For the tensile DIF curves, in the macro-level numerical model, the KCC material model cannot capture the aggregate interlocking that propagates the micro-cracking and energy dissipation beyond the localization zone [20][21][22].Therefore, the above tensile DIF curve (Equation ( 20)) with two branches is used in the numerical model.The tensile and compressive DIF curves of asphalt concrete used in numerical model are summarized in Figure 6.
DIF curve found that the initial segment of this curve matched the experimental results very well, while the numerical model results for strain rate larger than 100 s −1 seemed to overestimate the stress.This can be due to the "double counting" of the inertia effect in the numerical modelling when the strain rate exceeded 100 s −1 .Hence, in the current model, the second segment in compressive DIF (Equation ( 19)) is ignored when the strain rate exceeds 100 s −1 .Beyond this, the DIF is assumed to remain a constant value.For the tensile DIF curves, in the macro-level numerical model, the KCC material model cannot capture the aggregate interlocking that propagates the micro-cracking and energy dissipation beyond the localization zone [20][21][22].Therefore, the above tensile DIF curve (Equation ( 20)) with two branches is used in the numerical model.The tensile and compressive DIF curves of asphalt concrete used in numerical model are summarized in Figure 6.

Field Blast Test
In this section, the AC material is used in the multi-layer pavement.Figure 7  In this section, the AC material is used in the multi-layer pavement.Figure 7 shows the cross-sectional view of the multi-layer pavement slab including 100 mm-thick engineering cementitious composite (ECC) at the bottom, 100 mm-thick high strength concrete (HSC) in the middle and 75 mm-thick AC at the top.The geogrid (GST) was placed at the middle of the AC layer to reinforce the AC material.The numerical model for the multi-layer pavement system under blast load is developed based on the configuration of the full-scale field blast test [23].Selected key features of the field blast test and numerical model are presented below.to reinforce the AC material.The numerical model for the multi-layer pavement system under blast load is developed based on the configuration of the full-scale field blast test [23].Selected key features of the field blast test and numerical model are presented below.In the field blast test, a bomb with the equivalent of a 7.3 kg Trinitrotoluene (TNT) charge weight was selected for testing and placed at about 170 mm above the surface of the pavement.The charge weight was evaluated based on the typical terrorist weapon attack.This multi-layer pavement slab was cast at site with a dimension of 2800 mm × 2800 mm × 275 mm (width × depth × thickness).The pavement slab was anchored to the ground to simulate the practical boundary condition, and a vertical anchor was installed at each corner.Figure 8 presents the multi-layer pavement slab before the blast load.
Various instruments were installed onto the slab to measure its responses during blast loading.Figure 9 shows the instrumentation installed on the pavement slab.Four accelerometers were installed at the middle of the side of the slab to measure both vertical (V1 and V2 in Figure 9) and horizontal accelerations (H1 and H2 in Figure 9).The accelerometers were mounted onto steel In the field blast test, a bomb with the equivalent of a 7.3 kg Trinitrotoluene (TNT) charge weight was selected for testing and placed at about 170 mm above the surface of the pavement.The charge weight was evaluated based on the typical terrorist weapon attack.This multi-layer pavement slab was cast at site with a dimension of 2800 mm × 2800 mm × 275 mm (width × depth × thickness).The pavement slab was anchored to the ground to simulate the practical boundary condition, and a vertical anchor was installed at each corner.Figure 8 presents the multi-layer pavement slab before the blast load.
Various instruments were installed onto the slab to measure its responses during blast loading.Figure 9 shows the instrumentation installed on the pavement slab.Four accelerometers were installed at the middle of the side of the slab to measure both vertical (V1 and V2 in Figure 9) and horizontal accelerations (H1 and H2 in Figure 9).The accelerometers were mounted onto steel frames that were cast together with the slab.Three total pressure cells (TPC) (TPC1, TPC 2 and TPC3 in Figure 9) were buried in the soil just below the slab to measure the pressure transferred from the pavement slab.A numerical model is established using LSDYNA package [24].In this model, the slab and foundation soil are discretised in space with one point Gauss integration eight-node hexahedron Lagrange element.Only a quarter of the slab is modelled with symmetric boundary conditions (as shown in Figure 10).The geogrid is spatially discretised with the shell element, and it is assumed that the geogrid is fully bonded within the AC layer.The anchors on the pavement slab are simulated as the fixed points (fixed in vertical direction) in the corresponding position in the numerical model.The non-reflection boundary is applied on the side and bottom of the foundation soil to model semi-infinite space.Based on the mesh study on the numerical model, the 10-mm element size is adopted for the pavement slab, geogrid and soil mass.The blast pressure is extracted from AUTODYN and used for numerical model analysis.The detailed process of applying pressure to the pavement surface can be referred to elsewhere [20,25].Various instruments were installed onto the slab to measure its responses during blast loading.Figure 9 shows the instrumentation installed on the pavement slab.Four accelerometers were installed at the middle of the side of the slab to measure both vertical (V1 and V2 in Figure 9) and horizontal accelerations (H1 and H2 in Figure 9).The accelerometers were mounted onto steel frames that were cast together with the slab.Three total pressure cells (TPC) (TPC1, TPC 2 and TPC3 in Figure 9) were buried in the soil just below the slab to measure the pressure transferred from the pavement slab.A numerical model is established using LSDYNA package [24].In this model, the slab and foundation soil are discretised in space with one point Gauss integration eight-node hexahedron Lagrange element.Only a quarter of the slab is modelled with symmetric boundary conditions (as shown in Figure 10).The geogrid is spatially discretised with the shell element, and it is assumed that the geogrid is fully bonded within the AC layer.The anchors on the pavement slab are simulated as the fixed points (fixed in vertical direction) in the corresponding position in the numerical model.The non-reflection boundary is applied on the side and bottom of the foundation soil to model semi-infinite space.Based on the mesh study on the numerical model, the 10-mm element size is adopted for the pavement slab, geogrid and soil mass.The blast pressure is extracted from AUTODYN and used for numerical model analysis.The detailed process of applying pressure to the pavement surface can be referred to elsewhere [20,25].In the numerical model, AC, ECC and HSC are all grouped as concrete-like materials and modelled by the KCC model [20,23].The basic parameters for these three materials are listed in Table 3.The parameters for AC material in the KCC model can be found in Section 3.While the process of the determination of HSC and ECC parameters for KCC model is the same as that mentioned in Section 3, the detailed parameters can be referred to Wu [20].The Drucker-Prager model and plastic-kinematic model [17] are employed to model the foundation soil and geogrid, In the numerical model, AC, ECC and HSC are all grouped as concrete-like materials and modelled by the KCC model [20,23].The basic parameters for these three materials are listed in Table 3.The parameters for AC material in the KCC model can be found in Section 3.While the process of the determination of HSC and ECC parameters for KCC model is the same as that mentioned in Section 3, the detailed parameters can be referred to Wu [20].The Drucker-Prager model and plastic-kinematic model [17] are employed to model the foundation soil and geogrid, respectively.The AUTOMATIC_SURFACE_TO_SURFACE contact algorithm is employed to model the interaction between pavement slab and soil.The TIEBREAK contact algorithm in LSDYNA is used to simulate the interface behaviour between HSC and AC layers.The parameters for the foundation soil, geogrid and interface property are in Tables 4-6, respectively.The results of the numerical modelling of the multi-layer pavement under blast loading, with the incorporation of the above-mentioned material models, are summarized and compared with the blast test results.In the numerical results, the fringe level in the damage contour is the value for the scaled damage indicator δ, which is defined to describe the damage level of the material [14,20,23].A scaled damage indicator δ is related to the effective plastic strain λ in the material: (i) at the yield surface, λ = 0, leading to δ = 0; (ii) at the maximum strength surface, λ = λ m , leading to δ = 1; and (iii) at the residual strength surface, λ = λ r >> λ m , leading to δ = 1.99 ≈ 2. Thus, the δ value moving from 0 to 1 to 2 indicates that the failure surface migrates from the yield surface to the maximum strength surface and to the residual strength surface, respectively, as the material being stressed.In this study, when the residual strength of material reduces to 20% of its peak strength, the material seems to suffer severe failure.The plastic strain corresponded to that residual strength used to calculate the delta value.Both the laboratory and field test for AC material indicated that when δ was greater than 1.8, the material would be severely damaged.
The damage situation for the multi-layer pavement slab in field blast test is shown in Figure 11. Figure 11a shows that the blast pressure destroyed the upper half of the AC layer above the GST reinforcement.It is also noted that only the centre of the GST piece was burned off during the blast event.Figure 11b shows the resulting damage on the HSC layer after removing the top layer of asphalt.From this figure, it could be seen that the crater was very shallow and did not punch through the whole layer, and a crater of around 700 mm in diameter and a depth of 10 mm was formed on the HSC layer.The damage situation for the multi-layer pavement slab in field blast test is shown in Figure 11. Figure 11a shows that the blast pressure destroyed the upper half of the AC layer above the GST reinforcement.It is also noted that only the centre of the GST piece was burned off during the blast event.Figure 11b shows the resulting damage on the HSC layer after removing the top layer of asphalt.From this figure, it could be seen that the crater was very shallow and did not punch through the whole layer, and a crater of around 700 mm in diameter and a depth of 10 mm was formed on the HSC layer.The damage contour for the AC layer in numerical model is shown in Figure 12a.Comparing Figures 12a and Figure 11a, it is observed that the damage pattern in the numerical model is symmetrical, while that in the field measurement was skewed.This is because the bomb in the field was not placed at the centre of the slab, and one side of the AC layer was more severely damaged than the other.Shear cracking near the anchor point was observed in the numerical model, which The damage contour for the AC layer in numerical model is shown in Figure 12a.Comparing Figures 11a and 12a, it is observed that the damage pattern in the numerical model is symmetrical, while that in the field measurement was skewed.This is because the bomb in the field was not placed at the centre of the slab, and one side of the AC layer was more severely damaged than the other.Shear cracking near the anchor point was observed in the numerical model, which was similar to the experimental observations in the field test.It could be concluded that the basic failure pattern given by the numerical model agrees well with the results obtained from the field-testing.Figure 12b shows the damage pattern for the HSC layer.Comparing Figure 11b with Figure 12b, the damage pattern for HSC is very consistent between field measurement and numerical results.The diameter of the crater was about 750 mm in the numerical model, which is very close to that of the blast test result.As shown in Figure 12b, shear cracks are also observed near the anchor points.Based on the damage pattern in the field blast test, the crater on the top face of the HSC is shown to be shallow and with a thickness of less than 10 mm.However, after cracking occurred at the bottom of the HSC layer, the numerical model shows that the bottom of the HSC layer has experienced severe cracking.This might be due to the combination of the bending of the HSC layer under the blast load and the reflection of the stress wave at the bottom interface.In the numerical model, the interface between HSC and ECC is assumed to be fully bonded.However, ECC is more flexible than HSC, and thus, it would cause tensile stress at the bottom of the HSC layer when deformed together.The compression stress wave from the top face would also travel within the HSC layer and reflect as a tension stress at the interface, which could cause spalling.Based on the damage pattern in the numerical model, the HSC layer might be considered having failed, while the field observation suggests that HSC may have partially failed.The results from both the blast testing and the numerical simulation are comparable.The maximum difference of vertical acceleration between the blast testing and the numerical model is about 10%, and the numerical model predicted slightly higher in the vertical acceleration than that of the blast test.The pressure values in the corresponding points in the numerical model are compared with pressures obtained from the blast test, as summarized in Table 8.The layout of the total pressure cell in the blast could be referred to Figure 9.The pressure values from the numerical simulation are shown to be close to that from the blast test for TPC2; while for TPC3, it has a 20% discrepancy with the numerical simulation considering the inherent variation in the blast test.TPC1 was damaged during the blast test, and hence, no pressure reading was recorded from it.The numerical model predicts that the pressure might be as high as 13 MPa at that point, which is far beyond the In the field blast test, the four accelerometers were installed at the mid-side of pavement slab (as shown in Figure 9).These accelerometers were used to measure the vertical and horizontal acceleration of the pavement slab subjected to blast loading.For the horizontal acceleration, due to the centre of the charge being closer to one side of the pavement slab; there were two different horizontal acceleration readings; while in the numerical model, it was assumed that the explosive occurred in the centre of the pavement slab.Thus, in this section, only the vertical acceleration from the field blast test was compared with that of the numerical model.In the numerical model, the raw nodal acceleration contained considerable numerical noise.The ELEMENT_SEATBELT_ACCELEROMETER could be used to eliminate numerical noise and obtain more accurate node acceleration.The vertical acceleration from the blast testing is compared with that of the numerical model as shown in Table 7.The results from both the blast testing and the numerical simulation are comparable.The maximum difference of vertical acceleration between the blast testing and the numerical model is about 10%, and the numerical model predicted slightly higher in the vertical acceleration than that of the blast test.The pressure values in the corresponding points in the numerical model are compared with pressures obtained from the blast test, as summarized in Table 8.The layout of the total pressure cell in the blast could be referred to Figure 9.The pressure values from the numerical simulation are shown to be close to that from the blast test for TPC2; while for TPC3, it has a 20% discrepancy with the numerical simulation considering the inherent variation in the blast test.TPC1 was damaged during the blast test, and hence, no pressure reading was recorded from it.The numerical model predicts that the pressure might be as high as 13 MPa at that point, which is far beyond the maximum measurement capacity of the pressure cell installed.That can explain why TPC1 was destroyed due to the overwhelming blast loading.The multi-layer pavement slab was subjected to an 1181-kg drop weight impact.The cross-section of the multi-layer pavement is the same as that given in Figure 7.The basic mechanical properties of the, HSC, ECC and AC are determined according to ASTM standards, and the results are summarized in Table 9.The drop weight was a cylindrical projectile with a hemispheric head (100 mm in diameter and made with high strength steel), and the pavement slab was subjected to two times of impact from the same drop height of 1.5 m.During the test, the multi-layer pavement slab was placed on the top of compacted soil/sand in a steel strong box.Directly below the slab was the geocell (MiraCell MC-100) which was filled with compacted soil/sand.This was to enhance the strength of the soil/sand layer and to provide a solid sub-base as to simulate the practical condition.The setup for the multi-layer pavement slab is given in Figure 13.Various instruments were also installed to monitor the response of the pavement during the drop weight test.Figure 14 shows the positioning of the potentiometers.A photodiode system was used to trigger the data acquisition system during the test.It consists of two photodiodes and two laser sources placed 100 mm vertically apart.The data acquisition system would be triggered when the falling projectile crosses the top laser emitter.Impact velocity could be determined using the time interval that the projectile took to cross the second laser emitter.positioning of the potentiometers.A photodiode system was used to trigger the data acquisition system during the test.It consists of two photodiodes and two laser sources placed 100 mm vertically apart.The data acquisition system would be triggered when the falling projectile crosses the top laser emitter.Impact velocity could be determined using the time interval that the projectile took to cross the second laser emitter.In the numerical model, the AUTOMATIC_SURFACE_TO_SURFACE contact algorithm is employed to model the interaction between pavement slab and soil.The contact algorithm TIEBREAK is also used to simulate the interface behaviour between the AC and HSC layer, which is the same as that in the numerical model for blast loading as mentioned in Section 4.1.2.As the impactor might penetrate the AC layer, the erosion method is employed in the simulation, in which the maximum plastic failure strain of 0.2 is used to delete the distorted element once the actual strain exceeds this predefined failure strain.The multi-layer pavement slab, the drop weight head and soil mass are discretised in space with one point Gauss integration eight-node hexahedron elements.The geogrid is simulated with the four-node Belytschko-Tsay shell element that allows  no bending resistance.In addition, the PLASTIC_KINEMATIC material model is employed to describe the bi-linear behaviour of the geogrid under tensile loading [23].The parameters used for the geogrid are the same as given in Table 4.The Drucker-Prager model is used to simulate soil mass.In the laboratory drop weight test, the upper soil layers were compacted and reinforced with geocell material, which would enhance the strength of the soil, and the lower layer has no reinforcement.Hence, in the numerical model, it is necessary to consider the function of the geocell material.From the laboratory test on the geocell-encased sand [26], it is observed that the geocell confinement did not change the friction angle of soil while significant cohesion occurred in the granular soil, which indicated that for the geocell-reinforced sand layer, the strength and stiffness behaviour of the soil would be enhanced.However, in the numerical model, it is difficult to model and mesh the geocell material due to its complex geometry.Hence, it would be preferable to use the composite model to consider the enhancement of the shear strength and stiffness of the geocell-reinforced sand layer.Madhavi et al. [27] purposed an empirical equation to calculate Young's modulus of the geocell-reinforced sand using the secant tensile modulus of the geocell material and Young's modulus parameter of the unreinforced sand, which could be expressed as: in which   is the Young's modulus of the geocell-reinforced sand, M is the secant modulus of the In the numerical model, the AUTOMATIC_SURFACE_TO_SURFACE contact algorithm is employed to model the interaction between pavement slab and soil.The contact algorithm TIEBREAK is also used to simulate the interface behaviour between the AC and HSC layer, which is the same as that in the numerical model for blast loading as mentioned in Section 4.1.2.As the impactor might penetrate the AC layer, the erosion method is employed in the simulation, in which the maximum plastic failure strain of 0.2 is used to delete the distorted element once the actual strain exceeds this predefined failure strain.The multi-layer pavement slab, the drop weight head and soil mass are discretised in space with one point Gauss integration eight-node hexahedron elements.The geogrid is simulated with the four-node Belytschko-Tsay shell element that allows no bending resistance.In addition, the PLASTIC_KINEMATIC material model is employed to describe the bi-linear behaviour of the geogrid under tensile loading [23].The parameters used for the geogrid are the same as given in Table 4.
The Drucker-Prager model is used to simulate soil mass.In the laboratory drop weight test, the upper soil layers were compacted and reinforced with geocell material, which would enhance the strength of the soil, and the lower layer has no reinforcement.Hence, in the numerical model, it is necessary to consider the function of the geocell material.From the laboratory test on the geocell-encased sand [26], it is observed that the geocell confinement did not change the friction angle of soil while significant cohesion occurred in the granular soil, which indicated that for the geocell-reinforced sand layer, the strength and stiffness behaviour of the soil would be enhanced.However, in the numerical model, it is difficult to model and mesh the geocell material due to its complex geometry.Hence, it would be preferable to use the composite model to consider the enhancement of the shear strength and stiffness of the geocell-reinforced sand layer.Madhavi et al. [27] purposed an empirical equation to calculate Young's modulus of the geocell-reinforced sand using the secant tensile modulus of the geocell material and Young's modulus parameter of the unreinforced sand, which could be expressed as: (21) in which E r is the Young's modulus of the geocell-reinforced sand, M is the secant modulus of the geocell material at axial strain 2.5% in kN/m and σ 3 is the confining pressure from the geocell in kPa.k u is the dimensionless modulus parameter of the unreinforced sand, which is a modulus number in the hyperbolic model developed by Duncan and Chang [28].The confining pressure σ 3 could be calculated as: where D 0 is the initial diameter of the geocell and ε a is the axial strain of the geocell at failure; the induced cohesion in the geocell-reinforced sand is then related to the increase in the confining pressure σ 3 : in which c r is the enhanced cohesion and k p is the coefficient of passive earth pressure.
In the current study, the geocell MC-100 from Polyfelt is used.The geocell dimension is a rhombus with two diagonal lengths of 203 mm and 244 mm.Thus, the equivalent diameter is calculated as about 177.5 mm.The secant modulus M for the geocell is obtained as 278 kN/m from the tensile test [20].Additionally, the value of ε a is taken as 4.8%.The modulus parameter k u for unreinforced sand in the current study is taken as 727 MPa according to the curve fitting from the triaxial test, and hence, the confining pressure, enhanced cohesion and Young's modulus for geocell-reinforced sand could be calculated based on Equations ( 21)- (23).The parameters for the unreinforced sand and geocell reinforced sand are summarized in Table 10.During the impact test, the deformation of the drop weight head was negligible compared to the deformation of the pavement slab.Hence, the drop weight head is modelled with a rigid body in the current study.For the configuration of the drop weight head, the simple cylindrical shape is modelled instead of modelling the head with weight mass (as shown in Figure 15).The simple cylindrical head has a diameter of 100 mm with a length of 1292 mm.The total mass for the simple cylindrical head is about 1181 kg, from which the density of the drop head would be obtained.The properties for the drop head are listed in Table 11.The convergence study is conducted, and it was found that a 5-mm element size gave a stable response, which is therefore applied for the simulation.The numerical model of the multi-layer pavement slab under drop weight impact load is given in Figure 16.
Friction angle, φ (°) 40 40 During the impact test, the deformation of the drop weight head was negligible compared to the deformation of the pavement slab.Hence, the drop weight head is modelled with a rigid body in the current study.For the configuration of the drop weight head, the simple cylindrical shape is modelled instead of modelling the head with weight mass (as shown in Figure 15).The simple cylindrical head has a diameter of 100 mm with a length of 1292 mm.The total mass for the simple cylindrical head is about 1181 kg, from which the density of the drop head would be obtained.The properties for the drop head are listed in Table 11.The convergence study is conducted, and it was found that a 5-mm element size gave a stable response, which is therefore applied for the simulation.The numerical model of the multi-layer pavement slab under drop weight impact load is given in Figure 16.For the 1.5-m drop weight impact, the drop weight head is assigned with 5.02 m/s for the first impact.For the second impact, due to the penetration of the AC layer in the first impact, the distance between the laser diode system and the face of the pavement slab is increased, while in the numerical model, the impact head is just placed at the position right before reaching the surface.Hence, the initial velocity in numerical model for the second impact is determined as the sum of the experimental recorded velocity and velocity caused by gravity acceleration.Thus, the velocity is calculated as 5.06 m/s for the second impact with the gravity acceleration of 9.8 m/s 2 .
After the first impact, the fully restarted method in LSDYNA is used to conduct the second impact simulation.At the beginning of the second impact, the stress and residual velocity within the pavement slab and the velocity from the first impact are set to be zero.The damage factor and plastic strain is retained in order to check the accumulated damaged behaviour after the second  For the 1.5-m drop weight impact, the drop weight head is assigned with 5.02 m/s for the first impact.For the second impact, due to the penetration of the AC layer in the first impact, the distance between the laser diode system and the face of the pavement slab is increased, while in the numerical model, the impact head is just placed at the position right before reaching the surface.Hence, the initial velocity in numerical model for the second impact is determined as the sum of the experimental recorded velocity and velocity caused by gravity acceleration.Thus, the velocity is calculated as 5.06 m/s for the second impact with the gravity acceleration of 9.8 m/s 2 .
After the first impact, the fully restarted method in LSDYNA is used to conduct the second impact simulation.At the beginning of the second impact, the stress and residual velocity within the pavement slab and the velocity from the first impact are set to be zero.The damage factor and plastic strain is retained in order to check the accumulated damaged behaviour after the second impact.The simulation of the second impact is carried out when the downward velocity of impactor reached zero, as it is very time consuming to continue to simulate the vibration of pavement slab after impact.It should be noticed that the pavement slab stopped rebounding, which might bring numerical errors (i.e., energy unbalance).

Numerical Result
The damaged situations for the multi-layer slab under the first drop weight impact for the experiment and numerical model are presented in Figure 17.It is found from the numerical simulation that the AC layer is penetrated through, and the drop weight head is impeded by the HSC layer due to the high compressive strength.Such findings from the numerical model are consistent with the observation in the physical test.It is observed that no severe and moderate damage happened in the ECC layer.Additionally, the integrity of both HSC and ECC layers is kept after the first impact.The damage situations for the multi-layer pavement slab under the second drop weight impact and the numerical model are given in Figure 18.The simulation results indicate that the HSC layer stopped the impactor; however, severe damage happens at the top surface and in the middle at the side of the HSC layer.Severe damage occurs at the rear surface of the HSC layer, as well.For the ECC layer, it is shown that the cracking occurs and propagates from the centre at the top surface.However, the damage area is smaller than that of the HSC layer.The cracks are also found in the middle at the side of the ECC layer, and similar to the experiment, the crack does not propagate through the thickness of the layer.Figure 19 compares the numerical results with the experimental results about the rear surface of ECC layer after the second impact for the multi-layer slab.It is found that the damaged contour and cracking pattern were similar.Such results further indicate the reasonable close agreement between the numerical analysis and the experimental results using the KCC model for AC material.The damage situations for the multi-layer pavement slab under the second drop weight impact and the numerical model are given in Figure 18.The simulation results indicate that the HSC layer stopped the impactor; however, severe damage happens at the top surface and in the middle at the side of the HSC layer.Severe damage occurs at the rear surface of the HSC layer, as well.For the ECC layer, it is shown that the cracking occurs and propagates from the centre at the top surface.However, the damage area is smaller than that of the HSC layer.The cracks are also found in the middle at the side of the ECC layer, and similar to the experiment, the crack does not propagate through the thickness of the layer.Figure 19 compares the numerical results with the experimental results about the rear surface of ECC layer after the second impact for the multi-layer slab.It is found that the damaged contour and cracking pattern were similar.Such results further indicate the reasonable close agreement between the numerical analysis and the experimental results using the KCC model for AC material.
middle at the side of the ECC layer, and similar to the experiment, the crack does not propagate through the thickness of the layer.Figure 19 compares the numerical results with the experimental results about the rear surface of ECC layer after the second impact for the multi-layer slab.It is found that the damaged contour and cracking pattern were similar.Such results further indicate the reasonable close agreement between the numerical analysis and the experimental results using the KCC model for AC material.Table 12 compares the numerical results with the experimental data for the multi-layer slab under two times of the drop weight impact.After the first impact, the settlement from the numerical results is close to that obtained from the experiment at P2.A big deviation is found in between the numerical results and experiment results on the settlement at P3.This might be attributed to the dislocation of the potentiometer at P3 upon impact, which could lead to a misreporting of the settlement at P3.In addition, the settlement at P1 is lower in the numerical results compared with the experimental results.This could be attributed to the erosion process in the numerical model or, more specifically, the energy released phenomena.At the second impact, since no erosion technology is employed for the HSC layer, the numerical results are more reasonable.In the impact simulation, the HSC layer suffered deformation instead of being penetrated, and hence, no erosion technique was adopted to trigger the energy released phenomena.Thus, the settlement profile of the pavement slab decreased with the increase of the radial distance from the centre, which indicated that the whole pavement slab bent during the impact.In general, the numerical models provide reasonable estimations on the settlements of the multi-layer slab under impacts, in particular under the second impact load, though discrepancies exist at certain regions.

Conclusions
In this paper, the formulations of strength surfaces and damage factor in the KCC model are discussed and evaluated.It is found that the KCC model would reproduce the real behaviour of concrete-like material under severe loadings, such as blast and impact loadings.Then, the key parameters controlling the behaviour of AC material under dynamic loadings are discussed, and the parameters of AC materials using the KCC model are determined by laboratory static and dynamic tests.The proposed damage factor in this paper results in smoother descending than the original damage factor and also a higher failure strain that could improve the simulation of the behaviour of AC material.Furthermore, the numerical modelling of multi-layer pavement, in which the AC material served as the surface layer, under blast and impact loading is conducted using the KCC model.These results are compared with that of the field blast test and the laboratory drop weight impact tests, and it is found that the KCC model for the AC material can reproduce the real material behaviour under blast and impact loading.Table 12 compares the numerical results with the experimental data for the multi-layer slab under two times of the drop weight impact.After the first impact, the settlement from the numerical results is close to that obtained from the experiment at P2.A big deviation is found in between the numerical results and experiment results on the settlement at P3.This might be attributed to the dislocation of the potentiometer at P3 upon impact, which could lead to a misreporting of the settlement at P3.In addition, the settlement at P1 is lower in the numerical results compared with the experimental results.This could be attributed to the erosion process in the numerical model or, more specifically, the energy released phenomena.At the second impact, since no erosion technology is employed for the HSC layer, the numerical results are more reasonable.In the impact simulation, the HSC layer suffered deformation instead of being penetrated, and hence, no erosion technique was adopted to trigger the energy released phenomena.Thus, the settlement profile of the pavement slab decreased with the increase of the radial distance from the centre, which indicated that the whole pavement slab bent during the impact.In general, the numerical models provide reasonable estimations on the settlements of the multi-layer slab under impacts, in particular under the second impact load, though discrepancies exist at certain regions.

Conclusions
In this paper, the formulations of strength surfaces and damage factor in the KCC model are discussed and evaluated.It is found that the KCC model would reproduce the real behaviour of concrete-like material under severe loadings, such as blast and impact loadings.Then, the key parameters controlling the behaviour of AC material under dynamic loadings are discussed, and the parameters of AC materials using the KCC model are determined by laboratory static and dynamic tests.The proposed damage factor in this paper results in smoother descending than the original damage factor and also a higher failure strain that could improve the simulation of the behaviour of AC material.Furthermore, the numerical modelling of multi-layer pavement, in which the AC material served as the surface layer, under blast and impact loading is conducted using the KCC model.These results are compared with that of the field blast test and the laboratory drop weight impact tests, and it is found that the KCC model for the AC material can reproduce the real material behaviour under blast and impact loading.

Figure 2 .
Figure 2. Determination of parameters in KCC from experimental data.

Figure 2 .
Figure 2. Determination of parameters in KCC from experimental data.

23 Figure 3 .
Figure 3. Validation of the failure surface from the scaling law.

Figure 3 .
Figure 3. Validation of the failure surface from the scaling law.

Figure 4 .
Figure 4. Damage factor used for AC material.

Figure 4 .
Figure 4. Damage factor used for AC material.

Figure 4 .
Figure 4. Damage factor used for AC material.

Figure 6 .
Figure 6.Tensile and compressive dynamic increase factor (DIF) curve used in the numerical model for asphalt concrete.

Figure 6 .
Figure 6.Tensile and compressive dynamic increase factor (DIF) curve used in the numerical model for asphalt concrete.

Figure 7 .
Figure 7. Cross-section of the multi-layer pavement slab.HSC, high strength concrete; ECC, engineering cementitious composite.

Figure 7 .
Figure 7. Cross-section of the multi-layer pavement slab.HSC, high strength concrete; ECC, engineering cementitious composite.

Figure 8 .
Figure 8. Plan view of the multi-layer pavement slab before the blast event.

Figure 8 .
Figure 8. Plan view of the multi-layer pavement slab before the blast event.

23 Figure 9 .
Figure 9. Layout of instrumentation for field blast test.TPC, total pressure cell.

Figure 9 .
Figure 9. Layout of instrumentation for field blast test.TPC, total pressure cell.

Figure 9 .
Layout of instrumentation for field blast test.TPC, total pressure cell.

Figure 10 .
Figure 10.Finite element model of the multi-layer pavement slab.

Figure 10 .
Figure 10.Finite element model of the multi-layer pavement slab.

Figure 11 .
Figure 11.Damage of the multi-layer pavement after blast loading: (a) damage pattern of the AC layer in field blast test; (b) damage pattern of the HSC layer in field blast test.

Figure 11 .
Figure 11.Damage of the multi-layer pavement after blast loading: (a) damage pattern of the AC layer in field blast test; (b) damage pattern of the HSC layer in field blast test.

Figure 12 .
Figure 12.Damage contour for the AC and HSC layer in the multi-layers pavement: (a) damage contour of the AC layer in the numerical simulation; (b) damage contour of the HSC layer in the numerical simulation.

Figure 12 .
Figure 12.Damage contour for the AC and HSC layer in the multi-layers pavement: (a) damage contour of the AC layer in the numerical simulation; (b) damage contour of the HSC layer in the numerical simulation.

Figure 13 .
Figure 13.Setup for the impact test.

Figure 13 .
Figure 13.Setup for the impact test.

Figure 16 .
Figure 16.Numerical model for a multi-layer pavement slab under drop weight impact (quarter model).

Figure 16 .
Figure 16.Numerical model for a multi-layer pavement slab under drop weight impact (quarter model).

Figure 17 .
Figure 17.Damage of the multi-layer slab after the first impact: (a) damage pattern of slab in the laboratory test; (b) damage pattern of slab in the numerical simulation.

Figure 17 .
Figure 17.Damage of the multi-layer slab after the first impact: (a) damage pattern of slab in the laboratory test; (b) damage pattern of slab in the numerical simulation.

Figure 18 .
Figure 18.Damage of multi-layer slab after the second impact: (a) damage pattern of slab in the laboratory test; (b) damage pattern of slab in the numerical simulation.

Figure 18 .
Figure 18.Damage of multi-layer slab after the second impact: (a) damage pattern of slab in the laboratory test; (b) damage pattern of slab in the numerical simulation.

Figure 19 .
Figure 19.ECC bottom after the second impact: (a) damage pattern of the bottom of ECC layer in the laboratory test; (b) damage pattern of the bottom of ECC layer in the numerical simulation.

Figure 19 .
Figure 19.ECC bottom after the second impact: (a) damage pattern of the bottom of ECC layer in the laboratory test; (b) damage pattern of the bottom of ECC layer in the numerical simulation.

Table 1 .
Parameters for asphalt concrete (AC) material with f c = 4.6 MPa.

Table 3 .
Basic properties of materials in the multi-layer pavement for the blast load.

Table 4 .
Parameters for Geogrid MG-100 using the plastic-kinematic model.

Table 5 .
Material properties of the soil mass.

Table 6 .
Parameters for the interface simulation.

Table 7 .
Vertical acceleration of the multi-layer pavement slab.

Table 7 .
Vertical acceleration of the multi-layer pavement slab.

Table 8 .
Peak reading for the total pressure cell.

Table 9 .
Basic properties of the materials in the multi-layer pavement in the drop weight test.

Table 9 .
Basic properties of the materials in the multi-layer pavement in the drop weight test.

Table 10 .
Material properties of the foundation soil.

Table 11 .
Properties of the drop weight head.

Table 12 .
Summary of the settlements of the multi-layer slab from both numerical and experiment results.

Table 12 .
Summary of the settlements of the multi-layer slab from both numerical and experiment results.