Nonlocal Response in Electrolytic Cells: A Generalized Poisson–Nernst–Planck Model with Memory Effects

TL;DR

This paper introduces a new model for electrolytic systems that explains how memory effects influence ion movement and electrical impedance.

Contribution

The novel contribution is a generalized Poisson–Nernst–Planck model with memory effects to describe nonlocal ionic transport and anomalous diffusion.

Findings

The model shows that memory effects can cause non-Debye relaxation and fractional-like scaling in impedance.

The model accurately fits experimental impedance data from NH4Cl-glycerol solutions.

Memory-driven transport is shown to govern transitions between normal and anomalous diffusion regimes.

Abstract

We present an extension
of the standard Poisson–Nernst–Planck
model by incorporating temporal memory effects to describe the spectroscopy
impedance response in electrolytic systems. This model yields a modified
current-density relation in which the ionic flux depends nonlocally
on the applied electric field. The resulting electrical impedance
may exhibit non-Debye relaxation and fractional-like scaling at low
frequencies, providing a basis for anomalous diffusion in confined
electrolytes. We analyze impedance spectroscopy data from NH4Cl-glycerol solutions for various concentrations to validate the
model. The comparison demonstrates how the memory kernel governs the
transition between normal and anomalous diffusion regimes, enabling
accurate fits to experimental data. These results evidence the relevance
of memory-driven transport in complex fluids and suggest a pathway
to unify standard and fractional impedance models.

Figures

Real and imaginary parts of the impedance for different values of γ by considering μ­(iω) = 1 + (iωτ)γ–1/[1 + (iωτ)γ]. We consider, for simplicity, D = 2 × 10–9 m/s2, τ = 10 s, d = 10–3 m, S=π10−3m2 , ε = 90 ε 0 (ε 0 = 8.85 × 10–12 F/m), and λ = 9.97 × 10–7 m.

Nyquist diagram and electric conductivity for different values of γ by considering μ­(iω) = 1 + (iωτ)γ–1/[1 + (iωτ)γ]. We consider, for simplicity, D = 2 × 10–9 m/s2, τ = 10 s, d = 10–3 m, S=π10−3m2 , ε = 90 ε 0 (ε 0 = 8.85 × 10–12 F/m), and λ = 9.97 × 10–7 m.

Real and imaginary parts of the impedance for different values of γ by considering μ­(iω) = (iωτ)γ/[1 + (iωτ)γ]. We consider, for simplicity, D = 2 × 10–9 m/s2, τ = 10 s, d = 10–3 m, S=π10−3m2 , ε = 90 ε 0 (ε 0 = 8.85 × 10–12 F/m), k = 1 m/s, and λ = 9.97 × 10–7 m.

Nyquist diagram and electric conductivity for different values of γ by considering μ­(iω) = (iωτ)γ/[1 + (iωτ)γ]. We consider, for simplicity, D = 2 × 10–9 m/s2, τ = 10 s, d = 10–3 m, S=π10−3m2 , ε = 90 ε 0 (ε 0 = 8.85 × 10–12 F/m), k = 1 m/s, and λ = 9.97 × 10–7 m.

Panel a) shows a schematic representation of the circular steel electrode used in our experiments. Panel b) shows a schematic representation of electrophoretic ion migration under an alternating current (AC) field. Positively charged ions (cations, red spheres) migrate toward the negative electrode (cathode), while negatively charged ions (anions, blue spheres) migrate toward the positive electrode (anode). The polarity of the electrodes switches periodically with the applied AC, causing the ions to oscillate.

Real and imaginary parts of the impedance for different salt concentrations in glycerin. The red-dotted lines represent the experimental data, and the black solid lines represent the model used to describe the experimental data, with μ­(iω) = (iωτ)γ/[1 + (iωτ)γ]. For S1 the parameters used were ε ≈ 46ε 0, τ ≈ 10 s, γ ≈ 0.73, and λ ≈ 4.76 × 10–7 m. For S2 the parameters used were ε ≈ 43ε 0, τ = 4.0 s, γ = 0.66, and λ ≈ 2.34 × 10–7 m. For S3 the parameters used were ε ≈ 41 ε 0, τ ≈ 7.8 × 10–2 s, γ ≈ 0.75, and λ ≈ 4.8 × 10–8 m. For these cases, d = 10–3 m, S=π10−3m2 , ε 0 = 8.85 × 10–12 F/m, k = 1.0 m/s, and D ≈ 2 × 10–9 m/s2.

