• Part I
• Part II
• Part III

This is the fourth part of the review of “Ionic Transport in Nano-Porous Clays with Consideration of Electrostatic Effects” (Tournassat and Steefel, 2015) (referred to as TS15 in the following). For background and context please check the first part. This part covers the two sections “Constitutive equations for diffusion in bulk, diffuse layer, and interlayer water” and “Relative contributions of concentration, activity coefficient and diffusion potential gradients to total flux”.

“Constitutive equations for diffusion in bulk, diffuse layer, and interlayer water”

This section presents a mathematical formulation of ion diffusion in bentonite,¹ based on the material descriptions in the earlier sections. As we have previously noted, these descriptions are fundamentally flawed in several respects. In particular, compacted bentonite is presented as consisting of stacks (called “particles”), where it is supposed to make sense to differ between external and internal interface water. TS15 also mean that compacted bentonite (sometimes?) is supposed to contain a bulk water phase.

As I have commented on in earlier parts, the only reason I can see to provide this nonsensical material description is as an attempt to to motivate a macroscopic, multi-porous model of bentonite. Here, TS15 make this claim quite explicit, as they write

Still it is possible to define three porosity domains, or water domains, that can be handled separately: the bulk water, the diffuse layer water and the interlayer water, the properties for which can be each defined independently.

This is in essence what I have referred to as “the mainstream view” of bentonite. It is basically “possible” to define anything, but the real question is if provided definitions are relevant and useful. And, as we have already discussed in detail, there is no rationale for introducing these “porosity domains” when modeling water saturated, compacted bentonite.²

Here we will first comment on the conceptual aspects of the provided mathematical description. Thereafter, we will delve into the mathematical formulations, as I’m quite convinced that these are not correct. Unfortunately, this latter part will be quite burdened with equations and notation, but for the motivated reader I think it may be worth going through.

Conceptual aspects

TS15 choose the Nernst-Planck description of ion diffusion, and begin by commenting that this is more rigorous than using Fick’s law. I certainly agree with that a general description of ion diffusion in bentonite requires treating electrostatic couplings between the various system components (TOT-layers, ions). I don’t think, however, that putting up a massively complex description of multi-component diffusion in “three porosity domains” is the appropriate starting point for including such couplings. Since we have every reason to believe that e.g. no bulk water phase is present, I mean that this type of treatment only leads us astray from understanding the actual processes involved (we will return to this aspect in later parts of the review).

Also, as the “Fickian” aspect was the focus of the earlier section on diffusion, a reader of TS15 could here understandably get the impression that a Nernst-Planck treatment will “fix” the “issues” addressed there. But, as we have already discussed in some detail, the shortcomings of the traditional sorption-diffusion model are not solved by including multi-component diffusion in a bulk water phase. They are solved by removing the bulk water phase.

Although the above quotation states that the various “porosity domains” can be handled separately, and that their properties can be defined independently, this is not what is done in TS15. Rather, the treatment of any “porosity domain” assumes equilibrium with the corresponding bulk water phase. The entire description in TS15 is thus fully centered around the bulk water phase.

TS15 insist on treating their model quantities as functions of two spatial coordinates (\(x\) and \(y\)), in what they refer to as a “pseudo 2-D Cartesian system” (I don’t fully understand what that means). Diffusive flux is only assumed to take place in the \(x\)-direction, while the “\(y\)”-dimension is used for stacking the different “porosity domains”. The description can be schematically illustrated like this

Here we have for illustrative purposes discretized the various components in \(x\)- and \(y\)-directions. The bulk water domain is colored blue, the “interlayer” domain pink, and the “diffuse layer” domain green. For a given \(x\)-position, the “diffuse layer” and the “interlayer” domains are assumed to always be in equilibrium with the corresponding bulk water phase. TS15 nowhere consider the length scale in the \(y\)-direction (is it therefore the coordinate system is referred to as “pseudo 2-D”?), which in practice makes the model a collection of 1-dimensional domains that are in equilibrium locally. Note that even though diffusion only is accounted for in the \(x\)-direction, transport occurs also in the \(y\)-direction, as a consequence of equilibration between the “porosity domains”.

