Flow through a two-scale porous medium is here investigated by a unique comparison between simulations performed with computational fluid dynamics and the boundary element method with microparticle image velocimetry in model geometries.

1. Introduction

In general, prediction of porous media flow is straight-forward, as it is governed by Darcy’s law. However, when it comes to two-scale porous media, as in advanced composites manufacturing, paper making and drying of iron ore pellets, the detailed flow becomes complex [1].

In this letter a unique qualitative comparison of three methods to determine the flow dynamics in a dual-scale medium is carried out. The methods are Computational Fluid Dynamics (CFD), the Boundary Element Method (BEM), and microparticle image velocimetry (𝜇PIV). It is also suggested how these methods can interplay to produce the best results.

Nordlund et al. [2] have shown that CFD successfully can be used to calculate the permeability of model cells. The cells can then be connected to form a network for which the total permeability can be computed [3]. In order to do this the value of the permeability within a fibre bundle was approximated using the formulas suggested by Gebart [4] for flow along and perpendicular to regular fibre arrays according to

𝐾//=8𝐶(1𝑓)3𝑓2𝑅2,𝐾(1)=𝐶𝑓max𝑓15/2𝑅2,(2) where 𝑓 is the bundle fibre volume fraction, and 𝑅 is the fibre radius. The constant 𝐶 and the maximum 𝑓, 𝑓max are dependent on the arrangement of the fibres.

Due to its efficiency in very complex geometries (in the BEM only boundaries need to be meshed), the BEM can resolve every fibre and thus model the detailed microscale flow within the bundles [5, 6]. While in this approach the small scale can be perfectly modelled, it can be difficult to implement 3D effects on the fibre bundle scale as done with CFD. In this paper we apply both methods to the problem at hand.

With 𝜇PIV flows can be visualised, and quantitative measurements of instantaneous velocity fields can be performed [7, 8]. Hence, 𝜇PIV will here be applied to study the flow in a two-scale model geometry consisting of an array of parallel fibres placed in a cavity.

2. CFD Modelling

A computational model and the corresponding structural grid of an experimental flow cell were created using the software ANSYS CFX-11 and ICEM-CFD11. The geometry consists of an inlet pipe with 𝑅=0.6 mm leading to a 1.6 mm wide slit connected to the main channel, being a rectangular box with dimensions 5.3×7×7 mm3. The main channel has a 4×4 mm2 array of fibres, whose length takes up the entire channel depth and whose 𝑓40%.

The fibre array was modelled as a porous domain with permeability 𝐾 given by (2) with the constant 𝐶 and 𝑓max tuned for a hexagonal fibre arrangement. The flow rate at the inlet was set to 5.56𝑒–9 m3/s, and the outlet had a constant average static pressure. No slip conditions were applied on all walls. A second-order scheme was used, and measures were taken to ensure iterative and grid convergence as done in [2].

It is observed that the surrounding flow affects the velocity inside the bundle, resulting in a boundary layer flow within the porous medium with higher velocity near the edges. Figure 1 shows a contour plot of the velocity field. From additional simulations it can be observed that with increased permeability within the bundles the difference in velocity between the bundle flow and the open channel flow is decreasing, and the boundary layer becomes thicker.

3. BEM Modelling

