Next Article in Journal
E-Bayesian Estimation of Chen Distribution Based on Type-I Censoring Scheme
Next Article in Special Issue
Electric Double Layer and Orientational Ordering of Water Dipoles in Narrow Channels within a Modified Langevin Poisson-Boltzmann Model
Previous Article in Journal
Intermittent Information-Driven Multi-Agent Area-Restricted Search
Previous Article in Special Issue
Molecular Mean-Field Theory of Ionic Solutions: A Poisson-Nernst-Planck-Bikerman Model
 
 
Font Type:
Arial Georgia Verdana
Font Size:
Aa Aa Aa
Line Spacing:
Column Width:
Background:
Article

Review and Modification of Entropy Modeling for Steric Effects in the Poisson-Boltzmann Equation

Department of Applied Mathematics, Feng Chia University, Taichung 40724, Taiwan
Entropy 2020, 22(6), 632; https://doi.org/10.3390/e22060632
Submission received: 11 May 2020 / Revised: 3 June 2020 / Accepted: 4 June 2020 / Published: 8 June 2020

Abstract

:
The classical Poisson-Boltzmann model can only work when ion concentrations are very dilute, which often does not match the experimental conditions. Researchers have been working on the modification of the model to include the steric effect of ions, which is non-negligible when the ion concentrations are not dilute. Generally the steric effect was modeled to correct the Helmholtz free energy either through its internal energy or entropy, and an overview is given here. The Bikerman model, based on adding solvent entropy to the free energy through the concept of volume exclusion, is a rather popular steric-effect model nowadays. However, ion sizes are treated as identical in the Bikerman model, making an extension of the Bikerman model to include specific ion sizes desirable. Directly replacing the ions of non-specific size by specific ones in the model seems natural and has been accepted by many researchers in this field. However, this straightforward modification does not have a free energy formula to support it. Here modifications of the Bikerman model to include specific ion sizes have been developed iteratively, and such a model is achieved with a guarantee that: (1) it can approach Boltzmann distribution at diluteness; (2) it can reach saturation limit as the reciprocal of specific ion size under extreme electrostatic conditions; (3) its entropy can be derived by mean-field lattice gas model.

1. Introduction

One of the major limitations of the Poisson-Boltzmann (PB) and Poisson-Nernst-Planck (PNP) models is the assumption of point-like ions without considering their sizes. These models based on mean field theories work well for dilute electrolytes, but break down when the concentration is high and ions are crowded in it. A high concentration would generally cause steric repulsions and additional electrostatic correlations among ions, that cannot be described by classical PB/PNP models [1]. For example, the concentration of counter-ions, predicted by PB, can be unrealistically high near the electrode surface, when the electrode voltage is large. Another example occurs at the selectivity filter in a potassium channel, where potassium ions are strongly attracted into this extremely narrow filter by the strong negative charges of oxygens on the backbone of the filter. Employing classical PB/PNP models would overestimate the density of potassium inside the filter and give incorrect channel current predictions. Therefore, many researchers have worked on the modification of PB/PNP to include the steric effect of ions.
Steric effect has long been approached in modeling by modifying either the internal energy or entropy in the Helmholtz free energy. Through internal energy, the steric effect has been featured as excess hard-sphere energy either by density functional theory (DFT) [2,3] or Lennard-Jones potential [4]. These energies were all formulated using non-local potentials and cause the resultant modified PB/PNP to produce a series of complicated integro-differential equations, which are hard to compute in higher dimensions. For practical implementations, localization of hard-sphere potential and simplifying integro-differential equations into pure differential equations has been conducted in [5,6,7,8] for DFT and [9] for Lennard-Jones potential.
Through an entropy approach, Bikerman modified the classical Boltzmann distribution by adjusting bulk and local ion concentrations via the excluded volume concept [10]. Borukhov et al. [11] rigorously derived the same formula independently by adding solvent entropy through excluded volume into the Helmholtz free energy. Although the localized hard-sphere model-based DFT [5,6,7,8] also captures this solvent entropy as one of the terms accounting for excess hard-sphere chemical potential, the Bikerman model [10,11] has been a more popular steric model due to its easiness of application and qualitatively good agreement with experiments [12,13,14,15].
In order to obtain the potential and further derive a neat modified PB equation from free energy, Borukhov et al. [11] treated all ions as having identical size, which has been long criticized for neglecting specific ion sizes. Researchers have tried to address this shortcoming with specific ion sizes, and many of them simply extended original Bikerman model by replacing the identical ion size with specific ones without any rigorous justification. Although the resultant model has a better agreement with experiments than the original Bikerman model [16,17,18], it does not have a Helmholtz free energy to support it. Here modifications of the Bikerman model to include specific ion sizes have been developed iteratively in Section 4, Section 5 and Section 6, preceded by derivations of classical PB in Section 2, and the original Bikerman model in Section 3. Finally, in the Discussion and Conclusions section a specific-ion-size Bikerman model is presented with a guarantee that: (1) it can approach Boltzmann distribution at diluteness; (2) it can reach the saturation limit as the reciprocal of specific ion size under extreme electrostatic conditions; (3) its entropy can be derived by a mean-field lattice gas model.

2. Classical Poisson-Boltzmann Model

