Bayesian inversion algorithm for estimating local variations in permeability and porosity of reinforcements using experimental data

A novel Regularising Ensemble Kalman filter Algorithm based on the Bayesian paradigm was applied to RTM processes to estimate local porosity and permeability of fibrous reinforcements using measured values of local resin pressure and flow front positions during resin injection. The algorithm allows to detect locations of defects in the preform. It was tested in virtual experiments with two geometries, a two-dimensional rectangular preform and a more complex 3D shape, as well as in laboratory experiments. In both the virtual and laboratory experiments, it was demonstrated that the proposed methodology is able to successfully discover defects and estimate local porosity and permeability with good accuracy. The algorithm also provides confidence intervals for the predictions and estimations of defect probabilities, which are valuable for analysis of the process.


Introduction
Resin Transfer Moulding (RTM) is a cost-effective and versatile process for the manufacture of components from composite materials. In RTM, a preform from dry fibrous reinforcement is placed in a stiff tool. It is then impregnated with a liquid polymer system, typically a thermoset resin, upon application of a flow-driving pressure gradient. Once the reinforcement is fully impregnated, the composite is left in the tool until cure is complete. The component is then removed from the tool and finished. However, it is a well-known problem that fibrous reinforcements may show substantial local variability in fibre orientation and/or fibre volume fraction. This results in local variability in the reinforcement permeability [1][2][3], which affects resin flow in RTM processes [4][5][6], and in variability in mechanical properties of the finished component [7,8].
The resin flow in RTM is determined by the properties of the preform (which are related to fibre arrangement and fibre volume fraction) and of the injected resin, the geometry of the mould and the process parameters. The local properties of a fibre preform can be expressed in terms of two parameters, porosity, ( ), and permeability tensor, ( ), where can be any point within the preform. It is to be noted that ( ) depends on ( ), but also on a geometrical factor related to the fibre arrangement. The resin flow through this domain is characterised by the resin pressure, ( , ), at position and time and by the moving boundary (in other words, flow front), ( ), between resin and air. * Corresponding author. A transient one-phase resin flow through a preform defined by parameters ( ) and ( ) can be modelled by partial differential equations (PDEs) [5,[9][10][11][12] in a domain , which can be one-, two-or threedimensional. The boundary of the domain is = ∪ ∪ , where is the inlet, represents impermeable mould walls, and is the outlet. The domain is initially filled with the preform and air at a pressure 0 . This medium is injected with resin of viscosity through the inlet boundary at a pressure . The resin moves through occupying a time-dependent sub-domain * ( ) ⊆ , which is bounded by the moving boundary ( ), i.e. flow front, and the appropriate parts of . Flow in the fully resin-saturated part of the preform is frequently modelled by combining the continuity equation, with Darcy's law, where ( , ) is the phase-averaged flow velocity and is the resin viscosity. Eqs. (1) and (2) are accompanied by the following boundary and initial conditions ( , ) = 0 , ∈ ( ), > 0, ( , 0) = 0 , ∈ , (7) and the equation for the normal velocity ( , ) in the point on the moving boundary ( ): Here ( ) and ( , ) are the unit outer normals to the corresponding boundaries. The above model describes one-phase flow (only resin), but it can be extended to two-phase (resin and air) flow.
For the purpose of this paper, it is convenient to view a resin flow model as a map from parameters of interest (here, porosity, , and permeability, ) to the pressure, , and the boundary, : In the case of the above one-phase model, the porosity ( ) and permeability ( ) are independent input parameters. Eq. (10) means that given these input parameters the moving boundary problem described by Eqs. (1) to (9) can be solved to find the pair ( , ). In the language of Uncertainty Quantification, Eq. (10) is called a forward map, as it evaluates observables at given values of the parameters. The forward map can be a model expressed as a PDE problem (e.g. the problem given by Eqs. (1) to (9)) or its numerical approximation. In this paper, the models for transient two-phase flow in porous media implemented in Ansys Fluent© [13] are used as forward maps (see details in Section 3).
The material-inherent local variability of porosity, ( ), and permeability, ( ), as well as the presence of defects in the reinforcement structure, such as gaps resulting in race-tracking, have a detrimental effect on the robustness of RTM processes. Different approaches have been developed to characterise variability and race-tracking [14][15][16][17], as well as model their effect on the RTM process [4,14,[18][19][20]. While an effective global permeability of a reinforcement can be estimated using experimental techniques [21,22], the local variability of the permeability is more difficult to detect and estimate. Imaging techniques, such as micro-CT, can be used to scan small sections of the reinforcement or of a finished component. The local fibre arrangement can be determined with high accuracy, and the local permeability can be approximated numerically based on this information [23]. However, this method is not feasible for larger components. Alternatively, local permeability and defects can be estimated using in-process measurements combined with numerical modelling.
Unlike the forward problem outlined above, solving the inverse problem aims at finding parameters (here, porosity, , and permeability, ) based on given observations (here, pressure, , and/or location of the moving boundary, ). To this end, pressure sensors can be used at locations , and sensors to detect the position of the flow front, , so that ( ) is the position of the flow front measured at the th sensor at time . Pressure and flow front position at all sensors are recorded at times . All acquired data of pressure and flow front position are combined in the following variable The forward model given in Eq. (10) maps porosity and permeability onto predictions of pressure and flow front at the sensor locations and specified observation times, which is formally written as The main idea of inversion methods is to minimise the difference between the data obtained using estimated parameters in a mathematical model,  , and the data observed in an experiment. An iterative estimation of the parameters will yield gradually improving approximations of the experimental observations. Various inversion methods have been used for estimation of local permeability using in-process measurements. Comas-Cardona et al. [15] estimated the local permeability of a random chopped glass fibre mat by minimising the error between the predicted and observed resin flow front positions. The initial guess for the permeability was based on optical measurement of the local porosity. A similar approach was used by Caglar et al. [17] to estimate the permeability of preforms in the presence of race-tracking. Experimental flow front observations together with numerical modelling allowed intentionally induced racetracking to be detected. However, the predicted extent of race-tracking, in particular its length, deviated from the experiment because of a lack of data at the beginning of the experiment, where the flow front propagates too fast for tracking. Both studies used a MATLAB nonlinear optimisation solver to minimise the error between predicted and observed flow front positions.
A more efficient inversion method based on an Ensemble Squareroot Filter was used by Matsuzaki et al. [24,25]. It was demonstrated that this method can detect defects in complex 3D and thick structures fibre structures using flow front observations as input data. However, the detection of outlines of the defects was not very accurate. An approach based on Machine Learning, used by Gonzalez and Fernandez-Leon [26], showed that it can detect defects in RTM processes if they are similar to those used for training of a neural network. However, the neural network predictions become unreliable when type, shape, or total number of defects change.
Some studies focus on detecting race-tracking only, as it is one of the most common and severe issues in RTM processes. It is typically assumed that race-tracking affects flow along the entire lengths of reinforcement edges. This assumption implies that only a finite number of combinations of race-tracking regions can exist in a mould. It is possible to detect one of the race-tracking scenarios by simulating all possible scenarios, i.e. all combinations of race-tracking/no race-tracking along the edges, and performing a flow pattern recognition using experimental data [14,27,28]. An alternative assumption is that race-tracking affects only part of an edge. A similar methodology as before can be followed, though the total number of all possible combinations will be higher [29]. This approach can be efficient in detecting race-tracking but, because it relies on parametrised simulations, the estimations can be incorrect in cases where experiments severely deviate from the simulated scenarios.
This paper presents a methodology to predict local porosity and permeability of a preform based on Bayesian inference. In-process data, pressure values and resin arrival times at defined positions, are collected during resin injection as input data for the inversion algorithm. The inversion uses a novel Regularising Ensemble Kalman filter Algorithm (REnKA) [12] and ideas of Bayesian level set methods for geometric inverse problems [30,31]. The methodology is based on the recently developed infinite-dimensional Bayesian theory [32] and algorithms for data assimilation in problems with complex geometries and local non-uniformity [12,30,33].
Unlike in previous studies, the use of Bayesian inference makes it possible to detect defects in local porosity and permeability based on in-process data without making assumptions on the shape and positions of the defects. The presented methodology is validated using virtual experiments as well as lab experiments. The detection of defects can be used to create realistic RTM simulations and to characterise material variability.

Introduction
Bayesian inversion algorithms in the context of RTM are based on measurement of fluid (resin) pressure at several locations in the preform and/or of flow front positions at several moments in time during the process. A Bayesian inversion algorithm starts with an initial guess, a prior distribution in the Bayesian Statistics terminology, of local porosity and permeability. This distribution is characterised by an ensemble of samples of these material properties weighted with their probabilities. The prior distributions are based on available prior knowledge of porosity and permeability, e.g. their target values according to the design and expected level of variability, but otherwise they are random. Parametrisation of a prior distribution, which can be used to capture complex defect distributions, is described in Section 2.2.
Once an initial ensemble is generated, the RTM process is simulated for every sample of porosity and permeability fields to obtain pressure values and flow front positions at specified locations and moments in time using a mathematical/computational model,  . Based on these simulation results, the Bayesian inversion algorithm iteratively evolves the samples, updating the distribution of local porosity and permeability, so that the calculated pressure values and flow front positions become consistent with the measurements. Details of the REnKA algorithm are described in Section 2.3. Using the updated samples from the calculated distribution (the posterior distribution in the Bayesian Statistics terminology), estimates of local porosity and permeability can be computed together with their standard deviations. As described in Section 2.4, knowledge of the standard deviations is important to provide a measure of uncertainty (error) of the obtained estimates.

Parametrisation of local variabilities
The selection of the prior in the Bayesian algorithm is crucial for accurate estimation of physical properties of the material. Choosing priors is particularly challenging in cases where the local properties are determined by complex geometrical features of the component and/or show discontinuities arising from material defects. This challenge is addressed in this work by applying three levels of parametrisation of the physical properties: (i) level-set parametrisation to separate the regions with and without defects, (ii) random field parametrisation to account for material variability, and (iii) parametrisation of length scales to account for defects and variabilities of different sizes. Combined within the modern Bayesian methodology [30], this three-level parametrisation allows to infer the geometry of regions with higher or lower porosity and permeability than the nominal designed values.
The idea behind the level-set parametrisation in a 2D domain is shown in Fig. 1a for the porosity field. An underlying random field, ( ), often called the level-set function, is generated to parametrise porosity and permeability via truncation at some prescribed level, . Fig. 1b shows the truncated level set function which gives a piecewise constant porosity distribution. The region with low porosity (Region 1) corresponds to values of the level-set function greater than , while higher porosity values (Region 2) corresponds to values below . Heterogeneity within each region can be incorporated by assigning additional random fields with different mean values in each of those regions (Fig. 1c). The combined random fields correspond to the second level of parametrisation introduced above. Third-level parametrisation can be used to introduce variable correlation lengths of the random fields. The mathematical description of the parametrisations is given in Appendix A. For simplicity, the functions and scalar parameters comprising the parametrisation can be written in a single function ( ), and the resulting parametrisation of porosity and permeability can be denoted as where the parametrisation map, , is a relation defined for functions rather than only for the values of these functions. The inference problem is posed in terms of the Bayesian calibration of a parameter-to-output map, denoted as ( ), which maps all parameters in defining porosity and permeability ( ( ), ( )) into a variable, . This consists of values of pressure and flow front recorded at a set of specific times during the resin injection as explained below. More specifically, where  is the map described in the Introduction. The aim is to solve the following inverse problem: for given measurements, possibly corrupted by noise, of  defined in Eq. (14), find ( ). These measurements, denoted by  , are assumed to satisfy where is a vector of random noise. It is assumed that follows a Gaussian distribution with zero mean and covariance , which can be expressed as ∼ (0, ). Note that Eq. (15) simply states that, in the absence of modelling errors, the empirical measurements,  , can be obtained from the predictions of the RTM model by accounting for an additive random error in those predictions. Once from Eq. (15) has been found, Eq. (13) can be used to find and . It should be noted that and are assumed to be independent of each other in all examples discussed below. This assumption accounts for the fact that, in reality, there is no unique relation between porosity and permeability. Typically, the permeability is assumed to depend on the porosity and a geometry factor. The geometry factor, which is often assumed to be constant, depends on the cross-sectional shape, arrangement and orientation of the fibres. All are to some extent stochastic, and fibre arrangement and orientation also depend on the inter-fibre angle (affected by shear) and the material porosity (affected by shear and compaction). Anisotropy of the permeability, , can also be implemented by introducing an additional parametrisation of the tensor components. This additional parametrisation does not alter the algorithm proposed here, but it adds more variables which need to be estimated.
The parametrisation of porosity and permeability presented cannot describe all possible cases. However, other parametrisations can be easily combined with REnKA to accommodate more complex scenarios. For example, multiple level-sets [34] can be used to identify several types of defects with substantially different material properties. Since the truth is unknown in a general case, choosing an accurate parametrisation a priori can be a challenging task.

The Bayesian approach
The solution of the inverse problem is not unique. Different may be consistent with the observed data,  , i.e. satisfy Eq. (15). Therefore, a deterministic approach which produces a single estimate of can overlook a possible range of admissible solutions. Here, a probabilistic framework is considered that allows to compute a distribution for and thus facilitates the quantification of uncertainties associated with the estimated parameters. The infinite-dimensional Bayesian framework [32] is adopted, in which is a random function with a specified prior probability distribution, P( ). The prior is based on probabilistic knowledge of the material properties of the preform, before any data are acquired (i.e. before resin injection). For example, the prior may incorporate design values of porosity and permeability as well as information on regions where defects such as race tracking [9] are likely to be present.
In the Bayesian framework, the solution to an inverse problem is the conditional (posterior) distribution of the unknown for given  . The posterior, denoted by P( | ), can be expressed via Bayes' rule as follows [32]: where P( | ) is the probability of the observed measurements,  , given a particular realisation of the unknown, ( ). The term P( ) in Eq. (16) denotes the probability of  . This term is a normalisation constant defined by Taking into account Eq. (15) and the assumption that is a Gaussian, it follows that  | ∼ (( ), ), and thus Eq. (16) becomes Given a prior, Eq. (18) entirely defines the posterior, P( | ), provided that the normalisation constant, , defined in Eq. (19) can be computed. However, due to the nonlinearity of the parameter-to-output map, , which appears in Eq. (19), this normalisation constant can generally not be computed analytically. Hence, the resulting posterior distribution, P( | ), cannot be expressed in a closed form. Sampling methods need to be applied for the approximation of the Bayesian posterior [33]. Here, REnKA is used for this purpose [12,35]. The underlying idea of the method is to compute a sequence of intermediate distributions between the prior and the posterior: Each of these distributions is approximated via an ensemble of particles, { ( ) } =1 (note that ( ) ∼ P ( )). The transition between two consecutive distributions is performed by using Bayes' rule under Gaussian assumptions on the underlying parameter , which, in turn, determines porosity and permeability via the parametrisation introduced in Section 2.2 and Appendix A. Algorithm 1 gives all steps needed to approximately solve Eq. (18) by constructing the sequence in Eq. (20). The number of intermediate distributions, , is computed adaptively [35] in MATLAB©. The parameter in Algorithm 1 controls convergence. The algorithm is designed so that once +1 = 1, the algorithm is converged in the sense that the distribution P +1 ( ) is a Gaussian approximation of the sought posterior P( | ).
It is to be noted that the map  ( , ) used in Algorithm 1 is an RTM process simulation, which is used in a black-box fashion. The total computational cost of REnKA is given by = × RTM simulations, where as defined above is the total number of iterations. As REnKA scales with respect to the number of particles, , its computational execution time can be substantially reduced via the use of parallel computing. The cost of updating the samples of the ensemble is negligible when compared to the cost of RTM simulations.

Post-processing the results of the inversion algorithm
The ensemble of particles obtained via REnKA (Algorithm 1) can be used to compute an approximation of the posterior distributions of porosity and permeability as described above. For the permeability, the posterior mean and standard deviation can be computed as: where the ( ) ( ) are computed from ( ) using Eq. (13).
In analogy to Eq. (21), the corresponding posterior mean, ( ), and standard deviation, ( ), of the porosity can be computed. The dependence on in Eq. (21) indicates that statistical measures of the unknown properties are functions of the position within the preform. Incorporating local variability in those properties is crucial to identify locations of material defects.
In addition, it is possible to evaluate the probability of a material point to be in a defect region: for a fixed , the probability that the level-set function ( ) from the top-level parametrisation is above a threshold value is: This probability gives the confidence level for whether a defect is present or not at the point .
Finally, the posterior ensemble of porosity and permeability, , is used to approximate the posterior distribution of measured process data via where ( ) ∼ (0, ). The ensemble { ( ) } =1 can be used to assess whether the posterior distribution of material properties is consistent with experimental data. More importantly, confidence intervals for the posterior can be computed from this calculated distribution.
Note that the computational cost of setting up the prior and of post-processing the results is negligible in comparison with running REnKA.

Simulation setup
Algorithm 1 was validated in virtual experiments using computer generated input data. RTM processes were simulated in ANSYS Flu-ent©. Through-thickness flow was neglected, which made it possible to use a 2D implementation of the numerical algorithm. It was assumed that the permeability of the reinforcement does not change with time (no compaction or relaxation, no change of properties with saturation) but varies locally. The permeability was assumed to be uniform on each cell of the mesh used in simulations but varied between cells. The local distribution was generated in MATLAB© and imported into ANSYS Fluent© via a user-defined function. The viscosity of the fluid was assumed to be constant at a value of 0.1 Pa s, i.e. the process was assumed to be isothermal with no resin cure present. No-slip boundary conditions were imposed at the walls of the cavity. The settings of the ANSYS Fluent© solver are summarised in Appendix C.
True porosity, † ( ), and true permeability, † ( ), are the prescribed properties of a preform in virtual experiments. Data, , from virtual sensors are recorded during simulations of flow through a preform with these properties, i.e.  =  ( † , † ). At any time, data consisting of fluid pressures and flow front positions are superimposed with Gaussian random noise, , i.e.  =  + . It is assumed that the error covariance matrix, , is diagonal with elements equal to the variance of each of the measurements. The generated data,  , are used as the input for Algorithm 1.

Rectangular 2D component
The aim of this virtual experiment is to identify defects in a flat rectangular preform. A schematic drawing of a rectangular injection tool is shown in Fig. 2(a). Six pressure sensors and seven equally spaced linear flow sensors were placed in the tool. A constant pressure of 0.4⋅10 5 Pa was set at the three inlets and 0 Pa at the outlet. Mesh convergence studies indicated that the mesh consisting of 1346 cells shown in Fig. 2(b) is suitable to obtain accurate results.
It is assumed that the reinforcement contains two circular defects of radius = 0.025 m with centres at coordinates (0.042 m, 0.108 m) and (0.055 m, 0.205 m), respectively. The porosity and permeability values are lower in the defect regions than in the background regions (outside of defect regions). Both porosity and permeability are represented using Gaussian random fields (see details in Appendix A) to account for variability within the defects and the background regions. The true porosity, † ( ), and true permeability, † ( ), for this virtual experiment are shown in Fig. 2(c). It is to be noted that the permeability is assumed to be isotropic in this example. Hence, it can be represented as a scalar. Parameters of the Gaussian random field generation are given in Appendix B.
The true values, † ( ) and † ( ), are used to simulate the RTM process and generate virtual (noise-free) data, . Standard deviations of 475 Pa and 1.4 × 10 −3 m for pressure and flow measurements, respectively, were used to generate the corresponding components of . These values correspond to 0.5% of the maximum pressure and 1.5% of the maximum flow front distance. The (squares of) these values are used in the diagonal of the error covariance matrix, . Virtual measurements are used to infer † ( ) and † ( ) within the Bayesian approach discussed in Section 2.3.
A prior for , i.e. an initial guess of porosity and permeability, is selected to be un-informative, which in particular means that porosity and permeability can take values from a wide range. Selection of the prior does not rely on information about defects (neither on their sizes nor locations). Details of the selected values for the parameters are given in Appendix B. An initial ensemble { ( ) 0 } =1 with = 350 samples, which is generated from the prior , is used to compute the corresponding ensemble of porosity and permeability using Eq. (13). Several samples from the prior ensembles of porosity are shown in Fig. 3. Permeability samples have the same structure as porosity samples, but no functional dependence between permeability and porosity values is imposed. The samples from the priors show substantial variability between them in terms of location and shape of defects and values of porosity and permeability , which is the result of selecting a wide range of parameters of the prior . The mean and standard deviation of both distributions are almost uniform over the entire domain. The probability of the presence of a defect, as defined by Eq. (22), for this prior is close to zero everywhere in the domain.
Reconstructing the true material properties, the REnKA algorithm converged after 7 iterations using the selected prior. The computational time for a single sample was about 90 s on a desktop PC (Intel Xeon E5-1660 v4 @ 3.2 GHz) using Ansys Fluent© 19.1, which brings the total cost of the problem to about 61 CPU-hours. The total time required can be reduced significantly by running samples in parallel (e.g. the total time will be about 36 min if executed on 100 CPUs). Fig. 4 shows several samples from the posterior ensemble. Each of these samples contains defects of similar shape as well as some variability outside of these defects which varies from sample to sample. Fig. 5 shows a comparison between the true porosity and permeability and the inferred mean values, the standard deviation of the porosity and permeability, and the local defect probability as computed by Eq. (22). This comparison shows that, although the selected prior contained no knowledge about mean porosity values or correlation length of defects, the position and shape of the defects are predicted with high accuracy as can be seen from both mean estimates and defect probability plots. It also can be seen that the correlation length in porosity and permeability   field is restored and that it is possible to use these posterior prediction to restore a relationship between porosity and permeability. Another observation is that the variability of the estimates is significant around the defects and near some of the edges, which is a result of uncertainty in predictions.
The accuracy of the inversion can be estimated by comparing the input data with data from the posterior after the algorithm converged. Prior and posterior values of resin pressure at one of the pressure sensors (the bottom-left sensor on Fig. 2(a)) are shown in Fig. 6 along with prior and posterior predictions for flow front positions as measured by the flow sensor on the left in Fig. 2(a). It can be seen that the values from the prior ensemble have wide uncertainty around the experimental results, which shows that the selected prior covers the range of expected results. The resin pressure and flow front position from the posterior distributions have very narrow confidence intervals and have good fit to the results based on the input data.
From the good quality of the inversion (i.e. accurate estimates of pressure and flow front measurements), it is possible to conclude that the algorithm can produce accurate estimations of properties and defects using an un-informative (wide) prior together with sufficiently informative pressure and flow front data. Therefore, a prior has no substantial effect on the posterior estimates in a data-rich setting, provided that the prior reflects possible variety of properties encountered in the physical problem. However, selection of a prior in a less datarich setting becomes more important and would require employing additional techniques such as Bayesian model selection [36].
Apart from the prior selection, the number, type (pressure or flow) and accuracy of sensors also play an important role in the quality of the produced estimates. Numerical experiments conducted in a similar   setting [12] showed that lower number of sensors can be used to detect large defects but the uncertainty of the predictions will be higher than shown in the experiment above. However, other features such as correlation lengths are often not recovered precisely, unless the number of sensors is sufficiently high and the measurement error is relatively low. It also should be noted that flow front sensors can be replaced by pressure sensors with no reduction in the inversion quality [12].
The number of particles, , was chosen on the basis that it provided robust estimates of the unknown properties. Specifically, additional simulations using more particles for the initial ensemble produced similar posterior statistics but did not improve the quality of the inversion. While the optimal choice for the number of particles used within REnKA is still an open problem, existing work confirms that ∈ [200, 400] is a good compromise between computational cost and accuracy [12,35]

3D component
An adapted geometry of a 3D component was used to validate REnKA for a more complex 3D problem. The component, shown in    is set to 10 5 Pa. The pressure at the vent, located at the wide end of the component, is set to 0 Pa.
The preform in this virtual experiment contains three defects, where porosity and permeability are higher than the corresponding properties of the background. It was assumed that the true porosity, † ( ), has values of 0.5 in the background and of 0.7 in the defects, as shown in Fig. 7(right). For simplicity, there is no local variability in porosity and permeability within the background or defect regions. The true permeability, † ( ), was again assumed to be isotropic. It was calculated assuming that the dependence on porosity can be expressed as where the parameters were selected as = 10 −10 m 2 and = −1.4. However, this relationship is assumed only between the true porosity and true permeability and is not used within the inversion algorithm. As gaps may form between the preform and the female part of the RTM tool along the bends in the component, effects similar to racetracking may occur during mould filling. These defects are modelled explicitly as regions with higher porosity and permeability.
Virtual measurements,  , are used to infer † ( ) and † ( ) within the Bayesian approach discussed in Section 2.3. The virtual pressure measurements are superimposed with 1% noise. For the inversion, a 2D three-level parametrisation of porosity and permeability fields is employed as described in Appendix A.4. The parametrisation assumes that porosity and permeability for each region are fixed for each sample, but these values vary between samples. The generated 2D fields are mapped onto the 3D geometry preserving the in-plane distances between points.
An initial ensemble of = 400 particles is used to compute the corresponding ensembles of porosity and permeability. The computational time for a single sample was about 700 s on the desktop PC as specified above, which brings the total cost of the problem to about 855 CPU-hours. The computational time can be reduced significantly by running the analysis in parallel (e.g. the total time will be about 8.5 h if executed on 100 CPUs). Selected samples of the prior ensembles of porosity are shown in Fig. 8. Permeability samples follow same the pattern, but no functional dependence between permeability and porosity values is imposed. This choice of prior does not assume any knowledge of the positions of defects and has a substantial local variability.
Using the described setting, REnKA converged after 11 iterations. The mean is computed from the ensembles of porosity and permeability to identify locations of the defects. It is shown in Fig. 9 along with the true porosity, † ( ), and permeability, † ( ). It can be seen that the inferred porosity and permeability match the true porosity and permeability in the preform in terms of position and shape of the defects. However, the mean values of porosity and permeability in the background and defect regions are not exactly equal to the true porosity   and permeability values. Yet, the difference between these values is less than 5%. The calculated standard deviations do not show large variability anywhere but around the intentionally induced defects. The defect probability, shown in Fig. 10, also shows high values only for the regions where the defects were placed intentionally. This accurate detection of complex defects in the 3D preform is enabled by the large number of sensors employed in this virtual experiment. It is expected that smaller numbers of sensors will yield higher standard deviations of porosity and permeability, and hence a defect probability map with less sharp boundaries between regions with and without defects.

Lab experiments
A rectangular mould with a geometry and a configuration of pressure sensors identical to that used in Section 3.2 (three inlet gates and one outlet) was used in the experiments. The mould consists of a steel bottom, a steel spacer frame of 2 mm thickness, and a transparent PMMA top. Each of the three inlets is connected to a separate pressurised fluid tank fitted with an electro-pneumatic pressure regulator. Six pressure transducers are mounted in the bottom plate of the mould. A digital video camera was used to record the mould filling process and obtain flow front positions. The experimental setup is shown in Fig. 11.
A continuous glass fibre random mat with an areal weight of (259 ± 15) g/m 2 was used in the experiments. Preforms consisted of 7 layers of the reinforcement compacted to a porosity of approximately 0.71. Defects were created by placing circular patches of 4 additional layers of the reinforcement in the middle of the preforms, between the second and third layers, as shown in Fig. 12. Size and positions of the defects are the same as described in Section 3.2. The injected fluid was engine oil with a viscosity of 0.106 Pa s at the test temperature. The injection pressure at all three gates was set to 0.4⋅10 5 Pa.
Measurements from the pressure transducers were collected at a rate of 10 s −1 during the experiments. These measurements were resampled to create a data vector with length of about 200 readings. The video of the mould filling was analysed in MATLAB© to emulate data from seven linear flow front sensors. The collected experimental data,  , were used to infer preform porosity and permeability using the Bayesian approach. The algorithm was initialised with the same initial ensemble as in Section 3.2. Gaussian Process (GP) Regression was used to estimate the measurement error via the MATLAB© toolbox GPStuff [37]. Assuming that sensors are independent of each other, a GP was fitted to the experimental data from each sensor. Using the GP built-in routine for maximum likelihood estimation, the variance of each GP is inferred. These variances are used to construct the diagonal of the measurement error covariance .
REnKA converged after 8 iterations. The posterior ensemble of unknown parameters was used to infer the porosity and permeability distributions, the mean and standard deviation of which are shown in Fig. 13. The defect probability map computed using Eq. (22) is also shown in Fig. 13. The porosity and permeability distributions indicate the presence of defects which were intentionally induced. Additional variability in porosity and permeability (in the background) reflects material-inherent variability of the reinforcement but also possible residual errors in the algorithm. Apart from the two large defects in the porosity and permeability, the restored mean distributions of porosity and permeability show a number of smaller defects. However, the corresponding calculated standard deviations show higher values near these smaller defects, indicating higher uncertainty of results at these locations. The probability for most of these smaller defects is below 0.9, which is lower than the probability of the two large defects. These smaller defects can be attributed to uncertainty of the experiment or local variability as well as to the lack of experimental data, noise in the measurements or modelling errors.
The posterior values of pressure and flow front positions, which show how well the inversion algorithm approximated the input data, show an excellent fit to the data in Fig. 14.

Conclusions
The Bayesian inversion algorithm, REnKA, has been applied in this study to reconstruct local porosity and permeability values using data obtained during the resin injection process in RTM. The algorithm requires a smaller number of samples than a straightforward Monte Carlo algorithm. It can perform computations on each of the samples in parallel, which makes this algorithm faster and more scalable than the inversion methods based on non-linear optimisation.    Fig. 2(a); Right: Posterior predictive equal-tailed 95% confidence intervals for flow front position and experimental data for the two left flow sensors in Fig. 2(a).
The ensemble of samples makes it possible to estimate statistical characteristics of the material parameters, such as mean and standard deviation. Similarly, the probability for the presence of defects can be computed. This is not possible in deterministic inversion algorithms including those based on Machine Learning.
The algorithm also employs a novel three-level parametrisation, which is important in modelling complex random fields. In particular, the parametrisation enabled to describe random fields with two (possibly unconnected) regions with different properties, characterised, e.g. by different mean values and different length scales.
The algorithm was tested using both simulated and experimental in-process data. It was demonstrated that the algorithm can restore the defects with high precision in terms of shape and position, as well as values of porosity and permeability. The variability of the calculated values was shown to be high near the defect interfaces, and, in some cases, near the edges of the domain, where in-process data could not be obtained. It was shown that the algorithm is capable of detecting defects such as race-tracking in 3D components for which the two length scales differ by an order of magnitude.
The algorithm was validated using experimental data obtained from injection experiments with two intentionally induced defects. Both of the defects were detected with a good precision in terms of shape and position. The uncertainty in the calculations was higher than for the virtual experiments owing to the presence of experimental errors as well as possibly limited precision of the algorithm and accuracy of the resin flow model.
This study shows that in-process data can be successfully used to infer local porosity and permeability distributions. The inversion algorithm can be used either for characterisation of material or for defect detection to estimate the manufacturing process and product quality.
Additional studies on finding optimum numbers and positions of sensors for engineering problems are required based on the theory of optimal experimental design. Applicability of the algorithm in less datarich settings also depends on appropriate selection of the prior which can be constructed using more complex parametrisation than presented here or generated using a narrow range of parameters. Speeding up the algorithm is another area for future work. It can be made faster by having a quicker RTM solver (e.g. by making use of Model Order Reduction techniques) and better selection of priors (e.g. by quickly informing them by collected data). Further testing of the algorithm in virtual and lab settings is also of substantial interest with the ultimate objective to use it in an industrial environment or in active control systems.

Declaration of competing interest
The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.

A.2. Second level parametrisation: random fields
The second level of parametrisation defines variability within the two regions defined by the level-set parametrisation. For this purpose, each ( ) is parametrised in terms of random fields (RFs) [38,39]. For simplicity, the 2D case is considered here, but the formulation can be extended to 3D RFs.
For a fixed , the function ( ) is found by solving the following fractional stochastic PDE: where ( ) is a generalised random function (e.g., Gaussian white noise), is the identity operator, and is the gamma function. The value of log defines the mean value of , is a parameter that controls the smoothness (or roughness) of ( ), is an amplitude scale, and ,1 and ,2 are positive scalars which describe intrinsic lengthscales in the two orthogonal directions. These length-scales characterise correlation lengths of ( ). For simplicity, it is assumed here that the value of 1 is specified along the horizontal direction. Anisotropic variability along a specific direction can also be incorporated via suitable modifications of Eq. (A.6) [38].
Computing via solving Eq. (A.6) requires the following set of inputs parameters: where the subindex is included in the parameters to indicate that each ( ) has its own parametrisation. Therefore, Eq. (A.6) can be written in operator form as ↦  2 ( ) = .
(  [39], and a covariance operator given by an auto correlation function of the Matern class [38,40,41]: where is the modified Bessel function of the second kind of order , and While there are other approaches (e.g. Karhunen-Loeve expansion) to characterise Gaussian RFs, the advantages of the formulation in Eq. (A.6) is that it allows to consider the case of spatially variable length scales, which is essential to identify defects within complex geometries.

A.3. Third level parametrisation: non-constant length scales
The parametrisation in Eq. (A.9) assumes constant length scales, ,1 and ,2 , for the level-set function, , which parametrises the interface between two regions with different physical properties. However, in order to infer/capture defects with various geometries, the function ( ) in Eqs. (A.2) to (A.3) needs to have spatially varying length-scale. This is crucial in the case where the defect regions display significant changes throughout the preform. For example, regions of race tracking are typically long and thin [9], while defects due to manufacturing errors in laying up preforms are of arbitrary shape which can have the same characteristic lengths in both directions. Hence, in ( ) =  2 ( ) =  2 ( , , , ,1 , ,2 , ), (A.11) ,1 and ,2 are no longer assumed to be positive scalars but exponentials of RFs: (A.14)

A.5. Numerical implementation
The main challenge for the numerical implementation of the parametrisation maps, ( ), from the previous subsections lies in computing RFs by solving Eq. (A.6). This challenge is addressed by implementing the techniques proposed by Lindgren et al. [39] in MATLAB©. A regular grid of 100 × 100 cells is used to discretise these RFs which, in turn, define the porosity and permeability ( ( ) , ( ) ) = ( ( ) ). These properties are then interpolated on the mesh used in ANSYS Fluent© (see Section 3.2).
Eq. (A.6), which needs to be solved for each particle at every iteration of the REnKA algorithm, is a relatively simple linear equation that can be solved in a few seconds with standard computer resources. Hence the cost is negligible compared to the cost of running the flow simulations. As stated earlier, the equation is solved using an in-house MATLAB© script, but existing PDE solvers can be used as well. For the sake of simplicity, the examples presented here used 2D random fields or 2D random fields mapped into a 3D mesh (Section 3.3). However, REnKA can be applied to directly infer 3D random fields using the corresponding 3D version of Eq. (A.6) and more sophisticated PDE solvers (e.g. FEniCS [42]) which can work with arbitrary geometries.

A.6. Probability to find a defect
The probability of a point to be in the defective region is characterised via Eq. (A.4) (recall that is fixed a priori). For a fixed ∈ * , the probability for the value of the level-set function ( ) to be above the threshold that defines the defective region is: This probability can be expressed as where is the probability density of the random variable ( ) for the fixed . From the prior/posterior ensemble for , the ensemble ( ) ( ) can be constructed. This ensemble is a sample distribution of ( ). Therefore, using Eq. (A.4) again: While the proposed Bayesian inversion algorithm can be used to infer all the parameters encapsulated in , for simplicity, some of these parameter are kept constant. The amplitude scales, , were selected according to the changes (i.e., an order of magnitude) expected in the RFs for porosity and permeability within each region. Similarly, the smoothness parameters are chosen so that the samples of these Table 1 Fixed parameters for 2D problems. Unit of 1 and 2 is log(m 2 ). fields display realistic degrees of variability as shown in Fig. 3. For example, it is expected that the level-set function, which defines the interface between two region, is sufficiently smooth while the porosity and permeability fields exhibit more irregularities and variations. This is ensured by choosing = 3.0 and = 1.0, respectively. Values of the constant parameters are given in Table 1.
Since parameters and are constants and do not need to be inferred, the focus is on the parameters , ,1 , ,2 , and ( ). It is assumed that, under the prior, the unknown parameters are independent random variables/functions (for all ∈ ), and so the joint prior can be written as where P( ), P( ,1 ), P( ,2 ) and P( ) are the corresponding prior distributions. Additional prior knowledge of these parameters, such as correlations between them, can be also incorporated within the proposed framework. For all ∈ , are chosen to be Gaussian white noise (informally denoted as P( ) = (0, )). As discussed earlier, this choice implies that all RFs [39], { } ∈ , employed for the parametrisation are, under the prior, Gaussian with Matern covariance. The mean of each field is log . The prior distributions for the are selected as follows This selection of distributions for the 1 , 2 , 1 and 2 ensures that the average values for the RFs are within a reasonable range for the corresponding regions of low/high porosity and permeability.
The parameter is set to be equal to 1, which amounts to selecting a zero mean field (i.e., log = 0) throughout the inversion process. This is possible because the level-set function is an artefact that is employed to define two regions in the preform. For this experiment, the threshold in Eq. (A.1) is chosen to be equal to 0.5.
For the distribution of the length scales, the following priors were used

B.2. 3D component
For the experiment of Section 3.3, the simplified 3-level parametrisation described in Appendix A.4 was used. The parametrisation in Eq. (A.20) was simplified by fixing some of these parameters as displayed in Table 2, and the inversion was conducted on All elements in are assumed to be independent, so that the prior can be written as