Impact Responses and Wave Dissipation Investigation of a Composite Sandwich Shell Reinforced by Multilayer Negative Poisson’s Ratio Viscoelastic Polymer Material Honeycomb

This analysis investigated the impact wave response and propagation on a composite sandwich shell when subjected to a low-velocity external shock, considering hygrothermal effects. The sandwich shell was crafted using face layers composed of functional gradient metal–ceramic matrix material and a core layer reinforced with negative Poisson’s honeycomb. The honeycomb layer consisted of a combination of viscoelastic polymer material and elastic material. The equivalent parameters for the functional gradient material in the face layers were determined using the Mori–Tanaka and Voigt models, and the parameters for the negative Poisson’s ratio honeycomb reinforcement core layer were obtained through Gibson’s unit cell model. Parameters relevant to a low-velocity impact were derived using a modified Hertz contact law. The internal deformations, strains, and stress of the composite sandwich shell were described based on the higher-order shear deformation theory. The dynamic equilibrium equations were established using Hamilton’s principle, and the Galerkin method along with the Newmark direct integration scheme was employed to calculate the shell’s response to impact. The validity of the analysis was confirmed through a comparison with published literature. This investigation showed that a multilayer negative Poisson’s ratio viscoelastic polymer material honeycomb-cored structure can dissipate impact wave energy swiftly and suppress shock effectively.


