text:
# Vacuum polarization by a magnetic flux of special rectangular form

## 1  Introduction

Since the classical work  of H.B.G. Casimir, where the energy of vacuum polarized by two conducting planes was calculated, a number of similar problems was investigated for various external conditions (boundary conditions, background potentials etc.), see, for example, the recent review  or the books  . No general rule for the dependence of the vacuum energy on the background properties has been found so far. In particular, it is unknown how to forecast the sign of the energy.

An interesting problem is the calculation of the vacuum energy of a spinor field in the background of a magnetic field. This problem was first considered as early as in 1936 by W. Heisenberg and H. Euler  who were interested in the effective action in the background of a homogeneous magnetic field.

Presently, the interest is shifted to string like configurations. In  the contribution of the fermionic ground state energy to the stability of electro weak strings was addressed. In  the gluonic ground state energy in the background of a center-of-group vortex in QCD had been considered.

A very special inhomogeneous reflectionless magnetic background is for example the model for domain wall . In this case, the equation for the vacuum fluctuations was easy to solve analytically. However, the final expression for the effective action involves an additional integration over a momentum, thus bringing no real simplifications.

In QED, for the background of a flux tube with constant magnetic field inside the ground state energy was calculated in . Here it turned out to be negative remaining of course much smaller than the classical energy of the background field.

An interesting approach is that used in  where the issue of the sign and of bounds on fermionic determinants in a magnetic background had been considered.

The remarkably simple case of a magnetic background field concentrated on a cylindrical shell (i.e., a delta function shaped profile) was considered in . Here different signs of the vacuum energy turned out to be possible in dependence on the parameters (radius and flux).

In general, the issue of stability of the strings is of interest. In electro weak theory, and in QED in particular, the coupling is small. Hence the vacuum energy as being a one loop correction to the classical energy is suppressed by this coupling. While in QED a magnetic string is intrinsically unstable, in electro weak theory there are unstable and stable configuration known . Here quantum corrections may become important for the stability. A question of special interest is whether a strong or singular background may have a quantum vacuum energy comparable to the classical one.

For the calculation of ground state energy, it is necessary to subtract the ultraviolet divergences. This procedure is in general well known and we follow . Roughly speaking one has to subtract the contribution of the first few heat kernel coefficients ($`a_0`$ through $`a_2`$ in the given case of (3+1) dimensions). After that one can remove the intermediate regularization and one is left with finite expressions. However, in order to obtain these finite expressions in a form suitable for numerical calculations one has to go one more step. As described in , see also below in section 2, one has to add and to subtract of a certain part of the asymptotic expansion of the integrand.

In the present paper we generalize the analysis done in  to a rectangular shaped magnetic background field. So we consider the vacuum energy of a spinor field in QED for a background as given by Eqs. (2.1, 2.5). This problem is technically more involved and allows progress in two directions. First, it allows to refine the mathematical and numerical tools for such problems and, second, it allows to address the question how the vacuum energy behaves for an increasingly singular background (making the rectangle narrower). So this model interpolates to some extend between the flux tube with homogeneous field inside in  and the delta shaped one in .

The basic principles of this procedure are described briefly in Sec. 2. Namely it starts with the well known zeta functional regularization. The regularized ground state energy is represented as a zeta function of a hamiltonian spectrum and treated in termini of heat kernel expansion. The representation of regularized ground state energy as an integral of the logarithmic derivative of the Jost function for spinor wave scattering problem on the external magnetic background is obtained. The Jost function is obtained from the exact solutions of Dirac equation, derived in Sec. 3. The explicit form of the exact and asymptotic Jost function is considered in Sec. 4. In Sec.5 the representation useful for further numerical evaluations is derived. Sec. 6 is devoted to the calculation of the heat kernel coefficient $`a_{5/2}`$ and the Sec. 7 contains some numerical evaluations of the ground state energy. The divergent part is identified as that part of the corresponding heat kernel expansion which does not vanish for large $`m`$ (mass of the quantum spinor field). After the subtraction of this divergent part the remaining analytical expression must be transformed in order to lift the regularization. A part of the uniform asymptotic expansion of the Jost function is used for this procedure. Finally the analytical expression for the ground state energy has been evaluated numerically.

Throughout the paper $`\mathrm{}=c=1`$ is used.

## 2 The renormalized ground state energy

Consider a spinor field on a magnetic background $`\stackrel{}{B}`$ of the form

$$\stackrel{}{B}=\frac{\mathrm{\Phi }}{2\pi }h(r)\stackrel{}{e_z}$$
(2.1)
where $`\mathrm{\Phi }`$ is constant defining the magnetic flux, $`\stackrel{}{e_z}`$ is a unit vector in the cylindrical coordinate system ($`r,\varphi ,z`$). The Lagrangian is

$$=\frac{1}{4}F_{\mu \nu }F^{\mu \nu }+\overline{\psi }\left[i\gamma ^\mu (_\mu +eA_\mu )m\right]\psi .$$
(2.2)
The potential $`\stackrel{}{A}(r)`$

$$\stackrel{}{A}=\frac{\mathrm{\Phi }}{2\pi }\frac{a(r)}{r}\stackrel{}{e_\phi }$$
(2.3)
possesses cylindrical symmetry and the radial part $`a(r)`$ is taken to be

