text:
# Surface Brightness Measurements for APM Galaxies

## 1 Introduction

Surface brightness (SB) is one of the fundamental parameters describing a galaxy. It plays an important role in many diverse aspects of extragalactic astronomy and cosmology, from identifying the whole family of galaxy populations to modelling their different distributions and motions in the Universe. These aspects are generally related to the formation and evolution of galaxies, as well as to the nature of large scale structures in the Universe.

Unfortunately, accurate measurements of the SB of galaxies rely critically on high quality observations with perfect sky conditions and adequate telescope time, so it is not easy to get high quality SB data for galaxies, especially for low surface brightness galaxies (LSBG). Large, homogeneous data samples are even more difficult to construct, and so uniform SB data for galaxies is very scarce compared to other parameters, such as position, colour, redshift and morphology. Hitherto relatively few sets of SB measurements have been available; the total number of such galaxies is only a few thousand, and the samples have either concentrated on LSBGs or been limited only to bright galaxies (e.g. de Jong & van der Kruit 1994; de Jong 1996; Jansen et al. 2000; Heraudeau & Simien 1996; Lauberts & Valentijn 1989; Impey et al. 1996; Morshidi-Esslinger et al. 1999, hereafter MDS).

On the other hand it is not difficult to make a rough estimate of a galaxy’s SB using parameters related to both the luminosity and size of a galaxy. These two kinds of basic parameters, such as magnitude and effective radius of galaxies, are common in many catalogues. One can consequently obtain a larger sample of SB data of galaxies at the expense of some accuracy for individual objects. The statistical analysis of a large number of galaxies will not sensitively depend on the accuracy of the SB measurement for each individual galaxy because reasonable errors in SB can be allowed for in the analysis.

The APM catalogue includes both the magnitude and isophotal size for each galaxy so these can be used to estimate the SB of galaxies. Based on this strategy, some work has already used SB estimates from the APM survey. Côté et al. (1999) used the mean isophotal SB of APM galaxies to select the candidates for an HI survey of LSBGs. More recently Cross et al. (2001) calculated the effective SB of galaxies in the 2dF Galaxy Redshift Survey and estimated their bivariate brightness distribution, as well as the number and luminosity density of galaxies.

In this paper, we consider a variety of SB estimates, and present a method for determining the central surface brightness, $`\mu _0`$, of galaxies in the APM Galaxy Survey. We apply a field correction across each survey plate, and use a matching procedure over the whole survey area to yield a homogeneous SB sample of over 2 million galaxies with a wide magnitude range over a large volume of space. The results will be used in future papers to study the distribution of galaxies as a function of surface brightness, in particular their clustering properties. We describe the basic data in Section 2. The principles behind the method are described in Section 3, together with an introduction of the raw APM parameters that relate to the SB of galaxies. The field correction and matching process are briefly outlined in Section 4, where we also discuss the uncertainty of our SB estimates and compare to other SB measurements. In Section 5 we consider the distributions of $`\mu _0`$ and the SB profiles of bright APM galaxies. Finally we summarize our results in Section 6.

## 2 The APM Catalogue

The APM Galaxy Survey is hitherto the largest galaxy survey used for the investigation of the large-scale distribution of galaxies. It contains over 3 million galaxy positions and apparent magnitudes, compiled by scanning 269 UK Schmidt Telescope (UKST) survey plates over the southern sky with the Automatic Plate Measuring (APM) machine. A subset of 2.4 million galaxies in 185 UKST fields south of declination $`20^{}`$, compiled from photographic data of a consistently high quality and representing the first stage of the APM survey, was described by Maddox et al. (1990a,b). This subset covered the south galactic cap and is defined by the boundaries of UKST plates whose centres have galactic latitude between $`b=30^{}`$ and $`40^{}`$ (the exact limit depending on star number densities, and therefore on longitude), as well as declination. It is this subset of 185 plates, representing an area of 4300 square degrees, that we use in this paper.

The APM identified images using a connected-pixel algorithm, selecting groups of 16 or more pixels having intensities above a fixed threshold optical density. With a mean pixel size of 0.51 arcsec (the plate scale varies slightly from plate to plate and across individual plates), this corresponds to a minimum area of 4.3 arcsec2. The detection threshold $`t`$ was set for each plate to be about twice the rms noise in the sky, leading to a mean threshold of 12 APM density units, corresponding to $`25b_J\text{mag arcsec}^2`$ above the sky background. The threshold varies within each plate due to vignetting and emulsion effects and also varies from plate to plate, because of the different photometric zero-points for each plate. The standard deviation in the threshold due to the variation in photometric zero points from plate to plate is $`0.3b_J`$ mag arcsec-2.

The images were classified as stars, galaxies, blended images and noise based on comparisons with the median stellar profile as a function of magnitude, as described in detail by Maddox et al. (1990a). In our present analysis we consider only images classified as galaxies. The total magnitude was estimated from the isophotal magnitude and area of each detected image, assuming a Gaussian light profile. This process was optimized for faint images near the magnitude limit where galaxy profiles are significantly affected by seeing. Different plates have different limiting magnitudes, and so to ensure a uniform limit, the final sample was further constrained by the magnitude limit of the poorest plate, 20.5 bJ. Rejecting images fainter than 20.5 bJ in total magnitude produced the final APM Catalogue that we consider in this paper.

## 3 Measurement of Surface Brightness

### 3.1 APM data

The APM Catalogue includes measurements of the isophotal flux and the image area within the detection isophote for each image. These parameters are used here to estimate the SB of each galaxy assuming an exponential light profile. The APM parameters are calculated using emulsion density $`D`$, rather than incident flux $`f`$ measurements. Below the emulsion saturation $`D_{\mathrm{max}}`$, it is fairly well approximated by a linear relationship,

