Category Archives: Diffusion

The problem with geometric factors

5 Replies

Many papers on diffusion in compacted bentonite and claystone1 assume a direct connection between diffusivities measured in clay (\(D_\mathrm{clay}\)), and the geometrical configuration of the pore space. It seems close to mandatory to include
an expression of this type \begin{equation} D_\mathrm{clay} = \frac{1}{G}\cdot D_\mathrm{w}\tag{1} \end{equation}

where \(D_\mathrm{w}\) is the corresponding diffusion coefficient in bulk water (usually taken as the limiting value at infinite dilution), and \(G\) is a factor postulated to have pure geometric meaning. The geometric multiplier is written in a variety of ways, sometimes in terms of other factors which individually are supposed to convey specific geometric information, e.g. \(1/G = \delta/\tau^2\), where \(\delta\) is the “constrictivity” and \(\tau^2\) is the “tortuosity factor”2.

However, the above described procedure is futile:

  1. A geometric factor can not be deduced a priori, because the pore space geometry of any realistic system is way too complicated. So-called tortuosities have only been calculated in some trivial cases, relating e.g. to highly symmetric networks of identical cylinders, which has basically no resemblance to the chaotic smectite-lasagna expected to be found in bentonite. Indeed, many papers on these issues admit that the geometric factor in practice functions as a “fudge factor” (to be “optimized”).
  2. Even if a geometric factor could be deduced from only considering the pore space configuration, it is not expected to relate to \(D_\mathrm{clay}\) in accordance with eq. 1. If we imagine removing all influence of pore geometry on the diffusivity in the clay, what remains is diffusion in an environment heavily influenced by the presence of exchangeable cations. Thus, a correct geometric factor should not multiply \(D_\mathrm{w}\), but the diffusion coefficient measured on the short time and length scale in the clay. Such measurements exist — at least for water diffusion — and show values different from \(D_\mathrm{w}\).

Because it is practically impossible to independently deduce a geometric factor, it is de facto defined by the ratio \(D_\mathrm{clay}/D_\mathrm{w}\) in many studies. I dare to say that in any work within the bentonite research community, a geometric factor has no other meaning than \(D_\mathrm{clay}/D_\mathrm{w}\) (or some function of this parameter, depending on definitions).

This means — at best — that the use of a geometric factor does nothing for a model description, but only replaces the parameter \(D_\mathrm{clay}\) with the parameter \(D_\mathrm{clay}/D_\mathrm{w}\); it simply corresponds to the awkward choice of specifying diffusion coefficients in units of \(D_\mathrm{w}\).

Often, however, the procedure has worse consequences, because eq. 1 (as used in practice) represents an unjustified claim; it is not correct to use the value of \(D_\mathrm{clay}/D_\mathrm{w}\) to make specific statements about the pore space geometry, while ignoring other possible mechanisms contributing to the value of \(D_\mathrm{clay}\). If diffusion instead is described without invoking a geometric factor, \(D_\mathrm{clay}\) itself functions as as a model parameter, and no implicit claims are made about why it has the value it has.

There are many examples in the literature where this confusion shows up. In e.g. Appelo and Wersin (2007), it is concluded that tortuosity for “anions” is larger than for water, while tortuosity for “cations” is smaller than for water. This conclusion is based on fitting of geometric factors without independent information on the pore space geometry:

[…] the pore may become filled entirely with a diffuse double layer when it narrows sufficiently. This constricts the passage of anions, and since the anions must circumnavigate the obstacle, they have greater tortuosity than tritium. This explains that a model with a tortuosity factor for iodide that is 1.4 times higher than for tritium better matches the data.

[…] in interlayers and pore constrictions the cations pass in relatively larger amounts than in free porewater, and consequently, they have smaller tortuosity than tritium.

Similar descriptions are given in Altmann et al. (2012) (in a section describing “advances” in process understanding with the FUNMIG project)

At low ionic strengths, cations are mostly located in the diffuse layer. This layer is of relatively low volume and connected throughout the porosity — cations need less time to explore this small volume and therefore to pass from one side of the media to the other.[…]

Anions are excluded from the diffuse layer and, therefore,also have a smaller volume to explore than neutral tracers. They however must pass through small pores in order to explore the total porosity, and these pores act as electrostatic throats […], forming an energy barrier to anion movement, reducing the probability of passage (resulting in reduced \(D_e\)) and increasing the time needed to explore all of the porosity.

