Slide electrification charge can exist without residual liquid film: KPFM measurements near receding three-phase contact lines

Abstract

Spontaneous charge transport near the solid–gas–liquid three-phase contact line (TPCL) is a key phenomenon for interface-driven power generation and related applications, yet its details—especially near a receding TPCL—remain poorly understood. In this study, we experimentally investigated the charge distribution induced by the recession of the TPCL. Specifically, we mapped the charged regions left by glycerol–water droplets on SiO₂ surfaces using Kelvin probe force microscopy (KPFM). It was found that the charge distribution at solid–liquid interfaces near the TPCL is determined by the Debye length of the solution used, which is consistent with a theoretical model proposed in prior studies. However, nanodroplets or residual films, which have been proposed as necessary for the stabilization of charges at the solid–gas interface, were not observed, implying the existence of a different charge stabilization mechanism. Furthermore, we found that the charge distribution unpredictably exists on the solid–gas interface side, possibly induced by charge relaxation. This deviation between the KPFM results and theoretical prediction was resolved by newly introducing a charge-relaxation term to the model. These results deepen our understanding of charge-transfer mechanisms at TPCLs and provide valuable insights for interface-chemistry applications such as energy harvesting.

---

Introduction

Contact electrification at solid–liquid interfaces is widespread and has found applications in diverse fields, such as wettability control and energy harvesting.^1–11 In particular, spontaneous charge transfer near the solid–gas–liquid three-phase contact line (TPCL) has attracted increasing attention over the past decade because of its potential for droplet and flow manipulation,^12–16 the management of electrostatic effects in microfluidic and semiconductor liquid processes,^17–19 and kV-class energy harvesting driven by flows, evaporation, and droplet motion.^20–29 Therefore, a deeper understanding of this phenomenon is crucial for advancing sustainable technologies.

In recent years, it has become clear that charge transfer associated with a receding TPCL plays an important role in these phenomena. A receding TPCL continuously exposes fresh solid and collapses the nearby electrical double layer (EDL); this repeated cycle transports counterions along the moving line and leaves a net charge behind it. As a result, the polarity and magnitude of the generated potential depend on the retreat speed, wetting hysteresis, and Debye screening at the interface—rather than on bulk flow alone. Consistent with this view, studies involving sliding droplets and natural evaporation on low-permittivity and porous substrates have reported variations in wetting dynamics, charge separation, and electromotive force induced by TPCL motion.^{14,16,21,25,28,30} Together, these findings suggest that the receding contact line may act as the microscopic engine that couples interfacial mechanics to charge generation, governing droplet dynamics, slide electrification, and hydrovoltaic effects. However, the detailed mechanism of charge transport in its immediate vicinity remains unclear.

Recently, a model has been proposed to predict the surface charge density on the solid–gas surface left behind by a sliding droplet. The model also provides an analytic expression for changes in the solid–liquid interfacial charge that accounts for EDL distortion near the contact line by introducing the parameters contact angle (θ), Péclet number (Pe), and Debye length.^30–33 While this framework captures experimental trends, it has been reported to overestimate measured values due to atomistic-scale effects.³³ To reconcile this gap, an effective charge-transfer coefficient has been introduced based on hypotheses such as charge regulation by nanodroplets, residual liquid films, water layers or hydration-shell-like structures. However, their existence is still an assumption and has not yet been experimentally confirmed, largely due to the scarcity of direct nanoscale measurements near a receding contact line. Therefore, further model validation based on spatially resolved measurements is necessary.¹¹

In this study, we induced the recession of the TPCL of glycerol–water droplets on SiO₂ substrates (and, for comparison, on SiN substrates under selected conditions) to expose the contact-electrified area. We subsequently measured the surface potential and charge-density profile over that area using Kelvin probe force microscopy (KPFM). By varying the glycerol–water composition, we tuned the bulk Debye length (λ_bulk) and examined its effect on interfacial charging. By combining these results with theoretical analysis and a newly introduced charge-relaxation model on the solid–gas interface, we obtained a more quantitative picture of solid–liquid interfacial charging near the TPCL.

Experimental

