Investigation on Dynamic Propagation Characteristics of In-Plane Cracks in PVB Laminated Glass Plates

Polyvinyl butyral (PVB) laminated glass has been widely used as an important component of mechanical and construction materials. Cracks on PVB laminated glass are rich in impact information, which contribute to its impact resistance design. In this paper, a three-dimensional (3D) numerical simulation model describing PVB laminated glass under impact loading is firstly established and validated qualitatively and quantitatively compared with the corresponding experimental results recorded by the high-speed photography system. In the meantime, the extended finite element method (XFEM) is introduced to analyze the crack propagation mechanism of laminated glass based on dynamic stress intensity factors (DSIFs) and propagations of stress waves. Parametric studies are then carried out to investigate the influence of five critical parameters, that is, plate dimension, crack length, impact energy, glass properties, and PVB properties, on crack propagation characteristics of laminated glass. Results show that the interaction between crack tip and stress waves as well as the propagations of stress waves corresponds to the fluctuations of DSIFs at crack tip. Both the structure and material variables are proven to play a very important role in glass cracking DSIFs and thus govern the crack propagation behavior. Results may provide fundamental explanation to the basic crack propagation mechanism on radial cracks in PVB laminated glass under impact loading conditions, thus to instruct its impact design improvement.


Introduction
Due to its excellent transparency, formability, and impact resistance, polyvinyl butyral (PVB) laminated glass becomes great potential transparent protective material and has been widely used in both automotive and construction industry.The impact of foreign objects on laminated glass can generate microcracks, which will gradually extend and eventually lead to the failure of laminated glass.Regardless of whether fracture failure occurs, the fragility of glass layer and the good adhesion of PVB interlayer ensure the integrity of sandwiched structure and make the structure have certain energy absorption characteristics and toughness, which effectively improve the security of laminated structure.In order to utilize laminated glass as transparent protective material more effectively, its cracking propagation behavior under impact conditions needs to be comprehensively studied.However, dynamic failure characteristics of laminated glass are very complex; therefore studies on its dynamic cracking propagation characteristics become increasingly necessary.
In recent years, investigations on mechanical properties of PVB laminated glass have been thoroughly conducted, ranging from quasi-static mechanical properties [1][2][3][4] to dynamic impact response [4][5][6][7][8][9][10].In general, the impact test in accordance with its impact velocity can be divided into three categories, that is, low-speed drop-weight test [11][12][13], high-speed impact test with air gun [3,14,15], and explosion shock wave impact test [6,9,10,16], where low-speed drop-weight test is one of the most widely used testing methods to evaluate the material critical failure conditions and its fracture mechanism.Lee et al. [17] and Park and Chen [14], respectively, studied the static and dynamic crack propagation behavior in layered glass cross adhesive interlayers.Failure wave was introduced to account for their observations.Zhang et al. [18] and Hu et al. [19] experimentally and numerically studied failure conditions and crack morphologies in brittle plate with a thin polymer backing under various impact velocities and concluded the determining factors to predict its dynamic fracture patterns, that is, elastic properties, fracture energy, and accurate stress wave tracing model.Meanwhile, Vandenberghe et al. [20] investigated crack patterns in brittle plate upon controlled transverse impact using a global scaling law to predict the number of radial cracks.Also, the coexistence of radial and circumferential cracks was observed.As for a further step, Xu et al. [21] and Chen et al. [22] conducted out-ofplane dynamic loading test on PVB laminated glass based on drop-weight combined with high-speed photography system.Velocity and acceleration time history for both radial and circular cracks were calculated to understand the crack initiation and propagation mechanism.It is found that radial cracks appear earlier than circular ones with higher propagation speed.Furthermore, Chen et al. [23,24] investigated in-plane crack propagation behavior in PVB laminated glass plate via captured in-plane and interplane cracking process and found that the supported glass layer would always initiate before the loaded layer with the same final radial crack morphologies.The quantity effect of radial cracks on its propagation behavior was also thoroughly investigated.
Numerical simulation attracts more and more attentions to investigate the fracture mechanism of PVB laminated glass for its simplicity and efficiency [5,8,[25][26][27][28][29][30][31][32].Bennison et al. [33] studied the stress distribution and fracture sequence in PVB laminated glass upon biaxial flexure loading using finite element method (FEM) and found the relationship between crack initiation and maximum stress location.Meanwhile, Xu and Li [34] numerically investigated the crack evolutions on PVB laminated windshield plate subjected to pedestrian head.It was found that the plastic deformation, radial cracks, and circular cracks were mainly caused by shear stress, compressive stress, and tensile stress, respectively.Furthermore, a series of newly developed numerical simulation methods are increasingly applied to researches on fracture characteristics of PVB laminated glass for their unique advantages in cracking simulation.Xu and Zang [35], Zang et al. [36], and Gao and Zang [37] investigated the impact fracture mechanism of PVB laminated glass on the basis of discrete element method (DEM).Chen et al. [38] numerically investigated the cracking behavior of PVB laminated glass beam under lowvelocity impact in the framework of cohesive zone modeling with the glass-ply cracking modeled by extrinsic cohesive model.The propagations of stress waves during cracking process were also primarily analyzed.Based on extended finite element method (XFEM), Xu et al. [39] and Xu et al. [40] established PVB laminated glass plate model subjected to impact loading validated by both radial crack velocity and crack morphology and deeply investigated its radial and circular crack propagation mechanism.
However, as a key parameter to indicate the cracking driving force and interpret the crack propagation behavior in cracking mechanism analysis, dynamic stress intensity factor (DSIF) cannot be easily obtained via experimental methods.The propagations of stress waves during cracking can also have significant effect on cracking response of PVB laminated glass under impact loading.Besides, seldom numerical studies are available in exploring the cracking behavior of PVB laminated glass during out-plane impact loading from the perspective of DSIFs and stress waves.It is therefore necessary to investigate the dynamic crack propagation behavior of PVB laminated glass on the basis of DSIFs and stress waves.
Hence, the numerical model of PVB laminated glass plate with a layered structure is established in this study.The numerical model in simulations with propagation is then validated via the comparison between simulation and experimental results of radial crack propagation process, which are recorded by the high-speed photography system.Further, the time history of DSIFs is calculated in simulations with fixed crack, and the interaction between stress waves and crack tip is also introduced to interpret the relationship between crack propagation behavior and its driving force.Finally, parametric studies on five critical parameters, that is, plate dimension, crack length, impact energy, glass properties, and PVB properties, are conducted based on XFEM, which reveals the laws of crack propagation behavior on PVB laminated glass plates.

Experimental and Numerical Details
2.1.Experiment Setup.The drop-weight device combined with high-speed photography experimental system [21][22][23] is sketched in Figure 1.The 150 mm × 200 mm rectangular specimen of laminated glass is used in the experiment, which consists of a 0.76 mm PVB interlayer sandwiched by two 2 mm soda-lime glass sheets and is fabricated with the same processing technology of industrial automotive windshield.The drop-weight device is used to apply vertical load via free fall with different height.Meanwhile, the impact direction convertor with hemispherical end transfers the vertical load to horizontal load at the specimen center, that is, a load direction perpendicular to specimen plane.The high-speed photography system is adopted to record the real time crack propagation process on glass layer of specimen.The radial crack morphology is recorded on the film and crack propagation speed is calculated [22].

Computational Model.
In this study, the three-dimensional (3D) finite element (FE) model consisting of an impactor, a prefabricated crack, and a PVB laminated glass plate is shown in Figure 2. One-quarter of the laminated plate and impactor are modeled to take advantage of the structure symmetry.The cylinder impactor has a mass of 2 kg via equivalent density with the same mass of experimental impactor.The 2 mm thick crack with the length of a (shown in Figure 2) is placed at the center of laminated plate in the glass layer of supported side and the angle between crack and structure symmetry line is 45 degrees.In order to eliminate the influence of structural asymmetry, square laminated plate model with dimensions of 200 mm × 200 mm on plane direction is introduced to study its crack propagation characteristics.The laminated plate model is meshed with 2 mm × 2 mm (in plane direction) 3D eightnode reduced integration solid (C3D8R) elements.Since stress concentration exists and the stress gradients are very large at crack tip, the crack tip area needs to be meshed with sufficient fine elements to achieve stress and strain accurately.For this purposes, four mesh sizes at crack tip area in plane direction, that is, 0.35 mm × 0.35 mm, 0.7 mm × 0.7 mm, 1 mm × 1 mm, and 2 mm × 2 mm, are applied to achieve a stable convergence of the computational model.The results are shown in Figure 3, where the mode I DSIFs of node 1 (node 1 to node 9 at crack tip are defined in Figure 2) with fixed crack length  = 0.02 m under the impact speed of 5 m/s are illustrated.As can be seen, all simulations show the same oscillation frequencies and the mesh refinement has only effect on oscillation amplitudes.Therefore, the configuration of 0.35 mm × 0.35 mm can be used for further simulations sufficiently.In longitudinal direction, the laminated structure of specimen is modeled with 20 element layers, wherein the outside 8 element layers represent glass layers and the inner 4 element layers stand for the PVB interlayer (shown in Figure 2).Therefore, each layer of the laminated structure has more than three integration points to guarantee the accuracy of calculation and analysis in longitudinal direction.
Clamped boundary conditions of the laminated plate model and frictionless surface-to-surface contact between impactor and specimen are adopted to simulate the real experimental conditions.As typical viscoelastic material, PVB interlayer is simulated with a linear elastic model with a shear modulus  =  0 and bulk modulus , combined with a viscoelastic model.The elastic modulus  PVB and Poisson's ratio ] are described as [42] where  and  0 are bulk modulus and short time shear modulus of PVB, respectively.The generalized Maxwell model using Prony series with one term is introduced to simulate the viscoelasticity of PVB [42]: where () is the relaxation modulus which is the function of time ,  ∞ is long time shear modulus, and  is decay factor.As typical brittle material, glass layers are defined by linear elastic model.Reference parameters of glass layer, PVB interlayer, and impactor are listed in Table 1.
The computational model with the same dimensions of experimental system is set up using a cylinder impactor and a laminated plate model [39] to verify the effectiveness of the  abovementioned model.The experimental system model is then verified in simulations with propagation via qualitative method (radial crack morphology) and quantitative method (crack propagation velocity) shown in Figures 4 and 5, respectively.As can be seen from Figure 4, a total of 22 radial cracks exist in experimental images and most of them have propagated to the border of the specimen shown in Figure 4(a).In comparison, 20 radial cracks in all have propagated close to the border in simulation images shown in Figure 4(b).Thus, both the crack propagation process and the morphology of cracks in simulation agree well with the experimental results.Figure 5 shows the comparison of crack velocity time history curves between experiment and simulation results.As can be observed, the cracks propagation time history in simulation is recorded from the increment of the original crack, for initial flaws are preembedded in the model.Another difference between the experimental and simulation results is that the gradient of the descending stage of the velocity curves in simulation is smaller than that in experiment.The responsible reason should be the rubber pad used between the mental clamp and glass sheet to avoid possible damage during screwing in experiment.The rebound velocity in simulation is higher in simulation where no cushion pad is provided.So the peak value as well as the time history of the simulation results is in good agreement with the experimental results.The computational results are in good agreement with the corresponding experimental results from both qualitative perspective and quantitative perspective.Therefore, a reliable FE model with fairly high quality is established and can be adopted for further studies on crack behavior of PVB laminated glass.

Fundamentals of XFEM.
The extended finite element method XFEM, which is a numerical technique based on generalized finite element method and the partition of unity method, extends the piecewise polynomial function space of conventional finite element methods with extra enrichment functions.XFEM does not need to be updated to track the crack path and can obtain accurate solution with coarse mesh [43,44].The displacement interpolation of XFEM is described as below [43]: where  Γ is the number of nodes in elements cut by crack;  Λ is the number of nodes in elements with crack tip;   (x) are the conventional shape functions; u  are the nodal DOFs for conventional shape functions; a  and b   are the nodal DOFs associated with the Heaviside function (x), with the value above the crack and −1 below the crack.The crack tip enrichment functions   (x) are expressed as [43] where x is an integration point and (r, ) denote the coordinate values from a polar coordinate system located at the crack tip.Details of XFEM can be found in [43][44][45].

Calculation of Dynamic Stress Intensity
where Dyn ΙΙΙ ]  is SIF vector and  is the prelogarithmic energy factor matrix, which is diagonal for homogeneous and isotropic materials.Therefore the equation can be simplified to where  * =  for plane stress and  * = /(1 − ] 2 ) for plane strain.
The interaction integral method [46] is introduced to calculate DSIFs for cracks under mixed-mode loading, instead of according to a known J-integral.The J-integral can be described as With regard to an auxiliary pure mode I crack tip field, the DSIF  Dyn Ι can be defined as Then overlay the auxiliary field on the actual field yields Advances in Materials Science and Engineering the interaction integral for mode I field can be expressed as For mode II and mode III, similar results can be obtained: can be neglected within the scope of our discussion impact velocity [48,49].Therefore, the rate sensitivity of  Dyn ΙC for glass layer can be ignored in the analysis of laminated glass fracture mechanism, that is,  Dyn ΙC ≈  ΙC ≈ 0.7 MP√m.