And the same type of argument — that different diffusive behavior of cations and anions is caused by differences in “pathways” — is found in Charlet et al. (2017)

The different pathways of anions and cations in these regions lead to different effective diffusion rates (Van Loon et al., 2007; Mazurek et al., 2009)

But effective diffusion rates — which can basically be made to vanish e.g. for chloride, as has been demonstrated in bentonite — do not vary because of changes in “pathways”, but because of changes in boundary conditions, due to the presence of interfaces to external solutions.

To investigate the behavior of \(D_\mathrm{clay}/D_\mathrm{w}\) one must of course consider real — not effective — diffusion coefficients (unfortunately, actual diffusion coefficients are usually referred to as “apparent” in the bentonite research community). Real diffusion coefficients reflect the mobility of the diffusing species, and can be measured e.g. in non-steady-state closed-cell tests, which are not influenced by interfaces to external solutions3. A very good appreciation for how such diffusion coefficients behave is given by the vast amount of data measured by Kozaki et al. in Na-montmorillonite (Kunipia-F). Some of this data (sources: 1, 2, 3, 4, 5) is plotted below, in terms of \(D_\mathrm{clay}/D_\mathrm{w}\) (using a linear \(y\)-axis, for a change).

Kozaki's diffusion data.

This plot clearly shows that if \(D_\mathrm{clay}/D_\mathrm{w}\) is interpreted as having pure geometrical meaning, the conclusion must be that chloride and sodium have basically identical “tortuosity”, while the “tortuosity” of the cations calcium and cesium is considerably larger.4 Thus, the closed-cell data disproves any claim that “tortuosity” for anions is larger than for cations. Moreover, from this data it is quite a stretch to state that “tortuosity” for anions is larger than for water, and it is plain wrong to state that “tortuosity” for cations is smaller than for for water. (For the record, I’m strongly convinced that both calcium and cesium occupy the same pore volume as sodium in this material.)

I should emphasize that I am not primarily arguing for a better geometrical interpretation of \(D_\mathrm{clay}/D_\mathrm{w}\), but rather that this whole procedure (using eq. 1) should be avoided. This is not to say that pore geometry doesn’t influence clay diffusivity — most reasonably it does — but analysis based on eq. 1 will be flawed, because of the unjustified model assumptions.

When it comes to geometrical understanding of diffusion in bentonite, I think it would be more rewarding to focus on the quite large variation of reported diffusivities in seemingly identical systems. As an example, the figure below shows the diffusivity of cesium in very similar, pure Na-montmorillonite, measured in four different closed-cell studies: apart from Kozaki et al.’s data (blue), also Sato et al. (1992) (green), and Tachi and Yotsuji (2014) (orange).

Cs diffusion data from Kozaki et al. Sato et al., and Tachi and Yotsuji

In aggregation, this figure demonstrates almost an order of magnitude variability for the diffusion coefficient of cesium, even in very pure systems. If this kind of variability cannot be explained, it is of course not very relevant to develop models describing (smaller) differences in diffusivity between different types of ions in different types of bentonites.

If you claim to understand the variation seen in this plot, rather than simply finding an appropriate factor to blame (particle size?), here is a challenge: Give a prescription for how to prepare a Kunipia-F/P sample in order to minimize/maximize the resulting diffusivity.5

Footnotes

[1] In the following, the term “bentonite” is used as a shorthand for “compacted bentonite and claystone”.

[2] The nomenclature here is a mess, but I leave it to people who really want to use these concepts to actually define them in detail.

[3] They can also be deduced from the transient behavior in “steady-state” tests, if the transport capacity of the confining filters is not limiting the process.

[4] Note that “tortuosity”, as used e.g. by Appelo and Wersin (2007), is related to the inverse of \(D_\mathrm{clay}/D_\mathrm{w}\).

[5] The diffusion direction should still be along the axis of sample compaction, of course.

You cannot add fluxes from different domains in multiporosity models

2 Replies

Most models of compacted bentonite and claystone1 assume that the material contains different domains with distinctly different properties. In particular, it is common to assume that the system contains both an aqueous domain (a bulk water solution) and one or several domains associated with the solid surfaces.2 Models of this type is here called multiporosity models.

