Analysis of multiple impact tests ’ damage to three - dimensional four - directional braided composites

: This article was designed with a plurality of impact tests of three - dimensional four - directional braided composites, and the impact response of specimen impacted by a circular punch was studied. Ultrasonic C - scanning was used to detect the internal damage area to study the damage propagation process under multiple impact loads. The ﬁ nite - element software ABAQUS was used to model the meso - structure of three - dimensional four - directional braided composites. Based on material characteristics, the three - dimensional Hashin damage criterion was used for the ﬁ ber bundle, and the maximum stress criterion was used for the matrix to judge the material damage. Combined with test and simulation results, the failure mode and damage evolution process of the specimen under multiple impact loads were studied. The results showed that the impact resistance of the three - dimen sional woven composite material is a ﬀ ected by the braided angle. The larger the braided angle of the specimen, the better the impact resistance. The damaged area of the large braided angle material expanded to the periphery, and the damaged area of the small braided angle material was primarily developed in the longitudinal direction. The failure modes of the specimen during the impact process were primarily a longitudinal tensile failure of ﬁ ber bun dles, transverse tensile failure and transverse compression failure of ﬁ ber bundles and matrix.


Introduction
The three-dimensional braided composite material first uses three-dimensional braiding technology to weave fiber bundles into a preform. It then uses the preform as a reinforcement to impregnate and solidify the composite material. Due to the unique spatial structure of 3D woven composites, it has better mechanical properties than traditional laminated composites, such as high specific strength, high specific stiffness, good impact toughness [1][2][3][4][5], considerable delamination resistance, and integral formability. At present, three-dimensional braided composite materials have been widely used in aerospace, automobile manufacturing, sports products, and other fields [6][7][8][9][10]. The service environment of braided composites is relatively special, and it will frequently encounter friction, knock and impact with objects. Due to the particularity of the braided composite structure, although the structure surface is almost completely without damage, different types of damage will appear more or less inside the material, and it is difficult to detect. With the extension of service time, internal damage will continue to accumulate and gradually expand, and the degree of damage will gradually deepen. In order to ensure the reliability and safety range of braided composites in service, judge the best service life and damage degree within the bearing range of braided materials, and maximize the utilization rate of materials, it is necessary to study the evolution of multiple impact damage of braided composites.
Gao et al. [6] performed shock compression tests on composite materials through separate Hopkinson compression bars and recorded them with high-speed cameras. At the same time, they conducted a meso-level mechanical simulation model to analyze the specific performance of the material's failure process and failure form in each shock wave. The results have shown that the plastic deformation primarily occurs in the first shock wave during the multiple shocks, and the stiffness deduction is mainly in the second shock wave. The failure is related to the damage to the yarn; at the same time, the specimen with a large braided angle exhibits ductile fracture, and the final failure result is related to the injury of the interface and the damage to the matrix.
Shi et al. [11] compared the fracture behavior of three-dimensional four-directional braided composites under impact loading. The progressive damage and fracture behaviors, such as crack opening displacement and crack propagation, were characterized by a drop-weight impact testing machine combined with high-speed photography. Deformation profiles and strain fields were evaluated using two-dimensional digital image correlation methods. A finite-element model was established to reveal the process of impact damage, crack extension, stress distribution, and energy absorption from the microstructure level. The results revealed that the braided preform is the primary energy absorption site. Yang et al. [9] conducted an experimental study on the deformation and failure mechanism of three-dimensional braided carbon fiber-reinforced composites under impact load. The damage initiation and propagation, full-field displacement, and strain field of 3D woven carbon fiber composites under ballistic impact were simultaneously studied using a dual-camera system. In addition, the energy absorption mechanism of 3D woven carbon fiber composites was also studied. Typical deflection curves and damage patterns were obtained. The results have shown that the absorbed energy has an approximately linear relationship with the incident velocity.
Wu et al. [12] tested the in-plane dynamic compression properties of tubular 3D four-directional and 3D five-way braided composites. It has been found that stress-strain curves, compressive deformation, and failure modes are very sensitive to strain rate. Zhou et al. [13] studied the transverse shock compression properties of three-dimensional tubular braided composites with different layers. The peak load and energy absorption increased with the number of braided layers; the failure modes of different braided layers were different. Pan et al. [14] studied the progressive damage failure process of three-dimensional textile structural composites under high-speed impact loading. Li et al. [15] used computed tomography technology combined with a finiteelement analysis method to reveal the impact shear fracture mechanism of 3D braided composites. Wan et al. [16] established a multi-scale finite-element model to analyze the dynamic mechanical properties of 3D braided hybrid (carbon fiber/aramid/epoxy) composites. Wehrkamp-Richter et al. [17] studied the damage and expansion mechanism of triaxial braided composites under multiaxial tensile stress. Kubota et al. [18] studied the thermal-mechanical coupling properties of composite materials. Heating caused strain mismatch between the fiber and the matrix to create tensile and compressive stresses. Less porous material, tensile stress leads to transverse cracks. Shi et al. [19] studied the matrix cracks and delamination phenomena of composite orthogonal laminates by numerical simulation.
Miravete et al. [20] used a representative volumetric unit cell model and finite-element method to study the progressive damage of 3D braided composites. Zhang et al. [21] developed an element-based finite-element model to evaluate the mechanical properties of 3D braided composites under different biaxial tensile loads. To predict the initial damage and evolution of the braided yarn and matrix in a cell, an anisotropic failure model was proposed based on the Murakami failure criterion. The damage evolution of 3D braided composites under typical biaxial tensile loading was simulated, and the damage mechanism was revealed during the simulation. Zako et al. [22] established an anisotropic progressive damage model based on damage theory and used two damage factors to reduce the stiffness of fiber bundles. The initial failure state was 0, and the complete failure state was 1.
With the deepening and expansion of the research on 3D woven composite materials, the research direction has gradually extended from the basic properties of the material to the overall performance of the structure. Still, in terms of dynamics, shock loads in the form of shock waves are mainly used, and the research on the lowspeed impact is still to be developed.