This description is exactly what we have investigated in the blog post on why multi-porous models cannot be taken seriously. To summarize what was said there, without properly defining the length scales, it makes no sense to “short-cut” the model in the \(y\)-direction (to assume equilibrium for all domains at the same \(x\) is in a sense equivalent to assuming infinitely high mobility of all components in the \(y\)-direction). And even if we assume that such an assumption is valid — which would mean that we consider a thin strip of stacked parallel domains, where the extension in \(y\) is negligible in comparison to the extension in \(x\) — the resulting model has really nothing to do with actual bentonite. As we concluded in the multi-porosity blog post, the only way to make sense of this type of description is as a set of macroscopic continua that are assumed to be locally in equilibrium. How this equilibrium is supposed to be maintained has never been suggested by any proponent of this description. Note that this description (in particular the existence of a bulk water phase in equilibrium) disqualifies the model for describing swelling and swelling pressure.

Incorrect application of the Nernst-Planck framework

While the presented model makes little sense conceptually, TS15 also fail in applying the Nernst-Planck framework. The problem arises, as far as I can see, from that they don’t fully recognize the role of the electric potential.

As we now begin scrutinizing the details of the formulation, we will suppress the variables \(x\) and \(y\) in order to, hopefully, make the equations a little more readable. It should be understood that any quantity is evaluated for some specific value of \(x\), and that all “porosity domains” are supposed to be in equilibrium at the same value of \(x\).

The electro-chemical potential

In most standard thermodynamic text books we learn that the chemical potential governs the equilibrium associated with mass transfer. Just as e.g. pressure and temperature (which govern mechanical and thermal equilibrium, respectively), the chemical potential is defined by a specific derivative of a thermodynamic potential, e.g.³

\begin{equation} \bar{\mu} = \left ( \frac{\partial G}{\partial n} \right )_{T,p} \tag{1} \end{equation}

where \(G\) is the Gibbs free energy, \(n\) number of moles, \(p\) pressure, and \(T\) temperature. The corresponding mass flux is generally written

\begin{equation} j = -\frac{cD}{RT}\nabla \bar{\mu} \tag{2} \end{equation}

where \(c\) is concentration, \(D\) the diffusion coefficient⁴, and \(RT\) the usual absolute temperature factor. Here, and in the following, we use the symbol \(\nabla\), which denotes the general gradient operator, but since the model is effectively one-dimensional, it can simply be seen as a neat way of writing \(\partial/\partial x\).

For charged species, it is common to refer to the quantity defined in eq. 1 as the electro-chemical potential, and write it as composed of an “ordinary” and a purely “electrical” part

\begin{equation} \bar{\mu} = \mu + zF\Psi \tag{3} \end{equation}

where \(z\) is the charge number of the considered species, \(F\) is the Faraday constant and \(\Psi\) is the electric potential. The “ordinary” chemical potential \(\mu\) (without bar) is, perhaps a bit confusingly, also often referred to as the chemical potential. I will here continue to refer to this part as “ordinary”. The “ordinary” chemical potential is furthermore conventionally expressed in terms of a reference potential (\(\mu^0\)) and an activity \(a\)

\begin{equation} \mu = \mu^0 + RT\ln a \tag{4} \end{equation}

A lot can be said about the decomposition in eq. 3, but it is clear that singling out an electric potential term is useful in e.g. electrochemistry or for describing charged clay. It should, however, be kept in mind that mass transfer is fundamentally governed by gradients in \(\bar{\mu}\); always keeping eqs. 1 and 2 in mind will avoid us from making mistakes, because the mass transfer rate relates to the “total” (i.e. electrochemical) potential, and for charge neutral species the description reduces to gradients in the “ordinary” chemical potential.

Nernst-Planck flux

Combining eqs. 3 and 4 gives

\begin{equation} \bar{\mu} = \mu^0+RT \ln a + zF\Psi \tag{5} \end{equation}

with the corresponding flux (eq. 2)

\begin{equation} j = -cD\nabla \ln a -\frac{cDzF}{RT} \nabla \Psi \end{equation}

Expressing the activity in terms of an activity coefficient, \(a = \gamma c\), the flux can also be written (TS15 are quite fond of including activity coefficients explicitly)

\begin{equation} j = -D\nabla c -cD\nabla \ln \gamma -\frac{cDzF}{RT} \nabla \Psi \tag{6} \end{equation}

Considering an arbitrary set of diffusing charged species (using the index \(i\)), and utilizing that the electric current is zero, lead to an expression for the electric potential gradient

\begin{equation} \nabla \Psi = -\frac{RT}{F}\frac{\sum z_iD_i \left ( \nabla c_i + c_i \nabla \ln \gamma_i \right )}{\sum z_i^2 D_i c_i} \tag{7} \end{equation}

Misunderstanding the electric potential