To examine solid–liquid charging around the TPCL, we compared the charging behavior near a receding contact line using droplets of three glycerol–water solutions with glycerol weight fractions of 93, 86, and 60 wt% and assessed the λ_bulk dependence of interfacial charging. All experiments were conducted at room temperature (T ≈ 293 K and RH ≈ 30%) under ambient conditions. First, glycerol and water were mixed to prepare each solution, and 3 µL of the resulting glycerol–water solution was dispensed onto Si substrates coated with 300 nm of thermally grown SiO₂ (Fig. 1a and c). The substrate surface was first measured briefly prior to TPCL recession. The TPCL was then retracted by gently aspirating the droplet with a pipette, thereby exposing the previously wetted area (Fig. 1d). After the TPCL had receded, KPFM mapping was initiated within 5 min over the exposed region, including the former TPCL position (Fig. 1b and e). The acquisition of a single KPFM image took approximately 15 min. Here, the Si base substrate was grounded through the sample stage and used as the electrical reference. The measured potentials were converted into surface charge density σ_KPFM using the parallel-plate approximation for an order-of-magnitude estimation (see SI Note 1). We also measured the contact angle θ for each droplet.

Fig. 1 (a) and (b) Optical images acquired immediately before and after TPCL recession, respectively. Each scale bar is 200 µm. (c)–(e) Schematic showing the droplet being receded on the substrate (TPCL recession) and the exposed surface being measured by KPFM.

Results and discussion

Fig. 2 shows a distinct boundary between regions of positive and negative surface charges in the KPFM surface-charge maps (Fig. 2a–c). Additional data for an intermediate composition (80 wt%) are provided in SI Note 2 (Fig. S3). Away from the boundary, the previously wetted sides were consistently charged strongly. For convenience, we hereafter refer to the charged side as the “wetted region” and the opposite side across the boundary as the “dry region.”

Fig. 2 KPFM experimental results for glycerol–water solution droplets on SiO₂ substrates (93, 86, and 60 wt%) and a SiN substrate. Results for 80 wt% glycerol are shown in SI Note 2 (Fig. S3). For the SiO₂ cases, the calculated bulk Debye lengths are λ_bulk ≈ 1020 nm (93 wt%), 850 nm (86 wt%), and 520 nm (60 wt%); λ_bulk was estimated as described in SI Note 6. Surface charge-density maps for (a) SiO₂/93 wt%, (b) SiO₂/86 wt%, (c) SiO₂/60 wt%, and (d) SiN/approx. 80 wt%. Darker contrast indicates more negative charge in (a)–(d). (e)–(h) Corresponding topography images acquired over the same areas as in (a)–(d), respectively. Each scale bar is 10 µm. (i)–(l) Line profiles along the red A–B lines in (a)–(d), respectively. The light-orange shaded band indicates the W_25–75 transition width. The original CPD/surface-potential maps and profiles are provided in SI Note 1.

KPFM mapping of the wetted region showed that, for all three glycerol–water compositions, a pronounced negative charge appeared relative to the dry region. The topography images (Fig. 2e–g) show no distinct structures across the boundary region, except for occasional contamination and minor feedback perturbations due to electrostatic attraction. We extracted line profiles and defined a coordinate axis from the negative side (A) to the positive side (B) relative to the TPCL. In the line profiles of surface charge density (Fig. 2i–k), Position = 0 was defined as the point of maximum surface-charge gradient. To quantify the steepness of the charge-density transition, we defined the 25–75% transition width, W_25–75, as the distance over which the surface charge density changes from 25% to 75% of its total variation across the analyzed window (−10 to +10 µm). We obtained W_25–75 = 2.41 µm (93 wt%), 1.65 µm (86 wt%), and 1.04 µm (60 wt%), indicating that a spatially localized charging transition forms near the TPCL and becomes relatively steeper (i.e., narrower) as the water content increases. Thus, the steepness of the localized charging transition near the TPCL depends on the droplet composition.

Meanwhile, surface potential and topography images measured near the TPCL without retraction (Fig. 3) show no marked variation in surface charge even within a few micrometers of the TPCL; the surface charge density was nearly uniform. Hence, the gentle decrease of the charging gradient for Position > 0 in Fig. 2i–k should be attributed to processes occurring after the TPCL receded. The profile can be treated as quasi-static during one image acquisition (see SI Note 11).

Fig. 3 (a) KPFM surface-potential image acquired within a few micrometres of TPCL. (b) Corresponding topography image over the same area as in (a). Each scale bar is 10 µm. The droplet is on the left-hand side of the image.

