Pressure drop in flow across ceramic foams—A numerical and experimental study

Author's Accepted Manuscript

Pressure drop in flow across ceramic foams - a numerical and experimental study W. Regulski, J. Szumbarski, Ł. Łaniewski-Wołłk, K. Gumowski, J. Skibiński, M. Wichrowski, T. Wejrzanowski

www.elsevier.com/locate/ces

PII: DOI: Reference:

S0009-2509(15)00451-0 http://dx.doi.org/10.1016/j.ces.2015.06.043 CES12446

To appear in:

Chemical Engineering Science

Received date: 19 February 2015 Revised date: 11 June 2015 Accepted date: 13 June 2015 Cite this article as: W. Regulski, J. Szumbarski, Ł. Łaniewski-Wołłk, K. Gumowski, J. Skibiński, M. Wichrowski, T. Wejrzanowski, Pressure drop in flow across ceramic foams - a numerical and experimental study, Chemical Engineering Science, http://dx.doi.org/10.1016/j.ces.2015.06.043 This is a PDF file of an unedited manuscript that has been accepted for publication. As a service to our customers we are providing this early version of the manuscript. The manuscript will undergo copyediting, typesetting, and review of the resulting galley proof before it is published in its final citable form. Please note that during the production process errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain.

Pressure drop in flow across ceramic foams - a numerical and experimental study W. Regulskia,∗, J. Szumbarskia , Ł. Łaniewski-Wołłka , K. Gumowskia , J. Skibińskib , M. Wichrowskic , T. Wejrzanowskib a

Institute of Aeronautics and Applied Mechanics. Warsaw University of Technology b Faculty of Materials Engineering. Warsaw University of Technology c Institute of Fundamental Technological Research (IPPT). Polish Academy of Sciences

5

Abstract The unique properties of ceramic foams make them well suited to a range of applications in science and engineering such as heat transfer, reaction catalysis, flow stabilization, and filtration. Consequently, a detailed understanding of the transport properties (i.e. permeability, pressure drop) of these foams is essential. This paper presents the results of both numerical and experimental investigations of the morphology and pressure drop in 10ppi (pores per inch), 20ppi and 30ppi ceramic foam specimens with porosity in the range of 75 to 79%. The numerical simulations were carried out using a GPU implementation of the three-dimensional, multiple-relaxation-time lattice Boltzmann method (MRT-LBM) on geometries of up to 360 million nodes in size. The experiments were undertaken using a water channel. Foam morphology (porosity and specific surface area) was studied on post-processed, computed tomography (CT) images, and the sensitivity of these results to CT image thresholding was also investigated. Comparison of the numerical and experimental data for pressure drop exhibited very good agreement. Additionally, the results of this study were verified against other researchers’ data and correlations, with varying outcomes. Keywords: ceramic foam, pressure drop, lattice Boltzmann method, Darcy-Forchheimer equation, anisotropy, specific surface area, pore-scale simulation 10

1. Introduction The industrial importance of materials with open porosity structures in the form of ceramic or metallic foams has grown in recent years. These materials exhibit specific properties such as high specific surface area, high porosity, low density, favourable mechanical, thermal and corrosion resistance. Thus they are well suited to serve as compact heat exchangers, reaction catalyst support, flow stabilizers or filters (Twigg

15

and Richardson, 2007). This results in an increased need for a priori knowledge of their hydrodynamic ∗

Corresponding author Email address: [email protected] (W. Regulski) Preprint submitted to Chemical Engineering Science

June 25, 2015

1

INTRODUCTION

2

properties. Thus, much effort has been devoted to link the foam structural parameters to their pressure drop. An extensive review of pressure drop correlations for foams was provided by Edouard et al. (2008). All presented formulae are based on the Darcy-Forchheimer relation, ρ μ Δp = U + U 2, ΔL k1 k2 20

where

Δp ΔL

(1)

is the average pressure gradient, U is the mean flow velocity, ρ is the fluid density, μ is the

fluid dynamic viscosity, and the coefficients k1 and k2 are called viscous and inertial permeabilities. These permeabilities should be functions of the foam’s geometry. The key geometrical parameters of the foam are its porosity, ε =

Vvoid Vtotal ,

which is the ratio of the void space, Vvoid , to total volume, Vtotal , occupied by

the foam, foam strut diameter, dstrut , mean pore diameter, dp , pore per linear inch number, ppi, and 25

specific surface area, ac =

S Vtotal ,

ratio of the surface in contact with the flow, S, to the whole volume

occupied by the specimen. Obviously, some of these parameters are correlated. The primary purpose of this study is to present results of morphology assessment and a series of simulations performed for typical ceramic foam structures. Three foams of nominal pore densities 10ppi, 20ppi and 30ppi and porosities in the range ε = 75%−79% are investigated. The obtained results (porosity, 30

specific surface area and pressure drop) are compared to data available in the literature. Pressure drops are compared to those of similar specimens reported in other works. Additionally, the appropriateness of pressure drop correlations proposed by other authors is explored. Special attention is paid to the post-processing of the computer tomography (CT) images and its influence on the obtained parameters. 1.1. Experimental data and relations for pressure drop

35

The Darcy-Forchheimer equation is a general relation pressure drop in porous media. Although obtained empirically by Forchheimer (1901), it can be shown that it emerges from a proper averaging of the Navier-Stokes equations (Whitaker, 1996) and it can actually be applied to various kinds of porous structures like packed beds of spheres (van der Sman, 2002), rings, or virtually any other structure (Rautenbach, 2009). In particular, for the porous bed composed of uniformly spaced spheres, the specific

40

Darcy-Forchheimer-type formula was proposed by Ergun (1952), (1 − ε)2 1−ε Δp = 150 3 2 μU + 1.75 3 ρU 2 , ΔL ε ds ε ds

(2)

1

INTRODUCTION

3

where ds is the sphere’s diameter. Foam-like structures are viewed as an inverse of bed-like structures where solid spheres are substituted by pores. Therefore, many researchers use the Ergun formula to fit their own experimental data. They usually replace sphere diameter, ds , with the eqivalent hydraulic diameter, dh , that should somehow be linked to the foam geometry. Additionally, varying numerical 45

constants preceeding the geometry-dependent coefficients are often used. This kind of general correlation was proposed by Gibson and Ashby (1999), (1 − ε)2 a2c (1 − ε) ac 2 Δp =α μU + β ρU , 3 ΔL ε ε3

(3)

where the characteristic dimension is the reciprocal of the specific surface area, ac . This formulation was used by Richardson et al. (2000), where a series of high-porosity (ε > 90%) alumina foams of ppi number in range from 10 to 65 were investigated. They reviewed a few specific relations between ac and 50

the mean pore diameter, dp , and chose the one proposed by Kozeny (1927) for packed beds of particles. The pore diameter was found from the image analysis of slices of the foam specimens. The coefficients α and β were fitted to the experimental data, and they became functions of porosity and mean pore diameter. Additionally, the influence of the surface roughness of the foam was investigated and some correlations for this effect were drawn as well. The relative surface roughness of a specially prepared

55

30ppi foam was investigated by means of the nitrogen sorption. The roughness factor was defined as the ratio of this measured area versus ac of a "smooth" foam geometry. The work of Moreira and Coury (2004) concerned a series of high-porosity foams of nominal linear pore density of 8ppi, 20ppi and 45ppi. The work explored pressure drop dependence on morphological parameters. The authors compared their data to previously existing correlations and proposed their own

60

formula. Additionally, the tortuosity of the medium was investigated. Incera Garrido et al. (2008) investigated experimentally a series of 10ppi, 20ppi, 30ppi and 45ppi ceramic foams with porosity in the range ε = 75% − 85%. The morphology of the foams, pressure drop as well as heat transfer coefficients were assessed. Additionally, another correlation for pressure drop was developed. The foam data presented in their work will serve as a benchmark to validate our results and

65

will be reviewed in detail in Section 3.1.2 and Section 3.3.2. Inayat et al. (2011b) focused on 10ppi, 20ppi and 30ppi foams of porosity about ε = 85% and drew their morphology correlations for α and β from Eq. 3. In their review, Edouard et al. (2008) remarked that it is always possible to come up with a set of

1

INTRODUCTION

4

coefficients to fit one’s experimental data and that many previously proposed correlations have only limited 70

application. In their other works (Lacroix et al., 2007, Huu et al., 2009) they proposed a relation that would be universal, namely the Ergun equation with the original numerical coefficients but a properly chosen equivalent particle diameter, dh , instead of the sphere diameter, ds . The equivalent particle diameter for a particular foam was obtained from the assumption that the equivalent packed bed of spheres and the foam have equal specific surface areas and equal porosities. The total volume of the representative unit

75

cell containing the sphere can be found from the sphere volume and the porosity of the packing,

Vtotal =

Vs . 1−ε

On the other hand, the specific surface area of the packing is ac =

(4) Ss Vtotal .

Thus it is possible to link

porosity, specific surface area and the diameter,

ds = 6

1−ε . ac

(5)

Since the foam has the same porosity and specific surface area, the relation in Eq. 5 holds for it as well. It also defines the foam’s equivalent hydraulic diameter, dh . The authors went one step further 80

and postulated the specific geometrical shape of a foam cell. First, a cell composed of three cylindrical struts was utilised (Lacroix et al., 2007), and this was later modified to a decahedron (Huu et al., 2009). Based on assumptions about the shape of the strut and detailed geometrical calculations, the equivalent hydraulic diameter, dh , and the foam strut diameter, dstrut , were linked in those works. In this work, however, we will exploit the presented formulation that involves only Eq. 5, without postulating any

85

specific shape of the foam. By combining Eq. 2 with Eq. 5 the following expression for the pressure gradient can be written, a2 ac Δp = 6.25 3c μU + 0.29 3 ρU 2 . ΔL ε ε

(6)

Another approach that aims to be universal was proposed by Dietrich (2012). He provided an extensive review of over 100 previously available pressure drop data sets and proposed his pressure-flow velocity relation in a non-dimensional form,

Hg = 110Re + 1.45Re2 ,

(7)

1

90

INTRODUCTION

5 d3

Δp h with non-dimensional pressure gradient (the Hagen number), Hg = ρ ΔL μ2 , as a function of the

Reynolds number, Re =

ρU dh εμ .

Eq. 7 translates to the following dimensional formula, ρU 2 μU Δp = 110 2 + 1.45 2 , ΔL εdh ε dh

(8)

In previous works, Dietrich et al. (2009) investigated a few available formulations for dh . They were based on semi-empirical correlations using three distinct approaches. First involved pore strut diameter and pore window size, next used only the specific surface area, while the third one assumed regular 95

geometrical models for foams. In our case, due to the lack of pore window size and strut diameters, we will only use the second idea and connect dh to specific surface area,

dh = 4

ε , ac

(9)

The relation in Eq. 9 comes from the so-called pipe-model (Richardson et al., 2000, Große et al., 2009). This model postulates that the foam has the same porosity and surface area as the set of hollow circular channels of diameter, dh , traversing the medium. Substitution of Eq. 9 into Eq. 8 yields the 100

following expression, a2 ac Δp = 6.88 3c μU + 0.36 3 ρU 2 , ΔL ε ε

(10)

As one can observe, the choice of the appropriate hydraulic diameter is somewhat arbitrary and different models were used. It should be noted that the relations for dh given by Eq.5 and Eq. 9 exhibit completely different behaviour. The former increases with increasing ε while the latter decreases. This issue was highlighted by Edouard et al. (2008) in their review, and it has a straightforward explanation. In 105