Cracking Propagation Direction.
The maximum tangential stress criterion is adopted to compute the cracking propagation direction, where  Dyn ΙΙΙ is not taken into consideration in these criteria.In the case of   / = 0 (or   = 0, where  and  are normal stress and tangential stress, resp., and r and  are polar coordinates at the crack tip in a plane perpendicular to the crack plane), the crack propagation angle relative to the crack plane can be obtained [46]: where θ = 0 indicates the crack propagation in a straight direction.θ > 0 corresponds to  Dyn ΙΙ < 0, while θ < 0 stands for  Dyn ΙΙ > 0.

Simulations with Fixed
Crack.According to our previous experimental studies [23], there are three types of cracks in PVB laminated plate under impact loading, that is, radial cracks in supported side, radial cracks in loaded side, and circumferential cracks in supported side, where supported side and loaded side are defined in Figure 2. Besides, radial cracks appear before circumferential cracks, while radial cracks in supported side are generated earlier than radial cracks in loaded side and radial cracks in both sides overlap completely.Therefore, the study on radial cracking propagation behavior in supported side is fairly important for investigating the fracture properties of PVB laminated glass.
Then, a series of numerical simulations on dynamic response of laminated glass with a fixed crack are conducted under different impact conditions using commercial finite element software ABAQUS.In our simulations, the crack is specified to be stationary and the material in crack tip should be in pristine condition to predict the crack propagation.Five structure and material variables, that is, plate dimension, crack length, impact energy, glass properties, and PVB properties, are introduced to conduct parametric studies in Section 4, and the detailed variable settings in parametric study are listed in Table 2 with other parameters kept at their reference values shown in Table 1.