A “standard” procedure in multiporosity models is to assume that the diffusive flux in any domain can be fully expressed using quantities that relates only to that same domain, while a “total” flux is given by summing the contributions from each domain. In e.g. the case of having one bulk water domain and one surface domain, the diffusive flux is typically written something like

$$ \label{eq1} j = j_{bulk} + j_{surface} = -D_e\cdot\nabla c -D_s \cdot\nabla c_s \tag{1} $$

were \(D_e\) and \(c\) denote, respectively, the diffusion coefficient and concentration in the bulk water domain, and the corresponding parameters in the surface domain are denoted \(D_s\) and \(c_s\). It is important to keep in mind that what is referred to as a “total” flux in these descriptions is the experimentally accessible, macroscopic flux (here simply labelled \(j\)).

This approach is basically adopted by all proponents of multiporosity models, and has been so for a long time, see e.g.:

There is, however, a big problem with this procedure:

Equation 1 does not make sense

The statement of eq. 1 — that the macroscopic, measurable flux has two independent contributions — is only true if the different domains are completely isolated from each other. This can e.g. be seen by making \(D_s\) vanishingly small in eq. 1; the flux is still non-zero due to the contribution from the completely unaffected \(j_{bulk}\)-term. In contrast, for non-isolated domains the expected behavior is that a macroscopic diffusive flux will vanish if the diffusivity of any domain is made to approach zero.

Consequently, a large part of the compacted bentonite research field seems to treat diffusion by assuming the different domains to be isolated “pipes”. But this is of course not how bentonite works! If you postulate multiple diffusive domains in the clay, you must certainly allow for diffusing species to be able to be transferred between them. And here is where it gets mysterious: the default choice is to also assume that the domains are non-isolated. In fact, the domains are typically assumed to be in equilibrium. In e.g. the case of having one bulk water and one surface domain, the equilibrium is usually expressed something like

$$
c_s = K\cdot c
\tag{2}
$$

where \(K\) is some sort of partition coefficient.

Thus, the “standard” treatment is to combine eq. 1 — even though it is only valid for isolated domains — with an equilibrium assumption (eq. 2) to arrive at expressions which typically look something like this

$$
j = \left (D_e + K\cdot D_s\right)\cdot \nabla c
\tag{3}
$$

Note how this equation has the appearance of an addition of diffusion coefficients, as the “total” flux now is related to a single concentration gradient (which in practice always is a bulk solution concentration).

Since eq. 3 results from combining two incompatible assumptions (isolation and equilibrium), it should come as no surprise that it is not a valid description for multiporosity models.

I cannot get my head around why this erroneous description has prevailed within the compacted clay research field for so long. The problem of deriving “effective” properties of heterogeneous systems appears in many research fields, and applies to many other quantities than diffusivity, e.g. heat conductivity, electrical conductance, elastic moduli, etc. It should therefore be clear also in the compacted clay world that obtaining the effective diffusivity for a system as a whole is considerably more involved than simply adding diffusivities from each contributing domain.

Some illustrations

To discuss the issue further, let’s make some illustrations. A molecule or ion in a porous system conducts random “jiggling” motion, exploring all regions to which it has access, as illustrated in this figure

This animation shows a simulated random walk at the microscopic scale, in a “made-up” generic porous system, consisting of solids (brown), a bulk water domain (blue), and a surface domain (pink). The green dots show the trace of the random walker; light green when the walker is located in the bulk water domain (with higher mobility), and darker green when it is located in the surface domain (with lower mobility). Note how the walker frequently switches between the domains — which is exactly what is required for the two domains to be in equilibrium.

As diffusivity is fundamentally related to the statistical properties of the random walk,3 it should be crystal clear that the macroscopic diffusion coefficient (reflecting mobility on length scales beyond what is shown in the animation) depends in a complicated way on the mobility of the involved domains, just as it depends on the detailed geometrical configuration. Oddly, the complex geometrical dependency is usually acknowledged within the bentonite research field (by means of e.g. formation and/or tortuosity factors), which makes it even more incomprehensible that the issue of multiple domain diffusivities is treated as it is.

In contrast, the figure below illustrates the microscopic behavior implied by eq. 1.