Materials and equipment
The yarn of the three-dimensional four-directional prefab is made of Toray T700-SC-12000-50B carbon fiber, the single-strand carbon fiber bundle has a fineness of 800tex, and the fiber volume content is 60%. The length of the flower pitch varies according to the braiding angle as shown in Figure 1(c) [23]. The surface weaving angle is the angle between the plane of the yarn on the surface of the preform and the growth direction of the flower pitch, as shown in Figure 1(a) [24], and the space angle between the main yarn in the preform and the growth direction of the flower pitch is shown in Figure 1(b). Also, the matrix uses TED-86 epoxy resin; the parameters of epoxy resin and fiber bundle are given in Table 1. The size of the specimen is 120 mm × 80 mm × 5 mm, and the error range is within ±0.5 mm. Some of the material parameters are given in Table 2.

Experimental methods
The test equipment mainly used in the impact test includes a model Instron-9250HV drop weight testing machine and positioning fixture and a model Sonoscan-D500 ultrasonic C scanner. The test equipment is shown in Figure 2.
The impact test of the three-dimensional four-directional braided composite material was carried out on a drop weight testing machine; the drop hammer weighed 7.76 kg. The specimen before testing is shown in Figure 3. The energy of 30 J required in this experiment can be satisfied by lifting the height of the drop hammer and then using its gravitational potential energy. It can be converted into impact energy.
The positioning fixture needs to be specially customized according to the size of the test piece. The positioning fixture and the pneumatic fixture can fix the test piece in the upper and lower parts, preventing the uncontrollable second impact once the first impact occurs.
In this impact test, a circular punch and an impact energy of 30 J were used to impact the braided board with a braiding angle of 15 and 30°for multiple times. When the specimen was penetrated, and the impact test was no longer performed. Three specimens under the same conditions were taken for testing, and typical specimens were taken for analysis.   In order to detect the damaged area inside the specimen by ultrasonic C-scanning, it is necessary to understand its working principle first. Point scanning is a process in which a beam of ultrasonic waves is focused on one point through a lens in the probe, and the energy is highly concentrated to find the subtle defects in the specimen. The specimen can be placed within the focus scanning range during the detection, and the effective focal length range of the detection lens can be changed to 5-50 mm. Ultrasonic C-scanning is the scanning of the plane perpendicular to the direction of ultrasonic emission. The synthesized image of the C-scan represents the superposition area of the scanned specimen. This device can draw the defect damage range, and the shape of the defect damage area can be simply and directly seen through imaging, but the normal depth of the plane where the defect is located cannot be given. Its imaging method is to scan the scanning area of points per unit area, and finally generate scanning images according to the different waveforms of all points in the scanning area. During scanning, the impacted specimens and the ultrasonic probe need to be placed in deionized water to scan the area containing the impact damage, the scanning range is set to 80 mm × 80 mm, and the shape and size of the damaged area inside the specimen are recorded to provide an experimental basis for analyzing the damage mode.