Analysis of the Time History of DSIFs
In the condition of dynamic fracture, DSIFs are functions of time and the in-depth analysis of the time history of DSIFs helps to fully understand the crack propagation characteristics.Figure 6(a) shows the velocity and acceleration time history of impactor with the initial velocity of 5 m/s, and the corresponding mode I, II, and III DSIFs of node 1 at crack tip are shown in Figure 6(b).It can be seen that both the acceleration curve and DSIF curves exhibit regular fluctuations in the same law before 3 × 10 −3 s, while the acceleration curve almost has no fluctuation and DSIF curves still have slight regular fluctuations after 3 × 10 −3 s.Since the impactor has been away from the plate after 3 × 10 −3 s, the existence of stress waves leads to the fluctuations of DSIF curves in the entire time course and the crest appears when crack tip interacts with stress waves, indicating that the interaction between crack tip and stress waves accounts for the DSIFs' changes.Therefore, the time of interaction between the crack tip and stress waves corresponds to the fluctuation crests of DSIFs.
Since | 6), the crack tip field is mixed-mode crack tip filed which mode I is in dominance [21,22], that is, the mixedmode DSIF  Dyn ≈  Dyn Ι , and mode I DSIF can be introduced to represent mixed-mode DSIF. Figure 7 shows the mode I DSIF time history of node 1 to node 8 at crack tip under the impact velocity of 5 m/s (Figure 7(b) shows the data before 5 × 10 −5 s).As can be seen, mode I DSIFs of node 1 to node 8 at crack tip begin to rise rapidly after a time lag which is defined as lag time, then fluctuate in the same regularity under the impact stress waves, and finally fluctuate slightly with the crest decreasing under the residual stress waves, so the entire process can be divided into three stages: lag stage, interaction stage, and residual stage.The propagation process of stress waves is introduced to explain this phenomenon (shown in Figure 8): the impact generates the compressive stress waves, which propagate from the loaded side to the supported side in the thickness direction.Then the compressive stress waves reflect from the supported side (free surface) into tension stress waves, which consist of shear wave, longitudinal wave, and Rayleigh wave.The tension stress waves further propagate to the boundaries and reflect to the crack tip (shear wave and longitudinal wave through glass layer, while Rayleigh wave through free surface of supported side), which results in the crests of DSIFs.Besides, mode I DSIFs of node 1 to node 8 decrease one by one for the weakening of reflected tension stress waves during propagation, while there is almost no difference in peak times of these nodes for the extremely high propagation velocity of stress waves.Therefore, the round-trip propagation of tension stress waves leads to the fluctuations of DSIFs.
The velocities of involved stress waves are expressed as follows [50]: where   ,   , and   stand for the velocities of longitudinal wave, shear wave, and Rayleigh wave, respectively, besides   = max{  ,   ,   }; , ], and  represent Young's modulus, Poisson's ratio, and density, respectively.The layer interfaces can cause large amounts of stress wave reflection and further dramatically delay the stress wave propagation.Similar delays have also been discovered in other experimental researches [14].Therefore, the lag time of DSIF, which is the shortest propagation time of compressive stress waves from loaded side to supported side, can be expressed as where  Glass ,  PVB , and  Interface represent propagation time of longitudinal wave through glass layers, PVB interlayer, and layer interfaces, respectively.Analogously, the interval of propagation time between node 1 and node 8 should be Furthermore, the time when DSIF increases to dynamic fracture toughness with the stress waves accumulation, that is,  Dyn ≈  Dyn ΙC ≈  ΙC ≈ 0.7 MP√m, is considered as the simulated propagation time.The simulated crack propagation times of node 1 to node 8 and the corresponding linear fitting results are illustrated in Figure 9; then the simulated propagation time interval between node 1 and node 8 can be calculated as 4.55 × 10 −7 s which is in good agreement with theoretical calculation.
According to the stress wave propagation theory, crack propagation firstly occurs at node 1, so parametric studies on crack propagation characteristics focus on mode I DSIF at node 1 in Section 4. in plane direction are considered to study the dimension effect on cracking propagation characteristics of glass layer.Figure 10 gives the effect of plate dimension on mode I DSIF of node 1.As can be seen, the interaction stage is getting longer and the peak DSIF would be smaller with the increase of plate size under the same impact energy due to the decrease of structure stiffness.The time that the stress wave propagation needed from the crack tip to the plate borders and return to the crack tip becomes longer due to the increase of plate edge length; therefore the time interval of DSIF peaks, that is, fluctuation cycle, increases with the increase of plate size which is fully consistent with the simulation results as shown in Figure 10.11(a)).Since the pre-crack length increases, the stress wave propagation path from impact point to crack tip increases, which leads to the increase of lag time, and the propagation path from crack tip to proximal border decreases, while it increases to distal border, which results in the changes of fluctuation frequency.In addition, the DSIF decreases with the increase of pre-crack length, due to the increase of propagation path from impact point to crack tip which leads to the increased stress wave weakening.

