text:
# Higher-order binding corrections to the Lamb shift of 2𝑃 states

## I Introduction

The evaluation of the one-loop self-energy of a bound electron is a long standing problem in Quantum Electrodynamics. There are mainly two approaches. The first, developed by P. Mohr , relies on a multidimensional numerical integral involving a partial wave expansion of the electron propagator in the Coulomb field. This approach is particularly useful for heavy hydrogen-like ions. The second approach is based on an expansion of the electron self-energy in powers of $`Z\alpha `$,

$$\delta E_{\mathrm{SE}}=\frac{\alpha }{\pi }(Z\alpha )^4mF,$$
(1)
where

$`F`$ $`=`$ $`A_{40}+A_{41}\mathrm{ln}\left[(Z\alpha )^2\right]+(Z\alpha )A_{50}+`$ (3)
$`(Z\alpha )^2\left(A_{60}+A_{61}\mathrm{ln}\left[(Z\alpha )^2\right]+A_{62}\mathrm{ln}^2\left[(Z\alpha )^2\right]+o(Z\alpha )\right).`$
The leading contribution as given by $`A_{41}`$ has been originally calculated by Bethe in . Many others have contributed to the evaluation of higher orders corrections, for details see an excellent review by Sapirstein and Yennie in . A very general analytical method has been introduced by Erickson and Yennie in . Erickson and Yennie were able to calculate all the coefficients in (3) except for $`A_{60}`$. The calculation of corrections of $`(Z\alpha )^2`$ relative order is a highly nontrivial task because the binding Coulomb field enters in a nonperturbative way, and there is no closed form expression for the Dirac-Coulomb propagator. Additionally, one-loop electron self-energy contributes to all orders in $`Z\alpha `$, and the separation of the $`(Z\alpha )^2`$ relative contribution involves hundreds of terms. A very efficient scheme of the calculation has been introduced in . It was calculated there the $`A_{60}`$ coefficient for the $`1S`$ and $`2S`$ states in hydrogen atom. The method was based on the division of the whole expression into two parts, $`E_L`$ and $`E_H`$, by introducing an artificial parameter $`ϵ`$ which is a cutoff in the photon frequency. In the high energy part $`E_H`$ one expands the electron propagator in powers of the Coulomb field and uses Feynman gauge. In the low energy part one uses Coulomb gauge and applies a multipole expansion. The most important ingredient of this method is the expansion in the parameter $`ϵ`$ after the expansion in $`Z\alpha `$ is performed (for details see the next section).

The calculation presented in this paper is a further development of this original method. In the low energy part we use a Foldy-Wouthuysen transformation. The transformation clearly identifies the leading order contribution and separates out all higher order terms. An additional advantage is that the nonrelativistic Schrödinger-Coulomb propagator can be used here. A closed-form expression of this propagator is known in coordinate and in momentum space (for details see ). This method is applied to the $`2P_{1/2}`$ and $`2P_{3/2}`$ states. All coefficients including $`A_{60}`$ are obtained. We recover all the previously known results, and the new results for $`A_{60}`$ are in agreement with those obtained from the extrapolation of P. Mohr’s data. Our results are relevant for single electron, small $`Z`$ systems (for example atomic hydrogen and $`\mathrm{He}^+`$), which are currently investigated with very high precision. New theoretical values for the Lamb shift of the $`2P_{1/2}`$ and $`2P_{3/2}`$ states and the fine structure summarize our calculations.

## II The $`ϵ`$fw-method

The self-interaction of the electron leads to a shift of the hydrogen energy levels. This shift at the one-loop level is given by

$$\delta E_{\mathrm{SE}}=ie^2\frac{d^4k}{(2\pi )^4}D_{\mu \nu }(k)\overline{\psi }|\gamma ^\mu \frac{1}{\overline{)}p\overline{)}km\gamma ^0V}\gamma _\nu |\psi \overline{\psi }|\delta m|\psi ,$$
(4)
where $`\delta m`$ refers to the mass counter term, and it is understood that the photon propagator $`D_{\mu \nu }`$ has to be regularized to prevent ultraviolet divergences. $`\overline{\psi }`$ is the Dirac adjoint $`\overline{\psi }=\psi ^+\gamma ^0`$.

For the $`\omega `$-integration ($`k_0\omega `$), the lower part of the Feynman integration contour $`C_F`$ is bent into the “right” half plane with $`\mathrm{}(\omega )>0`$ and divided into two parts, the low energy contour $`C_L`$ and the high energy contour $`C_H`$, see Fig. 1. The $`ϵ`$ parameter corresponds to the cut-off $`K`$ which was introduced by H. Bethe in his original evaluation of the low energy part of the electromagnetic shift of energy levels  (specifically, $`K=ϵm`$). The two contours are separated along the line $`\mathrm{}(\omega )=ϵm`$, where $`ϵ`$ is some arbitrary dimensionless parameter, which we assume to be smaller than unity. This method of $`\omega `$-integration has been described in detail in . The two integrations lead to the high and low energy parts $`E_L`$ and $`E_H`$, which are functions of the fine structure constant $`\alpha `$ and of the free parameter $`ϵ`$. Their sum, however,

$$\delta E_{\mathrm{SE}}(\alpha )=E_L(\alpha ,ϵ)+E_H(\alpha ,ϵ),$$
(5)
does not depend on $`ϵ`$. The most important step is the expansion in $`ϵ`$ after the expansion in $`\alpha `$. It eliminates, without actual calculations, many terms that vanish in the limit $`ϵ0`$. To be more specific, in expanding $`E_L`$ and $`E_H`$ in $`ϵ`$ we keep only finite terms (the $`ϵ^0`$-coefficients) and the terms which diverge as $`ϵ0`$. The divergent terms cancel out in the sum, the finite terms contribute to the Lamb shift. This cancelation of the divergent terms is an important cross-check of the calculation. One may use different gauges of the photon propagator for the two parts, because the gauge-dependent term vanishes in the limit $`ϵ0`$. For convenience, we use the Feynman gauge for the high and the Coulomb gauge for the low energy part.

In this work, the treatment of the low energy part is largely simplified by the introduction of a Foldy-Wouthuysen (fw) transformation. It enables one to clearly separate out the leading (nonrelativistic dipole) term, which gives the $`\alpha (Z\alpha )^4`$-contribution, from the relativistic corrections, which give terms in $`\alpha (Z\alpha )^6`$. An additional advantage is the fact that all contributions to the low energy part can be evaluated using the nonrelativistic Schrödinger-Coulomb-Green’s function, whose closed-form solution is well known . Terms which contribute to the Lamb shift up to $`\alpha (Z\alpha )^6`$ can be readily identified, and each of these can be calculated independently. In the low energy part we may expand in the photon momentum $`k`$. The terms which contribute to the Lamb shift in the order of $`\alpha (Z\alpha )^6`$ correspond to the “non-relativistic dipole” term (involving the non-relativistic propagator and wave function), the “non-relativistic quadrupole” term and the “relativistic dipole” term (which involves the relativistic corrections to the wave function and the Dirac-Coulomb propagator). The terms of higher order in $`k`$ vanishes in the limit $`ϵ0`$.

Calculations of the high-energy part are performed almost entirely with the computer algebra system Mathematica . Because of the presence of an infrared cut-off, one can expand the Dirac-Coulomb propagator in powers of the Coulomb potential. A subsequent expansion of the propagator in electron momenta is also performed. This leads finally to the calculation of matrix elements of operators containing $`V`$ and $`𝐩`$ on the P-states. Because $`P`$-wave functions vanish at the origin, all of the relevant matrix elements are finite up to the order of $`(Z\alpha )^6`$.

## III The high-energy part

In this part we use the Feynman gauge ($`D_{\mu \nu }(k)=g_{\mu \nu }/k^2`$). and the Pauli-Villars regularization for the photon propagator

$$\frac{1}{k^2}\frac{1}{k^2}\frac{1}{k^2M^2},$$
(6)
so that the following expression remains to be evaluated:

$$E_H=ie^2_{C_H}\frac{d^4k}{(2\pi )^4}\left[\frac{1}{k^2}\frac{1}{k^2M^2}\right]\overline{\psi }|\gamma ^\mu \frac{1}{\mathit{}\mathit{}m\gamma ^0V}\gamma _\mu |\psi \overline{\psi }|\delta m|\psi $$
(7)
We start by calculating the matrix element

$$\stackrel{~}{P}=\overline{\psi }|\gamma ^\mu \frac{1}{\overline{)}p\overline{)}km\gamma ^0V}\gamma _\mu |\psi $$
(8)
up to the order of $`(Z\alpha )^6`$. The first step in the evaluation of $`\stackrel{~}{P}`$ is the expansion of the matrix

$$M=\gamma ^\mu \frac{1}{\overline{)}p\overline{)}km\gamma ^0V}\gamma _\mu $$

in powers of the binding field. We denote the denominator of the free electron propagator by $`D`$ ($`D=\overline{)}p\overline{)}km`$). Realizing that the binding field $`V=(Z\alpha )^2m/\rho `$ carries two powers of $`(Z\alpha )`$ (with $`\rho =r/a_{\mathrm{Bohr}}`$), we expand the matrix $`M`$ up to $`V^3`$, which leads in turn to four matrices, denoted $`M_i`$,

$`M_0=\gamma ^\mu {\displaystyle \frac{1}{D}}\gamma _\mu ,M_1=\gamma ^\mu {\displaystyle \frac{1}{D}}\gamma ^0V{\displaystyle \frac{1}{D}}\gamma _\mu ,M_2=\gamma ^\mu {\displaystyle \frac{1}{D}}\gamma ^0V{\displaystyle \frac{1}{D}}\gamma ^0V{\displaystyle \frac{1}{D}}\gamma _\mu ,`$ (9)
$`M_3=\gamma ^\mu {\displaystyle \frac{1}{D}}\gamma ^0V{\displaystyle \frac{1}{D}}\gamma ^0V{\displaystyle \frac{1}{D}}\gamma ^0V{\displaystyle \frac{1}{D}}\gamma _\mu .`$ (10)
with $`M=M_0+M_1+M_2+M_3+O((Z\alpha )^7)`$. Defining $`\stackrel{~}{P}_i=\overline{\psi }|M_i|\psi `$, we write the element $`\stackrel{~}{P}`$ as the sum

$$\stackrel{~}{P}=\stackrel{~}{P}_0+\stackrel{~}{P}_1+\stackrel{~}{P}_2+\stackrel{~}{P}_3+O((Z\alpha )^7).$$
(11)
This expansion corresponds to a division of the initial expression into 0-,1-,2- and 3-vertex parts. We then expand each of the matrices $`M_i`$ into the standard 16 $`\mathrm{\Gamma }`$ matrices, which form a basis set of $`4\times 4`$ matrices.

$$M_i=\underset{\beta =0}{\overset{15}{}}c_{i,\beta }\mathrm{\Gamma }^\beta \text{where}c_{i,\beta }=\frac{1}{4}\mathrm{Tr}(\mathrm{\Gamma }_\beta M_i).$$
(12)
The expansion coefficients $`c_{i,\beta }`$ are rational functions of the binding field, the electron and photon energy and momenta. They can therefore be expanded in powers of $`\alpha `$, leaving none of the electron momentum operators in the denominator. Next, we evaluate the matrix elements of these operators with the relativistic (Dirac) wave function $`\psi `$. It is a property of $`P`$ states, which vanish at the origin, that up to order $`(Z\alpha )^6`$, all of the desired matrix elements are finite.

As an example, we describe here the evaluation of the three-vertex matrix element $`\stackrel{~}{P}_3=\overline{\psi }|M_3|\psi `$. It takes on the same values for both $`2P`$ states. Expanding $`M_3`$ into the 16 $`\mathrm{\Gamma }`$-matrices, we find that up to order $`(Z\alpha )^6`$, all expansion coefficients vanish except for the identity $`\mathrm{Id}`$ and $`\gamma ^0`$-matrices. The expansion coefficients are explicitly

$$c_{3,\mathrm{Id}}=16V^3\frac{k^2k^2\omega 4\omega +3\omega ^2+2}{\left(k^2+2\omega \omega ^2\right)^4}b_{3,\mathrm{Id}}V^3,$$
(13)
where $`k=|𝐤|`$ and for simplicity $`m=1`$, and

$$c_{3,\gamma ^0}=2V^3\frac{k^48k^2+12k^2\omega +16\omega 6k^2\omega ^212\omega ^2+4\omega ^3\omega ^48}{\left(k^2+2\omega \omega ^2\right)^4}b_{3,\gamma ^0}V^3.$$
(14)
So up to order $`(Z\alpha )^6`$, the two $`c`$-expansion coefficients are (except for their dependence on $`k`$ and $`\omega `$) functions of the binding field only. Thus, the matrix element $`\stackrel{~}{P}_3`$ is given by

$$\stackrel{~}{P}_3=b_{3,\mathrm{Id}}\overline{\psi }|V^3|\psi +b_{3,\gamma ^0}\overline{\psi }|\gamma ^0V^3|\psi .$$
(15)
The relevant matrix element of the wave function is

$$\overline{\psi }|V^3|\psi =\overline{\psi }|\gamma ^0V^3|\psi =\frac{1}{24}(Z\alpha )^6m^3+O((Z\alpha )^7).$$
(16)
where the first equality holds only in the order of $`(Z\alpha )^6`$. The above matrix elements take on the same values for the $`2P_{1/2}`$ and $`2P_{3/2}`$ states because the radial parts of both $`2P`$ states are the same in the non-relativistic limit.

For the other vertex parts, many more terms appear, and the matrix elements contribute in the lower order also. We give one example here, to be evaluated for the 1-vertex part,

$$\overline{\psi }|\gamma ^0𝐩\left(V𝐩\right)|\psi =\frac{5}{48}(Z\alpha )^4m^3\frac{283}{1152}\left(Z\alpha \right)^6m^3\text{for }2P_{1/2}.$$
(17)
and

$$\overline{\psi }|\gamma ^0𝐩\left(V𝐩\right)|\psi =\frac{5}{48}(Z\alpha )^4m^3\frac{71}{1152}\left(Z\alpha \right)^6m^3\text{for }2P_{3/2}.$$
(18)
For a more detailed review of the calculations see . Having calculated $`\stackrel{~}{P}`$, we subtract the mass-counter-term before integrating with respect to $`k`$ and $`\omega `$. The final $`k`$ and $`\omega `$ integration is performed in the following way. Those terms which appear to be ultraviolet divergent are regularized and integrated covariantly using Feynman parameter approach. The remaining terms are integrated with respect to $`k`$ by residual integration and with respect to $`\omega `$ by changing the integration variable to

