Investigating microstructure-property relationships of nonwovens by model-based virtual materials testing

Quantifying the relationship between geometric descriptors of microstructure and effective properties like permeability is essential for understanding and improving the behavior of porous materials. In this paper, we employ a previously developed stochastic model to investigate microstructure-property relationships of nonwo-vens. First, we show the capability of the model to generate a wide variety of realistic nonwovens by varying the model parameters. By computing various geometric descriptors, we investigate the relationship between model parameters and microstructure morphology and, in this way, assess the range of structures which may be described by our model. In a second step, we perform virtual materials testing based on the simulation of a wide range of nonwovens. For these 3D structures, we compute geometric descriptors and perform numerical simulations to obtain values for permeability as an effective material property. We then examine and quantify the relationship between microstructure morphology and permeability by fitting parametric regression formulas to the obtained data set, including but not limited to formulas from literature. We show that for structures which are captured by our model, predictive power may be improved by allowing for slightly more complex formulas.


Introduction
Fiber-based materials, especially nonwovens, play an important role in a wide variety of applications like gas-diffusion layers in fuel cell technology [1], filtration [2], printing paper [3], and hygiene products [4].In different applications, various material properties are crucial for the performance of the final product and thus, must be optimized during design and production of the fiber mats.Aside from properties of the fiber material, the geometric structure of the fiber systems is directly linked to many effective properties like wettability or diffusivity.Using modern tomographic imaging techniques, the 3D microstructure of existing nonwovens can be investigated and used as an input for numerical simulations to determine effective properties.However, while tomographic imaging and numerical simulations alone may provide some insight into advanced properties which may not be accessible to experimental investigation, they are still limited in the range of feasible structures: Changing the properties of the nonwovens necessitates changes to the production process and is thus costly and time-consuming.
By instead using a stochastic model for the 3D microstructure of nonwovens, a huge variety of virtual but realistic structures can be easily generated on the computer and still be used for numerical simulations.By this so-called virtual materials testing, relationships between the 3D geometry and effective properties can be investigated at a low cost and the development of new materials with improved properties can be accelerated [5][6][7].Moreover, quantitative functional relationships between microstructure and effective properties can be obtained from this kind of data, as has already been shown in the literature.For example, in [8] such relationships were determined on structures derived from artificial models with no direct reference to measured image data for really existing materials.On the other hand, for fiber-based materials, functional relationships based on extensive experimental studies have been derived in [9].However, the latter investigations focused on the solid volume fraction and fiber diameter as the main contributing factors.When performing investigations based on simulated data of fiber systems, a larger number of geometric descriptors can be calculated and a more comprehensive view of the relationship between microstructure and effective properties can be obtained.
The model-based approach described above is primarily motivated by the shortcomings of real-world methods.In order to investigate relationships between manufacturing, microstructure, and effective properties of porous materials through physical experimentation, it is necessary to produce a large amount of samples.This step alone would already require immense time and financial resources.Furthermore, manufacturing processes are often rigid and varying individual production parameters is not always trivial.While it is possible in theory to produce the necessary amount of samples, it would not be feasible.Even if the necessary amount of samples was available, determining all properties of interest through experimentation could prove problematic for many reasons, e.g., if the resolution of a sensor is not high enough or if a material sample is not stable enough to perform the experiment.
Based on the flexible stochastic model for the 3D structure of nonwovens which was developed in [10], we present a simulation study to investigate the relationship between geometric descriptors and effective properties of nonwovens, where we focus on the permeability as the effective property of interest.We lay the foundation for our study by drawing a large number of virtual yet realistic microstructures from the model proposed in [10].By considering a wide range of different parameter sets for the model, we manage to obtain realistic structures with significantly varying properties.On these structures, we compute geometric descriptors and, by numerical simulations, permeability.
First, we investigate the relationship between model parameters and geometric descriptors by systematically varying the values of selected parameters.Not surprisingly, varying parameters which directly control the porosity of the resulting structures leads to a huge variation in a wide range of geometric descriptors.Varying other model parameters which control the shape and arrangement of fibers, we are able to adjust some selected geometric descriptors.Joint variation of multiple model parameters leads to drastic differences in the resulting structures which we investigate in a second step.
The major part of this paper is dedicated to investigating the relationship between geometric descriptors and permeability as illustrated in Figure 1.In particular, we fit parametric regression formulas for predicting the permeability from geometric descriptors.We discuss the performance of various approaches, partly taken from the literature and partly developed for this study.Finally, we discuss the results and outline possibilities for further research.
Fig. 1 Overview of the framework for virtual materials testing.Virtual but realistic structures are simulated using the presented stochastic model.Then, geometric descriptors and permeability are computed for these structures.Functional relationships are established to predict permeability from geometric descriptors.

