On the machine learning approach for computational homogenization of short fiber-reinforced composites with elastoplastic behavior

Abstract

This study aims to develop machine learning (ML)-based homogenization models to efficiently predict the effective elastoplastic properties of short fiber-reinforced composites (SFRCs), reducing the reliance on computationally expensive micromechanical simulations. Sobol sampling is employed to generate training data by varying constitutive and microstructural parameters of representative volume element. Several ML models including artificial neural networks (ANN), support vector regression (SVR), random forest (RF), and extreme gradient boosting (XGB) are trained to predict homogenized stress-strain responses. A new approach is introduced that decomposes the stress-strain response into elastic and plastic components, allowing the ML models to learn these components separately and effectively. Additionally, the Taguchi method (L27 orthogonal array) is employed to minimize simulation runs and evaluate parameter sensitivity through ANOVA. The best-performing ML model is implemented in a finite element analysis (FEA) of an automotive component. Among all models, ANN demonstrated the highest accuracy in predicting the macroscopic elastoplastic response across a wide range of input parameters. Finally, the ANN-based elastoplastic material model is validated on a real-world macroscopic structure by implementing it into the finite element analysis of an automotive component. The results demonstrate that ML-based homogenization closely matches the traditional homogenization methods and also highlights its ability to effectively capture nonlinear material behavior.

Keywords

Short fiber reinforced composite machine learning representative volume element homogenization elastoplastic material

Introduction

In recent years, research on the prediction and modeling of the mechanical properties of composite materials has significantly increased. Due to their lightweight structure and high strength, composite materials are widely used in the automotive, aerospace, and industrial manufacturing sectors. However, the inherent anisotropy and nonlinear mechanical behavior arising from the interaction between fibers and the matrix pose significant challenges for accurately predicting their response under various loading conditions.^1–3 Therefore, computational homogenization techniques have been developed to estimate the homogenized response of composites. Homogenized response refers to the effective macroscopic behavior obtained by averaging the heterogeneous microstructural properties of the constituents (e.g., matrix and fibers). This response captures how the composite behaves as a uniform material under external loading, based on detailed micro-level interactions. Voyiadjis and Deliktas⁴ introduced a multilevel homogenization technique to obtain the effective properties of the heterogeneous crystalline structures. Representative Volume Element (RVE) analysis plays a crucial role in modeling the microstructural properties and their effects on the overall elastoplastic behavior of short fiber-reinforced composites (SFRC).^5,6 The RVE is defined as being large enough to statistically represent microstructural heterogeneity while remaining significantly smaller than the overall structure.⁷ These computational homogenization techniques facilitate the prediction of material performance without the need for costly and time-consuming experimental procedures.

The Mori–Tanaka method is a widely used homogenization technique for modeling the elastic and plastic properties of SFRCs.^8,9 This method accounts for fiber-matrix interactions, allowing for the prediction of the effective elastic and plastic properties of the material. Furthermore, the J2 plasticity model is commonly employed to capture the nonlinear behavior of elastoplastic short fiber-reinforced composites.^10–13 Additionally, various direct numerical simulation methods have been developed for computational homogenization. These include the full-field method, the mean-field method, and the Fast Fourier Transform (FFT) method. While finite element (FE)-based RVE analysis ensures high accuracy, it has disadvantages such as complex meshing requirements and high computational costs for intricate microstructures.^14,15 Alternatively, mean-field models enable faster computations by considering the average strain and stress in each microstructural phase.^16,17 On the other hand, the FFT method enhances computational efficiency by performing RVE analysis without requiring a stiffness matrix.^18,19 Digimat MF analysis is one of the most popular and effective tools for predicting material behavior by simulating fiber orientation, distribution, and fiber-matrix interactions. The flexible framework of Digimat MF provides a versatile approach to analyzing the impact of microstructural properties on macroscopic material responses.²⁰ This software plays a key role in computational homogenization, particularly for predicting the mechanical behavior of SFRC. However, the high computational cost of RVE-based simulations presents a significant challenge, especially when dealing with complex material systems and large datasets. To address this issue, machine learning (ML) has emerged in recent years as a complementary approach. ML algorithms, trained on data obtained from detailed RVE analyses, have been proposed to enable rapid predictions of effective material properties and mechanical responses. These approaches are often highlighted for their ability to reduce computational time.²¹

The emergence of data-driven approaches such as ML has opened new opportunities for developing reliable predictive material models. Unlike traditional methods that rely on empirical material modeling, this approach derives computational models directly from experimental data, bypassing the need for intermediate assumptions. With the growing availability of extensive datasets that characterize material behavior, there is a rising trend in leveraging data from analytical models, numerical simulations experiments, or hybrid sources. The capability to extract meaningful insights directly from data has attracted significant interest among engineers and researchers across diverse disciplines such as geotechnical engineering,²² micro/nanostructual systems,²³ structural engineering,²⁴ and constitutive modeling.²⁵ This has particularly become possible due to advancements in computational hardware technologies. Recently, Mentges et al.²⁶ combined finite element analysis and orientation averaging to generate the training data. This data was used to train the micromechanics-based ANN for estimating the elastic properties of short fiber composites. The results demonstrated that a well-trained ML model can rapidly predict the elastic properties of the material with high accuracy. Tariq et al.²⁷ implemented various boosting ML algorithms to investigate the linear elastic properties of continuous fiber reinforced composite. They determined that XGB can effectively and accurately predict the homogenized behavior of composites with coefficient of determination (R²) values of up to 0.999. Recently, Chen et al.²⁸ trained an ANN model to predict the hyperelastic behavior of short fiber-reinforced rubber composites. Their results showed that the model’s predictions closely matched those of FE simulations while drastically reducing computational expenses. Xu et al.²⁹ combined backpropagation neural networks with micromechanical models to investigate the elastic properties of multilayer braided composites. By training the neural network on data from microscale micromechanical simulations, they achieved highly accurate predictions of the composite properties. Similarly, Ghane et al.³⁰ developed a multiscale ML framework to estimate the elastic behavior of woven composites. They utilized high-fidelity FE simulation data from a meso-scale model of a plain-weave fabric composite to train an ANN model. The ANN model provided a good balance between accuracy and efficiency in predicting microscale and mesoscopic parameters. El Said³¹ used convolutional layers with Long Short-Term Memory (LSTM) layers to predict the stress response of laminated composites. By training the model on data from the RVE models database, he achieved accurate predictions of the nonlinear behavior of laminated composites. Le et al.³² developed a computational homogenization approach that employs neural networks to model nonlinear elastic heterogeneous materials. Their technique demonstrated greater efficiency compared to the traditional nested multilevel FE² method. Xue et al.³³ investigated the macroscopic elastic properties of composites materials with an ANN model trained on numerical results from fast Fourier transform method. They determined that ANN model can predict the elastic properties with high accuracy and improve the performance of existing homogenization models. Lourenço et al.³⁴ applied a range of physics-based constraints to a recurrent neural network model for learning elastoplastic behavior. Their study demonstrated that incorporating such constraints significantly enhances training efficiency, extrapolation capability, and predictive accuracy. Similar data-driven methodologies employing ML techniques to explore material behavior have been extensively studied, as evidenced by the research of Maharana et al.,³⁵ Zhao et al.,³⁶ Tariq et al.,³⁷ Santos et al.,³⁸ Mouloodi et al.,³⁹ Tariq et al.,⁴⁰ and Mahakur et al.⁴¹

This comprehensive literature review highlights extensive research on utilizing ML techniques to predict material behavior. The researchers agree that ML models offer significant computational advantages over traditional numerical methods like finite element analysis or homogenization techniques. While ANNs have been widely employed to predict the homogenized elastic properties of composite microstructures, their application to capturing the nonlinear elastoplastic behavior of composites, particularly SFRC, remains limited. Furthermore, the adoption of advanced ML techniques such as ensemble learning methods, has been relatively limited in this domain, with most studies focusing on elastic properties rather than the more complex nonlinear plastic responses. This research identifies a critical gap in the systematic evaluation and comparison of advanced ML algorithms for predicting the homogenized nonlinear elastoplastic behavior of SFRC, which is essential for understanding the mechanical performance in composite materials under realistic loading conditions.

Therefore, the present study aims to evaluate the effectiveness of various advanced ML techniques in predicting the elastoplastic behavior of SFRC materials through homogenization, while significantly reducing computational complexity compared to conventional methods. The predictive accuracy of these ML models is validated by performing finite element analysis on a real-world automotive component, in which the predicted material model is compared against conventional numerical methods. The ML approaches considered in this work include ensemble learning-based models such as Random Forest (RF) and Extreme Gradient Boosting (XGBoost), alongside Artificial Neural Networks (ANNs) and Support Vector Regression (SVR). These models are trained on a dataset generated using Sobol sampling methods that systematically vary nine key micromechanical parameters for capturing the homogenized material response of SFRCs. The predictive performance of each trained ML model is assessed using statistical metrics such as R², Mean Absolute Error (MAE), Mean Absolute Percentage Error (MAPE), and Root Mean Square Error (RMSE), along with visualization tools like scatter plots and bar charts. Furthermore, a combination of Taguchi analysis⁴² and ANOVA⁴³ is applied to quantify the influence of individual microstructural parameters on the homogenized material response. The predicted data from the most accurate ML model is then incorporated into finite element anlysis to simulate the behavior of a macroscopic automotive component, demonstrating its practical applicability and reliability in real-world engineering scenarios.

Mori–Tanaka method for short fiber-reinforced composites with an elastoplastic matrix

Modeling and simulating the behavior of short fiber-reinforced composites (SFRCs) with elastoplastic matrices requires the considering of both the matrix plasticity and fiber orientation. Homogenization methods are employed to predict the macroscopic mechanical response of composite materials. Homogenization involves analyzing the microstructural properties of a heterogeneous material to determine its homogenized macroscopic behavior.^44,45