Though the classical PB model is well known, we still derive the model here for review and comparison with its modified versions discussed later. Starting by stating the Helmholtz free energy, internal energy and entropy, we have:
F = U T S
U = [ ε 2 | ϕ | 2 + z p e p ϕ + z n e n ϕ + q ϕ + p W s o l , p + n W s o l , n ] d V ,
T S = k B T [ p   log p c 0 p + n   log n c 0 n ] d V ,
where F is Helmholtz free energy; U is internal energy; T is temperature; S is entropy; ϕ is electric potential. p, n denote cation/anion concentrations, and z p , z n denote their valence, respectively. e denotes elementary charge. q denotes permanent charge. Permittivity ε = ε 0 ε r with ε 0 being the permittivity for vacuum and ε r being the relative permittivity or dielectric constant. c 0 is some reference concentration such as bulk concentration of electrolyte. W s o l , p and W s o l , n denote the solvation energies for cations and anions, respectively. Although the traditional PB model generally does not include solvation energy in the expression, it is important when modeling some electrolyte systems involving hydration/dehydration of ions and is therefore it is explicitly included in the energy here. Based on the Born model, the solvation energies for cations and anions are:
W s o l , i = z i 2 e 2 8 π ε 0 r i ( 1 ε r ( x ) 1 ) ,   i = p , n .
Differentiation of F with respect to ϕ gives the Poisson equation:
· ( ε ( x ) ϕ ) = z p e p + z n e n + q
By doing the differentiation of F with respect to p and n, we obtain the chemical potentials for p and n, respectively:
F p = μ p = z p e ϕ + k B T   log p c 0 + W s o l , p ,
F n = μ n = z n e ϕ + k B T   log n c 0 + W s o l , n
At equilibrium, the chemical potential is uniform everywhere and therefore the local chemical potential must be equal to its bulk value, which is usually known:
μ p = μ p , b ,   μ n = μ n , b ,
with the bulk chemical potential for cations and anions:
μ p , b = k B T   log ( p b c 0 ) + W s o l , p , b , μ n , b = k B T   log ( n b c 0 ) + W s o l , n , b ,
where the subscript b denotes the bulk situation. Equations (8) and (9) can be solved solve for p and n:
p = p b e β z p e ϕ e β Δ W s o l , p = p b e β E p ,   n = n b e β z n e ϕ e β Δ W s o l , n = n b e β E n ,
where β = 1 / k B T , E p = z p e ϕ + Δ W s o l , p , E n = z n e ϕ + Δ W s o l , n , and:
Δ W s o l , i = W s o l , i W s o l , i , b = z i 2 e 2 8 π ε 0 r i ( 1 ε r ( x ) 1 ε r , b ) ,   i = p , n .
From (10), as ϕ , we obtain p ,   n 0 . Likewise, as ϕ , we obtain p 0 ,   n . These unrealistic infinite concentrations for p and n are mainly because ions are treated as particles without size in the classical PB model. This pitfall has motivated modifications of the classical PB/PNP model to account for the finite-size effect, or so-called steric effect, of ions. In reality, the limit of p should be at most 1 / v p , where v p is the particle volume of p . This can be derived by considering a volume V fully occupied by cation p only, with the number of cation particles being N p , and then:
p m a x = N p V = N p N p v p = 1 v p .
Likewise, the limit of n is at most 1 / v n , where v n is the particle volume of n .
Substituting (10) into (4), and we obtain the classical PB equation:
· ( ε ( x ) ϕ ) = z p e p b e β E p + z n e n b e β E n + q .
For z:z electrolyte without considering solvation energy, the equation above reduces to:
· ( ε ( x ) ϕ ) = 2 z e c b s i n h ( β z e ϕ ) q ,
with c b = p b = n b .

3. Bikerman Model

