text:
# Chandra Limits on X-ray Emission Associated with the Supermassive Black Holes in Three Giant Elliptical Galaxies

## 1 Introduction

Knowledge about the frequency and masses of supermassive black holes (SMBH) in galactic nuclei, and the characteristics of their associated radiation, deepen our understanding of galaxy formation and clarify its connection with the active galactic nuclei (AGN) phenomenon (e.g., Blandford 1999). While the intense energetics and rapid variability in the most luminous active galaxies form the empirical foundation of the SMBH/AGN paradigm, such extreme objects are rare – especially in the local universe.

Nevertheless, there is now a consensus that most, if not all, galaxies host SMBH. This is based, in part, on the high frequency of low-level nuclear activity in nearby galaxies (Ho, Filipenko, & Sargent, 1997), and on statistical arguments connected to the cosmic X-ray background (e.g., Fabian & Iwasawa 1999). However, the strongest evidence for the ubiquity of SMBH comes from dynamical studies. Magorrian et al. (1998) found that, of 32 galaxies that could be fit with an axisymmetric dynamical model, 26 (29) require the presence of a central massive ($`10^810^{10}`$ M) dark object at 95% (68%) confidence. More unambiguous and accurate measurements based on Hubble Space Telescope (HST) stellar or ionized gas spectroscopy confirms the presence of SMBH while deriving systematically lower SMBH masses (Gebhardt et al., 2000a; Merritt & Ferrarese, 2001a).

Although kinematically-studied luminous AGN are found to host SMBH (Gebhardt et al., 2000b; McLure & Dunlop, 2000), the converse is not true: activity in most SMBH host galaxies is modest or absent. This can be quantified by considering the ratio of the nuclear luminosity to the Eddington luminosity ($`L_{\mathrm{edd}}=\mathrm{1.25\; 10}^{38}M_{\mathrm{SMBH}}`$, where the black hole mass $`M_{\mathrm{SMBH}}`$ is measured in solar units). The most powerful AGN have ratios on the order of unity, as expected from the theory of thin disks accreting at their maximal rate. However, this ratio is less than $`10^7`$ and $`10^4`$, respectively, in the two most securely measured SMBH – SgrA at the Galactic Center (Quataert, Narayan, & Reid, 1999) and NGC 4258 (Gammie, Narayan, & Blandford, 1999). If these SMBH accrete at the Bondi (1952) rate, they must do so in a fundamentally different – and much less radiatively efficient – way than in luminous AGN. The theory of advection-dominated accretion flows (ADAF; e.g., Narayan & Yi 1994), provides a predictive formalism for explaining the observed quiescence in the nuclei of SgrA, NGC 4258, and other systems (Fabian & Rees, 1995; Mahadevan, 1997).

Elliptical galaxies are the sites of the most luminous quasars, and also host the most extreme cases of quiescent SMBH. Although they contain nuclear black holes with the highest estimated masses and are observed to contain an ample supply of fuel in the form of X-ray emitting hot interstellar gas, they often emit less than $`10^6L_{\mathrm{edd}}`$ (Loewenstein et al., 1998; Di Matteo et al., 2000). Radio/X-ray observations of some systems are inconsistent with predictions of the simplest models of an ADAF accreting at the Bondi rate, but can be accommodated by variations that include outflow (Blandford & Begelman, 1999; Quataert & Narayan, 1999; Di Matteo et al., 1999) or convection (Narayan, Igumenshchev, & Abramowicz, 2000; Quataert & Gruzinov, 2000).

Prior to the launch of the Chandra X-ray Observatory, X-ray searches for nuclear point sources in elliptical galaxies were complicated in the imaging domain by centrally concentrated hot gas that dominates the soft X-ray flux; and, in the spectral domain by an extended hard component consisting of X-ray binaries. Thus the previous sensitivity to SMBH-associated X-ray emission, for even the nearest giant ellipticals, was $`10^{40}`$ erg s-1 (Roberts & Warwick, 2000; Sulkanen & Bregman, 2001), as neither the ASCA SIS/GIS nor the ROSAT PSPC/HRI had the necessary combination of sensitivity and spatial and spectral resolution to probe more deeply. Allen, Di Matteo, & Fabian (2000) claim to detect a very hard nuclear X-ray component in six elliptical galaxies with characteristics consistent with ADAF/outflow models, given the observed radio emission (Di Matteo et al., 1999). However, these results are based on a non-unique decomposition of ASCA spectra.

Because of the presence of bright stellar cores and dust, the best approach for constraining the level of SMBH activity is via high angular resolution X-ray observations (e.g., George et al. 2001): for a typical AGN X-ray/optical ratio, the Chandra detection threshold corresponds to a B magnitude $`24`$, while elliptical galaxy central surface brightnesses range from $`1116`$ mag arcsec-2. This is further demonstrated by deep Chandra surveys in combination with optical follow-up. Low-luminosity AGN, most not optically detected, dominate the X-ray background; and, nonthermal emission is found to be common – at any given time $`10`$% of all bright galaxies host $`L_X10^{41}`$ erg s-1 AGN (Mushotzky et al., 2000; Barger et al., 2001). There is also a high frequency of nuclear X-ray point sources associated with local AGN, though often at levels that overlap in luminosity with X-ray binaries (Ho et al., 2001). Chandra has detected unresolved nuclear emission associated with quiescent SMBH in M31 at $`5\times 10^9L_{\mathrm{edd}}`$ (Garcia et al., 2000), the Galactic Center at $`10^{11}L_{\mathrm{edd}}`$ (Baganoff et al., 2001), and NGC 6166 at $`10^7L_{\mathrm{edd}}`$ (Di Matteo et al., 2001).

In this paper we use Chandra imaging data to place luminosity limits on any nuclear point sources in the elliptical galaxies NGC 1399, NGC 4472, and NGC 4636. These limits raise new questions about the radiative efficiency and rate of accretion onto supermassive black holes. Following Merritt & Ferrarese (2001b), we adopt distances of 20.5, 16.7, and 15.0 Mpc for NGC 1399, NGC 4472, and NGC 4636, respectively.

## 2 Data Analysis and Results

### 2.1 Observations, Data Reduction and Analysis

NGC 1399, NGC 4472, and NGC 4636 were observed with the Chandra ACIS-S detector for 40–60 ksec during Cycle 1. Results on the diffuse hot gas emission and the X-ray binary population, and a more detailed data reduction description, are included elsewhere (e.g, Angelini, Loewenstein, & Mushotzky 2001). The latest revised standard pipeline processing was used for NGC 4472 and NGC 4636, while an earlier version in combination with additional bad pixel removal, gain map correction, and quantum efficiency map correction was used for NGC 1399. Final useful exposure times are 56, 40, and 52 ksec for NGC 1399, NGC 4472, and NGC 4636, respectively. NGC 1399 was approximately centered on the S3 chip; NGC 4636 and NGC 4472 were positioned at the standard ACIS-S aimpoint. 1.25–7 keV images of the central $`50^{\prime \prime }\times 50^{\prime \prime }`$ are shown for all three galaxies in Figures 1a, 1b, and 1c. The positions of the optical nuclei from archival HST NICMOS (NGC 4472, NGC 4636) or WFPC2 (NGC 1399) images are marked. The brightest point sources in the images are easily visible and indicated in the images; their counting rates are in the 0.005–0.008 cts s-1 range. The bright X-ray core in NGC 4472 is resolved and has a thermal spectrum. Only in the case of NGC 4636 is there a possible coincidence of the nucleus with a Chandra X-ray source (just to the southwest of the nuclear position indicated in Figure 1c), with a 0.3–7 keV count rate $`0.004`$ cts s-1 (see below). 0.3–7.0 keV surface brightness profiles are shown in Figure 2. For illustration purposes, the NGC 4636 profile is centered on the near-nuclear source. No excess emission from unresolved nuclear emission is apparent for the other two galaxies.

The position of the nucleus in NGC 1399 is identified based on registration of globular cluster X-ray sources with their WFPC2 counterparts (Angelini et al. 2001); the resulting HST and Chandra centroids are consistent. Similarly, the nucleus of NGC 4636 is located based on registration with a bright NICMOS source and with the USNO catalog (Monet et al., 1998), and is found not to be consistent with the position of the near-nuclear source (Figure 1c). For NGC 4472, we adopt the unregistered HST centroid that is consistent with the peak in the soft thermal emission seen in the broadband Chandra image (Figure 1b). $`3\sigma `$ upper limits to any nuclear X-ray point source emission are derived following Gehrels (1986). Counts are accumulated in inner (“source”) and outer (“background”) concentric annuli with radii that are varied to obtain the tightest upper limit. The Chandra calibration database enclosed energy function at the appropriate off-axis angle is used, and an exposure map correction is applied for NGC 1399, which was centered off-axis. The assumed position of the SMBH is allowed to vary in a $`3^{\prime \prime }\times 3^{\prime \prime }`$ box (corresponding to a very conservative adopted Chandra/HST relative systematic pointing uncertainty) centered on the assumed nuclear positions and the highest upper limit in this box is adopted. These boxes encompass the peaks of the diffuse emission for NGC 4472 and NGC 1399. The diffuse hot gas emission is implicitly assumed to be locally flat on average when treated as “background”, consistent with the flat central surface brightness profiles (Figure 2). Luminosity upper limits are derived both in the 0.3–7 keV and 2–10 keV energy bands assuming a slope 1.5 power-law spectrum – the limits increase by $`2`$ for slope 0.5. Due to reduced dilution from the diffuse hot gas emission, the 2–10 keV count rate limits are 3–5 times lower than in the softer band; however, because of the declining effective area the energy flux upper limits are comparable. The 0.3–7 keV detection limit varies from a few $`10^{37}`$ erg s-1 in regions of low surface brightness hot ISM emission to a few $`10^{38}`$ erg s-1 in the core of NGC 1399.

The 2–10 keV limits – all derived within $`1^{\prime \prime }`$ of the optical nucleus – are displayed in Table 1; these correspond to 28, 23, and 16 counts for NGC 1399, NGC 4472, and NGC 4636, respectively – comparable to the brightest sources in the subimages shown in Figure 1. Table 1 also includes less conservative, but perhaps more realistic, limits calculated by fixing the positions at the assumed nuclear positions described above.

The inner gas density profile is estimated by fitting the point-source-subtracted surface brightness profiles with $`\beta `$-models,

$$\mathrm{\Sigma }=\mathrm{\Sigma }_o\left(1+\frac{r^2}{a^2}\right)^{(3\beta 1/2)}$$
(1)
out to some maximum radius, followed by analytical deprojection with conversions obtained from spectral fits to the inner regions (best-fit temperatures $`kT=0.8`$, 0.8, and 0.6 keV for NGC 1399, NGC 4472, and NGC 4636, respectively). The quality of fit declines as the maximum radius is increased due to deviations from the $`\beta `$-model, but the profiles are sufficiently concentrated that the central density estimates – 4.6, 3.2, and $`1.1\times 10^{25}`$ gm cm-3 for NGC 1399, NGC 4472, and NGC 4636, respectively – are robust.

### 2.2 Accretion Rates and Eddington Ratio Limits

Individual black hole masses for the three galaxies studied here have been estimated by Magorrian et al. (1998) from two-integral dynamical modeling of ground-based kinematical data. A more accurate value can be derived using the correlation with bulge velocity dispersion, $`\sigma `$, derived for galaxies with limits based on HST spectroscopy. Therefore, we adopt SMBH masses from Merritt & Ferrarese (2001b) (Table 1); these are 3–5 times lower than in Magorrian et al. (1998). The resulting 2–10 keV luminosities in units of the Eddington luminosity and the Bondi (1952) accretion rates for adiabatic spherical inflow are displayed in Table 1. Our derived upper limits on $`L_X/L_{\mathrm{edd}}`$ are $`7.3\times 10^9`$, $`1.5\times 10^8`$, and $`2.8\times 10^8`$ for NGC 1399, NGC 4472, and NGC 4636, respectively. The corresponding efficiencies, assuming accretion at the Bondi rate $`\dot{M}_{\mathrm{Bondi}}`$, are $`4.1\times 10^6`$, $`2.4\times 10^5`$, and $`6.2\times 10^4`$.

