On Mechanical Properties of Polymer Composite Simulated Ice Structures for Flight Testing

1 Introduction

When an aircraft flies at subsonic speeds within cloud layers below 8000 m, ice may accrete on the leading edge of the wing. Ice formation is primarily caused by supercooled droplets, whose diameters typically range from 5 to 50 μm. Upon impact with the aircraft surface, these droplets freeze, leading to ice deposits [1,2,3]. The accumulation of ice on the wing surface adversely affects the aircraft’s aerodynamic characteristics, increasing drag coefficient and reducing lift coefficient, thereby posing a significant threat to flight safety [4]. Consequently, it is essential to assess aircraft flight performance under icing conditions and to demonstrate safe operation under both continuous and intermittent maximum icing conditions.

Flight testing is a key component of compliance verification for icing conditions. Its main objective is to confirm that the aircraft can maintain adequate flight quality after ice accretion and that the degraded performance remains within acceptable limits. Icing certification flight tests are generally classified into three categories. The first category is dry-air flight tests with anti-icing systems activated, which verify the normal operation of anti-icing equipment without compromising flight quality under clear atmospheric conditions and establish the operational envelope of anti-icing systems. The second category refers to dry-air flight tests with simulated ice shapes. Ice replicas are mounted on the aircraft, and tests are conducted under stable dry-air conditions to assess aerodynamic performance and flight handling characteristics. This approach greatly reduces the reliance on flight tests conducted in unpredictable natural icing environments. The third category is natural icing flight tests, which validate the actual location and geometric features of in-flight ice accretion. Such tests also confirm that the performance degradation trends caused by natural ice are consistent with those induced by simulated ice.

Among all test items, dry-air flight tests with simulated ice shapes play a decisive role in the entire airworthiness verification framework. The implementation workflow consists of four sequential steps: ice shape prediction based on numerical simulations and icing wind tunnel tests, aerodynamic assessment of predicted ice configurations, onboard installation of ice specimens while preserving original geometric features, and full-scale flight validation. To understand the development of such simulated ice testing methods, it is instructive to review the historical and technical evolution of icing research.

In 1978, NASA and the FAA held a seminar at the Lewis Research Center, restarting extensive icing research. With the rapid development of computational fluid dynamics (CFD), the focus shifted toward numerical simulation of icing phenomena and experimental validation. Potapczuk and Gerhart [5] used the Reynolds-averaged Navier–Stokes (RANS) method to numerically compute the two-dimensional (2D) aerodynamic characteristics of two ice shapes. Kwon and Sankar [6] were the first to perform numerical analysis of the aerodynamic characteristics of a three-dimensional (3D) iced NACA 0012 airfoil at 4° and 8° angles of attack (AOAs). Good agreement between numerical predictions and experimental data was achieved at 4° AOA, whereas obvious discrepancies emerged near the leading edge at 8° AOA due to severe local flow separation. Lee and Bragg [7] adopted quarter-circle profiles to simulate horn ice on a NACA 23012 airfoil and performed parametric tests on ice height and chordwise position. The results revealed that ice located at 12% chord caused the most severe lift loss, and this conclusion was also applicable to the NLF 0414 airfoil. Papadakis et al. [8] reported that when the ice angle was perpendicular to the airfoil chord line, the aerodynamic degradation was the most severe.

To improve icing prediction capability, a series of professional icing simulation software has been developed worldwide, including LEWICE [9] from the United States, ONERA [10] from France and FENSAP-ICE [11] from Canada. Existing studies have confirmed that flow separation bubbles downstream of horn ice were the primary cause of aerodynamic performance deterioration [12]. Upper-surface horn ice mainly affected lift coefficient and stall angle, while lower-surface horn ice dominates performance degradation at low and negative AOAs [13]. Moreover, ice surface roughness had a limited impact on overall aerodynamic characteristics, as verified by multiple published works [14].

