RTMsim - A Julia module for filling simulations in Resin Transfer Moulding with the Finite Area Method

RTMsim is a robust, easy-to-use and simple-to-extend software tool, implemented in Julia (Bezanson et al., 2017) for resin transfer moulding (RTM) filling simulations. A shell mesh, injection pressure, resin viscosity and the parameters describing the preform are required input. The software was validated with results from literature (Advani et al., 1994) and compared with results from well-established RTM filling simulation tools (Gion Barandun et al., 2012),(ANSYS, 2021).

: Schematic of a thin-shell RTM geometry for manufacturing a flat plate with linear (a) and radial (b) flow.
In the past, numerous models have been implemented in different software packages to perform filling simulations for RTM. Those can be divided into three groups: • General purpose CFD software packages, such as ANSYS Fluent or OpenFOAM • Commercially available software packages which are tailored for the simulation of the RTM process, such as PAM-RTM, RTM-Worx or LIMS • Easy-to-use simulation tools such as myRTM All packages describe the flow on a macroscopic level. The first group models the flow through the porous cavity using volume-averaged Navier-Stokes equations (Grössing et al., 2015). The second group makes use of some assumptions on flow physics and boundary conditions and at first solves a Laplace equation for the pressure inside the filling fluid (resin) and in a second step calculates the flow velocity field to propagate the flow front, see e.g. (RTM-Worx Documentation, Overview, User Manual and Recent Publications, 2021) or (Publications about LIMS, 2021). In (Grössing et al., 2015) it has been shown that the first and second group render very similar results. myRTM from the third group is easy-to-use and can predict the filling pattern properly but neither predict the filling time correctly nor consider orthotropic preform permeability, see (Gion Barandun et al., 2012). Solving conservation laws for fluid flow as in the first group requires a volume mesh of the cavity and consequently the solution is often time-consuming. The second and third group can be solved on a shell mesh where the thickness of the cavity is a property of the cell (similar to porosity and permeability) and slip boundary conditions at the top and bottom walls of the cavity are assumed.
Based on the analysis of the existing software tools for RTM filling simulations the following functional requirements for a new software tool were derived: • The simulation model shall give correct results for filling pattern and filling time.
• The simulation tool takes only composite-manufacturing related inputs and the simulation shall be robust and independent of under-relaxation coefficients or time step.
• The simulation tool takes a shell model of the geometry as input and the locationdependent properties are assigned directly on the shell elements.
• New modules can be implemented by either adding equations of the same type or modifying existing equations.
RTMsim ia a new software tool for RTM filling simulations which fulfills these requirements: Several test cases were used for successfully validating the implemented model. The simulations run robustly and independent of numerics related input such as under-relaxation coefficients or time step. The porous cavity is fully described by a mesh file with triangular cells on the part's mid-surface and cell set definitions (for specifying the location of pressure injection ports and vents and for specifying regions with different properties than the main preform by assigning different thickness, porosity and permeability values). Additional equations (e.g. for modeling the degree-of-cure) can be added with equations of the same type and modifications of existing equations (e.g. for variable cavity thickness as needed for vacuum assisted resin infusion simulations (Strong, 2008)) are possible with cell properties which vary over simulation time.

Modeling
From a flow physics point of view the filling process is described by an incompressible resin which is injected under pressure into a cavity which initially is either evacuated or filled with air.
Since the simulation tool shall run robust and independent of numerics-related input, a one phase compressible model is chosen. The governing equations shown in (Grössing et al., 2015) are modified. Conservation laws for a compressible continuity equation, momentum equations in a local cell coordinate system, an adiabatic law as equation of state and a volume-of-fluid equation to track the resin flow front must be solved. The local cell coordinate system is used to describe the flow in the coordinate system of the orthotropic in-plane permeability (permeability 1 in first principal direction and permeability 2 = ⋅ 1 in second principal direction).
The physical quantities in the governing equations are functions of space x and time where x ∈ Ω ⊂ ℝ 3 and Ω is the fluid body also called the fluid cavity. The governing equations are • is the mass density of the fluid, • u is the superficial velocity which is related to the physical velocity via the porosity given by the stationary fibrous media, • uu is the dyadic product, • is the pressure, • − K −1 u is the pressure loss from flow through a porous medium according to Darcy's law with dynamic viscosity and permeability tensor K, • is the continuous fraction function with 0 ≤ ≤ 1, = 0 if empty and = 1 if completely filled, • is the adiabatic index, and • = / is a constant determined from pressure and mass density reference values.
Since the filling is a result of the pressure difference between injection pressure and initial cavity pressure, gauge pressures are used for the simulation. The initial pressure in the cavity is set equal to zero and the injection pressure is set equal to the pressure difference. This normalization results in a filling which is independent of the pressure level and only depends on the pressure difference.
In order to numerically solve the flow model the computational domain (time and space) is discretized. The temporal domain (i.e. the simulation time) is split into a finite number of time steps where values for the physical quantities on the spatial domain are calculated in a time-marching manner. The spatial domain (i.e. the flow volume) is split into wedge-type cells.
Since the walls are assumed to be slip boundaries only one cell is used through the thickness and the cells are bounded by three control surfaces only. Therefore, the spatial domain is defined on the part's mid-surface and the local thickness of the flow volume is a property of the cell. The mid-surface model can be curved and cells can have edges where more than two cells are connected to each other such as for handling T-junctions. The mid-surface model is divided into a finite number of triangular control areas. The control areas cover the spatial domain completely without overlapping. The cell-centered method where the nodes are placed at the centroid of the control volume is used.
Since a shell mesh is used in the simulation tool, the conservation laws must be solved on a shell mesh using a generalization of the so-called finite area method (Tukovic & Jasak, 2012). The governing equations are discretized for every cell of the shell mesh. Consequently, the velocity vector has only two non-zero components and described in a local coordinate system. The local coordinate system is the result of a projection of a reference vector. The reference vector is chosen such that the local -and -directions correspond with the first and second principal permeability directions. For the evaluation of the numerical flux functions and the pressure gradient, the neighbouring triangular cells are rotated about the common edge to lie in the plane of the considered cell and a coordinate transformation for the velocity components must be performed for every neighbouring cell. Figure 2 shows two neighbouring cells before and after rotation and with their local coordinate systems.  Figure 3 shows the shell mesh which is used for two out of three test cases, which are reported in the subsequent section.