## 3 The Nature of the Accretion Flow

The conservative upper limits we derive are 1–2 orders of magnitude lower than the luminosities estimated by Allen et al. (2000) using ASCA spectral decomposition. Direct Chandra images (Angelini et al. 2001) indicate that most of the hard X-ray emission spatially unresolved by ASCA originates in X-ray binaries. Since upper limits on the nuclear X-ray emission and estimates of the Bondi accretion rate have each decreased by $`1030`$ from Di Matteo et al. (2000) (the latter because of lower black hole mass estimates)111For NGC 4636 our Bondi accretion rate is $`3000`$ times smaller than that of Di Matteo et al. (2000). A factor of $`10`$ is due to our smaller black hole mass; the remainder is primarily due to a lower inferred central gas density that may be underestimated since the accretion radius is unresolved in this case., a re-examination of the nature of the accretion flow is in order.

For our estimated black hole masses and Bondi accretion rates, we have computed spectra using the techniques described in Quataert & Narayan (1999); we aim to test the basic ADAF model and so do not include outflows or convection. The resulting predicted 2–10 keV luminosities, $`L_{\mathrm{ADAF}}`$, are $`2\times 10^{41}`$, $`10^{40}`$, and $`10^{36}`$ ergs s-1 for NGC 1399, NGC 4472, and NGC 4636, respectively (Table 1). Our observational upper limits are thus significantly below the ADAF predictions for NGC 1399 and NGC 4472, but consistent for NGC 4636. The predicted X-ray luminosities in the standard ADAF models for the cases of NGC 1399 and NGC 4472 are dominated by Compton scattered synchrotron emission that originates close to the black hole. The luminosity is primarily sensitive to the gas density (i.e., accretion rate) and electron temperature. A factor of $`10`$ decrease in the accretion rate below the Bondi value can reduce the predicted luminosity to the observational limits.

With the revised black hole masses and accretion rates we find that the radio flux predicted by the ADAF model is less than the observed flux for NGC 1399 and NGC 4636, but is in excess for NGC 4472; thus the radio flux problem highlighted by Di Matteo et al. (1999) is somewhat ameliorated. Finally, a preliminary HST measurement of the nuclear FUV emission in NGC 1399 (0’Connell et al., in preparation), – $`\nu L_\nu \mathrm{7\; 10}^{38}`$ erg s-1 – also conflicts with pure-inflow ADAF models, although dust may be an issue at this wavelength.

Our observations require a modification of the standard ADAF/Bondi model such that the gas density near the SMBH is reduced by at least a factor of $`10`$. This could be due to the role of strong outflows or convection in an ADAF. Alternatively, the central regions of these galaxies may not be in steady state due to heating by supernovae or by episodic accretion onto the central SMBH (e.g., Ciotti & Ostriker 1997). In addition the cooling time at the accretion radius $`R_{\mathrm{Bondi}}`$ (Table 1) is $`<10^7`$ yr for all three galaxies; so, inflow with mass deposition could be important (Brighenti & Mathews, 1999).

## 4 Conclusions

We have used Chandra observations to place new and unprecedented upper limits on the X-ray emission associated with $`10^810^9`$ M nuclear black holes in three giant elliptical galaxies. The black holes clearly are not fuel-starved: all have central hot gas densities exceeding $`10^{25}`$ gm cm-3. Adopting masses from the $`M_{\mathrm{SMBH}}\sigma `$ correlation, and assuming that mass is inflowing at the Bondi (1952) rate, we derive $`3\sigma `$ 2–10 keV luminosity upper limits of 7.3, 15 and $`28\times 10^9L_{\mathrm{Edd}}`$ and radiative efficiencies in this band of 4.1, 24, and $`620\times 10^6`$ for NGC 1399, NGC 4472, and NGC 4636, respectively. Unless the black hole masses are considerably lower than predicted by the $`M_{\mathrm{SMBH}}\sigma `$ correlation, these observations are inconsistent with the simplest ADAF models in NGC 1399 and NGC 4472, and imply accretion onto the SMBH at below the Bondi rate.

If the decline in the accretion rate is due to outflow, this raises the question of the ultimate fate of this gas as it encounters the dense hot interstellar medium in the core. Many questions about stability, intermittency, and regulation of gas flows in the inner regions of ellipticals remain to be addressed. The Chandra isophotes in the cores of NGC 1399 and NGC 4472 appear round and undisturbed, although there is an apparent spiral-like structure in the NGC 4636 diffuse emission.

The more stringent limits required to further constrain the nature of the central accretion flow will be difficult. Our observations indicate $`>0.01`$ super-Eddington X-ray binaries per square-arcsecond near the nucleus; so that detections below $`10^{39}`$ erg s-1 cannot be unambiguously identified as AGN. Although the hot gas emissivity precipitously declines at energies above 2 keV, so does the Chandra sensitivity, restricting the detection limit in gas-rich systems. Tighter limits can be obtained on gas poor systems that, although of interest, may be fuel-starved and would thus shed little light on the accretion process. It would be especially useful to spectrally decompose the inner $`1^{\prime \prime }`$ to place a limit on the fraction of non-ISM emission. This would appear as a hard excess and would require a very deep exposure due to the declining Chandra effective area mentioned above.

We are grateful to Dave Davis and Una Hwang for assistance with data reduction, and to Wayne Landsman for sharing unpublished HST results. We made use of the HST archive, the NASA/IPAC Extragalactic Database (NED), and IDL function fitting routines developed by Craig Markwardt.


format: text

text:
# Tunneling decay in a magnetic field

## I Introduction

Magnetic field can very strongly change the tunneling rate of charged particles. This change, in turn, strongly depends on properties of the system, as in the well-known effect of giant hopping magnetoresistance in solids . Therefore tunneling in a magnetic field has been broadly used as a sensitive and revealing probe of electron systems in solids, including quantum Hall systems , two-layer heterostructures away from the quantum Hall region , and correlated electron layers on the surface of liquid helium .

Correlated two-dimensional (2D) electron systems are currently attracting much interest . The possibility to extract information about electron correlations and dynamics through tunneling in a magnetic field is one of the motivations of the present work. The major motivation, however, comes from the fact that tunneling in a magnetic field is an interesting and in many respects unusual theoretical problem, even in the single-particle formulation. Existing results, although often highly nontrivial, are limited to the cases where the potential has either a special form  (e.g., linear  or parabolic ), or a part of the potential or the magnetic field are in some sense weak .

The problem of tunneling has two parts. One is to find the tail of the wave function of the localized intrawell state under the potential barrier $`U(𝐫)`$ and behind it, and the other is to find the escape probability $`W`$. In the magnetic field, $`W`$ differs exponentially from the probability to reach the boundary of the classically allowed range $`U(𝐫)=E`$, where $`E`$ is the energy of the particle. This is because, as it tunnels, the particle is accelerated by the Lorentz force, and it comes out from the barrier with a finite velocity $`𝐯`$. The standard argument that the exit point is the turning point $`𝐯=\mathrm{𝟎}`$ relies on time-reversal symmetry (see below) and does not apply in the presence of a magnetic field.

A simple potential $`U(𝐫)`$ and the wave function $`\psi (𝐫)`$ of the metastable state in this potential are sketched in Figs. 1, 2. The wave function decays away from the potential well. At some point $`𝐫`$, on the background of the decaying tail there emerges a propagating small-amplitude wave packet, which corresponds to the escaped particle. As a result, in a part of the classically allowed region $`U(𝐫)<E`$ the function $`\psi (𝐫)`$ is determined by this wave packet, whereas in the other part of this region $`\psi `$ is determined by a different branch of the tail of the localized state. The boundary between these areas is sharp, and the slope of the wave function changes on this boundary nearly discontinuously.

An important part of the WKB formulation of the tunneling escape problem in a magnetic field was found (b) in the analysis of decay for a special model of an atomic system \[see Eq. (10) below\]. In a general case, both the tail of a metastable state and the outgoing wave packet can be found using the approach briefly outlined in our Letter .

In the WKB approximation the wave function is sought in the form

$$\psi (𝐫)=D(𝐫)\mathrm{exp}[iS(𝐫)](\mathrm{}=1).$$
(1)
Here, $`S(𝐫)`$ is the classical action, and $`D`$ is the prefactor. In the classically allowed range, (1) describes a wave propagating with a real momentum $`𝐩=\mathbf{}S`$. On the other hand, in the classically forbidden range the wave function decays. For the ground intrawell state, the decay of $`\psi (𝐫)`$ is not accompanied by oscillations in the absence of a magnetic field. Then the action $`S(𝐫)`$ is purely imaginary under the barrier, and $`|\mathbf{}S|`$ is the decrement of the wave function.

Both behind and inside the barrier, the action can be obtained from the Hamilton equations of motion

$`\dot{S}=𝐩\dot{𝐫},\dot{𝐫}=H/𝐩,\dot{𝐩}=H/𝐫,`$ (2)
where $`H`$ is the Hamiltonian of the system,

$$H=\frac{1}{2m}\left[𝐩+\frac{e}{c}𝐀(𝐫)\right]^2+U(𝐫),$$
(3)
($`𝐀(𝐫)`$ is the vector potential of the magnetic field).

In the standard approach to tunneling decay, which applies for $`B=0`$ , one looks for the purely imaginary action $`S`$ under the barrier. It is calculated by changing to imaginary time and momentum in Eqs. (2), but keeping coordinates real,

$$tit,𝐩i𝐩,𝐫𝐫,U(𝐫)U(𝐫).$$
(4)
Eqs.(2) then take the form of equations of classical motion in an inverted potential $`U(𝐫)`$, with energy $`EU(𝐫)`$. The imaginary-time trajectory goes from the turning point inside the potential well to the turning point on the boundary of the classically allowed region, where it matches the appropriate classical trajectory of the escaped particle behind the barrier, with real $`𝐫,𝐩,t`$. The velocities on the both trajectories at the matching point can coincide only if $`𝐯=\mathrm{𝟎}`$.

In the presence of a magnetic field, because of broken time-reversal symmetry, the replacement (4) may not be performed. It would lead to a complex Hamiltonian, which makes no sense and indicates that a more general approach is required. The action $`S(𝐫)`$ is complex under the barrier for real $`𝐫`$. This complexity has also an important counterpart in the instanton formulation of the problem of tunneling decay in a magnetic field, see below. We note that a complex action arises also in other cases, like barrier penetration for oblique incidence  and scattering by a complex potential (as in the case of an absorbing medium) . The method discussed below can be applied to these problems as well.

In this paper we consider a single-particle tunneling decay in a three-dimensional potential of a general form, for arbitrary magnetic fields. We illustrate the approach using a toy model of a correlated 2D electron system. We show that the tunneling exponent $`S(𝐫)`$ and the escape rate in a magnetic field can be found from dynamical equations (2) by analytically continuing these equations to complex phase space and time. The initial conditions for the trajectories are determined by the analytical continuation of the usually known intrawell wave function. The resulting set of complex trajectories has singularities, caustics, on which there occurs branching of the complex action $`S(𝐫)`$ and the tails of the decaying and propagating waves are matched. Careful analysis allows us to find the complete semiclassical wave function and reveal the singular features of $`\psi (𝐫)`$ related to the branching of $`S`$.