The Mori–Tanaka method is a well-established and extensively validated in the literature and has long been recognized as one of the most reliable micromechanical models for predicting composite properties. Many studies reports that Mori–Tanaka predictions closely match experimental measurements for a variety of composites. For example, Katouzian and Vlase⁴⁶ applied the Mori–Tanaka model to a carbon-fiber/epoxy laminate and noted that the “results obtained experimentally and their comparison with the theoretically obtained values show a good agreement.” Similarly, Stránský et al.⁴⁷ used Mori–Tanaka to predict the thermal conductivity of composite and stated that the model’s results were “corroborated by available experimental data.” More recently, Yan et al.⁴⁸ derived elastic-modulus predictions for concrete with recycled aggregate. Their Mori–Tanaka-based model produced values “in good agreement with the experimental results” and was notably closer to measurements. These examples, among many others, demonstrate that Mori–Tanaka homogenization reliably captures the true behavior of composites and align very well with laboratory measurements.

In this context, the Digimat MF software provides a powerful tool for homogenizing short fiber reinforced elastoplastic composites.⁴⁹ By considering microstructural parameters such as fiber orientation, distribution, and mechanical properties, Digimat MF can simulate the macroscopic elastoplastic behavior of the material. The software utilizes advanced homogenization methods, such as Mori–Tanaka and Eshelby, to accurately model the effects of fiber-matrix interactions and plastic deformation processes.^50,51 Additionally, the models used to characterize elastoplastic behavior contribute to predicting the material’s response under various loading conditions.⁵² This study focuses on the homogenization processes conducted using Digimat MF to ensure an accurate simulation of the mechanical properties of short fiber reinforced elastoplastic matrix composites.

As shown in Figure 1, the aim of homogenization is to calculate the instantaneous effective stiffness tensor (C_eff) at the macro level based on the mechanical properties of the microstructural components: It can be mathematically represented as

\dot{σ} = C_{e f f} : \dot{ε}

(1)

Figure 1.

Homogenization process of short fiber-reinforced composites.

In this case, C_eff stands for the effective elastic stiffness tensor that needs to be ascertained, $\dot{σ}$ represents the macroscopic stress rate tensor, and $\dot{ε}$ for the macroscopic strain rate tensor. By embedding the inhomogeneities in a medium with effective thermomechanical properties, the Mori–Tanaka (MT) technique makes the assumption that the interactions between different types of inhomogeneities with finite volume concentrations are roughly taken into account. Therefore, the MT method uses equation (2) to determine the effective characteristics of composites⁵³:

C_{e f f} = C_{m} + v_{f} (C_{f} - C_{m}) A_{f}

(2)

Here, v_f is the fiber volume fraction, C_m is the matrix’s stiffness tensor, and C_f is the fiber’s stiffness tensor. The instantaneous stress concentration factor is defined as

A_{f} = {[I + S {(C_{m})}^{- 1} (C_{f} - C_{m})]}^{- 1}

(3)

In this case, S stands for the Eshelby tensor and I for the identity tensor. The shape, aspect ratios, and orientation of the fibers employed in composite materials all affect the Eshelby tensor:

S_{i j k l} = f (v_{m}, a s p e c t r a t i o s, o r i e n t a t i o n o f f i b e r s)

(4)

The Eshelby tensor for a general anisotropic material is calculated by the surface integral as shown in equation (5):

S_{i j k l} = \frac{1}{8 π} C_{m_{m n k l}} \int_{- 1}^{+ 1} d ζ_{3} \int_{0}^{2 π} {G_{i m j n} (\bar{ξ}) + G_{j m i n} (\bar{ξ})} d ω

(5)

where

G_{i j k l} (\bar{ξ}) = {\bar{ξ}}_{k} {\bar{ξ}}_{l} N_{i j} (\bar{ξ}) / D (\bar{ξ})

{\bar{ξ}}_{i} = \frac{ζ_{i}}{a_{i}}

ζ_{1} = {(1 - ζ_{3}^{2})}^{1 / 2} \cos w, ζ_{2} = {(1 - ζ_{3}^{2})}^{1 / 2} \sin w

ζ_{3 =} ζ_{3}

D (\bar{ξ}) = ϵ_{m n l} K_{m l} K_{n 2} K_{l 3}

N_{i j} (\bar{ξ}) = \frac{1}{2} ϵ_{i k l} ϵ_{j m n} K_{k m} K_{\ln}

K_{i k} = C_{m_{i j k l}} {\bar{ξ}}_{j} {\bar{ξ}}_{l}

(6)

The permutation tensor is represented by $ϵ_{i k l}$ and for a spherical inclusion the aspect ratios are equal (a₁ = a₂ = a₃). To account for fiber orientation effects on the mechanical behavior, the orientation distribution tensor Θ is introduced. By combining Θ with the Eshelby tensor, the orientation-averaged Eshelby tensor is computed as

{\bar{S}}_{i j k l} = \int S_{i j k l} (n) \cdot f (n) d Ω

(7)

The Eshelby tensor for a certain fiber orientation is denoted by S_ijkl (n). The probability density of fibers being orientated in direction n is described by the orientation distribution function (ODF), f(n). The differential solid angle across the unit sphere is denoted by dΩ. Another way to express the orientation-averaged Eshelby tensor is with the orientation distribution tensor:

{\bar{S}}_{i j k l} = 〈 S_{i j k l} (n) 〉 = \int S_{i j k l} n_{i} n_{j} n_{k} n_{l} f (n) d Ω

(8)

The tensor of orientation distribution, Θ, can be expressed as

Θ_{i j k l} = 〈 n_{i} n_{j} n_{k} n_{l} 〉 = \int n_{i} n_{j} n_{k} n_{l} f (n) d Ω

(9)

The unit vector n in this case indicates the fiber’s tiny direction. Software tools like Digimat evaluate this integration effectively using numerical methods and symmetry conditions. Advani and Tucker proposed an approach that uses the orientation distribution tensor O (or, more generically, a_ij and a_ijkl) to model the fiber orientation. The average of the fiber orientation for randomly oriented fibers is represented by the second-order orientation tensor, a_ij, in this instance, and a_ij = O. According to the Mori–Tanaka technique, the distribution of fiber orientation has a direct impact on the stiffness matrix’s constituent parts, which in turn effects the B coefficients that are employed in the Advani–Tucker formulation. Therefore, O is included in the stiffness tensor (C_ijkl) calculation for randomly oriented fibers.⁵⁴

C_{i j k l} = B_{1} a_{i j k l} + B_{2} (a_{i j} δ_{k l} + a_{k l} δ_{i j}) + B_{3} (a_{i k} δ_{j l} + a_{i l} δ_{j k} + a_{j l} δ_{i k} + a_{j k} δ_{i l}) + B_{4} (δ_{i j} δ_{k l}) + B_{5} (δ_{i k} δ_{j l} + δ_{i l} δ_{j k})

(10)

The fourth-order and second-order fiber orientation tensor components are denoted by the letters $a_{i j k l}$ and $a_{i j}$ , respectively. The components of the second-order unit tensor are denoted by $δ_{i j}$ . The components of a transversely isotropic unidirectional composite’s stiffness matrix are associated with the coefficients B.⁵⁵

In order to simulate short fiber-reinforced composites in a straightforward but accurate way, it is assumed that the matrix phase behaves elastoplastically and the fiber phase stays elastic throughout. According to this assumption, the stress state and internal state variables (such as damage, kinematic hardening, or isotropic hardening) determine the averaged concentration factors of the matrix phase. Different plasticity models have been proposed in the literature to describe plastic deformation in composite materials. Among these, the von Mises J2 plasticity model, which is based on the second invariant of the deviatoric stress tensor and assumes yield occurs when the distortion energy reaches a critical value, has been widely used to model the behavior of short fiber-reinforced polymers.^56–58 This model can reasonably approximate the macroscopic yield behavior and strain hardening response of polymers using isotropic or kinematic hardening rules. While alternative constitutive models, such as those based on shear transformation zones for post-yield behavior in polymers^59,60 or viscoplastic models for capturing brittle-to-ductile transitions in polycrystalline materials,⁶¹ could be considered in other contexts, they primarily account for large deformations and time-dependent effects. Since this study does not involve such effects, the J2 plasticity model with nonlinear isotropic hardening was chosen to characterize the matrix behavior, ensuring a balance between accuracy and computational efficiency. Therefore, the J2 plasticity model with nonlinear isotropic hardening is adopted, and the incremental elastoplastic behavior of the matrix phase is defined by the following equation (11):

\dot{σ} = C : (\dot{ε} - {\dot{ε}}^{p})

{\dot{ε}}^{p} = \dot{λ} \frac{\partial f}{\partial σ}

\dot{R} = [k (1 + m p) - m (R - R_{\infty})] \dot{p}

\dot{p} = \dot{λ}

\dot{f} = \frac{\partial f}{\partial σ} \dot{σ} + \frac{\partial f}{\partial R} \dot{R} = 0

f = J_{2} (σ) - R (p) - σ_{y p}

(11)

The formula is $R (p) = k p + R_{\infty} [1 - \exp (- m p)]$ . Both linear and exponential laws characterize the nonlinear isotropic hardening. The material yield strength is denoted by $σ_{y p}$ . The constants k, m, and $R_{\infty}$ depend on the material. The above-mentioned theoretical framework offers a fundamental method for precisely simulating and enhancing the mechanical properties of composites reinforced with short fibers. However, the accurate identification of characteristics including fiber orientation, elastic modulus, yield qualities, and hardening behavior is crucial to the validity and practicality of these theories.

Design of computational RVE experiments

Selection of input features

The behavior of glass fiber reinforced polypropylene (PPGF) composite under loading conditions is dependent on various material and structural properties that define an RVE model. Amongs these properties, elasticity modulus (E_m) and Poisson’s ratio ( $ν_{m}$ ) of the matrix, elasticity modulus (E_f) and Poisson’s ratio ( $ν_{f}$ ) of the fiber, fiber volume fraction ( $V_{f}$ ), length to diameter or fiber aspect ratio ( $L_{f} / d_{f}$ ) and fiber orientation tensor components ( $a_{11}$ , $a_{22}$ , $a_{33}$ ) are crucial in understanding the nonlinear response of composites at microstructural level. These properties/parameters can be utilized to study the resulting homogenous response of an RVE with computationally efficient ML techniques.