The next approach is to use BEM to model 2D Stokes flow through dual-scale porous media now characterised by their interbundle porosity 𝜑𝑖, the porosity of the intrabundle space 𝜑𝑡=1𝑓, the number of fibres in the bundle 𝑁𝑓, R , and the bundle axis ratio 𝜆. In addition, the mean nearest fibre spacing 𝛿1 is included, to account for the degree of uniformity in the fibre placement. It has been established by Chen and Papathanasiou [5, 6] that as the degree of heterogeneity in the fibre distribution increases, 𝛿1 decreases. These fibre distributions can be generated by a Monte Carlo process. The minimum achievable value of 𝛿1 is zero, corresponding to the situation of all fibres touching each other, while the maximum 𝛿1 equals the interfibre spacing of a hexagonal array. The linear system of equations is derived from the discretization of the boundary integral equation [9] using quadratic shape functions for the geometry, the velocities, and the tractions. Symmetry boundary conditions are applied on the two horizontal boundaries of the unit cell. Unidirectional flow and constant total pressure drop along the two vertical boundaries are furthermore set; see Figure 2. On the surface of each fibre, no-slip conditions are applied. This results in a linear system of equations, 𝐀𝐱=𝐛 in which the matrix A is full and nonsymmetric. The system is solved with a parallel LU solver from the ScaLapack library [10], and an in-house parallel code is used for the formation of the matrix A.

There is a boundary layer flow formed at the border of the bundle in which the velocity is higher than the Darcy velocity for the CFD results as well as the ones from BEM (cf. Figures 1 and 2). This flow is largely affected by 𝛿1 as exemplified in Figure 2 for two contrasting 𝛿1 and for 𝜙𝑖=0.5, 𝜙𝑡=0.5, 𝜆=2 and 𝑁𝑓=500. More fluid penetrates into the nonuniform fibre bundle due to the large pores along its perimeter. Inside the nonuniform fibre bundle, the regions of faster flow also extend further into the interior of the bundle. When comparing a number of such results with predictions of the Phelan-Wise model [11], based on homogenous fibre bundles, there is also generally a difference. Experimental measurement of the permeability of systems similar to Figure 2 and back-calculation of 𝐾𝑡 by fitting the results in models derived in [2, 11] could give better conformity. Hence there is a need for new experimental methods such as 𝜇PIV.

4. μPIV Measurements

Experiments were conducted on the flow cell of Figure 3 with dimensions and key geometrical/flow features as in the CFD model where borosilicate glass fibre rods with 𝑅=150𝜇m form the porous domain.

Experiments were carried out on a rectangular fibre array with relatively low 𝑓 and a much denser hexagonally arranged array with 𝑓40%. Paraffin oil was used as fluid based on its refraction index properties seeded with Rhodamine B fluorescent particles from microparticles GmbH 10.2±0.17𝜇m. A KDS Model 100 Series syringe pump was used with the volume flow rate set to 20 mL/h. All measurements were carried out in the middle of the channel where the out-of-plane component of the flow is small.

There is limited movement between the rows of fibre although a circular motion can be observed; see Figure 4. This is followed by an increase in velocity in the passage between fibres due to a decreasing cross-section area; see Figure 5. Notice that even small geometrical deviations result in noticeable differences in the velocity field in the channels between the fibres. This is also the case for the denser array; see Figure 6.

So far all the measurements were taken close to the objective. As the measurement plane was moved deeper into the channel, a decrease in the optical quality of the pictures was noted.

A way to improve the quality is to average over the entire time-series and then subtracting that average from every image in the series. This is exemplified in Figure 7 showing the velocity field obtained in the centre of the channel. The 𝜇PIV technique can hence be used to visualise the flow within an array of fibres, even far away from the edge of the cell being the nearest to the objective.

An investigation was finally made in which the flow around the fibre bundle was captured with good results although the velocity near the fibres is much lower than the bulk flow, and this creates experimental problems.

5. Conclusions

The Boundary Element Method is an excellent approach in capturing the microscale details of the flow The results obtained could easily be used as input to the porous data used in CFD. The results from the CFD and BEM are in qualitative agreement with 𝜇PIV. The 𝜇PIV in itself is a very promising technique for experimental observations as well as quantitative derivations of the detailed flow and the permeability.


This work was partly financed by NANOFUN-POLY European network of excellence through SWEREA-SICOMP AB and the US National Science Foundation, IREE program. Dr. Papathanasiou also acknowledges support of the EU_People Programme, award PIRG-01-GA-2007-208341.