Introduction
Laminate structures find extensive applications in various engineering fields, including aerospace, aviation, military equipment, transportation, and mechanical engineering.This is primarily due to their exceptional characteristics, such as high specific strength and stiffness, low density, and remarkable design adaptability.As a result, the investigation of their mechanical behavior and structural design has garnered increasing interest among researchers and scholars [1][2][3][4][5].Galos [1] presented a review on fiber reinforced polymer matrix composites formed using thin-ply laminates and Garg and Chalak [2] performed a critical review of the available literature for the prediction of the behavior of laminated composites and sandwich structures under hygrothermal conditions.
In the fields of materials science and solid mechanics, Poisson's ratio serves as a fundamental indicator of the Poisson effect, which refers to the deformation experienced by a material in directions perpendicular to the specific direction of loading.It quantifies the relationship between transverse strain and axial strain.More precisely, Poisson's ratio is defined as the negative ratio of transverse strain to axial strain.In simpler terms, it measures the extent of transversal elongation relative to axial compression when subjected to small changes.For isotropic materials, Poisson's ratio has been proven to fall within the range of −1 to 0.5, based on thermodynamic considerations of strain energy.This theoretical framework provides a comprehensive understanding of the mechanical behavior exhibited by such materials [6].However, for most solid materials, Poisson's ratio typically falls within the range of 0 to 0.5.However, it is worth noting that in the case of soft materials like rubber, where the bulk modulus significantly surpasses the shear modulus, Poisson's ratio tends to approach a value near 0.5 [7].In the case of open-cell viscoelastic polymer material foams, Poisson's ratio tends to approach zero due to the collapse of the foam cells under compression.This unique behavior can be attributed to the specific structure and mechanical properties of these foams.
It is worth noting that many common solids found in nature typically exhibit Poisson's ratios within the range of 0.2 to 0.3.This range is observed in a wide variety of traditional mechanical materials, highlighting their characteristic response to external forces.However, it is also important to acknowledge that certain artificial materials and microstructures have the potential to provide effective equivalent Poisson's ratios that span the range between −1 and 1.These exceptional cases show the remarkable versatility and engineered capabilities of these materials in tailoring their mechanical properties to meet specific requirements [8].These materials offer remarkable mechanical properties that are advantageous in various engineering applications.They demonstrate exceptional energy absorption capabilities and resistance to fractures, shear deformation, and indentation.To circumvent the need for a cumbersome expression, Evans et al. [9] proposed that materials with a negative Poisson's ratio be referred to as "auxetic materials" or simply "auxetics." Negative Poisson's ratio materials can be broadly classified into four categories: natural negative Poisson's ratio materials, cellular negative Poisson's ratio materials, metallic negative Poisson's ratio materials, and multi-material negative Poisson's ratio composites.The initial reports of materials with a negative Poisson's ratio date back to the early 20th century, primarily observed in the realms of chemistry and biology.Notably, certain body-centered and face-centered cubic crystal minerals, as well as animal tissues, exhibited these intriguing characteristics [10].In a study conducted by Lakes, the implications of negative Poisson's ratios in component design under stress were investigated [11].Their findings revealed that when the Poisson's ratio assumes a negative value, stress concentration factors are reduced under certain conditions.This discovery sheds light on the potential benefits of incorporating negative Poisson's ratio materials in the design process, as they can aid in minimizing stress concentration effects.Baughman et al. have astutely highlighted that negative Poisson's ratios are a prevalent characteristic among cubic metals when oriented in specific directions [12].Their research suggested that approximately 69% of cubic elemental metals exhibit a negative Poisson's ratio when subjected to stretching along the (110) direction.Moreover, they have proposed several initial applications for such materials, including the use of negative Poisson's ratio metals as electrodes in piezoelectric strain sensors and as vanes in aircraft gas turbine engines.These innovative applications demonstrate the potential utility of negative Poisson's ratio materials in various technological domains.
Cellular negative Poisson's ratio materials are classified as a type of porous material, which endows them with low density, exceptional mechanical energy absorption, noise reduction capabilities, and robust damping properties.In practical applications, common two-dimensional concave structures, including concave hexagonal honeycomb models, concave triangle models (also known as the double arrow model), and star models, have a negative Poisson's ratio.It is important to note that the Poisson's ratio characteristic of these concave structures is intricately linked to their unique concave angles.However, it is essential to recognize that the presence of internal concave-convex shapes alone is not enough to guarantee a negative Poisson's ratio in the structure.Additional factors must be considered to ascertain the presence of a negative Poisson's ratio.In the case of a star-shaped two-dimensional system, it is worth noting that auxetic angles below approximately 20 degrees do not result in a negative Poisson's ratio.Moreover, an increase in the auxetic angle will lead to a reduction in the effective pure shear modulus, as highlighted in a comprehensive study [8].To address these considerations, an analysis incorporating a concave hexagonal honeycomb structure can be employed to induce a negative Poisson's ratio.In addition, to enhance mechanical energy dissipation and absorption, the introduction of viscoelastic damping materials has proven instrumental in this endeavor.Through this combined approach, the desired mechanical properties are amplified, offering potential benefits in a range of applications.
Functionally graded materials and/or structures exhibit a gradual variation in composition and structure throughout their volume, leading to consequential alterations in their properties.This unique characteristic allows for the intentional design and fabrication of materials and structures tailored to specific functions and applications.The concept of functionally graded materials and structures was initially introduced in the field of aerospace engineering, where the focus was on developing thermal barriers capable of withstanding exceedingly high temperatures.By strategically manipulating their composition and structure, these materials and structures can effectively serve as thermal barriers, shielding against extreme thermal conditions [13].In engineering applications, the distribution of material gradients can be classified as either continuous or discontinuous.Functional gradient materials, exhibiting superior mechanical behaviors, are abundantly found in nature.For instance, plants showcase a distribution of fibers in bamboo and fruit peels, while animals exhibit a micro-porous distribution in bones and certain tissues.These natural occurrences serve as inspiration, highlighting the extensive presence of functional gradient materials and their potential for enhancing mechanical performance in practical engineering applications.
In the realm of practical engineering, the pursuit of specific mechanical behaviors has prompted the proposal and design of bioinspired functional gradient materials and microstructures.These innovative approaches aim to create system that can effectively absorb external loads, dissipate vibration energy within the structures, and resist external impacts.To fulfill these requirements, Zhou and Wang proposed a functional gradient structure with high stiffness in the outer layers and flexible inner layers.This ingenious design enables the structure to achieve its intended purpose of load absorption, energy dissipation, and impact resistance [14].Furthermore, extensive investigations have been conducted on fiber-reinforced resin matrix materials with varying distribution configurations of functional gradient models.The insightful analysis conducted in these studies has unequivocally demonstrated that the implementation of a specified design for the functional gradient distribution of fibers can yield a remarkable increase in the buckling load capacity of structures.This finding underscores the significant impact that a welldesigned functional gradient distribution can have on enhancing the structural integrity and load-bearing capabilities of the material.The potential for enhancing the buckling load through tailored functional gradient fiber distribution holds great promise in the field of engineering [14].
In the past few years, graphene platelets have emerged as a novel reinforcement material for enhancing the mechanical behavior of structures.Particularly, the concept of functionally graded designs in graphene-reinforced composite structures has garnered increasing interest from researchers and scientists.In this regard, Zhao et al. [15] have presented a comprehensive and critical review, shedding light on the advancements and potential of functionally graded graphene-reinforced composite structures.This review serves as a valuable resource for researchers and scientists who seek to delve deeper into this exciting area of study.Meanwhile, a monograph covering functional gradients inspired by biological material in terms of design principles, functions, and applications has been provided by Liu et al. [16].FGM structures are extensively applied in mechanical engineering, such as in composite beams [17,18], plates [14,19], and shells [20,21].
As we know, there is a conflict between a material's strength and toughness [22].Ceramic matrix composite material with a functional gradient design is regarded as an efficient method to modify composites, which can achieve a material with extra high strength and toughness under high temperature, and have excellent microwave transparent and broadband wave transmission ability at the same time [23,24].Sofiye [25] presented an exhaustive review of the literature on the vibration and buckling of functionally graded materials, especially for mechanical behavior investigations related to functionally graded and laminated sandwich conical shells.By applying the higher-order shear deformation theory and smeared stiffeners technique, assuming Coriolis acceleration and centrifugal forces, Aris and Ahmadi [20] derived a governing equation of a cross-stiffened, rotating, truncated conical shell with functional gradient materials and investigated the influences of the cone angle, volume fraction, and thermal conditions on the critical frequency.Assuming the thermo-mechanical characteristic was distributed by the temperature and geometric position of the structure, Bagheri et al. [26,27] constructed a functionally graded material assembled as a spherical conical shell and presented a series analysis under rapid heating through the generalized differential quadrature method and the Crank-Nicolson method.Except for thermal impedance, functionally graded materials show excellent behavior of sound transmission loss by a proper design.By applying the space harmonic approach and virtual work principle, Fu et al. [28] established a corrugated core functionally graded material sandwich plate filled with porous material and analytically investigated the sound transmission loss; meanwhile, their analysis was validated through published studies.Truncated conical shells are fundamental structures that are applied in engineering to increase the mechanical performance and functional grade, and are used to design stiffeners with different configurations.Considering complex loads, such as axial compressive loads and external uniform pressure, Dung and Chan [29] provided an analytical investigation of the buckling behaviors of truncated conical shells reinforced with orthogonal stiffeners.
As recalled in the above literature, viscoelastic polymer materials are widely utilized in engineering for reducing and minimizing the sound radiation and vibration level by dissipation and absorbing mechanical energy propagation in composite structures.As all of us know, service conditions are the crucial effects that influence the parameters of viscoelastic polymer materials, such as Young's modulus, damping ratio, and Poison's ratio.To accurately evaluate the dynamic behavior of composite structures, a number of researchers have investigated the thermal and moisture behaviors of viscoelastic polymer materials [30,31].In order to eliminate power machinery vibration and noise, Qu et al. [32] prepared a novel semi-interpenetrating viscoelastic material through linear polyimide and crosslinked epoxy, and investigated the damping, thermal, and mechanical performances of this new viscoelastic composite.By using the representative volume element mode with a random distribution of carbon fibers, Yuan et al. [33] proposed a temperature-dependent micromechanical model to predict the temperature-dependent mechanical behaviors of this composite, such as the transverse tensile strength and transverse compressive strength, failure behaviors, and parametric sensitivity.Moreover, the moisture and humidity of the environment influence the mechanical parameters of viscoelastic composites.Under a complex external load, the mechanisms of moisture intrusion, water penetration, and diffusion in viscoelastic composites has attracted a lot of investigation [34,35].Ye and Zhang [35] established a thermodynamically consistent moisture diffusion model for coupling the moisture diffusion effect and viscoelastic response of the multiphase viscoelastic composites, developed a two-constituent phase-field fracture model for describing the hygroscopic swelling behaviors at the multiphase interfacial decohesion of viscoelastic polymer material, and proposed a crack filter theory for characterizing the moisture flux propagation at the interface between the multiphases.In fact, the hardening and softening behaviors of viscoelastic polymer materials are influenced by other factors, such as the strain rate, frequency, duration, and external load amplitude [36][37][38][39][40]. Shitikova and Krusser [39] presented a comprehensive review on the historical development and formulation of viscoelastic materials' models in application and investigation, which mainly covers the classical models and fractional model.Barba et al. [36] provided experimental and modeling insights into the mechanical properties of PEEK under different strain rates and temperatures.
To achieve a balance between high strength and toughness while enhancing temperature resistance, a functional gradient design method was employed for a ceramic matrix with a metal-filled composite.Simultaneously, a composite reinforced layer with a negative Poisson's ratio was created using a viscoelastic polymer material and elastic honeycomb reinforcement.Investigations in the present work are outlined as follows: the problem is provided, and then mathematical modeling, derivation processing, equivalent performance, and governing equation construction are proposed in Section 2; a determination of the response and wave propagation behaviors of the composite laminated shells is presented in detail in Section 3; numerical results and a qualitative discussion are provided in Section 4, and a conclusion is given in Section 5.

Mathematical Modeling
In this paper, we present the design of a novel functionally graded composite sandwich doubly curved shell.The core of this structure features a multilayer honeycomb with a negative Poisson's ratio, as illustrated in Figure 1.The face-sheets of the shell are composed of functionally graded materials, combining ceramic and metal layers.The honeycomb core consists of three layers, with the middle layer being a viscoelastic polymer material and the inner and outer layers being solid metal.To investigate the behavior of the composite sandwich shell, it is subjected to impact by a sphere with a low velocity V 0 in hygrothermal environments.The structure is defined by a global coordinate system (x, y, z) with the mid-plane set at z = 0.The principal curvature radiuses at the mid-plane are denoted as R x and R y along the x and y axes, respectively.The lengths of the curved edge in the x and y directions are represented by L x and L y , respectively.Both the top and bottom face-sheets exhibit symmetry about the mid-plane, with an equivalent thickness denoted as h f.The core layer of the structure possesses a thickness of h c .Each layer is bonded together well, and there is no sliding during the impact.
Materials 2024, 17, 233 5 of 30 formulation of viscoelastic materials' models in application and investigation, which mainly covers the classical models and fractional model.Barba et al. [36] provided experimental and modeling insights into the mechanical properties of PEEK under different strain rates and temperatures.
To achieve a balance between high strength and toughness while enhancing temperature resistance, a functional gradient design method was employed for a ceramic matrix with a metal-filled composite.Simultaneously, a composite reinforced layer with a negative Poisson's ratio was created using a viscoelastic polymer material and elastic honeycomb reinforcement.Investigations in the present work are outlined as follows: the problem is provided, and then mathematical modeling, derivation processing, equivalent performance, and governing equation construction are proposed in Section 2; a determination of the response and wave propagation behaviors of the composite laminated shells is presented in detail in Section 3; numerical results and a qualitative discussion are provided in Section 4, and a conclusion is given in Section 5.