Methods
In the following, we will outline the methods applied throughout the present paper.In particular, we will give a short summary of the model proposed in [10] and discuss how we choose model parameters for virtual materials testing.Moreover, we describe the geometric descriptors which are used for the investigation of microstructure-property relationships and show how permeability is numerically computed on the simulated structures.Finally, we discuss how parametric regression formulas are fitted for the prediction of permeability from geometric descriptors.

Stochastic Model for Nonwoven Fiber Systems
For the present study, we will employ a spatial stochastic model developed in previous works [10], which was designed to represent nonwoven fiber systems.It is based on assuming fibers with a circular cross section of constant diameter d > 0 and by representing the centerlines of single fibers by polygonal tracks with constant segment length.This is, the centerline of a fiber is considered as a sequence of (random) points . . ., P −1 , P 0 , P 1 , P 2 , . . .∈ R 3 with |P n+1 − P n | = c for all n ∈ Z and some constant c > 0. Furthermore, it is assumed that {P n , n ∈ Z} forms a stationary, reversible thirdorder Markov chain [11].The transition function, which determines the distribution of P n conditional on P n−1 , P n−2 , P n−3 is constructed by modeling the joint distribution of the random vector (P 1 , . . ., P 4 ).Due to stationarity, this fully determines the transition function of the given third-order Markov chain.
For modeling the joint distribution of (P 1 , . . ., P 4 ), the z-coordinates and the (x, y)-coordinates of P i = (X i , Y i , Z i ) for i = 1, . . ., 4 are independently modeled.More precisely, Z 4 is assumed to be independent of Z 1 conditional on Z 2 , Z 3 , i.e., {Z n , n ∈ Z} forms a second-order Markov chain and, for simplicity of notation, the joint distribution of (Z 1 , Z 2 , Z 3 ) is considered.This, in turn, is equivalent to modeling the distribution of (Z 1 , Z 2 − Z 1 , Z 3 − Z 2 ).In the model proposed in [10], a copula approach is used to model this distribution.It consists of modeling the marginal distributions of Z 1 and Z 2 − Z 1 (which has the same distribution as Z 3 − Z 2 ) by generalized normal distributions.Then, in a second step, the correlation structure is modeled by a pair-copula approach, using Clayton copulas for the joint distribution of (Z 3 − Z 2 , Z 1 ) and the joint distribution of (Z 2 − Z 1 , Z 1 ), and a Student's t copula for the joint distribution of ( For modeling the sequence of random vectors (X 1 , Y 1 ), . . ., (X 4 , Y 4 ), the angle A between (X 3 , Y 3 ) − (X 2 , Y 2 ) and (X 2 , Y 2 ) − (X 1 , Y 1 ) and the angle B between (X 4 , Y 4 ) − (X 3 , Y 3 ) and (X 3 , Y 3 ) − (X 2 , Y 2 ) are considered.They determine the values of (X 1 , Y 1 ), . . ., (X 4 , Y 4 ) up to rigid transformations and the joint distribution of (A, B) is sufficient to describe the transition function of (X 1 , Y 1 ), . . ., (X 4 , Y 4 ).The joint distribution of (A, B) is again modeled by a copula approach where the (identical) marginal distributions of A and B are modeled by a generalized normal distribution and the correlation structure of (A, B) is modeled using a Student's t copula, see [10].Note that for example, the curl of fibers is largely determined by the distribution and correlation of the angles A and B. Higher values of the scale parameter α A generally correspond to stronger curl, whereas higher values of ρ (A,B) lead to less changes in direction.
Table 1 gives an overview of the 14 parameters of the model which will be varied for the present study.They will be denoted by P = (p 1 , . . ., p 14 ) in the following.Additional parameters are the extent of the structure in x-and y-direction which are fixed to 5 mm, and the length of individual segments of the polygonal tracks, fixed to 50 µm.