As stated earlier, the Bikerman model [10] has been a popular steric-effect model due to its easiness of application and qualitatively good agreement with experimental data. It modifies the free energy of the classical PB (1)–(3) by adding a solvent entropy term. This term also partially represents the excessive energy accounting for overcrowding of ions and solvent molecules in localized hard-sphere models based on DFT [5,6,7,8]. The free energy in the Bikerman model treats all species of ions and solvent molecules with an identical size, and is stated as follows:
F = U T S ,
U = [ ε 2 | ϕ | 2 + z p e p ϕ + z n e n ϕ + q ϕ + p W s o l , p + n W s o l , n ] d V ,
T S = k B T [ p   log ( p v ) p + n   log ( n v ) n + w   log ( w v ) w ] d V ,
where w is concentration of solvent (such as water); v is the universal particle volume. Why are all solute and solvent particles treated as having the same size? Why are specific sizes of ions and solvent molecules not used here? This was not explained in the original model [10,11], and the justification of using identical size for all species particles will be addressed later in Section 5.
If we assume that, besides occupation of ions, the rest of space is occupied by solvent molecules (which can be taken as water here). Then:
w = N w V = N w v V v = V N p v N n v V v = 1 v ( 1 p v n v ) ,
where V is the whole volume of electrolyte; N p ,   N n ,   N w are number of cation, anion and solvent particles in an electrolyte with volume V, respectively. Equation (18) can then be rewritten as:
w v + p v + n v = 1 ,
which simply means the sum of volume fractions of water, cation and anion is one. Note that here we assume that, besides occupation of ions, the rest of space is occupied by water.
Substituting (19) into (17) we can obtain:
T S = k B T [ p   log ( p v ) p + n   log ( n v ) n + 1 v ( 1 p v n v )   log ( 1 p v n v ) 1 v ( 1 p v n v ) ] d V .
Differentiation of F with respect to ϕ again gives the Poisson equation:
· ( ε ( x ) ϕ ) = z p e p + z n e n + q .
By doing the derivation of F with respect to p and n, we obtain the chemical potentials for p and n, respectively:
F p = μ p = z p e ϕ + k B T [ l o g ( p v ) l o g ( 1 p v n v ) ] + W s o l , p ,
F n = μ n = z n e ϕ + k B T [ l o g ( n v ) l o g ( 1 p v n v ) ] + W s o l , p .
At equilibrium, the chemical potential is uniform everywhere and therefore the local chemical potential must be equal to its bulk value, which is usually known:
μ p = μ p , b ,         μ n = μ n , b ,
with:
μ p , b = k B T [ log ( p b v ) log ( 1 p b v n b v ) ] + W s o l , p , b ,
μ n , b = k B T [ log ( n b v ) log ( 1 p b v n b v ) ] + W s o l , n , b .
By substituting (22), (23), (25) and (26) into (24), we can relate the local ion-to-solvent volume fraction ratios (denoted as γ i ,   i = p , n . ) to their counterparts in bulk solution in a Boltzmann manner for p and n, respectively:
γ p = p v 1 p v n v = p b v 1 p b v n b v e β E p ,
γ n = n v 1 p v n v = n b v 1 p b v n b v e β E n .
Summation of (27) and (28) gives the solute-to-solvent volume fraction ratio as:
p v + n v 1 p v n v = γ p + γ n ,
and we can further obtain the solute volume fraction:
p v + n v = γ p + γ n 1 + γ p + γ n .
From (29) and (19), we know then volume fraction for p, n and w, respectively.
p v = γ p 1 + γ p + γ n ,   n v = γ n 1 + γ p + γ n ,   w v = 1 1 + γ p + γ n ,  
which further gives p and n in terms of their bulk values p b , n b , identical particle size v, and local energy E p , E n :
p = p b e β E p ( 1 p b v n b v ) + p b v e β E p + n b v e β E n ,     n = n b e β E n ( 1 p b v n b v ) + p b v e β E p + n b v e β E n
Equations (27) and (28) can be re-arranged to obtain:
p = p b e β ( E p + S t r c ) ,     n = n b e β ( E n + S t r c ) ,
with ionic steric potential S t r c expressed as:
S t r c = k B T   log ( 1 p b v n b v 1 p v n v ) .
The steric potential S t r c , first described in [16,17,18,19], characterizes the crowding of ions and their finite-size effect by a bulk-to-local water fraction ratio. Larger local ion concentrations would have a larger steric potential.
Also, by letting 1 p v n v = w v , and 1 p b v n b v = w b v , Equations (27) and (28) can be simplified as:
p v w v = p b v w b v e β E p = γ p ,
n v w v = n b v w b v e β E n = γ n .
Therefore, by (30), we obtain:
z p p v + z n n v = z p p b v w b v e β E p + z n n b v w b v e β E n 1 + p b v w b v e β E p + n b v w b v e β E n = z p p b v e β E p + z n n b v e β E n w b v + p b v e β E p + n b v e β E n .
Since p b n b = z n z p due to electric neutrality in bulk conditions, therefore:
p b v = z n μ z p z n v , n b v = z p μ z p z n v ,
where μ = p b + n b . Also:
w b v = 1 p b v n b v = 1 μ v .
then (36) becomes:
z p p + z n n = z p z n μ ( e β E p + e β E n ) ( 1 μ v ) ( z p z n ) + μ v ( z p e β E n z n e β E p ) .
Substituting (39) into (21), we obtain the Bikerman-PB equation:
· ( ε ( x ) ϕ ) = z p z n e μ ( e β E p + e β E n ) ( 1 μ v ) ( z p z n ) + μ v ( z p e β E n z n e β E p ) + q .
For z:z electrolyte without considering the solvation energy, the equation above becomes:
· ( ε ( x ) ϕ ) = 2 z e c 0 s i n h ( β z e ϕ ) 1 + 2 r s i n h 2 ( β z e ϕ 2 ) q ,
as shown in [11] with c 0 = p b = n b , r = μ v = 2 c 0 v .
Two important criteria need to be checked for all modified PB/PNP models accounting for steric effects:
CRITERION I: When ion concentrations p and n are dilute, will they follow the classical Boltzmann distribution?
CRITERION II: As ϕ , will p and n approach their saturation limits 1 v p and 1 v n , respectively?
For CRITERION I when p and n are dilute here, it means their volume fractions are negligible, and therefore 1 p v n v 1 , and 1 p b v n b v 1 . Steric potential term S t r c then vanishes, and by (32) p = p b e β E p , and n = n b e β E n , which follows the Boltzmann distribution.
For CRITERION II, let us consider ϕ first, and ϕ can be derived similarly. As ϕ , γ p , and γ n 0 . Therefore, p v 1 , and n v 0 by (30), which further means p 1 v p = 1 v , and n 0 . Likewise, as ϕ , we can get n 1 v n = 1 v , and p 0 .

4. The Bikerman Model with Specific Ion Sizes

