# Dual-scale flow simulation for LCM processes

William Han, Quentin Govignon, Arthur Cantarel, Fabrice Schmidt

► To cite this version:

William Han, Quentin Govignon, Arthur Cantarel, Fabrice Schmidt. Dual-scale flow simulation for LCM processes. ICCM23 - 23rd international conference on composite materials, Jul 2023, Belfast, United Kingdom. hal-04198689

HAL Id: hal-04198689

https://imt-mines-albi.hal.science/hal-04198689v1

Submitted on 1 Oct 2023

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 CC BY 4.0 - Attribution - International License
# DUAL-SCALE FLOW SIMULATION FOR LCM PROCESSES

W. Han<fnref n="1" />, Q. Govignon<fnref n="2" />, A. Cantarel<fnref n="3" /> and F. Schmidt<fnref n="4" />

**Keywords:** Numerical simulation, dual-scale flow, Liquid Composite Moulding, realistic geometry

## ABSTRACT

This study aims to propose a simulation method for liquid composite moulding processes. The focus is done on injection in a continuous fibre preform. In this case, the injection is realised on a dual-scale porous media, as the flow regime of the resin differs depending on if it goes between fibres inside tows or between tows. To handle the two scales during a resin filling simulation, Brinkman's equation is used with a saturation conservation equation in OpenFOAM®, the open-source finite volume method toolbox. It allows to simulate the flow in both the region between tows and inside tows at meso-scale. To demonstrate the capability of the flow simulation method, a simulation is realised on a realistic textile geometry of the resin flow through a fibrous preform. In a future work the aim is to combine this simulation method with the effects of the kinetics of resins that can be used for the resin transfer moulding processes, such as a reactive ε-caprolactam resin to synthesize PA6.

## 1 INTRODUCTION

Liquid composite moulding (LCM) processes are often used as they permit manufacturing composite materials with a complex geometry with moderate material and energy costs. In these processes, it is crucial to understand how the reinforcement is filled. Indeed, the choice of processing conditions, and the resin cure during and after the injection can affect the final quality and mechanical properties of the manufactured composite.

For flow in a single-scale porosity medium, using Darcy's law is enough to simulate the flow at macro-scale. This happens when impregnating a medium with randomly oriented fibres. With a continuous fibre preform, fibre tows are periodically and geometrically arranged. Such a medium is said to be a dual-scale porosity medium, as the resin flows both between fibres inside tows, and between tows themselves. However, for dual-scale mediums, Imbert et al. showed that a different approach is needed because of the slow time-scale the resin impregnate the tows [1].

One solution is to study the flow at meso-scale, where the tows can be modelled as a porous media, while the flow in channels between tows will follow Stokes' law. This approach has been proposed before in works realized by Blais et al. and Li et al. [2,3]. They used the Beavers-Joseph-Saffman between the porous domain and the non-porous domain in order to couple Darcy's flow with Stokes's flow.

In this work, we propose an approach based on Brinkmann's equation (2), which allow to use the same equation in the whole media. The equation is equivalent to Stokes' equation outside the tows, while the Darcian term inside Brinkmann's equation directs the flow inside the tows. By coupling it with a saturation conservation equation (3), the filling of the dual-scale media can be simulated, even with 3D geometries. To showcase this, a real textile geometry obtained by micro-tomography [4] is incorporated in a filling simulation including both the intra-tow flow and the inter-tow flow.