$$D=\{\begin{array}{ccc}\alpha f+\beta \hfill & \text{for}& f<(D_{\mathrm{max}}\beta )/\alpha \hfill \\ D_{\mathrm{max}}\hfill & \text{for}& f>(D_{\mathrm{max}}\beta )/\alpha \hfill \end{array}$$
(1)
where $`\alpha `$ and $`\beta `$ are constants for any particular part of the plate (Maddox et al. 1990a,b, Cawson et al. 1987), and the value of $`D_{\mathrm{max}}`$ is found to vary by $`2\%`$ across the field (Maddox et al. 1990b). The incident flux $`f`$, is the sum of the sky flux $`f_{\mathrm{sky}}`$, and the object flux $`f_{\mathrm{obj}}`$. So, the APM measures the sky-subtracted density $`D_{\mathrm{obj}}`$ for particular pixels of objects, as

$$\begin{array}{ccc}\hfill D_{\mathrm{obj}}& =& DD_{\mathrm{sky}}\hfill \\ & =& \alpha f+\beta (\alpha f_{\mathrm{sky}}+\beta )\hfill \\ & =& \alpha f_{\mathrm{obj}}\hfill \end{array}$$
(2)
where $`D_{\mathrm{sky}}=\alpha f_{\mathrm{sky}}+\beta `$ is the local sky value which was obtained from a pre-scan and smoothing process (Maddox, 1988). Obviously, the upper limit of $`D_{\mathrm{obj}}`$ is, $`D_{\mathrm{sat}}=D_{\mathrm{max}}D_{\mathrm{sky}}`$, for any particular position in a plate. The saturation density $`D_{\mathrm{sat}}`$ corresponds typically to a magnitude surface brightness $`22.3b_J\text{mag}\text{arcsec}^2`$ above the sky, but this figure varies with the sky background level. Note that this saturation is quite severe: for a typical disk galaxy with a central SB of $`21.4b_J\text{mag}\text{arcsec}^2`$, the profile is saturated inside 0.8 disk scale-lengths.

For each image, the APM measures the relative brightness by summing $`D_{\mathrm{obj}}`$ for pixels above the threshold $`t`$ ( $`25`$ mag arcsec-2 ),

$$I_{\mathrm{iso}}=\underset{D_{\mathrm{obj}}t}{}D_{\mathrm{obj}}.$$
(3)
The threshold $`t`$ is constant across the field of each plate, but varies from plate to plate according to the fluctuation of the sky density for individual plates (Maddox et al. 1990a,b). When working with APM data it is convenient to define an ‘isophotal APM magnitude’, $`m_{\mathrm{iso}}`$, of an image as

$$m_{\mathrm{iso}}=+\mathrm{\hspace{0.25em}2.5}\mathrm{log}_{10}I_{\mathrm{iso}}$$
(4)
This APM parameter is related to the actual isophotal magnitude, $`b_{\mathrm{iso}}`$, which is given by

$$b_{\mathrm{iso}}=Z2.5\mathrm{log}_{10}I_{\mathrm{iso}}$$
(5)
where $`Z`$ is the magnitude zero-point of the photographic plate. Also, the isophotal size of the image is given by

$$A_{\mathrm{iso}}=\underset{D_{\mathrm{obj}}t}{}1$$
(6)
$`I_{\mathrm{iso}}`$, $`t`$ and $`A_{\mathrm{iso}}`$ are the three main APM parameters that we use to measure the SB profile of galaxies. In former work on APM galaxy survey, Maddox et al. (1990b) used these parameters together with the assumption of Gaussian galaxy profiles to calculate a total APM magnitude $`m_{\mathrm{tot}}`$ which was then calibrated to $`b_J`$ using CCD magnitudes. For faint galaxies near the magnitude limit of the APM survey this is a reasonable approximation, because the galaxy images are small and the profile shape is dominated by the seeing disc. But for brighter galaxies an exponential profile is a better approximation (see Section 5). In the following subsections, we mainly use these three fundamental APM parameters to investigate the SB of galaxies.

Additionally, there are two other APM parameters that are related to the profile of galaxies. Both of them are density weighted parameters. One is $`\sigma ^2`$, defined as

$$\sigma ^2=\frac{_{D_{\mathrm{obj}}t}(x^2+y^2)D_{\mathrm{obj}}(x,y)}{I_{\mathrm{iso}}}$$
(7)
where, $`D_{\mathrm{obj}}(x,y)`$ is the density at location $`x`$ and $`y`$ relative to the image centroid. Theoretically, $`\sigma ^2`$, together with $`A_{\mathrm{iso}}`$ can constrain the SB profile of the image, and judge whether it is best fit by a Gaussian, or exponential, or even an $`r^{1/4}`$ profile. Unfortunately, in practice this is not possible, because $`\sigma ^2`$ changes by only about $`10\%`$ between these three profiles, and this is the same order as the observational error in $`\sigma ^2`$ (see Section 5). Additionally, $`\sigma ^2`$ can be seriously affected by other factors, such as the ellipticity and substructure of the image, and any other asymmetry of the profile which will tend to increase $`\sigma ^2`$. This means that the measurement of $`\sigma ^2`$ is not as robust as $`A_{\mathrm{iso}}`$ or $`I_{\mathrm{iso}}`$, and so we do not use $`\sigma ^2`$ in the measurement of SB in this paper. In Section 5.2, we consider the distribution of $`\sigma ^2`$ for bright galaxies which have a large image and hence a smaller error in $`\sigma ^2`$, and show that an exponential profile is actually a reasonable approximation to the SB profile of most galaxies.

Another useful parameter is the ellipticity of the image, which is derived from density-weighted second moments. This is can be easily converted to the apparent axial ratio of the APM image (Maddox et al. 1990a; Lambas et al. 1992),

$$e=\left(\frac{b}{a}\right)_{\mathrm{apparent}}.$$
(8)
On the assumption that the areal profiles of an APM image have similar shapes, then $`e`$ is independent of other profile parameters, and will not affect the measurement of SB. However, $`e`$ will affect the estimation of $`\sigma ^2`$ as discussed in Section 5.2.

### 3.2 Modelling the SB profile of galaxies

If the image of a galaxy is represented by a circularly-symmetric pure exponential surface brightness profile,

$$p(r)=p_0\mathrm{exp}(r/r_0).$$
(9)
The central peak value $`p_0`$, and the characteristic size $`r_0`$ fully describe the image. In this case, $`r_0`$ is the exponential scale length. The total luminosity of the image, $`I_{\mathrm{tot}}`$, is simply the integral of $`p(r)`$. The integration can be divided into three parts which contribute to $`I_{\mathrm{tot}}`$,

$$I_{\mathrm{tot}}=I_{\mathrm{sat}}+I_{\mathrm{iso}}+I_{\mathrm{field}},$$
(10)
where $`I_{\mathrm{sat}}`$, $`I_{\mathrm{iso}}`$ and $`I_{\mathrm{field}}`$ represent the flux lost due to saturation, the observed flux and the unobserved flux beyond the isophotal size of the galaxy respectively. These can be written as:

$$I_{\mathrm{sat}}=2\pi _0^{r_{\mathrm{sat}}}p(r)r𝑑r2\pi r_{\mathrm{sat}}^2D_{\mathrm{sat}};$$
(11)
$$I_{\mathrm{iso}}=2\pi _{r_{\mathrm{sat}}}^{r_{\mathrm{iso}}}p(r)r𝑑r+2\pi r_{\mathrm{sat}}^2D_{\mathrm{sat}};$$
(12)
and

$$I_{\mathrm{field}}=2\pi _{r_{\mathrm{iso}}}^{\mathrm{}}p(r)r𝑑r,$$
(13)
where

$$r_{\mathrm{iso}}=\sqrt{A_{\mathrm{iso}}/\pi }=r_0\mathrm{ln}(p_0/t)$$
(14)
denotes the radius corresponding to the isophotal area for circular images and

$$r_{\mathrm{sat}}=r_0\mathrm{ln}(p_0/D_{\mathrm{sat}})$$
(15)
denotes the radius corresponding to the saturated area for circular images. These quantities are illustrated schematically in Figure 1. Of these three parts, $`I_{\mathrm{iso}}`$ is the only directly observed quantity, and is given by $`_{D_{\mathrm{obj}}t}D_{\mathrm{obj}}`$ (eq. 3). The measured values of $`I_{\mathrm{iso}}`$ and $`A_{\mathrm{iso}}`$ can be used together with eq.(14) and eq.(15) to find the parameters $`p_0`$ and $`r_0`$. The resulting equations are not soluble analytically, but are simple to solve numerically. The actual procedure that we adopt is as follows: take an initial value of ($`p_0^{}`$, $`r_0^{}`$); calculate the ($`A_{\mathrm{iso}}^{}`$, $`I_{\mathrm{iso}}^{}`$) from eq.(12), eq.(14) and eq.(15). Then compare with the observational APM parameters ($`A_{\mathrm{iso}}`$, $`I_{\mathrm{iso}}`$), and calculate a feedback of ($`\mathrm{\Delta }p_0^{}`$, $`\mathrm{\Delta }r_0^{}`$). Several iterations will lead to the final stable result of ($`p_0`$, $`r_0`$). The analysis can be adapted to non-circular images by assuming concentric elliptical isophotes and using the APM ellipticity.

The effects of seeing can also be included in this iterative process, but for most galaxies the correction is less than 0.1 mag. Our estimates of $`p_0`$ are based on $`I_{\mathrm{iso}}`$ and $`A_{\mathrm{iso}}`$, and for images larger than about 100 pixels these measurements are not sensitive to the seeing.

### 3.3 Surface Brightness

The pair of parameters ($`p_0`$, $`r_0`$) entirely describes the exponential profile of circularly-symmetric images, so they can be used to infer the value of any other characteristic SB parameter for such a galaxy. In particular we can define the central SB in units of magnitudes per unit area,

$$\mu _0=Z2.5\mathrm{log}_{10}(p_0)$$
(16)
where $`Z`$ is the zero-point for that particular UKST plate. We can also define an effective SB,

$$\mu _\mathrm{e}=\mu _0+C$$
(17)
where the constant C depends on the exact definition of the effective SB. It can also be useful to define an effective radius, $`r_\mathrm{e}`$, which is the radius containing half the light of the galaxy. This is related to the scale size by $`r_\mathrm{e}=1.678r_0`$ for a circularly-symmetric image with an exponential light profile. For the SB at $`r_\mathrm{e}`$, the constant is $`C=1.822`$, and for the mean SB within $`r_\mathrm{e}`$, $`C=1.124`$. Also the average SB within the isophotal area (for the APM scans of UKST blue plates this isophote is $`25`$ bJ mag arcsec-2) can be written as

$$\mu _{\mathrm{iso}}=Z2.5\mathrm{log}_{10}\left(\frac{I_{\mathrm{iso}}}{A_{\mathrm{iso}}}\right)$$
(18)
where $`A_{\mathrm{iso}}`$ and $`I_{\mathrm{iso}}`$ are the corresponding isophotal size and brightness, derived directly from eq.(9). Additionally, one can also calculate the ‘total APM magnitude’ from $`I_{\mathrm{tot}}`$,

$$m_{\mathrm{tot}}=+\mathrm{\hspace{0.25em}2.5}\mathrm{log}_{10}(I_{\mathrm{tot}})=+\mathrm{\hspace{0.25em}2.5}\mathrm{log}_{10}(2\pi p_0r_0^2),$$
(19)
for circular images. This is related to the conventional total magnitude of the image, which is given by

$$b_{\mathrm{tot}}=Z2.5\mathrm{log}_{10}(I_{\mathrm{tot}})=Z2.5\mathrm{log}_{10}(2\pi p_0r_0^2)$$
(20)
where $`Z`$ is the magnitude zero-point.

The relationship between these three SB measures can be seen in Figure 2. In this Figure, we fix the value of $`r_0`$ to illustrate the dependencies between these SB parameters, giving $`r_0`$ a value of 2.15 arcsec, chosen to be fairly typical of the APM sample. Though the different SB measures are all derived from the same set of ($`p_0`$, $`r_0`$) data, we see in Figure 2a that some galaxies can show much larger differences between $`\mu _{\mathrm{iso}}`$ and either $`\mu _0`$ or $`\mu _\mathrm{e}`$ than other galaxies. Figure 2b presents the isophotal magnitude, labelled $`m_{\mathrm{iso}}`$, as a function of the total magnitude for the same value of $`r_0`$. This is shown for both the standard magnitude and for the APM magnitude parameters of equations (4) and (19). $`m_{\mathrm{iso}}`$ can be significantly fainter than $`m_{\mathrm{tot}}`$, at both high and low SB. At high SB (bright total magnitude), the flux $`I_{\mathrm{sat}}`$ is lost because of the emulsion saturation; at low SB (faint total magnitude) the flux $`I_{\mathrm{field}}`$ is missed because of the isophotal limit. This also provides a rough estimate of how much $`\mu _0`$ will be underestimated if these two parts of missing flux are not accounted for. The effect of saturation and isophotal threshold on the magnitude is shown explicitly in Figure 3(a). As mentioned in Section 3.1, bright galaxies may be saturated over a large fraction of the image. This can be seen quantitatively in Figure 3(b), which shows the ratio of isophotal radius, $`r_{\mathrm{iso}}`$, and saturation radius, $`r_{\mathrm{sat}}`$, to the scale length, $`r_0`$ as a function of total magnitude.

As described above, this approach only provides a rough measurement of the SB of galaxies. It does not use the real profile of different types of galaxies, which may be different from an exponential disk. It also neglects any internal structures, such as arms, bars or a central bulge. Also, although we have rejected images that are most likely to be merged pairs, the automated image classification is not perfect, and there is a residual contamination of about 5-10% of merged images, whose profiles will not be well represented by our simple model.

Although for any individual galaxy our SB measurement is unlikely to be very accurate, it is helpful to constrain the shape of any given profile and it does allow us to make general comparisons of galaxy SB within the whole sample. Theoretically, this approach is suitable for any other profile which is specified by only two parameters, such as a Gaussian profile,

$$p(r)=p_0\mathrm{exp}\left((r/r_0)^2\right)$$
(21)
or an $`r^{1/4}`$ law profile,

$$p(r)=p_0\mathrm{exp}\left((r/r_0)^{1/4}\right).$$
(22)
Without any other effective observational constraints on the profile shape an exponential profile is the most reasonable choice since it is a good representation for a majority of galaxies. As discussed in Section 5.2, we find that $`\sigma ^2`$ as defined in eq.(7) is not very effective in distinguishing between different profiles. Additionally, the exponential profile has been a popular choice in earlier work, so it allows us to compare with other results. So, in this paper, we continue our analysis and discussions under the assumption that galaxies have an exponential SB profile.

Note that these estimates are in terms of the “APM” SB, which is based on APM density measurements. They are directly related to true SB only if there is a linear relationship between the real flux and the APM density (eq.1). Clearly the emulsion saturation violates this assumption, and CCD observations show that it is not exactly linear below saturation (Cawson et al. 1987). Our fitting technique approximates the saturation as a simple maximum APM density, and ignores any non-linearity below the saturation density. This non-linearity will not strongly influence the ordering of galaxies as a function of their SB, but more accurate measurements would have to allow for this factor. Furthermore, if the intrinsic SB is required for a particular analysis instead of the apparent SB, then the K-correction, cosmological dimming, and the correction of the apparent ellipticity of disk galaxies should also be taken into account.

## 4 Field correction and overall plate matching

To provide a uniform sample of SB measurements over the whole APM galaxy survey area, we make two further systematic corrections. First, we correct for the field response function to make the different parts of each UKST plate consistent from the central region to the edge. Then, we make an overall photometric zero-point adjustment to each plate, to give a uniform sample over the whole survey area.

### 4.1 Field response corrections

The variation in the response function across the field of each plate is caused by two main factors. One is the geometrical vignetting of the optical system of the telescope (UKST), that reduces the effective transmission off its axis. A theoretical analysis (Dawe & Metcalfe 1982) shows the transmission is very flat for images closer than $`2.5^{}`$ from the center of a plate, but quickly falls to $`80\%`$ of the central value at $`4^{}`$ off axis. This reduction is same for each plate. Another loss in sensitivity is caused by the differences in the chemical desensitization at different parts of the UKST survey plates. For most of the survey plates, the original plate holder has an air gap between the plate and the filter, and the greater air circulation near the plate edges caused greater desensitization of the hypersensitized III-aJ emulsion. This effect causes a multiplicative reduction in the measured flux, which is quite similar to the geometrical vignetting, but the pattern actually varies for different plates, because the local weather conditions during individual observations were significantly different.

Both of these effects can be considered as a change in the relation between APM densities and true flux. If two objects have same flux, $`f_{\mathrm{obj}}`$, but are located at different positions on a plate, then the geometrical vignetting and the different desensitization lead to different observed APM densities $`D_{\mathrm{obj}}`$ via eq.(2). The constant $`\beta `$ is essentially the fog level of the emulsion, which is constant, so the loss in sensitivity means that the slope $`\alpha `$ is not a constant across the plate. It is a function of position in the field, being reduced by the geometrical vignetting and emulsion desensitization from centre to the edge. Thus, both these effects can be considered as variations in $`\alpha `$.

Fortunately, $`\alpha `$ is only a function of position in the field, and does not depend on the specific properties of an image. Therefore, the correction process of Maddox et al. (1990b), which was devised to correct the total magnitude of images, can be simply applied here to correct the SB of galaxies. According to eq. (8) of Maddox et al. (1990b), the correction to the total APM magnitude is

$$\mathrm{\Delta }m_1=2.5\mathrm{log}_{10}\left(\frac{\alpha _1}{\alpha _0}\right)$$
(23)
where the subscript 0 denotes the centre of a plate and the subscript 1 denotes the particular position of an object. According to the relationship between the central SB and the total APM magnitude of eqs. (16) and (19), we can see that $`\mu _0`$ will need the same correction as $`m_{\mathrm{tot}}`$:

$$\mathrm{\Delta }\mu _1=2.5\mathrm{log}_{10}\left(\frac{\alpha _1}{\alpha _0}\right).$$
(24)
The value of $`\mathrm{\Delta }\mu _1`$ is determined individually for each galaxy by estimating the APM sky density in the surrounding 10 arc minutes. The apparent variation in sky density is converted to an effective change in $`\alpha `$ using an empirical relation which was chosen to ensure that the number density of galaxies is uniform when averaged over many UKS plates. This correction will be generally applicable for other SB measurements, except for the mean isophotal surface brightness, $`\mu _{\mathrm{iso}}`$. Since the threshold $`t`$ is fixed in terms of observed APM density, the actual threshold isophote varies over each field and a more complex correction would be required to ensure that $`\mu _{\mathrm{iso}}`$ is uniform.

Most survey galaxies are located within the central $`3^{}`$ radius of a UKST plate, and so the correction is very small, typically $`0.1`$ magnitude. For galaxies outside this region the correction can increases up to 0.4 mag. See figure 4 for the pattern of averaged $`\mathrm{}\mu _1`$ of field correction. See the figures of Maddox et. al (1990b) for further details.

Figures 5(a) to (b) show that the field correction slightly reduces the rms difference between the SB measurements for those galaxies imaged on two independent plates. F409 and F471 are a pair of contiguous plates with quite a large overlapping area, and there is a total of 7960 identified galaxies with measurements from both the plates. For quite a few of these galaxies, the distance from the centre of one of the plates can be large, which means that uncorrected values of $`\mathrm{\Delta }\mu _0`$ can be quite large. Figure 5 (a) plots the $`\mu _0`$ before field correction and (b) after it. The rms scatter in these two plots is reduced from 0.27 to 0.25 mag. This improvement suggests that the rms $`\mathrm{\Delta }\mu _0`$ from the field correction is $`0.1`$.

### 4.2 Matching calibrations in plate overlaps over the whole survey

Variations in the emulsion sensitivity and observing conditions mean that each UKST survey plate has a different photometric zero-point. To produce a uniform set of galaxy SB estimates, we need to adjust the zero-point of each plate to force a match between galaxies in overlapping regions. In the APM survey paper II, Maddox et al. (1990b) described a scheme to match the APM total magnitudes over the survey. In this paper we apply a similar procedure to match our SB estimates.

The matching procedure uses the fact that each plate covers an area of sky that overlaps with its neighbours, so there are many galaxies that have independent measurements from two neighbouring plates. For a particular overlap the mean difference between the two measurements gives an estimate of the difference between the zero-points of the two plates. Taking all overlaps together we can find the set of plate zero-points which is most consistent with these plate-to-plate offsets.

The main difference for this paper is that a galaxy SB profile is determined by two parameters ($`p_0`$, $`r_0`$) as described in Section 3, and so both parameters must be matched. In principle, we could choose to match any two independent parameters derived from these two, but for simplicity we choose $`\mu _0`$ and $`b_{\mathrm{tot}}`$. They are related to $`p_0`$ and $`r_0`$ through eqs. (16) and (20), while the zero-point in eq. (16) represents the overall corrected uniform zero-point. Since we have assumed exponential profiles in this analysis, the magnitude zero-points are not the same as the Gaussian total magnitudes that were applicable to galaxies near the plate limit, as used in previous work on the APM galaxy survey.

Images were matched in the overlaps of every pair of plates with centres closer than $`7^{}.`$ For the 269 plates in the south galactic cap, this gives 809 distinct overlaps with an average of 1405 galaxies brighter than $`b_J=20.5`$ in each overlap. To allow for non-linearities which introduce a non-unit slope in the APM magnitude and SB measurements we used a linear fit to $`(\mu _i\mu _j)`$ vs $`(\mu _i+\mu _j)/2`$. Before matching, the variance in zero-point for the different UKST survey plates is very large, as can be seen in the offset of $`\mu _0`$ for galaxies in the overlapping plate pairs. The distribution of offsets $`\overline{\mu }_1\overline{\mu }_2`$ for all 809 overlaps, is shown in Figure 6, and the rms scatter is 0.43 mag. After matching, this value is reduced to 0.08 mag. As shown in Figure 7, similar reduction in scatter is seen for $`m_{\mathrm{tot}}`$, where the rms before matching is 0.31 mag, and after matching is 0.05 mag.

For different plates, the zero-point corrections for $`\mu _0`$ and $`m_{\mathrm{tot}}`$ are well correlated, but not exactly the same. The scatter about the correlation is caused by the different seeing on different plates, which leads to a different apparent $`r_0`$. The matching process has effectively given $`r_0`$ for galaxies on each individual plate a systematic correction to a common “zero-point”.

### 4.3 Error estimation

For an individual galaxy, the deviation of the actual SB profile from an exponential profile is likely to be the main source of error in the SB estimate, but this effect is beyond the scope of this paper. If we assume a given profile, then we can discuss the “intrinsic” error in the APM SB caused by the observational errors in ($`I_{\mathrm{iso}}`$, $`A_{\mathrm{iso}}`$).

The errors in ($`I_{\mathrm{iso}}`$, $`A_{\mathrm{iso}}`$) come from a complex combination of sources which we consider in two main categories. First is the observing process, which includes a mixture of weather conditions and instrumental effects, including the residual calibration differences between the survey plates that are not corrected by the field correction function and overall plate matching. This category also includes photon noise. Second are the errors from the APM scan itself. These errors can be estimated directly from repeat scans of selected plates. For $`I_{\mathrm{iso}}`$ we find this source contributes only about 0.04 mag uncertainty, and for $`A_{\mathrm{iso}}`$ we find about 5%. These errors could be analytically transformed to $`\mu _0`$ following a standard propagation of errors procedure, but even if the resulting uncertainty in $`\mu _0`$ is slightly larger, it will be negligible compared to the observational errors.

The overlapping areas of neighbouring plates provide a direct opportunity to estimate the combined errors in our SB estimates. Figure 5(c) shows the residual errors for the overlapping field of F409 and F471, after applying the field correction and the overall matching correction, and it can be seen that there is still an rms scatter of about 0.25 mag between the two independent measurements for each galaxy. This implies the error in each $`\mu _0`$ will be $`0.18`$ mag. We have calculated the rms scatter for each one of the 809 plate overlaps, independently, and show the frequency distribution of these rms values in Figure 8. For most of these overlaps, the rms of $`\mathrm{\Delta }\mu _0`$ is less than 0.3 mag, and the mode is at 0.21 mag, corresponding to a typical uncertainty of 0.15 mag for an individual SB measurement.

If we consider all of the plate overlaps as a whole, the rms $`\mathrm{\Delta }\mu _0`$ over all 1,282,600 galaxy pairs is 0.23 mag, so the overall uncertainty for an individual SB measurement has an rms of 0.16. This value includes both the random error on each SB measurement, and the rms residual difference between the zero-points for all plates after applying the matching corrections. The difference between this rms, and the rms within each overlap suggests that the residual error in the zero-point matching has an rms of 0.07 mag per plate, consistent with the overlap residuals.

Additionally, we have considered the error as a function of $`\mu _0`$. Table 1 lists the rms difference between the $`\mu _0`$ measurements in overlaps for a series of $`\mu _0`$ ranges. Galaxies with $`\mu _0`$ about 21.0 $``$ 22.0 mag have the smallest error at 0.16 mag. At both low and high surface brightness the error increases up to $``$ 0.23. This tendency could also be seen in Figure 5(c). The reason for the increased error is that $`\mu _0`$ depends on the profile of a galaxy (see Section 5.2), and when we estimate the total intensity of the galaxy, $`I_{\mathrm{tot}}`$ (eqs. (10) to (13)), only $`I_{\mathrm{iso}}`$ is directly constrained by the observational data, while $`I_{\mathrm{field}}`$ and $`I_{\mathrm{sat}}`$, the other two contributions to $`I_{\mathrm{tot}}`$, involve extrapolations that are more extreme at low or high SB.

So in summary, the overall “intrinsic” uncertainty in $`\mu _0`$ including measurement errors, and zero-point errors is typically 0.2 mag. Comparing galaxies measured on one plate, the relative uncertainty is about 0.16 (also see the next subsection).

### 4.4 Galaxy surface density and mean $`\mu _o`$

After applying the corrections to the individual plate data, we should be able to define galaxy catalogues selected with uniform SB criteria. We can make a simple test of the uniformity of the sample by calculating the mean SB for galaxies on each plate. Cosmic structure will introduce some intrinsic variation in the mean SB in different directions on the sky, but this is likely to be small; indeed the variation in galaxy numbers from plate to plate at $`b_J=20`$ is only $`12\%`$ (Maddox, Efstathiou & Sutherland 1996). Figure 9 shows a plot of the mean $`\mu _0`$ for each plate in the survey, before and after applying the matching corrections. Using all fields, the rms scatter before correction is 0.25 mag, and after correction it reduces to 0.095 mag. However, it is clear that fields centred with declination $`>20^{}`$ have a larger scatter, even after correction. This is partly due to the fact that the scanned plates with declination $`>20^{}`$ are a rather heterogeneous mix of original plates and non-survey plate copies, and partly due to the significantly large galactic extinction in this part of sky (see Schlegel, Finkbeiner & Davis 1998). In principle it should be possible to correct these problems, but for simplicity we have simply excluded the fields with declination $`>20^{}`$ in the subsequent analysis.

A further simple check on the uniformity is to compare the average number of galaxies per plate for SB selected samples. Figure 10 shows the distribution of galaxies on the sky for three subsamples selected on SB from high to low SB respectively. The ranges in $`\mu _0`$ were chosen to give samples with roughly the same number of galaxies in each subsample. Note that the high and low SB samples are selected from the tails of the SB distribution (see Section 5.1) and so are most sensitive to errors in the plate matching.

For the 185 plates there are no plates which show a clearly discrepant number of galaxies, suggesting that our matching procedure does indeed produce a uniform set of SB measurements. This improvement can also be seen in figure 11, which plots the numbers of galaxies in each plate for these three subsamples before and after matching correction. After matching, the fractional rms scatters in the number of galaxies on each plate are 0.27, 0.13, 0.17 respectively for the three samples. These are close to what is expected for randomly selected $`5^{}\times 5^{}`$ patches of sky.

The higher rms for the high SB sample reflects the higher clustering amplitude of high SB galaxies, which is clearly apparent in Figure 10. We plan to quantify these differences in in clustering in a future paper.

### 4.5 Checking the APM results and photometric calibration with deeper CCD data

It is important to test the accuracy of the APM surface photometry, particularly the photometric calibration, through a comparison with well-calibrated higher-resolution CCD data. To this end, new CCD surface photometry has been performed on faint ($`b_J18`$ to 22 mag) galaxies in the region around the Hubble Deep Field South (HDF-S), which lies in the area of the available APM data (at galactic coordinates ($`l`$,$`b`$) = ($`328^{}`$, $`49^{}`$)). Much well-calibrated public imaging data are available around the HDF-S, and ground-based HDF-S CCD observations cover a sufficient area to provide the numbers of galaxies required for these tests. The Goddard Space Flight Center HDF-S public data (Teplitz et al. 2001; Palunas et al. 2000) were selected for the task because they combine a large image area, standard broadband photometric passbands and good photometric calibrations. The Goddard HDF-S data were recorded with the Big Throughput Camera on the 4-metre Blanco telescope at the Cerro-Tololo Inter-American Observatory (Wittman et al. 1998) and cover a 47 arcmin by 46 arcmin area. An independent measurement found that the seeing of the blue data was 1.7 arcsec in the central 15 arcmin by 15 arcmin region; although not good, this is comparable to or better than the photographic APM data. Pixels were 0.43 arcsec in width.

We measured the SB profiles of galaxies in this region through the following steps. Firstly, images were detected in the B-band CCD data using the SExtractor package (Bertin 1998; Bertin & Arnouts 1996). This provided image centroid positions, ellipticities and orientations. For each galaxy image, mean surface brightnesses were computed in a series of concentric elliptical annuli defined by these SExtractor image parameters. Then the surface brightness was plotted against semi-major axis for each image and an exponential light profile fitted, defining the range in semi-major axis for these fits by visual inspection. Care was taken to avoid fitting to the core of the image where seeing effects were significant. This provided a central surface brightness and exponential scale length for each selected galaxy. Finally, the surface brightness results were transformed from the instrumental blue magnitude system to the photographic $`b_J`$ system using the Goddard HDF-S transformation and the standard B to $`b_J`$ equation of Blair & Gilmore (1982). This used colours measured within identical elliptical apertures in the B and V data. The overall uncertainty in the individual central surface brightness measurements is about 0.1 mag.

An independent check of these results was made using other public ground-based HDF-S data. SExtractor total magnitudes were computed using the Anglo-Australian Observatory public CCD observations for the region immediately around the WFPC2 field. These were transformed to the $`b_J`$ band and compared with the total magnitudes under the fitted exponential profiles derived from the Goddard HDF-S data. Six galaxies with Goddard photometry were identified in the AAO field. The difference in the total magnitudes was $`+0.032\pm 0.030`$ mag, confirming the validity of the results obtained from the Goddard public data. The surface photometry from the Goddard data can therefore be used to check the APM photographic results.

The HDF-S CCD surface photometry gave central surface brightnesses $`\mu _0`$ and scale lengths $`r_0`$ for a large sample of galaxies with $`b_J=17.5`$ to 22 for which the exponential profile was judged to be an adequate fit. Of this sample, 24 objects appeared in the APM Galaxy Catalogue. Of these 24, two galaxies were rejected from further consideration because their angular sizes were so small (isophotal areas smaller than 100 pixels) that seeing was likely to affect appreciably the results from the APM data. This left 22 galaxies in the region around the HDF-S that had surface photometry both from the APM data and from the CCD observations. These galaxies were used to test the reliability of the surface brightness results derived from the APM data.

Figure 12 shows the relationship between the CCD $`\mu _0`$ and our APM measurements. The rms scatter is 0.29 mag arcsec-2, and the difference in their zero-points is 0.26 mag arcsec-2. This dispersion is larger than expected from simply adding the uncertainty of 0.18 in the APM measurements in quadrature with 0.1 from the CCD measurements. The extra scatter is probably due to the fact that we have not corrected the non-linear correspondence between flux and the emulsion density of the survey plates. Note also that nearly all of the galaxies in this sample have $`\mu _0<22`$ mag arcsec-2, which is brighter than the saturation level of the UKST emulsion and thus the APM measurements depend on extrapolated galaxy profiles.

The scatter could be reduced by restricting the range in radius used to fit the CCD data so that it matches more closely the range available in the APM data. However, this would not be a true representation of the uncertainties in the APM SB estimates. We want to compare the APM results with results from CCD surface photometry that have been obtained using standard surface photometry methods. Note that there will also be systematic biases correlated with morphological type. Since we have assumed an exponential profile, a galaxy with an $`r^{1/4}`$ profile will not be well represented by our model (see Section 5.2).

In summary, although the APM measurements have relatively large uncertainties, they are adequate to classify galaxies into broad surface brightness categories.

### 4.6 Comparison with other SB measurements of UKST survey plates

Morshidi-Esslinger et al. (1999, hereafter MDS) made a survey of a substantial number of UKST survey plates in the south galactic cap in an attempt to identify LSBGs with large angular size. They used the APM to measure profiles of nearly 3000 LSBGs and estimated the central surface brightness by extrapolating the best fit exponential profile. Their SB estimates are based on the same original plates as our data, and they use the same scanning machine. They also make the same assumption of exponential galaxy profiles, so the MDS survey is not a completely independent data set. Nevertheless it provides a useful comparison for our estimates.

There are two main differences between our measurements and those of MDS. First, they fitted exponential profiles to up to 8 isophotal areas $`A_0,A_1,\mathrm{},A_7`$ – see the details in Maddox et al. (1990a) and MDS (1999) – to solve for the ($`p_0`$, $`r_0`$) parameters, whereas we used only $`A_0`$ to calculate these values. The second difference is in the calibration procedure. MDS’s calibration of SB depends on zero-points for total magnitudes from other works for fields where this is available, and uses the magnitudes of stars in plate overlaps for fields with no external calibration. Our calibration depends mainly on matching $`\mu _0`$ directly for the galaxies plate overlaps, with an overall calibration from external photometry.

Figure 13 compares our estimates of $`\mu _0`$ with 2185 MDS galaxies that can be identified in our sample. The unidentified galaxies are mainly classified as ‘merged’ images in the APM galaxy survey, and thus not included in our catalogue. There is a very good correlation between the two estimates, but the scatter is rather large, with an rms of 0.35 mag. To investigate the cause of this dispersion, we have examined data from three fields (F157, F358 and F482) having more than 65 MDS galaxies in each, and have plotted them with different symbols in Figure 14. It is clear that the rms scatter for each plate is very small, only 0.16, 0.15 and 0.14 respectively. These scatters are slightly smaller than the typical rms uncertainties estimated from the APM plate overlaps. So, even if we assume that the MDS measurements are correct, and all of the scatter is due to the approximation involved in our simple method of estimating $`\mu _0`$, the errors introduced are smaller than the observational errors associated with measuring different UKST plates.

These three plates were deliberately chosen to lie at the extremes of the distribution plotted in Figure 13, and it is clear that the zero-points of these three plates are very different. The zero-point for field F358, is almost 1 mag away from other two plates. If our zero-point was in error by 1 mag, it would be clearly visible as a discrepant point in Figure 9. Also the area covered by field F358 would show up as a large overdensity in Figure 10(a), and an underdensity in Figure 10(c). None of these tests suggests there is any significant error in our zero-point for this plate, and so we believe the MDS zero-point is in error. Assuming that the rms error on our measurements is 0.16 mag as discussed in Section 4.3, this suggests that the MDS zero-points have an rms of 0.31 mag, and so F358 would be a 3 sigma deviation, which is not unexpected in a sample of around 100 plates.

MDS estimate that the uncertainty in their calibration is about 0.1 mag for plates with external calibration, and 0.2 for plates where they have used overlaps. There are two possible sources for the larger errors that we find. First, they do not apply any sort of vignetting or desensitization correction. Since their zero-points rely on edge matching and small-area CCD calibration data, this will introduce an error in the plate zero-points with an rms of order 0.1-0.2 mag. Second, MDS use the magnitudes of faint stars in the plate overlaps to estimate the zero-points on plates with no external calibration. Even faint stars tend to have a higher SB than galaxies, and this is likely to introduce a further bias. These two effects together might account for the observed 0.31 mag rms difference. This comparison shows that calibrating the zero-points of photographic data is a difficult problem that is still without a satisfactory solution.

Impey et al. (1996, hereafter ISIB) scanned 24 UKST survey plates to search for LSBGs with $`\mu _0>22`$ mag arcsec-2. Some of these plates are covered by the APM data used in this paper, and we find 202 galaxies in common with their data. Plots comparing the our magnitudes and surface brightnesses with the ISIB measurements are shown in Figure 15. There is a good correlation between the ISIB estimates and ours, but there are several differences in their analysis and definition of parameters which preclude a quantitative comparison with their data.

## 5 Discussion

### 5.1 The Distribution of Surface Brightnesses

Figure 16 shows a histogram of the frequency distribution of $`\mu _0`$. The distribution is centred on $`\mu _0=22.2`$ mag arcsec-2, and extends to $`\mu _020`$ mag arcsec-2 at the high SB end, and to $`\mu _024.5`$ mag arcsec-2 at the low SB end. The star-galaxy separation criteria (Maddox et al. 1990a) mean that some very high SB images will be rejected from the galaxy sample because their profile is so similar to stars. However, the rapid drop in the apparent number of high SB galaxies will be due to their small isophotal areas and a real drop in the space density of high SB galaxies (as shown by Cross et al. 2001). The lack of very low SB galaxies is less clearly a reflection of the true distribution of galaxy SBs, since a galaxy must have at least 4 square arcseconds brighter than 25 mag arcsec-2 to be selected in the APM galaxy sample. The distribution is further constrained by the 20.5 mag limit on the Gaussian profile magnitudes which will also bias against low SB galaxies. Cross et al. (2001) have used the 2dF Galaxy Redshift Survey to estimate the distribution of intrinsic surface brightnesses, and conclude that, although the APM selection criteria do exclude low SB galaxies, the missing galaxies would not contribute significantly to the overall luminosity density. We will consider how our selection criteria affect the observed distribution in a future paper.

It is clear from Figure 16 that the $`\mu _0`$ subsamples plotted in Figure 10(a) and 10(c) represent the galaxies at the extremes of SB in the APM sample. This means that the distribution of galaxies in the sub-samples is very sensitive to any residual errors in the uniformity of the $`\mu _0`$ measurements, and so the apparent uniformity of the galaxies in Figure 10(a) and 10(c) show that our corrections have given us a reliable and uniform set of $`\mu _0`$ measurements.

### 5.2 The profiles of APM bright galaxies

In this Section, we use galaxies from the APM Bright Galaxy Catalogue (Loveday 1996; hereafter APMBGC) to investigate how the profile of a galaxy depends on its morphological type, and test if the assumption of an exponential profile is valid. The galaxies in the APMBGC have $`b_J<16.44`$ and have a relatively large angular size ($`r_{\mathrm{iso}}`$) and large $`\sigma ^2`$. Although $`\sigma ^2`$ is not very sensitive to changes in the profile of galaxies, and the observational error is quite large, this parameter does contain useful information about the SB profile. Furthermore, the APMBGC galaxies have been morphologically classified by visual inspection of the UKS plates and we have analyzed the distribution of $`\sigma ^2`$ separately for both early and late type galaxies.

We calculated the ($`p_0`$, $`r_0`$) of each galaxy, assuming an exponential SB profile as in Section 3. Since we also know the apparent ellipticity $`e`$ for each galaxy (from the density weighted second moments), and the values of the detection threshold and saturation of the plate, we can simulate the observed APM density for each pixel. Any other observational parameters can then be estimated from this exponential model image. In particular, we have estimated the value of $`\sigma ^2`$, and compare this to the actual observed value. Obviously, if the original image is exactly an exponential profile, and there are no other observational errors, then the observed value, $`\sigma _{\mathrm{obs}}`$ will be exactly equal to that from the exponential model $`\sigma _{\mathrm{exp}}`$, so the ratio $`\sigma _{\mathrm{obs}}/\sigma _{\mathrm{exp}}`$ should be exactly unity. To test the sensitivity of $`\sigma `$ to variations in the SB profile, we created simulated images that follow the $`r^{1/4}`$ law profile (eq. 22), and analysed them using the assumption of an exponential profile. In this case we found that $`\sigma _{\mathrm{obs}}/\sigma _{\mathrm{exp}}0.9`$, with a small variation according to their real central surface brightness and apparent ellipticity. The same analysis for Gaussian profiles gave a ratio of $`\sigma _{\mathrm{obs}}/\sigma _{\mathrm{exp}}1.1`$.

In Figure 17, we plot the frequency distribution of $`\sigma _{\mathrm{obs}}/\sigma _{\mathrm{exp}}`$ for early and late galaxies from the APMBGC. There is a large dispersion in $`\sigma _{\mathrm{obs}}/\sigma _{\mathrm{exp}}`$ , even for only late type galaxies. Nevertheless, the peak of the distribution is very close to 1.0, which implies that the exponential profile is a reasonable approximation for most galaxies (77% of the sample has $`0.9<\sigma _{\mathrm{obs}}/\sigma _{\mathrm{exp}}<1.1`$). For early type galaxies (E and S0), the dispersion is larger than for late type galaxies, and most of them have smaller $`\sigma _{\mathrm{obs}}/\sigma _{\mathrm{exp}}`$ values that are more consistent with an $`r^{1/4}`$ law profile. This shows that, although the APM measurements span a rather limited range of SB, the $`\sigma `$ parameter can distinguish between early and late-type galaxies.

For fainter APM galaxies, the observational uncertainty on $`\sigma ^2`$ will be much larger, and the small differences between the $`r^{1/4}`$, exponential and Gaussian assumptions mean that $`\sigma ^2`$ is not effective in constraining the SB profile.

## 6 Summary

We have described a simple method to estimate the central surface brightness, $`\mu _0`$, of galaxies, allowing for saturation and isophotal effects, and applied it to all the galaxies in the APM galaxy survey. We have corrected the $`\mu _0`$ estimates for the effects of vignetting and differential desensitization over the field of each UKS plate. We have also used the plate overlap areas to match the individual plate measurements, and so provide homogeneous $`\mu _0`$ estimates for all the APM galaxies. Internal and external checks of the uncertainties in $`\mu _0`$ suggest that the uncertainty for an individual measurement is about 0.16 mag.

For very bright APM galaxies, the second radial moment of the SB profile, $`\sigma ^2`$, shows that late type galaxies are consistent with an exponential profile, whereas early type galaxies have a higher mean SB, consistent with an $`r^{1/4}`$ law.

This uniform set of SB estimates for the APM galaxies will be used in later papers to select samples of high and low SB galaxies and analyze differences between their distributions.

## Acknowledgements

Z.S. would like to thank all members of the astronomy group in Nottingham for their hospitality during his stay in Nottingham. Z.S. acknowledges support from grant NKBRSF G19990754 and the National Nature Sciences Foundation of China. This work was also partly supported by PPARC. We thank the referee for useful comments on the original manuscript. The comparison with CCD observations used the NASA Goddard Space Flight Center HDF-S public data recorded by J.P. Gardner and P. Palunas at the Cerro Tololo Inter-American Observatory, a division of the National Optical Astronomy Observatories. Further comparisons used the Anglo-Australian Observatory Hubble Deep Field-South public data, recorded by K. Glazebrook, R. Abraham and C. Tinney on the Anglo-Australian Telescope.


format: text

text:
# Negative Specific Heat of a Magnetically Self-Confined Plasma Torus

## Abstract

It is shown that the thermodynamic maximum entropy principle predicts negative specific heat for a stationary magnetically self-confined current-carrying plasma torus. Implications for the magnetic self-confinement of fusion plasma are considered.

CLASSIFICATION:
PHYSICAL SCIENCES
Applied Mathematics

©2002 The authors. Reproduction of this article, in its entirety, for non-commercial purposes is permitted.

The goal of the controlled thermonuclear fusion program is to make the energy source that powers our Sun available to human society. Deep in the Sun’s interior, favorable conditions for the quasi-stationary nuclear burning of the solar plasma prevail as a result of the immense gravitational self-forces that keep this huge accumulation of matter together. Since gravitational self-confinement is not operative at the reactor and laboratory scale, alternate means of confinement have to be employed to achieve sufficiently high plasma densities and temperatures in a reactor. In the perhaps most prominent stationary fusion reactor scheme, the tokamak, strong electric ring currents are induced in an electrically neutral plasma to achieve axisymmetric magnetic self-confinement in a rotationally invariant toroidal vessel $`𝒯`$. In a torus with sufficiently large major axis, such a magnetic self-confinement mimics gravitational self-confinement on account of the Biot-Savart law, according to which in a system of parallel current filaments all filaments attract each other magnetically with the same force law as would be the case gravitationally in a system of parallel mass filaments. What makes the magnetic forces more attractive (in the double sense of this phrase) than gravity for laboratory purposes is their very much bigger coupling constant ($`10^{40}v^2/c^2`$ for two electrons moving on parallel trajectories with speed $`v`$ as measured in the laboratory; note that the even stronger ($`v^2/c^2`$ is replaced by 1) electrostatic repulsive forces between the same two electrons are very effectively screened in a neutral plasma and are traditionally neglected to a good approximation). Of course, the analogy does not extend to all aspects of plasma self-confinement. In particular, the solenoidal-vectorial character of stationary current densities necessitates the toroidal topology of magnetic self-confinement, whereas gravity not only allows but manifestly prefers spherical confinement over toroidal. Unfortunately, axisymmetric toroidal magnetic self-confinement is not known for its stability either. Although major efforts are devoted to the stabilization of the plasma configuration, a vast reservoir of instabilities capable of destroying the confinement has dramatically slowed down the development of a operating tokamak fusion reactor.

Matters are not exactly helped by the fact that our theoretical understanding of the physics on the various space-time scales that govern magnetic plasma confinement is still quite incomplete. In particular, while the solenoidal character of the magnetic induction together with the axisymmetry and stationarity of the law of momentum balance tell us that the poloidal magnetic flux function $`\mathrm{\Psi }`$ and the toroidal current density $`j`$ must satisfy some local functional relation $`(\mathrm{\Psi },j,r)=0`$ (in which $`r`$ is the cylindrical distance from the axis of symmetry) which turns Ampere’s law into some in general nonlinear elliptic partial differential equation for $`\mathrm{\Psi }`$ (known in the fusion literature as the Grad–Shafranov equation), the actual relation $``$ is not fixed (except for the explicit $`r`$-dependence in $``$). Some information about $``$ should be contained in the law of energy balance between current drive (through an applied toroidal emf and other means), ohmic heating and, ultimately, radiation losses. Unfortunately only the so-called classical and neo-classical transport coefficients have been computed in some detail  whereas the small scale turbulent dissipation mechanisms in a tokamak plasma remain a largely challenging open problem. In this situation, theoreticians have been forced to rely on fair judgment and good taste when guessing some additional principle(s) that would effectively complete the characterization of the stationary magnetically self-confined plasma torus in a tokamak.

A sizeable fraction of the literature employs a linear approximation of $``$ to get $`j\mathrm{\Psi }`$, rendering the Grad–Shafranov equation linear, and some improved-accuracy modeling uses a third-order polynomial approximation of $``$ . Subsequently the filter of linearized dynamical stability analysis, based mostly on macroscopic magnetofluid theory and mesoscopic kinetic theory, is applied to sort out unstable configurations. While this approach has met with a certain limited success, one does not learn what the approximations are approximate to. Over the years a number of plasma theorists  have argued that an equilibrium thermodynamics-inspired maximum entropy principle with a few global dynamical constraints111Since *any* probability density maximizes the entropy relative to itself, a stationary plasma torus is necessarily a maximum entropy configuration under *some* constraints. What makes the maximum entropy proposal non-empty is the insistence on only a few global natural dynamical constraints. should give answers close to the truth. In essence the various formulations in  give for $``$ the answer $`j\mathrm{exp}(\mathrm{\Psi })`$, which leads to a nonlinear Grad–Shafranov equation that may have more than one solution $`\mathrm{\Psi }`$, depending on the domain $`𝒯`$, the boundary conditions for $`\mathrm{\Psi }`$, and the values of the physical parameters of the problem. In addition, this approach provides a global stability criterion within the class of axisymmetric states satisfying the same dynamical constraints. Only those solutions that maximize the relevant relative entropy functional will be globally stable.

Since in non-equilibrium statistical mechanics the maximum entropy principle has not acquired a status anywhere near as fundamental as in equilibrium statistical mechanics, it is mandatory to register some arguments in its favor for the case at hand. Thus, the relation $`j\mathrm{exp}(\mathrm{\Psi })`$ has been shown to be almost universally singled out also by a truly dissipative Fokker–Planck approach to stationary magnetically confined plasma . The perhaps most compelling reason to give it serious considerations, however, is the successful application of the maximum entropy approach to the physically distinct but mathematically quite similar problems of stationary planar incompressible flows,222The first qualitative predictions based on statistical mechanics of the Hamiltonian system of $`N`$ point vortices were made in . The quantitative evaluation began with ; impressive agreement with simulated flows is reported in . Its mathematical rigorous foundations are by now almost complete, the latest word being ; see  for a review. More recently a formulation based directly on continuum vorticity has gained much ground; see  for a state-of-the-arts report. where the vorticity plays a rôle closely analogous to the current density in the plasma torus, and the strongly magnetized pure electron plasma in a circular cylinder,333In the guiding center approximation the dynamics of this plasma system is identical to that of $`N`$ point vortices . Statistical mechanics in the corotating frame predicts that at high enough effective energies the nonlinear $`m=1`$ diocotron mode has higher entropy than any other configuration with the same energy and angular momentum , in accordance with remarkable real experiments . where the charge density plays that rôle. In this spirit, we have conducted a thorough investigation of the thermodynamic-type maximum entropy approach to the magnetically self-confined stationary plasma torus .

A most curious finding of the study  is that the gravity-inspired toroidal magnetic plasma self-confinement scheme inherits from the stars their gravo-thermal *negative specific heat*.444Ref. , pp. 60-63, explains why a homogeneous piece of “everyday matter” must have positive specific heat. See p. 62 of the same reference for why those arguments do not rule out negative specific heat in an isolated gravitating system. Indeed, the virial and the equipartition theorems imply that in a spherical equilibrium system the energy is distributed -2:1 between gravitational and kinetic. A decrease in total energy $`E`$ of a gravitational equilibrium gas ball will increase its thermal energy. Such a system grows hotter while losing energy through, say radiation. Negative specific heat in self-gravitational perfect gases is evaluated quantitatively already in  and is further discussed in ; however, none of these configurations with negative specific heat is thermodynamically stable though some are metastable. Thermodynamically stable self-gravitating configurations with negative specific heat can occur when the Newtonian $`1/r`$ singularity is stabilized either as in quantum mechanics  or in classical hard balls systems . This result is a little surprising, for it follows from what we said earlier that the plasma torus should actually more closely mimic a cylindrical caricature of a star, and the specific heat of a ‘cylindrical maximum entropy star’  and its plasma physical clone, the cylindrical Bennett pinch , is non-negative! The existence of a maximum entropy plasma torus with negative specific heat is therefore a truly nontrivial fact. The purpose of this note is to point out some potentially important consequences of our finding for plasma physical applications.555Recently, the existence and importance of negative specific heat was also reported for the diocotron mode of the guiding center plasma  alias point vortex gas, and for certain vorticity structures in geostrophic flows . However, very different from the gravo-thermal type negative specific heat that we report here to be a characteristic also of the magnetically self-confined plasma torus, the negative specific heat of these quasi-particle systems does *not* couple to the thermal motion of the underlying physical particle systems, which is evident from the fact that these quasi-particle systems also exhibit negative temperature . To pave the way for the discussion we first describe the model and our results.

The model

Since the whole problem is rotationally invariant, we work with conventional cylindrical coordinates $`r,\theta ,z`$. The magnetic induction field decomposes accordingly as $`𝐁=𝐁_T+𝐁_P`$, where $`𝐁_T𝐞_\theta `$ and $`𝐁_P𝐞_\theta `$. In an actual tokamak, the toroidal component $`𝐁_T`$ and a part, $`𝐁_0`$, of the poloidal component are externally generated harmonic fields that serve the purpose of azimuthal and radial stabilization. The total poloidal induction $`𝐁_P=\mathrm{\Psi }\times \theta `$ is the sum of $`𝐁_0`$ and a component which is generated by the electric plasma current density vector, $`j𝐞_\theta `$, via the toroidal Ampére’s law

$$r\left(r^2\mathrm{\Psi }\right)=4\pi c^1j.$$
(1)
In a similar manner one decomposes the electric field, $`𝐄=𝐄_T+𝐄_P`$, where $`𝐄_T`$ is driving the plasma current while the poloidal part is determined by Coulomb’s law $`𝐄=4\pi \rho `$, where $`\rho `$ is the electric charge density. Usually the so-called quasi-neutrality approximation is invoked, which determines the poloidal electric field in leading order through a singular perturbative approach to Coulomb’s law. To keep matters as simple as possible, we consider an ‘electron-positron’ plasma, which is totally charge symmetric with regard to the particle species so that $`\rho `$ vanishes exactly. In that case the poloidal electric field vanishes identically, too, while the toroidal one is implicitly contained in the electric plasma current $`I=Nq\omega /2\pi `$. Here $`N`$ is total number of plasma particles, $`q`$ the elementary charge, and $`\omega `$ the mean absolute angular frequency of a species. This settles the electromagnetic part of the model, and we turn to the statistical mechanics part.

While most works in  are formulated at the macroscopic magnetofluid level, they can in essence be recovered from the statistical mechanics approach of Kiessling et al. in , which begins with the Hamiltonian $`N`$ particle formulation and takes the kinetic limit. At this kinetic level, the plasma particles are described by distribution functions on single-particle phase space. We seek those distribution functions which maximize the familiar Boltzmann entropy functional under the constraints that the two separating integrals of motion, particle number and energy, take prescribed values $`N`$ and $`E`$, and given that the plasma carries a prescribed electric current $`I`$. Since the current is not a separating integral, one has to resort to a ruse and prescribe each species’ canonical angular momentum, which in the axisymmetric kinetic limit is a separating integral, and subsequently pass from this microcanonical angular momenta ensemble to its canonical convex dual, characterized by prescribed $`\omega `$, viz. $`I`$. The solutions of the corresponding Euler–Lagrange equations for this variational principle are also stationary solutions to the axisymmetric kinetic equations of Vlasov. Over velocity space the resulting distribution functions are simply rigidly rotating Maxwellians with temperature $`T=(k_\mathrm{B}\beta )^1>0`$, tied to the energy constraint, and angular frequencies $`\pm \omega (I)`$, microcanonically tied to the angular momentum constraints but canonically prescribed. This allows one to explicitly integrate over the velocity space to retain only the effective macroscopic entropy principle for the total number density of plasma particles, $`n(𝐱)=\overline{n}(r,z)`$, which in our charge symmetric plasma is just twice the value of the common space-dependent Boltzmann factor of each species’ distribution function. In effect this renders the entropy a functional of $`n`$,

$$S(n)=k_\mathrm{B}n(𝐱)\mathrm{ln}\left(\lambda _{\mathrm{dB}}^3n(𝐱)/2\right)𝑑𝒯+\frac{5}{2}Nk_\mathrm{B},$$
(2)
where $`\lambda _{\mathrm{dB}}=h/\sqrt{2\pi mk_\mathrm{B}T}`$ is the thermal de Broglie wave length. This entropy functional has to be maximized under the constraints that $`n0`$ is axisymmetric, $`n𝑑𝒯=N`$, and that the effective energy functional666The negative sign in front of the magnetic energy in $`W`$ is due to the canonical constraint of prescribed electric current; cf. the negative sign in front of the centrifugal contribution to the kinetic energy of a rotating thermal system in the co-rotating frame, see Landau-Lifshitz (op. cit., pp. 71-73). Incidentally, those very centrifugal contributions to $`W`$ are negligible in our plasma and have been omitted. Moreover, the toroidal field $`𝐁_T`$ does not show since we consider only axisymmetric configurations with toroidal current density.

$$W(n)=\frac{1}{2}n(𝐱)K(𝐱,𝐱^{})n(𝐱^{})𝑑𝒯𝑑𝒯^{}+\frac{3}{2}Nk_\mathrm{B}T$$
(3)
takes a prescribed value, say $`E`$. Here, $`K(𝐱,𝐱^{})=\left(2\pi I/cN\right)^2G(𝐱,𝐱^{})`$, and $`G`$ is the Green’s function for $`\left(r^2\right)`$ in $`𝒯`$ for boundary conditions detailed below. With the help of a variant of Moser’s corollary  of the Trudinger–Moser inequality it can be shown that given $`I`$, the entropy functional $`S(n)`$ takes its finite maximum on the set of non-negative axisymmetric densities $`n(𝐱)=\overline{n}(r,z)`$ satisfying $`n(𝐱)𝑑𝒯=N>0`$ and $`W(n)=E`$, and the maximizer is a regular solution of the Euler–Lagrange equation.777These are quite nontrivial facts. In particular, all this is not true if we relax the condition of axisymmetry. We remark that more than one maximizing density function $`n`$ might exist, and in addition the nonlinear Euler–Lagrange equation may have other types of solutions. We call a solution *$`S`$ stable* if it is a *global* maximizer of the entropy (for the given constraints),  *$`S`$ meta-stable* if it is merely a *local* maximizer, and unstable otherwise. Of course, only those $`S`$ stable solutions which exhibit magnetic self-confinement are of interest.

Results

Explicitly carrying out the variations and converting the Euler–Lagrange equation for $`n`$ into an equation for $`\mathrm{\Psi }`$, using $`\mathrm{\Psi }(𝐱)=c^1G(𝐱,𝐱^{})j(𝐱^{})/r^{}𝑑𝒯^{}`$ and $`j(𝐱)=qn(𝐱)\omega r`$, we obtain Pfirsch’s  nonlinear Grad–Shafranov equation

$$\left(r^2\mathrm{\Psi }\right)=8\pi ^2c^1I\frac{e^{\beta \omega q\mathrm{\Psi }/c}}{e^{\beta \omega q\mathrm{\Psi }/c}𝑑𝒯},$$
(4)
which is to be solved in the torus $`𝒯`$ for the boundary conditions encoded in $`G`$. Solving (4) is in general only possible numerically on a computer. However, some explicit analytical control is available if one simplifies the actual laboratory geometry somewhat and considers a torus $`𝒯`$ with rectangular cross section $`\{r,\theta ,z|r_i<r<r_o;\theta fixed;0<z<H\}.`$ The poloidal flux function $`\mathrm{\Psi }(𝐱)=\overline{\mathrm{\Psi }}(r,z)`$ is assumed to satisfy periodic conditions at the $`z`$ boundary and to be constant at $`r_i`$ and $`r_o`$, so that the radial component of $`𝐁_P`$ vanishes at the inner and outer boundaries of $`𝒯`$. In this setting the harmonic poloidal part is simply a homogeneous $`𝐁_0𝐞_z`$, which we choose so that $`\overline{\mathrm{\Psi }}(r_i,z)=\overline{\mathrm{\Psi }}(r_o,z)`$. By gauge freedom we can now even set $`\overline{\mathrm{\Psi }}(r_i,z)=0`$. Beside the desired self-confined configurations, these boundary conditions allow also unconfined ones, namely Pfirsch’s toroidal sheet pinch , given by the following $`𝐞_z`$-invariant solution of (4),

$$\mathrm{\Psi }_{\mathrm{Pf}}(r)=\frac{2c}{\beta \omega q}\mathrm{ln}\frac{\mathrm{cosh}(\kappa ^2[2r^2r_o^2r_i^2]/2)}{\mathrm{cosh}(\kappa ^2[r_o^2r_i^2]/2)},$$
(5)
with $`\kappa (0,\mathrm{})`$ a parameter and $`\beta (\kappa ^2)`$ given by

$$\beta =4c^2q^2N^1\omega ^2H\kappa ^2\mathrm{tanh}(\kappa ^2[r_o^2r_i^2]/2).$$
(6)
While these solutions do not describe a magnetically self-confined plasma torus, they serve as our jumping off point for the numerical computations of the confined configurations. Our strategy, which was also contemplated by K. Schindler,888Private communication. is to look for ringlike bifurcations from the toroidal sheet solution (5). At an infinite sequence of discrete values $`E_1>E_2>\mathrm{}`$, with $`E_kE_{\mathrm{}}=\pi ^2I^2(r_o^2r_i^2)/2Hc^2`$, other solutions bifurcate off of the sheet pinch sequence, breaking its $`z`$ invariance. The bifurcation points are determined by setting $`\mathrm{\Psi }(𝐱)=\mathrm{\Psi }_{\mathrm{Pf}}(r)+ϵ\psi (r,z)+O(ϵ^2)`$, with $`\psi (r_{i,a},z)=0`$ and $`\psi (r,z+H)=\psi (r,z)`$, and expanding (4) to first order in $`ϵ`$, giving the linearized problem

$$\left(r^2\psi \right)+V\psi V\psi =0,$$
(7)
where $`\psi =\psi (r,z)V(r)𝑑𝒯/V(r)𝑑𝒯`$, and

$$V(r)=8\kappa ^4\mathrm{sech}^2(\kappa ^2[2r^2r_o^2r_i^2]/2).$$
(8)
By Fredholm’s alternative , the solution of (7) is trivial except for certain discrete values of $`\kappa `$ at which the bifurcations occur. We have proved  that for our $`𝒯`$ all bifurcations off of the sheet pinch are due to modes $`\psi _k`$, $`k=1,2,\mathrm{}`$, that satisfy $`\psi _k=0`$. The first mode is of the form $`\psi _1(r,z)=R(r)\mathrm{cos}(2\pi [zz_0]/H)`$, with $`z_0`$ arbitrary, and with $`R(r)`$ satisfying

$$r\left(r^1R^{}\right)^{}+(2\pi /H)^2R+r^2VR=0$$
(9)
for $`R(r_i)=R(r_o)=0`$. With realistic domain dimensions $`r_i=1`$, $`r_o=\sqrt{2}`$, and $`H=2`$, a standard Runge-Kutta solver finds the unique nontrivial solution at $`\kappa =\kappa _1=1.62`$, giving $`E_1=2.72W_{}`$, with $`W_{}=2\pi ^2r_iI^2/25c^2`$. Numerical solutions of (4), with $`r_i=1`$, $`r_o=\sqrt{2}`$, and $`H=2`$, were then computed with a well-tested bifurcation code , based on a continuation method . Our code reproduced the analytical sheet pinch solution and its first bifurcation point in excellent agreement with our independently obtained semi-analytical results. We then numerically followed the first bifurcating branch that emerges from $`\psi _1`$ to nonlinear amplitudes, where it develops into a toroidal ring pinch with a double X magnetic structure similar to the double X structure in the PDX-PBX Tokamak experiment in Princeton.999For design information, see http://www.pppl.gov/oview/pages/pbxm\_design.html. Retrospectively, this vindicates our choice of boundary conditions.

(FIG. 1)

Our primary interest is in the energy-entropy diagram. Shown in Fig. 2 is $`\mathrm{}S(n)`$ versus $`W(n)`$, where $`\mathrm{}S(n)=S(n)Nk_\mathrm{B}\mathrm{ln}\sqrt{(4\pi e/N)^5(I/\mathrm{}c)^6m^3r_i^3}`$, with $`e`$ the Euler number, and with the density function $`n`$ running along the computed bifurcation sequences of ring and sheet pinch. At sufficiently high effective energies, Pfirsch’s sheet pinch is the unique solution of (4) for the stipulated boundary conditions, hence maximizing entropy. Numerically it appears to be the case for all $`W(n)>E_1=2.72W_{}`$, see Fig. 2. For all $`W(n)<E_1`$ down to $`W(n)=0.5W_{}`$ where we terminated the computation, the ring pinch has higher entropy than the toroidal sheet pinch at the same effective energy. By asymptotic analysis we found that also for $`W(n)\mathrm{}`$, and by continuity for $`W(n)W_{}`$, the maximum entropy configuration consists of a highly concentrated ring pinch which, in rescaled coordinates centered at the density maximum, converges to Bennett’s cylindrical pinch  as $`W(n)\mathrm{}`$. On the basis of this evidence we surmise that the ring pinch has maximum entropy for all $`W(n)<E_1`$, implying its $`S`$ stability in the class of rotationally invariant plasma with effective energy $`W(n)<E_1`$ and current $`I`$. We remark that the first bifurcation off of the toroidal sheet pinch into the $`S`$ stable toroidal ring pinch branch is then a symmetry-breaking second-order phase transition.

It remains to determine the specific heat of the configurations, which we recall is negatively proportional to the second derivative of $`S`$ with respect to $`E`$. Thus we inspect the curvature of the graphs of the entropy as function of energy for the various solutions, given in Fig. 2. The graph representing the sheet pinch is concave. However, the graph for the ring pinch is manifestly convex over the whole computed range of energies $`0.5W_{}<W(n)<E_1`$. We were also able to prove the convexity analytically to second order in perturbation theory away from the bifurcation point. This confirms what we have announced earlier: the specific heat of the ring pinch is negative!

(FIG. 2)

Discussion

At last, we discuss the potential implications that our finding of gravo-thermal type negative specific heat has for the problem of toroidal magnetic self-confinement of plasma. In a gravitationally bound plasma, negative specific heat on one hand aids the ignition of nuclear burning in a proto-star by heating it up when it loses energy by radiation, but it is also responsible for some more spectacular instabilities once the nuclear burning expires, like the onset of the red giant structure . It would be intriguing enough if the negative specific heat of a magnetically self-confined plasma torus should be confirmed to aid the ignition of nuclear burning in a tokamak. For this to be so, one would have to be able to hold $`N`$ and $`I`$ fixed and secure the toroidal invariance (which is what one wants to achieve anyhow), while $`E`$ would have to decrease (the plasma radiation would seem to help in this respect) *slow enough* so that one would essentially evolve along the ring pinch branch in Fig.2 to the left, thereby heating up the plasma while pinching it more strongly. This would not seem unwelcome. For now, however, energy leakage by radiation is a serious problem, while at the same time the *emf* current drive leads to yet uncontrolled ohmic heating of the plasma. In this case where $`E`$ is allowed to fluctuate too widely, the negative specific heat will have a very unwanted effect on the confinement. This can be illustrated by considering the temperature rather than energy $`E`$ to be controlled by the competition of ohmic heating and radiation (still assuming $`N`$ and $`I`$ fixed, and toroidal invariance). In that case the canonical ensemble determines the stability. But microcanonical and canonical ensembles are not equivalent when the microcanonical one exhibits states with negative specific heat , and sure enough, none of the computed ring pinches with negative specific heat minimizes the free energy functional

$$F(\mathrm{\Psi })=\frac{1}{8\pi }r^2|\mathrm{\Psi }|^2𝑑𝒯+N\beta ^1\mathrm{ln}\left(\frac{2e}{N\lambda _{\mathrm{dB}}^3}e^{\beta \omega q\mathrm{\Psi }/c}𝑑𝒯\right).$$
(10)
Actually, $`F`$ is unbounded below for these $`\beta `$, $`N`$ and $`I`$ values, any minimizing sequence concentrating on a singular ring current; cf.  for a good discussion of the translation-invariant analog. Of course, a real plasma would not get anywhere near such a singular ring current configuration, for a highly concentrated plasma ring is known to be susceptable to magnetofluid dynamical instabilities that destroy the axisymmetry.

Conclusion

To summarize, the $`S`$ stable ring pinches have negative specific heat of the gravo-thermal type and will therefore be stable if and only if $`N`$, $`E`$, and $`I`$ are essentially fixed and the toroidal invariance is secured, in which case the negative specific heat may aid the ignition of thermonuclear burning. The “if” part is good news; the bad news is the “only if” part.

MK gratefully acknowledges financial support through a DFG Habilitations-Fellowship in early, and through NSF Grants # DMS-9623220 and # DMS-0103808 in later stages of this project. TN gratefully acknowledges support by an ESA Research Fellowship in early, and by a PPARC Advanced Fellowship in later stages of this project. Sincere thanks go to S. Goldstein for so many illuminating conversations.


text:
# Decay process accelerated by tunneling in its very early stage

## Abstract

We examine a fast decay process that arises in the transition period between the Gaussian and exponential decay processes in quantum decay systems. It is usually expected that the decay is decelerated by a confinement potential barrier. However, we find a case where the decay in the transition period is accelerated by tunneling through a confinement potential barrier. We show that the acceleration gives rise to an appreciable effect on the time evolution of the nonescape probability of the decay system.

Quantum decay processes such as the nuclear alpha decay are usually well described by means of the exponential decay law; see, for example, ref:Gamow ; ref:We-Wi ; ref:Br-Wi . Theoretically, however, deviations from the exponential law are expected in the beginning and also toward the end of a decay process ref:review . It is understood that the quantum decay process in general proceeds through three different stages; initial, intermediate and final. The initial stage is characterized by the Gaussian law, the intermediate stage by the exponential law, and the final stage by the power law. The decay speeds in the initial and final stages are smaller than that in the intermediate stage. In particular, the slow decay process in the Gaussian period leads to the possibility of the quantum Zeno effect ref:review , in which the decay process is decelerated by repeated measurements.

Many years ago, in his one-dimensional model analysis of a decaying system, Winter found that the speed of the decay process exhibits an irregular behavior in the transition period between the Gaussian and exponential periods ref:Win . Very recently Dicus et al. reexamined the same system, which consists of a particle which is initially confined in a region and leaks out by tunneling through a delta-function potential barrier ref:DICUS . In the irregular decay process, the decay speed can be larger than that in the exponential region. Such fast decay in the transition period is interesting in the sense that it may give rise to acceleration of the decay process by repeated measurements, which is the so-called anti- (inverse-) Zeno effect ref:IQZE ; ref:Fisc .

The purpose of this report is to investigate the details of the fast decay process in the transition period between the Gaussian and exponential periods footnote:0 . In particular, we explore how the decay process in the transition period depends on the strength of the confinement potential barrier with a finite potential width. It is usually expected that the decay is suppressed as the confinement becomes stronger because the potential barrier suppresses the time evolution of the wave function. However, we will illustrate a case where the decay in the transition period is accelerated by tunneling through a confinement potential barrier. Further, we show that the acceleration gives rise to an appreciable effect on the time evolution of the nonescape probability of the decay system.

We consider a model in one dimension with a potential that consists of an infinite repulsive wall at $`x=0`$ and a repulsive square barrier at $`1<x<1+w`$,

$`V(x)=\{\begin{array}{cc}0& \mathrm{for}0<x<1,\hfill \\ h& \mathrm{for}1<x<1+w,\hfill \\ 0& \mathrm{for}1+w<x,\hfill \end{array}`$ (4)
where $`h`$ and $`w`$ are the height and width of the potential, respectively.

We assume that a particle of mass $`m`$ is initially confined within the potential barrier and it leaks out in time. The time-evolution of the system is determined by the time-dependent Schrödinger equation

$`i{\displaystyle \frac{\varphi (x,t)}{t}}=H\varphi (x,t),H={\displaystyle \frac{^2}{x^2}}+V(x)`$ (5)
where the units are such that $`\mathrm{}=1`$ and $`2m=1`$. We solve Eq. (5) numerically using the implicit solution methodref:NC , with the unitary time-evolution operator $`U(1iH\delta t/2)/(1+iH\delta t/2)`$, where $`\delta t`$ is the time mesh. For the range of $`x`$, we take . In order to suppress reflection of the wave function at the boundary ($`x=500`$) of the model space, we assume an imaginary potential at $`x490`$.

For the initial wave function of the particle which is confined in the region of $`0x1`$, we assume

$$\varphi (x,0)=\sqrt{2}\mathrm{sin}\pi x.$$
(6)
It is understood that $`\varphi (x,0)=0`$ for $`x>1`$. The wave function leaks out in time by tunneling through the potential barrier $`V(x)`$. The energy expectation value $`<H>`$ at $`t=0`$ is $`\pi ^2`$. In numerical illustrations we choose the height of the potential barrier such that $`h>\pi ^2`$. In this sense we regard the process as a tunneling decay process.

We define probability $`P(a,t)`$ for the particle being in the interval $`0xa`$ at $`t`$ by

$`P(a,t)={\displaystyle _0^a}𝑑x|\varphi (x,t)|^2.`$ (7)
If we set $`a=1+w`$, $`P(1+w,t)`$ represents the probability that the particle is confined by the potential barrier. For $`w`$ we consider a few different values. In order to compare the results for different values of $`w`$, we set $`a`$ such that $`a>1+w_{max}`$. Throughout this paper, we take $`a=4`$. The $`P(a,t)`$ is a kind of the nonescape probability, which is the probability that the particle has not escaped from the potential by time $`t`$ ref:F-T ; ref:GC-M-M .

Next we introduce function $`g(a,t)`$ defined by

$`g(a,t)={\displaystyle \frac{dP(a,t)/dt}{P(a,t)}}.`$ (8)
If probability $`P(a,t)`$ decays exponentially, that is, if $`P(a,t)e^{\gamma t}`$, then $`g(a,t)`$ is independent of time,

$`g(a,t)=\gamma .`$
If $`P(a,t)`$ obeys the Gaussian decay law, that is, if $`Pe^{t^2/\tau }`$, $`g(a,t)`$ is proportional to $`t`$:

$`g(a,t)=2t/\tau .`$ (9)
Thus, from the $`t`$ dependence of $`g(a,t)`$ we can see how well $`P(a,t)`$ obeys the exponential law or the Gaussian law.

Figure 1 shows the nonescape probability $`P(a,t)`$, where the potential height and width are taken as $`h=10`$ and $`w=0.6`$, respectively. In the initial period $`0t0.3`$, the decay is extremely slow. The period should correspond to the Gaussian decay process. In the period of $`t2`$, the system is subject to the exponential decay law. Between the two stages, there is a period ($`0.3t2`$) in which the decay process is neither Gaussian nor exponential. This is the transition period. In the period the decay of the nonescape probability seems to be much faster than that in the exponential period.

The function $`g(4,t)`$ calculated from $`P(4,t)`$ of Fig. 1 is shown by the solid line in Fig. 2. In the period of $`t2`$, $`g(4,t)`$ is almost constant, which means that the decay process is exponential. In the period of $`0t0.3`$, $`g(4,t)`$ is not exactly proportional to $`t`$, i.e., the decay process in the initial stage slightly deviates from the Gaussian decay law. However, the decay speed in this period is still smaller than that in the exponential period. Corresponding to the rapid decay process seen in $`P(t)`$ of Fig. 1, the maximum of $`|g(4,t)|`$ is obtained in the transition period ($`t0.6`$). The $`g(4,t)`$ starts from zero at $`t=0`$. Thus, the quantum Zeno effect is possible when we repeat measurements with a sufficiently small time interval. On the other hand, the anti-Zeno effect is possible by repeated measurements only when the net decay rate of the fast decay process in the transition period is large compared to that of the slow decay process in the Gaussian period (see footnote footnote:1 ).

In order to see how the decay process in the transition period depends on the strength of the potential barrier, we examine $`g(4,t)`$ for various potential heights and widths. The dashed and dotted-dashed lines in Fig. 2 show the $`g(4,t)`$ for $`h=20`$ and $`h=30`$, respectively with a fixed width $`w=0.6`$. The decay process in the Gaussian period $`0t0.3`$ does not depend strongly on the potential height. This is because in the initial stage the higher energy components of the initial wave function relative to the potential height contribute mainly to the decay of the system. On the other hand, the decay speed in the exponential period becomes much smaller as the confinement becomes stronger. The fast decay in the transition period depends strongly on the potential height. It tends to be suppressed as the confinement becomes stronger.

In Fig. 3, we show $`g(4,t)`$ for various potential widths and a fixed height $`h=15`$. The decay process in the Gaussian period $`0t0.3`$ does not seem to depend strongly on the potential width. On the other hand, the decay rate of the exponential period becomes much smaller as the potential width becomes broader. This is due to the increase of the confinement strength. The decay speed in the transition period becomes larger as the potential width becomes narrower.

As we have shown, the speed of the fast decay process becomes smaller as the potential barrier becomes stronger. Thus, we might guess that the fastest decay will be obtained in the decay process with no potential barrier. However, as we will show in the following, potential barriers with appropriate widths and heights can accelerate the decay process.

In Fig. 4, the solid line shows $`g(4,t)`$ calculated with no potential barrier. The dashed, dot-dashed and dotted lines exhibit $`g(4,t)`$ calculated with the potential barriers with $`h=10,20`$ and $`30`$ and fixed width $`w=0.2`$, respectively. It should be noted that $`w=0.2`$ is much thinner than those used in Figs. 2 and 3, and therefore, in this case, the confinement is very weak compared with that in Figs. 2 and 3. In the Gaussian period, the potential dependence of the decay speed is not very appreciable. However, in the transition region, the decay speed becomes faster as the potential height becomes lower. For $`h=10`$, the maximum decay speed at $`t1.8`$ exceeds that for no potential barrier at $`t0.7`$.

As shown in Fig. 5, such an acceleration of the decay speed by tunneling gives rise to an appreciable difference in the time evolution of $`P(4,t)`$. The solid and dashed lines are the $`P(4,t)`$ for no potential barrier and for the potential barrier with $`w=0.2`$ and $`h=10`$, respectively. The nonescape probability for $`h=10`$ becomes smaller than that for no potential barrier at $`t1.5`$. At this time the residual nonescape probability is still about ten percent. In this sense, the effect of this acceleration cannot be ignored. On the other hand, for $`h=20`$, at the time region in which the nonescape probability becomes smaller than that for no potential barrier, the residual nonescape probability is negligibly small.

The fluctuations in the behaviors of $`g(4,t)`$ indicate that the decay processes are still in the transition period from the Gaussian to exponential period. However, as shown in Fig. 5, the decay process has been almost completed before $`t=4`$. Therefore, even if the decay process proceeds to the exponential period eventually, the exponential decay has no importance in this case. The stability of quantum system is usually characterized by the magnitude of the imaginary part of the pole that gives the inverse of the lifetime in the exponential period ref:review ; ref:Win ; ref:DICUS . The result that we have shown implies that such a pole analysis may not be effective for highly unstable systems ref:ModOpt .

In Fig. 4, we have seen that the acceleration by tunneling can be obtained for $`h10`$. Next, we investigate the acceleration with a fixed height. We examine the time evolution of $`g(4,t)`$ for $`h=10`$ footnote:2 . The dashed, dot-dashed and dotted lines in Fig. 6 exhibit $`g(4,t)`$ calculated for the potential barriers with $`w=0.2,0.4`$ and $`0.6`$ and fixed height $`h=10`$, respectively. One can see that the decay rate becomes larger for thinner potential widths. For $`w=0.2`$, the maximum decay speed exceeds that for no potential barrier. Thus, one sees that the acceleration of the decay speed by tunneling can be obtained when the strength of the confinement by the potential barrier is sufficiently weak. In our illustrations, for $`h10`$ and $`w0.2`$, the accelerations are remarkable.

Finally, we mention that the $`g(4,t)`$ represented by dotted line in Fig. 4 ($`w=0.2`$ and $`h=10`$) takes positive values around $`t3.5`$. Recall that the probability current $`j(a,t)`$ is related to the nonescape probability $`P(a,t)`$ by $`j(a,t)=dP(a,t)/dt=g(a,t)P(a,t)`$. This means that, if $`g(a,t)`$ is positive, $`j(a,t)`$ is negative. However, in ref:Win ; ref:NTD  the negative currents were obtained at very late time region after the exponential decay period. Our result implies that in a highly unstable quantum system the negative current can occur even at the initial stage after the Gaussian period.

We have investigated the fast decay process in the transition period between the Gaussian and exponential decay processes. In most cases of the tunneling process, the decay speed becomes smaller as the potential barrier becomes stronger. As a special case, we have found that the fast decay process can be remarkably accelerated by tunneling through potential barriers with appropriately small widths and heights. A detailed analysis of the acceleration of the fast decay process by tunneling is a future project.

Acknowledgment

We would like to thank Professor Y. Nogami for useful discussions and comments. This work was supported in part by the Minstry of Education, Culture, Science, Sports and Technology of Japan.