The shortcoming of the Bikerman model is the usage of a universal particle size, denoted by v, for cations, anions and solvents. Using specific ion and solvent sizes would be closer to reality. Taking NaCl solution as an example, the spherical diameters for Cl, Na+ and water are D C l = 3.62 , D N a + = 2.04 , and D w = 2.08 , and then the particle volume ratio is v N a + : v C l : v w = 1 : 5.59 : 1.06 , in which using universal particle volume would be far from reality in the case of high ion concentrations. In appearance, it seems, and many researchers did, we can just simply modify the model to include specific ion sizes by changing p v and p b v to p v p and p b v p ; similarly, n v and n b v to n v n and n b v n for (22) to (41). With this straightforward extension, we obtain p, n as:
p = p b e β E p ( 1 p b v p n b v n ) + p b v p e β E p + n b v n e β E n ,
n = n b e β E n ( 1 p b v p n b v n ) + p b v p e β E p + n b v n e β E n
and the specific-ion-size Bikerman-PB equation:
· ( ε ( x ) ϕ ) = e ( z p p b e β E p + z n n b e β E n ) 1 ( p b v p + n b v n ) + ( p b + n b ) ( z p v n e β E n z n v p e β E p ) / ( z p z n ) + q
For z:z electrolyte without considering solvation energy, the equation above becomes:
· ( ε ( x ) ϕ ) = 2 z e c 0 s i n h ( β z e ϕ ) 1 c 0 ( v p + v n ) + c 0 ( v n e β z e ϕ + v p e β z e ϕ ) q .
Let us denote (42) as the specific-ion-size Bikerman model 1 (SISBM1) for convenience of notation. However, we can not find an energy functional like (15)–(17) to support this naive extension, which means chemical potentials (22) and (23) with universal particle size replaced by specific ion sizes cannot be derived. The correct specific-ion-size energy functional and chemical potentials should be derived as follows:
F = U T S ,  
U = [ ε 2 | ϕ | 2 + z p e p ϕ + z n e n ϕ + q ϕ + p W s o l , p + n W s o l , n ] d V ,  
T S = k B T [ p   log ( p v p ) p + n   log ( n v n ) n + w   log ( w v w ) w ] d V .
By w v w = 1 p v p n v n , (47) can be rewritten as:
T S = k B T [ p   log ( p v p ) p + n   log ( n v n ) n + 1 v w ( 1 p v p n v n )   log ( 1 p v p n v n ) 1 v w ( 1 p v p n v n ) ] d V .
Differentiation of F with respect to ϕ again gives the Poisson equation:
· ( ε ( x ) ϕ ) = z p e p + z n e n + q .
By doing the differentiation of F with respect to p and n, we can obtain the chemical potentials for cations and anions, respectively:
μ p = z p e ϕ + k B T [ log ( p v p ) log ( 1 p v p n v n ) k p ] + W s o l , p
μ n = z n e ϕ + k B T [ log ( n v n ) log ( 1 p v p n v n ) k n ] + W s o l , n
where k p = v p v w , k n = v n v w .
At equilibrium, the chemical potential is uniform everywhere and therefore the local chemical potential must be equal to its bulk value, which is usually known:
μ p = μ p , b ,     μ n = μ n , b ,
with:
μ p , b = k B T [ log ( p b v p ) log ( 1 p b v p n b v n ) k p ] + W s o l , p , b ,  
μ n , b = k B T [ log ( n b v n ) log ( 1 p b v p n b v n ) k n ] + W s o l , n , b .
To solve p and n from (52), there is no closed form solution like (31) for p and n due to the nonlinearity, unless some simplified case such as k p = k n = 1 is considered, which is actually reduced to the original Bikerman model with v p = v n = v w = v . Like (32), p and n at most can be expressed as:
p = p b e β ( E p + k p S t r c ) ,     n = n b e β ( E n + k n S t r c ) ,
with the steric potential S t r c being modified from (33) to include specific ion sizes:
S t r c = k B T   log ( 1 p b v p n b v n 1 p v p n v n )
Let us denote (55) as the specific-ion-size Bikerman model 2 (SISBM2) for convenience. Note that a similar model was also obtained in [18,19] without a rigorous derivation.
Again, we need to check criteria I and II for this specific-ion-size model. For CRITERION I, when p and n are dilute, it again means 1 p v p n v n 1 ,  and 1 p b v p n b v n 1 . Therefore S t r c 0 by (56), and then p = p b e β E p , and n = n b e β E n by (55), which follows a classical Boltzmann distribution. For CRITERION II, as ϕ in (50), k B T [ log ( p v p ) log ( 1 p v p n v n ) k p ] should approach + for μ p to be finite. This can only be achieved by n 0 , and p ( 1 v p ) (saturation). Applying the same reasoning for (51), as ϕ , k B T [ log ( n v n ) log ( 1 p v p n v n ) k n ] should approach + for μ n to be finite. Then p 0 , and n ( 1 v n ) (saturation). This specific-ion-size model seems correct and reasonable so far, but actually there is a pitfall. That is its entropy formula (48) cannot be derived by the traditional mean-field lattice gas model. This will be explained in the next section.

5. Mixing Entropy Derivation Based on the Mean-Field Lattice Gas Model