$$u=\frac{\sqrt{2m\omega \omega ^2}+i\omega }{\sqrt{2m\omega \omega ^2}i\omega }.$$
(19)
This integration procedure is described in details in . The final results for the high-energy-part are (for the definition of $`F`$ see Eq. (1))

$$F_H(2P_{1/2})=\frac{1}{6}+(Z\alpha )^2\left[\frac{4177}{21600}\frac{103}{180}\mathrm{ln}(2)\frac{103}{180}\mathrm{ln}(ϵ)\frac{2}{9ϵ}\right]$$
(20)
and

$$F_H(2P_{3/2})=\frac{1}{12}+(Z\alpha )^2\left[\frac{6577}{21600}\frac{29}{90}\mathrm{ln}(2)\frac{29}{90}\mathrm{ln}(ϵ)\frac{2}{9ϵ}\right].$$
(21)
## IV The low energy part

In this part we are dealing with low energy virtual photons, therefore we treat the binding field non-pertubatively. Choosing the Coulomb gauge for the photon propagator, one finds that only the spatial elements of this propagator contribute. The $`\omega `$-integration along $`C_L`$ is performed first, which leads to the following expression for $`E_L`$,

$$E_L=e^2_{|𝐤|<ϵ}\frac{d^3k}{(2\pi )^3\mathrm{\hspace{0.17em}2}|𝐤|}\delta ^{T,ij}\psi |\alpha ^ie^{i𝐤𝐫}\frac{1}{H_D(E_\psi \omega )}\alpha ^je^{i𝐤𝐫}|\psi (\omega |𝐤|).$$
(22)
$`H_D`$ denotes the Dirac-Coulomb-Hamiltonian, $`\delta ^T`$ is the transverse delta function, and $`\alpha ^i`$ refers to the Dirac $`\alpha `$-matrices. In the matrix element

$$P^{ij}=\psi |\alpha ^ie^{i𝐤𝐫}\frac{1}{H_D(E_\psi \omega )}\alpha ^je^{i𝐤𝐫}|\psi $$
(23)
we introduce a unitary Foldy-Wouthuysen transformation $`U`$,

$$P^{ij}=U\psi |(U\alpha ^ie^{i𝐤𝐫}U^+)\frac{1}{U(H_D(E_\psi \omega ))U^+}(U\alpha ^je^{i𝐤𝐫}U^+)|U\psi .$$
(24)
The lower components of the Foldy-Wouthuysen transformed Dirac wave function $`\psi `$ vanish up to $`(Z\alpha )^2`$, so that we may approximate $`|U\psi `$ by

$$|U\psi =|\varphi +|\delta \varphi \text{with}\varphi |\delta \varphi =0,$$
(25)
where $`|\varphi `$ is the nonrelativistic (Schrödinger-Pauli) wave function, and $`|\delta \varphi `$ is the relativistic correction.

We define an operator acting on the spinors as even if it does not mix upper and lower components of spinors, and we call the odd operator odd if it mixes upper and lower components. The Foldy-Wouthuysen Hamiltonian consists of even operators only. For the upper left $`2\times 2`$ submatrix of this Hamiltonian, we have the result

$$H_{\mathrm{FW}}=U(H_D(E_\psi \omega ))U^+=m+H_S+\delta H,$$
(26)
where $`H_S`$ refers to the Schrödinger Hamiltonian, and $`\delta H`$ is is the relativistic correction,

$$\delta H=\frac{\left(𝐩\right)^4}{8m^3}+\frac{\pi \alpha }{2m^2}\delta (𝐫)+\frac{\alpha }{4m^2r^3}𝝈𝐋$$
(27)
Now we turn to the calculation of the Foldy-Wouhuysen transform of the operators $`\alpha ^i\mathrm{exp}\left(𝐤𝐫\right)`$. The expression $`U\alpha ^i\mathrm{exp}\left(i𝐤𝐫\right)U^+`$ is to be calculated. Assuming that $`\omega =|𝐤|`$ is of the order $`O((Z\alpha )^2)`$, we may expand the expression $`U\alpha ^ie^{i𝐤𝐫}U^+`$ in powers of $`(Z\alpha )`$. The result of the calculation is

$`U\alpha ^ie^{i𝐤𝐫}U^+`$ $`=`$ $`\alpha ^i\left(1+i\left(𝐤𝐫\right){\displaystyle \frac{1}{2}}\left(𝐤𝐫\right)^2\right){\displaystyle \frac{1}{2m^2}}p^i\left(𝜶𝒑\right)`$ (31)
$`+\gamma ^0{\displaystyle \frac{p^i}{m}}\left(1+i\left(𝐤𝐫\right){\displaystyle \frac{1}{2}}\left(𝐤𝐫\right)^2\right)`$
$`\gamma ^0{\displaystyle \frac{1}{2m^3}}p^i𝐩^2{\displaystyle \frac{1}{2m^2}}{\displaystyle \frac{\alpha }{r^3}}\left(𝐫\times 𝚺\right)^i`$
$`+{\displaystyle \frac{1}{2m}}\gamma ^0\left(𝐤𝐫\right)\left(𝐤\times 𝚺\right)^i{\displaystyle \frac{i}{2m}}\gamma ^0\left(𝐤\times 𝚺\right)^i.`$
In the limit $`ϵ0`$ the odd operators in the above expression do not contribute to the self energy in $`(Z\alpha )^2`$ relative order, so one can neglect the odd operators. It can be shown easily that also the last term in the above expression (proportional to $`𝐤\times 𝚺`$) does not contribute to the Lamb shift in $`(Z\alpha )^2`$ relative order for $`ϵ0`$.

Because we can ignore odd operators, and because the lower components of the Foldy-Wouthuysen transformed wave function vanish, we keep only the upper left $`2\times 2`$ submatrix of Eq. (31), and we write $`U\alpha ^ie^{i𝐤𝐫}U^+`$ as

$`U\alpha ^ie^{i𝐤𝐫}U^+`$ $``$ $`{\displaystyle \frac{p^i}{m}}\left(1+i\left(𝐤𝐫\right){\displaystyle \frac{1}{2}}\left(𝐤𝐫\right)^2\right)`$ (34)
$`{\displaystyle \frac{1}{2m^3}}p^i𝐩^2{\displaystyle \frac{1}{2m^2}}{\displaystyle \frac{\alpha }{r^3}}\left(𝐫\times 𝝈\right)^i`$
$`+{\displaystyle \frac{1}{2m}}\left(𝐤𝐫\right)\left(𝐤\times 𝝈\right)^i,`$
This can be rewritten as

$$U\alpha ^ie^{i𝐤𝐫}U^+=\frac{p^i}{m}e^{i𝐤𝐫}+\delta y^i,$$
(35)
where $`\delta y^i`$ is of order $`(Z\alpha )^3`$. It is understood that the term $`\frac{p^i}{m}e^{i𝐤𝐫}`$ is also expanded up to the order $`(Z\alpha )^3`$. Denoting by $`E`$ the Schrödinger energy ($`E=(Z\alpha )^2m/8`$ for 2P states) and by $`\delta E`$ the first relativistic correction to $`E`$, we can thus write the matrix element $`P^{ij}`$ as

$$P^{ij}=\varphi +\delta \varphi |\left[\frac{p^i}{m}e^{i𝐤𝐫}+\delta y^i\right]\frac{1}{H_S(E\omega )+\delta H\delta E}\left[\frac{p^j}{m}e^{i𝐤𝐫}+\delta y^j\right]|\varphi +\delta \varphi .$$
(36)
In this expression, the leading term and the (first) relativistic corrections can be readily identified. Spurious lower order terms are not present in Eq. (36). By expansion of the denominator $`H_S(E\omega )+\delta H\delta E`$ in powers of $`\alpha `$, the whole expression can be written in a form which involves only the Schrödinger-Coulomb-Green’s function

$$G(E\omega )=\frac{1}{H_S(E\omega )},$$
(37)
whose closed-form expression in coordinate space is given in Eq. (40). We now define the dimensionless quantity

$$P=\frac{m}{2}\delta ^{T,ij}P^{ij}.$$
(38)
Using the symmetry of the $`P`$-wave functions and Eq. (36), we easily see that $`P`$ can be written, up to $`(Z\alpha )^2`$, as the sum of the contributions (39, 47, 48, 49, 50, 51). The leading contribution (the “non-relativistic dipole”) is given by

$$P_{\mathrm{nd}}=\frac{1}{3m}\varphi |p^i\frac{1}{H_S(E\omega )}p^i|\varphi .$$
(39)
The evaluation of this matrix element is described here as an example. For the Schrödinger-Coulomb propagator, we use the following coordinate-space representation ,

$$G(𝐫_1,𝐫_2,E\omega )=\underset{l,m}{}g_l(𝐫_1,𝐫_2,\nu )Y_{l,m}\left(\widehat{𝐫}_1\right)Y_{l,m}^{}\left(\widehat{𝐫}_2\right),$$
(40)
with $`E\omega \alpha ^2m/(2\nu ^2)`$.

$$g_l(r_1,r_2,\nu )=\frac{4m}{a\nu }\left(\frac{2r_1}{a\nu }\right)^l\left(\frac{2r_1}{a\nu }\right)^le^{(r_1+r_2)/(a\nu )}\underset{k=0}{\overset{\mathrm{}}{}}\frac{L_k^{2l+1}\left(\frac{2r_1}{a\nu }\right)L_k^{2l+1}\left(\frac{2r_2}{a\nu }\right)}{(k+1)_{2l+1}(l+1+k\nu )},$$
(41)
where $`a=a_{\mathrm{Bohr}}=1/(\alpha m)`$, and $`(k)_c`$ is the Pochhammer symbol. The evaluation of eq. (39) proceeds in the following steps: The angular integration is performed first. Secondly, the remaining integrals over $`r_1`$ and $`r_2`$ are evaluated using the formula (see e.g. ),

$`{\displaystyle _0^{\mathrm{}}}𝑑te^{st}t^{\gamma 1}L_n^\mu (t)={\displaystyle \frac{\mathrm{\Gamma }(\gamma )\mathrm{\Gamma }(n+\mu +1)}{n!\mathrm{\Gamma }(\mu +1)}}s^\gamma {}_{2}{}^{}F_{1}^{}(n,\gamma ,1+\mu ;{\displaystyle \frac{1}{s}}).`$ (42)
The following formula is useful for carrying out the summation with respect to $`k`$ ,

$$\underset{n=0}{\overset{\mathrm{}}{}}\frac{\mathrm{\Gamma }(n+\lambda )}{n!}s^n{}_{2}{}^{}F_{1}^{}(n,b;c;z)=\mathrm{\Gamma }(\lambda )(1s)^\lambda {}_{2}{}^{}F_{1}^{}(\lambda ,b;c;\frac{sz}{1s}).$$
(43)
The summations lead to hypergeometric functions in the result,

$`P_{\mathrm{nd}}(t)`$ $`=`$ $`{\displaystyle \frac{2t^2\left(36t3t^2+12t^3+29t^4+122t^5413t^6\right)}{9\left(1t\right)^5\left(1+t\right)^3}}+`$ (45)
$`{\displaystyle \frac{256t^7\left(3+11t^2\right)}{9\left(1t\right)^5\left(1+t\right)^5}}{}_{2}{}^{}F_{1}^{}(1,2t;12t;\left({\displaystyle \frac{1t}{1+t}}\right)^2)`$
where

$$t\frac{\sqrt{2mE}}{\sqrt{2m(E\omega )}}=\frac{1}{2}\nu .$$
(46)
In this expression the terms that gives divergent in $`ϵ`$ terms are separated out of the hypergeometric function, so the could be easily integrated out. The other contributions to $`P`$ (for definition fo $`P`$ see eq. (38)) are

* the non-relativistic quadrupole,

$$P_{\mathrm{nq}}=\frac{1}{3m}\varphi |p^ie^{i𝐤𝐫}\frac{1}{H_S(E\omega )}p^ie^{i𝐤𝐫}|\varphi P_{\mathrm{nd}},$$
(47)
* the corrections to the current $`\alpha ^i`$ from the Foldy-Wouthuysen transformation,

$$P_{\delta y}=\delta ^{T,ij}\varphi |\delta y^i\frac{1}{H_S(E\omega )}p^je^{i𝐤𝐫}|\varphi ,$$
(48)
* the contribution due to the relativistic Hamiltonian,

$$P_{\delta H}=\frac{1}{3m}\varphi |p^i\frac{1}{H_S(E\omega )}\delta H\frac{1}{H_S(E\omega )}p^i|\varphi ,$$
(49)
* the contribution due to the relativistic correction to the energy,

$$P_{\delta E}=\frac{1}{3m}\varphi |p^i\frac{1}{H_S(E\omega )}\delta E\frac{1}{H_S(E\omega )}p^i|\varphi ,$$
(50)
* and due to the relativistic correction to the wave function,

$$P_{\delta \varphi }=\frac{2}{3m}\delta \varphi |p^i\frac{1}{H_S(E\omega )}p^i|\varphi .$$
(51)

For almost all of the matrix elements we use the coordinate-space representation of the Schrödinger-Coulomb propagator given in Eq. (40). There are two exceptions: For the non-relativististic quadrupole, we use Schwinger’s momentum space representation and carry out the calculation in momentum space. A rather involved contribution is

$$P_{\delta H}=\frac{1}{3m}\varphi |p^iG(E\omega )\left[\frac{\left(𝐩\right)^4}{8m^3}+\frac{\pi \alpha }{2m^2}\delta (𝐫)+\frac{\alpha }{4m^2r^3}𝝈𝐋\right]G(E\omega )p^i|\varphi .$$
(52)
where $`G(E\omega )=1/(H_S(E\omega ))`$. The form of $`\delta H`$ implies a natural separation of $`P_{\delta H}`$ into three terms,

$$P_{\delta H}=P_{p^4}+P_\delta +P_{LS}.$$
(53)
For $`P_\delta `$,