Electrical conductivity and Nyquist diagram for different salt concentrations in glycerin. The red-dotted lines represent the experimental data, and the black solid lines represent the model used to describe the experimental data, with μ­(iω) = (iωτ)γ/[1 + (iωτ)γ]. The parameter values are the same as Figure .

Full text

Introduction

1

Anomalous diffusion is a common phenomenon in complex systems. ?−? ? It deviates significantly from normal Brownian motion, which is characterized by Markovian processes and, consequently, by a linear time dependence for the mean-squared displacement. Anomalous diffusion exhibits a nonlinear time dependence of the mean-square displacement, thereby characterizing subdiffusion and superdiffusion. These behaviors have been observed in different contexts, and many approaches have been considered to analyze them. ?−? ? ? ? For example, subdiffusive behavior has been consistently observed in live cell imaging of proteins and organelles due to macromolecular crowding and viscoelasticity of the cytoplasm. ?,? Anomalous diffusion models describe the non-Gaussian transport of neurotransmitters and signaling molecules in dendritic spines and synaptic clefts. ?,? Fractional dynamics capture the aging and memory effects in complex fluids, polymer networks, and colloidal glasses. ?,? Environmental science also employs these models to describe the subdiffusive migration of contaminants in fractured rocks and soils, where traditional Fickian laws fail. ?,? Furthermore, in finance and human mobility, Lévy-like superdiffusion and heavy-tailed processes have been used to model asset price fluctuations and urban movement patterns. ?,? These applications are increasingly supported by models based on fractional partial differential equations, ?−? ? ? generalized Langevin equations, ?−? ? and continuous-time random walks (CTRWs), ?,? providing powerful tools to incorporate memory, heterogeneity, and nonlocality in both space and time. Such mathematical frameworks are essential to accurately describe the underlying dynamics observed in complex systems ranging from crowded intracellular environments,? to ionic and subnuclear transport,? to pollutant migration in heterogeneous media,? and even to anomalous pattern recognition and classification in experimental data.? Complex systems, such as confined electrolytes, viscoelastic media, biological membranes, polymer gels, and ionic liquids, often exhibit nonlocal behaviors over time, manifested by slow relaxation, hysteresis, and anomalous transport. ?−? ? On the other hand, the Poisson–Nernst–Planck (PNP) model gives a fundamental description of ionic transport in electrolytic media and is widely used in various contexts, including electrochemistry, cell biology, and materials physics. ?,? However, this model assumes that ionic fluxes respond instantaneously to applied electric fields and concentration gradients, implying a local relationship between force and flux. To capture such effects that are not suitably described in all frequency ranges, it is necessary to incorporate other mechanisms for the dynamics of ions, for example, by using time-fractional derivatives in continuity equations, ?,? generalized Langevin or fractional PNP formulations. ?,?,?

Here, we analyze an extension of the Poisson–Nernst–Planck (PNP) model by modifying the drift term to include a nonlocal contribution, which can be related to the memory effects and noninstantaneous processes of relaxation associated with the electric field within the medium. To accomplish this task, we have structured our papers as follows. In Sec. 2, we extend the PNP model to incorporate anomalous behavior and in its framework we obtain closed expressions for the electrical impedance in the ac small–signal limit. In section, we compare the results of section with the experimental data obtained for the electrolytic cells prepared by dissolving ammonium chloride (NH_4_Cl) in glycerol (Sigma-Aldrich). In section, we present our conclusions and discussions.

Extension of the PNP-Model

2

Before starting our discussion about the PNP-model, we consider the following scenario: the dynamics of an ion with charge q and velocity v(t), under the influence of a time-dependent electric field E(t), is described by the generalized Langevin equation:

$$[eqn]$$

