text:
# Field emission theory beyond WKB - the full image problem

## Introduction

The classic theory of electron field emission from a cold metal surface due to Fowler and Nordheim (FN)  was one of the earliest applications of wave mechanics since its foundation in the early 1920’s. With increased interest in field emission processes, especially for use in display screens , it is worth examining closely the foundations of this theory. As well as metal-vacuum interfaces, proposed applications of the theory includes graphite-vacuum, graphite-insulator-vacuum, carbon-nanotube-vacuum and Schottky junctions . Since our approach can be generalized to other cases of interfaces between dissimilar media, it may have application to the systems discussed by Mott and Gurney  and Herring . These other cases might include vacuum breakdown and electrical forming processes . In section 5 we show that our new analysis removes discrepancies of many orders of magnitude between theory and experiment.

In the FN theory, the role of the Fermi surface, Fermi-Dirac statistics and in particular the unique property of quantum mechanical tunnelling were key concepts. The essentially equilibrium distribution of electrons is modeled by a one dimensional ($`x`$ dependent only) Hamiltonian:

$$H_0=\frac{p_x^2}{2m}+V(x),$$
(1)
where all energies are measured with the same reference as the effective potential $`V(x)`$. Then the current density $`J`$ of the emitted electrons is given by an integral over the energy $`W`$:

$$J=e_{W_a}^{\mathrm{}}P(W)𝑑W,$$
(2)
where $`P(W)dW`$ gives the number of electrons within $`dW`$ that emerge form the metal per second per unit area, with $`W_a`$ an appropriate lower energy reference. The quantity $`P(W)`$ is easily expressed as the product of a quantum mechanical tunnelling transmission coefficient $`T(W)`$ and a Fermi type supply function $`N(W)`$:

$$N(W)=\frac{4\pi mkT}{h^3}\mathrm{ln}\left\{1+\mathrm{exp}\left[(\frac{Wϵ_F}{kT})\right]\right\},$$
(3)
in which $`ϵ_F`$ is the Fermi energy. In view of the complexity of the Schrodinger equation, the transmission problem was originally treated within a WKB approximation or variants of it :

$$T(W)=\mathrm{exp}\left[_{x_1}^{x_2}\left(\frac{8m}{\mathrm{}^2}(V(x)W)\right)^{1/2}𝑑x\right],$$
(4)
in which $`x_1`$ and $`x_2`$ are the classical turning points of the potential $`V(x)`$. It is a well known result however that the validity of the WKB approximation involves a criterion, see section 1, that can be derived from the requirement that the variation of the deBroglie wavelength of the electron must be less than the dimension of the region of tunnelling . This restriction on the validity of the FN theory was unfortunately not discussed in many of the earlier papers, including the classic paper of . The limitations of the WKB approach are well known, however, in the treatment of thermionic emission; a common technique for improving the calculation is to make a parabolic approximation to the top of of the barrier. Nevertheless, FN theory has been applied in regimes where it is not valid. In section 2 we shall examine the image-free problem exactly 111Since image effects are less important for silicon based Schottky junctions due to the large dielectric constant $`ϵ_s11.9`$, these results are actually very relevant to the latter, which however we shall not take up here in this paper., whose solution in one dimension is amenable to analytical treatment in terms of Airy type functions. Some preliminary assessment of the FN-WKB can already be made for this case, which is a preamble to the full image problem which we shall treat in section 3 by various analytical methods. Interestingly the tunnelling problem close to the barrier top appears to be amenable to analytical treatment as well and has not been published to the best of our knowledge. We shall note that the failure of the WKB criteria coincides quite closely to this regime, hence the availability of analytical results are extremely welcomed. Next we shall note that the role of screening has also been omitted in the original FN theory. We shall introduce the additional screening potentials in section 4 and thereby propose that this is one mechanism for a significant barrier reduction, in addition to negative electron affinity attributed to diamond like carbon systems. In section 5 we shall apply our theory to the analysis of experimental data for carbon field emitters. Regardless of reasonable assumptions of point tip field enhancement effects, we shall see that the FN theory does not even get close to the order of magnitude agreement with the data by comparison with ours. We shall conclude in section 6 with a summary of our results, postponing further solutions of the problem involving screening to a future paper.

## 1 WKB Validity criteria

Following  it is straightforward to derive the criterion for the validity of WKB:

$$p^3>>m\mathrm{}eF,$$
(5)
where $`F`$ is the external field and we shall ignore the electron image for the moment. Here $`p`$ is the electron momentum so that:

$$(2m(\varphi eFx))^{3/2}>>m\mathrm{}eF,$$
(6)
where $`\varphi `$ is the work function. Assuming a carrier effective mass equal to the free electron mass and inserting numerical values appropriate to a field in V/$`\mu `$m and a work function in eV, (which shall be our units throughout this paper) gives

$$\zeta (F,\varphi )=0.046\frac{F^{1/3}}{(\varphi )^{1/2}}<<1,$$
(7)
in which the tunnelling length $`x`$ has been set to zero, representing the best case for WKB, the very top of the barrier ($`x`$ is typically of the order of nm). Typical numerical values for $`\zeta `$ are given table 1, where we see that the WKB approximation begins to be questionable for fields above about $`10\mathrm{V}/\mu \mathrm{m}`$ for typical work function values and is certainly seriously questionable in the range close to dielectric breakdown fields which are typically $`10^3\mathrm{V}/\mu \mathrm{m}`$ . A smaller effective mass such as for p-Si where ($`m^{}=0.6`$ for holes) would increase $`\zeta `$ further. Later on in the paper we shall show further anomalies with the FN theory when we are close to this limit of validity. In order to proceed and be able to assess the approximation further we shall have to resort to a full solution of the Schrodinger equation which we shall begin in the next section by first ignoring the image potential.

## 2 Exact solution without electron images

We shall start with the image free tunnelling problem given by the Schrodinger equation:

$$\psi ^{\prime \prime }(x)+\frac{2m}{\mathrm{}^2}[WV(x)]\psi (x)=0,$$
(8)
where $`W`$ is the energy eigenvalue of the Hamiltonian $`H_0`$ with the potential given by :

