1. Introduction
Composite materials, such as carbon-epoxy laminates, are increasingly favored in structural applications due to their high strength-to-weight ratio and customizable properties. Designing with these materials requires not only defining the overall geometry but also determining the optimal laminate layup for its construction (Reference Voelkl and WartzackVoelkl & Wartzack, 2018). Laminated composites are widely used in aerospace, automotive, and civil engineering for components where in-plane stresses dominate. Thus, laminate layup design is typically guided by in-plane shear and normal stress levels calculated using Classical Laminate Theory (CLT), which assumes negligible or inconsequential out-of-plane stresses. While this approach ensures adequate strength to resist critical in-plane loads, it cannot account for weaknesses or failures in the through-thickness or out-of-plane directions (Dong et al., Reference Dong, Pister and Taylor1962; Reissner & Stavsky, Reference Reissner and Stavsky1961). The out-of-plane strength of composite laminates is generally at least an order of magnitude lower than their in-plane strength (Reference Wowk, Marsden and ThibaudeauWowk et al., 2020). In certain scenarios, such as at component edges (Reference Yen and HwuYen & Hwu, 1993), cut-outs (Reference Tian, Yang and WangTian et al., 2016), or thickness variations caused by ply drop-offs (Reference Kim, Rhee and ChoKim et al., 2008), out-of-plane strength becomes a critical design factor. In these situations, out-of-plane stresses are often highly localized, and concentrated near laminate edges at the interfaces between adjacent plies. CLT is incapable of predicting these out-of-plane edge stresses, forcing designers to rely on computationally expensive 3D finite element (FE) models for accurate analysis when such stresses are critical (Reference Mittelstedt, Becker, Kappel and KharghaniMittelstedt et al., 2022). Although research has demonstrated that stacking sequence significantly affects the generation of out-of-plane interlaminar stresses, the time and resources required to create and run detailed composite FE models often prevent their use during the laminate design process (Reference HerakovichHerakovich, 2012; Reference Wowk, Marsden and ThibaudeauWowk et al., 2020). Instead, edge effects due to interlaminar stresses are often left unaddressed early on and are only tackled through approximate methods or testing further down the design process (Reference Wowk, Marsden and ThibaudeauWowk et al., 2020). By this stage, laminate stacking sequences can no longer be altered, and less efficient solutions must be implemented to address unforeseen failure modes.
2. State of the art
2.1. Interlaminar stresses
The presence of free edges in laminated material systems introduces an additional level of complexity to the design process. The state of stress in the vicinity of free edges is three-dimensional, with nonzero through-thickness stresses. The through-thickness stresses include the interlaminar normal stress σ z and two interlaminar shear stresses τ xz and τ yz as depicted in Figure 1, which includes a qualitative stress profile of the interlaminar shear stress τ xz along the zy-plane. These through-thickness stresses are called interlaminar stresses.
Finite-width laminated coupon under axial load
Stress-free boundaries typically present the most severe interlaminar stresses. The interlaminar stresses in the boundary layer are critical in structural applications because they can lead to delamination-type failures at loads significantly below those causing in-plane failure (Reference HerakovichHerakovich, 1998). These stresses should be considered in the design of laminated structures with free edges, including those around holes and cutouts (Reference Hajikazemi and Van PaepegemHajikazemi & Van Paepegem, 2018).
Analyzing interlaminar stresses presents a significant challenge because it requires solving a complex elasticity boundary-value problem, for which no exact solution currently exists (Reference Mittelstedt, Becker, Kappel and KharghaniMittelstedt et al., 2022). Over the years, researchers have developed various approximate methods to address this issue, broadly categorized into analytical and numerical approaches. Comprehensive literature reviews by Kant and Swaminathan Reference Kant and Swaminathan(2000), Mittelstedt and Becker Reference Mittelstedt and Becker(2007) and Mittelstedt et al. Reference Mittelstedt, Becker, Kappel and Kharghani(2022) outline the evolution of these methods and their contributions to understanding interlaminar stresses and boundary-layer phenomena in laminated composites. Among the numerical approaches available, FE analyses using fully three-dimensional continuum elements – referred to here as solid FE – are generally the most accurate for capturing the complex stress fields that arise near free edges. However, solid FE simulations are computationally intensive and time-consuming, making them impractical for iterative design processes or for analyzing full-scale structures during early design stages (Reference Wowk, Marsden and ThibaudeauWowk et al., 2020).
2.2. Surrogate modeling and machine learning for laminated composites
Machine learning (ML) offers a compelling way to reduce or replace the high computational overhead of FE simulations for composite analysis. It is important to note that the following examples represent only a subset of the many innovative applications of ML in this field. For instance, Bessa et al. Reference Bessa, Bostanabad, Liu, Hu, Apley, Brinson, Chen and Liu(2017) introduced a data-driven framework that compresses the dimensionality of microstructures, enabling rapid prediction of stress-strain responses. Liang et al. Reference Liang, Liu, Martin and Sun(2018) used deep learning to replicate full-field stress predictions at near-FE fidelity yet dramatically lower computation times. Similarly, Yan et al. Reference Yan, Zou, Ilkhani and Jones(2020) developed a multiscale surrogate model to handle progressive damage analyses, while Ammasai Sengodan Reference Ammasai Sengodan(2021) employed deep neural networks to predict composite properties from reduced-dimensional microstructure descriptions. Recent work by Wang et al. Reference Wang, Zhang, Xuan, Fan, Fu, Xue and Yao(2024) demonstrated how such ML surrogates can also capture damage evolution in laminates. To the authors’ knowledge, no published methods exist that specifically address free-edge interlaminar stress via an ML surrogate model.
3. Research problem and research goal
Predicting interlaminar stresses is a critical challenge in the design and analysis of laminated structural components. These stresses – caused by mismatches in the elastic properties of adjacent layers – are highly concentrated near free edges and can initiate failure modes such as delamination and matrix cracking (Reference Asur Vijaya Kumar, Dean, Reinoso and PaggiAsur Vijaya Kumar et al., 2021). Shell elements, commonly used in structural simulations based on CLT, are computationally efficient and effective for capturing global structural behavior, such as in-plane stresses and deformations. However, they are inherently inadequate for resolving through-thickness stress distributions, including interlaminar stresses. To overcome this limitation, designers must use computationally expensive solid elements combined with very fine meshing to accurately calculate through-thickness stresses (Reference Wowk, Marsden and ThibaudeauWowk et al., 2020). This approach, while accurate, is resource-intensive and unsuitable for simulating full-scale structures or for iterative analysis during early design stages.
This leads to the central research question of this study: How can interlaminar stresses in composites be predicted more efficiently? Specifically, the study aims to develop a method that allows designers to use efficient shell elements while still accounting for interlaminar stresses. ML techniques, including Gaussian process regression (GPR) and artificial neural networks (ANNs), will be used to achieve this. The following sections address this research question in detail.
4. Methods and procedures
4.1. Overview
The basic procedure of the new method is presented in this chapter. It builds upon the approach by Marian et al. Reference Marian, Mursak, Bartz, Profito, Rosenkranz and Wartzack(2023) but is specifically modified for the usage with laminated structures. The method integrates FE simulations with ML techniques, leveraging high-fidelity FE outputs to train ML models. These models approximate the complex relationships governing interlaminar stresses, allowing designers to account for through-thickness stresses while using computationally efficient shell elements. This approach enables accurate predictions of interlaminar stresses even during early design stages, facilitating rapid iterations without relying on computationally expensive 3D solid element simulations. The overall workflow is depicted in Figure 2, which outlines the key steps: defining simulation parameters, performing FE simulations, and training and testing different machine learning models.
Flowchart of the ML method proposed in the present contribution (based on (Reference Marian, Mursak, Bartz, Profito, Rosenkranz and WartzackMarian et al., 2023))
4.2. Finite element modeling
FE simulations were performed to generate the training data required for the ML models. A rectangular, laminated, and symmetric plate was subjected to uniaxial tension. One end of the plate was fixed in all directions, while a 0.16 mm longitudinal displacement was applied at the opposite end, resulting in an axial strain of 0.01. The simulation setup is illustrated in Figure 3a. The model dimensions were 16 × 8 mm, with a total thickness of 1.2 mm and a 4-ply stacking sequence of [θ 1/θ 2]s. The two ply angles, θ 1 and θ 2, were varied between -90° and 90° to evaluate interlaminar stresses for different configurations. Only discrete values in 0.1° increments were allowed. Latin hypercube sampling was employed to efficiently sample the design space, yielding 8100 design points – a number equivalent to that obtained from a uniform grid sampling with 2° increments.
FE modeling details. (a) Model geometry, boundary, and loading conditions (b) Mesh configuration with 10 elements through the thickness of each ply and 16 elements in the edge region
To achieve highly accurate predictions of interlaminar stresses, the plate was discretized using fully integrated solid elements with quadratic shape functions, resulting in 33,600 elements. Since interlaminar stresses are most pronounced near free edges, mesh refinement was applied specifically in these regions. This approach was guided by the observations of Pipes and Pagano Reference Pipes and Pagano(1970), who demonstrated that edge effects are typically confined to a region within half the ply thickness (t/2 = 0.15 mm). Based on this finding, the coupon geometry was subdivided into a central region and two edge regions, as shown in Figure 3b. The mesh was designed to include ten elements through the thickness of each ply, ensuring sufficient resolution to accurately capture the variations in interlaminar stresses near the free edges.
Simulations were performed using an implicit solution algorithm with automatic switching to explicit solving and adaptive time stepping, utilizing a preconditioned conjugate gradient solver for efficient numerical convergence.
4.3. Material properties
The material properties used for the simulations are provided in Table 1. In FE representations of orthotropic laminates, Poisson’s ratio v and shear modulus G can be defined in two common ways: The first approach assumes identical values for v and G in all three material directions (Pipes & Pagano, 1970). The second approach assigns distinct values for the y-z direction (Reference Narendra and LagaceNarendra & Lagace, 1994), reflecting its typically higher Poisson’s ratio and lower shear modulus compared to the x-y and x-z directions.
Material properties of the laminate used in this study, as reported by Narendra & Lagace Reference Narendra and Lagace(1994)
Using identical values for v and G across all directions simplifies the model but leads to significant inaccuracies in predicting the magnitude of interlaminar stresses. Specifically, the maximum normal stress σz is underestimated by a factor of three, while τxz predictions remain unaffected (Reference Wowk, Marsden and ThibaudeauWowk et al., 2020). To ensure accurate predictions of interlaminar stresses, this study applied distinct values for v and G in the y-z direction. Although the half-ply-thickness guideline by Pipes and Pagano Reference Pipes and Pagano(1970) was based on simplified material parameters, it is well established that the spatial extent of the boundary-layer region near free edges remains roughly half a ply thickness (Reference Mittelstedt, Becker, Kappel and KharghaniMittelstedt et al., 2022). Hence, the approximate half-ply-thickness zone near free edges remains a robust guideline even under more refined material definitions.
4.4. Machine learning models
Machine learning methods – including Gaussian process regression and artificial neural networks – were employed to predict interlaminar stresses using MATLAB´s Machine Learning and Deep Learning Toolbox. Following the general workflow described by Marian et al. Reference Marian, Mursak, Bartz, Profito, Rosenkranz and Wartzack(2023), all input and output data was normalized to values between -1 and 1 using the MATLAB built-in function mapminmax. The ML models were then developed using training data to provide accurate predictions for unseen test data. To achieve this, the dataset was divided into 85% for training and 15% for testing. The coefficient of determination R 2 was used to evaluate the prediction quality due to its effectiveness in assessing model performance (Reference Chicco, Warrens and JurmanChicco et al., 2021).
The prediction quality of the ML models was evaluated after training with the original input data from the FE simulations. The goal of the study is to enable accurate predictions of interlaminar stresses at the free edges of laminated composites. The inputs to the ML models are the in-plane stresses σx and σy, taken from the center of the coupon. These values are easily obtainable from standard shell FE simulations, making the approach practical and efficient. The ML models then predict the interlaminar stress components ( σz, τxz and τyz) at three critical points near the free edges, as defined in Figure 4. The critical points were selected at the interfaces between adjacent plies, where delamination is most likely to initiate. This placement ensures the predictions focus on the most structurally significant regions of the laminate. The method relies on only a small number of input values ( σx and σy at these points), which corresponds to the stress data typically available from conventional shell FE simulations. This minimal input requirement makes the approach highly adaptable and computationally efficient.
Inputs and outputs of the ML models from the FE simulations
4.4.1. Gaussian process regression
Gaussian process regression is a ML method used to model complex relationships in data by treating functions as probabilistic distributions. In simple terms, GPR assumes that any finite set of function values (e.g., as f(x 1), f(x 2),…,f(x n)) follows a Gaussian distribution, making it well-suited for capturing patterns in data (Reference Rasmussen and WilliamsRasmussen & Williams, 2008).
A GPR model is defined by two key components: the mean function and the covariance function (or kernel). The mean function represents the average behavior of the system, and it is often set to zero for simplicity. The kernel, on the other hand, defines how the model generalizes from the training data to make predictions for new, unseen inputs (Reference DuvenaudDuvenaud, 2014). Mathematically, this is expressed as in Equation (1):
$$f(x)~GP[m(x),k(x,x^\prime ,\theta _k )]$$
Here, θ k represents the hyperparameters of the kernel, which control the model´s flexibility and smoothness.
The present study builds on the GPR implementation strategy described in Marian et al. Reference Marian, Mursak, Bartz, Profito, Rosenkranz and Wartzack(2023). MATLAB´s Machine Learning and Deep Learning Toolbox was employed to develop the GPR model, combining the Gaussian process with explicit base functions to map the inputs into a higher-dimensional space, allowing for more accurate predictions. MATLAB´s fitrgp function was used to train the model by estimating key parameters, including the kernel hyperparameters θ k, the noise variance σ 2, and the coefficients of the base functions β b (The MathWorks, Inc., 2024).
4.4.2. Artificial neural networks
Artificial neural networks are powerful tools for modeling complex, nonlinear systems, gaining widespread attention due to advancements in computational power. They have been successfully applied across diverse fields, including engineering, to solve challenges such as estimating manufacturing costs (Reference MumaliMumali, 2022), evaluating energy system efficiencies (Reference Ding, Wang and BiDing et al., 2011), forecasting environmental phenomena (Reference Dabbakuti and GDabbakuti & G, 2019), and simulating physical processes (Reference Noori and KalinNoori & Kalin, 2016). At their core, ANNs serve as nonlinear approximators capable of modeling intricate relationships within data (Reference Ferrari and StengelFerrari & Stengel, 2005).
Adapting the methodology described by Marian et al. Reference Marian, Mursak, Bartz, Profito, Rosenkranz and Wartzack(2023) in this study, ANNs were used to predict interlaminar stresses at the free edges of laminated composites. Multi-hidden-layer feedforward backpropagation networks were trained using MATLAB’s nnstart app. The training process involved iteratively adjusting neuron weights, starting from the output layer and propagating backwards through the network (Reference Lämmel and CleveLämmel & Cleve, 2020). To reduce the risk of overfitting, 15% of the training data was allocated for validation. The optimal ANN architecture was determined by testing configurations with zero to three hidden layers and varying the number of neurons per layer between 10 and 15.
5. Results and discussion
This chapter presents the results of predicting interlaminar stresses at free edges in laminated composites using ML techniques. An ANN and GPR model were employed to approximate the stress components σ z, τ xz, and τ yz at critical points along the laminate edge. The performance of these ML models is evaluated using coefficients of determination and scatter plots comparing predicted and calculated values for unknown test data. Table 2 summarizes the coefficients of determination for the predicted interlaminar stresses at points along the interface between layers, where delamination is most likely to initiate (see Section 4.4). Across all stress components, both models achieved R 2 values close to 1, demonstrating high predictive accuracy. However, the GPR model outperformed the ANN in terms of R 2 for most cases, particularly τ yz.
Coefficients of determination of GPR and ANN predictions for interlaminar stresses against testing data after training
The scatter plots in Figures 5, 6, and 7 provide a detailed comparison of predicted versus calculated values for σ z, τ xz, and τ yz, respectively. Each figure contains six subplots: three for the ANN and three for the GPR model, corresponding to the three critical points along the laminate edge. The scatter plots reveal key differences in model performance:
Predicted versus calculated values (testing data) of the interlaminar normal stress σ z at the points 1, 2, and 3 for the GPR and ANN
Predicted versus calculated values (testing data) of the interlaminar shear stress τ xz at the points 1,2, and 3 for the GPR and ANN at Point 2
Predicted versus calculated values (testing data) of the interlaminar shear stress τ yz at the points 1, 2, and 3 for the GPR and ANN
-
Interlaminar normal stress σ z – Figure 5: Predictions from the GPR model align more closely with the diagonal line, indicating superior predictive accuracy. The ANN results show slightly more scatter, particularly in regions of high stress.
-
Interlaminar shear stress τ xz – Figure 6: The GPR model again demonstrates better alignment with the diagonal line. Predictions from the ANN exhibit noticeable scatter across the range, especially at the mid-range and extreme stress values.
-
Interlaminar shear stress τ yz – Figure 7: The ANN performs noticeably worse than the GPR model for τ yz, exhibiting greater scatter in its predictions. However, since the magnitude of τ yz is significantly lower than σ z and τ xz, this discrepancy is less critical.
The visual analysis provided by the scatter plots complements the R 2 metrics and highlights GPR model´s consistent edge in predictive accuracy. The ANN results, while slightly less accurate, remain robust, particularly given their faster training times and suitability for handling large datasets. In summary, GPR performed better overall, with higher R 2 values and superior performance in scatter plots, making it the preferred choice when predictive accuracy is paramount. Notably, GPR also offers uncertainty quantification, providing additional value in scenarios where confidence in predictions is critical.
A detailed analysis of the training dataset revealed that the interlaminar stress values vary smoothly with changes in ply angles. In other words, minor variations in ply orientation result in gradual, continuous changes in the target stress outputs rather than abrupt shifts. This inherent smoothness makes the data particularly amenable to GPR, which excels in interpolating smoothly across homogeneous datasets. By contrast, although ANNs are capable of capturing complex nonlinearities, the relatively smooth behavior of the interlaminar stresses likely favored the performance of GPR in this application.
In addition to accuracy, computational efficiency was evaluated. Each of the 8100 solid FE simulations, performed on a standard workstation (13th Gen Intel® Core™ i5 13600), took about three minutes of CPU time, totaling roughly 405 hours. However, parallelization on four cores reduced the wall-clock time to approximately 100 hours. Training the GPR and ANN models then required only about one minute per model, and subsequent inference – predicting interlaminar stresses for a new configuration – was effectively instantaneous. As a result, the initial cost of generating the simulation data is offset when multiple design iterations or large-scale optimization loops are anticipated. For instance, if a designer needs to explore many stacking sequences or run repeated sensitivity studies, the surrogate provides near-immediate stress predictions, which is especially advantageous in early-to-mid design phases where laminate architectures may still change frequently.
Beyond iterative design scenarios, this approach also benefits the simulation of large or complex composite parts. In such cases, efficient shell elements can be used to capture in-plane behavior, and the pretrained ML model can then calculate localized interlaminar stresses near free edges by drawing on normal stresses extracted from the shell solution. In this way, the method offers both (1) high-speed predictions for multiple configurations and (2) a practical means to handle big or intricate structures where 3D solid FE analysis might otherwise be prohibitively expensive.
6. Summary
In this contribution, an approach for predicting interlaminar stresses in laminated structures using ML techniques was presented. Currently, accounting for these stresses requires the use of computationally expensive FE models, which limits their application, especially during early design stages. To address this challenge, FE simulations were used to generate training data for ML models, specifically GPR and ANNs. These models were trained to predict interlaminar stresses at free edges based on readily available in-plane stress data from shell FE simulation.
Results demonstrate that both ML models achieve high predictive accuracy while significantly reducing computational effort. This method enables efficient early-stage design optimization for laminated composites, addressing critical edge effects that would otherwise require costly 3D FE analysis. The study’s findings pave the way for integrating ML techniques into composite design workflows, offering potential for improved structural efficiency and shorter design cycles.
While the study focuses on uniaxial tension, the versatility of ML models makes them applicable to more complex loading conditions, such as bending, torsion, and combined thermal-mechanical effects. Future work could refine stress predictions by incorporating these additional loading conditions and further simplifying the input requirements.
Acknowledgement
The presented research work is funded by the Deutsche Forschungsgemeinschaft (DFG, German Research Foundation) – 508331526. The authors thank the German Research Foundation for the financial support.
Author contributions
Conceptualization, M.G.; methodology, M.G.; software, M.G.; validation, M.G.; formal analysis, M.G.; investigation, M.G.; resources, M.G.; data curation, M.G.; writing—original draft preparation, M.G.; writing—review and editing, T.D., S.W. and D.K.; visualization, T.D. and M.G.; supervision, S.W. and D.K.; project administration, S.W. and D.K.; funding acquisition, S.W. and D.K.. All authors have read and agreed to the published version of the manuscript.
Open access