In this section, we would like to derive the entropy in (20) by the traditional mean-field lattice gas model. Consider the entropy for an aqueous electrolyte system:
T S = k B T   log   W ,
where W is the number of microstates at equilibrium which possess a maximum number of microstates. Mixing entropy in electrolyte studies macrostates through spherical particles’ (solute and solvent) occupation of identical cubic sites is based on the mean-field lattice gas model. The necessity of using identical cubic sites provides a combinatorial basis when computing the maximum number of microstates. The most probable distribution of all solute (ions) and solvent particles, reaching maximum number of microstates for each species, is that each identical cubic site generally would be at most occupied by one solute/solvent particle as depicted in Figure 1a. This is based on the concept that the size of each species’ particle is infinitesimal or finite but dilute. When the actual size for each species’ particle is considered and an aqueous electrolyte is extremely concentrated as depicted in Figure 1b, the most probable distribution above may not be available. The situation in Figure 1b will be addressed in the next section.
The entropy based on the most probable distribution of K-species solute (ions) and solvent (treated as K + 1 -th species) particles, under dilute situation, over a total of N = j = 1 K + 1 N j available identical sites in a system is:
W = j = 1 K + 1 W j = ( N N 1 ) ( N N 1 N 2 ) ( N N 1 N 2 N K 1 N K ) ( N N 1 N 2 N K N K + 1 ) = N ! N 1 ! ( N N 1 ) ! ( N N 1 ) ! N 2 ! ( N N 1 N 2 ) ! ( N j = 1 K 1 N j ) ! N K ! ( N j = 1 K N j ) ! ( N j = 1 K N j ) ! N K + 1 ! = N ! ( j = 1 K N j ! ) N K + 1 ! ,
where N j ,   j = 1 , , K , is the particle number of j-species ion. N K + 1 is the particle number of solvent, so the entropy becomes:
T S = k B T   log N ! ( j = 1 K N j ! ) N K + 1 ! .
Using the Stirling formula log M ! M log M M with M 1 , we can rewrite the entropy as:
T S = k B T   [ N log N N j = 1 K N j log N j + j = 1 K N j N K + 1 log N K + 1 + N K + 1 ] = k B T [ N log N j = 1 K N j log N j N K + 1 log N K + 1 ] = k B T [ N log N j = 1 K N j log N j ( N j = 1 K N j ) log ( N j = 1 K N j ) ] = k B T [ N log N N j = 1 K N j j = 1 K N j log N j N j = 1 K N j ]
using the following relations:
V = N v s ,   o r   N V = 1 v s ,
N j V = c j ,
N j N = N j v s N v s = N j v s V = c j v s ,
where c j is the concentration of j-species particle; V is the volume of system; v s is the volume of an identical cubic site that composes the volume of system. It is naturally requested that v s m a x 1 j K + 1 v j , where v j is the particle volume of j-species particle. Usually v s = m a x 1 j K + 1 v j in aqueous electrolyte system, where solute and solvent particles are generally crowded.
Applying (61)–(63) to (60), the entropy density can be expressed as:
T S V = k B T [ 1 v s log 1 1 j = 1 K c j v s j = 1 K c j log c j v s 1 j = 1 K c j v s ] .
For a binary electrolyte, (64) can be expressed as:
T S = k B T [ 1 v s log 1 1 p v s n v s p   log p v s 1 p v s n v s n   log n v s 1 p v s n v s ] d V ,
or:
T S = k B T [ p   log ( p v s ) + n   log ( n v s ) + 1 v s ( 1 p v s n v s ) log ( 1 p v s n v s ) ] d V ,
which, without loss of generality, can be augmented as:
T S = k B T [ p   log ( p v s ) p + n   log ( n v s ) n + 1 v s ( 1 p v s n v s ) log ( 1 p v s n v s ) 1 v s ( 1 p v s n v s ) ] d V .  
Equation (67) is exactly the same as (20) with:
v = v s = m a x { v p , v n , v w } .
This means the universal particle volume v in the Bikerman model is actually the volume of an identical occupation site v s , which is limited from below by the largest particle size among all solute and solvent particles. The original Bikerman model has long suffered criticism for assuming all ions have the same size instead of using specific ion sizes in the model. The above reasoning explains why specific ion size information is left out mainly due to the need for all cubic sites to be identical in order to support the combinatorial basis demanded by the mean-field lattice gas model. Actually, information of specific ion sizes is still carried but only implicitly as shown in (68). Researchers may prefer to use SISBM2 as illustrated in Section 4, but actually its entropy formula (48) cannot be derived by the mean-field lattice gas model described above. Note that usually solute and solvent particles are treated as spheres in modeling. If a p , a n , a w are diameters for p, n, and w, respectively and their maximum is a n for example, then v = v s = a n 3 not 4 π 3 ( a n 2 ) 3 since the identical occupation is cubic. This is why a 3 , instead of 4 π 3 ( a 2 ) 3 , used in [10,11].
In CRITERION II described above, as ϕ , p and n should approach their saturation limits 1 v p and 1 v n , respectively. Here, this would be changed to approach 1 v s instead of 1 v p and 1 v n respectively, although approaching 1 v p and 1 v n sounds more physically correct. This paradoxical conclusion is from entropy rigorously derived by the traditional mean-field lattice gas model based on combinatorics requiring identical occupation sites. Can this be fixed to resume the limit approach to 1 v p and 1 v n and still holding the ground of combinatorics at the same time? An attempt at this is discussed in the next section.

6. Entropy Fixing for Electrolytes under Extreme Concentration Conditions