$$P_\delta =\frac{1}{3m}\varphi |p^iG(E\omega )\left[\frac{\pi \alpha }{2m^2}\delta (𝐫)\right]G(E\omega )p^i|\varphi ,$$
(54)
which involves the zitterbewegungs-term (proportional to the $`\delta `$-function), we use a coordinate-space representation of the Schrödinger-Coulomb propagator involving Whittaker functions (this representation is also to be found in ). The result for $`P_\delta (t)`$ is

$$P_\delta (t)=\frac{\alpha ^2}{27}\frac{t^4\left(3+4t+7t^28tF_2(t)\right)^2}{\left(t^21\right)^4}$$
(55)
where

$$F_2(t)={}_{2}{}^{}F_{1}^{}(1,2t,12t,\frac{t1}{t+1}).$$
(56)
Both terms $`P_{p^4}`$ and $`P_{LS}`$,

$`P_{p^4}`$ $`=`$ $`{\displaystyle \frac{1}{3m}}\varphi |p^iG(E\omega )\left[{\displaystyle \frac{\left(𝐩\right)^4}{8m^3}}\right]G(E\omega )p^i|\varphi `$ (57)
$`P_{LS}`$ $`=`$ $`{\displaystyle \frac{1}{3m}}\varphi |p^iG(E\omega )\left[{\displaystyle \frac{\alpha }{4m^2r^3}}𝝈𝐋\right]G(E\omega )p^i|\varphi ,`$ (58)
involve two propagators $`G(E\omega )`$. We use the Schrödinger equation and the identity

$$[H_S(E\omega ),\frac{1}{r}\frac{}{r}r]=\frac{𝐋^2}{mr^3}\frac{Z\alpha }{r^2}.$$
(59)
to rewrite them to the form that contain only one propagator with modified parameters. Namely, to the desired order in $`(Z\alpha )`$, the expression with two propagators can be replaced by an expression with just one propagator, in which an $`(Z\alpha )^2`$-correction is added to the angular momentum parameter $`l`$ or to the fine structure constant $`\alpha `$ in the radial part of the Schrödinger-Coulomb propagator as given in Eq. (40). For the $`P_{p^4}`$ and $`P_{LS}`$ contributions, many more terms appear in the calculation, and derivatives of the hypergeometric functions with respect to parameters have to be evaluated. The result consists of terms involving elementary functions and hypergeometric functions only, and other terms which involve slightly more complex functions. Some of the summations give rise to the Lerch transcendent $`\mathrm{\Phi }`$. Summations of the form

$$\underset{k=0}{\overset{\mathrm{}}{}}k^n\xi ^k\frac{}{b}{}_{2}{}^{}F_{1}^{}(k,b,c,z).$$
(60)
can be evaluated with the help of Eq. (43), for more details see . Although we do not describe the calculations in detail, we stress that the summation with respect to the $`k`$-index is the decisive point in the calculation. In general, a sensible use of contiguous relations is necessary to simplify the result of any of the summations. Symbolic procedures were written to accomplish this. Through the compartmentalization of the calculation achieved by the Foldy-Wouthuysen transformation, it has been possible to keep the average length of intermediate expressions below 1000 terms.

The contribution to $`E_L`$ due to the $`\delta E_{\mathrm{SE}}`$ is given by

$$E_L=\frac{2\alpha }{\pi m}_0^ϵ𝑑\omega \omega P(\omega ).$$
(61)
Changing the integration variable to $`t`$, we have

$$F=\frac{1}{2}_{t_ϵ}^1𝑑t\frac{1t^2}{t^5}P(t).$$
(62)
The $`P`$-terms are integrated with respect to $`t`$ by the following procedure. Terms which give a divergence for $`ϵ0`$ are extracted from the integrand. The extraction can be achieved by a suitable expansion in the argument of the hypergeometric function(s) which appear in $`P(t)`$. The extracted terms consist of elementary functions of $`t`$ only, so they can be integrated analytically. After integration, the terms are first expanded in $`(Z\alpha )`$ up to $`(Z\alpha )^2`$, then in $`ϵ`$ up to $`ϵ^0`$. The remaining part, which involves hypergeometric functions, is integrated numerically with respect to $`t`$ by the Gaussian method.

The $`t`$-integration leads to $`F`$-terms which we name according to the $`P`$-terms $`F_{\mathrm{nd}}`$, $`F_{\mathrm{nq}}`$, $`F_{\delta y}`$, $`F_{\delta H}`$, $`F_{\delta E}`$ and $`F_{\delta \varphi }`$. The $`F_{\mathrm{nd}}`$-term, which is the same for both $`2P`$-states, is given by

$$F_{\mathrm{nd}}=\frac{4}{3}\mathrm{ln}k_0(2P)+\frac{2}{9}\frac{(Z\alpha )^2}{ϵ}.$$
(63)
We have recovered the first 9 digits of the Bethe logarithm with our (Gaussian) integration procedure (the value for the Bethe logarithm given in  is $`\mathrm{ln}k_0(2P)=0.0300167089(3)`$). The $`F_{\mathrm{nd}}`$-term has, for $`ϵ0`$, a divergence of $`+2/9(Z\alpha )^2/ϵ`$, which cancels the corresponding divergence in the high energy part. All other $`F`$-terms produce logarithmic divergences in $`(Z\alpha )^2\mathrm{ln}(ϵ)`$ (see Table I). The results for the low-energy parts of the $`2P`$-states are

$$F_L(2P_{1/2})=\frac{4}{3}\mathrm{ln}k_0(2P)+\left(Z\alpha \right)^2\left[0.79565(1)+\frac{103}{180}\mathrm{ln}\left((Z\alpha )^2\right)+\frac{103}{180}\mathrm{ln}\left(ϵ\right)+\frac{2}{9ϵ}\right]$$
(64)
and

$$F_L(2P_{3/2})=\frac{4}{3}\mathrm{ln}k_0(2P)+\left(Z\alpha \right)^2\left[0.58452(1)+\frac{29}{90}\mathrm{ln}\left((Z\alpha )^2\right)+\frac{29}{90}\mathrm{ln}\left(ϵ\right)+\frac{2}{9ϵ}\right].$$
(65)
The divergence in $`1/ϵ`$ and in $`\mathrm{ln}(ϵ)`$ cancels out when the low- and high-energy-parts are added. The results for the $`F`$-factors (sum of low-energy-part and high-energy-part) are:

$$F(2P_{1/2})=\frac{1}{12}\frac{4}{3}\mathrm{ln}k_0(2P)+\left(Z\alpha \right)^2\left[0.99891(1)+\frac{103}{180}\mathrm{ln}\left((Z\alpha )^2\right)\right]$$
(66)
for the $`2P_{1/2}`$-state and

$$F(2P_{3/2})=\frac{1}{6}\frac{4}{3}\mathrm{ln}k_0(2P)+\left(Z\alpha \right)^2\left[0.50337(1)+\frac{29}{90}\mathrm{ln}\left((Z\alpha )^2\right)\right]$$
(67)
for the $`2P_{3/2}`$-state. The $`A_{60}`$ coefficients are given by

$$A_{60}(2P_{1/2})=0.99891(1)$$
(68)
and

$$A_{60}(2P_{3/2})=0.50337(1).$$
(69)
The last digit is the cumulated inaccuracy of the numerical integrations. The values for the $`A_{40}`$ and $`A_{61}`$ coefficients are in agreement with known results .

These results can be compared to those obtained by P. Mohr  by extrapolation of his numerical data for higher Z,

$$G_{SE}(2)=0.96(4),G_{SE}(1)=0.98(4)\text{for }2P_{1/2},$$
(70)
and

$$G_{SE}(2)=0.46(2),G_{SE}(1)=0.48(2)\text{for }2P_{3/2},$$
(71)
where the function $`G_{\mathrm{SE}}(Z)`$ for $`2P`$-states is defined by

$$F=A_{40}+(Z\alpha )^2\left[A_{61}\mathrm{ln}\left((Z\alpha )^2\right)+G_{\mathrm{SE}}(Z)\right].$$
(72)
Because $`G_{\mathrm{SE}}(Z=0)=A_{60}`$, these values are clearly in very good agreement with the results of our analytical calculation. Using P. Mohr’s numerical data , we have obtained the following estimates for higher order terms summarized by $`G_{\mathrm{SE},7}`$

$$F=A_{40}+(Z\alpha )^2\left[A_{60}+A_{61}\mathrm{ln}\left((Z\alpha )^2\right)+(Z\alpha )G_{\mathrm{SE},7}(Z)\right],$$
(73)
$$G_{\mathrm{SE},7}(2P_{1/2},Z=1)=3.1(5)\text{and}G_{\mathrm{SE},7}(2P_{3/2},Z=1)=2.3(5).$$
(74)
One of the most important aspects of rather lengthy calculations such as those presented here is to avoid errors. The result has been checked in many ways. Except for checking the values of the terms divergent in $`ϵ`$, it was also checked the value of each $`P`$-contribution as $`\omega 0`$. It can be shown easily that the sum of all contributions to the matrix element $`P`$ in the low-energy part must vanish in the limit $`\omega 0`$. Care must be taken when checking the sum, because after the Foldy-Wouthuysen transformation, hidden terms are introduced which do not contribute to the Lamb shift, but contribute in the limit $`\omega 0`$. The hidden terms originate from the odd operators in Eq. (31). Taking into account these terms, the sum vanishes for both states.

## V Other contributions to the Lamb shift

For the Lamb shift $``$, we use the implicit definition

$$E=m_r\left[f(n,j)1\right]\frac{m_r^2}{2(m+m_N)}\left[f(n,j)1\right]^2++E_{\mathrm{hfs}},$$
(75)
where $`E`$ is the energy level of the two-body-system and $`f(n,j)`$ is the dimensionless Dirac energy, $`m`$ is the electron mass, $`m_r`$ is the reduced mass of the system and $`m_N`$ is the nuclear mass.

For the final evaluation of the Lamb shift the following contributions are added:

1. One-loop self energy. The coefficients are presented in this work. For the determination of the Lamb shift the reduced mass dependence of the terms has to be restored. The relevant formulae are given in . For example, the $`A_{60}`$ have a reduced mass dependence of $`(m_r/m)^3`$. We use Eq. (74) to estimate the theoretical uncertainty from the one–loop contribution.
2. Vacuum polarization correction. It enters for $`P`$-states in higher order (for the formulae see , p. 570).
3. Two-loop contributions due to the anomalous magnetic moment . It is given in analogy to the one-loop contribution as

$$\delta E_{2\mathrm{loop}}=\left(\frac{\alpha }{\pi }\right)^2m\frac{(Z\alpha )^4}{n^3}\left[B_{40}+\mathrm{}\right]$$
(76)
where the $`B`$-coefficients are labeled in analogy to the $`A`$-coefficients for the one-loop self energy. The $`B_{40}`$ coefficient is due to the anomalous magnetic moment of the electron. It is given as

$$B_{40}=\frac{C_{jl}}{2(2l+1)}\left[\frac{197}{72}+\frac{\pi ^2}{6}\pi ^2\mathrm{ln}2+\frac{3}{2}\zeta (3)\right]\left(\frac{m_r}{m}\right)^2,$$
(77)
where $`C_{jl}=2(jl)/(j+1/2)`$.
4. Two loop contributions in higher order. Recently, the logarithmic term

$$B_{62}=\left[\frac{4}{27}\frac{n^21}{n^2}\mathrm{ln}^2\left((Z\alpha )^2\right)\right]\left(\frac{m_r}{m}\right)^3,$$
(78)
has been calculated in . The $`B_{62}`$ term, which is enlarged by the logarithm, probably dominates the contributions to the two-loop self energy in higher order. So the result may also be used to estimate the theoretical uncertainty of the two–loop contribution, coming mainly from the unknown $`B_{61}`$ coefficient. It is taken to be half the contribution from $`B_{62}`$.
5. Three-loop self energy as given by the anomalous magnetic moment .

$$\delta E_{3\mathrm{loop}}=\left(\frac{\alpha }{\pi }\right)^3m\frac{(Z\alpha )^4}{n^3}\left[C_{40}+\mathrm{}\right]$$
(79)
where

$$C_{40}=\left[2\frac{C_{jl}}{2(2l+1)}\mathrm{\hspace{0.17em}1.17611}(1)\right]\left(\frac{m_r}{m}\right)^2.$$
(80)
6. The additional reduced mass dependence of order $`(m_r/m_N)^2(Z\alpha )^4`$ , which we will refer to as the $`(Z\alpha )^4`$ recoil contribution,

$$\delta E_{\mathrm{rec},4}=\frac{(Z\alpha )^4}{2n^3}\frac{m_r^3}{m_N^2}\left(\frac{1}{j+1/2}\frac{1}{l+1/2}\right)\left(1\delta _{l0}\right),$$
(81)
7. The Salpeter correction (relativistic recoil) in order $`(Z\alpha )^5`$ as given in . The formula is for $`P`$-states

$$\delta E_{\mathrm{rec},5}=\frac{m_r^3}{mm_N}\frac{(Z\alpha )^5}{\pi n^3}\left(\frac{8}{3}\mathrm{ln}k_0(n)\frac{7}{3}\frac{1}{l(l+1)(2l+1)}\right).$$
(82)
8. Relativistic recoil corrections in the order of $`(Z\alpha )^6m_r/m_N`$,

$$\delta E_{\mathrm{rec},6}=\frac{m^2}{m_N}(Z\alpha )^6\left[\frac{1}{2}\varphi |\frac{𝐋^2}{r^4}|\varphi \right].$$
(83)
The formula for P-states has been first calculated in . This general form has been obtained by us.

The above contributions are listed in table II for the $`2P`$ states.

## VI Results and Conclusions

The new theoretical values for the Lamb shifts of the $`2P_{1/2}`$ and $`2P_{3/2}`$ states are

$$(2P_{1/2})=12835.99(8)\mathrm{kHz}$$
(84)
and

$$(2P_{3/2})=12517.46(8)\mathrm{kHz}.$$
(85)
From the values of the $`2P`$ Lamb shifts, the fine structure can be determined. It turns out that the limiting factor in the uncertatinty is the experimental value of the fine structure constant $`\alpha `$. Using a value of  (1987)

$$\alpha ^1=137.0359895(61)\text{(44 ppb)},$$
(86)
the fine structure can be determined as

$$E(2P_{3/2})E(2P_{1/2})=10969043(1)\mathrm{kHz}.$$
(87)
With the most recent and most precise value of $`\alpha `$ available  (1995),