In Sec. II and Appendix A we provide a simple model which catches basic physics of tunneling from correlated 2D electron layers. In Sec. III we consider the tunneling exponent and formulate the boundary value problem for tunneling trajectories in a magnetic field in complex phase space. In Sec. IV we discuss matching of different semiclassical solutions across the caustic of the set of the tunneling trajectories. We show that a switching surface (one of the anti-Stokes manifolds) starts at the caustic. The wave function has an observable singular feature at this surface, which is a sharp change of the slope of $`\mathrm{ln}|\psi (𝐫)|`$. In Sec. V we provide explicit results for two simple exactly solvable models of physical interest, which also illustrate general features of tunneling in a magnetic field. In Sec. VI we discuss the path-integral formulation of the problem of tunneling decay in a magnetic field. Sec. VII contains concluding remarks.

## II A model of the tunneling barrier for a correlated 2D electron system

One of the most interesting systems where tunneling in a magnetic field has been investigated experimentally is a correlated 2D electron system. Here, electrons are localized in the $`z`$-direction in a metastable 1D potential well $`U_0(z)`$. The intrawell electron motion is quantized in the $`z`$-direction, and electrons can tunnel from the well into extended states. Many 2D systems of current interest are strongly correlated: electrons are far away from each other, exchange is weak, and there is at least short-range order in the $`(x,y)`$-plane. The tunneling electron can be then identified and “labeled”. Its tunneling motion is accompanied by the motion of other electrons. The many-electron dynamics of a correlated system can be described in terms of in-plane electron vibrations, and the corresponding Hamiltonian is given in Appendix A assuming that the electrons form a Wigner crystal. Here we will made a further simplification and think of an electron as tunneling in a static potential created by all other electrons. As we showed earlier, this is a good approximation for the tunneling problem.

The static potential from surrounding electrons can be assumed to be parabolic with respect to the in-plane coordinates $`x,y`$. If the characteristic width $`L`$ of the tunneling barrier is less than the interelectron distance, the overall potential is a sum of the parabolic in-plane part and $`U_0(z)`$,

$$U(𝐫)=\frac{m\omega _0^2}{2}(x^2+y^2)+U_0(z).$$
(5)
The form of $`U_0(z)`$ depends on the system. Inside the well $`U_0`$ is often singular, like in the case of electrons on helium where $`U_0`$ includes the image potential. The potential barrier itself can be close to a square barrier, as in the case of unbiased semiconductor heterostructures, or can be nearly linear, as in the presence of strong enough bias voltage, with

$$U_0(z)=\frac{\gamma ^2}{2m}\left(1\frac{z}{L}\right)$$
(6)
inside the barrier. The potential (6) describes, in particular, the barrier for a correlated 2DES on a helium surface. This system was experimentally investigated in Ref. , and showed an unexpected dependence of the tunneling rate on $`B`$ and electron density that we recently addressed .

## III The tunneling exponent

For a smooth tunneling barrier $`U(𝐫)`$, the underbarrier wave function $`\psi (𝐫)`$ (1) can be obtained from the tunneling trajectories (2). The initial conditions for these trajectories are determined by the tail of the intrawell wave function. They can be obtained even if the potential $`U(𝐫)`$ is singular within the well, as in the case of an image potential or a stepwise potential in a semiconductor heterostructure.

To obtain the initial conditions for the trajectories (2), we can take a surface $`\mathrm{\Sigma }`$ close to the well and yet in the range where $`U(𝐫)`$ is already smooth. The wave function $`\psi (𝐫)`$ on $`\mathrm{\Sigma }`$ is presumably known from the solution of the Schrödinger equation inside the well and is semiclassical. Only the exponent of this wave function is needed to find the initial conditions for (2), which take the form

$$𝐫(0)=𝐫|_\mathrm{\Sigma },𝐩(0)=i\left[\mathbf{}\mathrm{ln}\psi (𝐫)\right]_\mathrm{\Sigma },$$
(7)
with the action $`S(0)=i\left[\mathrm{ln}\psi (𝐫)\right]_\mathrm{\Sigma }`$. Only the lowest-order terms in $`\mathrm{}`$ should be kept in the expressions for $`𝐩(0),S(0)`$. The final result should be independent of the choice of $`\mathrm{\Sigma }`$.

The trajectories (2) with the initial conditions (7) form a two-parameter set, in the case of 3D tunneling. The two parameters are the initial coordinates on the surface $`\mathrm{\Sigma }`$. We can choose curvilinear coordinates $`(x_1,x_2,x_3)`$ so that $`x_3|_\mathrm{\Sigma }=0`$, and respectively $`x_3(0)=0`$. The trajectories are then parametrized by $`x_1(0),x_2(0)`$.

To illustrate these arguments we consider the initial conditions for an electron with the potential (5), (6), which tunnels from a 2D layer. Inside the metastable potential well the electron motion separates into a quantized motion in the normal to the layer $`z`$-direction and in-plane vibrations. To slightly simplify the analysis, we will neglect the effect of a magnetic field on the intrawell wave function, but not on the tail of $`\psi (𝐫)`$ deep under the barrier where the effect will have accumulated. This is a good approximation for not too strong fields provided the characteristic intrawell localization length $`1/\gamma `$ is small compared to the tunneling length $`L`$.

It is convenient then to choose the surface $`\mathrm{\Sigma }`$ as a plane $`z=`$ const close to the well, but behind the intrawell turning point. We set $`z=x_3=0`$ on $`\mathrm{\Sigma }`$ and choose $`x_1=x,x_2=y`$. If we set the energy of the out-of-plane motion $`E=0`$ and assume that the electron is in the ground intrawell state, we obtain from (5)

$`z(0)=0,p_z(0)=i\gamma ,S(0)=im\omega _0[x^2(0)+y^2(0)]/2,`$ (8)
$`p_x(0)=im\omega _0x(0),p_y(0)=im\omega _0y(0).`$ (9)
To find the trajectory which arrives at a given real $`𝐫`$ deep under the barrier, it may be necessary, particularly for $`B0`$, to start with a complex $`𝐫(0)`$. The corresponding values of $`𝐩(0)`$ can be found by analytically continuing the intrawell wave function to complex $`𝐫`$. The whole trajectory will then go in complex phase space (including the configurational space) and also in complex time $`t`$. The energy of the trajectory is given by the energy of the intrawell state from which the system tunnels. It remains real.

The rate of tunneling decay is determined by $`\mathrm{Im}`$ $`S`$ at the point where the particle emerges from the barrier as a semiclassical wave packet. This wave packet propagates along a real classical trajectory $`𝐫_{\mathrm{cl}}(t)`$, which is a real-time solution of the Hamiltonian equations (2). The underbarrier trajectory for tunneling escape should coalesce with this classical trajectory. Therefore at some time $`t`$ it should have both real coordinate and velocity,

$$\mathrm{Im}𝐫(t)=\mathrm{Im}𝐩(t)=0.$$
(10)
Eqs. (10) determine the complex starting point of the trajectory for tunneling escape $`𝐫(0)`$ \[i.e., the complex $`x_{1,2}(0)`$, since $`x_3(0)=0`$\] and also the imaginary part of the duration of motion along this trajectory Im $`t`$. The real part of $`t`$ remains undetermined: a change in Re $`t`$ in (10) results just in a shift of the particle along the classical trajectory $`𝐫_{\mathrm{cl}}(t)`$, see Fig. 3. Such a shift does not change Im $`S`$, since $`𝐩=\mathbf{}S`$ is real along $`𝐫_{\mathrm{cl}}(t)`$. We note that the number of equations (10) is equal to the number of variables Re $`x_{1,2}(0)`$, Im $`x_{1,2}(0)`$, and Im $`t`$ \[the value $`x_3(0)=0`$ is fixed on $`\mathrm{\Sigma }`$\], with account taken of energy conservation. The conditions (10) were first given (b) for a $`\delta `$-shape potential well and a linear tunneling barrier, but only the condition Im $`𝐫(t)=0`$ was used.

In the absence of a magnetic field, we can choose the surface $`\mathrm{\Sigma }`$ such that the momentum $`𝐩|_\mathrm{\Sigma }`$ is imaginary for real $`𝐫`$, i.e. the decay of the localized wave function is not accompanied by oscillations \[cf. (8)\]. The equations of motion (2) can then be solved in purely imaginary time, with real $`𝐫(t)`$ and imaginary $`𝐩(t)`$. The escape trajectory ends at the turning point $`𝐩=\mathrm{𝟎}`$, even for a multidimensional system .

The tunneling exponent $`R`$ is given by the value of Im $`S`$ at any point on the trajectory $`𝐫_{\mathrm{cl}}`$,

$$R=2\mathrm{Im}S(𝐫_{\mathrm{cl}}).$$
(11)
For a physically meaningful solution, Im $`S`$ should have a parabolic minimum at $`𝐫_{\mathrm{cl}}`$ as a function of the coordinates transverse to the trajectory. The outgoing wave packet will then be Gaussian near the maximum.

From (10), even in the presence of a magnetic field the tunneling exponent can be obtained by solving the equations of motion (2) in imaginary time, with complex $`𝐫`$. However, such solution does not give the wave function for real $`𝐫`$ between the well and the classical trajectory $`𝐫_{\mathrm{cl}}`$. Neither does it tell us where the particle shows up on the classical trajectory.

## IV Branching of the action and its observable consequences

The complete WKB solution of the tunneling problem can be obtained and the wave function $`\psi (𝐫)`$ can be found if one takes into account that the action $`S`$ as given by Eqs. (2) is a multivalued function of $`𝐫`$, even though it is a single-valued function of $`t`$ and $`x_{1,2}(0)`$. This means that several trajectories (2) with different $`t`$ and $`x_{1,2}(0)`$ can go through one and the same point $`𝐫`$. However, except for the points on the switching surface (see below), only one of the branches of the action $`S(𝐫)`$ contributes to the wave function $`\psi (𝐫)`$.

### A Caustics in a magnetic field

In multidimensional systems, branching of the semiclassical action generally occurs on caustics, or envelopes of the Hamiltonian trajectories , see Fig. 4. Caustics are multidimensional counterparts of turning points familiar from the analysis of tunneling in 1D systems . The prefactor $`D(𝐫)`$ in the WKB wave function (1) diverges at a caustic. In the case of 1D semiclassical motion along the $`z`$-axis we have $`Dp_z^{1/2}`$, and $`D\mathrm{}`$ at the turning points $`z_t`$, which are given by the condition $`p_z=0`$. The action is branching at turning points, $`S(z)S(z_t)(zz_t)^{3/2}`$. Its behavior near caustics in a multidimensional system is very similar (see below), with $`zz_t`$ corresponding to the distance from the caustic.

Since neighboring Hamiltonian trajectories (2), (7) touch each other on a caustic, the one-to-one correspondence between the coordinates $`x_1`$, $`x_2`$, $`x_3`$ on the trajectory and the parameters $`t`$, $`x_1(0)`$, $`x_2(0)`$ breaks. Therefore the equation for a caustic has the form

$$J(𝐫)=0,J(𝐫)=\frac{(x_1,x_2,x_3)}{(x_1(0),x_2(0),t)}.$$
(12)
The Jacobian $`J(𝐫)`$ can be related in a standard way to the prefactor $`D(𝐫)`$, which in turn is determined by the first-order (in $`\mathrm{}`$) correction to the action $`iS^{(1)}`$, $`D(𝐫)=\mathrm{exp}[S^{(1)}(𝐫)]`$. The equation for $`S^{(1)}(𝐫)`$ can be obtained from the Schrödinger equation by seeking the wave function in the form $`\psi =\mathrm{exp}(iS)`$, with $`S=S^{(0)}iS^{(1)}`$. This gives $`2𝐯S^{(1)}=\mathbf{}𝐯`$, where $`m𝐯=\mathbf{}S^{(0)}+(e/c)𝐀`$. The vector $`𝐯`$ gives the velocity on the Hamiltonian trajectory (2). Taking into account that $`𝐯S^{(1)}(𝐫)dS^{(1)}/dt`$ and that $`\mathbf{}𝐯=d\mathrm{ln}J(𝐫)/dt`$, where the time derivatives are taken along the trajectory, we obtain

