# Poisson–Boltzmann equation

The Poisson–Boltzmann equation describes the distribution of the electric potential in solution in the direction normal to a charged surface. This distribution is important to determine how the electrostatic interactions will affect the molecules in solution.
It is expressed as a differential equation of the electric potential $\psi$, which depends on the solvent permitivity $\varepsilon$, the solution temperature $T$, and the mean concentration of each ion species $c_{i}^{0}$:
    $\nabla ^{2}\psi =-{\frac {1}{\varepsilon }}\sum _{i}c_{i}^{0}q_{i}\exp \left({\frac {-q_{i}\psi (x,y,z)}{k_{B}T}}\right)$
The Poisson–Boltzmann equation is derived via mean-field assumptions. From the Poisson–Boltzmann equation many other equations have been derived with a number of different assumptions.

## Origins
### Background and derivation

The Poisson–Boltzmann equation describes a model proposed independently by Louis Georges Gouy and David Leonard Chapman in 1910 and 1913, respectively. In the Gouy-Chapman model, a charged solid comes into contact with an ionic solution, creating a layer of surface charges and counter-ions or double layer. Due to thermal motion of ions, the layer of counter-ions is a diffuse layer and is more extended than a single molecular layer, as previously proposed by Hermann Helmholtz in the Helmholtz model. The Stern Layer model goes a step further and takes into account the finite ion size.
| Theory | Important characteristics | Assumptions |
| Helmholtz | Surface charge neutralized by a molecular layer of counter-ions; surface charge potential linearly dissipated from surface to counter-ions to satisfy charge | Thermal motion, ion diffusion, adsorption onto the surface, solvent/surface interactions considered negligible |
| Gouy-Chapman | Thermal motion of ions accounted for; ions behave as point charges | Finite ion size ignored; uniformly-charged surface; non-Coulombic interactions ignored |
| Stern | Finite ion size and hydration sphere considered; some ions are specifically adsorbed by the surface in the plane, known as the Stern layer | Stern layer is thin compared to particle size; fluid velocity = 0 in Stern layer |
The Gouy–Chapman model explains the capacitance-like qualities of the electric double layer. A simple planar case with a negatively charged surface can be seen in the figure below. As expected, the concentration of counter-ions is higher near the surface than in the bulk solution.
The Poisson–Boltzmann equation describes the electrochemical potential of ions in the diffuse layer. The three-dimensional potential distribution can be described by the Poisson equation $\nabla ^{2}\psi ={\frac {\partial ^{2}\psi }{\partial x^{2}}}+{\frac {\partial ^{2}\psi }{\partial y^{2}}}+{\frac {\partial ^{2}\psi }{\partial z^{2}}}=-{\frac {\rho _{e}}{\varepsilon }},$ where
* $\rho _{e}$ is the local electric charge density in C/m$^3$,
* $\varepsilon$ is the permittivity of the solvent,
* ψ is the electric potential.
The freedom of movement of ions in solution can be accounted for by Boltzmann statistics. The Boltzmann equation is used to calculate the local ion density such that $c_{i}=c_{i}^{0}\cdot \exp \left({\frac {-W_{i}}{k_{\mathrm {B} }T}}\right),$ where
* $c_{i}^{0}$ is the ion concentration at the bulk,
* $W_{i}$ is the work required to move an ion closer to the surface from an infinitely far distance,
* $k_{\mathrm {B} }$ is the Boltzmann constant,
* $T$ is the temperature in kelvins.
The equation for local ion density can be substituted into the Poisson equation under the assumptions that the work being done is only electric work, and that the concentration of salt is much higher than the concentration of ions. The electric work to bring an ion of charge $q_{i}$ to a surface with potential ψ can be represented by $W_{i}=q_{i}\psi$. These work equations can be substituted into the Boltzmann equation, producing an expression for the concentration of each ion species $c_{i}=c_{i}^{0}\exp \left(-{\frac {q_{i}\psi (x,y,z)}{k_{B}T}}\right)$.
Substituting this Boltzmann relation into the local electric charge density expression, the following expression can be obtained $\rho _{e}=\sum _{i}c_{i}q_{i}=\sum _{i}c_{i}^{0}q_{i}\exp \left({\frac {-q_{i}\psi (x,y,z)}{k_{B}T}}\right).$
Finally the charge density can be substituted into the Poisson equation to produce the Poisson–Boltzmann equation:
$\nabla ^{2}\psi =-{\frac {1}{\varepsilon }}\sum _{i}c_{i}^{0}q_{i}\exp \left({\frac {-q_{i}\psi (x,y,z)}{k_{B}T}}\right)$
When distance is measured as multiples of Bjerrum length $l_{b}$ and potential is measured in multiples of $k_{B}T/e$ then the equation can be rearranged to dimensionless form $\nabla ^{2}\psi =2c_{0}(l_{b})^{3}\sinh \left(\psi \right).$