Here is shown one random walker in the bulk domain (green trace) and one in the surface domain (red trace). Since no exchange occurs between the domains, they are not in equilibrium (the concept does not apply). The lack of exchange between the domains is of course completely non-physical, but it is only under this assumption that the system is described by two independent fluxes (eq. 1). It may also be noted that the two types of walkers experience completely different geometrical configurations (i.e. formation factors), none of which corresponds to the geometrical configuration in the case of non-isolated domains.

A more reasonable formulation

If you insist on having a multiporosity model, the correct way to go about is to assume equilibrium within a so called representative elementary volume (REV), and to formulate the model on the macroscopic — averaged — level. On this level of description, it is legitimate to speak of both a bulk water and a surface concentration at any given point in the model (they represent REV averages), but there is of course only a single flux, as well as a single value of the chemical potential at any given point. The diffusive flux can be written quite generally as

\begin{equation} j = -\phi
D_\mathrm{macr.}\cdot\frac{\bar{c}}{RT}\nabla\mu
\end{equation}

where \(\bar{c}\) is the macroscopic, averaged, concentration, \(\phi\) is the porosity, and \(\mu\) is the (electro-) chemical potential. \(D_\mathrm{macr.}\) is the macroscopic diffusion coefficient, containing statistical information on how the diffusing substance traverses the porous system on length scales larger than the REV (in particular, it contains averaged information on the mobility in all involved domains, as well as their detailed geometrical configuration).

The value of \(\mu\) can be obtained from considering the microscopic scale. Since the system by assumption is in equilibrium within a REV, \(\mu\) has the same value everywhere within this volume. On the other hand, the (electro-) chemical potential may also be expressed e.g. as

\begin{equation}
\mu(\xi ) = \mu_0 + RT\cdot\ln a(\xi ) + zF\psi(\xi ) = \mathrm{const}
\end{equation}

where \(\xi\) is a microscopic variable denoting the position within the REV, \(\mu_0\) is a reference potential, \(a(\xi )\) is the chemical activity and \(\psi(\xi )\) the electrostatic potential. In the bulk water domain, the electrostatic potential is zero by definition (ignoring possible couplings between different types of charged species), and \(\mu\) can here be evaluated as (where the index “b” is used to denote bulk)

\begin{equation}
\mu_b = \mu_0 + RT\cdot\ln a_b
\tag{4}
\end{equation}

But this is by definition the value of the (electro-)chemical potential everywhere within the REV.

Thus, going back to the macroscopic scale and using eq. 4 for \(\mu\), we have for the flux

\begin{equation}
j = -\phi D_\mathrm{macr.}\cdot\frac{\bar{c}}{a_b}\nabla a_b
\end{equation}

Ignoring also possible gradients in activity coefficients, the flux can be simplified as

\begin{equation}
j = -\phi D_\mathrm{macr.}\cdot\frac{\bar{c}}{c}\nabla c
\end{equation}

Where, as before, \(c\) denotes a bulk solution concentration. For the case of having one bulk domain, and one surface domain, \(\bar{c}\) can be written

\begin{equation}
\bar{c} = \frac{\phi_b}{\phi}c + \frac{\phi_s}{\phi}c_s
\end{equation}

where \(\phi_b\) and \(\phi_s\) are the porosities for the bulk and surface domains, respectively (\(\phi = \phi_b + \phi_s\)). Using also eq. 2, the flux can finally be written

\begin{equation}
j = -D_\mathrm{macr.}\left (\phi_b+\phi_s\cdot K\right )\nabla c
\tag{5}
\end{equation}

This expression can be compared to the result from the “standard” procedure, eq. 3. Although eq. 3 is derived using incompatible, non-physical assumptions, it is similar to eq. 5 , but involves two (\(D_e\) and \(D_s\)) rather than a single diffusion parameter (\(D_\mathrm{macr.}\)). If these two parameters are used merely as fitting parameters (which they basically always are), eq. 3 can always be made to resemble eq. 5. In this sense, eq. 3 is not wrong. Rather, it is not even wrong.

Footnotes

[1] In the following I will only say bentonite, but I mean any type of clay system with significant swelling/ion exchange properties.

[2] I strongly oppose the way bulk water phases are included in such models, but this post is not a critique of that model choice. Here, I instead point out that virtually all published descriptions of multiporosity models handle the flux incorrectly.

[3] Specifically, the average of the square of the displacement of the walker is proportional to time, with the diffusion coefficient as proportionality constant (apart from a dimensional factor).