$$D(𝐫)=\mathrm{const}\times [J(𝐫)]^{1/2}.$$
(13)
It follows from Eqs. (12), (13) that the prefactor $`D`$ diverges on caustics, and the WKB approximation does not apply there.

There are both formal and physical distinctions between caustics for tunneling trajectories with and without a magnetic field. For $`B0`$, the trajectories are complex, and Eq. (12) specifies a complex surface. This surface intersects the real 3D space along a line. In contrast, for $`B=0`$, because of time-reversal symmetry, tunneling trajectories lie in real configuration space. In this symmetric case Eq. (12) specifies a surface  rather than a line in the 3D space. As we show below, this distinction leads to observable consequences.

### B Local analysis near caustics

The analysis of the wave function and branching of the action $`S(𝐫)`$ at complex caustics in the magnetic field is similar to that for caustics in real space, including turning points in the 1D case . Near a caustic, it is convenient to change to the variables $`x^{},y^{}`$, and $`z^{}`$ which are locally parallel and perpendicular to the caustic surface, respectively (we set $`z^{}=0`$ on the caustic). Since a caustic is an envelope of the Hamiltonian trajectories (2), the normal to the caustic component of the velocity is $`v_z^{}=0`$ for $`z^{}=0`$. However, for $`B0`$ the normal component of the momentum is not equal to zero. Therefore the wave function near a point $`𝐫_{\mathrm{caust}}`$ on the caustic can be sought in the form

$$\psi (𝐫_{\mathrm{caust}}+𝐫^{})=e^{i𝐩_{\mathrm{caust}}𝐫^{}}\varphi (z^{};𝐫_{\mathrm{caust}}),$$
(14)
where $`𝐩_{\mathrm{caust}}`$ is the momentum at the point $`𝐫_{\mathrm{caust}}`$ (note that $`𝐫_{\mathrm{caust}}`$ is a 2D complex vector, $`z_{\mathrm{caust}}^{}0`$). We assume that the dependence of $`𝐩_{\mathrm{caust}}`$ on $`𝐫_{\mathrm{caust}}`$ (i.e., along the caustic) is smooth, and the dependence of $`\varphi `$ on $`𝐫_{\mathrm{caust}}`$ is much more smooth than on $`z^{}`$. Generally, $`𝐩_{\mathrm{caust}}`$ is complex even where the caustic goes through real space, Im $`𝐫_{\mathrm{caust}}=0`$. Therefore the classical trajectory of the escaped particle does not go through the caustic, in contrast to the case of zero magnetic field, cf. Fig. 4.

Because $`v_z^{}=0`$, the equation for $`\varphi (z^{};𝐫_{\mathrm{caust}})`$, which follows from the 3D Schrödinger equation with a magnetic field, coincides with the 1D Schrödinger equation near a turning point

$$\left[\frac{1}{2m}\frac{d^2}{dz_{}^{}{}_{}{}^{2}}+U_z^{}^{}(𝐫_{\mathrm{caust}})z^{}\right]\varphi (z^{};𝐫_{\mathrm{caust}})=0$$
(15)
\[here, $`U_z^{}^{}U/z^{}`$\]. The boundary conditions to this equation are discussed below.

The function $`\varphi `$ is singled-valued. It is given by a linear combination of the Airy functions . For comparatively large $`|z^{}|`$ (but still close to the caustic) it becomes a linear combination of the functions

$$w_{1,2}=(z^{})^{1/4}\mathrm{exp}\left[i\alpha z_{}^{}{}_{}{}^{3/2}\right],$$
(16)
with $`\alpha \alpha (𝐫_{\mathrm{caust}})=(2/3)[2mU_z^{}^{}]^{1/2}`$. To make the functions $`w_{1,2}`$ uniquely defined, we have to make a cut on the complex $`z^{}`$-plane. Our choice of the cut is shown in Fig. 5.

With account taken of Eq. (14), we find that the action near the caustic is

$`S(x^{},y^{},z^{})S(x^{},y^{},0)+(p_{\mathrm{caust}})_z^{}z^{}+\alpha z_{}^{}{}_{}{}^{3/2},`$ (17)
with an appropriately chosen branch of $`z_{}^{}{}_{}{}^{3/2}`$.

Another way to understand Eq. (17) is based on the analysis of the set of the Hamiltonian trajectories (2). Because the caustic is an envelope of the trajectories and $`v_z^{}=0`$ on the caustic, $`z^{}`$ is quadratic in the increments $`\delta x_{1,2}(0),\delta t`$ of the parameters of the set. Therefore $`\delta x_{1,2}(0),\delta t`$ are nonanalytic in $`z^{}`$, as is also the action $`S`$. Taking into account cubic terms in $`\delta x_{1,2}(0),\delta t`$ we obtain (17). The coefficients in $`S`$ can be expressed in terms of the derivatives of $`S,𝐫`$ calculated along the trajectories over $`x_{1,2}(0),t`$.

### C Choosing the action branches

Eq. (15) describes how the WKB solutions, which correspond to different branches of the action $`S`$, connect on the caustic. Of interest for the problem of tunneling escape is the caustic where there are connected the tails of the intrawell wave function and the outgoing wave packet for the escaped particle. From (12), on this caustic there is a point through which there goes the Hamiltonian trajectory (2), (7) for escape, with the initial coordinates $`x_{1,2}(0)`$ given by the condition (10) of arriving, ultimately, at the classical outgoing trajectory $`𝐫_{\mathrm{cl}}`$, cf. Fig. 4.

In general, a caustic can be thought of as a “mirror”, which partly reflects the wave packet. The boundary condition to Eq. (15) for the tunneling escape problem is the “radiation condition”. Very far from the potential well the solution of the full Schrödinger equation is a semiclassical wave packet moving in real space away from the well. Respectively, the wave function far from the well has a form $`\psi (𝐫)\mathrm{exp}[iS(𝐫)]`$. This solution has to be continued to the caustic, which means that there is a range of directions in the complex $`z^{}`$-plane not too close to the caustic ($`|\alpha z_{}^{}{}_{}{}^{3/2}|\mathrm{}`$) where the wave function is described by only one exponential $`\mathrm{exp}[iS(𝐫)]`$. Away from this range, the wave function is a combination of two waves. This is again similar to a 1D problem, where behind the turning point $`z_t`$ the wave function for real $`z`$ is a propagating wave, whereas before $`z_t`$ there is a wave with an exponentially decaying amplitude coming from the intrawell state and the wave reflected back to this state.

To connect the tails of the wave functions near the caustic in our problem one can use the explicit solution of Eq. (15) with the radiation boundary condition, as in the 1D case .

An alternative way to see how the boundary condition works, which also allows us to reveal the specific feature of tunneling in a magnetic field, is to follow the transformation of the wave function $`\varphi `$ for $`|\alpha z_{}^{}{}_{}{}^{3/2}|\mathrm{}`$ as the argument of $`z^{}`$ varies by $`2\pi `$. This analysis is based on the notion of the Stokes and anti-Stokes lines . We count $`\mathrm{arg}z^{}`$ off from $`\mathrm{arg}\alpha ^{2/3}`$. Then, for the choice of the cut in Fig. 5 and with the functions $`w_{1,2}`$ given by Eq. (16), the Stokes lines are the rays $`\mathrm{arg}z^{}=(2n+1)\pi /3`$ with $`n=0,1,2`$. On the Stokes lines $`\mathrm{Re}z_{}^{}{}_{}{}^{3/2}=0`$ and the ratio $`|w_2/w_1|`$ is extremal (maximal or minimal). The anti-Stokes lines are the rays $`\mathrm{arg}z^{}=2n\pi /3`$ with $`n=0,1,2`$, where $`|w_2/w_1|=1`$.

From the radiation boundary condition, there is a range $`\mathrm{\Delta }`$ of $`\mathrm{arg}z^{}`$ where $`\varphi (z^{})`$ is given by only one of the functions $`w_{1,2}(z^{})`$, not by a superposition of $`w_1`$ and $`w_2`$. This condition is physically meaningful provided the corresponding $`w_i`$ is exponentially small compared to $`w_{3i}`$ in a part of the range $`\mathrm{\Delta }`$. \[Otherwise the condition $`\varphi (z^{})=\mathrm{const}\times w_i`$ is not a limitation, in the WKB approximation.\] It follows from the analysis below that $`\mathrm{\Delta }`$ includes one of the rays $`\mathrm{arg}z^{}=0`$ or $`\pi `$, along which $`\varphi (z^{})`$ is oscillating as $`\mathrm{exp}(i\alpha z_{}^{}{}_{}{}^{3/2})`$ or $`\mathrm{exp}(i\alpha z_{}^{}{}_{}{}^{3/2})`$.

For concreteness, we will assume that $`\mathrm{\Delta }`$ contains the semi-axis $`\mathrm{arg}z^{}=0`$, and that between $`\mathrm{arg}z^{}=0`$ and the cut in Fig. 5 $`\varphi (z^{})=Cw_2(z^{})\mathrm{exp}[i|\alpha |z_{}^{}{}_{}{}^{3/2}]`$, where $`C`$ is a constant (it is not incorporated into Fig. 5).

Since the function $`\varphi (z^{})`$ is single-valued, if we cross the cut in Fig. 5 by incrementing $`\mathrm{arg}z^{}`$, $`\varphi (z^{})`$ becomes equal to $`iCw_1(z^{})`$. It remains exponentially small as $`\mathrm{arg}z^{}`$ grows up to $`2\pi /3`$, including the Stokes line $`\mathrm{arg}z^{}=\pi /3`$. Then behind the anti-Stokes line at $`\mathrm{arg}z^{}=2\pi /3`$, the function $`w_1`$ becomes exponentially big. It is important that, on the Stokes line $`\mathrm{arg}z^{}=\pi `$, one has to take into account the admixture to $`\varphi `$ of an exponentially small term $`w_2(z^{})`$. This can be seen from the explicit solution for $`\varphi `$. The need to incorporate this term can also be understood by noticing that, when we increment $`\mathrm{arg}z^{}`$ by $`2\pi `$, we have to recover the original asymptotic form of $`\varphi `$. This latter argument explains the coefficient at $`w_2`$ in Fig. 5.

On the anti-Stokes line $`\mathrm{arg}z^{}=4\pi /3`$, the values of $`|w_1|`$ and $`|w_2|`$ become equal to each other, and $`\varphi (z^{})`$ is primarily determined by $`w_2`$ for larger $`\mathrm{arg}z^{}`$. After $`\mathrm{arg}z^{}`$ crosses the Stokes line $`5\pi /3`$, the exponentially small term $`w_1`$ in $`\varphi (z^{})`$ disappears, according to the Stokes prescription .

The wave function which connects to the outgoing wave on the caustic is that of the metastable state. Since the asymptotic behavior of $`\varphi (z^{})`$ for large $`|z^{}|`$ is given by $`w_2(z^{})`$, the asymptotic behavior of the wave function of the metastable state within the range $`2\pi /3<\mathrm{arg}z^{}<4\pi /3`$ is given by $`w_1(z^{})`$. We note that, for the radiation boundary condition, the switching between the wave functions occurs only on one of the anti-Stokes lines.

### D Switching between the branches of the wave function