At present, research on icing has primarily focused on developing more accurate icing prediction models and aerodynamic performance calculation methods. Min and Yee [15] introduced leading-edge surface roughness and convective coefficients to compute four ice shapes. The results outperformed LEWICE predictions and resolved the issue of excessive ice-horn growth. Suo et al. [16] proposed an icing prediction model based on a geometry-constrained neural network. With flight velocity, ambient temperature, liquid water content, median volumetric diameter, and icing time as inputs and ice thickness as the output, the model incorporated geometric constraints into the loss function to improve prediction accuracy. The results showed that the geometric constraints effectively reduced fluctuations in ice shape prediction and enhance accuracy. The trained model could be used for rapid assessment of wing icing conditions. Zhao et al. [17] established a deep operator network and conducted training and testing on the NACA 0012 iced airfoil. The predicted results were compared with those of the traditional direct CNN network. Unlike the traditional CNN, which can only predict aerodynamic coefficients under fixed flow conditions consistent with the training data, this network can flexibly predict aerodynamic coefficients under different flow conditions. Bornhoft et al. [18] performed aerodynamic analyses on realistic ice shapes obtained by laser scanning and on simplified ice shapes. The results indicated that realistic ice shapes yielded better computational accuracy and that introducing roughness features on simplified ice shapes was necessary.

While the above studies have advanced ice prediction and aerodynamic performance calculation capabilities, there is relatively little content regarding the actual simulation of ice type installation and flight test verification. Flemming et al. [19] defined the critical icing conditions for the S-92A helicopter to meet airworthiness requirements. Through CFD calculations, the critical ice shapes for various components and the associated aerodynamic data were obtained. For the vertical and horizontal tails, the simulated ice shapes were fabricated by wire-cut foam core wrapped with glass fiber and installed on the designated surfaces. The outer layer of the ice shapes was coated with 40-grit aluminum oxide, while 40-grit sandpaper was applied to the main and tail rotors. Flight tests showed that the ice shapes on the tail had a relatively minor impact on aerodynamic characteristics, whereas rotor roughness caused a 70% loss of power. Ratvasky et al. [20] conducted flight tests with simulated ice shapes on a commercial jet aircraft. Three ice shapes were used in the tests. The actual flight data agreed well with the simulation results; however, a stall dive tendency was observed when the flaps were extended. Deiler et al. [21] performed flight tests on an aircraft after installing a supercooled large droplet (SLD) ice shape. The results indicated a significant degradation in aircraft performance after ice accretion, but the SLD ice shape was not more critical compared with other ice shapes. In 2016, Boeing filed a patent for the installation of simulated ice shapes applicable to actual airworthiness certification processes [22]. The patent employed 3D-printed simulated ice shapes attached to the designated surface using double-sided tape, and linear components were used to remove the simulated ice shapes. Li et al. [23] performed CFD simulations on the icing configuration of a civil aircraft wing and proposed a simulated ice shape structure. The mechanical properties of the structure were obtained through finite element analysis, verifying its feasibility for actual flight test applications.

To fill the above-mentioned gaps, this paper presents an integrated study that couples aerodynamic analysis, structural design, finite element simulation, and experimental validation for a simulated ice shape intended for dry-air flight tests. Unlike most existing works that focus solely on ice prediction or aerodynamic degradation, the present work systematically quantifies the aerodynamic penalty of a realistic glaze ice shape on a civil aircraft wing at climb/descent and cruise phases, and extracts the corresponding surface loads as design inputs. A composite sandwich structure with a 3D-printed hollow core and extended GFRP bonding tabs is proposed, and its mechanical response under the peak aerodynamic load is evaluated using cohesive-zone and Tsai-Wu failure criteria. The predicted failure mode and load are validated through quasi-static tests on a 250 mm hand-lay-up specimen. This paper is based on the necessity of simulated ice for airworthiness certification. It adopts mature processing techniques and material systems for structural design, focusing on the mechanical properties of the designed structure. This provides a verifiable method for the engineering design and airworthiness certification of the simulated ice, thereby directly supporting airworthiness flight tests.

