Modeling of space charge dynamics in polyethylene under AC stress based on bipolar charge transport model

The bipolar charge transport model is based on bipolar charge injection, transport as well as recombination of opposite charges. The charge carrier dynamics can be described by the following equations [2, 9]:

$$\begin{eqnarray}&&\displaystyle \frac{\partial {n}_{e,h}\left(x,t\right)}{\partial t}+\displaystyle \frac{1}{e}\displaystyle \frac{\partial {j}_{e,h}\left(x,t\right)}{\partial x}={s}_{e,h}\left(x,t\right)\,{\rm{continuity}}\end{eqnarray} \tag{ 1 }$$

$$\begin{eqnarray}&&{j}_{e,h}\left(x,t\right)=e{\mu }_{e,h}E\left(x,t\right){n}_{e,h}^{f}\left(x,t\right)\,{\rm{transport}}\end{eqnarray} \tag{ 2 }$$

$$\begin{eqnarray}&&\displaystyle \frac{\partial E\left(x,t\right)}{\partial x}=\displaystyle \frac{e\left({n}_{h}-{n}_{e}\right)}{\varepsilon }\,{\rm{Poisson}}\end{eqnarray} \tag{ 3 }$$

where n is the number densities, j is the current density, e is the electron charge, s is the source term resulting from recombination processes, μ is the mobility, E is the electric field, ε is the permittivity (2.3 for polyethylene). The subscript 'e' and 'h' refer to the electrons and holes, respectively. The superscript 'f' refers to the mobile carriers. By assuming that there is an exponential distribution of the trap density and the temperature T is much less than the shape parameter T_0(e,h) , the density of mobile carriers can be approximated by

$$\begin{eqnarray}&&{n}_{e,h}^{f}=\displaystyle \frac{{n}_{e,h}}{1+a},a=\displaystyle \frac{{T}_{0\left(e,h\right)}}{T}\end{eqnarray} \tag{ 4 }$$

The source term can be expressed as

$$\begin{eqnarray}&&{s}_{e,h}={S}_{1}{n}_{e}^{f}\left({n}_{h}-{n}_{h}^{f}\right)+{S}_{2}\left({n}_{e}-{n}_{e}^{f}\right){n}_{h}^{f}\end{eqnarray} \tag{ 5 }$$

where S₁ and S₂ are the recombination coefficients for mobile electrons-trapped holes and trapped electrons-mobile holes, respectively.

Based on the hopping mechanism of charge transport in polyethylene, the mobility has the form

$$\begin{eqnarray}&&{\mu }_{e,h}=\displaystyle \frac{2\nu d}{E}\exp \left(-\displaystyle \frac{{{\rm{\Delta }}}_{f\left(e,h\right)}}{{k}_{B}T}\right)\sinh \left(\displaystyle \frac{eE{d}_{\left(e,h\right)}}{2{k}_{B}T}\right)\end{eqnarray} \tag{ 6 }$$

$$\begin{eqnarray}&&{{\rm{\Delta }}}_{f\left(e,h\right)}=-{k}_{B}{T}_{0\left(e,h\right)}\,\mathrm{ln}\left[\displaystyle \frac{{n}_{\left(e,h\right)}}{{{N}^{^{\prime} }}_{e,h}{k}_{B}{T}_{0\left(e,h\right)}}+\exp \left(-\displaystyle \frac{{{\rm{\Delta }}}_{m\left(e,h\right)}}{{k}_{B}{T}_{0\left(e,h\right)}}\right)\right]\end{eqnarray} \tag{ 7 }$$

Where ν is the frequency of escape (6.2 × 10¹² s⁻¹ at room temperature), Δ_f is the upper filled level, Δ_m is the maximum trap depth, N' is the pre-exponential factor of the trap density distribution, T₀ is the shape parameter, d_(e,h) is the average distance between electron or hole traps.

An injection flux based on Schottky law is of the form

$$\begin{eqnarray}&&{j}_{e}\left(0,t\right)=A{T}^{2}\exp \left(-\displaystyle \frac{e{w}_{ei}}{{k}_{B}T}\right)\exp \left(\displaystyle \frac{e}{{k}_{B}T}\sqrt{\displaystyle \frac{eE\left(0,t\right)}{4\pi \varepsilon }}\right)\end{eqnarray} \tag{ 8 }$$

$$\begin{eqnarray}&&{j}_{h}\left(d,t\right)=A{T}^{2}\exp \left(-\displaystyle \frac{e{w}_{hi}}{{k}_{B}T}\right)\exp \left(\displaystyle \frac{e}{{k}_{B}T}\sqrt{\displaystyle \frac{eE\left(d,t\right)}{4\pi \varepsilon }}\right)\end{eqnarray} \tag{ 9 }$$

for electrons at cathode and holes at anode. A = 1.2 × 10⁶ Am⁻²K⁻² is the Richardson constant, k_B is the Boltzmann constant, w_ei and w_hi are injection barrier for electrons and holes, d is the dielectric thickness.

The extraction fluxes are of the form

$$\begin{eqnarray}&&{j}_{e}\left(d,t\right)=e{\mu }_{e}E\left(d,t\right){n}_{e\mu }\left(d,t\right)\end{eqnarray} \tag{ 10 }$$

$$\begin{eqnarray}&&{j}_{h}\left(0,t\right)=e{\mu }_{h}E\left(0,t\right){n}_{h\mu }\left(0,t\right)\end{eqnarray} \tag{ 11 }$$

for electrons at anode and holes at cathode.

A sinusoidal, triangular or square alternating voltage is applied to the anode (x = d), and at the cathode (x = 0) the potential is set to zero. The conduction current density is