We first consider the origin of the charges induced by contact electrification. When an SiN substrate was used instead of SiO₂, the wetted region exhibited positive, rather than negative, charging (Fig. 2d, h and l). This polarity reversal can be explained by the difference in the points of zero charge (PZC) of the two surfaces: SiO₂, which is rich in silanol groups, has a PZC around pH 2.5–3.5, whereas SiN, which is rich in amino groups, has a PZC around pH 6.5.^34,35 In our glycerol–water solutions, ions arise mainly from water autodissociation, carbonic-acid dissociation, and the weak-acid dissociation of glycerol; therefore, the solutions are expected to be weakly acidic, with pH close to 6. Consistent with this estimate, the pH of the glycerol–water solutions used in our experiments was measured to be approximately 6.0 using a pH meter. Because the liquid pH is above the PZC of SiO₂ and slightly below that of SiN, the SiO₂ surface is expected to be negatively charged, whereas the SiN surface can be positively charged depending on the solution conditions (see SI Notes 4 and 5). These results suggest that dissociation of surface functional groups upon solid–liquid contact leads to charge separation between the solid surface and the liquid phase.

We conclude that the “wetted region” was already dry during the measurements, based on the low noise in the topography images and the stable imaging (Fig. 2(e)–(h)). The absence of scan distortions, which are typically observed on insufficiently dried, well-wetted surfaces, also supports this interpretation (see SI Note 3). In addition, the presence of a continuous nanometer-scale wetting layer is unlikely under the present conditions, as supported by additional AFM measurements (see SI Note 3). Therefore, the observed negative surface charge is attributed to contact electrification and remains stable in the absence of hydration-shell-like structures, residual liquid films, or adsorbed water layers proposed in prior studies,^14,31,33 although this does not exclude the presence of molecularly localized water molecules on the substrate surface.

Next, to quantitatively explain the trend of the surface charge profile near the boundary, we applied a framework that combines the effects of wedge-induced EDL distortion,³² recession-driven interfacial flow,³³ and, newly introduced here, charge relaxation at the solid–gas interface.

Wedge-induced EDL distortion and flow effect

First, we consider an absolute coordinate system along the solid surface, where the TPCL is located at x = 0 and the droplet occupies the region x ≤ 0. We introduce the dimensionless coordinate X = x/λ_bulk, where λ_bulk is the Debye length of the bulk solution. The λ_bulk was calculated from ionic strength and static permittivity of solution. Then, we discuss the effective surface charge distribution σ_CA(x) near the contact line, accounting for wedge-induced distortion of the EDL. Using the charge difference Δσ, the EDL-distortion function f(X), the far-field solid–liquid surface charge density σ_SL,0, and the far-field solid–gas surface charge density σ_SG,0, the surface charge distribution σ_CA(x) is expressed as follows:where

Δσ = σ_SL,0 − σ_SG,0 (2)

andwith c = 1 + 0.36α, , , and (see SI Note 8 for further details).

Next, we extend eqn (1) by incorporating the recession-induced flow to account for the EDL-stretching effect near the contact line. σ_Pe(x) is expressed in terms of the Péclet number, denoted by Pe(x), as follows:

where
with D_ion,eff denoting the ionic diffusion coefficient and u(|x|) the near-contact-line flow velocity. The effects of changes in viscosity and other properties depending on the glycerol composition are also taken into consideration. For the contact angle, we used θ = 66°, a representative value based on the measurements on SiO₂ across all three samples. The values of λ_bulk were estimated by calculation of the solution composition and were 1020, 850, and 520 nm for 93, 86, and 60 wt%, respectively. The diffusion coefficient D_ion,eff was estimated from the solution. u(|x|) is on the order of 10⁻⁵ m s⁻¹ near the contact line and decays away from it. The calculations based on eqn (1) and (4) correspond to the green dotted and light-blue dashed curves, respectively, in Fig. 4 (see SI Notes 6, 7, 9, and 10 for further details).