2 Aerodynamic Analysis of the Iced Airfoil

2.1 Computational Model

This study focuses on a civil aircraft airfoil with a chord length of 2.36 m. Based on experimental data [24] from the CIRA-IWT icing wind tunnel, the ice shape shown in Fig. 1 was obtained. The x-direction and y-direction represent the freestream direction and the wing normal direction, respectively.

Figure 1: Ice shape on the airfoil section.

The aerodynamic performance of this ice shape is first analyzed to evaluate its impact on the clean wing performance. Aerodynamic loads on the ice surface are then extracted to provide input for subsequent mechanical analyses.

Considering the operational scenarios of the civil aircraft, the climb/descent phase and the cruise phases are selected for computation. The k-ω SST turbulence model [25] is employed, which combines the standard k-ε [26] and k-ω [27] models and offers high accuracy and reliability in flow field simulations. Compared to other turbulence models, the k-ω SST turbulence model is more effective in simulating turbulence at high angles of attack and offers higher accuracy in aerodynamic performance evaluation [28]. The computational domain size is 100 m, greater than 25 times the chord length, satisfying the pressure far-field boundary condition. The specific computational inputs are listed in Table 1. Airfoil deformation is neglected, and the no-slip wall condition is considered.

2.2 Grid Independence Verification

Grid generation is completed in ICEM software with a O-block topology scheme. As shown in Fig. 2, Mesh refinement is implemented locally around ice horns, where the first-layer grid size on the airfoil was 0.003 mm to achieve a y+ value of 1, with a mesh growth rate of 1.1. The same O-block meshing strategy was applied to the clean airfoil for comparison.

Figure 2: Grid of airfoil with ice shape: (a) total grid; (b) refined grid of the leading edge.

CFD predictions are sensitive to grid density, necessitating a grid convergence verification. Following the International Towing Tank Conference (ITTC) recommended specifications, three grid sets are refined proportionally with a scaling factor of 2 under unchanged mesh topology, and verification is conducted on the iced airfoil at 0.2 Ma and 0° angle of attack (AOA). Table 2 summarizes the computed lift coefficient (CL) and drag coefficient (CD). The second grid provides sufficient accuracy while saving computational resources; therefore, it was used for all subsequent computations.

2.3 Comparison of Lift and Drag Coefficients

Leading-edge ice alters the original airfoil geometry and triggers premature boundary-layer separation, preventing full flow attachment over the airfoil surface and inducing prominent aerodynamic degradation.

Fig. 3 shows the lift and drag coefficients versus AOA for the clean and iced airfoils at 0.2 Ma and 0.6 Ma. For the clean airfoil, an increase in Mach number raises the lift coefficient at a given AOA but also causes an earlier stall angle. Higher incoming flow velocity accelerates suction-side flow and intensifies near-surface low-pressure regions to improve lift generation; meanwhile, compressibility-induced shock waves emerge on the upper airfoil, accelerating boundary-layer detachment and reducing stall angle.

Conversely, geometric distortion from leading-edge ice disrupts flow acceleration and suppresses shock formation, further facilitating flow separation. Consequently, the iced airfoil suffers substantial lift coefficient loss and pronounced drag coefficient increment. Compared with the clean airfoil, the maximum lift coefficient of the iced airfoil decreases from 1.57 to 0.65 at 0.2 Ma, and the stall angle reduces from 13° to 5.5°. At 0.6 Ma, the maximum lift coefficient drops from 1.08 to 0.59, and the stall angle reduces from 6.5° to 4.5°. Such dramatic aerodynamic degradation quantitatively demonstrates the detrimental influence of leading-edge icing.

Figure 3: Variations of aerodynamic coefficients at different AOA and Ma: (a) C_L vs. AoA; (b) C_D vs. AoA.

2.4 Flow Field Distributions

To understand the underlying flow mechanisms causing these aerodynamic penalties, the pressure coefficient distributions and velocity fields are examined. Fig. 4 presents the pressure coefficient (CP) distributions for the clean and iced airfoils at different AOA under 0.2 Ma. The presence of the leading-edge ice shape significantly affects the pressure coefficient distribution, with more dramatic pressure variations on the leeward side of the ice horn. At 0.6 Ma, the shock wave at the leading edge disappears due to the influence of the ice horn on the flow characteristics, and a pronounced adverse pressure gradient forms (Fig. 5).

Figure 4: Pressure coefficient distributions at 0.2 Ma.

Figure 5: Pressure coefficient distributions at 0.6 Ma.

Fig. 6 illustrates velocity distributions and near-wall streamline distributions at different AOA under 0.2 Ma and 0.6 Ma. The horn ice disrupts the flow field at the leading edge, and distinct separation vortices appear behind the ice horn. As the angle of attack increases, the separation vortex region behind the lower ice horn gradually shrinks, while that behind the upper ice horn expands. At the stall angle, separation vortices cover the entire upper surface, causing complete flow detachment and a sharp drop in lift.

Figure 6: Velocity distributions and streamlines: (a) 0.2 Ma; (b) 0.6 Ma.

2.5 Comparison of Aerodynamic Loads

By extracting the aerodynamic forces along the chord and normal directions of the ice surface, it is possible to effectively compare the aerodynamic loads experienced at different angles of attack. Fig. 7 shows the computed aerodynamic loads per unit length on the ice surface. As the AOA increases, the normal load increases, while the chordwise load varies only slightly. Normal load peaks at stall angle with a maximum magnitude of 1.89 kN/m and declines sharply afterwards due to large-scale flow separation. In addition, an increase of Mach number leads to a noticeable rise in aerodynamic loads, imposing higher requirements on the structural strength of the ice shape. These loads serve as the input for the subsequent structural analysis of the simulated ice.

Figure 7: Aerodynamic loads per unit length on the ice surface.

3 Finite Element Simulation of the Simulated Ice Shape

Based on the above analysis, the computed profile was extruded by 250 mm to create a mechanical simulation model for subsequent load application and strength assessment.

3.1 Structural Design of Simulated Ice

The design of the simulated ice structure follows four principles: (a) adoption of commercially available conventional raw materials; (b) feasibility of high-precision mechanical machining; (c) moderate structural elasticity to accommodate minor elastic deformation of the host wing during flight; (d) simplified structural layout and manufacturing workflow for on-site engineering implementation.

Based on these principles, a composite sandwich structure consisting of a glass-fiber reinforced plastic (GFRP) skin with a UV-curable resin 3D-printed core was adopted, as shown in Fig. 8. To reduce weight, the 3D-printed core is hollow with a wall thickness of 1 mm.

Figure 8: Schematic diagram of the proposed design.

Key advantages of the proposed design are summarized below: (a) aluminum foil tape is inserted between the ice specimen and airfoil to avoid direct contact and facilitate convenient post-test disassembly and surface cleaning; (b) 3D-printing guarantees high forming precision, low weight and satisfactory surface hardness of the inner core; (c) GFRP wraps the core and extends 50 mm chordwise to improve interfacial bonding strength; (d) assembly can be completed under ambient dry conditions without rigorous environmental constraints.

To meet the requirements of processing accuracy and material applicability, UV-curable resin was selected for processing. The UV light curing processing method is technologically mature and has no special process requirements. The UV-curable resin system adopts a cationic photopolymerization epoxy system. Since it does not require a heating source, it can undergo polymerization through UV irradiation. Traditional thermosetting epoxy requires a temperature increase of over 60°C, which can cause slight melting of ice corners and loss of corresponding geometric features. By using segmented and thin-layer UV exposure, the heat release and shrinkage stress during polymerization can be reduced; glass fiber fillers further inhibit part warping and dimensional deviation, making it suitable for high-precision and complex curved surface ice model printing [29].