Mathematical Modeling
In this paper, we present the design of a novel functionally graded composite sandwich doubly curved shell.The core of this structure features a multilayer honeycomb with a negative Poisson's ratio, as illustrated in Figure 1.The face-sheets of the shell are composed of functionally graded materials, combining ceramic and metal layers.The honeycomb core consists of three layers, with the middle layer being a viscoelastic polymer material and the inner and outer layers being solid metal.To investigate the behavior of the composite sandwich shell, it is subjected to impact by a sphere with a low velocity V0 in hygrothermal environments.The structure is defined by a global coordinate system (x, y, z) with the mid-plane set at z = 0.The principal curvature radiuses at the mid-plane are denoted as Rx and Ry along the x and y axes, respectively.The lengths of the curved edge in the x and y directions are represented by Lx and Ly, respectively.Both the top and bottom face-sheets exhibit symmetry about the mid-plane, with an equivalent thickness denoted as hf.The core layer of the structure possesses a thickness of hc.Each layer is bonded together well, and there is no sliding during the impact.Considering the remarkable properties exhibited by ceramics, such as high strength, exceptional hardness, and excellent temperature resistance, the face-sheets of the structure are fabricated using functionally graded materials that transition from ceramic to metal Considering the remarkable properties exhibited by ceramics, such as high strength, exceptional hardness, and excellent temperature resistance, the face-sheets of the structure are fabricated using functionally graded materials that transition from ceramic to metal from the exterior to the interior.The composition of the face-sheets is assumed to vary continuously throughout the thickness of the structure, following a power law distribution.The volume fraction of ceramic within the face-sheets can be mathematically expressed as: The volume fraction of metal can be deduced using a complementary approach as: The power law index, denoted as p, and the total thickness of the structure, represented by h (where h = 2h f + h c ), play crucial roles in determining the volume fraction of ceramic along the z direction.Figure 2 presents the volume fraction of ceramic for various values of the power law index p.As one moves from the inner surface to the outer surface of the face-sheets, the volume fraction of ceramic, V c , increases continuously from 0 to 1. Notably, it is evident that a decrease in the value of p leads to an increase in the volume fraction of ceramic.
from the exterior to the interior.The composition of the face-sheets is assumed to continuously throughout the thickness of the structure, following a power law distr tion.The volume fraction of ceramic within the face-sheets can be mathematically pressed as: The volume fraction of metal can be deduced using a complementary approach ( ) ( ) The power law index, denoted as p, and the total thickness of the structure, re sented by h (where h = 2hf + hc), play crucial roles in determining the volume fractio ceramic along the z direction.Figure 2 presents the volume fraction of ceramic for var values of the power law index p.As one moves from the inner surface to the outer sur of the face-sheets, the volume fraction of ceramic, Vc, increases continuously from 0 Notably, it is evident that a decrease in the value of p leads to an increase in the vol fraction of ceramic.

Equivalent Mechanical Parameters
Before proceeding with further derivation, it is necessary to determine the equiva mechanical properties of the sandwich shell.The face-sheets, being composed of funct ally graded materials, can have their equivalent mechanical properties obtained using law of mixture.This enables us to calculate the overall properties of the face-sheet considering the distribution and composition of the constituent materials.

( ) (
) ( ) where P(z) represents the equivalent mechanical properties, such as Young's modu shear modulus, mass density, Poisson's ratio, and more; and the subscripts c and m de the ceramic and metal constituents, respectively.Considering the influence of temp ture variation on the mechanical properties of both ceramic and metal, the temperat

Equivalent Mechanical Parameters
Before proceeding with further derivation, it is necessary to determine the equivalent mechanical properties of the sandwich shell.The face-sheets, being composed of functionally graded materials, can have their equivalent mechanical properties obtained using the law of mixture.This enables us to calculate the overall properties of the face-sheets by considering the distribution and composition of the constituent materials.
where P(z) represents the equivalent mechanical properties, such as Young's modulus, shear modulus, mass density, Poisson's ratio, and more; and the subscripts c and m denote the ceramic and metal constituents, respectively.Considering the influence of temperature vari-ation on the mechanical properties of both ceramic and metal, the temperature-dependence of the material properties can be expressed using the following relationships [41]: where P 0 , P 1 , P 2 , and P 3 are the temperature-dependent coefficients.For materials Si 3 N 4 and SUS304, the temperature-dependent coefficients are listed in Table 1.
To obtain the equivalent modulus of the core layer, the multilayer honeycomb walls should be firstly replaced with a single equivalent layer, as shown in Figure 3, and then the analytical expressions derived by Gibson et al. [43] can be used.Based on the classical lamination theory, the relationship between the forces and strains of multilayer structures can be expressed as: To obtain the equivalent modulus of the core layer, the multilayer honeycomb wa should be firstly replaced with a single equivalent layer, as shown in Figure 3, and the the analytical expressions derived by Gibson et al. [43] can be used.Based on the classic lamination theory, the relationship between the forces and strains of multilayer structur can be expressed as:   The extensional stiffness, coupling stiffness, and bending stiffness matrices, denoted as A, B, and D, respectively, play a crucial role in characterizing the mechanical behavior of the structure.The expressions for these matrices are as follows:

Material
Equation ( 5) can be reformulated as follows: In accordance with the concept of Young's modulus, let us assume that the structure experiences solely normal strain with no shear strains, denoted as κ = 0.Under this assumption, the expression simplifies to the following: Writing Equation ( 9) in its explicit form as follows: Due to the symmetry of matrices A, B, and D, the resulting matrix P is also symmetric.To obtain the equivalent axial modulus, the multilayer structure must be placed in a unidirectional tensile state.To achieve this, an additional shear force is applied to counteract any shear deformation.To ensure γ 0 xy = 0, the needed shear force: Substituting Equation (11) into Equation ( 10), according to the stress-strain relationship, the in-plane axial moduli can be expressed as: E yy = 1 P 22 − P 2 26 P 66 t where t represents the total thickness of the multilayer structure.The equivalent shear modulus can be expressed as follows when the structure is subjected to a pure shear state: Then we can calculate the equivalent mechanical properties of the multilayer honeycomb core based on the expressions [44]: Materials 2024, 17, 233 9 of 29 where θ is the inclined angle, and η 1 and η 3 are two dimensionless parameters related to the geometry of the honeycomb,