where m represents the mass of the ion under the action of an external electric field E(t), Γ­(t) is a kernel function representing the friction with memory effects and ξ­(t) is the random force with ⟨ξ­(t)⟩ = 0 and ⟨ξ­(t) ξ­(t′)⟩ = k _ B _ TΓ(|t – t′|) . For this, we consider the ensemble average, which yields

$$[eqn]$$

This equation can be solved using the Laplace transform $[eqn]$ , which, in Laplace space, results in the following:

$$[eqn]$$

with

$$[eqn]$$

Applying the inverse Laplace transform:

$$[eqn]$$

From this equation, we can derive the current density associated with the motion of these ions, i.e., $[eqn]$ , where $[eqn]$ represents the ion concentration, which implies

$$[eqn]$$

Equation shows that in the case of ions subjected to a medium whose dynamic aspect is governed by a generalized Langevin equation, the current density has a nonlocal dependence on the applied electric field, which, for the particular case μ­(t) ∝ δ­(t) (where δ­(t) is the Dirac’s delta) leads to obtaining $[eqn]$ . Thus, different choices of μ­(t) correspond to different friction kernels Γ­(t), capturing the properties of the medium such as viscoelasticity of the host medium and structural heterogeneity.

On the other hand, for the PNP model, we notice that the current density is given by

$$[eqn]$$

where D ± is the diffusion coefficient related to the positive and negative ions. The first term corresponds to the diffusion of the ions, and the second term corresponds to the interaction of the ions with the electric field determined by the Poisson equation.

From eqs and ?, one can observe that in the first case, the medium influences the current density through μ­(t), which defines the relaxation process. It can be connected to hydrodynamic effects with memory and correlated with stochastic noise inducing nonlocal temporal responses in the system.? In the second case, the properties of the medium can be directly related to the diffusion coefficients of each ion in the bulk, and the term associated with the electric field is not influenced by the medium, as in the previous case, taking into account memory and hydrodynamic effects. These two models can be combined to result in the current density:

$$[eqn]$$

which incorporates the effects of the eq into eq. Note that we are using the limit of integrations from −∞ to t (following ref ?) since we are interested in the solution in a steady-state process as a result of a periodic applied potential (see below). The memory kernel introduces a delayed response that reflects the medium’s intrinsic microscopic features. This nonlocal coupling can be related to viscoelastic rearrangements, correlated thermal fluctuations, solvent reorganization, and heterogeneous local mobility. As a result, the ionic flux is no longer driven solely by the instantaneous electric field but by its entire weighted history, giving rise to non-Debye relaxation in the impedance and to subdiffusive transport in the low-frequency limit.

To obtain the linear response for a system subjected to a periodic external potential in this scenario, we consider eq combined with the equations

$$[eqn]$$

and

$$[eqn]$$

and appropriate boundary conditions. First, let us consider an electrolyte confined in a sample of thickness d separated by two electrodes located at z = −d/2 and z = d/2. The sample is subjected to an alternating electric field and the system is governed by eq, the continuity equation, eq, and Poisson’s equation, eq. Considering the ac small–signal limit, we have $[eqn]$ , with $[eqn]$ (where $[eqn]$ ) and Φ­(z,t) = ϕ­(z)­e^iω t ^ (with Φ­(±d/2, t) = ± (V 0/2)­e^iωt ^). For perfect blocking electrodes, i.e., $[eqn]$ , it is possible to show that the electrical impedance is given by

$$[eqn]$$

where $[eqn]$ is the surface area of the electrode, μ­(iω) is the kernel connected with the relaxation process of the system, where μ­(iω) = ∫0 ^∞^ dt′ e^–iωt′^μ­(t′), $[eqn]$ , D + = D – = D, and $[eqn]$ is the Debye length. Equation recovers the standard case presented in ref ? for the PNP model with perfect blocking electrodes for μ­(iω) = 1, μ­(t) = δ­(t), and in the asymptotic limit of low frequency, we have

$$[eqn]$$

With respect to μ­(iω), it represents the physical process for ionic drift in a given medium. Hence, it is a phenomenological parameter chosen to represent the interaction between the ion and the medium. For example, for 1/μ­(iω) ≈ (iωτ)^γ^ (τ is a relaxation time) in the limit ω → 0 with 0 < γ < 1, we can approximate the previous result to