3.2 Finite Element Model

Simulations were performed using Abaqus/Standard. The airfoil was modeled as an analytical rigid body since the airfoil deformation during actual flight was small and negligible. C3D8R elements were used to mesh the composite material and the 3D-printed core. The finite element model is shown in the Fig. 9a, comprising 422,550 elements with a mesh quality greater than 0.9. In the model, three layers of zero thickness cohesive elements [30] were also inserted at the inner surface of the outer glass fiber layer (the blue outer adhesive layer in Fig. 9b), the inner surface of the core (the yellow middle adhesive layer in Fig. 9b), and the inner surface of the inner glass fiber layer (the red inner adhesive layer in Fig. 9b). By adding this adhesive layer, simulations can be conducted on the surface delamination between the composite material and the core material (the blue outer adhesive layer and the yellow middle adhesive layer), as well as between the structure and the airfoil (the red inner adhesive layer).

Figure 9: Finite element model: (a) total elements; (b) cohesive elements.

To evaluate the overall structural strength, damage criteria were defined. For interface failure, a traction-separation constitutive law was adopted for the cohesive elements. The elastic constitutive equation is expressed as: t=tntstt=Knn Kss Kttδnδsδt(1) where t and δ are traction stress and separation displacement vector consisting of the normal and two shear components. Under the mixed mode loading, damage initiation and evolution depend on damage mode I, II and III.

A quadratic stress criterion was used to determine damage initiation: tntn02+tsts02+tttt02=1(2) where t n 0 , t s 0 and t t 0 represent the peak values of the nominal stress when the deformation is either purely normal to the interface or purely in the first or the second shear direction, respectively. The McCauley operator is used here to show that the pure compressive stress has no contribution to the damage of cohesive elements. A linear degradation equation defines the damage evolution process: D=δmfδmmax−δm0δmmaxδmf−δm0(3) where δ m 0 , δ m f and δ m max are the initial failure effective displacement, ultimate failure equivalent displacement, and maximum effective displacement during the loading process. The effective displacement can be calculated by the following expression:

ABAQUS intrinsic variable SDEG quantifies cohesive-layer damage: damage initiates when SDEG > 0 and full debonding (complete load loss) occurs at SDEG = 1. The Tsai-Wu failure criterion [31] is used to predict intralaminar failure of GFRP and coded into a user subroutine. Composite remains intact when SDV1 < 1 and initiates damage when SDV1 = 1, where SDV1 denotes the failure index defined by FTasi−Wu=F1σ1+F2σ2+F11σ12+F22σ22+F66σ122+2F12σ1σ2(5) F1=1XT−XC,F11=1XTXC,F2=1YT−YC,F22=1YTYC,F66=1S122,F12=−F11F222(6) where XT, XC, YT, YC and S12 are tensile and compressive strength in principal direction 1, tensile and compressive strength in principal direction 2, and in-plane shear strength, respectively.

The UV-curable resin (9111A) has a density of 1100 kg/m3, an elastic modulus of 3.64 GPa, and a Poisson’s ratio of 0.33. The thermal deformation temperature of this material is 54°C, and the glass transition temperature is 80.8°C. Both of these are higher than the flight environment, ensuring that the material will not deform during the flight test. It undergoes curing under 355 nm ultraviolet light and can be processed on a 3D printer (ZRapid iSLA660), offering great processing convenience.

Table 3 lists key parameters of the epoxy-based GFRP to be used in finite element analysis, including elastic modulus (E1, E2, E3) in three principal directions, Poisson’s ratio (v12, v13, v23), shear modulus (G12, G13, G23), in-plane shear strength (S12), as well as directional tensile strengths (XT, YT) and compressive strengths (XC, YC).