Switching between the WKB wave functions of the localized state and the outgoing wave packet is an important observable consequence of the analysis in the previous subsection. It is due to branching of the WKB action. The switching manifold starts on the caustic and is given by the condition Im $`S_1(𝐫)=\mathrm{Im}S_2(𝐫)`$, where $`S_{1,2}`$ are the actions for the corresponding WKB branches. On the opposite sides of the switching manifold one of the WKB wave functions is exponentially bigger than the other.

In the presence of a magnetic field, caustics go through real space along the lines given by the condition $`J(𝐫)=0`$, Im $`𝐫=0`$, with $`𝐫𝐫(x_{1,2}(0),t)`$ being a point on the Hamiltonian trajectory (2), (7). The switching manifold in real space is a surface which starts from the caustic line and goes away from it in one direction. Although the wave function is continuous on this surface, the derivative of its logarithm sharply changes from $`\mathbf{}S_1`$ to $`\mathbf{}S_2`$.

The exit point $`𝐫`$ in the configuration space where the escaped particle emerges from under the barrier is determined by the intersection of the classical escape trajectory $`𝐫_{\mathrm{cl}}(t)`$ and the switching surface. This point can be found from the global analysis of the WKB wave function. It does not lie on the caustic, nor is it given by the condition $`𝐩=0`$ or equivalently, $`U(𝐫)=E`$. In a 2D system the caustic pierces real space at a point, and the switching surface becomes a line. An example of the caustic, the switching line, and the exit point for a 2D system is shown in Fig. 3.

In the absence of a magnetic field, as we mentioned before, the caustic for the tunneling trajectories is a surface rather than a line in real configuration space. The switching surface in real space should then coincide with the caustic. The exit point is the turning point $`𝐩=\mathrm{𝟎}`$. The whole topology thus differs qualitatively from that for $`B0`$.

## V Tunneling from a correlated 2D electron system

We now apply the general approach to tunneling from a correlated 2D electron system transverse to a magnetic field and illustrate the occurrence of the singularities discussed above. We will use the simple but nontrivial model of the electron system discussed in Sec. II and described by Eqs. (5), (6), and its generalization to the case where the in-plane symmetry is broken. We will assume that the magnetic field is parallel to the electron layer, and choose the $`y`$-axes along the field $`𝐁`$.

### A A model with in-plane symmetry

For an electron in the potential (5), (6), classical motion along the $`𝐁\widehat{𝐲}`$ axis is decoupled from the motion in the $`(x,z)`$-plane. The WKB tunneling problem then becomes two-dimensional, with complex classical trajectories (2) lying in this plane. The Hamiltonian equations are linear, and we can find the trajectories explicitly. We can also explicitly find the tunneling exponent and analyze  its dependence on the two dimensionless parameters $`\omega _0\tau _0`$ and $`\omega _c\tau _0`$, where $`\tau _0`$ is the tunneling time in the absence of the magnetic field, $`\tau _0=2mL/\gamma `$. The expression for the tunneling exponent also follows as a limiting case from the result of the next subsection. Here we will discuss the structure of the WKB action.

The symmetry of the potential $`U(x,y,z)=U(\pm x,\pm y,z)`$ in (5) gives rise to a specific symmetry of the set of the Hamiltonian trajectories,

$$tt^{},xx^{},zz^{},SS^{},$$
(18)
and of the singularities of this set. In particular, the caustic where there are connected the outgoing wave packet and the intrawell wave function goes through real plane $`(x,z)`$ at a point $`x_{\mathrm{caust}}=0`$ on the symmetry axis. The $`z`$-coordinate of this point depends on the magnetic field, with $`z_{\mathrm{caust}}=L`$ for $`B=0`$. The form of the action for $`zz_{\mathrm{caust}}`$ near the caustic in the symmetry plane $`x=0`$ is shown in Fig. 6. As seen from the inset, the slope of the action is $`\mathrm{Im}S/z>0`$ at $`z_{\mathrm{caust}}`$, in contrast to the 1D case where the slope is equal to zero at the turning point $`z_t`$. We note that the branches 1 and 2 are formed by the trajectories that go through real space at times -Im $`t`$ being, respectively, smaller and larger than -Im $`t_{\mathrm{caust}}`$ for the same trajectories \[$`t_{\mathrm{caust}}`$ is given by Eq. (12)\].

Using the explicit form of the trajectories (2) with the initial conditions (8), one can find the complex caustic $`z^{}(x,z)=0`$ near $`x=0,z=z_{\mathrm{caust}}`$. It has the form $`zz_{\mathrm{caust}}=iC^{}x`$ with real $`C^{}/B>0`$, cf. Fig. 4. It is seen from Fig. 6 that the singular parts of Im $`S_{1,2}`$ behave as $`(z_{\mathrm{caust}}z)^{3/2}`$ near the caustic. Therefore we can choose the coordinate $`z^{}`$ that gives the distance from the caustic as $`zz_{\mathrm{caust}}iC^{}x`$. The branching behavior near the caustic is then described by Fig. 5. The range $`\pi \mathrm{arg}z^{}<2\pi `$ corresponds to real $`z`$ and real positive $`x`$.

Close to the caustic, only one branch of the action describes the wave function for negative $`x`$ and real $`z`$ (the upper half of the complex-$`z^{}`$ plane). For positive $`x`$, we should keep both branches, and depending on $`x,z`$, the WKB wave function is given by the branch with the smaller Im $`S`$. Near the switching line where Im $`S_1=\mathrm{Im}S_2`$ the total wave function is given by a linear combination of the two WKB solutions.

In Fig. 6 for $`zz_{\mathrm{caust}}`$, the action branch 1 describes the tail of the intrawell wave function, and the branch 2 corresponds to the wave “reflected” from the caustic. The branch 2 is nonmonotonic in $`z`$ for $`x=0`$, with a minimum at $`z_m<z_{\mathrm{caust}}`$. By symmetry (18), the momentum component $`p_x`$ is real for $`x=0`$, whereas $`p_z=0`$ at $`z_m`$. Therefore the point $`z=z_m`$, $`x=0`$ belongs to the classical escape trajectory shown in Fig. 3(b). Again by symmetry, this is also the point where the escape trajectory comes closest to the well at $`z=0`$. However, this is not the exit point for the tunneling particle in the configuration space. Indeed, the wave function at this point is determined by the branch 1, because Im $`S_1<\mathrm{Im}S_2`$.

Cross-sections of the action surfaces by planes $`z=`$ const are shown in Fig. 7. For $`zz_{\mathrm{caust}}`$, both branches Im $`S_{1,2}`$ are symmetrical in $`x`$. However, the branch 2 is nonmonotonic in $`x`$ for $`z>z_m`$. It has a local maximum at $`x=0`$ and two symmetrical minima. These minima lie on the classical trajectory shown in Fig. 3(b). For $`z=z_m`$, the maximum and the minima merge together. We note that at this point Im $`S_2x^4`$ near the minimum.

Behind the caustic in real space, $`z>z_{\mathrm{caust}}`$, the action $`S(x,z)`$ on one of the two branches is equal to $`S^{}(x,z)`$ on the other branch (cf. Fig. 7). At their minima with respect to $`x`$, the values of Im $`S(x,z)`$ are independent of $`z`$. These minima lie on the escape classical trajectory, as seen from the comparison with Im $`S`$ for $`z<z_{\mathrm{caust}}`$. The branches of the action describe the wave packets incident on the barrier from large $`z`$ and the emitted wave packet. Only the latter is physically meaningful for the problem of tunneling escape.

As discussed above, switching between the action branches occurs for positive $`x`$ where the branches of Im $`S`$ cross each other. The resulting action, which determines the WKB wave function, is shown in Fig. 7 by solid lines. The switching line thus obtained coincides near the caustic with the anti-Stokes switching line discussed in the previous section.

The escaped particle can be observed as a semiclassical object if $`\mathrm{Im}S_2(𝐫_{\mathrm{cl}})<\mathrm{Im}S_1(𝐫_{\mathrm{cl}})`$. It “shows up” at the point where the classical escape trajectory intersects the switching line, cf. Fig. 3. The exit point is located for $`x>0`$ even though the potential (5) is symmetric. This is a consequence of the symmetry-breaking by a magnetic field.

### B A non-symmetric model

The problem considered in the previous section for the quadratic in $`x,y`$ potential $`U(𝐫)`$ (5) can be solved differently. The trick is  to make a canonical transformation to the new coordinate $`p_x`$ and the conjugate momentum $`x`$. The kinetic energy then becomes $`[m^2\omega _0^2x^2+p_z^2]/2m`$ and is independent of the new coordinates and the magnetic field. The time-reversal symmetry is thus “restored”, and the problem is mapped onto the standard problem of tunneling in the 2D potential $`U_0(z)+(p_x+m\omega _cz)^2/2m`$.

The general method discussed in this paper is not limited to potentials with these special properties. In this subsection we illustrate how this method works where variables do not separate. To this end, we consider a 2D problem of tunneling transverse to the field $`𝐁\widehat{𝐲}`$ in the potential

$$U(x,z)=\frac{1}{2}m\omega _0^2x^2+\mu xz+\frac{\gamma ^2}{2m}\left(1\frac{z}{L}\right)(z>0),$$
(19)
which differs from the potential discussed earlier by the term $`\mu xz`$. In the problem of tunneling from a 2D electron system, this term mimics the dependence of the tunneling potential on the displacements of neighboring electrons, see Appendix A.

The term $`\mu xz`$ breaks the symmetry of the Hamiltonian trajectories (18). However, the Hamiltonian equations (2) are still linear, and we explicitly solved them. The results for Im $`S`$ and the classical escape trajectory obtained using the initial conditions (8) are shown in Figs. 8, 9. Because of broken symmetry, the branching point of the action in real space, where the caustic of the set of trajectories (2) goes through real plane $`(x,z)`$, lies at $`x_{\mathrm{caust}}0`$. It is marked by the cross on Fig. 8(d). Its position depends on $`\mu `$ and other parameters of the system. Similarly, the time where the caustic crosses the real space has both real and imaginary part, in contrast to the case $`\mu =0`$ where it was purely imaginary, see Fig. 3.

For $`\mu 0`$, the surfaces Im $`S(x,z)`$ become asymmetric. The general structure of the solution, however, remains the same as in the case $`\mu =0`$. This can be seen by comparing the cross-sections of the action in Figs. 7 and 8. In both figures, the cross-sections in (a), (b), and (c) refer to the planes $`z<z_{\mathrm{caust}},z=z_{\mathrm{caust}}`$, and $`z>z_{\mathrm{caust}}`$, respectively. As in the symmetrical case, the branches 1 and 2 for $`\mu 0`$ are formed by the trajectories in complex time $`t`$, with -Im $`t`$ being, respectively, smaller and bigger than -Im $`t_{\mathrm{caust}}`$. At the point $`x_{\mathrm{caust}},z_{\mathrm{caust}}`$, the branches 1 and 2 touch each other.

For $`z<z_{\mathrm{caust}}`$ there occurs switching between the branches. It can be analyzed in the same way as for $`\mu =0`$. The WKB wave function is determined by the branch in Fig. 8 shown with the solid line. The switching line starts from the branching point $`x_{\mathrm{caust}},z_{\mathrm{caust}}`$ and goes in the direction of positive $`x`$.

The branch 2 has two minima as a function of $`x`$ for $`z_m<z<z_{\mathrm{caust}}`$, which are asymmetric for $`\mu 0`$. However, their depths, i.e. the minimal values of Im $`S`$ on this branch, remain equal to each other and are the same in all cross-sections $`z=`$ const. These minima lie on the classical trajectory along which the electron escapes. At $`z=z_m`$ they merge together, and Im $`S`$ becomes quartic in $`x`$ near the minimum. The value $`z_m`$ shows how close the escape trajectory comes to the localized intrawell state.