### Related theories

The Poisson–Boltzmann equation can take many forms throughout various scientific fields. In biophysics and certain surface chemistry applications, it is known simply as the Poisson–Boltzmann equation. It is also known in electrochemistry as Gouy-Chapman theory; in solution chemistry as Debye–Huckel theory; in colloid chemistry as Derjaguin–Landau–Verwey–Overbeek (DLVO) theory. Only minor modifications are necessary to apply the Poisson–Boltzmann equation to various interfacial models, making it a highly useful tool in determining electrostatic potential at surfaces.

## Solving analytically

Because the Poisson–Boltzmann equation is a partial differential of the second order, it is commonly solved numerically; however, with certain geometries, it can be solved analytically.

### Geometries

The geometry that most easily facilitates this is a planar surface. In the case of an infinitely extended planar surface, there are two dimensions in which the potential cannot change because of symmetry. Assuming these dimensions are the y and z dimensions, only the x dimension is left. Below is the Poisson–Boltzmann equation solved analytically in terms of a second order derivative with respect to x.
${\frac {d^{2}\psi }{dx^{2}}}={\frac {c_{0}e}{\varepsilon }}\cdot \left[e^{\frac {e\psi (x)}{k_{\mathrm {B} }T}}-e^{\frac {-e\psi (x)}{k_{\mathrm {B} }T}}\right]$
Analytical solutions have also been found for axial and spherical cases in a particular study. The equation is in the form of a logarithm of a power series and it is as follows: ${\frac {d^{2}\psi }{dr^{2}}}+{\frac {L}{r}}{\frac {d\psi }{dr}}=e^{\psi }-\delta e^{-\psi }$
It uses a dimensionless potential $\psi ={\frac {e\Phi }{kT}}$ and the lengths are measured in units of the Debye electron radius in the region of zero potential $R_{eD}={\sqrt {\frac {kT}{4\pi e^{2}n_{e0}}}}$ (where $n_{e0}$ denotes the number density of negative ions in the zero potential region). For the spherical case, L=2, the axial case, L=1, and the planar case, L=0.

### Low-potential vs high-potential cases

When using the Poisson–Boltzmann equation, it is important to determine if the specific case is low or high potential. The high-potential case becomes more complex so if applicable, use the low-potential equation. In the low-potential condition, the linearized version of the Poisson–Boltzmann equation (shown below) is valid, and it is commonly used as it is more simple and spans a wide variety of cases. $\psi =\psi _{0}e^{-\mathrm {K} x}$

#### Low-potential case conditions

Strictly, low potential means that $e\left\vert \psi \right\vert \ll k_{\mathrm {B} }T$; however, the results that the equations yields are valid for a wider range of potentials, from 50–80mV. Nevertheless, at room temperature, $\psi \leq \mathrm {25mV}$ and that is generally the standard. Some boundary conditions that apply in low potential cases are that: at the surface, the potential must be equal to the surface potential and at large distances from the surface the potential approaches a zero value. This distance decay length is yielded by the Debye length $\lambda _{D}$ equation.
$\mathrm {K} ={\sqrt {\frac {2c_{0}e^{2}}{\varepsilon k_{\mathrm {B} }T}}}$
$\lambda _{D}=\mathrm {K} ^{-1}$
As salt concentration increases, the Debye length decreases due to the ions in solution screening the surface charge. A special instance of this equation is for the case of $25^{\circ }C$ water with a monovalent salt. The Debye length equation is then:
$\lambda _{D}={\frac {\mathrm {0.304nm} }{\sqrt {c_{0}}}}$
where $c_{0}$ is the salt concentration in mol/L. These equations all require 1:1 salt concentration cases, but if ions that have higher valence are present, the following case is used.
$\mathrm {K} ={\sqrt {{\frac {e^{2}}{\varepsilon k_{\mathrm {B} }T}}\sum c_{i}{Z_{i}}^{2}}}$

## High-potential case