Simulated Structures
By systematic variation of the parameters P = (p 1 , . . ., p 14 ) of the model described in Section 2.1, we obtain a wide range of realistic, yet novel, nonwoven structures.The range of chosen parameters is based on two sets of model parameters P 1 = (p 1 1 , . . ., p 1 14 ), P 2 = (p 2 1 , . . ., p 2 14 ) which were obtained in [10] by fitting the model to measured data of two different nonwovens, where some further bounds for the model parameters must be taken into account, see Table 1.
parameter of Student's t copula modeling (A, B) 1 These are modeled via a generalized normal distribution with location parameter fixed to 0.
2 Z 1 is modeled via a generalized normal distribution.The location and scale parameters are chosen such that P(0 < Z 1 < zmax) = 0.9998 and E(Z 1 ) = zmax/2.
3 For Student's t copula, a higher value of parameter ρ corresponds to a stronger correlation. 4For the Clayton copula, a higher value of parameter α corresponds to a weaker correlation.
For the different aspects of the current study, we consider three different scenarios creating the following sets of structures: Scenario I The first set is obtained by keeping all parameters fixed with the exception of one single model parameter which is systematically varied.Each of the 14 parameters was varied in 10 steps between specifically chosen minimum and maximum values, see Table 1.By this method, we created two subsets of structures, for which all (but the varied) parameters were chosen to equal P 1 or P 2 , respectively.In total, this resulted in 280 = 14 • 10 • 2 cases, for each of which 3 structures were simulated, resulting in 840 structures.
Scenario II The second set of structures was created by randomly varying all model parameters at the same time.This means that, except for the fiber diameter which was kept fixed, each entry p i of the parameter vector P = (p 1 , . . ., p 14 ) was independently drawn from a specific (truncated) normal distribution which was chosen such that the mean value was given by µ = (p 1 i + p 2 i )/2, and the standard deviation by σ = |p 1 i − p 2 i |/Φ 0.75 , where Φ 0.75 denotes the 0.75-quantile of the standard normal distribution.Thus, the value of p i belongs to the interval (p 1 i , p 2 i ) with probability 0.5.The normal distributions were truncated to observe the parameter bounds shown in Table 1.For each of 500 parameter vectors chosen in this way, 3 structures were simulated, resulting in 1500 structures.
Scenario III The third set of structures is a small dataset showcasing the effect of fiber diameter d and fiber density D L on permeability.Note that the porosity ε introduced in Section 2.3 below can be approximately expressed by for some constant a > 0. Thus different combinations of fiber density and fiber diameter can lead to the same value of porosity but hugely different values of permeability.
To investigate this behavior, all model parameters were kept fixed, except for the fiber diameter and fiber density.Seven different target values for porosity ε between 0.9 and 1 were chosen and for each target value, 5 different combinations of fiber diameter and fiber density leading to the specific value of porosity, were considered.For each combination, 3 structures were simulated, resulting in 105 structures.