Fig. 4 Measured and modeled surface charge-density profiles near the receding TPCL on SiO₂. The black curve shows the measured profile. The green dotted curve represents the model including only wedge-induced EDL distortion (eqn (1)), the light-blue dashed curve represents the model including both wedge-induced distortion and the flow effect (eqn (4)), and the red solid curve represents the full model including wedge-induced distortion, the flow effect, and charge relaxation (eqn (5) and (6)). Panels (a)–(c) correspond to 93 wt%, 86 wt%, and 60 wt% glycerol, respectively, using λ_bulk = 1020 nm, 850 nm, and 520 nm. To compare with the model coordinate x (contact line at x = 0), the measured profiles are plotted after applying a constant horizontal shift Δx to the experimental coordinate x_exp → x_exp + Δx, where Δx = −1.712 µm (93 wt%), −0.919 µm (86 wt%), and −0.752 µm (60 wt%) determined by minimizing the root mean squared error (RMSE) over −10 ≤ x ≤ 10 µm; the model x-coordinate is left unchanged. Model curves were computed using the composition-dependent D_ion,eff and the fitted relaxation length ℓ_sm: D_ion,eff = 2.93 × 10⁻¹¹ m² s⁻¹ and ℓ_sm = 1.568 µm (93 wt%); D_ion,eff = 6.23 × 10⁻¹¹ m² s⁻¹ and ℓ_sm = 1.852 µm (86 wt%) and D_ion,eff = 2.88 × 10⁻¹⁰ m² s⁻¹ and ℓ_sm = 0.608 µm (60 wt%). Results for 80 wt% glycerol are provided in SI Note 12 (Fig. S9).

In this framework, note that X is normalized by λ_bulk, which implies that a shorter λ_bulk confines the wedge-induced EDL distortion to the immediate vicinity of the contact line, whereas a longer λ_bulk leads to appreciable EDL distortion farther from the contact line.

Consistent with the experiments, the framework predicts a steeper gradient for lower glycerol fractions (i.e., shorter Debye lengths), because the suppression of solid–liquid contact charging is confined closer to the contact line. However, the calculated σ(x) remains too steep and exhibits a discontinuity at x = 0 (green dotted and light-blue dashed curves in Fig. 4). To account for this remaining mismatch, we introduce a charge-relaxation term around x = 0.

Charge relaxation

To more satisfactorily explain the experimental findings, we introduce an exponential surface–charge decay factor with a fitting parameter ℓ_sm, motivated by the fact that local charging on SiO₂ generally relaxes into the surroundings.^36,37 This relaxation term is introduced as an effective representation of boundary charge smoothing rather than a single uniquely identified microscopic mechanism. Using these parameters, we analyzed the relationship between the solution characteristics and the charging state. Accordingly, the surface–charge distribution σ_rel(x) around the contact line is given as follows:

σ_rel(x) = σ_Pe(x) + δσ(x) (5)

whereHere, denotes the left-hand limit of the derivative of σ_Pe(x) (eqn (4)) evaluated at x = 0. The quantity Q_tr represents the total transferred charge associated with the fitting parameter ℓ_sm. The fitted ℓ_sm values were 1.568, 1.852, and 0.608 µm for 93, 86, and 60 wt% glycerol, respectively. See SI Note 11 for a detailed derivation.

As shown by the red solid curves in Fig. 4, the calculations based on eqn (5) and (6) reproduce the overall trend of the measured charging profiles (the black solid curves) with comparatively high accuracy. By combining high-resolution KPFM mapping with a physics-based model, we obtain quantitative agreement with the observed gradient profiles. Specifically, the composition-dependent change in the charging gradient can be explained consistently by three factors: wedge-induced EDL distortion associated with the contact-angle geometry, a flow-induced effect near the TPCL, and charge relaxation at the solid–gas interface after dewetting. Importantly, the relaxation term enables a more realistic estimate of the charging state in the vicinity of a receding contact line, where the surface transitions from the wetted to the dry region. The remaining mismatch between the model and the experimental data may originate from uncertainties in the near-contact-line flow field, the actual dissolved CO₂ concentration, and small variations in the contact angle.

Conclusions

Our KPFM measurements reveal that slide-electrification charge can persist on the solid–gas surface even when no clear evidence of nanodroplets or a residual liquid film is observed. This finding suggests that charge stabilization does not necessarily require residual wetting structures and may involve an alternative stabilization pathway for charges left behind by a receding contact line.

Focusing on the immediate vicinity of the solid–gas–liquid TPCL, we experimentally show that three factors—Debye length, contact angle, and flow—jointly control the charging near the TPCL. When the solution composition is adjusted to vary λ_bulk, the charging gradient becomes steeper as λ_bulk decreases, in quantitative agreement with a macroscopic theory that incorporates wedge geometry through the contact angle θ and flow through the Péclet number Pe.

Beyond the solid–liquid side, we further reconcile the experimentally observed charge distribution extending into the solid–gas region by newly introducing an exponential charge relaxation on the solid–gas interface, which is absent in prior models. This extension enables us to capture the continuous decay of charge at the microscale and provides a more realistic description of charge evolution after TPCL recession.