The range for each RVE parameter is chosen based on the minimum and maximum values observed in commonly used PPGF composite configurations. This selection ensures that the range encompasses various types of PPGF composites, specifically those with

V_{f}

ranging from 0.10% to 0.60%, covering materials from PPGF-10 to PPGF-60. The elasticity modulus of matrix and fiber is considered within the ranges 1300

\leq

E_m

\leq

1800 and 70,000

\leq

E_f

\leq

95,000. Whileas the Poisson’s ratio of matrix and fiber ranges between 0.30

\leq

ν_{m} \leq

0.45 and 0.15

\leq

ν_{f} \leq

0.26, respectively. The fiber aspect ratio vary between 10

\leq

L_{f} / d_{f} \leq

30, where its low value indicates shorter and thicker fibers, and high value indicates longer and thinner fibers. The summary of the training range of each parameter is presented in Table 1.

Table 1.

Training range of input parameters defined with their respective upper and lower limits.

	E_m (MPa)	$ν_{m}$	E_f (MPa)	$ν_{m}$	$V_{f}$ (%)	$L_{f} / d_{f}$	a ₁₁	$a_{22}$	$a_{33}$
Min	1300	0.30	70,000	0.15	0.10	10	0	0	0
Max	1800	0.45	95,000	0.26	0.60	30	1	1	1

Features sampling technique for input parameters

The quantity and quality of sample data are crucial for reliable conclusions drawn from simulations, machine learning, and experimental design. A regular grid of sample points can lead to overlap in high-dimensional spaces. To avoid this, a well-balanced and uniform sampling technique is needed. This study uses the Sobol sampling⁶² algorithm, a low-discrepancy quasi-random method that ensures uniform coverage of high-dimensional feature spaces. Sobol sampling is effective for sensitivity analysis and optimization, offering better uniformity than random methods while remaining computationally efficient.^63,64

The Sobol sampling method was used to generate a well-distributed design of experiments, ensuring comprehensive coverage of the design space. A base-2 quasi-random sequence with an exponent of 11 was employed, resulting in 2¹¹ = 2048 sample points. Each sample represents a unique combination of input parameters, that is, $(E_{m} {, ν}_{m}, E_{f}, ν_{f}, V_{f}, L_{f} / d_{f}, a_{11}, a_{22}, a_{33})$ . During the sample generation process, the sum of the orientation tensor components $(a_{11}, a_{22}, a_{33})$ must be equal to 1. To ensure this constraint, an adjustment is applied such that if their sum exceeds 1, each component is scaled down proportionally by dividing it by the total sum. Conversely, if the sum is less than 1, the components are increased proportionally to reach the required value. The pseudocode for this adjustment is provided in Figure 2. This approach ensures that the sampled data includes diverse orientation states, such as nearly uniform 3D random distribution (≈0.333, ≈0.333, ≈0.333), planar orientations (≈0.5, ≈0.5, ≈0), and unidirectional alignments along each principal direction (≈1, ≈0, ≈0), (≈0, ≈1, ≈0) and (≈0, ≈0, ≈1).

Figure 2.

Ajustments of the components of orientation tensor.

A discrepancy of 6.15 × 10⁻⁵ was calculated to assess how evenly the generated samples are distributed across the parameter space. This low-discrepancy value indicates minimal clustering or gaps, ensuring an efficient exploration of the design space. Additionally, Principal Component Analysis (PCA) is used in this study to reduce the dimensionality of high-dimensional data and visualize its distribution. PCA is an unsupervised technique that transforms correlated variables into a new set of uncorrelated variables through an orthogonal transformation. This allows the data to be represented in a lower-dimensional space while preserving as much variance as possible.⁶⁵ Here, PCA reduces the original 9-dimensional dataset to 2 dimensions, making it easier to assess the spread and uniformity of the Sobol-sampled data. This step ensures that the dataset is well-distributed before training ML models. A scatterplot of the dataset after PCA transformation is illustrated in Figure 3. The first principal component represents the most variance, while the second represents the next highest variance. The even spread of points in the plot indicates that the dataset is well-sampled across the entire feature ranges.

Figure 3.

Dimensionality reduction of input features using PCA.

Generation of target data and preprocessing

The effective elastoplastic stress-strain curve of each PPGF RVE sample is determined using Mori–Tanaka homogenization theory implemented in DIGIMAT software. The material expert (MX) module available in DIGIMAT contains variety of materials calibrated with J2 plasticity models. The database of these materials are usually uploaded by different reliable industries. Nonetheless, each RVE analysis accounted for both microstructural characteristics and material properties of the fiber and matrix in all DOE samples. A 3% loading strain was applied to each RVE to evaluate its mechanical response under deformation. To efficiently handle multiple simulations, a Python script was written to automate the workflow. This script facilitated key tasks, including extracting input parameters from the DOE, generating RVEs, performing homogenization analyses, and storing the computed stress-strain curves as effective material response. Automating this process significantly reduced manual effort, improved consistency, and ensured a systematic evaluation of the homogenized mechanical behavior across different samples.

The homogenized stress-strain curve derived from the RVE analysis for a specific material sample incorporates plasticity in it. Therefore, to enhance the training of ML models, each stress-strain curve is decomposed into two distinct components, that is, the linear elastic region and the nonlinear plastic region. The linear elastic portion is characterized by a single scalar value of elasticity modulus (E_H), which represents the homogenized material’s stiffness. Beyond the elastic limit, the nonlinear plastic region is more complex and requires a detailed representation. To achieve this, the plastic part of the stress-strain curve is interpolated at 11 evenly spaced points between 0% and 3% plastic strain, using an interval of 0.03%, as illustrated in Figure 4. This approach ensures that the ML model is trained to predict both the elasticity modulus and the stress values at these 11 interpolation points, effectively capturing the material’s transition from elasticity to plasticity. By decomposing the response in this way, the ML model gains a clearer understanding of both stiffness and nonlinear behavior, leading to improved accuracy and generalization when predicting the homogenized mechanical properties of new samples.

Figure 4.

Decomposition of nonlinear stress-strain curve and comparison of interpolated points with the simulated curve.

In this study, the Pearson correlation coefficient (PCC) was used to measure the linear relationship between input parameters and the homogenized elasticity modulus (E_H). PCC values range from −1 to 1, where a value above 0.5 indicates a significant positive correlation, meaning both variables increase together. A value below −0.5 represents a strong negative correlation, where one variable increases as the other decreases. A PCC of 0 means no linear correlation, while values between −0.5 and 0.5 suggest a weaker or less noticeable relationship. The PCC is mathematically expressed as follows:

PCC = \frac{\sum_{i = 1}^{n} (x_{i} - \hat{x}) (y_{i} - \hat{y})}{\sqrt{\sum_{i = 1}^{n} {(x_{i} - \hat{x})}^{2} \sqrt{\sum_{i = 1}^{n} {(y_{i -} \hat{y})}^{2}}}}

(12)

where x_i and y_i are individual sample points from two parameters. The

\hat{x}

and

\hat{y}

are the means of the respective parameters. The PCC heatmap in Figure 5 provides insight into how different RVE properties influence the homogenized elasticity modulus (E_H). The fiber volume fraction (V_f) exhibits the strongest positive correlation with E_H at 0.7, indicating that increasing the fiber content significantly enhances the composite’s stiffness. Additionally, the a₁₁ component of the orientation tensor shows the second-highest positive correlation (0.57) with E_H, suggesting that a higher alignment of fibers in the primary direction contributes to greater stiffness. In contrast, the a₂₂ and a₃₃ components of the orientation tensor have negative correlations of −0.29 with E_H, implying that fiber orientations in these directions tend to reduce the overall stiffness of the material. This figure also reveals that most RVE properties, that are considered independent input parameters, exhibit correlation values close to zero with each other. This indicates that there are no significant linear relationships among these parameters, ensuring that each variable contributes uniquely to the dataset.

Figure 5.

Heatmap for the PCC between input parameters and homogenized elasticity modulus.

In machine learning, several data normalization or scaling methods are commonly used to improve model performance and training stability. These include min-max scaling, which transforms values to a range between 0 and 1, and standardization, which re-centers data around the mean with unit variance. Other techniques involve scaling data to lie within −1 and 1 or normalizing each sample’s feature values so they sum to 1. These methods help address issues like varying feature magnitudes and skewed distributions, ultimately ensuring more effective training and generalization. In this work, the dataset was scaled in a range of 0 to 1 using the min-max scaling formulation represented as,

x (s c a l e d) = \frac{x - x (\min)}{x (\max) - x (\min)}

(13)

Here, x is the actual value, and x (max) and x (min) are the maximum and minimum values of a feature, respectively. After applying normalization, the dataset was split into training and testing sets to ensure effective model learning and evaluation. Approximately 80% of the data that consists of 1638 samples, was allocated for training, while the remaining 20%, comprising 409 samples, was set aside for testing. Model performance is assessed on the previously unseen data (i.e., testing set) to ensure that the predictions are not biased by the training data. In addition, model generalization and overfitting are monitored using evaluation metrics such as R², MAE, MAPE, and RMSE on both the train and the test set.

Theoretical background of various ML models

In this study, four ML models are evaluated to determine the most effective algorithm for predicting the homogenized elastoplastic stress-strain response of PPGF composite materials. The selected models include Random Forest (RF), Extreme Gradient Boosting (XGB), Artificial Neural Network (ANN), and Support Vector Regression (SVR). A brief overview of each model is provided in the following sections.

Random forest (RF)

RF is an ensemble learning method that combines multiple decision trees to create a more robust predictive model. Instead of relying on a single decision tree, RF builds many trees using different random subsets of the dataset, then averages their predictions. This approach, known as bootstrap aggregating or bagging, helps reduce overfitting and makes the model more resilient to noise and outliers. The RF prediction can be expressed as

