# Mesoscopic finite element modeling of biaxial non-crimp fabric including representative stitch pattern

Ruochen Zheng, Naïm Naouar, Julien Colmars, Auriane Platzer, Bastian

Schaefer, Fabrice Morestin, Luise Kaerger, Philippe Boisse

► To cite this version:

Ruochen Zheng, Naïm Naouar, Julien Colmars, Auriane Platzer, Bastian Schaefer, et al., Mesoscopic finite element modeling of biaxial non-crimp fabric including representative stitch pattern. *Composite Structures*, 2024, 339, 10.1016/j.compstruct.2024.118126. hal-04707620

HAL Id: hal-04707620

https://hal.science/hal-04707620v1

Submitted on 10 Oct 2024

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
<figure><img src="image_2.png" /></figure>

# Mesoscopic finite element modeling of biaxial non-crimp fabric including representative stitch pattern

Ruochen Zheng<sup>a</sup>, Naim Naouar<sup>a</sup>, Julien Colmars<sup>a,*</sup>, Auriane Platzer<sup>a</sup>, Bastian Schäfer<sup>b</sup>, Fabrice Morestin<sup>a</sup>, Luise Kärger<sup>b</sup>, Philippe Boisse<sup>a</sup>

## ARTICLE INFO

### Keywords:

- Non-crimp fabric

- Mesoscopic modeling

- Forming simulation

- Finite element modeling

## ABSTRACT

A comprehensive understanding of the mechanical properties of carbon fiber reinforcements is necessary to accurately simulate forming processes. The parametrization of macroscopic models requires extensive experimental characterization for different deformation modes. In order to reduce experimental effort, this work proposes a mesoscopic model of bi-axial non-crimp fabric (Bi-NCF) that can be used for virtual material characterization. The fiber yarns are considered as a continuous medium and modeled with solid elements, and their mechanical behavior is described by a hyperelastic constitutive law dedicated to anisotropic fibrous media. The stitches are modeled with 2-node beam elements. This model is developed based on precise geometry measurements obtained from X-ray tomography results. An additional stitch pre-tension step is incorporated in the model to improve contact and simulate the tension induced by the sewing process. Experimental Picture Frame Tests are used to validate the model at mesoscale, while experimental and numerical Bias-extension Tests demonstrate its potential to simulate larger scales and effectively predict local defects.

## 1. Introduction

Continuous fiber composites are commonly used in aeronautical, automotive and marine applications for their high strength and light weight. Among the possible reinforcements, bi-axial non crimp fabrics (Bi-NCF) are an interesting solution. The fibers in these fabrics are arranged in two main directions and are bound together by stitches that can be tailored to meet the specific needs of each application. As the fiber yarns remain straight along their direction, Bi-NCFs provide improved mechanical properties as compared to traditional woven composites [1].

For these reasons, Bi-NCFs are used in various industries [2,3]. For instance, in aircraft components, such as fuselage panels and wing structures, to reduce weight and improve fuel efficiency [1,4,5]. Automotive applications include weight reduction of structural parts [2,6,7]. NCFs are increasingly used in renewable energy facilities, specifically in the manufacturing of wind turbine blades, taking advantage of their high strength-to-weight ratio and fatigue resistance [1].

Nevertheless, the shaping of these fabrics is a major technical challenge. The deformation mechanisms involved are specific and there are few reliable approaches for the simulation of forming [8–12]. In particular, defects appear during dry-forming processes, leading to

"resin-rich" areas during injection, which weakens the final composite parts. Defects include gaps between the fiber yarns or wrinkling [13]. Several studies have thus focused on characterizing the mechanical behavior of these materials [14–17].

Numerical models are useful tools to predict the mechanical response of the fabric during forming processes. Finite element models at the mesoscopic level are the most able to capture defects during forming, especially the gaps between yarns. However this requires high computational costs. At this scale, the fiber yarns are often considered as homogeneous and continuous solids, made of the compact assembly of carbon fibers, while stitches are modeled by slender structural elements. Some authors proposed modeling for the interactions (frictional contact, slippage, etc.) between the fiber yarns and the stitches [18–20]. In these models, the unidirectional fiber layers are modeled by 3D elements and the stitches are modeled by bar elements. In general, mesoscopic models mainly focus on unit cell simulations [12,21–24], giving insightful analysis of the local mechanisms inducing defects during forming, or material variability [25]. However, these reduced models can hardly scale to the forming simulation of complete parts, mostly limited by computational costs. For this reason, there is very little simulation of NCF at the mesoscopic scale. Some examples can be

<footnote><sup>a</sup> Univ Lyon, INSA Lyon, CNRS, LaMCoS, UMR5259, Villeurbanne, F-69621, France</footnote>

<footnote><sup>b</sup> Karlsruhe Institute of Technology (KIT), Institute of Vehicle System Technology (FAST), Lightweight Design (LB), Karlsruhe 76131, Germany</footnote>

<footnote>* Corresponding author.
E-mail address: julien.colmars@insa-lyon.fr (J. Colmars).</footnote>
<figure><img src="image_3.png" /><figcaption>Fig. 1. Bi-NCF (a) Warp yarn view with a zig-zag stitch pattern (b) Weft yarn view with a looped stitch pattern.</figcaption></figure>

cited: in [19], the geometry of yarns was overly simplified, which did not allow to access mesoscopic information when performing forming simulations. Recently [26] proposed a micro-meso model of NCF based on multi-filament method. This method has the ability to capture meso-scale deformations. One limitation of this method is that the behavior of yarns is derived from the kinematic interaction of hundreds of beam elements (fibers); thus the mechanical parameters of this model are calibrated on complex simulations, and have no direct link with the yarn properties that can be measured in mesoscopic characterization tests. The objective of our study is to propose a numerical approach based on geometric and mechanical properties measured at mesoscopic scale.

In this study, the finite element modeling of NCF is constructed from a realistic geometrical description of the yarns and stitches. In recent years, two main approaches have been adopted to develop mesoscopic models of carbon fiber fabrics: an ideal geometry is generated from software tools such as TexGen [27] and WiseTex [28,29]; or a more realistic finite element mesh can be directly obtained from the results of X-ray tomography [30–32]. The former softwares are valuable tools for efficiently and automatically creating geometries, based on a generalized description of the internal structure of the fabric (width of the yarns, woven and/or stitching pattern, etc.). The latter approach based on tomography offers a higher level of accuracy as it gives access to smaller details and spatial variation of the geometrical parameters [33]. However, tomography images of Bi-NCF are difficult to directly turn into a geometrical or finite element model as the stitches are indistinguishable among warp and weft fiber yarns. To alleviate for this lack of information, stitches structure is often modeled a posteriori around of the fiber yarn, by making generalized assumptions about the internal unit cell structure [34–36]. In the present study, X-ray tomography was carried out to measure geometrical parameters of both warp and weft yarns while the surrounding stitches were described explicitly from the knitting pattern. The geometrical model is constructed on a parameterized geometry, and then fitted on tomography measurements. This method keeps the advantage of the flexibility of geometric modeling while getting as close as possible to tomography measurements.

This article proposes an integrative approach to mesoscopic modeling. The objective is to provide an FEM model powered by a realistic geometry measured by tomography measurements and mechanical properties characterized at mesoscopic level. A complete finite element model of a Bi-NCF at mesoscopic scale is developed within the commercial finite element software Abaqus-Explicit. One of the main goal of this model is to capture local details such as gaps between yarns, while performing simulations of relatively large samples, such as the

ones classically used in Bias Extension Tests (BET). The rest of the paper details the implementation of the numerical model and proposes a detailed validation of the model through experimental-numerical comparisons.

The paper is organized as follows. In Section 3, we present the geometrical modeling: a realistic geometry is generated from X-ray tomography measurements of the fiber yarns and a geometrical model of the knitting pattern for the stitches. The resulting finite element model is given in Section 4, where we describe the mechanical modeling approach based on specific constitutive modeling of the fiber yarns and the stitches respectively. The constitutive models are calibrated with original experimental tests also detailed in Section 4. Next, in Section 5, the mesoscopic model is validated on a unit cell simulation compared to a dedicated Picture Frame Test (PFT). Finally, the ability of the proposed mesoscopic model to scale to macroscopic simulation is demonstrated in Section 6 on the numerical-experimental comparison of a BET.