Impact Forces Determination
For analyzing the dynamic responses of this composite laminated shell, it is imperative to determine the transient impact load acting on the composite laminate shells.Considering the indentation period and contact area between the impacting body and the composite laminate, this contacting load can be regarded as a force varying over time.This force can be expressed as: where F c denotes the impact force between the indenter and the composite laminate shell, δ denotes the Dirac delta function, and (x c , y c ) defines the coordinates of the indentation center.For our investigation, we have selected the central point of the sandwich, denoted by (L x /2, L y /2).The impact force is determined through the modified Hertz contact law.In this analysis, two stages are defined, and the contact force can be obtained by applying the following equation [45]: Loading stage Unloading stage (27) The modified Hertz contact stiffness, denoted as k c , plays a crucial role in this analysis.By utilizing this stiffness, we can determine the indentation α as follows: In the above equation, w(t) expresses the displacement at the contact point in the transient analysis, and w 0 expresses the deflection of the composite sandwich shell just at the point where the impact contacts the impactor, as shown in Figure 4.
sis.By utilizing this stiffness, we can determine the indentation α as follows: ( ) ( ) In the above equation,⎯w(t) expresses the displacement at the contact point in the transient analysis, and w0 expresses the deflection of the composite sandwich shell just at the point where the impact contacts the impactor, as shown in Figure 4.In Equation (27), αmax and F max c express the maximum indentation and contacting force in the loading procedure, and α0 denotes the permanent indentation generated by the impact load.In the present analysis, permanent damage to the composite sandwich shell under a low-velocity impact is not taken into consideration, and thus it is assumed that α0 is equal to zero.According to Ref. [46], the modified Hertz contact stiffness kc can be yielded according the following equation: where Eimp is Young' modulus, νimp denotes Poisson's ratio, Rimp means the impactor's radius, and ν22 and E22 express the Poisson's ratio and Young' modulus of the top and bottom layers of the functional gradient sandwich structure, respectively.

Constitutive Equations and In-Plane Variables Determination
According to the higher-order shear deformation theory [47], the displacement of the functional gradient sandwich shell can be defined as: where u0, v0, and w0 are the displacement variables at the mid-plane of the laminated sandwich shell; ϕx and ϕy express the rotations of the normal to mid-plane about the y and x α w In Equation (27), α max and F max c express the maximum indentation and contacting force in the loading procedure, and α 0 denotes the permanent indentation generated by the impact load.In the present analysis, permanent damage to the composite sandwich shell under a low-velocity impact is not taken into consideration, and thus it is assumed that α 0 is equal to zero.According to Ref. [46], the modified Hertz contact stiffness k c can be yielded according the following equation: where E imp is Young' modulus, ν imp denotes Poisson's ratio, R imp means the impactor's radius, and ν 22 and E 22 express the Poisson's ratio and Young' modulus of the top and bottom layers of the functional gradient sandwich structure, respectively.

Constitutive Equations and In-Plane Variables Determination
According to the higher-order shear deformation theory [47], the displacement of the functional gradient sandwich shell can be defined as: where u 0 , v 0, and w 0 are the displacement variables at the mid-plane of the laminated sandwich shell; ϕ x and ϕ y express the rotations of the normal to mid-plane about the y and x axes; and function f (z) determines the distribution of the transverse shear strains along the shell thickness.A lot of different functions can be applied.In this work, the cubic function was selected as follows: where h is the total thickness of the functional gradient laminate structure, h = 2h f + h c .According to the strain displacement relationships of the shell structure: All nonzero strain components have the form: where the strains in-plane can be expressed as: And meanwhile, the transverse shear strains can be denoted as: The constitutive relationship for the nth layer of these functional gradient composite sandwich shell structures can be expressed as: In the above equations, [Q] (n) expresses the matrices of stiffness, {σ} (n) is the stress, and {ε} (n) denotes the strain vectors of the composite sandwich shell at the layer nth, respectively.The variables in Equation ( 37) can be expressed as: The coefficients of the equivalent stiffness are dependent on the Young's modulus and Poisson's ratio of the materials.For the functionally graded material face-sheets, the coefficients Q ij are also related to the space coordinates, which are: For the negative Poisson's ratio core layer of the composite structure, the coefficients Q ij are: (40) In Equation (38), two parameters, ε HT xx and ε HT yy are strains defined to describe the hygrothermal effects, and those two variables can be determined according to the following equations: In Equations ( 41) and ( 42), α 11 and α 22 are coefficients of thermal, and β 11 and β 22 are moisture expansions in 11 and 22 directions, respectively; T means temperature and H means moisture.In the present research, three different distribution patterns of temperature and moisture along the thickness of the structure are considered.

Dynamic Equilibrium Equations
For this functional gradient composite sandwich shell, according to Hamilton's principle, the dynamic equilibrium equations can be achieved as: where δU means virtual strain energy, δV the virtual work done by external forces, and δK describes the virtual kinetic energy.
The virtual kinetic energy of this functional gradient composite sandwich shell can be achieved as: By substituting Equation (31) into Equation (47), we can determine the virtual kinetic energy as: where the mass inertias can be defined as: Furthermore, the virtual strain energy of the functional gradient composite sandwich shell can be achieved as: Furthermore, the forces and moments variables of the sandwich can be obtained as: Moreover, the above resultants can be also arranged as: The stiffness coefficients of the functional gradient composite sandwich shell can be determined by: Meanwhile, the δV produced by external loads for this composite sandwich can be calculated as: where, hygrothermal forces N HT xx and N HT yy , can be determined according to [48]: Substituting Equations ( 47), (52), and (58) into Equation (46), and then applying the integration by parts, the dynamic equilibrium equations are expressed as: For the impactor, the equation of motion can be expressed as: in which the impactor's mass is defined by m 0 .

Solution Process
In order to perform a numerical analysis for the dynamic equations at an arbitrary boundary condition, the Galerkin method was introduced, and the Newmark direct integration was applied.The displacements of the functional gradient sandwich shell can be discretized as a sum of a series of functions of the following form according to reference [49]: where U mn (t), V mn (t), W mn (t), Φ mn (t), and Ψ mn (t) are amplitudes that are required to be defined, and X m (x) and Y n (y) are applied to satisfy the equation for an arbitrary boundary condition at defined edges.These spatial basis functions are provided in Table 2.Moreover, the following initial boundary conditions are achieved: U mn (0) = .
The initial displacement and velocity of the impactor can be expressed as: The composite sandwich shell is assumed to have simply supported (S) or clamped (C) edges, or combinations of them, and they are given as: Simply supported (S): Clamped (C): Substituting Equations ( 66)-(70) into Equations ( 60)-( 64), the motion of the composite shell can be expressed as: where [M], [C], [K], {X}, and {F} are the matrices of mass, damping, and stiffness, and the vectors of deformation and equivalent external load, respectively.In the presented analysis, the symmetric laminate sandwich shell is constructed; thus, these matrices are all symmetric.The mass matrices and stiffness matrices are provided in Appendix A. The Rayleigh's proportional damping model is employed in this analysis, incorporating a damping matrix that consists of both mass-and stiffness-proportional terms.This damping matrix can be expressed as follows: where α D and β D are Rayleigh's damping parameters, and those two parameters can be determined according to the viscoelastic material's damping ratio ξ and the two fundamental frequencies ω i and ω j of the composite laminated shell, as shown in the following equation as: where the damping coefficients α D and β D are associated with the Rayleigh's proportional damping model and can be determined based on the damping ratio ξ and the two natural frequencies of the structure ω i and ω j .The relationship between these parameters can be expressed as: In the presented investigation, the Young's modulus of the viscoelastic polymer material is assumed to be complex, as: where i is the imaginary unit, η VEM is the viscoelastic polymer material's loss factor, and E expresses the storage modulus.Then, equations of the free vibration can be derived as an eigenvalue problem, which can be expressed as: where ω* expresses the complex frequency.Furthermore, the natural frequency and damping ratio of the viscoelastic composite sandwich can be achieved as: For solving the impact response of the composite laminate shell structure, the Newmark direct integration algorithm is utilized, which discretizes differential Equation (75) into a time domain.Check Ref. [42] for more detail about the preformation of the Newmark direct integration algorithm.