the first model, the characteristic dimension refers to the object that is solid (the sphere) so its diameter must decrease with the increase of ε while in the second it refers to size of the void space itself (the pipe diameter). Eventually both models yield relations for pressure that differ only in terms of numerical constants. It is also important to emphasise that, quite often, the pressure drop in porous structures is investigated

110

by using idealized tesselations (based on the so-called representative unit cell, RUC) or random packings in the form of the Voronoi cells (Wejrzanowski et al., 2013), rather than real geometries. The RUC may have a shape of a cube (Fourie and Du Plessis, 2002), Kelvin tetrakaidecahedron (Inayat et al., 2011a,c)

1

INTRODUCTION

6

or the Weaire-Phelan cell (Große et al., 2009). In particular, the Weaire-Phelan structure seemed to be very promising when representing the real foam-like structures. The average number of pore-to115

pore connections in this structure is 13.5 (in the Kelvin cell it is 12) while real foams exhibit 13.7 connections according to Matzke (1946) who studied soap bubbles. Detailed comparison of the WeairePhelan conglomerate to real structures revealed, however, that some of their properties significantly differ from the ones of real foams (Lautensack and Sych, 2008). For this reason, structural models based on Kelvin and Weaire-Phelan tesselations were also reported to be inaccurate (Große et al., 2009,

120

Incera Garrido et al., 2008). 1.2. Flow simulations in real foam structures The literature on the simulation of real, foam-like, high-porosity structures is not very rich. The analysis of real porous structures requires reconstruction of the geometry obtained from computer tomography scans, which can be intractable. This reconstruction results in a voxelised image with each cube

125

belonging to a fluid or solid space. Thus the reconstructed fluid-solid interface appears as a "staircase" geometry. Moreover, in the case of the classical computational fluid dynamics (CFD) tools such as the finite element method (FEM) or the finite volume method (FVM), additional preprocessing of the surface is required in order to generate the computational mesh. Some simulations performed using FVM were reported by Petrasch et al. (2008) for 10ppi ceramic

130

foam for solar receiver applications. The pressure drop in linear and inertial flow regimes for a 12.5million-cell tetrahedral mesh was investigated. FVM from commercial code ANSYS CFX was used. The Nusselt number variation with flow velocity and foam morphological characteristics was studied as well. Habisreuther et al. (2008) investgated real CT-obtained and artificial porous structures by means of a commercial FVM code. Skibiński et al. (2012) investigated pressure drop data of a 15mm×15mm×15mm

135

sample cut from a 10ppi foam geometry, using FVM. Additionally, he compared this data with pressure drops of a set of randomly generated Laguere-Voronoi tesselations of apparently similar characteristics to the original foam, and he reported significant discrepancies. Another case concerned research aimed at an electromagnetic improvement of the alloy filtration (Kennedy et al., 2013), using the COMSOL module based on FEM. Only 2D flow simulations were carried out, however. Ranut et al. (2014) performed a

140

full 3D simulation of combined heat and fluid flow across 10ppi, 20ppi and 30ppi metallic foams. Again, the simulations were carried out using a commercial FEM ANSYS software with extensive preprocessing of the geometry and mesh generation in the ANSYS ICEM CFD package. Lastly, simulations by de

1

INTRODUCTION

7

Carvalho et al. (2014) concerned flow in a metallic foam for aero-engine separators. The CT-reconstructed real geometry and artificial Weaire-Phelan structures were studied. Simulations were performed on a 145

tetrahedral mesh with use of the ANSYS Fluent package. The realizable k − ε turbulence model with enhanced wall treatment (EWT) was utilized. Alternatives to conventional CFD are available for the simulation of foams. The lattice Boltzmann method (LBM), which has received increasing attention over the past two decades and is used in this work, is particularily attractive (Succi, 2001). It is a method based on the mesoscopic kinetic description

150

of fluid flow. The algorithm consists of two steps, namely streaming and collision. These are performed on particle populations moving on a Cartesian grid. The LBM can be formulated from a range of collision models (BGK, MRT, ELBM) and particle velocity sets (D2Q9, D3Q19, D3Q27), these issues are covered in detail in Section 2.2.1. LBM is straightforward to implement and parallelize. Because of the use of the Cartesian grid, there is no need for geometry-fitted mesh generation. The geometry, after

155

designation of solid and liquid voxels, can be implemented directly in the "staircase" form. LBM has been shown to reproduce detailed flow characterictics of confined flows (Regulski and Szumbarski, 2012) and it has received significant attention both in single- and multi-component flow across porous media, see an extensive review by Liu et al. (2015). However, few studies refer to simulations of flows in real foam-like structures.

160

The first D3Q19 BGK-LBM simulation of flow across a SiC foam was presented by Bernsdorf et al. (1999). The size of the simulation domain was Lx × Ly × Lz = 100 × 149 × 149 and the flow regime was linear with a coarse specimen resolution of δx = 0.5mm. Still, the pore-scale flow was resolved - pore throats had the size of at least several voxels. Graf von den Schulenburg et al. (2007) investigated the flow in a compressed polyurethane foam with different levels of compression using a D3Q19 BGK-LBM

165

model on a 128 × 128 × 47 node geometry. They verified it against Magnetic Resonance Imaging (MRI) measurements of flow field in the experiment and a good agreement was found. Work by Magnico (2009) concerned flows across a CT-obtained Ni-Cr foam. Two cases were analyzed with FVM and LBM. The first case was the realistic geometry and second one was the same geometry after a transformation that introduced anisotropy. The tensorial form of the Darcy-Forchheimer relation was proposed there. In the

170

work of Gerbaux et al. (2010) D3Q19 MRT-LBM and FVM simulations as well as experiments were used to investigate three specimens of metallic foams. In this work, the domain size for LBM simulations was of order 250 × 120 × 120 but the flow regime was only linear. Very good agreement between both simulation

2

MATERIALS AND METHODS

8

methods were reported. Experiment and simulations agreed in two cases, whereas they significantly diverged in the third one, with this issue attributed to improper geometry reconstruction (large difference 175

of porosities between the real structure and the CT-reconstructed image was reported). Beugre et al. (2010) performed a complete simulation of the flow across a highly-porous (ε = 0.955) nickel foam using D3Q19 MRT-LBM formulation. The size of the computational domain was 4003 nodes. They provided experimental results as well and reported good agreement. Two-dimensional flow simulations using slices of the geometry obtained from 3D imaging were also reported. A ceramic foam filter was investigated by

180

Laé et al. (2013) for metal filtration applications with D2Q9 model. Finally, the LBM simulation of fluid flow with heat transfer was performed by Prasianakis et al. (2013) in the porous gas diffusion layer of a fuel cell. The authors analysed the permeability and the relative effective diffusivity using the state-ofthe-art D3Q27 thermal LBM model. The size of the computational domain was 200 × 100 × 180 nodes. Yet another recent example (Ettrich et al., 2014) presented an LBM simulation of a ceramic foam used

185

as catalyst support, with the specific feature of this work being the use of a special diffusive boundary condition for solid-liquid interface. 2. Materials and methods 2.1. Ceramic foams In this work, three foam specimens were investigated experimantally and numerically. These were

190

Al2 O3 VUKOPOR A ceramic foam filters with 10ppi, 20ppi and 30ppi nominal pore densities. The frontal images of the filters are given in Fig. 1. Each specimen was a cube of 50mm × 50mm × 50mm size. The nominal porosities were not provided by the supplier. Inspection of the filters confirmed that these structures were polyurethane-matrix based (manufactured by the so-called "replica-technique") with hollow spaces in the struts left after burning out the polymer template (Richardson et al., 2000).

195

Polyurethane template based filters are very common. Filters of this kind were also investigated by other groups (Incera Garrido et al., 2008, Große et al., 2009, Inayat et al., 2011a). The hollow strut spaces will be reffered to as the "internal porosity" or "internal void space", see Fig. 1 for details. The presence of internal porosity significantly influences some of the foam’s characteristics (see Sec. 2.1.3). The porosity of foams was calculated on the post-processed image after removal of the internal void

200

space in the struts. Obviously, only this porosity has any practical meaning when one wants to draw the correlations for pressure drop. The measurements of the real stuctures’ morphological characteristics

2

MATERIALS AND METHODS

(a) 10ppi

9

(b) 20ppi

(c) 30ppi

Figure 1: Frontal images of Al2 O3 filters used in the investigation. Each specimen has dimensions 50mm × 50mm × 50mm.

were not performed. Therefore neither total porosity nor mean strut diameter and window diameter are available. 2.1.1. Reconstruction of CT-images 205

Image acquisition of the foams was conducted using the high resolution SkyScan X-Radia XCT-400 tomograph. The samples were radiated with a directional microfocus GE|phoenix|x-ray xs|240D tube at an acceleration voltage of 150kV and a current of 50A. The radiography images were reconstructed using the XCT-Reconstruction software. This process resulted in a set of 256-level grayscale bitmap-format tomograms. In order to recover the actual geometry of the specimen, this data set had to be filtered

210

and binarised i.e. each pixel had to be prescribed with a value 0 or 1 (black or white) denoting void or solid phase. Both processes were performed using SkyScan cTAN software. First, a median filter was used to remove the noise. This filtering procedure has an additional advantage of smoothing the solidliquid boundaries. Then, an appropriate threshold level for binarisation was chosen. Finally, the internal porosity was closed. The last process, however, was only partially successful. There existed regions that

215

should be regarded as internal porosity, but they still had contact with external void space. This resulted either from imperfections in the reconstruction process or the existence of some narrow connections with the void space in reality (see Fig. 2 for details). An additional algorithm that would completely close those spaces had to be used, this issue is addressed in the next section. The presence of the internal porosity as well as the effect of using a median filter is clearly visible in Fig. 3.

220

It is important to note that there exists no standard procedure for the appropriate choice of the threshold for the gray-scale CT-image. Große et al. (2008) postulated that the threshold should be taken

2

MATERIALS AND METHODS

10

2 1 3

Figure 2: Classificaiton of porosity present in the post-processed image. (1) "open porosity" or "external void space" - void region that is filled with fluid and where the flow takes place, only this region is relevant from the hydrodynamic point of view, (2) "internal/closed porosity" - void space in the foam’s struts that has no connection with the external void space, (3) strut porosity that should be regarded as internal porosity (practically no flow occures here) but which topologically belongs to external void space either due to imperfect CT-scan reconstruction or by being connected by narrow channels in the real structure. The last kind of porosity is most difficult to remove from the post-processed image.

(a) Reconstructed

(b) Binarised, non-filtered

(c) Binarised, median-filtered

Figure 3: CT-image of 10ppi foam specimen at various stages of reconstruction: (a) 256-level grayscale image, (b) binarised image, (c) binarised image with median filter (used before the binarisation). Internal void spaces are clearly visible. They are artefacts of binarisation because the threshold value is low so that some volume that is definitely solid is not captured in the resulting picture, causing an exaggeration of the internal void space.

at the level when the rate of change of the specific surface area with respect to the change of threshold value,

dac dtr ,

reaches minimum. They, however, did not support this claim by any evidence. Therefore, in our

case, the thresholding was performed at the level that seemed to be the most appropriate representation 225

of the structure based on the ’naked-eye’ comparison with the grayscale image. Nevertheless, the influence of the threshold level on the structure morphology and pressure drop is discussed in Sections 3.1.3 and 3.2.1.

2

MATERIALS AND METHODS

11

2.1.2. Removing internal porosity in the CT-scans The foam structure contains closed porosity in the foam struts that has to be removed in order to 230