The high-potential case is referred to as the “full one-dimensional case”. In order to obtain the equation, the general solution to the Poisson–Boltzmann equation is used and the case of low potentials is dropped. The equation is solved with a dimensionless parameter $y\equiv {\frac {e\psi }{k_{B}T}}$, which is not to be confused with the spatial coordinate symbol, y. Employing several trigonometric identities and the boundary conditions that at large distances from the surface, the dimensionless potential and its derivative are zero, the high potential equation is revealed.
$e^{-\mathrm {K} x}={\frac {(e^{y/2}-1)(e^{y_{0}/2}+1)}{(e^{y/2}+1)(e^{y_{0}/2}-1)}}$
This equation solved for $e^{y/2}$ is shown below.
$e^{y/2}={\frac {e^{y_{0}/2}+1+(e^{y_{0}/2}-1)\cdot e^{-\mathrm {K} x}}{e^{y_{0}/2}+1-(e^{y_{0}/2}-1)\cdot e^{-\mathrm {K} x}}}$
In order to obtain a more useful equation that facilitates graphing high potential distributions, take the natural logarithm of both sides and solve for the dimensionless potential, y.
$y=2\ln {\frac {e^{y_{0}/2}+1+(e^{y_{0}/2}-1)\cdot e^{-\mathrm {K} x}}{e^{y_{0}/2}+1-(e^{y_{0}/2}-1)\cdot e^{-\mathrm {K} x}}}$
Knowing that $y\equiv {\frac {e\psi }{k_{B}T}}$, substitute this for y in the previous equation and solve for $\psi$. The following equation is rendered.
$\psi ={\frac {2k_{B}T}{e}}\cdot \ln {\frac {e^{y_{0}/2}+1+(e^{y_{0}/2}-1)\cdot e^{-\mathrm {K} x}}{e^{y_{0}/2}+1-(e^{y_{0}/2}-1)\cdot e^{-\mathrm {K} x}}}$
$y_{0}={\frac {e\psi _{0}}{k_{B}T}}$

### Conditions

In low potential cases, the high potential equation may be used and will still yield accurate results. As the potential rises, the low potential, linear case overestimates the potential as a function of distance from the surface. This overestimation is visible at distances less than half the Debye length, where the decay is steeper than exponential decay. The following figure employs the linearized equation and the high potential graphing equation derived above. It is a potential-versus-distance graph for varying surface potentials of 50, 100, 150, and 200 mV. The equations employed in this figure assume an 80mM NaCl solution.

## General applications

The Poisson–Boltzmann equation can be applied in a variety of fields mainly as a modeling tool to make approximations for applications such as charged biomolecular interactions, dynamics of electrons in semiconductors or plasma, etc. Most applications of this equation are used as models to gain further insight on electrostatics.

### Physiological applications

The Poisson–Boltzmann equation can be applied to biomolecular systems. One example is the binding of electrolytes to biomolecules in a solution. This process is dependent upon the electrostatic field generated by the molecule, the electrostatic potential on the surface of the molecule, as well as the electrostatic free energy.
The linearized Poisson–Boltzmann equation can be used to calculate the electrostatic potential and free energy of highly charged molecules such as tRNA in an ionic solution with different number of bound ions at varying physiological ionic strengths. It is shown that electrostatic potential depends on the charge of the molecule, while the electrostatic free energy takes into account the net charge of the system.
Another example of utilizing the Poisson–Boltzmann equation is the determination of an electric potential profile at points perpendicular to the phospholipid bilayer of an erythrocyte. This takes into account both the glycocalyx and spectrin layers of the erythrocyte membrane. This information is useful for many reasons including the study of the mechanical stability of the erythrocyte membrane.

#### Electrostatic free energy