For the bulk water phase, TS15 indeed provide an expression for the flux that is essentially the same as eq. 6 (their eq. 37), and which they refer to as the Nernst-Planck equation. They claim, however, that the electrochemical potential in this case lack an electric potential term (my emphasis)^5,6

In absence of an external electric potential, the electrochemical potential in the bulk water can be expressed as (Ben-Yaakov 1981; Lasaga 1981) \begin{equation} \bar{\mu}_\mathrm{bulk} = \mu^0+RT \ln a_\mathrm{bulk} \end{equation}

But even without an externally applied electric field,⁷ a zero bulk electric potential cannot be assumed, of course, if the goal is to treat individual ion mobilities; as just shown, the gradient in electric potential that appears in eq. 6 is a result of a corresponding term in the electrochemical potential (eq. 5). Oddly, TS15 seem to treat the electric potential term in the flux as a quantity unrelated to the electrochemical potential, giving it a separate symbol, \(^\mathrm{b}\Psi_\mathrm{diff}\), and writing

\(^\mathrm{b}\Psi_\mathrm{diff}\) is the diffusion potential that arises because of the diffusion of charged species at different rates.

It may be natural for a reader at this point to simply assume that TS15 have missed writing out the term \(zF^\mathrm{b}\Psi_\mathrm{diff}\) when stating the electrochemical potential. But this seems to be a genuine misunderstanding rather than a mistake/typo, because the pattern repeats in the derivation of the flux in the other “porosity domains”.

For e.g. the “diffuse layer”,⁸ TS15 recognize the presence of an electric potential in the expression for the electrochemical potential, writing it (this is more awkwardly expressed in eq. 42 in TS15)

\begin{equation} \bar{\mu}_\mathrm{DL} = \mu^0 + RT\ln a_\mathrm{DL} + zF\Psi_\mathrm{DL} \tag{8} \end{equation}

where index “DL” refers to quantities in the “diffuse layer”.

However, the corresponding flux expression contains a different potential, labelled \(^\mathrm{DL}\Psi_\mathrm{diff}\) (eq. 46 in TS15)

\begin{equation} j_\mathrm{DL} = -\frac{c_\mathrm{DL}D_\mathrm{DL}}{RT}\nabla \bar{\mu}_\mathrm{DL} -\frac{c_\mathrm{DL}D_\mathrm{DL}zF}{RT} \nabla ^{\mathrm{DL}}\Psi_\mathrm{diff} \tag{9} \end{equation}

TS15 don’t further comment what \(^{\mathrm{DL}}\Psi_\mathrm{diff}\) is supposed to represent, but it must reasonably be understood as “the diffusion potential that arises because of the diffusion of charged species at different rates”, in analogy with what was claimed for the bulk water phase. Note that when eq. 8 is combined with eq. 9, the flux expression contains two different electric potential gradients! (TS15 never address this oddity.)

It is thus quite clear that TS15 misunderstand the function of the electric potential in the Nernst-Planck framework. When presenting the expression for the “diffuse layer” flux (eq. 9), they also refer to Appelo and Wersin (2007), who, in turn, express the misconception explicitly⁹

The gradient of the electrical potential [in the expression for the flux] originates from different transport velocities of ions, which creates charge and an associated potential. This electrical potential may differ from the one used in [the expression for the electro-chemical potential], which comes from a charged surface and is fixed, without inducing electrical current.

I cannot understand this passage in any other way than that Appelo and Wersin (2007) are under the impression that different electric potentials can simultaneously act independently in a given point. And it seems like TS15 are under some similar impression.

This ignorance leads to more errors in the description of the “diffuse layer” in TS15. We should remember that the promoted model requires the “diffuse layer” and bulk water domains to be in equilibrium (for the same coordinate value \(x\)). When TS15 express this condition, i.e. \(\bar{\mu}_\mathrm{DL} = \bar{\mu}_\mathrm{bulk}\), they again leave out the electric potential in the bulk water (eq. 42 in TS 15)

\begin{equation} \mu^0 + RT\ln a_\mathrm{DL} + zF\Psi_\mathrm{DL} = \mu^0 + RT\ln a_\mathrm{bulk}\;\:\;\;\;\;\;\mathrm{(WRONG)} \tag{10} \end{equation}

eq. 10 can be rewritten

\begin{equation} a_\mathrm{DL} = a_\mathrm{bulk}\cdot e^{-\frac{zF}{RT}\Psi_\mathrm{DL}} \;\:\;\;\;\;\;\mathrm{(WRONG)} \tag{11} \end{equation}

TS15 utilize a simplified version of eq. 11, expressed in terms of concentrations rather than activities, by assuming identical activity coefficients in the two domains¹⁰