Results and discussion
The load-time and peak load-time curves are shown in Figure 4. Comparing the curve changes in Figure 4(a) and (b), it can be seen that the peaks of the load-time curves of the first two shocks fluctuate significantly, indicating that the two specimens have the initial damage that mainly occurred in the first two impact processes. With an increase in the number of impacts, the initial damage expanded and the damage degree continued to increase, and the stiffness of the specimen decreased more significantly with the number of impacts. The peak contact force shows a different downward trend with an increase in the number of impacts.
As shown in Figure 4(c), the peak load decreased slightly for both specimens when the number of impacts was 5-20, and the peak load decreased gradually after 20 impacts. After 20 impacts, the peak load drop gradually became larger, the peak load drop was greater for the 15°specimen than for the 30°specimen, and the impact resistance of both specimens decreased. The specimens with 15°and 30°braiding angles were penetrated after the 28th and 44th impacts, respectively, at  which time the specimens completely lost their impact resistance. Through the comprehensive comparison of the two braided angle specimens, it can be seen that the impact response of the two specimens is not significantly different when the number of impacts is small. When the number of impacts is large, the specimen with a braid angle of 30°showed good impact resistance. The specimen with a braiding angle of 15°was more sensitive to the cumulative damage of impact.
The surface damage of the test piece is shown in Figures 5 and 6. After 1 impact, a circular crater appeared on the surface of the specimen with a braiding angle of 15°, and the resin peeled off the back of the specimen; after 5 impacts, cracks appeared on the surface of the specimen, matrix cracking and fiber pull-out appeared on the back, and the cracks extended along the braiding direction; after 20 impacts, the crater at the impact location and the back of the specimen showed obvious fiber bundle fracture, and the cracking length of the back of the specimen extended to 60 mm and then basically did not change; the number of impacts continued to increase to 28, and the fiber bundle at the impact location gradually fractured until the specimen was penetrated. For the specimen with a braiding angle of 30°, the damage extension requires more impact energy accumulation. The cracking direction of the matrix on the back of the specimen is more likely to extend at the impact position, and the shape is an irregular circle. This is because the properties of the woven material are greatly affected by the preform. The larger the braiding angle, the better the   transverse mechanical properties of the composite material, so the directionality of damage expansion is not prominent. Figures 7 and 8 show the results of ultrasonic C-scanning of the specimens, which show that the internal damage of the specimens with different braiding angles is different. The internal damaged area with a braiding angle of 15°mainly expands in the longitudinal direction, and its lateral expansion is not apparent. The internal damaged area of the samples with a braiding angle of 30°was in the shape of a lung lobe. The damaged area expanded longitudinally and laterally, and the damaged area was more significant than that of the specimen with a 15°angle. The reason is that the space angle between the fiber bundle in the prefabricated specimen with different braiding angles and the longitudinal angle is different, so the stress of the fiber bundle is different. Also, the transverse width of the fiber bundles affected by the impact load will be different. The fiber bundle axis of the 30°specimen has a larger angle with its longitudinal direction, so the damaged area in the specimen extends more obviously along the transverse direction than that of the 15°specimen.