Overall, we establish a microscale method to probe contact electrification near a receding contact line and, accordingly, develop a model that predicts the interfacial charge distribution more accurately than previous approaches. These findings deepen the physical understanding of charge transfer at TPCLs and offer actionable guidelines for applications such as droplet manipulation, mitigation of electrostatic damage, and interface-driven energy harvesting.

Author contributions

Yuki Ishihara: conceptualization, methodology, investigation, data curation, formal analysis, visualization, writing – original draft, writing – review and editing, and funding acquisition. Hideaki Teshima: supervision, project administration, conceptualization, methodology, resources, validation, writing – original draft, writing – review and editing, and funding acquisition. Qin-Yi Li: validation, writing – review and editing, and funding acquisition. Koji Takahashi: supervision, resources, validation, writing – review and editing, and funding acquisition. All authors discussed the results and approved the final version of the manuscript.

Conflicts of interest

There are no conflicts to declare.

Data availability

Supplementary information (SI): experimental setup, additional KPFM data, and supporting analysis. See DOI: https://doi.org/10.1039/d6cp00454g.

Acknowledgements

This work was supported by JST PRESTO grant no. JPMJPR23O8; JSPS KAKENHI grant nos. JP22KK0249, JP24K00822, JP24H00293, and JP25K22071; and a Grant-in-Aid for JSPS Research Fellow no. JP25KJ1916.

References