The mechanical properties of the adhesive interlayer, including cohesive strength, fracture energy, and stiffness parameters, are detailed in Table 4. The equivalent thickness of the adhesive layer is 0.1 mm.

As demonstrated in Section 2, aerodynamic loads reached maximum magnitude at stall conditions, and normal load dominates structural failure risk during practical flight certification tests. The axial load, when viewed from a macroscopic perspective, exerts a pressure load on the simulated ice model. However, it does not cause the ice model to fall off, so it is not taken into account. Since the simulation model is a uniform cross-section model, the aerodynamic loads applied to each section are consistent. To facilitate comparison with the subsequent test conditions, all nodes on the upper surface of simulated ice are kinematically coupled to a single reference point, and displacement-controlled normal loading is applied incrementally until global structural collapse. The structural safety factor is defined as the ratio between ultimate simulated failure load and peak operational aerodynamic load.

3.3 Numerical Simulation Results

The load-displacement curve shown in Fig. 10 indicates that the ultimate failure load of the 250 mm-span specimen is 2.20 kN, corresponding to a unit-length failure load of 8.80 kN/m, far beyond the maximum applied aerodynamic normal load of 1.89 kN/m. The resultant structural safety factor is calculated as 4.65, which fully satisfies flight-test structural integrity requirements.

Figure 10: The load-displacement curve.

Fig. 11 shows the adhesive damage cloud at displacement loading of 0.028 mm. Initial interfacial damage emerges at the leading edge of inner and middle adhesive layers, accompanied by minor localized damage at outer adhesive bonding junctions.

Figure 11: Adhesive damage cloud at displacement loading of 0.028 mm: (a) middle adhesive; (b) outer adhesive; (c) inner adhesive.

Fig. 12 shows the adhesive damage cloud at the peak load. With load progressively increasing to the peak load, damage propagates along all adhesive edges, while chordwise-extended outer adhesive layer remains intact and sustains residual structural load.

Figure 12: Adhesive damage cloud at peak load: (a) middle adhesive; (b) outer adhesive; (c) inner adhesive.

Fig. 13 shows the adhesive damage cloud under further loading after the peak load. It is seen that damage in the middle and outer adhesive layers expands from the core junction outward, reducing the effective bonding area and causing a continuous decrease in load-carrying capacity, eventually leading to debonding.

Figure 13: Adhesive damage cloud under further loading after the peak load: (a) middle adhesive; (b) outer adhesive; (c) inner adhesive.

At the final loading stage, the outer and middle adhesive layers are completely damaged, whereas the inner adhesive layer is still not fully deboned; nevertheless, the structure loses its load-carrying capacity (Fig. 14). The analysis shows that the extended GFRP strips on both sides play a crucial role in the structural load bearing capacity, and their integrity directly affects the bonding performance.

Figure 14: Adhesive damage cloud at the final loading stage: (a) middle adhesive; (b) outer adhesive; (c) inner adhesive.

In addition to adhesive damage, the integrity of the GFRP skins and the 3D-printed core was also evaluated. Fig. 15 shows the distribution of Tsai-Wu damage factor of GFRP at the peak load. Fig. 16 illustrates the von Mises stress distribution of the 3D printed sandwich core at the peak load. The Tsai-Wu damage factor of the GFRP is only 4.75 × 10−4, indicating a sufficient safety margin; the GFRP main body does not fail. The maximum von Mises equivalent stress of hollow 3D-printed core is 14.07 MPa, below material ultimate strength without structural strength.

Figure 15: The distribution of Tsai-Wu damage factor of GFRP at the peak load.

Figure 16: The von Mises stress distribution of the 3D printed sandwich core at the peak load.

4 Experimental Analysis