Numerical simulations 3.1 Finite simulation model
The internal damage to the braided composite material is caused by the impact. Due to the unique structure of the braided preform, it is difficult to detect and accurately locate the internal damage of the material by experimental methods. Therefore, in this section, the 15°specimen is used as a typical example to simulate the multiple spherical impacts through ABAQUS CAE commercial software. Study the failure modes of materials during impact. First, the 3D four-directional braided composite material is modeled by three-dimensional software UG and then imported into ABAQUS for analysis to obtain the single impact simulation result of the composite plate. Then, by controlling the speed and displacement of the punch, a second impact is repeated at a position, and so on, the specimen is impacted multiple times. The fiber is regarded as a transversely isotropic material. The damage forms of the yarn on the microscopic scale include  longitudinal stretching, longitudinal compression, transverse stretching, and transverse compression. The matrix is an isotropic material, and the main failure modes are a tensile failure and compression fracture. According to the characteristics of the material components, the fiber bundle adopts the three-dimensional Hashin initial damage criterion, and the matrix resin adopts the maximum stress criterion to judge the initial damage of the matrix. Then the damaged material is subjected to stiffness degradation to simulate the damage process of the material. Combined with experiments to verify the effectiveness of the simulation method and results, the damage failure mode and damage evolution process of the specimen are analyzed.

Geometric model
UG software has powerful geometric modeling functions, which can basically meet all modeling requirements. In addition, the models in UG can be easily imported into other software for analysis and processing. The molding process of the three-dimensional four-directional preform is relatively complicated. ABAQUS can adapt to the modeling of simple structures, and its preprocessing function is slightly lacking for geometric models of complex shapes. Therefore, this part uses UG software to model the yarn geometry. Then, it is imported into ABAQUS in the form of a solid geometric model for some pre-and postprocessing. Since the internal structure of the threedimensional four-directional braided composite material is very complex, it is complicated to restore the braided structure and yarn direction. This part mainly explores the failure mode of the fiber bundle yarn, so some idealized assumptions of the model are required without affecting the primary purpose: (1) The weaving process is uniform, and the weaving structure is complete and non-porous. To model the trajectory of a knot length as an example, firstly establish the xyz coordinate axis, the quarter of the knot length is h, and the starting coordinates are (0, b, 0) through the first coordinate transformation the coordinates are (0, b, h), continue the next coordinate transformation, the coordinates of the second become (a, b, 2h), then the coordinates of the third become (a, 2b, 3h) finally the fourth. The coordinates of the fourth transformation are (2a, 2, 4h), thus completing the construction of a knot length. The sample curves generated by the four coordinate transformations are swept after midpoint capture to generate yarns, and finally the yarns are arrayed according to the yarn weaving pattern.
The molding process of the three-dimensional fourdirectional preform is completed on the braiding machine, as shown in Figure 9(a) [24], which is a schematic diagram of the braided preform. First, the arrangement of the yarn carriers of the chassis needs to be initialized as shown in Figure 9(b). Then the yarn movement is driven by the cross movement of the yarn carriers, thereby forming the space of the yarns to interweave with each other, and then through the tightening process to tighten the fiber bundles wrapped together, every four steps during the operation of the machine is a machine cycle, so it is called the four-step method. The movement rule of the yarn carrier in a process is: first, the whole row of odd-numbered rows moves one column to the left, the entire row of even-numbered rows moves one column to the right, and the moving directions of two adjacent rows are opposite. The whole column moves up one row, and the two adjacent rows move in opposite directions. Then, the entire row of odd-numbered rows moves one column to the right, and the whole row of even-numbered rows moves one column to the left. Finally, the whole column of odd-numbered columns is moved up by one row, and the whole column of even-numbered columns is pushed down by one row; the schematic diagram of the weaving process is shown in Figure 9(c)-(f). The yarn arrangement of the preform is such that a yarn carrier carries only one yarn from beginning to end. The yarn arrangement and total number are shown in the formula (2.1).
where N is the total number of yarns, N m is the number of rows of the prefab, and N n is the number of columns of the prefab.
After the yarn has passed N s steps, the yarn carrier returns to the initial position.
where G is the greatest common divisor of N m and N n . According to the row material parameters, the diameter of the yarn is d, and the length of the flower pitch is h which is the length of the preform grows in the longitudinal direction in each machine cycle. The relationship between each braiding angle and the length of the flower pitch is as follows:

3)
According to the four-step weaving process, the weaving process is introduced by taking the 15°braided material as an example. Take a complete machine cycle process to track the movement track of the yarn carrier. Since the yarn carrier at a specific position has the same law, only a few typical yarn carriers can be tracked to obtain the same rules. Therefore only a few typical yarn carries can be tracked to obtain the coordinates of all the same regular yarns at the end of each step. After each machine cycle, the preform should go through a tightening process to make the structure tight. If the coordinate points are directly connected, the resulting yarn is a folded line, which does not conform to the actual yarn direction. Therefore, take the two adjacent coordinates, the midpoint is used as the control point of the yarn model, and the trajectory of the yarn is fitted by a cubic spline curve through the control point. The sweeping trajectory of the available yarn is shown in Figure 10  The yarn section is drawn according to the diameter of the yarn, and the braided preform is obtained as shown in Figure 10(b) and (c).
The yarn model of UG is then imported into the commercial finite-element software ABAQUS. First, the yarns in each position are assembled to generate the overall model of the prefab. According to the size of the prefab, the model of the matrix is established through the preprocessing module of the simulation software. The size is the same as that of the prefab, and the height is the length of three flower pitches. Then, the yarn and the matrix are assembled to obtain a three-dimensional four-directional braided composite material plate model.

Impact model
From the characteristics of the impact problem, it can be seen that the impact simulation calculation is a transient response problem, and the dynamic display solver should be used for the simulation calculation. The density of a given material is required for a dynamic response problem. In order to reduce the amount of calculation and thus shorten the calculation time, the size of the model was reduced in the length and width directions. The weaving arrangement of the prefab was taken as four rows and eight columns, and the length of three flower pitches was taken in the longitudinal direction. The overall size of the model was 20 mm × 6 mm × 25.2 mm. Combined with the experiment, a spherical punch was selected for simulations, the punch was set as a rigid body, and the diameter of the punch was 6 mm. In order to make the punch reach the same mass as the punch in the test, a mass point of 7.26 kg was applied to the punch. Since the yarn is a transversely isotropic material, this needs to define the main direction of the material separately, specifying the longitudinal direction of the fiber bundle; the matrix is an isotropic material, and the local coordinate system is consistent with fiber yarn. After assigning properties to the material, the local coordinate system to specify the direction of the material is selected. The material properties and mechanical properties of the fiber bundle yarn and matrix are given in Table 1. Figure 11 shows the prefabricated body model, positioning fixture model, three-dimensional four-directional braided composite flat plate model, and impact overall model.

Meshing and boundary loading
First, we mesh the model. The mesh size of the fiber bundle mesh is 0.3 mm. Because the shape of the yarn is a curved cylindrical with a circular cross-section, the yarn can be meshed from the bottom up. The shape is a hexahedral element, the mesh type is the first-order element C3D8R, and the cell seed size is 0.2 mm. Although the matrix is an isotropic material, the post-processing analysis still requires a local coordinate system, so the three directions of the global coordinates are used as the one direction of the matrix. The sweep method was adopted for the matrix mesh. The stacking direction of the custom element is selected. Along the growth direction of the flower pitch, the mesh type of the matrix is a C3D8R solid element, and the element size is 0.4 mm. The division results of fiber bundles are shown in Figure 12  Considering the reinforcement relationship between the fiber bundle and the matrix, the combination of the fiber bundle and the matrix is selected as the embedded form, because the embedded relationship does not require the relative geometric position of the nodes and allows the parts of the component to be embedded in each other. The punch and the specimen are selected in the general contact form. The normal direction is a hard contact, and the tangential direction is the penalty function form to simulate the non-smooth contact between the punch and the specimen. The friction coefficient is 0.25. According to the contact range between the punch and the specimen, a general contact pair is set. Considering that the punch still has a contact relationship with the surface of the internal unit after the unit is deleted in the impact process, all the surfaces of the unit that may be in contact with the punch are selected as the surface of the specimen. According to the clamping state of the specimen during the impact test, the degree of freedom of the direction of the specimen is limited by the model of the impact positioning fixture. After adjusting the spatial position of the fixture, all the degrees of freedom of the positioning fixture are limited. A small amount of preloading pressure is applied to the surface of the positioning fixture to simulate the pressing action of the pneumatic fixture to ensure that the specimen does not move in the direction of the speed of the punch. In the initial analysis step, the punch has only a translational degree of freedom along the one direction. The initial velocity field is applied to the punch, and the direction is along the negative direction of the x-axis. According to the kinetic energy theorem, the velocity at the beginning of contact between the punch and the specimen is calculated to get V = 2806.1 mm/s. The displacement load constraints of the model are shown in Figure 12(d).
After the first impact is completed, in order to avoid the problem of missing element material properties caused by element distortion and element deletion, this paper proposes a method to apply multiple impact loads by controlling the displacement of the punch and the initial contact velocity. After the first impact is completed, the distance from the punch to the specimen is obtained from the visualization module, and the punch is translated near the impact surface of the specimen by step 2. After the translation is completed, the pinch is given a velocity along the negative direction of the x-axis in step 3. Both analysis steps will generate displacement, and it is necessary to ensure that punch cannot contact the specimen after the two analysis steps are completed. After that, the punch can obtain the required initial velocity, and then the impact simulation is carried out on the specimen. Since the purpose of the translation and speed assignment of the punch is to obtain a certain initial speed before the punch contacts the specimen, it is not necessary to pay attention to the translation and acceleration process of the punch, which can greatly shorten the time of these two analysis steps and improve the calculation efficiency.