$$[eqn]$$

where $[eqn]$ and C _ bulk _ = εS/(2λ). This behavior has also been obtained in refs ?−? ? by considering the fractional approach, i.e., the fractional time derivative in bulk or related to the surface effects. It should be noted that this behavior for the electrical impedance connected to the memory kernel μ­(t) exhibits features analogous to those of the constant phase elements (CPE) in the asymptotic limit of low frequency, commonly used in electrical spectroscopy impedance to account for distributed relaxation times and surface heterogeneities. Furthermore, according to the developments described in refs ? and ?

$$[eqn]$$

which implies a connection between the electrical conductivity, σ­(ω), and the mean square displacement. This result implies that ⟨(z – ⟨z⟩)^2^⟩ ∼ t ^γ^ for Re­[σ­(iω) ] ∼ ω^γ^ (0 < γ < 1); i.e., the underlying process is subdiffusive. It can be related to a random walk process asymptotically characterized by a long-tailed waiting time distribution. It is worth mentioning that Re­[σ­(iω)] = constant implies ⟨(z – ⟨z⟩)^2^⟩ ∼ t (usual diffusion), which may be connected to a Poisson distribution for the waiting time. These arguments can be analyzed from the continuous time random walk approach, for example, by following the results presented in ref.,? which allows us to connect the second moment with the mean square displacement for a separable probability density function, ?,? i.e., ψ­(z, t) = λ­(z) φ­(t), as follows: ⟨(z – ⟨z⟩)^2^⟩ ∝ φ̃(s)/[s(1 – φ̃(s))], where φ̃(s) is the waiting time distribution in the Laplace space. We can connect this result with the electrical conductivity, yielding: σ­(iω) ∝ iωφ̃(iω)/[1 – φ̃(iω)], which for a long-tailed distribution, i.e., φ̃(iω) ∼ 1 – (iω)^γ^, implies in Re­[σ­(iω) ] ∼ ω^1−γ^ in the asymptotic limit of low frequency and, consequently, in a subdiffusion. For the case characterized by the eq, we obtain σ­(ω) ≈ σ_0_(ωτ)^1−γ/2^ with $[eqn]$ , which allows us to show that ⟨(z – ⟨z⟩)^2^⟩ ∝ t ^γ/2^, where 0 < γ < 1, corresponding to a subdiffusion. This discussion demonstrates that the parameter γ associated with the kernel, beyond promoting the impedance connection with CPE elements, can capture a broad spectrum of time scales associated with structural heterogeneity and viscoelasticity, reflecting ion motion.

Figure shows the behavior of the real and imaginary parts of the impedance when μ­(iω) = a + (iωτ)^γ–1^/[a + (iωτ)^γ^] for different values of γ (0 < γ ≤ 1). We also show the standard PNP result, as indicated by the dashed black curve. For γ = 1.0, and a = 1 we obtain μ­(iω) = 1 + (iωτ)/[1 + (iωτ)] and τ represents a relaxation time. This kernel represents a process with two contributions: one localized in time, described by the parameter a, and the second, which represents an exponential decay process (Debye relaxation) with decay time τ.? This decay is associated with non-Markovian behavior in the drift process, characterized by a memory effect. This kind of kernel is particularly important in the low-frequency regime. We also observe that the real part of the impedance behaves similarly to what is observed for the PNP model in the presence of adsorption and desorption.? The solid red curve is calculated for γ = 0.4 and a = 1, so μ­(iω) = 1 + (iωτ)^−0.6^/[1

• (iωτ)^0.6^]. In this case, the kernel replaces exponential decay by a nonexponential decay process. In Figure, we present the Nyquist diagram and the electrical conductivity, using the same kernels as in Figure.

Real and imaginary parts of the impedance for different values of γ by considering μ­(iω) = 1 + (iωτ)γ–1/[1 + (iωτ)γ]. We consider, for simplicity, D = 2 × 10–9 m/s2, τ = 10 s, d = 10–3 m, S=π10−3m2 , ε = 90 ε 0 (ε 0 = 8.85 × 10–12 F/m), and λ = 9.97 × 10–7 m.