Here we hope to construct a steric PB model with entropy able to be derived by the mean-field lattice gas model, and at the same time showing physically correct saturation limits for ions as ϕ . The mean-field lattice gas model is fixed here such that each identical cubic site is allowed to be occupied by more than one solute particle of the same species as illustrated in Figure 1b. Although this kind of distribution is no more a most probable distribution as stated earlier, it allows more efficient packing when space is extremely limited and size among species varies largely. Again, we consider the entropy for an aqueous electrolyte system:
T S = k B T   log W ,
Let N ˜ j ,   j = 1 , , K + 1 , be the particle number of species j and N j ,   j = 1 , , K + 1 , the number of identical sites occupied by j-species particles with N ˜ j N j . This means that in an extremely concentrated situation an identical site can be occupied by more than one particle of the same species. If an identical cubic site, on average, can allow r j j-species particles to occupy it, we can then relate N ˜ j and N j by N ˜ j = N j r j , or equivalently v s = r j v j . Again, v j is the particle volume of species j. v s is the volume of an identical cubic site with v s = m a x 1 j K + 1 v j . The entropy based on the most probable distribution of all ‘grouped’ species particles over a total of N = j = 1 K + 1 N j available identical sites in a system is:
W = j = 1 K + 1 W j = N ! N 1 ! ( N N 1 ) ! ( N N 1 ) ! N 2 ! ( N N 1 N 2 ) ! ( N j = 1 K 1 N j ) ! N K ! ( N j = 1 K N j ) ! ( N j = 1 K N j ) ! N K + 1 ! = N ! ( j = 1 K N j ! ) N K + 1 ! .
Note that, after all the ions (in group) are distributed, there are N j = 1 K N j = N K + 1 sites that will be filled by solvent molecules in group, so the entropy becomes:
T S = k B T   log N ! ( j = 1 K N j ! ) N K + 1 ! .
Using the Stirling formula log M ! M log M M with M 1 , we can rewrite the entropy as:
T S = k B T   [ N log N N j = 1 K N j log N j + j = 1 K N j ( N j = 1 K N j ) log ( N j = 1 K N j ) + ( N j = 1 K N j ) ] = k B T [ N log N N j = 1 K N j j = 1 K N j log N j N j = 1 K N j ] ,
Using the following relations:
V = N v s ,   or   N V = 1 v s ,
N j V = N ˜ j r j   V = c j r j ,
N j N = N j v s N v s = N j v s V = c j r j r j v j = c j v j ,
where c j is the concentration of species j; V is the volume of the system.
The entropy per unit volume can be expressed as:
T S V = k B T [ 1 v s log 1 1 j = 1 K c j v j j = 1 K c j r j log c j v j 1 j = 1 K c j v j ] .
Compared with (64), specific ion sizes can now appear explicitly in the entropy formula (76). For a binary electrolyte:
T S = k B T [ 1 v s log 1 1 p v p n v n p r p log p v p 1 p v p n v n n r n log n v n 1 p v p n v n ] d V ,
or:
T S = k B T [ p r p log ( p v p ) + n r n log ( n v n ) + 1 v s ( 1 p v p n v n ) log ( 1 p v p n v n ) ] d V ,
which, without loss of generality, can be augmented as:
T S = k B T [ p r p log ( p v p ) p r p + n r n log ( n v n ) n r n + 1 v s ( 1 p v p n v n ) log ( 1 p v p n v n ) 1 v s ( 1 p v p n v n ) ] d V .  
by:
μ p = δ f δ p = z p e ϕ + W s o l , p + k B T r p   log p v p 1 p v p n v n ,
μ n = δ f δ n = z n e ϕ + W s o l , n + k B T r n   log n v n 1 p v p n v n ,
Again, at equilibrium, the chemical potential is uniform everywhere and therefore the local chemical potential must be equal to its bulk value, which is usually known:
μ p = μ p , b ,     μ n = μ n , b .
Usually bulk solutions are dilute and chemical potentials under that condition can be formulated following (80) and (81) with r p = r n = 1 :
μ p , b = k B T   log p b v p 1 p b v p n b v n + W s o l , p , b ,
μ n , b = k B T   log n b v n 1 p b v p n b v n + W s o l , n , b .
By denoting γ p = e μ p , b z p e ϕ W s o l , p k B T / r p   ,   γ n = e μ n , b z n e ϕ W s o l , n k B T / r n ,
(82) forms
p v p 1 p v p n v n = γ p ,     n v n 1 p v p n v n = γ n ,
and can solve for p and n:
p v p = γ p 1 + γ p + γ n ,   n v n = γ n 1 + γ p + γ n ,
or:
p v p = ( p b v p 1 p b v p n b v n ) r p e β p E p 1 + ( p b v p 1 p b v p n b v n ) r p e β p E p + ( n b v n 1 p b v p n b v n ) r n e β n E n
n v n = ( n b v n 1 p b v p n b v n ) r n e β n E n 1 + ( p b v p 1 p b v p n b v n ) r p e β p E p + ( n b v n 1 p b v p n b v n ) r n e β n E n ,
where E p = z p e ϕ + Δ W s o l , p , E n = z n e ϕ + Δ W s o l , n , β p = ( k B T / r p ) 1 , β n = ( k B T / r n ) 1 . Let us denote (87), (88) as the specific-ion-size Bikerman model 3 (SISBM3) for convenience.
Again, we need to check this new model with criteria I and II. For CRITERION II, we can easily deduce from (87) ϕ ,   n 0 , p v p 1 , p 1 v p (saturation). Similarly, from (88), ϕ , p 0 , n v n 1 , n 1 v n (saturation). There is no constraint like p , n 1 v s as ϕ any more, and entropy here can be derived by mean-field lattice gas model.
For CRITERION I, p and n will not approach a Boltzmann distribution p b e β E p and n b e β E n at diluteness unless r p = r n = 1 . This violation of the Boltzmann distribution at the dilution limit is because we allow multiple ions of the same species to occupy an identical cubic site. This failure and a possible cure will be discussed in next section.

7. Discussion and Conclusions

