Mechanical response under low-velocity impact of composite/Aluminum interfaces in fibre metal laminates. Parameter identification via finite element simulations, data interpolation and minimization procedures

Fibre Metal Laminates (FMLs) are hybrid composites obtained by alternating thin metallic sheets, typically aluminium alloys, with fibre-reinforced polymer plies. Their adoption in the aeronautical field is motivated by the ability to combine high specific strength with improved damage tolerance. In service, however, these structures are frequently exposed to accidental impacts. While high-energy events may cause immediate perforation, low-energy events are often more insidious: they may not leave visible marks on the outer surface, but can initiate delamination at the composite/aluminium interfaces, reducing stiffness and compromising structural integrity over time. In order to reproduce such conditions in a controlled environment, low-velocity impact (LVI) tests have been standardised by international regulations. These tests, performed with instrumented drop-weight towers, are widely used to simulate realistic scenarios such as tool drops during maintenance, runway debris impacts, or hail strikes, and to assess the residual performance of composite and hybrid laminates after impact. This thesis addressed this problem by focusing on the identification of cohesive parameters that govern delamination in a particular type of laminates called GLARE (GLAss REinforced laminate) under LVI. The work originated from an experimental campaign carried out by Marco Scaglia at Politecnico di Milano. In that study, some GLARE-type specimens were produced at Politecnico using aluminium sheets and glass-fibre prepregs supplied by industrial partners, with the laminates cured in-house under controlled conditions. The subsequent impact tests were conducted at an external accredited facility (TEC Eurolab, Modena), using an instrumented drop-weight machine in accordance with the relevant standards. The campaign provided force–time and displacement–time histories of the impactor, together with qualitative information on delamination patterns. While those data were valuable, their interpretation in terms of interface properties remained unresolved. This thesis therefore continued that work with the objective of reconstructing the experiments numerically and extracting the underlying cohesive parameters through inverse analysis. The experimental setup was reconstructed in Abaqus/Explicit with a finite element (FE) model that included the entire system—laminate, impactor, and boundary conditions. Aluminium layers were described by the Johnson–Cook plasticity law; glass-fibre plies by a Hashin-type damage formulation; and composite/aluminium interfaces by a bilinear traction–separation cohesive law with Benzeggagh–Kenane mixed-mode fracture evolution. Given the complexity of the model and the number of elements, in order to reduce the computational burden of a single simulation, the entire FE model was reduced to a quarter by exploiting geometric and load symmetry: a choice that preserved the global response while decreasing the number of elements by a factor of four and lowering the wall-clock time of a single forward simulation from roughly four hours (full model) to about one hour (quarter model). Before any identification with inverse analysis, an extensive campaign of forward analyses was performed: in each forward run, the interface parameters of the bilinear laws in mixed mode behaviour—initial stiffnesses (Kn, Ks, Kt), nominal strengths (Tn, Ts, Tt), and fracture energies (Gn, Gs, Gt)—were systematically varied and the resulting force–time and displacement–time histories were examined in terms of peak force, maximum displacement, contact duration, and absorbed energy. This progressive exploration clarified parameter roles and produced numerical trends that qualitatively matched the experiments, thereby defining meaningful search intervals for the subsequent identification phase. Building on this basis, the calibration targeted both a single-parameter identification scenario (Gn, the Mode I fracture energy) and a two-parameter identification scenario (Gn, Tn, with Tn the nominal strength along Mode I direction): Gn was selected for its primary influence on energy absorption and displacement levels, while Tn was introduced because of its influence on the initial stiffness and peak response. The interface parameters were refined through inverse analyses using two complementary strategies. (i) A classical approach, where Abaqus and MATLAB were directly coupled through the Trust Region Reflective algorithm. At each iteration, the algorithm updated the unknown parameters, ran a full FE direct simulation, and evaluated the cost by least-squares comparison of force–time and displacement–time histories against the target data. The minimization proceeded by progressively adjusting the parameters to reduce the cost, until convergence was reached when no further reduction was possible. (ii) An interpolation-based approach, where a set of cohesive parameter combinations was selected a priori as nodes in the cost–parameter space. For each node, a forward simulation was performed and the corresponding cost was computed by least-squares comparison with the target curves (pseudo-experimental or experimental). These discrete values were then interpolated to reconstruct a continuous cost surface, on which MATLAB’s fminsearch function was applied to locate the minimum and the associated parameters. Validation of the methodology was first carried out with pseudo-experimental datasets, both without noise and with controlled noise levels, to ensure robustness of the procedure. Finally, the inverse analysis was applied to the real experimental data. In this case, displacement–time histories proved to be the most reliable source of information due to limitations in force measurement during the drop-weight tests. Across pseudo-experimental and experimental datasets, the results indicated that Gn exerted a markedly stronger influence on the overall response, whereas Tn was less influential; accordingly, its identification was consistently less accurate than that of Gn. In computational terms, the classical approach required between fifteen and more than thirty forward simulations per identification, often pushing the total runtime beyond one day even after model reduction; by contrast, the interpolation-based approach reduced the effort to a few hours—the time required to compute costs at the preselected interpolation nodes via forward simulations—at the price of an accuracy that depended on node density and could overlook secondary minima if the grid was too coarse. Possible extensions include the joint identification of three parameters—for instance adding the initial stiffness Kn in Mode I direction—and the use of denser or adaptively refined interpolation grids to enhance accuracy while retaining manageable computational costs.

