NPP for thin membranes
Nernst-Planck-Poisson Model for the description of electrodiffusion in thin membranes. It allows obtaining the potential and concentration profiles inside the membrane.
The Nernst-Planck equation relates the flux of ions to diffusion and migration processes and the Poisson equation describes the electric field occurring due to the movement of ions in the system. Software uses Chang-Jaffe boundary conditions for concentrations and Dirichlet boundary conditions for potential.
NPP model - basic equations
The NPP model is an initial-boundary problem. Below the solution based on Dirichlet boundary conditions for concentrations and potential is presented.
Mass-balance equation
No reactions are considered, therefore the mass-balance equation takes the form:
$$
\frac{\partial c_i}{\partial t} = - \frac{\partial J_i}{\partial x}
$$
where `c_i` represents the concentration of the `i`-th ion and `J_i` denotes its flux; `x` is space and `t` represents time.
Nernst-Planck equation
The Nernst-Planck equation for flux includes the description of diffusion and migration processes (as two components of the sum):
$$
J_i= - D_i \frac{\partial c_i}{\partial x} - \frac{F}{RT} z_i D_i c_i \frac{\partial \varphi}{\partial x}
$$
where `D_i` represents the diffusion coefficient of the i-th ion and `z_i` denotes its charge; `\varphi` represents the electrical potential, `F` is the Faraday constant, `R` is the gas constant and `T` represents temperature.
Poisson equation
The Poisson equation relates the potential with the occuring deviation from electroneutrality:
$$
\frac{\partial^2 \varphi}{\partial x^2} = - \frac{F}{\varepsilon} \sum_i z_i c_i
$$
where `\varepsilon` is the permittivity of the medium.
Initial and boundary conditions
Constant initial concentrations of all ions are assumed:
`c_i (x,0) = c_{i,ST} = const`
Dirichlet boundary conditions for potential are used, i.e. constant values of potential are assumed at the boundaries:
$$
\left\{ \begin{array}{l}
\varphi(0,t) = \varphi_L = const \\
\varphi(l,t) = \varphi_R = const \\
\end{array} \right.
$$
Dirichlet boundary conditions can be used also for concentrations:
$$
\left\{ \begin{array}{l}
c_i(0,t) = c_{i,L} \cdot k_{in} / k_{out} = const \\
c_i(l,t) = c_{i,R} \cdot k_{in} / k_{out} = const
\end{array} \right.
$$
Alternatively Chang-Jaffe boundary conditions can be applied:
$$
\left\{ \begin{array}{l}
J_i(0,t) = k_{in} \cdot c_{i,L} - k_{out} \cdot c_i(0,t)\\
J_i(l,t) = k_{out} \cdot c_i(l,t) - k_{in} \cdot c_{i,R}
\end{array} \right.
$$
where `k_{i n}` and `k_{out}` are the heterogeneous rate constants of `i`-th ion, entering and leaving the membrane, respectively.
NPP model - numerical solution
Numerical solution of the above equations is based on the Method of Lines (MoL) and on the explicit Euler method.
Method of Lines for concentrations
We divide the considered segment into `n+1` sub-segments, each having the length `h`, to obtain the discrete points `x_0` to `x_{n+1}`. Let's then define `c_{i,k} \equiv c_i(x_k)`, and define the additional points ` x_{k-1/2} ` and ` x_{k+1/2}` in the distance `h/2` from point `x_k`. In these points we will be establishing the fluxes (`J_{i,k-1/2}` and `J_{i,k+1/2}`).
Mass balance equation for the `i`-th component in `k`-th point of space is given as:
$$
\frac{dc_{i,k}}{dt} = - \frac{J_{i,k+1/2}-J_{i,k-1/2}}{h}
$$
Nernst-Planck equation for fluxes `J_{i,k-1/2}` and `J_{i,k+1/2}` takes the form:
$$
\left\{ \begin{array}{l}
J_{i,k-1/2}= - D_i \frac{d c_{i,k-1/2}}{dx} - \frac{F}{RT} z_i D_i c_{i,k-1/2} \frac{d \varphi_{k-1/2}}{dx} \\
J_{i,k+1/2}= - D_i \frac{d c_{i,k+1/2}}{dx} - \frac{F}{RT} z_i D_i c_{i,k+1/2} \frac{d \varphi_{k+1/2}}{dx}
\end{array} \right.
$$
Now the space derivatives, occuring in Nernst-Planck equation, can be substituted with finite differences:
$$
\left\{ \begin{array}{l}
J_{i,k-1/2}= - D_i \frac{c_{i,k}-c_{i,k-1}}{h} - \frac{F}{RT} z_i D_i \frac{c_{i,k}+c_{i,k-1}}{2} \frac{\varphi_k -\varphi_{k-1}}{h} \\
J_{i,k+1/2}= - D_i \frac{c_{i,k+1}-c_{i,k}}{h} - \frac{F}{RT} z_i D_i \frac{c_{i,k+1}+c_{i,k}}{2} \frac{\varphi_{k+1} -\varphi_k}{h}
\end{array} \right.
$$
After substituting `J_{i,k-1/2}` and `J_{i,k+1/2}` into the mass balance equation, one obtains:
$$
\frac{dc_{i,k}}{dt} = \frac{D_i}{h^2} \left( c_{i,k+1}-2c_{i,k}+c_{i,k-1} \right) + \frac{F z_i D_i}{2RT h^2} \left( (c_{i,k+1}+c_{i,k})(\varphi_{k+1} -\varphi_k) - (c_{i,k}+c_{i,k-1})(\varphi_k -\varphi_{k-1}) \right).
$$
Euler Method
If `c_{i,k}` denotes the concentration at time `t_j` and `c\text{_}n e w_{i,k}` denotes the the concentration at time `t_{j+1}` then:
$$
c\text{_}new_{i,k} = c_{i,k} + dt \cdot \frac{dc_{i,k}}{dt}
$$
After inserting the time derivative, one obtains:
$$
c\text{_}new_{i,k} = c_{i,k} + \alpha_i \left( c_{i,k+1}-2c_{i,k}+c_{i,k-1} \right) + \beta_i \bigl( (c_{i,k+1}+c_{i,k})(\varphi_{k+1} -\varphi_k) - (c_{i,k}+c_{i,k-1})(\varphi_k -\varphi_{k-1}) \bigr)
$$
where `\alpha_i = \frac{dt \cdot D_i}{h^2}` and `\beta_i = \frac{dt \cdot F z_i D_i}{2RT h^2}`
Method of Lines for Poisson equation
After substituting the space derivative with the finite difference, one obtains:
$$
\frac{\varphi_{k+1}-2\varphi_k+\varphi_{k-1}}{h^2} = - \frac{F}{\varepsilon} \sum_i z_i c_{i,k}
$$
After then multiplying by `h^2`, the equation of the `i`-th element is obtained:
$$
\varphi_{k-1} - 2\varphi_k + \varphi_{k+1} = - \gamma \sum_i z_i c_{i,k}
$$
where `\gamma = \frac{h^2 F}{\varepsilon}`
Knowing that ` \varphi_0 = \varphi_L = const` and `\varphi_{n+1} =\varphi_R =const`, we can substitute `k=1` and `k=n`, obtaining:
$$ -2\varphi_1 + \varphi_2 = - \varphi_L - \gamma \sum_i z_i c_{i,1} $$
$$ \varphi_{n-1} - 2\varphi_{n} = - \varphi_R - \gamma \sum_i z_i c_{i,n} $$
The set of these three equations represent the tri-diagonal matrix:
$$
\left[ \begin{array}{ccccccccc|c}
-2 & 1 & 0 & 0 & \cdots & 0 & 0 & 0 & 0 & - \gamma \sum z_i c_{i,1} - \varphi_L \\
1 & -2 & 1 & 0 & \cdots & 0 & 0 & 0 & 0 & - \gamma \sum z_i c_{i,2} \\
0 & 1 & -2 & 1 & \cdots & 0 & 0 & 0 & 0 & - \gamma \sum z_i c_{i,3} \\
0 & 0 & 1 & -2 & \cdots & 0 & 0 & 0 & 0 & - \gamma \sum z_i c_{i,4} \\
\vdots & \vdots & \vdots & \vdots & \ddots & \vdots & \vdots & \vdots & \vdots\\
0 & 0 & 0 & 0 & \cdots & -2 & 1 & 0 & 0 & - \gamma \sum z_i c_{i,n-3} \\
0 & 0 & 0 & 0 & \cdots & 1 & -2 & 1 & 0 & - \gamma \sum z_i c_{i,n-2} \\
0 & 0 & 0 & 0 & \cdots & 0 & 1 & -2 & 1 & - \gamma \sum z_i c_{i,n-1} \\
0 & 0 & 0 & 0 & \cdots & 0 & 0 & 1 & -2 & - \gamma \sum z_i c_{i,n} - \varphi_R \\
\end{array} \right]
$$
