On the response of re-entrant auxetic structure under dynamic crushing

. The present work focuses on the dynamic crushing response of 2D re-entrant auxetic honeycomb. The aim of this work is to propose analytic formulations to predict the crushing strength of the considered structure by including design parameters (specific geometrical ratios and/or material properties). The model for dynamic impact is based on shock waves propagation analogy in a rigid, perfectly plastic, locking (R-P-P-L) material model. The formulations depend on the geometric and the material characteristics of the auxetic but also on the impact velocity. Finite Element simulations of the impact of a plate on an auxetic structure were carried out using the RADIOSS ﬁ explicit solver. The impact velocities varies from 35 m / s to 200 m / s . The numerical simulations presented in this study show good accordance between analytical and Finite Element results.


Introduction
Structured materials, known for their multifunctional application [1], are remarkable for their particular stressstrain curves: after an elastic part where the curve is quasi-linear, a plateau stress interval is reached where the stress is almost constant and finally the curves rise abruptly at densification strain [2].During the past few decades [3,4], a new class of structured materials have emerged called ≪ auxetic ≫ (from Greek ≪ auxetos ≫: that may be increased).This kind of structure is characterized by a negative Poisson's ratio which gives it, compared to conventional, cellular structure a higher shear modulus, better indentation resistance and energy absorption ( [5][6][7]).
The present work focuses on one of these improved properties: the energy absorption capacity of auxetic structures (AS).Among the large range of structures having a negative Poisson's ratio, the re-entrant honeycomb is one of the most investigated geometry.Experimental samples are relatively easy to develop using 3D printing technologies [8].In addition, this geometry can be implemented in FE models without significant difficulties.That is why it has been chosen for this study.
* e-mail: adeline.petit@etu.univ-orleans.fr In-plane auxetic crushing test is the benchmark analysis to determine quasi-static and dynamic crushing strength of AS.This strategy is adopted by Ruan et al. [9], Hu et al. [10,11], and Hou et al. [12] to develop their models.In 1993, Reid et al. [13], on the basis that wood behave like cellular structure, proposed a shock model for wood using the rigid, perfectly-plastic, locking (R-P-P-L) simplification.This model of the dynamic crushing of cellular structure is validated in 2003 by Ruan et al. [9] thanks to numerical in-plane dynamic crushing of honeycombs.
Based on energy conservation principle and shockwave analogy, a theoretical (analytical) formulation of stress will be established.Finite Element simulations have confirmed the analytical model.
The objective is to propose a general methodology to assess the strength of an AS loaded by the impact of a plate.

Geometry of the structure
The complete auxetic structure, which will be investigated in compression along the Y -axis in this article, is an auxetic "re-entrant honeycomb" represented in Figure 1.The structure is composed of cells, as the one shown in Figure 2 with the notations used for the geometrical data.Note that S 0 is the total projected area on the X-Z plane, as shown in Figure 1.
A unit-cell geometry is characterised by its neutral axis (black dotted line in Fig. 2).The horizontal edge is denoted H and the bevel edge is denoted L. The angle formed by edges H and L is denoted θ.The cellwall thickness is denoted t thick and b is the out-of-plane thickness.
The edges dimensions h and ℓ are defined as follow : The initial cell height is defined by: At the initial time (t = 0), the total height in the Y -direction for n cell rows (layer) is: Two characteristic ratios are defined:

Relative density and densification strain
The relative density of a cellular structure is defined as the ratio of the cellular structure mass density ρ c to the mass density ρ s of the material of which it is constituted.
For the re-entrant honeycomb, it is defined as follow : where V c and V s are the cellular structure and material volumes, respectively.The relative density can be written in terms of the non-dimensional ratio τ and η : (8) The strain will be defined by the ratio of the height variation D 0 − D to the initial height D 0 , D being the height at time t: At densification the strain is ε DS and the height is D s , while the travelled distance of the upper wall is ε DS D 0 .
At this characteristic strain, the cells are all collapsed so that all cell-walls are in contact and the compressive stress is almost the compressive strength of the material itself.So, the height is D s = (4 n + 1) t thick and the densification strain ε DS is therefore expressed as: 3 Analytic prediction of stress