<footnote n="1">¹ Mines Paris, Université PSL, Centre des Matériaux (MAT), CNRS UMR 7633, BP 87 F-91003
Evry, France, [william.han@minesparis.psl.eu](mailto:william.han@minesparis.psl.eu), [https://www.mat.minesparis.psl.eu](https://www.mat.minesparis.psl.eu)</footnote>

<footnote n="1">¹,²,³,⁴ Institut Clément Ader (ICA), Université de Toulouse, CNRS, IMT Mines Albi, INSA, ISAE-SUPAERO, UPS, France, [https://ica.cnrs.fr](https://ica.cnrs.fr)</footnote>

<footnote n="2">² IMT Mines Albi, Campus Jarlard, F-81013 Albi, France, [quentin.govignon@mines-albi.fr](mailto:quentin.govignon@mines-albi.fr)</footnote>

<footnote n="3">³ IUT de Tarbes, 1 rue Lautréamont, F-65016 Tarbes, France, [arthur.cantarel@iut-tarbes.fr](mailto:arthur.cantarel@iut-tarbes.fr)</footnote>

<footnote n="4">⁴ IMT Mines Albi, Campus Jarlard, F-81013 Albi, France, [fabrice.schmidt@mines-albi.fr](mailto:fabrice.schmidt@mines-albi.fr)</footnote>
## 2 THEORY

### 2.1 Governing equations

The continuity equation (1) and the Navier-Stokes-Brinkmann equation (2) are used to compute the averaged velocity $\bar{u}$ and pressure $p$.

$$ \nabla \cdot \bar{u} = 0 \qquad (1) $$

$$ \frac{\rho}{\epsilon} \left[ \frac{\partial \bar{u}}{\partial t} + \nabla \cdot \left( \frac{\bar{u} \cdot \bar{u}}{\epsilon} \right) - \nu \Delta \bar{u} \right] = - \nabla p + \eta K^{-1} \bar{u} \qquad (2) $$

$\epsilon$ is the average porosity, $\rho$ the density, $\nu$ and $\eta$ are respectively the kinematic and the dynamic viscosity, and $K$ is the permeability tensor.

The saturation $S_r$ of the resin is calculated using the conservation equation (3). It is a modified form of the Volume of Fluid (VOF) method adapted to porous mediums [5].

$$ \frac{\partial \epsilon S_r}{\partial t} + \nabla \cdot (\bar{u} S_r) + \nabla \cdot \left( (\bar{u}_r - \bar{u}_g) \epsilon S_r S_g \right) = 0 \qquad (3) $$

The terms $\bar{u}_r$ and $\bar{u}_g$ respectively refer to the resin and the gas averaged velocities and are used to calculate an interface compression term [6], while $S_g = 1 - S_r$ is the gas saturation. Surface tension forces are not included in this work.

### 2.2 Method

The finite volume method is used with CFD toolbox OpenFOAM® 8.x to discretize the governing equations in a hexahedral mesh. The Crank-Nicolson temporal discretization scheme is used with a constant time step $\delta t = 10^{-5}$ s. In equation (2), the viscous term is calculated with the central differencing method, while the convection term is calculated with the “limitedLinearV” scheme used with maximum emphasis on convergence. The Darcy term is solved explicitly. The divergence terms in equation (3) are calculated with a van Leer discretization scheme. Those are second order spatial discretization schemes commonly used in OpenFOAM® that satisfies TVD (Total Variation Diminishing) conditions for stability and convergence. The Crank-Nicolson temporal discretization scheme is used. The time step $\delta t$ is chosen adaptively with equation (4) using the velocity and the characteristic length $\delta l$ of the most critical element. It calculates the Courant number $Co$ which is kept below 0.05 in order to ensure accuracy and convergence of the calculations.

$$ Co = \frac{|\bar{u}|\delta t}{\delta l} < 0.05 \qquad (4) $$

## 3 RESIN FLOW SIMULATION WITH A REALISTIC GEOMETRY

### 3.1 Geometry

The geometry (Figure 1) is a 3D representation of a four layers stack E-glass plain weave specimen, voxelated by Wijaya et al. [4] using an X-ray micro tomography (µCT). It has a global fibre volume fraction of 50%.

### 3.2 Permeability of the tows

The fibre volume fraction (FVF) was also obtained for each tow's voxel by Wijaya [7] which allowed longitudinal and transverse permeabilities to be approximated using Gebart theory [8]. To simplify implementation in OpenFOAM®, a unique fibre volume fraction and permeability tensor has been set respectively for warp tows and weft tows. For the FVF, the arithmetic average was calculated, giving $V_{f,warp} = 0.602$ and $V_{f,weft} = 0.667$. For the permeabilities, due to the non-periodic nature of the voxels' permeabilities values, the geometric average has been chosen, as it is recommended for averaging in media with random permeability distribution in reservoir engineering [9]. Thus, equation (5) where $h_i$ and $k_i$ are respectively the thickness and the permeability of a voxel has been
used to calculate the warp permeability tensor $\bar{K}_{tows,warp}$ and weft permeability tensor $\bar{K}_{tows,weft}$ (equations (6) and (7)).

$$ \bar{K}_g = \exp \left[ \frac{\sum_{i=1}^{n} (h_i \ln(k_i))}{\sum_{i=1}^{n} h_i} \right] \qquad (5) $$

$$ \bar{K}_{tows,warp} = \begin{bmatrix} 3.69 \cdot 10^{-12} & 0 & 0 \\ 0 & 7.68 \cdot 10^{-13} & 0 \\ 0 & 0 & 7.68 \cdot 10^{-13} \end{bmatrix} \text{m}^2 \qquad (6) $$

$$ \bar{K}_{tows,weft} = \begin{bmatrix} 1.78 \cdot 10^{-12} & 0 & 0 \\ 0 & 3.59 \cdot 10^{-13} & 0 \\ 0 & 0 & 3.59 \cdot 10^{-13} \end{bmatrix} \text{m}^2 \qquad (7) $$

<figure><img src="image_2.png" /><figcaption>Figure 1: 3D tow geometries of the textile: a/ weft tows (orange), b/ warp tows (green), c/ both tows and dimensions.</figcaption></figure>

### 3.3 Mesh and boundary conditions

The boundary conditions are described in Figure 2. The upper and lower surfaces represent the mould. They have a no-slip boundary condition, and their normal pressure gradient is set to zero. The transverse surface either have inlet and outlet boundary conditions or symmetry conditions, to enforce the periodicity of the textile specimen.

<figure><img src="image_3.png" /><figcaption>Figure 2. Simulation geometry and boundary conditions.</figcaption></figure>
The flow is directed by a pressure gradient $\Delta P = 100$ Pa between the inlet and the outlet. The resin enters through the inlet ($S_r = 1$ at the inlet) with a dynamic viscosity of $\eta = 2$ mPa.s, which is the viscosity of the $\varepsilon$-caprolactam reactive mix used to synthesize PA6 at 433 K [10]. The mesh has 154,283 orthogonal hexahedral elements, and a ratio of 1/8 between the element height and its length has been chosen for good geometry description [7].

Initially, the domain is filled by gas, which was set using the viscosity and density of air at 433 K ($\eta_g = 2.42 \cdot 10^{-5}$ Pa s, $\rho_g = 0.815$ kg m$^{-3}$). The simulation was conducted for a filling duration of 200 seconds.

## 3.4 Results and discussion

The results of the filling simulation are shown in Figure 3, with the tows represented in light grey and the resin in red.

<figure><img src="image_4.png" /><figcaption>Figure 3. Filling simulation at different times. The tows are plotted in light grey, while the resin is represented in red.</figcaption></figure>

At $t = 20$ s, the resin has already reached the other end of the domain. However, as the global resin saturation plot over time in Figure 4 shows, the saturation gets slower and slower as the injection progresses, barely reaching 0.94 at $t = 200$ s.

<figure><img src="image_5.png" /><figcaption>Figure 4. Global resin saturation of the domain during the simulation.</figcaption></figure>
To better observe the impregnation of the tows, a cross-section of the domain following the flow direction has been selected in Figure 5, and the saturation at different times of the injection is shown in Figure 6.

<figure><img src="image_6.png" /><figcaption>Figure 5. Position and cross-sectional observation of the textile geometry along warp tows (in green) and perpendicularly to weft tows (in orange). The yellow arrow indicates where the pressure curves in Figure 7 are sampled.</figcaption></figure>

<figure><img src="image_7.png" /><figcaption>Figure 6. Resin saturation at different times for the cross section indicated in Figure 5.</figcaption></figure>
At the start ($t = 0.2$ s), in this cross-section, the resin has only started to impregnate the tows. Since the warp tow longitudinal permeability ($3.69 \cdot 10^{-12}$ m²) is 10 times higher than the weft tow transverse permeability ($3.59 \cdot 10^{-13}$), the front has progressed further in the warp tows. At $t = 2$ s, while the first row of tows is still being impregnated, resin from the other part of the domain has arrived through the inter-tow channels: we can observe that this indeed happens in Figure 3. At $t = 20$ s, all the inter-tow channels are filled. We can see between the state at $t = 20$ s and $t = 100$ s how the air gets trapped as the resin circles the most difficult to reach parts of the tows. Between $t = 100$ s and $t = 200$ s, we can see that some air bubbles have fully formed, which will be only eliminated when they are pushed out of the domain.

In Figure 7, the pressure curves are drawn along the yellow arrow indicated in Figure 5 at different times. At $t = 0.02$ s and $t = 0.2$ s, the domain is mostly unsaturated, and thus the pressure drop is mostly located at the inlet. However, at $t = 2$ s, the first row of inter-tow channels is filled, which creates a pressure plateau around 50 Pa because of the unsaturated tows between the inlet and the resin front. At $t = 20$ s, four of these plateaus can be observed (counting the one at 0 Pa) which corresponds to all the tow spaces that are unsaturated between the inlet and the outlet. As these spaces slowly saturates, the pressure slowly reaches a linear profile, but since there is still trapped air, the plateaus did not disappear at $t = 200$ s.

<figure><img src="image_8.png" /><figcaption>Figure 7. Pressure profile along the arrow indicated in Figure 5.</figcaption></figure>

# 4 CONCLUSIONS AND PERSPECTIVES

This work presents a simulation based on Brinkmann's equation with a saturation conservation equation to simulate the resin injection during LCM processes at meso-scale. It has been realized on a realistic tow geometry. The progress of the resin saturation and the pressure profile has been studied. The formation of air bubbles is observed, and conforming to Kim et al. observations [11], a non linear pressure profile is observed, especially when the inter-tow channels are partially saturated. The high-resolution of the simulation has permitted to observe the piecewise linearity of the curve due to the air pockets in conditions where capillary forces are negligible before viscous forces.

This simulation can therefore be used as a base for the full process simulation of a reactive thermoplastic resin. To realize this, resin cure equations and rheokinetics for the synthesis of the PA6 will be coupled to the saturation equations in future works.

To improve the accuracy of the simulation, it will also be interesting to account for capillary forces that arises at the interface between channels and tows, and to devise an experimental confrontation.

# ACKNOWLEDGEMENTS

Financial support from the Occitanie Region is acknowledged. The authors wish to thank Wilsen Wijaya and the University of Auckland for graciously providing the textile geometry. The authors also acknowledge Olivier de Almeida for his help in the project.
REFERENCES

[1] M. Imbert, S. Comas-Cardona, E. Abisset-Chavanne and D. Prono, Experimental investigation of intra-tow fluid storage mechanisms in dual-scale fiber reinforcements, *Composites Part A: Applied Science and Manufacturing, **107*, 2018, pp. 70–82.
(doi: 10.1016/j.compositesa.2017.12.015).

[2] M. Blais, N. Moulin, P.-J. Liotier and S. Drapier, Resin infusion-based processes simulation : coupled Stokes-Darcy flows in orthotropic preforms undergoing finite strain, *International Journal of Material Forming, **10*, 2017, pp. 43–54. (doi: 10.1007/s12289-015-1259-2).

[3] C. Li, A. Cantarel and X. Gong, A two phases coupling free and porous flow method at multi-scales of LCM process, in: *Proceedings of the 2019 International Conference on Composite Materials*, Melbourne, Australia, 2019.

[4] W. Wijaya, M.A. Ali, R. Umer, K.A. Khan, P.A. Kelly and S. Bickerton, An automatic methodology to CT-scans of 2D woven textile fabrics to structured finite element and voxel meshes, *Composites Part A: Applied Science and Manufacturing, **125*, 2019, 105561. (doi: 10.1016/j.compositesa.2019.105561).

[5] F.J. Carrillo, I.C. Bourg and C. Soulaine, Multiphase flow modeling in multiscale porous media: *An open-source micro-continuum approach*, Journal of Computational Physics: X, **8**, 2020, 100073. (doi: 10.1016/j.jcpx.2020.100073).

[6] K.E. Wardle and H.G. Weller, Hybrid Multiphase CFD Solver for Coupled Dispersed/Segregated Flows in Liquid-Liquid Extraction, *International Journal of Chemical Engineering, **2013***, 2013, pp. 1–13. (doi: 10.1155/2013/128936).

[7] W. Wijaya, Permeability of 2D Woven Composite Textile Reinforcements: Textile Geometry and Compaction, and Flow Modelling, PhD Thesis, The University of Auckland, 2020.

[8] B.R. Gebart, Permeability of Unidirectional Reinforcements for RTM, *Journal of Composite Materials*, **26**, 1992, pp. 1100–1133, (doi: 10.1177/002199839202600802).

[9] T. Ahmed, *Fundamentals of Rock Properties*, Reservoir Engineering Handbook, Elsevier, 2019. (doi: 10.1016/B978-0-12-813649-2.00001-3).

[10] W. Han, Composites manufacturing by injection of reactive thermoplastic resin in a fibrous preform, PhD Thesis, IMT Mines Albi - University of Toulouse, 2022.

[11] S.H. Kim, J.W. Jung, M.X. Li, S.W. Choi, W.I. Lee, and C.H. Park, Unsaturated flow behavior in double-scale porous reinforcement for liquid composite molding processes, *Journal of Reinforced Plastics and Composites*, **36**, 2017, pp. 85–97, (doi: 10.1177/0731684416671422).