To validate the finite element predictions, quasi-static mechanical tests were conducted on a fabricated specimen. An aluminum fixture with the corresponding cross-section was machined as a base to simulate the wing leading-edge geometry. A simulated ice core with a wall thickness of 1 mm and a length of 250 mm was fabricated using 3D-printing. The simulated ice core was installed on the aluminum fixture surface using a hand lay-up composite molding process. The test components are shown in Fig. 17.

Figure 17: Test components: (a) 3D-printed ice core; (b) aluminum fixture.

A canvas bag was adhered to the upper surface of the simulated ice for uniform load via hydraulic actuator loading. The failure mode of simulated ice when the structure lost its load-carrying capacity is shown in Fig. 18. The inner surface of the core separated from the inner GFRP layer, and the upper extended part of the outer GFRP layer deboned, causing the structure to lose load-carrying capacity. The test failure mode is consistent with the simulated one. The test failure load of the 250-mm test section is 2.42 kN, while the finite element simulated one is 2.20 kN. The error is 9.10%, indicating that the simulation can accurately predict the structural load bearing capacity.

Figure 18: The failure mode of simulated ice.

5 Conclusions

Leading-edge ice accretion severely deteriorates airfoil aerodynamic performance and threatens in-flight safety. Therefore, dry-air flight certification with substitute simulated ice is indispensable for aircraft airworthiness approval. Based on mature processing methods and materials, this paper proposes a polymer composite simulated ice structure that meets the requirements of airworthiness certification, and conducts corresponding CFD aerodynamic analysis, finite element simulation analysis, and experimental verification to prove the reliability of the designed structure. Through the above analyses, the following conclusions are mainly drawn:

• (a)Aerodynamic performance degradation: The glaze ice at the leading edge induces disordered flow around the ice shape, affecting the entire flow field and causing earlier boundary layer separation and massive vortex formation. This leads to premature stall, posing a serious threat to flight safety. At 0.2 Ma, the maximum lift coefficient drops by approximately 58.6% (from 1.57 to 0.65), and the stall angle reduces from 13° to 5.5°. At 0.6 Ma, the maximum lift coefficient drops by approximately 45.4% (from 1.08 to 0.59), and the stall angle reduces from 6.5° to 4.5°. The maximum normal aerodynamic unit-length load acting on the icing surface reaches 1.89 kN/m at stall. This value serves as the critical design load for simulated ice structural design.
• (b)Structural design and analysis: A composite sandwich structure was designed, with the 3D-printed core providing the aerodynamic shape and inner/outer GFRP layers extending 50 mm chordwise to enhance structural strength. The axial load exerts a pressure load on the simulated ice model which does not cause the structure to fall off, so it is not taken into account. Finite element analysis shows that the safety factor is 4.65 against the maximum aerodynamic normal unit-length load (1.89 kN/m), meeting the safety requirements for actual flight tests. The structural failure mode is debonding of the core adhesive layer and the outer adhesive layer, leading to loss of load-carrying capacity, while the GFRP main body and the core remain undamaged.
• (c)Experimental validation: An aluminum fixture of the wing leading edge and a 3D-printed core were manufactured, and the simulated ice was installed using a hand lay-up molding process. A canvas bag was adhered to the designated area for loading. The error between test failure load and the finite element simulated one is 9.10%, indicating good prediction of structural load bearing capacity. Debonding between the core and the outer GFRP adhesive layer was also observed in the experiment, validating the simulation results.

Overall speaking, the designed simulated ice meets the structural requirements for dry-air flight tests in certification campaigns. This paper analyzes the mechanical properties of the designed polymer composite sandwich structures based on the evolution process of interface detachment and the stress of each component. The results are verified through experiments. For actual flight test applications, further research should focus on three-dimensional wing aerodynamic interference effects and multi-condition durability assessment of simulated ice. Since this work only focuses on one typical ice shape, future work should expand the research scope to various icing configurations and complex loads such as vibration under diverse atmospheric icing environments to enrich the application boundary of the proposed simulated ice.