The classical trajectory becomes observable in configuration space once it crosses the switching line. The shape of the trajectory and the exit point for several values of the asymmetry parameter $`\mu `$ are shown in Fig. 9. The outgoing wave packet is Gaussian near the maximum (Im $`S`$ is parabolic near the corresponding minimum).

By solving the Hamiltonian equations (2) for the potential (19), we obtained the following expression for the tunneling exponent (11),

$$R=2\gamma L[\tau _{\mathrm{rd}}+\nu _0\kappa (\tau _{\mathrm{rd}})]$$
(20)
Here, $`\tau _{\mathrm{rd}}`$ is the imaginary part of the time to reach the classical escape trajectory \[see Fig. 3(a)\] in the units of that same time for $`B=\mu =0`$, which is given by $`\tau _0=2mL/\gamma `$. Along with the function $`\kappa `$ in (20), $`\tau _{\mathrm{rd}}`$ can be found from the equation

$`\kappa (\tau _{\mathrm{rd}})`$ $`{\displaystyle \frac{\nu _0^2(\tau _{\mathrm{rd}}\mathrm{cos}\nu _{_{}}\tau _{\mathrm{rd}}\nu _{_{}}^{}{}_{}{}^{1}\mathrm{sin}\nu _{_{}}\tau _{\mathrm{rd}})+\nu _{_{}}^{}{}_{}{}^{2}\mathrm{cos}\nu _{_{}}\tau _{\mathrm{rd}}\nu _{_{}}(1\nu _0+\nu _0\tau _{\mathrm{rd}})\mathrm{sin}\nu _{_{}}\tau _{\mathrm{rd}}}{\stackrel{~}{\mu }^2(\nu _0\mathrm{cos}\nu _{_{}}\tau _{\mathrm{rd}}\nu _{_{}}\mathrm{sin}\nu _{_{}}\tau _{\mathrm{rd}})}}`$ (22)
$`={\displaystyle \frac{(\nu _0^2\tau _{\mathrm{rd}}\nu _{+}^{}{}_{}{}^{2})\mathrm{cosh}\nu _+\tau _{\mathrm{rd}}+\nu _+(1\nu _0+\nu _0\tau _{\mathrm{rd}}\nu _0^2/\nu _{+}^{}{}_{}{}^{2})\mathrm{sinh}\nu _+\tau _{\mathrm{rd}}}{\stackrel{~}{\mu }^2(\nu _0\mathrm{cosh}\nu _+\tau _{\mathrm{rd}}+\nu _+\mathrm{sinh}\nu _+\tau _{\mathrm{rd}})}}.`$
Here, $`\stackrel{~}{\mu }=\mu \tau _0^2/m`$ is the dimensionless asymmetry parameter. The dimensionless frequencies $`\nu _0=\omega _0\tau _0`$, $`\nu _c=\omega _c\tau _0`$, and $`\nu =(\nu _c^2+\nu _0^2)^{1/2}`$ characterize the motion under the barrier, $`\nu _{_\pm }^{}{}_{}{}^{2}=\pm \nu ^2/2+\sqrt{\nu ^4/4+\stackrel{~}{\mu }^2}`$.

The tunneling exponent $`R`$ depends on the interrelation between the in-plane electron dynamics, which is characterized by the frequency $`\omega _0`$, the cyclotron frequency $`\omega _c`$, the tunneling time $`\tau _0`$, and the asymmetry-induced mixing of in-plane and out-of-plane motions $`\stackrel{~}{\mu }`$. The dependence of $`R`$ on the asymmetry parameter and on the magnetic field is shown in Fig. 10. For $`\omega _c=\stackrel{~}{\mu }=0`$, we have $`\tau _{\mathrm{rd}}=1`$, and $`R=4\gamma L/3`$. The exponent $`R`$ increases with the magnetic field. A qualitative result of the in-plane confinement is that it eliminates the localization in the linear potential $`U_0(z)`$ due to the magnetic barrier, i.e. the divergence of $`R`$ for $`\omega _c\tau _0>1`$. The effect of the magnetic field on $`R`$ becomes small for strong confinement, $`\omega _0\tau _01`$ and $`\omega _0\omega _c`$.

The asymmetry results in lowering of the tunneling barrier for $`B=0`$ and the corresponding increase of the tunneling rate. This can be qualitatively understood, since a displacement in the $`x`$-direction with $`\mu x<0`$ increases the effective force in the $`z`$-direction which pulls the electron away from the layer. For small asymmetry, $`\stackrel{~}{\mu }1`$, the correction to $`R`$ is quadratic in $`\stackrel{~}{\mu }`$. In the limit of a thin and high tunneling barrier for $`x=0`$ or soft in-plane vibrations, where $`\omega _0\tau _01`$, tunneling is most likely to occur in the barrier, which is adiabatically prepared by the in-plane displacement $`x`$ (the “completely adjusted” barrier). The expression for $`R`$ takes a form

$$R=2\gamma L(3\nu _0/\stackrel{~}{\mu }^2)^{1/3}(\omega _0\tau _01),$$
(23)
it depends on the in-plane frequency $`\omega _0`$ and on $`\mu `$ nonanalytically. The role of the asymmetry increases with the magnetic field, as seen from Fig. 10.

In terms of comparison with the currently available experimental data on tunneling from a correlated many-electron system , of utmost interest is the situation where the asymmetry is small. The observed dependence of $`R`$ on $`B`$ did not show the divergence expected for an unconfined electron. The simple model (5) provides a qualitative explanation of the experiment . An excellent quantitative agreement, without adjustable parameters, was achieved by incorporating the curvature of the potential $`U_0(z)`$ due to electron correlations.

## VI The path-integral formulation in a magnetic field

In the absence of a magnetic field, the problem of tunneling decay is often considered using the instanton technique . This technique applies if the potential well is parabolic near the minimum and thermalization inside the well occurs much faster than escape from the well (in the case of 2D electron systems, both conditions are often violated ). Because the Schrödinger equation for metastable states has to be solved with the radiation boundary condition, the energies of these state acquire small imaginary parts, and so does the partition function $`Z`$. The escape rate $`W`$ for finite temperatures is simply related to Im $`Z`$,

$$W2T\mathrm{Im}Z/\mathrm{Re}Z$$
(24)
(we have set $`k_B=1`$).

The partition function is given by the integral over periodic paths $`𝐫(\tau )`$ in imaginary time ,

$$Z=_{𝐫(0)=𝐫(\beta )}𝒟𝐫(\tau )\mathrm{exp}\left[S_E[𝐫(\tau )]\right],$$
(25)
where $`\beta =T^1`$ and $`S_E`$ is the Euclidean action (the action in imaginary time). It is real for $`B=0`$ and for real trajectories $`𝐫(\tau )`$.

The general expressions (24), (25) should also apply in the presence of a magnetic field. However, the Euclidean action for an electron

$$S_E=_0^\beta 𝑑\tau \left[\frac{m}{2}\left(\frac{d𝐫}{d\tau }\right)^2+U(𝐫)+i\frac{e}{c}𝐀(𝐫)\frac{d𝐫}{d\tau }\right]$$
(26)
is now complex. Therefore the standard way  of evaluating the escape rate has to be revised, except for special symmetric cases like the one discussed in Sec. V A, where one can change to new varibales in which $`S_E`$ becomes real . Although we will often refer to the action functional of the form (26), much of the results below apply also to a more general retarded Euclidean action, which is of interest for systems coupled to a bath.

In the spirit of the WKB approximation, the path integral (25) will be evaluated by the steepest descent method. The extremal paths $`𝐫(\tau )`$ of the functional $`S_E`$ satisfy the equation

$$m\frac{d^2𝐫}{d\tau ^2}=U(𝐫)+i\frac{e}{c}\left[\frac{d𝐫}{d\tau }\times 𝐁\right].$$
(27)
The equation of motion (27) has to be solved with the periodic boundary condition $`𝐫(0)=𝐫(\beta )`$. Note that the sign of the potential has been inverted compared to the case of classical motion in real time.

For low temperatures and for the potential $`U(𝐫)`$, which is parabolic near its intrawell minimum $`𝐫_{\mathrm{well}}`$, Eq. (27) has a solution $`𝐫(\tau )=𝐫_{\mathrm{well}}`$, with $`S_E=0`$. It gives the real part of the partition function, see below. As in the case $`B=0`$, the imaginary part of $`Z`$ is determined by another solution of (27), which is of the bounce type. This solution, $`𝐫_\mathrm{b}(\tau )`$, starts near $`𝐫_{\mathrm{well}}`$, slides downhill in the inverted potential $`U(𝐫)`$, and in time $`\beta `$ comes back. For $`B=0`$ the corresponding path is a symmetric real trajectory, $`𝐫_\mathrm{b}(\tau )=𝐫_\mathrm{b}(\beta \tau )`$, which bounces off the turning point $`\dot{𝐫}_\mathrm{b}(\beta /2)=0`$.

For $`B0`$, because of broken time-reversal symmetry, the path $`𝐫_\mathrm{b}(\tau )`$ is complex, and the velocity along this path does not become equal to zero. The path is not symmetrical in time, because the replacement $`\tau \tau `$ changes Eq. (27). However, if we simultaneously change $`ii`$, the equation remains unchanged. Therefore the bounce-type path has the symmetry

$$𝐫_\mathrm{b}(\tau )=𝐫_\mathrm{b}^{}(\beta \tau ).$$
(28)
An immediate and very important consequence of Eq. (28) is that the value of $`S_E(𝐫_\mathrm{b})`$ for the bounce-type path is real. This value gives the tunneling exponent, see below.

### A The eigenvalue problem

The prefactor in $`Z`$ can be found by integrating over the tubes of paths around the extremal paths. It can be done by expanding $`S_E`$ in deviations from the extremal paths to the second order, and then expanding $`𝐫(\tau )𝐫_{\mathrm{well}}`$ and $`𝐫(\tau )𝐫_\mathrm{b}(\tau )`$ in the eigenfunctions $`𝝍_n(\tau )`$ of the appropriate eigenvalue problem,

$`\widehat{}𝝍_n{\displaystyle _0^\beta }𝑑\tau ^{}\widehat{𝐅}(\tau ,\tau ^{})𝝍_n(\tau ^{})=\lambda _n𝝍_n(\tau ),`$ (29)
$`\widehat{F}_{ij}(\tau ,\tau ^{})=\delta ^2S_E/\delta r_i(\tau )\delta r_j(\tau ^{}).`$ (30)
Here, the derivatives of the action are calculated on the corresponding extremal trajectory $`𝐫_{\mathrm{well}}`$ and $`𝐫_\mathrm{b}(\tau )`$, and periodic boundary conditions are assumed. The operator $`\widehat{F}`$ is simplified for a non-retarded action (26), $`\widehat{F}_{ij}(\tau ,\tau ^{})=\delta (\tau \tau ^{})\widehat{\mathrm{f}}_{ij}(\tau )`$.

For $`B=0`$, the operator $`\widehat{}`$ is Hermitian, with $`\widehat{\mathrm{f}}_{ij}=m\delta _{ij}(d^2/d\tau ^2)+^2U/r_ir_j`$, if the action is given by (26). Therefore the functions $`𝝍_n(\tau )`$ form complete and orthogonal sets for each extreme trajectory (27), and the eigenvalues $`\lambda _n`$ are real.

For $`B0`$, the operator $`\widehat{}`$ becomes non-Hermitian. For example, in the case of a uniform magnetic field in (26), $`\widehat{\mathrm{f}}_{kl}(\tau )`$ acquires an extra term $`i(e/c)ϵ_{klj}B_j(d/d\tau )`$ ($`ϵ_{klj}`$ is the Levi-Civita symbol). Therefore some of the eigenvalues $`\lambda _n`$ become complex. The eigenvectors $`𝝍_n`$ with different $`n`$ are orthogonal not to each other, but to the eigenvectors $`\mathit{\varphi }_n`$ of the Hermitian conjugate operator,

