format: text

text:
# Optical properties of self-assembled quantum wires for application in infra-red detection

## I Introduction

Quantum-well infra-red photodetectors (QWIP’s) have been extensively studied in recent years. The main mechamism used in QWIPs is the inter-subband optical transition, because the wavelengths for these transitions in typical III-V quantum wells can be tailored to match the desired operating wavelength (1-20 $`\mu m`$) for infra-red (IR) detection. Due to its narrow band absorption, QWIP’s are complementary to the traditional HgCdTe detectors, which utilize the inter-band absorption, and therfeore are applicable only for broad-band absorption. The main drawback of QWIP’s is the lack of normal-incidence capability, unless some processing is made to create diffraction gartings on the surface, which tends to reduce the responsivity of the material to the incident radiation. Because electrons in quantum wells have translational invariance (within the effective-mass model) in the plane normal to the growth axis, the electron inter-subband transitions for normal-incident radiation is zero (or very small even if the coupling with other bands is considered). One way to break the translational invariance is to introduce the surface diffraction grating as commonly adopted in many QWIP’s fabricated today. A better (and less expensive) way to break the in-plane translational invariance is to utilize the strain-induced lateral modulation provided in self-assembled nano-structure materials. These nano-structures inculde quantum dots and quantum wires. Because the lateral modelation is formed via self-assembly, the fabrication of this type of materials will be much more efficient once the optimized growth parameters are known. Hence, it will be cost effective to use them for device fabrications.

Self-assembled III-V QWR’s via the strain-induced lateral-layer ordering (SILO) process have attracted a great deal of attention recently. \[-\] The self-assembly process occurs during the growth of short-period superlattices (SPS) \[e.g. (GaAs)2/(InAs)2.25\] along the  direction on InP substrate. The excess fractional InAs layer leads to stripe-like islands during the initial MBE growth. The presence of stripes combined with strain leads to natural phase separation as additional layers of GaAs or InAs are deposited and the structure becomes laterally modulated in terms of In/Ga composition. A self-assembled QWR heterostructure can then be created by sandwiching the laterally modulated layer between barrier materials such as Al0.24Ga0.24In0.52As (quarternary), Al0.48In0.52As (ternary), or InP (binary).\[4-6\] It was found that different barrier materials can lead to different degree of lateral composition modulation, and the period of lateral modulation ranges from 100Å to 300 Å depending on the growth time and temperature.

In this paper, we explore the usefulness of InGaAs quantum wires (QWR’s) grown by the strain-induced lateral ordering (SILO) process for IR detection. Our theoretical modeling inculdes the effects of realistic band structures and microscopic strain distributions by combining the effective bond-orbital model (EBOM) with the valence-force-field (VFF) model. One of the major parameters for the IR detectors is the absorption quantum efficiency which is directly related to the absorption coefficient by $`\eta =1e^{\alpha l}`$ where $`\alpha `$ is the absorption coefficient and $`l`$ is the sample length. Thus, to have a realistic accessment of the materials for device application, we need to perform detailed calculations of the absorption coefficient, taking into account the excitonic and band structure effects. Both inter-subband and inter-band transitions are examined systematically for a number of structure parameters (within the experimentally feasible range) chosen to give the desired effect for IR detection.

It is found that the wavelengths for the inter-subband transitions of InGaAs self-assembled QWR’s range from 10 to 20 $`\mu m`$, while the inter-band transitions are around 1.5 $`\mu m`$. Thus, the material provides simultaneneous IR detection at two contrasting wavelengths, something desirable for application in multi-colored IR vedio camera.

Several structure models with varying degrees of alloy mixing for lateral modulation are considered. For the inter-band absorption, the excitonic effect is important, since it gives rise a large shift in transition energy and substantial enhancement of the absorption spectrum. To study the excitonic effect on the absorption spectrum for both discrete and contunuum states, we use a large set of basis functions with a finite-mesh sampling in the k-space and diaginalize the exciton Hamiltonian directly. Emphasis is put on the analysis of line shapes of various peak structures arising from discrete excitonic states of one pair of subbands coupled with the excitonic (discrete and continuum) states associated with other pairs of subbands. We find that the excitonic effect enhances the first absorption peak around 1.5 times and shifts the peak position by 20-30meV.

## II Theoretical model

The QWR structures considered here consist of 8 pairs of (GaAs)2(InAs)2.25 short-period superlattices (SPS) sandwiched between Al0.24Ga0.24In0.52As barriers. The SPS structure prior to strain induced lateral ordering (SILO) is depicted in Fig. 1. With lateral ordering, the structure is modeled by a periodic modulation of alloy composition in layers with fractional monolayer of (In or Ga) in the SPS structure. In layers 7 and 9 (starting from the bottom as layer 1), we have