Geometric Descriptors
In the following, we introduce the geometric descriptors on which further analysis of microstructure-property relationships is based.Most descriptors are similar to those used in [12], where a more detailed discussion can be found.Recall that for the computation of the descriptors we may use the representation of a structure as set of polygonal tracks as well as voxel image.
Porosity Probably the most fundamental geometric descriptor when it comes to any kind of porous material is the porosity ε ∈ [0, 1], or the solid volume fraction 1−ε.We estimate porosity from voxel images by simply counting the number of voxels belonging to the void space and dividing by the total number of voxels.
Tortuosity The notion of tortuosity is often considered when investigating transport phenomena.While there exists a wide variety of different specifications of tortuosity [13,14], we only consider the geodesic tortuosity τ ≥ 1 which is based purely on geometric information and does not involve numerical simulations.It is related to the windedness of shortest transportation paths.We compute the tortuosity by considering the voxelized structures and defining a transport direction, in our case, from z = 0 to z = z max .For each target voxel on the plane through z = z max , we compute the shortest path from the opposite side to that voxel through the pore space using Dijkstra's algorithm [15] on the voxel grid.This voxel-based tortuosity is then computed as the ratio of the length of the shortest path and the thickness z max of the structure.For the present study, the mean geodesic tortuosity µ(τ ) and its standard deviation σ(τ ) are considered which are estimated from the tortuosity values for all target voxels from a given structure.
Constrictivity A further directly transport-related property is the constrictivity β ∈ [0, 1] which measures bottleneck effects along transportation pathways.Again, we define the transport direction from z = 0 to z = z max .Then, the constrictivity β is defined as the ratio β = rmin rmax of two characteristic pore radii r min and r max , where r max is the maximum radius such that 50 % of pore space can be covered in (overlapping) spheres of that radius.These spheres are not allowed to penetrate into the solid phase, i.e., the fibers.Similarly, r min is the maximum radius of spheres which can intrude from the plane at z = 0 into the structure and fill 50 % of pore space.This concept may be viewed as a digital version of mercury intrusion porosimetry.
Specific surface area The specific surface area S is the (mean) area of the solidpore interface per unit volume.It is estimated by computing the total surface area within a given sampling window and dividing it by tthe window's volume.In the case of fiber systems, the specific surface area is directly related to the porosity, the fiber diameter and the curl and overlap of fibers.
Chord length distribution Chords are straight lines in a given direction fully contained in the pore space which start and end at the pore-solid interface or the boundary of the sampling window.The distribution of the length C of such segments is called chord length distribution.In this work, we consider chords in z-direction and estimate the chord length distribution by collecting all possible chords within the voxelized structure.Note that other definitions use other notions of "random chord" and essentially lead to a length-weighted distribution of chord lengths.In the following, we consider the mean chord length µ(C).
Mean spherical contact distance The average distance from a randomly selected point within the pore space to the nearest point within the solid phase (i.e., to the nearest fiber) is called mean spherical contact distance, In the following, it will be denoted by µ(H).We estimate it by randomly choosing a large number of points within the pore space and averaging over their distances to the nearest fiber, computed directly on the polygonal track data.
The geometric descriptors stated above are determined for all simulated structures of Scenarios I, II and III.The obtained results are stored in a database, which contains the values of model parameters used for the simulation of a given structure, along with the corresponding values computed for geometric descriptors and permeability, where the computation of permeability is explained in the next section.