extract hydrodynamically relevant foam characteristics such as the surface area in contact with the flow. An initial idea for the removal was to additionally threshold the foams after performing numerical flow simulations. The dead ends of strut porosity would be erased based on the criterion of very small velocity there. This approach, however, had two drawbacks. First, it was very memory-consuming because the new geometry had to be obtained after post-processing of simulation data. Furthermore, it incorrectly

235

included areas in the separation zones on the external surface of the struts. Therefore this approach was rejected. A different kind of algorithm was then implemented. This procedure relied only on two-dimensional images without any need of flow simulations at all. The consecutive steps of the correction algorithm effects are demonstrated in Fig. 4. Fig. 4a shows a fragment of a sample slice of the 20ppi foam after

240

post-processing of a grayscale CT-scan. One can clearly see that there exist internal areas with no contact with the surrounding fluid as well as many regions that seem to belong to the internal porosity (one can deduce that from a preprocessed grayscale CT-image) but due to the imperfect post-processing of the scan remain unclosed. Yet, some of these areas can be in actual contact with the external area, but this issue is insignificant due to the fact that the connecting throat is very narrow and there is little flow in

245

there. This procedure requires the adjustment of an appropriate thickening level and the maximal size of the removed spot. These parameters are chosen by method of trial and error. The effects of corrections with various thickening levels are shown in Fig. 5. The procedure is performed along three perpendicular directions of the foam specimen. It must be noted that this approach is totally heuristic. On one hand, it

250

is very efficient because it operates on 2D images. On the other hand, it has a few drawbacks. It results in closing the ’tips’ of large pores that emerge as small entities within the solid regions. Additionally, it generates some artefacts on the solid surface and sometimes creates false connections (like the top left corner in Fig. 5d). The latter is removed by using an additional procedure not described here. It is worth underlining, however, that it is quite easy to detect those artefacts during the inspection of the

255

three-dimensional geometry. Eventually, however, they have to be removed manually.

2

MATERIALS AND METHODS

12

(a) basic scan

(b) thickened scan

(c) reversed scan

(d) small solid regions removed

(e) thickened reversed scan

(f) post-processed scan

Figure 4: Correction procedure shown on a fragment of slice of the 20ppi foam. (a) base scan, (b) morphological thickening of the solid part - the solid part starts to grow a specific number of pixels (shown in red), some of the open regions become closed, (c) reverse liquid and solid pixels with each other - obtain a negative image, (d) remove new ’hanging’ solid spots below specified size (those are in fact former inner porosities), (e) thicken new solid material - that effectively goes back to initial shape, (f) reverse solid with liquid to come back to positive image which is the final state.

2.1.3. Measurement of the specific surface area The specific surface area of the post-processed images had to be measured. Two distinct approaches were used. First, the Cauchy-Crofton theorem was explored (Petrasch et al., 2008, Große et al., 2009, Liu et al., 2010). It states that the surface area of the object can be calculated by generating a random set of 260

lines that cross the object and counting the number of intersections. More specifically, the following ratio is defined,

´ SΓ L nΓ (φ, ψ, h) dφdψdh , =´ Ω SΩ LΩ nΩ (φ, ψ, h) dφdψdh

(11)

in which SΓ and SΩ are the areas of two objects. The first, Γ, is the targeted geometry while the second, Ω, is the encapsulating ball. LΩ represents the set of all lines crossing Ω. These lines are parametrised with the set of variables (φ, ψ, h). The number of line-geometry cross sections nΩ (with the ball) and nΓ 265

(with the foam) are taken for each line. Obviously, for the ball nΩ ≡ 2 almost always. Eq. 11 is valid for the complete set LΩ which is infinite. Therefore a Monte-Carlo approximation is used and a random set of lines that cross that ball and internal geometry is generated, implying that the

2

MATERIALS AND METHODS

(a)

(c)

13

(b)

(d)

Figure 5: Correction of closed porosity in the 20ppi foam slices, showing 336x400 pixel images. Blue color represents the fluid pixels, green - solid, red - solid after correction. Thickening levels: (a) input image, (b) 4 pixels, (c) 6 pixels, (d) 8 pixels.

algorithm requires an appropriate amount of lines. We used a set of one and two million lines for each geometry and the resulting specific surface areas were very similar (discrepancy less than 0.1%). 270

Another approach to determine the specific surface area was also used. Lindblad (2005) demonstrated the method of the so-called local weighted configurations and proved its accuracy. The key idea of this technique was to come up with a proper surface for each possible configuration of 2 × 2 × 2 voxel cube consisting of either solid or void pixels. Results from both methods differ no more than by ±1% but only the results of the Cauchy-Crofton

2 275

MATERIALS AND METHODS

14

algorithm are presented, and they are listed in Table 1. 2.2. Simulations 2.2.1. Lattice Boltzmann method The simulations were performed using the lattice Boltzmann method (LBM), which is an approach to simulate fluid flow by repetetive streaming and collision of a certain set of discrete probability density

280

functions (often called ’populations’) fi (t, xj ) placed on a Cartesian grid (Chen and Doolen, 1998). The populations move (stream) with prescribed velocities. The velocity set is constructed in a way that this dynamic system possesses enough degrees of freedom to recover the governing equations of hydrodynamics. It usually consists of vectors that connect the neighboring nodes. The standard choice for two-dimensional problems is the D2Q9 lattice while for three-dimensional problems D3Q15, D3Q19 and D3Q27 sets are

285

used ("D" stands for dimension, "Q" for the number of velocity vectors). In our implementation, the set D3Q19 was used. The evolution equation, known as the lattice Boltzmann equation (LBE), relates the populations at neighboring grid nodes, fi (t + Δt, xj + ci Δt) − fi (t, xj ) = Ωi (f ) ,

(12)

where xj represents the point on the grid, ci is the lattice velocity and Ω (f ) is the collision operator. The macroscopic flow variables, namely the density and velocity, are the consecutive moments of the 290

discrete probability density functions, ρ=

N −1  i=0

fi ,

u=

−1 1 N ci fi . ρ i=0

(13)

The LBM has received ever more recognition in the CFD community due to its straightforward implementation and ease of handling complex geometries via the so-called bounce-back boundary condition (Pan et al., 2006). Originally, the LBM was derived from the lattice gas cellular automata (LGCA) which were formulated 295

to recover the behavior of gas systems (Higuera and Jiménez, 1989). Nowadays, LBM is viewed as a certain discrete form of the Boltzmann transport equation from the kinetic theory of gases. Using the Chapman-Enskog analysis it can be proved (Dellar, 2003, Viggen, 2009) that the LBE recovers the pseudo-incompressible Navier-Stokes equations in the limit of vanishing Mach number.

2

MATERIALS AND METHODS

15

The LBM algorithm has other advantages, namely there is no need to solve the Poisson equation for 300

pressure and the computational operations are local. However, the lack of a pressure solver translates to the lack of the incompressibility constraint. For this reason, simulations must be run with limited velocities and pressure in the fluid is recovered from the equation of state, p = ρc2s . The speed of sound on the standard lattice is cs =

√1 3

and the practical limit for flow velocity is |u| < 0.1 (that translates to

M a < 0.17). 305

Regarding the collsion operator, Ω (f ), its oldest and simplest form is the BGK collision model (Bhatnagar et al., 1954) of the form, 1 Ωi (f ) = − (fi − fieq ) , τ

(14)

where the population fieq represents the equilibrium distribution given by the truncated Maxwell-Boltzmann form,



fieq

= ρwi

u · ci (u · ci )2 u2 1+ 2 + − cs 2c4s 2c2s



,

(15)

and τ is the collision relaxation time that is strictly connected to the kinematic viscosity of the fluid, 

νLB = c2s τ −

310



1 , 2

(16)

and wi are the weights associated with specific lattice velocity vectors. It must be noted that the velocities and viscosity in the LBM simulation are non-dimensional and the physical values of velocity and viscosity follow from an appropriate rescaling. This issue is addressed in Section 2.2.1. MRT collision model . It can be shown that the BGK collision model used together with the bounce-back boundary condition results in spurious permeability-viscosity dependence (d’Humieres and Ginzburg,

315

2009). Additionally, the BGK collision operator is not the most stable collision operator - the stability drops with increasing τ (i.e. decreasing viscosity). In order to alleviate both effects, the MultipleRelaxation-Time model was proposed by d’Humieres et al. (2002). The key idea of the MRT model is to transform the set of populations, fi , into a set of moments using the transformation matrix, M (see Appendix A of d’Humieres et al. (2002) for the exact form of the matrix). The constructed moments

320

have various hydrodynamic and non-hydrodynamic interpretations: some are density, velocity and energy, while others lack relevant interpretation - they are the so-called ghost modes. Then the system is relaxed towards its equilibrium (also defined in the moment space) using a set of relaxation values,

2

MATERIALS AND METHODS

16

S = diag {s0 , s1, . . . , s18 }. After that, the distribution functions are recovered via the inverse transormation, M −1 . The MRT-lattice Boltzmann equation takes then the form,

fi (t + Δt, xj + ci Δt) − fi (t, xj ) = −(M −1 SM )i (f − f eq ) .

(17)

In our case, we use the so-called two-relaxation time model, where the relaxation rates are si =

325

ω,

i = 0, ..., 18 except s4 = s6 = s8 = 8 2−ω 8−ω . This formulation results in the fixed value of the so-called

magic number Λ =



1 2

−

1 ω1



1 2

−

1 ω2



=



1 2

−

1 ω



1 2

−

8−ω 8(2−ω)



=

3 16 .

This keeps the error of simulation

results completely viscosity-independent in the linear flow regime and decreases spurious dependence in other cases (Pan et al., 2006, d’Humieres and Ginzburg, 2009). It must be noted that there exist other 330

advanced collision models such as the Cascaded model (Geier et al., 2007) or the so-called Entropic LBM (Karlin et al., 2014). Both exhibit advantageous behaviour especially for high Reynolds numbers flows. Nevertheless, the MRT model performed satisfactorily in our case. Additionally, the mentioned advanced models have not been thoroughly analysed yet in combination with the bounce-back boundary condition. More complex boundary conditions dedicated to Entropic LBM have been proposed as well (Chikatamarla

335

and Karlin, 2013). Implementation of body force. There exist several methods to introduce the body-force term (e.g. gravity) into LBM simulations. An extensive review can be found in Huang et al. (2011). It is stated there that in the case of single-phase flows all methods give practically the same results. In this work, the formulation of Kupershtokh et al. (2009) known as the exact difference method (EDM) was used. The body-force

340

is introduced by adding appropriate difference of the equilibrium distribution functions. The LBE with force source term takes the form, fi (t + Δt, xj + ci Δt) − fi (t, xj ) = −(M −1 SM )i (f (u) − f eq (u)) + fieq (u + Δu) − fieq (u) .

(18)

Consecutive equilibrium distributions need to be calculated either at current velocity u or the one increased by Δu =

F Δt ρ .

Additionally, the actual flow velocity is obtained after introduction of a correction

term: ureal = u + 12 Δu. 345

Scaling from non-dimensional to physical units. Simulations are carried out in non-dimensional, or lattice, units. Non-dimensionality results from the normalisation of the populations and setting lattice cell size to unity. Thus there is a need to rescale the quantities from physical units to non-dimensonal and the

2

MATERIALS AND METHODS

17

other way round. The rescaling is performed based on the dynamic similarity of dimensional and nondimensional flows. Flows need to have the same Reynolds number,

Re = 350

ULB N UL = . ν νLB

(19)

