AugmentedPBSystem(grid, data)

Create MPB system with pressure. Ion conserving if data.conserveions==true. In that case, the solution is a sparse matrix.

 unknowns(sys, data)

Initialize and return unknown vector.

VoronoiFVM.unknowns(cell::AbstractAugmentedPBCell)

Initialize and return the unknown vector for a given cell.

Arguments

• cell::AbstractAugmentedPBCell: The cell configuration

Returns

• Initial unknown vector with appropriate structure and values

The unknowns include:

• Ion mole fractions (y_α for α = 1,...,N)
• Solvent mole fraction (y_0)
• Electric potential (φ)
• Pressure (p)
• Electric field strength (E)
• For ion-conserving cells: bulk ion concentrations at the domain center

Initial values are set to reasonable defaults (mole fractions ≈ 0.1 for ions, adjusted for solvent to maintain sum = 1).

 calc_cnum(sol,sys)

Obtain ion number densities from system

 calc_c0num(sol,sys)

Obtain solvent number density from system

 calc_cmol(sol,sys)

Obtain ion molarities (molar densities in mol/L) from system

calc_cmol(sol, cell::AbstractAugmentedPBCell)

Calculate ion concentrations in mol/L from the solution.

Arguments

• sol: Solution vector from the solver
• cell::AbstractAugmentedPBCell: The cell configuration

Returns

• Matrix of ion concentrations in mol/L (one row per species, one column per node)

 calc_c0mol(sol,sys)

Obtain solvent molarity (molar density in mol/L) from system

calc_c0mol(sol, cell::AbstractAugmentedPBCell)

Calculate solvent concentration in mol/L from the solution.

Arguments

• sol: Solution vector from the solver
• cell::AbstractAugmentedPBCell: The cell configuration

Returns

• Vector of solvent concentrations in mol/L (one value per node)

 calc_χ(sol,sys)

calc_χ(sol, cell::AbstractAugmentedPBCell)

Calculate electric susceptibility from the solution.

Arguments

• sol: Solution vector from the solver
• cell::AbstractAugmentedPBCell: The cell configuration

Returns

• Vector of susceptibility values (one value per node)

calc_spacecharge(sys, sol)

 ysum(sys,sol)

qsweep(sys)

Sweep over series of surface charges and calculate resulting potential difference.

struct AugmentedPBData

Data structure containing data for equilibrium calculations. All data including molarity in SI basic units

• z::Vector{Int64}: Ion charge numbers.

• N::Int64: Number of ionic species

• κ::Vector{Float64}: Ion solvation numbers

• molarity::Float64: Bulk molarity transformed to number density

• n_E::Vector{Float64}: Bulk ion number densities

• n_avg::Vector{Float64}: Average ion number densities

• q::Vector{Float64}: Surface charges

• φ::Float64: Applied potential

• n0_ref::Float64: Solvent molarity

• v0::Float64: Solvent molecular volume

• vu::Vector{Float64}: Unsolvated ion molecular volume

• χ0::Float64: Dielectric susceptibility of solvent

• χ::Any: Dielectric susceptibility

• χvar::Bool: Dielectric susceptibility model flag

• i0::Int64: Solvent molar fraction index

• iφ::Int64: Electric potential species index

• ip::Int64: Pressure species index

• iE::Int64: Field strength species index

• coffset::Int64: Offset of n_E in species list

• p_ref::Float64: Reference pressure

• pscale::Float64: Pressure scaling nparameter

• cscale::Float64: Concentration scaling parameter

• qscale::Float64: Charge scaling parameter

• Escale::Float64: Electric field scaling parameter

• E_ref::Float64: Reference voltage

• T::Float64: Temperature

• δ0::Any: Variable susceptibility parameter for solvent

• δ::Any: Variable susceptibility parameters

• kT::Float64: Temperature times Boltzmann constant

• e::Float64: Electron charge

• ε_0::Float64: Vacuum permittivity

• conserveions::Bool: Ion conservation flag

• nv::Vector{Float64}: node volumes

apply_voltage!(data, v)

Set working electrode voltage

apply_charge!(data,q)

Set surface charge in data

set_molarity!(data,M)

Set the molarity of the electrolyte and update depending data

set_molarity!(cell::AbstractAugmentedPBCell, M)

Set the bulk electrolyte molarity.

Arguments

• cell::AbstractAugmentedPBCell: The cell configuration
• M: Molarity in mol/L