$$x_{\text{In}}=\{\begin{array}{cc}x_m[1\mathrm{sin}(\pi y^{}/2b)]/2\hfill & \text{ for }y^{}<b\hfill \\ 0\hfill & \text{ for }b<y^{}<L/2b\hfill \\ x_m\{1+\mathrm{sin}[\pi (y^{}L/2)/2b]\}/2\hfill & \text{ for }L/2b<y^{}<L/2+b\hfill \\ x_m\hfill & \text{ for }L/2+b<y^{}<Lb\hfill \\ x_m\{1\mathrm{sin}[\pi (y^{}L)/2b]\}/2\hfill & \text{ for }y^{}>Lb,\hfill \end{array}$$
(1)
where $`x_m`$ is the maximum In composition in the layer, $`2b`$ denotes the width of lateral composition grading, and $`L`$ is the period of the lateral modulation in the  direction. The experimental feasible range of $`L`$ is between 100 Åand 300 Å. The length of $`L`$ is controled by the growth time and temperature.

In layers 3 and 13, we have

$$x_{\text{In}}=\{\begin{array}{cc}0\hfill & \text{ for }0<y^{}<5L/8b\hfill \\ x_m\{1+\mathrm{sin}[\pi (y^{}5L/8)/2b]\}/2\hfill & \text{ for }5L/8b<y^{}<5L/8+b,\hfill \\ x_m\hfill & \text{ for }5L/8+b<y^{}<7L/8b\hfill \\ x_m\{1\mathrm{sin}[\pi (y^{}7L/8)/2b]\}/2\hfill & \text{ for }7L/8b<y^{}<7L/8+b\hfill \\ 0\hfill & \text{ for }7L/8+b<y^{}<L.\hfill \end{array}$$
(2)
Similar equation for $`x_{\text{Ga}}`$ in layers 5 and 11 can be deduced from the above. By varying the parameters $`x_m`$ and $`b`$, we can get different degrees of lateral alloy mixing. Typically $`x_m`$ is between 0.6 and 1, and $`b`$ is between zero and $`15a_{[110]}`$ 62 Å.

A VFF model\[13-15\] is used to find the equilibrium atomic positions in the self-assembled QWR structure by minimizing the lattice energy. The strain tensor at each atomic (In or Ga) site is then obtained by calculating the local distortion of chemical bonds.

Once the microscopic strain distribution in the model structure is determined, the energy levels and wave-functions of self-assembled quantum wires are then calculated within the effective bond-orbital model (EBOM). Detailed description of this method can be found in Refs. ,-. EBOM used here is a tight-binding-like model in which two s- like conduction bands (including spin) and four valence bands with total angular momentum $`J=3/2`$ (due to spin-orbit coupling of p-like orbitals with the spinor). Thus, the present model is comparable to the six-band $`𝐤𝐩`$ model as adpoted in Ref. \[Ehe\]

To minimize the computing effort, we express the electron and hole states for the quantum wire structures in terms of eigen-states of a quantum well structure with different in-plane wave vectors. The quantum well consists of 8 pairs of (GaAs)2(InAs)2 short-period superlattice (SPS) plus two InAs monolayers (one inserted after the second pair of SPS and the other after the sixth pair of SPS), so the the total In/Ga composition ratio is consistent with the (GaAs)2(InAs)2.25 SPS. The whole stack of SPS’s is then sandwiched between two slabs of Al0.24Ga0.24In0.52 barriers. Let us denote the quantum well eigen-states as $`|n,k_1,k_2_{QW}`$ where $`n`$ labels the subband, $`k_1`$ denotes the wave vector along the wire ($`[1\overline{1}0]`$) direction and $`k_2`$ labels the wave vector in the  direction, which is perpendicular to the wire and the growth axis. Expanding the quantum well states in terms of bond-orbitals, we have

$`|n,k_1,k_2_{QW}={\displaystyle \frac{1}{\sqrt{L}}}{\displaystyle \underset{\alpha ,𝐑}{}}f_{n,k_1,k_2}(\alpha ,R_z)\mathrm{exp}(ik_2R_2)\mathrm{exp}(ik_1R_1)|u_\alpha (𝐑),`$

where $`L`$ is the sample length along the wire axis, $`f_{n,k_1,k_2}(\alpha ,R_z)`$ is the eigen-vector for the quantum well Hamiltonian and $`u_\alpha (𝐑)`$ denotes an $`\alpha `$-like bond orbital state at site $`𝐑`$ ($`\alpha =1,\mathrm{},6`$ for two s-like conduction-band and four $`J=3/2`$ valence-band orbitals). Here $`𝐑`$ runs over all lattice sites within the SPS layer (well region) and AlGaInAs layer (barrier region).

We then diagonalize the hamiltonian for the quantum wire (QWR) within a basis which consists of the quantum well states with $`k_2`$’s separated by reciprocal lattice vectors $`g_m=m(2\pi /a_{[110]})`$; $`m`$= integers. Namely,

$`|i,k_1,k_2={\displaystyle \underset{n,m}{}}C_{i,k_1}(n,k_2+g_m)|n,k_1,k_2+g_m_{QW}`$

where $`C_{i,k_1}(n,k_2+g_m)`$ is the eigen-vector for the quantum-wire hamiltonian matrix for the $`ith`$ QWR subband at wave vector $`(k_1,k_2)`$.

In terms of the bond orbitals, we can rewrite the QWR states as

$`|i,k_1,k_2={\displaystyle \underset{\alpha ,𝐑}{}}F_{i,k_1,k_2}(\alpha ,𝐑)|u_\alpha (𝐑)`$

where

$`F_{i,k_1,k_2}(\alpha ,𝐑)={\displaystyle \frac{1}{\sqrt{L}}}{\displaystyle \underset{n,m}{}}C_{i,k_1}(n,k_2+g_m)f_{n,k_1,k_2+g_m}(\alpha ,R_z)\mathrm{exp}[i(k_2+g_m)R_2]\mathrm{exp}(ik_1R_1)`$

is the QWR envelope function. For the laterally confined states, the dispersion of bands versus $`k_2`$ is negligible; thus, the $`k_2`$ dependence can be ignored.

The absorption coefficeient for inter-subband transitions between subbands $`i`$ and $`j`$ is given by

$$\alpha _{ij}(\mathrm{}\omega )=\frac{4\pi ^2e^2\mathrm{}}{n_rm_e^2cV\mathrm{}\omega }\underset{k_1,k_2}{}|i,k_1,k_2|\widehat{ϵ}𝐩|j,k_1,k_2|^2[f_j(k_1,k_2)f_i(k_1,k_2)]\delta (E_j(k_1,k_2)E_i(k_1,k_2)\mathrm{}\omega )$$
(3)
where $`n_r`$ is the refractive index of the QWR, $`V`$ is the volume of the QWR sample restricted within the SPS region, $`f_i(f_j)`$ is the Fermi-Dirac distribution function for subbnad $`i`$ ($`j`$). The optical matrix elements between QWR subband states are related to those between bond orbitals by

$`i,k_1,k_2|\widehat{ϵ}𝐩|j,k_1,k_2={\displaystyle \underset{\alpha ,\alpha ^{},\tau }{}}F_{i,k_1,k_2}^{}(\alpha ,𝐑)F_{j,k_1,k_2}(\alpha ^{},𝐑)u_\alpha (𝐑)|\widehat{ϵ}𝐩|u_\alpha ^{}(𝐑+\stackrel{}{\tau }),`$

where $`\stackrel{}{\tau }`$ runs over on-site or the 12 nearest-neighbor sites in the fcc lattice. The optical matrix elements between bond orbitals are related to the band parameters by requiring the optical matrix elements between bulk states near the zone ceneter to be identical to those obtained in the $`𝐤𝐩`$ theory. We obtain

$`langleu_s(𝐑)|p_\alpha u_\alpha ^{}(𝐑)=\sqrt{m_0E_p/2}\delta _{\alpha ,\alpha ^{}};\alpha =x,y,z,`$

$`langleu_s(𝐑)|p_\alpha u_s(𝐑+\stackrel{}{\tau })=(\mathrm{}/2\sqrt{2}a)(E_p/Egm_0/m_e^{})\tau _\alpha ;\alpha =x,y,z,`$

where $`\tau _\alpha `$ is the $`\alpha `$-th of the lattice vector $`\tau `$ in units of a/2, $`E_p`$ is the inter-band optical matrix element as defined in Ref. , and $`m_e^{}`$ is the electron effective mass.

Next, we study the inter-band transitions. For this case, the excitonic effect is important. Here we are only interested in the absorption spectrum near the band edge due to laterally confined states. Thus, the dispersion in the $`k_2`$ direction can be ignored. The exciton states with zero center-of-mass momentum can then be written as linear combinations of products of electron and holes states associated with the same $`k_1`$ (wave vector along the wire direction). We write the electron-hole product state for the $`i`$-th conduction subband and $`j`$-th valence subband as

$`|i,j;k_1_{ex}=|i,k_1|j,k_1`$

$`{\displaystyle \underset{\alpha ,\beta ,𝐑_e,𝐑_h}{}}F_{i,k_1}(\alpha ,𝐑_e)G_{j,k_1}(\beta ,𝐑_h)|u(\alpha ,𝐑_e)>|u(\beta ,𝐑_h)>.`$

The matrix elements of the exciton Hamiltonian within this basis is given by

$$i,j,k_1|H_{ex}|i^{},j^{},k_1^{}=[E_i(k_1)\delta _{i,i^{}}E_j(k_1)\delta _{j,j^{}}]\underset{𝐑_e,𝐑_h}{}_{ii^{}}^{}(𝐑_e)v(𝐑_e𝐑_h)𝒢_{jj^{}}(𝐑_2),$$
(4)
where $`v(𝐑_e,𝐑_h)=\frac{4\pi e^2}{ϵ(0)|𝐑_e𝐑_h|}`$ is the coulomb interaction between the electron and hole screened by the static dielectric constant $`ϵ(0)`$, and

$`_{ii^{}}(𝐑_e)={\displaystyle \underset{\alpha }{}}F_{i,k_1}^{}(\alpha ,𝐑_e)F_{i^{},k_1}(\alpha ,𝐑_e)`$

describes the charge density matrix for the electrons. Similarly,

$`𝒢_{jj^{}}(𝐑_h)={\displaystyle \underset{\beta }{}}G_{j,k_1}^{}(\beta ,𝐑_h)G_{j^{},k_1}(\beta ,𝐑_h)`$

describes the charge density matrix for the holes. In Eq. (x), we have adopted the approximation

$`u(\alpha ,𝐑_e)|u(\beta ,𝐑_h)|v|u(\alpha ^{},𝐑_e^{})|u(\beta ^{},𝐑_h^{})v(𝐑_e𝐑_h)\delta _{\alpha ,\alpha ^{}}\delta _{\beta ,\beta ^{}}\delta _{𝐑_e,𝐑_e^{}}\delta _{𝐑_h,𝐑_h^{}},`$

since the Coulomb potential is a smooth function over the distance of a lattice copnstant, except at the origin, and the bond orbitals are orthonormal to each other. At the origin ($`𝐑_e=𝐑_h`$), the potential is singular, and we replace it by an empirical constant which is adjusted so as to give the same exciton binding energy as obtained in the effective-mass theory for a bulk system. The results are actually insensitive to the on-site Coulomb potential parameter, since the Bohr radius of the exciton is much larger than the lattice constant.

After the diagonalization, we obtain the excitonic states as linear combinations of the electron-hole product states, and the inter-band absorption coefficient is computed according to

$$\alpha ^{ex}(\mathrm{}\omega )=\frac{4\pi ^2e^2\mathrm{}}{n_rm_e^2cA\mathrm{}\omega }\underset{n}{}|\underset{k_1,i,j}{}C_{ij}^n(k_1)\widehat{ϵ}𝐩_{ij}(k_1)|^2\delta (\mathrm{}\omega E_n^{ex}),$$
(5)
where $`A`$ is the cross-sectional area of the SPS region within the QWR unit cell (as depicted in Fig. 1). $`C_{ij}^n(k_1)`$ is the $`n`$-th eigen-vector obtained by diagonalizing the exciton Hamiltonian of Eq. (4). The momentum matrix elements $`𝐩_{ij}(k_1)`$ are given by

$$𝐩_{ij}(k_1)=\underset{\alpha ,\beta ,\stackrel{}{\tau }}{}F_{i,k_1}^{}(\alpha ,𝐑_e)G_{j,k_1}(\beta ,𝐑_h+\stackrel{}{\tau })u_\alpha (𝐑_e)|\widehat{ϵ}𝐩|u_\beta (𝐑_h+\stackrel{}{\tau }).$$
(6)
In the absence of band mixing, the conduction subband state reduces to a pure $`s`$-like state and the valence subband state reduces to a pure $`p`$-like state. In that case, only the element $`u_s(𝐑_e)|p_\alpha |u_{p_\alpha }(𝐑_e)=\sqrt{m_0E_p/2}`$ is needed,

In order to obtain a smooth absorption spectrum, we replace the $`\delta `$ function in Eq. (1) by a Lorentzian function with a half-width $`\mathrm{\Gamma }`$,

$$\delta (E_iE)\mathrm{\Gamma }/\{\pi [(E_iE)^2+\mathrm{\Gamma }^2]\}$$
(7)
$`\mathrm{\Gamma }`$ is energy width due to imhomogeneous broadening, which is taken to be 0.01 eV (??).

## III Results And Discussions

We have performed calculations of inter-subband and inter-band absorption spectra for the QWR structure depicted in Fig. 1 with varaying degree of alloy mixing and different lengths of period ($`L`$) in lateral modulation. We find that the inter-subband absorption spectra are sensitive to the length of period ($`L`$), but rather insensitive to the degree of of alloying mixing. Thus, we only present results for the case with moderate alloy mixing, which are characterized by parameters $`b=33\AA `$ and $`x_m=1.0`$. In all the calculations, the bottom layer atoms of QWR’s are bounded by the InP substrate, while the upper layer atoms and GaAS capping layer atoms are allowed to move freely. This strucure is corresponding to the unclamped struture as indicated in reference .

For different period length L of QWR’s, the strain distribtion profiles are qualitatively similar as shown in reference . As L decreases, the hydrostatic strain in rich In region (i.e. right half zone of QWR’s unit) increase, while it decreases in rich Ga region. The bi-axial strain has the opposite change with L. The variation of hydrostatic and bi-axial strains with deducing QWR’s period reflects in the potential profiles as the difference of CB and VB band eage increases, which can be seen in Figure 2.

It can be easily understood that the shear strains increase when L is decuded.

The potential profiles due to strain-induced lateral ordering seen by an electron in two QWR structures considered here ($`L=50a_{[110]}`$ and $`L=40a_{[110]}`$) are shown in Fig. 2. more discussions…

The conduction subband structures for the self-assembled QWRs with alloying mixing ($`x_m=1.0`$ and $`b=8a_{[110]}`$) for ($`L=50a_{[110]}`$ and $`L=40a_{[110]}`$) are shown in Fig. 3. All subband are grouped in pairs with a weak spin splitting (not resolved on the scale shown). For $`L=50a_{[110]}`$, the lowest three pairs of subbands are nearly dispersionless along the $`k_2`$ direction, indicating the effect of strong lateral confinement. The inter-subband transition between the first two pairs give rise to the dominant IR response at photon energy around 60 meV. For $`L=40a_{[110]}`$, only the lowest pair of subbands (CB1) is laterally confined (with a weak $`k_2`$ dispersion). The higher subbands corresponding to laterally unconfined states (but remain confined along the growth axis) and they have large dispersion versus $`k_2`$. We find three pairs of subbands (CB2-CB4) are closely spaced in energy (within 5 meV?). State orgion of degeneracry??

The valence subband structures for the self-assembled QWRs with alloying mixing ($`x_m=1.0`$ and $`b=8a_{[110]}`$) for ($`L=50a_{[110]}`$ and $`L=40a_{[110]}`$) are shown in Fig. 4. more discussions??

A. Inter-subband absorption

Inter-subband absorption spectrum is the most relavent quantity in determining the usefulness of self-assembled QWR’s for application in IR detection. Fig. 5 shows the calculated inter-subband absorption spectra of the self-assembled QWR structure (as depicted in Fig. 1) for three different lengths of period: $`L=72,50`$, and $`40a_{[110]}`$ (approximatley 300 Å, 200 Å, and 160 Å, respectively). In the cacluation, we assume that these QWR structures are n-type doped with linear carrier density around $`1.65\times 10^6cm^1`$ (which corresponds to a Fermi level around 25 meV above the condunction band minimum).

For comparison purposes, we show results for polarization vector along both the  (solid curves) and  directions (dashed curves). The results for \[1$`\overline{1}`$0\] polarization are zero due to the strict translational invariance imposed in our model calculation.

The peak positions for the inter-subband transition with normal incidence (with  polarization) are around 65 meV, 75 meV, and 110 meV for the three cases considered here. All these are within the desirable range of IR detection. As expected, the transition energy increases as the length of period decreases due to the increased degree of lateral confinement. However, the transition energy will saturate at around 110 meV as we further reduce the length of period, since the bound-to-continuum transition is already reached at $`L=40a_{[110]}`$.

The absorption strengths for the first two cases ($`L=72a_{[110]}`$ and $`L=50a_{[110]}`$) are reasonably strong (around 400 $`cm^1`$ and 200 $`cm^1`$, respectively). They both correspond to the bound-to-bound transitions. In contrast, the absorption strength for the third case is somewhat weaker (around 50 $`cm^1`$), since it corresponds to the bound-to-continuum transition. For comparison, the absorption strength for typical III-V QWIPs is around ??

The inter-subband absorption for the  polarization is peaked around ?? meV. The excited state involved in this transition is a quantum confined state due to the Al0.24Ga0.24In0.52 barriers. Thus, it has the same physical origin as the inter-subband transition used in typical QWIP structure. Although this peak is not useful for IR detection with normal incidence, it can be used as the second-color detection if one puts a diffraction grating on the surface as typically done in the fabrication of QWIPs.

B. Inter-band absorption

The inter-band optical transitions are important for the characterization of self-assembled QWR’s, since they are readily observable via the Photoluminescence (PL) or optical transimission experiment. For IR-detector application, they offer another absorption peak at mid IR wavelengths, which can be used together with the inter-subband transitions occured at far IR wavelengths for multi-colored detection. Thus, to understand the full capability of the self-assembled QWR material, we also need to analyze the inter-band absorption.

Fig. 6 shows the squared optical matrix elements versus $`k_2`$ for two self-assembled QWR’s considered in the previous section (with $`L=50`$ and $`40a_{[110]}`$). For the case with $`L=50a_{[110]}`$, the optical matrix elements for both  and \[1$`\overline{1}`$0\] polarizations are strong with a polarization ratio $`P_{[1\overline{1}0]}/P_{[110]}`$ around 2. This is similar to the case with $`L=72a_{[110]}`$ as reported in Ref. xx. For the case with $`L=40a_{[110]}`$, the optical matrix elements for both  and \[1$`\overline{1}`$0\] polarizations are weak. This is due to the fact that the electrons and hole are laterally confined in different regions in the QWR, as already indicated in the potential profile as shown in Fig. 2(b). Thus, the inter-band absorption for this case will be uninteresting.

Fig. 7 shows the inter-band absorption spectra for SILO QWR’s with $`L=72`$ and $`50a_{[110]}`$, including the excitonic effects. The PL properties of the $`L=72a_{[110]}`$ structure with alloying mixing characterized by $`x_m=0.1`$ and $`b=8a_{[110]}`$ has been studied in our previous paper. The QWR structure has a gap around 0.74 eV with a PL polarization ratio ($`P_{[1\overline{1}0]}/P_{[110]}`$) around 3.1. The absorption coefficient for this structure has a peak strength around 250 cm-1. The binding energy for the ground state exciton labeled 1-1 (derived primarily from the top valence subband and the lowest conduction subband) is around 20 meV. Thus, the peak position in the absorption spectrum shifts from 0.76 meV (without the excitonic effect) to 0.74 meV (with the excitonic effect). The exctionic effect also enhances the peak strength from 200 cm-1 to 250 cm-1. The other peak structures (labeled 2-2, 2-3,… etc.) are derived primarily from the transitions between the lower valence subbands to the higher conduction subbands).

For the QWR structure with $`L=50a_{[110]}`$, we obtain similar absorption spectrum with a peak strength around 400 $`cm^1`$(??). The exciton binding is around 40?? meV, and the excitonic enhancement factor of the first peak is around 1.15 (??), higer than the case with $`L=72a_{[110]}`$. This indicates that the case with $`L=50a_{[110]}`$ has stronger lateral confinement for electrons and holes, which leads larger exciton binding energy and stronger excition oscillator strength (due to the larger probability that the electron and hole appear at the same position). The secondary peaks due to excitonic states derived from higher subbands are also subtantially stronger than their counterparts in the $`L=72a_{[110]}`$ case.

## IV Summary and Discussions

We have studied the inter-subband and inter-band absorption spectra for self-assembled InGaAs quantum wires for consideration in IR-detector application. Detailed band structures, microscopic strain distributions, and excitonic effects all have been taken into account. A number of realistic structures grown via strain-induced lateral ordering process are examined. We find that the self-assembled InGaAs quantum wires are good candidate for multi-colored IR detector materials. They offer two groups of strong IR absorption peaks: one in the far-IR range with wavelengths covering 10 - 20 $`\mu m`$ (via the inter-subband transition), the other in the mid-IR range with wavelengths centered around 1.5 $`\mu m`$ (via the inter-band transition). Due the strain induced lateral modulation, the inter-subband transition is strong for normal incident light with polarization along the direction of lateral modulation (). This gives the self-assembled InGaAs quantum wires a distinct advantage over the quantum well systems for application in IR detection.

The inter-subband absorption is found to be sensitive to the length of period ($`L`$) of laterial modulation with the aborption peak position varying from 60 meV to 110 meV as the length of period is reduced from 300 Å to 160 Å. However, further reduction in the length of period does not shift the absorption peak very much, as the excited states become laterally unconfined.

For the inter-band transition, we find that the excitonsic effect enhances the absorption peak strength by about 10-20 %, and shift the peak position by about 20-40 meV for the structures considered. The reduction in the period length ($`L`$) leads to stronger lateral confinement, hence larger exciton binding and stronger absorprtion strength. As conclusion, this paper should give the experiment the realistic guidance in the growth of the IR detector and present the interesting physical thoughts for the theoretists and experimentists.

In conclusion, we successfully demonstrated that self-assembled quantum wires are promising IR-detector materials and we provided theoretical modeling for the optical characteristics for realistic QWR structures, which can be used to guide future fabrication of quantum wire infrared detectors.

Figure Captions

Fig. 1. Schematic sketch of the unit cell of the self-assembled quantum wire for the model structure considered. Each unit cell consists of 8 pairs of (2/2.25) GaAs/InAs short-period superlattices (SPS). In this structure, four pairs of (2/2.25) SPS (or 17 diatomic layers) form a period, and the period is repeated twice in the unit cell. Filled and open circles indicate Ga and In rows (each row extends infinitely along the $`[1\overline{1}0]`$ direction).

Fig. 2. Conduction band and valence band edges for self-assembled QWR structure depicted in Fig. 1 for (a) $`L=50a_{[110]}`$ and (b) $`L=40a_{[110]}`$. Dashed: without alloy mixing. Solid: with alloy mixing described by $`x_m=1.0`$ and $`b=8a_{[110]}`$.

Fig. 3. Conduction subband structure of self-assembled QWR for (a) $`L=50a_{[110]}`$ and (b) $`L=40a_{[110]}`$ with $`x_m=1.0`$ and $`b=8a_{[110]}`$.

Fig. 4. Valence subband structure of self-assembled QWR for (a) $`L=50a_{[110]}`$ and (b) $`L=40a_{[110]}`$ with $`x_m=1.0`$ and $`b=8a_{[110]}`$.

Fig. 5. Inter-subband absorption spectra of self-assembled QWR for (a) $`L=72a_{[110]}`$, (b) $`L=50a_{[110]}`$, and (c) $`L=40a_{[110]}`$ with $`x_m=1.0`$ and $`b=8a_{[110]}`$. Solid:  polarization, dashed :  polarization.

Fig. 6. Inter-band optical matrix elements squared versus $`k_1`$ of self-assembled QWR’s for (a) $`L=50a_{[110]}`$ and (b) $`L=40a_{[110]}`$ with $`x_m=1.0`$ and $`b=8a_{[110]}`$.

Fig. 7. Inter-band absorption spectra of self-assembled QWR’s for (a) $`L=72a_{[110]}`$ and (b) $`L=50a_{[110]}`$ with $`x_m=1.0`$ and $`b=8a_{[110]}`$. Solid:  polarization with excitonic effect. Dotted : \[1$`\overline{1}`$0\] polarization with excitonic effect. Dashed:  polarization without excitonic effect.


text:
# A note on symplectic singularities

Introduction

Let $`X`$ be an affine (or Stein) open subset of an algebraic variety of $`dimn`$ over $``$. $`X`$ is called a symplectic singularity if the smooth locus $`U`$ admits an everywhere non-degenerate holomorphic closed 2-form $`\omega `$ and, for an arbitrary resolution $`f:YX`$ with $`f^1(U)U`$, $`\omega `$ extends to a regular 2-form on $`Y`$. Let $`\mathrm{\Sigma }`$ be the singular locus. In this note we shall prove the following.

Theorem. Let $`X`$ be a symplectic singualrity. Then $`\mathrm{\Sigma }`$ has no codimension 3 irreducible components.

Corollary 1. Let $`X`$ be a symplectic singularity. Then $`X`$ is terminal if and only if $`\mathrm{Codim}(\mathrm{\Sigma }X)4`$ .

Proof of Corollary 1. Assume that $`X`$ has only terminal singularities. Then $`\mathrm{Codim}(\mathrm{\Sigma }X)3`$. By Theorem, $`\mathrm{\Sigma }`$ has no codimension 3 irreducible components. Therefore $`\mathrm{Codim}(\mathrm{\Sigma }X)4`$.

Conversely, assume $`\mathrm{Codim}(\mathrm{\Sigma }X)4`$. Take a resolution $`f:YX`$ of $`X`$ and denote by $`\{E_i\}`$ the $`f`$-exceptional divisors. Since $`X`$ has canonical singularities, we have $`K_Y=f^{}K_X+\mathrm{\Sigma }a_iE_i`$ with non-negative integers $`a_i`$. One can prove that $`a_i`$ are actually all positive by the same argument as §1 111Although it is not explicitly mentioned, the following fact is proved in §1: Let $`X`$ be a symplectic singularity and $`f:YX`$ a resolution of singularities. Assume that $`E_i`$ is an irreducible $`f`$-exceptional divisor such that $`a_i=0`$. Then $`dimf(E_i)=n2`$, where $`n=dimX`$..

We can define a symplectic variety $`Z`$ similarly as a compact normal Kaehler space $`Z`$ whose regular locus $`V`$ admits an everywhere non-degenerate holomorphic closed 2-form which extends to a regular 2-form on an arbitrary resolution $`\stackrel{~}{Z}`$ of $`Z`$. Then the following result is also a corollary of our theorem.

Corollary 2. Let $`Z`$ be a projective symplectic variety with terminal singularities. Then the Kuranishi space $`\mathrm{Def}(Z)`$ is smooth.

Proof of Corollary 2. We already know the result when $`Z`$ is a symplectic variety with $`\mathrm{Codim}(\mathrm{\Sigma }Z)4`$ by \[Na 2, Theorem (2.4)\]. Now the result follows from Corollary 1. Q.E.D.

In the remainder we put $`dimX=n`$. Our strategy is the following. Assume that $`\mathrm{\Sigma }`$ has an $`n3`$ dimensional irreducible component $`\mathrm{\Sigma }^{}`$. Let $`H:=H_1H_2\mathrm{}H_{n3}`$ be a complete intersection of $`(n3)`$ general hyperplane sections of $`X`$. $`H`$ intersects $`\mathrm{\Sigma }^{}`$ in a finite number of points. Pick a point $`p`$ from these points and we shall only consider a small open neighborhood of $`pX`$. For a resolution $`f:YX`$, $`\stackrel{~}{H}:=f^1(H)`$ is a resolution of $`H`$ by the theorem of Bertini. Since $`X`$ has canonical singularities, we have $`K_Y=f^{}K_X+\mathrm{\Sigma }a_iE_i`$ where $`E_i`$ are $`f`$-exceptional divisors and the coefficients $`a_i`$ are all non-negative. By the adjunction formula, $`K_{\stackrel{~}{H}}=f^{}K_H+\mathrm{\Sigma }a_i(E_i\stackrel{~}{H})`$. As $`F_i:=E_i\stackrel{~}{H}`$ are exceptional divisors of $`f|_{\stackrel{~}{H}}`$, we know that $`H`$ has also canonical singularities. By taking an (n-3)-parameter family of complete intersections $`H_s`$ with $`H_0=H`$, we can regard $`X`$, locally around $`p`$, as the total space $``$ of the (n-3)-parameter family. By definition, there is a map $`\pi :\mathrm{\Delta }^{n3}`$ such that each fiber is $`H_s`$.

If the 3-dimensional canonical singularity $`(H,p)`$ is not compound Du Val (cDV) singularity, then there is at least one $`f|_{\stackrel{~}{H}}`$-exceptional divisor $`F_i`$ over $`p`$ such that $`a_i=0`$ by \[Re, Theorem 2.2\]. Therefore, there is at least one $`f`$-exceptional divisor $`E_i`$ with $`a_i=0`$. This $`E_i`$ is mapped onto $`\mathrm{\Sigma }^{}`$ by $`f`$. We shall show that this contradicts the fact that $`X`$ is a symplectic singularity in §1.

Next consider the case where $`(H,p)`$ is a cDV point. $`(,p)`$ is a hypersurface singularity in $`(\mathrm{\Delta }^{n+1},0)`$ defined by $`g(x,y,z,t,s_1,s_2,\mathrm{},s_{n3})=0`$, where $`(x,y,z,t,s_1,s_2,\mathrm{},s_{n3})`$ are coordinates of $`\mathrm{\Delta }^{n+1}`$, and $`\mathrm{\Delta }^{n3}`$ is given by $`(x,y,z,t,s_1,s_2,\mathrm{},s_{n3})(s_1,s_2,\mathrm{},s_{n3})`$. Since $`(,p)`$ is an (n-3)-parameter deformation of an isolated cDV point $`(H,p)`$ and since $`(H,p)`$ is a one parameter deformation of a Du Val singularity $`S`$, $`(,p)`$ is an (n-2)-parameter deformation of a Du Val singularity $`S`$. We may assume that the n-2 parameters are given by $`t`$ and $`s_1,s_2,\mathrm{},s_{n3}`$. Let $`\mathrm{Def}(S)`$ be the semi-universal deformation space (=Kuranishi space) of $`S`$. Then we have a holomorphic map $`\varphi :\mathrm{\Delta }^{n2}\mathrm{Def}(S)`$ such that $`\varphi (0)=0\mathrm{Def}(S)`$. Let $`\omega `$ be a holomorphic 2-form defined on the regular part of $``$. By using the above description of $``$, we shall write $`\omega `$ rather explicitly in terms of $`dx`$, $`dy`$, $`dz`$, $`dt`$, $`ds_1`$, …, $`ds_{n3}`$. This explicit description of $`\omega `$ will tell us that $`^{n/2}\omega `$ is somewhere-vanishing (cf. §2).

Remark. Since a symplectic singularity is a Gorenstein canonical singularity, there is a closed subset $`\mathrm{\Sigma }_0`$ of $`\mathrm{\Sigma }`$ such that each point $`p\mathrm{\Sigma }\mathrm{\Sigma }_0`$ has an open neighborhood (in $`X`$) isomorphic to $`(R.D.P)\times (^{n2},0)`$ (\[Re\]). Then our theorem can be generalized to the following form:

Let $`X`$ be a symplectic singularity. Then $`\mathrm{Codim}(\mathrm{\Sigma }_0X)4`$.

The proof is almost the same as Theorem. When $`X`$ has a symplectic resolution, this result is proved in \[Na 2\] by a different method.

§1. non-cDV case

We use the same notation as Introduction. This section deals with the case $`(H,p)`$ is not a cDV point. As in Introduction, for a resolution $`f:YX`$, there is an $`f`$-exceptional divisor $`E_i`$ with $`a_i=0`$ which dominates $`\mathrm{\Sigma }^{}`$. We put $`\stackrel{~}{}:=f^1()`$. Denote by the same notation $`f`$ (resp. $`E_i`$) the restriction of $`f`$ to $`\stackrel{~}{}`$ (resp. the restriction of $`E_i`$ to $`\stackrel{~}{}`$). We shall derive a contradiction by assuming that there is an everywhere non-degenerate 2-form $`\omega `$ on the smooth locus of $``$.

Since $`a_i=0`$, $`^{n/2}f^{}\omega `$ does not vanish at general point $`qE_i`$. This means that $`f^{}\omega `$ is a non-degenerate 2-form at $`q`$. We can choose a system of local coordinates $`(x_1,x_2,\mathrm{},x_n)`$ at $`q\stackrel{~}{}`$ in such a way that $`x_1=0`$ is a defining equation of $`E_i`$ and $`f^{}\omega (q)H^0(\stackrel{~}{},\mathrm{\Omega }_\stackrel{~}{}^2k(q))`$ is of the following form

$$dx_1dx_2+dx_3dx_4+\mathrm{}+dx_{n1}dx_n.$$

Let $`\omega ^{}H^0(E_i,\mathrm{\Omega }_{E_i}^2)`$ be the restriction of $`f^{}\omega `$ to $`E_i`$. We see by the above expression of $`f^{}\omega (q)`$ that $`^{n/21}\omega ^{}0`$ at $`q`$.

Take a birational proper morphism $`\nu :𝒲\stackrel{~}{}`$ from a smooth variety $`𝒲`$ to $`\stackrel{~}{}`$ in such a way that $`F:=(f\nu )^1(\mathrm{\Sigma }^{})`$ is a simple normal crossing divisor of $`𝒲`$. Denote by $`\widehat{\mathrm{\Omega }}_F^i`$ (resp.$`\widehat{\mathrm{\Omega }}_{F/\mathrm{\Sigma }^{}}^i`$) the sheaf of i-forms of $`F`$ modulo torsion (resp. the sheaf of relative i-forms of $`F`$ over $`\mathrm{\Sigma }^{}`$ modulo torsion). Let $`\omega ^{\prime \prime }H^0(F,\widehat{\mathrm{\Omega }}_F^2)`$ be the restriction of $`(f\nu )^{}\omega `$ to $`F`$. We may assume that $`\nu `$ is an isomorphism at $`q`$, and $`\omega ^{\prime \prime }=\omega ^{}`$ at $`q`$. Therefore, $`^{n/21}\omega ^{\prime \prime }0`$ at $`q`$. We conclude that $`\omega ^{\prime \prime }`$ is not the pull-back of any 2-form of $`\mathrm{\Sigma }^{}`$ by the map $`F\mathrm{\Sigma }^{}`$ because $`dim\mathrm{\Sigma }^{}n3`$. Denote by $`F_x`$ the fiber of the map $`F\mathrm{\Sigma }^{}`$ over $`x\mathrm{\Sigma }^{}`$.

We consider two exact sequences

$$0\widehat{\mathrm{\Omega }}_F^2\widehat{\mathrm{\Omega }}_{F/\mathrm{\Sigma }^{}}^20,$$

$$0(f\nu )^{}\mathrm{\Omega }_\mathrm{\Sigma }^{}^2(f\nu )^{}\mathrm{\Omega }_\mathrm{\Sigma }^{}^1\widehat{\mathrm{\Omega }}_{F/\mathrm{\Sigma }^{}}^10.$$

We know that $`H^0(F_x,\widehat{\mathrm{\Omega }}_{F_x}^2)0`$ or $`H^0(F_x,\widehat{\mathrm{\Omega }}_{F_x}^1)0`$ for a general point $`x\mathrm{\Sigma }^{}`$ by taking global sections of these exact sequences because $`\omega ^{\prime \prime }`$ is not the pull-back of any 2-form of $`\mathrm{\Sigma }^{}`$. By the mixed Hodge structure on $`H^1(F_x)`$ or $`H^2(F_x)`$, $`H^1(F_x,𝒪_{F_x})0`$ if $`H^0(F_x,\widehat{\mathrm{\Omega }}_{F_x}^1)0`$ and $`H^2(F_x,𝒪_{F_x})0`$ if $`H^0(F_x,\widehat{\mathrm{\Omega }}_{F_x}^2)0`$.

On the other hand, let $`_x`$ be a general complete intersection of n-3 general hyperplanes passing through $`x`$. Put $`𝒲_x:=(f\nu )^1(_x)`$. Then the map $`(f\nu )|_{𝒲_x}:𝒲_x_x`$ is a resolution of singularities, and $`F_x`$ is the exceptional locus of this resolution. Since $`_x`$ has only rational singularities, $`R^i((f\nu )|_{𝒲_x})_{}𝒪_{𝒲_x}=0`$ for all $`i>0`$. We deduce from this fact that $`H^i(F_x,𝒪_{F_x})=0`$ for all $`i>0`$. This is a contradiction.

§2. cDV case

We use the same notation as Introduction. This section deals with the case $`(H,p)`$ is a cDV point. Let $`S`$ be a general hyperplane section of $`H`$ passing through $`p`$. If $`S`$ is of type $`A_n`$ (resp. $`D_n`$, $`E_n`$), then we call $`(H,p)`$ is of type $`cA_n`$ (resp. $`cD_n`$, $`cE_n`$). As in Introduction, $`(,p)`$ is a hypersurface singularity defined by $`g(x,y,z,t,s_1,\mathrm{},s_{n3})=0`$ in $`(\mathrm{\Delta }^{n+1},0)`$, where $`(x,y,z,t,s_1,s_2,\mathrm{},s_{n3})`$ are coordinates of $`\mathrm{\Delta }^{n+1}`$. The projection $`\mathrm{\Delta }^{n+1}\mathrm{\Delta }^{n2}`$ ($`(x,y,z,t,s_1,s_2,\mathrm{},s_{n3})(t,s_1,s_2,\mathrm{},s_{n3})`$ induces a Du Val fibration $`\mathrm{\Delta }^{n2}`$. $`S`$ is the central fiber of this fibration. Moreover, the projection $`\mathrm{\Delta }^{n+1}\mathrm{\Delta }^{n3}`$ ($`(x,y,z,t,s_1,s_2,\mathrm{},s_{n3})(s_1,s_2,\mathrm{},s_{n3})`$) induces a cDV fibration $`\mathrm{\Delta }^{n3}`$. The central fiber of this fibration is $`H`$. We may assume that each fiber $`H_s`$ of this cDV fibration has a unique singular point at $`(0,0,0,0,s_1,\mathrm{}s_{n3})`$ and its type as a cDV point is independent of $`(s_1,\mathrm{},s_{n3})`$ (cf. Introduction). The Du Val fibration $`\mathrm{\Delta }^{n2}`$ induces a holomorphic map $`\varphi :\mathrm{\Delta }^{n2}\mathrm{Def}(S)`$ such that $`\varphi (0)=0\mathrm{Def}(S)`$. $``$ is isomorphic to the family of Du Val singularities over $`\mathrm{\Delta }^{n2}`$ obtained by pulling back the semi-universal family over $`\mathrm{Def}(S)`$ by $`\varphi `$. Therefore $`g(x,y,z,t,s_1,s_2,\mathrm{},s_{n3})`$ can be written in the following way according to the type of the Du Val singularity $`S`$ (cf. \[Lo\]):

$`A_n`$: $`x^2+y^2+z^{n+1}+\varphi _1(t,s)z^{n1}+\mathrm{}+\varphi _{n1}(t,s)z+\varphi _n(t,s)`$

$`D_n(n4)`$: $`x^2+y^2z+z^{n1}+\varphi _1(t,s)z^{n2}+\mathrm{}+\varphi _{n1}(t,s)+\varphi _n(t,s)y`$

$`E_6`$: $`x^2+y^3+z^4+\varphi _1(t,s)z^2+\mathrm{}+\varphi _3(t,s)+y(\varphi _4(t,s)z^2+\mathrm{}+\varphi _6(t,s))`$

$`E_7`$: $`x^2+y^3+yz^3+y(\varphi _1(t,s)z+\varphi _2(t,s))+\varphi _3(t,s)z^4+\mathrm{}+\varphi _7(t,s)`$

$`E_8`$: $`x^2+y^3+z^5+y(\varphi _1(t,s)z^3+\mathrm{}+\varphi _4(t,s))+\varphi _5(t,s)z^3+\mathrm{}+\varphi _8(t,s)`$

We abbreviated $`s`$ for $`s_1,\mathrm{},s_{n3}`$. The $`\varphi _i(t,s)`$’s satisfy the following properties, which will be used below:

(a) $`\varphi _i(0,s)0`$;

(b) $`\varphi _i(t,s)t^2\{t,s_1,\mathrm{},s_{n3}\}`$ except the following cases

$`\varphi _1(t,s)`$ in the case $`(D_n)`$;

$`\varphi _4(t,s)`$ in the case $`(E_6)`$;

$`\varphi _3(t,s)`$ in the case $`(E_7)`$;

$`\varphi _1(t,s)`$ in the case $`(E_8)`$.

(a) is satisfied because the fibers of the Du Val fibration $`\mathrm{\Delta }^{n2}`$ over $`(0,s_1,\mathrm{},s_{n3})\mathrm{\Delta }^{n2}`$ are Du Val singularities of the same type. We can check (b) in the following way. Assume that some other $`\varphi _i(t,s)`$’s are not contained in $`t^2\{t,s_1,\mathrm{},s_{n3}\}`$. Let $`H_s`$ be the fiber of the cDV fibration $`\mathrm{\Delta }^{n3}`$ over $`s=(s_1,\mathrm{},s_{n3})`$. Take $`s`$ general. Then the type of a cDV point $`H_s`$ does not coincide with the type of a Du Val singularity $`S`$. This contradicts our assumption.

Let $`^0`$ and $`S^0`$ be the smooth locus of $``$ and $`S`$ respectively. Denote by $`p`$ the Du Val fibration $`\mathrm{\Delta }^{n2}`$. We write $`p^0`$ for the restriction of $`p`$ to $`^0`$. We have two exact sequences

$$0\mathrm{\Omega }_^0^2\mathrm{\Omega }_{^0/\mathrm{\Delta }^{n2}}^20,$$

$$0(p^0)^{}\mathrm{\Omega }_{\mathrm{\Delta }^{n2}}^2\mathrm{\Omega }_{^0/\mathrm{\Delta }^{n2}}^1(p^0)^{}\mathrm{\Omega }_{\mathrm{\Delta }^{n2}}^10.$$

Let $`\omega `$ be a holomorphic 2-form defined on $`^0`$, that is, $`\omega H^0(^0,\mathrm{\Omega }_^0^2)`$. Let $`\omega |_{S^0}H^0(S^0,\mathrm{\Omega }_^0^2|_{S^0})`$ be the image of $`\omega `$ by the natural map $`H^0(^0,\mathrm{\Omega }_^0^2)H^0(S^0,\mathrm{\Omega }_^0^2|_{S^0})`$. We shall write $`\omega |_{S^0}`$ explicitly and compute $`^{n/2}\omega |_{S^0}`$.

By restricting above exact sequences to $`S^0`$ and taking global sections, we have the exact sequences

$$0H^0(S^0,|_{S^0})H^0(S^0,\mathrm{\Omega }_^0^2|_{S^0})\stackrel{\alpha }{}H^0(S^0,\mathrm{\Omega }_{S^0}^2),$$

$$0H^0(S^0,(p^0)^{}\mathrm{\Omega }_{\mathrm{\Delta }^{n2}}^2|_{S^0})H^0(S^0,F|_{S^0})H^0(S^0,\mathrm{\Omega }_{S^0}^1(p^0)^{}\mathrm{\Omega }_{\mathrm{\Delta }^{n2}}^1).$$

As the first step we shall construct an explicit lift $`\omega _1H^0(S^0,\mathrm{\Omega }_^0^2|_{S^0})`$ of $`\alpha (\omega |_{S^0})H^0(S^0,\mathrm{\Omega }_{S^0}^2)`$.

Step 1. Consider the sheaf $`\mathrm{\Omega }_{^0/\mathrm{\Delta }^{n3}}^2`$ of relative 2-forms with respect to the cDV fibration $`\mathrm{\Delta }^{n3}`$. By a generator $`dydzdt/(g/x)`$ of the relative dualizing sheaf $`\omega _{/\mathrm{\Delta }^{n3}}`$, the relative tangent sheaf $`\mathrm{\Theta }_{^0/\mathrm{\Delta }^{n3}}`$ and $`\mathrm{\Omega }_{^0/\mathrm{\Delta }^{n3}}^2`$ are identified. Thus we have $`H^0(^0,\mathrm{\Theta }_{^0/\mathrm{\Delta }^{n3}})H^0(^0,\mathrm{\Omega }_{^0/\mathrm{\Delta }^{n3}}^2)`$. An element of $`H^0(^0,\mathrm{\Theta }_{^0/\mathrm{\Delta }^{n3}})`$ can be expressed as

$$f_1/x+f_2/y+f_3/z+f_4/t$$

with holomorphic functions $`f_i`$ on $``$ which satisfy $`f_1g/x+f_2g/y+f_3g/z+f_4g/t=0`$ on $``$. Therefore an element of $`H^0(^0,\mathrm{\Omega }_{^0/\mathrm{\Delta }^{n3}}^2)`$ can be expressed as

$$f_2dzdt/(g/x)f_3dydt/(g/x)+f_4dydz/(g/x).$$

Lemma (2.1). $`f_i`$ vanish at the origin for $`i=1,2,3,4`$.

Proof. We restrict the equation $`f_1g/x+f_2g/y+f_3g/z+f_4g/t=0`$ on $``$ to $`H:=\{s_1=\mathrm{}=s_{n3}=0\}`$. Note that $`g(x,y,z,t,s)/x|_{s_1=\mathrm{}=s_{n3}=0}=g(x,y,z,t,0)/x`$, and similar facts hold for the partial derivative of $`y`$, $`z`$ and $`t`$. $`H`$ is an isolated hypersurface singularity of $`^4(x,y,z,t)`$ defined by $`g(x,y,z,t,0)=0`$. Assume that some $`f_i`$ is unit. Then, applying Theorem of Frobenius, we can write on $`^4`$

$$f_1/x+f_2/y+f_3/z+f_4/t=(unit)/x^{}$$

by a suitable coordinates change of $`(x,y,z,t)`$ to $`(x^{},y^{},z^{},t^{})`$. Therefore, in $`\{x^{},y^{},z^{},t^{}\}`$

$$g/x^{}+hg=0$$

with some $`h\{x^{},y^{},z^{},t^{}\}`$. Put $`g^{}:=gexp(h𝑑x^{})`$. Then we have $`g^{}/x^{}=0`$. Since $`g^{}=0`$ is a defining equation of $`H`$ in $`^4(x^{},y^{},z^{},t^{})`$, this means that $`H`$ has non-isolated singularities or $`H`$ is smooth, which is a contradiction. Q.E.D.

In particular, the image of $`\omega `$ by the map $`H^0(^0,\mathrm{\Omega }_^0^2)H^0(^0,\mathrm{\Omega }_{^0/\mathrm{\Delta }^{n3}}^2)`$ has such an expression. We take $`f_1,\mathrm{},f_4`$ in such a way that the above 2-form is this image. Since $`\mathrm{\Omega }_^0^2\mathrm{\Omega }_{^0/\mathrm{\Delta }^{n2}}^2`$ factors through $`\mathrm{\Omega }_{^0/\mathrm{\Delta }^{n3}}^2`$, $`f_4dydz/(g/x)H^0(^0,\mathrm{\Omega }_{^0/\mathrm{\Delta }^{n2}}^2)`$ coincides with the image of $`\omega `$ by the map $`H^0(^0,\mathrm{\Omega }_^0^2)H^0(^0,\mathrm{\Omega }_{^0/\mathrm{\Delta }^{n2}}^2)`$. Therefore $`\alpha (\omega |_{S^0})=\overline{f}_4dydz/(g/x)H^0(S^0,\mathrm{\Omega }_{S^0}^2)`$, where $`\overline{f}_4`$ is the restriction of $`f_4`$ to $`S^0`$.

Now the 2-form $`\overline{f}_2dzdt/(g/x)\overline{f}_3dydt/(g/x)+\overline{f}_4dydz/(g/x)`$ can be considered as an element of $`H^0(S^0,\mathrm{\Omega }_^0^2|_{S^0})`$ because $`(g/s_i)|_{t=0}0`$ by the property (a) above222This is the reason why we consider $`\omega |_{S^0}`$ instead of $`\omega `$. In general, $`f_2dzdt/(g/x)f_3dydt/(g/x)+f_4dydz/(g/x)`$ is not an element of $`H^0(^0,\mathrm{\Omega }_^0^2)`$.. So we put

$$\omega _1:=\overline{f}_2dzdt/(g/x)\overline{f}_3dydt/(g/x)+\overline{f}_4dydz/(g/x).$$

Step 2. Put $`\omega _2:=\omega |_{S^0}\omega _1`$. By the construction of $`\omega _1`$, $`\omega _2H^0(S^0,|_{S^0})`$. We shall write $`\omega _2`$ explicitly in terms of $`dx`$, $`dy`$, $`dz`$, $`dt`$, $`ds_1`$, …, $`ds_{n3}`$.

By a generator $`dydz/(g/x)`$ of the dualizing sheaf $`\omega _S`$, $`\mathrm{\Theta }_{S^0}`$ and $`\mathrm{\Omega }_{S^0}^1`$ are identified. Therefore $`H^0(S^0,\mathrm{\Theta }_{S^0})H^0(S^0,\mathrm{\Omega }_{S^0}^1)`$. An element of $`H^0(S^0,\mathrm{\Theta }_{S^0})`$ can be expressed as

$$h_1/x+h_2/y+h_3/z$$

with holomorphic functions $`h_i`$ on $`S`$ which satisfy $`h_1g/x+h_2g/y+h_3g/z=0`$ on $`S`$. Therefore, an element of $`H^0(S^0,\mathrm{\Omega }_{S^0}^1)`$ can be expressed as

$$h_2dz/(g/x)h_3dy/(g/x).$$

Lemma (2.2) $`h_i`$ vanish at the origin for $`i=1,2,3`$.

Proof. If some $`h_i`$ is a unit, then $`S`$ should have non-isolated singularities or $`S`$ should be smooth, by the same argument as the proof of Lemma (2.1). Q.E.D.

The image of $`\omega _2`$ by the map $`H^0(S^0,F|_{S^0})H^0(S^0,\mathrm{\Omega }_{S^0}^1(p^0)^{}\mathrm{\Omega }_{\mathrm{\Delta }^{n2}}^1)`$ has the expression

$$\mathrm{\Sigma }_{1jn3}\{h_2^{(j)}dz/(g/x)h_3^{(j)}dy/(g/x)\}ds_j+\{h_2^{}dz/(g/x)h_3^{}dy/(g/x)\}dt$$

for some $`h_i^{(j)}`$ and $`h_i^{}`$ which satify $`h_1^{(j)}g/x+h_2^{(j)}g/y+h_3^{(j)}g/z=0`$ and $`h_1^{}g/x+h_2^{}g/y+h_3^{}g/z=0`$. Let $`D`$ be the divisor of $`S`$ defined by $`g/x=0`$. Then the above 2-form can be thought of as an element of $`H^0(S^0,|_{S^0}(D))`$. Now let us consider $`\omega _2`$ as an element of $`H^0(S^0,|_{S^0}(D))`$ by the injection $`H^0(S^0,|_{S^0})H^0(S^0,|_{S^0}(D))`$. We put

$$\omega _3:=\omega _2\mathrm{\Sigma }_{1jn3}(h_2^{(j)}dz/(g/x)h_3^{(j)}dy/(g/x))ds_j\{h_2^{}dz/(g/x)h_3^{}dy/(g/x\}dt.$$

By the exact sequence

$$0H^0(S^0,(p^0)^{}\mathrm{\Omega }_{\mathrm{\Delta }^{n2}}^2|_{S^0}(D))H^0(S^0,F|_{S^0}(D))H^0(S^0,\mathrm{\Omega }_{S^0}^1(D)(p^0)^{}\mathrm{\Omega }_{\mathrm{\Delta }^{n2}}^1)$$

we have $`\omega _3H^0(S^0,(p^0)^{}\mathrm{\Omega }_{\mathrm{\Delta }^{n2}}^2|_{S^0}(D))`$. We write

$$\omega _3=\mathrm{\Sigma }_{1jn3}k^{(j)}ds_jdt/(g/x)+\mathrm{\Sigma }_{1l,mn3}k^{(l,m)}ds_lds_m/g/x)$$

with holomorphic functions $`h^{(j)}`$, $`h^{(l,m)}`$ on $`S^0`$. We shall check the condition that $`\omega _2H^0(S^0,|_{S^0})`$. In particular, $`\omega _2\mathrm{\Gamma }(S^0\{g/y=0\},|_{S^0})`$ and $`\omega _2\mathrm{\Gamma }(S^0\{g/z=0\},|_{S^0})`$. In order to check these conditions, we substitute $`dy`$ by

$$1/(g/y)\{g/xdx+g/zdz+g/tdt\}$$

or substitute $`dz`$ by

$$1/(g/z)\{g/xdx+g/ydy+g/tdt\}.$$

Note that the terms of $`ds_j`$ do not appear because $`(g/s_i)|_{t=0}0`$ by the property (a).

By using the equations $`h_1^{(j)}g/x+h_2^{(j)}g/y+h_3^{(j)}g/z=0`$ and $`h_1^{}g/x+h_2^{}g/y+h_3^{}g/z=0`$, we finally have the following two expressions of $`\omega _2`$:

(y): $`\omega _2=\mathrm{\Sigma }_{1jn3}\{h_3^{(j)}dx/(g/y)h_1^{(j)}dz/(g/y)\}ds_j`$

$`+\{h_3^{}dx/(g/y)h_1^{}dz/(g/y)\}dt`$

$`+\mathrm{\Sigma }_{1jn3}\{k^{(j)}/(g/x)h_3^{(j)}(g/t)/(g/x)(g/y)\}ds_jdt`$

$`+\mathrm{\Sigma }_{1l,mn3}k^{(l,m)}ds_lds_m/(g/x)`$.

(z): $`\omega _2=\mathrm{\Sigma }_{1jn3}\{h_2^{(j)}dx/(g/z)+h_1^{(j)}dy/(g/z)\}ds_j`$

$`+\{h_2^{}dx/(g/z)+h_1^{}dy/(g/z)\}dt`$

$`+\mathrm{\Sigma }_{1jn3}\{k^{(j)}/(g/x)+h_2^{(j)}(g/t)/(g/x)(g/z)\}ds_jdt`$

$`+\mathrm{\Sigma }_{1l,mn3}k^{(l,m)}ds_lds_m/(g/x)`$.

Lemma (2.3). $`k^{(j)}(x)+m^2𝒪_{S,0}`$ $`(1jn3)`$, and $`k^{(l,m)}x𝒪_{S,0}`$ $`(1l,mn3)`$.

Proof. $`k^{(l,m)}/(g/x)`$ should be regular on $`S^0`$. This implies that $`k^{(l,m)}x𝒪_{S,0}`$ because $`g/x=2x`$. We shall prove that $`k^{(j)}m𝒪_{S,0}`$ by using property (b) of $`g(x,y,z,s,t)`$ in each case $`A_n`$, $`D_n`$, $`E_6`$, $`E_7`$ or $`E_8`$.

$`(A_n)`$: Since $`g/t=0`$, $`k^{(j)}x𝒪_{S,0}`$ by the same reason as $`h^{(l,m)}`$ are so.

$`(D_n)`$: $`\omega |_{S^0}`$ should be regular on $`S\{g/z=0\}`$. In the expression (z), we denote by $`G^{(j)}`$ the coefficients of $`ds_jdt`$. $`G^{(j)}`$ is regular on $`S\{g/z=0\}`$. In our case, $`g/x=2x`$, $`g/z=y^2+(n1)z^{n2}`$ and $`g/t=cz^{n2}`$ (c: constant). Since it is easily checked that $`\{g/x=0\}`$ and $`\{g/z=0\}`$ have no common irreducible components on $`S`$, this implies that $`k^{(j)}(g/z)+h_2^{(j)}(g/t)(g/x)𝒪_{S,0}`$.

More explicitly, we have $`\{y^2+(n1)z^{n2}\}k^{(j)}+cz^{n2}h_2^{(j)}x𝒪_{S,0}`$. We see at first that $`k^{(j)}m𝒪_{S,0}`$ because $`y^2`$ does not appear in the quadratic part of $`g(x,y,z,0,0)=x^2+y^2z+z^{n1}`$.

Next the linear part of $`k^{(j)}`$ does not contain $`y`$ because any element of the form $`\alpha y^3`$ \+ (other terms), with $`\alpha ^{}`$, is not contained in the ideal $`(x,z^{n2},x^2+y^2z+z^{n1})`$ of $`\{x,y,z\}`$.

Finally the linear part of $`k^{(j)}`$ does not contain $`z`$ because $`z`$ does not appear in the linear part of $`h_2^{(j)}`$ (by an explicit calculation using the equation $`h_1^{(j)}g/x+h_2^{(j)}g/y+h_3^{(j)}g/z=0`$ in $`𝒪_{S,0}`$) and therefore any element of the form $`\alpha \{y^2z+(n1)z^{n1}\}`$ \+ (other terms), with $`\alpha ^{}`$, is not contained in the ideal $`(x^2+y^2z+z^{n1},x,h_2^{(j)}z^{n2})`$ of $`\{x,y,z\}`$.

$`(E_6)`$: $`\omega |_{S^0}`$ should be regular on $`S\{g/y=0\}`$. In the expression (y), we denote by $`F^{(j)}`$ the coefficients of $`ds_jdt`$. $`F^{(j)}`$ is regular on $`S\{g/y=0\}`$. In our case, $`g/x=2x`$, $`g/y=3y^2`$ and $`g/t=cyz^2`$ (c: constant). Since it is easily checked that $`\{g/x=0\}`$ and $`\{g/y=0\}`$ have no common irreducible components on $`S`$, this implies that $`k^{(j)}(g/y)h_3^{(j)}(g/t)(g/x)𝒪_{S,0}`$. More explicitly, we have $`3yk^{(j)}cz^2h_3^{(j)}x𝒪_{S,0}`$. We have $`k^{(j)}m𝒪_{S,0}`$ because $`g(x,y,z,0,0)`$ has no linear terms. Moreover, $`y`$ and $`z`$ do not appear in the linear part of $`k^{(j)}`$ because any element of the form $`\alpha y^2`$ \+ (other terms) or of the form $`\alpha yz`$ \+ (other terms), $`(\alpha ^{})`$, is not contained in the ideal $`(x,x^2+y^3+z^4,z^2)`$ of $`\{x,y,z\}`$.

$`(E_7)`$: $`\omega |_{S^0}`$ should be regular on $`S\{g/y=0\}`$. In the expression (y), we denote by $`F^{(j)}`$ the coefficients of $`ds_jdt`$. $`F^{(j)}`$ is regular on $`S\{g/y=0\}`$. In our case, $`g/x=2x`$, $`g/y=3y^2+z^3`$ and $`g/t=cz^4`$ (c: constant). Since it is easily checked that $`\{g/x=0\}`$ and $`\{g/y=0\}`$ have no common irreducible components on $`S`$, this implies that $`k^{(j)}(g/y)h_3^{(j)}(g/t)(g/x)𝒪_{S,0}`$. More explicitly, we have $`(3y^2+z^3)k^{(j)}cz^4h_3^{(j)}x𝒪_{S,0}`$. We have $`k^{(j)}m𝒪_{S,0}`$ because the term $`y^2`$ is not contained in the quadratic part of $`g(x,y,z,0,0)`$. Moreover, $`y`$ and $`z`$ do not appear in the linear part of $`k^{(j)}`$ because any element of the form $`\alpha (3y^3+yz^3)`$ \+ (other terms) or of the form $`\alpha y^2z`$ \+ (other terms), $`(\alpha ^{})`$, is not contained in the ideal $`(x,h_3^{(j)}z^4,x^2+y^3+yz^3)`$ of $`\{x,y,z\}`$.

$`(E_8)`$: $`\omega |_{S^0}`$ should be regular on $`S\{g/y=0\}`$. In the expression (y), we denote by $`F^{(j)}`$ the coefficients of $`ds_jdt`$. $`F^{(j)}`$ is regular on $`S\{g/y=0\}`$. In our case, $`g/x=2x`$, $`g/y=3y^2`$ and $`g/t=cyz^3`$ (c: constant). Since it is easily checked that $`\{g/x=0\}`$ and $`\{g/y=0\}`$ have no common irreducible components on $`S`$, this implies that $`k^{(j)}(g/y)h_3^{(j)}(g/t)(g/x)𝒪_{S,0}`$. More explicitly, we have $`3yk^{(j)}cz^3h_3^{(j)}x𝒪_{S,0}`$. We have $`k^{(j)}m𝒪_{S,0}`$ because $`g(x,y,z,0,0)`$ has no linear terms. Moreover, $`y`$ and $`z`$ do not appear in the linear part of $`k^{(j)}`$ because any element of the form $`\alpha y^2`$ \+ (other terms) or of the form $`\alpha yz`$ \+ (other terms), $`(\alpha ^{})`$, is not contaied in the ideal $`(x,z^3,x^2+y^3+z^5)`$ of $`\{x,y,z\}`$.

Step 3. We shall write $`\omega |_{S^0}`$ explicitly and calculate $`^{n/2}\omega |_{S^0}`$. Summing up $`\omega _1`$ and $`\omega _2`$, we have

$`\omega |_{S^0}=\{\overline{f}_2+h_2^{}\}dzdt/(g/x)\{\overline{f}_3+h_3^{}\}dydt/(g/x)`$

$`+\overline{f}_4dydz/(g/x)`$

$`+\mathrm{\Sigma }_{1jn3}\{h_2^{(j)}dz/(g/x)h_3^{(j)}dy/(g/x)\}ds_j`$

$`+\mathrm{\Sigma }_{1jn3}k^{(j)}ds_jdt/(g/x)+\mathrm{\Sigma }_{1l,mn3}k_{}^{}{}_{}{}^{(l,m)}ds_lds_m.`$

By Lemmas (2.1), (2.2) and (2.3) the functions $`\overline{f}_i`$, $`h_i^{(j)}`$, $`h_i^{}`$ are all contained in $`m𝒪_{S,0}`$, and $`k^{(j)}(x)+m^2𝒪_{S,0}`$. By Lemma (2.3), $`k_{}^{}{}_{}{}^{(l,m)}:=k^{(l,m)}/(g/x)`$ is regular.

An explicit calculation shows that

$$^{n/2}\omega |_{S^0}=Bdydzdtds_1\mathrm{}ds_{n3}$$

where $`B`$ is of the following form

$`(\overline{f}_2+h_2^{})/(g/x)\mathrm{\Sigma }_{1jn3}\{h_3^{(j)}/(g/x)`$(regular function)}

$`+(\overline{f}_3+h_3^{})/(g/x)\mathrm{\Sigma }_{1jn3}\{h_2^{(j)}/(g/x)`$(regular function)}