$$\alpha ^1=137.03599944(57)\text{(4.2 ppb)},$$
(88)
we obtain a value of

$$E(2P_{3/2})E(2P_{1/2})=10969041.52(9)(8)\mathrm{kHz},$$
(89)
where the first error originates from the uncertainty in $`\alpha `$ and the second from the uncertainty in the Lamb shift difference. Our result for the fine structure disagrees with that used by Hagley and Pipkin in  for the determination of $`L(2S2P_{1/2})`$. Therefore their result of $`L(2S2P_{1/2})=1057839(12)`$ is to be modified and according to our calculation it should be

$$L(2S2P_{1/2})=1057842(12)\mathrm{kHz}.$$
(90)
Precise theoretical predictions for $`P`$-states could be used to compare two different kind of measurements of Lamb shifts in the hydrogen. One is the classic $`2S_{1/2}`$-$`2P_{1/2}`$ Lamb shift measured by several groups , , , and the second is the combined Lamb shift $`(4S2S)\frac{1}{4}(2S1S)`$ as measured by the Hänsch group (for a review see ). The experimental value of 2S Lamb shift can be extracted from E($`2S`$-$`2P_{1/2}`$) having the precise value for $`2P_{1/2}`$ Lamb shift, and can also be determined from the combined Lamb shift through the formula

$$(2S)=\frac{8}{7}\left[\left((4S)\frac{5}{4}(2S)+(1S)\right)_{\mathrm{exp}}\left((4S)\frac{17}{8}(2S)+(1S)\right)_{\mathrm{theo}}\right],$$
(91)
where the subscript exp denotes experimental, and the subscript theo denotes theoretical values. This theo combination has the property that terms scaling $`1/n^3`$ cancel out, which means that almost all QED effects do not contribute, and therefore the quantity can be precisely determined. Such a comparison of completely different experimental techniques is an interesting and valuable test of high precision experiments.

The method of calculation presented in this paper could be directly applied for the evaluation of Lamb shifts and the fine structure in two electron systems, for example in helium or positronium. It was a purpose of this method to use only a Schrödinger-Coulomb propagator, and relativistic effects are incorporated through the Foldy-Wouthuysen transformation. This method clearly separates out the lower and the higher order terms, and expresses the energy shift through the matrix elements of nonrelativistic operators.

## Acknowledgments

This work was done while one of us (K. P.) was a guest scientist at the Max-Planck-Institute for Quantum Optics. The authors would like to thank T. W. Hänsch for hospitaliy, encouragement and stimulation. We are very grateful to P. Mohr for supplying extrapolation data of his 1992 calculations, and to M. Weitz and A. Weis for carefully reading the manuscript. (U. J.) would also like to thank H. Kalf for helpful discussions with respect to the treatment of hypergeometric functions.
Note added (2000): The analytic results for higher-order binding corrections to the Lamb shift of $`2P_{1/2}`$ and $`2P_{3/2}`$–states (in particular, the $`A_{60}`$–coefficient) have recently been confirmed by an improved numerical calculation in the range of low nuclear charge numbers $`Z=1`$$`5`$. For details see the e-print physics/0009090.


format: text

text:
# The mass of a halo

## 1 Introduction

One of the most fundamental predictions of a theory of structure formation is the number density of objects of a given mass, the mass function, at a given redshift. Accurate mass functions are used in a number of areas in cosmology; in studies of galaxy formation, in measures of volumes (e.g. galaxy lensing) and in attempts to infer the normalization of the power spectrum and the density parameter from the abundance of rich clusters. In the latter case the mass function is the point of contact allowing us to bootstrap the excellent statistics of large N-body simulations with observable properties of clusters normalized for example by hydrodynamic simulations of smaller volumes. In this way we can obtain reliable estimates of e.g. the number of clusters as a function of temperature and redshift, which can in turn be used to constrain the matter density, the normalization of the power spectrum and the statistics of the initial density field.

Several different definitions for the “mass” of a halo are in common use, each having different advantages. For example there are at least 3 different definitions of mass for the mass-temperature relation as computed from hydrodynamical simulations of galaxy clusters and all are different from the definition of mass commonly used in the mass function which is itself different from the mass usually employed in analytic studies based on the Press-Schechter (PS (1974); hereafter PS) theory. Observational data have improved to the point where it is important to distinguish between these different definitions of mass, lest we bias our theoretical predictions. We give an example of how such a conversion can be made, at least approximately, for a certain class of mass estimator.

## 2 Halo mass definitions

In this section we describe some of the definitions of the mass of a dark matter halo in common use. This list is obviously not exhaustive, but it is representative. We begin by recalling some background about the spherical top-hat collapse model (see e.g. Peebles PPC (1993); Peacock Peacock (1999); Liddle & Lyth LidLyt (2000), and references therein) from which much of the language in this field has been borrowed.

The spherical top-hat ansatz describes the formation of a collapsed object by solving for the evolution of a sphere of uniform overdensity $`\delta `$ in a smooth background of density $`\overline{\rho }`$. By Birkhoff’s theorem the overdense region evolves as a positively curved Friedman universe whose expansion rate is initially matched to that of the background. The overdensity at first expands but because it is overdense the expansion slows (relative to the background) and eventually halts before the region begins to recollapse. Technically the collapse proceeds to a singularity but it is assumed in a “real” object virialization occurs at twice111There is a small correction to this in the presence of a cosmological constant which contributes a $`\mathrm{\Lambda }r^2`$ potential. the turn-around time, resulting in a sphere of half the turn-around radius. In an Einstein-de Sitter model the overdensity (relative to the critical density) at virialization is $`\mathrm{\Delta }_\mathrm{c}=18\pi ^2178`$. We shall always use $`\mathrm{\Delta }_\mathrm{c}`$ to indicate the overdensity relative to critical of a virialized halo, which will be lower for smaller $`\mathrm{\Omega }_\mathrm{M}`$. Note that some authors use a different convention in which $`\mathrm{\Delta }_\mathrm{c}`$ is specified relative to the background matter density – our $`\mathrm{\Delta }_\mathrm{c}`$ is $`\mathrm{\Omega }_\mathrm{M}`$ times theirs. The linear theory extrapolation of this overdensity is normally denoted $`\delta _\mathrm{c}`$ and is $`(3/20)(12\pi )^{2/3}1.686`$ in an Einstein-de Sitter model. This overdensity is often used as a threshold parameter in PS theory and its extensions and has a very weak cosmology dependence which is often neglected.

We now turn to some of the definitions of mass in use in the literature. Let us suppose that we have identified some object, a problem to which we shall return in §4, and have chosen a fiducial ‘center’ about which to take spherically averaged profiles. We define $`M_\mathrm{\Delta }`$ as the mass contained within a radius $`r_\mathrm{\Delta }`$ inside of which the mean interior density is $`\mathrm{\Delta }`$ times the critical density

$$_0^{r_\mathrm{\Delta }}r^2𝑑r\rho (r)=\frac{\mathrm{\Delta }}{3}\rho _{\mathrm{crit}}r_\mathrm{\Delta }^3.$$
(1)
The ‘virial mass’ from the spherical top-hat collapse model would then be simply $`M_{\mathrm{\Delta }_\mathrm{c}}`$. We shall refer to this mass as $`M_{\mathrm{th}\mathrm{vir}}`$.

An alternative interpretation of the spherical collapse model is that the virial radius corresponds to that point within which the material is virialized and external to which the mass is still collapsing onto the object. Some simulations suggest that this occurs at $`\mathrm{\Delta }=200`$ more or less independent of cosmology, and so a common mass estimator is $`M_{200}`$ which is approximately the virial mass if $`\mathrm{\Omega }_\mathrm{M}=1`$. Authors also use $`M_{500}`$ or $`M_{1000}`$ in some applications.

The final mass estimator we shall consider is the one used by Jenkins et al. (JFWCCEY (2000)) in their analysis of simulations performed by the Virgo consortium. The mass they assign to a halo is simply the sum of the masses of the particles identified as members of the halo by their halo finder. This mass estimator is the least easily interpreted theoretically, though it is well defined algorithmically, but has the very convenient feature that the mass function for this estimator is independent of cosmology!

Unfortunately there is no unique algorithmic definition of a dark matter halo, even within a 3D simulation itself. To define the mass of an object we must first specify the object in question! Although other halo finders are in common use, we shall deal exclusively with halos found using the Friends-of-Friends (Davis et al. DEFW (1985)) algorithm, hereafter called FOF. Our reason for this choice is that it is the group finder employed by the Virgo consortium in deriving the mass functions from their very large N-body simulations (Jenkins et al. JFWCCEY (2000)) and in terms of which they find a universal mass function. The FOF algorithm has one free parameter, $`b`$, the linking length in units of the mean inter-particle spacing. Commonly used values of $`b`$ are 0.1, 0.15 and 0.2, although other choices exist (e.g. Gardner Gar (2000) chooses instead $`b^3=\mathrm{\Omega }_\mathrm{M}\mathrm{\Delta }_\mathrm{c}/3`$). Jenkins et al. (JFWCCEY (2000)) find that the mass function is universal if they take $`b=0.2`$, independent of the cosmology under consideration. Jenkins et al. (JFWCCEY (2000)) claim that in the limit of a large number of particles per halo, FOF finds all particles within an iso-density contour $`b^3`$ times the mean matter density or $`\mathrm{\Omega }_\mathrm{M}/b^3`$ times critical density. We shall return to this point in §4.

In the next sections we describe how these different masses are related.

## 3 An approximate conversion

It would be advantageous if we could convert, at least approximately, between the different definitions of mass in common use. In this section we make a first attempt in this direction. We assume that all halos are spherical and have a universal profile, for example222There is some controversy in the literature about the slope of the density profile as $`r0`$. Luckily departures from the NFW form are expected to occur only at very small radii, containing a vanishingly small amount of the total mass, and thus do not concern us here. the NFW form (Navarro, Frenk & White NFW (1996)):

$$\rho (r)x^1\left(1+x\right)^2,$$
(2)
where $`x=r/r_\mathrm{s}`$ and $`r_\mathrm{s}`$ is a scale radius usually specified in terms of the concentration parameter $`cr_{200}/r_\mathrm{s}`$. Navarro et al. (NFW (1996)) refer to $`r_{200}`$ and $`M_{200}`$ throughout as the ‘virial radius’ and ‘virial mass’ respectively. Again, N-body simulations have shown that the concentration parameter is a weak function of virial mass, taking the value $`c5`$ for masses characteristic of clusters which shall be our focus here.

With a given profile it is straightforward to relate the various mass definitions as shown in Fig. 1. Note that the masses can differ significantly in this model and that the cosmology dependence is quite strong. How well does this crude model predict the relationship between different mass definitions in practice?

## 4 Simulations

The model of the previous section made several approximations. To compare these different mass estimates with each other more carefully we have used some N-body simulations, originally run for another purpose. We simulated the Ostriker & Steinhardt (OstSte (1995)) concordance model, which has $`\mathrm{\Omega }_\mathrm{m}=0.3`$, $`\mathrm{\Omega }_\mathrm{\Lambda }=0.7`$, $`H_0=100h\mathrm{kms}^1\mathrm{Mpc}^1`$ with $`h=0.67`$, $`\mathrm{\Omega }_\mathrm{B}=0.04`$, $`n=1`$ and $`\sigma _8=0.9`$ (corresponding to $`\delta _H=5.02\times 10^5`$). For this cosmology $`\mathrm{\Delta }_\mathrm{c}101`$ from the top-hat collapse model.

We have used two ‘high’ resolution and two lower resolution N-body simulations. The first was a $`256^3`$ particle, dark matter only simulation in a $`100^3h^1`$Mpc box. The simulation was run on 16 processors of the Origin2000 at NCSA with the TreePM-SPH code (White et al. TreePM (2000)) operating in collisionless (dark matter only) mode. The gravitational force softening was of a spline form (e.g. Hernquist & Katz HerKat (1989)), with a “Plummer-equivalent” softening length of $`15h^1\mathrm{kpc}`$ comoving. The simulation was evolved from $`z=100`$ until the present and took a total of 1300 (normalized) CPU hours. A second $`256^3`$ particle simulation was run from $`z=70`$ in a $`256h^1`$Mpc box with a $`35h^1`$kpc softening. This run took a total of 880 (normalized) CPU hours, also on 16 processors. To provide additional statistics on the high-mass end of the mass function we additionally ran several smaller ($`128^3`$ particles) simulations of the same model in boxes of side $`200h^1`$Mpc.

From the $`z=0`$ outputs we generated group catalogues using the FOF algorithm with $`b=0.2`$. For each halo we defined the center for our spherical averages as the particle with the minimum potential energy. This corresponds closely to the most bound particle and the density peak for a halo in all but the most disturbed systems, and is more robust than the center of mass. With the group catalog and centers so defined it is straightforward to calculate each of the estimators described above. Throughout we shall use $`M_{200}`$ as our “base” mass to which the others are compared, since this fits best into the philosophy of §3.

We show in Fig. 2 that the different mass estimators are indeed highly correlated as one would expect (see also Table 1). The ratio of the FOF mass to $`M_{200}`$ is seen to span a reasonable range and be offset from unity by a non-negligible amount. Fig. 3 shows a histogram of the ratio of $`M_{\mathrm{fof0}.2}`$ and $`M_{\mathrm{th}\mathrm{vir}}`$ to $`M_{200}`$ for all of the halos above a given number of particles.

Let us focus first on $`M_{\mathrm{th}\mathrm{vir}}`$. The shape of the histogram is roughly the same regardless of particle number indicating that the distribution is not strongly affected by finite particle numbers in the halos chosen. The top-hat virial mass and $`M_{200}`$ are fairly tightly correlated, as one might expect since they are both spherically averaged statistics and the density contrasts are not too different. The approximation of §3 provides a reasonable estimate of the ratio of masses, presumably because of the spherical averaging being performed.

The situation with regards $`M_{\mathrm{fof0}.2}`$ is more complicated. Firstly the scatter is much larger, since the groups found by FOF can be quite irregularly shaped (see Fig. 4). The actual ratio of $`M_{\mathrm{fof0}.2}`$ to $`M_{200}`$ is larger than predicted by our model. By eye the halo finder looks to have merged close groups. It is possible that with better mass resolution this would occur to a lesser extent, although the lack of a clear trend with particle number in Fig. 3 argues that convergence would be quite slow.