Damage criterion
The different structures of composite materials make their failure modes more complex, so different damage failure criteria need to be defined according to different types of composite materials. Since the damage laws in ABAQUS are mainly aimed at plane stress states, it is necessary to introduce three-dimensional damage criteria through the VUMAT subroutine to judge the damage to the matrix and fiber bundle materials. Then the yarn is a carbon fiber-reinforced transversely isotropic material, and the failure modes of the fiber bundle are mainly the extension of the fiber bundle in the axial direction, the compression of the fiber bundle in the axial direction, the transverse stretching in the cross-sectional direction of the fiber bundle, and the transverse compression in the cross-sectional direction of the fiber bundle, taking the growth direction of the braided flower pitch as the longitudinal direction. The three-dimensional Hashin initial damage criterion with strong applicability is used for yarns [25]. The matrix failure mode is mainly matrix tension and matrix compression, the matrix is an isotropic material, and the maximum stress failure criterion is used for the matrix [25]. The formula for each failure criterion is as follows: Longitudinal stretching of fiber bundles:  Fiber bundle lateral compression: 1. When the stress is greater than or equal to the strength, it is judged as failure.
Where Matrix stretch failure: Matrix compression failure: 0.

Damage evolution
When the impact test is not performed, the fiber beam yarn and matrix material of the three-dimensional fourdirectional composite are not damaged, and the rigidity of the fibers and the substrate is the material's inherent stiffness. After the impact test starts, when the stress state of the fiber unit satisfies the 3D Hashin initial criterion and the matrix unit satisfies the maximum stress failure criterion, damage to the matrix and fiber bundle materials begins to occur and stiffness degradation of the damaged unit occurs. In this paper, the damage-containing stiffness matrix of the material is recalculated by setting damage variables of each type, and stiffness degradation of the material is performed through this method. When the fiber bundle unit satisfies the initial damage criterion, the terms of the stiffness matrix are reduced accordingly. The stiffness matrix changes from a non-destructive state to a damaged state, and the unit is destroyed when the damage state variable reaches 1.
The coefficients that account for the stiffness matrix that considers damage are as follows. Initial stress-strain relationship: After satisfying the initial damage criterion, the stress-strain relationship is: where σ is the stress, C 0 is the non-destructive stiffness matrix, C d is a matrix containing damage, and ε is the strain.
The horizontal isotropic materials of the yarn fiber bundle, C 0 and C d relationships are shown in formula (2.12). The relationship between the damage stiffness matrix and the non-damaged stiffness matrix is the same as the form of the fiber bundle but the material is isotropic material, which is more concise in form. where C ij d are the coefficients of the stiffness matrix after damage, C ij 0 are the coefficients of the non-destructive stiffness matrix; S t 2 and S c 2 are the influence factors of residual shear stiffness; G ij d is the shear modulus after damage, and G ij 0 is the initial shear modulus before damage.
The realization process of damage criterion and damage evolution is as follows. First, the computer reads the initial parameters and strength values of the non-destructive material, calculates the non-destructive stiffness matrix of the material, calculates the initial strain and state variables, and updates the stress. The damage to the unit is judged by the state variables. When the material is damaged, the coefficients of the      material stiffness matrix are reduced, the stiffness matrix is recalculated, and then the strain is updated, the state variables are recalculated, and the stress is updated again. A damage judgment is made, and the process is repeated until the material fails completely.