$`+\overline{f}_4/(g/x)\mathrm{\Sigma }_{1jn3}\{k^{(j)}/(g/x)`$(regular function)}

$`+\mathrm{\Sigma }_{1p,q,rn3}h_2^{(p)}/(g/x)h_3^{(q)}/(g/x)k^{(r)}/(g/x)`$(regular function).

When $`S`$ is of type $`A_n`$, then $`k^{(p)}x𝒪_{S,0}`$ by the proof of Lemma (2.3). Hence $`B=C/x^2`$ for a suitable $`Cm^2`$, where $`m`$ is the maximal ideal of $`𝒪_{S,0}`$. If $`^{n/2}\omega |_{S^0}`$ is nowhere-vanishing (i.e. $`^{n/2}\omega |_{S^0}`$ is a generator of $`\omega _{}|_S`$), then we should have $`B=(unit)/x`$. This means that $`C=(unit)x`$, hence $`xm^2`$, which is absurd.

When $`S`$ is one of other types. Then $`B=C/x^3`$ for a suitable $`Cxm^2+m^4`$. If $`\omega |_{S^0}`$ is nowhere-vanishing, then we should have $`B=(unit)/x`$. This means that $`C=(unit)x^2`$, hence $`x^2xm^2+m^4`$. However this is absurd because $`g(x,y,z,0,0)`$ always contains a cubic term which cannot be divided by $`x`$.