\begin{equation} c_\mathrm{DL} = c_\mathrm{bulk}\cdot e^{-\frac{zF}{RT}\Psi_\mathrm{DL}} \;\:\;\;\;\;\;\mathrm{(WRONG)} \tag{12} \end{equation}

Note that the exponential in eqs. 11 and 12 actually should contain the electric potential difference between “diffuse layer” and bulk (see below).

As TS15 have not included any electric potential in the bulk water phase, they continue by incorrectly substituting \(RT\nabla \ln a_\mathrm{bulk}\) for \(\nabla\bar{\mu}_\mathrm{DL}\) in eq. 9 (i.e. they use the incorrect relation in eq. 10), giving (TS15 eq. 47)

\begin{equation} j_\mathrm{DL} = -c_\mathrm{DL}D_\mathrm{DL}\nabla \ln a_\mathrm{bulk} -\frac{c_\mathrm{DL}D_\mathrm{DL}zF}{RT} \nabla ^{\mathrm{DL}}\Psi_\mathrm{diff} \;\;\;\;\mathrm{(WRONG)} \tag{13} \end{equation}

Note that this additional error “solves” the earlier pointed out problem of having two electric potential gradients.

By utilizing the requirement of zero electric current, eq. 13 gives

\begin{equation} \nabla ^{\mathrm{DL}}\Psi_\mathrm{diff} = -\frac{RT}{F}\frac{\sum z_iD_{\mathrm{DL},i}c_{\mathrm{DL},i} \nabla \ln a_{\mathrm{bulk},i} }{\sum z_i^2 D_{\mathrm{DL},i} c_{\mathrm{DL},i}} \;\:\;\;\;\;\;\mathrm{(WRONG)} \tag{14} \end{equation}

By substituting eq. 12 into this expression, we end up with the formula for the gradient of the mysterious potential \(^{\mathrm{DL}}\Psi_\mathrm{diff}\) (TS15 eq. 48)

\begin{equation} \nabla ^{\mathrm{DL}}\Psi_\mathrm{diff} = \end{equation} \begin{equation} -\frac{RT}{F}\frac{\sum z_iD_{\mathrm{DL},i} e^{-\frac{zF}{RT}\Psi_\mathrm{DL}} \left ( \nabla c_{\mathrm{bulk},i} + c_{\mathrm{bulk},i}\nabla \ln \gamma_{\mathrm{bulk},i} \right )} {\sum z_i^2 D_{\mathrm{DL},i} e^{-\frac{zF}{RT}\Psi_\mathrm{DL}}c_{\mathrm{bulk},i}} \;\mathrm{(WRONG)} \tag{15} \end{equation}

At face value, eq. 15 is a quite weirdly looking equation, as it relates two electric potentials — \(^{\mathrm{DL}}\Psi_\mathrm{diff}\) and \(\Psi_\mathrm{DL}\) — that both are supposed to be associated with the “diffuse layer”. But, as we will see below, there is actually a way to make some sense of eq. 15, by completely reinterpreting what these potentials represent.

A “correct” formulation

Most of the errors pointed out above are corrected by including the electric potential in the bulk water and writing the condition for equilibrium as (compare eq. 10)

\begin{equation} \mu^0 + RT\ln a_\mathrm{DL} + zF\Psi_\mathrm{DL} = \mu^0 + RT\ln a_\mathrm{bulk} + zF\Psi_\mathrm{bulk} \tag{16} \end{equation}

Writing the electric potential difference between “diffuse layer” and bulk water as¹¹

\begin{equation} \Psi^\star \equiv \Psi_\mathrm{DL} – \Psi_\mathrm{bulk} \tag{17} \end{equation}

eq. 16 can be rewritten

\begin{equation} a_\mathrm{DL} = a_\mathrm{bulk}\cdot e^{-\frac{zF}{RT}\Psi^\star} \tag{18} \end{equation}

Note that, when correctly derived, eq. 18 naturally contains the difference in electric potential between “diffuse layer” and bulk.

The flux in the “diffuse layer” is (eq. 6)

\begin{equation} j_\mathrm{DL} = -c_\mathrm{DL} D_\mathrm{DL} \nabla \ln a_\mathrm{DL} – \frac{c_\mathrm{DL} D_\mathrm{DL} zF}{RT} \nabla \Psi_\mathrm{DL} \tag{19} \end{equation}

But if we now plug in eq. 18 in eq. 19 we of course get