In Eq. 19, U and ν denote velocity and kinematic viscosity in physical units, L is the physical size of some characteristic dimension of the system (e.g. size of the foam) while ULB and νLB are respective quantities in the lattice (i.e. non-dimensional) units. N is the characteristic dimensions given in voxels (it can be imposed as the scan resolution in the case of this work). The expression for physical size of the LBM cell is straightforward,

δx = 355

L . N

(20)

Lattice resolution, δx, together with the physical duration of one simulation step, δt, link the nondimensional and physical velocities, U = ULB

δx . δt

(21)

Using Eq. 19 - Eq. 21 one can show that the simulation time step is given by,

δt =

νLB (δx)2 . ν

(22)

The physical time step is fixed to the lattice resolution. This results from the fact that the postulated velocity set on the lattice is scaled together with the lattice resolution (populations must travel from node 360

to node in one time step). The resulting physical velocity of the flow is thus inversely proportional to the image resolution and lattice viscosity,

U = ULB

δx ν = ULB . δt δx · νLB

(23)

From this relation one sees the main limitation of LBM, namely the fact that the physical velocity is tied to the parameters of the simulation. In order to obtain bigger physical velocity in the simulation, we must either increase lattice velocity, or increase the resolution of specimen (decrease δx) or decrease the 365

lattice viscosity. The first is limited by the incompressibility constraint, the second necessitates significant memory consumption, whilst the third results in the loss of stability of the computation at some level of

2

MATERIALS AND METHODS

18

νLB . 2.2.2. Solver On the contemporary computational hardware the LBM would not be efficient without proper paral370

lelisation. Since this computational method relies on local operations and the computational data is very well-structured (the foam geometry resides in a prism), the parallelisation on a GPU-type architecture seems natural. A GPU processor can be considered as a vectorized unit that is capable of performing numerical operations faster than a CPU, provided that those operations are alike and that memory access patterns are regular. The code was written using NVIDIA-CUDA language that provides extensions to

375

the C language which enable the user to create functions (called kernels) that are executed on a GPU. Due to limited operating memory resources of a GPU (up to 6GB RAM on NVIDIA Tesla M2090 in our case) and a huge memory demand of LBM (compared to FEM of FVM), the parallelisation had to also cover the multi-GPU communication level. This was performed using the Message Passing Interface library. An efficient GPU implementation has to consider issues such as hierachy of memory types in the GPU

380

and memory access latency. In order to consider the above remarks, a special implementation technique was utilized. About 70% of the code was generated in an automatic way using a specially developed tools based on the R-language (Łaniewski-Wołłk, 2014). The calculations were performed on up to 24 NVIDIA Tesla M2090 GPUs divided into 3 racks. Code tests revealed very good scalability even for up to 48 GPUs (Łaniewski-Wołłk and Rokicki, 2015).

385

2.2.3. Simulation setup Pressure drop measurements can be carried out in various ways. Some approaches were analyzed by Narváez Salazar and Harting (2010). The simplest idea is to postulate a velocity field at the inlet and pressure value at the outlet. Then, in the course of the LBM simulation, the pressure at the inlet adjusts itself to the proper value. However, this configuration raises a few problems. First, it requires a certain

390

free space upstream and downstream of the specimen because, in the case of a porous structure, the velocity field at its inlet and pressure field at the outlet are not uniform. This results in an increased memory demand. Secondly, this configuration also introduces an additional error into the simulation because of the compressible nature of LBM - the pressure in front of the specimen is bigger, thus the density is bigger as well. While flowing across the porous structure, the fluid expands and increases its

395

velocity (obviously mass flow rate is kept constant). Therefore the question arises: Which velocity - inlet

2

MATERIALS AND METHODS

19

Figure 6: Boundary conditions used in the flow simulation: red - periodic inlet/outlet, gray - the channel walls (visualised only at the bottom and the back of the sample). Driving force was directed along the horizontal direction.

or outlet one is the "real" flow velocity? This effect becomes more pronounced when the pressure drop increases. The approach used in this work employed periodic boundary conditions on the inlet and outlet together with no-slip (bounce-back) boundary conditions at the sides of the sample (see Fig. 6). Thus, a channel 400

of a rectangular cross-section was created. The driving force was aligned with the channel. It is important to highlight that there is no free space at the inlet and the outlet. Obviously in the case of a real foam geometry, the foam shape at the inlet is non-conformal with the outlet shape, thus the velocity profile in this area is disturbed. Our tests revealed that the difference between the pressure drop measured in such a configuration and one with velocity inlet-pressure outlet is insignificant. A similar configuration

405

was used for example by Prasianakis et al. (2013). Actually, we noticed that the relative difference in both configurations dropped as we increased the average velocity of flows. Just to illustrate, we analyzed two times compressed 20ppi foam geometry. For the average flow velocity of Uin = 0.002 ms the relative difference reached in pressure drop was 5% while for Uin = 0.02 ms discrepancy dropped to less than 1%. The smaller discrepancy for flows with higher velocity is attributed to the increasing contribution of the

410

inertial drag occuring within the structure. The influence of the presence of walls is investigated in Section 3.2.2, because it has more significant effect on the pressure drop. 2.3. Experimental setup 2.3.1. Construction of the experimental stand

415

In order to investigate the flow resistance of the porous structures experimentally, a dedicated apparatus was constructed. It is presented in Fig. 7. The measurements were taken automatically - the

2

MATERIALS AND METHODS

20

Figure 7: Experimental apparatus consisting of a plexiglass channel, piping and auxiliary equipment. Everything is mounted on an aluminum frame. 9a

7

8

9b

10 5

11

13 14

6

12 1a

1b 2

3

4

Figure 8: Scheme of the experimental setup. The elements of the set-up: (1) Pumps. Maximum power Pmax = 90W each with three gears available, (2) water filter to protect the flowmeter, (3) flowmeter with the impulse output, (4) gauge system, (5) differential manometer Kobold DOM A20 with analog output, (6) water container, (7) honeycomb flow straightener with eyes of diameter 5mm and a net at its inlet with holes of size 0.75mm, (8) space for investigated structure, (9) auxiliary pneumatic containers compensating for rapid pressure variations, (10) main flow channel, (11) analog-digital converter National Instruments USB6009, (12) impulse counter, (13) computer, (14) inverter. The connections are marked with different colors: black - the plumbing, red - electric measurement circuit, brown - circuit powering the pumps, blue auxiliary water cables, green - USB cables.

control software adjusted the forcing (in the form of the power of the pump) and the response was measured (pressure drop and volumetric flow rate). The measurement output consisted of the set of points in the pressure-velocity space as well as the coefficients of the least-square fitted parabola representing the 420

Darcy-Forchheimer curve. The channel had a square cross-section of side length L = 50mm. The flow was in horizontal direction. The schematic view of the installation is shown in Fig. 8. Siemens Sitrans P DSII differential manometer was used. It permitted measurements with up to Δp = 6 000P a, with measurement error less than 0.075%. It was equipped with a display to check the

425

parameters. Data acquisition was realized with use of the current of intensity of 4− 20mA. The flowmeter

3

RESULTS AND DISCUSSION

21

Figure 9: Measurement space with a specimen. Tightening screws and plate are visible.

model was Kobold DOM A20. It was capable of measuring the flow rate in the range of 1 − 40l/min with error not greater than 0.5%. This translated to a velocity range of 0.0066 − 0.2667m/s. In our case, the velocity magnitude was limited to 0.2m/s. The main channel was made from transparent plexiglas. The total length of the chanel was 30L = 430

1.5m. A honeycomb flow straightener was placed 8L downstream from the inflow. After the next 8L the place for the investigated specimen was prepared. The measurement space was designed so that the channel cross-section was kept almost constant. The specimen was fixed in place by an additional elastic board (see Fig.9). 2.3.2. Measurement procedure

435

The control and data acquisition was realized in the LabView software. The created software enabled the control of the pumps and registration of the volumetric flow rate and pressure difference between the front and the back of the specimen. The measurement of pressure drop lasted one second. The impulse number (representing the flow rate) was gathered at this time. This procedure was then repeated and the average of two such measurements was taken as the representative value. Then, the control software

440

changed the pump power in order to change the flow rate.

3

RESULTS AND DISCUSSION

22 2

2

No.

ppi

Lx × Ly × Lz [−]

resolution [μm]

Lx × Ly × Lz [mm]

εpre [−]

εpost [−]

apre [ m ] m3

apost [ m ] m3

1. 2. 3.

10 20 30

672x672x800 672x672x800 672x672x800

45.252 45.252 25

30.41x30.41x36.20 30.41x30.41x36.20 16.8x16.8x20

0.805 0.768 0.797

0.777 0.751 0.787

1039.7 1337.9 1666.9

721.5 1075.0 1361.7

Table 1: Parameters of foam images before and after post-processing. Values of specific surface area ac are based on the Cauchy-Crofton method with two million generated lines. The discrepancy between measurements with Lindblad and Cauchy-Crofton methods are approx. 1%. Similar discrepancies as well as absolute accuracy of 1% are reported in Liu et al. (2010). Therefore the error of the measurement is at the level of 2%. The number of significant digits for presented numbers (both for  and ac ) follows the convention used in other works (eg. Incera Garrido et al. (2008).)

(a) 10ppi

(b) 20ppi

(c) 30ppi

Figure 10: Rendered foam specimens used in the simulations. Note that 30ppi was scanned with higher resolution.

3. Results and discussion 3.1. Foam morphological parameters 3.1.1. Parameters of reconstructed images Complete foam data is given in Table 1. Data presented there does not refer to the entire foam because 445

the scanning process could not cover the entire foam specimen. The reasons are the following: Firstly, serious disturbances occurred on the borders of the measured specimen. Secondly, the total scanned volume was related to the scanning resolution. We were able to read about 1000 pixels in each direction. After removing external layers full of disturbances and cropping the geometry to multiplicity of 32 (due to GPU memory alignment constraints), we arrived at final specimen sizes of 672 × 672 × 800 voxels. 10ppi

450

and 20ppi foam scans were performed with resolution of about 45μm while the 30ppi foam was scanned at 25μm. The last choice was dictated by the fact that we did not want to risk the quality of the image, as the scanner may miss the thinnest struts and the foam topology would be distorted. One can observe that the corrections of void space did not affect the porosity much, but had a tremendous impact on the specific surface area (20%-30% reduction). Rendered foam specimens are shown in Fig. 10.

3

RESULTS AND DISCUSSION

23 2

No.

ppi

εnom

εtot

εopen

ac [ m ] m3

Porosity measurement

1. 2. 3. 4. 5. 6. 7. 8. 9. 10. 11. 12.

10 20 30 30 10 20 30 20 20 10 30 30

0.8 0.8 0.8 0.85 0.8 0.8 0.8 0.8 0.8 0.8 0.8 0.85

0.818 0.804 0.816 0.852 0.808 0.802 0.806 0.741-0.781 0.721-0.775 0.794-0.815 0.816-0.808 0.870-0.842

0.752 0.751 0.766 0.807 0.765 0.748 0.752 0.736-0.757 0.718-0.745 0.773 0.758 0.793

675.4 1187 1437.8 1422 664 1204 1402 1187-1244 1213-1268 664 1474 1520

Hg-intrusion Hg-intrusion Hg-intrusion Hg-intrusion Hg-intrusion Hg-intrusion Hg-intrusion CT-post-processing CT-post-processing Gravimetry, Hg-intrusion Gravimetry, Hg-intrusion Gravimetry, Hg-intrusion

Source Incera Garrido et al. (2008)