Results and discussion
Combined with Figure 4 and 13, whether it is the impact process or the peak value of the impact load, the result of the finite-element simulation and the experiment are relatively similar. In the finite-element simulation, the unit is gradually damaged during the impact of the spherical punch on the specimen. Since the contact between the punch and the specimen is spherical, the surface of the specimen is prone to unit deletion, so the load-time curve has an obvious sawtooth shape, especially the brief drop in load was most pronounced during the first impact.
Combined with the experimental results of multiple impacts, the damage evolution analysis was carried out by taking the multiple impact process of a 15°specimen as a typical example. As shown in Figure 14(a), the first impact process is divided into the following stages. The punch starts to contact the specimen to the linear elastic stage, and the surface and interior of the specimen are not damaged. Then the load first decreased slightly and afterward increased to the maximum value of the contact force, indicating that the specimen was slightly damaged. In the next stage, the contact force fluctuated. The impact velocity of the punch head is obviously changed by the blocking action of the specimen, so the load fluctuates obviously. In the final stage, the contact force begins to decline and the punch begins to rebound. The subsequent multiple impact process is different from the first impact process. As the specimen has been damaged after the first impact, damage expansion and damage accumulation occur in the multiple impact process, so the stiffness of the specimen does not change by stages as in the first impact. The impact process lasts for a short time when the last impact occurs. The reason is that when the units at the impact position are damaged, the damaged units will be deleted and lose their blocking effect. A comparison of load peaks between the test and simulation is shown in Figure 14(d). Since the boundary conditions of the finite-element model and the experiment cannot be guaranteed to be completely consistent, there is some error in the test and simulation results, but the error is less than 20%, which verifies the validity and reliability of the simulation method process. Figure 15 shows the expansion process of longitudinal tensile damage along the thickness direction of a three-dimensional four-directional braided flat plate during the first impact. SDV1 represents the state variable of tensile damage of the material along the longitudinal direction. When SDV1 is greater than or equal to 1, tensile damage occurs along the longitudinal direction of the material. The tensile damage of fiber bundles mainly spreads from the surface of the impacted side of the specimen to the inside because the fiber bundles near the punch side are stretched first. As the impact process progresses, the inner fibers are gradually bent and deformed by the extrusion of the outer fibers and the matrix, followed by different degrees of tensile damage to the fiber bundles. In addition, because the fiber bundles in the preform are intertwined with each other, the fiber tensile damage spreads unevenly. The tensile damage of the matrix is mainly concentrated in the contact point between the specimen and the punch. The reason is that the toughness of the epoxy resin matrix is poor. Therefore, when the depth of impact embedded in the specimen is shallow, the matrix only has a small amount of tensile failure at the contact position. Figure 16 shows that the direction of movement of the impact is perpendicular to the longitudinal direction of the matrix, so the compression damage in the matrix is not obvious. The angle between the axis of the fiber bundle and the longitudinal direction of the substrate is approximately 15°, and the punch will impact laterally with respect to the fiber bundle. Therefore, during the pressing process of the punch, a small amount of longitudinal compression damage to the fiber bundle will appear on the inside of the bent fiber bundle. Figure 17 shows the simulation results of transverse tensile damage of the three-dimensional four-directional braided composite flat plate. The punch first squeezes the matrix along the impact direction, and the epoxy resin has poor toughness, so a small amount of resin on the impact side surface is crushed and peeled off. The damage to the fiber bundles on the side away from the punch is gradually prominent, because the impact energy is not enough to fracture the fiber bundles in the entire thickness direction. However, the longitudinal and transverse tensile mechanical properties of the fiber bundles are quite different. The fiber bundles on the back of the specimen are obviously subjected to transverse stretching and downward extrusion, so there is obvious transverse tensile damage to the interior and back of the specimen. Figure 18 shows the simulation results of transverse compression damage during the impact of the braided plate. It can be seen from the figure that the lateral compression damage of the specimen mainly occurs near the impact surface of the specimen, mainly due to the crushing of the matrix and the fiber bundle. The transverse mechanical properties of the matrix depend to a large extent on the mechanical properties of the matrix. In addition, because the two sides of the specimen are limited by the positioning fixtures and cannot move laterally, the lateral compression damage of the fiber bundle is more likely to expand laterally, while the damage expansion of the matrix is small, which only occurs in the two parts where the punch is embedded in the specimen.
After the first impact is completed, the specimen is subjected to multiple impact simulations, and the damage state of the central layer section of the specimen along the thickness direction is selected for analysis. Figures 19-22 show the results of various damage extensions of the woven composite panel at the end of the 1st, 5th, 10th, and 20th impacts. Observing the extension process of longitudinal tensile damage, it can be seen that longitudinal tensile damage of fiber bundles appeared in the specimen after the first impact, but the damaged area was small. After the fifth impact, the longitudinal tensile damage expanded significantly and gradually extended to the longitudinal direction, and the tensile damage of the back fiber bundles began to appear. The reason is that with the accumulation of impact damage, the fiber bundles on the impact surface fail during the cumulative impact process, the back fiber bundles are stretched, and the longitudinal tensile stress of the fiber bundles increases. When the tensile stress meets the failure criterion, damage occurs. Observing the expansion process of longitudinal compression damage, it can be seen that the compression damage of the fiber bundle along the longitudinal direction is less, and the compression damage of the matrix along the longitudinal direction is mainly concentrated in the contact part between the punch and the specimen. Also, with an increase in the number of impacts, the damaged area of the matrix does not expand significantly.
The overall transverse tensile damage and transverse compression damage of the specimen are more obvious in the central section, and the damage speed is faster. This is because the longitudinal tensile strength of the fiber bundles is high. During the impact process, the fiber bundles on both sides of the punch are compressed, and the fiber bundles under the punch are pulled until the fiber bundles and the matrix under the punch are transversely cracked. The vertical cracks appear on the lower surface of the specimen as shown in Figure 23. With the accumulation of the number of impacts, the cracks gradually expand and extend in the longitudinal direction. Combined with the ultrasonic C-scanning results of the impact specimen in the experiment, as shown in Figure 6, the damaged area of the fiber bundle and the matrix as a whole is in good agreement with the shape of the damaged area inside the ultrasonic C-scan specimen.