To investigate whether FOF has merged adjacent halos and whether a FOF halo can be considered as all particles above a given density threshold (as assumed in §3) we have examined the particles in the two groups of Fig. 4. For each particle we calculate its distance from the halo center and its density, defined from the distance ($`r_{\mathrm{neig}}`$) to the 32nd nearest neighbour333We use all of the particles in the simulation, not just group members, in computing $`r_{\mathrm{neig}}`$. as

$$\rho \frac{3M}{4\pi }r_{\mathrm{neigh}}^3$$
(3)
where $`M`$ is the enclosed mass. This density estimate is somewhat noisy, but will be sufficiently accurate for our purposes.

The results are shown in Fig. 5 along with the appropriate NFW profile (assuming $`c=5`$) for comparison. Note that the substructure in the halos is readily apparent as density ‘spikes’ as functions of radius and can be matched to ‘sub-halos’ in Fig. 4.

It is clear from Fig. 5 that FOF is (at least approximately) linking all particles above a density $`\mathrm{\Omega }_\mathrm{M}b^3`$ of critical. However the chosen contrast is low enough that this includes particles well outside $`r_{200}`$ (or the virial radius) and the presence of these ‘nearby’ particles is biasing the masses above $`M_{200}`$ and contributing to the extra scatter. This problem could be mitigated to some extent by reducing $`b`$, but it is not clear whether the mass function would remain universal if a different (cosmology independent) $`b`$ were chosen.

Finally we show in Fig. 6 the different mass functions that would be derived from these halos using the 3 different mass definitions. Note that our results, obtained with a completely independent N-body code and analyzed with a completely independent set of software, agree well with the universal mass function of Jenkins et al. (JFWCCEY (2000)). However, especially at the high mass end, the number density is quite sensitive to the mass definition used. Intriguingly the fitting formula quoted by Jenkins et al. (JFWCCEY (2000)) provides a good fit to the mass function if we use $`M_{\mathrm{th}\mathrm{vir}}`$ as our definition of mass, an even better fit than if we use $`M_{\mathrm{fof0}.2}`$.

## 5 An example: the temperature function

An example of where this mass uncertainty may be important is in using the observed temperature function of local X-ray clusters to constrain the normalization of the present day power spectrum. Although Pen (Pen (1998)) has argued for a direct prediction of the temperature function from simulations, most people use a hybrid approach in obtaining this constraint. First a mass function, usually calibrated by simulations such as those of Jenkins et al. (JFWCCEY (2000)), is assumed and then a mass temperature relation from hydrodynamic simulations is used to predict the temperature function.

Hydrodynamic simulations show good agreement for the total mass and X-ray temperature properties of clusters (e.g. Frenk et al. Frenketal (1999)). The mass-temperature relation from these simulations generally follows the virial relation quite well. Ignoring a small correction from the $`\mathrm{\Lambda }r^2`$ potential, for an object virialized at a redshift $`z`$

$$M(T,z)\left[T/(1+z)\right]^{3/2}$$
(4)
and the proportionality constant is then fixed by the simulations. There still remains some disagreement over the precise coefficient and its dependence on cosmology (see e.g. Henry (Henry (2000)) for a recent compilation) but in addition to this different authors choose different conventions for $`M`$. For reasonable $`\mathrm{\Omega }_\mathrm{M}`$, the difference in mass definitions is about the same size as the difference in proportionality constants between groups and the $`1\sigma `$ scatter in $`M`$ at fixed $`T`$ found by the different groups. Unfortunately this alone does not account for all of the discrepancy.

To take only two relations as examples, Evrard, Metzler & Navarro (EMN (1996)) define the mass in Eq. 4 as $`M_{200}`$. The $`M`$ of Bryan & Norman (BryNor (1998)), on the other hand, is the ‘virial’ mass in the sense of the spherical top-hat model. In both cases there is ample justification for the definition chosen, and the authors are explicit about their choice. However, there is an $`\mathrm{\Omega }_\mathrm{M}`$ dependent ratio between these two definitions which cannot be neglected if the temperature function is to be predicted to the factor of 2 level!

## 6 Conclusions

We have shown that different definitions of the “mass” of a halo exist, and have different strengths and weaknesses. It is important to be consistent when combining relations which use different definitions of mass, and we have given an approximate method for converting between some commonly used mass estimators. Mass estimates based on spherical averages are much more tightly correlated with each other than with the mass obtained simply by summing the particles in the group, and can be quite well estimated by assuming a ‘universal’ spherical profile (§3). The hydrodynamical simulations which calibrate observables as a function of cluster mass typically use such a spherically averaged mass definition.

Unfortunately these definitions are not very tightly correlated with the particle based mass used in the universal mass function of halos reported by Jenkins et al. (JFWCCEY (2000)). We have argued that this is because, with a linking length of $`b=0.2`$, FOF is merging neighbouring halos. Such a problem would be mitigated by reducing $`b`$, but it has not been demonstrated that the mass function is universal for $`b0.2`$.

In the cosmological model we have simulated, the fitting form of Jenkins et al. (JFWCCEY (2000)) provides a good match to the mass function of our halos, strengthening the claim of those authors that it is universal. However we find the best fit when the mass estimator used is the top-hat virial mass, rather than the FOF mass. While we have not investigated other cosmological models, we expect that the profiles (from the point of view of estimating masses) will not be very cosmology dependent. It would be interesting to see whether the Jenkins et al. (JFWCCEY (2000)) mass function is ‘universal’ if one uses the (cosmology dependent) top-hat virial mass of the halos.

I would like to thank V. Springel for the use of his FOF group finder, J. Mohr for useful conversations on this issue, and C. Metzler for encouraging me to write it up. M. White was supported by the US National Science Foundation and a Sloan Fellowship. Parts of this work were done on the Origin2000 system at the National Center for Supercomputing Applications, University of Illinois, Urbana-Champaign.


text:
# On the ordering of energy levels in homogeneous magnetic fields

## 1 Introduction

We consider the non-relativistic quantum mechanical theory of a single particle in a homogeneous magnetic field and in an external potential. The study of such systems is currently of interest in the context of theories of atoms in strong magnetic fields. In the way of describing the atoms as an assembly of electrons in an effective potential, an essential question is where the electrons are located. In the case of a very strong magnetic field this problem is connected with the comparison of individual energy levels, as the atom is not to be described with a completely semiclassical theory. See \[LSY94, BSY00, HS00\]. In all these studies the superharmonicity of the potential off the nucleus turned out to be the right property to deduce theorems on the localization of the electrons and on level ordering. It is precisely this property, superharmonicity off a certain axis, which is also used in our extension of the present knowledge presented in this paper.

The potential is assumed to be invariant under rotations around the $`z`$-axis, the magnetic field in $`z`$-direction, so that the eigenvalue $`m`$ of $`L_z`$, the $`z`$-component of the angular momentum, is a “good quantum number”. We may ignore spin and intrinsic magnetic moment, since this degree of freedom decouples in this case from the spatial behavior. Assuming appropriate gauging and appropriate units, either with a positive charge of the particle and the magnetic field pointing into the positive $`z`$-direction, or with a negative charge of the particle and the magnetic field pointing into the negative $`z`$-direction, the Hamiltonian acts in $`^2(^3,d^3𝐱)`$ as

$$H=\mathrm{\Delta }BL_z+\frac{B^2}{4}r_{}^2+V(r_{},z),$$
(1)
where $`B>0`$ denotes the absolute value of the strength of the magnetic field, $`r_{}=\sqrt{x^2+y^2}`$, $`𝐱=(x,y,z)`$. See \[AHS78\] for details on the definition of $`H`$.

The Hamiltonian can be considered as a direct sum of $`H_m`$, where each $`H_m`$ acts in the subspace of eigenfunctions of $`L_z`$ with $`L_z=m`$. We compare the spectra of the operators $`H_m`$, denoting

$$E_m=\mathrm{inf}\mathrm{spec}H_m.$$
(2)
Since each $`H_m`$ is unitarily equivalent to $`H_m+2mB`$, we can restrict the investigation to non-negative $`m`$. (Note that, for $`V=0`$, $`E_m=B`$ for all $`m0`$.)

If $`V(r_{},z)=V_{}(r_{})+V_z(z)`$, the three-dimensional system splits into a one-dimensional $`z`$-dependent and a two-dimensional $`r_{}`$-dependent system, and all the level comparison theorems concerning comparison with a two-dimensional oscillator can be applied, with obvious modifications. See \[BGM85\], the Subsections 4.7 and 6.4 of \[B85\], and \[B91\] for a simpler proof; see also Theorem B in \[GS95\].

But here we are interested in general cases without a splitting into $`r_{}`$-dependent and $`z`$-dependent systems. What is known up to now in these general situations is:

1. There are examples, where $`inf\{E_m\}`$ is attained as a minimum, but not at $`m=0`$, \[LC77, AHS78\].
2. For all $`m0`$ we have $`E_mE_{m+1}+B`$, with strict inequality if the $`E_m`$’s are eigenvalues.
3. If $`E_m`$ and $`E_{m+1}`$ are eigenvalues, $`V(r_{},z)`$ not constant in $`r_{}`$, then $`V/r_{}0`$ implies $`E_{m+1}>E_m`$, and $`V/r_{}0`$ implies $`E_{m+1}<E_m`$, as shown by Grosse and Stubbe in \[GS95\].
4. If $`\mathrm{inf}\mathrm{spec}H`$ is an eigenvalue, $`\mathrm{\Delta }V0`$ (in the sense of distributions) and $`\mathrm{\Delta }V0`$, both at $`r_{}>0`$, the corresponding eigenvector is unique and has angular momentum $`m=0`$, \[BS00\].

The general bound (ii), which holds also for higher eigenvalues, is shown in Lemma 13 in the Appendix, where also a short version of the proof of (iii) is given.

###### 1 Remark (On proving (iv)).

In \[BS00\] we considered $`V=\frac{1}{|𝐱|}+\frac{1}{|𝐱|}\rho `$, with $`\rho `$ a nice axially symmetric positive function. But the proof can be extended, without any change, to the more general case, where negative electric charges are located on the $`z`$-axis, and where $`\rho `$ is any repulsive axially symmetric distribution which is not zero everywhere away from the $`z`$-axis. Actually, this theorem will be a corollary of the central mathematical result of this paper, Theorem 4 (at least if $`\mathrm{inf}\mathrm{spec}H`$ is a discrete eigenvalue).

Posing conditions on the potential $`V`$ in Section 2, we try to allow for a wide variety of physical applications. One of special importance is dealing with the energies of a charged particle in the electric field generated by attractive charges situated on the $`z`$-axis and by an axially symmetric repulsive charge-distribution. It is

###### 2 THEOREM (Monotonicity and concavity for finite charges).

Let

$$V(r_{},z)=_{\mathrm{}}^{\mathrm{}}\frac{1}{|𝐱(0,0,z^{})|}\sigma (z^{})𝑑z^{}+_^3\frac{1}{|𝐱𝐱^{}|}\rho (r_{}^{},z^{})d^3𝐱^{},$$
(3)
with $`\sigma (z)dz`$ a positive finite Borel measure on $``$, and $`\rho (r_{},z)d^3𝐱`$ a non-negative finite Borel measure on $`^3`$. In these cases the sequence $`E_m`$ is non-decreasing, concave, and $`lim_m\mathrm{}E_m=B`$. Moreover, if $`\rho (r_{},z)d^3𝐱<\sigma (z)𝑑z`$, the sequence $`E_m`$ is strictly monotone increasing and strictly concave.

The proof will be given at the end of Section 2.

Another physical application deals with the entire lowest band of energies in an infinite periodic system. It will be stated in Theorem 10.

## 2 Ordering in external potentials

Consider the “annihilation operator”

$$a=\frac{1}{2}\left(\sqrt{\frac{2}{B}}\left(\frac{}{x}i\frac{}{y}\right)+\sqrt{\frac{B}{2}}\left(xiy\right)\right),$$
(4)
and its adjoint, the “creation operator” $`a^{}`$. These operators obey the commutation relation

$$[a,a^{}]=1,$$
(5)
and the operator $`a`$ lowers the angular momentum,

$$[L_z,a]=a.$$
(6)
In absence of external potentials both $`a`$ and $`a^{}`$ commute with the Hamiltonian

$$H_B=\mathrm{\Delta }BL_z+\frac{B^2}{4}r_{}^2.$$
(7)
Given any multiplication operator $`V`$, a straightforward calculation shows the formal equation

$$[a^{},[a,V]]=[[V,a^{}],a]=\frac{1}{2B}\left(\frac{^2V}{x^2}\frac{^2V}{y^2}\right).$$
(8)
It holds in the sense of distributions, if $`V_{\mathrm{loc}}^1`$, and this condition will be fulfilled as a consequence of the following

Technical condition on $`V`$:

For fixed angular momentum $`L_z=m`$, $`V`$ is relatively bounded with respect to $`H_B`$, with some bound less than $`1`$; i.e., for some $`b_1<1`$ and $`b_2`$, depending on $`m`$,

$$V\psi b_1H_B\psi +b_2\psi $$
(9)
for all $`\psi `$ with $`L_z\psi =m\psi `$. Note that this condition implies in particular that $`V_{\mathrm{loc}}^2`$.

The condition (9) guarantees that each $`H_m`$ is semibounded and essentially self-adjoint on $`𝒞_0^{\mathrm{}}`$, and $`H`$ is well defined by $`H=_mH_m`$, without having to assume the semiboundedness of $`H`$.

The following proposition provides a large class of examples of potentials that fulfill the technical condition stated above, and is applicable for most practical purposes.

###### 3 PROPOSITION.

Assume that $`V(𝐱)`$ is locally $`^2`$, and there exists a constant, which we denote as $`V_{2,\mathrm{loc}}`$, such that for all $`n`$ and some $`R>0`$

$$_n^{n+1}𝑑z_{r_{}<R}𝑑x𝑑y|V(𝐱)|^2V_{2,\mathrm{loc}}^2.$$
(10)
Moreover, we assume that there exist constants $`C`$ and $`C_V<B^2/4`$, such that

$$|V(r_{},z)|C+C_Vr_{}^2\mathrm{for}\mathrm{all}r_{}>R.$$
(11)
Then $`V`$ fulfills the condition (9).

