Effective Diffusivity in Porous Materials

This example introduces the concept of effective diffusivity in porous media by comparing the transport through an artificial porous structure described in a detailed model with a simplified homogeneous porous media approach using effective transport properties.

The exercise consists of two parts. The first part describes how to create the model with a detailed geometry. The second part shows how to define a homogeneous model for porous media using an effective diffusivity calculated using the detailed model from the first part.

Introduction

Transport through porous structures is usually treated using simplified homogeneous models with effective transport properties. This is in most cases a necessity, since the typical dimensions of the pores and particles making up the porous structure are several orders of magnitude smaller than the size of the domain that is to be modeled.

However, it might be interesting to investigate the assumptions and simplifications done when homogenizing a porous structure by comparing a homogeneous model with a model defined using the detailed structure.

The artificial porous structure used in this example is shown in Figure 1 below.

Figure 1: Artificial porous structure. The domain colored in red is accessible for diffusion.

Model Definition

The model equation in the modeled domain shown in Figure 1 is the time-dependent equation

where c denotes concentration (mol/m3 using SI units) and D the diffusion coefficient (m2/s) of the solute.

The boundary conditions are of three different types. A concentration boundary condition applies at the left vertical boundary in Figure 1. It is expressed as

where c0 is a given concentration.

The right vertical boundary in Figure 1 is set according to

where km is the mass transfer coefficient (m/s), and c1 is the concentration in a bulk solution outside of the porous structure.

All other boundaries are insulating boundaries according to

The initial condition is given by a bell-shaped profile along the x-axis with its maximum at x = 0 and a corresponding value of c = c0:

Assume a gaseous solution with a solute content of 3 mol/m3 at the concentration boundary. The diffusion coefficient is set to 1·10−5 m2/s.

The second part of this exercise uses a homogenized 1D model geometry with effective transport properties and an average porosity. The model equation then becomes:

where ε denotes the average porosity and Deff is the effective diffusivity. These properties are calculated from the results of the detailed structure; see the next section. At the boundaries, the concentration and flux conditions described above apply.

Results and Discussion

The simulations are run for t = 0 to 0.1 s, when the simulation reaches steady state. Figure 2 below shows the concentration profile after 0.05 s in the porous structure. Already at this stage the concentration has almost reached steady state, which is visible by the nearly linear concentration profile across the structure.

Figure 2: Concentration profile in the modeled artificial porous structure at t = 0.05 s.

When modeling porous media, the exact concentration in the pore structure is not the most important issue because the description of the structure is homogenized and not detailed as in Figure 2. The most interesting issue is then the description of the flux. To calculate the average flux, integrate over the flux boundary and divide by its length, L0, which yields the following expression:

Figure 3 shows the value of this integral as a function of time. If you let the process reach steady-state, the average flux becomes 8.051·10−3 mol/(m2/s). Considering the almost linear profile across the structure, it is natural to replace the porous structure with a 1D homogenized structure along the x-axis. It is then possible to calculate the effective diffusivity according to the following:

where cout is the average concentration (mol/m3) at the flux boundary, and L1 is the length of the geometry along the x-axis. The average concentration is obtained by integrating according to the expression below:

This gives an average concentration cout = 1.63·10−3 mol/m3. Using L1 = 8·10−4 m, the effective diffusivity is:

which yields a value for the effective diffusivity of 2.15·10−6 m2/s compared to the “free” diffusivity of 1·10−5 m2/s. The effective and “free” diffusivities are usually related according to the equation

where ε is the porosity of the structure and τ is the tortuosity, which is a measure of the actual length per unit effective length a molecule has to diffuse in a porous structure. To calculate the porosity of the modeled structure, you integrate the value 1 over the structure and then divide this by the length and width of the structure:

resulting in a value of 0.383. The value of τ can then be calculated to 1.62. In addition, the tortuosity is usually expressed as a power of the porosity, resulting in an expression for the effective diffusivity according to

If you use the calculated values for porosity and effective diffusivity, the value for p is 1.60. The experimental values for p for porous structures used for transport in catalysts, soils, and other porous structures is usually in the range 1.5–2.

Using the value of the effective diffusivity, a simple homogenized 1D model provides the possibility to compare the value of the flux using a homogenized model to the value using the detailed 2D structure. Figure 3 shows that there is an excellent agreement between the model using a detailed geometry and the homogenized model. The difference in the time-dependent flux is hardly visible between the two cases in the graph.

Figure 3: Average flux at the flux boundary in the detailed 2D model (solid blue line) and the 1D homogenized approximation (dashed green line).

Notes About the COMSOL Implementation

Both models described above are straightforward to define in COMSOL Multiphysics. One feature that is of great use in this example is the ability to define integration coupling operators to generate the values of the integrals needed to evaluate the results. The definition of these integrals is described in detail in the step-by-step instructions below.

COMSOL_Multiphysics/Diffusion/effective_diffusivity

Modeling Instructions

From the menu, choose .

New

In the window, click

Model Wizard

In the window, click

In the tree, select .

Click .

Click

In the tree, select .

Click

Insert the geometry sequence to focus on the simulation. If you want to know how to build the geometry, you find instructions in the appendix.

Geometry 1

In the toolbar, click

Browse to the model’s Application Libraries folder and double-click the file effective_diffusivity_geom_sequence.mph.

In the toolbar, click

Global Definitions

Parameters 1

In the window, under click .

In the window for , locate the section.

In the table, enter the following settings:


Name	Expression	Value	Description
D2	1e-5[m^2/s]	1E-5 m²/s	Diffusion coefficient
c_max	3[mol/m^3]	3 mol/m³	Peak initial concentration
k_f	5[m/s]	5 m/s	Mass transfer coefficient
a	1000	1000	Dimensionless constant