Validation
Case 1: the validation of dimensionless natural frequencies of simply supported Al/Al 2 O 3 functionally graded plates and shells.The current analysis model can be applied to research the dynamic characteristics of a plate structure if the radiuses R x = R y = ∞.In this case, the geometric parameter of the structure is L x = L y = 10h and the material properties of the structure are shown in Table 3.As shown in Table 4, the dimensionless frequencies of the structures obtained by the present method are almost the same as those obtained in Ref. [51].Hence, the effectiveness of the present model can be validated.predicted by FSDT is larger than that by predicted by HSDT in the present work, and thus the frequencies calculated based on FSDT will be larger.The comparison of the results shows that the contact force obtained through the present method is close to the results in Ref. [55] during the loading stage and in Refs.[54,56].during the unloading stage.In the early stage, there is minimal discrepancy in the displacements of the impactor between these results.However, as the unloading stage progresses, deviations gradually increase due to the accumulation of the differences in the contact force over time.Overall, these comparisons indicate that the solutions obtained by the present method are in excellent agreement with those reported in the previous literature.The slight differences observed can be attributed to the use of different theories, such as the first-order shear deformation plate theory in Ref. [56] and the classical plate theory in Ref. [54], as well as variations in numerical solution procedures.
In the Galerkin method, the number of terms of the spatial function is finite and must be truncated.With an increase in the number of selected spatial function terms, the calculation results gradually converge to the exact solution, but meanwhile, the amount of calculation also increases.Therefore, on the premise of ensuring the convergence of the results, the appropriate number of truncated items should be selected to achieve a balance between computational accuracy and computational efficiency.Figure 6 presents the convergence analysis of the Galerkin method used in the present study.As the truncation The comparison of the results shows that the contact force obtained through the present method is close to the results in Ref. [55] during the loading stage and in Refs.[54,56].during the unloading stage.In the early stage, there is minimal discrepancy in the displacements of the impactor between these results.However, as the unloading stage progresses, deviations gradually increase due to the accumulation of the differences in the contact force over time.Overall, these comparisons indicate that the solutions obtained by the present method are in excellent agreement with those reported in the previous literature.The slight differences observed can be attributed to the use of different theories, such as the first-order shear deformation plate theory in Ref. [56] and the classical plate theory in Ref. [54], as well as variations in numerical solution procedures.
In the Galerkin method, the number of terms of the spatial function is finite and must be truncated.With an increase in the number of selected spatial function terms, the calculation results gradually converge to the exact solution, but meanwhile, the amount of calculation also increases.Therefore, on the premise of ensuring the convergence of the results, the appropriate number of truncated items should be selected to achieve a balance between computational accuracy and computational efficiency.Figure 6 presents the convergence analysis of the Galerkin method used in the present study.As the truncation numbers (m and n) increase, the difference in calculation results decreases.It can be observed that convergence is achieved when the truncation numbers reach m = 11 and n = 11.Therefore, for the remainder of the paper, we selected these two truncated coefficients to ensure accurate and reliable results.
contact force over time.Overall, these comparisons indicate that the solutions obtained by the present method are in excellent agreement with those reported in the previous literature.The slight differences observed can be attributed to the use of different theories, such as the first-order shear deformation plate theory in Ref. [56] and the classical plate theory in Ref. [54], as well as variations in numerical solution procedures.
In the Galerkin method, the number of terms of the spatial function is finite and must be truncated.With an increase in the number of selected spatial function terms, the calculation results gradually converge to the exact solution, but meanwhile, the amount of calculation also increases.Therefore, on the premise of ensuring the convergence of the results, the appropriate number of truncated items should be selected to achieve a balance between computational accuracy and computational efficiency.Figure 6 presents the convergence analysis of the Galerkin method used in the present study.As the truncation numbers (m and n) increase, the difference in calculation results decreases.It can be observed that convergence is achieved when the truncation numbers reach m = 11 and n = 11.Therefore, for the remainder of the paper, we selected these two truncated coefficients to ensure accurate and reliable results.