Dietrich et al. (2009) Große et al. (2008) Große et al. (2009)

Table 2: Parameters of similar foams found in the literature. Both open and closed porosity were determined by Hg-intrusion except Große et al. (2008) where CT-analysis was used. Specific surface area was determined using the Cauchy-Crofton formula (Eq. 11) after reconstruction of the CT-images.

455

3.1.2. Comparison with other data Similar foam specimens were investigated by other groups. Each of the following cases concerns Al2 O3 ceramic sponges and the selected specimen data is presented in Table 2. Große et al. (2008) provided extensive experimental data for a series of 20ppi foams of nominal porosity (i.e. given by the supplier) of εnom = 80%. Data was obtained by two teams after post-processing of CT-images. Discrepancies are

460

attributed to different preprocessing of the specimen and insufficient resolution. They are considerable in the case of porosity. In their other paper, Große et al. (2009), the morphological data of series of 10ppi, 20ppi and 30ppi foams were given. In that case, Hg-intrusion porosimetry was used to obtain both open and total porosity (penetration abilities of mercury are dependent on external pressure applied). Additionally, the gravimetric method was utilised, which uses the weight of the sample, its external

465

dimensions and material density. The latter was obtained by helium pycnometry. Their third work, Dietrich et al. (2009), focused on investigation of the series of foams of the same kind as before. It is not clear if it concerned exactly the same specimens or ones very similar to them. Some data presented there are exactly the same (ac ), while other parameters are very close to each other. Incera Garrido et al. (2008) investigated the same ppi series. Foams of various porosities were taken into account. Two last

470

works also provided extensive pressure drop data that will be compared againts our results in Section 3.3.2. One can immediately conclude that foam specimen no. 2, 10 and 12 match porosities of our foam 2

while the 20ppi and 30ppi specimens have closely (±0.5%), but the 10ppi foam has lower ac by 100 m m3 bigger ac . On the other hand, presented results also demonstrate that the range of discrepancies when 475

using experimental methods can be significant (porosity errors of up to 5% for specimen no. 9).

3

RESULTS AND DISCUSSION

24 2

No.

Threshold case

ε

ac [ m ] m3

1. 2. 3.

low medium high

0.7914 0.7769 0.7658

696.4 721.5 746.3

20000

Table 3: Parameters of the 10ppi foam obtained at different thresholding levels.

10000 0

5000

Δp P a ΔL [ m ]

15000

ε = 79.1% ε = 77.7% ε = 76.6% DF-fit

0.00

0.05

0.10

0.15

0.20

Uavg [ m s]

Figure 11: Pressure drop for 10ppi foam obtained from LBM simulations with the same base CT-image with different thresholding levels. DF-fit denotes Darcy-Forchheimer curves that were least-square fitted to the data from the LBM simulations.

3.1.3. Influence of threshold value on foam parameters The influence of the thresholding level on the foam parameters was investigated in the case of the 10ppi foam. Thresholding was performed after appropriate filtering processes and closing of internal void spaces. Threshold levels are referred to as "low", "medium" and "high". Obtained data is shown in Tab. 3. 480

Porosity varies relatively about 3.3% while the specific surface area ac about 7%. Obviously it affects the predicted pressure drop when typical Darcy-Forchheimer-type correlations are used and it also distorts input geometries for the numerical simulations. The scale of the potential error in the latter case will be investigated in Section 3.2.1.

3

RESULTS AND DISCUSSION

25

3.2. Simulation results 485

3.2.1. Influence of the threshold on the pressure drop Fig. 11 shows the Darcy-Forchheimer curves obtained from simulations of the 10ppi specimen with different threshold levels. From Section 3.1.3 we know that the relative differences in porosity and specific surface area are 3.3% and 7% respectively. The simulations showed that these differences translate to the relative error between lowest- and highest-porosity foam at a level of more than 20%. This error was

490

m measured for the following velocities: Uinlet = 0.001 ms , 0.1 m s , 0.2 s . The obtained results emphasize the

significance of careful image analysis and reconstruction. 3.2.2. Influence of channel walls The measurement of pressure drop for the foam specimen confined in the channel is obviously affected by the presence of channel walls. These generate additional drag in two ways. Firstly, they introduce 495

another surface in contact with the liquid, and secondly, their presence creates an additional boundary layer so that more fluid passes through central part of the foam rather than the region closer to the sides of the channel. The narrower the channel is, the greater the relative influence on the obtained pressure drop. As discussed earlier, one of the key differences between the experimental and simulation foam configuration was the fact that the CT-obtained geometries for simulations did not encompass the entire

500

physical specimen. Thus, we would expect that the presence of channel walls introduces relatively more drag to the specimen in the simulated case than in the experimental one. The degree to which the foam confinement increased the pressure drop was investigated for each foam seperately for multiple values of driving force. The reference case was the configuration where no side walls were used. Instead, periodic boundary conditions on the sides were implemented. Since the flow was

505

driven by the same external force field for both configurations, the difference is observed in the average velocity. The discrepancy is defined as, =

open − U wall Uavg avg , wall Uavg

(24)

wall and U open are the mean flow velocities with and without the channel walls, respectively. in which Uavg avg

The error values are presented in Table 4. The relative difference in each case is at the level of several percent. It is the biggest in the case of the 30ppi foam probably because it is the smallest physical 510

specimen. Interestingly, this difference is almost constant across all values of driving force. Nevertheless,

3

RESULTS AND DISCUSSION No. 1. 2. 3.

ppi 10 20 30

26

avg. rel. difference 4.41% 7.00% 8.89%

rel. difference standard deviation σ 0.40% 0.13% 0.17%

Table 4: Mean velocity difference (given by Eq. 24) between unconfined and confined foam configuration. For the whole range of driving force, in each case the difference is practically constant. The simulations were performed for non-dimensional −7 −6 −6 force values of g = 4 · 10 , 1 · 10 , 3 · 10 , 5 · 10−6 for all specimens and additional g = 7 · 10−6 in the case of 20ppi foam. That translates to all but the smallest and biggest pressure drops reported in Fig. 12 for the 10ppi and 30ppi foams and all but the smallest one for the 20ppi foam. Direction x y z average

10ppi k1 [m2 ] × 10−9 k2 [m] × 10−3 122.86 1.73 171.83 2.51 165.16 2.81 149.90 2.25

20ppi k1 [m2 ] × 10−9 k2 [m] × 10−3 81.37 1.41 77.64 1.20 71.48 0.95 76.61 1.16

30ppi k1 [m2 ] × 10−9 k2 [m] × 10−3 45.28 0.81 49.60 0.93 48.14 0.83 47.61 0.85

Table 5: Darcy-Forchheimer coefficients for the foams based on the LBM simulations. Water parameters at 20◦ C: μ = kg 1.002 · 10−3 P a · s, ρ = 998.3 m 3 . The average value of the coefficient is the harmonic average because the coefficients appear in the denominator of the DF-formulae (Eq. 1). The accuracy of the presented coefficients follows the approach presented in other works on this subject (Incera Garrido et al., 2008).

we observe that it is not as significant as the error introduced by improper image preprocessing. 3.2.3. Simulation results for all three foams All three foams were investigated with use of the LBM solver and all three perpendicular flow directions were taken into account. The simulation results together with the Darcy-Forchheimer curves obtained for 515

them are presented in Fig. 12. The Darcy-Forchheimer curves were obtained by the least-square fit. In Section 3.3.1, the pressure drops is compared against experimental values. The latter were obtained for velocities up to U = 0.2 m s . In the simulations, however, such speeds were impossible to achieve due to the speed and stability limitations of the LBM discussed in Section 2.2.1. Flow across the 30ppi foam, due to a better scan resolution, resulted in the highest physical velocity.

520

One can clearly observe the anisotropy of the foam structure. This is a common issue that was previously reported in other works (Große et al., 2008, 2009) and can be explained by the features of the formation process of the foam specimen (Scheffler and Colombo, 2006). In contrast, in the case of metallic foams it was reported that the pressure drop is lower in one direction (Ranut et al., 2014). The coefficients k1 and k2 are provided in Table 5.

RESULTS AND DISCUSSION

40000 30000

Sim. x Sim. y Sim. z DF-fit

10000

10000

Δp P a ΔL [ m ]

15000

Sim. x Sim. y Sim. z DF-fit

0

0

5000

Δp P a ΔL [ m ]

27

20000

20000

3

0.00

0.05

0.10

0.15

0.20

0.00

0.05

0.10

Uavg [ m s]

Uavg [ m s]

(a) 10ppi

(b) 20ppi

0.15

0.20

20000 0

10000

Δp P a ΔL [ m ]

30000

40000

Sim. x Sim. y Sim. z DF-fit

0.00

0.05

0.10

0.15

0.20

Uavg [ m s]

(c) 30ppi Figure 12: Simulation results of pressure gradients for 10ppi, 20ppi and 30ppi foams. Sim. represent values obtained in the LBM simulations, DF-fit stands for the Darcy-Forchheimer curves that were fitted to this data. 525

3.3. Comparison of simulation results to other data 3.3.1. Comparison to experimental results Data shown in Fig. 13 shows very good agreement between results obtained from numerical simulations and experiments for 10ppi and 20ppi foams. The relative discrepancy is defined in the following way,

3

RESULTS AND DISCUSSION

28

=

ΔpSIM − ΔpEXP . ΔpEXP

(25)

In the case of velocities Uavg > 0.025 ms for 10ppi foam, the difference did not exceed 10% for x530

direction and 15% in two other directions. For these latter directions there was a tendency for the discrepancy to drop with the increase of the flow velocity ( < 10% for Uavg > 0.1 m s ). In the case of 20ppi foam, the agreement was even better. For the x-direction the discrepancy grew from = 3% for m 0.05 ms > Uavg > 0.02 m s to = 8% at Uavg = 0.18 s . The y- and z-directions also showed the deviation

well below = 10% with strong tendency for convergence for higher velocities ( = 4% for Uavg = 0.1 m s 535

m and = 1% for Uavg = 0.18 m s ). Velocities smaller than Uavg = 0.025 s were neglected in the estimation

because the experimental results showed very high oscillation. These oscillations occurred because the pump driving the flow in the experiment was working with minimal power and its output was very unstable in that regime. It should be noted that in the case of 10ppi and 20ppi foams a very good agreement was observed 540

despite three important differences between the experimental and numerical configuration. Firstly, numerical data for 10ppi and 20ppi foams were obtained for flow velocities only up to 50% of the maximal flow speed in the experimental case. Furthermore, the boundary conditions were different (as discussed in Section 2.2.3). Finally, the foam geometries used in the simulations did not represent the complete volume of the real foams used in the experiment. Thus, these results confirm that the inlet/outlet effects

545

are negligible. Additionally, the specimen size turned out to be big enough to diminish the influence of the presence of the channel walls and can be treated as a representative volume of the whole foam. In the case of the 30ppi foam, the agreement was not as good as in the previous cases, with results for only one flow direction agreeing quite well. For the x-direction the discrepancy was in the range of = 10% while for the two other directions deviation exceeded = 20%. These results suggest that the sample

550

size was not representative for the whole specimen. Another difference between the simulation and the experiment that is more significant than in the case of 10ppi and 20ppi specimen was the influence of the channel walls, because a relatively smaller portion of the specimen is confined in the channel in the case of the 30ppi foam. Additionally, inspection of the geometry of this foam revealed significant heterogeneity within the porosity distribution of the structure. This heterogeneity can be attributed to the phenomena

555

occuring during the forming process and it was also reported by Große et al. (2009). Maximal velocities obtained in this simulation reached 60% of the maximal values reached in the experimental case.