If we wish to obtain a model for electrolytes such that: (1) it can approach a Boltzmann distribution at diluteness; (2) it can reach the saturation limit as the reciprocal of specific ion size under extreme electrostatic conditions; (3) its entropy can be derived by a mean-field lattice gas model. The only options here is SISBM3 with r p = r n = 1 , since SISBM2 satisfies (1) and (2) but not (3). How can we justify r p = r n = 1 for SISBM3 here? Interpreting all ion sizes as being about the same is certainly not acceptable. Remember SISBM3 is designed for extremely high ion concentrations motivated by the more efficient packing shown in Figure 1b. Actually for situations that would give rise to extremely high ion concentrations and make the steric effect not negligible, such as the Stern layer in the electric double layer (EDL) of a charged wall (discussed next) and the selectivity filter of a K channel [20], there would be ‘locally’ one species only, which is the counter-ion of the local electrostatic environment, since co-ions (and even water) would be totally expelled. Taking the K channel selectivity filter as an example, its extreme narrowness and the strong negative oxygen charges inside it would definitely justify only one species being inside the selectivity filter, which is definitely potassium. This implies r p to be 1 locally in the filter, and we can justify r n = 1 inside the filter as well since anions would be extremely dilute there due to strong electrostatic repulsion. For the rest of the K channel where ions are at most moderately concentrated, the steric effect is much less significant, and basically the original Bikerman model would be appropriate for it. Since SISBM3 with r p = r n = 1 would be a very good approximation of the original Bikerman model under mild ion concentrations, we can use SISBM3 with r p = r n = 1 globally for the whole K channel then. Notice, under r p = r n = 1 , SISBM3 is actually same as SISBM1, but with a rigorous derivation now. This model has been useful and proven to fit the experimental data quite well [16,17,18]. Although here we just discussed steric-effect modifications for PB, modifications for PNP can be likewise derived.
Here we compare SISBM1 and PB by computing ion distributions in a 1D charged wall problem. Many researchers have used this physical model to investigate the surface differential capacitance of electrodes adjacent to electrolyte solutions [12,13,14,15]. Here (14) for PB and (44) for SISBM1 were used to calculate the ion distributions of a binary KCl electrolyte solution without considering the solvation energy and permanent charges. The associated boundary conditions are ϕ ( 0 ) = V w a l l , and ϕ ( ) = 0 . The bulk concentration of KCl as x is set to be c 0 = 100 mM, and dielectric constant is set to 80 for the whole domain ( 0 , ) . The Debye length, featuring the order of thickness of EDL, is λ D = ε 0 ε r k B T c 0 e 2 = 13.78 . The simulation result is shown in Figure 2 with Figure 2a being the distributions of [ K + ] for SISBM1 and PB under V w a l l = 0.1 V and 2 V. Figure 2b is the counterpart plot of Figure 2a for [ C l ] . In Figure 2a, the [ K + ] distributions for SISBM1 and PB are very close to each other and almost indistinguishable in the graph at a weak wall voltage V w a l l = 0.1 V. When the wall voltage increases to V w a l l = 2 V, [ K + ] calculated by SISBM1 reaches its saturation limit 1 / v K = 7.90 × 10 4 mM right adjacent to the wall, but [ K + ] unrealistically increases beyond the saturation limit when computed by the PB model. The main effect of SISBM1 is to offer a saturation limit for counter-ions (K here) when electrostatic attraction from electrode is strong enough, while it is very close to the result of PB when the electrostatic attraction is weak. In Figure 2b, [ C l ] distributions calculated by SISBM1 and PB are very close to each other for both strong and weak wall voltages due to the diluteness caused by electrostatic repulsion to the co-ion (Cl here) of the electrode. Note that, corresponding to a saturation layer of [ K + ] adjacent to wall at V w a l l = 2 V (see Figure 2a), [ C l ] almost vanishes at that layer as well (see Figure 2b). This implies a total exclusion of Cl over there due to the saturation of K, and justifies the locally one-species argument above. If we use the original Bikerman model (41), in which ion sizes are universal, a similar saturating phenomenon for counter-ion concentration can still be obtained. However, specific ion sizes are particularly desired when electrolyte solutions are ternary, like a mixture of KCl and NaCl solutions since K and Na have different sizes, which would saturate at different limiting concentrations. These would be otherwise indistinguishable if using the original Bikerman model.
Above we assume the rest of space after the occupation of ions is exclusively occupied by solvent particles such as water. [9,10,11,12] have suggested that the rest of space should be occupied by solvent or void, so the K + 1 species in (58) and (70) should be interpreted as solvent or void. This may make more sense. Taking the selectivity filter of a K channel as an example, more and more evidences have shown the selectivity filter of a K channel is exclusively occupied by potassium and voids, and water is not allowed there due to the strong solvation energy barrier. [9,10,11,12] even explicitly separate water and voids as two species in their modeling. However, that means the species transport equation (Nernst-Planck equation) of water needs to be modeled explicitly when constructing a PNP type model. This water equation is generally hard to model due to its physical complexity, especially for its electrostatic behavior since water is a dipole.

Funding

This research was funded by Ministry of Science and Technology of Taiwan under Grant Nos. MOST 108-2115-M-035 -002 -MY2 and MOST 108-2218-E-035-005.

Acknowledgments

The author would like to thank Robert Eisenberg and Jinn-Liang Liu for inspiring suggestions and discussions at Hsinchu Taiwan.

Conflicts of Interest

The author declares no conflict of interest.