$`{\displaystyle _0^\beta }𝑑\tau ^{}\widehat{𝐅}^{}(\tau ,\tau ^{})\mathit{\varphi }_n(\tau ^{})=\lambda _n^{}\mathit{\varphi }_n(\tau ).`$
Taking into account the symmetry (28) of the extremal trajectories, we find that

$$\widehat{𝐅}^{}(\tau ,\tau ^{})=\widehat{𝐅}(\beta \tau ,\beta \tau ^{})=\widehat{𝐅}^{}(\tau ,\tau ^{}).$$
(31)
The energy spectra for several complex Hamiltonians with the symmetry, which is similar to (31) (and was called the $`𝒫𝒯`$-symmetry), were investigated earlier numerically and using the WKB approximation.

The symmetry (31) has several consequences. First, it shows that $`𝝍_n(\tau )=\alpha _n\mathit{\varphi }_n^{}(\tau )`$, where $`\alpha _n`$ is a constant. This means that, with proper normalization, the orthogonality relation becomes

$$_0^\beta 𝑑\tau 𝝍_m(\tau )𝝍_n(\tau )=\delta _{mn}$$
(32)
(here, we assumed that the eigenvalues are nondegenerate; for degenerate eigenvalues, the condition (32) can be satisfied by choosing appropriate linear combinations of the eigenfunctions with same $`\lambda _n`$).

It also follows from Eq. (31) that, if $`𝝍_n(\tau )`$ is an eigenfunction of (29) with an eigenvalue $`\lambda _n`$, then $`𝝍_n^{}(\beta \tau )`$ is also an eigenfunction of the same boundary value problem, but with the eigenvalue $`\lambda _n^{}`$. This means that a part of the eigenvalues $`\lambda _n`$ in (29) are real, whereas another part are pairs of complex conjugate numbers.

Pairs of complex conjugate eigenvalues are formed in the following way. For $`B=0`$, all eigenvalues are real. With increasing $`B`$ some eigenvalues approach each other, pairwise, while still remaining real. Eventually they merge, and for larger $`B`$ become complex conjugate, as described in Appendix B. For 1D Schrödinger-type equation with complex Hamiltonians such behavior with the varying control parameter was indeed observed numerically .

#### 1 Eigenvalues near the potential well

As an illustration, we consider the eigenvalue problem near $`𝐫_{\mathrm{well}}`$ for the action functional (26). Here, Eq. (29) becomes linear, and the eigenfunctions $`𝝍_n(\tau )`$ can be sought in the form of linear combinations of $`\mathrm{exp}(\pm i\omega _n\tau )`$, with $`\omega _n=2\pi n/\beta `$. The eigenvalues are obtained from the equation

$$det\left[(m\omega _n^2\lambda _{n\nu })\delta _{kl}+m\mathrm{\Omega }_{kl}^2\frac{e}{c}\omega _nϵ_{klj}B_j\right]=0,$$
(33)
where $`m\mathrm{\Omega }_{kl}^2=[^2U/r_kr_l]_{𝐫_{\mathrm{well}}}`$, and $`𝐁`$ is the magnetic field at $`𝐫_{\mathrm{well}}`$. The subscript $`\nu `$ enumerates the eigenvalues $`\lambda `$ for a given Matsubara frequency.

If, for example, $`𝐁`$ is pointing along a principal axes of the tensor $`\mathrm{\Omega }_{kl}^2`$ (say, the axes 1), then we have $`\lambda _{n\mathrm{\hspace{0.17em}1}}/m=\omega _n^2+\mathrm{\Omega }_1^2`$, and

$`m^1\lambda _{n\mathrm{\hspace{0.17em}2},3}=`$ $`\omega _n^2+{\displaystyle \frac{1}{2}}(\mathrm{\Omega }_2^2+\mathrm{\Omega }_3^2)`$ (35)
$`\pm {\displaystyle \frac{1}{2}}\left[(\mathrm{\Omega }_2^2\mathrm{\Omega }_3^2)^24\omega _c^2\omega _n^2\right]^{1/2}`$
where $`\mathrm{\Omega }_\nu ^2>0`$ are the principal values of the tensor $`\mathrm{\Omega }_{kl}^2`$. Clearly, the eigenvalues $`\lambda _{n\mathrm{\hspace{0.17em}2},3}`$ are complex conjugate pairs, for large enough $`\omega _n^2\omega _c^2`$.

Eq. (35) shows explicitly also how pairs of complex eigenvalues emerge with varying magnetic field as a result of merging of adjacent real eigenvalues, as discussed for the general case in Appendix B.

#### 2 Eigenvalues for the bounce trajectory

A specific feature of the eigenvalue problem (29) for the bounce trajectory $`𝐫_b(\tau )`$ at low temperatures is that one of the eigenvalues is $`\lambda _1=0`$. It corresponds to the eigenfunction $`𝝍_1(\tau )\dot{𝐫}_\mathrm{b}(\tau )`$. For $`B=0`$, the vector function $`𝝍_1(\tau )\dot{𝐫}_\mathrm{b}(\tau )`$ has one zero for all components. Therefore, as can be shown using standard arguments, it is the eigenfunction of the first excited state of the multicomponent Schrödinger-type equation (29) (except for a nongeneric case where the components of $`𝝍`$ separate). Since the eigenvalue problem (29) is Hermitian for $`B=0`$, all eigenvalues $`\lambda _n`$ with $`n2`$ are positive, and the eigenvalue of the ground state is negative, $`\lambda _0<0`$. .

We are not aware of the oscillation theorem for non-Hermitian problems. However, since $`𝝍_1(\tau )\dot{𝐫}_\mathrm{b}(\tau )`$ is an eigenfunction with zero eigenvalue even for $`B0`$, as $`B`$ increases from zero, the eigenvalue $`\lambda _1`$ does not merge with other real eigenvalues to form a pair of complex conjugate eigenvalues. Therefore, pairs of complex conjugate eigenvalues will be only formed from $`\lambda _n`$ that were positive for $`B=0`$. The negative root $`\lambda _0`$ will remain real and negative. In principle, as a result of coalescence of complex conjugate eigenvalues, there may emerge pairs of negative real eigenvalues. However, the total number of negative real eigenvalues will be odd.

### B The prefactor

We are now in a position to calculate the prefactor in the partition function $`Z`$. The standard step is to expand the deviation $`\delta 𝐫(\tau )`$ of the integration path in (25) from the extreme trajectory $`𝐫_{\mathrm{well}}`$ or $`𝐫_\mathrm{b}`$ in terms of the eigenfunctions $`𝝍_n`$ of the corresponding eigenvalue problem, $`\delta 𝐫(\tau )=c_n𝝍_n(\tau )`$. With account taken of the orthogonality condition (32), the increment $`\delta S_E`$ of the Euclidean action related to the deviation of the trajectory $`\delta 𝐫`$ then becomes $`\delta S_E=\lambda _nc_n^2/2`$.

The above expansion assumes that the set $`\{𝝍_n\}`$ is complete. The completeness is known for $`B=0`$, where the eigenvalue problem (29) is Hermitian. As $`B`$ changes, the number of states does not change. From the orthogonality condition (32), none of the wave functions becomes a linear combination of other wave functions. This makes us believe that the functions $`𝝍_n`$ form a complete set even for $`B0`$ and justifies the above expansion.

The path integral (25) can be obtained as a limit $`N\mathrm{}`$ of integrals over $`d𝐫(\tau _k)`$ at discretized instants of time $`\tau _k=k\mathrm{\Delta }\tau ,\mathrm{\Delta }\tau =\beta /N`$. In the standard way, we change to integration over $`dc_n`$. Because of the orthogonality relation (32), the determinant $`det[𝝍_n(\tau _k)]`$ of the transformation of variables is real and is equal to $`\pm (\mathrm{\Delta }\tau )^{N/2}`$. Integration of $`\mathrm{exp}(\delta S_E)`$ over $`dc_n`$ gives $`\mathrm{const}\times _n\lambda _n^{1/2}`$.

Let us now consider the contribution to the partition function $`Z_{\mathrm{well}}`$ from trajectories close to the potential minimum $`𝐫_{\mathrm{well}}`$. The corresponding eigenvalues $`\lambda _n^{(\mathrm{well})}`$ are either positive or belong to complex conjugate pairs, cf. (35). Therefore $`Z_{\mathrm{well}}`$ is real. Since $`S_E[𝐫_{\mathrm{well}}]=0`$, there is no exponentially small factor in $`Z_{\mathrm{well}}`$. This term gives the partition function for low-lying intrawell excitations in the presence of the magnetic field.

In evaluating the contribution $`Z_\mathrm{b}`$ from paths close to the bounce trajectory, special care has to be taken of the eigenvalue $`\lambda _1^{(\mathrm{b})}=0`$. A standard analysis  shows that integration over $`dc_1`$ gives the factor $`\beta `$ in $`Z_\mathrm{b}`$. The positive and complex conjugate eigenvalues $`\lambda _n^{(\mathrm{b})}`$ give a real positive factor in $`Z_\mathrm{b}`$, whereas the negative eigenvalue $`\lambda _0^{(\mathrm{b})}`$ (or an odd number of negative eigenvalues) make $`Z_\mathrm{b}`$ purely imaginary. In addition, $`Z_\mathrm{b}`$ contains the exponential factor $`\mathrm{exp}\{S_E[𝐫_\mathrm{b}(\tau )]\}`$. Overall, this gives the tunneling rate (24)

$$W2T|Z_\mathrm{b}|/Z_{\mathrm{well}}\mathrm{exp}\{S_E[𝐫_\mathrm{b}(\tau )]\}.$$
(36)
Eq. (36) shows that the instanton technique can be used in the presence of a magnetic field in spite of the fact that the field breaks time-reversal symmetry. The actual calculation is in many respects different from that for $`B=0`$. In particular, the bounce trajectory is complex. The eigenvalues which determine the prefactor should be found from a non-Hermitian eigenvalue problem. They may be complex, in which case they form pairs of complex-conjugate numbers. The bounce trajectory touches the real escape trajectory at a real point $`𝐫_\mathrm{b}(\beta /2)`$, with a finite real velocity $`\dot{𝐫}_\mathrm{b}(\beta /2)`$. However, from our general WKB analysis of tunneling, it follows that this is not the point where the particle “shows up” as a semiclassical wave packet.

## VII conclusions

In conclusion, the problem of tunneling in a magnetic field can be solved in the semiclassical limit by analyzing the Hamiltonian trajectories of the particle in complex phase space and time. The boundary conditions are determined by the intrawell wave function and its analytic continuation. This approach allows one to find both the tunneling exponent and the tail of the wave function of the localized state. It does not require to consider either a part of the potential or the magnetic field as a perturbation, and it can be applied to a three-dimensional potential of a general form.

The escape rate is generally exponentially smaller than the probability for a particle to reach the boundary of the classically accessible range $`U(𝐫)=E`$. The escaped particle “shows up” from the tunneling barrier with finite velocity and behind the surface $`U(𝐫)=E`$. The connection of the decaying and propagating waves occurs on caustics of the set of the Hamiltonian trajectories, where the action is branching. The caustics are complex surfaces in 3D space. In the presence of a magnetic field, they go through real space along lines (instead of surfaces, for $`B=0`$).

An interesting feature of tunneling in a magnetic field is the occurrence of a switching surface, where there merge different WKB branches of the wave function. The slope of the logarithm of the wave function sharply changes at the switching surface, from its value on one of the branches to that on the other branch. The escaped particle first shows up as a propagating semiclassical wave packet on the switching surface. It happens where the classical escape trajectory crosses the switching surface.