###### Proof.

Note that for operators $`A0`$ and $`D^{}=D`$ the inequality

$`(A+D^2)^2`$ $`=`$ $`A^2+D^4+2DAD[[D,A],D]`$ (12)
$``$ $`A^2+D^4[[D,A],D]`$
holds. Using this with $`A=\mathrm{\Delta }`$ and $`D=Br_{}/2`$ we get, for any $`\psi `$,

$$\mathrm{\Delta }\psi ^2+\left(\frac{B}{2}\right)^4r_{}^2\psi ^2\left(\mathrm{\Delta }+\frac{B^2}{4}r_{}^2\right)\psi ^2+B^2\psi ^2.$$
(13)
Denoting $`V_1(𝐱)=V(𝐱)\theta (Rr_{})`$ and $`V_2=VV_1`$, we can estimate

$`V\psi `$ $``$ $`V_1\psi +V_2\psi `$ (14)
$``$ $`V_1\psi +{\displaystyle \frac{4C_V}{B^2}}H_B\psi +\left({\displaystyle \frac{4C_V}{B}}(|m|+1)+C\right)\psi ,`$
where we used $`L_z\psi =m\psi `$. Now $`V_1`$ is relatively $`\mathrm{\Delta }`$-bounded with arbitrary small bound (\[RS78\], Thm. XIII.96), and by (13) the same holds with $`\mathrm{\Delta }`$ replaced by $`H_B`$. ∎

For $`m0`$ let now $`\psi _m`$ denote the ground state wave function of $`H_m`$, if it exists. Writing

$$\psi _m(𝐱)=e^{im\phi }r_{}^mf(r_{},z),$$
(15)
where $`(r_{},\phi )`$ denote polar coordinates for $`(x,y)`$, we see that $`f`$ is a ground state for

$$\stackrel{~}{H}_m=\frac{^2}{r_{}^2}\frac{2m+1}{r_{}}\frac{}{r_{}}\frac{^2}{z^2}+\frac{B^2}{4}r_{}^2mB+V(r_{},z)$$
(16)
on $`^2(^3,r_{}^{2m}d^3𝐱)`$. If $`\stackrel{~}{H}_m`$ has a ground state, it is unique and strictly positive. So $`\psi _m`$ is unique. Moreover, $`f`$ is a bounded Hölder continuous function, so $`\psi _m`$ behaves like $`r_{}^m`$ for small $`r_{}`$; i.e. $`|\psi _m(𝐱)|C|r_{}|^m`$.

We now have the necessary prerequisites to prove

###### 4 THEOREM (Pseudoconcavity in angular momentum).

If $`V(r_{},z)`$ fulfills the technical condition stated above, and if the Laplacian of the potential off the z-axis is non-positive, the ground state energies $`E_m`$ of the restricted Hamiltonians $`H_m`$ obey the following inequalities, if $`m1`$:

$$E_m\mathrm{min}\{E_{m1},E_{m+1}\},$$
(17)
$$E_m\frac{1}{2}(E_{m1}+E_{m+1})\mathrm{if}E_mE_{m1},$$
(18)
$$E_m\frac{1}{2m+1}(mE_{m1}+(m+1)E_{m+1})\mathrm{if}E_mE_{m+1}.$$
(19)
If $`\mathrm{\Delta }V`$ is not vanishing everywhere off the z-axis, and if $`E_m`$ is a discrete eigenvalue, the inequalities (17), (18) and (19) are strict.

###### Proof.

We add $`H_BE_m`$ to $`V`$ in (8), use $`[H_B,a]=0`$, note that $`H=H_B+V`$, take the expectation value with $`\phi 𝒞_0^{\mathrm{}}`$ and expand the double commutator. We get, using moreover $`^2/x^2+^2/y^2=\mathrm{\Delta }^2/z^2`$,

$`\phi |a^{}(HE_m)a+a(HE_m)a^{}|\phi ={\displaystyle \frac{1}{2B}}\phi |\mathrm{\Delta }V|\phi {\displaystyle \frac{1}{2B}}\phi |{\displaystyle \frac{^2V}{z^2}}|\phi `$
$`+\phi |a^{}a(HE_m)+(HE_m)aa^{}|\phi .`$ (20)
Now assume that $`E_m`$ is a discrete eigenvalue and let $`\psi _m`$ be the corresponding ground state. Since $`H`$ is a closed operator and $`𝒞_0^{\mathrm{}}`$ is a core of $`H`$, there exist sequences of wave functions $`\phi _k𝒞_0^{\mathrm{}}`$ with $`\phi _k=1`$ and angular momentum $`m`$, converging in norm to $`\psi _m`$ such that also $`H\phi _kH\psi _m=E_m\psi _m`$ in norm. Since

$$2Ba^{}a=H_B+\frac{^2}{z^2}+2BL_zB,$$
(21)
condition (9) on $`V`$ guarantees that the operators $`a^{}a`$ and $`aa^{}`$ are bounded relative to $`H_m`$ on $`_m`$, the subspace of $``$ where $`L_z=m`$. It follows that, as $`k\mathrm{}`$, also $`a^{}a\phi _k`$ and $`aa^{}\phi _k`$ converge and are bounded in norm, and the second line in (2), with $`\phi =\phi _k`$, converges to $`0`$.

Set $`𝐚=(0,0,a)`$. We claim that it is no restriction to assume that the functions

$$f_k(a)|\phi _k(𝐱+𝐚)|^2V(𝐱)d^3𝐱$$
(22)
have their minimum at $`a=0`$. To see this, let $`\stackrel{~}{\phi }_k(𝐱)=\phi _k(𝐱𝐚_k)`$, where $`𝐚_k=(0,0,a_k)`$ is chosen such that $`f_k(a_k)f_k(a)`$ for all $`a`$. This is possible since by assumption $`E_m`$ is a discrete eigenvalue, so one can not attain any expectation value in the gap above $`E_m`$ by shifting $`\phi _k`$ to infinity. Now

$$\stackrel{~}{\phi }_k|H|\stackrel{~}{\phi }_k\phi _k|H|\phi _kE_m,$$
(23)
so there exists a subsequence, again denoted by $`\stackrel{~}{\phi }_k`$, such that $`\stackrel{~}{\phi }_k\psi _m`$ weakly in $`^2`$. Since also the norms converge, there is even strong convergence. In particular, $`a_k0`$ as $`k\mathrm{}`$. Therefore, since $`H_B`$ is translation invariant in $`z`$-direction, and since $`V`$ is relatively bounded,

$$H(\stackrel{~}{\phi }_k\psi _m)0.$$
(24)
Now denote $`\stackrel{~}{\phi }_k`$ by $`\phi _k`$, which proves our claim. With $`\phi _k`$ chosen as above, we have

$$|\phi _k|^2\frac{^2V}{z^2}=\frac{^2}{a^2}f_k(a)|_{a=0}0.$$
(25)
Concerning $`\phi _k|\mathrm{\Delta }V|\phi _k`$, the Hölder continuity of $`\psi _m`$ (\[LL97\], Theorem 11.7) and the fact that $`\psi _m(r_{}=0,z)=0`$ if $`m1`$ allows us to integrate $`|\psi _m|^2\mathrm{\Delta }V`$ and apply Fatou’s Lemma to conclude that

$$\underset{k\mathrm{}}{lim\; sup}\frac{1}{2B}\phi _k|\mathrm{\Delta }V|\phi _k\frac{1}{2B}\psi _m|\mathrm{\Delta }V|\psi _m,$$
(26)
where we used that $`\phi _k(𝐱)\psi _m(𝐱)`$ pointwise. (The convergence is even uniform in $`𝐱`$, since $`\mathrm{\Delta }\phi _k`$ is uniformly bounded by (13).) If $`\mathrm{\Delta }V<0`$ somewhere away from the $`z`$-axis, the fact that $`|\psi _m|>0`$ away from $`r_{}=0`$ implies that the right hand side of (26) is strictly negative. Using this and (25) in (2) we can conclude that

$$\underset{k\mathrm{}}{lim\; sup}\phi _k|a^{}(HE_m)a+a(HE_m)a^{}|\phi _k<0.$$
(27)
Now, since $`a^{}|\phi _k`$ has angular momentum $`L_z=m+1`$, and $`a|\phi _k`$ has angular momentum $`L_z=m1`$, this implies that the spectrum of either $`H_{m+1}`$ or $`H_{m1}`$ must reach below $`E_m`$, that is

$$E_m>\mathrm{min}\{E_{m1},E_{m+1}\}\mathrm{if}m1.$$
(28)
To prove the relations (18) and (19) in their strict form we use the commutation relation (5) and observe that $`a^{}\phi _k^2=a\phi _k^2+1`$. So the inequality (27) tells us moreover that

$$\underset{k\mathrm{}}{lim\; sup}\left((E_{m1}E_m)a\phi _k^2+(E_{m+1}E_m)(a\phi _k^2+1)\right)<0.$$
(29)
In the case $`E_mE_{m1}`$, we note that the strict form of (18) holds trivially, if $`E_m>E_{m+1}`$. Otherwise, if $`E_mE_{m+1}`$, (29) gives

$$\underset{k\mathrm{}}{lim\; sup}a\phi _k^2(E_{m1}+E_{m+1}2E_m)<E_mE_{m+1}0.$$
(30)
In the case $`E_mE_{m+1}`$ we use $`a\phi _k^2m>0`$, which is a consequence of $`a^{}aL_z`$. Using $`1/m1/a\phi _k^2`$, the inequality (29) gives

$$\underset{k\mathrm{}}{lim\; sup}a\phi _k^2\left((E_{m1}E_m)+(E_{m+1}E_m)(1+\frac{1}{m})\right)<0.$$
(31)
So our assertions on the strict inequalities follow.

If $`H_m`$ does not have a discrete eigenvalue at the bottom of its spectrum, or if $`\mathrm{\Delta }V`$ happens to be zero everywhere off the $`z`$-axis, add $`\epsilon W`$ to the potential, with the superharmonic function $`W(r_{},z)=(z^2r_{}^2)`$, and $`\epsilon `$ small enough (more precisely, $`0<\epsilon <B^2(1b_1)/4`$, with $`b_1`$ as in (9)). The Hamiltonians $`H_m(\epsilon )=H_m+\epsilon W`$ do have discrete ground states for each $`m`$ (\[RS78\], Theorem XIII.67), and (28), as well as (30) and (31) hold for the corresponding energies $`E_m(\epsilon )`$. Strictly speaking, the potentials $`V+\epsilon W`$ are not relatively bounded with respect to $`H_B`$ on $`_m`$, but it is not difficult to see, using boundedness relative to $`H_B+\epsilon z^2`$ instead, that the conclusions above remain valid. Now taking the limit $`\epsilon 0`$ leads to the inequalities (17), (18) and (19). ∎

###### 5 Remark (Change of the conditions).

Theorem 4 and the following corollaries hold also when the condition $`\mathrm{\Delta }V0`$ off the $`z`$-axis is replaced by $`^2V/x^2+^2V/y^20`$ off the $`z`$-axis. The proof is essentially the same as above, but without any involvement of $`^2V/z^2`$, and discreteness of the ground state is not needed.

Considering the entire sequence $`E_m`$, the relation (17) obviously implies that an increase of $`E_m`$ in $`m`$ can turn into a decrease, but not vice versa. We state this and the global forms of the other relations as

###### 6 COROLLARY (Global pseudoconcavity).

Under the same conditions on $`V`$ as above, there is an $`M`$, possibly $`0`$ or $`\mathrm{}`$, such that the sequence $`\{E_m\}`$ is strictly increasing and concave for $`0mM`$ and non-increasing for $`mM`$. Moreover, the ground state energies of two neighboring values of angular momentum, $`\mathrm{},\mathrm{}+1`$, determine upper bounds on all of $`\{E_m\}`$ by a tangential sequence: For all $`\mathrm{}0`$ and $`m0`$,

$$E_mE_{\mathrm{}}+(m\mathrm{})(E_{\mathrm{}+1}E_{\mathrm{}})\mathrm{if}E_{\mathrm{}}E_{\mathrm{}+1},$$
(32)
$$E_mE_{\mathrm{}}+(S_mS_{\mathrm{}})(\mathrm{}+1)(E_{\mathrm{}+1}E_{\mathrm{}})\mathrm{if}E_{\mathrm{}}E_{\mathrm{}+1},$$
(33)
where

$$S_m=\underset{\mu =1}{\overset{m}{}}\frac{1}{\mu }.$$
(34)
These inequalities are strict, if some $`E_\mu `$ for $`m\mu \mathrm{}`$ or $`\mathrm{}+1\mu m`$ is a discrete eigenvalue, if $`m\{\mathrm{},\mathrm{}+1\}`$ and if the Laplacian of the potential is not identically zero off the $`z`$-axis.

###### Proof.

For $`E_m`$ in the increasing part of the sequence $`\{E_m\}`$, the inequality (32) is the standard statement that a tangent line lies above a concave sequence. Since the linear tangential sequence is increasing, it is trivial to extend this inequality to the $`E_m`$ in the decreasing part.

To prove (33), consider $`\mathrm{\Delta }_m=(E_mE_{m1})`$ and observe that (19) can be written as a monotonicity relation for $`m\mathrm{\Delta }_m`$,

$$(m+1)\mathrm{\Delta }_{m+1}m\mathrm{\Delta }_m\mathrm{if}E_mE_{m+1}.$$
(35)
Starting with $`\mathrm{\Delta }_{\mathrm{}+1}`$, which is not positive, this means $`\mathrm{\Delta }_\mu \frac{\mathrm{}+1}{\mu }\mathrm{\Delta }_{\mathrm{}+1}`$ if $`\mu \mathrm{}+1`$, and $`\mathrm{\Delta }_\mu \frac{\mathrm{}+1}{\mu }\mathrm{\Delta }_{\mathrm{}+1}`$ if $`\mu \mathrm{}`$. The extension of the latter inequality to the increasing part of $`\{E_\mu \}`$ is trivial. Summing these inequalities for $`\mathrm{\Delta }_\mu `$, either for $`m+1\mu \mathrm{}`$ or $`\mathrm{}+1\mu m`$, gives (33). Moreover, any strict inequality for one of the $`\mathrm{\Delta }_\mu `$ gives (33) as a strict inequality. ∎

