Source Code
Warning
If the Newton solver diverges for any of the solutions, you may try to reduce the relaxation parameter.
- Eq04.create_refined_mesh(refinement_style: str, number_cells: int) Mesh
Creates a one-dimensional mesh with a refined region at the left boundary
- Parameters:
refinement_style (str) – How the mesh should be refined. Options are ‘log’, ‘hard_log’, ‘hard_hard_log’
number_cells (int) – Number of cells in the mesh
- Returns:
One-dimensional mesh, ready for use in FEniCSx
- Return type:
Mesh
- EqN.solve_System_Neq(phi_left: float, phi_right: float, p_right: float, z_alpha: list, y_R: list, K: float | str, Lambda2: float, a2: float, number_cells: int, solvation: float = 0, PoissonBoltzmann: bool = False, relax_param: float = None, x0: float = 0, x1: float = 1, refinement_style: str = 'uniform', return_type: str = 'Vector', rtol: float = 1e-08, max_iter: float = 500)
Solve the dimensionless system of equations presented in: Numerical Treatment of a Thermodynamically Consistent Electrolyte Model, B.Sc. Thesis Habscheid 2024
- System of equations:
λ²Δ φ =−L²n^F
a²∇p=−n^F∇ φ
div(J_α)=0 α∈ {1,…,N−1}
with φ the electric potential, p the pressure, n^F the total free charge density, J_α the diffusion fluxes of species α, λ² a dimensionless parameter, L²=1, a² a dimensionless parameter, N the number of species, and α the species index.
! If the Newton solver diverges, you may try to reduce the relaxation parameter.
- Parameters:
phi_left (float) – Value of φ at the left boundary
phi_right (float) – Value of φ at the right boundary
p_right (float) – Value of p at the right boundary
z_alpha (list) – Charge numbers for species α = 1,…,N-1
y_R (list) – Atomic fractions at right boundary for species α = 1,…,N-1
K (float | str) – Dimensioness bulk modulus of the electrolyte. If ‘incompressible’, the system is solved for an incompressible electrolyte
Lambda2 (float) – Dimensionless parameter
a2 (float) – Dimensionless parameter
number_cells (int) – Number of cells in the mesh
solvation (float, optional) – solvation number, not implemented yet, by default 0
PoissonBoltzmann (bool, optional) – Solve classical Nernst-Planck model with the use of the Poisson-Boltzmann formulation if True, else solve the presented model by Dreyer, Guhlke, Müller, Not implemented yet, by default False
relax_param (float, optional) – Relaxation parameter for the Newton solver xₙ₊₁ = γ xₙ f(xₙ)/f’(xₙ) with γ the relaxation parameter , by default None -> Determined automatically
x0 (float, optional) – Left boundary of the domain, by default 0
x1 (float, optional) – Right boundary of the domain, by default 1
refinement_style (str, optional) – Specify for refinement towards zero Options are ‘uniform’, ‘log’, ‘hard_log’, ‘hard_hard_log’ by default ‘uniform’
return_type (str, optional) – ‘Vector’ or ‘Scalar’ (not implemented yet, should be implemented in a later version), ‘Scalar’ returns dolfinx.fem type and ‘Vector’ numpy arrays of the solution, by default ‘Vector’
rtol (float, optional) – Relative tolerance for Newton solver, by default 1e-8
max_iter (float, optional) – Maximum number of Newton iterations, by default 500
- Returns:
Returns atomic fractions for species A and C, electric potential, pressure, and the mesh If return_type is ‘Vector’, the solution is returned as numpy arrays Only return_type ‘Vector’ is implemented yet
- Return type:
y_A, y_C, phi, p, msh
- Eq04.solve_System_4eq(phi_left: float, phi_right: float, p_right: float, z_A: float, z_C: float, y_A_R: float, y_C_R: float, K: float | str, Lambda2: float, a2: float, number_cells: int, solvation: float = 0, PoissonBoltzmann: bool = False, relax_param: float = None, x0: float = 0, x1: float = 1, refinement_style: str = 'uniform', return_type: str = 'Scalar', rtol: float = 1e-08, max_iter: float = 500)
Solve the dimensionless system of equations presented in: Numerical Treatment of a Thermodynamically Consistent Electrolyte Model, B.Sc. Thesis Habscheid 2024
- System of equations:
λ²Δ φ =−L²n^F
a²∇p=−n^F∇ φ
div(J_A)=0
div(J_C)=0
with φ the electric potential, p the pressure, n^F the total free charge density, J_α the diffusion fluxes of species α, λ² a dimensionless parameter, L²=1, a² a dimensionless parameter, N the number of species, and α the species index.
! If the Newton solver diverges, you may try to reduce the relaxation parameter.
- Parameters:
phi_left (float) – Value of φ at the left boundary
phi_right (float) – Value of φ at the right boundary
p_right (float) – Value of p at the right boundary
z_A (float) – Charge number of species A
z_C (float) – Charge number of species C
y_A_R (float) – Atomic fractions of species A at right boundary
y_C_R (float) – Atomic fractions of species C at right boundary
K (float | str) – Dimensioness bulk modulus of the electrolyte. If ‘incompressible’, the system is solved for an incompressible electrolyte
Lambda2 (float) – Dimensionless parameter
a2 (float) – Dimensionless parameter
number_cells (int) – Number of cells in the mesh
solvation (float, optional) – solvation number, by default 0
PoissonBoltzmann (bool, optional) – Solve classical Nernst-Planck model with the use of the Poisson-Boltzmann formulation if True, else solve the presented model by Dreyer, Guhlke, Müller, by default False
relax_param (float, optional) – Relaxation parameter for the Newton solver xₙ₊₁ = γ xₙ f(xₙ)/f’(xₙ) with γ the relaxation parameter , by default None -> Determined automatically
x0 (float, optional) – Left boundary of the domain, by default 0
x1 (float, optional) – Right boundary of the domain, by default 1
refinement_style (str, optional) – Specify for refinement towards zero Options are ‘uniform’, ‘log’, ‘hard_log’, ‘hard_hard_log’ by default ‘uniform’
return_type (str, optional) – ‘Vector’ or ‘Scalar’, ‘Scalar’ returns dolfinx.fem type and ‘Vector’ numpy arrays of the solution, by default ‘Scalar’
rtol (float, optional) – Relative tolerance for Newton solver, by default 1e-8
max_iter (float, optional) – Maximum number of Newton iterations, by default 500
- Returns:
Returns atomic fractions for species A and C, electric potential, pressure, and the mesh If return_type is ‘Vector’, the solution is returned as numpy arrays
- Return type:
y_A, y_C, phi, p, msh
- Eq02.solve_System_2eq(phi_left: float, phi_right: float, p_right: float, z_A: float, z_C: float, y_A_R: float, y_C_R: float, K: float | str, Lambda2: float, a2: float, number_cells: int, solvation: float = 0, relax_param: float = None, x0: float = 0, x1: float = 1, refinement_style: str = 'uniform', return_type: str = 'Scalar', rtol: float = 1e-08, max_iter: float = 500)
Solve the simplified dimensionless system of equations presented in: Numerical Treatment of a Thermodynamically Consistent Electrolyte Model, B.Sc. Thesis Habscheid 2024
- System of equations:
λ²Δ φ =−L²n^F
a²∇p=−n^F∇ φ
- with the space charge
n^F = z_A y_A(φ, p) + z_C y_C(φ ,p)
- if the mixture is compressible:
y_alpha = C_alpha * (K+p−1)^(−κ+1)a²K exp(−z_α φ)
- if the mixture is incompressible:
y_alpha = D_alpha * exp(−(κ+1)a²p−z_α φ)
with φ the electric potential, p the pressure, n^F the total free charge density, J_α the diffusion fluxes of species α, λ² a dimensionless parameter, L²=1, a² a dimensionless parameter, N the number of species, and α the species index.
! If the Newton solver diverges, you may try to reduce the relaxation parameter.
- Parameters:
phi_left (float) – Value of φ at the left boundary
phi_right (float) – Value of φ at the right boundary
p_right (float) – Value of p at the right boundary
z_A (float) – Charge number of species A
z_C (float) – Charge number of species C
y_A_R (float) – Atomic fractions of species A at right boundary
y_C_R (float) – Atomic fractions of species C at right boundary
K (float | str) – Dimensioness bulk modulus of the electrolyte. If ‘incompressible’, the system is solved for an incompressible electrolyte
Lambda2 (float) – Dimensionless parameter
a2 (float) – Dimensionless parameter
number_cells (int) – Number of cells in the mesh
solvation (float, optional) – solvation number, by default 0
relax_param (float, optional) – Relaxation parameter for the Newton solver xₙ₊₁ = γ xₙ f(xₙ)/f’(xₙ) with γ the relaxation parameter , by default None -> Determined automatically
x0 (float, optional) – Left boundary of the domain, by default 0
x1 (float, optional) – Right boundary of the domain, by default 1
refinement_style (str, optional) – Specify for refinement towards zero Options are ‘uniform’, ‘log’, ‘hard_log’, ‘hard_hard_log’ by default ‘uniform’
return_type (str, optional) – ‘Vector’ or ‘Scalar’, ‘Scalar’ returns dolfinx.fem type and ‘Vector’ numpy arrays of the solution, by default ‘Scalar’
rtol (float, optional) – Relative tolerance for Newton solver, by default 1e-8
max_iter (float, optional) – Maximum number of Newton iterations, by default 500
- Returns:
Returns atomic fractions for species A and C, electric potential, pressure, and the mesh If return_type is ‘Vector’, the solution is returned as numpy arrays
- Return type:
y_A, y_C, phi, p, msh
- ElectrolyticDiode.ElectrolyticDiode(Bias_type: str, phi_bias: float, g_phi: float, z_A: float, z_C: float, y_A_bath: float, y_C_bath: float, K: float | str, Lambda2: float, a2: float, number_cells: list, solvation: float = 0, PoissonBoltzmann: bool = False, relax_param: float = None, Lx: float = 2, Ly: float = 10, rtol: float = 1e-08, max_iter: float = 500, return_type: str = 'Vector')
Solves the system of equations for the example of an electric diode
Solve the dimensionless system of equations presented in: Numerical Treatment of a Thermodynamically Consistent Electrolyte Model, B.Sc. Thesis Habscheid 2024
! If the Newton solver diverges, you may try to reduce the relaxation parameter.
- Parameters:
Bias_type (str) – ForwardBias, NoBias, BackwardBias
phi_bias (float) – Bias in φ
g_phi (float) – Neumann boundary condition for φ
z_A (float) – Charge number of species A
z_C (float) – Charge number of species C
y_A_bath (float) – Bath concentration of anions
y_C_bath (float) – Bath concentration of cations
K (float | str) – Dimensioness bulk modulus of the electrolyte. If ‘incompressible’, the system is solved for an incompressible electrolyte ! Only implemented for incompressible mixtures, yet
Lambda2 (float) – Dimensionless parameter
a2 (float) – Dimensionless parameter
number_cells (int) – Number of cells in the mesh
solvation (float, optional) – solvation number, by default 0
PoissonBoltzmann (bool, optional) – Solve classical Nernst-Planck model with the use of the Poisson-Boltzmann formulation if True, else solve the presented model by Dreyer, Guhlke, Müller, by default False
relax_param (float, optional) – Relaxation parameter for the Newton solver xₙ₊₁ = γ xₙ f(xₙ)/f’(xₙ) with γ the relaxation parameter , by default None -> Determined automatically
Lx (float, optional) – Length of domain in x-direction, by default 2
Ly (float, optional) – Length of domain in y-direction, by default 10
rtol (float, optional) – Relative tolerance for Newton solver, by default 1e-8
max_iter (float, optional) – Maximum number of Newton iterations, by default 500
return_type (str, optional) – Type of return value, by default ‘Vector’ ‘Vector’ -> Returns the solution as a numpy array for triangle elements ‘Extended’ -> Returns the solution as both, a numpy array for triangle elements and the fenicsx function u
- Returns:
Returns atomic fractions for species A and C, electric potential, pressure, and the mesh as numpy arrays for triangle elements. x_vector is a list of both dimensions [x, y] If return_type is ‘Extended’, also returns the fenicsx function u
- Return type:
y_A_vector, y_C_vector, phi_vector, p_vector, x_vector
Jan Habscheid Jan.Habscheid@rwth-aachen.de
- Helper_DoubleLayerCapacity.C_dl(Q_DL: ndarray, Phi_pot: ndarray) ndarray
Double Layer Capacity
- Parameters:
Q_DL (np.ndarray) – Charge of the system
Phi_pot (np.ndarray) – Electric potential values
- Returns:
Double Layer Capacity
- Return type:
np.ndarray
- Helper_DoubleLayerCapacity.Phi_pot_center(Phi_pot: ndarray) ndarray
Returns vector with the center of the electric potential values.
- Parameters:
Phi_pot (np.ndarray) – Input vector with electric potential values.
- Returns:
Vector with the center of the electric potential values.
- Return type:
np.ndarray
- Helper_DoubleLayerCapacity.Q_DL_dim_ana(y_A_R: float, y_C_R: float, y_N_R: float, z_A: float, z_C: float, phi_L: float, phi_R: float, p_R: float, K: str | float, Lambda2: float, a2: float, nR_m: float, e0: float, LR: float, solvation: float) float
Calculates charge of the system using the analytical method in dimensionless units
Q = sgn(φᴸ −φᴿ)λ√(2/(κ+1))(√(Λᴸ)−√(Λᴿ))
with Λ = ln(∑_α D_α exp(−z_α φ + Eₚ a²)
- Parameters:
y_A_R (float) – Value of anion fraction at the right boundary
y_C_R (float) – Value of cation fraction at the right boundary
y_N_R (float) – Value of neutral fraction at the right boundary
z_A (float) – Charge number of anions
z_C (float) – Charge number of cations
phi_L (float) – Value of electric potential at the left boundary
phi_R (float) – Value of electric potential at the right
p_R (float) – Value of pressure at the right boundary
K (str | float) – Bulk modulus, use ‘incompressible’ for an incompress
Lambda2 (float) – Dimensionless parameter
a2 (float) – Dimensionless parameter
nR_m (float) – Reference number density in 1/m^3
e0 (float) – Dielectric constant
LR (float) – Reference length in m
solvation (float) – Solvation number
- Returns:
Charge of the system in µAs/cm³
- Return type:
float
- Helper_DoubleLayerCapacity.Q_DL_dimless_ana(y_A_R: float, y_C_R: float, y_N_R: float, z_A: float, z_C: float, phi_L: float, phi_R: float, p_R: float, K: str | float, Lambda2: float, a2: float, solvation: float) float
Calculates charge of the system using the analytical method in dimensionless units
Q = sgn(φᴸ −φᴿ)λ√(2/(κ+1))(√(Λᴸ)−√(Λᴿ))
with Λ = ln(∑_α D_α exp(−z_α φ + Eₚ a²)
- Parameters:
y_A_R (float) – Value of anion fraction at the right boundary
y_C_R (float) – Value of cation fraction at the right boundary
y_N_R (float) – Value of neutral fraction at the right boundary
z_A (float) – Charge number of anions
z_C (float) – Charge number of cations
phi_L (float) – Electric potential at the left boundary [V]
phi_R (float) – Electric potential at the right boundary
p_R (float) – Pressure at the right boundary
K (str | float) – Bulk modulus, use ‘incompressible’ for an incompressible mixture
Lambda2 (float) – Dimensionless parameter
a2 (float) – Dimensionless parameter
solvation (float) – Solvation number
- Returns:
Charge of the system in dimensionless units
- Return type:
float
- Helper_DoubleLayerCapacity.Q_num_(y_A: ndarray, y_C: ndarray, n: ndarray, x: ndarray, z_A: float = -1.0, z_C: float = 1.0) float
Calculates the charge of the system
Q = ∫_Ω n^F dΩ
- Parameters:
y_A (np.ndarray) – Anion fraction
y_C (np.ndarray) – Cation fraction
n (np.ndarray) – Total number density
x (np.ndarray) – Spatial discretization
z_A (float, optional) – Charge number of anions, by default -1.0
z_C (float, optional) – Charge number of cations, by default 1.0
- Returns:
Charge of the system
- Return type:
float
- Helper_DoubleLayerCapacity.Q_num_dim(y_A: ndarray, y_C: ndarray, n: ndarray, x: ndarray, z_A: float, z_C: float, nR_m: float, e0: float, LR: float) float
Calculates the charge of the system
Q = ∫_Ω n^F dΩ
- Parameters:
y_A (np.ndarray) – Anion fraction
y_C (np.ndarray) – Cation fraction
n (np.ndarray) – Total number density
x (np.ndarray) – Spatial discretization
z_A (float, optional) – Charge number of anions, by default -1.0
z_C (float, optional) – Charge number of cations, by default 1.0
nR_m (float) – Reference number density in 1/m^3
e0 (float) – Dielectric constant
LR (float) – Reference length in m
- Returns:
Charge of the system
- Return type:
float
- Helper_DoubleLayerCapacity.dx(Phi_pot: ndarray) float
Returns the difference between the first two electric potential values.
Assumes that the electric potential values are equally spaced.
- Parameters:
Phi_pot (np.ndarray) – Input vector with electric potential values.
- Returns:
Difference between the first two electric potential values.
- Return type:
float
- Helper_DoubleLayerCapacity.n(p: ndarray, K: str | float) ndarray
Calculates the total number density
- Parameters:
p (np.ndarray) – Pressure
K (str | float) – Bulk modulus, use ‘incompressible’ for an incompressible mixture
- Returns:
Total number density
- Return type:
np.ndarray