The Poisson–Boltzmann equation can also be used to calculate the electrostatic free energy for hypothetically charging a sphere using the following charging integral: $\Delta G^{\text{el}}=\int ^{\tau }qU(\tau ')\,d\tau '$ where $\tau q$ is the final charge on the sphere
The electrostatic free energy can also be expressed by taking the process of the charging system. The following expression utilizes chemical potential of solute molecules and implements the Poisson-Boltzmann Equation with the Euler-Lagrange functional: $\Delta G^{\text{el}}=\int _{V}\left(kT\sum _{i}c_{i}^{\infty }\left[1-\exp \left({\frac {-z_{i}qU}{kT}}\right)\right]+p^{f}U-{\frac {-\varepsilon ({\boldsymbol {\nabla }}U)^{2}}{8\pi }}\right)dV$
Note that the free energy is independent of the charging pathway [5c].
The above expression can be rewritten into separate free energy terms based on different contributions to the total free energy $\Delta G^{\text{el}}=\Delta G^{\text{ef}}+\Delta G^{\text{em}}+\Delta G^{\text{mob}}+\Delta G^{\text{solv}}$ where
* Electrostatic fixed charges = $\Delta G^{\text{ef}}=\int _{V}{\frac {p^{f}U}{2}}dV$
* Electrostatic mobile charges = $\Delta G^{\text{em}}=\int _{V}{\frac {\sum _{i}c_{i}z_{i}qU}{2}}dV$
* Entropic free energy of mixing of mobile species = $\Delta G^{\text{mob}}=kT\int _{V}c_{i}\ln {\frac {c_{i}}{c_{i}^{\infty }}}dV$
* Entropic free energy of mixing of solvent = $\Delta G^{\text{solv}}=kT\int _{V}\sum _{i}c_{i}^{\infty }\left[1-\exp \left({\frac {-z_{i}qU}{kT}}\right)\right]dV$
Finally, by combining the last three term the following equation representing the outer space contribution to the free energy density integral $\Delta G^{\text{out}}=\Delta G^{\text{em}}+\Delta G^{\text{mob}}+\Delta G^{\text{solv}}$
These equations can act as simple geometry models for biological systems such as proteins, nucleic acids, and membranes. This involves the equations being solved with simple boundary conditions such as constant surface potential. These approximations are useful in fields such as colloid chemistry.

### Materials science

An analytical solution to the Poisson–Boltzmann equation can be used to describe an electron-electron interaction in a metal-insulator semiconductor (MIS). This can be used to describe both time and position dependence of dissipative systems such as a mesoscopic system. This is done by solving the Poisson–Boltzmann equation analytically in the three-dimensional case. Solving this results in expressions of the distribution function for the Boltzmann equation and self-consistent average potential for the Poisson equation. These expressions are useful for analyzing quantum transport in a mesoscopic system. In metal-insulator semiconductor tunneling junctions, the electrons can build up close to the interface between layers and as a result the quantum transport of the system will be affected by the electron-electron interactions. Certain transport properties such as electric current and electronic density can be known by solving for self-consistent Coulombic average potential from the electron-electron interactions, which is related to electronic distribution. Therefore, it is essential to analytically solve the Poisson–Boltzmann equation in order to obtain the analytical quantities in the MIS tunneling junctions. Applying the following analytical solution of the Poisson–Boltzmann equation (see section 2) to MIS tunneling junctions, the following expression can be formed to express electronic transport quantities such as electronic density and electric current $f_{1}f^{0}-f_{0}+{\frac {eE_{z}\tau _{0}}{m}}{\frac {\partial f_{0}}{\partial v_{z}}}\left(1-e^{\frac {-\tau }{\tau _{0}}}\right)-\int _{0}^{t}{\frac {e}{m}}e{^{\frac {t-\tau '}{\tau _{0}}}}\nabla \rho [r-v(t-t')]\times {\frac {\partial f_{0}}{\partial v}}dt'$
Applying the equation above to the MIS tunneling junction, electronic transport can be analyzed along the z-axis, which is referenced perpendicular to the plane of the layers. An n-type junction is chosen in this case with a bias V applied along the z-axis. The self-consistent average potential of the system can be found using $\rho \rho _{1}+\rho _{2}$ where
* $\rho _{1}\approx {\frac {aE_{z}}{2\lambda _{D1}}}e^{-\lambda _{D1}z}$ and
* $\rho _{2}\approx {\frac {ne{\sqrt {\pi }}G(i\lambda _{D1})e^{{\frac {-t}{\tau _{0}}}-\lambda _{D1}z}}{3{\sqrt {3}}\varepsilon \lambda _{D1}}}\left(1-e^{1-{\sqrt {\frac {2ne^{2}t^{2}}{m\varepsilon }}}}\right)$
λ is called the Debye length.
The electronic density and electric current can be found by manipulation to equation 16 above as functions of position z. These electronic transport quantities can be used to help understand various transport properties in the system.

## Limitations

As with any approximate model, the Poisson–Boltzmann equation is an approximation rather than an exact representation. Several assumptions were made to approximate the potential of the diffuse layer. The finite size of the ions was considered negligible and ions were treated as individual point charges, where ions were assumed to interact with the average electrostatic field of all their neighbors rather than each neighbor individually. In addition, non-Coulombic interactions were not considered and certain interactions were unaccounted for, such as the overlap of ion hydration spheres in an aqueous system. The permittivity of the solvent was assumed to be constant, resulting in a rough approximation as polar molecules are prevented from freely moving when they encounter the strong electric field at the solid surface.
Though the model faces certain limitations, it describes electric double layers very well. The errors resulting from the previously mentioned assumptions cancel each other for the most part. Accounting for non-Coulombic interactions increases the ion concentration at the surface and leads to a reduced surface potential. On the other hand, including the finite size of the ions causes the opposite effect. The Poisson–Boltzmann equation is most appropriate for approximating the electrostatic potential at the surface for aqueous solutions of univalent salts at concentrations smaller than 0.2 M and potentials not exceeding 50–80 mV.
In the limit of strong electrostatic interactions, a strong coupling theory is more applicable than the weak coupling assumed in deriving the Poisson-Boltzmann theory.