References

  1. Bazant, M.Z.; Storey, B.D.; Kornyshev, A. Double layer in ionic liquids: Overscreening versus crowding. Phys. Rev. Lett. 2012, 109, 149903. [Google Scholar] [CrossRef]
  2. Rosenfeld, Y. Free-energy model for the inhomogeneous hard-sphere fluid mixture and density-functional theory of freezing. Phys. Rev. Lett. 1989, 63, 980–983. [Google Scholar] [CrossRef]
  3. Rosenfeld, Y. Free energy model for the inhomogeneous fluid mixtures: Yukawa-charged hard spheres, general interactions, and plasmas. J. Chem. Phys. 1993, 98, 8126–8148. [Google Scholar] [CrossRef]
  4. Eisenberg, B.; Hyon, Y.; Liu, C.J. Energy variational analysis of ions in water and channels: Field theory for primitive models of complex ionic fluids. J. Chem. Phys. 2010, 133, 104104. [Google Scholar] [CrossRef] [Green Version]
  5. Ji, S.; Liu, W. Poisson-Nernst-Planck systems for ion flow with Density Functional Theory for hard-sphere potential: I-V relations and critical potentials. Part I: Analysis. J. Dyn. Diff. Equat. 2012, 24, 955–983. [Google Scholar] [CrossRef]
  6. Lin, G.; Liu, W.; Yi, Y.; Zhang, M. Poisson-Nernst-Planck systems for ion flow with a local hard-sphere potential for ion size effects. SIAM J. Appl. Dyn. Syst. 2013, 12, 1613–1648. [Google Scholar] [CrossRef]
  7. Liu, W.; Tu, X.; Zhang, M. Poisson–Nernst–Planck systems for ion flow with density functional theory for hard-sphere potential: I–V relations and critical potentials. Part II: Numerics. J. Dyn. Diff. Equat. 2012, 24, 985–1004. [Google Scholar] [CrossRef]
  8. Qiao, Y.; Lu, B.; Chen, M. A local approximation of fundamental measure theory incorporated into three dimensional Poisson-Nernst-Planck equations to account for hard sphere repulsion among ions. J. Stat. Phys. 2016, 163, 156–174. [Google Scholar] [CrossRef] [Green Version]
  9. Horng, T.-L.; Lin, T.-C.; Liu, C.; Eisenberg, B. PNP equations with steric effects: A model of ion flow through channels. J. Phys. Chem. B 2012, 116, 11422–11441. [Google Scholar] [CrossRef] [PubMed]
  10. Bikerman, J.J. Structure and capacity of electrical double layer. Philos. Mag. 1942, 33, 384–397. [Google Scholar] [CrossRef]
  11. Borukhov, I.; Andelman, D.; Orland, H. Steric effects in electrolytes: A modified Poisson-Boltzmann equation. Phys. Rev. Lett. 1997, 79, 435–438. [Google Scholar] [CrossRef] [Green Version]
  12. Gavish, N.; Promislow, K. Systematic interpretation of differential capacitance data. Phys. Rev. E 2015, 92, 012321. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  13. López-García, J.J.; Horno, J.; Grosse, C. Differential capacitance of the diffuse double layer at electrode- electrolyte interfaces considering ions as dielectric spheres: Part I. Binary electrolyte solutions. J. Colloid Interface Sci. 2017, 496, 531–539. [Google Scholar] [CrossRef] [PubMed]
  14. Nakayama, Y.; Andelman, D. Differential capacitance of the electric double layer: The interplay between ion finite size and dielectric decrement. J. Chem. Phys. 2015, 142, 044706. [Google Scholar] [CrossRef] [Green Version]
  15. Zhao, H. Influence of nonelectrostatic ion-ion interactions on double-layer capacitance. Phys. Rev. E 2012, 86, 051502. [Google Scholar] [CrossRef]
  16. Liu, J.-L.; Eisenberg, B. Poisson-Nernst-Planck-Fermi theory for modeling biological ion channels. J. Chem. Phys. 2014, 141, 22D532. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  17. Liu, J.-L.; Eisenberg, B. Poisson-Fermi model of single ion activities in aqueous solutions. Chem. Phys. Lett. 2015, 637, 1–6. [Google Scholar] [CrossRef]
  18. Liu, J.-L.; Eisenberg, B. Poisson-Fermi modeling of ion activities in aqueous single and mixed electrolyte solutions at variable temperature. J. Chem. Phys. 2018, 148, 054501. [Google Scholar] [CrossRef]
  19. Liu, J.-L.; Xie, D.; Eisenberg, B. Poisson-Fermi Formulation of Nonlocal Electrostatics in Electrolyte Solutions. Comput. Math. Biophys. 2017, 5, 116–124. [Google Scholar] [CrossRef] [Green Version]
  20. Song, Z.; Cao, X.; Horng, T.-L.; Huang, H. Selectivity of the KcsA potassium channel: Analysis and computation. Phys. Rev. E 2019, 100, 022406. [Google Scholar] [CrossRef] [PubMed] [Green Version]
Figure 1. (a) Moderately concentrate situation with each solute/solvent particle only occupy one identical site. (b) Extremely concentrate situation with each identical site can be allowed to be occupied by multiple solute/solvent particles of the same species in order to increase packing efficiency in space.
Figure 1. (a) Moderately concentrate situation with each solute/solvent particle only occupy one identical site. (b) Extremely concentrate situation with each identical site can be allowed to be occupied by multiple solute/solvent particles of the same species in order to increase packing efficiency in space.
Entropy 22 00632 g001
Figure 2. (a) K distributions in charged wall problem computed by SISBM1 and PB models under different wall voltages. (b) Cl distributions with conditions same as (a).
Figure 2. (a) K distributions in charged wall problem computed by SISBM1 and PB models under different wall voltages. (b) Cl distributions with conditions same as (a).
Entropy 22 00632 g002

Share and Cite

MDPI and ACS Style

Horng, T.-L. Review and Modification of Entropy Modeling for Steric Effects in the Poisson-Boltzmann Equation. Entropy 2020, 22, 632. https://doi.org/10.3390/e22060632

AMA Style

Horng T-L. Review and Modification of Entropy Modeling for Steric Effects in the Poisson-Boltzmann Equation. Entropy. 2020; 22(6):632. https://doi.org/10.3390/e22060632

Chicago/Turabian Style

Horng, Tzyy-Leng. 2020. "Review and Modification of Entropy Modeling for Steric Effects in the Poisson-Boltzmann Equation" Entropy 22, no. 6: 632. https://doi.org/10.3390/e22060632

APA Style

Horng, T. -L. (2020). Review and Modification of Entropy Modeling for Steric Effects in the Poisson-Boltzmann Equation. Entropy, 22(6), 632. https://doi.org/10.3390/e22060632

Note that from the first issue of 2016, this journal uses article numbers instead of page numbers. See further details here.

Article Metrics

Back to TopTop