Nyquist diagram and electric conductivity for different values of γ by considering μ­(iω) = 1 + (iωτ)γ–1/[1 + (iωτ)γ]. We consider, for simplicity, D = 2 × 10–9 m/s2, τ = 10 s, d = 10–3 m, S=π10−3m2 , ε = 90 ε 0 (ε 0 = 8.85 × 10–12 F/m), and λ = 9.97 × 10–7 m.

Now, let us consider the boundary condition,

$$[eqn]$$

which corresponds to the Chang–Jaffe boundary condition. ?−? ? ? The Chang–Jaffe boundary condition, eq, is often used to represent specific ion adsorption.? In eq, k represents the adsorption rate and $[eqn]$ is the concentration at thermodynamic equilibrium. In this case, it is possible to obtain the electrical impedance in the asymptotic limit of a small-ac and show that it is given by

$$[eqn]$$

with $[eqn]$ . Figure illustrates the real and imaginary parts of the impedance. Figure illustrates the Nyquist diagram and the electrical conductivity obtained from eq.

Real and imaginary parts of the impedance for different values of γ by considering μ­(iω) = (iωτ)γ/[1 + (iωτ)γ]. We consider, for simplicity, D = 2 × 10–9 m/s2, τ = 10 s, d = 10–3 m, S=π10−3m2 , ε = 90 ε 0 (ε 0 = 8.85 × 10–12 F/m), k = 1 m/s, and λ = 9.97 × 10–7 m.

Nyquist diagram and electric conductivity for different values of γ by considering μ­(iω) = (iωτ)γ/[1 + (iωτ)γ]. We consider, for simplicity, D = 2 × 10–9 m/s2, τ = 10 s, d = 10–3 m, S=π10−3m2 , ε = 90 ε 0 (ε 0 = 8.85 × 10–12 F/m), k = 1 m/s, and λ = 9.97 × 10–7 m.

These figures show that the nonlocal effects introduced in the drift term have a significant influence on the system behavior in the asymptotic limit of low frequency, where the diffusion process has a pronounced impact. Similarly to the previous case, the behavior exhibited in Figures and ? can be attributed to anomalous diffusion.

Experimental Data and Impedance

3

We prepared samples by dissolving ammonium chloride (NH_4_Cl) in glycerol (Sigma-Aldrich) to compare our method with experimental data. Three different concentrations (all by weight %) were tested: (i) S_1_ was made with 0.0007% of NH_4_Cl; (ii) S_2_ containing 0.0048% of NH_4_Cl; and (iii) S_3_ with 0.1% of NH_4_Cl. Electrochemical impedance spectroscopy measurements were performed to compare the effect of the salt. Measurements were taken using a Hioki IM3533 LCR meter at room temperature, with a frequency range of 0.01 Hz to 200 kHz, and a small voltage of 25 mV. A circular sanded stainless steel electrode with a radius of 10 mm and a spacing of 1 mm between surfaces (see Figurea,b), was immersed in a container of 10 mL filled with the selected sample.

Panel a) shows a schematic representation of the circular steel electrode used in our experiments. Panel b) shows a schematic representation of electrophoretic ion migration under an alternating current (AC) field. Positively charged ions (cations, red spheres) migrate toward the negative electrode (cathode), while negatively charged ions (anions, blue spheres) migrate toward the positive electrode (anode). The polarity of the electrodes switches periodically with the applied AC, causing the ions to oscillate.