\begin{equation} j_\mathrm{DL} = -c_\mathrm{DL} D_\mathrm{DL} \nabla \ln a_\mathrm{bulk} + \frac{c_\mathrm{DL} D_\mathrm{DL} zF}{RT} \nabla \Psi^\star – \frac{c_\mathrm{DL} D_\mathrm{DL} zF}{RT} \nabla \Psi_\mathrm{DL}, \end{equation} which can be simplified to \begin{equation} j_\mathrm{DL} = -c_\mathrm{DL} D_\mathrm{DL} \nabla \ln a_\mathrm{bulk} – \frac{c_\mathrm{DL} D_\mathrm{DL} Fz}{RT} \nabla \Psi_\mathrm{bulk}, \tag{20} \end{equation} and, by identifying the electro-chemical potential in the bulk \begin{equation} j_\mathrm{DL} = -\frac{c_\mathrm{DL} D_\mathrm{DL}}{RT} \left ( RT \nabla \ln a_\mathrm{bulk} + Fz \nabla \Psi_\mathrm{bulk} \right ) = -\frac{c_\mathrm{DL} D_\mathrm{DL}}{RT} \nabla \bar{\mu}_\mathrm{bulk} \end{equation}

This whole “derivation” leads back to the rather trivial result that the flux in the diffuse layer is given by eq. 2, which we could have written down from the start! (because the model assumes \(\bar{\mu}_\mathrm{bulk} = \bar{\mu}_\mathrm{DL}\); eq. 16)

As TS15 have established the expression for the gradient of the electrochemical potential in the bulk water phase (which is implicit in their eq. 40), there should strictly be no need to consider a new expression for the same quantity in any other phase. Rather, they could simply have used the bulk water expression in all “porosity domains”, as a consequence of the assumption that these are all supposed to be in equilibrium. In a sense, this is actually what is done in TS15 — mainly by chance! — by establishing eq. 13 (their eq. 47).

Comparing with eq. 20, we see that the incorrect eq. 13 can be “saved” by reinterpreting \(^{\mathrm{DL}}\Psi_\mathrm{diff}\) as \(\Psi_\mathrm{bulk}\). Similarly, as TS15 assume the bulk electric potential to be zero, eq. 15 can be “saved” by also reinterpreting \(\Psi_\mathrm{DL}\) as \(\Psi^\star\) in that expression.¹² I find this quite hilarious: By making several errors in its derivation, eq. 15 is in a sense a correct expression for the electric potential gradient in the bulk water — a potential that TS15 has put identically equal to zero.

But even if the total flux in the “diffuse layer” is correctly given by combining eqs. 12, 13 and 15 (and by completely ignoring what TS15 mean \(\Psi_\mathrm{DL}\) and \(^{\mathrm{DL}}\Psi_\mathrm{diff}\) represent), TS15 continue by defining the separate terms in eq. 13 as contributions from the “concentration gradient”, and the “diffusion potential”. As we will explore next, this interpretation fails miserably.

“Relative contributions of concentration, activity coefficient and diffusion potential gradients to total flux”

According to TS15, the “concentration gradient” and the “diffusion potential” contributions to the “diffuse layer” flux are given by, respectively (TS15 eq. 50 and below)

\begin{equation} j_\mathrm{conc,TS15} = -D_\mathrm{DL}A\nabla c_\mathrm{bulk} \tag{21} \end{equation}

and

\begin{equation} j_\mathrm{E,TS15} = zD_\mathrm{DL}c_\mathrm{bulk}A\frac{\sum z_iD_{\mathrm{DL},i}A_i\left ( \nabla c_{\mathrm{bulk},i} + c_{\mathrm{bulk},i} \nabla \ln \gamma_{\mathrm{bulk},i} \right )} {\sum z_i^2 D_{\mathrm{DL},i} A_i c_{\mathrm{bulk},i}} \tag{22} \end{equation}

Here we use the index “conc” for the “concentration gradient” contribution, and “E” for the “diffusion potential” contribution. \(A\) is referred to as a “DL enrichment factor”, and is essentially defined as the concentration ratio \(c_\mathrm{DL}/c_\mathrm{bulk}\). Using the incorrect relation in eq. 12, TS15 write these as \(A = e^{-\frac{zF}{RT}\Psi_\mathrm{DL}}\), but, as we see from eq. 18, they are really given by¹³ (we continue assuming identical activity coefficients in the two domains)

\begin{equation} A = e^{-\frac{zF}{RT}\Psi^\star} \tag{23} \end{equation}