## 2. Materials

The selected material is a 0°/90° Bi-NCF with a tricot stitch pattern, produced by American manufacturer Zoltek from their PX35 continuous carbon fiber tow. This NCF consists of two layers of non-woven fibrous plies laid in perpendicular directions. The polyester stitches are used to pierce the two plies at specific positions between the fiber yarns and to bond the plies together. The warp direction of fiber yarns, with a so-called *zig-zag* stitch pattern, is shown in Fig. 1(a), and the weft direction, with a so-called *looped* stitch pattern, is shown in Fig. 1(b). Throughout the article, both faces of the material are then referred to as *zig-zag* and *looped* face. The specifications of the material, as given by the supplier, are listed in Table 1.

## 3. Geometrical modeling

### 3.1. Tomography

X-ray micro-Computed Tomography (XRµCT) was conducted on a new microCT platform developed by the company RX Solutions in collaboration with LaMCoS and Mateis laboratories. The Dual Tomograph for High Energy (DTHE) is designed around one rotational axis and two 300 kV Xray lines. The data acquisition system recorded a total of 1120 projections, uniformly distributed across a 360° rotation around the sample's vertical axis. To obtain the 3D reconstruction, the recorded projections were processed using Xact software. A 18 mm × 18 mm ×
<figure><img src="image_4.png" /><figcaption>Fig. 2. Tomography 3D view of the sample, with (a) warp yarn side, and (b) weft yarn side.</figcaption></figure>

<table><caption>**Table 1**

Characteristics of the material Zoltek™ PX35 multi-directional fabrics.</caption><thead><tr><th colspan="2">General specifications</th></tr></thead><tbody><tr><td>Fiber</td><td>Carbon, 50 k filaments</td></tr><tr><td>Fiber diameter</td><td>7.2 µm</td></tr><tr><td>Fiber orientation</td><td>0°/90°</td></tr><tr><td>Stitch</td><td>PES</td></tr><tr><td colspan="2">Areal density g m<sup>-2</sup></td></tr><tr><td>Fiber yarn</td><td>300 g m<sup>-2</sup></td></tr><tr><td>Stitch</td><td>8 g m<sup>-2</sup></td></tr><tr><td>Fabric</td><td>608 g m<sup>-2</sup></td></tr></tbody></table>

2.25 mm volume sample of Bi-NCF was analyzed with a voxel resolution of 15 µm. The resulting 3D view of the sample is shown in Fig. 2.

Several local geometrical parameters were then measured, as illustrated in Fig. 3:

**gap** The gap is defined as an opening between fiber yarns when the knitting needle pierces the two fibrous plies during the knitting process. It is easily detected on the tomography images ((X, Y) plane) as the width of the black area between the fiber yarns, as shown in Fig. 3(a).

**width of a yarn** As above mentioned, a yarn is defined as a fiber bundle, bound together by the stitch. The width of a yarn is thus defined as the distance between the edges of the fiber bundle, as shown in Fig. 3(a).

**thickness of a yarn** The thickness of a yarn is measured on the cross-sectional view of the sample ((X, Z) plane), as shown in Fig. 3(b).

The measurements were performed using Fiji software [37] on different views and areas of the sample. Ten positions were selected for each parameter, and the average value was calculated and reported in Table 2. Note that during the forming process, the thickness and width of both warp and weft yarns will change due to the stretching of the stitch, leading to a larger gap, causing a future “resin-rich” area. The careful measurement of these geometrical parameters in the initial configuration of the fabric thus improves the prediction of local defects occurring during forming, as compared to an idealized geometry generated from the manufacturer’s data sheet only. However, the proposed approach is based on the hypothesis that a unit representative cell can be found on the NCF geometry.

<table><caption>**Table 2**

Geometrical parameters (average values) measured on tomography images of a Bi-NCF sample.</caption><thead><tr><th></th><th>Width<br>W<sub>tomo</sub></th><th>Thickness<br>H<sub>tomo</sub></th><th>Gap<br>G<sub>tomo</sub></th></tr></thead><tbody><tr><td>Warp yarn</td><td>4.59 mm (±0.13)</td><td>0.38 mm (±0.06)</td><td>0.39 mm (±0.10)</td></tr><tr><td>Weft yarn</td><td>2.08 mm (±0.08)</td><td>0.38 mm (±0.04)</td><td>0.40 mm (±0.04)</td></tr></tbody></table>

## 3.2. Geometrical model of the yarn

Unlike woven reinforcements, the weft and warp yarns of NCF remain straight along the fiber direction, neglecting the minor local compression of the yarns caused by the stitch during the knitting process. Thus, only their cross-section needs to be defined. Several model of the cross-section have included circular, flat, elliptical, biconvex or rectangular shapes [34,38,39].

In this study, the cross-section geometry was first approximated by an elliptical shape in the (X, Z) plane, given by the following equation:

$$ \left(\frac{X}{W_{\text{ini}}/2}\right)^2 + \left(\frac{Z}{H_{\text{ini}}/2}\right)^2 = 1, \qquad (1) $$

where $H_{ini}$ is the initial thickness and $W_{ini}$ is the initial width of the fiber yarn. The ellipses are slightly truncated at the ends which enables to mesh with 8-nodes hexahedral elements; this provides a first idealized geometry of the yarn sections (Fig. 4(a)).

To enhance the fidelity of the final model and ensure that it captures the tightening of fibers during the knitting process, a pre-tension step is introduced in Section 4.3.2. This step compresses the sides of the ellipse, transforming the cross-section into a semi-ellipse, as shown in Fig. 4(b). The initial input dimensions of the yarn, $H_{ini}$ and $W_{ini}$, are adjusted by trial-and-error so that the pre-tensioned dimensions, H and W in Fig. 4(b), match the respective tomography measurements, $H_{tomo}$ and $W_{tomo}$, given in Table 2. Additionally, the gap between each fiber yarn, denoted by $G_{ini}$ in the initial configuration, is assumed to be uniformly distributed across the fabric and remains constant along the direction of the fibers. However after the pre-tension step the gap, denoted G, can exhibit variability along the yarn due to non-homogeneous tightening of stitches, as observed in tomography measurements and discussed in previous studies [25].

## 3.3. Knitting pattern

The explicit description of the stitches is essential due to the existence of various stitch patterns resulting from different NCF types and
<figure><img src="image_5.png" /><figcaption>**Fig. 3.** Measurements of geometrical parameters on tomography results: (a) orthogonal view of a warp yarn ply in the (X, Y) plane, (b) cross-sectional view of the fabric in the (X, Z) plane.</figcaption></figure>

<figure><img src="image_6.png" /><figcaption>**Fig. 4.** Cross-sections of the Bi-NCF geometrical model (a) before and (b) after the stitch pre-tension step.</figcaption></figure>

manufacturing processes. These stitch patterns significantly influence the mechanical properties of the entire NCF [40]. Several studies have provided mathematical descriptions of different stitch patterns [41–43], and the use of a trace fiber technique for stitching path analysis was introduced in [34].