Effects of Impact Energy.
Two aspects, that is, impact velocity and impact mass, are introduced to investigate the influence of impact energy on cracking propagation characteristics of glass layer.The mode I DSIFs of node 1 with different impact velocities and impact masses are given in Figures 12(a) and 12(b), respectively.As shown in Figure 12, both the increase of impact velocity and impact mass lead to the increase of peak mode I DSIF at all stages and almost have no effect on fluctuation frequency of DSIF curve.Besides, the interaction time between impactor and plate increases with the increase of impact mass but is not affected by impact velocity.Since the acceleration of impactor is mainly affected by impact velocity, the rebound time of impactor increases with the impact mass and almost has no change with impact speed.In summary, the increase of impact energy results in the increase of stress wave accumulation, that is, peak DSIF, while it has little influence on stress wave propagation velocity, that is, fluctuation frequency.

Effects of Glass Properties.
In static case, stress intensity factors have nothing to do with material's Young's modulus and density.However, under impact loading conditions, the stress wave accumulation at crack tip leads to the crack propagation; thus stress wave propagation speed plays a critical role in crack propagation characteristics.As can be seen from ( 14), longitudinal wave, shear wave, and Rayleigh wave velocities are all functions of Young's modulus and density.Besides, stress wave velocities are proportional to Young's modulus and inversely proportional to density.
Figures 13(a) and 13(b) show the effects of glass Young's modulus  Glass and density  Glass on mode I DSIF at node 1.As can be observed, larger Young's modulus corresponds to the shorter interaction stage and larger peak DSIF owing to the increase of structure stiffness, while density changes have little effect on the length of interaction stage and peak DSIF with structure stiffness unchanged.Besides, the lag time increases and the peak times are postponed, for the decrease of Young's modulus or the increase of density leads to the decrease of stress wave velocities.Furthermore, the stress wave propagation time becomes longer with the decrease of Young's modulus or the increase of density; therefore the fluctuation frequency increases with the increase of Young's modulus or the decrease of density, which agrees well with the simulation results shown in Figure 13.