$$\begin{eqnarray}&&{j}_{{\rm{c}}}(t)=\displaystyle {\int }_{0}^{d}\left({j}_{e}\left(x,t\right)+{j}_{h}\left(x,t\right)\right)dx/d\end{eqnarray} \tag{ 12 }$$

Note that the displacement current is not discussed here.

As the spaced average combination rate is proportional to the electroluminescence intensity, the electroluminescence intensity can be defined by

$$\begin{eqnarray}&&{\rm{EL}}=\displaystyle {\int }_{0}^{d}\left({s}_{e}\left(x,t\right)+{s}_{h}\left(x,t\right)\right)dx/d\end{eqnarray} \tag{ 13 }$$

Sets of parameters appear suitable to describe space charge dynamics under DC stress [2]. Almost the same parameters are adopted for AC conditions except for the lower initial charge density and the mobility. Parameters used in the model are shown in table 1.

The domain of calculation is divided into 2000 elements, of vary lengths for 100 μm thick polyethylene thin film: 50 elements within a distance of 10 nm from the electrode, 90 elements at the distance of 10 nm ∼ 100 nm, 180 elements at the distance of 100 nm ∼ 1 μm, 360 elements at the distance of 1 μm ∼ 10 μm, and 640 elements at the distance of 10 μm ∼ 90 μm from the cathode. It is noted that it is necessary to use variable space discretization, because space charge layer thickness is much thinner compared with those under DC stress. The time step Δt is chosen to satisfy the Courant-Friedrichs-Lewy condition whose typical value is 2.0 × 10⁻⁵ s. As the continuity equation (1) is the convection dominated convection–diffusion equation, accurate numerical algorithm without under- and overshoot is vital for reproducing charge dynamics in polyethylene. Three-order flux limiter is adopted to achieve better accuracy than first-order upwind, as well as positivity without introducing numerical diffusion or oscillation, which is described in detail in [2].

Square, triangular or sinusoidal voltage is applied to the anode, as shown in figure 1. Numerical simulations are performed for a stress time of 5 s. The amplitude and frequency of all three waveforms are 6 kV and 50 Hz, respectively. Current density under sinusoidal waveform is shown in figure 2. It can be seen that the phase difference between applied voltage and current density is 180°. The reason is that when the applied voltage is positive, the direction of electron drift is along the positive x direction, and current density is negative. The maximum current density is about 1.1 × 10⁻²⁰ A m⁻², which disagrees with the results of Baudoin et al [9]. The maximum current density of Baudoin et al is about 1.7 × 10⁻¹ A m⁻². When the initial charge density is 0.5 Cm⁻³ and 60 kV mm⁻¹ DC electric field is applied, the current density is 1.3 × 10⁻⁹ A m⁻², which agrees with the experimental results [2]. It can be inferred that the current density of Baudoin et al far exceeds the actual values.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Current densities under square, triangular or sinusoidal voltages are shown in figure 3. It can be seen that the frequency and the maximum of three current densities are the same, and the phase differences between current densities and applied voltages are all 180°.

Zoom In Zoom Out Reset image size

Electroluminescence intensity under sinusoidal waveform is shown in figure 4. It can be seen that electroluminescence intensity reaches its peak when the applied voltage reaches its positive and negative peak, while electroluminescence intensity reaches its minimum value when the applied voltage reaches zero. It is noted that there is no much difference between the maximum and minimum of the electroluminescence intensity.

Zoom In Zoom Out Reset image size

Electroluminescence intensities under square, triangular or sinusoidal voltages are shown in figure 5. It can be seen that there is much difference between electroluminescence intensities under different waveforms. Electroluminescence intensity under square voltage is almost straight, and electroluminescence intensity under sinusoidal voltage is almost 100 Hz sinusoidal waveform. However, electroluminescence intensity under triangular voltage appears the bathtub-shaped. Electroluminescence intensity decreases relatively rapidly, and increases much more slowly, which agrees the results of Baudoin et al [9].

Zoom In Zoom Out Reset image size

Space charge density under sinusoidal waveform is shown in figure 6. It can be seen that the distribution of total electrons and mobile electrons is similar, and the distribution of total holes and mobile holes is similar. The reason is that mobile carrier densities are approximated by formula (4). Electron and hole densities reach the maximum at t = 5 ms and t = 15 ms. The maximum electron and hole densities are 1.4 × 10⁸ m⁻³ and 1.1 × 10⁹ m⁻³, respectively. This is because that electron injection barrier is larger than hole injection barrier as shown in table 1, so more holes are injected from the anode. The charge layer thickness is about 0.18 ∼ 0.30 nm from the electrode, which agrees with the results of Baudoin et al in magnitude. It is hereby certified that variable space discretization, being tightened close to the electrodes is essential to achieve enough spatial resolution of space charge dynamics.

Zoom In Zoom Out Reset image size

Space charge density under triangular waveform is shown in figure 7. It can be seen from figure 7 that electron and hole densities reach the maximum at t = 0 ms, t = 10 ms and t = 20 ms. The maximum electron and hole densities under triangular voltage are the same as those under sinusoidal voltage. However, high carrier densities under triangular voltage last shorter than those under sinusoidal voltage.

Zoom In Zoom Out Reset image size

Space charge density under square waveform is shown in figure 8. It can be seen from figure 8 that carrier distribution under square voltage is quite different from that under sinusoidal and triangular voltage. The carrier distribution at the same distance from the electrodes is approximately the same. When the voltage polarity reversal at the anode occurs, the carrier distribution changes little. The reason is that the characteristic time to redistribute charge carriers is the order of several seconds. When square voltage with 50 Hz frequency is applied, the time interval of 10 ms is not enough to redistribute charge carriers. The charge layer thickness under square voltage is a little larger than that under sinusoidal or triangular voltage.

Zoom In Zoom Out Reset image size