We remark that the asymptotic behavior of the decreasing tangential sequence, $`S_m\mathrm{ln}m`$, is optimal. This is demonstrated on the example with the infinite charged tube in Section 3.

An immediate consequence of this logarithmic divergence is

###### 7 COROLLARY (Either no decrease or infinite decrease).

The
sequence $`\{E_m\}`$ is either not decreasing at all or it is decreasing to $`\mathrm{}`$ as $`m\mathrm{}`$.

Moreover, we have the generalization of the theorem on the ground state stated in \[BS00\].

###### 8 COROLLARY (The ground state has zero angular momentum).

If
$`V`$ fulfills the conditions as in Theorem 4, and if $`\mathrm{inf}\mathrm{spec}H`$ is a discrete eigenvalue, the corresponding eigenvector is unique and has zero angular momentum.

###### Proof.

Corollaries 6 and 7 imply that either $`E_0<E_m`$ or $`E_0=E_m`$ for all $`m`$. Since $`\mathrm{inf}\mathrm{spec}H`$ is a discrete eigenvalue by assumption, it cannot be infinitely degenerate, so our assertion is proved. ∎

###### 9 Remark (The hydrogen atom).

Corollary 8 applies also to the hydrogen atom in a constant magnetic field, and provides a simple proof of $`L_z=0`$ in the ground state; the third proof following \[AHS81\] and \[GS95\].

We now have the tools to prove Theorem 2.

###### Proof of Theorem 2..

The potential $`V`$ is relatively bounded with respect to $`\mathrm{\Delta }`$ with relative bound $`0`$ on the entire Hilbert space. This implies that $`V`$ fulfills the technical condition which is necessary to apply Theorem 4 and Corollary 7, and moreover that the Hamiltonian is bounded below. So the sequence $`E_m`$ is non-decreasing and concave, and $`limE_m=B`$. Since the asymptotic behavior of the potential, as $`r_{}\mathrm{}`$, is $`V(𝐱)(ZC)/|𝐱|`$, where $`Z=\sigma 𝑑z`$, $`C=\rho d^3𝐱`$, it is easy, if $`C<Z`$, to construct trial functions for large $`m`$, such that the expectation value of $`H_m`$ is lower than the edge of its essential spectrum, which is at $`B`$. Such trial wave functions are of the form $`\mathrm{\Phi }_m(x,y)\epsilon \phi (\epsilon z)`$, where $`\mathrm{\Phi }_m`$ is the wave function in $`^2(^2)`$ with $`L_z=m`$ in the lowest Landau band, and $`\epsilon `$ has to be small enough. So each $`E_m`$ is a discrete eigenvalue, and increase and concavity are strict. ∎

## 3 Physical applications and extensions

The mean field model of a positively charged ion
Consider

$$V=\frac{1}{|𝐱|}+\rho \frac{1}{|𝐱|}$$
(36)
with $`\rho =\rho (r_{},z)0`$ and $`C=\rho <1`$. We can apply Theorem 2, and infer that the $`E_m`$’s are monotone increasing eigenvalues.

The mean field model of a negatively charged ion
The potential is as in (36), but now with $`C>1`$. Since $`V(𝐱)(C1)/|𝐱|`$ as $`r_{}\mathrm{}`$, and since the angular momentum barrier shields the nucleus, there are no bound states for large $`m`$. The energies $`E_m`$ may be increasing up to some finite $`M`$, and will all be equal to $`B`$ for $`mM`$.

Attractive point charge in a repulsive homogeneously charged hollow tube
Consider the potential

$$V=\frac{1}{|𝐱|}\tau \theta (r_{}R)\mathrm{ln}(r_{}/R)$$
(37)
where $`\tau `$ and $`R`$ are positive constants. The second part of $`V`$ is the “renormalized” potential of an infinite hollow tube with radius $`R`$. It is unbounded below as $`r_{}\mathrm{}`$, and this implies $`E_m\mathrm{}`$ as $`m\mathrm{}`$. On the other hand, the ground states with angular momentum $`m`$ are localized near $`r_{}\sqrt{2m/B}`$, and will not really “feel” the hollow tube, if they are inside. The $`E_m`$ will increase like the corresponding states of the free hydrogen atom, at least up to $`mR^2B/2`$. So there is an increase of ground state energies followed by a decrease. Moreover, since the $`E_m`$ are varying continuously as functions of the parameters $`\tau `$ and $`R`$, there will be cases, where the maximum of the $`E_m`$ is attained twice, at two neighboring angular momenta simultaneously. The asymptotics of $`\{E_m\}`$ as $`m\mathrm{}`$ is, to leading order,

$$E_m\frac{\tau }{2}\mathrm{ln}m.$$
(38)
An infinite periodic chain of atoms
Let $`𝐚`$ be the vector $`(0,0,a)`$, defining the periodicity in $`z`$-direction. Consider

$$V(r_{},z)=\underset{n,\mathrm{ren}}{}\left(\frac{1}{|𝐱n𝐚|}+\rho \frac{1}{|𝐱n𝐚|}\right),$$
(39)
where the axially symmetric positive charge density $`\rho `$ is localized in the elementary slice $`0z<a`$, and $`\rho `$ is assumed to be finite. If $`\rho 1`$, the infinite sum has to be coupled with a renormalization: Let $`D=1\rho `$.

$$\underset{n,\mathrm{ren}}{}\mathrm{}\underset{N\mathrm{}}{lim}\left(\underset{N}{\overset{N}{}}\mathrm{}+2D\mathrm{ln}N\right).$$
(40)
(Convergence is guaranteed at least for $`\rho `$ with compact support.) This potential is e.g. of interest for the study of chains of atoms in strong magnetic fields within the DM theory of \[LSY94\], as will be discussed in \[JRY00\]. Finiteness of $`\rho `$ guarantees that V can be split into the sum of two parts, one of them uniformly locally $`^2`$, the other diverging only logarithmically as $`r_{}\mathrm{}`$. So Proposition 3 applies and $`V`$ fulfills the technical condition.

Since the Hamiltonian commutes with a group of discrete translations in $`z`$-direction, the Hilbert space can be decomposed into the direct integral of $`_\alpha `$, $`0\alpha <2\pi `$, such that the Hamiltonian $`H`$ is represented as a direct integral of $`H_\alpha `$, where the domain of definition of $`H_\alpha `$ contains the Bloch-functions, continuous functions with the property

$$\psi (x,y,z+a)=e^{i\alpha }\psi (x,y,z),$$
(41)
normalized by $`\psi ^2=_^2𝑑x𝑑y_0^a|\psi |^2𝑑z`$. The ground state energies $`E_m`$ are attained in the space $`_0`$, but we can extend our results to all of the lowest band, to the $`E_m(\alpha )`$, the ground state energies of $`H_m`$ restricted to $`_\alpha `$.

###### 10 THEOREM (Pseudoconcavity of Bloch function energies).

Consider potentials which are symmetric under a group of discrete $`z`$-translations, i.e. $`V(r_{},z+a)=V(r_{},z)`$ for some $`a>0`$, and fulfill the conditions of Prop. 3. Let $`H_\alpha `$ and $`H_{\alpha ,m}`$ be the restrictions of $`H`$ and $`H_m`$ onto the space of Bloch functions, and denote $`E_m(\alpha )=\mathrm{inf}\mathrm{spec}H_{\alpha ,m}`$. Theorem 4 and the Corollaries 6, 7 and 8 hold also for each of the $`E_m(\alpha )`$ instead of the $`E_m`$. The pseudoconcavity is moreover strict if $`\mathrm{\Delta }V`$ is not identically zero off the $`z`$-axis.

###### Proof.

One can mimic the proof of Theorem 4. All the operations we made there commute with the $`z`$-translation. Since each $`H_{\alpha ,m}`$ has discrete spectrum, as we show in Proposition 12 in the Appendix, the inequalities (17) – (19) hold in strict form if $`\mathrm{\Delta }V0`$. ∎

## Appendix

###### 11 LEMMA.

Let $`\mathrm{\Delta }`$ be the Laplacian with Neumann boundary conditions on

$$𝒢_n=\{𝐱,nz<n+1\}.$$
(42)
Then, for $`p>3/2`$ and for some constant depending on $`p`$,

$$\mathrm{\Delta }+V(𝐱)1C(p)V_{}_p^{2p/(2p3)},$$
(43)
where $`V_{}=min\{V,0\}`$, and the norm is for the restriction of $`V_{}`$ onto $`𝒢_n`$.

###### Proof.

By the Sobolev inequality on $`𝒢_n`$ there is a constant $`C>0`$ such that

$$f_2^2+f_2^2Cf_6^2$$
(44)
for all $`f^1(𝒢_n)`$ (cf. \[A75\], Thm. 5.4). Therefore we can estimate, for $`\psi _2=1`$,

$$\psi |\mathrm{\Delta }+V|\psi 1+C(R)\psi _6^2V_{}_p\psi _6^{3/p},$$
(45)
where we used Hölders inequality, with $`p^1+q^1=1`$, twice:

$`|\psi |V_{}|\psi |`$ $``$ $`V_{}_p\psi ^2_q,`$ (46)
$`\psi ^2_q`$ $``$ $`\psi _6^{3/p}\psi _2^{23/p}.`$ (47)
Optimizing (45) with respect to $`\psi _6`$ yields the desired result. ∎

###### 12 PROPOSITION (Discrete spectrum).

Each operator $`H_{\alpha ,m}`$ defined in Theorem 10 has discrete spectrum.

###### Proof.

We split the potential $`V`$ into

$`V_{\mathrm{loc}}`$ $`=`$ $`\theta (Rr_{})V`$
$`\mathrm{and}V_{\mathrm{}}`$ $`=`$ $`\theta (r_{}R)V.`$
Then we split the Hamiltonian $`H`$ into

$`H_{\mathrm{loc}}`$ $`=`$ $`\delta H_B+V_{\mathrm{loc}}`$
$`\mathrm{and}H_{\mathrm{}}`$ $`=`$ $`(1\delta )H_B+V_{\mathrm{}},`$
where $`\delta =\frac{1}{2}(14C_V/B^2)`$, $`C_V`$ as defined in Prop. 3. We do the same splittings for all the $`H_{\alpha ,m}`$. We use that $`\mathrm{inf}\mathrm{spec}H_{\alpha ,m,\mathrm{loc}}\mathrm{inf}\mathrm{spec}H_{N,m,\mathrm{loc}}`$, where $`N`$ denotes the Neumann boundary conditions, use Lemma 11 and the condition (10) to find a lower bound to $`H_{\alpha ,m,\mathrm{loc}}`$ by a constant operator,

$$H_{\alpha ,m,\mathrm{loc}}\mathrm{const}.(1+V_{2,\mathrm{loc}}^4)C(V).$$
(48)
Now we add $`H_{\alpha ,m,\mathrm{loc}}`$ and $`H_{\alpha ,m,\mathrm{}}`$ and use the condition (11):

$$H_{\alpha ,m}C(V)(1\delta )(\mathrm{\Delta }+Bm)C+\frac{1}{2}(\frac{B^2}{4}C_V)r_{}^2.$$
(49)
The operator on the right side obviously has discrete spectrum, and the min-max principle implies that the same is true for $`H_{\alpha ,m}`$. ∎

We add a short version of the proof of the Grosse-Stubbe inequality. We need the following

###### 13 LEMMA (General bound on the decrease).

Let $`E_{m,n}`$ be the $`n`$’th eigenvalue (counted from below) of $`H_m`$, or the edge of its essential spectrum, if there are less than $`n`$ discrete eigenvalues. Then, if $`m0`$,

$$E_{m,n}E_{m+1,n}+B.$$
(50)
###### Proof.

We map each subspace $`_m`$ unitarily onto $`^2(_+\times ,2\pi r_{}dr_{}dz)`$, by writing, for $`\psi _m`$,

$$\psi (𝐱)=e^{im\phi }\chi (r_{},z),$$
(51)
and mapping $`\psi \chi `$.

The Hamiltonian $`H_m`$ is then unitarily equivalent to

$$\widehat{H}_m=\frac{^2}{r_{}^2}\frac{1}{r_{}}\frac{}{r_{}}\frac{^2}{z^2}+\frac{B^2}{4}r_{}^2+\frac{m^2}{r_{}^2}mB+V(r_{},z),$$
(52)
and the right side of (50) is the $`n`$’th eigenvalue of

$$\widehat{H}_{m+1}+B=\widehat{H}_m+\frac{2m+1}{r_{}^2}.$$
(53)
This means that $`\widehat{H}_{m+1}+B\widehat{H}_m`$, so (50) follows with the help of the min-max principle. ∎

###### Proof of the Grosse-Stubbe inequality.

The following formula can be interpreted as calculating the acceleration of the particle, subtracting the effect of the Lorentz force, and then taking a linear combination of the $`x`$\- and $`y`$-components, denoted as $`x_+=x+iy`$:

$$[H,[H,x_+]]+2B[H,x_+]=2\frac{x_+}{r_{}}\frac{V}{r_{}}.$$
(54)
Now we take the matrix elements between $`\psi _{m+1}`$ and $`\psi _m`$,

$$\left((E_{m+1}E_m)^2+2B(E_{m+1}E_m)\right)\psi _{m+1}|x_+|\psi _m=2\psi _{m+1}|\frac{x_+}{r_{}}\frac{V}{r_{}}|\psi _m,$$
(55)
and choose the phases appropriately, so that we can write $`\psi _m=(x_+/r_{})^m\chi _m`$ with $`\chi _m0`$. We get

$$(E_{m+1}E_m)(E_{m+1}E_m+2B)\chi _{m+1}|r_{}|\chi _m=2\chi _{m+1}|\frac{V}{r_{}}|\chi _m.$$
(56)
The matrix element on the left side is positive, and by Lemma 13 we know that $`E_{m+1}E_m+2BB>0`$. So the sign of $`\frac{V}{r_{}}`$, which determines the sign of the right side, determines also the sign of $`(E_{m+1}E_m)`$ on the left side. ∎

The Grosse-Stubbe inequality can be extended in a weakened form to cases where the $`E_m`$ are not necessarily eigenvalues. To see this, we approximate $`V(r_{},z)`$ by $`V_R(r_{},z)V(\mathrm{min}\{r_{},R\},z)+z^2/R`$. With this potential the ground state energies $`E_{R,m}`$ have to be eigenvalues, and are strictly ordered. Then we consider the limit $`R\mathrm{}`$. Since $`E_{R,m}E_m`$, it follows that the monotonicities of $`V`$ in $`r_{}`$ imply monotonicities of the $`E_m`$ in $`m`$, which don’t have to be strict if they are not eigenvalues.