TS15 also define a third contribution, related to the gradient of the bulk water activity coefficient. Here we will not further discuss this contribution, as it does not give any additional insight. Moreover, since TS15 anyway derive their model under the unjustified assumption that activity coefficients in the “diffuse layer” and the bulk water are identical, I cannot see the use of including their spatial variation in the description.¹⁴ (TS15 spend a couple of pages on activity coefficient models that we will ignore.)

Examples

To explore the various couplings in the presented model, TS15 apply the Nernst-Planck framework in three examples. We can see immediately from the presented graphs that their partitioning of the total flux in “concentration gradient” and “diffusion potential” contributions makes no sense.

“Example 2” imposes constant concentration gradients in the bulk water of NaCl and corresponding \(^{22}\mathrm{Na}^+\) and \(^{36}\mathrm{Cl}^-\) tracers; the NaCl concentration drops from 0.1 M to 0.001 M, and the tracer concentrations drop from 10^-9 M to 10^-11 M (domain length is 10 mm).

The corresponding sodium and chloride tracer concentrations in the “diffuse layer” look like this¹⁵

These profiles make sense: bulk water ionic strength decreases with distance, but so do the tracer concentrations. For the process of accumulating \(^{22}\mathrm{Na}^+\) in the diffuse layer, these two effects oppose each other, resulting in a quite flat profile. We thus expect the corresponding “concentration gradient” contribution to the flux to be quite moderate, and to fall off with distance (as the profile flattens with distances). The corresponding flux graph presented in TS15, however, looks completely different¹⁶

This plot makes no sense: The “concentration gradient” contribution is seen to increase quite dramatically with distance, rather than falling off. The value of this contribution is also orders of magnitude too large, given the imposed sodium diffusion coefficient of 1.33⋅10^-10 m²/s. Moreover, the “concentration gradient” contribution is “compensated” by an equally nonsensical “diffusion potential” contribution. Note, for instance, that the “diffusion potential” contribution is negative, which implies that the corresponding electric field is supposed to be directed towards higher concentrations. This can certainly not be the case, as the electric potential gradient is caused by the negative ion having higher mobility than the positive ion (chloride diffusivity is set to 2.03⋅10^-10 m²/s).

In “example 3”, the tracer concentrations in the bulk is set to a constant value (1⋅10^-9 M), while the same concentration gradient as in “example 2” is maintained for the main NaCl electrolyte (from 0.1 M to 0.001 M). We thereby expect the corresponding \(^{22}\mathrm{Na}^+\) concentration in the “diffuse layer” to strongly increase with distance, which is also what is presented in TS15

while the corresponding flux plot looks like this¹⁶

This plot is almost comically absurd. According to TS15, the highly skewed concentration profile above is supposed to give no (zero, nil, 0) contribution to the flux (we see from eq. 21 that this is a consequence of that this “contribution” is directly proportional to the concentration gradient in the bulk). Instead, the huge flux is supposed to be caused entirely by an electric field that has the wrong direction! I can’t even really begin to imagine how these two plots have ended up next to each other in a peer-reviewed published article.

Note that the flux associated with a concentration gradient is what we may reasonably call a “Fickian” contribution. If TS15 mean (and they do) that these examples demonstrate how ion diffusion in bentonite works, we can understand the focus on the “Fickian” aspect at the beginning of the article (covered here). But the only reasonable response to these outlandish results is that they demonstrate that the definitions of eqs. 21 and 22 simply make no sense.

The real concentration gradient and electric field contributions

The only reasonable way to define “concentration gradient” and “diffusion potential” contributions to the “diffuse layer” flux is as the two terms in eq. 19, respectively. To rewrite these, we utilize eq. 16 (or 18), giving for the “concentration gradient” contribution (we continue ignoring activity coefficients)

\begin{equation} j_\mathrm{conc, corr} = -c_\mathrm{DL} D_\mathrm{DL} \nabla \ln c_\mathrm{DL} = \end{equation} \begin{equation} -c_\mathrm{DL} D_\mathrm{DL} \nabla \ln c_\mathrm{bulk} + \frac{c_\mathrm{DL} D_\mathrm{DL}zF}{RT} \nabla \Psi^\star = \end{equation} \begin{equation} -D_\mathrm{DL} A \nabla c_\mathrm{bulk} + \frac{A c_\mathrm{bulk} D_\mathrm{DL}zF}{RT} \nabla \Psi^\star = j_\mathrm{conc, TS15} + j^\star \tag{24} \end{equation}

where we have defined

