# Space Charge Accumulation and Decay in Dielectric Materials with Dual Discrete Traps
## 2 Dynamic Responses of Charge Trapping and De-trapping
### Instantaneous Charge Injection
First, the charge trapping term and the charge de-trapping term in the charge trapping/de-trapping Equation (3) is split into two parts. Secondly, we just solve the charge trapping term in the differential equation of free charge density, and the de-trapping terms in the differential equations of trapped charge densities by the NSFD method. \[\frac{\mathrm{d}q_{free}(t)}{\mathrm{d}t}\ =\ -q_{free}(t)\sum_{i\ =\ 1}^{n}P_{Ti}\!\!\left(1-\frac{q_{trapi}(t)}{Q_{Ti}}\right) \tag{17}\]  
\[\frac{\mathrm{d}q_{trapi}(t)}{\mathrm{d}t}\ =\ -P_{Di}q_{trapi}(t),(i\ =\ 1,2) \tag{18}\]  
It is assumed that \(\Delta t\) is the time step size, and \([t_{k}=k\Delta t]_{k\geq 0}\) is a sequence of equally-spaced time points. Then, the differential equations can be transformed into difference equations. Adopting the NSFD scheme, we can obtain the numerical form of Equations (17) and (18) as follows:  
\[q_{free}(t_{k+1})\ =\ q_{free}(t_{k})\exp\!\left[-\sum_{i\ =\ 1}^{n}P_{Ti}\!\! \left(1-q_{trapi}(t_{k})/Q_{Ti}\right)\Delta t\right] \tag{19}\]  
\[q_{trapi}(t_{k+1})\ =\ q_{trapi}(t_{k})\exp(-P_{Di}\Delta t),\ i\ =\ 1,2 \tag{20}\]  
The numerical forms are unconditionally stable regardless of \(\Delta t\), which is superior to the finite different method (FDM). According to [PERSON] et al. [29], the Lambert \(W\) function is used in Equations (19) and (20) (the NSFD scheme) to increase the accuracy of numerical solutions. Thirdly, the charge de-trapping term in the differential equation of free charge density and the charge trapping terms in the differential equation of trapped charge densities can be computed explicitly using the trapped charge densities and free charge densities at \(t_{k+1}\) obtained by the NSFD scheme. Then, we can obtain numerical solutions of free charge density \(q_{free}(t_{k+1})\) and trapped charge densities \(q_{trapi}(t_{k+1})\) at \(t_{k+1}\). Figure 2a,b shows the analytical and numerical results of free charge density and trapped charge densities as a function of time with the total charge densities of 2 cm\({}^{-3}\) and 10 cm\({}^{-3}\), respectively. When the total charge density is relatively low (e.g., 2 cm\({}^{-3}\)), the approximations, \(q_{trapi}(t)/Q_{Ti}\approx 0\), are valid, and the analytical results are in good agreement with the numerical results. Nevertheless, when the total charge density is relatively high (e.g., 10 cm\({}^{-3}\)), the approximations, \(q_{trapi}(t)/Q_{Ti}\approx 0\), are not always valid. Therefore, a small difference occurs between the analytical results and the numerical results in the time range from \(10^{3}\) to \(10^{4}\) s. However, the shape of the curves is nearly the same. The free charge density decays exponentially with time both in Figure 2a,b. The trapped charge density in trap level 2, \(q_{trap2}(t)\), increases with time during the charge trapping/de-trapping process. Whereas, the trapped charge density in trap level 1, \(q_{trap1}(t)\), increases firstly and then decreases with time. Since the trapping probability of trap center is proportional to the trap density, the trapping probability of trap level 1 will be five times larger than that of trap level 2. Consequently, more free charges are captured by trap level 1 at the beginning. In addition, the time needed to reach steady state for charges in trap level 2 is prolonged in the presence of trap level 1.