text:
# CMB Signatures of a Primordial Magnetic Field

## I Introduction

The presence of a magnetic field in the early universe influences the evolution of metric perturbations, and as a result, might leave observable traces in the cosmic microwave background (CMB). During the past few years the gravitational effects of a homogeneous magnetic field were presented by several authors; for a detailed review see grasso  and references therein. More realistic is the case of a stochastic magnetic field, because any causal generation mechanism results in stochastic fields. Some estimates of CMB temperature power spectra from density perturbations induced by a stochastic magnetic field are given in lemoine  and numerical simulations are presented in koh . The vorticity induced by a stochastic magnetic field was studied in SB  and CMB temperature fluctuations were obtained in the case of a single $`k`$ mode. The gravitational waves generated by tangled magnetic fields and the resulting tensor CMB temperature fluctuations spectra are given in DFK .

Here we outline a comprehensive analytic calculation of temperature and polarization power spectra due to all types of metric perturbations arising from a primordial stochastic magnetic field, which we assume to be statistically homogeneous and isotropic and described by a power law. Detailed results will be presented elsewhere.

## II Metric Perturbations Induced by Magnetic Fields

The energy density of the magnetic field is treated as a first order perturbation, which thus does not alter the background cosmology. Prior to decoupling, the conductivity of the primordial plasma can be assumed to be infinite, which in the comoving frame implies the “frozen-in” condition $`𝐄=𝐯\times 𝐁`$, where $`𝐯`$ is the plasma velocity and $`𝐄`$ is the electric field induced by plasma motions. Infinite conductivity leads to a vanishing electric field in linear perturbations theory ($`v1`$) and permits a split of the spatial structure and time dependence of the magnetic field. As the universe expands, magnetic field lines are simply conformally diluted due to flux conservation: $`𝐁(\eta ,𝐱)=𝐁(\eta _0,𝐱)/a^2`$, where $`a`$ is the scale factor.

We assume that magnetic field is distributed homogeneously and isotropically with a power law power spectrum DFK :

$$B_i(𝐤)B_j^{}(𝐤^{})=(2\pi )^5\frac{\lambda ^{n+3}}{\mathrm{\Gamma }(\frac{n+3}{2})}B_\lambda ^2(\delta _{ij}\widehat{k}_i\widehat{k}_j)k^n\delta (𝐤k^{}),k<k_D.$$
(1)
The spectrum vanishes for all scales smaller than a damping scale $`k>k_D`$; we require the spectral index $`n>3`$ so as not to overproduce superhorizon coherent fields. $`B_\lambda `$ is the rms magnetic field strength today smoothed over a co-moving length scale $`\lambda `$.

The electromagnetic stress-energy tensor may be geometrically decomposed into scalar (density), vector (vorticity) and tensor (gravitational wave) perturbation modes: $`\tau _{ij}^{(B)}=\mathrm{\Pi }_{ij}^{(S)}+\mathrm{\Pi }_{ij}^{(V)}+\mathrm{\Pi }_{ij}^{(T)}`$. From the tensor $`\mathrm{\Pi }_{ij}^{(S)}`$ we can construct scalar $`\mathrm{\Pi }^{(S)}`$ which describes the density perturbations, and from the tensor $`\mathrm{\Pi }_{ij}^{(V)}`$ a vector $`\mathrm{\Pi }_i^{(V)}`$ which describes the vorticity perturbations. Then $`\mathrm{\Pi }^{(S)}`$ is proportional to the magnetic field energy density $`|\rho _B(𝐤,\eta _0)|\frac{3}{2}|\mathrm{\Pi }^{(S)}(𝐤,\eta _0)|`$, and $`\mathrm{\Pi }_i^{(V)}`$ is proportional to the vortical piece of the Lorentz force, $`L_i^{(V)}=k\mathrm{\Pi }_i^{(V)}`$. The isotropic spectra $`\mathrm{\Pi }^{(S),(V),(T)}`$ appear as source terms for scalar, vector and tensor perturbations (for details see DFK , KKM1 ).

To obtain the evolution equations for perturbations, we must derive the scalar, vector and tensor magnetic field correlation functions for $`\mathrm{\Pi }^{(S)}`$, $`\mathrm{\Pi }_i^{(V)}`$ and $`\mathrm{\Pi }_{ij}^{(T)}`$ as described in DFK  and obtain the corresponding spectra as a functions of $`n`$, $`B_\lambda `$ and $`k_D`$. Terms induced by the magnetic field (e.g. Lorentz force, magnetic field anisotropic stresses) appear only in the dynamical equations describing baryons. We use notational conventions of Ref. cmbform1  and employ the Newtonian gauge; $`a`$ is the scale factor, $`\eta `$ is conformal time, and 0 subscripts refer to the present time.

Scalar perturbations: Since we assume that the background spacetime is unperturbed and linear order perturbations arise only from the magnetic field, the stress energy of the magnetic field is not compensated by anisotropic stress of the plasma and we have non-zero initial conditions for gravitational potentials. Such initial conditions might arise from an inflation-like process; the compensated case, arising from causal processes like phase transitions, has more complicated dynamical effects, but our estimates should still be of the correct order of magnitude. We limit the spectral index of the magnetic field to $`n<2`$ DFK ; the case $`n>2`$ corresponds to causally generated fields, which also must be compensated.

The scalar constraint equations from the Einstein equations are two Poisson equations,

$$k^2\mathrm{\Phi }=4\pi Ga^2[\rho _f\delta _f+\rho _B],k^2(\mathrm{\Psi }+\mathrm{\Phi })=12\pi Ga^2\mathrm{\Pi }^{(S)},$$
(2)
where $`\rho _f\delta _f`$ is the total fluid perturbation. The dynamical Einstein equations determine the evolution of the photon and baryon densities (see also Adams ):

$$\dot{\delta }_\gamma =\frac{4}{3}(kv_\gamma ^{(S)}+3\dot{\mathrm{\Phi }}),\dot{\delta }_{\mathrm{cdm}}=(kv_{\mathrm{cdm}}^{(S)}+3\dot{\mathrm{\Phi }}),\dot{\delta }_b=(kv_b^{(S)}+3\dot{\mathrm{\Phi }})3\frac{\dot{a}}{a}\frac{\mathrm{\Pi }^{(S)}}{\rho _b}.$$
(3)
Using the zero-order approximation for tight coupling regime, $`v_bv_\gamma v^{(S)}`$, we obtain the following equation of motion KKM1  from momentum conservation equations for baryons and photons:

$$\dot{v}^{(S)}\frac{k}{4}\delta _\gamma k\mathrm{\Psi }+\frac{k\mathrm{\Pi }^{(S)}}{\rho _\gamma +p_\gamma }=0.$$
(4)
Using the equations of motion, it is possible to show that the two terms on the right side of the first Poisson equation are roughly equal. An approximate solution for the gravitational potential is then

$$\mathrm{\Phi }(\eta )\frac{12\pi G\mathrm{\Pi }^{(S)}(\eta _0,k)}{k^2a^2}.$$
(5)
Vector perturbations: The vector Einstein equation which describes the evolution of the vector potential sourced by a stochastic magnetic field is

$`\dot{V}_i+2{\displaystyle \frac{\dot{a}}{a}}V_i=16\pi Ga^2\mathrm{\Pi }_i^{(V)}(\eta ,𝐤)/k.`$ (6)
Another constraint equation relates the vector potential to the vorticity $`𝛀`$:

$$k^2V_i=16\pi Ga^2(\rho +p)\mathrm{\Omega }_i.$$
(7)
Introducing the Lorentz force term into the baryon Euler equation, neglecting baryon inertia, and solving the baryon momentum conservation equation in the tight-coupling approximation, we obtain an approximate solution for the vorticity,

$$𝛀(\eta ,𝐤)\frac{k𝚷^{(V)}(\eta _0,𝐤)\eta }{\rho _{\gamma 0}+p_{\gamma 0}}.$$
(8)
Tensor perturbations: The tensor Einstein equation that describes the evolution of gravitational waves sourced by a stochastic magnetic field is lemoine ; DFK

$`\ddot{H}_{ij}+2{\displaystyle \frac{\dot{a}}{a}}\dot{H}_{ij}+k^2H_{ij}=8\pi Ga^2\mathrm{\Pi }_{ij}^{(T)}(\eta ,𝐤).`$ (9)
The tensor case has no additional constraint equations. As in the case of vector perturbations, we neglect the tensor anisotropic stress of the plasma, which is in general negligible. A Green function technique gives an approximate solution for $`\eta >\eta _{\mathrm{eq}}`$ DFK

$`\dot{H}(\eta ,k)4\pi G\eta _0^2z_{\mathrm{eq}}\mathrm{ln}\left({\displaystyle \frac{z_{\mathrm{in}}}{z_{\mathrm{eq}}}}\right)k\mathrm{\Pi }^{(T)}(\eta _0,k){\displaystyle \frac{j_2(k\eta )}{k\eta }},`$ (10)
where $`H(\eta ,k)`$ and $`\mathrm{\Pi }^{(T)}(\eta _0,k)`$ are isotropic correlation functions of the tensors $`H_{ij}`$ and $`\mathrm{\Pi }_{ij}^{(T)}`$ respectively, $`z_{\mathrm{in}}`$ is the initial redshift at which the magnetic field is generated, and $`z_{\mathrm{eq}}`$ is the redshift at matter-radiation equality.

## III CMB Anisotropies

The solutions for the metric perturbations presented above can be used to compute the power spectra of temperature and polarization anisotropies in the CMB. The contributions from scalar, vector, and tensor perturbations are all uncorrelated. For all three cases, the temperature fluctuations arise simply from the Sachs-Wolfe effect, while the polarization fluctuations are sourced by these temperature fluctuations. For the scalar case, the temperature perturbation is just the familiar $`\mathrm{\Phi }/3`$ at decoupling. The vector and tensor cases are slightly less familiar; expressions are given in Ref. Durrer  for vector perturbations and Ref. DFK  for tensor perturbations.

We have derived analytic approximations to the angular power spectra of temperature and polarization for all three types of metric perturbation, neglecting small-scale damping effects due to photon diffusion or magnetic damping. As examples, for scalar perturbations with $`n<3/2`$, the temperature power spectrum is

$$\mathrm{}^2C_{\mathrm{}}^{TT}\frac{(2\pi )^5}{3}\frac{G^2\eta _0^4}{a_{\mathrm{dec}}^4}\frac{2^{2n+1}n}{\left(n+3\right)\left(2n+3\right)}\frac{\mathrm{\Gamma }\left(2n\right)}{\mathrm{\Gamma }^2\left(n+\frac{1}{2}\right)\mathrm{\Gamma }^2\left(\frac{n+3}{2}\right)}\left(\frac{\lambda }{\eta _0}\right)^{2n+6}\mathrm{}^{2n+2}B_\lambda ^4,$$
(11)
while the power spectrum for the polarization-temperature cross correlation is

$$\mathrm{}^2C_{\mathrm{}}^{TE}\frac{8\pi ^4}{9}\frac{G\eta _0^2}{\rho _{\gamma 0}+p_{\gamma 0}}\frac{2^{2n1}n}{\left(n+3\right)\left(2n+3\right)}\frac{\mathrm{\Gamma }\left(2n2\right)}{\mathrm{\Gamma }^2\left(n\frac{1}{2}\right)\mathrm{\Gamma }^2\left(\frac{n+3}{2}\right)}\left(\frac{\lambda \mathrm{}}{\eta _0}\right)^{2n+6}\left(\frac{\eta _0}{\eta _{\mathrm{dec}}}\right)^2B_\lambda ^4.$$
(12)
Complete expressions for all cases will be presented elsewhere; only the temperature power spectrum for tensors DFK  and a partial calculation of the temperature power spectrum for vectors SB  have been worked out previously.

The total power spectra for $`n=2`$ are shown in the Figure. For the smoothing scale we take $`\lambda =0.1`$ Mpc, which corresponds to the maximum comoving length scale of galaxies (as in DFK ); for the magnetic field damping scale we take the Alfven cut-off scale SB ; DFK , which is roughly $`k_D4.5h`$ Mpc<sup>-1</sup>. Photon diffusion and magnetic damping will give an exponential cutoff to this fluctuation power spectra above a characteristic value of $`\mathrm{}`$. Note that for $`n>2`$, the polarization power spectra are actually larger than the temperature power spectrum for $`l>50`$.

The magnetic field spectral index is confined to the range $`3<n<2`$. As shown in DFK ; KKM1 , no spectral index will produce a scale-invariant CMB power spectrum ($`\mathrm{}^2C_{\mathrm{}}`$ constant) for scalar temperature and polarization anisotropies while $`n3`$ results in a scale-invariant cross correlation spectrum $`\mathrm{}^2C_{\mathrm{}}^{TE}`$. The vector perturbations cannot produce scale-invariant spectra for allowed values of $`n`$; for tensor perturbations both, the temperature and polarization spectra $`\mathrm{}^2C_{\mathrm{}}`$ are scale-invariant for $`n3`$. Scalar perturbations are always subdominant if compared to vector and the dominant tensor perturbations.

The CMB power spectra due to a magnetic field vary in amplitude like the square of the energy density perturbations, or as $`B^4`$. It will thus be impractical to obtain upper limits on magnetic field strengths which are significantly more stringent than around $`10^9`$ Gauss; limits improve as $`n`$ decreases. On the other hand, if primordial fields are present at this level, the combined signature in the various microwave background power spectra will give an unmistakable signal. We also note that the vector and tensor perturbations generated by magnetic fields are one of the few cosmological sources of “B” polarization kkt , along with primordial tensor perturbations and gravitational lensing of the microwave background. Stochastic magnetic fields also will induce Faraday rotation in the microwave background polarization, providing an additional “B” polarization signal; we are currently considering this question.

Acknowledgments: This work has been supported by the COBASE program of the U.S. National Research Council and by the NASA Theory Program through grant NAG5-7015. T.K. acknowledges a grant from the Swiss NSF for participation in CAPP2000. A.K. is a Cotrell Scholar of the Research Corporation.