\begin{equation} j^\star \equiv \frac{A c_\mathrm{bulk} D_\mathrm{DL}zF}{RT} \nabla \Psi^\star. \tag{25} \end{equation}

In the same manner, the correct “diffusion potential” contribution is

\begin{equation} j_\mathrm{E, corr} = – \frac{c_\mathrm{DL} D_\mathrm{DL} zF}{RT} \nabla \Psi_\mathrm{DL} = \end{equation} \begin{equation} – \frac{D_\mathrm{DL}Ac_\mathrm{bulk} zF}{RT} \nabla \Psi_\mathrm{bulk} – \frac{D_\mathrm{DL}Ac_\mathrm{bulk} zF}{RT} \nabla \Psi^\star = \end{equation} \begin{equation} D_\mathrm{DL}Ac_\mathrm{bulk} z \frac{\sum z_iD_iA_i\left ( \nabla c_{\mathrm{bulk},i} + \nabla c_{\mathrm{bulk},i} \ln \gamma_{\mathrm{bulk},i} \right )} {\sum z_i^2 D_i A_i c_{\mathrm{bulk},i}} – \frac{D_\mathrm{DL}Ac_\mathrm{bulk} zF}{RT} \nabla \Psi^\star = \end{equation} \begin{equation} j_\mathrm{E, TS15} – j^\star \tag{26} \end{equation}

where we have utilized that \(\nabla \Psi_\mathrm{bulk}\) is actually what is expressed in eq. 15 (where \(\Psi_\mathrm{DL}\) should be replaced by \(\Psi^\star\)).

We note that, to compensate the nonsensical expressions given in TS15, we should add the term \(j^\star\) (eq. 25) to the “concentration concentration” contribution (eq. 21), and subtract the same term from the “diffusion potential” contribution (eq. 22). Making these corrections gives the following components of the tracer fluxes in “example 2”

This is an infinitely more reasonable situation than what is depicted in TS15. Although the sodium flux has a non-negligible contribution from the electric field, the larger contribution is still from the concentration gradient (and none of these are gigantic terms that cancel). The concentration contribution also falls off with distance, in accordance with the shape of the concentration profile.

For chloride, the field contribution to the flux is negligible, i.e. this flux is essentially fully governed by the concentration gradient. The electric field contributions for both ions are also seen to have the correct signs: the electric field is directed from high to low concentration, and mainly functions to boost the sodium transport, in order to “keep up” with the faster chloride ions.

For “example 3” we get the following picture

The corrected \(^{22}\mathrm{Na}^+\) flux is essentially fully due to the concentration gradient, in absolute contrast to what is concluded in TS15, who mean that this flux is completely governed by an electric field in the wrong direction. Also the \(^{36}\mathrm{Cl}^-\) transport is basically solely governed by the concentration gradient, rather than by an incorrectly directed electric field (as stated in TS15). In conclusion, most of the “diffuse layer” diffusion in these examples can actually be classified as “Fickian”.

We may also investigate the electric potential profile in the “diffuse layer” in both of these examples (this is the same in the two cases, as the main electrolyte distribution does not change)

Here we have chosen the reference \(\Psi_\mathrm{DL}(0) = 0\). The total potential drop is only about 1 mV. Such a relatively small drop is reasonable because the denominator in the Nernst-Planck expression for the electric potential gradient (eq. 7) will always be large due to the ever-present counter-ions in the “diffuse layer”. The electric potential gradient — and thus the corresponding electric potential drop — is therefore suppressed. Physically, this means that since many (equally charged) charge carriers are always present, smaller potential differences are required to cancel electric currents caused by differences in mobility (a “diffuse layer” is a quite good conductor).

Even worse problems?

Even though some sense can be made out of the derived expression for the flux in the “diffuse layer” domain — by completely reinterpreting the electric potentials involved — it seems as the overall model is too constrained. Specifically, for an imposed set of concentration profiles in one domain it is not possible, as far as I can see, to simultaneously have zero current in all domains, while also maintaining (Donnan) equilibrium. As this blog post is already quite massive, I will elaborate on this point in the next part of the review.

Summary

Here is an attempt to sum up the main messages of this blog post.

• Conceptually, the clay model presented in TS15 is exactly what was discussed in the blog post on multi-porous models, and the same issues that are identified there are present here. In particular, no attention is paid to length scales (perhaps that is why TS15 call the coordinate system “pseudo-2D”…), and no mechanism whatsoever is suggested for how the different diffusing domains are supposed to maintain equilibrium.
  