Computation of Permeability
Mathematically, the permeability κ is a tensor and in general depends on the direction of the flow.For nonwoven, fibers tend to be oriented strongly anisotropically nearly in a plane, and only the through direction perpendicular to that plane is measured.By our convention for modelled structures, the through-direction is the z-direction.With this in mind, a scalar permeability κ is computed for the through-direction (or z-direction) based on appropriately rearranging Darcy's law Here, Q is the flow rate, A the area perpendicular to the flow direction over which the flow rate occurs, i.e., in the x-y plane, ū is the (macroscopic) flow velocity, ∆P is the pressure drop, L is the thickness of the nonwoven, i.e., length in the z-direction, and µ is the dynamic viscosity of the fluid, that contains the kinematic viscosity of the fluid and its density.The permeability κ identifies the proportionality between the pressure drop and flow velocity.
Darcy's crucial observation was that the permeability is determined by the geometry of the pore space alone and thus is a constant for porous materials.The viscosity µ and media thickness L are fluid and material parameters, respectively, that must be fixed and known in order to compute κ.Darcy's Law, see Eq. ( 2), is valid under some assumptions on the representativeness of the material and parameters of the flow, i.e., that the fluid is incompressible, Newtonian and flowing rather slowly, corresponding to a moderate pressure drop.In real experiments, u and ∆P are measured, and the proportionality only holds for small pressure drops and velocities, but independently of the fluid viscosity µ and experimental value ∆P .For faster flows and higher pressure drops, there exists a generalization of Darcy's law called Forchheimer's law.The latter uses two material constants which can also be determined by computer simulations, but is not considered here.
In the computational experiments, the permeability for the nonwoven models is determined using the FlowDict module of the GeoDict software [16,17].For a given ∆P , the approach is to solve a large linear system of equations resulting from a discretization of the Stokes equations, see Eqs, (3) -( 5), on the microstructure where the pressure p and the x-, y-, and z-components of the (local, interstitial) velocity vector u at each pore space voxel of the binary image are the unknowns.Then the local velocities get averaged to find the macroscopic or superficial ū [18].For the sake of concreteness the simulations assume that the pressure drop is 0.02 Pa and that the dynamic viscosity µ = 1.834 • 10 −5 i.e., that the pore space is filled with air at 20 • C. Note some symbols conflict with other definitions in this paper, but are chosen in this section to be consistent with existing literature.
Each simulated structure is initially given as a set of polygonal tracks (fibers) within a fixed bounding cuboid W = [0, x max , y max , z max ] and equipped with a fiber diameter d.However, for the numerical flow simulations, a conversion into (binary) voxel images is done by the software.The voxel length h must be chosen small enough to resolve the fiber diameter, and so that the results of the flow simulation do not depend on it as would be the case for too large voxel sizes.For highly porous fibrous materials, as in the case under consideration here, voxel length h = d/10 is sufficient.These binary images also form the base for the computation of some of the geometric descriptors.
Finally, to compute the permeability in the z-direction we add an inlet of the length l in and an outlet of the length l out to the domain to allow for the flow to be uniform far in front and behind the nonwoven, and thus to be able to use periodic boundary conditions for the velocity in all three spatial directions as well as periodic boundary conditions for the pressure.
Under the assumption of slow, viscous, steady state incompressible flow through the nonwoven, the Reynolds number is zero, and the time-derivative as well as the inertial term in the Navier-Stokes equations can neglected.This means we can use the steady state incompressible Stokes equations with periodic boundary conditions on the domain and no-slip boundary conditions, see Eq. ( 5), on the fiber surfaces where G is the volume that is occupied by fibers, Ω is the computational box [0, and ∂G is the surface of G. u is the periodic velocity vector with components, u, v and w. f = (0, 0, c) is a constant body force vector where c is derived from the desired pressure drop ∆P and the length of the computational domain in the z-direction, and p is the periodic pressure.That pressure can be made physical by adding the linear function that is independent of x and y and assumes the values ∆P at the inlet boundary and 0 at the outlet boundary.Figure 2 illustrates this behavior.It clarifies also the fact that L is only the thickness of the nonwoven media and does not include the inlet and outlet, as longer or slightly shorter inlets and outlets would not change the behavior of the flow.
Fig. 2 The physical pressure resulting from a simulation with inlet and outlet.There is no pressure loss in the empty space before and behind the nonwoven, and the slight deviation from linear over the nonwoven is typical for highly porous and uniform porous materials.
The linear systems of equations resulting from the discretization of Eqs. ( 3) -( 5) are solved in [17] using the LIR approach from [19,20] with improved discretization of the no-slip boundary conditions (following [21]) and greater efficiency than were proposed in [18] .

Fitting of Parametric Regression Formulas
For the prediction of permeability from geometric descriptors, we consider different parametric regression formulas which will be discussed in Section 3.2 below.For fitting and evaluation of the prediction formulas, we use the data obtained for structures from Scenario II, where we split it into a subset for fitting the formulas and a (smaller) subset for evaluating their performance.Fitting of the formulas is performed by means of SciPy's [22] built-in methods for least-squares optimization using the Levenberg-Marquardt algorithm.Note that the values of permeability obtained for the structures from Scenario II are within a much tighter range than, e.g., those considered in [8].Nevertheless, they span two orders of magnitude, such that the behavior for large values of permeability would disproportionately influence the fit of the formulas.Thus, we apply a log-transform prior to fitting.
To assess the quality of the fit, we use the mean absolute percentage error (MAPE), which is given by and the coefficient of determination R 2 given by , where y 1 , . . ., y k is the ground truth data, ŷ1 , . . ., ŷk are the corresponding predictions and ȳ = 1 k k j=1 y j is the mean of the ground truth data.The values of MAPE and R 2 are computed on (not log-transformed) test data which has not been used for fitting.However, note that the least-squares fit essentially optimizes the R 2 -value for the log-transformed data.

Results
For each structure generated in Scenarios I, II and III as outlined in Section 2.2, we computed the geometric descriptors described above as well as the corresponding permeability.In the following, we present the results of the statistical analysis of this data.This includes the investigation of microstructure-property relationships, where we consider various parametric regression formulas for the prediction of permeability from a range of geometric descriptors and discuss their performance.
Recall that for each specification of model parameters, 3 structures were generated which vary with respect to geometric descriptors and permeability.For structures of Scenario II, the coefficient of variation for the numerically computed values of permeability is equal to 60.10 when considering the full set of simulated structures.However, the mean coefficient of variation when considering each specification of model parameters individually is only equal to 3.06, indicating that structures which were simulated using the same specification of model parameters are rather similar.Fig. 3 Examples of simulated structures of Scenario II.The left structure exhibits a very low tortuosity, the middle and right structures have a low and high constrictivity, respectively.

Statistical Analysis of Simulated Fiber Systems
In a first step, we illustrate the range of feasible structures and investigate some basic relationships between model parameters and geometric descriptors resp.permeability.Structures of Scenario II, which were simulated using randomly chosen model parameters as discussed above, exhibit a wide range of values for many geometric descriptors, see Figure 3, which shows cutouts of 3 simulated structures of Scenario II.This can be made more precise by histograms of geometric descriptors and permeability computed for each of the simulated structures of Scenario II, see Figure 4.It can be clearly seen that none of these histograms can be modeled by a Gaussian distribution.Furthermore, we may use the data of Scenario III to investigate the relationship between fiber diameter, fiber density and permeability.Recall that the volume fraction of fibers (i.e., 1 − ε) is roughly a function of the squared fiber diameter while fiber density has a linear effect, see Eq. (1).As expected, increasing the fiber diameter while decreasing the fiber density to keep the porosity fixed has a positive effect on permeability, see Figure 5.
Finally, to investigate the influence of single model parameters on geometric descriptors, we consider structures of Scenario I. Figure 6 shows this kind of relationships for selected pairs of model parameters and descriptors.As expected, increasing 0.9 0.92 0.94 0.96 0.98 0.02 0.04 0.06 0.08 0.10 0.12 Fig. 5 Relationship between permeability, porosity, and fiber diameter.As expected, for given porosity, permeability is low when the structure is composed of a huge number of thin fibers as compared to a smaller number of thicker fibers.Each dot corresponds to a structure of Scenario III.
the fiber density changes many geometric descriptors, including the mean chord length.But other parameters also influence the morphology of the fiber systems.The value of β A which influences the curvature of fibers projected onto the x-y-plane has a strong influence on the constrictivity of the resulting structures.The dependence of geometric descriptors on other parameters is more complex and the effect of some model parameters on the morphology of the fiber systems may not be captured be the geometric descriptors stated in Section 2.3.Fig. 6 Relationships between selected pairs of model parameters and geometric descriptors, for structures of Scenario I.For each pair, we used the model parameters P 1 = (p 1 1 , . . ., p 1 14 ) (blue dots) and P 2 = (p 2 1 , . . ., p 2 14 ) (orange dots), which were obtained in [10] by fitting the model to measured data of two different nonwovens, and varied only one single parameter.Left: Constrictivity vs shape parameter of the distribution of A. Center: Mean spherical contact distance vs parameter α of the copula corresponding to (Z 2 − Z 1 , Z 1 ).Right: Mean chord length vs fiber density.
However, further information can be obtained when considering pairs of descriptors simultaneously.Figure 7 shows scatter plots of various pairs of descriptors obtained for structures of Scenario II.As expected, the porosity ε has a huge impact on the permeability, but other properties play an important role as well.

Quantitative Microstructure-Property Relationships
In this section we consider various parametric regression formulas which express permeability in term of geometric descriptors of the simulated fiber systems.In this context, we first analyze the correlation between different geometric descriptors and permeability.Figure 8 shows these correlations based on data of Scenario II.Note that  various geometric descriptors are strongly correlated permeability.As expected, porosity and permeability are positively correlated whereas mean geodesic tortuosity and permeability are negatively correlated.Moreover, most geometric descriptors are highly correlated with porosity, i.e., suggesting that a large amount of information about the morphology of the fiber systems is contained in a single geometric descriptor.However, besides porosity, the other geometric descriptors considered in this paper deliver further (refined) morphological information.Based on the data of Scenario II, we fitted different regression formulas to predict the permeability κ from geometric descriptors.In the following, we will discuss these formulas and their goodness of fit.Note that some formulas which were proposed in the literature, specifically for fibrous porous media, predict the quantity κ/d 2 instead of κ [9].However, as the fiber diameter d was not varied for the structures of Scenario II, we will ignore the fiber diameter and fit formulas to directly predict κ.
While some prediction formulas are designed observing physical units, others are solely targeted to obtain the best numerical values at the cost of a limited interpretability.As a baseline formula which does not preserve meaningful units, we consider the which represents one of the most simple relationships between porosity and permeability.Based on the data of Scenario II, we obtain the fitted parameters c 1 = 3.85 × 10 −8 and c 2 = 61.07.Additionally incorporating constrictivity β and specific surface are S, the following formula was introduced in [23], which preserves the physical units: where fitting led to the parameters c 1 = 5.03 × 10 −12 , c 2 = −24.88,c 3 = 0.13 and c 4 = 59.51.Furthermore, a slightly simplified version of Eq. ( 7) given by was considered in [24].Then, the fitted parameters are c 1 = 4.96 × 10 −12 , c 2 = −24.17and c 3 = 59.91.Deviating from pure power-law type formulas as described above, in [8] still another formula was introduced which is based on the same descriptors as Eq. ( 7): Again ignoring physical units, we consider the following slight modification of Eq. ( 7) which allows for a flexible choice of the exponent of S: Finally, we consider a power-law type formula incorporating a wider range of geometric descriptors to investigate the theoretically achievable predictive power of these types of relationships, where fitting the formula A visual impression of the predictive power of Eqs. ( 6) -( 13) can be obtained from Figure 9. Visually, all prediction formulas perform reasonably well which is corroborated by the values of R 2 and MAPE shown in Table 2.However, there is a clear trend for most prediction formulas showing a better fit for lower values of permeability.As the fits were performed on log-transformed data, this is probably due to the fact that only few structures exhibit the higher values of permeability while a huge amount of data is present for lower values.Especially, when considering Eq. ( 13), it becomes clear that the prediction can be improved significantly when considering further geometric descriptors as compared to, e.g., Eq. ( 6) which is solely based on the porosity ε.
Fig. 9 Permeability predicted by the parametric formulas given in Eqs. ( 6) -( 13), against the ground-truth permeability computed by the methods described in Section 2.4.Each dot corresponds to one structure from the test dataset of Scenario II.

Conclusion
For the efficient development of improved functional materials, understanding the relationship between easy to manipulate geometric descriptors of their 3D morphology and desired material properties is important.Traditional experimental setups for the investigation of these relationships are expensive and time-consuming and, furthermore, some important aspects may not be accurately measured.By using a stochastic microstructure model combined with numerical simulations, virtual materials testing may be performed to investigate these relationships.Compared to real experiments, our approach based on modeling and simulation is cheap and can be used to generate a huge database of virtual but realistic structures for which all the desired geometric descriptors and effective properties can be easily computed.From this data, a thorough understanding of the relationship between 3D morphology and effective properties of nonwovens can be obtained.Using this knowledge, the production process may be optimized by targeting solely geometric descriptors which are easier to manipulate than the effective material properties.
In this paper, we employed a stochastic model for nonwoven fiber materials to generate a large set of different structures and computed various geometric descriptors for each simulated structure.Furthermore, through numerical simulations, we computed permeability.To exemplify the approach of virtual materials testing, we then deduced parametric regression formulas to describe the relationship between 3D morphology and permeability.We used different types of formulas found in the literature as well as some custom-designed formulas and were able to show that permeability does not only depend on porosity, but very much on other geometric dscriptors as well.
Our approach enables material scientists and designers to make fast predictions about a material's effective properties based purely on easy-to-measure geometric descriptors.Using the prediction formulas considered in this paper, it might also be easier to design or modify a material in order to achieve a desired permeability, as it becomes clear how individual descriptors factor into the permeability of the material.For real-world production processes this means that instead of finding the right production parameters to obtain the desired material properties, it is only necessary to optimize the process to get 3D morphologies with the right geoemtric descriptors.If the relationship between production parameters and geometric descriptors of the material is known to the manufacturer, the entire pipeline of the material's design can be done virtually.
The methods shown in this paper can easily be extended to other material properties like wettability.Further work could use the large amount of virtual samples as training data for a machine learning based approach that could establish a fast and more precise, yet less interpretable link between geometric descriptors and effective properties of fiber-based materials.

Fig. 4
Fig. 4 Histograms of geometric descriptors and permeability computed for structures of Scenario II.

Fig. 7
Fig. 7 Relationship between selected pairs of geometric properties and permeability (color of dots) for structures of Scenario II.Left: Porosity vs constrictivity.Center: Mean chord length vs specific surface area.Right: Standard deviation of tortuosity vs mean chord length.

Fig. 8
Fig. 8 Correlation between geometric descriptors and permeability, based on data of Scenario II.

Table 1
Overview of the model parameters which will be varied for the present study.The given limits are chosen based on the parameter sets P 1 and P 2 .If applicable, a meaningful unit is given for each parameter.
A (0.3, 1.5) shape parameter of the distribution 1 modeling A and B β Z 1 (1, 30) shape parameter of the distribution 2 modeling Z 1

Table 2
Performance metrics of the fitted prediction formulas computed on the test dataset.R 2 log denotes the R 2 on the log scale which was essentially optimized during the fit.