Department of Mathematics, Graduate school of Sience, Osaka University, Toyonaka 560, Japan


text:
# 1 Left plot–The smooth spectral shape (for 𝛾=5, large angle scattering) were the plateau structures slighlty appear near the top side of the spectrum. The plateau like stucture (due to consecutive cycling), and the spectral index is in consistency with similar computational work presented by Baring (1999). Right plot–The spectral shape for the same gamma but for pitch angle diffusion. The same plateau-like spectrum structure hardly appears.

Particle Spectra Of Ultra-Relativistic Shock Waves111Some of these results were presented in the 17th European Cosmic Ray Symposium, July 2000, Lodz Poland

Athina Meli and John Quenby

Astrophysics Group, Blackett Laboratory Imperial College of Science,Technology and Medicine, London, UK

## Abstract

Numerical Monte Carlo simulations of the diffusive shock acceleration in the test particle limit are investigated. We simulate high relativistic flow astrophysical plasmas for upstream $`\gamma `$ $`5`$ and up to $`\gamma `$ $`1000`$. These gamma values are relevant to the diffusive shock acceleration models of Active Galactic Nuclei (AGN) and Gamma Ray Bursts (GRBs). The spectral shape for the above high relativistic processes will be calculated. This work presents distinctive spectral shapes with plateau-like characteristics due to the relativistic effects and the spectral indexes are found to be as flat as the flow becomes relativistic; the latter is in a very good agreement with previous theoretical and numerical works concerning relativistic plasma flows.

Introduction

It is known that multiwavelength observations from both Gamma Ray Bursts (GRB) (Meszaros & Rees, 1993) and a certain number of Active Galactic Nuclei (AGN), indicate as well as invoke, the presence of relativistic and ultrarelativistic shock waves. In these shocks the mechanism of diffusive shock acceleration at collisionless plasma shocks takes place (Peacock 1981) and explains the appearance of the non-thermal particle populations from these sources. Specifically, authors such as Vietri (1995) and Waxman (1995) have pointed out that these extremely energetic events with hyper-relativistic flows ($`\mathrm{\Gamma }`$ $`1001000`$) in GRBs, could efficiently accelerate protons which could be the source of the Ultra High Energy Cosmic Rays (UHECR). It is an important fact that the analytical models of Vietri and results on UHECR, depend critically on past (Lieu & Quenby, 1989) pioneering large angle scattering numerical and analytical simulations of particle diffusion shock acceleration for $`\mathrm{\Gamma }`$ flows $`3`$ where they reported a spectral flattening when the flow becomes relativistic (Kirk & Schneider, 1987 first noticed a spectral flattening). For $`U_1`$=0.96c, $`U_2`$=0.32c (that is compression ratio = 3) they find the spectral index to be flattened to around 1.2 in excellent agreement with theoretical and numerical work of Kirk and Schneider (1987a,b). The numerical results obtained by our current Monte Carlo simulations in the regime of diffusive shock acceleration for $`\mathrm{\Gamma }`$$`510^3`$, are consistent again with the most plausible predictions regarding the flatteness of the spectrum. The first order Fermi acceleration theory for parallel non-relativistic flows, yields a spectral index $``$2 for the differential spectrum dn/dE$``$$`E^\alpha `$, where $`\alpha =r+2/r+1`$ and $`r=u_1/u_2`$. In the relativistic and ultra-relativistic limits for diffusive acceleration many studies investigated the issue of the spectrum index’s behaviour in relation to flow velocities and magnetic field irreqularities present in the relativistic plasma flow and do find a considerable spectral flattening as the plasma flow becomes relativistic; but the spectral shape at the source is not yet clearly understood. The aim of this work is to address the problem of the spectral shape at source for both large angle scattering and pitch angle diffusion models in order to investigate the variations between the two. Monte Carlo numerical simulations have been used in the test-particle diffusive limit employing highly rleativistic flows. We compare our results with similar studies for $`\gamma `$$``$5 (Baring, 1999) and find a good agreement concerning the spectral index and the spectral shape, for both large angle scattering and pitch angle diffusion mechanisms.

Monte Carlo Method

In our Monte Carlo code, the particle splitting mechanism is used in order to improve the statistics of the number of the remaining particles in the acceleration process. The injection of energetic seed particles ($``$$`10^7`$), with weight equal to 1.0, is far upstream. Particles are allowed to scatter randomly in the respective fluid frames towards the shock. A guiding centre approximation is used and the particle propagation is followed in a 1D geometry. Each particle trajectory is been followed across the shock according to known ’jump conditions’ and made to leave the system from the moment that it escapes far downstream at the spatial boundary or if it reaches a well defined maximum energy. When a number of particles escapes through the above mentioned boundaries, we replace them with new ones -in order to have again the same initial number of particles- but in such a way that as their weight is decreased, the number of particles remains almost constant throughout the simulation. Lorentz transformations are used between frames on shock passage. The results then are recorded just downstream, in the shock frame. For these simulations, we assume the shocks to be ”parallel”, either because that is the field configuration, or turbulence removes ”reflection” at the interface. We keep a compression ratio of 4 and consider both isotropic and pitch angle diffusion in the respective plasma rest frames. The probability that a particle will move a distance z along the field lines at pitch angle $`\theta `$ before a scattering is given by the following expression:

$$Prob\left(z\right)exp\left(z/\lambda \left|cos\theta \right|\right)$$
(1)
The new pitch angle $`\theta ^{}`$, for pitch angle diffusion, is calculated by the simple trigonometric formula:

$$cos\theta ^{}=cos\theta \sqrt{1sin^2\delta \theta }+sin\delta \theta \sqrt{1cos^2\theta }cos\varphi $$
(2)
where $`\varphi `$(0,2$`\pi `$) is the azimuth angle with respect to the original momentum direction. We note here that in the pitch angle diffusion case we keep the scattering within the 1/$`\gamma _{up}`$, where $`\gamma _{up}`$ is the upstream gamma flow measured in the shock frame, in order to directly compare our results with theoretical investigations and predictions.

Results

The primary purpose of this work is to investigate the raw spectral shape and consequently the behaviour of particles being accelerated in ultra-relativistic shocks. These simulations are expected to be relevant within the context of a highly relativistic shock propagation which could be found in GRBs. All results have been recorded in the downstream side, in the shock frame. In the figures shown, there are results for isotropic scattering and pitch angle diffusion for a various number of upstream $`\gamma `$. In figures 1, 2, 3 and 4 we can observe the very characteristic spectral shapes for large angle scattering and pitch angle diffusion. A distinctive feature is, the plateau-like structure for both large angle scattering and pitch angle diffusion mechanisms, at the left side and right side plots respectively, for Lorentz $`\gamma `$=5-990. It can be seen how the relativistic effects at each discrete acceleration cycle is evident, for both mechanisms.

The striking effect though is that the spectra do develop a smoother power-law shape in parts which scales with the energy gain of the particles through the acceleration process. By observing closely the spectral shapes for large angle scattering and pitch angle diffusion it is seen that in the case of pitch angle dissusion the plateau-like parts of the spectrun -which correspond to the particle’s shock cycling- is steeper than the plateau structures of the large angle scattering mechanism. This implies that in pitch angle diffusion case, there is an indication showing that the probability of convection downstream of the shock drops much more with increasing energy than in the case of large angle scattering diffusion. This is a fact that has been shown theoretically in the relativistic regime, where for parallel shocks the large angle scattering gives flatter spectra than the pitch angle diffusion. So, one could infer that at relativistic flow velocities is difficult to assume or take for granted that the spectral shape of the accelerated particles follows a smooth power-law shape although there is a characteristic monotonic increase in energy which is of course a signature of the 1st order Fermi acceleration. We can also understand the steepening behaviour partly because our calculations show that about $`90\%`$ of all particles in the simulation box are lost downstream each cycle. We would like to note here that for the $`\gamma =5`$ for large angle scattering our spectral index has almost the same value (1.48) found in a similar work presented by Baring (1999) for $`\gamma `$5. There is under way a computational approach and a more detailed study on parallel relativistic shocks and on highly relativistic oblique shocks which could be presented in a following paper.

Conclusions

A Monte Carlo numerical investigation has been reported in the test particle limit of parallel diffusive shock acceleration. Very high gamma flow astrophysical plasmas have been used, from $`\gamma _{up}`$ $`5`$ up to $`\gamma _{up}`$ $`1000`$ which are relevant to AGNs but mostly to highly relativistic astrophysical supersonic plasma flows in GRBs. Very distinctive spectral shapes (plateau-like structures) have been demonstrated for both large angle scattering and pitch angle diffusion and the spectral indexes have been found to be in a good agreement with theoretial predictions, exhibiting a spectral flattening in the relativistic regime. These results reveal categorically a problem of the spectral shape of particles accelerated at ultra-relativistic shocks and provide fresh evidence regarding the behaviour of particles regarding the energy gain at every consecutive shock cycle concerning the large angle scattering and pitch angle diffusion mechanisms. This crucial issue, needs further investigation as many other parameters need to be included within the simulation codes in order for the models to be as realistic as possible.

References

Baring, M.G, 1999, (ICRC Salt Lake City, Proceedings)
Kirk,J.G., & Schneider, P., 1987a, ApJ 315, 425
Kirk, J.G, $`\&`$ Schneider, P., 1987, AJ 322, 256-265
Kirk,J.G., & Webb, G.M., 1988, ApJ 331, 336
Lieu, R., & Quenby, J.J., 1990, ApJ 350, 692
Meszaros,P., & Rees, M.J., 1993, ApJ 405, 278
Peacock, J.A, 1981, MNRAS 196, 135-152
Quenby, J.J., & Lieu, R.,1989, Nature 342, 654


text:
# Complex Diffusion Monte-Carlo method for the systems with complex wave function: test by the simulation of 2𝐷 electron in uniform magnetic field

## I Introduction

Quantum-mechanical wave function is an essentially complex in many cases. There are well-known examples – electron in the external uniform magnetic field when vector potential has a central-symmetric form and system of anyons. The complexity of wave function does not give a possibility for the simulations of such kind of systems by using well-known Green Function Monte-Carlo method (see review ). This method essentially demand the reality of the system’s wave function which considered as probability weight in during stochastic process.

There had been undertaken the several attempts  to construct a Monte-Carlo method for the simulation of quantum systems with complex wave function. Authors of the works  suggested to take as a probability weight the module of the complex distribution function. All quantities are calculated in  by averaging over this complex distribution function. The phase of this function $`\alpha `$ is taken accordingly $`\alpha =\alpha _i`$, where $`\alpha _i`$ is the phase at the $`i`$-time step of the stochastic process.

Another quantum Monte Carlo method using algorithm without branching for the simulating of complex problems was developed in . The main trouble of this approach was an increasing of the statistical error as a function of the whole time of the simulation.

The basic difficulty of the numerical simulation of the fermions is essentially the same as for the systems with the complex wave function. Their wave function can change the sign and therefore can not be used as the probability weight in the simulation process.

For the simulation of the continuum (not on the lattice) fermionic systems it were developed widely used fixed node Monte-Carlo method ( see also reviews ) and recently proposed constrained path Monte-Carlo method . In these both methods it was assumed the restriction on the random walks connected with the uncertainties in the space localization of the wave function node surfaces. The comparison of these methods was done in .

Our method is very close to fixed phase Diffusion Monte-Carlo method developed in , which was applied to two-dimensional electrons in magnetic field. Laughlin’s wave function  was used as a trial one. In the framework of Diffusion Monte-Carlo (DMC) method it was calculated only modulus of the system’s wave function. The phase of system’s wave function was not calculated but considered as fixed and equal to the phase of Laughlin wave function.

Here it is proposed new Complex Diffusion Monte-Carlo (CDMC) method including also the simulation of the phase factor of the wave function.

CDMC are tested first by the calculations of the ground state wave function and energy of two-dimensional electron placed into the external uniform magnetic field. There is an excellent agreement between simulations’ results and an analytical ones. Namely, it is reproduced the ground state energy and his degeneracy on the orbital quantum number $`m`$ and also the simulated wave function is in exact correspondence with analytical prediction. The description of CDMC method is given at the section II of the paper. In the section III we represent the description of the simulation algorithm. Analytical expressions for the algorithm quantities are given at section IV. Section V represent the discussion of the simulation results.

## II Complex Diffusion Monte-Carlo method for the quantum systems with complex wave function

Hamiltonian of the electron placed into external uniform magnetic field with vector potential $`\stackrel{}{A}={\displaystyle \frac{1}{2}}[\stackrel{}{H}\stackrel{}{r}]`$, where $`\stackrel{}{H}`$ \- magnetic field, $`\stackrel{}{r}`$ \- electron’s radius vector, and in external potential $`V(r)`$ has a form:

$$\widehat{H}=\frac{1}{2M}(\stackrel{}{p}+\frac{|e|}{c}\stackrel{}{A})^2+V(r).$$
(1)
Here $`M`$ \- mass and $`e=|e|`$ charge of electron, $`\stackrel{}{p}=i\mathrm{}\stackrel{}{}`$, where $`\stackrel{}{}=\stackrel{}{i}{\displaystyle \frac{}{x}}+\stackrel{}{j}{\displaystyle \frac{}{y}}`$, $`c`$ \- light velocity. By taking in to account a relation $`\stackrel{}{p}\stackrel{}{A}\stackrel{}{A}\stackrel{}{p}=div\stackrel{}{A}=0`$ for the such kind of vector potential one can rewrite Hamiltonian (1) in the following form:

$$\widehat{H}=\frac{\stackrel{}{p}^2}{2M}+\frac{|e|}{Mc}\stackrel{}{A}\stackrel{}{p}+\frac{|e|^2}{Mc^2}\stackrel{}{A}^2+V(r).$$
(2)
Let us introduce a complex distribution function

$$f(\stackrel{}{r},t)=\mathrm{\Psi }_T^{}(\stackrel{}{r})\mathrm{\Psi }(\stackrel{}{r},t).$$
(3)
Here $`\mathrm{\Psi }_T^{}(\stackrel{}{r})`$ is a complex conjugated trial wave function of electron in magnetic field. Wave function $`\mathrm{\Psi }(\stackrel{}{r},t)`$ satisfies a Schrodinger equation with imaginary time ( expressed in $`\mathrm{}`$ units )

$$\frac{\mathrm{\Psi }(\stackrel{}{r},t)}{t}=(\widehat{H}E_T)\mathrm{\Psi }(\stackrel{}{r},t).$$
(4)
It is necessary to choose a trial energy $`E_T`$ in such manner  that at $`t\mathrm{}`$ $`\mathrm{\Psi }(\stackrel{}{r},t)\mathrm{\Psi }_0(\stackrel{}{r})`$ \- an exact stationary wave function of ground state of Hamiltonian (1).

So, by accounting (4), one can write the equation for the distribution function $`ff(\stackrel{}{r},t)`$:

$$\frac{f}{t}=\frac{\stackrel{}{p}^{2}}{2M}f\frac{\stackrel{}{p}}{M}(f\stackrel{}{F}_Q(\stackrel{}{r}))+\frac{|e|}{Mc}\stackrel{}{A}\stackrel{}{p}f+(E_L(\stackrel{}{r})E_T)f,$$
(5)
where

$$\stackrel{}{F}_Q(\stackrel{}{r})=\mathrm{\Psi }_T^1(\stackrel{}{r})\stackrel{}{p}\mathrm{\Psi }_T^{}(\stackrel{}{r}),$$
(6)
$$E_L(\stackrel{}{r})=\mathrm{\Psi }_T^1(\stackrel{}{r})\widehat{H}^{}\mathrm{\Psi }_T^{}(\stackrel{}{r}),$$
(7)
$$\widehat{H}^{}=\frac{1}{2M}(\stackrel{}{p}\frac{|e|}{c}\stackrel{}{A})^2+V(r).$$
(8)
For the $`\stackrel{}{p}=i\mathrm{}\stackrel{}{}`$ and $`D={\displaystyle \frac{\mathrm{}^2}{2M}}`$, and by renaming

$$\stackrel{}{F}_Q(\stackrel{}{r})=2\mathrm{\Psi }_T^1(\stackrel{}{r})\stackrel{}{}\mathrm{\Psi }_T^{}(\stackrel{}{r}),$$
(9)
we have an equation:

$$\frac{f}{t}=D\mathrm{\Delta }f+D\stackrel{}{}(f\stackrel{}{F}_Q(\stackrel{}{r}))\frac{i\mathrm{}|e|}{Mc}\stackrel{}{A}\stackrel{}{}f+(E_L(\stackrel{}{r})E_T)f.$$
(10)
Here $`\mathrm{\Delta }=\stackrel{}{}^2`$.

In general case when $`\mathrm{\Psi }_T^{}(\stackrel{}{r})`$ is a complex wave function the quantity $`\stackrel{}{F}_Q(\stackrel{}{r})`$ has a form:

$$\stackrel{}{F}_Q(\stackrel{}{r})=Re\stackrel{}{F}_Q(\stackrel{}{r})+iIm\stackrel{}{F}_Q(\stackrel{}{r}).$$
(11)
So, a final view of an equation for the distribution function $`f`$ is

$`\begin{array}{c}{\displaystyle \frac{f}{t}}=D\mathrm{\Delta }f+D\stackrel{}{}(fRe\stackrel{}{F}_Q(\stackrel{}{r}))+i[\stackrel{}{}(DfIm\stackrel{}{F}_Q(\stackrel{}{r}))\\ {\displaystyle \frac{\mathrm{}|e|}{Mc}}\stackrel{}{A}\stackrel{}{}f]+(E_L(\stackrel{}{r})E_T)f.\end{array}`$ (14)
Following , we assume: when the time step of integration of equation (14) $`\tau 0`$ the function $`\stackrel{}{F}_Q(\stackrel{}{r})`$ remains constant, i.e. we will suppose that $`\stackrel{}{F}_Q(\stackrel{}{r})\stackrel{}{F}_Q(\stackrel{}{r}^{})`$ at $`\tau 0`$, where vector $`\stackrel{}{r}`$ corresponds to time point $`t+\tau `$ and $`\stackrel{}{r}^{}`$ to point $`t`$.

Let us introduce a new quantity

$$\stackrel{}{A}_Q(\stackrel{}{r},\stackrel{}{r}^{})=\frac{1}{2}Im\stackrel{}{F}_Q(\stackrel{}{r}^{})\frac{\mathrm{}|e|}{2DMc}\stackrel{}{A}(\stackrel{}{r}).$$
(15)
Then at $`\tau 0`$ the Green function of equation (14) has a form:

$$G(\stackrel{}{r},\stackrel{}{r}^{};\tau )=G_1(\stackrel{}{r},\stackrel{}{r}^{};\tau )\mathrm{exp}\left[\tau (E_L(\stackrel{}{r})E_T)\right]\mathrm{exp}\left[i\stackrel{}{A}_Q(\stackrel{}{r},\stackrel{}{r}^{})(\stackrel{}{r}\stackrel{}{r}^{}D\tau Re\stackrel{}{F}_Q(\stackrel{}{r}^{}))\right],$$
(16)
where

$$G_1(\stackrel{}{r},\stackrel{}{r}^{};\tau )=\frac{\mathrm{exp}[D\tau \stackrel{}{A}_Q^2(\stackrel{}{r},\stackrel{}{r}^{})]}{4\pi D\tau }\mathrm{exp}\left[\frac{(\stackrel{}{r}\stackrel{}{r}^{}D\tau Re\stackrel{}{F}_Q(\stackrel{}{r}^{}))^2}{4D\tau }\right].$$
(17)
One can see that the Green function $`G(\stackrel{}{r},\stackrel{}{r}^{};\tau )`$ (16) and distribution function $`f`$ given by (3) are the complex functions. These both quantities are related by usual integral equation:

$$f(\stackrel{}{r},t+\tau )=𝑑\stackrel{}{r}^{}G(\stackrel{}{r},\stackrel{}{r}^{};\tau )f(\stackrel{}{r}^{},t).$$
(18)
From equation (18) it is followed that the modulus and the phase of the distribution function at consequent time point are determined by the modulus and the phase of the Green function and of ones of the distribution function at previous time point of integration of diffusion equation (14).

We have to note that for the general case of complex wave function $`\mathrm{\Psi }_T^{}(\stackrel{}{r})`$ the energy $`E_L(\stackrel{}{r})`$ is also a complex, so real and imaginary part of $`E_L(\stackrel{}{r})`$ contribute to ones of the Green function. Therefore the last two exponents in (16) has a form:

$`\begin{array}{c}\mathrm{exp}[\tau (ReE_L(\stackrel{}{r})E_T)]\times \\ \times \mathrm{exp}\left[i\stackrel{}{A}_Q(\stackrel{}{r},\stackrel{}{r}^{})(\stackrel{}{r}\stackrel{}{r}^{}D\tau Re\stackrel{}{F}_Q(\stackrel{}{r}^{}))i\tau ImE_L(\stackrel{}{r})\right].\end{array}`$ (21)
## III  The description of the simulation algorithm

1) Let us to have the $`N_c`$ of the initial configurations $`\stackrel{}{r}`$, i.e. a set of initial systems in which the electron position has random and uniform distribution. One can choose the configurations of $`\stackrel{}{r}`$ and with the distribution function $`f(\stackrel{}{r},0)=|\mathrm{\Psi }_T(\stackrel{}{r})|^2`$, because at $`t=0`$ (3) is real.

The choice of the boundary condition depends on the problem. For the one electron in the external magnetic field a boundary is periodic, i.e., if electron cross a boundary from one side of simulated cell, it enters into cell from opposite side. For example, for the system of bosons in 2D parabolic well and in external magnetic field the boundaries must be free because parabolic well determines itself the boundary of space distribution of the particles.

2) It is performed the quantum drift and diffusion of the particle from $`k`$ \- th configuration, for example, in accordance with the formula

$$\stackrel{}{r}_k=\stackrel{}{r}_k^{}+D\tau Re\stackrel{}{F}_Q(\stackrel{}{r}_k^{})+\chi .$$
(22)
Here $`\chi `$ \- is a gaussian random number having mean value zero and dispersion $`2\sqrt{D\tau }`$.

3) The transition into new space point in this configuration is accepted with probability

$$P(\stackrel{}{r}^{}\stackrel{}{r},\tau )\mathrm{min}(1,W(\stackrel{}{r},\stackrel{}{r}^{})),$$
(23)
where

$$W(\stackrel{}{r},\stackrel{}{r}^{})=\frac{|\mathrm{\Psi }_T(\stackrel{}{r})|^2G_1(\stackrel{}{r}^{},\stackrel{}{r};\tau )}{|\mathrm{\Psi }_T(\stackrel{}{r}^{})|^2G_1(\stackrel{}{r},\stackrel{}{r}^{};\tau )}.$$

Here $`G_1(\stackrel{}{r},\stackrel{}{r}^{};\tau )`$ is given by (17).

If transition of electron is accepted then in accordance with (21) it has new phase, if no then electron keep his old phase.

4) After changing of the electron position from $`k`$ \- configuration into new space point, are calculated $`ReE_L(\stackrel{}{r}_k)`$, $`ImE_L(\stackrel{}{r}_k)`$ and other quantities of interest.

5) By using of the first exponential factor in (21) it is calculated the multiplicity $`M_k`$ (the branching probability) for the configuration $`k`$ accordingly

$$M_k=\mathrm{exp}[\tau (ReE_L(\stackrel{}{r}_k)E_T)]$$
(24)
If $`M_k`$ is not integer, we add an uniformly distributed random number between 0 and 1 to it and take $`M_k`$ equal to nearest integer.

6) If $`M_k0`$ then $`M_k`$ copies of new $`k`$-th configuration place in the list of the new $`N`$ configurations that one is the initial at the next step $`\tau `$ of integration of the diffusion equation. If $`M_k=0`$ then there is no $`k`$-th configuration in the list of new $`N`$ configuration for the next time step $`\tau `$.

All quantities of interest as $`ReE_L(\stackrel{}{r}_k)`$, $`ImE_L(\stackrel{}{r}_k)`$ and so on are multiplied by the factor $`M_k`$ for the calculating of mean values of these ones.

7) It is repeated steps 2)-6) of algorithm until all $`N_c`$ configurations will not overlooked and electrons on these configurations will not simulated on the displacement and the having a new phase.

8) It are calculated a mean energy and others mean quantities on the $`N`$ number of configurations, got in the point 6) of algorithm, at this time step $`\tau `$ in accordance with formula (29) (see below, where it is necessary to change $`M`$ by $`N`$ and to take into account that the phase of electron is equal to the phase of the Green function $`\alpha _G`$ and at the calculation of a mean energy that the energy $`E_L(\stackrel{}{r}_k)`$ has a real $`ReE_L(\stackrel{}{r}_k)`$ and an imaginary $`ImE_L(\stackrel{}{r}_k)`$ parts.)

9) It are repeated steps 1)-8) an integer number of time steps $`\tau `$ of integration of diffusion equation. After that it are determined the mean values of quantities $`Re\overline{E}`$, $`Im\overline{E}`$ and others on this integer number too. It is redetermined the new value of $`E_T`$ in accordance with $`(E_T)_{new}=[(E_T)_{old}+Re\overline{E}]/2`$ accordingly with assumption that $`Im\overline{E}<<Re\overline{E}`$. An integer number of time steps $`\tau `$ represents an one time block $`\mathrm{\Delta }t`$.

10) Every time block $`\mathrm{\Delta }t`$ has the $`N_c`$ of initial configurations. The $`N_c`$ of initial configurations are filled randomly by configurations of just ended time block. The random choice of configurations consists of two steps:
a) a random choice of a number of time step $`\tau `$ in every time block $`\mathrm{\Delta }t`$;
b) a random choice of configuration from $`N`$ configurations at this fixed time step $`\tau `$.
In this manner filled list of $`N_c`$ configurations will be initial one for the next time block $`\mathrm{\Delta }t`$.

11) The repeating of big number of time blocks $`\mathrm{\Delta }t`$ decreases essentially a correlation between configurations in neighbor time blocks and provides a right calculation of mean quantities.

Let us discuss in detail the calculation of the mean quantities. In general case (see ) the mean value of the some quantity $`F(\stackrel{}{R}(t))`$ with the complex distribution function $`f(\stackrel{}{R},t)`$ at the time $`t`$ of the running process is

$$<F(t)>=\frac{{\displaystyle \underset{i=1}{\overset{M}{}}}\mathrm{exp}[i\alpha (\stackrel{}{R}_i(t)]F(\stackrel{}{R}_i(t))}{{\displaystyle \underset{i=1}{\overset{M}{}}}\mathrm{exp}[i\alpha (\stackrel{}{R}_i(t)]}.$$
(25)
Here $`\stackrel{}{R}_i(t)`$ is the coordinates $`\stackrel{}{r}_1,\stackrel{}{r}_2,\mathrm{},\stackrel{}{r}_N`$, $`N`$ \- number of particles of system, $`M`$ \- number of configurations at time point $`t`$. The particles coordinates $`\stackrel{}{r}_1,\stackrel{}{r}_2,\mathrm{},\stackrel{}{r}_N`$ in (21) are weighted with probability $`|f(\stackrel{}{R}_i,t)|`$, and the quantity $`\alpha (\stackrel{}{R}_i(t))`$ is a phase of the distribution function $`f(\stackrel{}{R}_i,t)`$.

As it is clear from the integral expression (18), the phase of the distribution function at consequent time moment determines through the phase of the Green function and the one of the distribution function at previous time moment. So, we have:

$$\alpha (\stackrel{}{R}_i(t+\tau ))=\alpha _G(\stackrel{}{R}_i(t+\tau ),\stackrel{}{R}_j(t))+\alpha (\stackrel{}{R}_j(t)).$$
(26)
Here $`\alpha _G(\stackrel{}{R}_i(t+\tau ),\stackrel{}{R}_j(t))`$ is the phase of the Green function (16) (with account of (21)), and index $`j`$ shows that $`\stackrel{}{R}_j`$ taken from configurations at time moment $`t`$. By substituting (26) into (25) we have

$`\begin{array}{c}<F(t+\tau )>={\displaystyle \frac{{\displaystyle \underset{i=1}{\overset{M}{}}}e^{i\alpha _G(\stackrel{}{R}_i(t+\tau ),\stackrel{}{R}_j(t))+i\alpha (\stackrel{}{R}_j(t))}F(\stackrel{}{R}_i(t+\tau ))}{{\displaystyle \underset{i=1}{\overset{M}{}}}e^{i\alpha _G(\stackrel{}{R}_i(t+\tau ),\stackrel{}{R}_j(t))+i\alpha (\stackrel{}{R}_j(t))}}}=\\ ={\displaystyle \frac{{\displaystyle \underset{i=1}{\overset{M}{}}}e^{i\alpha _G(\stackrel{}{R}_i(t+\tau ),\stackrel{}{R}_j(t))}F(\stackrel{}{R}_i(t+\tau ))}{{\displaystyle \underset{i=1}{\overset{M}{}}}e^{i\alpha _G(\stackrel{}{R}_i(t+\tau ),\stackrel{}{R}_j(t))}}}.\end{array}`$ (29)
The mean quantity $`<F(t+\tau )>`$ in (29) is determined only by the phase $`\alpha _G`$ of the Green function, because all quantities under sum in the numerator and the denominator are weighted with the probability $`|f(\stackrel{}{R}_i(t+\tau ),t+\tau )|`$, i.e. at time moment $`t+\tau `$ with new configurations $`\stackrel{}{R}_i(t+\tau )`$.

Next, let us consider the expression for the phase of the Green function (21). We have

$$\alpha _G=\stackrel{}{A}_Q(\stackrel{}{r},\stackrel{}{r}^{})(\stackrel{}{r}\stackrel{}{r}^{}D\tau Re\stackrel{}{F}_Q(\stackrel{}{r}^{}))\tau ImE_L(\stackrel{}{r}).$$
(30)
From (6) and (7) it is seen that the expressions $`Re\stackrel{}{F}_Q(\stackrel{}{r}^{}))`$ and $`E_L(\stackrel{}{r})`$ have no $`\tau `$ dependence. From (22) one can show that at $`\tau 0`$ $`|\stackrel{}{r}||\stackrel{}{r}^{}|`$ as $`\tau ^{1/2}`$, because $`\chi `$ has $`\tau ^{1/2}`$ dependence. At the same limit $`\stackrel{}{A}_Q(\stackrel{}{r},\stackrel{}{r}^{})\stackrel{}{A}_Q(\stackrel{}{r}^{},\stackrel{}{r}^{})`$ (see the expression (15) for the $`\stackrel{}{A}_Q(\stackrel{}{r},\stackrel{}{r}^{})`$), i.e. there is no $`\tau `$ dependence, also. So, we have $`\alpha _G\tau ^{1/2}`$.

## IV Analytical expressions for the algorithm quantities

We take trial wave function of electron in the form:

$$\mathrm{\Psi }_T^{}(\stackrel{}{r})=\frac{C}{a_H}\left(\frac{\alpha x+i\beta y}{a_H}\right)^m\mathrm{exp}\left(\gamma \frac{(x^2+y^2)}{4a_H}\right).$$
(31)
Here $`C`$ is normalization constant; $`a_H=(\mathrm{}/Mw_H)^{1/2}`$ is magnetic length, where $`w_H=|e|H/Mc`$ \- cyclotron frequency; $`x`$ and $`y`$ the two spatial components of electron’s coordinate; $`m`$ is orbital quantum number; $`\alpha `$, $`\beta `$ and $`\gamma `$ are arbitrary numerical constants.

At $`\alpha =\beta =\gamma =1`$ (31) gives an exact analytical expression of electron ground state wave function. We deform exact electron wave function, i.e. take as the input $`\alpha `$, $`\beta `$ and $`\gamma `$ not equal to $`1`$, and then study a relaxation of the trial function to the exact one via simulation process.

As length unit we take $`a_H`$, as energy unit - $`\mathrm{}w_H/2`$, and as time unit - $`2/\mathrm{}w_H.`$

By substituting of the wave function (31) into (9), we have

$`\begin{array}{c}Re\stackrel{}{F}_Q(\stackrel{}{r})={\displaystyle \frac{1}{\alpha ^2x^2+\beta ^2y^2}}[\stackrel{}{i}(2m\alpha ^2x\gamma \alpha ^2x^3\gamma \beta ^2xy^2)+\hfill \\ +\stackrel{}{j}(2m\beta ^2y\gamma \alpha ^2yx^2\gamma \beta ^2y^3)],\hfill \\ Im\stackrel{}{F}_Q(\stackrel{}{r})={\displaystyle \frac{2m\alpha \beta }{\alpha ^2x^2+\beta ^2y^2}}[y\stackrel{}{i}+x\stackrel{}{j}].\hfill \end{array}`$ (35)
Here and below $`\stackrel{}{i}`$ and $`\stackrel{}{j}`$ the unit vectors in $`x`$ and $`y`$ directions correspondingly.

For the $`\stackrel{}{A}_Q(\stackrel{}{r},\stackrel{}{r}^{})`$ (15) gives

$`\begin{array}{c}\stackrel{}{A}_Q(\stackrel{}{r},\stackrel{}{r}^{})={\displaystyle \frac{1}{\alpha ^2x^2+\beta ^2y^2}}[\stackrel{}{i}(2m\alpha \beta y^{}+y(\alpha ^2x^2+\beta ^2y^2))+\hfill \\ +\stackrel{}{j}(2m\alpha \beta x^{}x(\alpha ^2x^2+\beta ^2y^2))].\hfill \end{array}`$ (38)
The wave function (31) and eq. (7) (here $`V(r)=0`$) gives

$$ReE_L(\stackrel{}{r})=ReE_L^{(1)}(\stackrel{}{r})+ReE_L^{(2)}(\stackrel{}{r}),$$
(39)
where

$$ReE_L^{(1)}(\stackrel{}{r})=\gamma (m+1)+(x^2+y^2)\left[\frac{1\gamma ^2}{4}\frac{m\alpha \beta }{\alpha x^2+\beta ^2y^2}\right],$$

$$ReE_L^{(2)}(\stackrel{}{r})=\frac{m(m1)(\beta ^2\alpha ^2)(\alpha x^2\beta ^2y^2)}{(\alpha ^2x^2+\beta ^2y^2)^2};$$

and also

$$ImE_L(\stackrel{}{r})=\frac{mxy(\beta ^2\alpha ^2)}{\alpha ^2x^2+\beta ^2y^2}\left[1\frac{(m1)2\alpha \beta }{\alpha ^2x^2+\beta ^2y^2}\right].$$
(40)
By using formulas (35)-(40), it is easy to write the analytical expressions for other quantities.

## V The discussion of the simulation results.

The simulation of the ground state of one electron in magnetic field is a simple problem and was performed on the usual personal computer.

Figs.1-4 are presented the result of the calculations of the real $`Re\overline{E}`$ and imaginary part $`Im\overline{E}`$ of the ground state energy for different orbital quantum numbers $`m`$ and different starting trial functions defined by the parameters $`\alpha ,\beta ,\gamma `$ as a functions of the number $`T`$ of time blocks $`\mathrm{\Delta }t`$ (see point 9) of algorithm) ( a full number of time blocks are chosen ten in each running).

In every time block $`\mathrm{\Delta }t`$ initial number of configurations $`N_c`$ was chosen equal 1000 and number of time steps $`\tau `$ equal 100. We everywhere suggested $`\tau `$ equal 0.01. Typical mean number of the population number $`N`$ ( see point 6) of the algorithm) was around 2500 at initial blocks $`\mathrm{\Delta }t`$ and around 1000 at final blocks $`\mathrm{\Delta }t`$.

Fig.1 present $`Re\overline{E}`$ and Fig.2 present $`1000Im\overline{E}`$ for $`m=0,4,8`$ with parameters $`\alpha =\beta =\gamma =1`$ as a functions of $`T.`$

It is seen from Figs.1,2 $`Re\overline{E}`$ is approaching fast to the exact value of ground state energy and reproduce very well an orbital quantum number $`m`$ degeneracy, while the imaginary part of energy $`Im\overline{E}10^3Re\overline{E}.`$

The additional calculations shows, that the increasing of the statistics, i.e. the increasing of $`N_c`$ and a number of steps $`\tau `$ in each block $`\mathrm{\Delta }t`$, decreases $`Im\overline{E}`$ more.

The Figs.3,4 represent the $`Re\overline{E}`$ and $`Im\overline{E}`$ for the different choices of the $`m=0`$ trial function. From these Figs. we can see that the deformation of wave function in some limits does not alter of the relaxation of the running process to exact final state.

Fig.5 present the initial uniform spatial distribution of electron before simulation, while Fig.6 - the final spatial distribution of electron after simulation for $`m=13`$.

In this simulation the trial wave function $`\mathrm{\Psi }_T^{}`$ is taken with parameters $`\alpha =\beta =\gamma =1`$, $`m=13`$ and $`N_c=500`$. Since the trial function coincide with exact wave function we have to expect fast convergence of initial distribution to exact one if the algorithm is correct. In accordance with analytical solution (see ) the spatial distribution for one electron in this case must have a ring like form with mean radius $`(2m+1)^{1/2}`$ and wide 1 ( in magnetic length $`a_H`$ units). As it is seen from Fig.6, there is an excellent agreement between simulation and analytical results.

We suppose also further tests of the method by simulations of anyons .