Effects of PVB Properties.
As longitudinal boundaries of cracks in glass layers, PVB interlayer has an important influence on DSIFs of glass layer.Figures 14(a) and 14(b) show the effects of PVB Young's modulus  PVB and density  PVB on mode I DSIF curves at node 1.As can be seen, the interaction stage and peak DSIF decrease with the increase of PVB Young's modulus due to the increase of structure stiffness, while density changes have little effect on the length of interaction stage and peak DSIF with structure stiffness unchanged.Moreover, the fluctuation frequency increases with the increase PVB Young's modulus or the decrease of PVB density.According to (14), the stress wave propagation velocities in glass layer are far higher than that in PVB interlayer; thus the existence of PVB interlayer has an inhibiting effect on the stress wave propagation in glass layer, which leads to the increase of stress wave propagation time in glass layer.
The viscosity  in the one-term generalized Maxwell model can be expressed as As seen, shear modulus ratio  ∞ / 0 and delay factor  are two main viscosity factors, which are inversely proportional to viscosity.Figures 14(c  of viscosity, which explains the viscosity effect.However, the stress wave propagation time is far less than the PVB relaxation time and only a small amount of relaxation occurs during the impact.Therefore, the influence of delay factor on DSIF is less obvious than shear modulus ratio.

Concluding Remarks
In this paper, the crack propagation characteristics of PVB laminated glass under low-speed impact are investigated in the framework of XFEM.A PVB laminated glass 3D model is proposed, which has been qualitatively and quantitatively validated via comparing with the experimental results captured by high-speed photography system.The verified numerical model with fixed crack is then used to systematically investigate the crack propagation behavior.The time history of brittle DSIFs as well as the propagations of stress waves is derived to interpret the crack propagation mechanism.
Further, five critical structure and material variables, that is, plate dimension, crack length, impact energy, glass properties, and PVB properties, are introduced to conduct parametric studies.It is found that the interaction between the crack tip and stress waves and the round-trip propagations of tension stress waves lead to the fluctuations of DSIFs, which corresponds to the influence of plate dimension and crack length.The impact energy and glass properties result in the changes of DSIF peak values and fluctuation frequencies, respectively.As longitudinal boundaries of cracks in glass

Figure 1 :
Figure 1: The schematic diagram of experimental system combing drop-weight test device with high-speed photography system.

Figure 2 :K 2 )Figure 3 :
Figure 2: The overview of PVB laminated glass computational model subjected to impact loading.

