4.2.5 The Second Iteration of the Forward Resolvent Series

The structure of the algorithm can be seen by derivation of the second iteration term of the Neumann series for (4.24) using the forward formulation to obtain the ensemble Monte Carlo algorithm. This algorithm does not give the distribution function value at a given point, but instead the relative number of carriers in a region $\Omega$ with volume $\Delta\vec{k}$ around point $\vec{k}$ in the quasi-momentum space^4.7.

The relative number of carriers in $\Omega$ determined by the small perturbation $f_{1}(\vec{k},t)$ is defined as:

$$(f_{1})_{\Omega}=\int_{\Delta_{\vec{k}}}f_{1}(\vec{k},t)\,d\vec{k}=\int f_{1}(\vec{k},t)\theta_{\Omega}(\vec{k})\,d\vec{k},$$ (4.37)

where the domain indicator $\theta_{\Omega}(\vec{k})$ has bee introduced as a function which evaluates to one if $\vec{k}\in\Omega$ and zero otherwise. In order to use the phase space volume conservation law described in Section 2.1.3, an additional integral over the real space is added and (4.38) is rewritten as

$$(f_{1})_{\Omega}=\frac{1}{V_{cr}}\int\,d\vec{k}\int\,d\vec{r}f_{1}(\vec{k},t)\theta_{\Omega}(\vec{k}),$$ (4.38)

where $V_{cr}$ is a crystal volume. From (4.37) the second iteration term is obtained:

$$(f_{1}^{(2)})_{\Omega}=\frac{1}{V_{cr}} \int_{0}^{t}\,dt_{1}\int_... ...exp\biggl(-\int_{0}^{t_{2}}\widetilde{\lambda}[\vec{K}_{2}(y)]\,dy\biggr)\times$$ (4.39)

$$\times\widetilde{S}(\vec{k}_{2},\vec{K}_{1}(t_{2}))\exp\biggl(-\i... ..._{t_{1}}^{t}\widetilde{\lambda}[\vec{K}(y)]\,dy\biggr)\theta_{\Omega}(\vec{k}).$$

Considering $\vec{r}$ as $\vec{R}(t)$ which is the real space trajectory corresponding to the quasi-momentum trajectory (4.31) it follows from the phase space volume conservation^4.8 that $d\vec{k}d\vec{r}=d\vec{K}(t_{1})d\vec{R}(t_{1})$. Denoting $\vec{K}(t_{1})$ with $\vec{k}_{1}^{a}$ (4.31) can be rewritten in a forward initialization:

$$\vec{K}(t)=\vec{k}_{1}^{a}+\frac{q}{\hbar}\vec{E}_{s}(t-t_{1}).$$ (4.40)

In the same manner combining $d\vec{k}_{1}$ with $d\vec{R}(t_{1})$ gives $d\vec{k}_{1}d\vec{R}(t_{1})=d\vec{K}_{1}(t_{2})d\vec{R}_{1}(t_{2})$. Denoting $\vec{K}_{1}(t_{2})$ with $\vec{k}_{2}^{a}$ (4.31) can be rewritten in a forward initialization:

$$\vec{K}_{1}(t_{1})=\vec{k}_{2}^{a}+\frac{q}{\hbar}\vec{E}_{s}(t_{1}-t_{2}).$$ (4.41)

Finally $d\vec{k}_{2}d\vec{R}_{1}(t_{2})=d\vec{K}_{2}(0)d\vec{R}_{2}(0)$. Denoting $\vec{K}_{2}(0)$ as $\vec{k}_{i}$ one obtains for the corresponding forward initialization:

$$\vec{K}_{2}(t_{2})=\vec{k}_{i}+\frac{q}{\hbar}\vec{E}_{s}t_{2}.$$ (4.42)

Integrating out $R_{2}(0)$ cancels the crystal volume $V_{cr}$ and the final expression is:

$$(f_{1}^{(2)})_{\Omega}=\int_{0}^{t}\,dt_{1}\int_{0}^{t_{1}}\,dt_{... ...exp\biggl(-\int_{0}^{t_{2}}\widetilde{\lambda}[\vec{K}_{2}(y)]\,dy\biggr)\times$$ (4.43)

$$\times\widetilde{S}(\vec{K}_{2}(t_{2}),\vec{k}_{2}^{a})\exp\biggl... ...{1}}^{t}\widetilde{\lambda}[\vec{K}(y)]\,dy\biggr) \theta_{\Omega}(\vec{K}(t)).$$

To express the time integrals in a forward way the following identity is used:

$$\int_{0}^{t}\,dt_{1}\int_{0}^{t_{1}}\,dt_{2}=\int_{0}^{t}\,dt_{2}\int_{t_{2}}^{t}\,dt_{1},$$ (4.44)

which is schematically shown in Fig. 4.3.

Figure 4.3: The same integration area can be covered either vertically a) or horizontally b).

Now the second iteration term has the form:

$$(f_{1}^{(2)})_{\Omega}= \int_{0}^{t}\,dt_{2}\int_{t_{2}}^{t}\,dt_{1}\int\,d\vec{k}_{2}^{a}\int\,d\vec{k}_{1}^{a}\int\,d\vec{k}_{i}\{G(\vec{k}_{i})\}\times$$

$$\times\biggl\{\exp\biggl(-\int_{0}^{t_{2}}\widetilde{\lambda}[\ve... ...t_{2}),\vec{k}_{2}^{a}]}{\widetilde{\lambda}[\vec{K}_{2}(t_{2})]}\biggr\}\times$$ (4.45)

$$\times\biggl\{\exp\biggl(-\int_{t_{2}}^{t_{1}}\widetilde{\lambda}... ...t_{1}),\vec{k}_{1}^{a}]}{\widetilde{\lambda}[\vec{K}_{1}(t_{1})]}\biggr\}\times$$

$$\times\exp\biggl(-\int_{t_{1}}^{t}\widetilde{\lambda}[\vec{K}(y)]\,dy\biggr)\theta_{\Omega}(\vec{K}(t)),$$

where $\vec{k}^{a}$ stands for an after-scattering wave vector and $\vec{k}_{i}$ denotes an initial wave vector. The quantity $\widetilde{S}[\vec{k},\vec{k}^{'}]/\widetilde{\lambda}[\vec{k}]$ represents a normalized after-scattering distribution. As can be seen from (4.13) and (4.14) it is normalized to unity. It follows from (4.47), that during Monte Carlo simulation a particle trajectory is constructed in terms of new quantities $\widetilde{S}$ and $\widetilde{\lambda}$.

S. Smirnov: