text:
# The 𝑧=2.51 Extremely Red Submillimeter Galaxy SMM J04431+0210

## 1. Introduction

Since the discovery of the submillimeter population of galaxies (Smail, Ivison, & Blain 1997; Hughes et al. 1998; Barger et al. 1998; Eales et al. 1999), extensive follow-up observations across the entire electromagnetic spectrum ranging from X-ray to radio wavelengths have been obtained (see review by Blain et al. ). Despite these observational efforts, fundamental questions remain about the population as a whole, foremost among these are the redshift distribution and dominant energy source heating the dust. Detailed studies of individual sub-mm galaxies are crucial in answering these questions, however, progress has been slow since many sub-mm sources are extremely faint in their observed optical bands. Several of the original potential bright optical counterparts are now suspected to be mis-identifications based on deep near-infrared observations (Smail et al. 1999; Frayer et al. 2000; Dunlop et al. 2003) and redshift estimates based on radio data (Smail et al. 2000; Webb et al. 2003).

The current data show that the sub-mm population has a mixture of AGN and starburst characteristics with properties which are roughly consistent with those of the local population of ultraluminous ($`L>10^{12}L_{}`$) infrared galaxies (ULIRGs). The molecular CO line data for the sub-mm galaxies show the presence of sufficient molecular gas to fuel high levels of star-formation activity (Frayer et al. 1998, 1999). However, the relative importance of AGN and starburst activity in powering the high luminosities of the sub-mm population is still open to debate, as is the case with local ULIRGs (e.g., Sanders & Mirabel 1996; Genzel et al. 1998; Veilleux, Kim, & Sanders 1999). In order to understand the nature of the sub-mm population, we have been carrying out multi-wavelength observations of individual systems in the SCUBA Cluster Lens Survey (Smail et al. 2002). This survey consists of 15 background sub-mm sources identified from sensitive sub-mm mapping of seven massive, lensing clusters. The advantage of this sample is that the amplification of the background sources allows for deeper source frame observations.

In this paper we report a redshift of $`z=2.5092\pm 0.0008`$ for the sub-mm galaxy SMM J04431+0210 (N4), obtained from near-infrared observations with the NIRSPEC instrument on the Keck II telescope. The galaxy N4 at $`04^\mathrm{h}43^\mathrm{m}07\stackrel{\mathrm{s}}{\mathrm{.}}1`$, $`+02\mathrm{°}10\mathrm{}25\mathrm{}`$ (J2000) was previously identified as the likely sub-mm counterpart based on its extremely red I$``$K color (Smail et al. 1999). An increasing number of sub-mm sources are now believed to be extremely red objects (EROs \[$`RK>6`$\]) due to their high dust content (Smail et al. 1999; Gear et al. 2000; Lutz et al. 2001; Ivison et al. 2001; Wehner, Barger, & Kneib 2002). Since N4 is the brightest K-band (2.2$`\mu `$m) ERO counterpart in the Cluster Lens Survey, it is an ideal candidate for testing whether redshifts can be obtained for the SCUBA population in the near-infrared. In addition, N4 has a 850$`\mu `$m flux density of only 1.6 mJy, after correcting for lensing, making it one of the intrinsically weakest sub-mm sources yet identified. Studies of weak sub-mm sources such as N4 are crucial to our understanding of the population as a whole since the majority of the sub-mm background arises from sources with S(850$`\mu `$m)$`<2`$ mJy (Blain et al. 1999a; Cowie, Barger, & Kneib 2002). We also report on follow-up CO observations with the OVRO Millimeter Array to search for molecular gas, as well as optical spectroscopy with the ESI instrument on the Keck II telescope to search for additional rest-frame optical and ultraviolet (UV) emission lines. The observations presented in this paper are summarized in Table 1.

A cosmology of H$`{}_{o}{}^{}=70\mathrm{km}\mathrm{s}^1\mathrm{Mpc}^1`$, $`\mathrm{\Omega }_\mathrm{M}=0.3`$, and $`\mathrm{\Omega }_\mathrm{\Lambda }=0.7`$ is assumed throughout this paper. We adopt the lensing amplification factor of $`f=4.4`$ which was derived by Smail et al. (1999) and includes both the amplification due to the foreground cluster potential and the spiral galaxy near the position of N4.

## 2. OBSERVATIONS AND DATA REDUCTION

### 2.1. Spectroscopy

N4 was observed with NIRSPEC (Near Infrared Spectrograph) on Keck II111The W. M. Keck Observatory is operated as a scientific partnership among the California Institute of Technology, the University of California, and the National Aeronautics and Space Administration. The Keck Observatory was made possible by the generous financial support of the W. M. Keck Foundation. We are most fortunate to have the opportunity to conduct observations from the summit of Mauna Kea which has a very significant cultural role within the indigenous Hawaiian community. on 2000 December 26. NIRSPEC is a spectrograph operating over 0.95 to 5.4$`\mu `$m (McLean et al. 1998). Observations were obtained in low-resolution mode (R$``$1500) with a $`0\stackrel{}{\mathrm{.}}76`$ (4 pixel) wide slit. The seeing was typically 0$`\stackrel{}{\mathrm{.}}5`$ during the spectroscopic observations. The 42$`\mathrm{}`$ length slit was placed at a position angle of $`14^{}`$ to contain both N4 and the core of the bright nearby galaxy $`3\stackrel{}{\mathrm{.}}1`$ from N4 (Fig. 1). We used the nspec7 filter to search the entire K-band and obtained 130 minutes of useful on-source integration. H$`\alpha `$ emission was detected at $`z=2.51`$ (Sec. 3.1).

In 2002 August we followed up the K-band detection of H$`\alpha `$ by searching for the \[OIII\] 5007, 4959 doublet and H$`\beta `$ line in H-band using NIRSPEC. We used the same slit setup as done previously for K-band. The nspec6 filter was used since it is optimized on the redward side of H-band at the wavelengths of the redshifted emission lines.

The NIRSPEC data were combined in pairs at two different positions along the slit. This provided a zeroth-order background subtraction as well as removing the dark and bias levels from the data. Pixel-to-pixel gain variations were removed by dividing by a normalized flat produced from halogen lamp images. The spectra were transformed to produce a linear spatial and wavelength scale using standard tasks in the IRAF LONGSLIT package.222IRAF is the Image Reduction and Analysis Facility sofware package written and supported by the National Optical Astronomy Observatories. Observations of a star nodded along the slit provided spatial calibration, while neon and argon arc lamp lines as well as sky lines provided wavelength calibration. A linear background was subtracted from the frames as a function of wavelength to remove time variations in the sky emission. The spectra of bright calibration stars were used to correct for variations in the detector response as a function of wavelength.

Given its proximity, residual light from the nearby galaxy is present in the extracted spectrum of N4 after standard data reduction. In order to obtain the true continuum level and spectral slope for N4, the wings of the galaxy light were subtracted from N4 as a function of wavelength. There is a natural trade-off between the signal-to-noise (S/N) ratio of the final spectrum and the wavelength resolution over which the subtraction is done. We smoothed over two resolution elements in wavelength for the galaxy subtraction, which left a slight residual narrow Paschen-$`\alpha `$ emission line from the foreground galaxy ($`z=0.18`$) at $`2.1\mu `$m. The resulting continuum level of the K-band spectrum was calibrated by averaging over the emission from 2 to 2.4$`\mu `$m and placed on the K-band magnitude scale derived from the NIRC imaging (Sec. 2.2). The absolute flux level of the H-band spectrum was calculated using the continuum levels of the bright nearby galaxy and N4 at 2.0$`\mu `$m, where the H- and K-band spectra have overlapping wavelength coverage.

In 2002 August we also obtained follow-up spectroscopy observations with the ESI (Echellette Spectrograph and Imager) instrument (Sheinis et al. 2002) on Keck II to search for additional emission lines in the rest-frame ultraviolet and optical bands. The ESI data were taken in the high-resolution echelle mode which provided wavelength coverage from 3900Å to 1100Å at a resolution of approximately $`11\mathrm{km}\mathrm{s}^1`$. A slit with a width of $`1\stackrel{}{\mathrm{.}}25`$ was used. The ESI optical spectra were calibrated using the spectrophotometric standard Feige 110 (Massey et al. 1988).

### 2.2. Near-Infrared Imaging

We observed N4 using the Near Infrared Camera (NIRC) on Keck I in 1999 October in K-band and 2001 August in J and K-band. NIRC employs a $`256\times 256`$ pixel InSb detector with a pixel scale of $`0\stackrel{}{\mathrm{.}}15`$ (Matthews & Soifer 1994). The standard $`K`$-band filter was used instead of the bluer $`K_\mathrm{s}`$ filter since the object is extremely red. We observed using a random dithered pattern with about $`10\mathrm{}`$ between adjacent dither positions. Sixty second integrations were taken at each position. A total of 22 and 74 minutes of on-source data were obtained in the K and J filters, respectively. The seeing-disks of the stars observed throughout the exposures varied from about $`0\stackrel{}{\mathrm{.}}5`$ to $`0\stackrel{}{\mathrm{.}}7`$ (FWHM).

Sets of dark frames from each night were subtracted from each exposure to remove the dark current as well as the bias level. The dark-subtracted exposures were divided by a normalized skyflat. Frames were sky-subtracted using temporally–adjacent images. Individual exposures were aligned to the nearest pixel using common objects in the frames. The data were placed on the Vega-magnitude scale from observations of a set of near-infrared standard stars (Persson et al. 1998) taken at a range of air masses to correct for atmospheric extinction. Based on the dispersion in the zero-points derived throughout the nights, the uncertainty of the derived magnitude scale is estimated to be better than 0.06 magnitude.

### 2.3. Millimeter Interferometry

Based on the redshift derived from the NIRSPEC observations, N4 was observed using the Owens Valley Millimeter Array333The Owens Valley Millimeter Array is a radio telescope facility operated by the California Institute of Technology and is supported by NSF grant AST 9981546. to search for CO emission between 2001 January and April. A total of 23 hours of usable integration time on source was obtained in two configurations of six 10.4m telescopes. The array configurations yielded a $`6\mathrm{}\times 5\mathrm{}`$ synthesized beam adopting natural-weighting. We used a digital correlator configured with $`112\times 4`$ MHz channels centered on 98.5174 GHz to search for CO(3$``$2) emission at the corresponding H$`\alpha `$ line redshift of $`z=2.510`$. The spectral-line bandwidth corresponds to 1400$`\mathrm{km}\mathrm{s}^1`$ ($`\mathrm{\Delta }z=\pm 0.005`$), which is sufficiently large to observe the entire H$`\alpha `$ line width, including the redshift uncertainty (Sec. 3.1). In addition to the CO line data, 3mm continuum data were recorded with a bandwidth of 2 GHz. The nearby quasar J0433+053 was observed every 25 minutes for gain and phase calibration. Absolute flux calibration was determined from observations of the standard calibrator 3C 273 whose flux history is monitored by observations of Uranus and Neptune. The absolute flux calibration uncertainty for the data is approximately 15%. The data were reduced using the OVRO MMA software (Scoville et al. 1983) and standard tasks in AIPS.444AIPS is the Astronomical Image Processing System sofware package written and supported by the National Radio Astronomy Observatory.

## 3. RESULTS

The identification of N4 as the sub-mm counterpart is fairly secure based on previous studies. As discussed by Smail et al. (1999), the nearby $`z=0.18`$ spiral galaxy can be reliably ruled out based on its radio/sub-mm flux density ratio and its 450$`\mu `$m limit. The random probability of finding a bright ERO similarly as red ($`RK>6`$) as N4 within the $`3\mathrm{}`$ positional uncertainty of the SCUBA data is less than 0.1%, based on the observed ERO surface densities for $`K<20`$ (Thompson et al. 1999; Smith et al. 2002). If we relax the sub-mm positional uncertainty and account for lensing, the likelihood of a chance coincidence with an ERO is still less than 1%, and no other candidate red objects are nearby. The identification of N4 as the sub-mm counterpart has recently been confirmed by the detection of CO emission from N4 at the H$`\alpha `$ redshift with the IRAM PdB interferometer (R. Neri et al. 2003, in preparation). The OVRO CO upper-limit presented in this paper is consistent with the more sensitive PdB observations.

### 3.1. Spectroscopy

An emission-line complex is clearly detected in the K-band spectrum of N4 (Fig. 2), which is identified as H$`\alpha +`$\[NII\] at $`z=2.51`$. We also detect evidence for a weak emission line at 1.756$`\mu `$m, which corresponds to \[OIII\] 5007Å at the H$`\alpha `$ redshift (Fig. 3). Fits to the lines and continuum were made using the SPECFIT package in IRAF (Table 2). In fitting the H$`\alpha +`$\[NII\] line complex, we have allowed the widths and fluxes of the H$`\alpha `$ and \[NII\]6583Å lines to be free, while fixing the \[NII\]6548Å line to have a flux of 30% that of \[NII\]6583Å. Both \[NII\] lines were constrained to have the same FWHM. The continuum was fitted with a linear function. Initially the position of the \[NII\]6583Å line was fixed with respect to H$`\alpha `$, but this requirement was later relaxed to derive the redshift and to confirm the identification. From the fit to the total K-band spectrum, we derive a redshift of $`z=2.5092\pm 0.0008`$ based on the H$`\alpha `$ and NII 6583 lines. The redshift uncertainty includes the error from the fit and the systematic uncertainty associated with the derived wavelength scale. The intrinsic H$`\alpha `$ line width is $`520\pm 40\mathrm{km}\mathrm{s}^1`$, corrected for instrumental resolution, while the \[NII\] lines are slightly narrower, having an intrinsic FWHM of $`440\pm 60\mathrm{km}\mathrm{s}^1`$. The \[NII\] 6583/H$`\alpha `$ line flux ratio is found to be $`0.47\pm 0.06`$. No broad-line H$`\alpha `$ component is readily visible in the data. The 1$`\sigma `$ upper-limit on a $`5000\mathrm{km}\mathrm{s}^1`$ broad H$`\alpha `$ line is about 10% of the H$`\alpha `$ narrow-line flux.

The H-band data were fitted adopting the same techniques as used for the K-band data. The \[OIII\] line was fitted to derive the wavelength position, line width, and line flux. In deriving the H$`\beta `$ limit, we assumed the H$`\alpha `$ redshift and linewidth. The observed limit on the \[OIII\]5007/H$`\beta `$ flux ratio is greater than 1.5 ($`3\sigma `$).

The observed Balmer decrement of H$`\alpha `$/H$`\beta >6`$ ($`3\sigma `$) implies significant extinction. Assuming the standard intrinsic H$`\alpha `$/H$`\beta `$ value of 2.85 (Case-B Balmer recombination decrement for T$`=10^4`$ K and N$`{}_{e}{}^{}=10^4`$ cm-3 \[Osterbrock 1989\]) and the Galactic reddening function parameterized by Lequeux et al. (1979), we derive an extinction at H$`\alpha `$ of greater than $`1.6`$ mag and a color excess of E(B-V) $`>0.7`$ mag. The high-level of extinction implied by the H$`\alpha `$/H$`\beta `$ ratio is consistent with the non-detection of Ly$`\alpha `$ emission from the ESI observations (Table 2).

Figure 4 shows the two-dimensional image of the K-band spectrum. The H$`\alpha `$ emission is spatially resolved, while the \[NII\] emission is consistent with a compact object. The size of H$`\alpha `$ emission region is about 1.1 arcsec, after deconvolution with a Gaussian representing the seeing point spread function. Corrected for lensing, the intrinsic diameter of the detected H$`\alpha `$ emission regions is about 2 kpc. If the observed velocity gradient is due to rotation, a rough estimate of the dynamical mass can be made based on the H$`\alpha `$ line width. The H$`\alpha `$ emission is spread over 750$`\mathrm{km}\mathrm{s}^1`$ (FWZI), and using a radius of 1.0 kpc we calculate a dynamical mass of $`3\times 10^{10}\mathrm{sin}(i)^2M_{}`$, where the inclination ($`i`$) is unknown. The estimated dynamical mass of N4 is comparable with masses found in the central regions of normal galaxies in the local universe.

Figure 5 shows spectra extracted from the different regions of N4. The \[NII\] 6583/H$`\alpha `$ line flux ratio varies measurably as a function of distance from the nucleus of N4. The off-nuclear data have a \[NII\] 6583/H$`\alpha `$ line flux ratio of $`0.16\pm 0.06`$, consistent with star-forming regions, while the nuclear region has a much higher ratio of $`0.66\pm 0.09`$, indicative of an AGN or a LINER (Osterbrock 1989). Therefore, the data for N4 are consistent with a narrow-line AGN surrounded by an extended region of active star-formation (Baldwin, Phillips, & Terlevich 1981; Veilleux & Osterbrock 1987) or a powerful starburst with large-scale nuclear shocks driven by a nuclear superwind. (Armus, Heckman, & Miley 1989; Heckman, Armus, & Miley 1990) .

### 3.2. Imaging

The sensitive HST R-band (F702W filter) and Keck LRIS I-band images of N4 have been previously discussed by Smail et al. (1999). We compare these optical images with the deep J- and K-band data presented in this paper. In all bands we subtracted the nearby galaxy using the IRAF task ELLIPSE, which enabled source masking while successfully removing the wings of the galaxy emission. The images were aligned using common sources in the frames. For optimization in K-band, we adopted a $`2\mathrm{}`$ diameter aperture for photometry. Although the S/N ratios from the aperture measurements in the other bands are lower than for K-band (Table 3), the emission peak was detected in each of the four bands with a S/N ratio greater than five. The errors on the measurements represent the rms of multiple aperture measurements on the sky around N4, avoiding the galaxy residuals from the central disk. We measure the colors of $`RK=6.34\pm 0.27`$ mag and $`JK=3.15\pm 0.29`$ mag for N4. The $`RK`$ measurement presented here is slightly bluer than the original $`RK>6.7`$ color quoted by Smail et al. (1999), which is attributed to the difference in the adopted aperture size and methods of galaxy subtraction.

We have resolved N4 in both the R and K-band images (Fig. 6). After deconvolution of the $`0\stackrel{}{\mathrm{.}}6`$ stellar PSF, approximated by a Gaussian, we derive a FWHM size of $`1\stackrel{}{\mathrm{.}}1\pm 0\stackrel{}{\mathrm{.}}2`$ for the K-band emission, consistent with the H$`\alpha `$ source size (Sec. 3.1). The K-band image is extended in the north-south direction, along the expected direction of lensing amplification, while the R-band image appears extended east-west. The K-band peak is marginally offset along the line of R-band elongation, but its derived position is coincident with the R-band peak within the errors. Offsets between the K-band and R-band emission regions are not unexpected given the high level of extinction within N4 and given previous results for sub-mm galaxies (Ivison et al. 2001).

### 3.3. Dust and Molecular Gas

The OVRO CO observations yield a 3$`\sigma `$ upper-limit of 4.5 mJy/beam, averaging over the H$`\alpha `$ line-width of 520$`\mathrm{km}\mathrm{s}^1`$ with a natural-weighted beam size of $`6\mathrm{}\times 5\mathrm{}`$. The limit on the integrated CO(3$``$2) flux density of $`S(CO)<2.5\mathrm{Jy}\mathrm{km}\mathrm{s}^1`$ corresponds to a CO luminosity of $`L^{}(\mathrm{CO})<1.9\times 10^{10}\mathrm{K}\mathrm{km}\mathrm{s}^1\mathrm{pc}^2`$, correcting for lensing. The CO to H2 conversion factor is still uncertain for the sub-mm population. Adopting a conversion factor of $`L^{}(\mathrm{CO})/\mathrm{M}(\mathrm{gas})=2\mathrm{K}\mathrm{km}\mathrm{s}^1\mathrm{pc}^2(\mathrm{M}_{}\mathrm{yr}^1)^1`$, which is appropriate for low redshift ULIRGs (Solomon et al. 1997; Scoville, Yun, & Bryant 1999), the observations imply M(gas)$`<4\times 10^{10}M_{}`$ for N4. The limit on the gas mass is comparable to the total dynamical mass estimated from the H$`\alpha `$ line-width (Sec. 3.1). The CO(3$``$2) luminosity limit is slightly lower than the previous CO detections for sub-mm selected galaxies (Frayer et al. 1998, 1999), after converting to the same cosmology, but is well within expectations given the lower $`850\mu `$m flux density for N4.

We failed to detect 3mm continuum emission from N4, achieving a 3$`\sigma `$ upper-limit of S(100 GHz)$`<1.5`$ mJy. The 3mm continuum limit is consistent with the previous SCUBA measurements of N4 (Smail et al. 1999). The spectral energy distribution (SED) for N4 is consistent with those found for previous sub-mm galaxies (Ivison et al. 2000) and local ULIRGs. Based on its redshift and the limit on its 850$`\mu `$m/1.4GHz flux ratio, we estimate a dust temperature of T$`>40`$ K from the Carilli & Yun (1999) relationship plotted as a function of temperature by Blain et al. (2002). The dust temperature is consistent with both cool and warm infrared SEDs found for ULIRGs (Soifer et al. 1989).

## 4. DISCUSSION

As mentioned in Section 3, the rest-frame optical emission-line flux ratios of N4 imply a composite source. The relatively large nuclear \[NII\]/H$`\alpha `$ and \[OIII\]/H$`\beta `$ ratios are suggestive of a non-thermal ionizing source, although we cannot rule out shock excitation from a wind. The spatially-resolved H$`\alpha `$ line emission and the lower off-nuclear \[NII\]/H$`\alpha `$ ratio suggest a change to hot stars as the dominant source of ionizing photons outside the nucleus.

The observed \[NII\] 6583/H$`\alpha `$ and \[OIII\] 5007/H$`\beta `$ line ratios for N4 are consistent with the average values found in low redshift ULIRGs with L$`{}_{\mathrm{ir}}{}^{}>2\times 10^{12}L_{}`$, albeit there is a large range of ratios found for local ULIRGs (Veilleux et al. 1999). The limits on the \[OI\] 6300 and \[SII\] $`\lambda \lambda `$6716,6731 lines do not provide any additional constraints favoring one classification over another. Unfortunately, the \[OI\] 6300 line for N4 is redshifted to the wavelength of the Pa$`\alpha `$ from the foreground galaxy. Although there is tentative evidence for \[OI\] 6300 emission even after careful galaxy subtraction, no robust measurement could be derived for this line. The H$`\alpha `$ equilvalent width (EQW) and relatively narrow line width suggest that N4 does not contain a broad-line Type-I AGN, or that the broad-line region is heavily obscured. A higher S/N spectrum would be required to place more stringent limits on the H$`\alpha `$ broad-line flux from N4.