In our case, the initial configuration of tricot stitch path is computed to ensure contact with the elliptical cross section of fiber yarns [2] (accounting for stitch radius). After the pre-tension step, one can observe that the stitch knitting points are not perfectly aligned, as evidenced by the stitch piercing positions in Fig. 3(a).
<figure><img src="image_7.png" /><figcaption>Fig. 5. Tricot stitch model of Bi-NCF. (a) One loop of stitch (red) on a 10 mm × 10 mm model with a zoomed-in view of its loop connection position. (b) Sectional view ((Y, Z plane) of the stitch model in the initial configuration (before pre-tension step).</figcaption></figure>

<figure><img src="image_8.png" /><figcaption>Fig. 6. Geometrical model of a 30 mm × 30 mm Bi-NCF periodic sample, generated from tomography measurements and a tricot knitting pattern. (a) Zig-zag face with warp yarns. (b) Looped face with weft yarns.</figcaption></figure>

Fig. 5 illustrates the periodic stitch pattern surrounding one warp yarn and piercing through four weft yarns (red path).

### 3.4. Geometrical model of bi-NCF

The geometrical description of both yarns and stitches were implemented in a Matlab script to automatically generate a numerical sample of Bi-NCF in Abaqus. The resulting complete geometrical model of the Bi-NCF is shown in Fig. 6.

## 4. Mechanical modeling

In this section, the choice of finite elements and constitutive law are presented for the constituents (yarn and stitches). Note that the entire mechanical modeling process is conducted within Abaqus-Explicit.

### 4.1. Yarn modeling

#### 4.1.1. Finite elements for yarns

Based on prior studies [44–47], the mesoscopic scale considers each fiber yarn as a homogeneous and continuous medium. Therefore, solid hexahedral elements were chosen to effectively approximate the semi-elliptical cross-section of each fiber yarn and facilitate the extrusion

of the finite element mesh along the fiber direction. Abaqus C3D8 elements were used to ensure a sufficient number of integration points in the cross-section. A minimum of two elements in thickness direction of the cross section was empirically determined as a suitable compromise between accuracy and computational cost.

#### 4.1.2. Constitutive modeling for yarns

Fiber yarns consist in a compact assembly of numerous fibers (up to 50 k in our case) that are linked with a binder; here the yarns are treated as a continuous medium at the mesoscale [48]. Consequently, a three-dimensional constitutive equation is derived to account for the specific behavior of the fibers. The homogenized material demonstrates a preferred direction along the fiber orientation, while the fibers are distributed isotropically within the cross-section. Thus, the constitutive model exhibits transverse isotropy. Considering the occurrence of large deformations during the forming of Bi-NCF, a hyperelastic law is employed. Following the rationale presented by Charmetant et al. in [47], the mechanical response of a yarn can be decomposed into multiple deformation modes, as illustrated in Fig. 7: fiber elongation, cross-sectional compaction, transverse shear (cross-sectional distortion), and longitudinal shear along the fiber direction.

In this work, the constitutive equation proposed by Charmetant et al. [47] for the fiber yarns was adopted. The hyperelastic strain energy potential relies on four physically based invariants ($I_{\text{elong}}$, $I_{\text{comp}}$,
<figure><img src="image_9.png" /><figcaption>Fig. 7. Four deformation modes of fiber yarn (a) elongation (b) cross-sectional compaction (c) transverse shear (d) longitudinal shear, from [47].</figcaption></figure>

$I_{\text{dis}}$, $I_{\text{sh}}$) that represent the four deformations modes. These invariants are defined as functions of the classical invariants $I_1$ to $I_5$ of the right Cauchy–Green strain tensor $\underline{C}$, expressed as [49]:

$$I_1 = \text{tr} \underline{C}, \qquad (2a)$$

$$I_2 = \frac{1}{2} (\text{tr} \underline{C}^2 - \text{tr} \underline{C}^2), \qquad (2b)$$

$$I_3 = \det \underline{C}, \qquad (2c)$$

$$I_4 = \underline{C} : \underline{M}, \qquad (2d)$$

$$I_5 = \underline{C}^2 : \underline{M}. \qquad (2e)$$

Here, $\underline{M}$ is a structural tensor defined by a unit vector representing the preferred direction $\underline{M}$ (fiber direction):

$$\underline{M} = \underline{M} \otimes \underline{M}. \qquad (3)$$

Based on the classical invariants Eq. (2), Charmetant et al. defined a set of four physical invariants [47]:

$$I_{\text{elong}} = \frac{1}{2} \ln I_4 \qquad (4a)$$

$$I_{\text{comp}} = \frac{1}{4} \ln \left( \frac{I_3}{I_4} \right) \qquad (4b)$$

$$I_{\text{sh}} = \sqrt{\frac{I_5}{I_4^2} - 1} \qquad (4c)$$

$$I_{\text{dist}} = \frac{1}{2} \ln \left( \frac{I_1 I_4 - I_5}{2 \sqrt{I_3 I_4}} + \sqrt{\left( \frac{I_1 I_4 - I_5}{2 \sqrt{I_3 I_4}} \right)^2 - 1} \right) \qquad (4d)$$

The total strain energy density $w$ is written as a function of these physical invariants, with separated contributions for each mode:

$$\begin{aligned} w &= w(I_{\text{elong}}, I_{\text{comp}}, I_{\text{dist}}, I_{\text{sh}}) \\ &= w_{\text{elong}}(I_{\text{elong}}) + w_{\text{comp}}(I_{\text{comp}}) + w_{\text{dist}}(I_{\text{dist}}) + w_{\text{sh}}(I_{\text{sh}}) \end{aligned} \qquad (5)$$

The strain energy density functions associated with elongation ($w_{\text{elong}}$), transverse ($w_{\text{dist}}$) and longitudinal ($w_{\text{sh}}$) shear are defined as standard quadratic elastic potentials [47]:

$$w_{\text{elong}} = \frac{1}{2} K_{\text{elong}} I_{\text{elong}}^2, \qquad (6a)$$

$$w_{\text{dist}} = \frac{1}{2} K_{\text{dist}} I_{\text{dist}}^2, \qquad (6b)$$

$$w_{\text{sh}} = \frac{1}{2} K_{\text{sh}} I_{\text{sh}}^2. \qquad (6c)$$

For the cross-sectional compaction strain energy density $w_{\text{comp}}$, Charmetant et al. proposed a power-based function to represent the microscopic nonlinear phenomena that occur during transverse compaction of the yarn (e.g., rearrangement and bending of the fibers) [47]:

$$w_{\text{comp}} = \begin{cases} K_{\text{comp}} |I_{\text{comp}}|^p & \text{if } I_{\text{comp}} \le 0, \\ 0 & \text{otherwise.} \end{cases} \qquad (6d)$$

The second Piola–Kirchhoff stress tensor is then computed using the classical formula:

$$\underline{S} = 2 \frac{\partial w}{\partial \underline{C}} = \sum_{\text{mode}} 2 \frac{\partial w}{\partial I_{\text{mode}}} \frac{\partial I_{\text{mode}}}{\partial \underline{C}} \qquad (7)$$

<caption>Table 3</caption>

Designation, value, physical interpretation, and calibration method of the material parameters in the Charmetant constitutive model used in this study (see Eq. (6)).

<table><thead><tr><th>Symbol</th><th>Value</th><th>Physical interpretation</th><th>Calibration method</th></tr></thead><tbody><tr><td>K<sub>elong</sub></td><td>108 500 MPa</td><td>Tensile stiffness along the fibers</td><td></td></tr><tr><td>K<sub>dist</sub></td><td>0.6 MPa</td><td>Transversal shear stiffness</td><td>values from [50], Fabric</td></tr><tr><td>K<sub>comp</sub></td><td>1.25 MPa</td><td>Cross-sectional compaction stiffness</td><td>B</td></tr><tr><td>p</td><td>1.84</td><td>Cross-sectional compaction exponent</td><td></td></tr><tr><td>K<sub>sh</sub></td><td>2.2 MPa</td><td>Longitudinal shear stiffness</td><td>experimental test (see Figs. 8 and 9)</td></tr></tbody></table>

In this study, the constitutive model, including Eqs. (2), (4), (6) and (7), was implemented in a VUMAT within Abaqus-Explicit.

### 4.1.3. Characterization

Referring to the definition of the strain energy densities in Eq. (6), the constitutive model employed requires the calibration of five material parameters, as presented in Table 3. Most of these parameters were previously estimated for carbon fiber yarn in [50], where two types of yarns were tested and the parameters from the type with lower crimp were chosen as an estimate. The tensile stiffness, $K_{\text{elong}}$, and the compaction parameters, $K_{\text{comp}}$ and $p$ were calibrated using tension and compaction tests on a single yarn. The transverse shear stiffness, $K_{\text{dist}}$, was estimated from fabric compression tests and an inverse parameter fit of the thickness–compression force curve. Since these parameters were determined for carbon fiber yarns, they are applicable to our Bi-NCF case. However, in [50], the longitudinal shear stiffness, $K_{\text{sh}}$, was calibrated based on tensile tests conducted on the woven fabric, taking advantage of the presence of crimp in the yarns.

Here an novel experimental protocol is proposed for the determination of longitudinal shear stiffness. The test, performed on a single carbon fiber yarn, is depicted in Fig. 8. This test can be seen as a variant of the shear test presented in [51]. Our way of conducting this test was based on the inextensibility of carbon fibers rather than on the use of specific jaws. A single yarn was extracted from the Bi-NCF fabric and then cut at both ends, leaving only a minimal number of carbon fiber filaments, giving a Z-shape configuration to the specimen (Fig. 8(a)). The boundary conditions were imposed on the external fibers which transmitted them to the sample by cohesion; this further avoided the problems of boundary conditions in the jaws while retaining the hypothesis of a homogeneous shear in the specimen. The test was conducted using a Lloyd electromechanical traction machine with a 1 kN capacity and a 10 N sensor. To enhance clamping and prevent slippage, adhesive tape was applied to secure the carbon fiber filaments at each end, thereby providing a larger surface area for the clamping tool.

The resulting middle shear zone had a length of 100 mm, with the carbon fiber orientation strictly parallel to the loading direction. A slow stretching speed was employed, with the top tape experiencing a speed of 0.3 mm s⁻¹ while the bottom tape remained fixed.

The test results are provided in Fig. 9 where we show the force measured by the load cell, $F$, plotted against the cross-head displacement,
<figure><img src="image_10.png" /><figcaption>Fig. 8. Longitudinal shear test on a single yarn: (a) Schematic illustration (b) Initial configuration (c) Sample after decohesion of the sample.</figcaption></figure>

<figure><img src="image_11.png" /><figcaption>Fig. 9. Longitudinal shear test on a single yarn: force vs. cross-head displacement.</figcaption></figure>

D. Initially, the shear force on the single fiber yarn exhibits a linear increase with displacement. However, as the testing machine continues to stretch the yarn, the shear force gradually decreases until it reaches zero, indicating progressive decohesion. To estimate the longitudinal shear stiffness, we focused solely on the stretching portion of the curve. By knowing all the other parameters of the constitutive model listed in Table 3, an inverse approach can be employed to identify the longitudinal shear stiffness through numerical simulations of same single yarn test. The resulting value of $K_{sh} = 2.2$ MPa is then added to Table 3.

## 4.2. Stitch modeling

### 4.2.1. Constitutive modeling and finite elements

Stitches are slender elements represented by beams. In addition, as the bending stiffness of the filament is much lower than its tensile stiffness, a stress resultant approach is adopted, utilizing the general

<table><caption>Table 4
Designation, value, and calibration method of the material parameters for the stitch.</caption><thead><tr><td>Symbol</td><td>Value</td><td>Designation</td><td>Calibration method</td></tr></thead><tbody><tr><td>E</td><td>350 MPa</td><td>Young's modulus</td><td>Tensile test (see Fig. 10)</td></tr><tr><td>G<sub>s</sub></td><td>0.000 56 N mm<sup>2</sup></td><td>Bending stiffness</td><td>Bending test (see Fig. 11, see Eq. (9))</td></tr><tr><td>I</td><td>0.000 001 6 mm<sup>4</sup></td><td>Moment of inertia for bending</td><td>Tensile and bending test (see Eq. (10))</td></tr><tr><td>r</td><td>0.1 mm</td><td>Cross-section radius</td><td>Tomography measurements</td></tr></tbody></table>

section beam elements provided by Abaqus/Explicit. Both tensile and bending responses are assumed to be linear-elastic. The cross-section is circular, with a radius of $r = 0.1$ mm, as measured by tomography; this parameter is mainly used by contact algorithm. All contacts in the model, including contact between stitches and yarns and contact between yarns, are modeled by Coulomb friction law with a friction coefficient of 0.2; the General Contact algorithm of Abaqus/Explicit is used in all simulations.

### 4.2.2. Characterization

All material parameters for the stitch were calibrated through original experimental tests, as also referenced in Table 4.

The test depicted in Fig. 10 was again conducted using the Lloyd electromechanical traction machine with a 1 kN capacity and a 100 N sensor. A piece of the stitch was extracted from the Bi-NCF, and both ends were wrapped around two small metal pieces multiple times to prevent slipping during stretching. The metal pieces were then clamped to the machine. The resulting stretching length was $L = 90$ mm. The stitch was then stretched until failure at the speed of 90 mm/min. The resulting force, $F$, was measured and plotted against the stretched length (cross-head displacement), $dL$, in Fig. 10(c). As can be seen, the behavior is reasonably represented by the following linear constitutive equation, yielding the Young's modulus value listed in Table 4:

$$F = E A \frac{dL}{L}, \qquad (8)$$

with $A = \pi r^2$ the cross-sectional area.
<figure><img src="image_12.png" /><figcaption>Fig. 10. Tensile test on a piece of stitch (a) Schematic (b) Experimental test initial configuration (c) Force vs. cross-head displacement.</figcaption></figure>

<figure><img src="image_13.png" /><figcaption>Fig. 11. Single stitch cantilever test.</figcaption></figure>

Secondly, the bending stiffness was determined using a cantilever bending test, a well-established method utilized in previous studies [17, 52–54]. However, as the stitch exhibits very low bending stiffness, careful sample preparation is essential. Separating the stitch causes the filament to naturally curl, making it unsuitable for the cantilever test. Although this curling effect does not impact the tensile test, it significantly affects the bending tests. To address this, the separated stitch was straightened over a period of more than three months under small dead weights. Subsequently, the straightest portions were selected as samples for the bending test. As shown in Fig. 11, the test configuration involved clamping one end of the stitch (point B) while the other end bent toward point A under its own weight, forming an angle of inclination, $\theta$, with the horizontal plane. After the test, the stitch was cut, and the overhanging length $l$ between point A and B was measured. The bending stiffness of the stitch, $G_s$, was then determined using the formula provided in Eq. (9) [53], and the resulting value is included in Table 4:

$$G_s = \frac{l^3 \cos\left(\frac{\theta}{2}\right)}{8 \tan \theta} \times \omega, \qquad (9)$$

where $\omega = 4.33 \times 10^{-7} \text{ N mm}^{-1}$ is the weight per unit length of the stitch, calculated from the density of polyester (1.38 g cm$^{-3}$) and the sectional area of the stitch (0.0314 mm$^2$).

Finally, the moment of inertia for bending, $I$, is a user-defined value based on the tensile and bending tests results. It is reported in Table 4 and can be calculated as

$$I = \frac{G_s}{E}. \qquad (10)$$

### 4.3. Finite element assembly of bi-NCF

In previous sections, we have discussed the finite element models for the yarn and the stitch, where each constituent is represented by a specific type of finite element to capture its unique behavior within the carbon fiber fabric. Now, leveraging the geometrical model outlined in Section 3, we proceed to assemble these into a unit cell. Next section presents the methodology for constructing the unit cell and introduces an original approach for ensuring adequate stitch pre-tension to secure the fiber yarns before conducting any forming simulations.

#### 4.3.1. Unit cell

The choice of the unit cell is not unique: for instance, the different types of representative unit cells of woven fabrics are presented in [46]. In the case of NCF, the choice of the elementary cell must possibly adapt to the stitch pattern. The unit cell shown in Fig. 12 was obtained by cutting the periodic structure at half yarn, which is a big advantage for the application of boundary and periodic conditions, as detailed next. It contains two warp yarns and four weft yarns and one entire stitch loop in the center. The model length is 10 mm × 10 mm. The yarns are meshed with 4488 solid elements. The stitches are meshed with 481 beam elements.

In order to apply periodic boundary conditions, the displacement field $\phi(\underline{X})$ of the structure can be decomposed into two components (see Eq. (11)): the macroscopic displacement field $\phi_m(\underline{X})$, which is known and serves as the input for the simulation, and the periodic local displacement $\underline{w}(\underline{X})$, which is unknown and is determined by the mechanical properties [46,55]. In order to ensure periodicity, the local displacements of two opposite sides must be equal. For instance, $\underline{X}^+$ and $\underline{X}^-$ are two opposite set of nodes showed in Fig. 12: periodicity condition is thus given by Eq. (12).

$$\underline{\phi}(\underline{X}) = \underline{\phi}_m(\underline{X}) + \underline{w}(\underline{X}) \qquad (11)$$

$$\underline{w}(\underline{X}^+) = \underline{w}(\underline{X}^-) \qquad (12)$$

By combining Eq. (11) and Eq. (12), the boundary conditions of Eq. (13) are applied in Abaqus for each pair of opposite nodes.

$$\underline{\phi}(\underline{X}^+) - \underline{\phi}(\underline{X}^-) = \underline{\phi}_m(\underline{X}^+) - \underline{\phi}_m(\underline{X}^-) \qquad (13)$$

#### 4.3.2. Stitch-yarn interaction: pre-tension step

During the manufacturing of Bi-NCF, the fiber yarns are tightly bound together by tension stitches through a knitting process. This
<figure><img src="image_14.png" /><figcaption>Fig. 12. Finite element representative unit cell of Bi-NCF (a) upper side, warp yarn with zig-zag stitch (b) lower side, weft yarn with looped stitch.</figcaption></figure>

tension imparts a higher frictional resistance to fiber sliding [56] and leads to gaps between adjacent yarns [2]. While lower knitting tension reduces distortions in fiber yarns and minimizes gaps, it can also loosen the loops of stitches, affecting the interaction between the yarns and stitches. Considering the significant slippage between yarns and stitches in Bi-NCF compared to woven fabric [57], this interaction plays a crucial role in the final results of mesoscopic modeling. Later in this section, we introduce a method to achieve the required pre-tension in the stitch, ensuring that it securely holds the fiber yarns together before any forming simulations take place.

The unit cell described above is assembled based solely on the geometrical considerations given in Section 3. Therefore, the complex interaction between the yarns and the stitches still needs to be incorporated into the mechanical model. Specifically, a contact model is added to simulate the interaction throughout the forming process. To address this challenge in mesoscopic simulation of NCF, various approaches have been proposed and developed in the literature. Among them, an additional artificial bar element was used in [19] to connect the nodes of fiber yarns with their adjacent nodes on the stitch, acting as a substitute for elastic contact. In [57], at macroscopic level, an anisotropic Coulomb friction law was introduced, while a viscous component was added in [58] to enhance numerical stability. The approach developed in this study draws inspiration from previous works where thermal shrinkage is used to mimic the knitting-process tension in yarns [59] or stitches [26]. In our work, prior to any forming simulation on the mechanical model, a pre-tension step was conducted to artificially shrink the stitch under thermo-elastic loading.

To achieve a homogeneous tension in the stitch, an initial homogeneous thermal strain is applied using linear thermo-elastic framework, with

$$ \sigma = E(\epsilon - \epsilon^{th}), \qquad (14) $$

where $(\epsilon, \sigma)$ are the strain and stress in the stitch beam element, respectively, $E$ is the Young's modulus of the stitch (calibrated in Section 4.2.2), and $\epsilon^{th}$ denotes the thermal strain.

At the end of the pre-tension step, the geometry of the yarn is modified asymmetrically: the width of the warp yarn has decreased, whereas the width of the weft yarn remains almost constant. This phenomenon can be attributed to the constraining effect of the stitch loops on the warp direction (see Figs. 4(b) and 5(a)). To maintain a balanced structure and ensure that the stitch neither becomes “too tight” nor “too loose”, various thermal strain parameters of the stitch were tested by trial-and-error, and a reasonable value was found to be $\epsilon^{th} = -0.15$. The initial geometrical parameters $H_{ini}$, $W_{ini}$ and $G_{ini}$ were simultaneously adjusted so that the pre-tensioned geometrical model matches the dimensions of the tomography measurements, as explained in the flow chart Fig. 13. In this figure, the whole process for FE mesoscopic modeling is summarized. The dotted-line arrows connects steps of the adjustment loop, which modifies the thermal strain $\epsilon^{th}$ and initial geometrical parameters in a trial-and-error process. All

ingredients of the final pre-tensioned finite element model are also represented: geometrical considerations in blue, mechanical and finite element parameters in red, pre-tension step in yellow.

## 5. Mesoscopic validation

In order to validate the unit-cell model, an experimental Picture Frame Test (PFT) was performed. Simulations and experiments are then compared in terms of in-plane shearing loading curves.

### 5.1. Experimental PFT

The Picture Frame Test (PFT) is commonly used for characterizing in-plane shear behavior [60]. It consists in fixing a cross-shaped sample inside a hinged frame made up of four rigid bars of the same length. The sample is positioned such that the direction of the fibers imperatively follows the edges of the frame. The frame is sheared in a universal testing machine. In this study, the NCF was embedded in a frame of 120 mm side. The set up is shown in Fig. 14.

It has been shown that $\pm 45^\circ$ Bi-NCF exhibits differences in force under different loading directions during shear deformation, resulting in positive and negative shear [61]. However, the $0^\circ/90^\circ$ Bi-NCF used here possesses a symmetric stitching structure, ensuring that one stitch segment undergoes tension while the other remains loose. As a result, it exhibits consistent force response regardless of the loading orientation ($+45^\circ$ or $-45^\circ$). Consequently, the PFT results do not depend on the initial orientation of the sample.

Yarns are fixed in the frame. A pure shear deformation is considered in the center square zone. It has been observed that the deviation between the theoretical shear angle and the local shear angle is small in this kind of tests [62,63]; this will be confirmed later. In this case, the shear angle $\gamma$ is given by Eq. (15) as a function of the machine displacement $d$ and the frame side length $L_{frame}$.

$$ \gamma = 45^\circ - 2 \arccos \left( \frac{\sqrt{2L_{frame} + d}}{2L_{frame}} \right) \qquad (15) $$

The net force is obtained by subtracting an empty frame force $F_f$ from the tested machine force $F_m$. Then the shear force $F_s$ can be calculated by this net force and the shear angle using Eq. (16).

$$ F_s = \frac{F_m - F_{frame}}{2 \cos \left( 45 - \frac{\gamma}{2} \right)} \qquad (16) $$

A normalized shear force $F_{nor}$, given in N mm$^{-1}$, is calculated with Eq. (17) based on an energy approach [64], in order to compare with other shearing tests.

$$ F_{nor} = F_s \cdot \frac{L_{frame}}{L_{fabric}^2} \qquad (17) $$
<figure><img src="image_15.png" /><figcaption>Fig. 13. Flow chart of FE mesoscopic modeling, including the adjustment loop of thermal strain and initial geometrical parameters, represented by dotted-line connectors.</figcaption></figure>

## 5.2. Experimental-numeric comparison on unit cell

Three PFTs have been conducted on the same material. The relation between shear angle $\gamma$ (calculated with Eq. (15)) and normalized force $F_{nor}$ (calculated with Eq. (17)) was fitted by a polynomial law. The curves $F_{nor}(\gamma)$ are shown in Fig. 15.

In order to reproduce the PFT, the FEM unit cell was submitted to a pure shear test. Due the presence of stitches and periodic boundary conditions, this can be done with few additional boundary conditions, showed in Fig. 14(b): DOFs are put at zero on the lower corner, and displacement is imposed on the top corner. As a reminder, the FEM unit cell is made of 4488 solid elements for yarns, and 481 beam elements for stitches. The calculation lasted 16 h on 2 CPU. The reported run-time corresponded to simulations without any mass-scaling. Time parameters were chosen to avoid any dynamical effect. The set of explicit parameters could be improved to achieve better numerical performances. The mesh sensitivity was conducted on the unit-cell, concluding that a minimum number of two elements in the yarn thickness was required in all simulations. Simulation results are shown in Fig. 16. Yarns are well embedded by the stitches, even the

central ones, on which no boundary conditions were put (except periodic conditions). Significant slippages occur between yarns and stitches (see Fig. 16(a)). The normalized shear force was calculated from the finite element results, and added to Fig. 15. Experimental PFT and unit cell modeling are in very good agreement in terms of shear force. The deviation is higher at the beginning of the test, probably because residual fiber tension is difficult to avoid in experimental PFT. More detailed comparison will be performed in the next section, especially looking at gaps and stitch tension, but these results already demonstrate that the macroscopic behavior of NCF can be deduced from mesoscopic finite element modeling, with mesoscopic mechanical properties that have been characterized experimentally on yarns and stitches.

# 6. Macroscopic applications

In this section, the performance of the model will be illustrated on a larger scale, using a Bias Extension Test (BET). An experimental BET was carried out on Bi-NCF, and was then simulated by reproducing the exact shape and boundary conditions of the BET. The two are compared
<figure><img src="image_16.png" /><figcaption>Fig. 14. Picture Frame Test on Bi-NCF (a) Initial configuration on experiments (b) schematic representation of Bi-NCF inside the PFT.</figcaption></figure>

<figure><img src="image_17.png" /><figcaption>Fig. 15. Normalized shear force of experimental PFT and unit cell shear simulation. Error bars for three tests.</figcaption></figure>

in terms of shear force, and the gaps between yarns will be analyzed, in addition to the tension in the stitches.

## 6.1. Experimental BET

By comparing with PFT, the BET is a simpler test to perform since it eliminates the potential for spurious tension. Five BETs have been conducted on 240 mm × 80 mm samples of Bi-NCF. The samples were clamped in a universal testing machine (Lloyd 1 kN with a 100 N sensor), with the fibers initially positioned at ±45° to the loading direction.

The machine moved at a constant rate of 30 mm/min until it reached a total displacement of 60 mm. Global view of the setting is proposed in Fig. 17.

Fig. 17(b) illustrates the theoretical solution of BET when no slip-page occur at any fiber intersection [65]. This is mostly the case in woven reinforcements, for which shear angle can be calculated from

the machine displacement d by using Eq. (18):

$$ \gamma = \frac{\pi}{2} - 2 \cos^{-1} \left( \frac{D+d}{\sqrt{2}D} \right) \qquad (18) $$

where D = H - W represents the length of the central shear zone C, obtained from the height, H, and the width, W, of the sample.

In the case of NCF, slippage cannot be disregarded [57,58], and Eq. (18) is not necessarily verified. Therefore, to accurately measure the shear angle γ in the central zone C and visualize the gap between fiber yarns, a zoom-in camera was placed on the central zone. A set of images, captured every 2 mm displacement on the machine, was collected. Several yarns were marked with a white dotted line and their edges were subsequently tracked in Fiji software. This technique allowed to accurately measure the shear angles, widths, and any changes in gaps between the marked yarns during deformation.

In order to compare the results with PFT, a normalized shear force was calculated in BET using Eq. (19), which comes from the virtual energy theorem applied to a three zone sample. The calculation was done with the angle γ directly measured on the NCF. Note that Eq. (19) remains valid only in the case where zone A exhibits no shearing, and zone B exhibits half shear of zone C. The non-linearity of Eq. (19) [63] was solved by using an iterative algorithm.

$$ F_{\text{nor}}(\gamma) = \frac{1}{(2H - 3W) \cos \gamma} \times \left( \left( \frac{H}{W} - 1 \right) F_m \left( \cos \frac{\gamma}{2} - \sin \frac{\gamma}{2} \right) - W F_{\text{nor}} \left( \frac{\gamma}{2} \right) \cos \frac{\gamma}{2} \right) \qquad (19) $$

## 6.2. Numerical simulations

The possibilities of the model will now be demonstrated on a larger scale, i.e. on an entire BET specimen. The size of the model is 240 mm × 80 mm with 47 941 solid elements for fiber yarns and 51 249 beam elements for stitches. All the nodes on the bottom side are clamped, and a displacement of 60 mm is imposed on the top nodes. The computation was launched on a cluster that consolidates independent computers with a global management to utilize abundant resources. It took 80 h with 20 CPU cores and 64 GB RAM. The results are shown in Fig. 18; both sides of the Bi-NCF (warp and weft) are
<figure><img src="image_18.png" /><figcaption>Fig. 16. Deformed configuration of unit-cell under pure shear test (a) upper side, warp yarn with zig-zag stitch (b) lower side, weft yarn with looped stitch.</figcaption></figure>

<figure><img src="image_19.png" /><figcaption>Fig. 17. Bias Extension Test (a) Initial configuration with sample dimensions 240 mm × 80 mm (b) Theoretical solution of BET with zones A (no-shear), B (half shear) and C (maximum shear).</figcaption></figure>

presented in Figs. 18(a) and 18(b). For comparison, the warp view of the corresponding experimental test is given in Fig. 18(c).

## 6.3. *Experimental to numerical comparison*

First observations show that the mesoscopic model is able to capture the sliding between warp and weft yarns. This is particularly noticeable on the edges of the sample in the lower B zone (Fig. 18(b)). This phenomenon was also observed in experimental results (Fig. 18(d)) where unraveling of stitches occurred. Moreover, the model is able to capture the tension level in the stitches. Figs. 18(e) and 18(f) show a zoom-in view of the specimen around the central zone. The model faithfully reproduces the experimental results with half of the stitches being in tension, while the other half undergoes compression and buckling.

In order to deepen the analysis of the results, data of the BET were also converted into normalized shear force curves (see Fig. 21(d)). In order to do so, a specific attention was put into BET post-processing. As already mentioned large slippage occurs in BET between the two

plies. For this reason, it is not straightforward to observe shear angle. A specific post-processing of the simulation results is introduced here to reconstruct the shear angle value. Yarn orientation is extracted from FEM results (see Fig. 19(a)); then crossing points between warp and weft yarns are located in the mid-plane of the sample. Using yarn orientation map at these points, a shear angle field was calculated. The resulting shear angle map is shown in Fig. 19(b). Note that three shear zones are indeed distinguishable for BET on NCF.

This effective shear angle was compared to the theoretical shear angle that would be obtained without slippage (Eq. (18)). The curves shown in Fig. 20 demonstrate that numerical and experimental shear angles are in very good agreement, but also that both are far from the so-called theoretical value, indicating the occurrence of slippage within the Bi-NCF. The discrepancy tends to increase as the shear angle reaches 25°.

In Fig. 21, further experimental to numerical comparisons are shown. Figs. 21(a) to 21(c) focus on the geometrical aspects (width of yarns, gaps, and their sum called broadness). The yarn undergoes continuous compression due to the tension in the stitch, resulting in an
<figure><img src="image_20.png" /><figcaption>**Fig. 18.** Results of BET: Deformed configuration of simulated (a, warp yarn, and b, weft yarn) and experimental BET (c, warp yarn); (d) zoom-in view of the experimental edge revealing the occurrence of slips; zoom-in view of the simulated (e) and experimental (f) central shear zones with the presence of both under-tension and bent stitches.</figcaption></figure>

<figure><img src="image_21.png" /><figcaption>**Fig. 19.** Calculation of effective shear angle (in degree) in FEM results; (a) Yarn orientation of weft yarn (red) and warp yarn (blue); (b) Shear angle between the two directions, in degree.</figcaption></figure>

initial increase in the gap prior to a shear angle of 30°. The broadness, considered to be the combination of the gap and yarn width, helps validating the measurement of both gap and width. These results show the effectiveness of the model in describing geometrical evolution of the NCF during shear tests. Finally, Fig. 21(d) compares the results in terms of normalized shear force, demonstrating great consistency of the results, and once again validating the relevance of the model.

<figure><img src="image_22.png" /><figcaption>**Fig. 20.** Comparison of experimental and simulated effective shear angles in BET with theoretical shear angle in pure shear tests, without slippage. Error bars for five experiments.</figcaption></figure>

## 7. Conclusion

In this paper, a finite element model for Bi-NCF was developed. It combines 3D elements with hyperelastic transversely isotropic behavior
<figure><img src="image_23.png" /><figcaption>Fig. 21. Comparison between simulation results and experimental tests: (a) Width of yarn vs. shear angle; (b) Development of gap vs. shear angle; (c) Broadness vs. shear angle; (d) shear force vs. shear angle. Error bars for five experiments.</figcaption></figure>

for the yarns and uses beam elements with a stress resultant approach for the stitches. Objectives were twofold: firstly, to represent NCF at both the mesoscopic scale (representative unit cell) and at larger scale (specimen level); secondly, to detect mesoscopic defects, like gaps between yarns.

Our comprehensive approach was validated through comparison with picture frame and bias extension tests, effectively replicating observed macroscopic shear behavior. This demonstrates its potential for virtual material characterization and meso-to-macro homogenization, reducing the need for extensive experimental testing. Additionally, mesoscopic simulations provided valuable insights into mesoscopic defects like gaping, critical in NCF composite forming, showcasing the model's ability to capture local phenomena. Moreover, the model allows flexibility in adjusting fiber yarn and stitch geometry for analyzing various stitching patterns and NCF structures

Next stage of this work is to simulate the forming process using the presented approach. This remains feasible with the computing capabilities currently available. However, the performance of the calculations presented here may represent a limit to their application in the industrial context. Other improvements of the model could be proposed. For instance, the inclusion of bending rigidity of the fibers in the constitutive model could enhance model accuracy, and it would be important to carefully assess its impact. Additionally, the current model relies on the assumption that the fabric is periodic; different ways of introducing variability into the geometry could be proposed in future work.

## CRediT authorship contribution statement

**Ruochen Zheng:** Writing – review & editing, Writing – original draft, Validation, Software, Methodology, Investigation, Formal analysis, Data curation, Conceptualization. **Naim Naouar:** Writing – review & editing, Writing – original draft, Validation, Software, Methodology, Investigation, Conceptualization. **Julien Colmars:** Writing – review & editing, Writing – original draft, Validation, Supervision, Methodology, Investigation, Conceptualization. **Auriane Platzer:** Writing – review & editing, Writing – original draft, Visualization, Supervision. **Bastian Schäfer:** Writing – review & editing, Software, Resources. **Fabrice Morestin:** Supervision. **Luise Kärger:** Writing – review & editing, Project administration, Funding acquisition. **Philippe Boisse:** Supervision, Project administration, Funding acquisition, Conceptualization.

## Declaration of competing interest

The authors declare the following financial interests/personal relationships which may be considered as potential competing interests: Ruochen Zheng reports financial support was provided by French National Research Agency.

## Data availability

Data will be made available on request.
Acknowledgments

This work was supported by the French National Research Agency, France (ANR) and the Deutsche Forschungsgemeinschaft, Germany (DFG, German Research Foundation) in the collaborative project AME-COMP "Composite forming simulation for non-crimp fabrics based on generalized continuum approaches, Germany" [grant number ANR-19-CE06-0031-012]. A CC-BYpubliccopyrightlicense has been applied by the authors to the present document and will be applied to all subsequent versions up to the Author Accepted Manuscript arising from this submission, in accordance with the grant's open access conditions.

References

[1] Mouritz AP. 1 - introduction to aerospace materials. In: Introduction to aerospace materials. Woodhead Publishing; 2012, p. 1-14.

[2] Lomov S. 4 - understanding and modelling the effect of stitching on the geometry of non-crimp fabrics. In: Lomov SV, editor. Non-crimp fabric composites. Woodhead publishing series in composites science and engineering, Woodhead Publishing, 2011, p. 84-102.

[3] Kimpara I. Use of advanced composite materials in marine vehicles. Mar Struct 1991;4(2):117-27.

[4] Sihn S, Kim RY, Kawabe K, Tsai SW. Experimental studies of thin-ply laminated composites. Compos Sci Technol 2007;67(6):996-1008.

[5] Budwal N, Kasper K, Goering J, Ward C. Flexible low-cost tooling solutions for a one-shot resin infusion of a 3D woven and multi-textile preform. Procedia Manuf 2020;51:856-63.

[6] Advani S, Sozer E, Mishnaevsky Jr. L. Process modeling in composite manufacturing. Appl Mech Rev 2003;56.

[7] Feraboli P, Masini A. Development of carbon/epoxy structural components for a high performance vehicle. Composites B 2004;35(4):323-30.

[8] González A, Graciani E, París F. Prediction of in-plane stiffness properties of non-crimp fabric laminates by means of 3D finite element analysis. Compos Sci Technol 2008;68(1):121-31.

[9] Mikhaluk D, Truong T, Borovkov A, Lomov S, Verpoest I. Experimental observations and finite element modelling of damage initiation and evolution in carbon/epoxy non-crimp fabric composites. Eng Fract Mech 2008;75(9):2751-66.

[10] Edgren F, Mattsson D, Asp L, Varna J. Formation of damage and its effects on non-crimp fabric reinforced composites loaded in tension. Compos Sci Technol 2004;64(5):675-92.

[11] Mattsson D, Joffe R, Varna J. Damage in NCF composites under tension: Effect of layer stacking sequence. Eng Fract Mech 2008;75(9):2666-82.

[12] Drapier S, Wisnom M. Finite-element investigation of the compressive strength of non-crimp-fabric-based composites. Compos Sci Technol 1999;59(8):1287-97.

[13] Li X, Hallett S, Wisnom M. Modelling the effect of gaps and overlaps in automated fibre placement (AFP) manufactured laminates. Sci Eng Compos Mater 2015;22:115-29.

[14] Arnold S, Sutcliffe M, Oram W. Experimental measurement of wrinkle formation during draping of non-crimp fabric. Composites A 2016;82:159-69.

[15] Shen H, Wang P, Legrand X, Liu L. Characterisation and optimisation of wrinkling during the forming of tufted three-dimensional composite preforms. Composites A 2019;127:105651.

[16] Viisainen J, Hosseini A, Sutcliffe M. Experimental investigation, using 3D digital image correlation, into the effect of component geometry on the wrinkling behaviour and the wrinkling mechanisms of a biaxial NCF during preforming. Composites A 2021;142:106248.

[17] Ghazimoradi M, Trejo EA, Butcher C, Montesano J. Characterizing the macroscopic response and local deformation mechanisms of a unidirectional non-crimp fabric. Composites A 2022;156:106857.

[18] Pickett A, Creeg G, de Luca P. Simplified and advanced simulation methods for prediction of fabric draping. Revue Européenne des Éléments Finis 2005;14:677-91.

[19] Creeg G, Pickett A. Meso-modelling of non-crimp fabric composites for coupled drape and failure analysis. J Mater Sci 2006;41:6725-36.

[20] Sirtautas J, Pickett A, Lépicier P. A mesoscopic model for coupled drape-infusion simulation of biaxial non-crimp fabric. Composites B 2013;47:48-57.

[21] Drapier S, Wisnom M. A finite-element investigation of the interlaminar shear behaviour of non-crimp-fabric-based composites. Compos Sci Technol 1999;59(16):2351-62.

[22] Edgren F, Asp L. Approximate analytical constitutive model for non-crimp fabric composites. Composites A 2005;36(2):173-81.

[23] Hes H, Roth Y, Himmel N. Elastic constants estimation of stitched NCF CFRP laminates based on a finite element unit-cell model. Compos Sci Technol 2007;67(6):1081-95.

[24] Pierreux G, Hemelrijck DV, Massart TJ. Automated RVE computations for evaluation of microdamage initiation in structural stitched non-crimp fabric composites. J Compos Mater 2020;54(30):4751-71.

[25] Colin D, Bel S, Hans T, Hartmann M, Drechsler K. Virtual description of non-crimp fabrics at the scale of filaments including orientation variability in the fibrous layers. Appl Compos Mater 2020;27(4):337-55.

[26] Thompson AJ, El Said B, Belnoue JP-H, Hallett SR. Modelling process induced deformations in 0/90 non-crimp fabrics at the meso-scale. Compos Sci Technol 2018;168:104-10.

[27] Long A, Brown L. 8 - modelling the geometry of textile reinforcements for composites: TexGen. In: Boisse P, editor. Composite reinforcements for optimum performance. Woodhead publishing series in composites science and engineering. Woodhead Publishing; 2011, p. 239-64.

[28] Verpoest I, Lomov SV. Virtual textile composites software WiseTex: Integration with micro-mechanical, permeability and structural analysis. Compos Sci Technol 2005;65(15):2563-74, 20th Anniversary Special Issue.

[29] Lomov S, Verpoest I, Cichosz J, Hahn C, Ivanov D. Meso-level textile composites simulations: Open data exchange and scripting. J Compos Mater 2014;48:621-37.

[30] Naouar N, Vidal-Sallé E, Schneider J, Maire E, Boisse P. Meso-scale FE analyses of textile composite reinforcement deformation based on X-ray computed tomography. Compos Struct 2014;116:165-76.

[31] Naouar N, Vidal-Salle E, Schneider J, Maire E, Boisse P. 3D composite reinforcement meso f.e. analyses based on X-ray computed tomography. Compos Struct 2015;132:1094-104.

[32] Liu Y, Straumit I, Vasiukov D, Lomov SV, Panier S. Prediction of linear and non-linear behavior of 3D woven composite using mesoscopic voxel models reconstructed from X-ray micro-tomography. Compos Struct 2017;179:568-70.

[33] Naouar N, Vasiukov D, Park CH, Lomov SV, Boisse P. Meso-FE modelling of textile composites and X-ray tomography. J Mater Sci 2020;55:16969-89.

[34] Wu W, Li W. Parametric modeling based on the real geometry of glass fiber unidirectional non-crimp fabric. Textile Res J 2019;89:004051751882484.

[35] Du G, Ko F. Analysis of multiaxial warp-knit preforms for composite reinforcement. Compos Sci Technol 1996;56(3):253-60.

[36] Lomov S, Belov E, Bischoff T, Ghosh S, Truong Chi T, Verpoest I. Carbon composites based on multiaxial multiply stitched preforms. Part 1. Geometry of the preform. Composites A 2002;33(9):1171-83.

[37] Schindelin J, et al. Fiji: an open-source platform for biological-image analysis. Nature Methods 2012;9(7):676-82.

[38] Lomov S, Huysmans G, Luo Y, Parnas R, Prodromou A, Verpoest I, et al. Textile composites: modelling strategies. Composites A 2001;32(10):1379-94.

[39] Ishikawa T, Chou T-W. Nonlinear behavior of woven fabric composites. J Compos Mater 1983;17(5):399-413.

[40] Krieger H, Gries T, Stapleton SE. Shear and drape behavior of non-crimp fabrics based on stitching geometry. Int J Mater Form 2018.

[41] Robitaille F, Clayton BR, Long AC, Souter BJ, Rudd CD. Geometric modelling of industrial preforms: Warp-knitted textiles. Proc Inst Mech Eng L 2000;214(2):71-90.

[42] Lomov S, Verpoest I, Peeters T, Roosevelt D, Zako M. Nesting in textile laminates: geometrical modelling of the laminate. Compos Sci Technol 2003;63(7):993-1007.

[43] Jiang G, Gu L, Honglian C, Miao X, Zhang A, Gao Z. Geometric model for multi-axial warp-knitted fabric based on NURBS. Fibres Textiles Eastern Eur 2014;105:91-7.

[44] Peng X, Cao J. A dual homogenization and finite element approach for material characterization of textile composites. Composites B 2002;33(1):45-56.

[45] Lomov SV, Ivanov DS, Verpoest I, Zako M, Kurashiki T, Nakai H, et al. Meso-FE modelling of textile composites: Road map, data flow and algorithms. Compos Sci Technol 2007;67(9):1870-91.

[46] Badel P, Vidal-Sallé E, Boisse P. Computational determination of in-plane shear mechanical behaviour of textile composite reinforcements. Comput Mater Sci 2007;40(4):439-48.

[47] Charmetant A, Vidal-Sallé E, Boisse P. Hyperelastic modelling for mesoscopic analyses of composite reinforcements. Compos Sci Technol 2011;71(14):1623-31.

[48] Latil P, Orgéas L, Geindreau C, Dumont PJJ, Rolland du Roscoat S. Towards the 3D in situ characterisation of deformation micro-mechanisms within a compressed bundle of fibres. Compos Sci Technol 2011;71(4):480-8.

[49] Spencer AJM. Part III - theory of invariants. In: Eringen AC, editor. Mathematics. Academic Press; 1971, p. 239-353.

[50] Iwata A, Inoue T, Naouar N, Boisse P, Lomov SV. Coupled meso-macro simulation of woven fabric local deformation during draping. Composites A 2019;118:267-80.

[51] Wendling-Hivet A, Ferré MR, Allaoui S, Nunez R, Loison S, Hivet G. Study of the cohesion of carbon fiber yarns: in-plane shear behavior. Int J Mater Form 2017;10(5):671-83.

[52] FTPeirce BSc, FInstP, FTI. 26-The "handle" of cloth as a measurable quantity. J Textile Inst Trans 1930;21(9):T377-416.
[53] Boisse P, Colmars J, Hamila N, Naouar N, Steer Q. Bending and wrinkling of composite fiber preforms and prepregs. a review and new developments in the draping simulations. Composites B 2018;141:234-49.

[54] Bai R, Chen B, Colmars J, Boisse P. Physics-based evaluation of the drapability of textile composite reinforcements. Composites B 2022;242:110089.

[55] Badel P, Vidal-Sallé E, Maire E, Boisse P. Simulation and tomography analysis of textile composite reinforcement deformation at the mesoscopic scale. Compos Sci Technol 2008;68(12):2433-40.

[56] Kong H, Mouritz A, Paton R. Tensile extension properties and deformation mechanisms of multiaxial non-crimp fabrics. Compos Struct 2004;66(1):249-59.

[57] Bel S, Hamila N, Boisse P, Dumont F. Finite element model for NCF composite reinforcement preforming: Importance of inter-ply sliding. Composites A 2012;43(12):2269-77.

[58] ten Thije R, Akkerman R. Finite element simulation of draping with non-crimp fabrics. In: Verijenko V, Adali S, Morozov E, von Klemperer C, editors. Proceedings of the 15th ICCM conference (ICCM-15). Elsevier; 2005.

[59] Pham MQ, Wendt E, Häntzsche E, Gereke T, Cherif C. Numerical modeling of the mechanical behavior of textile structures on the meso-scale for forming process simulations of composite 3D preforms. Eng Rep 2020;4.

[60] Lomov S, Willems A, Verpoest I, Zhu Y, Barburski M, Stoilova T. Picture frame test of woven composite reinforcements with a full-field strain registration. Textile Res J 2006;76:243-52.

[61] Li L, Zhao Y, gia-nam Vuong H, Chen Y, Yang J, Duan Y. In-plane shear investigation of biaxial carbon non-crimp fabrics with experimental tests and finite element modeling. Mater Des 2014;63:757-65.

[62] Lomov S, Barburski M, Stoilova T, Verpoest I, Akkerman R, Loendersloot R, et al. Carbon composites based on multiaxial multiply stitched preforms. Part 3: Biaxial tension, picture frame and compression tests of the preforms. Composites A 2005;36(9):1188-206.

[63] Cao J, Akkerman R, Boisse P, Chen J, Cheng H, de Graaf E, et al. Characterization of mechanical behavior of woven fabrics: Experimental methods and benchmark results. Composites A 2008;39(6):1037-53.

[64] Peng X, Cao J, Chen J, Xue P, Lussier D, Liu L. Experimental and numerical analysis on normalization of picture frame tests for composite materials. Compos Sci Technol 2004;64(1):11-21.

[65] Boisse P, Hamila N, Guzman-Maldonado E, Madeo A, Hivet G, dell'Isola F. The bias-extension test for the analysis of in-plane shear properties of textile composite reinforcements and prepregs: a review. Int J Mater Form 2017;10(4):473-92.