I Fibre Metal Laminates (FML) sono materiali compositi ibridi ottenuti alternando sottili lamine metalliche, tipicamente leghe di alluminio, con strati in materiale fibrorinforzato. Il loro impiego in ambito aeronautico è motivato dalla capacità di combinare un’elevata resistenza specifica con una buona tolleranza al danno. In esercizio, tuttavia, queste strutture sono spesso esposte a impatti accidentali. Se gli eventi ad alta energia possono provocare una perforazione immediata, quelli a energia più bassa risultano più subdoli: possono non lasciare segni visibili in superficie, ma innescare una delaminazione all’interfaccia composito/alluminio, riducendo la rigidezza e compromettendo nel tempo l’integrità strutturale. Per riprodurre tali condizioni in modo controllato, le prove di impatto a bassa velocità (Low-Velocity Impact, LVI) sono state normate a livello internazionale. Queste prove, eseguite con torri di caduta a peso strumentate, sono largamente impiegate per simulare scenari realistici (caduta di attrezzi di manutenzione, detriti di pista, grandine) e per valutare le prestazioni residue di laminati compositi e ibridi dopo l’impatto. La presente tesi ha affrontato il tema dell’identificazione dei parametri coesivi che governano la delaminazione nei laminati GLARE (GLAss REinforced laminate) sottoposti a impatti a bassa velocità (LVI). Il lavoro si è inserito in continuità con una campagna sperimentale condotta da Marco Scaglia presso il Politecnico di Milano. In quello studio, alcuni provini GLARE sono stati realizzati al Politecnico utilizzando lamiere di alluminio e prepreg in fibra di vetro forniti da partner industriali, con ciclo di cura eseguito in condizioni controllate. Le prove di impatto sono state successivamente svolte presso un laboratorio esterno accreditato (TEC Eurolab, Modena), mediante una torre di caduta strumentata conforme alle normative di riferimento. La campagna sperimentale ha fornito dati quantitativi, sotto forma di curve forza–tempo e spostamento–tempo dell’impattatore, oltre a osservazioni qualitative sui meccanismi di delaminazione. Tuttavia, la traduzione diretta di tali risultati in termini di proprietà d’interfaccia era rimasta irrisolta. La tesi si è quindi posta l’obiettivo di proseguire quel lavoro, ricostruendo numericamente le prove sperimentali e identificando, tramite analisi inverse, i parametri coesivi sottostanti. L’assetto sperimentale è stato ricostruito in Abaqus/Explicit mediante un modello a elementi finiti (FE) che comprendeva l’intero sistema — laminato, impattatore e condizioni al contorno. Gli strati in alluminio sono stati descritti dalla legge plastica di Johnson–Cook; i compositi in fibra di vetro da una formulazione di danno di Hashin; le interfacce composito/alluminio da una legge coesiva bilineare trazione–separazione con evoluzione del danno in modo misto secondo Benzeggagh–Kenane. Vista la complessità del modello e il numero di elementi, per ridurre l’onere computazionale di una singola simulazione, l’intero modello FE è stato ridotto a un quarto sfruttando la simmetria geometrica e di carico: questa scelta ha preservato la risposta globale riducendo il numero di elementi di un fattore quattro e abbassando il tempo di una singola simulazione diretta da circa quattro ore (modello completo) a circa un’ora (modello a un quarto). Prima di qualsiasi identificazione tramite analisi inversa, è stata svolta un’ampia campagna di analisi dirette: in ciascun lancio, i parametri coesivi delle leggi bilineari in comportamento misto — rigidezze iniziali (Kn, Ks, Kt), resistenze nominali (Tn, Ts, Tt) ed energie di frattura (Gn, Gs, Gt) — sono stati variati sistematicamente, analizzando gli effetti su picco di forza, massimo spostamento, durata del contatto ed energia assorbita. Questo percorso progressivo ha chiarito il ruolo dei parametri e ha prodotto risultati numerici qualitativamente simili a quelli sperimentali, definendo intervalli di ricerca significativi per la successiva fase di identificazione. Su tali basi, la calibrazione si è concentrata su due scenari: identificazione di un singolo parametro (Gn, l’energia di frattura in Modo I) e identificazione di due parametri (Gn, Tn, con Tn la tensione massima lungo la direzione di Modo I). Il parametro Gn è stato scelto per il suo effetto principale sull’energia assorbita e sui livelli di spostamento, mentre Tn è stato introdotto per la sua influenza sulla rigidezza iniziale e sul valore massimo della risposta. L’identificazione dei parametri d’interfaccia è stata condotta tramite analisi inverse seguendo due strategie complementari. (i) Nell’approccio classico, Abaqus e MATLAB sono stati accoppiati direttamente tramite l’algoritmo Trust Region Reflective: a ogni iterazione i parametri incogniti venivano aggiornati, si avviava una simulazione diretta completa e il costo veniva valutato mediante un confronto ai minimi quadrati tra i dati forza–tempo e spostamento–tempo e i dati di riferimento. La minimizzazione procedeva per aggiustamenti successivi fino a quando ulteriori riduzioni del costo non risultavano più significative. (ii) Nell’approccio basato su interpolazione, è stato selezionato a priori un insieme di combinazioni dei parametri coesivi come nodi nello spazio costo–parametri. Per ciascun nodo è stata eseguita una simulazione diretta ed è stato calcolato il costo mediante un confronto ai minimi quadrati con le curve target (pseudo-sperimentali o sperimentali). I valori discreti così ottenuti sono stati poi interpolati per ricostruire una superficie continua del costo, sulla quale la funzione fminsearch di MATLAB ha individuato il minimo e i parametri corrispondenti. La metodologia è stata inizialmente validata con dataset pseudo-sperimentali, sia privi di rumore sia con rumore controllato, e quindi applicata ai dati reali; in quest’ultimo caso, i dati spostamento–tempo si sono rivelati più affidabili, alla luce delle note limitazioni nella misura della forza nelle prove a caduta. Nel complesso dei casi pseudo-sperimentali e sperimentali, i risultati hanno mostrato che Gn ha inciso in misura molto maggiore sulla risposta globale, mentre Tn si è rivelato meno influente; di conseguenza, la sua identificazione è risultata sistematicamente più incerta rispetto a quella di Gn. Dal punto di vista computazionale, l’approccio classico ha richiesto in media da quindici a oltre trenta simulazioni dirette per ogni identificazione, con tempi complessivi che spesso hanno superato le 24 ore anche dopo la riduzione del modello; al contrario, l’approccio con interpolazione ha ridotto l’impegno a poche ore — in pratica il tempo necessario a popolare la griglia di interpolazione valutando il costo nei nodi tramite simulazioni dirette — a fronte però di un’accuratezza dipendente dalla densità dei nodi e del rischio di non cogliere minimi secondari se la griglia era troppo diradata. I possibili sviluppi includono l’identificazione congiunta di tre parametri — ad esempio aggiungendo la rigidezza iniziale Kn in Modo I — e l’uso di griglie di interpolazione più dense o adattive, per incrementare la precisione mantenendo il costo computazionale entro limiti gestibili.