\hat{f} (x) = \frac{1}{n_{T}} \sum_{i = 1}^{n_{T}} f_{i} (x^{'})

(14)

where

\hat{f} (x)

represents the predictions from the total number of decision trees (n_T),

f_{i} (x^{'})

is the prediction of each individual tree and

x^{'}

signifies a random selection of features. Unlike individual decision trees, which can be sensitive to variations in data, RF implements pruning techniques to provide more stable predictions and prevent overfitting. Due to its ability to handle high-dimensional data and maintain strong performance on unseen samples, RF is widely used for regression tasks across various fields.

Extreme gradient boosting (XGB)

XGB is an advanced version of the traditional Gradient Boosting algorithm, designed for high efficiency, scalability, and predictive accuracy. It constructs an ensemble of decision trees by optimizing a loss function using gradient descent while incorporating both L1 (Lasso) and L2 (Ridge) regularization to enhance generalization and prevent overfitting. The algorithm iteratively builds trees, refining predictions by minimizing errors from previous iterations. A key advantage of XGB is its ability to reduce processing time through optimized computations, making it suitable for both regression and classification tasks. The inclusion of a regularization term in the objective function helps control model complexity, reducing the risk of overfitting and improving overall performance. Its objective function is represented as follows:

{o b j}^{(m)} = \sum_{i = 1}^{n} L_{X G B} (y_{i}, {\hat{y}}_{i}^{(m)}) + \sum_{t}^{m} Ω (f_{t})

(15)

here,

L_{X G B}

is the loss function that is dependent on actual value

y_{i}

and predicted value

{\hat{y}}_{i}^{(m)}

Ω (f_{t})

is a regularization term to prevent overfitting.^66,67 Due to its speed, robustness, and ability to handle large datasets, XGB has become one of the most widely used machine learning algorithms in various fields.

Artificial neural network (ANN)

ANNs are computational models inspired by the structure and functioning of the human brain that are designed to recognize patterns from complex data. ANNs consist of interconnected layers of artificial neurons, including an input layer that receives data, one or more hidden layers that process information, and an output layer that generates predictions. Each neuron processes incoming signals using weighted connections, applies an activation function, and passes the result to the next layer. Through iterative learning, often using backpropagation and optimization techniques like stochastic gradient descent, ANNs adjust their weights to minimize prediction errors. Due to their ability to model nonlinear relationships and extract meaningful patterns from raw data, ANNs are widely used in diverse fields such as image and speech recognition, financial forecasting, medical diagnosis, and engineering applications. ANNs have gained significant popularity with the rise of deep learning architectures, demonstrating state-of-the-art performance in tasks involving large-scale and high-dimensional datasets.⁶⁸

Support vector regression (SVR)

SVR is a machine learning technique designed for predicting continuous values in regression tasks. Unlike traditional regression methods that minimize the error between predicted and actual values, SVR aims to find a hyperplane that best fits the data where the margin between the hyperplane and the support vectors (data points closest to the hyperplane) is maximized. It uses kernel functions (e.g., linear, polynomial, or radial basis function) to map input data into higher-dimensional spaces, enabling the modeling of complex nonlinear relationships. The optimization problem in SVR can be formulated as

\min_{w, b, ϵ} \frac{1}{2} {‖ w ‖}^{2} s u b j e c t e d t o | y_{i} - (w^{T} x_{i} + b) | < ϵ

(16)

where w is the weight vector, b is the bias term,

ϵ

is the margin of tolerance, and y_i represents the true values. The ability of SVR to generalize well even with limited data makes it a popular choice for regression problems in various domains.

Description of the automobile component and its finite element modeling

The elastoplastic material model predicted through ML-based homogenization is validated on a real-world macroscopic structure by implementing it into the finite element analysis of an automotive component. Specifically, the validation is performed on the cover part of the top mount, which serves as a crucial link between the vehicle’s chassis and suspension strut (Figure 6). The top mount plays a vital role in enhancing ride comfort by absorbing road-induced impacts and dampening vibrations, preventing them from being transmitted to the vehicle. As a key component of the suspension system, it contributes to the overall stability and performance of the automobile.

Figure 6.

Cover part of the top mount used for validation model in this study.

The structural model of the cover component was imported into Abaqus software with its actual geometry and dimensions. In the finite element modeling process, the part was assigned the homogenized elastoplastic material properties obtained from the ML model. To simulate realistic boundary conditions, the bottom surface was fully constrained in all degrees of freedom (DOF 1, 2, 3). A vertical load of 16 kN was applied to the top face, while lateral movements were restricted in the remaining degrees of freedom (DOF 1, 2), as illustrated in Figure 7.

Figure 7.

Finite element model of the cover part with boundary and loading conditions.

To ensure an optimal balance between computational accuracy and efficiency, a mesh convergence study was carried out using 10-node quadratic tetrahedral elements (C3D10). As illustrated in Figure 8, both the maximum von Mises stress and displacement values exhibit high sensitivity to mesh refinement between 6 mm and 5 mm, with the maximum stress increasing by approximately 18%. This substantial variation highlights the inadequacy of coarser mesh in capturing the mechanical response with sufficient accuracy. Conversely, for mesh sizes finer than 5 mm (i.e., 4 mm, 3 mm, and 2 mm), the relative changes in stress and displacement drop below 1%, indicating convergence of the solution. This behavior is supported by Figure 8(b), which presents the percentage change between successive mesh refinements. Based on a convergence threshold of 1% relative change, a mesh size of 4 mm was selected as the optimal discretization, offering the suitable balance between numerical precision and computational cost.

Figure 8.

Mesh sensitivity analysis of the automobile component with (a) variation of maximum stress and displacement with different mesh sizes ranging from 7 mm to 2 mm and (b) relative percent change in maximum stress and displacement between successive mesh refinements.

Modeling results and discussions

Optimal model confgurations

The ML models were implemented in Python 3.11 and trained on a high-performance system equipped with a 12th Gen Intel® Core™ i9-12,900H processor (2.5 GHz) and 64 GB of RAM. A comprehensive hyperparameter tuning process was conducted for each ML model using an iterative trial-and-error approach. This method was selected due to its flexibility and ease of interpretability, allowing us to manually adjust hyperparameters based on model performance trends observed during validation. Additionally, given the small size of the training dataset and the limited number of tunable hyperparameters, this approach provides a computationally efficient alternative than resource-intensive optimization methods such as Bayesian optimization, particle swarm optimization, or grid/random search. To ensure fair comparisons and reproducibility, a seed value of 20,250,201 was used throughout the tuning process. The final set of optimized hyperparameters for each model is presented in Table 2, while any unspecified parameters were left at their default settings. This table presents the optimal hyperparameter values for the ML models used to predict the homogenized elasticity modulus and nonlinear plastic stress.

Table 2.

Optimal hyperparameter values of the ML models.

ML model	Hyperparameter	Hyperparameter range	Optimal value (elasticity modulus)	Optimal value (plastic stress)
Random forest	Maximum tree depth	10–300	22	210
	Number of estimators	10–100	30	15
	Minimum samples leaf	1–5	1	1
Extreme gradient boosting	Maximum tree depth	1–10	7	6
	Learning rate	0.01–0.5	0.2	0.07
	Number of estimators	100–500	400	380
Artificial neural network	Number of hidden layers	2–4	3	3
	Number of neurons/layer	16–64	64/64/32	64/64/64
	Activation function	Relu, Sigmoid, Tanh	Sigmoid	Sigmoid
	Learning rate	0.001–0.1	0.001	0.005
	Iterations	500–3000	1000	1000
Support vector regression	Kernel function	Linear, Poly, RBF	RBF	RBF
	Epsilon	0.001–0.1	0.005	0.006
	Penalty parameter (C)	1–5	3	5

Prediction performance evaluation of all trained models

The comprehensive performance evaluation of the implemented ML models for predicting the homogenized elasticity modulus are provided in Table 3. This table reveals that among all models, ANN performs exceptionally well in generalization by exhibiting minimal performance drop between training and testing, while maintaining high accuracy. On testing dataset, ANN achieved the highest R² value (0.9995), along with the lowest errors across MAE, MAPE, and RMSE metrics valuing at 56.89, 0.0094, and 80.62, respectively. This suggests that ANN not only learns the complex micomechanical relationships effectively but also maintains strong predictive accuracy on unseen data. XGB also demonstrates strong performance during training but a substantial drop in accuracy is observed in the testing phase, indicating potential overfitting. Similarly, SVR exhibits higher error values in testing and RF records the highest error values across all performance metrics. These results suggest that XGB, RF, and SVR struggle to maintain predictive accuracy and indicate their poor generalization ability to capture intricate relationships compared to ANN.

Table 3.

Performance evaluation of models predicting the homogenized elasticity modulus.

Models	Training dataset (80%)				Testing dataset (20%)
Models	R²	MAE	MAPE	RMSE	R²	MAE	MAPE	RMSE
RF	0.998	107.61	0.016	171.69	0.982	292.03	0.043	474.46
XGB	0.9999	28.89	0.0052	40.62	0.9892	226.79	0.0327	375.89
ANN	0.9995	54.08	0.0089	77.37	0.9995	56.89	0.0094	80.62
SVR	0.9993	79.46	0.0152	92.91	0.9979	105.07	0.0187	165.09

Similarly, the performance of the trained ML models for predicting the homogenized plastic stresses are provided in Table 4. The results demonstrate that all ANN performed exceptionally well in predicting the plastic stresses, with highest R² value of 0.9979 and lowest errors across all metrics on testing dataset. ANN model again demonstrated the excellent generalization capability by exhibiting minimal performance drop between training and testing phases. In contrast, XGB, RF and SVR showed higher errors compared to ANN on the testing dataset. Moreover, these models have poor generalization ability to capture nonlinear relationships compared to ANN. Overall, ANN emerge as the most reliable model in terms of demonstrating both high accuracy and generalization for predicting the homogenized behavior of composite materials.

Table 4.

Performance evaluation of models predicting the homogenized plastic stresses.

Models	Training dataset (80%)				Testing dataset (20%)
Models	R²	MAE	MAPE	RMSE	R²	MAE	MAPE	RMSE
RF	0.998	0.84	0.011	1.49	0.984	2.23	0.030	3.71
XGB	0.9995	0.42	0.0063	0.59	0.9931	1.47	0.0196	2.44
ANN	0.9981	0.77	0.0111	1.1	0.9979	0.81	0.0114	1.19
SVR	0.9984	0.92	0.0145	1.11	0.9951	1.31	0.0194	2.0

The learning curves for the ANN model presented in Figure 9(a) and (b) illustrate the training and validation loss evolution for predicting the elasticity modulus and nonlinear plastic stress, respectively. In both cases, the loss values exhibit a consistent decreasing trend across 1000 iterations, indicating effective learning and convergence of the model. The y-axis, plotted on a logarithmic scale of mean squared error (MSE) in MPa, shows that the errors decrease from an initial order of 10⁻² to approximately 10⁻⁵ by the end of the iterations, suggesting a substantial reduction in prediction error over time. Notably, the validation loss curves closely overlap with the training loss curves throughout the process, which is a positive sign of good generalization and an absence of overfitting. If any indications of overfitting had emerged, appropriate mitigation techniques such as early stopping, L2 regularization, or dropout would have been considered to enhance the model’s generalization performance. The presence of some noise in the curves, particularly in the validation loss, could indicate minor instability or sensitivity to the data, though this does not appear to undermine the overall trend. Overall, these curves suggest that the ANN models are well-trained and robust for predicting both elasticity modulus and nonlinear plastic stress.

Figure 9.

ANN learning curves for (a) elasticity modulus and (b) nonlinear plastic stress.

The scatter plots in the Figure 10 compare the predicted homogenized elasticity modulus E_M against the actual values for four trained ML models including RF, XGB, ANN, and SVR. In these plots, the 45° diagonal line (y = x) called as ideal fit line represents the perfect prediction accuracy if the data points lie very close to it with no scattering. Among the models, ANN exhibits the highest predictive accuracy as its data points align almost perfectly with the ideal fit line, with all predictions falling within the ±5% error bounds. The density of predicted points on the ideal fit line highlights the consistency and reliability of ANN model in predicting the E_M. SVR also performs well with predictions closely following the ideal fit line and only a few outliers slightly exceeding the 5% error margins, particularly at higher E_M values. In contrast, RF and XGB demonstrate greater dispersion with a significant number of data points deviating beyond the error bounds. The higher scattering in RF and XGB predictions highlights their limitations in capturing the underlying patterns in the data, leading to reduced predictive accuracy.

Figure 10.

Comparison of predicted versus actual homogenized elastic modulus (E_H) for different ML models: (a) RF, (b) XGB, (c) ANN, and (d) SVR. The black line denotes the ideal fit (y = x), while dashed lines represent ±5% error bounds to assess prediction accuracy. The predictions are based on the full set of RVE input parameters: $E_{m}$ , $ν_{m}$ , $E_{f}$ , $ν_{f}$ , $V_{f}$ , $L_{f} / d_{f}$ , $a_{11}$ , $a_{22}$ and $a_{33}$ .

Figure 11 presents plastic stress-strain curves predicted by the trained ML models alongside the actual curve derived from the Mori–Tanaka homogenization method. The microstructural parameters for obtaining these results are as follows: $E_{m}$ = 1519.2, $ν_{m}$ = 0.427, $E_{f}$ = 94,064.3, $ν_{f}$ = 0.212, $V_{f}$ = 0.35, $L_{f} / d_{f}$ = 15.7, $a_{11}$ = 0.36, $a_{22}$ = 0.32, $a_{33}$ = 0.32. Although it appears from this figure that ML the predicted stress-strain curves closely align with the actual curve, the zoomed section reveals the true and more accurate nature of results. It can be seen that there are subtle differences in the curves, with the ANN-predicted curve showing the closest match to the actual curve. This highlights the superior ability of ANNs to capture the intricate nonlinear relationships in the data. The SVR curve follows as the next closest, while the XGB and RF curves exhibit slightly greater deviations particularly at higher strain levels, indicating a descending order of accuracy (ANN > SVR > XGB > RF). Overall, the figure demonstrates that ML models can effectively predict the homogenized plastic behavior of the composite material, with some models (ANN in this case) achieving near-exact alignment with the physically derived Mori–Tanaka solution.

Figure 11.

Comparison of actual and ML predicted plastic stress-strain curves for short fiber-reinforced composites. An inset zoom highlights differences in the high-strain regime.

Computational performance of ANN model

In addition to accuracy and generalization, computational efficiency is a crucial factor in evaluating the practicality of homogenization methods. Therefore, the prediction runtime of the proposed ML-based homogenization approach with the conventional Mori–Tanaka mean-field method is compared. The ANN training was performed using a workstation equipped with an 12th Gen Intel® Core™ i9-12,900H processor, 64 GB RAM, and an NVIDIA RTX A2000 GPU. The complete training process of the ANN model using the Adam optimizer and network architecture presented in Table 2 required about 12 minutes to converge for 2048 RVE-based samples. Once the ANN model is trained, it can predict the homogenized material properties of SFRC in less than one second per sample, owing to their purely algebraic forward-pass computation. In contrast, the Mori–Tanaka method requires constructing an RVE, which involves defining the microstructural geometry, assigning material properties, and performing simulations. This process typically takes approximately 3–4 minutes per sample, depending on the complexity of the RVE and computational resources. Therefore, while the training phase of the ML model involves a one-time computational cost, it is orders of magnitude more efficient than running large number of high-fidelity simulations. This comparison highlights the significant runtime advantage of ML-based homogenization. The ability to obtain instantaneous predictions without the need for repeated simulations makes the ML approach particularly well-suited for parametric studies, uncertainty quantification, and real-time structural analysis workflows.

Taguchi and ANOVA method for the evaluating SFRC performance

This section employs the Taguchi method to analyze the behavior of a SFRC materials during homogenization and implements ANOVA technique to investiagte the signicance of various RVE parameters on the response of interest. The Taguchi method is a robust statistical approach that utilizes orthogonal arrays to systematically design experiments. This allows for the evaluation of multiple factors with a significantly reduced number of trials compared to traditional full factorial designs.⁶⁹ By strategically selecting parameter combinations, this method ensures a balanced and efficient assessment of factor effects while minimizing experimental costs and material waste.^70,71 Taguchi method is widely adopted in manufacturing industries to enhance the process optimization and quality improvement by identifying optimal conditions with minimal resources.

In this study, nine RVE parameters, that is, (

E_{m} {, ν}_{m}, E_{f}, ν_{f}, V_{f}, L_{f} / d_{f}, a_{11}, a_{22}, a_{33}

) were selected as factors. Each factor was assigned three levels, that is, minimum, center and the maximum value of parameters provided in Table 1. Specifically, the levels are set as follows:

E_{m}

= 1300/1550/1800,

ν_{m}

= 0.3/0.4/0.45,

E_{f}

= 70,000/82,000/95,000,

ν_{f}

= 0.15/0.2/0.26,

V_{f}

= 0.1/0.35/0.6,

L_{f} / d_{f}

= 10/20/30,

a_{11}

a_{22}

a_{33}

= 0/0.33/0.66/1. Typically, a full factorial design that considers all possible parameter combinations is calculated using L^H formula, where L represents the levels and F represents number of factors. This design would require 3⁹ = 19,638 trials which demanding extensive computational resources and time. To overcome this challenge, we employed the L₂₇ Taguchi orthogonal array in this work with nine factors with three levels each. This orthogonal array requires only 27 trials, while significantly reducing the computational costs. The generated L₂₇ orthogonal arrays are presented in Table 5.

Table 5.

Taguchi L₂₇ (3⁹) orthogonal arrays and the calculated homogenization response of RVE.

Trial	Factors (parameters)									RVE response
Trial	E_m (MPa)	$ν_{m}$	E_f (MPa)	$v_{f}$	$V_{f}$ (%)	$L_{f} / d_{f}$	a ₁₁	a ₂₂	a ₃₃	$E_{H}$ (MPa)	$σ_{y i e l d}$ (MPa)
T1	1300	0.3	70,000	0.15	0.1	10	1	0	0	3523	18.8
T2	1300	0.3	70,000	0.2	0.35	20	0	1	0	2513	13.5
T3	1300	0.3	70,000	0.26	0.6	30	0	0	1	4404	15.2
T4	1300	0.4	82,000	0.15	0.1	20	0.66	0.33	0	3520	20.7
T5	1300	0.4	82,000	0.2	0.35	30	0.66	0	0.33	12,280	53.6
T6	1300	0.4	82,000	0.26	0.6	10	0.33	0.66	0	8435.6	24
T7	1300	0.45	95,000	0.15	0.1	30	0.33	0	0.66	2763.6	18.2
T8	1300	0.45	95,000	0.2	0.35	10	0	0.66	0.33	3786.8	18.1
T9	1300	0.45	95,000	0.26	0.6	20	0	0.33	0.66	8190	29.4
T10	1550	0.3	82,000	0.15	0.35	10	0.33	0.33	0.33	4562.4	16.8
T11	1550	0.3	82,000	0.2	0.6	20	0.66	0.33	0	21,612	55.1
T12	1550	0.3	82,000	0.26	0.1	30	0.33	0.66	0	2927.3	17.1
T13	1550	0.4	95,000	0.15	0.35	20	0.33	0	0.66	6694.8	25.6
T14	1550	0.4	95,000	0.2	0.6	30	0	0.66	0.33	7699.6	26.6
T15	1550	0.4	95,000	0.26	0.1	10	0.66	0	0.33	3027.6	15.5
T16	1550	0.45	70,000	0.15	0.35	30	0	0.33	0.66	5252	25.2
T17	1550	0.45	70,000	0.2	0.6	10	0.66	0	0.33	14,741	33.7
T18	1550	0.45	70,000	0.26	0.1	20	0.33	0.66	0	2675.4	15.4
T19	1800	0.3	95,000	0.15	0.6	10	0	1	0	6088	14.8
T20	1800	0.3	95,000	0.2	0.1	20	1	0	0	7245.2	28.4
T21	1800	0.3	95,000	0.26	0.35	30	0	0	1	3480.4	13.3
T22	1800	0.4	70,000	0.15	0.6	20	0.33	0.33	0.33	10,998	27.9
T23	1800	0.4	70,000	0.2	0.1	30	0.66	0.33	0	4583.2	22.8
T24	1800	0.4	70,000	0.26	0.35	10	0	0.66	0.33	4270	17.0
T25	1800	0.45	82,000	0.15	0.6	30	0.33	0.66	0	14,037	38.1
T26	1800	0.45	82,000	0.2	0.1	10	0.66	0	0.33	3365.4	14.9
T27	1800	0.45	82,000	0.26	0.35	20	0	0.33	0.66	5695.2	22.2

These orthogonal arrays are utilized for conducting the homogenization analyses of short fiber reinforced composite materials. In the homogenization analysis, the elastic modulus (E_H) and yield stress ( $σ_{y i e l d}$ ) are selected as key response parameters to evaluate the macroscopic mechanical behavior of the composite system. Since SFRC materials exhibit anisotropic mechanical properties influenced by fiber orientation, volume fraction, and matrix stiffness, E_H helps quantify how these microstructural factors contribute to the overall stiffness of the homogenized composite. The yield stress ( $σ_{y i e l d}$ ) is another crucial parameter that is affected by the microstructural factors. A higher yield stress indicates improved load-bearing capacity and resistance to plastic deformation, which is essential for structural applications requiring high strength and durability. These homogenization results are utilized for carrying out the ANOVA (Analysis of Variance) in subsequent section.

A statistical tool like ANOVA is used to analyze the effect of each parameter in a Taguchi L27 orthogonal array on a response variable. This will help determine how each of the nine RVE parameters influences the elastic modulus ( $E_{H}$ ) and the yield stress ( $σ_{y i e l d}$ ) of the composite. ANOVA works by decomposing the total variability in the response into contributions from different parameter and the residual error. The key formulations for this analysis involves calculating the total sum of squares ( ${S S}_{t o t a l}$ ), sum of squares for each parameter ( ${S S}_{P}$ ), error sum of squares ( ${S S}_{e r r o r}$ ), degrees of freedom (DoF), mean square (MS) and F-value ( $F_{v}$ ), which are formulated as follows:

{S S}_{t o t a l} = \sum_{i = 1}^{T r i a l s} {(Y_{i} - \bar{Y})}^{2}

{S S}_{P} = \sum_{i = 1}^{L e v e l s} {n_{j} ({\bar{Y}}_{j} - \bar{Y})}^{2}

{S S}_{e r r o r} = {S S}_{t o t a l} - \sum_{i = 1}^{F a c t o r s} {S S}_{P}

{M S}_{P} = \frac{{S S}_{P}}{{D o F}_{P}}, {M S}_{e r r o r} = \frac{{S S}_{e r r o r}}{{D o F}_{e r r o r}}

F_{v} = \frac{{M S}_{P}}{{M S}_{e r r o r}}

where

Y_{i}

is the observed response at ith trial,

\bar{Y}

is the mean response,

n_{j}

is the number of trials at level j,

{\bar{Y}}_{j}

is the mean response at level j and DoF = Levels - 1 is the degree of freedom for the factor. In ANOVA analysis, the sum of squares represents the variation in the response attributed to a specific parameter. Its value indicates a stronger influence of that parameter on the response. The F-value is the ratio of the variance explained by the parameter to the unexplained variance (error). A higher F-value suggests that the parameter has a more significant effect relative to the residual error. The p-value quantifies the probability that the observed response is due to random chance. A lower p-value indicates statistical significance, meaning the parameter has a meaningful impact on the response variable. In this work, the ANOVA analysis was performed using a significance level of 5% with 95% confidence interval. This means if the p-value is less than 0.05, the parameter is considered to have a statistically significant effect on the response variable.^72,73

The ANOVA results in Table 6 provide valuable insights into the parameters influencing the homogenized elastic modulus (

E_{H}

) and the yield stress (

σ_{y i e l d}

) of the nonlinear elastoplastic response obtained through RVE homogenization based on the Taguchi orthogonal arrays. For

E_{H}

, the fiber volume fraction (

V_{f}

) is the most influential parameter with a very low p-value (0.000,038) and the highest sum of squares, indicating that it significantly governs the homogenized stiffness of the SFRC. Additionally, among the orientation tensor components, a₁₁ is highly significant with p-values of 0.001,571 followed by a₂₂, and a₃₃ whose p-values are 0.007,452, and 0.011,040, respectively. This highlights the critical role of fiber alignment in determining the overall elastic response. The fiber length-to-diameter ratio (

L_{f} / d_{f}

) also show significant effects on the homogenized stiffness with p-value of 0.029,340. In contrast, matrix elasticity modulus (E_m), fiber elasticity modulus (E_f ), and fiber Poisson’s ratio (

v_{f}

) do not exhibit a statistically significant influence on

E_{H}

Table 6.

ANOVA analysis results for the homogenized elasticity modulus ( $E_{H}$ ) and yield stress ( $σ_{y i e l d}$ ).

Parameters	DoF	ANOVA for $E_{H}$			ANOVA for $σ_{y i e l d}$
Parameters	DoF	Sum of squares	F-value	p-value	Sum of squares	F-value	p-value
E _m	2.0	1.29709e+07	3.13789	0.116763	58.32089	1.49810	0.296671
$ν_{m}$	2.0	3.05492e+07	7.39036	0.024070	122.4767	3.14609	0.116296
E _f	2.0	7.61329e+06	1.84178	0.237874	36.13145	0.92811	0.445461
$v_{f}$	2.0	3.77449e+06	0.91311	0.450605	51.30074	1.31777	0.335415
$V_{f}$	2.0	3.57770e+08	86.55083	0.000038	1445.143	37.12183	0.000418
$L_{f} / d_{f}$	2.0	2.78061e+07	6.726782	0.029340	491.3916	12.62252	0.007081
a ₁₁	3.0	1.24451e+08	20.07123	0.001571	1103.588	18.89881	0.001847
a ₂₂	3.0	6.83480e+07	11.02303	0.007452	331.1951	5.67167	0.034747
a ₃₃	3.0	5.82087e+07	9.38778	0.011040	523.0812	8.95769	0.012354
Residual	6.0	1.24009e+07			116.7892

For the $σ_{y i e l d}$ of the nonlinear elastoplastic response, $V_{f}$ again emerges as the dominant parameter that enhances the yield stress with F-value of 37.12 and p-value of 0.000,418. The fiber length-to-diameter ratio ( $L_{f} / d_{f}$ ) and fiber orientation (a₁₁, a₂₂, and a₃₃) are also significant with p-values of 0.007,081, 0.001,847, 0.034,747 and 0.012,354, respectively. This indicates that fiber dispersion and alignment significantly impact the material’s load-bearing capacity. Other parameters, including E_m, E_f, $ν_{m}$ , and $v_{f}$ does not show a statistically significant effect. This finding suggests that the yield stress of the homogenized material is primarily governed by geometric and volume-related fiber characteristics rather than intrinsic material properties.

Overall, these results align with the fundamental principles of composite material behavior, where fiber volume fraction, distribution and orientation dictate the overall mechanical performance, both in terms of stiffness and elastoplastic strength. It can be noted that the statistical significance of the fiber elasticity modulus (E_f) in the ANOVA results is low, although it is known that E_f physically contributes to composite stiffness. The apparent statistical insignificance arises primarily from the specific range of input parameters considered in this study. The variations in E_f were relatively narrow compared to those of other parameters such as fiber V_f and $a_{11}, a_{22}, a_{33}$ . As a result, the influence of E_f on the homogenized modulus (E_H) was not dominant compared to parameters with broader variability and stronger effects. This finding is further supported by the Pearson Correlation Coefficient (PCC) matrix presented in Figure 5, where the correlation between E_f and E_H is notably weak (PCC ≈0.06). This confirms that variations in E_f had only a marginal impact on the predicted E_H within the defined RVE dataset. This is consistent with literature reports where composite behavior is more sensitive to volume fraction, fiber distribution, and orientation than to moderate variations in fiber stiffness. In high-stiffness fiber regimes (e.g., carbon or glass fibers), the effect of further increasing E_f may reach a plateau in terms of its impact on E_H, particularly when V_f is held constant.

Finite element analysis of automobile cover part based on different homogenization approaches

Three trials from the Taguchi L27 arrays presented in Table 6 are selected based on performance diversity, parameter coverage, and engineering relevance. The homogenization of these trials is carried out using the traditional Mori–Tanaka-based mean-field homogenization and ML-based homogenization techniques. Subsequently, homogenized material data from both approaches is implemented into Abaqus software to carry out the finite element analyses of an automobile cover part. The finite element results as per both homogenization schemes are compred to illustrate the accuracy of ML-based homogenization.

The trials selected for further finite element analysis of the cover part are T1, T5 and T11. Trial T1 represents the lowest stiffness (E_H = 3523 MPa) and yield strength (σ_yield = 18.8 MPa), reflecting a baseline material configuration (V_f = 0.1, E_f = 70,000 MPa). This trial is useful for cost-driven applications with lower material requirements. Trial T5 with an intermediate fiber volume fraction (V_f = 0.35) balances moderate stiffness (E_H = 12,280.6 MPa) with high strength (σ_yield = 53.6 MPa), highlighting the impact of high aspect ratio (L_f/d_f = 30) and fiber orientation (a₁₁ = 0.66). This trial captures a typical performance scenario and helps illustrate the composite behavior when reinforcement is at an average level. Trial T11 exhibits the highest stiffness (E_H = 21,612 MPa) and yield strength (σ_yield = 55.1 MPa) due to high fiber volume fraction (V_f = 0.6), high fiber elasticity modulus (E_f = 82,000 MPa) and favorable fiber alignment (a₁₁ = 0.66). This trial represents the upper bound of the composite performance that is ideal for simulating high-stress regions. The nonlinear elastoplastic stress-strain data is computed for the selected trials (T1, T5, and T11) using the Mori–Tanaka (MT) mean-field homogenization and ML-based homogenization techniques. The obtained homogenized data is then used to calibrate the plasticity model in Abaqus. Since ANN performed best among all the ML models tested in this study, it is used for predicting the homogenization response of the SFRC composites. The macroscopic mechanical behavior of the material is obtained by carryıng out finite element analysis of the automotive component with material data obtained from MT homogenization and ANN-based homogenization techniques.

As illustrated in Figure 12, the stress contour patterns of the cover part using the ANN-homogenized material exhibit a high degree of similarity to those obtained with the MT approach. This strongly supports the validity of the ANN model in capturing the material response effectively of the composite. Notably, stress concentration regions, particularly in areas of geometric discontinuities align well between the two approaches. In Trial T1, the MT homogenized method yields a maximum stress of approximately 99.21 MPa concentrated at the component’s legs, while the ANN-based method produces a nearly similar stress peak. Similarly, for Trial T5, both methods show a maximum stress around 171 MPa at the same critical regions, with minimum stresses in the central body. While in Trial T11, a maximum stress of about 168 MPa is observed in both methods. The slight variations observed in the stress values are within an acceptable range. This confirms the robustness of the ANN-based homogenization in accurately predicting the macroscopic stress response of composite materials.

Figure 12.

Stress contour plots of the cover part characterized with MT homogenized material model and ANN-homogenized material model for the three trials (i.e., T1, T5, and T11) selected from Taguchi orthogonal arrays.

Figure 13 illustrates the stress-strain curves obtained from finite element analyses using Mori–Tanaka (MT) mean-field homogenization and the ANN-based homogenization technique for Trials T1, T5, and T11. The figure shows that for each trial, the ANN-based predictions exhibits a strong correlation with the MT-based results, capturing both the overall stiffness and nonlinear stress-strain behavior. Notably, Trial 1, which has the lowest fiber volume fraction (V_f = 0.1) and aspect ratio (L_f/d_f = 10), exhibits the lowest stress levels across the strain range. In contrast, Trials 5 and 11 display significantly stiffer responses due to higher fiber content and fiber aspect ratios. The consistency between the ANN and MT results across different material parameter configurations highlights the ANN model’s ability to generalize well and accurately replicate the complex homogenized material behavior. Overall, this comparison validates the ANN-based homogenization technique as an efficient and reliable alternative to traditional mean-field methods without compromising accuracy in predictions.

Figure 13.

Comparison of stress-strain curves obtained from the finite element analysis of three trials with Mori–Tanaka (MT) and ANN homogenization technique.

Conclusion

This study explores the potential of various advanced ML techniques including random forest (RF) and extreme gradient boosting (XGBoost), artificial neural networks (ANNs) and support vector regression (SVR) in predicting the homogenized nonlinear elastoplastic behavior of short fiber-reinforced (SFRC) composite materials. The training dataset was systematically generated using the Sobol sampling method, considering a wide range of microstructural parameters such as fiber and matrix elasticity modulus, Poisson’s ratio, fiber volume fraction, fiber aspect ratio, and fiber orientation tensor components. Additionally, Taguchi and ANOVA analyses were employed to quantify the influence of different microstructural parameters on the homogenized response of SFRC. The results revealed the most significant factors affecting the macroscopic elastoplastic behavior. The predictive performance of the ML models was evaluated using multiple statistical metrics, including R², MAE, MAPE, and RMSE, along with visualization tools such as scatter plots and bar charts.

Among the models investigated, ANN demonstrated the highest predictive accuracy, with most predictions falling within a ±5% error range. These results underscore the ANN model’s ability to replace conventional homogenization techniques, offering a computationally efficient alternative while maintaining reliability in capturing the complex nonlinear behavior of composite materials. The effectiveness of ANN-based homogenization was further validated by integrating it into a finite element analysis of an automotive component, where the ANN-predicted material model produced results consistent with conventional numerical techniques. This close alignment validates the ANN-based homogenization’s capability to replicate the physics-driven Mori–Tanaka approach and its practical applicability in real-world engineering problems.

Overall, this study underscores the advantages of ML-driven approaches in predicting the homogenized behavior of composite materials, offering an accurate alternative to traditional numerical methods. The proposed methodology can facilitate rapid material design and optimization for engineering applications. In principle, the ML models can be extended to other types of composites such as particulate-reinforced, layered, or textile composites provided that a sufficiently diverse and representative dataset is available for the new material system and input features are adapted to capture the relevant micromechanical descriptors of the new composite (e.g., particle shape descriptors instead of fiber anisotropy). Similarly, while the trained models currently reflect predictions based on numerical homogenization using Mori–Tanaka, they could be adapted to predict results from alternative homogenization frameworks (e.g., self-consistent schemes, FE-based analysis, and FFT-based methods). However, this would require retraining the models on data generated via those methods. Future research can extend this work by exploring physics-informed architectures to enhance predictive capabilities further, as well as incorporating uncertainty quantification to improve the reliability of ML-based material models. It is also intended to extend this study by incorporating experimental results to further support the model validation.

Footnotes

Acknowledgments

The authors would like to thank Scientific and Technological Research Council of Turkey (TÜBİTAK) for their support within the scope of the International Cornet project (Project No: 123N452) and the SME user committee member YAMAS Rubber & Metal Solutions for providing the part information used in this research. We also sincerely thank TÜBİTAK for the scholarship awarded to the first author under the 2244 - Industrial PhD Program (Project No: 119C089).

Declaration of conflicting interests

The authors declared no potential conflicts of interest with respect to the research, authorship, and/or publication of this article.

Funding

The authors disclosed receipt of the following financial support for the research, authorship, and/or publication of this article: This work was studied under funding of the TÜBİTAK-ARDEB 1071 International Cornet projects program (Project No: 123N452) and TÜBİTAK for the scholarship awarded to the first author under the 2244 - Industrial PhD Program (Project No: 119C089).

ORCID iDs

Ayşe Polat

Aiman Tariq

Babür Deliktaş

Murat Yazıcı

Data Availability Statement

Some or all data that support this study’s findings are available from the corresponding author upon reasonable request.*

References

Sun

, et al. Mixed modes crack propagation of orthogonal woven layer in carbon/aramid/epoxy laminates. Compos Appl Sci Manuf 2024; 177: 107929.

Hessman

Welschinger

Hornberger

, et al. On mean field homogenization schemes for short fiber reinforced composites: unified formulation, application and benchmark. Int J Solid Struct 2021; 230-231: 111141.

Müller

Brylka

Dillenberger

, et al. Homogenization of elastic properties of short-fiber reinforced composites based on measured microstructure data. J Compos Mater 2016; 50(3): 297–312.

Voyiadjis

Deliktas

. Modeling of strengthening and softening in inelastic nanocrystalline materials with reference to the triple junction and grain boundaries using strain gradient plasticity. Acta Mech 2010; 213(1): 3–26.

Elmoghazy

Abuelgasim

EMO

Osman

, et al. Effective mechanical properties evaluation of unidirectional and bidirectional composites using virtual domain approach at microscale. Arch Adv Eng Sci 2023; 1(1): 27–37.

Higuchi

Yokozeki

Nishida

, et al. High-fidelity computational micromechanics of composite materials using image-based periodic representative volume element. Compos Struct 2024; 328: 117726.

Voyiadjis

Deliktas

Aifantis

. Multiscale analysis of multiple damage mechanisms coupled with inelastic behavior of composite materials. J Eng Mech 2001; 127(7): 636–645.

Mori

Tanaka

. Average stress in matrix and average elastic energy of materials with misfitting inclusions. Acta Metall 1973; 21(5): 571–574.

Benveniste

. A new approach to the application of Mori-Tanaka’s theory in composite materials. Mech Mater 1987; 6(2): 147–157.

10.

Doghri

Ouaar

. Homogenization of two-phase elasto-plastic composite materials and structures: study of tangent operators, cyclic plasticity and numerical algorithms. Int J Solid Struct 2003; 40(7): 1681–1712.

11.

Fang

Liang

, et al. Investigation on the compressive properties of the three dimensional four-directional braided composites. Compos Struct 2011; 93(2): 392–405.

12.

Kammoun

Doghri

Brassart

, et al. Micromechanical modeling of the progressive failure in short glass–fiber reinforced thermoplastics–First Pseudo-Grain damage model. Compos Appl Sci Manuf 2015; 73: 166–175.

13.

Petrolo

Pagani

Trombini

, et al. Voids effect on micromechanical response of elastoplastic fiber-reinforced polymer composites using 1D higher-order theories. Mech Mater 2023; 184: 104747.

14.

Guedes

Kikuchi

. Preprocessing and postprocessing for materials based on the homogenization method with adaptive finite element methods. Comput Methods Appl Mech Eng 1990; 83(2): 143–198.

15.

Feyel

. A multilevel finite element method (FE2) to describe the response of highly non-linear structures using generalized continua. Comput Methods Appl Mech Eng 2003; 192(28-30): 3233–3244.

16.

Selmi

Doghri

Adam

. Micromechanical simulations of biaxial yield, hardening and plastic flow in short glass fiber reinforced polyamide. Int J Mech Sci 2011; 53(9): 696–706.

17.

Heidari-Rarani

Bashandeh-Khodaei-Naeini

Mirkhalaf

. Micromechanical modeling of the mechanical behavior of unidirectional composites–A comparative study. J Reinforc Plast Compos 2018; 37(16): 1051–1071.

18.

Michel

Suquet

. Nonuniform transformation field analysis. Int J Solid Struct 2003; 40(25): 6937–6955.

19.

Yvonnet

. The reduced model multiscale method (R3M) for the non-linear homogenization of hyperelastic media at finite strains. J Comput Phys 2007; 223(1): 341–368.

20.

Sharafi

Santare

Gerdes

, et al. A multiscale modeling approach of the fused filament fabrication process to predict the mechanical response of 3D printed parts. Addit Manuf 2022; 51: 102597.

21.

Kirchdoerfer

Ortiz

. Data-driven computational mechanics. Comput Methods Appl Mech Eng 2016; 304: 81–101.

22.

Tariq

Abualshar

Deliktas

, et al. ANN-based evaluation system for erosion resistant highway shoulder rocks. Int J Geoeng 2024; 15(1): 17.

23.

Tariq

Uzun

Deliktaş

, et al. Application of machine learning methodology for investigating the vibration behavior of functionally graded porous nanobeams. J Strain Anal Eng Des 2025; 60(2): 131–151.

24.

Polat

Tariq

Okay

, et al. Investigation of the critical buckling load of a column with linearly varying moment of inertia using analytical, numerical, and hybrid machine learning approaches. J Strain Anal Eng Des. 2025; 03093247251337987.

25.

Tariq

Deliktaş

. An inverse parameter identification in finite element problems using machine learning-aided optimization framework. Exp Mech 2025; 65: 325–349.

26.

Mentges

Dashtbozorg

Mirkhalaf

. A micromechanics-based artificial neural networks model for elastic properties of short fiber composites. Compos B Eng 2021; 213: 108736.

27.

Tariq

Polat

Deliktaş

. Boosting machine learning algorithms for predicting the macroscopic material behavior of continuous fiber reinforced composite. J Reinforc Plast Compos. 2024; 07316844241292694.

28.

Chen

Dong

, et al. A micromechanics-based artificial neural networks model for rapid prediction of mechanical response in short fiber reinforced rubber composites. Int J Solid Struct 2024; 305: 113093.

29.

You

. An integrated micromechanical model and BP neural network for predicting elastic modulus of 3-D multi-phase and multi-layer braided composite. Compos Struct 2015; 122: 308–315.

30.

Ghane

Fagerström

Mirkhalaf

. A multiscale deep learning model for elastic properties of woven composites. Int J Solid Struct 2023; 282: 112452.

31.

El Said

. Predicting the non-linear response of composite materials using deep recurrent convolutional neural networks. Int J Solid Struct 2023; 276: 112334.

32.

Yvonnet

Q-C

. Computational homogenization of nonlinear elastic materials using neural networks. Int J Numer Methods Eng 2015; 104(12): 1061–1084.

33.

Xue

Cao

Shao

, et al. Study of effective elastic properties of heterogeneous materials with an artificial neural network model. Mech Mater 2023; 179: 104597.

34.

Lourenço

Tariq

Georgieva

, et al. On the use of physics-based constraints and validation KPI for data-driven elastoplastic constitutive modelling. Comput Methods Appl Mech Eng 2025; 437: 117743.

35.

Kumar Maharana

Soni

Mitra

. A machine learning based prediction of elasto-plastic response of a short fiber reinforced polymer (SFRP) composite. Model Simulat Mater Sci Eng 2023; 31(7): 075001.

36.

Zhao

Chen

Jian

. A high-generalizability machine learning framework for analyzing the homogenized properties of short fiber-reinforced polymer composites. Polymers 2023; 15(19): 3962.

37.

Tariq

Uzun

Akpinar

, et al. Size dependent dynamics of a bi-directional functionally graded nanobeam via machine learning methods. Adv Nano Res 2025; 18(1): 33–52.

38.

Santos

Aguado

Corazza

, et al. Artificial neural network for aspect ratio prediction of lignocellulosic micro/nanofibers. Cellulose 2022; 29(10): 5609–5622.

39.

Mouloodi

Rahmanpanah

Burvill

, et al. Experimental, regression learner, numerical, and artificial neural network analyses on a complex composite structure subjected to compression loading. Mech Adv Mater Struct 2022; 29(17): 2437–2453.

40.

Tariq

Uzun

Deliktaş

, et al. Bending analysis of cantilever microbeams with three porosity distributions using physics-informed neural network and modified couple stress theory. Eng Appl Artif Intell 2025; 159: 111589.

41.

Mahakur

Paul

Bhowmik

, et al. Influence of surface modification on mechanical and tribology performance of jute filler polymer composites and prediction of the performance using artificial neural network. Polym Bull 2023; 80(11): 11953–11974.

42.

Karna

Sahai

. An overview on Taguchi method. Int J Eng Math Sci 2012; 1(1): 1–7.

43.

Gopalsamy

Mondal

Ghosh

. Taguchi method and ANOVA: an approach for process parameters optimization of hard machining while machining hardened steel. J Sci Ind Res 2009; 68: 686–695.

44.

Wetzel

Schneider

Sambale

, et al. RVE simulations of short fiber reinforced polyamide: direct and inverse matrix parameter identification in view of the semi-crystalline polymer structure. Compos Sci Technol 2025; 260: 110961.

45.

Praud

Schneider

Chatzigeorgiou

, et al. Microstructure generation and full-field multi-scale analyses for short fiber reinforced thermoplastics: application to PA66GF composites. Compos Struct 2024; 341: 118175.

46.

Katouzian

Vlase

. Mori–Tanaka formalism-based method used to estimate the viscoelastic parameters of laminated composites. Polymers 2020; 12(11): 2481.

47.

Stránský

Vorel

Zeman

, et al. Mori-Tanaka based estimates of effective thermal conductivity of various engineering materials. Micromachines 2011; 2(2): 129–149.

48.

Yan

Wang

Sun

, et al. Micromechanics-based prediction models and experimental validation on elastic modulus of recycled aggregate concrete. Sustainability 2021; 13(20): 11172.

49.

Schubert

. Novel theoretical self‐consistent mean‐field approach to describe the conductivity of carbon fiber filled Thermoplastics—Part I—Theory. Macromol Theory Simul 2018; 27(4): 1700104.

50.

Rashidinejad

Ahmadi

Hajikazemi

, et al. Modeling of geometric configuration and fiber interactions in short fiber reinforced composites via new modified Eshelby tensors and enhanced mean-field homogenization. Mech Mater 2021; 162: 104059.

51.

Doghri

Brassart

Adam

, et al. A second-moment incremental formulation for the mean-field homogenization of elasto-plastic composites. Int J Plast 2011; 27(3): 352–371.

52.

Polat

Deliktaş

Yazıcı

, et al. Computational modeling of weld-line impacts on mechanical behavior of fiber-reinforced thermoplastics. Eur J Mech Solid 2025; 109: 105485.

53.

Jain

. Micro and mesomechanics of fibre reinforced composites using mean field homogenization formulations: a review. Mater Today Commun 2019; 21: 100552.

54.

Advani

Tucker

III . The use of tensors to describe and predict fiber orientation in short fiber composites. J Rheol 1987; 31(8): 751–784.

55.

Gusev

. Finite element estimates of viscoelastic stiffness of short glass fiber reinforced composites. Compos Struct 2017; 171: 53–62.

56.

Orban

Manea

Aciu

, et al. Virtual manufacturing and mechanical properties of synthetic fiber-reinforced mortars. Procedia Manuf 2018; 22: 262–267.

57.

Isaincu

Dan

Ungureanu

, et al. Numerical investigation on the influence of fiber orientation mapping procedure to the mechanical response of short-fiber reinforced composites using Moldflow, Digimat and Ansys software. Mater Today Proc 2021; 45: 4304–4309.

58.

Friemann

Dashtbozorg

Fagerström

, et al. A micromechanics‐based recurrent neural networks model for path‐dependent cyclic deformation of short fiber composites. Int J Numer Methods Eng 2023; 124(10): 2292–2314.

59.

Voyiadjis

Samadi-Dooki

. Constitutive modeling of large inelastic deformation of amorphous polymers: free volume and shear transformation zone dynamics. J Appl Phys 2016; 119(22): 225104.

60.

Samadi-Dooki

Malekmotiei

Voyiadjis

. Characterizing shear transformation zones in polycarbonate using nanoindentation. Polymer 2016; 82: 238–245.

61.

Shojaei

Voyiadjis

Tan

. Viscoplastic constitutive theory for brittle to ductile damage in polycrystalline materials under dynamic loading. Int J Plast 2013; 48: 125–151.

62.

Sobol

. On the distribution of points in a cube and the approximate evaluation of integrals. USSR Comput Math Math Phys 1967; 7: 86–112.

63.

Renardy

Joslyn

Millar

, et al. To Sobol or not to Sobol? The effects of sampling schemes in systems biology applications. Math Biosci 2021; 337: 108593.

64.

Kucherenko

Albrecht

Saltelli

. Exploring multi-dimensional spaces: a comparison of Latin hypercube and quasi Monte Carlo sampling techniques. arXiv preprint arXiv:1505.02350. 2015.

65.

Jolliffe

Cadima

. Principal component analysis: a review and recent developments. Philos Trans A Math Phys Eng Sci 2016; 374(2065): 20150202.

66.

Tariq

Kadıoğlu

Uzun

, et al. Modeling the viscoelastic behavior of a FG nonlocal beam with deformable boundaries based on hybrid machine learning and semi-analytical approaches. Arch Appl Mech 2025; 95(4): 79.

67.

Chen

Guestrin

. Xgboost: a scalable tree boosting system. In: Proceedings of the 22nd ACM SIGKDD international conference on knowledge discovery and data mining, San Francisco, CA, USA, 13–17 August 2016, pp. 785–794.

68.

LeCun

Bengio

Hinton

. Deep learning. Nature 2015; 521(7553): 436–444.

69.

Nikzad

Heidari-Rarani

Mirkhalaf

. A novel Taguchi-based approach for optimizing neural network architectures: application to elastic short fiber composites. Compos Sci Technol 2025; 259: 110951.

70.

Ross

. Taguchi techniques for quality engineering: loss function, orthogonal experiments, parameter and tolerance design. New York: McGraw-Hill, 1996.

71.

Freddi

Salmon

. Introduction to the Taguchi method. In: Ceschin

Ruggeri

Sassanelli

Terzi

(eds) Design principles and methodologies: from conceptualization to first prototyping with examples and case studies. Berlin: Springer, 2019, pp. 159–180.

72.

Montgomery

. Design and analysis of experiments. 7th ed. Hoboken, NJ: Wiley, 2008.

73.

Chen

Carrera Uribe

Kwon

, et al. A comprehensive review of thermoelectric generation optimization by statistical approach: Taguchi method, analysis of variance (ANOVA), and response surface methodology (RSM). Renew Sustain Energy Rev 2022; 169: 112917.