1. H. W. Gibson, J. Am. Chem. Soc., 1975, 97(13), 3832–3833,  DOI:10.1021/ja00846a053.
2. L. S. McCarty and G. M. Whitesides, Angew. Chem., Int. Ed., 2008, 47(12), 2188–2207,  DOI:10.1002/anie.200701812.
3. A. Šutka, K. Mālnieks, L. Lapčinskis, P. Kaufelde, A. Linarts, A. Bērziņa, R. Zābels, V. Jurķāns, I. Gorņevs, J. Blūms and M. Knite, Energy Environ. Sci., 2019, 12, 2417–2421,  10.1039/C9EE01078E.
4. R. K. Pandey, H. Kakehashi, H. Nakanishi and S. Soh, J. Phys. Chem. C, 2018, 122, 16154–16160,  DOI:10.1021/acs.jpcc.8b04357.
5. H. T. Baytekin, A. Z. Patashinski, M. Branicki, B. Baytekin, S. Soh and B. A. Grzybowski, Science, 2011, 333(6040), 308–312,  DOI:10.1126/science.1201512.
6. R. Digilov, Langmuir, 2000, 16, 6719–6723.
7. L. E. Helseth, Langmuir, 2023, 39, 1826–1837.
8. Y. Sun, X. Huang and S. Soh, Angew. Chem., Int. Ed., 2016, 55, 9956–9960.
9. V. Artemov, L. Frank, R. Doronin, P. Stärk, A. Schlaich, A. Andreev, T. Leisner, A. Radenovic and A. Kiselev, J. Phys. Chem. Lett., 2023, 14(20), 4796–4802.
10. D. Choi, H. Lee, D. J. Im, I. S. Kang, G. Lim, D. S. Kim and K. H. Kang, Sci. Rep., 2013, 3, 2037,  DOI:10.1038/srep02037.
11. S. Lin, C. Zhang and T. Shao, Friction, 2025, 13(2), 9440968,  DOI:10.26599/FRICT.2025.9440968.
12. J. Hartmann, M. T. Schür and S. Hardt, Nat. Commun., 2022, 13, 289,  DOI:10.1038/s41467-021-27879-0.
13. N. Sinn, M. T. Schür and S. Hardt, Appl. Phys. Lett., 2019, 114, 213704,  DOI:10.1063/1.5091836.
14. A. Z. Stetten, D. S. Golovko, S. A. L. Weber and H.-J. Butt, Soft Matter, 2019, 15, 8667–8679,  10.1039/c9sm01348b.
15. J. Gon and C.-J. Kim, J. Microelectromech. Syst., 2008, 17(2), 257–264,  DOI:10.1109/JMEMS.2007.912698.
16. X. Li, P. Bista, A. Z. Stetten, H. Bonart, M. T. Schür, S. Hardt, F. Bodziony, H. Marschall, A. Saal, X. Deng, R. Berger, S. A. L. Weber and H.-J. Butt, Nat. Phys., 2022, 18, 713–719,  DOI:10.1038/s41567-022-01563-6.
17. S. Kumari, M. Keswani, S. Singh, M. Beck, E. Liebscher, L. Q. Toan and S. Raghavan, ECS Trans., 2011, 41(5), 93–99,  DOI:10.1149/1.3630831.
18. DISCO Corporation, Engineering R&D Division, Dicer Engineering Department, DISCO Technical Review, 2016, 1–3. Available at: https://www.disco.co.jp/eg/solution/technical_review/doc/TR16-02_Wafer%20ESD%20in%20dicing%20saws%20and%20the%20effect%20of%20the%20countermeasure_20160610.pdf.
19. A. Iqbal, M. S. Hossain and K. H. Bevan, Phys. Chem. Chem. Phys., 2016, 18, 29466–29477.
20. F. Yu, J. Li, Y. Jiang, L. Wang, X. Yang, Y. Yang, X. Li, K. Jiang, W. Lü and X. Sun, Adv. Sci., 2023, 10(30), 2302941,  DOI:10.1002/advs.202302941.
21. S. Fang, J. Li, Y. Xu, C. Shen and W. Guo, Joule, 2022, 6(3), 690–701,  DOI:10.1016/j.joule.2022.02.002.
22. Z. Zhang, X. Li, J. Yin, Y. Xu, W. Fei, M. Xue, Q. Wang, J. Zhou and W. Guo, Nat. Nanotechnol., 2018, 13, 1109–1119.
23. X. Li, G. Feng, Y. Chen, J. Li, J. Yin, W. Deng and W. Guo, Nano Res. Energy, 2024, 3, e9120110,  DOI:10.26599/NRE.2024.9120110.
24. S. Fang, H. Lu, W. Chu and W. Guo, Nano Res. Energy, 2023, 3, e9120108,  DOI:10.26599/NRE.2024.9120108.
25. G. Xue, Y. Xu, T. Ding, J. Li, J. Yin, W. Fei, Y. Cao, J. Yu, L. Yuan, L. Gong, J. Chen, S. Deng, J. Zhou and W. Guo, Nat. Nanotechnol., 2017, 12, 317–321.
26. L. Li, S. Feng, Y. Bai, X. Yang, M. Liu, M. Hao, S. Wang, Y. Wu, F. Sun, Z. Liu and T. Zhang, Nat. Commun., 2022, 13, 1043,  DOI:10.1038/s41467-022-28689-8.
27. B. Shao, Z. Song, X. Chen, Y. Wu, Y. Li, C. Song, F. Yang, T. Song, Y. Wang, S.-T. Lee and B. Sun, ACS Nano, 2021, 15, 7472–7481.
28. C. Hinduja, H.-J. Butt and R. Berger, Soft Matter, 2024, 20, 3349–3358.
29. X. Zhou, Y. Ji, Z. Ni, J. Garcia Lopez, K. Peneva, S. Jiang, N. Knorr, R. Berger, K. Koynov and H.-J. Butt, Adv. Mater., 2025, 37, 2420263,  DOI:10.1002/adma.202420263.
30. S. Das and S. K. Mitra, Phys. Rev. E: Stat., Nonlinear, Soft Matter Phys., 2013, 88, 033021,  DOI:10.1103/PhysRevE.88.033021.
31. A. D. Ratschow, H.-J. Butt, S. Hardt and S. A. L. Weber, Soft Matter, 2025, 21, 1251–1262.
32. A. Dorr and S. Hardt, Phys. Rev. E: Stat., Nonlinear, Soft Matter Phys., 2012, 86, 022601,  DOI:10.1103/PhysRevE.86.022601.
33. A. D. Ratschow, L. S. Bauer, P. Bista, S. A. L. Weber, H.-J. Butt and S. Hardt, Phys. Rev. Lett., 2024, 132, 224002,  DOI:10.1103/PhysRevLett.132.224002.
34. A. E. Kaloyeros, F. A. Jove, J. Goff and B. Arkles, ECS J. Solid State Sci. Technol., 2017, 6, 691–714.
35. S. Lin, L. Xu, A. C. Wang and Z. L. Wang, Nat. Commun., 2020, 11(1), 399,  DOI:10.1038/s41467-019-14278-9.
36. W. Olthuis and P. Bergveld, IEEE Trans. Electr. Insul., 1992, 27(4), 691–697,  DOI:10.1109/14.155784.
37. X. Bai, A. Riet, S. Xu, D. J. Lacks and H. Wang, J. Phys. Chem. C, 2021, 125, 11677–11686,  DOI:10.1021/acs.jpcc.1c02272.

---