Stress enhancement through the shock front
Now, an impact on the AS is considered.It corresponds to impact the top surface (y = D 0 , Fig. 1) with a plate of mass M , while the other face (y = 0) is fixed.The velocity of the plate is v 0 when it impacts the AS at time t 0 = 0, and v for t > 0.
When the AS is dynamically crushed, the compaction starts from the impacted face and propagates like a shock-wave through the structure, involving a strong discontinuity (or "front") between the undeformed and deformed regions.The main steps of the shock-wave analogy applied to cellular material have been presented by [13].As a similar approach is adopted in this study, these steps are briefly recalled in Figure 3.The AS is idealized with the R-P-P-L model.Therefore, the stress is raised instantaneously to the static plateau stress σ stat at the initial time.The σ stat expression can be obtained by considering the plastic collapse of the cells network [14] : with σ y the yield stress of the material.Figures 3c and d show the deformed states at times t and t + dt.Behind the front, the stress is raised to the dynamic value σ and the structure is compressed up to densification strain ε DS .The structure is supposed to keep the same area S 0 and to deform only in the uniaxial direction.
In Figure 3, the shock wave propagates from proximal to distal end with the wave speed v f .From kinematic analysis, this velocity in lagrangian coordinates is expressed as: where u is the displacement of the plate and λ is the deformed length of the crushed AS.As the current velocity is v = du/dt, the front velocity relative to the impactor is : The conservation of mass gives : where ρ DS is the AS density when the structure is crushed up to densification.
During an infinetely small duration dt, the element of length dD (hashed in Fig. 3c) passes from an undeformed to a crushed state (hashed in Fig. 3d) and its lengths is dλ.The mass included within this volume hashed in Figure 3d, S 0 dλ, is expressed in terms of the cellular density ρ c : During this time interval, the volume element S 0 dλ considered above is submitted to the stress σ stat on its right face and to the stress σ on its left face.The velocity of the element increases from 0 to v during this duration.Thus the conservation of momentum of this element reads: From equation ( 16) it can be deduced that: This formulation is the same as the one proposed by Reid [13] to assess the stress inside the crushed part (shaded part in Fig. 3).After the impact, the crushing velocity v and the reaction forces of the AS vary constantly with time due to the dynamic coupling between the plate and the AS at the interface.It is proposed to evaluate this evolution with respect to time by using the energy conservation approach.

Analytical determination of displacement, velocity and stress
The conservation of energy of the whole system during the period t to t + dt (Figs.3c and 3d) reads: The left hand side of equation ( 18) includes a part of the structure which is set in motion by the impactor.This part shares a fraction of the system kinetic energy with the impactor.A fraction of the initial kinetic energy is transferred to the structure in the form of internal (strain) energy.Equation ( 18) is expanded with the help of equations ( 13), ( 12), ( 14) for dλ, dD and ρ DS respectively.In addition the second order term ρ c S 0 du dv/ε DS is neglected in the left hand side of equation (18).Finally, the non-linear differential equation describing the plate decelleration is: Where the constants B and C are defined by: Equation ( 19) can be written in term of the displacement u only by using the second and first time derivatives of u: Note that equation ( 22) holds only if u, du/dt are all both non constant functions of time.
Solving equation ( 22) with the initial conditions du/dt = v 0 and u = 0 at t = 0 gives : The time derivative of of u(t) gives the impactor velocity: Substituting equation (24) into equation ( 17) gives the formulation of dynamic stress with respect to time and impact conditions.When u(t) and v(t) are known the stress σ is directly calculated by equation (17).In addition, the following physical quantities of interest can also be determined : -Stop time for the rigid plate.According to equation (24), if the time t equals the velocity of the plate is theoretically zero.-Time for densification.At the particular time value t DS the AS begins to densified and the corresponding strain is ε DS given by equation (10).By equating the crushed length to ε DS D 0 , the time for  densification is the the appropriate (physical) time root t DS of the following equation:

Finite element model
A finite element model was developed to simulate the impact of the plate on the AS represented in Figure 1 and to evaluate the validity of the previous formulations (Sect.3.2).The commercial explicit non-linear FE software RADIOSS ® (Altair HyperWorks) [15] is used to implement this model.The model is composed of a top falling plate, an auxetic structure supported by a fixed plate at y = 0.The auxetic is made of six cells in height Y -direction and five cells in horizontal X-direction.The characteristics dimensions for the AS and the plates are given in Tables 1 and 2.An elastic-perfectly plastic aluminium is used for the auxetic material.A perfectly-elastic material is used for the two plates.The material properties are given in Tables 3. In the simulation, the crushing plate is supposed elastic with an elastic Young modulus equal to that of steel.The mass density is adjusted so that the crushing plate (with its imposed dimensions) has a mass equal to 4 g.This gives the plate a density of ρ s = 1.11 × 10 4 kg/m 3 .
In this simulation, Poisson's ratio is not a critical parameter, as the behaviour of the plates is assumed to be elastic.As its value can vary between 0.28 and 0.34, we have chosen 0.33 to represent an average value for the steel.Poisson's ratio also has no influence on the perfectly plastic elastic behaviour of the auxetic structure material.Thus, this material property is not a critical parameter in the simulation.For the sake of simplicity, we have set the same value of 0.33 for plates and aluminium.
The mesh size sensitivity has been analyzed with a coarse mesh and four gradually refined meshes (from 2 to 5 elements in the thickness, corresponding to sizes from t thick /2 to t thick /5).Variations of the impactor velocity results v for these meshes are presented in Figure 4.It is found that the order of the deviation is of 1% between the one and five elements in thickness meshes.If the element sizes are equal or lower than t thick /2, the results differ of less than 0.5%.Therefore, a mesh with two elements in the thickness (0.15 mm size) was retained.
The mesh of the three parts is composed of underintegrated 8-node hexaedral solid elements with physical hourglass stabilization.The plates are meshed with elements of 1 mm length.The auxetic structure is meshed with three elements in Z-direction.The fixed plate is constraint in all directions, the top plate is able to move only in Y -translation direction and, to avoid out-of-plane displacement, the auxetic is fixed in Z-translation direction.A penalty method general contact is applied to model all contact interactions at the plate/auxetic and auxetic/auxetic interfaces.The contact friction is not taken into account.