RESULTS AND DISCUSSION

20000

3

10000

x y z x y z

0

5000

Δp P a ΔL [ m ]

15000

Exp. Exp. Exp. Sim. Sim. Sim.

29

0.00

0.05

0.10

0.15

0.20

Uavg [ m s]

Exp. Exp. Exp. Sim. Sim. Sim.

20000

Δp P a ΔL [ m ]

30000

40000

x y z x y z

x y z x y z

0

0

10000

10000

Δp P a ΔL [ m ]

30000

Exp. Exp. Exp. Sim. Sim. Sim.

20000

40000

(a) 10ppi

0.00

0.05

0.10

0.15

0.20

0.00

0.05

0.10

Uavg [ m s]

Uavg [ m s]

(b) 20ppi

(c) 30ppi

0.15

0.20

Figure 13: Comparison between numerical and experimental data for 10ppi, 20ppi and 30ppi foams. Exp. stands for experimental results, Sim. denotes the Darcy-Forchheimer curves fitted to the results from the LBM simulations.

3.3.2. Comparison to other results In this section data from the LBM simulations is compared with the experimental data of pressure drop from other research groups (Incera Garrido et al., 2008, Dietrich et al., 2009). Darcy-Forchheimer curves

3 560

RESULTS AND DISCUSSION

30

for foams most similar to the foams used in our research are plotted in Fig. 14. The Darcy-Forchheimer curves corresponding to lowest and highest permeability from our simulations are also presented. The values of the coefficients used for plotting are given in Table 6. In the case of the 10ppi foam shown in Fig. 14a, one can observe that porosity values are very close to each other but the specific surface areas of the experimental specimens are smaller. Yet, the experimental

565

curves fit well into the range spanned by our data for higher velocity range. However, due to the smaller specific surface area in the experimental specimens, one would expect smaller pressure drop than in our case. Also in the case of the 20ppi foam (Fig. 14b) we have the specimens with characteristics very close to our structure. Both reference results possess higher specific surface areas. Here, however, only data

570

reported by Dietrich et al. (2009) fits well into our range, while data of Incera Garrido et al. (2008) presents a larger pressure drop in the investigated range. One can easily notice that Darcy-Forchheimer (DF) curves from our results have a greater growth rate and will eventually overtake both results, despite lower specific surface area of our specimen. The results for the 30ppi foam are shown in Fig. 14c. The specimen of Incera Garrido et al. (2008)

575

that has lower porosity and higher specific surface area (no. 3 in Table 2) exhibits a higher flow resistance, but it can be seen, that its DF curve has slower growth rate. The more porous structure investigated by this group (no. 4 in Table 2) initially remains within our range, but also exhibits smaller growth rate later on. Dietrich et al. (2009) result practically overlaps with the lower bound of our case. To summarise, most of the presented experimental data is in a fairly good agreement with our results.

580

Still, one should recognize the limitations of the results of the previous works. First, both groups report the anisotropy of ceramic sponges in their works. Neither of them, however, investigates its influence in their measurements. Our results shown in Section 3.2.3 demonstrate that the directional differences can be significant. In addition, we must remember that the porosity measurements in other groups were performed on real samples. In the work of Große et al. (2008) it was reported that variations of porosity

585

measurements for the same sample were on the level of 2-3% based on the preprocessing of acquired data. Finally, Incera Garrido et al. (2008) used tomography resolutions at the level of 50 − 86μm. From the authors’ experience 50μm is an upper resolution bound for foams of this kind. It seems that the insufficient CT-resoultion can also be a source of error, especially in case of specific surface area estimation.

4

CONCLUSIONS AND FUTURE WORK Source Dietrich et al. (2009) Incera Garrido et al. (2008) Incera Garrido et al. (2008)

10ppi k1 [m2 ] × 10−9 k2 [m] × 10−3 77 1.87 28.58 3.13 -

31 20ppi k1 [m2 ] × 10−9 k2 [m] × 10−3 54 1.14 9.17 1.67 -

30ppi k1 [m2 ] × 10−9 32 7.23 11.07

×10−3 0.98 1.45 1.89

Table 6: Darcy-Forchheimer coefficients for the foams from other groups.

3.3.3. Comparison with existing correlations 590

The results obtained in our LBM simulations can be compared with known correlations which relate pressure drop with specific surface area. For this we used two relations. The first is given by Eq. 10 by Dietrich (2012) and the second was based on the formulation proposed in Lacroix et al. (2007) and Huu et al. (2009) given by Eq. 6. The results are presented in Fig. 15. Due to the structure anisotropy, data for all three directions is

595

drawn. One can notice that in all cases the correlations of both groups fit well to our results for 10ppi and 20ppi. Obviously, the correlation by Dietrich gives a higher pressure drop estimate. In the case of the 30ppi foam, Dietrich’s correlation fits just below the lowest DF curve while the data of Huu, Lacroix and co-workers (designated in Fig. 15 as ’HLE’) does not fit at all. Although, for the reasons discussed earlier, geometry of this specimen is not considered as reliable as data for 10ppi and 20ppi foams, we can

600

still discuss the applicability of the pressure drop correlations in this case because accurate parameters of the digital 30ppi specimen are available. It should be stressed that verification of correlations agains simulation results should in principle be prone to lower discrepancy than their verification to experimental data.The correlations involve parameters that were obtained by investigating the digital images of foams and exactly these images were used in the

605

simulations. Thus the inaccuracies are only related to determination of the structures’ parameters based on their voxelized representation. 4. Conclusions and future work This work presents the results of numerical and experimental investigations of pressure drop in flow across three specimens of ceramic foams. Very good agreement for 10ppi and 20ppi foams is observed.

610

Results for the 30ppi foam, probably due to the smaller and heterogeneous (thus not representative) physical volume investigated, do not show such a good agreement. Still, they can be considered satisfactory. Results confirm the capability of the lattice Boltzmann method to resolve the flow field in porous media at pore-scale level.

CONCLUSIONS AND FUTURE WORK

32

TW ε = 75.1% ac = 1075 m−1 G’08 (2) ε = 75.1% ac = 1187 m−1 D’09 (6) ε = 74.8% ac = 1204 m−1

10000 0

50000

0.00

0.05

0.10

0.15

0.20

0.00

0.05

0.10

Uavg [ m s]

Uavg [ m s]

(a) 10ppi

(b) 20ppi

0.15

0.20

30000

TW ε = 78.7% ac = 1361.7 m−1 G’08 (3) ε = 76.6% ac = 1437.8 m−1 G’08 (4) ε = 80.7% ac = 1422 m−1 D’09 (5) ε = 75.2% ac = 1402 m−1

0

10000

Δp P a ΔL [ m ]

20000

Δp P a ΔL [ m ]

10000 0

5000

Δp P a ΔL [ m ]

30000

40000

TW ε = 77.7% ac = 721.5 m−1 G’08 (1) ε = 75.2% ac = 675.4 m−1 D’09 (5) ε = 76.5% ac = 664 m−1

20000

4

0.00

0.05

0.10

0.15

0.20

Uavg [ m s]

(c) 30ppi Figure 14: Comparison of the LBM simulation results with the experimental results of other groups estimated for water kg parameters at 20◦ C: μ = 1.002 · 10−3 P a · s, ρ = 998.3 m 3 . The numbers in the brackets refer to consecutive specimens from Table 2, abreviations: TW - this work, D’09 - Dietrich et al. (2009), G’08 - Incera Garrido et al. (2008).

CONCLUSIONS AND FUTURE WORK

20000

4

10000

15000

Sim. x Sim. y Sim. z Eq.(10) (Dietrich) Eq.(6) (HLE)

0

5000

Δp P a ΔL [ m ]

33

0.00

0.05

0.10

0.15

0.20

Uavg [ m s]

40000

Sim. x Sim. y Sim. z Eq.(10) (Dietrich) Eq.(6) (HLE)

30000 Δp P a ΔL [ m ]

20000

0

0

10000

10000

Δp P a ΔL [ m ]

30000

Sim. x Sim. y Sim. z Eq.(10) (Dietrich) Eq.(6) (HLE)

20000

40000

(a) 10ppi

0.00

0.05

0.10

0.15

0.20

0.00

0.05

0.10

Uavg [ m s]

Uavg [ m s]

(b) 20ppi

(c) 30ppi

0.15

0.20

Figure 15: Comparison of pressure drop from LBM simulations against other groups’ correlations. Eq. 10 is based on Dietrich (2012) while Eq. 6 is based on the idea of Lacroix et al. (2007) and Huu et al. (2009) (referred to as HLE in the legend).

The influence of different factors on pressure drop was investigated. Those factors included image 615

thresholding levels and the presence of channel walls. The importance of proper image segmentation

4

CONCLUSIONS AND FUTURE WORK

34

was demonstrated. The obtained results were also compared against other researchers’ data for similar structures as well as to selected correlations for pressure drop. Good agreement was found in comparison with similar structures from the other groups. Also, the pressure drop correlations of Huu et al. (2009) and Dietrich (2012) agreed very well with presented data in most cases. In all cases, foam anisotropy was 620

detected. The obtained results imply further research areas. First, the analysis of pressure drop with respect to the CT-scan resolution and total specimen size (especially for the 30ppi foam) will be performed. The ultimate goal of our work is to define the pressure-drop correlations that will be based purely on the geometrical characteristics of these type of foams and image analysis is a key issue in this case. As soon

625

as this capacity is developed, the detailed internal characteristics will be investigated. Additionally, the real geometry will be manipulated in order to modify the degree of anisotropy and come up with data for more general correlations. Nomenclature Latin letters

630

ac - specific surface area



m2 m3



ci - lattice velocity vector associated with i-th distribution function (−) cs - lattice speed of sound (−) dh - equivalent hydraulic diameter (m) ds - sphere diameter (m) 635

dstrut - foam strut diameter (m) dp - mean pore diameter (m) fi - particle distribution function in the LBM simulations (−) f - vector of all distribution functions in the LBM simulations (−) F - force vector in the LBM simulations (−)

640

h - 3rd parameter used in random line generation for the Cauchy-Crofton formula(−) Hg - Hagen number(−)

k1 - viscous permeability coefficient m2 k2 - inertial permeability coefficient (m) L - physical length (m)



4 645

CONCLUSIONS AND FUTURE WORK

35

LΩ - set of all lines crossing the encapsulationg ball Ω in the Cauchy-Crofton formula ΔL - foam length (m) M - MRT transformation matrix (−) M a - Mach number(−) N - size of the object in lattice units (−)

650

nΓ,Ω - number of intersections of generated lines and object of interest used in the Cauchy-Crofton formula (−) p - lattice pressure (−) Δp - pressure drop in flows across foams (P a) P - power (W )

655

Re - Reynolds number (−) si - relaxation rates in the MRT collision operator (−)

S - surface area m2



S - diagonal matrix of the MRT relaxation rates (−) t - lattice time (−) 660

tr - threshold level (−) δt - physical simulation time step (s) Δt - lattice time step (−) u - lattice fluid velocity vector (−) Δu - lattice fluid velocity increase due to body force (−)

665

U - mean flow velocity

m

s

ULB - lattice flow velocity to rescale(−)

V - volume m3



wi - weights for the equilibrium distribution function (−) xi - position of the lattice node (−) 670

δx - physical lattice resolution (m) Greek letters α - coefficient in the viscous term of the generalized Ergun equation (−) β - coefficient in the inertial term of the generalized Ergun equation (−) Γ - index of the internal geometry for the Cauchy-Crofton formula