The switching surfaces are observable via experiments in which the electron density is measured, although such experiments are extremely hard to do. For tunneling from 2D electron systems in heterostructures, one can think of NMR experiments with samples that contain delta-doped layers of nuclear spins (of an isotope that differs from that in the host material). If the nuclei are sufficiently far from the electron layer, they will detect the local electron density on the tail of the wave function and its variation with varying fields. One could also use light-scattering measurements.

Switching between branches of the WKB wave function for $`B0`$ is similar to the switching between different branches of the probability distribution in classical systems away from thermal equilibrium. Such systems lack time-reversal symmetry, as do also quantum systems in a magnetic field. The tail of the classical distribution is formed by infrequent fluctuations. Fluctuational paths to a given state from the equilibrium position (attractor) form a narrow tube centered at the most probable path. This path is given by a solution of the variational problem of finding the maximum of the logarithm of the probability distribution . In many cases of physical interest the corresponding Euler equations are similar to Eqs. (2). However, in contrast to the tunneling problem, classical optimal paths can be observed. Switching surfaces in the phase space of fluctuating nonequilibrium systems separate areas reached along topologically different optimal paths(a). They have been seen in analog simulations (b).

It follows from the results of this paper that, for potential wells which are parabolic near the minimum, even in the presence of a magnetic field one can still use the instanton technique in order to find the escape rate. However, the bounce trajectory, which gives the tunneling exponent, is now complex. Also in contrast to the $`B=0`$ case, the evaluation of the prefactor requires solving a non-Hermitian boundary value problem, which generally has pairs of complex conjugate eigenvalues.

The results for the model of tunneling from a strongly correlated 2D electron system illustrate the general conclusions about tunneling in a magnetic field. They show that the developed method allows us to find the tunneling rate and the wave function in a generic system which does not have any special symmetry. They also confirm the general conclusions about the structure of singularities related to the branching of the WKB wave function, and the occurrence of the switching line. We found that the tunneling rate in the magnetic field is highly sensitive to the in-plane electron dynamics and exponentially increases when electrons are more strongly confined in the plane. It also increases if electrons in the 2D layer can adjust to the tunneling electron and thus decrease the potential barrier.

This research was supported in part by the NSF through Grant No. PHY-0071059.

## A The many-electron Hamiltonian

A simple and important model which allows us to consider the effect of electron correlations on tunneling from a 2D electron system is the model of a Wigner crystal. In this model, the in-plane electron motion is small-amplitude vibrations about equilibrium positions. Because of strong correlations, exchange effects are not important, and the tunneling electron can be identified. Its tunneling motion is affected by the interaction with other electrons.

We will assume that the equilibrium in-plane position of the tunneling electron is at the origin. Then, in the presence of a magnetic field $`𝐁`$ parallel to the electron layer, the full Hamiltonian is of the form

$$H=\frac{p_z^2}{2m}+U_0(z)+H_v+H_B,$$
(A1)
with

$$H_v=\frac{1}{2}\underset{𝐤,j}{}\left[m^1𝐩_{𝐤j}𝐩_{𝐤j}+m\omega _{𝐤j}^2𝐮_{𝐤j}𝐮_{𝐤j}\right]$$
(A2)
and

$`H_B=`$ $`{\displaystyle \frac{1}{2}}m\omega _c^2z^2\omega _czN^{1/2}{\displaystyle \underset{𝐤,j}{}}[\widehat{𝐁}\times 𝐩_{𝐤j}]_z`$ (A4)
$`+U_{\mathrm{int}}(z,\{𝐮_{𝐤j}\}).`$
Here, $`𝐩_{𝐤j}`$, $`𝐮_{𝐤j}`$, and $`\omega _{𝐤j}`$ are the 2D momentum, displacement, and frequency of the WC phonon of branch $`j`$ ($`j=1,2`$) with a 2D wave vector $`𝐤`$. The in-plane momentum of the tunneling electron is $`N^{1/2}𝐩_{𝐤j}`$ for $`B=0`$. The term $`U_0(z)`$ describes the tunneling barrier \[cf. Eq. (6)\] for the electron at the origin provided all other electrons are at their in-plane lattice sites.

The term $`H_B`$ couples the out-of-plane tunneling motion to lattice vibrations. The problem of many-electron tunneling is thus mapped onto a familiar problem of a particle coupled to a bath of harmonic oscillators . A part of the coupling is due purely to the magnetic field. Another part comes from the term $`U_{\mathrm{int}}`$, which describes the change of the tunneling barrier because of electron vibrations. Its simplest form is given by the lowest-order term of the expansion of the electron energy,

$`U_{\mathrm{int}}(z,\{𝐮_{𝐤j}\})=z{\displaystyle _{𝐤,j}}𝐠_{𝐤j}𝐮_{𝐤j}`$ (A5)
where $`𝐠_{𝐤j}`$ are coupling constants \[for electrons on a thick helium film  the major term in $`U_{\mathrm{int}}`$ is $`z^2`$\]. The coupling (A5) leads to lowering of the tunneling barrier as a result of appropriate displacements of the electrons surrounding the tunneling electron.

The major effect on tunneling comes from high-frequency in-plane vibrations, which have large density of states . Therefore it is not unreasonable to use the Einstein model of the Wigner crystal, in which all vibrations have same frequency $`\omega _0`$. Then, except for the term $`U_{\mathrm{int}}`$, the Hamiltonian (A1) becomes a sum of Hamiltonians of confined in the plane noninteracting electrons. The Hamiltonian of the tunneling electron has the form (3), with the potential $`U(𝐫)`$ given by (5).

For $`B=0`$, the out-of-plane motion of the tunneling electron is coupled only to in-plane displacements of other electrons. In the Einstein model, it means that the out-of-plane motion is decoupled from the in-plane motion of the tunneling electron itself. Instead, it is coupled to an in-plane oscillator with the coordinate given by a (totally symmetric) linear combination of displacements of the other electrons. This maps the problem onto the problem discussed in Sec. V B, with the in-plane electron coordinate $`x`$ in Eq. (19) corresponding to the coordinate of this oscillator, and with $`\mu `$ being a linear combination of the (weighted) coefficients $`𝐠_{𝐤j}`$. Of course, for $`B0`$ this mapping no longer applies, and extra degrees of freedom have to be taken into account. Yet we expect that the model (19) catches much of the qualitative features of many-electron tunneling.

## B Emergence of complex eigenvalues

In this appendix we consider how, with the varying control parameter, two real eigenvalues merge and then become complex. Near this transition, the eigenvalues can be sought by perturbation theory. We start with a value of $`𝐁=𝐁_0`$ (we can also use another control parameter), where the given adjacent eigenvalues $`\lambda _n,\lambda _m`$ are close to each other and are real. For small $`|\delta 𝐁|=|𝐁𝐁_0|`$, the functional $`\widehat{}`$ is close to its value for $`𝐁_0`$, $`\widehat{}\widehat{}_0+\delta \widehat{}`$ in (29). To first order in $`\delta \widehat{}`$, the eigenvalues are given by the expressions $`[\lambda _n(𝐁_0)+\lambda _m(𝐁_0)]/2+\lambda _\pm `$, with

$`\lambda _\pm ={\displaystyle \frac{1}{2}}\left(\delta \widehat{}^{nn}+\delta \widehat{}^{mm}\right)\pm {\displaystyle \frac{1}{2}}\left[\left(\delta \widehat{}^{nn}\delta \widehat{}^{mm}\delta \lambda \right)^2+4\delta \widehat{}^{nm}\delta \widehat{}^{mn}\right]^{1/2},\delta \widehat{}^{nm}=𝝍_n|\delta \widehat{}|𝝍_m.(\mathrm{B1})`$

with the wave functions calculated for $`𝐁_0`$, and with $`\delta \lambda =\lambda _m\lambda _n`$ for $`𝐁=𝐁_0`$.

Because of the symmetry (31), the matrix elements of $`\delta \widehat{}`$ in (B1) are real. However, the product $`\delta \widehat{}^{nm}\delta \widehat{}^{mn}`$ does not have to be positive, and in fact we are interested in the case where it is negative. In this case, instead of level anticrossing, we have the dependence of the eigenvalues on the distance $`\delta \lambda `$ shown in Fig. 11. In the gap the eigenvalues are complex conjugate.

We note that the control parameter in the physical system is not $`\delta \lambda `$, and it may be more interesting to look at the eigenvalues as functions of the amplitude of $`\delta \widehat{}^{nm}`$. Their behavior is similar to what is shown in Fig. 11, if the diagonal and off-diagonal matrix elements depend on the control parameter in the same way. Otherwise, once the eigenvalues become complex with changing control parameter, they do not have to become real again, as is the case for the eigenvalues given by Eqs.(35) as functions of $`\omega _c`$. We note that there is also an opposite process of merging of complex conjugate eigenvalues, which is also described by Eq. (B1).


text:
# A search for nearby young stars among the flare stars

## 1. Introduction

According to the recent definition given by Gershberg et al. (1999), these stars are on the lower part of the MS and show activity similar to the sun (sporadic flares, dark spots, variable emission from the chromosphere and the corona, radio, X-ray and UV outbursts).

The UV Cet-type stars are relatively young (possibly zero-age or pre-MS stars) with ages $`100`$ Myrs or younger. They may be former members of recently dispersed T or OB associations or may be ejected stars from associations by a three body encounter (like run-away T Tauri stars). A local dispersed association does not appear unlikely: The translucent high-latitude cloud MBM 12 with ongoing star formation is located at $`65`$ pc (Hearty et al. 2000) and the TW Hya association at $`55`$ pc, which has dispersed its gas and dust (Webb et al. 1999) are already known. Hence, there are indeed young T associations within $`100`$ pc. A few isolated young nearby stars are also known, e.g. GJ 182 at 27 pc and 20 Myrs.

## 2. The observations and first results

We observe all these stars with FOCES and FEROS, the high resolution echelle spectrograph of the Calar Alto Observatory in Spain and at ESO La Silla Observatory in Chile, respectively. A first reduction of the spectra of some stars shows the existence of strong Lithium absorption line at 6708 Å. As an example, HIP 71631 shows strong Lithium absorption at 6707 Åwith $`\mathrm{W}_\lambda (\mathrm{Li})=0.213`$ Å, stronger than Calcium at 6718 Å. This star has been classified as G0.Ve. The parallax given by HIPPARCOS measurement is $`29.46\pm 0.61`$ mas which corresponds to a distance of $`33.9\pm 0.7`$ pc. The V-magnitude is $`7.613\pm 0.011`$ mag. According to the HR-diagramm and the D’Antona & Mazzitelli (1994) tracks and isochrones the star has $`1.1\pm 0.1\mathrm{M}_{}`$ and an age of $`5\pm 2.410^7`$ years.

## 3. Future work

The high resolution echelle spectra enables us to study the stellar properties of the stars in detail. Especially we can measure radial velocity, temperature, metallicity, gravity and rotational velocity ($`v\mathrm{sin}i`$) to distinguish between pre-MS-dwarfs and Post-MS-giants. Stars in the sample which appear young will be proposed for VLBI observations to spatially resolve possible coronal loops and they are prime targets for direct imaging of extra-solar planets in the IR.

## References

D’Antona F. & Mazzitelli I. (1994) ApJS 90, 467

Gershberg R.E., Katsova M.M., Lovkaya M.N., Terebizh A.V., Shakhovskaya N.I. (1999), A&AS 139, 555

Hearty T., Fernández M., Alcalá J.M., Covino E., Neuhäuser R. (2000), A&A, 357, 681

Webb, R. A., Zuckerman, B., Platais, I., Patience, J., White, R. J., Schwartz, M. J., McCarthy, C. (1999) ApJ, 512, 63