Figure 4 :
Figure 4: Qualitative verification: comparison of radial crack morphologies between simulation results and experimental results.
Factors.In linear fracture mechanics, for example, soda-lime glass, the dynamic stress intensity factors ( Dyn Ι ,  Dyn ΙΙ ,  Dyn ΙΙΙ ) are widely used to characterize the stress and displacement fields near the crack tip, where  Dyn Ι ,  Dyn ΙΙ , and  Dyn ΙΙΙ are mode I, II, and III DSIFs.DSIFs relevant to energy release rate (Jintegral) can be described as

Figure 5 :
Figure 5: Quantitative verification: comparison of crack propagation velocities between simulation results and experimental results.

K 2 )
(MPa•m 1/(b) Mode I, II, and III DSIFs of node 1 at crack tip

Figure 6 :
Figure 6: Simulation results with the impact velocity of 5 m/s and crack length of 0.02 m.

2 ) 2 )Figure 7 :Figure 8 :
Figure 7: The mode I DSIF time history of node 1 to node 8 at crack tip with the impact velocity of 5 m/s and crack length of 0.02 m.

4. 1 .Figure 9 :
Figure 9: The simulated crack propagation times of node 1 to node 8 and the corresponding linear fitting curve.

2 )Figure 10 :
Figure 10: Mode I DSIF curves at node 1 of different plate dimensions with the impact velocity of 5 m/s and crack length of 0.02 m.