4 675

CONCLUSIONS AND FUTURE WORK

36

- relative error (−) ε - porosity (−) Λ - magick number linking the relaxation rates in the TRT collision operator(−) μ - fluid dynamic viscosity (P a · s) ν - fluid kinematic viscosity

680



m2 s



νLB - fluid kinematic viscosity in lattice units(−) ρ - fluid density



kg m3



ρ - lattice fluid density (−) σ - standard deviation τ - lattice relaxation time (−) 685

φ - 1st parameter used in random line generation for the Cauchy-Crofton formula (−) ψ - 2nd parameter used in random line generation for the Cauchy-Crofton formula (−) ω - relaxation rate in the LBM collision operators (−) Ω - index of the bounding sphere for the Cauchy-Crofton formula Ω (f ) - LBM collision operator (−)

690

Acknowledgements The work was financially supported by the Polish Ministry of Science and Higher Education research grant no. N N507 273636. This work was also supported by the European Union in the framework of European Social Fund through the “Didactic Development Program of the Faculty of Power and Aeronautical Engineering of the Warsaw University of Technology”. Finally, this research was supported

695

in part by PL-Grid Infrastructure. The calculations were performed on ZEUS HPC unit at the Academic Computer Centre CYFRONET in Cracow. The author would also like to thank Dr. Christopher R. Leonardi for his proof-reading of the manuscript and his comments. References Bernsdorf, J., Delhopital, F., Brenner, G., and Durst, F.

700

porous media using the lattice Boltzmann method.

Prediction of pressure losses in

In Keil, F., Mackens, W., Voss, H., and

Werther, J., editors, Scientific Computing in Chemical Engineering II, pages 336–343. Springer Berlin Heidelberg, 1999.

ISBN 978-3-642-64295-1.

http://dx.doi.org/10.1007/978-3-642-60185-9_39.

doi:

10.1007/978-3-642-60185-9_39.

URL

4

CONCLUSIONS AND FUTURE WORK

37

Beugre, D., Calvo, S., Dethier, G., Crine, M., Toye, D., and Marchot, P. 705

flow simulations on a metallic foam.

Lattice Boltzmann 3D

Journal of Computational and Applied Mathematics, 234

(7):2128 – 2134, 2010. ISSN 0377-0427.

doi: http://dx.doi.org/10.1016/j.cam.2009.08.100.

URL

http://www.sciencedirect.com/science/article/pii/S0377042709005482. Fourth International Conference on Advanced Computational Methods in Engineering (ACOMEN 2008). Bhatnagar, P. L., Gross, E. P., and Krook, M. A model for collision processes in gases. i. small amplitude 710

processes in charged and neutral one-component systems. Phys. Rev., 94:511–525, May 1954. doi: 10.1103/PhysRev.94.511. URL http://link.aps.org/doi/10.1103/PhysRev.94.511. de Carvalho, T. P., Morvan, H. P., and Hargreaves, D.

Pore-level numerical simulation of open-

cell metal foams with application to aero engine separators. Expo 2014: 715

In Proceedings of ASME Turbo

Turbine Technical Conference and Exposition, Duesseldorf, Germany, 2014.

URL

http://proceedings.asmedigitalcollection.asme.org/proceeding.aspx?articleid=1908195. Chen, S. and Doolen,

G. D.

Lattice Boltzmann method for fluid flows.

view of Fluid Mechanics, 30(1):329–364, 1998.

doi:

Annual Re-

10.1146/annurev.fluid.30.1.329.

URL

http://dx.doi.org/10.1146/annurev.fluid.30.1.329. Chikatamarla, S. and Karlin, I. 720

lations:

Boundary conditions.

1925 – 1930, 2013.

Entropic lattice Boltzmann method for turbulent flow simuPhysica A: Statistical Mechanics and its Applications, 392(9):

ISSN 0378-4371.

doi: http://dx.doi.org/10.1016/j.physa.2012.12.034.

URL

http://www.sciencedirect.com/science/article/pii/S0378437113000113. Dellar, ple 725

ISSN

P.

J.

relaxation

Incompressible times.

0021-9991.

limits

Journal doi:

of

of

lattice

Computational

boltzmann Physics,

equations 190(2):351

using

multi-

370,

2003.

–

http://dx.doi.org/10.1016/S0021-9991(03)00279-1.

URL

http://www.sciencedirect.com/science/article/pii/S0021999103002791. Dietrich, B. Pressure drop correlation for ceramic and metal sponges. Chemical Engineering Science, 74(0):192 – 199, 2012.

ISSN 0009-2509.

doi: http://dx.doi.org/10.1016/j.ces.2012.02.047.

URL

http://www.sciencedirect.com/science/article/pii/S0009250912001431. 730

Dietrich, B., Schabel, W., Kind, M., and Martin, H. ramic sponges - Determining the hydraulic diameter.

Pressure drop measurements of ceChemical Engineering Science, 64(16):

4

CONCLUSIONS AND FUTURE WORK 3633 – 3640, 2009.

ISSN 0009-2509.

38 doi:

http://dx.doi.org/10.1016/j.ces.2009.05.005.

URL

http://www.sciencedirect.com/science/article/pii/S0009250909003133. Edouard, 735

D.,

Lacroix,

SOLID foam: 311,

2008.

M.,

Huu,

P.,

C.,

and Luck,

State-of-the art correlation. ISSN 1385-8947.

doi:

F.

Pressure drop modeling on

Chemical Engineering Journal,

144(2):299 –

http://dx.doi.org/10.1016/j.cej.2008.06.007.

URL

http://www.sciencedirect.com/science/article/pii/S1385894708003586. Ergun, S. Fluid flow through packed columns. Chemical Engineering Progress, 48(2):89–94, 1952. URL http://pubs.acs.org/doi/pdf/10.1021/ie50474a011. 740

Ettrich, J., August, A., and Nestler, B. Open cell metal foams: Measurement and numerical modelling of fluid flow and heat transfer. In Cellular Materials - CellMat 2014, Dresden, Germany, 2014. URL http://www.dgm.de/download/tg/1589/1589-84.pdf. Forchheimer, P. Wasserbewegung durch Boden. Zeitschrift des Vereines Deutscher Ingenieuer, 1901.

745

Fourie,

J.

lar

metallic

ISSN

G.

and foams.

0009-2509.

Du

Plessis, Chemical

doi:

J.

P.

Pressure

Engineering

Science,

drop

modelling

57(14):2781

–

in

cellu-

2789,

2002.

http://dx.doi.org/10.1016/S0009-2509(02)00166-5.

URL

http://www.sciencedirect.com/science/article/pii/S0009250902001665. Geier, M., Greiner, A., and Korvink, J. Properties of the cascaded lattice Boltzmann automaton. International Journal of Modern Physics C, 18(04):455–462, 2007. doi: 10.1142/S0129183107010681. URL 750

http://www.worldscientific.com/doi/abs/10.1142/S0129183107010681. Gerbaux, O., Buyens, F., Mourzenko, V. V., Memponteil, A., Vabre, A., Thovert, J., and Adler, P.

Transport properties of real metallic foams.

(1):155–165, 2010.

ISSN 0021-9797.

doi:

Journal of Colloid and Interface Science, 342

http://dx.doi.org/10.1016/j.jcis.2009.10.011.

URL

http://www.sciencedirect.com/science/article/pii/S0021979709013241. 755

Gibson, L. J. and Ashby, M. F. Cellular Solids: Structure and Properties. Cambridge University Press, 1999. Graf von den Schulenburg, D. A., Paterson-Beedle, M., Macaskie, L., Gladden, L., and Johns, M.

Flow through an evolving porous media-compressed foam.

Journal of Materials Sci-

4

CONCLUSIONS AND FUTURE WORK ence, 42(16):6541–6548,

760

2007.

ISSN 0022-2461.

39 doi:

10.1007/s10853-007-1523-z.

URL

http://dx.doi.org/10.1007/s10853-007-1523-z. Große, J., Dietrich, B., Martin, H., Kind, M., Vicente, J., and Hardy, E. H. Volume image analysis of ceramic sponges. Chemical Engineering & Technology, 31(2):307–314, 2008. ISSN 1521-4125. doi: 10.1002/ceat.200700403. URL http://dx.doi.org/10.1002/ceat.200700403. Große, J., Dietrich, B., Garrido, G. I., Habisreuther, P., Zarzalis, N., Martin, H., Kind, M., and Kraushaar-

765

Czarnetzki, B. Morphological characterization of ceramic sponges for applications in chemical engineering. Industrial & Engineering Chemistry Research, 48(23):10395–10401, 2009. doi: 10.1021/ie900651c. URL http://pubs.acs.org/doi/abs/10.1021/ie900651c. Habisreuther, P., Djordjevic, N., and Zarzalis, N. rent in porous inert structure.

770

Numeric simulation of the micro cur-

Chemie-Ingenieur-Technik,

80(3):327–341,

2008.

URL

http://www.scopus.com/inward/record.url?eid=2-s2.0-41449089264|&partnerID=40|&md5=bf96a508229fd Higuera, F. J. and Jiménez, J. Boltzmann approach to lattice gas simulations. EPL (Europhysics Letters), 9(7):663, 1989. URL http://stacks.iop.org/0295-5075/9/i=7/a=009. Huang, H., Krafczyk, M., and Lu, X. Forcing term in single-phase and Shan-Chen-type multiphase lattice Boltzmann models. Phys. Rev. E, 84:046710, Oct 2011. doi: 10.1103/PhysRevE.84.046710.

775

URL http://link.aps.org/doi/10.1103/PhysRevE.84.046710. d’Humieres, D. and Ginzburg, I. Viscosity independent numerical errors for lattice Boltzmann models: From recurrence equations to "magic" collision numbers. Computers & Mathematics with Applications, 58(5):823 – 840, 2009. ISSN 0898-1221. doi: http://dx.doi.org/10.1016/j.camwa.2009.02.008. URL http://www.sciencedirect.com/science/article/pii/S0898122109000893. Mesoscopic Methods

780

in Engineering and Science. d’Humieres, D., Ginzburg, I., Krafczyk, M., Lallemand, P., and Luo, L.-S. Multiple-relaxation-time lattice Boltzmann models in three dimensions. Phil. Trans. R. Soc. Lond. A, 360:437–451, 2002. URL http://royalsocietypublishing.org/content/roypta/360/1792/437.full.pdf. Huu, T. T., Lacroix, M., Huu, C. P., Schweich, D., and Edouard, D. Towards a more realistic mod-

785

eling of solid foam: Use of the pentagonal dodecahedron geometry. Chemical Engineering Science,

4

CONCLUSIONS AND FUTURE WORK

40

64(24):5131 – 5142, 2009. ISSN 0009-2509. doi: http://dx.doi.org/10.1016/j.ces.2009.08.028. URL http://www.sciencedirect.com/science/article/pii/S0009250909005661. Inayat, A., Freund, A., Zeiser, T., and Schwieger, W. ceramic foams: 790

Determining the specific surface area of

The tetrakaidecahedra model revisited.

1179 – 1188, 2011a.

ISSN 0009-2509.

Chemical Engineering Science, 66(6):

doi: http://dx.doi.org/10.1016/j.ces.2010.12.031.

URL

http://www.sciencedirect.com/science/article/pii/S000925091000744X. Inayat, A., Freund, H., Schwab, A., Zeiser, T., and Schwieger, W. Predicting the specific surface area and pressure drop of reticulated ceramic foams used as catalyst support.

Advanced Engi-

