# Refining the physical description of charge trapping and detrapping in a transport model for dielectrics using an optimization algorithm

Khaled Hallak, Fulbert Baudoin, Virginie Griseri, Florian Bugarin, Stephane Segonds, Séverine Le Roy, Gilbert Teyssedre

<figure><img src="image_1.png" /></figure>

## To cite this version:

Khaled Hallak, Fulbert Baudoin, Virginie Griseri, Florian Bugarin, Stephane Segonds, et al., Refining the physical description of charge trapping and detrapping in a transport model for dielectrics using an optimization algorithm. *COMPEL: The International Journal for Computation and Mathematics in Electrical and Electronic Engineering*, 2023, 10.1108/COMPEL-04-2023-0143. [hal-04297725](https://hal.science/hal-04297725v1)

HAL Id: hal-04297725

https://hal.science/hal-04297725v1

Submitted on 3 Sep 2024

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_2.png" /></figure>
Cite as: K. Hallak, V. Griseri, F. Baudoin, F. Bugarin, S. Segonds, S. Le Roy, G. Teyssedre, "Refining the physical description of charge trapping and detrapping in a transport model for dielectrics using an optimization algorithm", COMPEL - International J. Computation Mathematics Electrical Electronic Engg 42, 1923-1937, 2023.
DOI: 10.1108/COMPEL-04-2023-0143.

# Refining the physical description of charge trapping and detrapping in a transport model for dielectrics using an optimization algorithm

Khaled Hallak<fnref n="1" />, Fulbert Baudoin<fnref n="1" />, Virginie Griseri<fnref n="1" />, Florian Bugarin<fnref n="2" />, Stephane Segonds<fnref n="2" />, Séverine Le Roy<fnref n="1" /> and Gilbert Teyssedre<fnref n="1" />

## Abstract

**Purpose:** The purpose of this paper is to optimize and improve a Bipolar Charge Transport (BCT) model used to simulate charge dynamics in insulating polymer materials, specifically Low Density Polyethylene (LDPE).

**Design/methodology/approach:** An optimization algorithm is applied to optimize the BCT model by comparing the model outputs with experimental data obtained using two kinds of measurements: space charge distribution using the Pulsed Electro-Acoustic (PEA) method and current measurements in nonstationary conditions.

**Findings:** The study provides an optimal set of parameters that offers a good correlation between model outputs and several experiments conducted under varying applied fields. The study evaluates the quantity of charges remaining inside the dielectric even after 24 hours of short-circuit. Moreover, the effects of increasing the electric field on charge trapping and detrapping rates are addressed.

**Originality:** This work provides a unique set of optimal parameters that best match both current and charge density measurements for a BCT model in LDPE, and demonstrates the use of Trust Region Reflective Algorithm (TRRA) for parameter optimization. The study also attempts to evaluate the equations used to describe charge trapping and detrapping phenomena, providing a deeper understanding of the physics behind the model.

**Research limitations/implications:** This study only examined experiments with different applied electric fields, and thus the obtained parameters may not suit the experimental outputs if the experimental temperature varies. Further improvement may be achieved by introducing additional experiments or another source of measurements.

**Keywords:** dielectric material, space charge, bipolar charge transport model, parameters optimization, Langevin's recombination, trapping and detrapping phenomena

## 1. Introduction

The increased interest in energy transmission through High-Voltage Direct Current (HVDC) necessitate a deep understanding of space charge generation and transport in polyethylene-based materials, commonly employed as high-voltage cable insulation. The presence of a space charge within the dielectric produces electric field distortions and may shorten the life of an electric cable. In addition to the experimental observations currently available for LDPE, numerical models can significantly help to better understand the dynamic behavior of space charge inside this dielectric. For this purpose, a numerical model has been developed [1], aiming to describe the mechanisms of charge generation and transport in Low-Density Polyethylene (LDPE) under DC stress. The model is based on Poisson's equation and the conservation law of charges. The model features charge transport of electronic carriers with the help of a field and temperature-dependent mobility (hopping), trapping and de-trapping phenomena, charge injection at the electrodes, and recombination of carriers. These features are guided by the

<footnote n="1">¹LAPLACE (Laboratoire Plasma et Conversion d'Energie), UPS, INPT, CNRS, University of Toulouse, Toulouse, France</footnote>

<footnote n="2">²ICA (Institut Clément Ader), CNRS/INSA/ISAE/Mines Albi/UPS, University of Toulouse, Toulouse, France</footnote>
experimental behavior observed in this kind of material. Transport and trapping are considered at different energy levels; shallow traps for transport and deep traps for trapping. The proposed model, which is briefly described in the following section, requires a set of parameters, including injection barrier, trap depth for the hopping mechanism, deep trapping and de-trapping rates, and recombination coefficients. Most of these parameters cannot be estimated using experimental measurements. In other words, it is challenging to identify the most suitable parameter values that best fit the experimental measurements. An optimization approach known as the TRRA was introduced to estimate the values of the unknown parameters [2]. This technique is used in our study to minimize the difference between experimental and simulated data in order to achieve the optimal set of parameters that best suit experimental measurements. The net charge density is determined by the PEA [3] and the current density acquired from charging current measurements. On the other side, the established BCT model for LDPE under DC stress is used to estimate the simulated charge and current densities.

A previous study [4] showed that incorporating several experiments into the optimization process improves the method's convergence rate toward the global set of parameters. In our case, experiments are considered with varying electric fields while the temperature was fixed for all experiments. Furthermore, TRRA will be employed to understand the impact of the applied electric field on the trapping and detrapping rates.

This paper is composed of five sections. First, the used physical model is described and then the experimental approach on which the model is based is presented. The third part introduces the principle of TRRA and its implementation. In the main section, parameters will be computed using several experimental measurements for current and charge densities. Finally, the results will be analyzed and the model limitations and possible improvements will be discussed.

## 2. Bipolar charge transport model

### 2.1 Bases of the model

Consider a flat sample of insulating material of thickness D sandwiched between two electrodes. As shown in Fig 1, the injected positive and negative charge carriers (holes and electrons) flow towards the opposite electrode due to the influence of an applied electric field. The injected charge carriers can be captured in trapping states due to the presence of localized energy states in the band-gap of the dielectric material. The trapped charges can also escape from these traps, forming mobile ones. Moreover, charge carriers with opposite polarity may recombine.

<figure><img src="image_3.png" /><figcaption>**Fig 1.** Schematic representation of the charge transport mechanisms in LDPE. $w_e$ and $w_h$ are injection barriers for electrons and holes; $S_0, S_1, S_2$, and $S_3$ are the recombination coefficients, $B_e$ and $B_h$ are electron and hole trapping coefficients, $D_e$ and $D_h$ are electrons and holes de-trapping coefficients, $w_{mo_e}$ and $w_{mo_h}$ are depths for a single shallow electron and hole trapping level.</figcaption></figure>

The model is 1D, function of the depth in the dielectric. Space and time-dependent equations are considered to describe the charging behavior in this model, known as the Poisson's equation, the continuity equation, and the transport equation neglecting diffusion [5], and can be written as follows:

Transport Equation:
$$j_{e,h}(x,t) = \mu_{e,h}(x,t) \cdot E(x,t) \cdot n_{e,h}(x,t) \qquad (1)$$

Poisson's Equation:

$$\frac{\partial E(x, t)}{\partial x} = \frac{n(x, t)}{\varepsilon_0 \varepsilon_r} = \frac{n_{h\mu}(x, t) + n_{ht}(x, t) - n_{e\mu}(x, t) - n_{et}(x, t)}{\varepsilon_0 \varepsilon_r} \qquad (2)$$

Continuity Equation, of the form:

$$\frac{\partial n_{e,h}(x, t)}{\partial t} + \frac{\partial j_{e,h}(x, t)}{\partial x} = s_i(x, t) \qquad (3)$$

where the indexes e and h refer to the type of charge, electrons or holes, respectively; $j_{e,h}$ are the fluxes of carriers; $n_{e,h}$ is the charge density of electrons and holes (in C/m³); respectively; E is the electric field, n is the net charge density; $n_{h\mu}$ and $n_{ht}$ are the mobile and trapped holes density respectively; $n_{e\mu}$ and $n_{et}$ are the mobile and trapped electrons density respectively. $\varepsilon_0 = 8.85 \times 10^{-12}$ F/m is the vacuum permittivity and $\varepsilon_r = 2.3$ is the relative dielectric permittivity for LDPE. The term $s_i$ is the source term (e.g., $s_{e\mu}$, $s_{et}$, $s_{h\mu}$ and $s_{ht}$, respectively, for mobile electrons, trapped electrons, mobile holes, and trapped holes), which describes the microscopic processes inside the dielectric, other than transport phenomena, such as, recombination, trapping, and de-trapping. An example of the mobile electrons source term $s_{e\mu}$ is presented below, where most of variables are in terms of time t and position x and are omitted:

$$\frac{\partial n_{e\mu}}{\partial t} + \frac{\partial j_{e\mu}}{\partial x} = s_{e\mu} = -R_1 \frac{\mu_e}{\varepsilon_0 \varepsilon_r} n_{ht} n_{e\mu} - R_3 \frac{(\mu_e + \mu_h)}{\varepsilon_0 \varepsilon_r} n_{e\mu} n_{h\mu} - B_e n_{e\mu} \left(1 - \frac{n_{et}}{N_{0et}}\right) + D_e n_{et} \qquad (4)$$

Here, the recombination coefficients (of the form R.μ/ε) involve the reduction factors $R_i$ and the Langevin's recombination rate which is mobility dependent, $B_{e,h}$ are the trapping coefficients for electrons and holes, $D_{e,h}$ are the de-trapping rates for electrons and holes, as shown in Fig 1. $N_{0et}$ and $N_{0ht}$ are the trap densities for electrons and holes, respectively.

The recombination reduction pre-factors $R_i$ are introduced in order to weight this physical phenomenon in relation to the others while respecting Langevin form.

The deep trapping is described through a single trapping level of each kind of carrier by a trapping coefficient $B_{e,h}$, which expresses the trapping rate of charged species per unit of time (1/s). The intensity of deep trapping is defined as:

$$K_e = B_e \cdot n_{e\mu} \cdot \left(1 - \frac{n_{et}}{N_{0et}}\right) \qquad (5)$$

$$K_h = B_h \cdot n_{h\mu} \cdot \left(1 - \frac{n_{ht}}{N_{0ht}}\right) \qquad (6)$$

The detrapping process is described by the de-trapping coefficients $D_{e,h}$, which reflect the probability of de-trapping charged species per unit of time (1/s), as:

$$D_{tre} = D_e \cdot n_{et} \qquad (7)$$

$$D_{trh} = D_h \cdot n_{ht} \qquad (8)$$

In this model, the charge generation in the additive-free LDPE is considered only from the injection at the electrodes following the modified Schottky law. The injected current densities can be expressed as:
$$j_e(0, t) = AT^2 \exp\left(-\frac{qw_e}{k_B T}\right) \left[ \exp\left(\frac{q}{k_B T} \sqrt{\frac{qE(0, t)}{4\pi\varepsilon_0\varepsilon_r}}\right) - 1 \right] \quad (9)$$

$$j_h(D, t) = AT^2 \exp\left(-\frac{qw_h}{k_B T}\right) \left[ \exp\left(\frac{q}{k_B T} \sqrt{\frac{qE(D, t)}{4\pi\varepsilon_0\varepsilon_r}}\right) - 1 \right] \quad (10)$$

where $j_e(0, t)$ and $j_h(D, t)$ are the current density of electrons and holes at cathode and anode respectively in A/m². A = 1.2× 10⁶ A/m²K² is the Richardson constant; $w_e$ and $w_h$ are the injection barriers for electrons and holes, respectively, $k_B = 8.617 \times 10^{-5}$ eV/K stands for Boltzmann constant, and T is the temperature in K.

## 2.2 Field dependent mobility

Charge transport is described by a hopping mechanism in which carriers move from site to site by tunneling through a potential barrier. This mobility accounts for the possible trapping and de-trapping in shallow traps with a short-term residence. To introduce the hopping mechanism, mobilities $\mu_{e,h}$ in Equation (1) are function of the hopping barrier representing the shallow traps, $w_{mo_{e,h}}$:

$$\mu_{e,h}(x, t) = \frac{2va_{e,h}}{E(x, t)} \exp\left(-\frac{qw_{mo_{e,h}}}{k_B T}\right) \sinh\left(\frac{qE(x, t)a_{e,h}}{2k_B T}\right) \quad (11)$$

where $w_{mo_{e,h}}$ is the shallow trap depth. $a_{e,h}$ is the average spacing between shallow traps, for electrons and holes.

Some previous studies used the same value of inter-trap distance for electrons and holes, despite their different transport environments in polyethylene. For instance, hole transport occurs via electron vacancies in or closely associated with the valence band, and is confined to polymer chains, requiring inter-chain transfer via tunneling. This makes hole conduction more sensitive to the morphology of polyethylene [6, 7]. This description of hole transport shows the complexity of describing accurately the physical phenomena and hence of getting precise constants or parameters related to the polymer physics. The inter-trap distances are an example of this kind of complexity, with reported values ranging from a = 2.8 to 10 nm [6, 8-12]. This complexity addressed in this study by optimizing the inter-trap distance parameter, considering the asymmetry between electrons and holes.

## 2.3 Field depending deep trapping and detrapping rates

In the literature, most of the presented BCT models consider the trapping and detrapping rates independent of the electric field [1-5]. Nevertheless, some works have addressed experimentally and theoretically the relationship between trapping and detrapping phenomena and the electric field [12, 14]. As regards the temperature dependency of the detrapping rate, it is described by a thermally activated process.

In order to see if there is an impact of the applied electric field on the trapping rates, the trapping and detrapping coefficients have been optimized for each simulated applied electric field, and are noted $B_{e,i}$ and $B_{h,i}$ for electrons and holes, respectively, for i =20, 40 and 60 kV/mm. In the same way, the detrapping rates for electrons and holes, have also been optimized for each applied electric field.

## 2.4 Post treatment of simulated space charge profiles

Model outputs cannot be directly compared with the experimental data, especially the charge profiles at the electrodes [3]. Induced charges, constituted by image and capacitive charges are theoretically surface charges represented by a plan of charges.

However, because of the limited spatial resolution of any space charge measurement method, the experimental results display a Gaussian-shaped charge spread at the dielectric interfaces. For this reason, a numerical filter (Gaussian) has been applied to the data in order to lose a part of the resolution reached by the simulated profiles as in **Supplementary material S1**. The Gaussian coefficient ($\sigma$) representing the standard deviation of the gaussian function has also been optimized.

To account for the influence charges (capacitive and image charges), surface charge densities have been added at the electrodes according to the following equations:
Capacitive charge:

$$n_{\text{cap}} = \pm \varepsilon_0 \cdot \varepsilon_r \cdot \frac{V}{D} \qquad (12)$$

Image charge at each electrode:

$$n_{\text{img}}(0) = -\frac{1}{D} \cdot \sum_{i=1}^{i=d} n(x_i) \cdot (D - x_i) \cdot \Delta x_i \qquad (13) \\ n_{\text{img}}(D) = -\frac{1}{D} \cdot \sum_{i=1}^{i=d} n(x_i) \cdot x_i \cdot \Delta x_i$$

where $\Delta x_i$ are the size of cells, d is the number of cells.

Moreover, it is not obvious that the material is free of charges prior to any experiment. To account for this, an additional parameter, $\alpha$, is proposed to represent an initial charge density for electrons and holes, prior any experiment, with a zero net charge density everywhere inside the dielectric:

$$\alpha = n_{h\mu} = n_{ht} = n_{e\mu} = n_{et} \qquad (14)$$

In total, 23 different model parameters are included in the optimization process and displayed in **Table I**. Among them, 12 parameters account for trapping-detrapping processes and their field dependency. The lower and upper bound constraints for the parameters are specified in order to have only physically acceptable results.

# 3. Experimental approach

## 3.1 Sample preparation and experimental protocol

LDPE pellets with a melt flow rate of 0.8 g/10 min and no antioxidants were provided by Borealis. Films of approximately 150 µm thickness were obtained by press-molding and metallized with a 50 nm gold layer on both sides.

Two different experiments were carried out using the same LDPE sample. Experiment 1 aims to measure the charge density in LDPE at three levels of fields of 60, 40, and 20 kV/mm. The same protocol is used for Experiment 2 for measuring the current density. Both experiments consist of polarization for 3 hours and a short-circuiting period of 24 hours each. Space charge measurements were performed using the PEA method. The current density was obtained by the electrode current density measurements. For further information about the protocol and instruments used, see **Supplementary Material S1**.

Since the study is based on a dual approach (experiment/simulation), all the measurements were performed using the same sample and exactly on the same location. This was done to ensure that all experiments were performed using the same thickness and same material properties. According to the protocol presented in **Supplementary Material S1**, 6 experiments are performed; thus, 6 additional parameters will be added to the optimization procedure to estimate the initial charge density before each experiment. The parameters $\alpha_{\text{Exp1}}^{20}$, $\alpha_{\text{Exp1}}^{40}$, $\alpha_{\text{Exp1}}^{60}$ represents the initial charge density before starting space charge measurements for fields 20, 40 and 60 kV/mm, respectively. Besides, $\alpha_{\text{Exp2}}^{20}$, $\alpha_{\text{Exp2}}^{40}$ and $\alpha_{\text{Exp2}}^{60}$ represent the initial charge density before starting current measurements for the same fields, respectively.
<table><caption>**Table I.** Symbols, units, lower, and upper bounds of the unknown parameters.</caption><thead><tr><th>Parameters</th><th>Barrier height of injection</th><th>Trapping depth</th><th>Trapping coefficient</th><th>Detrapping coefficient</th><th>Reduction pre-factor</th><th>Inter-Trap distance</th><th>Gaussian coefficient</th><th>Initial charges</th></tr></thead><tbody><tr><td>Symbol</td><td>w<sub>e</sub> w<sub>h</sub></td><td>w<sub>mo<sub>e</sub></sub> w<sub>mo<sub>h</sub></sub></td><td>B<sub>e</sub> B<sub>h</sub></td><td>D<sub>e</sub> D<sub>h</sub></td><td>R<sub>1</sub> R<sub>2</sub> R<sub>3</sub></td><td>a<sub>e</sub> a<sub>h</sub></td><td>σ</td><td>α</td></tr><tr><td>Units</td><td>eV</td><td>eV</td><td>s<sup>-1</sup></td><td>s<sup>-1</sup></td><td>Unitless</td><td>nm</td><td>μm</td><td>C/m<sup>3</sup></td></tr><tr><td>Lower Bound (lb)</td><td>1</td><td>0.59</td><td>0</td><td>0</td><td>0</td><td>0</td><td>0</td><td>0</td></tr><tr><td>Upper Bound (ub)</td><td>1.3</td><td>0.72</td><td>10</td><td>10</td><td>1</td><td>20</td><td>10</td><td>20</td></tr></tbody></table>

## 3.2 Experimental results

The absolute current density with respect to time for the three applied fields is presented in Fig 2 in black. It is obvious that the current density increases with respect to the applied electric field. Besides, even after three hours of polarization, the experimental current density seems to be decreasing with respect to polarization time and does not attain a stationary state. Furthermore, all the experimental current density curves have nearly the same slope.

<figure><img src="image_4.png" /><figcaption>**Fig 2.** Current density over time for varied electric fields (20 kV/mm, 40 kV/mm, 60 kV/mm): Experimental (Black), simulated using P1 (Green), and simulated using P3 (Red) Curves. Top: 60 kV/mm, Middle: 40 kV/mm, Bottom: 20 kV/mm.</figcaption></figure>

The obtained space charge patterns are presented in Fig 3 (upper row). The amount of injected charge increases as the applied voltage is increased. Positive (red) and negative (blue) charges are detected immediately in the bulk after the highest fields (40 and 60 kV/mm) are applied. These charges are injected at the anode and cathode, respectively, and move through the dielectric towards the opposite electrode. A low amount of charges are observed inside the bulk with a field of 20 kV/mm. Negative charges appear dominant at short time since massively injected in a short period of time. Moreover, the mobility of electrons seems to be higher than that of holes at any electric field. After 1800 s, the positive charges dominate inside the bulk of the insulator, compared to negative charges. Eventually, after 2000 s the behavior of the insulator seems to be reaching an equilibrium state in terms of the space charge distribution, and amount of injected charge.
<figure><img src="image_5.png" /><figcaption>Fig 3. Space charge profiles as a function of the time and sample thickness at room temperature, for fields of 20, 40 and 60 kV/mm. Upper row: experimental charge density. Lower row: simulated charge density using P3.</figcaption></figure>

# 4. Parameters optimization using experimental measurements

## 4.1 Optimization using current density measurements

In this section, TRRA is employed using only the current density measurements under three different applied fields (20, 40 and 60 kV/mm). For further information about TRRA, a detailed explanation is presented in the **Supplementary Material S2**. The cost function to be minimized is defined as:

$$C_1(P) = C_1(P)_{20} + C_1(P)_{40} + C_1(P)_{60} \qquad (15)$$

where

$$C_1(P)_E = \frac{\|j_{\text{exp}}(t) - j_{\text{sim}}(t, P)\|_E^2}{\|j_{\text{exp}}(t)\|_E^2} \qquad (16)$$
$$\text{s.t.} \quad lb \le P \le ub$$

where:

* P are the unknown parameters presented in **Table I**,

* $j_{\text{exp}}(t)$ are the experimental current density (A/m²),

* $j_{\text{sim}}(t, P)$ are the simulated current density computed using the BCT model in terms of P,

* lb and ub are lower and upper bounds of the parameters respectively.

Once $C_1(P)$ is minimized, TRRA produces the optimal set of parameters $P_1$. The comparison between the experimental and the simulated current density results is displayed in **Fig 2**. It seems that $P_1$ produce a good match between experimental and simulated current density for the three fields under consideration. This match could be illustrated numerically by computing the relative error between experimental and simulated results using $P_1$, as shown in **SupTable 1 (S4)**.

The same set of parameters $P_1$, has been used to reproduce the experimental space charge data. An example of simulated results given for an applied electric field of 40 kV/mm is presented in the **Fig 4** significant difference could be noticed between experimental and simulated outcomes. Thus, $P_1$ is found as an optimal solution that only fits the experimental currents and contradicts the space charge experimental data. It is obvious that current alone
cannot lead to identify all the set of parameters. So, the optimal solution found here is probably an example among many solutions.

<figure><img src="image_6.png" /><figcaption>Fig 4. Simulated charge density pattern in C/m³ using parameters set P₁, for a field of 40 kV/mm</figcaption></figure>

It is worth noting that charge density data may offer more comprehensive insights compared to current measurements, as space charge measurements consider time and position, whereas current measurements solely account for time. For this reason, TRRA is employed using the charge density measurements considering the same varying electric fields. Thus, the cost function is defined as:

$$C_2(P) = C_2(P)_{20} + C_2(P)_{40} + C_2(P)_{60} \quad (17)$$

where

$$C_2(P)_E = \frac{\|n_{\text{exp}}(x, t) - n_{\text{sim}}(x, t, P)\|_E^2}{\|n_{\text{exp}}(x, t)\|_E^2} \quad (18)$$
$$\text{s.t.} \quad \text{lb} \le P \le \text{ub}$$

where

* $n_{\text{exp}}(x, t)$ is the experimental net charge density data (C/m³)

* $n_{\text{sim}}(x, t, P)$ is the simulated net charge density.

The results obtained from the optimization of $C_2(P)$ indicate that relying solely on space charge measurements is inadequate for accurately predicting the material's conduction behavior. For further details regarding the outcomes of the optimization using charge density only, please refer to the **Supplementary Material S3**.

## 4.2 Optimization using both current and charge density measurements

As mentioned previously, our main target is to find a set of parameters able to fit both current and charge experimental data simultaneously. For this reason, TRRA is used in this section to optimize a multi-objective function by combining both charge and current densities in $C_3(P)$. The cost function to be optimized is defined as:

$$C_3(P) = C_1(P) + C_2(P) \quad (19)$$

The optimized set of parameters $P_3$ produced by minimizing $C_3(P)$ is displayed in **Table II**. All the optimal parameters will be discussed in the following section. The simulated results obtained using $P_3$ are represented in **Fig 2 and Fig 3**.

**Fig 2** shows that $P_3$ can produce a good match between experimental and simulated current density, for all considered fields. It is worth noting that it is extremely difficult to find a single set of parameters that can match numerous current measurements simultaneously, because the experimental current measurements are not stable,
and the curve shape may change if the same experiment is repeated. This is one of the factors that influenced the correlation between experimental and simulated current.

<table><caption>**Table II.** Optimized parameters P<sub>3</sub> obtained by minimizing C<sub>3</sub> (P)</caption><thead><tr><th>Parameters</th><th>20 kV/mm</th><th>40 kV/mm</th><th>60 kV/mm</th></tr></thead><tbody><tr><td>w<sub>e</sub></td><td></td><td>1.184</td><td></td></tr><tr><td>w<sub>h</sub></td><td></td><td>1.200</td><td></td></tr><tr><td>w<sub>mo<sub>e</sub></sub></td><td></td><td>0.563</td><td></td></tr><tr><td>w<sub>mo<sub>h</sub></sub></td><td></td><td>0.602</td><td></td></tr><tr><td>B<sub>e</sub></td><td>0.6770</td><td>0.6772</td><td>0.6802</td></tr><tr><td>B<sub>h</sub></td><td>0.0125</td><td>0.0301</td><td>0.0963</td></tr><tr><td>D<sub>e</sub></td><td>0.0030</td><td>0.0034</td><td>0.0035</td></tr><tr><td>D<sub>h</sub></td><td>0.0003</td><td>0.0012</td><td>0.0092</td></tr><tr><td>R<sub>1</sub></td><td></td><td>1</td><td></td></tr><tr><td>R<sub>2</sub></td><td></td><td>1</td><td></td></tr><tr><td>R<sub>3</sub></td><td></td><td>1</td><td></td></tr><tr><td>a<sub>e</sub></td><td></td><td>2.404</td><td></td></tr><tr><td>a<sub>h</sub></td><td></td><td>1.103</td><td></td></tr><tr><td>σ</td><td></td><td>1.377</td><td></td></tr><tr><td>α<sub>Exp<sub>1</sub></sub></td><td>0.026</td><td>0.021</td><td>0.009</td></tr><tr><td>α<sub>Exp<sub>2</sub></sub></td><td>14.7210</td><td>5.4923</td><td>1.2520</td></tr></tbody></table>

**Fig 3** shows the simulated space charge cartography for all electric fields. The results indicate that P<sub>3</sub> is capable of producing a satisfactory match between simulated and experimental charge density for all electric fields. For 20 and 40 kV/mm, the simulated results did not perfectly fit the experimental measurements for the first one hour of polarization. However, after one hour, the simulated measurements were consistent with experimental measurements for all electric fields. The model reproduces the essential phenomena observed in the experiment, including the dominance of negative charges initially and positive charges until the end of the experiment. Using charge and current density measurements individually produced better results than combining them into a single cost function. This may be due to model limitations or experimental measurement accuracy. Nevertheless, the relative error displayed in the **Supplementary Material S4** reveals also that P<sub>3</sub> provided a good correlation between experimental and simulated current and charge densities.

# 5. Discussion

Section 5 showed that the optimized parameters P<sub>3</sub> are more reliable than P<sub>1</sub> and P<sub>2</sub> since P<sub>3</sub> can reproduce both current and charge densities. Thus, the discussion part will be focused on analyzing and interpreting the results provided only by P<sub>3</sub> (Table II).

## 5.1 Trapping and detrapping rates

In order to identify the relation between trapping and detrapping phenomena and the applied electric field, the trapping and detrapping coefficients are considered as dependent on the electric field with a fixed temperature for all experiments. **Fig 5** represents the optimal values of trapping and detrapping rates at fields ranging from 20 kV/mm to 60 kV/mm found by TRRA. Based on **Table II**, both trapping and detrapping rates are found to be increasing with the increase of the electric field. Moreover, the rise of the trapping and detrapping rates of electrons with respect to the electric field is small compared to the increase of the coefficients for holes. Besides, the trapping rates of electrons and holes are greater than the detrapping rates at any electric field. Also, the trapping rate of electrons is much higher than that of holes.
In previous studies, the temperature dependence of the detrapping rates $D_{e,h}$ is considered of the form:

$$D_{e,h} = v \cdot \exp\left(\frac{-w_{tr_{e,h}}}{k_B T}\right) \qquad (20)$$

where $v$ is the frequency of the phonons, $w_{tr_{e,h}}$ are the barrier heights to escape from deep traps. Therefore, the energy level of deep traps could be evaluated from the optimal values of the detrapping coefficient using Equation (20). The obtained values are $w_{tr_e} = 0.846$ eV and $w_{tr_h} = 0.903$ eV at room temperature. These values are found in same range as in the literature [1, 5, 13, 14].

<figure><img src="image_7.png" /><figcaption>**Fig 5.** Trapping and detrapping rates in s⁻¹ for electrons and holes as a function of applied field at 20°C</figcaption></figure>

## 5.2 Mobility of carriers

The mobility of electrons and holes obtained by the BCT model using $P_3$ can be evaluated by substituting $w_{moe}$, $w_{moh}$, $a_e$ and $a_h$ in the mobility equation (11), considering the applied electric field and temperature. Some values are given in **Supplementary Material S5**. The mobility of electrons, is found greater than the mobility of holes for any applied field, by over one decade (5 10⁻¹⁴ vs. 1.5 10⁻¹³ at 20°C and 40 kV/mm). These results on mobility values are consistent with literature reports [1, 3, 5, 12], where both vary in the same range [10⁻¹⁴, 10⁻¹² m²·V⁻¹·s⁻¹].

The estimated values of the average inter-trap distance between two shallow traps are of the order of nm, consistently with the values reported in the literature [7-13]. The intertrap distance for electrons, $a_e = 2.5$ nm is greater than that for holes, $a_h = 1.85$ nm. A consequence of that is a simulated mobility of electrons more dependent on the electric field compared to the mobility of holes.

## 5.3 Recombination coefficients

According to **Table II**, the optimal values of the recombination reduction pre-factors ($R_1 = R_2 = R_3 = 1$), prove that either Langevin's recombination form fully describes the recombination process with no reduction taking place, or that recombination is not a dominant process in the present protocol making optimization challenging for this parameter.

The space- and time-averaged recombination between different kinds of carriers with respect to the applied electric field is presented in **Fig 6**. The evaluation of space- and time-averaged recombination is presented in **Supplementary Material S5**. The results show that recombination between trapped electrons and mobile holes has the lowest rate compared to other recombination processes. Furthermore, the recombination between mobile electrons and mobile holes is equivalent to that between mobile electrons and trapped holes. In addition to the discussed field dependencies of the recombination, their temperature dependency also exists and follows the consideration of the mobility equation.

Moreover, the obtained recombination results are comparable with the literature [1, 14], ranging between 5 × 10⁻³ and 5 × 10⁻² for fields ranging between 20 and 60 kV/mm.
<figure><img src="image_8.png" /><figcaption>**Fig 6:** Space- and time-averaged recombination coefficients between different types of carriers.</figcaption></figure>

## 5.4 Initial charge conditions

The optimal values of the initial charge density parameters are displayed in **Table II**. It is worth noting that Experiment 1 started with the highest applied electric field (60 kV/mm). In the simulations related to the space charge measurement protocol, the initial charge density inside the bulk is found to be increasing after performing a new experiment on the same LDPE sample (i.e., $\alpha_{Exp1}^{60} < \alpha_{Exp1}^{40} < \alpha_{Exp1}^{20}$). The same holds for the current measurements (i.e., Experiment 2, **SupFig 4**), where $\alpha_{Exp2}^{60} < \alpha_{Exp2}^{40} < \alpha_{Exp2}^{20}$. The quantity of charges obtained during the current measurements are very high compared to the charges found initially during the space charge measurements (i.e., $\alpha_{Exp2}^{20} \gg \alpha_{Exp1}^{20}$).

This high charge density is however necessary to obtain always decreasing current densities as a function of time for each applied electric field. The fact that these charges are necessary to obtain a consistent comparison between experimental and simulated data may be due to the neglection of the polarization process [15], which is a function of the applied field, and is observable short time after voltage application in the current measurements. Further research is needed to address this aspect, but the main conclusion regarding the initial charge density is that it is a necessary parameter for obtaining consistent simulated current densities, especially at short times.

## 5.5 Barrier height for injection

The optimal parameters produced by TRRA displayed in **Table II**, show that the barrier height of electrons injection is lower than that of holes. In other words, the injection of electrons is much easier than holes injection. Besides, the influence of the barrier height could be experimentally observed at short time. Indeed, directly after applying the voltage, the injection process is the dominant phenomenon. Based on the experimental measurements, it is clear that the injection of electrons is found to be stronger than that of holes at the beginning of each experiment for any electric field. These results are consistent with the simulated results obtained by the model using $P_3$.

# 6. Conclusion

The BCT model was improved by using TRRA to identify optimal parameters that match both current and charge density measurements at varying electric fields. The optimization process involved three cases: one based solely on current density data, another based solely on charge density data, and a third case that combined both datasets. The first and second cases resulted in optimal parameters that effectively matched each dataset individually. However, when the third cost function was applied, the findings showed a satisfactory agreement between experimental and simulated results for both current and charge density data.

The optimal values of the recombination reduction pre-factors prove that the recombination process is either entirely characterized by Langevin's recombination form with no reduction taking place, or that it is not possible to optimize the recombination processes with the chosen experiments, as it is not a dominant process. The optimal values of the mobility indicate that the mobility of electrons is higher than that of holes, which is consistent with the literature.

The optimal values of trapping and detrapping rates increase with the electric field in case of electrons but seem constant in case of holes. Besides, the trapping rate is found to be greater than the detrapping rate for both carriers
at any field. Moreover, the trapping rate of electrons was found to be much greater than that of holes. Whereas the detrapping rates of electrons and holes are found to be approximately equal. Finally, the initial charges inside the bulk were observed to be rising after performing a new experiment on the same LDPE sample, even though the sample was discharged for 24 hours.

The obtained optimal parameters properly suited both current and charge densities at varying electric fields. As noted earlier, the major objective of this work is to find a unique set of parameters that match both current and charge using any experimental conditions. In this study, only experiments with different applied electric fields were examined. However, if the polarity of the electric field is reversed, the parameters may not suit the experimental outputs anymore. Thus, our technique might be further improved by introducing additional experiments with alternative experimental conditions (i.e., stress protocol with polarity reversal or varying temperatures) or adding another source of measurements (i.e., electroluminescence as a signature of recombination processes).

## References

[1] Le Roy, S., Baudoin, F., Griseri, V., Laurent, C., & Teyssedre, G. (2004). Space charge dynamics in LDPE under temperature and electric field stress. Journal of Physics D: Applied Physics, 37, 298-305.

[2] Hallak, K., Baudoin, F., Griseri, V., Bugarin, F., & Segonds, S. (2020). Numerical Optimization Applying Trust Region Algorithm to Optimize Parameters Related to Charge Transport Model in LDPE. IEEE Transactions on Dielectrics and Electrical Insulation, 27(6), 2048-2055. https://doi.org/10.1109/TDEI.2020.008876.

[3] Thomas, C., Teyssedre, G., & Laurent, C. (2011). Space charge accumulation and decay in polyethylene: An experimental study. Journal of Physics D: Applied Physics, 44, 015401. https://doi.org/10.1088/0022-3727/44/1/015401.

[4] Hallak, K., Baudoin, F., Griseri, V., Bugarin, F., & Segonds, S. (2021). A New Approach for Optimizing a Bipolar Charge Transport Model for Dielectric Materials: Theoretical Framework. IEEE Transactions on Dielectrics and Electrical Insulation, 28(3), 872-879. https://doi.org/10.1109/TDEI.2021.009450.

[5] Hoang, A. T., Serdyuk, Y., & Gubanski, S. M. (2016). Space charge in polymeric insulation materials. Polymers, 8, 1-16. https://doi.org/10.3390/polym8010005.

[6] Lewis, T. J., & Llewellyn, J. P. (2013). Measurement of space charge dynamics in crosslinked polyethylene using PEA. Journal of Applied Physics, 113, 223705. https://doi.org/10.1063/1.4811653.

[7] Chen, G., & Zhao, J. (2011). Space charge dynamics in LDPE under DC voltage. Journal of Applied Physics, 109, 212001. https://doi.org/10.1063/1.3596609.

[8] Chen, G., Li, G., & Li, S. (2016). Space charge and electrical conduction in cross-linked polyethylene. Applied Physics Letters, 109, 062901. https://doi.org/10.1063/1.496

[9] Dennison, J. R. & Brunson, J. (2008). Space charge dynamics in crosslinked polyethylene. IEEE Transactions on Plasma Science, 36, 2246-2252. DOI: 10.1109/TPS.2008.2004886

[10] Ziari, Z., Sahli, S., & Bellel, A. (2010). Space charge in LDPE under DC voltage. Materials Journal of Condensed Matter, 12, 218-222. DOI: 10.1088/0953-8984/12/3/307

[11] Nath, R., Kaura, T., & Perlman, M. M. (1990). Space charge measurements in polyethylene. IEEE Transactions on Dielectrics and Electrical Insulation, 25, 419-425. DOI: 10.1109/14.55654

[12] Kim, M., Kim, S. H., & Lee, S. H. (2020). Surface charge accumulation and space charge dynamics in LDPE under DC voltage. IEEE Access, 8, 200051-200062. https://doi.org/10.1109/ACCESS.2020.3035204.

[13] Boufayed, F., Teyssèdre, G., & Laurent, C. (2006). Space charge measurements in low-density polyethylene and its nanocomposites. Journal of Applied Physics, 100, 104-105. DOI: 10.1063/1.2210028

[14] Doedens, E., Jarvid, M., Guffond, R., & Serdyuk, Y. V. (2020). Space charge dynamics in cable insulation materials. Energies, 13, 812-815. DOI: 10.3390/en13040812

[15] Adamec, V., & Calderwood, J. H. (1983). Electrical breakdown of insulating materials: a review. Journal of Physics D: Applied Physics, 16, 203-211. DOI: 10.1088/0022-3727/16/2/003
# Supplementary Material

## S1. Protocol and experimental details

### S1.1. Current measurement

The cell used for current measurements consists of a flat sample of insulating material of thickness D that is sandwiched between two electrodes with 1 cm diameter that have been optically polished (see **SupFig 1**). A positive or negative DC voltage $V_0$ is applied to the upper electrode at zero time, while the other electrode is grounded. Keithley 6512 ammeter is used to measure the current $I(t)$ in 2s steps during the application of the electrical stress.

<figure><img src="image_9.png" /><figcaption>**SupFig 1.** Current measurement cell (Source: Created by authors)</figcaption></figure>

It is possible to estimate the total current density $J_c$ by:

$$J_c(t) = \frac{I(t)}{S} \qquad (9)$$

where S is the area of the surface electrode.

### S1.2. Pulsed Electro-Acoustic method (PEA): Charge density measurement

The schematic diagram of the PEA device is shown in **SupFig 2**. It is made up of the following elements:

* A pulse generator that provides an electrical signal, for the present experiments, with a frequency of 1 kHz, 8 ns width and amplitude 250 V.

* A high voltage DC source used to polarize the sample.

* A digital oscilloscope and a computer for data acquisition and digital processing, respectively.

The upper electrode of the PEA cell is made of semi-conducting polymer (carbon black doped polyethylene), with a diameter 10 mm and the lower electrode is made of aluminum. Below the lower electrode, there is a 9 µm thick PVDF (polyvinylidene fluoride) piezoelectric sensor. The PEA cell we use has a spatial resolution of about 10 µm and a sensitivity of 0.1 C/m³.

<figure><img src="image_10.png" /><figcaption>**SupFig 2.** The schematic diagram of the PEA device (Source: Created by authors)</figcaption></figure>
### S1.3. Image and capacitive charges

Induced charges, constituted by image and capacitive charges are theoretically surface charges represented by a plan of charge as in the left of **SupFig 3**.

However, because of the limited spatial resolution of any space charge measurement method, the experimental results display a Gaussian-shaped charge spread at the dielectric interfaces. For this reason, a numerical filter (Gaussian) has been applied to the simulation data in order to lose a part of the resolution reached by the simulated profiles as shown to the right of the **SupFig 3**.

<figure><img src="image_11.png" /><figcaption>**SupFig 3.** Profiles featuring influence charges on the electrodes. (left) without Gaussian filter, (right) with Gaussian filter
(Source: Created by authors)</figcaption></figure>

### S1.4. Experimental protocol and measurements

Two different experiments were carried out using the same LDPE sample. Experiment 1 aims to measure the charge density in LDPE at three levels of fields of 60, 40, and 20 kV/mm. The same protocol is used for Experiment 2 for measuring the current density. Both experiments consist of polarization for 3 hours and a short-circuiting period of 24 hours each, as shown in **SupFig 4**.

<figure><img src="image_12.png" /><figcaption>**SupFig 4.** Experimental protocol of polarization/short-circuit (150 µm thick LDPE film) (Source: Created by authors)</figcaption></figure>

## S2. Trust region reflective algorithm

TRRA is one of the most important numerical optimization methods in solving Nonlinear Least Squares (NLS) problems. The TRRA is classified as one of the derivative-based algorithms where it uses both gradient and Hessian information to find a good search direction. It is well known for its high efficiency and accuracy. Here it is used to optimize several BCT model parameters simultaneously. The methodology is as follows: assume that we need to find the value of a set of parameters vector of real variables $P^*$ such that it minimizes the cost function $C(P)$, where $C(P)$ describes the difference between experimental and simulated data. The problem to be solved is defined as:
$$ \min_{P \in \mathbb{R}^n} C(P) = \min_{P \in \mathbb{R}^n} \|f_{\text{exp}} - f_{\text{sim}}(P)\| \\ \text{s.t. } lb \le P \le ub \qquad (2) $$

where $f_{\text{exp}}$ represents the experimental measurements and $f_{\text{sim}}$ represents the simulated data in terms of the unknown parameters. The principle of the optimization technique is illustrated in the **SupFig 5**.

<figure><img src="image_13.png" /><figcaption>**SupFig 5.** Principle of the optimization technique. P is the set of unknown parameters, $f_{\text{exp}}$ are the experimental data, $f_{\text{sim}}$ are the simulated data (Source: Created by authors)</figcaption></figure>

The idea of TRRA is to approximate $C(P)$ with a quadratic function $m_k$, to indicate the behavior of function $C(P)$ in a neighborhood N, which is known as the trust region, around the current point x. Applying Taylor's theorem, $C(P_k)$ can be well approximated by the quadratic function $m_k(s)$, in the form:

$$ m_k(s) = C(P_k) + \nabla C(P_k)^\top \cdot s + \frac{1}{2} s^\top \nabla^2 C(P_k)^\top \cdot s \qquad (3) $$

where $\nabla C(P_k) \in \mathbb{R}^n$ is the gradient of $C(P)$ and $\nabla^2 C(P_k) \subset \mathbb{R}^{n \times n}$ the Hessian of $C(P)$ evaluated at $P_k$ and $\nabla$ is the operator of partial derivative. The sub-problem of the trust region method is to set a trial step $s$ by decreasing the area N. If $C(P + s) < C(P)$, the current point P is updated to be $P + s$; the step is considered effective, and the trust region can remain for the next step. If the inequality is not obeyed, the step is ineffective, and consequently P remains unchanged. Moreover, the region N, will be minimized for the next step and the trust region sub-problem will be solved again. Therefore, the challenging issue of the TRRA is to solve its sub-problem, defined as:

$$ \min_{s \in \mathbb{R}^n} m_k(s) \\ \text{s.t. } \|s\|_2 \le \Delta_k \qquad (4) $$

For example, there are directions along which $m_k(s)$ is unbounded from below. In this situation, $\|s\|$ is infinite. Thus, to guarantee the local convergence of the iteration, the search for a minimizer of $m_k$ should be limited to some region around $P_k$, called “Trust Region”. This trust region is usually a sphere defined by Euclidean norm $\|s\| \le \Delta_k$, where the scalar $\Delta_k$ is called the trust region radius.

In our situation, TRRA approach is to estimate the unknown variables to achieve the minimum squared difference between the experimental and simulated data of current and charge densities:

For current density data and a given field E:

$$ \min_{P \in \mathbb{R}^n} C_1(P)_E = \min_{P \in \mathbb{R}^n} \frac{\|j_{\text{exp}}(t) - j_{\text{sim}}(t, P)\|^2}{\|j_{\text{exp}}(t)\|^2} \qquad (5) \\ \text{s.t. } lb \le P \le ub $$

For charge density data:

$$ \min_{P \in \mathbb{R}^n} C_2(P)_E = \min_{P \in \mathbb{R}^n} \frac{\|n_{\text{exp}}(x, t) - n_{\text{sim}}(x, t, P)\|^2}{\|n_{\text{exp}}(x, t)\|^2} \qquad (6) \\ \text{s.t. } lb \le P \le ub $$
where:

* P are the unknown parameters

* j<sub>exp</sub>(t) are the experimental current density (A/m<sup>2</sup>)

* n<sub>exp</sub>(x, t) are the experimental data on net charge density (C/m<sup>3</sup>)

* j<sub>sim</sub>(t, P) are the simulated current density computed using the BCT model in terms of P.

* n<sub>sim</sub>(x, t, P) are the simulated charge densities computed using the BCT model in terms of P.

* C<sub>1</sub>(P) is the cost function that represents the relative error between j<sub>exp</sub>(t) and j<sub>sim</sub>(t, P). It is a unit-less scalar.

* C<sub>2</sub>(P) is the cost function that represents the relative error between n<sub>exp</sub>(x, t) and n<sub>sim</sub>(x, t, P).

* ||·|| is the Euclidean norm or 2-norm.

* lb and ub are lower and upper bounds of the parameters respectively.

## S3. Optimization using charge density data

In this section, TRRA is employed using the charge density measurements considering the same varying electric fields. Thus, the cost function is defined as:

$$C_2(P) = C_2(P)_{20} + C_2(P)_{40} + C_2(P)_{60} \qquad (7)$$

where

$$C_2(P)_E = \frac{\|n_{\text{exp}}(x, t) - n_{\text{sim}}(x, t, P)\|_E^2}{\|n_{\text{exp}}(x, t)\|_E^2} \qquad (8)$$
s. t. $lb \le P \le ub$

where

* n<sub>exp</sub>(x, t) is the experimental net charge density data (C/m<sup>3</sup>)

* n<sub>sim</sub>(x, t, P) is the simulated net charge density.

The optimal set of parameters found by minimizing C<sub>2</sub>(P) is defined as P<sub>2</sub> and the obtained results are displayed in **SupFig 6**. Based on this figure, the simulated patterns are comparable with the experimental ones for all applied fields. Furthermore, most of the properties highlighted in the experiment are reflected by the model, such as the fast injection of negative charges at short time and then the dominance of positive charges inside the bulk after 500 seconds. The relative error between experimental and simulated results using P<sub>2</sub> are presented in **SupTable 2** (S4).

<figure><img src="image_14.png" /><figcaption>**SupFig 6.** Simulated charge density (C/m<sup>3</sup>) profiles using P<sub>2</sub> as a function of the time at room temperature, for fields of 20, 40 and 60 kV/mm (Source: Created by authors)</figcaption></figure>
**SupFig 7** displays the comparison between the simulated current density using $P_2$, and the experimental current for a field of 40 kV/mm. A clear discrepancy is observed between the experimental and simulated current, with a difference by more than a decade and totally different trends versus time. Similar simulation results were obtained with the other fields. Therefore, relying solely on measuring space charge dynamics is insufficient to predict the material's conduction behavior.

<figure><img src="image_15.png" /><figcaption>**SupFig 7.** Experimental and simulated current density using the optimal parameters $P_2$, for applied field of 40 kV/mm (Source: Created by authors)</figcaption></figure>

## S4. Relative Error using $P_1$, $P_2$, $P_3$

The following table shows the relative error between experimental and simulated data using all the obtained parameters, $P_1$, $P_2$ and $P_3$.

* $C_1(P_1)_E$ is the relative error between experimental and simulated **current density** data for each electric field individually using $P_1$.

* $C_2(P_2)_E$ is the relative error between experimental and simulated **charge density** data for each electric field individually using $P_2$.

* $C_1(P_3)_E$ is the relative error between experimental and simulated **current density** data for each electric field individually using $P_3$.

* $C_2(P_3)_E$ is the relative error between experimental and simulated **charge density** data for each electric field individually using $P_3$.

<table><caption>**SupTable 2.** Relative error between experimental and simulated current density for each field individually using parameters sets $P_1$, $P_2$, $P_3$, with $C_1$ the cost function considering current data and $C_2$ the cost function considering charge data. (Source: Created by authors)</caption><thead><tr><th>Relative error</th><th>20 kV/mm</th><th>40 kV/mm</th><th>60 kV/mm</th></tr></thead><tbody><tr><td>$C_1(P_1)_E$</td><td>0.006</td><td>0.001</td><td>0.003</td></tr><tr><td>$C_2(P_2)_E$</td><td>0.021</td><td>0.030</td><td>0.041</td></tr><tr><td>$C_1(P_3)_E$</td><td>0.023</td><td>0.016</td><td>0.044</td></tr><tr><td>$C_2(P_3)_E$</td><td>0.323</td><td>0.499</td><td>0.070</td></tr></tbody></table>
# S5. Mobility and recombination results

## S5.1. Mobility

The mobility values of electrons and holes as functions of the field, obtained using P3, are presented **SupTable 3**. These values were obtained by substituting $w_{moe}$, $w_{moh}$, $a_e$ and $a_h$ in the mobility equation.

<table><caption>**SupTable 3.** The simulated mobility values of electrons and holes at 20°C in $m^2/Vs$. (Source: Created by authors)</caption><thead><tr><th>Electric field</th><th>μ<sub>e</sub></th><th>μ<sub>h</sub></th></tr></thead><tbody><tr><td>20 kV/mm</td><td>3.50 × 10<sup>-13</sup></td><td>1.35 × 10<sup>-14</sup></td></tr><tr><td>40 kV/mm</td><td>5.00 × 10<sup>-13</sup></td><td>1.50 × 10<sup>-14</sup></td></tr><tr><td>60 kV/mm</td><td>9.00 × 10<sup>-13</sup></td><td>1.80 × 10<sup>-14</sup></td></tr></tbody></table>

## S5.2. Recombination coefficients

The space- and time-averaged recombination coefficients could be evaluated by substituting $R_{i=1,...,3}$ in Equations 9, 10 and 11 in order to obtain apparent values and compare the relative importance according to types of charges involved. (See SupFig 8).

For recombination coefficient between mobile electrons (eµ) and trapped holes (ht):

$$S_{e\mu-ht}(t) = \frac{R_1}{Dt_p} \iint \frac{\mu_e(x, t)}{\varepsilon_0 \varepsilon_r} dxdt \qquad (9)$$

Where $t_p$ is the charging time (3 h). For recombination between mobile holes (hµ) and trapped electrons (et):

$$S_{et-h\mu} = \frac{R_2}{Dt_p} \iint \frac{\mu_h(x, t)}{\varepsilon_0 \varepsilon_r} dxdt \qquad (10)$$

For recombination between mobile electrons (eµ) and mobile holes (hµ):

$$S_{e\mu-h\mu} = \frac{R_3}{Dt_p} \iint \frac{(\mu_e(x, t) + \mu_h(x, t))}{\varepsilon_0 \varepsilon_r} dxdt \qquad (11)$$

<figure><img src="image_16.png" /><figcaption>**SupFig 8.** Space- and time-averaged recombination coefficients between different types of carriers. (Source: Created by authors)</figcaption></figure>