$$a(r)=\{\begin{array}{c}0,r<R_1,R_2\\ \frac{1}{\kappa }(\frac{r^2}{R_1^2}1)R_1rR_2\\ 1,r>R_1,R_2\end{array}$$
(2.4)
The profile function for the magnetic field is

$$h(r)=\frac{1}{r}\frac{}{r}a(r)=\{\begin{array}{c}0,r<R_1,R_2\\ \frac{2}{\kappa R_1^2}R_1rR_2\\ 0,r>R_1,R_2\end{array}$$
(2.5)
with $`\kappa =\frac{R_2^2R_1^2}{R_1^2}`$.The shape of the background can be interpreted geometrically as an infinitely long flux tube empty inside.

We use the zeta-function regularization for the vacuum energy $`E_0`$. Since the background is static, the following representation holds,

$$E_0(s)=\frac{1}{2}\mu ^{2s}\underset{(n)}{}\omega _{(n)}^{12s},$$
(2.6)
where

$$\omega _{(n)}=ϵ_{(n)}=\sqrt{\stackrel{}{p}^2+m^2}$$

are eigenfrequences resp. eigenvalues of energy for one particle states, (the spectrum of the corresponding hamiltonian, see the next section), $`\mu `$ is a parameter with the dimension of mass which is introduced to keep the correct dimension of energy in this expression.

For technical reasons we assume that the system is contained in a large but finite cylinder of radius $`R`$ in order to have discrete eigenvalues in the transversal directions. Because of the translational invariance along the z-axis we can separate the $`z`$-component of momentum,

$$E_0(s)=\frac{1}{2}\mu ^{2s}\underset{\mathrm{}}{\overset{\mathrm{}}{}}\frac{dp_z}{2\pi }\underset{(n)}{}(p_z^2+k_{(n)}^2+m^2)^{\frac{1}{2}s}.$$
(2.7)
Here is it necessary to remark that $`E_0`$ has the meaning of density per unit length. Integrating out $`p_z`$ we get

$$E_0(s)=\mu ^{2s}\frac{1}{4\sqrt{\pi }}4\frac{\mathrm{\Gamma }(s1)}{\mathrm{\Gamma }(s\frac{1}{2})}\underset{(n)}{}(k_(n)^2+m^2)^{1s}.$$
(2.8)
In this expression a factor 4 resulting from the summation over spin states and over the sign of the one particle energies appeared. The remaining sum over $`(n)`$ is over the eigenvalues $`k_(n)`$ of a two dimensional wave equation for one mode in the perpendicular plane. The further transformation of this expression is described more detailed in , , . The first step is based on the property of logarithmic derivative of eigenfunctions with respect to momentum k. The contour integral,

$$\underset{\gamma }{}\frac{dk}{2\pi i}(k^2+m^2)^{1s}_k\mathrm{ln}\varphi _l(kr),$$

is equal to the sum over the spectrum, $`\underset{(n)}{}(k_n^2+m^2)^{1s}`$, if the contour $`\gamma `$ encircles the real positive axis. Then (2.8) becomes

$$E_0(s)=\mu ^{2s}\frac{1}{\sqrt{\pi }}\frac{\mathrm{\Gamma }(s1)}{\mathrm{\Gamma }(s\frac{1}{2})}\underset{l=\mathrm{}}{\overset{\mathrm{}}{}}\underset{\gamma }{}\frac{dk}{2\pi i}(k^2+m^2)^{1s}_k\mathrm{ln}\varphi _l(kr).$$
(2.9)
Here $`\varphi _l(r)`$ are the eigenfunctions of hamiltonian $`\widehat{H}`$ labeled by the orbital momentum $`l`$.

In the present formalism instead of wave functions $`\varphi _l(r)`$ we use the corresponding Jost functions $`f_l(k)`$ which contain all necessary informations on the spectrum of scattering problem. For the space domain outside of cylindrical magnetic background, at $`r>R_2`$, the wave functions can be chosen as a linear composition of Hankel functions

$$\varphi _l(kr)=\frac{1}{2}[\overline{f_l}(k)H_l^{(1)}(kr)+f_l(k)H_l^{(2)}(kr)],$$
(2.10)
where $`\overline{f}_l(k)`$ and $`f_l(k)`$ are the corresponding Jost functions.

Now after the subtraction of contribution from empty Minkowski space (identified as the term divergent at $`R\mathrm{}`$, which is independent of the background and corresponds to the heat kernel coefficient $`a_0`$,(2.15) below), and rotation of the path towards the imaginary positive axis one obtains

$$E_0(s)=\frac{1}{\sqrt{\pi }}\frac{\mathrm{\Gamma }(s1)}{\mathrm{\Gamma }(s\frac{1}{2})}\frac{\mathrm{sin}\pi s}{\pi }\underset{l=\mathrm{}}{\overset{\mathrm{}}{}}\underset{m}{\overset{\mathrm{}}{}}𝑑k(k^2m^2)^{1s}_k\mathrm{ln}f_l(ik),$$
(2.11)
which is a very useful representation of the regularized ground state energy. A merit of this representation is the absence of oscillations of the integrand for large arguments111Note that bound states are also taken into account in $`E_0`$ (2.11). In the considered problem bound states appear if the flux is larger than one flux unit. Strictly speaking, these are zero modes located on the lower end of the continuous spectrum, i.e. at k=0 (this is known from )..

To discuss divergences in $`E_0`$ we need a more general setting. $`E_0(s)`$ can be expressed through the zeta-function of the associated differential operator $`\widehat{𝐏}`$

$$E_0(s)=\frac{\mu ^{2s}}{2}\zeta _{(\widehat{𝐏})}(s\frac{1}{2}).$$
(2.12)
$`\zeta _{(\widehat{𝐏})}(s)`$ admits an integral representation

$$\zeta _{\widehat{𝐏}}(s)=\underset{0}{\overset{\mathrm{}}{}}𝑑t\frac{t^{s1}}{\mathrm{\Gamma }(s)}K(t)$$
(2.13)
where the kernel of this integral representation $`K(t)=\underset{n}{}e^{t(k_{(n)}^2+m^2)}`$ is the so called heat kernel. The associated “heat kernel expansion” for the kernel $`K(t)`$ of this integral representation at $`t0`$ is

$$K(t)\frac{e^{tm^2}}{(4\pi t)^{3/2}}\underset{n0}{}a_nt^n.$$
(2.14)
In accordance to the heat kernel expansion, the coefficients $`a_n`$ occurred in (2.14) must be for our background :

$$a_0=V$$
(2.15)
is simple an (infinite) volume of configuration space; the part proportional to it has been even dropped before during the renormalization procedure as a contribution of the empty Minkowski space. All other coefficients up to $`a_2`$ for this background must be zero through dimensional reasons and requirement of gauge invariance, (for details see the paper ).

$$a_{1/2}=a_1=a_{3/2}=0.$$
(2.16)
The first nonzero coefficient is known to read

$$a_2=\frac{2}{3}F_{\mu \nu }^2$$
(2.17)
and only the term with $`n=2`$ contributes to the “divergent” part of energy , . It will be shown below that only the contribution proportional to $`a_2`$ in $`E_0(s)`$ contains a simple pole at $`s0`$ and in the case of pure magnetic background is proportional to the classical energy

$$E^{class}=\frac{1}{2}\stackrel{}{B}^2=\frac{8e^2}{3}a_2$$
(2.18)
In the paper  it has been found, that the regularized ground state energy of the magnetic flux tube of finite radius $`R`$ behaves as $`\frac{R^3}{m}`$ at $`R\mathrm{}`$ which corresponds to the contribution proportional to $`a_{5/2}`$. It holds in accordance to the observed fact that in the case of non-smooth backgrounds the heat kernel coefficients with half-integer number starting from some value corresponding to the smoothness-class of the background are different from zero. It occurs namely in our case, where the background potential $`\stackrel{}{A}(r)`$ is continous and the magnetic field $`h(r)`$ has a discontinuity. The coefficient $`a_{5/2}`$ is nonzero in our case, see below Sec.6. The next terms of this asymptotic may be $`\frac{R^4}{m^2}`$ corresponding to $`a_3`$, and $`a_4`$ delivering terms $`\frac{R^6}{m^4}`$.

We define the renormalized ground state energy as

$$E^{ren}=E_0E^{div},$$
(2.19)
where $`E^{div}`$ is obtained from the heat kernel expansion (2.14) and $`E^{ren}`$ fulfills the normalization condition

$$limE^{ren}=0atm\mathrm{}.$$
(2.20)
Inserting the heat kernel expansion (2.14) into the zeta function and using (2.12) we have

$$E^{div}=\frac{a_2}{32\pi ^2}(\frac{1}{s}2+ln\frac{4\mu ^2}{m^2}),$$
(2.21)
where $`a_2`$ is the only non-zero heat kernel coefficient $`a_n,0<n2`$, contributing to the divergent part of vacuum energy.

For the renormalized ground state energy the asymptotical dependence on powers of $`m`$ at $`m\mathrm{}`$ follows from Eq.(2.14) to have the form

$$E_0^{ren}\stackrel{}{_m\mathrm{}}\underset{n>2}{}\frac{e_n}{m^{2n}}$$
(2.22)
with some coefficients $`e_n`$.

In accordance with the interpretation of $`E^{ren}`$ it must vanish in the limit of $`m\mathrm{}`$ since it is the energy of vacuum fluctuations. Through the subtraction of terms containing all non-negative powers of $`m^2`$ (which are the terms of heat kernel expansion up to $`a_2`$) the condition (2.20) is satisfied automatically.

Some comments on the subtraction scheme are in order. This scheme is motivated by the physical assumption that the quantum fluctuations should vanish if the mass of the fluctuating field becomes large. The scheme had been used in a number of Casimir energy calculations. In  it had been shown to be equivalent to the so called “no tadpole” normalization condition which is common in field theory. It should be noticed that this scheme does not apply to massless fluctuating fields, for a discussion of this point see .

So $`E^{ren}`$ as given by Eq.(2.19)is now finite at $`s0,`$ but it is not possible to use this expression immediately since the integral in this limit does not exist; i.e. one cannot carry out the analytical continuation to $`s=0`$. In order to get a representation where this can be done one can make a trick , namely add and subtract a part $`\mathrm{ln}f_l^{asym}`$ of the uniform asymptotical expansion for $`\mathrm{ln}f_l`$ in (2.11) containing a minimal number of terms to provide the convergence of the remaining part after the subtraction. Thus we can separate the resulting expression into four terms

$$E^{ren}=(E_0E^{asym})+(E^{as}E^{div})$$
(2.23)
where

$$E^{asym}=2C_s\underset{l}{}\underset{m}{\overset{\mathrm{}}{}}𝑑k(k^2m^2)^{1s}\frac{}{k}\mathrm{ln}f_l^{asym}(ik)$$
(2.24)
and redefine two terms in brackets as

$`E^{ren}=E^f+E^{as}`$
$`E^f=E_0E^{asym}`$ (2.25)
$`E^{as}=E^{asym}E^{div},`$ (2.26)
thus $`E^f`$ is defined to be

$$E^f=\frac{1}{2\pi }\underset{l=\mathrm{}}{\overset{\mathrm{}}{}}\underset{m}{\overset{\mathrm{}}{}}𝑑k(k^2m^2)\frac{}{k}[lnf_l(ik)\mathrm{ln}f_l^{as}(ik)]$$
(2.27)
The “finite” part $`E^f`$ of the renormalized vacuum energy is not only well defined at $`s0`$, but also provides a representation well suited for numerical analysis. The analytical continuation for $`E^{as}`$ at this limit will be constructed below.

## 3  Solution of the Dirac equation

We consider a spinor quantum field $`\psi `$ in the background of the classical magnetic flux. We start with the Dirac equation for this field

$$\left\{i\gamma ^\mu \frac{}{x^\mu }me\gamma ^\mu A_\mu \right\}\psi =0$$
(3.1)
with the electromagnetic potential (2.3). The $`\psi `$ in (3.1) is a spinor quantum field of mass $`m`$ and charge $`e`$, interacting with $`\stackrel{}{A}`$ (2.3). The gamma matrices in our representation are chosen to be the same as in

$$\gamma ^0=\left\{\begin{array}{cc}\sigma _3& 0\\ 0& \sigma _3\end{array}\right\}\gamma ^1=\left\{\begin{array}{cc}i\sigma _2& 0\\ 0& i\sigma _2\end{array}\right\}\gamma ^2=\left\{\begin{array}{cc}i\sigma _1& 0\\ 0& i\sigma _1\end{array}\right\}\gamma ^3=\left\{\begin{array}{cc}0& 1\\ 1& 0\end{array}\right\}.$$
(3.2)
Now we follow the standard procedure and separate the variables. Using the ansatz

$$\varphi (r,\phi ,z)=e^{ip_0x^0}e^{ip_zz}\mathrm{\Psi }(r,\phi )$$
(3.3)
we obtain the equation for the 4- component spinor $`\mathrm{\Psi }=\left[\begin{array}{c}\varphi \\ \chi \end{array}\right]`$

$$\left\{\begin{array}{cc}p_0+\widehat{L}m\sigma _3& p_3\sigma _3\\ & \\ p_3\sigma _3& p_0+\widehat{L}+m\sigma _3\end{array}\right\}\mathrm{\Psi }=0,$$
(3.4)
where $`\widehat{L}=i\underset{i=1}{\overset{2}{}}\sigma _i(_i+ieA_i)`$, $`\varphi ,\chi `$ are the two-component spinors.

Respecting the translational invariance of the system it is sufficient to solve the equations only for $`p_3=0`$. Then one of the two decoupled equations $`(\varphi ,\chi )`$ reads

$$\left\{\begin{array}{cc}p_0m& \frac{}{r}\frac{l\delta a(r)}{r}\\ & \\ \frac{}{r}\frac{l+1\delta a(r)}{r}& p_0+m\end{array}\right\}\mathrm{\Phi }(r)=0$$
(3.5)
(here the standard ansatz

$$\mathrm{\Phi }=\left[\begin{array}{c}i\psi ^u(r)e^{i(l+1)\phi }\\ \\ \psi ^l(r)e^{il\phi }\end{array}\right]$$
(3.6)
has been used, $`l`$ is the orbital quantum number), the radial part denoted as

$$\mathrm{\Phi }(r)=\left\{\begin{array}{c}\psi ^u(r)\\ \\ \psi ^l(r)\end{array}\right\}$$
(3.7)
and

$$\delta =\frac{e\mathrm{\Phi }}{2\pi }.$$
(3.8)
In order to use later the symmetry properties of the Jost function we redefine the parameter $`l`$ in (3.6) (orbital number) as $`\nu `$ according to

$$\nu =\{\begin{array}{c}l+\frac{1}{2}forl=0,1,2,\mathrm{}\\ \\ l\frac{1}{2}forl=1,2,\mathrm{}\end{array}$$
(3.9)
For the given construction of potential (2.4) we get three equations and three types of solutions for 3 areas of space respectively
Domain I: $`r<R_1`$. The free wave equation (the magnetic flux is zero), the solutions are the Bessel functions $`J_{\nu \pm \frac{1}{2}}(kr)`$

$`\mathrm{\Phi }_I^{}(r)=\left[\begin{array}{c}\psi _I^u(r)\\ \\ \psi _I^l(r)\end{array}\right]={\displaystyle \frac{1}{\sqrt{p_0m}}}\left[\begin{array}{c}\sqrt{p_0+m}J_{\nu \frac{1}{2}}(kr)\\ \\ \sqrt{p_0m}J_{\nu +\frac{1}{2}}(kr)\end{array}\right],`$ (3.12)
$`\mathrm{\Phi }_I^+(r)=\left[\begin{array}{c}\psi _I^{u+}(r)\\ \\ \psi _I^{l+}(r)\end{array}\right]={\displaystyle \frac{1}{\sqrt{p_0m}}}\left[\begin{array}{c}\sqrt{p_0+m}J_{\nu +\frac{1}{2}}(kr)\\ \\ \sqrt{p_0m}J_{\nu \frac{1}{2}}(kr)\end{array}\right]`$ (3.15)
$`p_0=\sqrt{m^2k^2}`$
This solution of equation (3.5) in the domain I is chosen to be the so called regular solution which is defined as to coincide for $`r0`$ with the free solution.
Domain II: $`R_1<r<R_2`$ The equation with a homogeneous magnetic field
has the solution

$`\mathrm{\Phi }_{II}^{}(r)=C_r^{}\left[\begin{array}{c}\psi _{II.r}^u(r)\\ \\ \psi _{II.r}^l(r)\end{array}\right]+C_i^{}\left[\begin{array}{c}\psi _{II.i}^u(r)\\ \\ \psi _{II.i}^l(r)\end{array}\right]=`$ (3.18)
$`C_r^{}\left[\begin{array}{c}\frac{(p_0+m)}{2(\stackrel{~}{\alpha }+1)}r^{\stackrel{~}{\alpha }+1}e^{\frac{\beta r^2}{4}}{}_{1}{}^{}F_{1}^{}(1\frac{k^2}{2\beta },\stackrel{~}{\alpha }+2;\frac{\beta r^2}{2})\\ \\ r^{\stackrel{~}{\alpha }}e^{\frac{\beta r^2}{4}}{}_{1}{}^{}F_{1}^{}(\frac{k^2}{2\beta },\stackrel{~}{\alpha }+1;\frac{\beta r^2}{2})\end{array}\right]`$ (3.20)
(3.21)
$`+C_i^{}\left[\begin{array}{c}\frac{1}{p_0m}\frac{2\stackrel{~}{\alpha }}{r}(\frac{\beta r^2}{2})^{\stackrel{~}{\alpha }}e^{\frac{\beta r^2}{4}}{}_{1}{}^{}F_{1}^{}(\frac{k^2}{2\beta }\stackrel{~}{\alpha },\stackrel{~}{\alpha };\frac{\beta r^2}{2})\\ \\ (\frac{\beta r^2}{2})^{\stackrel{~}{\alpha }}e^{\frac{\beta r^2}{4}}{}_{1}{}^{}F_{1}^{}(\frac{k^2}{2\beta }\stackrel{~}{\alpha },1\stackrel{~}{\alpha };\frac{\beta r^2}{2})\end{array}\right].`$ (3.23)
$`\mathrm{\Phi }_{II}^+(r)`$ is the same, but $`\stackrel{~}{\alpha }`$ replaced by $`\alpha `$
here: $`\alpha =\nu \frac{1}{2}+\frac{\delta }{\kappa }`$ for $`l0`$
$`\stackrel{~}{\alpha }=\frac{\delta }{\kappa }\nu \frac{1}{2}`$ for $`l<0`$
$`\beta =\frac{2\delta }{\kappa R_1^2}`$
$`{}_{1}{}^{}F_{1}^{}`$ is a confluent hypergeometric function , .

The coefficients $`C_i,C_r`$ are some constants that will be irrelevant for expressing the Jost function. The indices u and l denote “upper” and “lower” components of spinor respectively. The lower index “i” or “r” corresponds to “regular” or “irregular” part of solutions dependent on the behaviour of the function if continued to $`r=0`$. If the external background vanishes ($`\delta 0`$), contributions of irregular parts to the solution disappear.
Domain III: $`r>R_2`$ The free wave equation outside of the magnetic flux has the solutions

$`\mathrm{\Phi }_{III}^{}(r)`$ $`=`$ $`{\displaystyle \frac{1}{2}}\overline{f}_\nu ^{}(k){\displaystyle \frac{1}{\sqrt{p_0m}}}\left[\begin{array}{c}\psi _{III.r}^u(r)\\ \\ \psi _{III.r}^l(r)\end{array}\right]+{\displaystyle \frac{1}{2}}f_\nu ^{}(k){\displaystyle \frac{1}{\sqrt{p_0m}}}\left[\begin{array}{c}\psi _{III.i}^u(r)\\ \\ \psi _{III.i}^l(r)\end{array}\right]`$ (3.26)
$`=`$ $`{\displaystyle \frac{1}{2}}\overline{f}_\nu ^{}(k)\left[\begin{array}{c}\hfill \eta H_{\nu \frac{1}{2}+\delta }^{(1)}(kr)\\ \\ \hfill H_{\nu +\frac{1}{2}+\delta }^{(1)}(kr)\end{array}\right]+{\displaystyle \frac{1}{2}}f_\nu ^{}(k)\left[\begin{array}{c}\hfill \eta H_{\nu \frac{1}{2}+\delta }^{(2)}(kr)\\ \\ \hfill H_{\nu +\frac{1}{2}+\delta }^{(2)}(kr)\end{array}\right],`$ (3.29)
$`\mathrm{\Phi }_{III}^+(r)`$ $`=`$ $`{\displaystyle \frac{1}{2}}\overline{f}_\nu ^+(k){\displaystyle \frac{1}{\sqrt{p_0m}}}\left[\begin{array}{c}\psi _{III.r}^{u+}(r)\\ \\ \psi _{III.r}^{l+}(r)\end{array}\right]+{\displaystyle \frac{1}{2}}f_\nu ^+(k){\displaystyle \frac{1}{\sqrt{p_0m}}}\left[\begin{array}{c}\psi _{III.i}^{u+}(r)\\ \\ \psi _{III.i}^{l+}(r)\end{array}\right]`$ (3.32)
$`=`$ $`{\displaystyle \frac{1}{2}}\overline{f}_\nu ^+(k)\left[\begin{array}{c}\hfill \eta H_{\nu +\frac{1}{2}\delta }^{(1)}(kr)\\ \\ \hfill H_{\nu \frac{1}{2}\delta }^{(1)}(kr)\end{array}\right]+{\displaystyle \frac{1}{2}}f_\nu ^+(k)\left[\begin{array}{c}\hfill \eta H_{\nu +\frac{1}{2}\delta }^{(2)}(kr)\\ \\ \hfill H_{\nu \frac{1}{2}\delta }^{(2)}(kr)\end{array}\right],`$ (3.35)
$`\eta `$ $`=`$ $`{\displaystyle \frac{\sqrt{p_0+m}}{\sqrt{p_0m}}}`$
and $`f_\nu (k)`$ is the desired Jost function in accordance to the definition ( 2.10). The functions $`f_l(k)`$ and $`\overline{f}_l(k)`$ are conjugated to each other because of the choice of the solution in the domain I to be regular (3).

## 4 The exact and asymptotic Jost function for the background

To obtain a Jost function we need to impose certain matching conditions for the spinor wave at the boundaries between I-II and II-III domains. It is known that for a continuous potential $`a(r)`$ as we consider here, the spinor wave function must be continuous, hence for its components it holds

$`\mathrm{\Phi }_{II}^{}(r)=\mathrm{\Phi }_I^{}(r)|_{r=R_1},`$
$`\mathrm{\Phi }_{III}^{}(r)=\mathrm{\Phi }_{II}^{}(r)|_{r=R_2}`$ (4.1)
Resolving these equations for $`f_\nu ^\pm `$ we obtain, that

$`f_\nu ^\pm (k)`$ $`=`$ $`2{\displaystyle \frac{A+B}{W_{III}W_{II}}}`$
$`A`$ $`=`$ $`(\psi _I^{u\pm }\psi _{II.i}^{l\pm }\psi _I^{l\pm }\psi _{II.i}^{u\pm })|_{R_1}(\psi _{II.r}^{u\pm }\psi _{III.r}^{l\pm }\psi _{II.r}^{l\pm }\psi _{III.r}^{u\pm })|_{R_2},`$
$`B`$ $`=`$ $`(\psi _I^{l\pm }\psi _{II.r}^{u\pm }\psi _I^{u\pm }\psi _{II.r}^{u\pm })|_{R_1}(\psi _{II.i}^{u\pm }\psi _{III.r}^{l\pm }\psi _{II.i}^{l\pm }\psi _{III.r}^{u\pm })|_{R_2}`$ (4.2)
The denominator of this expression can be written using the Wronskians of hypergeometric and Bessel functions as follows

$`W_{III}=(\psi _{III.i}^u\psi _{III.r}^l\psi _{III.i}^l\psi _{III.r}^u)|_{R_2}=({\displaystyle \frac{p_0+m}{p_0m}})^{\frac{1}{2}}{\displaystyle \frac{4i}{\pi kR_2}}`$
(4.3)
$`W_{II}=(\psi _{II.r}^u\psi _{II.i}^l\psi _{II.r}^l\psi _{II.i}^u)|_{R_1}={\displaystyle \frac{\alpha \beta }{p_0m}}({\displaystyle \frac{\beta R_1}{2}})^{\alpha 1}`$ (4.4)
To calculate $`E_0`$ we need this function with imaginary argument (2.11) so we need to replace $`k`$ by $`ik`$. As a result we obtain the new expression for $`f(ik)`$ that contains now modified Bessel functions $`I_{\nu \pm \frac{1}{2}}`$ instead of $`J_{\nu \pm \frac{1}{2}}`$ and modified Bessel functions $`K_{\nu \pm \frac{1}{2}\delta }`$ instead of $`H_{\nu \pm \frac{1}{2}\delta }^{(1,2)}`$.

$`f_\nu ^+(ik)`$ $`=`$ $`{\displaystyle \frac{\pi kR_1R_2e^{\delta /\kappa +\frac{i\pi \delta }{2}}}{\frac{1}{2}\stackrel{~}{\nu }}}\times `$
$`\{[{\displaystyle \frac{kR_1}{2\stackrel{~}{\nu }+1}}{}_{1}{}^{}F_{1}^{}(1+{\displaystyle \frac{k^2}{2\beta }},\stackrel{~}{\nu }+{\displaystyle \frac{3}{2}};{\displaystyle \frac{\beta R_1^2}{2}})I_{\nu \frac{1}{2}}(kR_1)`$
$`{}_{1}{}^{}F_{1}^{}({\displaystyle \frac{k^2}{2\beta }},\stackrel{~}{\nu }+{\displaystyle \frac{1}{2}};{\displaystyle \frac{\beta R_1^2}{2}})I_{\nu +\frac{1}{2}}(kR_1)]\times `$
$`[{\displaystyle \frac{2\stackrel{~}{\nu }+1}{R_2}}{}_{1}{}^{}F_{1}^{}({\displaystyle \frac{1}{2}}\stackrel{~}{\nu }+{\displaystyle \frac{k^2}{2\beta }},{\displaystyle \frac{1}{2}}\stackrel{~}{\nu };{\displaystyle \frac{\beta R_2^2}{2}})K_{\nu \frac{1}{2}\delta }(kR_2)`$
$`k{}_{1}{}^{}F_{1}^{}({\displaystyle \frac{1}{2}}\stackrel{~}{\nu }+{\displaystyle \frac{k^2}{2\beta }},{\displaystyle \frac{3}{2}}\stackrel{~}{\nu };{\displaystyle \frac{\beta R_2^2}{2}})K_{\nu +\frac{1}{2}\delta }(kR_2)]`$
$`[k{}_{1}{}^{}F_{1}^{}({\displaystyle \frac{1}{2}}\stackrel{~}{\nu }+{\displaystyle \frac{k^2}{2\beta }},{\displaystyle \frac{3}{2}}\stackrel{~}{\nu };{\displaystyle \frac{\beta R_1^2}{2}})I_{\nu +\frac{1}{2}}(kR_1)+`$
$`{\displaystyle \frac{2\stackrel{~}{\nu }1}{R_1}}{}_{1}{}^{}F_{1}^{}({\displaystyle \frac{1}{2}}\stackrel{~}{\nu }+{\displaystyle \frac{k^2}{2\beta }},{\displaystyle \frac{1}{2}}\stackrel{~}{\nu };{\displaystyle \frac{\beta R_1^2}{2}})I_{\nu \frac{1}{2}}(kR_1)]\times `$
$`[{\displaystyle \frac{kR_2}{2\stackrel{~}{\nu }+1}}{}_{1}{}^{}F_{1}^{}(1+{\displaystyle \frac{k^2}{2\beta }},\stackrel{~}{\nu }+{\displaystyle \frac{3}{2}};{\displaystyle \frac{\beta R_2^2}{2}})K_{\nu \frac{1}{2}\delta }(kR_2)+`$
$`{}_{1}{}^{}F_{1}^{}({\displaystyle \frac{k^2}{2\beta }},\stackrel{~}{\nu }+{\displaystyle \frac{1}{2}};{\displaystyle \frac{\beta R_2^2}{2}})K_{\nu +\frac{1}{2}\delta }(kR_2)]\},`$
where $`\stackrel{~}{\nu }=\nu +\frac{\delta R_1^2}{R_1^2R_2^2}`$, and

$`f_\nu ^{}(ik)`$ $`=`$ $`{\displaystyle \frac{\pi kR_1R_2e^{\delta /\kappa \frac{i\pi \delta }{2}}}{\frac{1}{2}+\overline{\nu }}}\times `$
$`\{[{\displaystyle \frac{kR_1}{2\overline{\nu }1}}{}_{1}{}^{}F_{1}^{}(\overline{\nu }+{\displaystyle \frac{1}{2}}+{\displaystyle \frac{k^2}{2\beta }},\overline{\nu }+{\displaystyle \frac{1}{2}};{\displaystyle \frac{\beta R_1^2}{2}})I_{\nu +\frac{1}{2}}(rR_1)`$
$`{}_{1}{}^{}F_{1}^{}(\overline{\nu }+{\displaystyle \frac{1}{2}}+{\displaystyle \frac{k^2}{2\beta }},\overline{\nu }+{\displaystyle \frac{3}{2}};{\displaystyle \frac{\beta R_1^2}{2}})I_{\nu \frac{1}{2}}(kR_1)]\times `$
$`[{\displaystyle \frac{2\overline{\nu }1}{R_2}}{}_{1}{}^{}F_{1}^{}(1+{\displaystyle \frac{k^2}{2\beta }},{\displaystyle \frac{3}{2}}\overline{\nu };{\displaystyle \frac{\beta R_2^2}{2}})K_{\nu +\frac{1}{2}+\delta }(kR_2)`$
$`k{}_{1}{}^{}F_{1}^{}({\displaystyle \frac{k^2}{2\beta }},{\displaystyle \frac{1}{2}}\overline{\nu };{\displaystyle \frac{\beta R_2^2}{2}})K_{\nu \frac{1}{2}+\delta }(kR_2)]`$
$`[k{}_{1}{}^{}F_{1}^{}({\displaystyle \frac{k^2}{2\beta }},{\displaystyle \frac{1}{2}}\overline{\nu };{\displaystyle \frac{\beta R_1^2}{2}})I_{\nu \frac{1}{2}}(kR_1)+`$
$`{\displaystyle \frac{2\overline{\nu }+1}{R_1}}{}_{1}{}^{}F_{1}^{}(1+{\displaystyle \frac{k^2}{2\beta }},{\displaystyle \frac{3}{2}}\overline{\nu };{\displaystyle \frac{\beta R_1^2}{2}})I_{\nu +\frac{1}{2}}(kR_1)]\times `$
$`[{\displaystyle \frac{kR_2}{2\overline{\nu }1}}{}_{1}{}^{}F_{1}^{}(\overline{\nu }+{\displaystyle \frac{1}{2}}+{\displaystyle \frac{k^2}{2\beta }},\overline{\nu }+{\displaystyle \frac{1}{2}};{\displaystyle \frac{\beta R_2^2}{2}})K_{\nu +\frac{1}{2}+\delta }(kR_2)+`$
$`{}_{1}{}^{}F_{1}^{}(\overline{\nu }+{\displaystyle \frac{1}{2}}+{\displaystyle \frac{k^2}{2\beta }},\overline{\nu }+{\displaystyle \frac{3}{2}};{\displaystyle \frac{\beta R_2^2}{2}})K_{\nu \frac{1}{2}+\delta }(kR_2)]\},`$
where $`\overline{\nu }=\nu \frac{\delta R_1^2}{R_1^2R_2^2}`$, and

Now we need to obtain the asymptotic Jost function and asymptotic part of energy. We use the representation of the asymptotic Jost function $`\mathrm{ln}f_{as}(ik)`$ in the following form (see App.)

$$lnf_{as}(ik)=\underset{n=1}{\overset{3}{}}\underset{j=1}{\overset{9}{}}\underset{0}{\overset{\mathrm{}}{}}\frac{dr}{r}X_{nj}\frac{t^j}{\nu ^n}$$
(4.6)
where $`t=(1+[\frac{kr}{v}]^2)^{\frac{1}{2}}`$ and $`X_{nj}`$ (given explicitly in App.A) are represented in terms of $`r,\delta ,a(r)`$ and their derivatives.

This expansion had been obtained in  by iterations of Lippmann-Schwinger equation up to the order $`\nu ^3`$. In general it is possible to obtain higher orders using this formalism. However as the calculations  showed,the complication of the involved expressions increases very fast. It is remarkable that this expression does not contain a term with power $`\nu ^2`$. In the finite part of the energy the corresponding term is canceled in the sum of terms corresponding to positive and negative orbital momenta $`\nu `$ as well. The absence of the power $`\nu ^2`$ is a succession of the zero heat kernel coefficient $`a_{3/2}`$. Also it is a non-trivial fact that both, the fourth power of the magnetic flux $`\delta `$ and the second one is present.

It can be checked numerically that $`lnf^{as}`$ (4.6) is indeed the uniform asymptotic expansion of logarithm of (4) for $`k,\nu `$\- large, $`k/\nu =z`$-fixed.

## 5 The finite and asymptotic parts of vacuum energy

The expression (4.6) for the uniform asymptotic of Jost function, substituted into the expression for $`E^{as}`$ (2.26) yields

$$E^{as}=2C_s\underset{v=1/2,3/2,\mathrm{}}{}\underset{m}{\overset{\mathrm{}}{}}𝑑k(k^2m^2)^{1s}\frac{}{k}\underset{0}{\overset{\mathrm{}}{}}\frac{dr}{r}\underset{n=1}{\overset{3}{}}\underset{j=1}{\overset{9}{}}\underset{0}{\overset{\mathrm{}}{}}X_{nj}\frac{t^j}{\nu ^n}E^{div}$$
(5.1)
with the constant $`C_s=\frac{1}{2\pi }[1+s(1+2ln2\mu )]`$. $`E^{div}`$ as defined in (2.21) is

$$E^{div}=\frac{\delta ^2}{6\pi (R_2^2R_1^2)}(\frac{1}{s}2+\mathrm{ln}\frac{4\mu ^2}{m^2})$$
(5.2)
with the coefficient $`a_2`$ according to (2.17) has been used in the form

$$a_2=\frac{2}{3}F_{\mu \nu }^2=\frac{8\pi }{3}\delta ^2\underset{0}{\overset{\mathrm{}}{}}𝑑rrh(r)^2=\frac{16\pi \delta ^2}{3}\frac{1}{R_2^2R_1^2}$$
(5.3)
Then the sum over $`\nu `$ can be now transformed into two integrals by means of Abel-Plana formula (9.1). The first one cancels the $`E^{div}`$ exactly. Then we integrate the second one over $`k`$ using the identities (9.2, 9.3). It gives the form

$$E^{as}=\frac{1}{\pi }m^2\underset{n=1}{\overset{3}{}}\underset{j=n}{\overset{3n}{}}\underset{0}{\overset{\mathrm{}}{}}\frac{dr}{r}X_{n,j}\mathrm{\Sigma }_{nj}(rm)$$
(5.4)
with

$`\mathrm{\Sigma }_{nj}(x)={\displaystyle \frac{\mathrm{\Gamma }(s+j/21)}{\mathrm{\Gamma }(j/2)}}{\displaystyle \frac{i}{x^j}}`$
$`\times `$ $`{\displaystyle \underset{0}{\overset{\mathrm{}}{}}}{\displaystyle \frac{d\nu }{1+exp(2\pi \nu )}}\left({\displaystyle \frac{(i\nu )^{jn}}{[1+(i\nu /x)^2]^{s+j/21}}}{\displaystyle \frac{(i\nu )^{jn}}{[1+(i\nu /x)^2]^{s+j/21}}}\right)`$
In order to obtain an analytical continuation in $`s=0`$ of each term of the sum we integrate it over $`\nu `$ by parts several times till the divergency at $`\nu =0`$ through the power $`s+j/21`$ abrogates, and after that we can perform the integration over $`r`$. Further we integrate by parts resulting in the relations $`raa^{}\frac{1}{2}a^2r_r`$ and $`r^2aa^{\prime \prime }a^2+\frac{1}{2}a^2r_r^2r`$ (that hold because of the continuity of the potential $`a(r)`$), and we obtain finally the form

$$E^{as}=\frac{4}{\pi }\underset{0}{\overset{\mathrm{}}{}}\frac{dr}{r^3}[\delta ^2a(r)^2g_1(rm)\delta ^2r^2a^{}(r)^2g_2(rm)+\delta ^4a(r)^4g_3(rm)]$$
(5.6)
with $`g_i(x)=\underset{x}{\overset{\mathrm{}}{}}𝑑\nu \sqrt{\nu ^2x^2}f_i(\nu )`$ that will be calculated numerically (see the next section), the functions $`f_i(\nu )`$ are shown explicitly in (9.7).

The finite part of the ground state energy $`E^f`$ can be finally represented as follows

$$E^f=\frac{1}{2\pi }\underset{\nu =1/2,3/2,\mathrm{}}{}\underset{m}{\overset{\mathrm{}}{}}𝑑k(k^2m^2)\frac{}{k}[lnf_\nu ^+(ik)+lnf_\nu ^{}(ik)2\mathrm{ln}f^{as}(ik)]$$
(5.7)
In order to use this form for numerical evaluations we integrate by parts and get

$$E^f=\frac{1}{\pi }\underset{\nu =1/2,3/2,\mathrm{}}{}\underset{m}{\overset{\mathrm{}}{}}𝑑kklnf_\nu ^{sub}(ik)$$
(5.8)
(the logarithmic expression in square brackets in (5.7) denoted as $`lnf^{sub}`$).

## 6 Higher orders of the uniform asymptotic expansion of the Jost function and the heat kernel coefficient $`a_{5/2}`$

The background considered in this paper has singular surfaces where the magnetic field jumps. The heat kernel expansion for the case of singularities concentrated at surfaces has been considered in \[, , , \]. Although the general analysis of  is valid for our background, an explicit expression for $`a_{5/2}`$ has not been calculated yet.

We can use our obtained Jost function (4) to calculate the coefficient $`a_{5/2}`$ in the heat kernel expansion (2.14).

Suppose we have obtained the value of $`E_0(s)`$ (2.11)in the point $`s=1/2`$. It follows from (2.12 \- 2.14) that

$$E_0(s)\frac{\mu ^{2s}}{2(4\pi )^{3/2}}\frac{m^{42s}}{\mathrm{\Gamma }(s\frac{1}{2})}\underset{n=0,\frac{1}{2},1,\mathrm{}}{}\frac{a_n}{m^{2n}}\mathrm{\Gamma }(s2+n)$$
(6.1)
and at the limit of $`s\frac{1}{2}`$ only the term containing $`a_{5/2}`$ remains to be nonzero in the sum.

$$E_0(\frac{1}{2})=\frac{\mu ^{2s}}{2(4\pi )^{3/2}}a_{5/2}$$
(6.2)
From the other hand we have (2.11, 2.24)

$$E_0(s)=C_s\underset{\nu =1/2,3/2,\mathrm{}}{}\underset{m}{\overset{\mathrm{}}{}}𝑑k(k^2m^2)^{1s}\frac{}{k}f_\nu (ik)=C_sh(s).$$
(6.3)
Here we substitute the exact Jost function by its uniform asymptotic represented in the form

$$f_\nu ^{ua}(ik)=\underset{n=1,3,4,5\mathrm{}}{}\frac{h_n(t_1,t_2)}{\nu ^n},$$
(6.4)
where the coefficients $`h_n`$ are functions of $`t_1=(1+(\frac{kR_1}{\nu })^2)^{1/2}`$, $`t_2=(1+(\frac{kR_2}{\nu })^2)^{1/2}`$, the power $`\nu ^2`$ is absent, as noticed above (4.6),

$$C_s=\frac{\mu ^{2s}\mathrm{\Gamma }(s1)}{4\sqrt{\pi }\mathrm{\Gamma }(s1/2)}\frac{4\mathrm{sin}(\pi s)}{\pi }$$
(6.5)
At the limit $`s\frac{1}{2}`$ it yields

$$E_0(\frac{1}{2})=\frac{4}{3\pi }Res_{s\frac{1}{2}}h(s)$$
(6.6)
and therefore we obtain for $`a_{5/2}`$

$$a_{5/2}=\frac{64\sqrt{\pi }}{3}Res_{s\frac{1}{2}}h(s).$$
(6.7)
To obtain the explicit form of $`Res_{s\frac{1}{2}}h(s)`$ we use the uniform asymptotic expansion of the Jost function (9.8). The terms $`h_n(t_1,t_2)`$ in 6.4 can be obtained either by iterations of Lippmann-Schwinger equation (see , ) or by using the explicit form of the Jost function as well. All the further terms up from $`h_4`$ are produced from the explicit form of the Jost function (4) because of complication of the first way for higher orders $`n`$ (see the remark to 4.6 in Sec.6). Namely, we obtain several higher orders $`1/\nu `$ of uniform asymptotic expansion for special functions $`I_\nu ,K_\nu ,_1F_1`$ which the exact Jost function (4) consists of (it can be done starting with the explicit form for two first orders and executing the recursive algorithm several times), then after substitution of each of functions $`I,K,_1F_1`$ by its corresponding uniform asymptotic expansion and separation of powers of $`\nu `$ we arrive at the form (9.8). The coefficients $`h_n(t_1,t_2),n=1,\mathrm{},4`$ are given in the Appendix (9.8).

If the function $`h_n(t_1,t_2)`$ is a polynomial over $`t_1,t_2`$, so we can consider some term $`t^j,(t=t_i,i=1,2)`$ of it. Notice, that for $`h_1(t_1,t_2),h_3(t_1,t_2)`$ it is not the case, but we can treat the terms of kind $`\frac{1}{1+t}`$ and $`\frac{1}{(1+t)^k}`$, $`k`$-integer, as an infinite sum of powers $`t`$.

Apart from the construction of $`t=t_1,t_2`$ (6.4) they are strictly positive and less than 1, therefore the series $`\frac{1}{1+t}=\underset{i=0}{\overset{\mathrm{}}{}}t^i`$ converges regular and uniform and respecting that we have only finite integer powers $`k`$ so that $`\left[\frac{1}{1+t}\right]^k=\left(\underset{i=0}{\overset{\mathrm{}}{}}t^i\right)^k`$ converges as well, the sum over $`i`$ can be interchanged with the one over $`\nu `$ (Princeheim’s Theorem) thus the following procedure is valid. Performing the sum

$$h(s)=\underset{\nu =1/2,3/2,\mathrm{}}{}\underset{m}{\overset{\mathrm{}}{}}𝑑k(k^2m^2)^{1s}\frac{}{k}\underset{n=1,3,4,5}{}h_{n,j}(R_1,R_2)\frac{t^j}{\nu ^n}$$
(6.8)
by meaning of (9.1) we obtain the sum of two parts

$`h(s)`$ $`=`$ $`{\displaystyle \underset{0}{\overset{\mathrm{}}{}}}𝑑\nu {\displaystyle \underset{m}{\overset{\mathrm{}}{}}}𝑑k(k^2m^2)^{1s}{\displaystyle \frac{}{k}}{\displaystyle \underset{n}{}}h_{n,j}(R_1,R_2){\displaystyle \frac{t^j}{\nu ^n}}+`$
$`{\displaystyle \underset{0}{\overset{\mathrm{}}{}}}{\displaystyle \frac{d\nu }{e^{2\pi \nu }+1}}{\displaystyle \frac{1}{i}}\left[{\displaystyle \underset{m}{\overset{\mathrm{}}{}}}𝑑k(k^2m^2)^{1s}{\displaystyle \frac{}{k}}{\displaystyle \underset{n}{}}h_{n,j}(R_1,R_2){\displaystyle \frac{t^j}{\nu ^n}}\right]_{\nu =i\nu }^{\nu =i\nu },`$
where $`t=t_{1,2}`$. The first summand of (6) gives for each power $`j`$ of $`t=\{t_1,t_2\}`$ (using the (9.2))

$$\frac{m^3}{2}\frac{\mathrm{\Gamma }(2s)\mathrm{\Gamma }(\frac{1+jn}{2})\mathrm{\Gamma }(s+\frac{n3}{2})}{(Rm)^{n1}\mathrm{\Gamma }(j/2)}$$
(6.10)
(where R denotes $`R_1`$ and $`R_2`$ respectively) at the limit $`s\frac{1}{2}`$ only the terms corresponding to $`n=2`$ and $`n=4`$ (which will be calculated below) have a pole.

For the second part of the Abel-Plana formula (6) we have

$`{\displaystyle \frac{1}{i}}{\displaystyle \underset{0}{\overset{\mathrm{}}{}}}{\displaystyle \frac{d\nu }{1+e^{2\pi \nu }}}{\displaystyle \frac{\mathrm{\Gamma }(2s)\mathrm{\Gamma }(s+j/21)}{jR^{42s}\mathrm{\Gamma }(j/2)}}`$
$`\times `$ $`\left[(i\nu )^{2+jn}[(mR)^2+(i\nu )^2]^{1j/2s}(i\nu )^{2+jn}[(mR)^2+(i\nu )^2]^{1j/2s}\right]`$
It can be seen using the Taylor expansions at $`s\frac{1}{2}`$ of Gamma functions entering the (6), that only the terms containing -1, 1, and 3 powers of $`t_i`$ can contribute to the residuum at $`s\frac{1}{2}`$. But the expression in square brackets in (6) can be performed as

$$\nu ^2(\nu ^2(mr)^2)[(i\nu )^{jn}e^{i\pi (1j/2s)}(i\nu )^{jn}e^{i\pi (1j/2s)}]$$
(6.12)
thus for $`jn`$ even the expression in square brackets produces (dropping the non sufficient coefficient $`\pm 2`$ or $`\pm 2i`$):
for $`j=1,1:\mathrm{sin}\pi (1j/2s)=\mathrm{cos}\pi s`$
for $`j=3:\mathrm{sin}\pi (1j/2s)=\mathrm{cos}\pi s`$
and for $`jn`$ odd respectively
for $`j=\pm 1:\mathrm{cos}\pi (1j/2s)=\pm \mathrm{sin}\pi s`$
for $`j=3:\mathrm{cos}\pi (1j/2s)=\mathrm{sin}\pi s`$
But for $`s\frac{1}{2}`$ these functions behave as
$`\mathrm{sin}\pi s1+\frac{\pi ^2}{2}(\frac{1}{2}+s)^2,\mathrm{cos}\pi s\pi (\frac{1}{2}+s)`$ and it means that only the contribution of terms with odd powers $`jn`$ of $`i\nu `$ could survive and these corresponds for the possible values of $`j`$ to $`\nu ^2,\nu ^4,\nu ^6,\mathrm{}`$. In fact the coefficient $`h_2(s)`$ at $`\nu ^2`$ is zero, and all other possible terms up from $`\nu ^4`$ does not contain any powers of $`t_i`$ lower than 4; one can see it for example in the explicit form of the uniform asymptotical expansion of special functions entering the $`\mathrm{ln}f_v(ik)`$ (4). Therefore the second summand of (6) does not produce any contribution to the $`Res_{s\frac{1}{2}}h(s)`$ Thus we have that only the contribution from the first summand of (6) remains, and since the term of $`n=2`$ is zero the searched residuum resulting from the term of $`n=4`$ is:

$$Res_{s\frac{1}{2}}h(s)=\frac{1}{2}.\frac{3\sqrt{\pi }}{4}\frac{\delta ^2}{4(R_1^2R_2^2)^2}(R_1+R_2)\underset{j=4,6}{}\frac{\mathrm{\Gamma }(\frac{1+jn}{2})}{j/2},$$
(6.13)
where the (6.10) and the explicit form of $`h_4(t_1,t_2)=\frac{\delta ^2}{4(R_1^2R_2^2)^2}[R_1^4(t_1^4t_1^6)+R_2^4(t_2^4t_2^6)]`$ have been used. Finally we have

$$Res_{s\frac{1}{2}}h(s)=\frac{15\pi }{128}\frac{\delta ^2(R_1+R_2)}{(R_1^2R_2^2)^2}$$
(6.14)
and therefore

$$a_{5/2}=\frac{5\pi ^{3/2}}{2}\frac{\delta ^2(R_1+R_2)}{(R_1^2R_2^2)^2}$$
(6.15)
This is the heat kernel coefficient $`a_{5/2}`$ for the configuration of the magnetic background field as given by Eqs.(2.3-2.5).

We can calculate the heat kernel coefficient $`a_{5/2}`$ in a more general situation when the magnetic field jumps on an arbitrary surface $`\mathrm{\Sigma }`$. The coeffcients $`a_n`$ for $`n=1/2,\mathrm{},2`$ can be read off rather general expressions of the paper . Let $`B^\pm `$ be values of the magnetic field on two sides of $`\mathrm{\Sigma }`$. According the analysis of , the coeffcient $`a_{5/2}`$ must be an integral over $`\mathrm{\Sigma }`$ of a local invarinat of canonical mass dimension 4, which is symmetric under the exchange of $`B^+`$ and $`B^{}`$ and which vanishes if $`B^+=B^{}`$ (i.e. when the singularity disappears). There is only one such invariant which gives rise to the following expresssion:

$$a_{5/2}=\xi \underset{\mathrm{\Sigma }}{}(\stackrel{}{B}^+\stackrel{}{B}^{})^2𝑑\mu (\mathrm{\Sigma }),$$
(6.16)
where the integration goes over the surface $`\mathrm{\Sigma }`$ and $`\stackrel{}{B}^\pm `$ are the values of the magnetic field on both sides of $`\mathrm{\Sigma }`$ in the given point. The yet undefined constant $`\xi `$ can be found using Eq.(6.15) which constitutes a special case of (6.16). Here the surface $`\mathrm{\Sigma }`$ consists of two circles in the $`(\stackrel{}{X},\stackrel{}{Y})`$-plane so that

$$\underset{\mathrm{\Sigma }}{}(\stackrel{}{B}^+\stackrel{}{B}^{})^2𝑑\mu (\mathrm{\Sigma })=2\pi (R_1+R_2)\stackrel{}{B}^2,$$
(6.17)
where the jump $`(\stackrel{}{B}^+\stackrel{}{B}^{})`$ is just the value of $`\stackrel{}{B}`$ at $`r[R_1,R_2]`$. With this, Eq.(6.16) takes the form

$$a_{5/2}=2\pi e^2\xi \stackrel{}{B}^2(R_1+R_2).$$
(6.18)
On the other side, from Eq.(2.1) we have

$$\mathrm{\Phi }=B\underset{r[R_1,R_2]}{}Bd^2x=\pi B(R_2^2R_1^2)$$
(6.19)
and with Eq.(3.8) from Eq.(6.15) it follows that

$$a_{5/2}=\frac{5}{8\sqrt{\pi }}e^2\mathrm{\Phi }^2\frac{R_1+R_2}{(R_2^2R_1^2)^2}.$$
(6.20)
comparing (6.18) with (6.20) we get

$$\xi =\frac{5}{16\pi ^{3/2}}.$$
(6.21)
## 7 Graphics and numerical results

The asymptotic part of the ground state energy as given by Eqs. (5.4) and (5.6) can be represented as sum of a part proportional to the second power of the flux and one proportional to the fourth power,

$$E^{as}=\frac{4}{\pi }[e_1(R_1,R_2)\delta ^2+e_2(R_1,R_2)\delta ^4].$$
(7.1)
The corresponding coefficients $`e_1(R_1,R_2)`$ and $`e_2(R_1,R_2)`$ can be calculated numerically without problems. They are shown as functions of $`R_2`$ for fixed $`R_1`$ in Fig.1 (multiplied by $`R_2^2`$).

The finite part of the ground state energy, $`E_f`$, is used in the form as given by Eq. (5.8). Let $`E_f(\nu ,k)`$ denote the function to be integrated and summed over in that expression. In Fig.2 it is shown as function of $`k`$ for several values of the orbital momentum $`\nu `$. In order the make the behavior better visible it is multiplied there by $`\nu ^4k^3`$. These functions are smooth for all values of $`k`$, starting from some finite values at $`k=0`$. For large $`k`$, the functions $`\nu ^4k^2E_f(\nu ,k)`$ shown here tend to a constant thus the integral over $`k`$ is convergent. All integrals have been truncated at $`k=1500`$. The error caused by this is quite small and does not change the results shown in the Table 1. These integrals we denote by $`E_f(\nu )`$. They are shown as function of $`\nu `$ in Fig.3 in a logarithmic scale. Again, we multiplied by a power of argument, here by $`\nu ^2`$, in order to make the behavior for large $`\nu `$ visible. It is seen that $`\nu ^2E_f(\nu )`$ tends to a constant so that the sum over $`\nu `$ is convergent. The sum is taken up to $`\nu =232.5`$ and again the remainder is small.

The calculations have been performed for several values of the parameters. The results are displayed in Table.1. The computations are performed with an adapted arithmetical precision. In intermediate steps compensations between sometimes very large quantities appeared. The precision was adapted accordingly. For example, for $`R_1=0.99,\nu =250.5,k=1600`$ as much as 1404 decimal positions have been necessary to get at least 16 digits precision of the integrand $`E_f(\nu ,k)`$. This was a factor causing large computation time.

Table 1. The numerical evaluations for several values of $`R_1`$ and $`\delta `$
$`R_1`$=0.0001
$`\delta `$ $`E^f`$ $`E^{as}`$ $`E^{ren}`$ $`E^{class}`$ $`E^{tot}`$    0.5 -0.0130152 0.00152886 -0.0114863 3.14159 3.13011  1. -0.052403 0.00509591 -0.0473071 12.5664 12.5191  3. -0.450093 -0.052012 -0.502105 113.097 112.595  6. -0.91953 -1.52936 -2.44889 452.389 449.94  10. 4.85954 -12.9483 -8.08871 1256.64 1248.55  15. 46.3333 -67.3661 -21.0328 2827.43 2806.4  21. 215.273 -261.526 -46.2534 5541.77 5495.52  30. 987.62 -1095.29 -107.666 11309.7 11202.1
$`R_1`$=0.9
$`\delta `$ $`E^f`$ $`E^{as}`$ $`E^{ren}`$ $`E^{class}`$ $`E^{tot}`$    0.5 -0.278679 -0.00105651 -0.279736 16.5347 16.255  1. -1.11431 -0.00451479 -1.11882 66.1388 65.02  3. -10.0633 -0.0683525 -10.1317 595.249 585.117  6. -40.4372 -0.647622 -41.0848 2381. 2339.92  10. -112.17 -4.2629 -116.433 6613.88 6497.45  20. -423.279 -63.2506 -486.53 26455.5 25969.  40. -1079.05 -992.187 -2071.24 105822. 103751.
$`R_1`$=0.95
$`\delta `$ $`E^f`$ $`E^{as}`$ $`E^{ren}`$ $`E^{class}`$ $`E^{tot}`$    0.5 -0.72166 -0.00212769 -0.723788 32.2215 31.4977  1. -2.86633 -0.00877711 -2.87511 128.886 126.011  3. -25.8293 -0.104564 -25.9339 1159.97 1134.04  10. -287.834 -4.39364 -292.228 12888.6 12596.4
$`R_1`$=0.99
$`\delta `$ $`E^f`$ $`E^{as}`$ $`E^{ren}`$ $`E^{class}`$ $`E^{tot}`$    0.5 -5.59211 -0.0106385 -5.60275 157.869 152.266  1. -23.2975 -0.0428038 -23.3403 631.476 608.136  3. -202.057 -0.409235 -202.466 5683.28 5480.82
$`R_1`$=0.997
$`\delta `$ $`E^f`$ $`E^{as}`$ $`E^{ren}`$ $`E^{class}`$ $`E^{tot}`$    1. -91.996 -0.142036 -92.138 2097.54 2005.4
$`R_1`$=0.999
$`\delta `$ $`E^f`$ $`E^{as}`$ $`E^{ren}`$ $`E^{class}`$ $`E^{tot}`$    1. -311.182 -0.425555 -311.608 6286.33 5974.72

## 8  Conclusions and Discussions

In the preceding sections the ground state energy for a spinor in the background of a rectangular shaped flux tube had been numerically calculated. The corresponding Jost function had been written down explicitly, also its asymptotic part. The numerical calculation required work with high arithmetic precision. The results are displayed mainly in Table 1. For small inner radius of the flux the results are close to them of  where the same problem for a flux tube with homogeneous magnetic field inside, which corresponds to $`R_1=0`$ here, was considered. Especially, it is seen that for large flux $`\delta `$ there is a compensation of the $`\delta ^4`$-contribution between the finite and the asymptotic parts of the ground state energy leaving a behaviour proportional to $`\delta ^2\mathrm{ln}\delta `$ as shown in Fig.4. Here the asymptotic part gave an essential contribution. The ground state energy remains negative, but numerically small. Only for very large flux it could overturn the corresponding classical energy, but these values of the flux are clearly unphysical.

For values of the inner radius close to the outer one, $`R_11`$ (where we have put $`R_2=1`$) the picture changes. Here the vacuum energy grows faster than the classical one. Generally, both diverge proportional to $`(1R_1)^1`$, the classical energy is equal to $`E_{\mathrm{class}}=\delta ^22\pi (1R_1)^1`$. The vacuum energy, multiplied by $`(1R_1)`$, is shown in Fig.5 in a logarithmic scale. It is negative and growing a bit faster than the classical one which would be a constant in this plot. Here the asymptotic part of the ground state energy becomes increasingly unimportant (see Table 1). The question whether the vacuum energy for sufficiently small $`(1R_1)`$ may become larger than the classical one cannot be answered by the numerical results obtained. The problem is that the computations become very time consuming because of the increasing precision which is required. Also, one has to take higher $`k`$ and $`\nu `$ into account. The weakening of the growth for $`R_1=0.997`$ and $`R_1=0.999`$ seen in Fig.5 may be caused just by this reason. Here one has to note that the integrand is for large $`k`$ and $`\nu `$ always negative (see Fig.2) so that dropping some part (as we did within the numerical procedure) diminishes the result. So, as a result, we cannot exclude from the given calculation that the vacuum energy grows for a strong background field faster than the classical energy.

Further work is necessary to better understand these questions. An improvement of the numerical procedure is certainly desirable. It could go along two lines. First, in the calculation of the Jost function the compensation of large exponentials should be avoided by taking them into account analytically. Second, in the compensation between the logarithm of the Jost function and its asymptotic expansion in the integrand of $`E_f`$ in Eq.(5.7,5.8) higher orders of the asymptotic expansion could be used. However, for this reason one would have to continue the procedure invented in  for this expansion using the Lippman-Schwinger equation or to find some equivalent procedure.

## Acknowledgments

I thank M. Bordag for advice.
I thank D.V. Vassilevich and K. Kirsten for valuable discussions and helpful remarks.
I thank the Graduate college Quantenfeldtheorie at the University of Leipzig for support and friendly environment.

## 9 Appendix

The sum over $`v`$ has been transformed to integrals using the Abel-Plana formula as follows:

$$\underset{l=0}{\overset{\mathrm{}}{}}(l+\frac{1}{2})=\underset{0}{\overset{\mathrm{}}{}}𝑑\nu f(\nu )+\underset{0}{\overset{\mathrm{}}{}}\frac{d\nu }{1+e^{2\pi \nu }}\frac{f(i\nu )f(i\nu )}{i}$$
(9.1)
The integrations over $`\nu `$ and $`k`$ can be done using identities:

$$\underset{0}{\overset{\mathrm{}}{}}𝑑\nu \underset{m}{\overset{\mathrm{}}{}}𝑑k(k^2m^2)^{1s}\frac{}{k}\frac{t^j}{\nu ^n}=\frac{m^{22s}}{2}\frac{\mathrm{\Gamma }(2s)\mathrm{\Gamma }(\frac{1+jn}{2})\mathrm{\Gamma }(s+\frac{n3}{2})}{(rm)^{n1}\mathrm{\Gamma }(\frac{j}{2})}$$
(9.2)
$$\underset{m}{\overset{\mathrm{}}{}}𝑑k(k^2m^2)^{1s}\frac{}{k}t^j=m^{22s}\frac{\mathrm{\Gamma }(2s)\mathrm{\Gamma }(s+\frac{j}{2}1)}{\mathrm{\Gamma }(\frac{j}{2})}\frac{(\frac{\nu }{rm})^j}{(1+(\frac{\nu }{rm})^2)^{s+\frac{j}{2}1}}$$
(9.3)
The expansion in powers of $`\nu `$ for logarithm of asymptotic Jost function can be obtained in the form (see)

$$lnf_{as}(ik)=\underset{n=1}{\overset{3}{}}\underset{j=1}{\overset{9}{}}\underset{0}{\overset{\mathrm{}}{}}\frac{dr}{r}X_{nj}\frac{t^j}{\nu ^n}$$
(9.4)
where
$`X_{1,1}=X_{1,3}=X_{2,6}=\frac{1}{2}(a\delta )^2`$
$`X_{2,2}=\frac{1}{4}\delta ^2(a^2raa^{})`$
$`X_{2,4}=\frac{1}{4}\delta ^2(3a^2+raa^{})`$
$`X_{3,3}=\frac{1}{4}\delta ^2(a^2raa^{}+\frac{1}{2}r^2aa^{\prime \prime }\frac{1}{2}\delta ^2a^4)`$
$`X_{3,5}=\frac{1}{8}\delta ^2(\frac{39}{2}a^2+7raa^{}r^2aa^{\prime \prime }+6\delta ^2a^4)`$
$`X_{3,7}=\frac{1}{8}\delta ^2(35a^25raa^{}5\delta ^2a^4)`$
$`X_{3,9}=\frac{35}{16}\delta ^2a^2`$

For the representation of $`E_{as}`$ as 5.6

$$E^{as}=\frac{4}{\pi }\underset{0}{\overset{\mathrm{}}{}}\frac{dr}{r^3}[\delta ^2a(r)^2g_1(rm)\delta ^2r^2a^{}(r)^2g_2(rm)+\delta ^4a(r)^4g_3(rm)]$$
(9.5)
here $`f_i`$ are

$`f_1(x)`$ $`=`$ $`{\displaystyle \frac{1}{2}}f_{1,1}(x){\displaystyle \frac{1}{2}}f_{1,3}(x)+{\displaystyle \frac{1}{4}}f_{3,3}{\displaystyle \frac{39}{16}}f_{3,5}(x)+{\displaystyle \frac{35}{8}}f_{3,7}(x){\displaystyle \frac{35}{16}}f_{3,9}(x)`$
$``$ $`{\displaystyle \frac{1}{2}}x_x({\displaystyle \frac{1}{4}}f_{3,3}(x)+{\displaystyle \frac{7}{8}}f_{3,5}(x){\displaystyle \frac{5}{8}}f_{3,7}(x))`$
$`+`$ $`{\displaystyle \frac{1}{2}}x_x^2({\displaystyle \frac{x}{8}}f_{3,3}(x){\displaystyle \frac{x}{8}}f_{3,5}(x)),`$
$`f_2(x)`$ $`=`$ $`{\displaystyle \frac{1}{8}}(f_{3,3}(x)f_{3,5}(x)),`$
$`f_3(x)`$ $`=`$ $`{\displaystyle \frac{1}{8}}(f_{3,3}(x)6f_{3,5}(x)+5f_{3,7}(x)).`$ (9.6)
with $`f_{i.j}`$ are

$`f_{1,1}(x)`$ $`=`$ $`{\displaystyle \frac{1}{1+e^{2\pi x}}}`$
$`f_{1,3}(x)`$ $`=`$ $`({\displaystyle \frac{}{1+e^{2\pi x}}})^{}`$
$`f_{3,3}(x)`$ $`=`$ $`({\displaystyle \frac{1}{x}}{\displaystyle \frac{1}{1+e^{2\pi x}}})^{}`$
$`f_{3,5}(x)`$ $`=`$ $`{\displaystyle \frac{1}{3}}({\displaystyle \frac{1}{x}}({\displaystyle \frac{x}{1+e^{2\pi x}}})^{})^{}`$
$`f_{3,7}(x)`$ $`=`$ $`{\displaystyle \frac{1}{15}}({\displaystyle \frac{1}{x}}({\displaystyle \frac{1}{x}}({\displaystyle \frac{x^3}{1+e^{2\pi x}}})^{})^{})^{}`$
$`f_{3,9}(x)`$ $`=`$ $`{\displaystyle \frac{1}{105}}({\displaystyle \frac{1}{x}}({\displaystyle \frac{1}{x}}({\displaystyle \frac{1}{x}}({\displaystyle \frac{x^5}{1+e^{2\pi x}}})^{})^{})^{})^{}`$ (9.7)
The asymtotic of the logarithmic Jost function can be obtained in the form: $`\mathrm{ln}f_\nu (ik)=\underset{n}{}h_n(t_1,t_2)`$,(see 6.4) with:

$`h_1`$ $`=`$ $`{\displaystyle \frac{4}{3}}\lambda ^2[R_{1}^{}{}_{}{}^{4}t_1(2+t_1)(1+t_2)^23R_{1}^{}{}_{}{}^{2}R_{2}^{}{}_{}{}^{2}(1+t_1)^2t_2(1+t_2)+R_{2}^{}{}_{}{}^{4}(1+t_1)^2t_2(1+2t_2)]`$
$`h_2`$ $`=`$ $`0`$
$`h_3`$ $`=`$ $`{\displaystyle \frac{1}{6}}\lambda ^2[2R_{1}^{}{}_{}{}^{4}t_{1}^{}{}_{}{}^{3}(1+t_2)^2+R_{1}^{}{}_{}{}^{2}R_{2}^{}{}_{}{}^{2}(1+t_1)^2t_{2}^{}{}_{}{}^{3}(1+t_2)^2(3t_{2}^{}{}_{}{}^{2}4)`$
$`R_{2}^{}{}_{}{}^{4}(1+t_1)^2t_{2}^{}{}_{}{}^{3}(3t_{2}^{}{}_{}{}^{4}+6t_{2}^{}{}_{}{}^{3}t_{2}^{}{}_{}{}^{2}8t_22)]+`$
$`{\displaystyle \frac{2}{15}}\lambda ^4[4R_{1}^{}{}_{}{}^{8}t_{1}^{}{}_{}{}^{3}(4+t_1)(1+t_2)^45R_{1}^{}{}_{}{}^{6}R_{2}^{}{}_{}{}^{2}(1+t_1)^4t_{2}^{}{}_{}{}^{3}(1+t_2)^4+`$
$`15R_{1}^{}{}_{}{}^{4}R_{2}^{}{}_{}{}^{4}(1+t_1)^4t_{2}^{}{}_{}{}^{3}(1+t_2)^2(t_{2}^{}{}_{}{}^{2}+2t_21)`$
$`5R_{1}^{}{}_{}{}^{2}R_{2}^{}{}_{}{}^{6}(1+t_1)^4t_{2}^{}{}_{}{}^{3}(3t_{2}^{}{}_{}{}^{4}+12t_{2}^{}{}_{}{}^{3}+6t_{2}^{}{}_{}{}^{2}4t_21)+`$
$`R_{2}^{}{}_{}{}^{8}(1+t_1)^4t_{2}^{}{}_{}{}^{3}(5t_{2}^{}{}_{}{}^{4}+20t_{2}^{}{}_{}{}^{3}4t_21)]`$
$`h_4`$ $`=`$ $`{\displaystyle \frac{\delta ^2[R_{1}^{}{}_{}{}^{4}t_{1}^{}{}_{}{}^{4}(1t_{1}^{}{}_{}{}^{2})+R_{2}^{}{}_{}{}^{4}t_{2}^{}{}_{}{}^{4}(1t_{2}^{}{}_{}{}^{2})]}{4(R_1^2R_2^2)^2}}`$
$`\lambda `$ $`=`$ $`{\displaystyle \frac{\delta }{(R_1^2R_2^2)(1+t_1)(1+t_2)}}`$ (9.8)


format: text

text:
# Comparing gas-phase and grain-catalyzed H₂ formation

## 1 Introduction

In local molecular clouds, molecular hydrogen ($`\mathrm{H}_2`$) forms primarily on the surface of dust grains: two hydrogen atoms are adsorbed onto the surface of the grain and react to form $`\mathrm{H}_2`$, which subsequently escapes back into the interstellar medium. However, $`\mathrm{H}_2`$ can also form in the gas-phase, primarily through the reactions

$`\mathrm{H}+\mathrm{e}`$ $``$ $`\mathrm{H}^{}+\gamma `$ (1)
$`\mathrm{H}^{}+\mathrm{H}`$ $``$ $`\mathrm{H}_2+\mathrm{e},`$ (2)
although some also forms via the slower reactions

$`\mathrm{H}+\mathrm{H}^+`$ $``$ $`\mathrm{H}_2^++\gamma `$ (3)
$`\mathrm{H}_2^++\mathrm{H}`$ $``$ $`\mathrm{H}_2+\mathrm{H}^+.`$ (4)
In dense gas, three-body reactions can also be important (Palla, Salpeter & Stahler, 1983), but these are ineffective at number densities $`n<10^8\mathrm{cm}^3`$.

Discussions of $`\mathrm{H}_2`$ formation have tended to concentrate on the role played by dust, with little attention given to the gas-phase reactions. However, as I show in section 4, in some circumstances these reactions can dominate the $`\mathrm{H}_2`$ formation rate.

In this paper, I briefly outline the chemistry of $`\mathrm{H}_2`$ formation and show that it is easy to identify a critical dust-to-gas ratio $`𝒟_{\mathrm{cr}}`$, above which grain-catalyzed formation dominates. I illustrate the method by applying it to various situations of astrophysical interest, and show that it can be a useful tool for estimating the importance of gas-phase $`\mathrm{H}_2`$ formation.

## 2 The formation of molecular hydrogen

### 2.1 Gas-phase formation

Most of the molecular hydrogen that forms in the gas-phase does so via the formation of an intermediate $`\mathrm{H}^{}`$ ion, as outlined in reactions 1 and 2 above. The first of these reactions occurs much more slowly than the second, and so the equilibrium abundance of $`\mathrm{H}^{}`$ is small and is rapidly reached. Thereafter, the $`\mathrm{H}_2`$ formation rate is determined by two factors: the rate at which $`\mathrm{H}^{}`$ forms, and the fraction of $`\mathrm{H}^{}`$ ions that survive to form $`\mathrm{H}_2`$. The latter quantity is determined by competition between $`\mathrm{H}_2`$ formation via reaction 2 and $`\mathrm{H}^{}`$ destruction by mutual neutralization with $`\mathrm{H}^+`$ ions

$$\mathrm{H}^{}+\mathrm{H}^+2\mathrm{H},$$
(5)
and by photodetachment by the incident radiation field

$$\mathrm{H}^{}+\gamma \mathrm{H}+\mathrm{e}.$$
(6)
Various other reactions also destroy $`\mathrm{H}^{}`$, but these are either significantly slower than those above, or become important only at high temperatures, in which case any $`\mathrm{H}_2`$ that does form will very quickly be collisionally dissociated. For more details, the reader is referred to the recent reviews of Abel et al. (1997), Galli & Palla (1998), Stancil et al. (1998) and Lepp et al. (2002).

If we assume, for simplicity, that $`\mathrm{H}^{}`$ has already reached its equilibrium abundance, then we can write the $`\mathrm{H}_2`$ formation rate as

$$R_{\mathrm{H}_2,\mathrm{H}^{}}=k_1n_\mathrm{e}n_\mathrm{H}\frac{k_2n_\mathrm{H}}{k_2n_\mathrm{H}+k_5n_{\mathrm{H}^+}+k_6}$$
(7)
where $`n_i`$ is the number density of species $`i`$, and where the rate coefficients $`k_i`$ for the various reactions are listed in table 1.

If $`\mathrm{H}_2`$ formation via reaction 2 occurs much faster than the destruction of $`\mathrm{H}^{}`$ by the other reactions, then this reduces to

$$R_{\mathrm{H}_2,\mathrm{H}^{}}k_1n_\mathrm{e}n_\mathrm{H};$$
(8)
in other words, the $`\mathrm{H}_2`$ formation rate is approximately the same as the $`\mathrm{H}^{}`$ formation rate.

If, on the other hand, mutual neutralization dominates over $`\mathrm{H}_2`$ formation or photodetachment as a means of removing $`\mathrm{H}^{}`$, then equation 7 becomes

$$R_{\mathrm{H}_2,\mathrm{H}^{}}k_1n_\mathrm{e}n_\mathrm{H}\frac{k_2}{k_5x},$$
(9)
where $`x=n_{\mathrm{H}^+}/n_\mathrm{H}`$ is the fractional ionization of hydrogen. As long as $`n_\mathrm{e}n_{\mathrm{H}^+}`$, this equation can be further simplified to

$$R_{\mathrm{H}_2,\mathrm{H}^{}}\frac{k_1k_2}{k_5}n_\mathrm{H}^2.$$
(10)
Comparing this equation with equation 8, we see that for a small fractional ionization $`R_{\mathrm{H}_2}x`$, but that once the fractional ionization becomes large enough that mutual neutralization dominates, $`R_{\mathrm{H}_2}`$ becomes independent of the ionization: although increases in $`x`$ still increase the $`\mathrm{H}^{}`$ formation rate, this is balanced by the increase in the mutual neutralization rate and consequent decrease in the fraction of $`\mathrm{H}^{}`$ ions surviving to form $`\mathrm{H}_2`$. This change in behaviour occurs for fractional ionizations near a critical value $`x_{\mathrm{cr}}`$, defined by

$$x_{\mathrm{cr}}=\frac{k_2}{k_5}.$$
(11)
The precise value of $`x_{\mathrm{cr}}`$ is somewhat uncertain, due to the significant uncertainty that remains in the determination of the mutual neutralization rate. In this paper, I have chosen to adopt the rate listed in Galli & Palla (1998), which is derived from the data of Moseley et al. (1970). This is a conservative choice, in that it gives the lowest value of $`x_{\mathrm{cr}}`$; other possibilities include the rates of Duley & Williams (1984), Dalgarno & Lepp (1987) and Croft et al. (1999), with the last-named being preferred by the most recent compilation (Lepp et al., 2002). For the temperature range of interest, the Galli & Palla rate gives us a value $`x_{\mathrm{cr}}5\times 10^3`$, with only a slight dependence on temperature. The alternative rates typically give values of $`x_{\mathrm{cr}}`$ that are factors of a few larger.

When the destruction rate of $`\mathrm{H}^{}`$ ions is dominated by photodetachment, we obtain another limiting case of equation 7

$$R_{\mathrm{H}_2,\mathrm{H}^{}}k_1n_\mathrm{e}n_\mathrm{H}\frac{k_2n_\mathrm{H}}{k_6}.$$
(12)
This can be written as the $`\mathrm{H}^{}`$ formation rate divided by a suppression factor $`f_{\mathrm{rad}}`$:

$$R_{\mathrm{H}_2,\mathrm{H}^{}}\frac{k_1n_\mathrm{e}n_\mathrm{H}}{f_{\mathrm{rad}}},$$
(13)
where

$$f_{\mathrm{rad}}=\frac{k_6}{k_2n_\mathrm{H}}.$$
(14)
Determination of the photodetachment rate, and hence $`f_{\mathrm{rad}}`$, requires knowledge of the incident radiation field. Provided that the opacity of the gas is low, we can write the photodetachment rate as

$$k_6=4\pi _{\nu _{\mathrm{th}}}^{\mathrm{}}\frac{\sigma _\nu J_\nu }{h\nu }d\nu ,$$
(15)
where $`J_\nu `$ is the mean specific intensity, $`\sigma _\nu `$ is the photodetachment cross-section (de Jong, 1972)

$$\sigma _\nu =7.928\times 10^5\frac{(\nu \nu _{\mathrm{th}})^{3/2}}{\nu ^3}\mathrm{cm}^2,$$
(16)
and where $`h\nu _{\mathrm{th}}=0.755\mathrm{eV}`$ is the energy threshold for $`\mathrm{H}^{}`$ photodetachment.

Evaluating equation 15 for the local interstellar radiation field, as estimated by Mathis et al. (1983), gives

$$k_6=4.2\times 10^7\mathrm{s}^1,$$
(17)
and so, locally,

$$f_{\mathrm{rad}}=3.2\times 10^2n_\mathrm{H}^1.$$
(18)
If the opacity of the gas is high, then this will overestimate the effects of radiation. However, at frequencies near the $`\mathrm{H}^{}`$ photodetachment threshold, the continuum opacity of interstellar gas is low and absorption is dominated by dust. Consequently, a high opacity implies a high dust content, in which case grain-catalyzed formation will dominate.

These limiting cases provide useful insight into the physics of gas-phase $`\mathrm{H}_2`$ formation, but in general we must use the full form of equation 7, which we can rewrite as

$$R_{\mathrm{H}_2,\mathrm{H}^{}}=\frac{k_1n_\mathrm{e}n_\mathrm{H}}{1+x/x_{\mathrm{cr}}+f_{\mathrm{rad}}}.$$
(19)
If we now turn to $`\mathrm{H}_2`$ formation via the $`\mathrm{H}_2^+`$ ion, we find that the basic chemistry is remarkably similar. $`\mathrm{H}_2^+`$ is created by the radiative association of $`\mathrm{H}`$ and $`\mathrm{H}^+`$ (reaction 3), and destroyed by $`\mathrm{H}_2`$ formation (reaction 4), dissociative recombination

$$\mathrm{H}_2^++\mathrm{e}2\mathrm{H},$$
(20)
and photodissociation

$$\mathrm{H}_2^++\gamma \mathrm{H}+\mathrm{H}^+.$$
(21)
As with $`\mathrm{H}^{}`$, the formation of the molecular ion is the limiting step, with subsequent reactions occurring orders of magnitude faster. If we again assume that the $`\mathrm{H}_2^+`$ abundance has reached equilibrium, then we can write the $`\mathrm{H}_2`$ formation rate as

$$R_{\mathrm{H}_2,\mathrm{H}_2^+}=k_3n_{\mathrm{H}^+}n_\mathrm{H}\frac{k_4n_\mathrm{H}}{k_4n_\mathrm{H}+k_7n_\mathrm{e}+k_8},$$
(22)
which has the same form as equation 7. Indeed, we can rewrite it as

$$R_{\mathrm{H}_2,\mathrm{H}_2^+}=\frac{k_3n_{\mathrm{H}^+}n_\mathrm{H}}{1+x/x_{\mathrm{cr}}+f_{\mathrm{rad}}},$$
(23)
only now

$$x_{\mathrm{cr}}=\frac{k_4}{k_7},$$
(24)
and

$$f_{\mathrm{rad}}=\frac{k_8}{k_4n_\mathrm{H}},$$
(25)
where the photodissociation rate, $`k_8`$, is calculated in a similar fashion to the photodetachment rate above. Evaluating these, we find that $`x_{\mathrm{cr}}`$ for $`\mathrm{H}_2^+`$ is typically an order of magnitude larger than for $`\mathrm{H}^{}`$, and that in the local ISM

$$f_{\mathrm{rad}}n_\mathrm{H}^1,$$
(26)
where I have again used the Mathis et al. (1983) radiation field, together with the $`\mathrm{H}_2^+`$ photodissociation cross-section from Shapiro & Kang (1987).

Comparing these values with those for $`\mathrm{H}^{}`$, we see that $`\mathrm{H}_2^+`$ is significantly more robust. However, it forms at a much slower rate (between two and three orders of magnitude, depending on temperature) and so in most cases $`\mathrm{H}^{}`$ dominates. Nevertheless, there are exceptions, as we will see in section 4.

Finally, a few other possible mechanisms have been suggested for gas-phase $`\mathrm{H}_2`$ formation. Latter & Black (1991) propose that $`\mathrm{H}_2`$ can form as a result of direct radiative association

$$\mathrm{H}+\mathrm{H}\mathrm{H}_2+\gamma ,$$
(27)
provided that one of the hydrogen atoms is in an excited electronic state. Rawlings et al. (1993) show that a more efficient mechanism is formation of $`\mathrm{H}_2^+`$ by associative ionization

$$\mathrm{H}+\mathrm{H}\mathrm{H}_2^++\mathrm{e}$$
(28)
with the $`\mathrm{H}_2^+`$ thereafter forming $`\mathrm{H}_2`$ by reaction 4 above. This mechanism again requires one of the hydrogen atoms to be in an excited atomic state. However, this requirement means that in general these reactions are not important, as the necessary population of excited atomic hydrogen is only found in a few unusual circumstances, such as in the universe immediately after recombination.

### 2.2 Grain-catalyzed $`\mathrm{H}_2`$ formation

Despite its importance in local interstellar chemistry, the rate of $`\mathrm{H}_2`$ formation on dust grains is still uncertain. In local molecular clouds, observations suggest a formation rate (Jura, 1975)

$$k_93\times 10^{17}\mathrm{cm}^3\mathrm{s}^1.$$
(29)
Observations of $`\mathrm{H}_2`$ in the LMC and SMC with the fuse satellite (Tumlinson et al., 2002) suggest a value that is an order of magnitude smaller, but this is consistent with the underlying rate per unit dust mass being the same, since the mean dust-to-gas ratio within these galaxies is significantly smaller than in the Milky Way (Issa et al., 1990).

Unfortunately, direct measurements of this kind can only give us information about $`\mathrm{H}_2`$ formation in physical conditions that are easily accessible to observations, and provide little basis on which to predict the $`\mathrm{H}_2`$ formation rate in different regimes. For this, we must turn to theory.

A large body of theoretical work exists on the subject of $`\mathrm{H}_2`$ formation on grains (see, for example, the review of Pirronello et al., 2000, and the many references therein), stretching back almost forty years to the pioneering work of Gould and Salpeter (1963). In a highly influential paper, Hollenbach et al. (1971) parameterized the $`\mathrm{H}_2`$ formation rate as

$$R=0.5\overline{v}_\mathrm{H}\sigma _\mathrm{d}Sf_an_\mathrm{H}n_\mathrm{d},$$
(30)
where $`n_\mathrm{d}`$ is the number density of dust grains, $`\sigma _\mathrm{d}`$ is their mean geometric cross-section, $`\overline{v}_\mathrm{H}`$ is the mean velocity of the hydrogen atoms striking the grains, $`S`$ is the sticking coefficient (the probability that a hydrogen atom striking the grain will stick to the surface) and $`f_a`$ is the fraction of adsorbed hydrogen atoms that actually form $`\mathrm{H}_2`$, rather than simply escaping back into the gas phase. They argued that for gas and grain temperatures typical of molecular clouds, both $`S`$ and $`f_a`$ should be of order unity.

Hollenbach & Mckee (1979) later used this prescription to derive an $`\mathrm{H}_2`$ formation rate for the local ISM that continues to be widely cited:

$$k_9=3\times 10^{17}\frac{T_2^{0.5}f_a}{y}\mathrm{cm}^3\mathrm{s}^1,$$
(31)
where

$$y=1+0.4\sqrt{T_2+T_{\mathrm{gr},2}}+0.2T_2+0.08T_2^2,$$
(32)
and where $`T_2`$ and $`T_{\mathrm{gr},2}`$ are the gas and grain temperatures in units of $`100\mathrm{K}`$. They argued that $`f_a`$ should be approximately constant and of order unity for grain temperatures below some critical value $`T_{\mathrm{cr}}`$, but that for $`T_{\mathrm{gr}}>T_{\mathrm{cr}}`$, it should fall off exponentially, with most of the hydrogen atoms evaporating from the grain surface before they have time to form $`\mathrm{H}_2`$. The value of $`T_{\mathrm{cr}}`$ has proved hard to determine precisely, but is of the order of $`100\mathrm{K}`$.

Although this rate has been widely adopted in the literature, recent experiments have cast doubt on its validity at high temperatures, and suggest that the $`\mathrm{H}_2`$ formation rate may be smaller than previously assumed (Pirronello et al. 1997a,b, 1999; Katz et al. 1999; Biham et al. 2001). However, since this conclusion is not entirely clear (Cazaux & Tielens, 2002) and their work is still ongoing, I have tentatively adopted the Hollenbach & Mckee rate below, with the proviso that the values of $`𝒟_{\mathrm{cr}}`$ that I derive may prove to be lower limits if the results of Pirronello et al. are borne out by future work.

## 3 Comparing the different modes of formation

We can combine equations 19 and 23 to write the total gas-phase $`\mathrm{H}_2`$ formation rate as

$$R_{\mathrm{H}_2,\mathrm{gas}}=R_{\mathrm{H}_2,\mathrm{H}^{}}+R_{\mathrm{H}_2,\mathrm{H}_2^+},$$
(33)
while the grain-catalyzed rate can be written as

$$R_{\mathrm{H}_2,\mathrm{dust}}=k_9n_{\mathrm{tot}}n_\mathrm{H}\left(\frac{𝒟}{𝒟_{\mathrm{MW}}}\right),$$
(34)
where $`n_{\mathrm{tot}}`$ is the total particle number density, $`𝒟`$ is the dust-to-gas ratio and $`𝒟_{\mathrm{MW}}`$ is its value in the local ISM.

Combining these equations, we can easily solve for the dust-to-gas ratio at which $`R_{\mathrm{H}_2,\mathrm{gas}}`$ and $`R_{\mathrm{H}_2,\mathrm{dust}}`$, which I denote as $`𝒟_{\mathrm{cr}}`$:

$$𝒟_{\mathrm{cr}}=\frac{R_{\mathrm{H}_2,\mathrm{gas}}}{k_9n_{\mathrm{tot}}n_\mathrm{H}}𝒟_{\mathrm{MW}}.$$
(35)
In the common case that $`R_{\mathrm{H}_2,\mathrm{H}^{}}R_{\mathrm{H}_2,\mathrm{H}_2^+}`$, this equation reduces to

$$𝒟_{\mathrm{cr}}=\frac{k_1}{k_9}\frac{x}{1+x/x_{\mathrm{cr}}+f_{\mathrm{rad}}}𝒟_{\mathrm{MW}},$$
(36)
where $`x_{\mathrm{cr}}`$ and $`f_{\mathrm{rad}}`$ are given by equations 11 and 14 respectively. A similar equation can be written in the much less common case that $`R_{\mathrm{H}_2,\mathrm{H}_2^+}R_{\mathrm{H}_2,\mathrm{H}^{}}`$.

In order to help illustrate the behaviour of these equations, I plot in figures 1 to 4 the value of $`𝒟_{\mathrm{cr}}`$ as a function of temperature for gas illuminated by the Mathis et al. (1983) radiation field in four different scenarios: low ionization, low density gas ($`x=10^4,n_\mathrm{H}=1\mathrm{cm}^3`$; figure 1), low ionization, high density gas ($`x=10^4,n_\mathrm{H}=10^3\mathrm{cm}^3`$; figure 2), high ionization, low density gas ($`x=10^2,n_\mathrm{H}=1\mathrm{cm}^3`$; figure 3) and high ionization, high density gas ($`x=10^2,n_\mathrm{H}=10^3\mathrm{cm}^3`$; figure 4). In each case, I adopt a fixed grain temperature $`T_{\mathrm{gr}}=20\mathrm{K}`$, although small changes in $`T_{\mathrm{gr}}`$ have little effect on the results provided that it remains less than $`T_{\mathrm{cr}}`$.

A striking feature of these plots is the strong temperature dependence of $`𝒟_{\mathrm{cr}}`$. At low temperatures, grain-catalyzed $`\mathrm{H}_2`$ formation is relatively efficient and very little dust is needed before grain catalysis dominates. Above a few hundred K, however, the efficiency of grain catalysis decreases significantly, while the efficiency of gas-phase $`\mathrm{H}_2`$ formation continues to grow. As a result, the required dust abundance rises sharply with increasing temperature.

## 4 Astrophysical examples

From the behaviour outlined in figures 14, it is clear that gas-phase $`\mathrm{H}_2`$ formation is at its most effective in warm, dense gas with a high fractional ionization. However, most of the molecular gas that we observe in our galaxy is in the form of molecular clouds with low temperatures ($`T20\mathrm{K}`$) and very low fractional ionizations ($`x10^7`$) and in these conditions grain catalyzed formation dominates by many orders of magnitude.

A more promising place to look for gas-phase $`\mathrm{H}_2`$ formation is in the so-called warm neutral medium (WNM). In models of the multiphase ISM that assume thermal pressure equilibrium between phases (Field et al., 1969; McKee & Ostriker, 1977; Wolfire et al., 1995), this is predicted to have a temperature of approximately $`8000\mathrm{K}`$, high enough to collisionally dissociate $`\mathrm{H}_2`$. However, recent observations (Heiles & Troland, 2002) and simulations that include the effects of turbulence (Gazol et al., 2001; Mac Low et al., 2001) suggest that much of this gas is actually at much lower temperatures; for instance, Heiles & Troland quote a temperature range of $`500<T<5000\mathrm{K}`$.

Taking representative values for the temperature and ionization of the WNM to be $`T=2000\mathrm{K}`$ and $`x=10^2`$ (Heiles, 2001), I find that

$$\frac{𝒟_{\mathrm{cr}}}{𝒟_{\mathrm{MW}}}=\frac{4.1}{2.7+f_{\mathrm{rad},\mathrm{H}^{}}}+\frac{0.047}{1.1+f_{\mathrm{rad},\mathrm{H}_2^+}}.$$
(37)
In other words, gas-phase formation would dominate if we could ignore the effects of the radiation field. In practice, this is not possible; at the densities characteristic of the WNM ($`n0.1\mathrm{cm}^3`$), we have

$`f_{\mathrm{rad},\mathrm{H}^{}}`$ $`=`$ $`3.9\times 10^3`$ (38)
$`f_{\mathrm{rad},\mathrm{H}_2^+}`$ $`=`$ $`1.0\times 10^2,`$ (39)
and equation 37 becomes

$$𝒟_{\mathrm{cr}}=5.2\times 10^3𝒟_{\mathrm{MW}}.$$
(40)
These two examples demonstrate that gas-phase $`\mathrm{H}_2`$ formation is unimportant in the bulk of the gas in the Milky Way: either the temperature and ionization are too low, as in molecular clouds, or the gas is too diffuse and $`\mathrm{H}_2`$ formation is suppressed by the photodissociation of $`\mathrm{H}^{}`$ and $`\mathrm{H}_2^+`$.

However, there are a few counterexamples. For instance, gas-phase $`\mathrm{H}_2`$ formation has long been known to play an important role in the chemistry of nova ejecta (Rawlings, 1988) and protostellar outflows (Glassgold et al., 1989), where the gas initially has little or no dust (although more generally forms later) and where the high gas densities help mitigate the effects of photodissociation. Gas-phase formation is also predicted to dominate the molecular chemistry of freely-expanding supernova remnants such as SN1987A (Culhane & McCray, 1995). Finally, Lepp & McCray (1983) suggest that gas-phase formation may dominate in X-ray dissociation regions (dense clouds illuminated by hard X-rays); modelling by Maloney et al. (1996) would appear to confirm this.

These are somewhat unusual conditions, however, and in general dust abundances significantly below the typical galactic value are required before gas-phase $`\mathrm{H}_2`$ formation becomes competitive with grain-catalyzed formation.

One place in which we might expect to find these low dust abundances is in the metal-poor gas within dwarf galaxies. Kamaya & Hirashita (2001) examine a well-studied example, the metal-poor dwarf IZw18, and show that provided that its neutral ISM is clumpy (with clumps densities $`n100\mathrm{cm}^3`$) and moderately ionized ($`x10^3`$), then gas-phase formation will dominate. They also show that the formation of $`\mathrm{H}_2`$ in this manner would not conflict with the upper limit on the $`\mathrm{H}_2`$ column density of IZw18 obtained by fuse (Vidal-Madjar et al., 2000), due to the small filling factor of the clumps.

Another place we might look for significant gas-phase $`\mathrm{H}_2`$ formation is in damped Lyman-$`\alpha `$ (DLA) systems, many of which have low dust abundances (see, eg Lopez et al., 2002). An interesting example is the absorber at $`z=3.025`$ in the spectrum of Q0347-3819 recently studied by Levshakov et al. (2002). This system has a temperature $`T400\mathrm{K}`$ (as inferred from the Doppler broadening of its many associated $`\mathrm{H}_2`$ and metal absorption lines), and a fractional ionization $`x2\times 10^5`$. If we assume that $`\mathrm{H}^{}`$ photodetachment is negligible, we find that for this system

$$D_{\mathrm{cr}}=4.6\times 10^4𝒟_{\mathrm{MW}}.$$
(41)
Comparing this with a measured dust-to-gas ratio of $`𝒟0.05𝒟_{\mathrm{MW}}`$, we see that gas-phase $`\mathrm{H}_2`$ formation contributes no more than about 1% of the total $`\mathrm{H}_2`$ in this system. Including the effects of radiation merely strengthens this conclusion.

There is no reason to suspect that this situation is particularly unusual; all damped Lyman-$`\alpha `$ systems by definition have large H i column densities, and consequently will have small fractional ionizations. We would therefore expect grain-catalysis to dominate in these systems.

Finally, gas-phase $`\mathrm{H}_2`$ formation has long been known to play an important role in the early stages of galaxy formation. In primordial gas, this is obvious: there is no dust, so any $`\mathrm{H}_2`$ that forms must form in the gas phase. A more interesting problem is determining the value of $`𝒟_{\mathrm{cr}}`$ for these systems; in other words, at what point does grain-catalyzed formation overtake gas-phase formation?

For the purposes of this discussion, I adopt the example of an $`\mathrm{H}_2`$-cooled protogalaxy with temperature $`T=1000\mathrm{K}`$ and fractional ionization $`x=2\times 10^4`$ (Tegmark et al., 1997). These values are appropriate for the first generation of star-forming protogalaxies, and while they may be underestimates for later generations, my analysis can easily be rescaled for higher values. For this example protogalaxy, I find that

$$\frac{𝒟_{\mathrm{cr}}}{𝒟_{\mathrm{MW}}}=\frac{2.1\times 10^2}{1+f_{\mathrm{rad},\mathrm{H}^{}}}+\frac{1.3\times 10^4}{1+f_{\mathrm{rad},\mathrm{H}_2^+}}.$$
(42)
Thus, if the radiation field is unimportant, $`𝒟_{\mathrm{cr}}0.02𝒟_{\mathrm{MW}}`$, comparable to the values seen in some metal-poor dwarf galaxies at the present-day (Lisenfeld & Ferrara, 1998).

How strong is the radiation field within a protogalaxy? There are potentially three main contributors to this field: the protogalaxy’s own stellar population, emission from neighbouring galaxies and radiation from the cosmological background produced by distant sources. It is simplest to consider these separately.

Much of the optical and near-infrared radiation responsible for destroying $`\mathrm{H}^{}`$ and $`\mathrm{H}_2^+`$ is produced by long-lived stars, and so the contribution of the protogalaxy’s stellar population depends as much on its star formation history as on its current star formation rate. This makes it very difficult to parameterize its effects in the general case; it is much easier to examine a simple example that will hopefully be broadly representative.

For the purposes of this example, I assume:

1. That the protogalaxy underwent an instantaneous (or near-instantaneous) starburst $`10^8\mathrm{yr}`$ ago, following which it has formed no more stars.
2. That the stars which did form are located in the centre of the protogalaxy, within a small enough region that I can approximate their emission as coming from a point source.
3. That the luminosity and spectral energy distribution of this stellar cluster are well described by the $`Z=0.05Z_{}`$ model of Leitherer et al. (1999).
4. That the protogalaxy itself is well-described by a truncated isothermal sphere density profile (Iliev & Shapiro, 2001).

All of these assumptions are debatable, but they do provide us a basis on which to estimate the effects of the stellar radiation field. Moreover, these assumptions are somewhat conservative, and tend to minimize the effectiveness of the stellar radiation. For instance, if we reduce the time since the starburst from $`10^8\mathrm{yr}`$ to $`10^7\mathrm{yr}`$, then the photodetachment rate increases by a factor of fifty. Similarly, if we assume continuous star-formation rather than an instantaneous starburst, then we obtain a similar (or slightly larger) photodetachment rate once the total mass of stars formed has reached a comparable level.

We could also criticize the adoption of the Leitherer et al. model, on the basis that it assumes a standard Salpeter IMF, while there is considerable evidence that the the primordial IMF is biased towards high masses (Larson, 1998). However, this again means that we will underestimate the photodetachment rate (although we will significantly overestimate the lifetime of the stellar population).

Returning to my example, the first two assumptions allow me to write the $`\mathrm{H}^{}`$ photodetachment rate at a distance $`R`$ from the stars as

$$k_6=\frac{1}{R^2}_{\nu _{\mathrm{th}}}^{\nu _0}\frac{\sigma _\nu L_\nu }{h\nu }d\nu ,$$
(43)
where $`L_\nu `$, the stellar luminosity per unit frequency, is given by the Leitherer et al. model. Using this value, I obtain

$$k_6=4.9\times 10^8F(M_{},R)\mathrm{s}^1,$$
(44)
where

$$F(M_{},R)=\left(\frac{M_{}}{10^6\mathrm{M}_{}}\right)\left(\frac{R}{1\mathrm{kpc}}\right)^2,$$
(45)
and where $`M_{}`$ is the mass of stars formed in the starburst. Similarly, we can write the $`\mathrm{H}_2^+`$ photodissociation rate as

$$k_8=9.4\times 10^{10}F(M_{},R)\mathrm{s}^1,$$
(46)
and from these rates calculate $`f_{\mathrm{rad},\mathrm{H}^{}}`$ and $`f_{\mathrm{rad},\mathrm{H}_2^+}`$

$`f_{\mathrm{rad},\mathrm{H}^{}}`$ $`=`$ $`40n_\mathrm{H}^1F(M_{},R)`$ (47)
$`f_{\mathrm{rad},\mathrm{H}_2^+}`$ $`=`$ $`1.5n_\mathrm{H}^1F(M_{},R).`$ (48)
By comparing these values and equation 42, we can see that formation via $`\mathrm{H}_2^+`$ contributes at most about 10% of the $`\mathrm{H}_2`$ produced in the gas phase, with the rest coming from $`\mathrm{H}^{}`$.

To evaluate these numbers, I use the fact that for a truncated isothermal sphere,

$$n_\mathrm{H}\left(\frac{R}{1\mathrm{kpc}}\right)^26.9\times 10^3\frac{\mathrm{\Omega }_b}{\mathrm{\Omega }_m}$$
(49)
in regions outside of the core. The final unknown, $`M_{}`$, can be written as

$$M_{}=1.5\times 10^7\epsilon _{}(1+z)^{3/2}\left(\frac{\mathrm{\Omega }_b}{\mathrm{\Omega }_m^{3/2}h}\right)\mathrm{M}_{},$$
(50)
where $`\epsilon _{}`$ is the star formation efficiency of the protogalaxy, $`z`$ is its redshift of formation and $`h`$ is the Hubble constant in units of $`100\mathrm{km}\mathrm{s}^1\mathrm{Mpc}^1`$.

For a protogalaxy that formed in a standard $`\mathrm{\Lambda }`$CDM cosmology ($`\mathrm{\Omega }_m=0.3`$, $`\mathrm{\Omega }_b=0.04`$, $`h=0.7`$) at a redshift $`z=10`$, and that formed stars with an efficiency $`\epsilon _{}=0.01`$, we find that

$$𝒟_{\mathrm{cr}}3.8\times 10^4𝒟_{\mathrm{MW}},$$
(51)
Thus, in this particular example, radiation from the existing stellar population reduces $`𝒟_{\mathrm{cr}}`$ by almost two orders of magnitude.

In view of the uncertainties involved in producing this estimate, it would be unwise to over-generalize. However, since my assumptions verge on the conservative side, it seems likely that in realistic protogalactic models we would see similar effects, and that gas-phase $`\mathrm{H}_2`$ formation will rapidly be overtaken by grain-catalyzed formation.

What about protogalaxies that have yet to form stars? In this case, there is no significant local contribution to the radiation field, which instead is produced by neighbouring sources and/or the cosmological background.

For neighbouring sources, we can reuse the above formalism, as long as we set $`R`$ to the distance to the extragalactic source. However, this is typically an order of magnitude or more greater than the size of a protogalaxy, implying that the effect of the radiation will be at least two orders of magnitude smaller than the effects discussed above. Consequently, radiation from protogalaxies of the size discussed here will have little or no effect on gas-phase $`\mathrm{H}_2`$ formation within their neighbours, unless their emitted flux is substantially larger than has been assumed here.

For the background, we again face the problem that any conclusions that we can draw are strictly limited by our poor knowledge of the star formation history, this time on a cosmological rather than protogalactic scale. The best that we can do is to determine how strong the background needs to be before it has a significant effect. Modeling the background below the Lyman limit as a power-law,

$$J_\nu =J_{21}\left(\frac{\nu _0}{\nu }\right)^\alpha ,$$
(52)
where $`J_{21}=10^{21}\mathrm{erg}\mathrm{s}^1\mathrm{cm}^2\mathrm{Hz}^1\mathrm{sr}^1`$ and where $`h\nu _0=13.6\mathrm{eV}`$, I find that for $`\alpha =1`$,

$$k_6=8.1\times 10^{10}J_{21}\mathrm{s}^1,$$
(53)
and hence

$$f_{\mathrm{rad},\mathrm{H}^{}}=0.6J_{21}n_\mathrm{H}^1.$$
(54)
The significance of the background varies with $`n_\mathrm{H}`$ and hence with position within the protogalaxy. For my example protogalaxy, formed at a redshift $`z=10`$, gas near the truncation radius has a density $`n_\mathrm{H}7\times 10^3\mathrm{cm}^3`$ and thus is affected for $`J_{21}0.01`$; on the other hand, gas in the central core has $`n_\mathrm{H}1\mathrm{cm}^3`$ and is only affected for $`J_{21}1`$.

## 5 Conclusions

The simplicity of the basic chemistry involved in gas-phase $`\mathrm{H}_2`$ formation means that it is easy to construct a fairly accurate expression for the formation rate in terms of only a few parameters: the temperature, density and fractional ionization of the gas, plus the strength of the radiation field near the $`\mathrm{H}^{}`$ and $`\mathrm{H}_2^+`$ photodissociation thresholds. Expressions for the rate of $`\mathrm{H}_2`$ formation via the $`\mathrm{H}^{}`$ and $`\mathrm{H}_2^+`$ ions are given by equations 19 & 23 respectively, and the total formation rate is simply the sum of these two values.

Using these expressions, together with an analytical expression for the grain-catalyzed $`\mathrm{H}_2`$ formation rate, one can solve for $`𝒟_{\mathrm{cr}}`$, the dust-to-gas ratio required for grain-catalyzed $`\mathrm{H}_2`$ formation to overtake gas-phase formation. The results demonstrate that, in principle, gas-phase $`\mathrm{H}_2`$ formation could be comparable to grain-catalyzed formation in galactic gas, particularly at high temperatures where the latter is inefficient. In practice, however, it is usually significantly slower, either because of a shortage of free electrons and protons (which reduces the formation rate of the intermediate ions) or because the incident radiation field destroys the ions before they have a chance to form $`\mathrm{H}_2`$.

Finally, in order to demonstrate the simplicity and potential utility of this approach, I have applied it to a number of astrophysical examples. Not surprisingly, I find that in most cases $`𝒟_{\mathrm{cr}}`$ is significantly less than the mean galactic value, often by several orders of magnitude. Nevertheless, there are counterexamples, such as X-ray photodissociation regions (Lepp & McCray, 1983; Maloney et al., 1996), nova ejecta (Rawlings, 1988), or the high-redshift protogalaxies analyzed in detail here. A common thread linking many of these exceptions seems to be the fact that they have dust-to-gas ratios (but not necessarily metallicities) significantly lower than the mean galactic value.

I would like to acknowledge useful comments on an earlier draft of this paper from Mordecai Mac-Low and Michael D. Smith. I would also like to thank the anonymous referee for a timely and useful report, and for bringing the work of Lepp & McCray to my attention. This work was supported by NSF grant AST99-85392.