$$V(x)=\{\begin{array}{cc}W_a\hfill & \text{if }x<0\text{ – region 1;}\hfill \\ V_0eFx\hfill & \text{if }x>0\text{ – region 2.}\hfill \end{array}$$
(9)
Here $`V_0`$ under normal circumstances is zero which yields the work function $`\varphi `$ as the barrier height and $`F`$ is as before the applied external field. For later convenience we shall recast Eq(8) into the more compact form:

$$\psi ^{\prime \prime }(x)+(ϵ+\alpha x)\psi (x)=0,$$
(10)
with

$`\alpha `$ $`=`$ $`{\displaystyle \frac{2meF}{\mathrm{}^2}}F\times 2.55989\times 10^7(\mu \mathrm{m})^3,`$
$`ϵ`$ $`=`$ $`{\displaystyle \frac{2m(WV_0)}{\mathrm{}^2}}(WV_0)\times 2.55989\times 10^7(\mu \mathrm{m})^2,`$ (11)
in which $`F`$ is measured in $`\mathrm{V}/\mu \mathrm{m}`$. To the best of the authors’ knowledge this elementary quantum exercise does not seem to have been fully documented in the literature, although variants of it have been scattered around . Following Landau and Lifshitz , we shall transform to the new variable $`\xi =(x+ϵ/\alpha )\alpha ^{1/3}`$. Then equation (10) in region 2 becomes the typical Airy’s equation 222The advantage of using the Airy function formulation is that it naturally takes care of the long range boundary condition posed by the electric field.:

$$\psi ^{\prime \prime }(\xi )+\xi \psi (\xi )=0.$$
(12)
However some care needs to be exercised in selecting the appropriate type of Airy functions, since the standard Ai function has asymptotic form :

$$\psi (\xi )=\mathrm{Ai}(\xi )\xi ^{1/4}\mathrm{sin}(\frac{2}{3}\xi ^{3/2}+\pi /4).$$
(13)
As this also contains a reflected wave, it is not appropriate to our boundary conditions. The correct approach is to revert to the original 1/3 fractional Bessel functions definition for the Airy functions and then construct appropriate linear combinations to yield the forward travelling wave solution only. This approach of course leads to a modified appropriate Hankel function which as we shall see has the correct asymptotic form required. It is this construction that seems to be rather rarely seen in the literature except perhaps in the WKB derivation of Moll  . With this proviso, we can write down the solution for the wave functions of equation (12) as:

$$\psi (x)=\{\begin{array}{cc}e^{ikx}+re^{ikx}\hfill & \text{in region 1;}\hfill \\ t[u_1(\xi )e^{i\pi /3}u_2(\xi )]\hfill & \text{in region 2,}\hfill \end{array}$$
(14)
where $`k=(2m(W+W_a)/\mathrm{}^2)^{1/2}`$ and the $`u_1`$ and $`u_2`$ are given by Bessel functions of order $`1/3`$, here defined as:

$$u_1(\xi )=\frac{1}{3}(\pi \xi )^{1/2}J_{1/3}(\frac{2}{3}\xi ^{3/2}),$$
(15)
and

$$u_2(\xi )=\frac{1}{3}(\pi \xi )^{1/2}J_{1/3}(\frac{2}{3}\xi ^{3/2}),$$
(16)
respectively. Note however that analytical continuation to negative $`\xi `$ requires some care: see Appendix A. It is then straightforward to obtain the reflection and transmission amplitudes by matching the wave functions and their derivatives at the interface $`x=0`$. We provide the results here in terms of $`H_{1/3}=u_1e^{i\pi /3}u_2`$ which is a slightly modified Hankel function

$`t`$ $`=`$ $`{\displaystyle \frac{2}{H_{1/3}(\xi _0)ik_0^{1/3}H_{1/3}^{}(\xi _0)}}`$
$`r`$ $`=`$ $`{\displaystyle \frac{H_{1/3}(\xi _0)+ik_0^{1/3}H_{1/3}^{}(\xi _0)}{H_{1/3}(\xi _0)ik_0^{1/3}H_{1/3}^{}(\xi _0)}},`$ (17)
where $`\xi _0=\xi (0)=\overline{ϵ}=ϵ/\alpha ^{2/3}`$ while $`k_0^{1/3}=\alpha ^{1/3}/k`$. An important check must be made to these formulas which requires that the transmission coefficient:

$$T=\frac{k_0^{1/3}}{4}|t|^2,$$
(18)
and reflection coefficient: $`R=|r|^2`$ satisfy the unitarity property $`T+R=1`$. This is easily checked to be the case when appropriate use is made of the Wronskian identity:

$$u_1u_2^{}u_2u_1^{}=\frac{1}{2\sqrt{3}}.$$
(19)
For sufficiently small fields, we can obtain an asymptotic expansion (where $`\xi _0\mathrm{}`$ and $`k_00`$) that yields:

$$T_{asymp}\frac{4k_0^{1/3}|\xi _0|^{1/2}}{(1+k_0^{2/3}|\xi _0|)}e^{\frac{4}{3}\xi _0^{3/2}},$$
(20)
the prefactor being a new result not obtainable within WKB. For convenience we set $`W_a=2\varphi `$ which for $`\varphi `$ of order $`1\mathrm{eV}`$ is adequate as a lower bound (at room temperatures), the exception being very narrow band semiconductors which we are not concerned with here. We shall at the moment ignore any possible barrier height lowering (i.e. image free case) so that $`V_0=0`$ is the vacuum level and consider tunnelling near the Fermi level, so that $`W=\varphi `$ and thus $`k_0^{2/3}|\xi _0|=1`$ . Then we see that the above prefactor differs from the WKB result by a mere factor of two i.e.:

$$T_{asymp}2e^{\frac{4}{3}\xi _0^{3/2}}.$$
(21)
That this asymptotic result is valid for only small fields is clearly seen now for otherwise if we naively allow the field to be large, then we could arrive at an unphysical result in which $`T>1`$ which is absurd. In fact this occurs for $`\varphi =1\mathrm{eV}`$ when $`F`$ reaches the order of $`9.9\mathrm{kV}/\mu \mathrm{m}`$. Hence the use of the WKB formula for large transmissions without knowledge of the prefactor in the transmission can be dangerous. Unfortunately this caution does not seem to have been exercised by Good and Müller , nor by Murphy and Good . The latter have derived widely used criteria  for the validity of FN theory on the basis of the WKB integrand properties equation (4) alone, which is thus a logically inadequate procedure. If we examine the classic FN formula for the emission current:

$$J=\frac{1.54\times 10^2F^2}{\varphi t^2(0.0379\times F^{1/2}/\varphi )}\mathrm{exp}(\frac{6830\times \varphi ^{3/2}}{F}\upsilon (0.0379\times F^{1/2}/\varphi ))\mathrm{A}\mathrm{cm}^2$$
(22)
in which $`t(x)`$ and $`\upsilon (x)`$ are given in terms of elliptic integrals, for $`\varphi =2\mathrm{eV}`$ say, the exponent practically vanishes at $`F3000\mathrm{V}/\mu \mathrm{m}`$. Hence the entire emission current comes from the supply function which is the $`F^2`$ factor and, as discussed above, this can be untenable. In fact it is well known from elementary quantum mechanics that the transmission coefficient is in general not unity at the top of the barrier .

Later on, in section 3, we shall see that an exact solution is obtainable for the full image theory at the top of the barrier, in which the transmission coefficient has a finite value less than unity, but nevertheless with significant tunnelling current contributions. As far as we are aware, this solution has not been noted in the field emission literature for the last 40 years. In figure 1 we provide two curves that compare the exact transmission coefficient equation (18) for the image free case against the asymptotic formula equation (21) close to the Fermi level for very large fields where the latter becomes problematic. Fortunately these fields are very large and being much greater than typical breakdown fields this failure is unimportant in practice. We note that for smaller fields figure 2 the WKB formula is quite good, in spite of its failure to satisfy our earlier criterion of section 1. The accuracy of the WKB for small fields is however due to us staying away from the top of the barrier. We now compare equation (21) with the exact result equation (18) as we approach the barrier, i.e. $`V_0\varphi `$. One mechanism for this to occur could be due to a large image term i.e. Schottky effect, but there are other screening effects as well as we shall see later. For simplicity we shall not include the image effect in the Schrodinger equation yet but merely consider the effect of approaching the barrier top in the image free case first. Here we can compare the results of three formulas,the exact equation (18), the asymptotic formula equation (20) and the standard WKB result in which the prefactor in equation (20) is missing: these are shown in figure 3, where we have set $`V_0=0.95\varphi `$. Here we see that for quite small fields $`F<100V/\mu \mathrm{m}`$, significant departures from WKB are found. These results motivate us to further include the Schottky image term in the Schrodinger equation which we take up in the next section.

## 3 Exact solution including electron images

We now turn to the quantum tunnelling problem including the effect of the image term as in the original FN theory. Then equation (10) becomes:

$$\psi ^{\prime \prime }(x)+\left[ϵ+\alpha x+\frac{\beta }{x}\right]\psi (x)=0,$$
(23)
where the image potential is characterized by the parameter $`\beta =me^2/(2\mathrm{}^2)`$. Transforming to the parameter $`\xi `$ as before we now have the equation:

$$\psi ^{\prime \prime }(\xi )+\left[\xi +\frac{\overline{\alpha }}{\xi \overline{ϵ}}\right]\psi (\xi )=0,$$
(24)
where $`\overline{\alpha }=\beta \alpha ^{1/3}`$, and $`\overline{ϵ}=ϵ\alpha ^{2/3}`$. Clearly the image-free case $`\overline{\alpha }=0`$ reduces to the usual Airy differential equation equation (12) and the above can be rewritten as:

$$(\xi \overline{ϵ})\psi ^{\prime \prime }(\xi )+\left[\xi (\xi \overline{ϵ})+\overline{\alpha }\right]\psi (\xi )=0.$$
(25)
We shall in fact further rewrite this as:

$$\psi ^{\prime \prime }(\xi )+\left[\frac{(\xi \lambda _1)(\xi \lambda _2)}{(\xi \overline{ϵ})}\right]\psi (\xi )=0.$$
(26)
in which the two roots of the quadratic are given by:

$$\lambda _1=\frac{1}{2}[\overline{ϵ}+\sqrt{\overline{ϵ}^24\overline{\alpha }}]\mathrm{and}\lambda _2=\frac{1}{2}[\overline{ϵ}\sqrt{\overline{ϵ}^24\overline{\alpha }}].$$
(27)
Let us first discuss the WKB approximation , which is obtained by approximating the wavefunction as:

$$\psi (\xi )=\frac{1}{\sqrt{p(\xi )}}\mathrm{exp}\left[i^\xi p(\xi ^{})𝑑\xi ^{}+\frac{i\pi }{4}\right]$$
(28)
where:

$$p(\xi )=\left[\frac{(\xi \lambda _1)(\xi \lambda _2)}{(\xi \overline{ϵ})}\right]^{1/2}$$
(29)
and the lower limit in the integral in equation (28) is the smaller of $`\lambda _1`$ or $`\lambda _2`$ which depends on the sign of $`\overline{ϵ}`$ i.e. above or below the barrier respectively. The evaluation of the integral in equation (28) is then a straightforward exercise in Elliptic integral reduction which we need not repeat here . Equation (26) itself does not appear to have a closed form solution in terms of known functions although several methods of solution are available depending on the parameters. We postpone these to the appendices as we do not intend to pursue detailed calculations here, but merely to demonstrate the regimes where WKB is in need of correction. For this purpose we can concentrate on the region close to the top of the barrier at $`\xi =\overline{\alpha }^{1/2}+\overline{ϵ}`$. Let us start at the barrier top itself for which WKB predicts $`T=1`$ since in this case $`\upsilon (x)=0`$ in equation (22). Then equation (26) simplifies to the form:

$$\psi ^{\prime \prime }(\xi )+\left[\frac{(\xi \overline{ϵ}/2)^2}{(\xi \overline{ϵ})}\right]\psi (\xi )=0.$$
(30)
In view of the magnitude of the parameter $`\overline{ϵ}`$ it suffices to obtain a solution in the limit when this parameter is large. For $`\varphi =1\mathrm{eV}`$ and $`100<F<1000\mathrm{V}/\mu \mathrm{m}`$, $`\overline{ϵ}`$ ranges from 3 to 14. In particular for the matching conditions we are only interested in the solution near to $`\xi _0=\overline{ϵ}`$ which corresponds to the origin $`x=0`$ so that we can transform to the variable: $`y=\xi \overline{ϵ}`$ which is small for our purpose $`y<<\overline{ϵ}`$ so that equation (30) can be approximated by:

$$y\psi ^{\prime \prime }(y)+(y+\overline{ϵ}/2)^2\psi (y)y\psi ^{\prime \prime }(y)+(\overline{ϵ}/2)^2\psi (y)=0.$$
(31)
The appropriate solution for the above is then given by the Hankel function of order $`1`$ which upon transforming back to the original variables $`\xi `$ is given by:

$$\psi (\xi )=e^{i\pi /4}\sqrt{\frac{\pi }{2}}\overline{ϵ}(\xi \overline{ϵ})^{1/2}H_1^{(1)}(\overline{ϵ}(\xi \overline{ϵ})^{1/2}).$$
(32)
We can now obtain the transmission coefficient as was done in the previous section by matching the wave function and its derivative at the origin. The result is now given by:

$$T=\frac{2\overline{ϵ}^2k_0^{1/3}}{|\psi (\xi _0)|^2+ϵ^2k_0^{1/3}+k_0^{2/3}|\psi ^{}(\xi _0)|^2}$$
(33)
where the derivative of $`\psi `$ is easily shown to be given in terms of the Hankel function of order $`0`$:

$$\psi ^{}(\xi )=e^{i\pi /4}\sqrt{\frac{\pi }{8}}\overline{ϵ}^2H_0^{(1)}(\overline{ϵ}(\xi \overline{ϵ})^{1/2})$$
(34)
As usual we can verify that this solution satisfies unitarity $`T+R=1`$ through the use of the appropriate Wronskian identities for the Hankel functions. That this is satisfied exactly in spite of our approximation in equation (31) means that we have picked out the essential features in the theory. Equation (33) gives zero transmission at the top of the barrier because of a weak logarithmic singularity (typical of the $`1/x`$ potential near the origin) in $`\psi ^{}`$. We may write

$$T\frac{4\pi }{\overline{ϵ}^{3/2}(\mathrm{ln}(\delta ^{1/2}\overline{ϵ}))^2}$$
(35)
where $`\delta `$ is a renormalization constant, which we can arbitrarily set to $`\delta <<\overline{ϵ}`$. Physically, we may think of this as correcting for the situation in which the electron and the image charge are both at the interface by limiting the distance of closest approach to atomic dimensions. Equation (35) corrects for the error made by the WKB approximation at the top of the barrier. As a convenient choice of $`\delta `$, we require that it is chosen to limit the closest approach to the origin as $`x=a=1`$Å 333By way of comparison, the typical tunnelling length $`\zeta =\frac{\varphi }{eF}`$ falls in the range of $`10^4`$Å to $`10`$Å for $`1<F<10^3`$ V/$`\mu `$m. so that by the definition of $`\xi `$:

$$\delta =\frac{10^4F}{\varphi }.$$
(36)
Before proceeding further we shall briefly mention the results when we move slightly away from the top of the barrier. In this case the two roots become:

$$\lambda _1\frac{\overline{ϵ}}{2}(1+\sigma )and\lambda _2\frac{\overline{ϵ}}{2}(1\sigma )$$
(37)
where $`\sigma =(14\overline{\alpha }/\overline{ϵ}^2)^{1/2}`$ and in the limit $`\overline{ϵ}`$ large reduces equation (31) to:

$$y\psi ^{\prime \prime }(y)+\gamma ^2\psi (y)=0.$$
(38)
where $`\gamma =\frac{1}{2}\overline{ϵ}(1\sigma ^2)^{1/2}`$. This poses no difficulties and the method outlined above can be adapted to the solution. We summarize the appropriate formulas without further ado:

$$T=\frac{8\overline{\gamma }^2k_0^{1/3}}{|\psi (\xi _0)|^2+4\gamma ^2k_0^{1/3}+k_0^{2/3}|\psi ^{}(\xi _0)|^2}$$
(39)
where:

$$\psi (\xi )=e^{i\pi /4}\sqrt{\frac{\pi }{2}}2\gamma (\xi \overline{ϵ})^{1/2}H_1^{(1)}(2\gamma (\xi \overline{ϵ})^{1/2})$$
(40)
and

$$\psi ^{}(\xi )=e^{i\pi /4}\sqrt{\frac{\pi }{8}}(2\gamma )^2H_0^{(1)}(2\gamma (\xi \overline{ϵ})^{1/2})$$
(41)
In figure 4 we compare the transmission coefficient equation (33) versus the FN formula equation (22) less the prefactors. We see that the FN formula is not only in error as we trespass the barrier giving a non physical value of $`T>1`$, for the region near the top of the barrier (around $`F=700V/\mu \mathrm{m}`$ here), it drops exponentially fast as opposed to the exact solution.

Using equation (39) we can estimate the current near the barrier top which occurs for $`F`$ 700 V/$`\mu `$ in this case, remembering that we have considered $`\varphi `$=1 eV and $`W_a=2\varphi `$. Hence for temperatures much less than $`ϵ_F/k_B`$, the current is obtained from the integral:

$$J=\frac{4\pi me}{h^3}_{W_a}^{ϵ_F}T(W)(ϵ_FW)𝑑W$$
(42)
which has to be evaluated numerically. We can however make an estimate using the approximate expression equation (35). We can set the Fermi level at the top of the barrier and integrate only over a small region of order $`\delta `$ near to it, since outside this region the contributions are exponentially small. We choose $`\delta `$ according to equation (36) as before and thus the $`\overline{ϵ}^{3/2}`$ pre-factor dominates in the transmission coefficient since the logarithmic term is slowly varying in the region of integration. As we take $`W_a/\varphi =2`$ which is large compared with $`\delta `$, we now arrive at the expression for the current density as:

$$J0.1674\varphi ^{3/2}F^3\mathrm{A}\mathrm{cm}^2$$
(43)
Note that the field dependence of the current is in this case $`F^3`$ which differs from the FN prediction of $`F^2`$ and could be experimentally discernable in view of its non-exponential dependence. We can now compare equation (43) with the standard FN formula equation (22) see figure 5.

Having shown that we can have a large (non-exponential) current in spite of a smaller field than in the image-free case, due to proximity to the barrier top, we shall next examine other mechanisms by which the barrier lowering can occur, whose magnitude is comparable with, or might exceed the Schottky image effect. The effect to consider is due to the dielectric screening of the applied field F which we shall take up next.

## 4 The Effects of Screening

Screening is an important property of metals or dielectrics, and the mechanism itself plays a crucial role in the previous section in the form of the image potential which is a screening effect. In recent years there has been a revival of interest in screening,  especially in semiconductors, where surface effects can significantly alter the properties of devices. On the surface the screening length, typically of the order of $`\lambda _{TF}`$ from Thomas-Fermi theory ranges between the order of a few Å for metals to a few $`\mu \mathrm{m}`$ for semiconductors or even higher to tenths of mm due to surface effects since it is the charge carrier density at the surface layer that matters. The (Laplace) screening problem for a flat surface is easily solved for an external field $`F`$. This modifies the external potential $`V(x)`$ in equation (9) by having a non-zero $`V_0`$ such that $`V_0=\lambda _{TF}eF`$. We shall neglect an exponential penetration field inside which in fact modifies the inner potential to

$`V(x)`$ $`=`$ $`W_a\lambda _{TF}eFe^{x/\lambda _{TF}}\mathrm{for}x<0\mathrm{region}1`$ (44)
$`=`$ $`\lambda _{TF}eFeFx\mathrm{for}x>0\mathrm{region}2`$
The latter leads to a 2D electron gas which requires a self-consistent theory to treat since this in turn affects our $`\lambda _{TF}`$ as already mentioned. Further the proper inclusion of screening must now include its effects on the image potentials which pick up Yukawa type terms. It suffices here to estimate the barrier lowering from our expression for $`V_0`$. Since $`F`$ is in $`\mathrm{V}/\mu \mathrm{m}`$ then a $`\lambda _{TF}1\mu \mathrm{m}`$ will completely wipe out a barrier of $`1\mathrm{eV}`$. For 3D this corresponds to a carrier density of order $`10^{25}\mathrm{m}^3`$ but for a 2D layer this can be down by a $`2/3`$ power i.e. $`10^{17}\mathrm{m}^2`$. If we include the effects due to a diminished electron affinity as in the case of diamond-like carbon bonds, then we see that the present theory provides a strong candidate as an explanation for the high current densities in carbon granular samples recently observed in materials for field emission displays , as we shall see in the next section.

## 5 Relation to Experiments

Comparison with experiments such as for carbon field emitters  can be made by assuming that there is a barrier lowering mechanism such as the dielectric screening as mentioned above. The tunnelling current will then be assumed to be given by equation (43) but we shall include a factor $`\sigma `$ to account for granularity or inefficient screening (and hence barrier lowering) which diminishes the current. This modifies equation (43) trivially to:

$$J0.1674\varphi ^{3/2}\sigma ^3F^3\mathrm{A}\mathrm{cm}^2$$
(45)
In figure 6 we have plotted the results which shows that for a granularity factor of $`\sigma =0.15`$ our results compare favorably with experiments  bearing in mind that our results are only estimates and not accurate numerical evaluations444Finite temperature effects should also be included in a more accurate calculation, providing essential improvements to thermionic field emission theory.. In particular the parameters as given in the figures have not been fine tuned in anyway. On the contrary without assuming a significant field enhancement due to pointed tips, the FN theory is out by numerous orders of magnitude. However even with the assumption of a large field enhancement factor $`\gamma =10`$ which is reasonable for ellipsoids , and a barrier lowering $`V_0=0.7\varphi `$, the FN curve figure 6 does not in any way resemble the experimental data, a feature we have found to occur quite commonly in analyzing data from the carbon field emitting samples. Note the vast order of magnitude differences and the fast drop in the slope of the FN curve which is shown more clearly in the inset to figure 6.

These results clearly show that our theory has captured some essential physics associated to tunnelling near the barrier top which may be relevant to the carbon field emitters.

## 6 Conclusion

In conclusion we have presented a closer look at the traditional FN theory of field emission and found that the latter is untenable for calculations of emission currents when significant barrier lowering occurs. A full quantum mechanical solution of the problem provides us with an estimate which, coupled with a proposed barrier lowering due to screening or negative electron affinity, constitutes a reasonable order of magnitude estimate for the high current densities in samples with embedded particles, not available via FN theory. Further work is necessary to provide more quantitative comparisons noting that the granularity of the systems would lead to less total current density although on the other hand some field enhancements due to shape effects  not considered here might partially compensate for this loss.

## Acknowledgements

We wish to acknowledge the support of the Engineering and Physical Sciences Research Council through grants GR/M71404/01 and GR/R97047/01.

## 7 References

## Appendix A Analytical continuation of Bessel Functions

Confusion arises with Bessel functions of order $`1/3`$ due to special care needed to treat the phase factors: some popular mathematical software uses incorrect continuation formulae. Firstly our convention is the same as Watson’s so that the appropriate continuation formula is given by:

$$I_\nu (z)=\{\begin{array}{cc}e^{i\nu \pi /2}J_\nu (ze^{i\pi /2})\hfill & \text{if }\pi <\mathrm{arg}(z)\pi /2\hfill \\ e^{3i\nu \pi /2}J_\nu (ze^{3i\pi /2})\hfill & \text{if }\pi /2<\mathrm{arg}(z)\pi \text{.}\hfill \end{array}$$
(46)
An alternative approach is to perform the continuation for negative $`z`$ as:

$`u_1`$ $`=`$ $`{\displaystyle \frac{\sqrt{\pi }\sqrt{x}J_{1/3}(\frac{2}{3}x^{\frac{3}{2}})}{3}}={\displaystyle \frac{i}{3}}\sqrt{\pi }\sqrt{|x|}J_{1/3}({\displaystyle \frac{2i}{3}}|x|^{\frac{3}{2}})`$ (47)
$`=`$ $`{\displaystyle \frac{i}{3}}e^{i\pi /6}\sqrt{\pi }\sqrt{|x|}I_{1/3}({\displaystyle \frac{2|x|^{\frac{3}{2}}}{3}}),`$
which is incorrect. The negative argument prevents the use of the first form of equation(46). Instead the branch cut properties requires the continuation formula :

$$J_\nu (e^{im\pi /2}z)=e^{im\pi /2}J_\nu (z)$$
(48)
which holds similarly for the $`I_\nu `$ functions. Then:

$$u_1=\frac{\sqrt{\pi }\sqrt{x}J_{1/3}(\frac{2}{3}x^{\frac{3}{2}})}{3}=\frac{i}{3}\sqrt{\pi }\sqrt{|x|}J_{1/3}(\frac{2i}{3}|x|^{\frac{3}{2}})=\frac{1}{3}\sqrt{\pi }\sqrt{|x|}I_{1/3}(\frac{2|x|^{\frac{3}{2}}}{3}),$$
(49)
without any complex phase factors. Similarly $`u_2`$ picks up a minus sign:

$$u_2=\frac{\sqrt{\pi }\sqrt{x}J_{1/3}(\frac{2}{3}x^{\frac{3}{2}})}{3}=\frac{i}{3}\sqrt{\pi }\sqrt{|x|}J_{1/3}(\frac{2i}{3}|x|^{\frac{3}{2}})=\frac{1}{3}\sqrt{\pi }\sqrt{|x|}I_{1/3}(\frac{2|x|^{\frac{3}{2}}}{3}).$$
(50)
Hence one has to manually continue the $`J_\nu `$ to the $`I_\nu `$ functions for negative arguments when if using the incorrect continuation, equation (47). We have found that Mathematica v4.1 uses equation (47).

## Appendix B Power series solution for the full image problem I

The generic form of the full image problem is given by the following second order differential equation

$$\psi ^{\prime \prime }+\left(ϵ+\alpha x+\frac{\beta }{x}\right)\psi =0.$$
(51)
As before a straightforward linear transformation: $`\xi =(x+\frac{ϵ}{\alpha })\alpha ^{\frac{1}{3}}`$ leads to the following form:

$$\psi ^{\prime \prime }+\xi \psi +\frac{\overline{\alpha }}{\xi \overline{ϵ}}\psi =0,$$
(52)
where $`\overline{\alpha }=\beta \alpha ^{\frac{1}{3}}`$, and $`\overline{ϵ}=ϵ\alpha ^{2/3}`$. Clearly the image free case $`\overline{\alpha }=0`$ reduces to the usual Airy differential equation. Equation (52) can be rewritten as:

$$(\xi \overline{ϵ})\psi ^{\prime \prime }+\xi (\xi \overline{ϵ})\psi +\overline{\alpha }\psi =0.$$
(53)
Unfortunately this equation does not belong to the standard hypergeometric class and cannot therefore be treated by the usual mathematical physics functions. However using the standard Frobenius method one can obtain the following power series solutions. Let

$$u_m(\xi )=\underset{n=0}{\overset{\mathrm{}}{}}a_n\xi ^{m+n},$$
(54)
then equation (54) leads to the following 5 term recursion relation for the determination of the coefficients $`a_n`$:

$`a_{n1}(m+n1)(m+n2)`$ $``$ $`a_n\overline{ϵ}(m+n)(m+n1)`$ (55)
$`+`$ $`a_{n4}\overline{ϵ}a_{n3}+\overline{\alpha }a_{n2}=0.`$
All the coefficients are well determined in terms of the coefficient $`a_0`$. The first term satisfies

$$\overline{ϵ}m(m1)a_0=0,$$
(56)
giving $`m=0orm=1`$ which shows as expected that we shall have two linearly independent solutions. Carrying on we have $`a_1=0`$ and

$$a_2=\frac{\overline{\alpha }}{\overline{ϵ}}\frac{a_0}{(m+2)(m+1)}.$$
(57)
Notice that $`a_1=0`$ regardless but $`a_2`$ is only finite for $`\overline{\alpha }0,`$ i.e. $`\beta 0`$. The following coefficients can be sequentially determined:

$`a_3`$ $`=`$ $`{\displaystyle \frac{a_2(m+1)}{\overline{ϵ}(m+3)}}{\displaystyle \frac{a_0}{(m+3)(m+2)}},`$ (58)
$`a_4`$ $`=`$ $`{\displaystyle \frac{a_3(m+2)}{\overline{ϵ}(m+4)}}+{\displaystyle \frac{a_0}{\overline{ϵ}(m+4)(m+3)}}+{\displaystyle \frac{\overline{\alpha }}{\overline{ϵ}}}{\displaystyle \frac{a_2}{(m+4)(m+3)}},`$ (59)
and so on. The general 5 term recursion relation can be re-written as

$`a_n`$ $`=`$ $`{\displaystyle \frac{1}{(m+n)(m+n1)}}`$ (60)
$`\left[{\displaystyle \frac{1}{\overline{ϵ}}}a_{n1}(m+n1)(m+n2)+{\displaystyle \frac{\overline{\alpha }}{\overline{ϵ}}}a_{n2}a_{n3}+{\displaystyle \frac{1}{\overline{ϵ}}}a_{n4}\right].`$
Unfortunately there is no way to simplify this any further. When $`\overline{\alpha }=0`$ (image free case) the above in fact reduces to a two term relation given by:

$$a_n=\frac{a_{n3}}{(m+n)(m+n1)}.$$
(61)
This immediately leads us to the two independent solutions, for $`m=0`$

$`u_0(\xi )`$ $`=`$ $`1+{\displaystyle \underset{n=1}{\overset{\mathrm{}}{}}}{\displaystyle \frac{(1)^n\mathrm{1.4.7}\mathrm{}(3n2)}{(3n)!}}\xi ^{3n}`$ (62)
$`=`$ $`{\displaystyle \frac{\mathrm{\Gamma }(2/3)}{3^{1/3}}}\xi ^{1/2}J_{1/3}({\displaystyle \frac{2}{3}}\xi ^{3/2}),`$
and for $`m=1`$

$`u_1(\xi )`$ $`=`$ $`\xi +{\displaystyle \underset{n=1}{\overset{\mathrm{}}{}}}{\displaystyle \frac{(1)^n\mathrm{2.5.8}\mathrm{}(3n1)}{(3n+1)!}}\xi ^{3n+1}`$ (63)
$`=`$ $`3^{1/3}\mathrm{\Gamma }(4/3)\xi ^{1/2}J_{1/3}({\displaystyle \frac{2}{3}}\xi ^{3/2}),`$
which are the series expansions for the Bessel functions of order $`1/3`$. The oscillatory nature of these functions means that the convergence is in general poor and some more sophisticated re-summation method must be used to obtain in particular the asymptotic properties of these functions. In the next appendix we shall show just one method as to how this can be achieved.

## Appendix C Power series solution for the full image problem II

An alternative solution for the full image problem which is an effective re-summation of the previous power series solution that will allow us to perform numerical calculations, will now be shown. The approach is to first use a Laplace transform on equation (51). Let

$$\psi (\xi )=_Ce^{\xi t}f(t)𝑑t,$$
(64)
over some contour $`C`$, then the resulting differential equation in $`\psi `$ becomes converted into a second-order differential equation in $`f`$:

$$f^{\prime \prime }(t^2\overline{ϵ})f^{}+(\overline{\alpha }\overline{ϵ}t^22t)f=0.$$
(65)
It is easy to show that the solution for $`\overline{\alpha }=0`$ reduces to the Airy integral where $`f(t)=e^{t^3/3}`$ so we shall factorize $`f`$ in the form: $`f=ug`$ where $`u=e^{t^3/3}`$. The equation for $`g`$ now simplifies tremendously:

$$g^{\prime \prime }+(t^2+\overline{ϵ})g^{}+\overline{\alpha }g=0.$$
(66)
Unfortunately there is still no closed form solution for this equation although it can be proved that solutions for $`g`$ must be an entire function, by which we can develop an expansion as an infinite series:

$$g=\underset{n=0}{\overset{\mathrm{}}{}}a_nt^{m+n}.$$
(67)
The boundary condition requires that when $`\overline{\alpha }0`$ then $`g1`$. With this we can show that $`m=0`$ is the only allowed index. Then we have the following 4 term recursion relation:

$$(n1)(na_n+\overline{ϵ}a_{n1})+(n3)a_{n3}+\overline{\alpha }a_{n2}=0.$$
(68)
The first few solutions are easily obtained: $`a_0=1`$, $`a_2=\frac{\overline{\alpha }}{2}`$$`a_3=\frac{\overline{ϵ}\overline{\alpha }}{6}`$$`a_4=\frac{\overline{\alpha }}{24}(\overline{\alpha }\overline{ϵ}^2)`$ and so on. As required only the first term survives in the limit $`\overline{\alpha }=0`$. Interestingly the 4 term relation equation (68) can actually be reduced to a three term one. Let $`\gamma _n=\frac{a_n}{a_{n1}}`$ then

$$(n1)(n\gamma _n+\overline{ϵ})\gamma _{n1}+\frac{(n3)}{\gamma _{n2}}+\overline{\alpha }=0,$$
(69)
which may be more convenient for computation. Finally we have now effectively re-summed the series to obtain

$$\psi (\xi )=\underset{n=0}{\overset{\mathrm{}}{}}a_n_Ct^ne^{\xi t+\frac{1}{3}t^3}𝑑t.$$
(70)
We now see that

$`{\displaystyle _C}t^ne^{\xi t+\frac{1}{3}t^3}𝑑t`$ $`=`$ $`{\displaystyle \frac{d^n}{d\xi ^n}}{\displaystyle _C}e^{\xi t+\frac{1}{3}t^3}𝑑t`$ (71)
$`=`$ $`{\displaystyle \frac{d^n}{d\xi ^n}}\left[\xi ^{1/2}Z_{\pm 1/3}({\displaystyle \frac{2}{3}}\xi ^{3/2})\right],`$
the last term following from the equivalence between the Airy functions and the Bessel functions of order $`1/3`$. In view of the recursion properties for the derivatives of the Bessel functions, we have now converted the solution effectively into an infinite series of Bessel functions. Only the first term survives in the image free case which is now trivial and all the oscillatory bits are now absorbed into the Bessel functions.


format: text

text:
# Why Hot Horizontal Branch Stars Can Appear Redder Than Red Giants Based on observations with the ESO/MPI 2.2m telescope, located at La Silla Observatory (Chile) and on observations with the NASA/ESA Hubble Space Telescope.

## 1 Introduction

Color-magnitude diagrams (CMD) of resolved stellar populations are the fundamental tool for testing predictions of stellar evolutionary theory and the capability of stellar models to finely reproduce the observations. Photometric studies are no longer limited to the classical $`BVI`$ bands, as high quality near-infrared and ultraviolet CMDs are becoming more and more available. Ultraviolet (UV) bands certainly cover the ideal wavelength region for the study of hot stellar populations. Indeed, UV observations help to remove the degeneracy between color and temperature for hot blue horizontal branch (HB) stars, allowing a more reliable measurement of atmospheric parameters (in particular temperature estimates).

Recently, the blue horizontal branch (BHB) population has been the subject of specific UV studies, and some interesting features were revealed only in ultraviolet CMDs: jumps (Grundahl et al. grun99 (1999), G99, Momany et al. moma02 (2002)), and peculiar horizontal branch extensions, like the blue-hook (D’Cruz et al. cruz00 (2000), Brown et al. brow01 (2001)). The ultimate goal of these studies is a better understanding of the physical mechanisms governing core-helium burning stars, and driving the horizontal branch morphology. To this end, our group has undertaken a long-term project to obtain wide-field multi-wavelength data of BHB globular clusters (GCs). Mainly, we used the WFPC2 on board of HST to map the central core of our target clusters (Piotto et al. piotto02 (2002)), and the wide field imager (WFI) at the ESO/MPI 2.2m for the coverage of the outer regions. Unfortunately, observations carried out through non-standard photometric systems sometimes generate puzzling and apparently un-physical properties in the observed CMDs.

In this paper we report on a curious feature in the observed HB morphology in the UV color-magnitude diagrams, and an explanation for it. Basically, in some of the $`V`$, ($`UB`$) diagrams of our database, we find that the blue HB stars hotter than the instability strip and cooler than the G99 $`u`$–jump ($`0.1BV0.2`$) have $`UB`$ colors that place them on the right hand side of the red giant branch (RGB).

According to stellar evolutionary theory blue HB stars should not “appear” cooler than their RGB progenitor. Indeed, HB stars have hotter effective temperature than RGB stars with the same metallicity. In a color-magnitude diagram, we do not expect any red incursion of BHB stars to the cooler side of the RGB (which almost corresponds to the Hayashi track): for stars in hydrostatic equilibrium and having convection, the Hayashi line separates a permitted region (on its left) from a forbidden one (on its right). Despite this theoretical constraint, we will show that the unexpected finding of the BHB red incursion is present in at least 6 $`V`$, ($`UB`$) CMDs of BHB clusters we have observed. We will also show that this peculiarity is not present in all the CMDs from broad band photometry of BHB clusters, as the effect depends on the metallicity.

As discussed in Section 3, it turns out that the apparent anomaly rises when “non-standard” $`U`$ filters are used. We decided to publish these color-magnitude diagrams, and discuss the origin of the apparently un-physical HB morphology, in order to call the attention on the importance of using appropriate filter sets when acquiring photometric data. Alternatively, and for a proper interpretation of the empirical data, this paper shows the importance of using exactly the photometric passbands of the observational system when transforming the theoretical tracks to the observational plane.

## 2 Observations and Data Reductions

We present data collected with the Hubble Space Telescope (HST), and the ESO/MPI 2.2m telescope. $`UBV`$ data of NGC7089 and NGC5139 were taken on June 2002 with the Wide-Field-Imager (WFI) at the ESO/MPI 2.2m, employing the “new” $`U`$ filter (Germany germ02 (2002)). A second $`UBV`$ data set of NGC7078 and NGC7099 images was obtained on July 2000, employing the “old” $`U`$ filter at the same telescope.

Basic reduction of the CCD mosaic was performed using the IRAF package MSCRED (Valdes 1998), while stellar photometry was done using DAOPHOT/ALLFRAME (Stetson 1994). Calibration to the standard $`UBV`$ system was obtained using a set of standard stars from Landolt (1992). HST data of NGC7099, NGC7078, and NGC6205 were taken in the F$`336`$W, F$`439`$W, and F$`555`$W bands. For these clusters we present only instrumental diagrams. Besides the HST and 2.2m data, we present $`UBV`$ published data of NGC1904 by Kravtsov et al. (krav97 (1997)), and NGC1841 by Alcaino et al. (alca96 (1996)), both collected at the ESO/NTT telescope. Lastly, we present composite $`UBV`$ catalogs of NGC4833 and NGC4590. The $`BV`$ photometry of these clusters were taken from Melbourne et al. (mel00 (2000)) and Walker (walker94 (1994)), while the $`U`$ photometry were obtained employing the WFI “new” $`U`$ filter at the 2.2m telescope.

## 3 The Blue Horizontal Branch Red Incursion

Figure 1 shows a zoom around the HB in the $`V`$, ($`BV`$) and $`V`$, ($`UB`$) diagrams of the 8 GCs, obtained with the 4 different telescopes (HST, NTT, CTIO and 2.2m ESO telescope). The clusters span a range in metallicity of $`0.7`$ dex, from the metal-poor regime to intermediate metallicity. For all clusters in the sample, the presence of the HB elbow and its incursion through the RGB is quite evident. Due to the larger sample of HB stars, it is really impressive in NGC7089, NGC7078, and NGC5139111$`I`$, ($`UB`$) diagrams of these clusters show even a greater extent of the red incursion. Note how the blue HB of NGC7099 lies almost entirely on the red side of the RGB in the $`V`$, ($`UB`$) diagram.

Figure 2 shows some of the conflicting results we faced in the course of our investigation. The $`V`$, ($`UB`$) diagrams of NGC7099, from data collected at two different telescopes, show completely different HB morphologies. In the case of the diagram obtained with HST, the BHB incursion is clearly visible, while in the other diagram (based on data collected with the 2.2m telescope) the BHB does not cross the cluster RGB at all. These differences are not to be attributed to the limit of the UV region reachable from the ground (at variance with the HST space observations), as Fig. 1 shows two CMDs (for NGC7089 and NGC5139), from data obtained at the ESO/MPI 2.2m telescope, clearly disclosing the HB red incursion. Neither the effect can be attributed to red leakage in the employed $`U`$ filters (the WFPC2 F336W and WFI new-$`U`$ filters are known to be affected by red-leaks). Red leakage is affecting mainly RGB stars, causing them to have smaller than expected $`U`$ magnitudes and therefore bluer $`UB`$ colors. Indeed, considering the F336W filter (the most affected by red leakage among our filters) we calculate that for an RGB star at the level of the HB, with \[Fe/H\]$`=1.8`$, the red leakage in the F336W band causes a blue-ward shift of only $`0.012`$ magnitude while the effect we discuss in this paper is of the order of few tenths of a magnitude.

The apparent contradictory results shown in Fig. 2 lead us to search for differences among the filters adopted to acquire the different data sets, in particular at the ESO/MPI 2.2m telescope. We found that the 2.2m data showing the BHB incursion were obtained using the two “new” $`U`$ and $`B`$ filters (Germany germ02 (2002)), while data not showing the incursion were obtained with the “old” $`U`$ and $`B`$ filters.

Prompted by this result, we have analyzed the normalized transmission curves of a variety of $`U`$ and $`B`$ filters (Fig. 3). We found that the “new” $`U`$ filter at the 2.2m, as well as the HST F336W, and the Strömgren $`u`$ filters, have a transmission curve with a peak around $`\lambda 3400\AA `$. On the other hand, the “old” $`U`$ filter at the 2.2m, and the $`U`$Buser band (Buser buser78 (1978)) of the synthetic colors (Castelli cast99 (1999)) that were originally used for the transformations from the theoretical to the observational plane, have a transmission curve peaking around $`\lambda 3650\AA `$. It is worth mentioning that the latter group of filters have transmission curves that are closer to the standard $`U`$ Johnson (although this does not necessarily mean that the observational band, which is a convolution of the filter, CCD, telescope optics, and atmosphere transmission curves is exactly the standard $`U`$ Johnson (john55 (1955)) band). Interestingly enough, the clusters observed with the latter group of filters do not show the BHB incursion, whereas clusters observed through filters whose transmission peak at $`\lambda 3400\AA `$ do. Also the $`B`$ filter can be important. For example, clusters observed with the new 2.2m $`U`$ filter and the old 2.2m $`B`$ filter (see Fig. 3) do not show the incursion. Clearly, it is the combination of non-standard $`U`$ and $`B`$ passband that creates the observed anomaly.

In Fig. 3 we show the normalized transmission curve of the $`UB`$ filters used at the 2.2m, the HST F336W and F439W filters, and the Strömgren $`u`$ and $`v`$ filters. For two different metallicities (\[Fe/H\]$`=1.8`$ and $`0.7`$), we also show synthetic spectra (Castelli & Kurucz cast03 (2003)) for BHB stars, and for RGB stars with surface gravity and effective temperature appropriate of red giants at the level of the HB luminosity. When comparing a BHB spectrum with the transmission curves of different $`U`$ and $`B`$ filters, it is clear that the occurrence of the BHB red incursion is related to the Balmer discontinuity: $`U`$ filters peaking around $`\lambda 3400\AA `$ do not encompass the Balmer jump. These filters collect less $`U`$ flux with respect to $`U`$ filters peaking at $`\lambda 3650\AA `$. Consequently, the $`U`$ magnitudes of BHB stars are fainter and therefore their ($`UB`$) colors will be redder. On the other hand, RGB stars have a shallower Balmer jump, hence these are less affected by differences in $`U`$ filters. Also important is the broadness of the new-$`B`$ filter at the 2.2m with respect to the old-$`B`$ filter. The blue side of the new-$`B`$ encompasses the Balmer jump, and is the cause of the anomalous morphology of the HB when the corresponding $`B`$ magnitudes are combined with new-$`U`$ filter at the 2.2m.

A deeper understanding of the BHB red incursion can be achieved when considering the dependence of the ($`UB`$) color on stellar parameters such as surface gravity, effective temperature, and chemical composition. The $`UBV`$ synthetic photometry from the ATLAS9 models (Castelli cast99 (1999), Fig. 11) shows that, for a given gravity, the ($`UB`$) index does not monotonically increase with decreasing temperature (see also Fig. 4). Instead, between $`67505500`$K, the ($`UB`$) index shows a minimum which increases with decreasing metallicity, and is almost independent from gravity. Moving towards higher temperatures ($`7000`$K), the ($`UB`$) index increases from the quoted minimum and reaches a maximum value, whose intensity depends on gravity, but not on metallicity. The lower is the gravity, the higher is the maximum. Moreover, the ($`UB`$) maximum moves from $`7500`$K for log $`g=2.0`$ to about $`9000`$K for log $`g=5.0`$. This means that in a metal-poor cluster, the ($`UB`$) index of a star with $`T_{\mathrm{eff}}8000`$K and log $`g=3.0`$ (i.e. an HB star) may be close to the maximum value and therefore exceeds the ($`UB`$) index of a RGB star with $`T_{\mathrm{eff}}55006000`$K.

This anomalous ($`UB`$) behavior is not encountered when the “standard” $`U`$-Buser filter or the old-$`U`$ filter at the 2.2m telescope are used, simply because the difference between the maximum and the minimum in the $`T_{\mathrm{eff}}`$($`UB`$) curve is not large enough to produce the red incursion. On the other hand, the difference between the minimum and maximum ($`UB`$) index becomes larger for $`U`$ filters with transmission curves shifted towards shorter wavelengths, to the point of producing the red incursion for the F336W filter on HST, Strömgren $`u`$ filter, and the new $`U`$ filter at the 2.2m telescope.

To better explain the red incursion observed for low gravity stars, we replaced in the synthetic photometry (Castelli, cast99 (1999)) the $`U`$ filter from Buser (buser78 (1978)) and the $`B`$ filter from Ažusienis & Straižys (azus69 (1969)) with the new $`U`$ and $`B`$ filters of the 2.2m photometric system. Figure 4 shows that, for typical gravities of HB stars (log $`g=3.0`$), there is a large difference in the $`T_{\mathrm{eff}}`$($`UB`$) curves for the two sets of filters; particularly at the ($`UB`$) maximum, where the difference is $`0.2`$ mag. Figure 4 also shows that while the difference between the ($`UB`$) index of an HB star (log $`g=3.0`$) and an RGB star (log $`g=2.0`$) for the classical filters is approximately the same for all the temperatures, for the new $`UB`$ filters at the 2.2m the ($`UB`$) color of an HB star is redder than the ($`UB`$) of an RGB star.

The effect we have just discussed is not new. It has been known as the “ultraviolet deficiency” of cool, BHB stars in the ($`UB`$, $`BV`$) two-color diagram (see Markov et al. mark01 (2001) and references therein). Analyzing data for 6 GCs, Markov et al. compared the position of BHB stars with respect to i) a Population I sequence; and ii) a zero-age horizontal branch (ZAHB) two-color line, and quantified the $`\mathrm{\Delta }(UB)`$ deficiency. In analogy with our results, they also found that this observable is not common among all BHB clusters. Markov et al. argued that there exists a connection between this observable and cluster parameters. In particular, they showed how a satisfactory overlapping of the observed and theoretical sequences is achieved when using low-gravity (log $`g2.0`$) atmospheric models.

As expected from the above discussion, using a $`U`$ filter-set showing the BHB red incursion (e.g. the F336W filter on HST or the “new” $`U`$ filter at the 2.2m telescope) not all BHB clusters show the red incursion. It is clearly present in a metal-poor cluster like NGC7078, whereas completely missing in an intermediate metallicity cluster like NGC6205 (see also Fig. 6). Due to the small sample of clusters, we can neither observationally quantify the dependence of the effect on the cluster metal content nor precisely estimate the cutoff metallicity (if any) beyond which no incursion is observed. The fact that the HB of NGC1904 (\[Fe/H\]$`1.6`$, $`E_{BV}=0.01`$) shows a red incursion, while NGC6205 (\[Fe/H\]$`1.5`$, $`E_{BV}=0.02`$) does not show any incursion could suggest that this effect might be present mainly in clusters with metallicities lower than \[Fe/H\]$`1.5`$. Inspection of the HST F555W, (F336W$``$F439W) diagram of NGC6388 (\[Fe/H\]$`0.60`$, showing no indication of the incursion) suggests that the effect is not present in the CMDs from broad band photometry of high metallicity BHB clusters. In view of the strong dependence of the effect on the filter transmission curve, only a photometrically homogeneous observational data set, covering a large metallicity range will allow us to characterize the phenomenon, and investigate whether this feature can be used to constrain other cluster parameters, like metallicity, distance, and reddening.

The suggestion for a metallicity threshold of the incursion relies on a set of images collected with filters which are supposed to reproduce the standard $`UBV`$ Johnson bands at 3 widely used telescopes: HST, NTT and the 2.2m. It must be noted, as pointed out by the referee, that the red incursion effect is well visible also in other bands, like the Strömgren intermediate ones. However, the metallicity threshold is not present in all the photometric systems. Indeed, theoretical isochrones and ZAHBs (from the Cassisi & Salaris (cassisi97 (1997)) library) show the incursion in the Strömgren MV, ($`uv`$) plane even at metallicities of $`Z=0.006`$, representative of metal rich GCs.

Figure 5 helps to understand the metallicity effect of the BHB red incursion. For the HST flight and the Strömgren photometric systems we have plotted the following three color differences as a function of metallicity: (1) $`\mathrm{\Delta }(HBRGB)`$ measures the extent of the incursion, i.e. the difference in color between the reddest HB point and the RGB color at the level of the reddest HB point; (2) $`\mathrm{\Delta }(HBHB_{Z=0.0001})`$ is the difference in color between the reddest HB point for the quoted metallicities and that for $`Z=0.0001`$; and lastly, (3) $`\mathrm{\Delta }(RGBRGB_{Z=0.0001})`$ is the difference between the RGB colors (at the level of the reddest HB point for the three indicated metallicities) and that for $`Z=0.0001`$. There are two important conclusions that can be drawn from Fig. 5. First of all, we confirm that, while the red incursion of the HB is present at least up to metallicities of $`Z=0.006`$ in the Strömgren photometric system, the phenomenon is expected to disappear for clusters more metal rich than $`Z=0.001`$ in the HST flight system, where the corresponding $`\mathrm{\Delta }(RGBHB)`$ becomes negative. The second important evidence disclosed by Fig. 5 is that the extension of the red incursion (and its disappearance for more metal rich clusters in the broad band CMDs) is due to the fact that the RGB becomes redder and redder at increasing metallicity, while the color of the reddest HB point remains more less the same at varying the cluster metallicity.

In Fig. 6, we also compare the observed WFPC2/HST F555W, (F336W$``$F439W) CMDs of NGC7078 and NGC6205 with the isochrones from Cassisi & Salaris (cassisi97 (1997)) transformed into the HST flight system. The data have not been calibrated, but the calibration to the HST flight system is simply a zero point translation, and does not affect the overall morphology of the CMD. It is worth noticing that the models reproduce very well the BHB “elbow” morphology, and, most importantly, reproduce the BHB red incursion when present. We have also compared the models with the CMDs from ground-based photometry for those clusters showing the BHB red incursion. In this case the models are not able to reproduce such a feature at all. This is a consequence of the fact that the $`U`$-filter used to derive the color-$`T_{\mathrm{eff}}`$ relations (adopted for transforming the models from the theoretical to the empirical plane) is the $`U`$-Buser filter, clearly inadequate to fit data obtained with “non-standard” $`U`$ filters (cf Fig. 3).

This occurrence strongly reinforces the need for appropriate color-$`T_{\mathrm{eff}}`$ relations, based on exactly the same filter passband used in the observations, in order to realize a meaningful comparison between theory and observations. Most importantly, the anomalous CMD morphology (produced by $`U`$ filters whose transmission is far from the standard ones) calls attention to another problem, often ignored by those in charge of the development and maintenance of imaging instruments: filter passbands must be as close as possible to the bands they pretend to reproduce. This simple rule makes the calibration to the standard photometric system more accurate, and avoids the production of apparently un-physical features in the observed CMDs.

## 4 Summary and final remarks

In the present work, we show the existence of a curious feature in the $`V`$, ($`UB`$) CMDs of some GCs. Quite surprisingly, in some observed CMDs, hot HB stars ($`6000T_{\mathrm{eff}}`$(K) $`10000`$) become redder than the much cooler RGB. An investigation of this effect has shown the following:

* the size and the presence itself of this feature is strongly dependent on the passband of the adopted $`U`$ filter. In particular, $`U`$ filters which do not encompass the Balmer jump break, do not show the BHB red incursion;
* also the $`B`$-filter passband concur to cause the anomaly;
* the occurrence of this feature is due to the peculiar dependence of the $`U`$ and $`B`$ filters on star effective temperature, gravity, and metallicity;
* clusters more metal-rich than \[Fe/H\]$`1.5`$ do not show any BHB red incursion in the broad band Johnson and HST flight system, while there is no metallicity threshold in Strömgren $`uvby`$ photometry for metallicities lower than $`0.006`$;
* evolutionary models using the correct filter passband, when transformed to the observational plane, finely reproduce the observed HB morphology.

We report this anomalous HB morphology to point out the importance of having observational photometric systems which accurately reproduce the standard ones. Unfortunately, the fact that the anomalous CMDs shown in this paper have been obtained from data acquired at some of the most widely used instruments for astronomical imaging shows how this problem is yet to be addressed.

###### Acknowledgements.

We warmly thank the referee, Frank Grundahl, for the useful comments which helped to improve the paper. In particular, we thank him for pointing out the presence of the effect described in this paper also in the CMDs from Strömgren photometry. We also thank G. Bono and V. Castellani for useful discussions on this matter. Y.M thanks S. Lucatello for useful discussions on stellar spectra. L.R.B., Y.M., S.O., G.P., and A.R.B. acknowledge the support of the MIUR (Cofin2001). S.C acknowledges financial support by MIUR (Cofin2002).


text:
# The CMS Integration Grid Testbed

## I Introduction

The CMS Integration Grid Testbed (IGT) was commissioned by USCMS in the Fall of 2002 in order to provide a stable platform for integration testing of Grid middleware with the existing CMS Monte Carlo production environment. CMS (Compact Muon Solenoid) is a high energy physics detector planned for the Large Hadron Collider (LHC) at CERN, just outside of Geneva, Switzerland. While CMS will not begin taking data until after 2007, hundreds of physicists around the world are taking part in Monte Carlo simulation studies of the detector and its potential for discovering physics.

The IGT is intended to comprise a small number of nodes at USCMS Tier-1 and Tier-2 centers from hardware that is assigned to USCMS Facilities. Within the IGT, focused development and integration prototyping takes place in support of a set of middleware products and CMS applications to be used in a production grid setting. Currently, the IGT actually comprises most of the available USCMS Tier-1 and Tier-2 resources. This is because

* This was the first time we have tried seriously to put Grid middleware into production on this scale in a continuously available fashion.
* There was as yet no documentation or level of support in place appropriate to a production setting

The goals of this first running of the IGT are to prepare for a turnover of newly integrated software onto the first Production Grid and to start writing documentation appropriate to a production setting, while participating in official CMS Monte Carlo production. The turnover is scheduled to happen in Spring 2003. At that time, the resources assigned to the IGT will mostly migrate to the production grid (PG).111From time to time, it will happen that large amounts of resources will be temporarily assigned to the IGT for scalability tests in the integration setting. In addition, CMS has a development grid testbed (DGT) in order to provide a less stringent environment for more speculative development and test deployments. With these three grids in place, USCMS is well positioned to deploy a Grid with a view to continuous upgrades and improved software quality control while maintaining a development environment conducive to new concepts and experiments in Grid middleware. This also provides an effective response to the declared intention of the LHC Computing Grid (LCG) to produce a continuously available 24x7 production grid for the LHC experiments.

## II Hardware Included in the IGT

Table 1 shows the available hardware on the IGT. Participating sites included the California Institute of Technology, Fermi National Accelerator Laboratory, University of California at San Diego, and the University of Florida at Gainesville as of November 1, 2002. A group from the LCG at CERN also joined the effort after November 15, 2002. While the University of Wisconsin at Madison Computer Science department did not participate directly in the IGT, they participated directly in the USCMS Development Grid Testbed222The DGT is a smaller clone of the IGT on which it is permitted to deploy untested software.(DGT) in support of the IGT efforts by troubleshooting problems and by doing important regression tests on the middleware.

In Spring 2003, most of this hardware will be turned over to a production grid. A small number of machines will be saved to do integration testing in support of production grid operations, which is the original purpose of the IGT in the first place.

## III Software Running on the IGT

In the Grid environment, many different layers of software have to be present. At the lowest level is the OS itself. The IGT ran a Linux platform with either RedHat 6.X based operating systems or with RedHat 7.X based operating systems333The use of RedHat 6.X based operating systems was required for CMS production with applications based on CMS ORCA which used Objectivity as an object persistency layer. This was not because of limitations within Objectivity, but was rather due to licensing issues.. The Grid middleware was distributed using the Pacman bb:pacman  based Virtual Data Toolkit (VDT) bb:vdt . The version which ran on the IGT was VDT 1.1.3. The next level of software included the Job Creation level, consisting of mainly CMS developed tools and some PPDG provided tools to wrap CMS specific jobs in Grid aware wrappers, and finally the CMS applications themselves.

Most Grid sites ran the PBS batch system or the Condor High Throughput Computing System, the Farm Batch System Next Generation (FBSNG) was run at Fermilab.bb:fbsng

### III.1 VDT Middleware

The VDT is composed of three types of grid software, although the lines between these types are sometimes blurry. These three types are:

* *Core Grid Software* This is the grid middleware, and it includes Globus, bb:globus  Condor-G, and Condor bb:condor .
* *Virtual Data Software* This is the software that is able to either compute or fetch data on demand, depending on whether it needs to be computed or is already available.bb:virtualdata
* *Utilities* This is a selection of software that provides smaller but still important utilities, such as as software for fetching certification authority revocation lists on a regular schedule.

The IGT made use primarily of the core Grid software and some of the utilities.

The software in Table 2 was included with VDT version 1.1.3, which ran on the IGT.

### III.2 CMS Specific Software

The CMS software was distributed to the remote sites using DAR bb:dar  before any jobs were submitted. In principle, the CMS software can be installed as part of the job, as described below in the discussion of MOP. However, the particular production request handled by the IGT was large enough (greater than 100,000 CPU hours on a single GHz CPU) to justify special pre-installation.

CMS Monte Carlo production consists of several steps bb:cms-prod-note . First, vector representations of simulated physics collision events are generated with the Pythia-based CMKIN application. Second, the responses of the CMS detector are simulated in the GEANT 3 based CMSIM application. A third CMS specific step is to re-format the CMSIM events into an Objectivity DB format with the writeHits application444This step will be combined with the previous step in future versions of the software.. Forth, the native detector signals must be mixed with noise simulations and with simulated by-products of nearly contemporaneous collisions called pileup (PU) in the writeDigis application. Pileup events are typically pre-processed and stored in locally resident files ready to be mixed with signal events in writeDigis, leading to I/O bound behavior when the number of pileup events per signal event is large. (For the $`10^{34}`$ luminosity expected in the LHC, this ratio is about 200:1 on average.) Finally, the last stage of production involves the creation of analysis objects (ntuples) to be analyzed by the physicists. Table 3 summarizes the average behaviors of the executables used in CMS Monte Carlo production. Actual production depends most critically on the size of each event at the CMKIN stage: the more particles produced per event translated directly into higher processing times.

CMS production usually proceeds by breaking up production requests into sets of 250 events each555250 events was found to be a reasonably large round number of events that lead to manageable filesizes and processing times. and processing each collection serially through all steps. For the IGT production during Fall 2002, there were two requests for events. The first request was for 1M events processed through all steps. The second request was for 500K events processed only through the CMSIM stage.

#### III.2.1 Job Creation: McRunjob

McRunjob bb:mcrunjob4  is a package containing Python scripts designed to aid in organizing large scale CMS and DZero production processing applications. While initially developed specifically for DZero offline executables, McRunjob was easily extended to the CMS experiment. There are two basic classes within McRunjob. The Configurator encapsulates all of the knowledge needed to run an application or perform some simple task and exposes only the metadata with a customizable interface. The Linker is a container class for Configurators responsible for instantiation, maintaining a list of Configurators, and handling communication among the Configurators and the user and other Configurators. The user interacts with the Linker, sending commands to attach and configure various Configurators. Conceptually, McRunjob provides a language by which the use can specify a workflow pattern abstractly and then the McRunjob implementation takes care to turn it into a set of submittable jobs. Depending on which modules are included, McRunjob can target Virtual Data Language of the Chimera system or jobs using the SAM system at DZero. For the IGT, McRunjob comes with a set of Configurators that target the CMS legacy Impala runtime scripts which were used in official Spring 2002 Monte Carlo production. McRunjob produced 4000 jobs for the million event full production request, while the balance of 2000 CMSIM-only jobs were produced by the Impala scripts.

### III.3 MOP

The jobs as produced by McRunjob and Impala for the IGT run were not specially “grid aware.” The connection to the grid was provided by a thin software layer called MOP (Monte Carlo Production). MOP basically represents existing McRunjob or Impala produced jobs as Directed Acyclic Graphs (DAGs). There were four generic types of DAG nodes that help accomplish this. Stage-in nodes were responsible for transporting the execution environment to the worker node. Run nodes were responsible for running the executables. Stage-out nodes were responsible for transporting results back to the submit site. Finally, clean-up nodes were responsible for removing any left over job state from the worker nodes.

During the IGT running, MOP was invoked to create DAG representations of each job at job submit time. Once a DAG was produced, MOP submitted the DAG to the DAGMan package of Condor. DAGMan usually runs DAG nodes using the Condor system; for the IGT Condor-G was used as a gateway to allow DAGMan to run DAG nodes on remote sites running Globus Job Managers. In turn, these Globus Job Managers are able to run the jobs using local batch queues. For the IGT, only Condor and the Farm Batch System of Fermilab were used as local batch interfaces.

Scheduling functionality was not implemented in the MOP system during the Fall 2002 IGT run. Rather, jobs were distributed by direct operator specification at job submission time. MOP was logically divided into the MOP master site and the MOP worker sites. Jobs were created and submitted from the MOP master site, all input files were staged in from the MOP master site, and all output was returned to the MOP master site. No replica catalogs were used during the production process itself, but resulting data products were registered in a replica catalog at the end of processing. During Fall 2002, Fermilab hosted the IGT MOP master site, while UW Madison hosted a MOP master site for the DGT.

### III.4 Virtual Organization

The Virtual Organization (VO) was implemented using a system to help us to organize the Grid mapfiles of Globus. We used the Caltech Virtual Organization Group Manager bb:vo  developed at Caltech to manage the VO users on the testbed. Group Manger stored the user information in an LDAP database and organized the users into different groups. A VO administrator could create/add groups and populate users. This info could either be uploaded by user certificate or through the LDAP certificate server if the Department of Energy Science Grid (DOESG) certificates were used. The VO tool was installed on a server machine at Fermilab. Each site (including Fermilab) used the EDG mkgridmap script to generate the gridmap-file entries based on the information stored in the LDAP database at Fermilab. For the IGT Fall 2002 run, we used both Globus certificates and DOESG certificates.

### III.5 Monitoring

In the Fall 2002 IGT run, monitoring systems were not used for automatic controls such as scheduling, but rather for creating human readable displays only. Monitoring was divided logically into two different concerns. Local monitoring consists of gathering useful information from a single cluster or grid site, while Grid-wide monitoring consists of transporting and integrating information from local sites. IGT-wide monitoring was accomplished using the MonaLisa monitoring tool developed at Caltech bb:monalisa . Local monitoring was accomplished using modules written for MonaLisa at some sites and Ganglia at other sites.

## IV Running a Production Grid Continuously for the First Time

The IGT ran CMS production by integrating CMS applications with Grid tools. Several problems were identified and fixed during this run. These include integration issues arising from non-grid CMS tools integrated with Grid tools, bottlenecks arising from operating system limitations (eg- default limits not set high enough), and bugs in middleware and application software. Every component of the software contributed to overall problem count in some way. However, we found that with the current level of functionality, we were able to operate the IGT with 1.0 FTE effort during quiescent times over and above normal system administration and up to 2.5 FTE during crises. A sampler of problems appears below.

* (Pre-IGT) During Spring 2002, the Globus 2.0 GASS Cache was found to not support the required level of performance for CMS production. The software was re-engineered in consultation with Condor developers and Globus developers over the summer of 2002, and released in Globus 2.2.
* The Impala environment requires many helper files to run a production job. It was found that many simultaneous globus-url-copy operations originating from the MOP master site when submitting many jobs would cause some globus-url-copy operations to hang. Globus-url-copy operations were wrapped in Fault Tolerant Shell (FTSH) scripts. FTSH contains semantics to time-out and retry shell commands.
* Keeping in view the need of automatic resubmission of failed jobs; MOP used the auto restart features of Condor\_G. There have been several cases of this when actual jobs were still running on remote site. While the jobs eventually finish, there was wasted CPU time.
* There were some instances when jobs failed due to application code problems, or some real reason like full disk partitions. With auto-restart option, it becomes difficult to find out actual cause of trouble, while failed jobs are restarted and resubmitted in an infinite loop. The actual problem remain hidden unless jobs are watched closely.
* Condor\_G running on the MOP Master site uses ’gahp\_server’ to handle its communication with processes running under Globus on remote worker sites, one thread per tracked process. With over 400 CPUs available to IGT at later stages of production, running two assignments to produce 1.5 Million events, we had to divide production over two physically separate MOP master machines, to avoid scaling limit of the number of gahp\_server threads.

## V Analysis Environment on the IGT: CLARENS

A need was identified to provide a simple universal access method or ”portal” to data and CPU resources as part of the end-user physics analysis process. Clarens bb:clarens  is a flexible web services layer accessible through SOAP or XML-RPC which is designed for high security, high throughput request processing. Server functionality is easily extensible through administrator- or user-installed server-side modules written in C/C++ or the Python scripting language. Individual service requests are handled by a multi-process server, providing crash protection and support for non-blocking long-lived requests.

Security is based on X509 certificates for authentication, and optional transport encryption using SSL/TLS. A full Virtual Organization (VO) implementation coupled with fine-grained access control lists (ACL) provide powerful and easy to administer security for server methods and files. The traditional mode of certificate to user mapping for execution of server-side jobs is also supported.

Server methods available include file repository access, VO/ACL administration, job execution, and a remote interface to the Sql2ROOT/SOCATS (Stl based Object Caching And Transport System) RDBMS-based analysis framework. Client access to server methods is through ROOT (command line, interpreted and compiled C++), Python (command line and script-based), and Java.

## VI SC2002 Demo

A continuously updated display of simulated event production being run on the IGT as part of the Fall 2002 production was demonstrated at the Fermilab/SLAC booth at Supercomputing 2002 (SC2002) bb:sc2002 . The display was constructed using the ROOT analysis package using data obtained from Clarens servers involved in the IGT production at Caltech, Fermilab, UCSD and UFL. The data used were ROOT files produced in the last step in the production process at each site. The Clarens ROOT client continuously monitored the servers for newly produced data files, which were incorporated into the analysis on the show floor.

A second demo at the Caltech Center for Advanced Computing Research (CACR) booth showed an interactive analysis of event data stored at many remote RDBMS analyzed via the Clarens interface to SOCATS. The results of the analyses performed at Caltech and the Startlight point of presence (POP) in Chicago were made available as ROOT files that were interactively accessed through the ROOT Clarens client.

## VII Conclusions

The IGT was a success in that it produced all of the required events and provided many useful insights into operating a grid in production mode. Also, many problems were uncovered with the software at all levels. Figure 1 shows the progress of IGT full ntuple production during Fall of 2002.

Despite the problems, the production was remarkably smooth and sustained for over two months. The two notable flat spots occur during the SC2002 conference and during the winter holidays, which reflect loss of manpower during those periods.

In order to better quantify efficiency, the IGT run period was divided into 12 periods of about five days each. The average daily production rate in each interval was compared to the theoretical maximum daily rate of 45 K events per day IGT-wide666This number was estimated by scaling the daily throughput observed on one machine to all machines weighted by rated CPU speeds.. The average efficiency was just under 40%.

This performance is not much worse than the conventional CMS Spring 2002 Production. The Spring 2002 production was more complicated in that it involved more events with pileup and involved a lot more file transfers. Also, it is hard to calculate efficiency of the Spring 2002 production because it is hard to determine when a site was unavailable due to problems or just idle for lack of a request.

The EDG stresstest bb:edg-stresstest  ran during Fall of 2002 also. The EDG stresstest involved more functionality than the IGT in that it used a resource broker which relied on MDS to supply it with timely information. In short, they found that with the current state of Grid middleware, more functionality led to more problems, (and more FTE expended to track them down.) Viewed for what it was, the EDG stresstest was an immense success and found problems that were complementary to those found in the contemporaneous IGT run. Finally, our experience with the IGT has led to a plethora of documentation bb:igt1 ,bb:igt2  that can be a good start to providing production level support.

###### Acknowledgements.

We would like to acknowledge the CMS Core Computing and Software group and the USCMS Software and Computing projects for supporting this effort. We would especially like to thank Veronique Lefebure and Tony Wildish of the CMS Production Team for their support and helpful discussions. This effort was also supported by the USCMS Development Grid Testbed, including additional testbed hardware and support from the University of Wisconsin at Madison.