The only X-ray data that currently exists for N4 are from ROSAT, and these data are not sensitive enough to constrain the fraction of the bolometric luminosity which may arise from an AGN. Given the unknown contribution from a possible AGN and the uncertainty in the level of extinction, the star-formation rate (SFR) is fairly uncertain for N4. We can derive upper and lower bounds on the SFR based on the H$`\alpha `$ and sub-mm emission. The off-nuclear H$`\alpha `$ line luminosity uncorrected for extinction provides a lower-limit on the SFR in N4. At the observed resolution, approximately 40% of the total H$`\alpha `$ emission is outside the nucleus. This off-nuclear H$`\alpha `$ emission corresponds to a line luminosity of L(H$`\alpha )=2\times 10^8L_{}`$, corrected for lensing. Adopting a standard relationship between the H$`\alpha `$ line luminosity and the SFR (Kennicutt 1983), we find that SFR in N4 is greater than $`6M_{}\mathrm{yr}^1`$, uncorrected for extinction. Unfortunately, the H$`\beta `$ limit is not sensitive enough to provide a useful estimate of the extinction level outside the nucleus. For comparison, the total H$`\alpha `$ line luminosity corrected for the extinction found from the limit on the Balmer decrement would imply a total SFR of $`>70M_{}\mathrm{yr}^1`$.

A upper-limit to the SFR can be derived by assuming that all the far-infrared (FIR) radiation is due to star-formation, neglecting any contribution due to a possible AGN. The FIR luminosity for N4 is derived assuming that its SED is somewhere between that of the infrared-warm ULIRG Mrk 231 and the infrared-cool ULIRG Arp 220. Placing these ULIRGs at the redshift of N4 and scaling to the de-lensed 850$`\mu `$m flux of N4, we estimate an intrinsic FIR luminosity of L(FIR)$`=`$2–5$`\times 10^{12}L_{}`$ for N4. This FIR luminosity corresponds to a SFR of massive stars of SFR(M$`>5M_{})`$200–500$`M_{}\mathrm{yr}^1`$, using the relationship given by Condon (1992). Including the presence of lower mass stars with the IMF adopted by Kennicutt (1983), we obtain a total SFR$``$800–2000$`M_{}\mathrm{yr}^1`$. Hence, from the H$`\alpha `$ and sub-mm emission, the star-formation rate of N4 is probably at least $`100M_{}\mathrm{yr}^1`$ and could be as large as $`2000M_{}\mathrm{yr}^1`$.

Despite the high level of obscuration inferred for the sub-mm galaxies, the recent Ly$`\alpha `$ detections of radio-selected sub-mm sources by Chapman et al. (2003) suggest that significant rest-frame UV light can escape from the sub-mm sources. For N4 we failed to detect Ly$`\alpha `$ emission. The sensitivity of the ESI observations would have been sufficient to detect many of the Ly$`\alpha `$ sources seen by Chapman et al. (2003). However, given the fact that the Ly$`\alpha `$/H$`\alpha `$ ratio can span three orders of magnitudes for ULIRGs ($`<0.01`$ – 10; HST STIS observations by Surace et al. 2003, in preparation), the observed limit of Ly$`\alpha `$/H$`\alpha <0.14`$ for N4 is not particularly constraining.

N4 has a composite spectrum consisting of a narrow-line AGN/LINER nucleus surrounded by a resolved starburst. Other well-studied sub-mm galaxies have shown a wide range of characteristics. The system SMM J02399-0126 at $`z=2.8`$ shows an AGN spectrum (Ivison et al. 1998) which has been classified as broad absorption line QSO (Vernet & Cimatti 2001), while SMM J14011+0252 at $`z=2.6`$ shows no evidence of an AGN and is thought to be dominated by star formation (Barger et al. 1999; Ivison et al. 2000). The most distant know sub-mm galaxy SMM J09431+4700 ($`z=3.3`$) contains a narrow-line Seyfert-1 AGN (Ledlow et al. 2002), and the sub-mm galaxy N2 850.4 ($`z=2.4`$) is a composite source (AGN+starburst) showing strong outflows and P-Cygni profiles indicative of stellar winds (Smail et al. 2003). All of these sub-mm galaxies are significantly more luminous than N4 and can be classified as hyperluminous sources ($`L10^{13}L_{}`$). It may not be surprising that most these of the hyperluminous sub-mm sources contain AGN given that this is seen at lower redshifts (e.g., Evans et al. 1998). The results for N4 suggest that even the lower luminosity sub-mm sources could also be partially powered by AGN.

Although there is evidence for AGN in many sub-mm sources, the bulk of the infrared light for the high-redshift sub-mm population is still thought to be dominated by star-formation (Blain et al. 1999b, 2002). The X-ray data support this conclusion (Barger et al. 2001; Alexander et al. 2003; Almaini et al. 2003). The deepest X-ray data from Alexander et al. (2003) suggest that even though a significant faction of the bright sub-mm sources contain an AGN, the AGN typically have low luminosities and contribute negligibly to the total bolometric luminosity of the population.

In the local universe, the apparent fraction of AGN dominated ULIRGs increases at luminosities above L$`{}_{\mathrm{ir}}{}^{}>2\times 10^{12}L_{}`$ (Veilleux et al. 1999). It is unclear if the high-redshift sub-mm population of galaxies will follow this trend. The relative importance of AGN and star-formation activity in the sub-mm population will be constrained with future sensitive X-ray observations, SIRTF observations that will distinguish between infrared warm (AGN) versus cool (starburst) sub-mm galaxies, and high-resolution mm and CO interferometric observations that will constrain gas masses and source sizes.

## 5. CONCLUSIONS

We report the redshift of $`z=2.51`$ for the counterpart of the sub-mm galaxy N4. We have detected the H$`\alpha `$, \[NII\], and \[OIII\] lines and present upper-limits for Ly$`\alpha `$, H$`\beta `$, and CO(3$``$2) lines. The line ratios and observed SED are consistent with the range of properties found in local ULIRGs. The data show that N4 is comprised of a nuclear region showing a LINER/Type-II AGN spectrum with an extended H$`\alpha `$ component, presumably arising from an extended starburst. It is still unclear whether N4 is predominately powered by AGN or star-forming activity. Future multi-wavelength observations should provide better constraints on the source of its immense luminosity.

We thank the staff at Keck Observatory and the Owens Valley Millimeter Array who have made these observations possible. We thank B. T. Soifer and E. Egami for obtaining an early H-band spectrum of N4. DTF and LA are supported by the Jet Propulsion Laboratory, California Institute of Technology, under contract with NASA. AWB acknowledges support from the NSF under grant AST-0205937. NAR acknowledges support from a NSF Graduate Research Fellowship, and IS acknowledges support from the Royal Society and the Leverhulme Trust.


format: text

text:
# Ballistic quantum transport at high energies and high magnetic fields

## I INTRODUCTION

Accurate simulations of ballistic transport through quantum dots have remained a computational challenge, despite the conceptional simplicity of the problem. This is in part due to the fact that many of the most interesting phenomena occur in a parameter regime of either high magnetic field $`B`$ or high Fermi energy $`E_F`$. The regime of strong magnetic field $`B`$, where the magnetic length (in a.u.) $`l_B=\sqrt{c/B}`$ is small compared to the linear dimension $`D`$ of the dot, $`l_BD`$, gives rise to the Quantum Hall effect,pran87  the Hofstadter butterfly,hofs76  and Aharonov-Bohm oscillations of transport coefficients.ashc76  The high energy domain, where the De Broglie wavelength $`\lambda _D=\sqrt{2E_F}`$ satisfies $`\lambda _DD`$, is of particular relevance for approaching the semiclassical limit of quantum transport and for investigations of quantum signatures of classical chaos.gutz91 ; stoe99 ; hell84  Both of these regimes pose considerable difficulties for a numerical treatment. Methods based on the expansion of the scattering wavefunction in plane or spherical waves become invalid at high fields since diamagnetic contributions are generally neglected.yang95  Methods employing a discretization on a grid are limited by the constraint that the magnetic flux per unit cell must be small, which, in turn, requires high grid densities for large $`B`$.mack83  The same requirement has to be met for high $`E_F`$, where many grid points are needed to accurately describe the continuum limit. This implies however a large number of inversions of high-dimensional matrices and therefore limits the applicability of these approaches for large $`B`$ and (or) large $`E_F`$.
In the current communication we propose an approach that allows accurate treatment of these regimes. We present an extension of the previouslyrott00  introduced Modular Recursive Green’s function Method (MRGM) to include an additional magnetic field perpendicular to the two-dimensional scattering surface. The underlying idea for our approach goes back to Sols et al.sols89  and to the widely used Recursive Green’s Function Method (RGM).mack83 ; ferr97  In the standard RGM the Green’s function is propagated through the scattering region from one transverse strip to the next by repeated solutions of a matrix Dyson equation. We show that the efficiency of this conventional discretization can be increased considerably by taking the symmetries of a scattering problem into account from the outset. Specifically, when the two-dimensional nonseparable open quantum dot can be built up out of simpler separable substructures (referred to in the following as modules), the Green’s functions for each of these modules can be calculated efficiently and virtually exactly. Calculation of the $`S`$ matrix and of the scattering wavefunction is thus reduced to “welding” together the modules by a very small number of recursions. Key to this approach are tight-binding grids which are symmetry-adapted for each module. This leads to the separability of the eigenfunctions in the modules and allows an efficient incorporation of boundary conditions. As a result, a much higher grid density can be easily handled, which, in turn, is prerequisite for treating short magnetic lengths $`l_B`$ and short wave lengths $`\lambda _D`$. Matrix Dyson equations have to be solved only for each junction between the modules. The total number of necessary recursions (i.e. high-dimensional matrix inversions) is thereby reduced to the number of modules needed to reconstruct the quantum dot.
The efficiency of the MRGM will be demonstrated by applying it to transport through a circular and a stadium shaped quantum dot. These systems are known as prototype structures for regular and chaotic dynamics and have been studied thoroughly in the literature.gutz91 ; stoe99 ; hell84 ; marc92  Concerning the theoretical approaches for the investigation of electron dynamics in quantum dots, considerable attention has been dedicated to reach higher energiessimo97 ; bies01 ; zozo96 ; akis97  and higher magnetic fields.zozo96 ; lent91 ; ji95 ; chris98 ; horn02  Especially for the study of transport through open stadium billiards, several different methods have been employed.naka92 ; naka94 ; wang94 ; yang95 ; akis97 ; nazm02  In the following we will present numerical results obtained by the MRGM which attain a parameter range, to our knowledge not yet explored by other approaches. For small $`\lambda _D`$ we investigate the localization of the scattering wavefunction near classical scattering trajectories. Characteristic differences in the dynamics of generically regular and chaotic systems will be highlighted. In the high magnetic field regime, which is governed by edge states, differences between the dynamics in different geometries disappear and are replaced by universal quasi-periodic conductance oscillations. At a critical magnetic field these oscillations break off and transport terminates entirely. In the regime where more than one edge state is excited in the dot, we find interference fluctuations which we analyze in terms of a multi-channel Fano interference model.kim01  The key to the understanding of the observed fluctuations is that inter-channel scattering between different edge states takes place only by diffractive scattering at the lead junctions.
This paper is organized as follows. In section II the method for inclusion of a magnetic field in the MRGM is presented. Section III is dedicated to a discussion of numerical results, illustrating the high magnetic field and high energy behaviour in quantum dots. The paper concludes with a short summary in section IV.

## II METHOD

We consider ballistic nanostructures with a constant electrostatic potential inside the two-dimensional cavity, impose hard-wall boundary conditions, and assume a constant magnetic field to be oriented perpendicular to the scattering plane. The shape of the quantum dot will be chosen to be either a circle or a stadium (see Figs. 5,6 below), which represent prototype systems for regular and chaotic classical dynamics, respectively. Two semi-infinite waveguides of width $`d`$ at different electrochemical potentials ($`\mu _1,\mu _2`$) are attached. The aperture of the leads is chosen to be very small $`d/D=d/\sqrt{A^{\mathrm{dot}}}=0.0935`$, where $`A^{\mathrm{dot}}=4+\pi `$ is the scaled area of all the cavities studied and $`D`$ is a characteristic linear dimension of the cavity. At asymptotic distances, i.e. far away from the quantum dot, scattering boundary conditions are imposed. The asymptotic scattering state can be factorized into a longitudinal flux-carrying plane wave and a transverse standing wave. The latter is a simple sine wave in the field-free case and a combination of Kummer functions when the magnetic field is turned on.peet88 ; schu90  In our local coordinate system the longitudinal (transverse) direction in the $`i`$-th lead is always denoted by $`x_i`$ $`(y_i)`$. The wavefunctions in the waveguides thus vanish at $`y_i=\pm d/2`$. Atomic units $`(\mathrm{}=|e|=m_{\mathrm{eff}}=1)`$ will be used from now on, unless explicitly stated otherwise.

### II.1 Brief review of the MRGM

In order to highlight the technical difficulties in incorporating a magnetic field we start by briefly reviewing the MRGM for the field-free case. Starting point is the observation that a large class of dot geometries with non-separable boundaries can be decomposed into separable two-dimensional substructures, referred to in the following as modules. For each of these modules the discretization of the corresponding tight-binding (tb) Hamiltonian can be performed on a symmetry-adapted grid. The grid for each module is chosen such that the eigenfunctions of the tb Hamiltonian

$$\widehat{H}^{\mathrm{tb}}=\underset{i}{}\epsilon _i|ii|+\underset{i,j}{}V_{i,j}|ij|$$
(1)
separate into two generalized coordinates. $`\widehat{H}^{\mathrm{tb}}`$ contains hopping potentials $`V_{i,j}`$ for nearest-neighbour coupling and site energies $`\epsilon _i`$. Both quantities are chosen such that the Schrödinger equation, $`\widehat{H}^{\mathrm{tb}}|\psi _m=E_m|\psi _m`$, converges towards the continuum Schrödinger equation in the limit of high grid point density. The most straightforward application of this approach refers to modules for which the boundaries are nodal lines of Cartesian ($`x,y`$) or polar coordinates ($`\varrho ,\phi `$). For these coordinate systems we have rott00  at $`B=0`$

$`\begin{array}{ccccccccc}\hfill V_{i,i\pm 1}^x& =& \frac{1}{2\mathrm{\Delta }x^2},\hfill & \hfill V_{j,j\pm 1}^y& =& \frac{1}{2\mathrm{\Delta }y^2},\hfill & \hfill \epsilon _i& =& \frac{1}{\mathrm{\Delta }x^2}+\frac{1}{\mathrm{\Delta }y^2},\hfill \\ \hfill V_{i,i\pm 1}^\varrho & =& \frac{\varrho _{i\pm 1/2}}{2\varrho _i\mathrm{\Delta }\varrho ^2},\hfill & \hfill V_{j,j\pm 1}^\phi & =& \frac{1}{2\varrho _i^2\mathrm{\Delta }\phi ^2},\hfill & \hfill \epsilon _i& =& \frac{1}{\mathrm{\Delta }\varrho ^2}+\frac{1}{\varrho _i^2\mathrm{\Delta }\phi ^2},\hfill \end{array}`$ (4)
with $`\varrho _i=|i1/2|\mathrm{\Delta }\varrho `$. For separable energy eigenfunctions of the general form $`|E_m=|E_k|E_{k,n}`$ the spectral representation of the retarded $`(+)`$ and advanced $`()`$ Green’s function $`G^\pm (𝐫,𝐫^{},B,E_F)`$ of the module is simply given by

$$G^\pm (𝐫,𝐫^{},B,E_F)=\underset{k}{}\alpha |E_kE_k|\alpha ^{}\underset{n}{}\frac{\beta |E_{kn}E_{kn}|\beta ^{}}{E_F\pm iϵE_{kn}}.$$
(5)
where $`(\alpha ,\beta )`$ stand for the (generalized) coordinates $`(x,y)`$ or $`(\varphi ,\rho )`$. The indices $`(k,n)`$ represent the quantum numbers of the separable eigenfunctions $`|E_k,|E_{k,m}`$ associated with the degrees of freedom $`\alpha `$ and $`\beta `$ respectively.
The Green’s functions of the separate modules are joined by solving a matrix Dyson equation,

$$G=G^0+G^0\overline{V}G,$$
(6)
where $`G^0`$ and $`G`$ denote Green’s functions of the disconnected and the connected modules, respectively. The matrix $`\overline{V}`$ denotes the hopping potential $`V`$ multiplied by the size of the unit cell $`\overline{V}=V\mathrm{\Delta }_R`$, which in a Cartesian (polar) grid is $`\mathrm{\Delta }_R=\mathrm{\Delta }x\mathrm{\Delta }y(=\varrho _i\mathrm{\Delta }\varrho \mathrm{\Delta }\phi )`$. The complete scattering structure can thus be assembled from the individual modules (much like a jigsaw puzzle). The number of necessary recursions \[i.e. solutions of (6)\] is (approximately) equal to the number of modules. The exact number depends on the number of link modules required for different grid structures. For example, in order to connect a half-circle with a rectangle we need one additional link module which is plugged in between \[see Ref. rott00,  for details\]. The key property of these link modules is their adaption to two grid symmetries \[see Fig. 3b\]. Mathematically speaking, the transition from a polar to a Cartesian grid requires a link module in order to preserve the hermiticity of the tb Hamiltonian at the junction. In the recursion the link module is connected to the Cartesian (polar) grid by means of the hopping potential $`\overline{V}^x(\overline{V}^\phi )`$, respectively.
Once the Green’s function $`G^+`$ for the combination of all modules is assembled, the transmission amplitudes $`t_{nm}`$ from entrance lead mode $`m`$ into exit lead mode $`n`$ can be calculated by projecting $`G^+`$ onto the transverse wavefunctions in the leads $`\chi _n(y_i)`$. With the corresponding longitudinal wave numbers $`k_{x_i,n}`$ we have (at zero magnetic field),

$$t_{nm}(E_F)=i\sqrt{k_{x_2,n}k_{x_1,m}}_{d/2}^{d/2}𝑑y_2_{d/2}^{d/2}𝑑y_1\chi _n^{}(y_2)G^+(y_2,y_1,E_F)\chi _m(y_1).$$
(7)
Together with the reflection amplitudes $`r_{nm}`$ (for which an analogous relation holds) the $`S`$-matrix is completely determined and satisfies the unitarity condition implied by current conservation,

$$\underset{n=1}{\overset{M}{}}(|t_{nm}|^2+|r_{nm}|^2)=1.$$
(8)
The integer $`M`$ denotes the number of open channels in the leads. According to the Landauer formula, the total conductance $`g`$ through the quantum dot is given by

$$g=\frac{1}{\pi }\underset{m,n=1}{\overset{M}{}}|t_{nm}|^2=\frac{1}{\pi }T^{\mathrm{tot}}\mathrm{with}T^{\mathrm{tot}}+R^{\mathrm{tot}}=M.$$
(9)
### II.2 Inclusion of the magnetic field

Incorporation of the magnetic field into the MRGM poses a number of complications. The solutions of these difficulties will be presented in this section. At the core of the problem is the preservation of separability of the Schrödinger equation. The usage of gauge transformations as well as of Dyson equations for decomposing non-separable structures into separable substructures plays a key role in accomplishing this goal. The field $`𝐁=(0,0,B)`$ enters the tb Hamiltonian (1) by means of a Peierls phase factor,ferr97 ; peie33

$$V_{𝐫,𝐫^{}}V_{𝐫,𝐫^{}}\mathrm{exp}\left[i/c_𝐫^𝐫^{}𝐀(𝐫)\mathrm{𝐝𝐫}\right],$$
(10)
with which the field-free hopping potential $`V_{𝐫,𝐫^{}}`$ is multiplied. The vector potential $`𝐀(𝐫)`$ satisfies $`\times 𝐀(𝐫)=𝐁`$. The Peierls phase will, of course, in most cases destroy the separability of the eigenfunctions of $`\widehat{H}^{\mathrm{tb}}`$. The difficulties can be, in part, circumvented by exploiting the gauge freedom of the vector potential, i.e.,

$$𝐀𝐀^{}=𝐀+\lambda ,$$
(11)
where $`\lambda (𝐫)`$ is a scalar function. By an appropriate choice of $`\lambda `$ the wavefunction remains separable on a given symmetry adapted grid. Specifically, to preserve separability we employ the Landau gauge for a Cartesian grid

$$𝐀=(By,0,0),$$
(12a)
and the “symmetric” or circular gauge for a polar grid
$$𝐀=B/2(y,x,0)=\rho \times 𝐁/2.$$
(12b)
The scalar gauge potential generating the gauge transformation from (12a) to (12b) is $`\lambda (x,y)=Bxy/2`$.
A major complication results from the fact, that in the presence of the magnetic field the separability on an unrestricted grid of a given symmetry does not imply the separability in the presence of boundary conditions of the same symmetry. We illustrate this problem with the help of one typical example, the semi-infinite quantum wire with lead width $`d`$ (Fig. 1). We impose hard-wall boundary conditions $`\psi (x,y=\pm d/2)=0`$ and consider first the infinite quantum wire along the $`x`$ direction. Because of the Cartesian boundary conditions, the symmetry adapted gauge is the Landau gauge $`𝐀=By\widehat{𝐱}`$. Consider, for notational simplicity, the Schrödinger equation in the continuum limit,

$$H\varphi (y,x)=\frac{1}{2}\left(𝐩+\frac{1}{c}𝐀\right)^2\varphi (x,y)=\frac{1}{2}\left(\frac{^2}{x^2}\frac{^2}{y^2}\frac{i2B}{c}y\frac{}{x}+\frac{B^2y^2}{c^2}\right)\varphi (x,y)=E_F\varphi (x,y).$$
(13)
Since the longitudinal momentum $`p_x=i/x`$ commutes with $`H`$, the separability of the wavefunction persists in the presence of the magnetic field: $`\varphi (x,y)=f_k(x)\chi (y)`$ with $`f_k(x)=e^{ikx}`$. If, however, one introduces an additional Cartesian boundary condition along the $`y`$-axis \[i.e. $`\psi (x=0,y)=0`$ for a semi-infinite lead\] the situation changes. In the absence of the magnetic field, $`B=0`$, the linear term in $`p_x`$ vanishes and thus the choice $`f(x)=\mathrm{sin}(kx)`$ \[i.e. a linear combination of $`f_{\pm k}(x)`$\] satisfies the boundary condition and preserves the separability, even though $`p_x`$ is no longer conserved in the semi-infinite lead. However, for $`B0`$ and the same boundary condition $`\psi (x=0,y)=0`$, the term linear in $`B`$ and $`p_x`$ destroys the separability. The wavefunction takes now the general form

$$\varphi (x,y)=\underset{m}{}e^{ik_mx}\underset{n}{}c_{mn}\chi _{mn}(y).$$
(14)
The breakdown of separability by the introduction of an additional boundary condition indicates that the Green’s function of confined modules will be more complex than for extended systems for the same symmetry adapted grid and the same gauge. Therefore, the program of the modular method of building-up extended complex structures by “welding together” smaller modules of higher symmetry will be executed in reverse: non-separable confined modules will be generated by “cutting in pieces” larger separable modules. Confining boundary conditions will be introduced rather than removed by the matrix Dyson equation. In the example above, the semi-infinite quantum wire is generated by cutting the infinite wire at the line $`x=0`$, thereby imposing the additional boundary condition. Just as connecting modules, so is disconnecting a given module equivalent to the application of a matrix Dyson equation,

$$G^E=G^C+G^C\overline{V}G^E.$$
(15)
In this context $`G^E`$ ($`G^C`$) is the Green’s function of the extended (confined) module and $`\overline{V}`$ is the hopping potential that connects the modules. Solving (15) in reversed mode (i.e. for $`G^C`$ rather than for $`G^E`$) amounts to dissecting the larger module.
Provided that the Green’s functions of all the necessary modules are available, we have to link them with each other to assemble the entire scattering geometry. However, in the presence of a magnetic field we have to take into account that the different modules will be calculated in different symmetry-adapted gauges. Joining different modules requires, therefore, in general a gauge transformation. For the Green’s function on the grid $`G(𝐫_i,𝐫_j^{})`$ this transformation is simplified by the fact that the matrix of gauge transformations

$$\left[\mathrm{\Lambda }(𝐫_j)\right]_{jk}=\mathrm{exp}\left[i\lambda (𝐫_j)/c\right]\delta _{jk}$$
(16)
is diagonal in the grid representation. Correspondingly, the transformation of both the hopping potential $`\overline{V}`$ and the Green’s function is local, i.e.

$$\overline{V}(𝐫_i,𝐫_j^{})\overline{V}^{}(𝐫_i,𝐫_j^{})=\mathrm{\Lambda }(𝐫_i)\overline{V}(𝐫_i,𝐫_j^{})\mathrm{\Lambda }^{}(𝐫_j^{})$$
(17)
$$G(𝐫_i,𝐫_j^{})G^{}(𝐫_i,𝐫_j^{})=\mathrm{\Lambda }(𝐫_i)G(𝐫_i,𝐫_j^{})\mathrm{\Lambda }^{}(𝐫_j^{}).$$

It is thus not necessary to transform gauges of different modules to one global gauge. Instead, it is sufficient to perform a local gauge transformation at the points of the junctions $`\{𝐫_i\}`$, such that the gauges of the two modules to be joined agree at these points.
Finally, in order to extract the $`S`$-matrix, i.e. the amplitudes $`t_{nm}`$ and $`r_{nm}`$, matrix elements of the current operator must be of gauge-invariant form. This requirement can be fulfilled by employing a double-sided gradient operator which is defined asbara89

$$f\stackrel{}{𝐃}g=f(𝐱)𝐃g(𝐱)g(𝐱)𝐃^{}f(𝐱)=g\stackrel{}{𝐃}f\mathrm{with}𝐃=\frac{i}{c}𝐀(𝐱).$$
(18)
With its help the transmission amplitudes can be evaluated asbara89 ; ando91 ; zozo96

$`t_{nm}(E_F,B)`$ $`=`$ $`{\displaystyle \frac{i}{4\sqrt{\theta _n\theta _m}}}{\displaystyle _{d/2}^{d/2}}𝑑y_2{\displaystyle _{d/2}^{d/2}}𝑑y_1^{}\chi _n^{}(y_2)e^{ik_nx_2}(\stackrel{}{𝐃}\widehat{𝐱}_2)^{}`$ (19)
$`G^+(𝐱_2,𝐱_1^{},E_F,B)(\stackrel{}{𝐃}{}_{}{}^{}\widehat{𝐱}_1^{})\chi _m(y_1^{})e^{ik_mx_1^{}}.`$
The unit vectors $`\widehat{𝐱}_n`$ are assumed to be pointing in outward direction of the $`n`$-th lead and $`\theta _m`$ denotes the outgoing particle flux carried by $`\chi _m(y_1^{})e^{ik_mx_1^{}}`$ through the lead cross section. Determination of transverse states $`\chi _m(y_i)`$ and of the corresponding longitudinal momentum $`k_m`$ as well as the normalization factors $`\theta _m`$ will be discussed below. For reflection amplitudes $`r_{nm}`$, a relation similar to (19) holds.bara89  From $`t_{nm}`$ and $`r_{nm}`$ the conductance can be calculated by means of the Landauer formula \[Eq. (9)\].

### II.3 Calculation of modules

This section is dedicated to the evaluation of the Green’s functions for those modules, which we need to assemble a circle and a stadium billiard: the half-infinite leads, the rectangle, the circle, and the half-circle. With the exception of the circle, for all these modules Eq. (5) is not applicable. This is due to the non-separability for confined geometries as discussed above. Moreover the spectrum in open structures like the semi-infinite lead is continuous rather than discrete. Unlike in the field-free case,econ79  the resulting integrals cannot be calculated analytically. However, both problems can be overcome by applying the matrix Dyson equation in a non-standard way.

#### II.3.1 Rectangular module

As illustrated above for the semi-infinite wave guide, the Dirichlet boundary condition for the confined structure of a rectangle with magnetic field is not separable, no matter which gauge is chosen. The separability can however be restored by imposing periodic boundary conditions on two opposing sides of the rectangle. Topologically, this corresponds to folding the rectangle to the surface of a cylinder (Fig. 2a). In this case we connect the first ($`P`$) and the last ($`Q`$) transverse grid slice of a rectangular grid by a hopping potential $`|V_{PQ}^x|=|V_{QP}^x|=\frac{1}{2\mathrm{\Delta }x^2}`$. The Green’s function of this “cylinder surface” $`(cs)`$ will be denoted by $`G^{\mathrm{cs}}`$ in the following. The calculation of the rectangle Green’s function $`G^\mathrm{r}`$ will be obtained out of $`G^{\mathrm{cs}}`$ by a Dyson equation used here in “reversed” mode, i.e. for disconnecting tb grids. This method for calculating the rectangular module may seem like a detour, but it is numerically more efficient than a strip-by-strip recursion. For completeness we mention that an alternative way to calculate $`G^\mathrm{r}`$ was proposed in Ref. skja94, .
The Green’s function for the cylinder surface $`G^{\mathrm{cs}}`$ can be constructed from separable eigenfunctions, $`|E_m=|E_k^x|E_{kn}^y`$, according to Eq. (5). Solving the tight-binding Schrödinger equation for the cylinder surface, we obtain for the longitudinal eigenstates $`x_j|E_k^x=(N_x\mathrm{\Delta }x)^{1/2}\mathrm{exp}(i2\pi kj/N_x)`$, which results in a tridiagonal, symmetric matrix-eigenproblem of size $`N_y\times N_y`$ for the transverse modes,numr86

$`E_{kn}y_l|E_{kn}^y={\displaystyle \frac{1}{\mathrm{\Delta }x^2}}\left[\mathrm{cos}\left({\displaystyle \frac{2\pi k}{N_x}}+{\displaystyle \frac{B}{c}}y_l\mathrm{\Delta }x\right)1\right]`$
$`\times y_l|E_{kn}^y{\displaystyle \frac{1}{2\mathrm{\Delta }y^2}}\left(y_{l1}|E_{kn}^y2y_l|E_{kn}^y+y_{l+1}|E_{kn}^y\right).`$ (20)
By “cutting the cylinder surface open” along a line of constant $`x`$, we obtain from $`G^{\mathrm{cs}}`$ the desired Green’s function $`G^\mathrm{r}`$ for the rectangle (Fig. 2a). We demonstrate this for the rectangle Green’s function $`G_{PX}^\mathrm{r}`$ from the first transverse slice $`P`$ to any other slice $`X`$. To determine $`G_{PX}^\mathrm{r}`$ we solve the following system of Dyson equations,

$`G_{PX}^\mathrm{r}`$ $`=`$ $`G_{PX}^{\mathrm{cs}}G_{PQ}^\mathrm{r}\overline{V}_{QP}G_{PX}^{\mathrm{cs}}G_{PP}^\mathrm{r}\overline{V}_{PQ}G_{QX}^{\mathrm{cs}}`$ (21)
$`G_{PQ}^{\mathrm{cs}}`$ $`=`$ $`G_{PQ}^\mathrm{r}+G_{PQ}^\mathrm{r}\overline{V}_{QP}G_{PQ}^{\mathrm{cs}}+G_{PP}^\mathrm{r}\overline{V}_{PQ}G_{QQ}^{\mathrm{cs}}`$ (22)
$`G_{PP}^{\mathrm{cs}}`$ $`=`$ $`G_{PP}^\mathrm{r}+G_{PQ}^\mathrm{r}\overline{V}_{QP}G_{PP}^{\mathrm{cs}}+G_{PP}^\mathrm{r}\overline{V}_{PQ}G_{QP}^{\mathrm{cs}},`$ (23)
where the first line is the “reversed” Dyson equation. The three unknowns in the above equations are the Green’s functions connecting the slices $`(P,X),(P,Q)`$ and $`(P,P),`$ $`G_{PX}^\mathrm{r},G_{PQ}^\mathrm{r},G_{PP}^\mathrm{r}`$. By solving these three equations, the unknowns can be uniquely determined.

#### II.3.2 Circle and half-circle

In symmetric gauge, $`𝐀=B/2(y,x,0)`$, the Dirichlet boundary value problem for the circle with magnetic field is separable, $`|E_m=|E_k^\phi |E_{kn}^\varrho `$. On a discrete tb lattice this statement remains true, provided that a circular grid is employed. With the eigenstates for the azimuthal degree of freedom, $`\phi _j|E_k^\phi =(N_\phi \mathrm{\Delta }\phi )^{1/2}\mathrm{exp}(i2\pi kj/N_\phi )`$ and radial eigenstates $`g_{kn}(\varrho _i)=\sqrt{\varrho _i}\times \varrho _i|E_{kn}^\varrho `$, the finite difference equation for the $`g_{kn}(\varrho _i)`$ results in a tridiagonal symmetric eigenproblem

$`E_{kn}g_{kn}(\varrho _i)`$ $`=`$ $`{\displaystyle \frac{1}{\varrho _i^2\mathrm{\Delta }\phi ^2}}\left[\mathrm{cos}\left({\displaystyle \frac{2k\pi }{N_\phi }}{\displaystyle \frac{\varrho _i^2B\mathrm{\Delta }\phi }{2c}}\right)1\right]g_{kn}(\varrho _i)`$ (24)
$`{\displaystyle \frac{1}{2\mathrm{\Delta }\varrho ^2}}\left[{\displaystyle \frac{\varrho _{i1/2}}{\sqrt{\varrho _{i1}}\sqrt{\varrho _i}}}g_{kn}(\varrho _{i1})2g_{kn}(\varrho _i)+{\displaystyle \frac{\varrho _{i+1/2}}{\sqrt{\varrho _i}\sqrt{\varrho _{i+1}}}}g_{kn}(\varrho _{i+1})\right].`$
The Green’s function for the circular module is then calculated by a straight-forward application of Eq. (2.3).
For the Green’s function of the half-circle we employ an analogous procedure as in the previous subsection: we dissect the circle Green’s function into half-circles by means of a “reversed” Dyson equation. We demonstrate this by way of the example depicted in Fig. 2b, where the “full circle” (fc) is split up into two “half-circles” (hc). The resulting two halves are almost identical, with the exception of the two additional radial grid slices, by which the right half-circle is larger. For assembling the stadium billiard we have to make sure that the tb grid of the half-circle module can be linked directly to a vertical grid \[see Fig. 3b\]. For this reason, only the left one of the two half-circles in Fig. 2b can be used for this purpose.
Consider as example the Green’s function $`G_{PX}^{\mathrm{hc}}`$ describing the propagation from the grid slice $`P`$ at the junction of the two half-circles to any radial grid slice $`X`$ situated on the “left half-circle” (see Fig. 2b). $`G_{PX}^{\mathrm{hc}}`$ is determined by the following system of Dyson equations

$`G_{PX}^{\mathrm{hc}}`$ $`=`$ $`G_{PX}^{\mathrm{fc}}G_{PP}^{\mathrm{hc}}\overline{V}_{PQ}G_{QX}^{\mathrm{fc}}`$ (25)
$`G_{PP}^{\mathrm{fc}}`$ $`=`$ $`G_{PP}^{\mathrm{hc}}+G_{PP}^{\mathrm{hc}}\overline{V}_{PQ}G_{QP}^{\mathrm{fc}}`$ (26)
which yields a unique solution for $`G_{PX}^{\mathrm{hc}}`$.

#### II.3.3 Semi-infinite lead

Because of its continuous spectrum, the Green’s function for the semi-infinite lead poses an additional challenge beyond that of the non-separability of the wavefunction discussed above. We therefore apply one further “trick” to bypass this problem. Our approach is based on the observation that adding a slice to a semi-infinite quantum wire leaves this wire (up to irrelevant phases) invariant (see Fig. 3a). We assume a semi-infinite lead with $`x[\mathrm{\Delta }x,\mathrm{})`$ and hard-wall boundary conditions at $`x=\mathrm{\Delta }x`$ and $`y=\pm d/2`$. To this object we add a slice consisting of just one transverse chain of tb grid points which we place at $`x=0`$. The system of Green’s functions for the propagation from the transverse chain at $`x=0`$ ($`P`$) back to itself $`(P)`$ or to the first transverse slice of the semi-infinite lead ($`Q`$) at $`\mathrm{\Delta }x`$ reads

$`G_{PP}=G_{PP}^0+G_{PP}^0\overline{V}_{PQ}G_{QP},`$ (27)
$`G_{QP}=G_{QQ}^0\overline{V}_{QP}G_{PP}.`$ (28)
Each multiplication involves a matrix product with a dimension equal to the number of transverse grid points. The key point is now that the system of Eqs. (27,28) can be closed through the invariance condition (Fig. 3a) for the semi-infinite lead, i.e. $`G_{PP}=G_{QQ}^0`$. In Landau gauge $`𝐀=(By,0,0)`$ the latter relation does not involve additional gauge phases since these are already contained in the hopping matrix element. We additionally note that an equivalent point of departure for the derivation of $`G_{PP}`$ is the Bloch condition for states in the lead.ando91 ; zozo96
Setting $`Z=G_{PP}\overline{V}_{QP}`$ and using the hermiticity condition $`\overline{V}_{QP}=\overline{V}_{PQ}^{}\overline{V}^{}`$, Eqs. (27,28) can be converted to a quadratic matrix equation

$$ZZ\overline{V}^1(G_{PP}^0)^1Z+\overline{V}^1\overline{V}^{}=0.$$
(29)
Solvents $`Z`$ of a quadratic matrix equation $`Q(Z)=0`$ can be constructed from the eigenpairs $`(\beta _i,\chi _i)`$ of the corresponding quadratic eigenvalue equation $`Q(\beta _i)\chi _i=0,i[1,\mathrm{},2N]`$ in the diagonal form,high00

$$Z=MBM^1\mathrm{with}M=[\chi _1,\mathrm{},\chi _N],B=\mathrm{diag}(\beta _i).$$
(30)
The quadratic eigenvalue equation is equivalent to a generalized eigenvalue problem $`A\stackrel{~}{\chi }=\beta C\stackrel{~}{\chi }`$ of twice the original dimension.numr86  Its $`2N`$ dimensional eigenvectors $`\stackrel{~}{\chi }=(\chi ,\beta \chi )`$ are solutions of the symmetric eigenproblem

$$\left(\begin{array}{cc}\overline{V}^{}& 0\\ 0& \overline{V}\end{array}\right)\left(\begin{array}{c}\chi \\ \beta \chi \end{array}\right)=\beta \left(\begin{array}{cc}(G_{PP}^0)^1& \overline{V}\\ \overline{V}& 0\end{array}\right)\left(\begin{array}{c}\chi \\ \beta \chi \end{array}\right),$$
(31)
where $`\left(G_{PP}^0\right)^1=E_F\widehat{H}_{1\mathrm{D}}^{\mathrm{tb}}`$ and $`\widehat{H}_{1\mathrm{D}}^{\mathrm{tb}}`$ is the Hamiltionian of the one-dimensional transverse tb strip at $`x=0`$. The Fermi energy $`E_F`$ and the magnetic field $`B`$ enter (31) as independent parameters at which the eigenstates $`\stackrel{~}{\chi }_m`$ and eigenvalues $`\beta _m`$ are evaluated. The longitudinal momenta of the lead states $`\xi _m(x,y)=\chi _m(y)e^{ik_mx}/\sqrt{\theta _m}`$ are related to the eigenvalues by the relation $`\beta =\mathrm{exp}(ik\mathrm{\Delta }x)`$. The orthonormalisation and the completeness relations of the $`2N`$ eigenvectors $`\stackrel{~}{\chi }_m`$ can be formulated in terms of matrix relations, for the generalized eigenproblem,

$$\frac{1}{\sqrt{\theta _m\theta _n}}\stackrel{~}{\chi }_m^TC\stackrel{~}{\chi }_n=2i\frac{k_m}{|k_m|}\delta _{mn}\mathrm{and}\underset{m}{\overset{2N}{}}\frac{\stackrel{~}{\chi }_m\stackrel{~}{\chi }_m^T}{\theta _m}=2i\frac{k_m}{|k_m|}C^1.$$
(32)
With this specific choice of normalization the norm factors $`\theta _m`$ are determined such that every propagating state carries unit flux. We note parenthetically that the quadratic eigenvalue equation could also be applied to the semi-infinite lead at zero $`B`$ field. However, in that case, the Green’s function for quantum wires can be calculated analyticallyrott00 ; econ79  by complex contour integration.

#### II.3.4 Scattering wave functions and efficiency of the MRGM

The MRGM is particularly well-suited to determine transport coefficients as the Green’s function is then required only at the junctions between the modules and does not have to be evaluated throughout the interior of the entire quantum dot. Since for the calculation of the scattering wavefunction the Green’s function throughout the entire scattering region is needed, this particular advantage cannot be made use of here. However, also in this case, the MRGM is still more efficient than the standard RGM, as will be explained below.
The wavefunction $`\psi (𝐱)`$ can be obtained at any point $`𝐱`$ by projecting the retarded Green’s function (by means of the operator $`\stackrel{}{𝐃}`$) on the incoming wave (in mode m),bara89 ; zozo96

$$\psi _m(𝐱)=\frac{1}{2\sqrt{\theta _m}}_{d/2}^{d/2}dy_1^{}G^+(𝐱,𝐱^{},E_F,B)(\stackrel{}{𝐃}{}_{}{}^{}\widehat{𝐱}_1^{})\chi _m(y_1^{})e^{ik_mx_1^{}}.$$
(33)
$`G^+`$ contains the solution of the Dyson equations for all linked modules. That the evaluation of Eq. (33) can be done very efficiently results from two properties: First the number of recursions (i.e., of matrix inversions) needed to obtain $`G^+`$ is given by the fixed number of modules required to build up the scattering structure. This number is independent of the De Broglie wavelength. The latter enters only in terms of the size of the matrices involved in the recursion, since with increasing $`E_F`$ (decreasing $`\lambda _D`$) more grid points are required to represent the continuum limit. Compared to the standard RGM the numerical effort is reduced since in that approach the number of recursions scales with the grid density, i.e. $`k_F`$. A second advantage of the MRGM concerns the incorporation of the boundary conditions. In the modular method the boundaries follow by construction the nodal lines of the symmetry-adapted coordinate system for the module. Due to this reason the convergence towards the continuum limit is enhanced as compared to the slice-by-slice recursion. The calculation of the transport coefficients as a function of the Fermi wavenumber $`k_F`$ (or Fermi energy $`E_F`$) is simplified by the fact that the solution of the eigenvalue problem ($`|E_m,E_m`$) entering the Green’s function for each module \[Eq. (5)\] is independent of $`E_F`$. For the evaluation of the Green’s function at different values of $`E_F`$ the eigenproblem $`\widehat{H}^{\mathrm{tb}}|E_m=E_m|E_m`$ therefore has to be solved only once. Unfortunately, this feature does not extend to the variation of the magnetic field since both $`|E_m`$ and $`E_m`$ are dependent on $`B`$. Because of this property a new solution of the tb eigenproblem is required for each value of the field. The most severe restriction of the MRGM is, however, that its applicability is limited to those scattering structures which can be assembled from or cut out of separable modules. Also random potentials and soft walls can only be included as long as they preserve the separability of each module. We mention at this point, that a “hybrid RGM” for dealing with arbitrary boundary geometries was presented in the literature.zozo96

## III Numerical results for high $`B`$ fields and large $`k_F`$

In this section we present first magnetoconductance results which were calculated within the MRGM at high magnetic fields $`B`$ and large Fermi wavenumbers $`k_F`$. As prototype shapes of the cavity we use the circle and the Bunimovich stadium and consider different geometries for the attached quantum wires.

### III.1 Accuracy checks

Several checks for the accuracy of the numerical results have been performed. Exact relationships for transport coefficients such as conservation of unitarity and the Onsager relations are fulfilled with an accuracy of better than $`10^{10}`$. The grid density is chosen such that the magnetic flux per unit cell is $`B\mathrm{\Delta }_r/c<0.01`$ (as in Ref. bara91, ). Moreover, the typical number of grid points per Fermi half-wavelength is greater than $`30`$. Only for very high energy calculations (Fig. 5) the relative grid density is lower. For low magnetic fields, we can compare our results for $`|t_{nm}(k_F)|^2`$ with previous methods. As an example we show in Fig. 4 a comparison for $`|t_{11}(k_F)|^2`$ with the calculation by Yang et al.,yang95  which is based on a wave function expansion in spherical waves. The agreement for the circle is very good although diamagnetic terms are neglected in the approach of Ref. yang95, . For the stadium, the differences between the two methods are somewhat larger. This is due to the fact that the expansion of the stadium wave functions in spherical waves leads to a unitarity deficiency (see Fig. 4). We can also reproduce previous results of Ref. yang95,  concerning weak-localization line shapes and statistical magnetoconductance properties in chaotic and regular cavities. These will not be treated again. Our focus will be on the high magnetic field and high energy regime, where other methods failed.

### III.2 Wavefunctions

The starting point for the analysis of the scattering states $`\psi (𝐱)`$ for ballistic transport through quantum dots is Eq. (33). Figures 5 and 6 display the resulting electron density $`|\psi (𝐫)|^2`$ in the scattering region. In Fig. 5 we consider the wavefunctions at very high $`k_F`$ for both the circle and the stadium billiard, which are prototypical structures for regular and chaotic dynamics respectively. Large $`k_F`$ corresponds to the regime where the convergence towards classical scattering trajectories is expected to emerge. Figures 5a and 5b illustrate the different dynamics for an injection at high ($`m=20`$) and at low mode number ($`m=1`$), respectively. Since high mode numbers correspond classically to a large injection angle, the wavefunction condenses around a pentagon-shaped whispering gallery trajectory. For low-mode injection, a small circle representing the centrifugal barrier (or caustic) is seen, as well as rays representing the asterisk orbits.ishi95  Figures 5c and 5d display scattering states for the stadium. At low magnetic fields, the dynamics is chaotic and a typical wavefunction features a quasi-random pattern with a modest density enhancement near classically unstable periodic orbits (not shown). For special configurations of $`k_F`$ and $`B`$ “scars” emerge in the scattering wavefunctions (Fig. 5c). By contrast, for high magnetic fields the classical motion becomes regular. In the present example the wavefunction condenses around a “bundle” of cyclotron orbits executing three bounces at the cavity wall before exiting by the entrance lead (Fig. 5d). There has been an extensive discussion in the literature as to the existence of scars in open quantum billiards.akis97  Our present results clearly underscore that scars, defined here as the condensation of the wavefunction near classical (not necessarily unstable) trajectories, clearly exist for large $`k_F`$. Figure 6 illustrates the formation of edge states at high fields. With increasing $`B`$ fewer edge states can be excited in the cavity. In Fig. 6c ($`B=68.5`$) three transverse edge states are present while in Fig. 6d ($`B=125`$) only a single edge state remains. For two edge states carrying flux across the quantum dot, interferences give rise to a stationary nodal pattern with a fixed number of antinodes along the boundary (see Fig. 6a,b). We note that we are not aware of any other method that has so far been capable of investigating scattering states of open structures in this high-magnetic field regime.

### III.3 Transport coefficients

The interference between different edge states gives rise to characteristic fluctuations in the transport coefficients. Figure 7 shows the high magnetic field regime of the transmission probability in the first mode $`m=n=1`$ for both circle and stadium. Different orientations of the exiting quantum wire were chosen (oriented $`90^{}`$ and $`180^{}`$ relative to the incoming lead). A few universal trends are easily discernible: above a certain critical value of the magnetic field (denoted by $`B_c^1`$), the strongly fluctuating transmission probability gives way to very regular oscillations in all four cases \[see insets of Fig. 7 for magnification\]. The threshold value $`B_c^1`$ and the magnetic field, at which transport is terminated (separately displayed in Fig. 8) are identical for all systems investigated. Below $`B_c^1`$ the transport signal displays “beats”, i.e. the Fourier transform of the signal is characterized by several frequencies. The “universality” (i.e. geometry independence) of these features is related to the fact that in the high magnetic field regime transport is controlled by edge states (as depicted in Fig. 6). These states play a very prominent role in the Quantum Hall effect and have been studied extensively.pran87 ; ferr97 ; davi98 ; datt95  At magnetic fields, where the magnetic length is smaller than the system dimensions, $`l_BD`$, they are the only states coupling to the quantum wire since bulk Landau states cannot be excited by the leads. The edge states shown in Figs. 6a to 6d correspond to the points in the transmission spectrum also labeled by (a) to (d) in Fig. 7. By comparison with the scattering wavefunctions (as in Fig. 6), we observe that in the magnetic field region $`B_c^n<B<B_c^{n1}`$ edge states have up to $`n1`$ transverse nodes in the direction perpendicular to the boundary. Furthermore, the number of longitudinal antinodes from entrance to exit lead (see the corresponding numbers in Fig. 6a,b) can be directly mapped onto successive maxima in Fig. 7 (see numbers there). The range of $`B`$ depicted in Fig. 7 corresponds to $`BB_c^2`$ at $`k_F=1.5\pi /d`$. The transmission spectrum becomes increasingly complex as $`B`$ is reduced or equivalently $`k_F`$ is increased (not shown).
To determine the positions of the values $`B_c^n`$ we consider the energy shift of Landau levels near the boundary. Bulk Landau levels are degenerate since their quantized energy $`E_n=(n+1/2)B/c`$ is independent of their positition in space. This degeneracy is lifted if a Landau state is placed in the vicinity of the cavity wall: with decreasing distance to this boundary the energy of the state increases. Therefore the energies of edge states associated with the quantum number $`n`$ lie above the asymptotic bulk value $`E_n`$. When the incoming electron is diffracted at the mouth of the entrance lead, only those edge states whose energy is below the Fermi energy can carry the flux. Due to the sharp edges at the junction between lead and quantum dot, all energetically accessible edge channels are populated. The blow-up of the scattering wavefunction near the lead mouth (Fig. 9) highlights the diffractive edge scattering. This is a contrast to smooth edges where states in the lead could cross the lead junction adiabatically, i.e. without changing their state of quantization.glaz90  With sharp lead junctions however all edge states with quantum number $`n^{}n`$ will participate in transport up to a magnetic field where $`E_n`$ touches the Fermi energy, $`E_n=E_F`$, i.e. at the critical magnetic fields $`B_c^n/c=E_F/(n+1/2)`$. These threshold values are indicated by the dot-dashed vertical bars in Fig. 7,8 for $`B_c^271.1`$, $`B_c^1118.4`$ and $`B_c^0355.3`$ for a lead width $`d=0.25`$ and $`k_F=1.5\pi /d`$. In our numerical data, both the position of these threshold values as well as their independence of the geometry are in excellent agreement with this prediction. The only exception is the critical magnetic field $`B_c^0`$. Its value ($`355.3`$) lies slightly above the point where the transmission spectrum ceases (at $`B351.8`$) (see Fig. 8). The reason for this deviation is the fact that the termination point of the spectrum is not determined by the magnetic field of the lowest bulk Landau level in the cavity (i.e. $`B_c^0`$), but by the highest field at which the leads still carry flux. In the leads, however, the magnetic length does not satisfy the condition $`l_Bd`$ (at $`k_F=1.5\pi /d:l_Bd/4.7`$). Contrary to the cavity, the lead wavefunctions still “feel” the constriction by the opposing wall. For this reason the threshold magnetic field values of the transverse lead states lie slightly below those of the bulk Landau levels, leading to a termination already below $`B_c^0`$.

### III.4 Multi-channel interferences

The regular oscillations above $`B_c^1`$ as well as the complex fluctuating pattern below $`B_c^1`$ can be explained by a multi-channel scattering description. This model can be viewed as a generalizationkim01 ; wirt02  of a single-channel picture.siva89 ; wees89  For this description to be applicable, the cavity of the dot has to have smooth boundaries and disorder must be absent. Under these circumstances the flux transported by edge states is conserved in the interior and changed only by diffractive scattering at lead junctions: At the junction, a fraction of the flux will exit through the lead while the remaining portion of the flux will continue to propagate along the boundary. The stationary scattering state can be viewed as the coherent superposition of repeated loops around the billiard. In order to translate this picture into an analytic expression we define amplitudes for transmission and reflection at the two lead junctions. We denote the amplitudes for transmission from transverse mode $`m`$ in the entrance lead to the edge state in the dot with quantum number $`i`$ by $`\stackrel{~}{t}_{mi}`$. The amplitudes $`\stackrel{~}{t}_{in}^{}`$ stand for transmission from edge state $`i`$ in the dot to the transverse mode $`n`$ in the exit lead. The amplitudes $`\stackrel{~}{r}_{ij}(\stackrel{~}{r}_{ij}^{})`$ describe edge state reflection at the entrance (exit) lead from mode $`i`$ to mode $`j`$. (The tilde signs serve to distinguish these amplitudes from the transport coefficients for the whole geometry $`t_{nm}`$ and $`r_{nm}`$.) We further define the following matrices

$`[\stackrel{~}{T}]_{ij}`$ $`=`$ $`\stackrel{~}{t}_{ij}e^{ik_jL_jiBA_j/c},[\stackrel{~}{T}^{}]_{ij}=\stackrel{~}{t}_{ij}^{}`$ (34)
$`[\stackrel{~}{R}^{}]_{ij}`$ $`=`$ $`\stackrel{~}{r}_{ij}^{}e^{ik_jL_j^{}iBA_j^{}/c},[\stackrel{~}{R}]_{ij}=\stackrel{~}{r}_{ij}e^{ik_jL_jiBA_j/c},`$ (35)
where $`L_j,A_j`$ ($`L_j^{},A_j^{}`$) denote the lengths $`L`$ and areas $`A`$ the edge state $`j`$ covers from entrance to exit (from exit to entrance) of the dot. The areas $`A`$ can be determined in gauge-invariant form, although the corresponding classical orbits are not necessarily periodic.wirt99  The transmission through the whole cavity $`t_{ji}=[T]_{ij}`$ is then written as a geometric series of matrices,

$`T`$ $`=`$ $`\stackrel{~}{T}(1+\stackrel{~}{R}^{}\stackrel{~}{R}(1+\stackrel{~}{R}^{}\stackrel{~}{R}(1+\mathrm{})))\stackrel{~}{T}^{}`$ (36)
$`=`$ $`\stackrel{~}{T}\left({\displaystyle \underset{i=0}{\overset{\mathrm{}}{}}}(\stackrel{~}{R}^{}\stackrel{~}{R})^i\right)\stackrel{~}{T}^{}`$
$`=`$ $`\stackrel{~}{T}(1\stackrel{~}{R}^{}\stackrel{~}{R})^1\stackrel{~}{T}^{}.`$
Equation (36) serves as a convenient starting point for the analysis of the transmission fluctuations. Consider first the regime $`B>B_c^1`$, where only the lowest transverse edge state is excited. In this case Eq. (36) reduces to its scalar versionsiva89 ; wees89

$$T^{\mathrm{tot}}=|t_{11}|^2=\frac{|\stackrel{~}{t}_{11}|^2|\stackrel{~}{t}_{11}^{}|^2}{12\mathrm{R}\mathrm{e}\left[\stackrel{~}{r}_{11}^{}\stackrel{~}{r}_{11}e^{i\gamma }\right]+|\stackrel{~}{r}_{11}^{}|^2|\stackrel{~}{r}_{11}|^2},$$
(37)
with $`\gamma =k_1(L_1+L_1^{})B(A_1+A_1^{})/c`$. As expected, the fluctuations of $`|t_{11}(B/c)|^2`$ are determined by an Aharonov-Bohm type phase $`\gamma `$. At fixed $`k_F`$, the oscillation period is $`\mathrm{\Delta }B=2\pi c/A_1^{\mathrm{tot}}`$. By $`A_1^{\mathrm{tot}}=A_1+A_1^{}`$ we denote the area which the edge state acquires in one revolution around the dot. Taking into account that the dynamically accessible area of the edge state is somewhat smaller than the geometric area, $`A_1^{\mathrm{tot}}<A^{\mathrm{dot}}=4+\pi `$ (see Fig. 6), the prediction for the oscillation period is in excellent agreement with our numerical findings. Equation (37) also explains why the oscillation period is increasing with increasing $`B`$ (see Fig. 8). This explanation makes use of the somewhat counterintuitive fact, that for increasing magnetic field skipping orbits with fixed quantum number $`n`$ have an increasing mean distance from the boundary.wees89  Consequently, a larger $`B`$ field implies a smaller enclosed area $`A_1^{\mathrm{tot}}`$ and therefore a higher oscillation period $`\mathrm{\Delta }B`$. Furthermore Eq. (37) accounts for the fact that for most structures the successive maxima of $`T^{\mathrm{tot}}`$ reach unity.siva89  (The small deviation from this rule of the stadium with $`90^{}`$ lead orientaion will be explained below). In addition to unitarity, $`[|\stackrel{~}{t}_{11}^{}|^2+|\stackrel{~}{r}_{11}^{}|^2=1]`$ we have for identical lead junctions $`\stackrel{~}{r}_{11}=\stackrel{~}{r}_{11}^{}`$. (We call two junctions identical if the local environment of their lead mouths is the same and their respective distance is larger than a few wavelenths.) Above $`B_c^1`$ scattering of an edge state at a junction is essentially a one-dimensional process, for which the probability for transmission from left to right has the same magnitude as vice versa. Identical lead junctions therefore also imply $`\stackrel{~}{t}_{11}^{}=\stackrel{~}{t}_{11}`$, provided that the two corners of the lead junction have the same shape. If and only if all of the three above conditions are fulfilled, Eq. (37) yields $`T^{\mathrm{tot}}=1`$ at the resonance condition $`\gamma =2\pi n,n`$. Since for the two circle geometries and for the $`180^{}`$-stadium the two lead junctions are identical, we indeed find in these cases that $`|t_{11}(B)|^2`$ periodically reaches unity. On the other hand, when the leads are attached to the stadium at an angle of $`90^{}`$, one lead is attached to the straight section while the other is attached to the semicircle. The local environment of the two lead mouths is in this case different (i.e. the lead junctions are not identical), for which reason our numerical results do not quite reach $`|t_{11}|^2=1`$, when the resonance condition is fulfilled for this geometry (see inset of Fig. 7d and Fig. 8). Another relation exists between the resonance condition and the behaviour of edge states. In the closed cavity an edge channel always encloses an integer number of flux quanta \[$`BA/(\varphi _0c)=m`$\]. Therefore, the resonance condition is met whenever the energy of an edge state in the closed cavity crosses the Fermi edge.wees89
One interesting feature of the transmission fluctuations in the single-channel regime of the circular dot ($`B>B_c^1`$) is their invariance with respect to the lead orientation. The numerical results for the transmission probabilities of the circle with $`180^{}`$ and $`90^{}`$ lead orientation differ only at the tenth (!) decimal digit. This fact, as well as the observation, that in the case of the stadium billard the two lead orientations give different results, can again be explained by Eq. (37). The important point to note is that the interference phase ($`\gamma =k_1L_1^{\mathrm{tot}}BA_1^{\mathrm{tot}}/c`$) does not change when changing the positions of the leads around the circle. Due to the rotational symmetry also the coefficients $`\stackrel{~}{t}_{11},\stackrel{~}{t}_{11}^{},\stackrel{~}{r}_{11}`$ and $`\stackrel{~}{r}_{11}^{}`$ are the same for different angles between the leads. The same is thus true for the total transmission $`T^{\mathrm{tot}}`$ through the circular dot. The only exception to this rule occurs when the two leads are in close proximity to within a few wavelengths. In this case the transmission probability changes as compared to the results for the $`180^{}`$ and $`90^{}`$ lead orientation (not shown).
The fluctuations in the regime $`B<B_c^1`$ displayed in Fig. 7 can be analyzed with the help of a multi-channel scattering description. In the interval $`B_c^2BB_c^1`$ two channels corresponding to two edge states are open in the cavity and one channel in each of the leads. From the entrance to the exit lead mouth the two edge channels acquire the phases $`e^{ik_1L_1iBA_1/c}`$ and $`e^{ik_2L_2iBA_2/c}`$ respectively. Their interference at the exit lead will therefore give an oscillatory contribution to the total transmission $`T^{\mathrm{tot}}(B)=|t_{11}|^2`$ of the form $`\mathrm{cos}^2[B(A_1A_2)/(2c)]`$. For closer analysis we need to evaluate Eq. (36) which involves the inversion of $`2\times 2`$ matrices. In the case of parallel lead orientation the corresponding expressions are simplified due to the fact that the phases acquired from entrance to exit lead and vice versa are the same ($`A_n=A_n^{}`$ and $`L_n=L_n^{}`$),

$`t_{11}`$ $`=`$ $`[e^{i\phi _1}\stackrel{~}{t}_{11}\stackrel{~}{t}_{11}^{}+e^{i\phi _2}\stackrel{~}{t}_{12}\stackrel{~}{t}_{21}^{}+e^{i(2\phi _1+\phi _2)}(\stackrel{~}{r}_{11}\stackrel{~}{t}_{12}\stackrel{~}{r}_{12}\stackrel{~}{t}_{11})(\stackrel{~}{r}_{21}^{}\stackrel{~}{t}_{11}^{}\stackrel{~}{r}_{11}^{}\stackrel{~}{t}_{21}^{})`$ (38)
$`+e^{i(\phi _1+2\phi _2)}(\stackrel{~}{r}_{21}\stackrel{~}{t}_{12}\stackrel{~}{r}_{22}\stackrel{~}{t}_{11})(\stackrel{~}{r}_{22}^{}\stackrel{~}{t}_{11}^{}\stackrel{~}{r}_{12}^{}\stackrel{~}{t}_{21}^{})]/[1e^{i2\phi _1}\stackrel{~}{r}_{11}\stackrel{~}{r}_{11}^{}e^{i2\phi _2}\stackrel{~}{r}_{22}\stackrel{~}{r}_{22}^{}`$
$`e^{i(\phi _1+\phi _2)}(\stackrel{~}{r}_{21}\stackrel{~}{r}_{12}^{}+\stackrel{~}{r}_{12}\stackrel{~}{r}_{21}^{})e^{i2(\phi _1+\phi _2)}(\stackrel{~}{r}_{11}\stackrel{~}{r}_{22}\stackrel{~}{r}_{12}\stackrel{~}{r}_{21})(\stackrel{~}{r}_{12}^{}\stackrel{~}{r}_{21}^{}\stackrel{~}{r}_{11}^{}\stackrel{~}{r}_{22}^{})],`$
with the abbreviated notation $`\phi _n=k_nL_nBA_n/c`$. In Fig. 10 we show one half-period of the beats in $`T^{\mathrm{tot}}(B)=|t_{11}(B)|^2`$ for \[$`\pi n<B(A_1A_2)/(2c)<\pi (n+1)`$\], as calculated with Eq. (38). The absolute square of the numerator (dashed line, $`N`$) and denominator (dotted line $`D`$) of Eq. (38) display very similar oscillations, both in frequency and amplitude. However, since $`T^{\mathrm{tot}}=N/D`$, a series of dips are superimposed on the term $`\mathrm{cos}^2[B(A_1A_2)/(2c)]`$ at the points where $`N`$ and $`D`$ have their common minima. To classify these dips (i.e. antiresonances) we make use of the fact,kim01  that the unitarity of Eq. (38) can be satisfied by mapping the transport coefficients at the lead junctions (which are assumed to be identical) onto six independent parameters $`(r,ϵ,\varphi ,\vartheta ,\varphi _1,\varphi _2)`$,

$`\stackrel{~}{t}_{11}=\stackrel{~}{t}_{11}^{}=\sqrt{(1r^2)ϵ}e^{i[(\varphi _1+\psi )/2+\vartheta ]},`$
$`\stackrel{~}{t}_{12}=\stackrel{~}{t}_{21}^{}=\sqrt{(1r^2)(1ϵ)}e^{i[(\varphi _2+\psi )/2+\vartheta ]},`$
$`\stackrel{~}{r}_{11}=\stackrel{~}{r}_{11}^{}=[(1ϵ)+ϵr]e^{i(\varphi _1+\vartheta )},`$ (39)
$`\stackrel{~}{r}_{12}=\stackrel{~}{r}_{12}^{}=\stackrel{~}{r}_{21}=\stackrel{~}{r}_{21}^{}=(1r)\sqrt{ϵ(1ϵ)}e^{i[(\varphi _1+\varphi _2)/2+\vartheta ]}`$
$`\stackrel{~}{r}_{22}=[ϵ+(1ϵ)r]e^{i(\varphi _2+\vartheta )}.`$
Out of this set of parameters, only two ($`ϵ,r`$) are physically relevant. The variable $`ϵ`$ represents the coupling of the incoming lead state to the edge state $`n=1`$ in the interior and $`r`$ is related to the reflection coefficient $`\stackrel{~}{r}_{11}^{\prime \prime }`$ of an incident channel at the lead mouth, $`\stackrel{~}{r}_{11}^{\prime \prime }=re^{i(\psi +\vartheta )}`$. Both quantities $`ϵ,|r|`$ are restricted to the interval $`[0,1]`$. With these terms the absolute square of $`t_{11}`$ \[Eq. (38)\] can be considerably simplified,kimfoot

$$T^{\mathrm{tot}}=|t_{11}|^2=\frac{(1r^2)^2}{|\alpha |^2|\beta |^2}\frac{\mathrm{sin}^2(\eta /2+\vartheta _0)}{\mathrm{sin}^2(\eta /2+\vartheta _0+\mathrm{\Delta })+\mathrm{\Gamma }_0^2},$$
(40)
with $`\varphi =(\phi _2+\varphi _2)(\phi _1+\varphi _1),\eta =(\phi _2+\varphi _2)+(\phi _1+\varphi _1),r^{}=(1ϵ)e^{i\varphi /2}+ϵe^{i\varphi /2},\delta =\mathrm{arg}(r^{}),\vartheta _0=\vartheta +\delta ,\alpha =1+re^{i(\eta +2\vartheta )},\beta =1+re^{2i\delta },\mathrm{\Delta }=\mathrm{arg}(\beta /\alpha )`$. The linewidth $`\mathrm{\Gamma }_0`$ is given by

$$\mathrm{\Gamma }_0=\left|\frac{1|r^{}\beta /\alpha |^2}{2r^{}\beta /\alpha }\right|.$$
(41)
In the generic case of $`r0`$, resonances occuring in Eq. (40) show a typical Fano profile of the formfano61

$$T^{\mathrm{tot}}|t^{\mathrm{bg}}|^2\frac{(B/cB_n/c)^2}{(B/cB_n/c+\mathrm{\Delta })^2+\mathrm{\Gamma }_0^2},$$
(42)
with $`t^{\mathrm{bg}}`$ being the coefficient for background scattering.nock94  The Fano resonances at $`B/c=B_n/c\mathrm{\Delta }`$ will have an asymmetric lineshape unless $`\mathrm{\Delta }=0`$ (i.e. $`r=0`$). This is however the case for the billiard systems we consider, since almost no reflection of incoming lead states takes place at the lead mouths, $`\stackrel{~}{r}_{11}^{\prime \prime }0`$, and therefore $`r0`$. Under this assumption Eq. (40) simplifies to

$$T^{\mathrm{tot}}\frac{\mathrm{sin}^2(\eta /2+\vartheta _0)}{\mathrm{sin}^2(\eta /2+\vartheta _0)+\mathrm{\Gamma }_0^2},$$
(43)
with linewidth $`\mathrm{\Gamma }_0=(1|r^{}|^2)/(2|r^{}|)`$. This equation describes symmetric resonance lineshapes which can be identified as window resonances (also called Breit-Wigner dips/antiresonances) of the form

$$T^{\mathrm{tot}}\frac{(B/cB_n/c)^2}{(B/cB_n/c)^2+\mathrm{\Gamma }_0^2}.$$
(44)
The physical picture resulting from this analysis is the following: In the magnetic field region $`B_c^2BB_c^1`$, where two edge states are present in the interior of the dot and one in each of the leads, the transmission probability shows large-scale oscillations intersected by sharp window resonances. The large-scale envelope function is given by $`1/(1+\mathrm{\Gamma }_0^2)`$. Its maxima perfectly match with the roughly estimated term $`\mathrm{cos}^2[B(A_1A_2)/(2c)]`$ from above and can therefore be identified with the numbered points in Fig. 7, each of which corresponds to an integer number of longitudinal antinodes in the wavefunction along the boundary (see Fig. 6). The antiresonances superimposed on these oscillations occur at magnetic fields $`B=B_n`$ (where $`\eta /2+\vartheta _0=n\pi ,n`$) and their linewidth is given by $`\mathrm{\Gamma }_0`$. As a result, resonances which are situated on maxima of the term $`1/(1+\mathrm{\Gamma }_0^2)`$ are sharper than at its mimima \[see numerical data in Fig. 7 for confirmation\]. For an increasing number of edge states populated in the cavity ($`B<B_c^2`$) our numerical results show that the density of antiresonances is rapidly growing. This behaviour finally leads to a resonance overlap for a large number of edge states, which is prerequisite for the onset of Ericson fluctuations (i.e. universal conductance fluctuations).
For completeness we remark that the above analysis for the $`B`$-dependence of $`T^{\mathrm{tot}}`$ can likewise be carried out with $`k_F`$ instead of $`B`$ as the variable parameter. We can similarly identify threshold values $`k_c^n`$ for $`k_F`$, below which a number of $`n`$ edge states survive. The numerical results for the transmission probablility $`T^{\mathrm{tot}}(k_F)`$ for the case of one or two participating edge states (not shown) can again be described by Eq. (37) or Eq. (40) respectively.

### III.5 Comparison with experiments

A series of experimentsmarc92 ; wees89 ; webb85 ; pede00  have been performed where Aharonov-Bohm oscillations (ABOs) similar to the ones discussed here have been observed in ballistic transport measurements. The origin of the ABOs in these experiments is however twofold: In Refs. marc92, ; wees89, , it is the presence of edge states in a quantum dot which gives rise to the observed oscillations. In Refs. webb85, ; pede00,  on the contrary, the investigated scattering devices have the form of a ring, to which the scattering wave function is confined. The latter setup thus gives rise to ABOs already at low magnetic fields and has therefore been more readily accessible to a theoretical description.pich97  However, to our knowledge, no quantitative description for magnetotransport through a quantum dot in the regime of only one or two participating edge states has yet become available. We therefore discuss in the following similarities and differences between the experimental data and our calculations in this field. One important observation is that the magnetic fields where these quasi-regular transmission fluctuations appear in the experiment are lower than in the present calculation. For example, in the experiment for circle and stadium shaped quantum dots in a GaAs/AlGaAs-heterostructure,marc92  the threshold magnetic field values would be (in SI-units)

$$B_c^n=\frac{2\pi h}{(2n+1)\lambda _F^2e}\mathrm{with}\lambda _F=\sqrt{\frac{2\pi }{n_s}}.$$
(45)
With a given sheet density of $`n_s=3.6\times 10^{11}\mathrm{cm}^2`$ in the interior of the dot, the threshold magnetic fields are given by $`B_c^23`$ Tesla and $`B_c^15`$ Tesla. However, in the experiment regular oscillations were already observed below 2 Tesla. At those field values we find highly irregular transmission fluctuations corresponding to a threshold magnetic field $`B_c^n`$ with $`n1`$, indicative of a high density of resonances and Ericson fluctuations. We expect the origin of this discrepancy to lie in the absence of sharp edges in the experiment and, hence, of diffractive edge scattering. In the experimental quantum dot, the edges should be fairly smooth, leading to near-adiabatic transitions to edge states at the entrance to the quantum dot. Therefore fewer edge channels are excited than by diffractive edge scattering, where all energetically accessible channels up to $`n`$ are populated. Our present results suggest that the observed transmission fluctuations are a direct measure of the sharpness of the edges at the lead mouth. Therefore, investigations of quantum dots with varying sharpness of edges would be desirable. Since these are, however, difficult to fabricate we point to a different experimental approach, which is based on the analogy between transport in the edge state regime and field-free transport through a rectangle where only few propagating modes participate. Such structures are accessible for microwave experiments.stoe99 ; blom02  The measured transmission through such a microwave device could provide a stringent test for the multi-channel interference model presented above.

## IV SUMMARY AND OUTLOOK

We have presented a new technique for calculating ballistic magnetotransport through open quantum dots. The Modular Recursive Green’s Function Method (MRGM) is an extension of the widely used standard recursive Green’s function technique and is based on the decomposition of non-separable scattering geometries into separable substructures (modules). An unprecedented energy and magnetic field range can thereby be explored with high accuracy. We applied the MRGM to transport coefficients and scattering wavefunctions in the two extreme cases of high magnetic fields and short wavelengths. For very small cyclotron radii we found periodic oscillations in the transmission spectrum and beating phenomena, which are restricted to well defined intervals (as a function of $`B`$ and $`k_F`$ likewise). These features could be explained by interferences between edge states, travelling along the boundary of the cavity. For these states scattering only takes place at the lead junctions, whose sharp edges play a crucial role for the dynamics of the system. For a detailed analysis a multi-channel interference model was employed. This model allows to classify the observed transmission fluctuations in the framework of Fano resonances. For only one edge state present in the circular dot transport is independent of the lead orientation provided that the lead mouths are identical and separated from each other. Future envisioned applications include the investigation of Andreev billiards,cser02  quantum Hamiltonian ratchets,scha01  fractal conductance fluctuations,ketz96 ; baec02  and shot noise.ober02  Furthermore, the MRGM also seems suitable to perform ab-initio calculations of the integer Quantum Hall effect.kosc01  The challenge is in this case the inclusion of a disorder potential which is compatible with the separability conditions.

###### Acknowledgements.

Helpful discussions with Prof. Langer, C. Stampfer and L. Wirtz are gratefully acknowledged. Many thanks are also due to W. Gansterer and X. Yang for their computer codes. This work was supported by the Austrian Science Foundation (FWF).