neering Materials, 13(11):990–995, 2011b. ISSN 1527-2648. doi: 10.1002/adem.201100038. 795

http://dx.doi.org/10.1002/adem.201100038. Inayat, W. ideal 2011c.

800

URL

A.,

Schwerdtfeger,

Periodic

open-cell

tetrakaidecahedra ISSN

J.,

Freund,

foams:

packing.

0009-2509.

doi:

H.,

Koerner,

C.,

Pressure

drop

Chemical

Engineering

Singer,

measurements Science,

R., and

and

Schwieger,

modeling

66(12):2758

of –

http://dx.doi.org/10.1016/j.ces.2011.03.031.

an

2763, URL

http://www.sciencedirect.com/science/article/pii/S0009250911001977. Incera Garrido, G., Patcas, F., S., L., and Kraushaar-Czarnetzki, B. Mass transfer and pressure drop in ceramic foams: A description for different pore sizes and porosities. Chemical Engineering Science, 63(21):5202 – 5217, 2008. ISSN 0009-2509. doi: http://dx.doi.org/10.1016/j.ces.2008.06.015. URL http://www.sciencedirect.com/science/article/pii/S0009250908003485.

805

Karlin, I. V., Bösch, F., and Chikatamarla, S. S. Gibbs’ principle for the lattice-kinetic theory of fluid dynamics.

Phys. Rev. E, 90:031302, Sep 2014.

doi: 10.1103/PhysRevE.90.031302.

URL

http://link.aps.org/doi/10.1103/PhysRevE.90.031302. Kennedy, M., Zhang, K., Fritzsch, R., Akhtar, S., Bakken, J. A., and Ragnhild, E. A. acterization of ceramic foam filters used for liquid metal filtration. 810

Transactions B, 44(3):671–690, 2013.

ISSN 1073-5615.

Char-

Metallurgical and Materials

doi: 10.1007/s11663-013-9799-7.

URL

http://dx.doi.org/10.1007/s11663-013-9799-7. Kozeny, J. Uber kapillare leitung des wassers im boden. Sitzungsber Akad. Wiss. Wien, pages 271–306, 1927.

4

CONCLUSIONS AND FUTURE WORK

Kupershtokh, 815

A.,

Medvedev,

tice Boltzmann method. 2009.

ISSN

D.,

41

and Karpov,

D.

On equations

of state in a lat-

Computers & Mathematics with Applications,

0898-1221.

doi:

58(5):965 – 974,

http://dx.doi.org/10.1016/j.camwa.2009.02.024.

URL

http://www.sciencedirect.com/science/article/pii/S0898122109001011. Lacroix, M., Nguyen, P., Schweich, D., Huu, C. P., Savin-Poncet, S., and Edouard, D. sure drop measurements and modeling on SiC foams. 820

3259 – 3267, 2007.

ISSN 0009-2509.

doi:

Pres-

Chemical Engineering Science, 62(12):

http://dx.doi.org/10.1016/j.ces.2007.03.027.

URL

http://www.sciencedirect.com/science/article/pii/S0009250907002618. Laé, E., Duval, H., Riviere, C., Brun, P. L., and Guillot, J.-B. Experimental and Numerical Study of Ceramic Foam Filtration, pages 285–290. John Wiley & Sons, Inc., 2013. ISBN 9781118647783. doi: 10.1002/9781118647783.ch34. URL http://dx.doi.org/10.1002/9781118647783.ch34. 825

Łaniewski-Wołłk, Ł. R-language-based packages for code generation and symbolic operations, 2014. URL https://c-cfd.meil.pw.edu.pl/software/r-packages/. Łaniewski-Wołłk, Ł. and Rokicki, J. Adjoint lattice Boltzmann method for topology optimisation. Computers and Mathematics with Applications, 2015. URL http://arxiv.org/abs/1501.04741. Lautensack, C. and Sych, T. A random Weaire Phelan foam. In Proceedings of the VIII International

830

Conference on Stereology and Image Analysis in Materials Science STERMAT, Zakopane, 2008. URL http://www.mathematik.uni-kl.de/uploads/tx_sibibtex/RandomWeairePhelan.pdf. Lindblad, J. Surface area estimation of digitized 3D objects using weighted local configurations. Image Vision Comput., 23(2):111–122, February 2005. ISSN 0262-8856. doi: 10.1016/j.imavis.2004.06.012. URL http://dx.doi.org/10.1016/j.imavis.2004.06.012.

835

Liu, H., Schmieschek, S., Leonardi, C. R., Jones, B. D., Kang, Q., Narváez, A., Williams, J. R., Valocchi, J., A., and Harting, J. Multiphase lattice Boltzmann simulations for porous media applications. submitted to publication, 2015. doi: http://arxiv.org/pdf/1404.7523.pdf. Liu, Y.-S., Yi, J., Zhang, H., Zheng, G.-Q., and Paul, J.-C. digitized 3D objects using quasi-Monte Carlo methods.

840

3909, 2010.

ISSN 0031-3203.

doi:

Surface area estimation of

Pattern Recognition, 43(11):3900–

http://dx.doi.org/10.1016/j.patcog.2010.06.002.

http://www.sciencedirect.com/science/article/pii/S0031320310002669.

URL

4

CONCLUSIONS AND FUTURE WORK

Magnico, and

P.

anisotropic

2009. 845

Analysis

ISSN

of

metallic

42

permeability

foams.

0009-2509.

and

Chemical doi:

effective

viscosity

Engineering

by

Science,

CFD

on

64(16):3564

isotropic –

3575,

http://dx.doi.org/10.1016/j.ces.2009.04.036.

URL

http://www.sciencedirect.com/science/article/pii/S0009250909002954. Matzke, E.

The three-dimensional shape of bubbles in foam - An analysis of role of sur-

face forces in three-dimensional cell shape determination.

Am. J. Bot., 33:58–80, 1946.

URL

http://www.jstor.org/stable/2437492?seq=1#page_scan_tab_contents. Moreira, E. A. and Coury, J. R. The influence of structural parameters on the permeability of ceramic 850

foams. Brazilian Journal of Chemical Engineering, 21:23 – 33, 01 2004. ISSN 0104-6632. Narváez Salazar, A. and Harting, J. Evaluation of pressure boundary conditions for permeability calculations using the lattice Boltzmann method. Advances in Applied Mathematics and Mechanics, 2:685–700, 2010. doi: 10.4208/aamm.10-10S11. URL http://mtp.phys.tue.nl/publications/2010/NH10. Pan,

855

C.,

for ISSN

Luo,

porous

L.-S., medium

and Miller, flow

0045-7930.

C. T.

simulation.

doi:

An evaluation of lattice Boltzmann schemes Computers

&

Fluids,

35(8-9):898–909,

http://dx.doi.org/10.1016/j.compfluid.2005.03.008.

http://www.sciencedirect.com/science/article/pii/S0045793005001520.

2006. URL

Proceedings

of

the First International Conference for Mesoscopic Methods in Engineering and Science. Petrasch, J., Meier, F., Friess, H., and Steinfeld, A. 860

Tomography based determination of

permeability, Dupuit-Forchheimer coefficient, and interfacial heat transfer coefficient in reticulate porous ceramics. 2008.

ISSN 0142-727X.

International Journal of Heat and Fluid Flow, 29(1):315 – 326, doi:

http://dx.doi.org/10.1016/j.ijheatfluidflow.2007.09.001.

URL

http://www.sciencedirect.com/science/article/pii/S0142727X07001282. Prasianakis, 865

N. I.,

Rosen,

T.,

Kang,

J.,

Eller,

J.,

Mantzaras,

J.,

and Buchi,

F. N.

Simulation of 3d porous media flows with application to polymer electrolyte fuel cells. Commun. Comput. Phys.,

13:851–866,

2013.

doi:

10.4208/cicp.341011.310112s.

URL

http://www.diva-portal.org/smash/record.jsf?pid=diva2%3A605731&dswid=8661. Ranut, P., Nobile, E., and Mancini, L. Microtomography-based CFD analysis of transport in open-

4

CONCLUSIONS AND FUTURE WORK cell aluminum metal foams.

870

43

Journal of Physics: Conference Series, 501(1):012021, 2014.

URL

http://stacks.iop.org/1742-6596/501/i=1/a=012021. Rautenbach,

C.

Modelling

sorption tower.

of

PhD thesis,

flow

through

porous

Stellenbosch:

packing

University

elements

of

of Stellenbosch,

a

CO2

2009.

abURL

http://scholar.sun.ac.za/handle/10019.1/2319. Regulski, W. and Szumbarski, J. Numerical simulation of confined flows past obstacles - the comparative 875

study of lattice Boltzmann and spectral element methods. Archives of Mechanics, 64:423–456, 2012. URL http://am.ippt.pan.pl/am/article/view/v64p423/pdf. Richardson,

J.

T.,

lyst supports: ISSN 880

Peng,

Y.,

and

pressure drop.

0926-860X.

doi:

Remue,

D.

Properties

of

Applied Catalysis A: General,

ceramic

204(1):19

foam

cata-

– 32,

2000.

http://dx.doi.org/10.1016/S0926-860X(00)00508-1.

URL

http://www.sciencedirect.com/science/article/pii/S0926860X00005081. Scheffler, M. and Colombo, P. Cellular ceramics: structure, manufacturing, properties and applications. John Wiley & Sons, 2006.

885

Skibiński,

J.,

tional

design

IX,

volume

Wejrzanowski, of 74,

the pages

flow 109

T.,

Szumbarski,

properties –

118,

of

J.,

and

foams.

2012.

Kurzydłowski, In

doi:

Advances

in

K.

Computa-

Fluid

Mechanics

http://10.2495/AFM120101.

URL

http://www.witpress.com/elibrary/wit-transactions-on-engineering-sciences/74/23555. Succi, S. The Lattice Boltzmann Equation: For Fluid Dynamics and Beyond. Numerical Mathematics and Scientific Computation. Clarendon Press, 2001. ISBN 9780198503989. Twigg, M. V. and Richardson, J. T. Fundamentals and applications of structured ceramic foam catalysts. 890

Industrial & Engineering Chemistry Research, 46(12):4166–4177, 2007. doi: 10.1021/ie061122o. URL http://dx.doi.org/10.1021/ie061122o. van

der

Sman,

Forchheimer ISSN 895

R.

Prediction

equation.

0260-8774.

of

Journal doi:

airflow of

Food

through

a

vented

Engineering,

box

by

the

Darcy-

55(1):49

–

57,

2002.

http://dx.doi.org/10.1016/S0260-8774(01)00241-2.

http://www.sciencedirect.com/science/article/pii/S0260877401002412.

URL

4

CONCLUSIONS AND FUTURE WORK

44

Viggen, E. M. The lattice Boltzmann method with applications in acoustics. Master Thesis, 2009. URL http://www.iet.ntnu.no/groups/akustikk/publications/pdfs/NTNU_MSc2009_Viggen.pdf. Wejrzanowski, T., Skibiński, J., Szumbarski, J., and Kurzydłowski, K. modeled by Laguerre-Voronoi tessellations. 900

221, 2013.

ISSN 0927-0256.

doi:

Structure of foams

Computational Materials Science,

67(0):216 –

http://dx.doi.org/10.1016/j.commatsci.2012.08.046.

URL

http://www.sciencedirect.com/science/article/pii/S0927025612005423. Whitaker, S. The Forchheimer equation: A theoretical development. Transport in Porous Media, 25(1):27– 61, 1996. ISSN 0169-3913. doi: 10.1007/BF00141261. URL http://dx.doi.org/10.1007/BF00141261.