Parametric Analysis
In this subsection, a detailed investigation is performed on several parameters that influence the response of the sandwich shell during impact.These parameters include the temperature, moisture, power law index, thickness of the honeycomb wall, and initial velocity of the impactor.Each parameter is studied thoroughly to understand its influence on the behavior of the sandwich shell under impact conditions.The effects of varying these parameters are analyzed and discussed in depth to provide a comprehensive understanding of their role in determining the impact responses of the sandwich shell.Unless otherwise stated, the geometric parameter of the sandwich shell is L x = L y = 0.8 m and the core-to-face layer thickness ratio is h c /h f = 10.The thickness ratio of the multilayer honeycomb is 1:8:1 and the geometric parameters of the honeycomb core are η 1 = 3, η 3 = 0.1, θ = −30 • .The radius and initial velocity of the impactor are R 0 = 15 mm and V 0 = 5 m/s.The parameters of the Si 3 N 4 /SUS304 functionally graded material are listed in Table 1.Its material properties are determined through references [14,[57][58][59]: The material of the face layers and honeycomb are the same: The material properties of VEM are as follows: The material properties of the impactor are:

Effect of Temperature
Figure 7a shows the changes in the natural frequencies of the composite sandwich shell with the temperature and the inner angle of the honeycomb under the uniform distribution pattern of temperature.The natural frequency of the structure decreases with an increase of temperature and the inclined angle.The main reason is that the elastic modulus of the material will decrease with an increase of temperature.Within the range of inclined angle illustrated in Figure 7a, the equivalent moduli of the honeycomb structure decreased with the increase of the inclined angle, resulting in a decrease of the overall stiffness of the structure, and the natural frequency decreased accordingly.Figure 8 illustrates the effect of temperature variations on the contact force and displacement of sandwich shells at the contact point, considering three different temperature distribution patterns for ΔT = 200 K and ΔT = 0 K.The peak value of the contact force will decrease with a temperature rise, but there are no significant differences among the contact force curves under the three temperature distribution patterns.In view of the characteristics of the temperature distribution in this study, the temperature at the top surface of the structure remains constant under the three distribution patterns, and the mechanical properties of both the outermost surface of the sandwich shell and the impactor are also the same, resulting in an unchanging contact stiffness during the low-velocity impact.Hence, the impact process shows minimal variation in both the contact force and contact duration time across the different temperature distribution patterns.Although there is no significant difference in the contact force, the displacement of the sandwich shell at the impact point exhibits noticeable variations.As the temperature increases, the displacement at the impact point increases dramatically and also shows distinct differences under different temperature distribution patterns.The maximum displacement amplitude occurs when the temperature is uniformly distributed, followed by Figure 7b shows the differences in the natural frequencies with the temperature under three temperature distribution patterns when the inner angle of the honeycomb is θ = −30 • .This figure indicates that the natural frequency of the structure experiences the most significant decrease under a uniform distribution of temperature and the least under a sinusoidal distribution of temperature.Combined with the expressions of the different temperature distribution patterns, the temperature change at each spatial position along the thickness direction of the structure is ∆T in the case of a uniform distribution.However, for the other two distribution patterns, the temperature changes along the thickness direction from 0 to ∆T in a certain gradient, and the total temperature change is small, resulting in less influence on the natural frequency of the structure.
Figure 8 illustrates the effect of temperature variations on the contact force and displacement of sandwich shells at the contact point, considering three different temperature distribution patterns for ∆T = 200 K and ∆T = 0 K.The peak value of the contact force will decrease with a temperature rise, but there are no significant differences among the contact force curves under the three temperature distribution patterns.In view of the characteristics of the temperature distribution in this study, the temperature at the top surface of the structure remains constant under the three distribution patterns, and the mechanical properties of both the outermost surface of the sandwich shell and the impactor are also the same, resulting in an unchanging contact stiffness during the low-velocity impact.Hence, the impact process shows minimal variation in both the contact force and contact duration time across the different temperature distribution patterns.
Although there is no significant difference in the contact force, the displacement of the sandwich shell at the impact point exhibits noticeable variations.As the temperature increases, the displacement at the impact point increases dramatically and also shows distinct differences under different temperature distribution patterns.The maximum displacement amplitude occurs when the temperature is uniformly distributed, followed by a linear distribution and a minimum when the temperature shows a sinusoidal distribution.In order to further describe the attenuation characteristics and decay time of the structure during impact, the peak point of the displacement time curve for the wave at the indentation point can be fitted by an exponential function in terms of [60]: Figure 8 illustrates the effect of temperature variations on the contact force and displacement of sandwich shells at the contact point, considering three different temperature distribution patterns for ΔT = 200 K and ΔT = 0 K.The peak value of the contact force will decrease with a temperature rise, but there are no significant differences among the contact force curves under the three temperature distribution patterns.In view of the characteristics of the temperature distribution in this study, the temperature at the top surface of the structure remains constant under the three distribution patterns, and the mechanical properties of both the outermost surface of the sandwich shell and the impactor are also the same, resulting in an unchanging contact stiffness during the low-velocity impact.Hence, the impact process shows minimal variation in both the contact force and contact duration time across the different temperature distribution patterns.Although there is no significant difference in the contact force, the displacement of the sandwich shell at the impact point exhibits noticeable variations.As the temperature increases, the displacement at the impact point increases dramatically and also shows distinct differences under different temperature distribution patterns.The maximum displacement amplitude occurs when the temperature is uniformly distributed, followed by a linear distribution and a minimum when the temperature shows a sinusoidal distribution.In order to further describe the attenuation characteristics and decay time of the structure during impact, the peak point of the displacement time curve for the wave at the indentation point can be fitted by an exponential function in terms of [60]: Then, the decay time is denoted as t 1 , representing the time required for the displacement amplitude to decrease to 1/e (approximately 36.8%) of its maximum value [60].This value can be determined by utilizing a fitted function that accurately describes the decay behavior of the displacement amplitude over time.By analyzing the fitted function, the specific time t 1 can be obtained, providing valuable insight into the rate of decay and the overall behavior of the system.
As shown from the data in Table 6, the displacement response at the impact point of the structure is nearly identical under both linear temperature distribution and sinusoidal distribution conditions, with only a 3.65% difference in displacement amplitude and a 2.41% difference in decay time.However, it is noteworthy that a uniform temperature distribution leads to a more severe deterioration of damping performance of the structure.The displacement amplitude increases by 22.84% compared with that under the condition of a linear distribution, and the decay time is prolonged by 19.83%.When moisture is uniformly distributed along the thickness direction, the natural frequency of the composite sandwich shell varies with moisture and the inclined angle of the honeycomb, as shown in Figure 9a.With an increase of the inclined angle, the natural frequency of the sandwich shell gradually decreases, especially in the range of −55 • to −80 • , where an increase in the inclined angle will lead to a more significant reduction in the natural frequency.Figure 10 illustrates the impact responses of a sandwich shell under different moisture distribution patterns, taking into account moisture variations of ΔH = 20% and ΔH = 0%.The figure depicts that the influence of moisture on the contact force and displacement at the impact point is minimal.However, there is a slight reduction in displacement at the impact point as the moisture content increases.This study examined three moisture distribution patterns, namely a uniform distribution, linear distribution, and sinusoidal distribution.The results indicate that the uniform distribution pattern leads to the greatest decrease in displacement, followed by the linear distribution pattern, while the sinusoidal distribution pattern causes the least decrease.This behavior is consistent with the influence of temperature, suggesting a similar relationship between these two parameters.The results of function fitting on the displacement time curves of the structure at the impact point under different moisture distributions are listed in Table 7.Compared with the case of reference moisture, when ΔH = 20% and the moisture is uniformly distributed, the displacement amplitude decreases the most, reaching 9.23%, and the decay time extends 35.62% correspondingly.It can be seen that although the moisture change has no significant effect on the stiffness of the structure, it will significantly affect the damping performance of the structure.Figure 9b illustrates the changes trend of the natural frequency of a structure with moisture in three moisture distribution patterns when the inner angle θ = −30 • .With an increase of moisture, the natural frequency of this composite structure slightly increases.On the one hand, the mechanical parameters of the viscoelastic polymer materials in the structure are sensitive to moisture changes; on the other hand, only the middle layer of the negative Poisson's ratio honeycomb is made of viscoelastic polymer material, and its composition proportion in the structure is not high, and thus it has little influence on the performance of the structure.With an increase of moisture, the elastic modulus of the viscoelastic polymer materials will decrease slightly, and the stiffness of the structure should have decreased.However, due to the negative Poisson's ratio effect of the honeycomb core, the tensile stiffness and bending stiffness of the structure will change, and finally the natural frequency will increase slightly.
Figure 10 illustrates the impact responses of a sandwich shell under different moisture distribution patterns, taking into account moisture variations of ∆H = 20% and ∆H = 0%.The figure depicts that the influence of moisture on the contact force and displacement at the impact point is minimal.However, there is a slight reduction in displacement at the impact point as the moisture content increases.This study examined three moisture distribution patterns, namely a uniform distribution, linear distribution, and sinusoidal distribution.The results indicate that the uniform distribution pattern leads to the greatest decrease in displacement, followed by the linear distribution pattern, while the sinusoidal distribution pattern causes the least decrease.This behavior is consistent with the influence of temperature, suggesting a similar relationship between these two parameters.Figure 10 illustrates the impact responses of a sandwich shell under different moisture distribution patterns, taking into account moisture variations of ΔH = 20% and ΔH = 0%.The figure depicts that the influence of moisture on the contact force and displacement at the impact point is minimal.However, there is a slight reduction in displacement at the impact point as the moisture content increases.This study examined three moisture distribution patterns, namely a uniform distribution, linear distribution, and sinusoidal distribution.The results indicate that the uniform distribution pattern leads to the greatest decrease in displacement, followed by the linear distribution pattern, while the sinusoidal distribution pattern causes the least decrease.This behavior is consistent with the influence of temperature, suggesting a similar relationship between these two parameters.The results of function fitting on the displacement time curves of the structure at the impact point under different moisture distributions are listed in Table 7.Compared with the case of reference moisture, when ΔH = 20% and the moisture is uniformly distributed, the displacement amplitude decreases the most, reaching 9.23%, and the decay time extends 35.62% correspondingly.It can be seen that although the moisture change has no significant effect on the stiffness of the structure, it will significantly affect the damping The results of function fitting on the displacement time curves of the structure at the impact point under different moisture distributions are listed in Table 7.Compared with the case of reference moisture, when ∆H = 20% and the moisture is uniformly distributed, the displacement amplitude decreases the most, reaching 9.23%, and the decay time extends 35.62% correspondingly.It can be seen that although the moisture change has no significant effect on the stiffness of the structure, it will significantly affect the damping performance of the structure.11 illustrates how the power index of the functionally graded face-sheets affects the contact force and displacement of a sandwich shell at the impact point.While the contact stiffness remains constant for different power law indices during the impact, variations in the structural stiffness and damping performance lead to different levels of indentation on the structure's surface, resulting in distinct peak values of the contact force.With an increase in the power law index, there is a slight elevation in the peak contact force, but the contact duration remains unaltered throughout the impact process.11 illustrates how the power index of the functionally graded face-sheets affects the contact force and displacement of a sandwich shell at the impact point.While the contact stiffness remains constant for different power law indices during the impact, variations in the structural stiffness and damping performance lead to different levels of indentation on the structure's surface, resulting in distinct peak values of the contact force.With an increase in the power law index, there is a slight elevation in the peak contact force, but the contact duration remains unaltered throughout the impact process.With an increase of the power law index, the volume fraction of ceramic material at the same position on the face-sheets will decrease, resulting in a reduction in both the stiffness and mass of the face-sheets.Consequently, the natural frequency of the structure will also decrease.As the stiffness of the core layer is lower than that of the face-sheets, the deformation primarily occurs within the core layer.Therefore, the displacement of the structure at the impact point will decrease when the stiffness of the face-sheets decreases.
Table 8 shows the results of function fitting of the displacement time curve for different power law indexes.As evidenced by data in the table, with an increase of the power law index, the displacement amplitude of the impact point decreases, and the fitting index B also decreases, resulting in prolongation of the vibration attenuation.It can be inferred that increasing the power law index will mitigate the deformation of the structure while slowing down the vibration attenuation.Hence, the selection of an appropriate power law index in structural design is crucial, as it not only prevents excessive deformation of the structure but also enhances its damping performance.This facilitates rapid dissipation of the vibration energy and minimizes the generation of sound radiation noise caused by vibration.With an increase of the power law index, the volume fraction of ceramic material at the same position on the face-sheets will decrease, resulting in a reduction in both the stiffness and mass of the face-sheets.Consequently, the natural frequency of the structure will also decrease.As the stiffness of the core layer is lower than that of the face-sheets, the deformation primarily occurs within the core layer.Therefore, the displacement of the structure at the impact point will decrease when the stiffness of the face-sheets decreases.
Table 8 shows the results of function fitting of the displacement time curve for different power law indexes.As evidenced by data in the table, with an increase of the power law index, the displacement amplitude of the impact point decreases, and the fitting index B also decreases, resulting in prolongation of the vibration attenuation.It can be inferred that increasing the power law index will mitigate the deformation of the structure while slowing down the vibration attenuation.Hence, the selection of an appropriate power law index in structural design is crucial, as it not only prevents excessive deformation of the structure but also enhances its damping performance.This facilitates rapid dissipation of the vibration energy and minimizes the generation of sound radiation noise caused by vibration.12 present the influence of the honeycomb wall thickness of the core layer on the contact force and displacement at the impact point under a low-velocity impact.When the thickness of the honeycomb walls increases, the equivalent elastic moduli of the core layer will increase, resulting in an increase in stiffness of the sandwich structure, and then the resistance to deformation of the structure will be enhanced.Additionally, there is a slight increase in the peak contact force and a decrease in the displacement amplitude at the impact point under a low-velocity impact.Figure 12 present the influence of the honeycomb wall thickness of the core layer on the contact force and displacement at the impact point under a low-velocity impact.When the thickness of the honeycomb walls increases, the equivalent elastic moduli of the core layer will increase, resulting in an increase in stiffness of the sandwich structure, and then the resistance to deformation of the structure will be enhanced.Additionally, there is a slight increase in the peak contact force and a decrease in the displacement amplitude at the impact point under a low-velocity impact.With the increase of honeycomb wall thickness, the volume fraction of the viscoelastic polymer material increases, the damping performance of sandwich shell structure is enhanced, and then the vibration attenuation time is also shortened, which is consistent with the fitting data in Table 9.Therefore, the stiffness and damping performance of the structure can be improved by increasing the thickness of the honeycomb wall, and the vibration of the structure can be effectively suppressed.13 depicts the time history curves of the contact force and displacement of the sandwich shell at the impact point for different initial velocities of the impactor.With an increase of the initial velocity of the impactor, the impact energy acting on the structure increases, leading to a significant increase in the contact force and a decrease in the contact duration time.When the initial velocity of the impactor is doubled, the peak value of the contact force increases by approximately 2.3 times.With the increase of honeycomb wall thickness, the volume fraction of the viscoelastic polymer material increases, the damping performance of sandwich shell structure is enhanced, and then the vibration attenuation time is also shortened, which is consistent with the fitting data in Table 9.Therefore, the stiffness and damping performance of the structure can be improved by increasing the thickness of the honeycomb wall, and the vibration of the structure can be effectively suppressed.An increase in the initial velocity of the impactor results in a significant rise in the displacement amplitude of the structure at the impact point.Table 10 provides function fitting results for different initial velocities of the impactor, showing that while the structure remains unchanged, there is a variation in the displacement amplitude.Interestingly, the fitting parameter B remains nearly identical across different velocities, indicating a consistent relationship between the impactor's initial velocity and the displacement amplitude.Furthermore, the vibration attenuation time remains constant, suggesting that it is independent of the initial velocity of the impactor.

Conclusions
This study focused on the design of a new functionally graded composite sandwich doubly curved shell with a multilayer negative Poisson's ratio honeycomb core.The dynamic responses of this structure when subjected to low-velocity impacts in hygrothermal environments were investigated.The dynamic equilibrium equations for the structure under different temperature and moisture distribution patterns were derived using higherorder shear deformation theory and Hamilton's principle.To numerically solve the dynamic problem, the Galerkin method and the Newmark direct integration scheme were employed to discretize the differential dynamic equations in the space and time domains, respectively.Subsequently, a parametric study was conducted to analyze the effects of various factors on the contact force and the displacement at the impact point of the structure.Based on this analysis, the following five conclusions can be obtained: (1) Increasing temperature will reduce the contact force during the impact process, and then the damping performance of the structure will deteriorate, resulting in an increase of displacement at the impact point and prolongation of the vibration decay time.There was no significant difference between the contact force under three different temperature distribution patterns, but the damping performance of the structure decreased most seriously when the temperature was uniformly distributed throughout the thickness.An increase in the initial velocity of the impactor results in a significant rise in the displacement amplitude of the structure at the impact point.Table 10 provides function fitting results for different initial velocities of the impactor, showing that while the structure remains unchanged, there is a variation in the displacement amplitude.Interestingly, the fitting parameter B remains nearly identical across different velocities, indicating a consistent relationship between the impactor's initial velocity and the displacement amplitude.Furthermore, the vibration attenuation time remains constant, suggesting that it is independent of the initial velocity of the impactor.

Conclusions
This study focused on the design of a new functionally graded composite sandwich doubly curved shell with a multilayer negative Poisson's ratio honeycomb core.The dynamic responses of this structure when subjected to low-velocity impacts in hygrothermal environments were investigated.The dynamic equilibrium equations for the structure under different temperature and moisture distribution patterns were derived using higherorder shear deformation theory and Hamilton's principle.To numerically solve the dynamic problem, the Galerkin method and the Newmark direct integration scheme were employed to discretize the differential dynamic equations in the space and time domains, respectively.Subsequently, a parametric study was conducted to analyze the effects of various factors on the contact force and the displacement at the impact point of the structure.Based on this analysis, the following five conclusions can be obtained: (1) Increasing temperature will reduce the contact force during the impact process, and then the damping performance of the structure will deteriorate, resulting in an increase of displacement at the impact point and prolongation of the vibration decay time.
There was no significant difference between the contact force under three different temperature distribution patterns, but the damping performance of the structure

Figure 1 .
Figure 1.Schematic diagram of the composite sandwich shell and its constituents: (a) the structure under a low-velocity impact, (b) three layers of the structure; and (c) constituents of the multilayer honeycomb core.

Figure 1 .
Figure 1.Schematic diagram of the composite sandwich shell and its constituents: (a) the structure under a low-velocity impact, (b) three layers of the structure; and (c) constituents of the multilayer honeycomb core.

Figure 2 .
Figure 2. The volume fraction of ceramic in face-sheets for different values of the power law in p.

Figure 2 .
Figure 2. The volume fraction of ceramic in face-sheets for different values of the power law index p.

Figure 3 .
Figure 3.The single-layer equivalence for the multilayer negative Poisson's ratio honeycomb wi viscoelastic polymer material.

Figure 3 .
Figure 3.The single-layer equivalence for the multilayer negative Poisson's ratio honeycomb with viscoelastic polymer material.

Figure 4 .
Figure 4.A schematic of the indentation process.

Figure 4 .
Figure 4.A schematic of the indentation process.

Case 3 :
In previous literature, there is a lack of available results regarding the impact responses of a composite sandwich shell with a multilayer negative Poisson's ratio honeycomb core.Therefore, in order to verify the accuracy and effectiveness of the current method, a degenerated model is presented.The example focuses on the low-velocity impact problem of a simply supported square isotropic plate and is discussed in detail in the references cited[53,54].The isotropic plate in consideration has dimensional parameters of a × b × h = 200 × 200 × 8 mm.It experiences an impact load from a spherical impactor with a radius of 10 mm and an initial velocity of 1 m/s.Both the plate and impactor share the same material properties, including a mass density of ρ = 7810 kg/m −3 , Young's modulus of E = 206.8GPa, and Poisson's ratio of ν = 0.3.Figure 5 presents a comparison of the time history of the contact force and displacement of the impactor throughout the impact duration.Materials 2024, 17, 233 18 of 30

Figure 5 .
Figure 5.Comparison of the low-velocity impact responses of a simply supported isotropic square plate: (a) contact force and (b) displacement of the impactor [54-56].

Figure 5 .
Figure 5.Comparison of the low-velocity impact responses of a simply supported isotropic square plate: (a) contact force and (b) displacement of the impactor [54-56].

Figure 6 .
Figure 6.Galerkin method convergence investigation: (a) contact force and (b) displacement at the impact point.

Figure 6 .
Figure 6.Galerkin method convergence investigation: (a) contact force and (b) displacement at the impact point.

Materials 2024, 17 , 233 20 of 30 Figure 7 .
Figure 7.The natural frequencies of the composite sandwich shell: (a) variation with the inclined angle and temperature and (b) variation with temperature for θ = −30° under three different distribution patterns.

Figure 8 .
Figure 8. Influence of different temperature distribution patterns on the impact responses of the composite structure: (a) contact force, and (b) displacements at the impact point.

Figure 7 .
Figure 7.The natural frequencies of the composite sandwich shell: (a) variation with the inclined angle and temperature and (b) variation with temperature for θ = −30 • under three different distribution patterns.

Figure 8 .
Figure 8. Influence of different temperature distribution patterns on the impact responses of the composite structure: (a) contact force, and (b) displacements at the impact point.

Figure 8 .
Figure 8. Influence of different temperature distribution patterns on the impact responses of the composite structure: (a) contact force, and (b) displacements at the impact point.

Figure 9 .
Figure 9.The natural frequencies of the composite sandwich shell: (a) variation with the inclined angle and moisture and (b) variation with moisture for θ = −30° under three different distribution patterns.

Figure 10 .
Figure 10.Effect of different moisture distribution patterns on the impact responses of the structure: (a) contact force and (b) displacement of the structure at the impact point.

Figure 9 .
Figure 9.The natural frequencies of the composite sandwich shell: (a) variation with the inclined angle and moisture and (b) variation with moisture for θ = −30 • under three different distribution patterns.

Figure 9 .
Figure 9.The natural frequencies of the composite sandwich shell: (a) variation with the inclined angle and moisture and (b) variation with moisture for θ = −30° under three different distribution patterns.

Figure 10 .
Figure 10.Effect of different moisture distribution patterns on the impact responses of the structure: (a) contact force and (b) displacement of the structure at the impact point.

Figure 10 .
Figure 10.Effect of different moisture distribution patterns on the impact responses of the structure: (a) contact force and (b) displacement of the structure at the impact point.

Figure 11 .
Figure 11.Effect of different power law indexes on the dynamics responses of the structure: (a) contact force and (b) displacement at the impact point.

Figure 11 .
Figure 11.Effect of different power law indexes on the dynamics responses of the structure: (a) contact force and (b) displacement at the impact point.

Figure 12 .
Figure 12.Effect of different thicknesses of the honeycomb wall on responses of the structure: (a) contact force and (b) displacement at the impact point.

Figure 12 .
Figure 12.Effect of different thicknesses of the honeycomb wall on responses of the structure: (a) contact force and (b) displacement at the impact point.

Figure 13 30 Figure 13 .
Figure13depicts the time history curves of the contact force and displacement of the sandwich shell at the impact point for different initial velocities of the impactor.With an increase of the initial velocity of the impactor, the impact energy acting on the structure

Figure 13 .
Figure 13.Effect of different initial velocities of the impactor on the responses of the structure: (a) contact force and (b) displacement at the impact point.

Table 5 .
Comparison of dimensionless natural frequencies of Al/Al 2 O 3 functionally graded structures with a negative Poisson's ratio honeycomb core.

Table 6 .
Exponential function parameters and decay times are determined for various temperature distribution patterns.

Table 7 .
Exponential function parameters and decay times obtained through curve fitting.

Table 7 .
Exponential function parameters and decay times obtained through curve fitting.

Table 8 .
Fitted exponential function parameters and decay times for different power law indexes.

Table 8 .
Fitted exponential function parameters and decay times for different power law indexes.

Table 9 .
Fitted exponential function parameters and decay times for different thicknesses of the honeycomb wall.

Table 9 .
Fitted exponential function parameters and decay times for different thicknesses of the honeycomb wall.

Table 10 .
Fitted exponential function parameters and decay times for different initial velocities of the impactor.

Table 10 .
Fitted exponential function parameters and decay times for different initial velocities of the impactor.