The dielectric behavior of different NH_4_Cl concentrations in glycerol was investigated using the extended PNP approach, presented in eq. The real part of the impedance is directly related to the system’s Ohmic resistance, which is influenced by the salt concentration, as evidenced by the plateau displacement observed in Figure. At lower frequencies, surface processes, adsorption, and desorption take place and are represented by parameters k. The diffusion behavior within the system can be analyzed by comparing the samples S_1_, S_2_, and S_3_. For S_1_, the parameters were ε ≈ 46 ε 0, τ ≈ 10 s, γ ≈ 0.73, and λ ≈ 4.76 × 10^–7^ m, while for S_2_, ε ≈ 43ε 0, τ = 4.0 s, γ = 0.66, and λ ≈ 2.34 × 10^–7^ m. For the highest concentration sample, S_3_, the parameters were ε ≈ 41ε 0, τ ≈ 7.8 × 10^–2^ s, γ ≈ 0.75, and λ ≈ 4.8 × 10^–8^ m. These variations in λ reflect the salt concentration increases, which are inversely proportional to the number of charges. Also noteworthy is the fact that the dielectric constant decreases with the concentration of added salt. This dielectric decrement is directly related to the suppression of the collective response of the glycerol molecules.? Note also that the experimental impedance data show that in glycerol-based solutions, the hydrogen-bonding network introduces delays in solvent reorganization, resulting in noninstantaneous transport. ?,?,?

Figure shows the electrical conductivity and Nyquist diagrams for the three NH_4_Cl concentrations. The theoretical model accurately reproduces the experimental behavior through the memory kernel μ­(iω), which accounts for the distribution of relaxation times of glycerol-based electrolytes. This point is evident in the Nyquist representation, where the experimental data exhibit semicirculus and noncircular behaviors in the low-frequency, which is signature of non-Debye relaxation. These effects significantly influence the diffusion process and give rise to a memory effect in the ionic response. In addition, the viscosity and local mobility heterogeneities contribute to distributed relaxation times, which can be represented by power-law-type kernels. ?,? These features have implications for the kernel that must capture this behavior and provide a phenomenological explanation for the observed results.

Real and imaginary parts of the impedance for different salt concentrations in glycerin. The red-dotted lines represent the experimental data, and the black solid lines represent the model used to describe the experimental data, with μ­(iω) = (iωτ)γ/[1 + (iωτ)γ]. For S1 the parameters used were ε ≈ 46ε 0, τ ≈ 10 s, γ ≈ 0.73, and λ ≈ 4.76 × 10–7 m. For S2 the parameters used were ε ≈ 43ε 0, τ = 4.0 s, γ = 0.66, and λ ≈ 2.34 × 10–7 m. For S3 the parameters used were ε ≈ 41 ε 0, τ ≈ 7.8 × 10–2 s, γ ≈ 0.75, and λ ≈ 4.8 × 10–8 m. For these cases, d = 10–3 m, S=π10−3m2 , ε 0 = 8.85 × 10–12 F/m, k = 1.0 m/s, and D ≈ 2 × 10–9 m/s2.

Electrical conductivity and Nyquist diagram for different salt concentrations in glycerin. The red-dotted lines represent the experimental data, and the black solid lines represent the model used to describe the experimental data, with μ­(iω) = (iωτ)γ/[1 + (iωτ)γ]. The parameter values are the same as Figure .

Discussions and Conclusions

4

We have extended the standard Poisson–Nernst–Planck (PNP) model to include temporal memory effects by means of a nonlocal relaxation kernel. By deriving a modified current-density relation from the generalized Langevin equation, we evidenced that ionic flux in electrolytic cells may depend nonlocally on the applied electric field, leading to a generalized impedance expression that accounts for non-Debye relaxation phenomena. The theoretical framework effectively reproduces the main features observed in experimental impedance spectroscopy (EIS) data obtained from the NH_4_Cl–glycerol solutions across a broad frequency range. Specifically, the model captures the transition between normal diffusion and anomalous subdiffusive regimes, as indicated by the scaling behavior of the memory kernel and the impedance response. This implies that memory-driven ionic transport may play an important role in confined electrolytic systems. Note also that the connection of the electrical conductivity with the mean square displacement and, consequently, with the waiting time distributions allows us to connect the diffusion of the ions with an anomalous scenario characterized by subdiffusion in the asymptotic limit of low frequency. The approach developed here bridges traditional diffusion models with more phenomenological descriptions that include nonlocal components, providing a rigorous yet flexible framework for understanding anomalous transport. While fractional differential operators are not directly employed in our formulation, the memory kernel can simulate fractional-like behaviors commonly associated with complex fluids.