Results and dicussion
We consider the crushing of the AS by a falling plate of initial velocity v 0 in Y -direction.The simulations are done for a 4 g mass plate with initial velocity v 0 : 35, 50 and 200 m s −1 .The dynamic response of the AS obtained through the FE simulation is expressed in terms of stress σ (left scale of the graph) and the instantaneous plate velocity v (right scale of the graph) vs. time are plotted in Figures 5b, 6b and 7b, respectively.Figures 5a, 6a and 7a give respectively the AS deformation at selected times.
The plate is supposed always in contact with the AS.The stress σ FE is calculated by the ratio of the contact force (crushing plate -AS) to the initial transverse (X-Z plane) area S 0 .The dynamic stress σ is calculated with equation (17) by using the instantaneous velocity equation (24).The constant plateau stress σ stat , equation (11), is also plotted.The times t stop and t DS are calculated with equation (25) and equation (26), respectively.
The calculated stress, equation ( 17), leads to accurate estimation of the AS average stress values for the different impact conditions.
In Figures 5b, 6b and 7b an initial peak stress can be observed for a short duration which causes an important absorption of energy.Hence, during this same time lapse, the FE velocity decreases sharply.This initial peak is not taken into account in the theory developed in Section 3.2.That is why FE and analytic velocity differ from each other at the beginning.
At 35 m s −1 , the FE simulation shows that the plate is completely halted 1.85 ms after impact.The time t stop , calculated with equation (25), is 2.16 ms, which indicates that the method presented in Section 3.2 overestimates the stop time for low initial velocity.After the period of initial peak stress, values of plate FE velocity and analytical velocity, determined with equation (24), are very close.However, once the velocity is lower than about 20 m s −1 , the two values come apart from each other.
It can be seen, in Figure 5a, that the deformed region differs from that of the shock-wave model used to established the formulation of σ.In fact, the deformed region observed in the FE results corresponds to deformation modes, which depend upon v 0 and the AS topology and geometry, (see for example Ruan et al. [9] for constant crushing velocity).
At 50 m s −1 and 200 m s −1 , the AS was not able to stop the plate before densification.When densification begins the plate velocity is not zero.The plate will continue to crush the structure while the stress increases up to the peak and rebound of the plate will occurs.The prediction of t DS with equation (26), is reasonably accurate regarding FE results.For both, the evaluation of analytic plate velocity are very close to the FE velocity.
The formulations established in Section 3.2 seem to be more accurate regarding FE results for medium (50 m s −1 ) and high initial impact velocity (200 m s −1 ).
The alloy properties chosen in this article are close to those of the AL3024-T351 alloy.Khan et al. [16] have shown that if the strain rate does not exceed 1500 s −1 , the behaviour of the material is not very sensitive to strain rate.As we have verified that in the simulations of the dynamic compression of the structure, the dε/dt rate remains below 100 s −1 , we have chosen not to take the strain rate into account in the behaviour of the material making up the AS.

Conclusion
In this paper, closed-form formulas have been proposed for predicting the crushing strength of architectured   auxetic lattices loaded by the impact of a plate.The approach is based on energy conservation principle and shock-waves propagation analogy in a rigid, perfectly plastic, locking material model.It enables the analysis of periodic collapse of the structure.We have extended formulation in order to include more design parameters (specific geometrical ratios and/or material properties).Our analysis enables to theoretically predict the dynamic crushing strength.The formulation depends on the geometric and the material characteristics of the auxetic but also on the impact velocity.
The results of analytical calculations are consistent with the results of finite-element simulations.It is then possible to adapt the geometrical and material parameters in order to optimize the energy absorption ability for facing specific load cases.In-depth studies are in progress to understand the influence of the deformation modes on energy absorption.

Fig. 3 .
Fig. 3. Schematics and parameters defining a shock propagation into cellular bar after a rigid mass impact.

Fig. 4 .
Fig. 4. Finite element impactor velocity results v for several mesh sizes (from one to five elements in the thickness t thick ) at initial velocity v0 = 200 m s −1 .

Table 1 .
Data for the AS and plates in the FE model (lengths are in millimetre).

Table 2 .
Dimensions of the plates.

Table 3 .
Material properties for the FE model.