Validation and verification
Different test cases are available, successfully validating the Julia implementation of the RTM filling model. Test cases 1 and 2 are so-called numerical permeameter experiments, see e.g. (Fauster et al., 2019). The third test case is a filling study of a typical composite part which is manufactured using RTM, see (G. Barandun et al., 2013). The testcases are: 1. Validation of the software for radial flow with isotropic in-plane permeablity: The simulated flow front position after 200 s is compared with the calculated flow front position from literature (Advani et al., 1994). The preform size is 0.6 × 0.6 × 0.003 m 3 with a central injection gate with 0.02 m diameter. The preform has porosity 0.70 and isotropic permeability 3.0 ⋅ 10 −10 m 2 . Resin with dynamic viscosity 0.06 Pas is injected with 35000 Pa pressure.
2. Verification of the software for radial flow with tilted orthotropic in-plane permeablity: The simulated tilted elliptical flow front is analysed and the calculated orthotropic permeablity is compared with the input (permeability 1 and 2 , angle of the first principal direction to the horizontal). The preform size is 0.6 × 0.6 × 0.003 m 3 with a central injection gate with 0.02 m diameter. The preform has porosity 0.70 and permeability 3.0 ⋅ 10 −10 m 2 in first principal direction 30 ∘ to the horizontal and permeability 1.5 ⋅ 10 −10 m 2 in second principal direction. Resin with dynamic viscosity 0.06 Pas is injected with 35000 Pa pressure.
3. The filling of a complex annulus filler-like part is simulated with different software tools. The bounding box of the part is approximately 0.4 × 0.15 × 0.1 m 3 . The flow front positions at different time instances simulated with Ansys Fluent and RTMsim are compared in case 3.1. The filling patterns for same filling degrees simulated with myRTM and RTMsim are compared in case 3.2. Figure 4, Figure 5 and Figure 6 show the binary filling fraction for cases 1, 2 and 3. For every cell the fraction function is either 0 or 1, depending on the filling fraction being ≤ or > than a threshold value. The threshold for the volume-of-fluid method is determined from the permeameter simulation with isotropic permeabilities (case 1). The circular and elliptical flow fronts for cases 1 and 2 are calculated from the simulated filling fraction contour plots with an image processing algorithm (Fauster et al., 2019). Simulation input for the porosity was 0.70 and for the isotropic permeability was 3.0 ⋅ 10 −10 m 2 . The flow front position is similar but the shape is no smooth circle for the coarse mesh.
For case 1 also a comparison with results from an analytical formula (Advani et al., 1994) is performed. The following table shows the flow front positions as radial distance from the center at different time instances calculated with an analytical formula and calculated from the simulated flow front with an image processing algorithm for a coarse and a fine mesh.
[ The values for the fine mesh agree well with the results from the analytical formula. The values for the coarse mesh show an error of ≤ 15% which decreases significantly in the course of time.  : Results for test case 3. The results for the new simulation tool with the fine mesh are in good agreement with the ANSYS Fluent results for both filling pattern and filling time. myRTM shows similar filling patterns for same filling degrees as the new simulation tool with the coarse mesh. This comparison shows that a coarse mesh can be used for filling simulations where only the flow front propagation is investigated. If the filling time is important a mesh refinement study must be performed.
A follow-up paper will discuss the physical model in detail and show the results of a validation with experimental data from a radial flow experiment with two patches with different in-plane permeability and porosity levels.

Outlook
The existing software packages for RTM mould filling simulations perform well for most scientific and engineering applications, but there are applications where a new adaptable software package is required. For example, a transient pressure drop in the cavity after closing a valve in the hose between pressure tank and test rig can be studied with RTMsim after some extensions. Furthermore, parts of the RTMsim package can be integrated in separate codes. For example, a customized lightweight version of the solver was already used as training environment for Reinforcement learning (RL) algorithms to keep the flow front straight in a rectangular flow with reinforcement patches and three inlet ports (Stieber et al., 2023).
The source code is prepared for the following extensions: • Import of mesh file in different format.
• Additional functionalities, e.g. adding temperature and degree-of-cure equations with variable resin viscosity or for VARI with variable porosity and permeability.
• Different methods for numerical differentiation and for the calculation of the numerical flux functions.

Software download
The current release can be downloaded from https://github.com/obertscheiderfhwn/RTMsim/releases/. Installation instruction and use cases are described in the included README.md file.