Conclusion
In this paper, the damage region expansion law of three-dimensional four-direction braided composites with 15°and 30°braided angles is studied during the  accumulation of impact number. The finite-element simulation software ABAQUS was used to simulate the composites, and the failure modes and damage propagation process of the preform and matrix in the impact process were analyzed.
(1) The main damage forms of the material during the accumulation of impact times are longitudinal cracking at the impact position and compression damage of fiber bundles and matrices. When the braided angle is small, the damage mainly expands in the longitudinal direction and its streamlined. When the braided angle is large, the damage area expands around, and the shape of the damage area is approximately circular. The mechanical properties of the braided plate with a larger braiding angle are more stable. (2) The braiding angle is the main factor affecting the impact resistance and axial compression performance of the three-dimensional four-directional braided composite material. The larger the braiding angle, the better the transverse impact resistance of the material, but the worse the longitudinal compression performance. High-performance materials can be obtained by adjusting the weaving angle according to the application of the material. (3) According to the finite-element simulation results, the damage modes of the specimen are mainly the longitudinal tension of the fiber bundle, the transverse tensile damage and the compression damage of the fiber bundle and the matrix; the damage first expands on the surface and then along the thickness until the specimen is penetrated.
Funding information: This research is supported by "National Natural Science Foundation of China" (No. 11972140).

Conflict of interest:
The authors declare that they have no conflicts of interest to report regarding the present study.
Data availability statement: The datasets generated during and/or analysed during the current study are available from the corresponding author on reasonable request.