E ective permeability of a porous carbon composite is computed using the direct simulation Monte Carlo technique. The microstructure of the carbon compos- ite is synthetically generated using an in-house solver. Permeabilities obtained using synthetic microstructures of the precursor (carbon bers) is compared to two independent experimental dataset and good agreement is observed. An approach to digitally infuse matrix into the precursor is proposed. Compari- son of the permeability for the full composite (carbon bers with matrix) with experimental data demonstrates good agreement indicating that representative microstructures generated digitally can be used to compute and predict e ective permeability of porous carbon composites. Simulations of gases penetrating the full composite are performed, and an intrinsic material permeability (K o) of 11:8661011m2, and a Klinkenberg constant (b) of 5:655108that is only dependent on the temperature and the molecular weight of the gaseous species is obtained. 1. Introduction Porous materials have many present and potential future applications for catalysis, adsorption and separation, life science, energy, environmental and high-end technologies [1]. E ective permeability is an important property of a Corresponding author Email address: saviopoovathingal@uky.edu (Savio J. Poovathingal) Preprint submitted to International Journal of Heat and Mass Transfer January 29, 2021 porous material as it dictates the ow of gases and liquids through the pores of 5 the material. E ective permeability is used as a design property for extracting shale gas in the energy sector [2] and for molecular sieves for ltration [3]. E ective permeability is also an important property for carbon composites that are used as thermal protection systems (TPS) to mitigate the heat load during the entry process of a space capsule. 10 Carbon composites used as TPS materials are low-density materials that consist of a network of carbon bers infused with a resin (also referred to as the matrix). These materials relieve heat load through a number of processes, in- cluding absorption processes such as thermochemical decomposition of the ller resin and sublimation, and heat rejection processes through modi cation of the 15 high-temperature boundary layer through transpiration of ablation products, and systematic sacri cial removal of hot material through ablation [4]. The composites used on space capsules are highly porous with the ber and resin occupying not more than 20% of the total volume. The high porosity allows gases to ow through the material. For example, the resin decomposes (pro- 20 cess called pyrolysis) and travels through the charred (completely pyrolyzed) network of bers. Similarly, boundary layer gases can penetrate the material causing in-depth gas-phase and gas-surface reactions [5]. Flow through TPS materials are computed by solving the combined solid-gas momentum transport equations. Historically, the momentum transport equations have been formu- 25 lated as either surface models, one-dimensional volume averaged models [6, 7], and more recently have been formulated as two/three-dimensional volume av- eraged models [8{16]. Such volume-averaged approaches [17] for momentum transport require model closure, in this case, the e ective permeability of the material. Since the entry process occur at high altitudes, the ow inside these 30 carbon composites occur in the non-continuum regime. The analysis of Stardust trajectory indicates that the ow is either in the transition or slip regime while it ows through the pores of the carbon composite [18]. Therefore, analytical relations that are derived either in the continuum or the free-molecular limit cannot be used as the e ective permeability for ight relevant conditions. 35 2 A widely used carbon composite on various space missions consists of a carbon ber network called FiberForm (manufactured by FiberMaterials Inc.) into which a phenolic resin is infused. FiberForm is manufactured by chopping bers and pressing them into billets which creates a random network of bers with a preferential direction normal to the pressing plane creating an orthotropic 40 material. The resin is infused into this network which attaches to the bers. E ective permeability for FiberForm and the full composite has been obtained empirically by Marschall and Milos [19, 20]. Marschall and Milos performed the experiments using a ow tube apparatus with air as the test gas and all experiments were performed at room temperature. Recently, Panerai et al. [5] 45 performed ow tube experiments to obtain e ective permeability of FiberForm over nine di erent samples with di erent densities at di erent temperatures using Ar as the test gas. They found signi cant variability in the permeability values for di erent samples. E ective permeability of FiberForm has been computationally obtained by 50 various groups [21{26]. Borner et al. [21], White et al. [22], and Jambunathan et al. [23] used direct simulation Monte Carlo (DSMC) as the ow solver while Ho et al. [26, 27] used Lattice Boltzmann as the ow solver. In most of these cases, the microstructure of FiberForm was imported from x-ray computed to- mography (XRCT) images of the FiberForm sample. Since XRCT images were 55 obtained only for a few samples, it limits the validation e orts with experi- ments and the ability to explore the material parameter space that in uences e ective permeability. Stern et al. [25] developed a synthetic microstructure generating code called Fibergen to reproduce the microstructure of FiberForm. E ective permeability of the synthetically generated microstructure was com- 60 puted through DSMC simulations using the molecular gas dynamics solver [24]. It was shown that the microstructures of FiberForm generated using Fibergen reproduced the experimental results of Marschall and Milos [19]. Since Fiber- gen is a computational code, several parameters of the code can be controlled, allowing a systematic study of the parameters that in uence e ective material 65 properties. A systematic study was performed by the authors by varying the 3 ber radius and elevation angles. It was found that a ber radius of 5 mand an elevation angle of 25resulted in the best agreement between the experiments and simulations. However, it was also observed that moderate variations in ber radius, angle, and the inherent randomness in ber arrangement combine 70 to give approximately 30% variation in the computed e ective permeability value for a given pressure condition. The comparisons of permeability for mi- crostructures generated using Fibergen in the previous study was limited to the data of Marschall and Milos [19]. The objective of this article is twofold. New DSMC simulations are per- 75 formed to compute e ective permeability using synthetic microstructures gen- erated from Fibergen to demonstrate the wider applicability of using synthetic microstructures. The e ective permeability computed from the DSMC sim- ulations are compared against the experimental dataset of Panerai et al. [5]. Second, a method to synthetically infuse resin is proposed and this approach 80 is validated against experimental data that is available for the full composite (FiberForm infused with resin) [20]. All previous computational work have been limited to computing the permeability of FiberForm. It is also shown that an expression developed for FiberForm [5] can be extended to the full composite and e ective permeability values are derived that can be used in computational 85 uid dynamics (CFD) or material response solvers. 2. Methods Synthetic microstructures are generated using the Fibergen code and the permeability calculations are performed using the DSMC technique. Both ap- proaches are brie y discussed in this section. A more detailed description of the 90 code layouts, eciency considerations, and coupling strategies can be found in previous publications [24, 25]. 2.1. Fibergen: A Synthetic Microstructure Generating Code Fibergen uses dimensions of the domain to be simulated (also referred to as box in this article), target bulk porosity, nominal ber orientation, nominal 95 4 ber radius, and speci ed variances to generate cylindrical bers. To repli- cate the FiberForm microstructure, the radius of a cylindrical ber is sampled from a Gaussian distribution with a mean of 5 mand standard deviation of 0:1m[28]. FiberForm is orthotropic because of the manufacturing process wherein bers are pressed into billets. This pressing a ects the ber angles 100 because they have a bias toward being oriented parallel to the pressing plane (also referred to as the through the thickness plane) [24]. Therefore, elevation angle of the ber is sampled from a uniform distribution of 25. The nominal ber radius, elevation angle, and the standard deviations are chosen based on the conclusions of the parametric study performed previously [24]. After the 105 rst cylindrical ber is generated, successive bers are added only if the new ber does not intersect or overlap with other bers in the domain. Previously, bers were generated until a target porosity was attained. Several details of the Fibergen code along with the algorithm, code layout and architecture have been described before and is available in Ref. 25. Two new features have been 110 added to Fibergen which are described below. First, bers can now be generated until a target sample density is attained. In this mode, Fibergen uses the density of a single carbon ber as an input and bers are generated in the same manner as described above until the density of the simulated box reaches a target value. Data from Panerai et al. [5] is 115 available for eight samples with di erent densities. Therefore, modifying the convergent criteria of Fibergen allows a direct comparison with experimental data of Panerai et al. [5]. A sample microstructure generated using Fibergen is shown in Fig. 1a. The choice of density for an individual ber is discussed in Sec. 3. Previously, Fibergen created surface triangles on the bers that can be 120 used as the surface mesh in the DSMC solver. The second new feature that is added is to allow the output of the microstructure in three-dimensional voxel format. In the voxel format, a three-dimensional box containing the bers are discretized to a desired resolution. Nodal points inside and outside the bers are marked as 1 and 0, respectively. This voxelized feature is leveraged to 125 incorporate resin into the microstructure. A sample microstructure generated 5 at a voxel resolution of 2 :5mwithout the resin is shown in Fig. 1b. To incorporate the resin, the bers are initially generated as described above. The box containing the bers is discretized to a desired resolution and the vol- ume occupied by the bers in each cell is computed using the cut-cell tech- 130 nique [29]. In the manufacturing of the composite, the resin is infused into the microstructure of FiberForm which bonds to the carbon bers in a random manner. Since the resin randomly infuses into the composite, a similar approach is taken here. A discretized cell within the three-dimensional box is chosen ran- domly. Cells that are completely inside the ber or that do not contain any 135 ber (no cut-volume) are ignored and a new cell is chosen until a cell with par- tially occupied ber is selected. The volume occupied by the ber inside the cell is known through the cut-cell technique [29] and the resin is added to the remaining volume. All corners of the chosen cell are marked as 1 which assumes that the entire remaining volume is occupied by the resin. Resin is added until 140 the density of the composite reaches the target value. The nal voxelized mesh is used as the input into the DSMC solver. A sample microstructure generated at a voxel resolution of 2 :5mwith the resin is shown in Fig. 1c. 2.2. Sparta: DSMC Solver The numerical investigation is performed using the DSMC technique. The 145 DSMC technique is a stochastic approach that simulates the Boltzmann equa- tion. A detailed description of the DSMC technique and models can be found in Refs. 30, 31. The relevance and advantage of using the DSMC technique to com- pute e ective properties for porous TPS materials can be found in Refs. 24, 25. The DSMC simulations are performed using the open-source DSMC solver 150 Stochastic Parallel Rare ed-gas Time-accurate Analyzer (SPARTA) [32{34], which uses a multi-level Cartesian mesh to track and collide particles. Sur- faces can be identi ed in two ways. A set of triangles representing the surfaces is directly read into Sparta or a set of voxels are converted to a closed sur- face with triangles using the marching cube algorithm [32, 35]. The triangles 155 are sorted (for marching cube, sorting is not required) and the ow volume in 6 each cell is computed using the cut-cell technique. Both approaches are used to import microstructures of FiberForm (see Sec. 3 for related discussion). The standard ray-tracing approach is used to move the simulation particles. The ow setup for the permeability simulations is shown in Fig. 1d. As seen 160 from Fig. 1d, the microstructure is incorporated into the center of the simulated domain. Symmetric boundary conditions are applied in the two span-wise di- rections (y and z). In the upstream boundary, pressure (P in) and temperature (T) are speci ed. The in ow velocity (u) is obtained using a simple zeroth order extrapolation, i.e. u =u j, wherejrefers to the value at the rst interior cell. 165 In the downstream boundary, the method of characteristics of Nance et al. [36] and Fang and Liou [37] is used. Speci cally, the out ow pressure (P out) is spec- i ed and the remaining macroscopic properties are computed using information from the interior solution. The relevant equations to obtain the boundary in- formation are available in Refs. [25, 36, 37]. These boundary conditions have 170 been validated in the previous study [24]. For the mean ow quantities from the interior cell, instantaneous mean ow quantities can be used or a sub-relaxation averaging can be applied [25]. We found that both approaches provide the same accuracy but the sub-relaxation technique is more stable and robust. The simulations are run until steady state and a permeability force F is 175 computed as described in Marschall and Milos [19] (Eqn. 1). F_mTRL AMP=K0(Pav+b); (1) In Eqn. 1, F is the permeability force, is the viscosity, _ mis the mass ow rate, T is the temperature, R is the universal gas constant, L is the length of the sample, A is the area, M is the molecular mass,  Pis the pressure di erence across the sample, K 0is the permeability in the continuum limit, and b is the 180 Klinkenberg constant. A power-law viscosity model is used to compute F, which is shown in Eqn. 2. =15p mkTref 2(52!)(72!)d2 refT Tref! (2) 7 In Eqn. 2, k is the Boltzmann constant, T refis a reference temperature of 273 K,drefis the reference diameter, and !is the variable hard sphere (VHS) exponent. The reference diameter ( dref) of Ar, N 2, and O 2are 4:171010, 185 4:171010, and 3:961010m, respectively. The VHS exponent of Ar, N 2, and O 2are 0.81, 0.77, and 0.74 respectively. We note here that viscosity is not directly used in DSMC simulations but rather through the speci cation of dref and the VHS exponent. The physical meaning and derivation of Eqn. 1 can be found in Refs. [5, 24]. 190 In the experiments [5, 19, 20], F is plotted as a function of P avgand a curve t is performed to obtain K oand b. K odepends on the geometric property of the material (arrangement of bers and resin) and b depends on the type of gas and ow temperature. As discussed in Sec. 3, in the validation step, we compare the value of F directly with experimental data. When new data is generated for the 195 full composite, a linear least square t is performed to obtain K oand b. This formulation for the e ective permeability follows the Klinkenberg formulation, which is given in Eqn. 3. K=Ko 1 +8c dp (3) In Eqn. 3, K ois the permeability constant, d pis the pore diameter, is the mean free path, and c is a proportionality constant. 200 3. Results 3.1. Permeability computations of FiberForm In all simulations discussed in this article, a nominal pressure di erence of 100 Pa is prescribed for the boundary conditions. For example, for an average pressure (P avg) of 500 Pa, P inand Poutis set to 550 and 450 Pa, respectively. 205 The ow is allowed to evolve based on the boundary conditions. The nal Pavgand Pare then computed by extracting the pressure before and after the sample. Therefore, the nal P avgand Pwill not be exactly equal to the speci ed values, although it will be close to the speci ed values. A convergence 8 (a) Surface mesh for a microstructure gener- ated as a set of triangles. (b) Voxel representation of a microstructure with- out resin at a resolution of 2.5 m. (c) Voxel representation of a microstructure with resin at a resolution of 2.5 m. (d) Flow setup for all simulations. Figure 1: Surface meshes and ow setup schematic for all simulations. In Figs. 1b and 1c, black voxels represent carbon bers, while the red voxels in Fig. 1c represents the resin. test is initially performed on the size of the microstructure that needs to be 210 simulated. This is referred to as the representative elementary volume (REV). Microstructures are generated with a xed target porosity of 0.9. Box size in the span-wise direction is varied between 125 300mand between 200500m in the ow direction (the through thickness direction). P avg, P, and F for all cases are shown in Table 1 along with the error which is computed in comparison 215 with a microstructure with a domain size of 500 300300m. For the simulations described in Table 1, Ar is used as the test gas. It is observed that the error is low once the REV is at least 200 200200m. It would 9 appear that the error is higher with a box size of 500 250250m compared to a box size of 500 200200m(entries 2 and 3 in Table 1) 220 but it is merely an artifact of the variability of FiberForm. A dominant feature of FiberForm microstructure is the variability in the material properties [5] as seen in experiments and our previous computations [24]. Fibergen attempts to capture this e ect by sampling over distributions in a stochastic manner. Changes in the overall domain size will result in small changes in how the 225 bers are arranged, which in turn in uences the e ective permeability. It was previously observed that moderate changes in the ber radii and elevation angles result in variations up to 30% [24]. Therefore, all simulations with low error are equally valid solutions with sucient accuracy. To balance REV size with computational load, an REV of 200 200200mis used for all the 230 simulations in this article. The REV is used to perform a convergence test on the DSMC simulations for the number of simulated particles ( Np) , cell resolution ( x), and timestep ( t). Based on the convergence test, Np= 50, x= 2, t= 20is used where is the mean free path and is the mean collision time. A convergence test is also performed on the number of time steps 235 simulated until steady state and the number of timesteps over which mean ow mean quantities are averaged. Convergence is obtained by running simulations for 100,000 timesteps until steady state and averaged for an additional 100,000 timesteps. The simulated domain in the x-direction is 600 mwith 200mon either side (equal to the REV length) and a convergence test is also performed 240 on the required domain length in the x-direction. Next, the e ect of sample density on the permeability force (F) is examined. As discussed in Sec. 2.1, Fibergen uses individual ber density as an input to generate a microstructure with a target sample density. FiberForm is manu- factured from rayon-based bers for which the individual ber density varies 245 from 14001700 kg m3[38]. For a given sample density, a higher individual ber density results in lower number of bers inserted into the REV resulting in higher porosity, thereby making the microstructure more permeable. The permeability force F obtained for a sample density of 181 and 157 kg m3, and 10 Table 1: Permeability force (F) obtained with various REV sizes. All error is computed with reference to an REV size of 500 300300mfor which the force, F = 1.378e-7 N. Microstructure size (inm) X x Y x ZPavg(Pa) P(Pa) F (N)Absolute Error (%) 500 x 125 x 125 533.5 77.06 1.69e-07 22.74 500 x 200 x 200 514.6 90.19 1.417e-07 2.4 500 x 250 x 250 509.74 93.41 1.343e-07 2.68 200 x 125 x 125 535.0 74.67 1.734e-07 25.64 200 x 200 x 200 515 89.6 1.45e-07 5.3 200 x 250 x 250 521.96 86.86 1.45e-07 5.3 200 x 300 x 300 507.14 95.18 1.404e-07 1.673 ber densities varying from 1400 1700 kg m3is shown in Table 2. The cor- 250 responding porosity is also shown in Table 2 for reference. For all simulations in Table 2, P inand Poutis set to 550 and 450 Pa, respectively and ow tem- perature is set to 297 K. Ar is used as the test gas for the simulations with a sample density of 181 kg m3while air is used for computing the permeability for the sample density of 157 kg m3. These test gases correspond to the ap- 255 propriate gases used in the experimental study. As expected, it is seen that for a given sample density, lower values of F are obtained when the ber density is lower (material is less permeable). The sample density of 157 kg m3corre- sponds to one of the experimental data of Marschall and Milos [19], while the sample density of 181 kg m3corresponds to the dataset of Panerai et al. [5]. 260 First, we observe that small changes in porosity result in huge variations in the permeability values, which reinforces the conclusion that e ective properties of FiberForm exhibit high variability. Panerai et al. used the average density of all the experimental samples and the reported open porosity of FiberForm to estimate the ber density as 1400 kg m3. We also see excellent agreement be- 265 tween the experimental data and the simulation results with an individual ber 11 density of 1400 kg m3. Speci cally, we are able to capture the permeability force for independent datasets of Marschall and Milos [19] and Panerai et al. [5] merely by specifying the sample density. In the dataset of Panerai et al., two samples (TT04 and TT09) with the same sample density of 181 kg m3result 270 in drastically di erent permeability values [5] (see Table 1 of Ref. [5]). The permeability force (F) of sample TT09 is very high and can be reproduced only with a porosity of 0.91. Such variability is expected in FiberForm and could also arise from defects and impurities [28]. To assess the robustness of Fibergen and the DSMC simulations, the comparison of the DSMC simulations with all of 275 the experimental data of Panerai et al. is shown in Table A.3 in the appendix. While the average pressure for most simulations correspond to 500 Pa, some cases are also simulated at 250 Pa to assess if the DSMC simulations capture the e ective permeability values at di erent pressures. The ow temperature corresponds to the appropriate experimental value (Table 1 of Ref. [5]). It is 280 observed that the experimental values are being reproduced over a wide range of sample densities and temperatures. There is a systematic trend in some samples where the error is higher at lower pressures and higher temperatures. Consider the sample at a density of 181 kg m3(entry 4 in Table A.3) where the error is 0.8% at 297 K, but at 723 K, it is 29% and 37% for an average pressure of 500 285 and 250 Pa, respectively. There are two possible reasons for this discrepancy. The experiments were performed at higher pressures (all greater than 1000 Pa) and the experimental t used to obtain F could lose accuracy at values outside this range. However, one would expect this error to be random and not exhibit a systematic trend based on pressure or temperature. A plausible reason is 290 the limitation of Klinkenberg relations in the rare ed regime. The Klinkenberg formulation to which the experimental data is t (Eqn. 3) accounts for correc- tion related to rst-order rarefaction e ects. Beskok et al. [39] have derived a second-order correction factor (f c) for the permeability, given by Eqn. 4. In Eqn. 4, Kn is the Knudsen number, and b are constants [40]. The Klinkenberg 295 formulation only accounts for the (1+ Kn ) term. All values predicted from the DSMC computations are higher than the experimental permeability force (es- 12 timated through Eqn. 1) when the error is high suggesting that a second-order correction term may be required for some of the simulated pressures and tem- peratures. New experiments are being planned and since pressures of 250 and 300 500 Pa are relevant for ight conditions [18], second-order rarefaction e ects will be examined in the future. Despite a few cases with errors greater than 30%, it is seen that the experimental data is being reproduced with good accuracy for a wide range of FiberForm samples. It is important to emphasize here that the data is being reproduced through physical values of ber density and sample 305 densities into Fibergen as an input and not adjustable parameters. This result provides further evidence that representative microstructures generated using Fibergen can be used to compute and predict e ective material properties of FiberForm. fc= [1 + Kn ] 1 +4Kn 1bKn (4) Table 2: Comparison of permeability force computed using DSMC simulations with the ex- perimental data of Marschall and Milos (157 kg m3) [19] and Panerai et al.(181 kg m3) [5]. Sample density (kg/m3)Fiber density (kg/m3)Porosity ( )Fexpt (N)Fdsmc (N) 1571600 0.902 1.33e-071.71e-07 1400 0.887 1.38e-07 1810.911 1.19e-072 1.96e-0732.01e-07 1700 0.893 1.58e-07 1600 0.887 1.45e-07 1400 0.87 1.22e-07 1Microstructure generated with higher porosity to compare with experimental data of TT09 sample from Panerai et al. 2Sample TT04 in Ref. 5. 3Sample TT09 in Ref. 5 In the simulations described above, the microstructure was read into Sparta 310 13 as a set of triangles (herein referred to as the stl format). We now compute e ective permeability by generating surfaces in the form of voxels. In the voxel representation, constructing a perfect cylinder would require in nite resolution because it attempts to reproduce the cylinders with a set of cubes. Therefore, a convergence test is performed to determine the resolution required to capture 315 the e ective permeability. A microstructure with an individual ber density of 1400 kg m3and sample density of 157 kg m3is generated with an REV size of 200200200m. The surface mesh is generated at a resolution of 5, 2.5, and 1.25 mrespectively. The permeability force computed at a nominal average pressure of 500 Pa at di erent voxel resolution is shown in Fig. 2. The 320 permeability force values computed by the voxel approach is lower than the stl format at all resolutions. This is expected because a voxelized cylinder diameter will always be greater than the true cylinder diameter which the stl format captures. A larger cylinder diameter implies that the pathway for the uid is more constricted making the material less permeable (hence the lower value of 325 F). A resolution of at least 1.25 mis required to reduce the error to 6% when voxels are used to represent the surface mesh. A similar resolution requirement was also required for lattice Boltzmann simulations [26]. A resolution of 1.25 m results in computationally expensive calculations compared to the stl format. In DSMC simulations, the required cell size for permeability computations is 330 O(). The mean free path ( ) for Ar at a pressure of 500 Pa and 300 K is 105m. For a total simulation domain of 600 200200m, a resolution of 1.25m(required voxel resolution) will result in 12 million cells and 600 million simulation particles. While these computations are certainly feasible, the calculations consume signi cant compute hours. In comparison, the stl 335 format requires only 192,000 cells (resolved to O()) with 10 million particles resulting in computations that are orders of magnitude faster allowing rapid computations using the stl format. 14 Figure 2: Permeability force computed at an average pressure of 500 Pa and 300 K with various voxel resolutions. The sample density is 157 kg m3and Ar is used as the test gas for the simulations. True value represents the force computed using the stl format. 3.2. Permeability computations of the full composite In this section, DSMC simulations to compute e ective permeability of the 340 full composite is discussed. Although expensive, the voxel format is required for the full composite because the phenolic resin that bonds to the bers does not form an ordered structure. The resin is infused into the FiberForm microstruc- ture digitally through the process described in Sec. 2.1. All microstructures are generated with a resolution of 1.25 mand the simulations are performed at 300 345 K with air as the test gas. The simulations of the full composite are compared to the data from Marschall and Milos [20]. The construction of the microstructure for the full composite requires the speci cation of the individual ber density, FiberForm (precursor) density, phenolic resin density, and the full composite density. Marschall and Milos only provide the full composite densities [20]. For 350 the synthetic microstructures, the individual ber density is xed at 1400 kg m3based on the results discussed in Sec. 3.1, and the density of phenolic resin 15 is set to 1200 kg m3[41]. To assess the e ect of FiberForm density, microstructures are generated with FiberForm and full composite densities of 149 and 217 kg m3, respectively. It 355 is observed that the computations over-predict the permeability force with the speci ed densities of 149 and 217 kg m3(red curve in Fig. 3). To assess the sensitivity, F is computed for a microstructure generated with a FiberForm den- sity of 149 kg m3and a full composite density of 210 kg m3(pink symbol in Fig. 3). It is observed that there is no change in F when the full composite 360 density is varied indicating that F is insensitive to small changes in the amount of resin infused into the microstructure. The volume fractions of resin added for a composite density of 210 and 217 kg m3are 0.07 and 0.09, respectively. The calculations are repeated with a sample FiberForm (precursor) density of 157 kg m3. Additionally, simulations are also performed to compute F for a 365 microstructure generated with the density of phenolic resin set to 1400 kg m3 (orange symbol in Fig. 3). A similar behavior is observed with a FiberForm density of 157 kg m3, where small changes in the amount of resin added to the microstructure does not in uence the e ective permeability of the composite. The FiberForm sample density for the full composite is not directly available, 370 and it was initially estimated as 149 and 157 kg m3based on the pure Fiber- Form sample data (Refs. [19, 20]). The DSMC simulations indicate that small changes in the amount of resin infused does not signi cantly alter the e ec- tive permeability. In contrast, small changes in the FiberForm density result in signi cant changes in the e ective permeability (see Sec. 3.1 and associated 375 discussion). The experimental data in Fig. 3 indicates variability with minor changes in the full composite density but the trends are not consistent. The per- meability force (F) of the composite with a density of 211 kg m3is lower than the composite with a density of 209 kg m3(Fig. 3). The inconsistency could arise from di erences in the density of FiberForm into which the resin is infused. 380 The permeability force F computed with a FiberForm density of 180 kg m3 and a full composite density of 210 kg m3lies within 5 % of the experimental data. The FiberForm density of 180 kg m3lies within the expected range of 16 FiberForm densities [5, 19] suggesting that the approach of infusing resin dis- cussed in this article does reproduce the microstructures of the full composite. 385 However, new targeted experiments with precursor and full composite densities are required for further validation. Figure 3: Comparison of DSMC simulations with experimental data of the full composite. The FiberForm, full composite, and resin density is denoted as f,c, andr, respectively. Finally, the full composite ( c= 210 kg m3) generated with a FiberForm density of 180 kg m3is used to compute the permeability force using Ar as the test gas to predict values for the Klinkenberg formulation. The approach 390 follows the work in Ref. [5]. A constant (b b) is rst computed from Eqn. 5 which is used to normalize the data. In Eqns. 5, refers to a reference state, which is considered to be Ar at 300 K. The normalized force (F) is then given by Eqn. 6. The normalized force (F) collapses to a single curve for both Ar and air as shown in Fig. 4. A linear curve t is performed using Fand Pavgusing 395 data from both Ar and air to obtain new K oand b. 17 bi b=i(T) (T)r  2RT Mir 2 M RT(5) F=Ko(Pavg+b) (6) The curve t provides an e ective permeability, K eff=Ko(P+b), where Ko= 11:8661011m2and b= 5:655108. This e ective permeability is valid in the pressure range of the simulations for both Ar and air. Since the normalization is based purely on temperature, viscosity relation, and molecular 400 weight, the derived relation could be used for other conditions. To assess the feasibility of the new relation for di erent conditions, the permeability for the full composite is computed with CO as the test gas at 300 K and 650 Pa resulting in Fdsmc = 6:4108N, while the predicted value from the new relation is 7:3108N (error of 11 %). The error is low in comparison to the expected 405 variability of the material ( 30%) suggesting that the new relation could be extended to other gases that have not been explicitly simulated in this article. We note here that the phenolic resin undergoes thermal decomposition ejecting a wide range of gaseous species [42]. Additional validation and re nement of this relation with di erent gases ejected during pyrolysis will be performed in 410 the future. 4. Conclusion The direct simulation Monte Carlo (DSMC) method was used to compute the e ective material permeability of porous carbon composites. A synthetic microstructure generating called Fibergen was used to generate representative 415 microstructures of the composite consisting of carbon bers and phenolic resin. The simulations were compared against the experimental dataset of Marschall and Milos and Panerai et al. [5]. It was observed that the microstructures of FiberForm generated using Fibergen reproduced the e ective permeability of both datasets with errors within 10 % at low temperatures and high pres- 420 sures. A systematic increase in the error was observed at high temperatures 18 Figure 4: Normalized force (F) and the corresponding curve t to obtain the intrinsic per- meability (K o) and the permeability constant b. and low pressures for some samples, which could arise from the limitation of the Klinkenberg formulation that was used to interpret the experimental and simulation results. The initial validation with FiberForm was extended to the full composite 425 through a novel infusion approach for the digital microstructures. Comparison of the e ective permeability with the experimental data indicated good agreement. A new relation for the e ective permeability was obtained through additional DSMC simulations, and it was demonstrated that the relation could be used for gaseous species that were not explicitly simulated to derive the relation. Avenues 430 for new targeted experiments have been identi ed and additional validation will be performed in the future. 5. Acknowledgements Financial support from NASA Kentucky EPSCoR under NASA grant num- ber 80NSSC19M0052 is acknowledged. Brendan Soto would like to acknowledge 435 19 support from NASA Kentucky Space Grant Fellowship under grant number 80NSSC20M0047. Cameron Brewer is a fellow with the Chellgren Center for Undergraduate Excellence. Insightful discussions with Dr. Alexandre Martin on the Klinkenberg formulation is deeply appreciated. We also thank the Uni- versity of Kentucky Center for Computational Sciences and Information Tech- 440 nology Services Research Computing for their support and use of the Lipscomb Compute Cluster and associated research computing resources.