• Mathematically (or perhaps physically), the presented Nernst-Planck flux expressions are incorrectly derived. The source of the error, as far as I can see, appears to be a misunderstanding of how electric potentials function.
  
• TS15 define “contributions” to the “diffuse layer” flux, claimed to be related to the concentration gradient and the “diffusion potential” (i.e. the electric field), respectively. It is, however, quite obvious that these “contributions” are completely nonsensical: highly skewed concentration profiles are claimed to not have any concentration gradient contributions, and several “diffusion potential” “contributions” have the electric field in the wrong direction. We have shown that these “contributions” can be corrected, where the correction term involves the gradient of the Donnan potential. With these corrections, fluxes in the provided examples must be interpreted completely differently (they’re basically “Fickian”).
  
• As far as I can see, the proposed model has even larger problems, related to the imposed Donnan equilibrium. We will address this issue in the next part.
  

Footnotes

[1] As I have commented in the earlier parts: TS15 are fond of using the general terms “clays” and “clay minerals”, while it is clear that the publication mainly focus on systems with substantial ion exchange capacity and swelling properties. Here we will continue to use the term “bentonite” for these systems, and ignore the frequent references in TS15 to more general terms.

[2] It is of course crucial to include a component that represents compartments where the exchangeable ions reside. This is done in the TS15 model by both the “diffuse layer water” and the “interlayer water” domains. But the distinction made between these domains is based on the flawed “stack” concept.

[3] This equation assumes a single component. The formulation of the Nernst-Planck framework naturally involves several different charged species. When several species are involved, we will indicate this with an index \(i\) in the equations.

[4] In some of their equations, TS15 use (electrical) mobility, \(u\), rather than diffusivity, \(D\). These quantities are related via the Einstein relation \(D = uRT/(F|z|)\). I don’t see the point in involving \(u\), as it typically makes expressions even more cluttered, and since we here ultimately are interested in diffusion coefficients anyway.

[5] In order to not cause too much confusion, and to try to simplify a bit, I use slightly different mathematical notation than what is actually used in the quotation. In particular, I use the notation \(\bar{\mu}\) for the electro-chemical potential, while TS15 don’t use a bar (\(\mu\)). I also try to avoid the index \(i\) as much as possible.

[6] Fun fact: this statement is nowhere found in neither (Ben-Yaakov, 1981) nor (Lasaga, 1981) (at least I can’t find any).

[7] I whined about electrostatics being poorly understood in the bentonite research field in an earlier part of this review, but here is more fuel for my argument. The statement “absence of an [external] potential” has no physical meaning, as we are free to choose the reference point (the absolute value of a potential has no physical meaning). What TS15 must mean in the quote is “the absence of an external electric field”. The electric field relates to the potential as \(E = -\nabla \Psi\). Thus, all gradients of electric potentials that occur in this text are synonymous with electric fields (electric fields drive electric currents).

[8] This post focus almost entirely on the “diffuse layer” domain, but a similar analysis can be made for the “interlayer” domain. This is left as an exercise for the reader.

[9] It should of course also rather read “…which creates a charge separation and an associated potential gradient.”, or simply “…which induces an electric field.” (showing that this part of the sentence is redundant). See also footnote 7.

[10] TS15 write cryptically that equating the activity coefficients (and the reference potentials) in bulk and “diffuse layer” is assumed “by following the [Modified Gouy-Chapman] model”. But I don’t see why this model has to be alluded to here, these assumptions can just be made.

[11] Yes, this is a Donnan potential. We will discuss this more in the next part part of the review.

[12] Again, this is related to Donnan equilibrium between the bulk and “diffuse layer” domains, that we will discuss further in the next part.

[13] This is \(f_D^{-z} \), where \(f_D\) is the Donnan factor.

[14] Rather, I would argue for that the activity coefficients in a “diffuse layer” domain will be quite insensitive to the imposed external (bulk) concentration, for details see Birgersson (2017).

[15] In producing these graphs we have used the Donnan equilibrium framework to calculate the “diffuse layer” concentrations. These are given from eq. 23, where \(\Psi^\star\) is calculated from

\begin{equation} f_D = e^\frac{F\Psi^\star}{RT} = – \frac{q}{2c_\mathrm{bulk}} +
\sqrt{\frac{q^2}{4c_\mathrm{bulk}^2} + 1} \end{equation}

where \(q\) is a measure of the structural charge in the “diffuse layer”, in the examples set to \(q\) = 0.33 M.

[16] Note that I have not included activity coefficient gradients when producing the plots in this section. They may therefore differ slightly from the published plots. This does not in any way influence the conclusions drawn here.