Figure 11
shows the mode I DSIF curves at node 1 with different fixed crack lengths (0.01 m, 0.02 m, and 0.04 m) under the impact velocity of 5 m/s.As can be seen, with the increase of fixed crack length, the lag time increases proportionally (Figure11(b)) and the peak times are postponed with fluctuation frequency (reciprocal of fluctuation cycle) changing slightly (Figure

2 ) 2 )Figure 11 :
Figure 11: Mode I DSIF curves at node 1 with different fixed crack lengths under the impact velocity of 5 m/s.

Figure 12 :
Figure 12: Mode I DSIF curves at node 1 with crack length of 0.02 m under different impact energy.

Figure 13 :Figure 14 :
Figure 13: Mode I DSIF at node 1 of different glass properties with the impact velocity of 5 m/s and crack length of 0.02 m.

Table 1 :
Basic parameters for numerical model. IC ,  IIC ,  IIIC  IC = 10 J/m 2 ,  IIC =  IIIC = 50 J/m 2 Shear modulus:  0 ,  ∞  0 = 0.33 GPa,  ∞ = 0.69 × 10 −3 GPa Decay factor  12.6 s −1 ) 2.5.Crack Propagation Criteria.In static case, the open mode crack propagation criteria are generally expressed as  Ι ≥  ΙC , where  Ι and  ΙC (for glass:  ΙC ≈ 0.7 MP√m) are static stress intensity factor and static fracture toughness, respectively.Similarly, the dynamic open mode crack propagation criteria are defined as  Dyn ΙC is determined through dynamic structural analysis.Since the cracks are generated in linear elastic glass layer, the influence of loading rate on