Updates the bulk ion concentrations and related parameters in the cell data.

   L_Debye(data)

\[L_{Debye}=\sqrt{ \frac{(1+χ)ε_0k_BT}{e^2n_E}}\]

dlcap0(data)

Double layer capacitance at $φ=0$

\[C_{dl,0}=\sqrt{\frac{2(1+χ) ε_0e^2 n_E}{k_BT}}\]

 capscalc(sys, molarities)

Calculate double layer capacitances using qsweep results.

This provides an "inverse" method to calculate these capacitances. Usually one calculates charges dependent on voltages, here we calculate voltages dependent on charges.

struct DerivedData

Struct holding some derived data

• v::Vector{Float64}: Effective ion volumes

• y_E::Vector: Bulk ion mole fractions

• y0_E::Any: Bulk solvent mole fraction

DerivedData(augmentedpbdata)

Calculate bulk mole fractions from incompressibiltiy

DerivedData(augmentedpbdata, n_E)

Calculate bulk mole fractions from incompressibiltiy:

\[\begin{aligned} \sum\limits_αv_αn_α^E&=1\\ n_0^E&=\frac1{v_0}\left(1-∑\limits_{α>0}v_αn_α^E\right)\\ n^E&=\frac1{v_0}\left(1-∑\limits_{α>0}v_αn_α^E\right)+ ∑\limits_{α>0}n_α^E\\ &=\frac1{v_0}\left(1-∑\limits_{α>0}(v_α-v_0)n_α^E\right)\\ &=\frac1{v_0}\left(1-∑\limits_{α>0}((1+ κ_α)v_0-v_0)n_α^E\right)\\ &=\frac1{v_0}\left(1-∑\limits_{α>0}κ_αv_0n_α^E\right)\\ &=\frac1{v_0}-∑\limits_{α>0}κ_αn_α^E\\ y_α^E&=\frac{n_α^E}{n^E} \end{aligned}\]

makeδ(v, χ, T)

Calculate δ parameter for susceptibility models

W(x)

$W(x)= 3\frac{\mathcal L(x)}{x}= \frac{\coth(x)-\frac1x}{x}$

$\mathcal L(x)$ is called Langevin function.

 Λ(x)

$\Lambda(x)=\ln\left( \frac{\sinh(x)}x \right)$ This is the antiderivative of the Langevin function $\mathcal L(x)$.

y_α(φ,p,E,α,data, ddata)

Ion molar fractions

Equilibrium expression for mole fractions ($α≥0$) (16)

\[y_α(φ,p, E)=y_α^E\exp\left(\frac{-z_αe}{k_BT}(φ- φ^E)-\frac{v_α}{k_BT}(p-p^E) + \Lambda\left(\frac{\delta_\alpha |E|}{k_BT}\right)¸\right)\]

y0(p, E, data, ddata)

Solvent molar fraction

ysum(u, data)

spacecharge(u, data)

Solvated ion volumes:

\[\begin{aligned} q(φ,p)&=e∑\limits_α z_αn_α = ne∑\limits_α z_αy_α\\ &=e\frac{∑\limits_α z_αy_α(\phi,p)}{∑\limits_α v_α y_α(\phi,p)}\\ v_α&=v^u_α+κ_αv_0 \end{aligned}\]

susceptibility(u, data)

Susceptibility model

c_num!(c,φ,p, data)

Calculate number concentrations at discretization node

\[\begin{aligned} n&=\sum_{i=0}^N y_α v_α\\ n_α&=ny_α \end{aligned}\]

c0_num(c,φ,p, data)

Calculate number concentration of solvent at discretization node

reaction!(f,u,node, data)

Callback which runs in every grid point.

• Calculate space charge density and add this to the poisson equation

poisson_and_p_flux!(f, u, edge, data)

Runs on every grid edge. Calculate fluxes for the Poisson and the pressure equations.

bcondition!(y, u, bnode, data)

Boundary condition callback. The Dirichlet condition for the pressure in the mid of the domain ensures uniqueness of the pressure equation.

 ionconservation!(f, u, sys, data)

"Generic callback" which shall ensure the ion conservation constraint. This method runs over the full grid, and its sparsity pattern is automatically detected. It is called if ion conservation is required. In this case, $n^E_lpha$ are additional unknowns scaled by cscale (default: $N_A$) which are attached to the mid of the domain (node i3), and additional equations need to be assembled. These are N-1 ion conservation constraints and the an electroneutrality constraint.