text:
# Quadratic tranverse anisotropy term due to dislocations in Mn12-Ac crystals directly observed by EPR spectroscopy

## Abstract

High-Sensitivity Electron Paramagnetic Resonance experiments have been carried out in fresh and stressed Mn12-Acetate single crystals for frequencies ranging from 40 GHz up to 110 GHz. The high number of crystal dislocations formed in the stressing process introduces a $`E(S_x^2S_y^2)`$ transverse anisotropy term in the spin hamiltonian. From the behaviour of the resonant absorptions on the applied transverse magnetic field we have obtained an average value for $`E`$ = 22 mK, corresponding to a concentration of dislocations per unit cell of $`c`$ = 10-3.

Since J. Friedman et al. found the stepwise magnetic hysteresis on Mn12-Acetate molecular cluster, and interpreted it in terms of resonant quantum tunneling , a huge number of experimental measurements have been carried out on this compound (see also references herein), showing its quantum behaviour under many different experimental techniques. Mn12-Ac is a spin 10 molecular cluster which forms macroscopic single crystals in which all the molecules are identical and they take the same spacial orientation respect to the crystallographic axes. The magnetic structure of the spin $`S`$ = 10 of each molecule is represented, in first approximation, by the Hamiltonian $`H=DS_z^2g\mu _B\text{S}\text{H}`$, where $`D`$ is the uniaxial anisotropy constant, g is the Lande factor, $`\mu _B`$ is the Bohr magneton and H is the applied magnetic field. In the absence of external field, opposite sign levels of the 2$`S`$+1 $`S_z`$ components of the spin $`S`$ are degenerated in energy two to two in the both sides of the anisotropy barrier $`U`$ = D$`S^2`$. When a longitudinal magnetic field, $`H_z`$, is applied, the spin levels in the direction of the field become more stable and the system relaxes towards the new equilibrium state. It has been firmly established experimentally that this relaxation can occur via quantum tunneling. However, in the absence of transverse components of external magnetic fields, the uniaxial Hamiltonian cannot explain the quantum relaxation. Hamiltonian terms not commuting with $`S_z`$ are needed in order to explain the observed quantum phenomena. Due to the tetragonal symmetry of the Mn12-Ac crystals, transverse anisotropy terms which are quadratic on the spin operator are not allowed. The lowest order transverse terms are those proportional to $`S_+^4+S_{}^4`$ , otherwise, these terms imply that quantum tunneling can not occur at all the resonant fields, $`H_z=nD`$, contrary to the experimental observations in Mn12-Ac. To explain the presence of quantum tunneling at any resonant field there have been several theoretical approaches, including the transverse component of the magnetic field due to dipolar intermolecule interactions or hyperfine interaction with the nuclei spin of the manganese atoms . However, from the magnetic relaxation experiments , it seems to be necessary stronger magnetic fields than dipolar or hyperfine fields in order to explain the high value of the quantum relaxation that Mn12-Ac exhibits. Recently, Chudnovsky and Garanin  have suggested a new theoretical approach which explains the quantum behaviour of Mn12-Ac in terms of dislocations existing in the crystals. They propose that crystal dislocations introduce quadratic terms on $`S_x`$ and $`S_y`$, producing tunneling in a lower order of perturbation theory than the transverse-field. Unpublished experimental results seem to support this new theory . In this paper we will show a new experimental approach which indicates that dislocations formed in a strongly stressed single crystal introduce the $`E(S_x^2S_y^2)`$ term suggested by Chudnovsky and Garanin, where $`E`$ is the transverse anisotropy constant.

The Mn12-Ac organometallic cluster forms a molecular crystal of tetragonal symmetry with the lattice parameters $`a`$ = 1.732 nm and $`c`$ = 1.239 nm . The unit cell contains two Mn12O12 molecules surrounded by four water molecules and two acetic acid molecules. In the crystallization process point defects usually appear in a low number along the whole crystal. It has been shown experimentally that dislocations can be created in a Mn12-Ac single crystal by rapid thermal cycles, in which the high change of temperature in a low scale of time produces radial and tangential tensions between the core and the surface of the crystal . It is easy to visualize that tensions and pressure forces produced by mechanical distortion of the crystal may generate the same kind of dislocations in a single crystal. This is exactly what we have done with a Mn12-Ac single crystal. When a dislocation is created inside the crystal it has produced a disorder of the molecules in the vicinity of the dislocation. As the number of dislocations in the crystal grows the disorder extends over the whole crystal converting it in a mosaic crystal. The mosaicity is related to the number of dislocations existing in the crystal and can be determined experimentally by analyzing the width of the X-ray difraction peaks .

In our experiments we have used two single crystals: a) fresh single crystal and b) strongly stressed single crystal. To distort the crystal we put it in a glue by one of the extremities. Then a gradually increasing force was applied to the other extreme, perpendicular to the longer length, until the the crystal fractured. Then, the longer resulting part of the crystal was adequately cleaned and tested by X-ray diffraction to be sure that it was a single crystal. Both samples were characterized by X-ray analysis before making the EPR experiments. To obtain a fine check of the variation of the mosaicity between fresh and stressed crystals we have used low angle difraction peaks in order to minimize the widening associated to the lack of monochromaticity of the $`K_{\alpha 1}`$ and $`K_{\alpha 2}`$ X-ray emission lines of Mo element. We used a four-circle single-crystal X-ray diffractometer (Enraf-Nonius CAD4, MoKa radiation) in the characterisation. Reflections ($`\pm `$ 2,$`\pm `$ 2,$`\pm `$ 2) were studied for both fresh and stressed single crystals. Peak intensities of the stressed crystal were normalised according to the measured volume of the studied fragment. Figure 1 shows $`\mathrm{\Delta }\omega \mathrm{\Delta }(2\theta )`$ plot of (2,2,-2) reflection for the fresh (left) and stressed crystal (right). An enlargement of the peak width along the $`\omega `$ direction is clearly observed although there is not a significant change along the 2$`\theta `$ direction (see inset of figure 1). Assuming that the distance between dislocations is inversely related to the widening of the reflection peak in the $`\omega `$ direction  we can conclude that the stressed crystal has a larger mosaicity than the fresh one. As the $`\omega `$ peak-width in the stressed crystal is a factor about 2 higher than in the fresh crystal, we can conclude that the number of dislocations increases about one order of magnitude in the distortion process.

The high frequency resonance experiments have been carried out using the AB millimeter wave vector network analyzer (MVNA) . The base frequency obtained from this source (range 8 \- 18 GHz) is multiplied by Q, V and W Schottkey’s diodes to obtain frequency range used in our experiment (37 - 109 GHz). The sample, a single Mn12-Ac crystal, is placed on the bottom of the cylindrical resonant cavity, halfway between its axis and perimeter. The applied dc magnetic field is parallel to the cavity axis and approximately perpendicular to the easy (c) axis of the crystal. The experiment frequencies are TE0np ($`n`$, $`p=`$ 1, 2, 3, …) which are the resonant frequencies for the cavity used. Resonance Q-factor varies from 20000 at TE011 mode (41.6 GHz) to a few thousand at higher frequencies. Due to the high sensitivity at resonance, there is an increase over conventional EPR of almost 3 order of magnitude, allowing for the detaction of the absorption peaks suggested theoretically by the diagonalization of the corresponding spin Hamiltonian.

The spin Hamiltonian used to explain the experimental data obtained in the last years is

$$=DS_z^2BS_z^4+C(S_+^4+S_{}^4)g\mu _B𝐇𝐒$$
(1)
where $`D`$, $`B`$, and $`C`$ parameters have been experimentally obtained by EPR, Neutron spectroscopy and magnetic relaxation experiments . H is the applied magnetic field. The component of the applied magnetic field parallel to the easy axis direction, called longitudinal component, $`H_{}`$ = $`H`$ $`\mathrm{sin}\theta `$, changes the barrier height between the two classical spin orientations, while the component of the field on the hard plane, transverse component, $`H_{}`$ = $`H`$ $`\mathrm{cos}\theta `$, affects the overlapping of the respective wave functions, which determines the quantum splitting of the degenerate spin states. In the absence of longitudinal field, the quantum splitting of the different $`m`$ spin levels, $`\mathrm{\Delta }_m`$, and consequently the rate of resonant tunneling between the spin levels depend on the magnitude of the transverse component $`H_{}`$. The evolution of the spin levels in the two wells, as a function of the intensity of the transverse component of the field, can be deduced from the diagonalization of the spin Hamiltonian of equation (1). In the transmission specta of our EPR experiments we can detect absorptions peaks corresponding to the absorption of radiation of frequency $`f`$ by the transitions effectuated between the spin levels of the Hamiltonian with energy difference equal to $`hf`$. From the field position of these peaks we can extract the behaviour of the quantum splittings on the transverse magnetic field. The results obtained with the fresh crystal are plotted in figure 2 (solid circles). The experiments have been done at frequencies ranging from 50 GHz up to 110 GHz and fields up to 9 T, placing the crystal with the easy $`z`$-axis perpendicular to the field direction. We have found the best fitting of our data for $`D`$ = 555 mK, $`B`$ = 1.3 mK, and $`C`$ = 2.2$`\times 10^2`$ mK (solid lines in figure 2), in good agreement with the values given in references . From this figure, it is clearly observed that the dependence of the quantum splittings ($`\mathrm{\Delta }_{10}`$ first right line, $`\mathrm{\Delta }_9`$ second, $`\mathrm{\Delta }_8`$ third, and so on) on the transverse field is matching perfectly with the theoretical calculation. The lines that appear at high frequencies and round going down and up with the field correspond to the energy difference between different quantum splittings. However, there are four points appearing at low frequencies between ground splitting and first excited splitting lines which do not have a conventional explanation. These peaks appear also in the experiments done by Hill et al.  and the authors can not give any explanation. The presence of these peaks constitute an incognita for us too.

In figure 3 we show the EPR absorption spectra recorded at $`f`$ = 67 GHz for both fresh (A) and stressed (B) Mn12-Ac single crystals. The labelling used in the figure $`\alpha _{m,\phi }`$ refers to the $`\mathrm{\Delta }_m`$ splitting absorptions with the field applied perpendicular to the easy axis with an angle, $`\phi `$, with respect to the $`x`$ magnetic axis. For a fresh crystal, represented by the Hamiltonian (1), there is a symmetry between any direction perpendicular to the easy axis. The peaks observed for fresh crystal correspond to the resonances with the splittings $`\mathrm{\Delta }_9`$, $`\mathrm{\Delta }_8`$, $`\mathrm{\Delta }_7`$, and so on, as has been shown in figure 2. However, for the stressed crystal, each $`\alpha _m`$ absorption peak appears doubled. This phenomena can be explained by the addition of a $`E(S_x^2S_y^2)`$ term to the Hamiltonian (1), as Chusdnovsky and Garanin suggest as an effect of dislocations . This term introduces the hardest anisotropy along the $`x`$-axis, while $`y`$-axis remains as a medium anisotropy axis. For this reason, different directions of the applied field give different values of each quantum splitting, $`\mathrm{\Delta }_{m,\phi }`$. This is the same behaviour observed in the powder sample of Fe8 molecular clusters . The angular dependence of $`\mathrm{\Delta }`$ ($`\phi `$) at a fixed value of the transverse component of the field is not monotonic. Because of the shape of the function $`\mathrm{\Delta }`$ ($`\phi `$), for a sample with hard axis oriented at random, that is with not preference for any angle $`\phi `$, there are two values of $`\mathrm{\Delta }`$ for which the density of states has a peak. These are the values of the splitting corresponding to $`\phi =`$ 0 and $`\phi =\frac{\pi }{2}`$ (see references  and ). In the absence of dissipation, the contribution of each Mn12 molecule to the imaginary part of the susceptibility is proportional to $`\delta (\omega \frac{\mathrm{\Delta }[\phi ,H_{}]}{\mathrm{}})`$. However, the total imaginary part of the susceptibility is,

$$\chi ^{\prime \prime }_0^\pi g(\phi )\delta (\omega \frac{\mathrm{\Delta }[\phi ,H_{}]}{\mathrm{}})𝑑\phi ,$$
(2)
where $`g(\phi )`$ is the distribution of molecules on $`\phi `$. For a fresh crystal, having no significant number of dislocations, there is an equivalence between $`x`$ and $`y`$-axes as the Hamiltonian (1) has no preference for any transverse direction. Due to this, the splitting does not depend on $`\varphi `$ and the amplitude $`A`$ of the absorption of electromagnetic radiation must have only one peak corresponding to $`\mathrm{\Delta }[H_{}]`$ = h$`f`$, as it can be seen in figures 2 and 3A. On the contrary, in a strongly stressed single crystal of Mn12 the dislocations are randomly affecting the magnetic structure of the molecular clusters, introducing the term $`E(S_x^2S_y^2)`$ in different manner for each molecule . That is, the fact that the effect of the dislocations existing in a single crystal depends on its direction along the crystal and on the distance to a given molecule may be considered as a random generation of a $`x`$-hard-axis for each molecule. For this reason, a Mn12-Ac single crystal with a large mosaicity can be aproximated as a powder sample with $`x`$-axis of the molecules oriented at random. In this case, Eq. (2) can be rewritten as

$$\chi ^{\prime \prime }_0^{\mathrm{}}\delta (\omega \frac{\mathrm{\Delta }[\phi ,H_{}]}{\mathrm{}})\left(\frac{d\mathrm{\Delta }}{d\phi }\right)^1𝑑\mathrm{\Delta }=\left(\frac{d\mathrm{\Delta }}{d\phi }\right)^1|_{\mathrm{\Delta }=\mathrm{}\omega }.$$
(3)
Therefore, there are two field values, solutions of the equations $`\mathrm{\Delta }[0,H_1]`$ = h$`f`$ and $`\mathrm{\Delta }[\frac{\pi }{2},H_2]`$ = h$`f`$, at which the amplitude of the absorption is maximal. Due to this, the two doubled peaks of each splitting observed in the stressed single crystal (figure 3) can be attributed to the two mean orientations of the splitting $`\mathrm{\Delta }`$ on the angle $`\phi `$, $`\mathrm{\Delta }_{m,0}`$ and $`\mathrm{\Delta }_{m,\pi /2}`$. As the distance between the two doubled peaks basically depends on the parameter $`E`$ we can extract a first value of $`E`$ 20 mK.

A more precise analysis of the effect of dislocations on the spin Hamiltonian of Mn12-Ac can be achieved by studying the behaviour of the quantum splittings on the transverse magnetic field. Thus, EPR experiments have been done in the stressed crystal for many frequencies ranging from 40 Ghz up to 110 GHz. In figure 4, it has been plotted the field position of the EPR absorption peaks found at the experiment frequencies $`f`$ (open circles). The EPR data have been fitted by using the magnetic level structure resulting from the diagonalization of the Hamiltonian (1), adding the term $`E(S_x^2S_y^2)`$ attributed to the effect of the dislocations. The results of our fitting procedure are shown in figure 4. Black lines correspond to $`\phi `$ = 0 and blue lines correspond to $`\phi `$ = $`\pi `$/2. The values of the Hamiltonian parameters used in our fitting procedure are $`D`$ = 675 mK, $`B`$ = 0.9 mK, $`C`$ = 1.8$`\times `$10-2 mK, and $`E`$ = 22 mK. Comparing these values with the values extracted from the fitting of the fresh crystal absorption peaks of figure 2 one can conclude that dislocations introduce a general variation of the spin Hamiltonian, apart for the introduction of the transverse anisotropy term. This effect of dislocations on the magnetic structure of the molecules on the vicinity is expected in the theoretical model of Chudnovsky and Garanin . As the authors observe, dislocations may introduce other effects, for example, transverse magnetic fields. This effect has not been considered in our analysis because its complexity. However, we can extract quantitative information of the number of dislocations existing in the stressed single crystal by analyzing the distribution of the generated transverse anisotropies. From the theoretical model , one can extract the distribution of the logarithm of the transverse anisotropy, $`\mathrm{ln}(E/2D)`$, as a function of the concentration of dislocations per unit cell, $`c`$. This distribution function has a maximum at a different value of $`E/2D`$ depending on $`c`$. If we assume that the observed absorption peaks correspond to this maximum in the distribution, we can extract an approximated value of $`c`$ from the resulting Hamiltonian parameters from the fit of the experimentally obtained radiation spectra. We have obtained $`\mathrm{ln}(E/2D)`$ = -4.1. This corresponds to a concentration of dislocations per unit cell of $`c`$ 10-3, in good agreement with the theoretical estimation , and with the experimental results  for a similar samples. Using the comparison X-ray analysis of the mosaicity for both fresh and stressed crystals of figure 1, one can observe that the number of dislocations is increased by almost one order of magnitude than the fresh sample. We conclude that in a fresh crystal the number of dislocations per unit cell is approximately $`c`$ 10-4.

Through high-sensitivity EPR experiments carried out in a strongly stressed single crystal of Mn12-Acetate we have directly obtained the magnitude of the transverse anisotropy term $`E(S_x^2S_y^2)`$, with $`E`$ = 22 mK, associated to dislocations existing in the crystal. We have estimated the value of the concentration of dislocations per unit cell in both fresh and stressed crystals. It may be also possible that the the combined effect of a heavy X-ray irradiation dose and subsequent thermal stressing treatment create a large number of new defective sites and extends the original ones. These sites will necessarily have a lower symmetry and could lead to new EPR absorption peaks, with an $`E`$-term, such as determined here. A clear effect of lattice defects on magnetization tunneling has been detected recently . We also note that an EPR line-broadening effect of naturally present defects in Mn12 and Fe8 single crystals has been reported recently . Additional investigations are thus needed in order to clearly understand the origin of the newly found EPR peaks in the present work.

This work was supported by MEC Grant number PB-96-0169 and EC Grant number IST-1929. Work at NHMFL/FSU was supported through a contractual agreement between the NSF through Grant No. NSF-DMR-95-27035 and the State of Florida.

FIGURE CAPTIONS

Figure 1: $`\mathrm{\Delta }\omega `$-$`\mathrm{\Delta }(2\theta )`$ plot of both fresh (left) and stressed (right) Mn12-Ac single crystals for (2,2,-2) reflection. Increase of imperfections is evidenced by broadening along $`\omega `$ direction.

Figure 2: Resonant peaks from fresh single crystal of Mn12-Ac for different frequencies, ranging from 50 GHz up to 110 GHz, as a function of the magnetic field applied perpendicularly to the magnetic easy axis direction of the crystal. The solid lines are the fitting result of the diagonalization of the spin Hamiltonian of eq. (1).

Figure 3: EPR spectra recorded at 67 GHz for both fresh (A) and stressed (B) single crystals of Mn12-Ac. The experiments were done at $`T`$ = 10 K. The peaks are associated to the different quantum splittings through the next nomenclature: $`\mathrm{\Delta }_{m,\phi }`$.

Figure 4: EPR peaks for the stressed single crystal of Mn12-Ac. The data are fitted by adding the term $`E(S_x^2S_y^2)`$ to the Hamiltonian of eq. (1). The fitting lines represent the field behaviour of the quantum splittings for $`\phi `$ = 0 (black lines) and $`\phi `$ = $`\pi `$/2 (blue lines).


format: text

text:
# LABORATORY MODELLING OF THE EARTH RADIATION BELT

## Abstract

Method of the laboratory modelling of the Earth radiation belt is presented. Method can be used for the estimation of consequences of global energetic and communication projects realizations.

The radiation belt of the Earth is the inner part of the magnetosphere, in which the geomagnetic field hold charged particles with kinetic energy from 10 KeV to 100 MeV. This belt play the essential role for the sun radiation regime for the Earth and for the electromagnetic waves propagation.

As it is well known, even small alteration of solar activity influent essential upon the biological and ecological balance at the Earth for all levels - from viruses and micro-organisms to biological societies and ecological systems at all. The most bright indicator of the solar activity alteration are changes of the pathogenic organisms activity as a result of a displacement of biological equilibrium. But the other phenomena, which are not connected with pathogenic biological objects are very important for the life on the Earth and ecological balance too, although they can not be so obviously observed.

In this century the realisation of some global projects will be very actual for the supply of the Earth with the electrical energy to be converted from the solar energy and for the creation of global communication systems. The realisation of these projects cannot be neutral for the ionosphere and for the radiation belt state.

For example, for the Glazer Project it will be necessary to create a system of powerful solar energy converters, which should send powerful SHF beams through the radiation belt. But even small alterations in the radiation belt will change the solar irradiation regime of the Earth the same way, as for alteration of solar activity by constant magnetosphere.

Thus it is very dangerous to realise global energetic and communication projects without the all-round investigation of consequences of these realisations.

That is why it is necessary in particular to have a possibility of the radiation belt laboratory modelling and its experimental investigation as a complete system under the influence of different factors.

Besides that this task is actual also immediate for the solving of problems of radiation propagation and communication too.

This problem was solved in  and  for a general case, i.e. for all kinds of planets, which have magnetical fields. The laboratory Modelling of the Earth radiation belt is in fact a particular case of this solution. Therefore no additional patent descriptions for any details will be required.

In the Earth radiation belt the charged particles realised a motion, which can be considered as a sum of three components: 1 - quick cycle motion due to Lorenz force in the geomagnetic field; 2 - the ”longitude-drift” of these cycles along the Equator in the crossed magnetic and gravitational fields; and 3 - the oscillation of particles in directions, normal to the cycles planes and to the direction of drift.

Change of movement directions in these oscillations occurs due to geometrical phenomenon of so-called specular reflection.

The main problem for the realisation of radiation belt laboratory modelling is that it is on principle impossible to create the 3-dimensional radial gravitation field in a laboratory device.

This problem can be solved through the replacing of gravitation field with the electrostatical field, because the mathematical description of these two fields are similar. For the realisation of the drift-phenomenon it is principally necessary the existence of two crossed fields. It is not important for this process, what kind of field will be used together with the magnetic field, gravitational field or the electrostatical one.

Finally it is proposed to use for this modelling the Penning Trap  or our physically similar traps with other construction of electrodes -, which can make this modelling more effective.

By doing this way, the cycle motion in Earth magnetic field corresponds to the ”cyclotron motion” in the Penning Trap, the drift along the Equator corresponds to the ”magnetron motion”, and the interpolar oscillation corresponds to the ”axial motion”.

For the realisation of laboratory experiments with this model it is necessary to provide the possibility of the ceaseless change of field configuration during the experiment. It is achieved through the special construction of the trap ring electrode. This electrode in the form of rotation-hyperboloid of one sheet is collected from straight plates, which lies in hyperbolic surface. This technical solution can be used not only for the aims of the problem under consideration, but in constructions of some antenna systems with smooth rearrangeable geometry too.


text:
# References

UCRHEP-T326

TIFR/TH/01-48

November 2001

Anomalous Neutrino Interaction, Muon g-2,
and Atomic Parity Nonconservation
Ernest Maa and D. P. Royb
a Physics Department, University of California, Riverside, California 92521, USA
b Tata Institute of Fundamental Research, Mumbai (Bombay) 400005, India

## Abstract

We propose a simple unified description of two recent precision measurements which suggest new physics beyond the Standard Model of particle interactions, i.e. the deviation of $`\mathrm{sin}^2\theta _W`$ in deep inelastic neutrino-nucleon scattering and that of the anomalous magnetic moment of the muon. Our proposal is also consistent with a third precision measurement, i.e. that of parity nonconservation in atomic Cesium, which agrees with the Standard Model.

The minimal Standard Model (SM) of particle interactions is consistent with all present experimental data with only a few possible exceptions. One such is a recent measurement  of the electroweak parameter $`\mathrm{sin}^2\theta _W`$ from $`\nu _\mu `$ and $`\overline{\nu }_\mu `$ interactions with nucleons, which claims a three-standard-deviation departure from the SM prediction. Another is the measurement  of the anomalous magnetic moment of the muon, which originally claimed a value higher than the SM prediction by 2.6 standard deviations , but is now revised down to only 1.6$`\sigma `$ after a theoretical sign error has been corrected . A third important constraint comes from the measurement  of parity nonconservation in atomic Cesium, which was thought to be in disagreement with the SM, but subsequent improved theoretical calculations  have shown it to be in good agreement. In addition, the phenomonena of neutrino oscillations are now well-established  which suggest strongly that neutrinos have mass and mix with one another.

In this paper we propose a simple unified description of all the above effects by extending the SM to include the gauge symmetry $`L_\mu L_\tau `$ . The relevance of this symmetry to the muon $`g2`$ value and neutrino mass has been discussed by us in a previous paper . Here we focus on how it can also explain the NuTeV result  and its other possible experimental consequences.

Our model assumes the anomaly-free gauge symmetry $`U(1)_X`$ with gauge boson $`X`$ which couples to $`(\nu _\mu ,\mu )_L`$, $`\mu _R`$ with charge $`+1`$ and to $`(\nu _\tau ,\tau )_L`$, $`\tau _R`$ with charge $`1`$, but not to any other fermion. This means that it has the contribution

$$\mathrm{\Delta }a_\mu =\frac{g_X^2m_\mu ^2}{12\pi ^2M_X^2}$$
(1)
to the muon anomalous magnetic moment. It also contributes to $`\nu _\mu `$ and $`\overline{\nu }_\mu `$ interactions, but since $`X`$ does not couple to quarks, the NuTeV result  is only affected if $`X`$ mixes with the $`Z`$ boson of the SM. This also applies to atomic parity nonconservation.

In our previous paper , we assume for simplicity that $`XZ`$ mixing is zero by the imposition of an interchange symmetry in the Higgs sector, but we also mention that this symmetry cannot be maintained for the entire theory, so that a small deviation is to be expected. This small deviation (corresponding to a mixing angle of order $`10^3`$) turns out to be just what is needed to explain the NuTeV result, as shown below.

The Higgs sector of our model consists of three doublets: $`\mathrm{\Phi }=(\varphi ^+,\varphi ^0)`$ with charge 0 and $`\eta _{1,2}=(\eta _{1,2}^+,\eta _{1,2}^0)`$ with charge $`\pm 1`$ under $`U(1)_X`$. The mass matrix spanning $`X`$ and $`Z`$ is then given by

$$_{XZ}^2=\left[\begin{array}{cc}2g_X^2(v_1^2+v_2^2)& g_Xg_Z(v_1^2v_2^2)\\ g_Xg_Z(v_1^2v_2^2)& (g_Z^2/2)(v_0^2+v_1^2+v_2^2)\end{array}\right],$$
(2)
where $`v_0\varphi ^0`$ and $`v_{1,2}^2\eta _{1,2}^0`$ with $`v_0^2+v_1^2+v_2^2=(2\sqrt{2}G_F)^1`$. Assuming that $`v_1v_2`$ so that the $`XZ`$ mixing is small, we then have

$$M_Z^2\frac{1}{2}g_Z^2(v_0^2+2v_1^2),M_X^24g_X^2v_1^2,$$
(3)
with the $`XZ`$ mixing angle given by

$$\mathrm{sin}\theta \frac{g_Xg_X(v_1^2v_2^2)}{M_X^2M_Z^2}.$$
(4)
The effective $`\nu _\mu `$ and $`\overline{\nu }_\mu `$ interactions with quarks has the same structure as the SM, but the effective strength is changed from $`g_Z^2/M_Z^2`$ to

$`g_Z^2\left({\displaystyle \frac{\mathrm{cos}^2\theta }{M_Z^2}}+{\displaystyle \frac{\mathrm{sin}^2\theta }{M_X^2}}\right)2g_Xg_Z\mathrm{sin}\theta \mathrm{cos}\theta \left({\displaystyle \frac{1}{M_Z^2}}{\displaystyle \frac{1}{M_X^2}}\right)`$
$`{\displaystyle \frac{g_Z^2}{M_Z^2}}\left[1+{\displaystyle \frac{2g_X}{g_Z}}\left({\displaystyle \frac{M_Z^2}{M_X^2}}1\right)\mathrm{sin}\theta \right]{\displaystyle \frac{g_Z^2}{M_Z^2}}\rho _\mu .`$ (5)
Note that the factor of 2 in the $`\mathrm{sin}\theta `$ term comes from the fact that $`X`$ couples to $`\nu _\mu `$ with strength 1 whereas $`Z`$ couples to $`\nu _\mu `$ with strength 1/2 $`(=I_3)`$.

In the NuTeV analysis, if $`\rho _\mu =1`$ is assumed, then $`\mathrm{sin}^2\theta _W=0.2277\pm 0.0013\pm 0.0009`$, which deviates from the SM prediction of $`0.2227\pm 0.00037`$ by approximately $`3\sigma `$. On the other hand, if a simultaneous fit to both $`\rho _\mu `$ and $`\mathrm{sin}^2\theta _W`$ is made, they obtain

$$\rho _\mu =0.9983\pm 0.0040,\mathrm{sin}^2\theta _W=0.2265\pm 0.0031,$$
(6)
with a correlation coefficient of 0.85 between the two parameters. They then suggest that one but not both of them may be consistent with SM expectations. Here we choose to consider the deviation of the NuTeV result as being due to $`\rho _\mu `$.

The NuTeV analysis also makes a two-parameter fit in terms of the isoscalar combinations of the effective neutral-current quark couplings, resulting in

$$(g_L^{eff})^2=0.3005\pm 0.0014,(g_R^{eff})^2=0.0310\pm 0.0011,$$
(7)
with a negligibly small correlation coefficient, whereas the SM predictions are

$$(g_L^{eff})_{SM}^2=0.3042,(g_R^{eff})_{SM}^2=0.0301.$$
(8)
Now if we take for example $`\rho _\mu =0.9962`$, then the above two values become $`(g_L^{eff})^2=0.3019`$ and $`(g_R^{eff})^2=0.0299`$, placing them both within $`1\sigma `$ of the experimental measurements.

In atomic parity nonconservation, because $`X`$ does not couple to electrons, we have

$$\rho _e=\mathrm{cos}^2\theta +\mathrm{sin}^2\theta \left(\frac{M_Z^2}{M_X^2}\right)1$$
(9)
to a very good approximation. Thus there should be no deviation from the SM, in agreement with experiment.

From Eq. (5) we obtain

$$\mathrm{sin}\theta =(\rho _\mu 1)\left(\frac{g_Z}{2g_X}\right)\left(\frac{M_X^2}{M_Z^2M_X^2}\right),$$
(10)
which is of order $`10^3`$ for $`\rho _\mu =0.9962`$. This will affect precision data at the $`Z`$ resonance in the following way. First, the observed resonance is of course the physical $`Z`$ boson which has a small $`X`$ component. However, since $`X`$ does not couple to electrons, the production of $`Z`$ is only suppressed by $`\mathrm{cos}^2\theta `$ which is indistinguishable from 1. The decay of $`Z`$ to most fermions is also unaffected because the suppression factor is again just $`\mathrm{cos}^2\theta `$. The exceptions are $`Z\mu ^+\mu ^{},\overline{\nu }_\mu \nu _\mu ,\tau ^+\tau ^{},\overline{\nu }_\tau \nu _\tau `$. Their effective couplings are

$`\mu `$ $`:`$ $`g_V={\displaystyle \frac{1}{2}}+2\mathrm{sin}^2\theta _W2\left({\displaystyle \frac{g_X}{g_Z}}\right)\mathrm{sin}\theta ,g_A={\displaystyle \frac{1}{2}},`$ (11)
$`\nu _\mu `$ $`:`$ $`g_V={\displaystyle \frac{1}{2}}2\left({\displaystyle \frac{g_X}{g_Z}}\right)\mathrm{sin}\theta ,g_A={\displaystyle \frac{1}{2}}2\left({\displaystyle \frac{g_X}{g_Z}}\right)\mathrm{sin}\theta ,`$ (12)
$`\tau `$ $`:`$ $`g_V={\displaystyle \frac{1}{2}}+2\mathrm{sin}^2\theta _W+2\left({\displaystyle \frac{g_X}{g_Z}}\right)\mathrm{sin}\theta ,g_A={\displaystyle \frac{1}{2}},`$ (13)
$`\nu _\tau `$ $`:`$ $`g_V={\displaystyle \frac{1}{2}}+2\left({\displaystyle \frac{g_X}{g_Z}}\right)\mathrm{sin}\theta ,g_A={\displaystyle \frac{1}{2}}+2\left({\displaystyle \frac{g_X}{g_Z}}\right)\mathrm{sin}\theta .`$ (14)
Precision measurements of $`Z`$ couplings at LEP-I give

$$g_V^\mu =0.0359\pm 0.0033,g_V^\tau =0.0366\pm 0.0014,$$
(15)
where the smaller error on $`g_V^\tau `$ is due to the use of $`\tau `$ polarization along with the forward-backward asymmetry. Thus

$$g_V^\tau g_V^\mu =4(g_X/g_Z)\mathrm{sin}\theta =0.0007\pm 0.0036,$$
(16)
adding the two errors in quadrature. Consider now Eq. (10) with the more conservative choice

$$\rho _\mu =0.9976$$
(17)
which is within $`1.6\sigma `$ of the NuTeV measurement of $`(g_L^{eff})^2`$. Comparing it to Eq. (16), we then obtain the following $`2\sigma `$ bounds on $`M_X`$:

$$M_X<72\mathrm{GeV}\mathrm{or}M_X>178\mathrm{GeV}.$$
(18)
A lower bound on $`M_X`$ as a function of $`g_X`$ is also available from LEP-I data on $`Z`$ decay into the 4-muon final state via $`Z\mu ^+\mu ^{}X`$ . For example, if $`g_X=0.2`$, then $`M_X>58`$ GeV. Furthermore, Eq. (3) requires

$$g_X>\frac{g_ZM_X}{2M_Z}.$$
(19)
In Figure 1 we show the above lower limit on $`g_X`$ as well as the $`2\sigma `$ upper limits on $`g_X`$ as functions of $`M_X`$ from $`Z\mu ^+\mu ^{}X`$ decay and the difference of the $`Ze^+e^{}`$ and $`Z\mu ^+\mu ^{}`$ partial widths as the result of the $`X`$ radiative contribution. Details are provided in Ref. . The $`Z`$ decay limit essentially rules out $`M_X<60`$ GeV. The analogous process $`e^+e^{}\mu ^+\mu ^{}X`$ at LEP-II does not improve this bound, as already shown . Thus we conclude that $`M_X`$ between 60 and 72 GeV is still allowed, but perhaps $`M_X>178`$ GeV is more likely.

Going back to Eq. (1) for the muon $`g2`$ discrepancy, we note that there is a theoretical lower bound  of $`1.56\times 10^9`$ in this model, whereas the corrected  range of the experimental discrepancy is $`2.65\pm 1.65\times 10^9`$. This is entirely consistent with the low $`M_X`$ solution, while in the case of the high $`M_X`$ solution, the maximum deviation we get is $`2.7\times 10^9`$. In either case, the $`X`$ boson signal will be too small to be observable at the Fermilab Tevatron, but will be clearly visible at the CERN LHC  via the associated production processes $`u\overline{u}(d\overline{d})\mu \mu X`$ and $`u\overline{d}(d\overline{u})\mu \nu X`$. At a future muon collider, $`X`$ would be copiously produced, especially if it turns out to be light.

To obtain naturally small Majorana neutrino masses, we may add one heavy neutral fermion singlet $`N_R`$ with $`U(1)_X`$ charge 0 as in our previous paper, but then an extra charged scalar boson $`\zeta ^+`$ with charge +1 is needed there to get a second neutrino mass term, i.e. $`\nu _e\nu _\tau `$, radiatively. A possible alternative is to add two $`N_R`$’s. One is assumed to couple only to a linear combination of $`(\nu _\mu \eta _2^0\mu _L\eta _2^+)`$ and $`(\nu _\tau \eta _1^0\tau _L\eta _1^+)`$, and the other to $`(\nu _e\varphi ^0e_L\varphi ^+)`$ as well. Using the canonical seesaw mechanism , this structure allows for the appearance of two massive neutrinos: one is predominantly a mixture of $`\nu _\mu `$ and $`\nu _\tau `$, the other is a linear combination of $`\nu _e`$ and the orthogonal $`\nu _\mu \nu _\tau `$ mixture. This may then lead to a consistent pattern of neutrino masses and mixing for explaining the present atmospheric  and solar  neutrino data.

The interchange symmetry $`\eta _1\eta _2`$ in the Higgs sector allows us to assume $`v_1=v_2`$, but this cannot be maintained for the entire theory. If we try to extend this to the gauge sector, then $`\mu \tau `$ is implied. Hence $`m_\mu m_\tau `$ in the Yukawa sector would break this symmetry. However, the size of this breaking is only of order $`(m_\tau ^2m_\mu ^2)/v_0^2`$ which is smaller than what we require for $`\mathrm{sin}\theta `$. In other words, $`XZ`$ mixing of order $`10^3`$ is a very reasonable value.

In conclusion we have shown in this paper how the gauge symmetry $`L_\mu L_\tau `$ (as realized specifically by us in a previous paper ) explains naturally the recent NuTeV result  on the possible deviation from the Standard Model in $`\nu _\mu `$ and $`\overline{\nu }_\mu `$ scattering with nucleons. Our proposal also explains the possible discrepancy in the recent measurement  of the anomalous magnetic moment of the muon. It further explains why there is no deviation from the Standard Model in atomic parity nonconservation . Our model is constrained by the precision measurements of $`Z\mu ^+\mu ^{}`$ and $`Z\tau ^+\tau ^{}`$, from which we predict that the new gauge boson $`X`$ is likely to have a mass between 60 and 72 GeV, or be heavier than 178 GeV. As such, our model is verifiable experimentally in the future at the LHC.

We thank S. N. Ganguli and A. Gurtu for discussions on the LEP data. This work was supported in part by the U. S. Department of Energy under Grant No. DE-FG03-94ER40837.


text:
# Precision Determination of the Muon 𝑔_𝜇-2 and 𝛼̄_{Q.E.D.}(𝑀²_𝑍) 11footnote 1To appear in the Proceedings of the Budapest 2001 EPS HEP Conference

## 1 Introduction

We present a summary of the results of two recent precise calculations of the muon anomalous magnetic moment ($`g_\mu 2`$ and the electromagnetic coupling on the $`Z`$ ($`\overline{\alpha }_{\mathrm{Q}.\mathrm{E}.\mathrm{D}.}(M_Z^2)`$. The main sources of uncertainty are underlined.

## 2 Calculation of the muon $`g_\mu 2`$

The appearance of a new, very precise measurement of the muon magnetic moment  has triggered the interest in theoretical calculations of this quantity. Particularly, because the experimental figure (we give the result for the anomaly, averaged with older determinations )

$$a_\mu (\text{exp.})=\frac{g_\mu 2}{2}=(\mathrm{116\hspace{0.17em}592\hspace{0.17em}030}\pm 150)\times 10^{11}$$
(1)
lies slightly above theoretical evaluations based on the standard model, as much as $`2.6\sigma `$ in some cases.

It should be noted that all modern theoretical determinations  are compatible among themselves within errors (of order $`100\times 10^{11}`$) and that, with few exceptions, they are also compatible with the experimental result, at the level of $`1.5\sigma `$ or less. Because of this, it is our feeling that a new, complete evaluation would be welcome since, in fact, there exists as yet no calculation that takes fully into account all theoretical constraints and all the new experimental data. These data allow an improved evaluation of the low energy hadronic contributions to $`a_\mu `$, both directly from $`e^+e^{}`$ annihilations (in the $`\rho `$ region  and around the $`\omega `$ and $`\varphi `$ resonances ) and, indirectly, from $`\tau `$ decays  and, also indirectly, from measurements of the pion form factor in the spacelike region . We have also used perturbative QCD in the region $`s2GeV^2`$ (away from quark thresholds). In particular, the recent LEP $`\tau \nu _\tau +\mathrm{hadrons}`$, and BES $`e^+e^{}\mathrm{hadrons}`$ data justify the QCD result for $`2GeV^2s3^2GeV^2`$. Moreover, the BES  data, covering $`e^+e^{}`$ annihilations in the vicinity of $`\overline{c}c`$ threshold, permit a reliable evaluation of the corresponding hadronic pieces. Last but not least, the next order radiative corrections have been taken into account.

In the next sections we concentrate in the two most delicate pieces, from the point of view of the final uncertainty: the pion form factor and the hadronic light-by-light scattering contributions.

### 2.1 The pion form factor

To obtain $`F_\pi (s)`$ we fit the recent Novosibirsk data  on $`e^+e^{}\pi ^+\pi ^{}`$ and the tau decay data of Aleph and Opal . We take into account, at least partially, isospin breaking effects by allowing different masses and widths for the $`\rho ^0`$, $`\rho ^+`$ resonances. Moreover, and to get a good grip in the low energy region where data are inexistent or very poor, we also fit $`F_\pi (s)`$ at spacelike $`s`$ . This is possible in our approach because we use an expression for $`F_\pi `$ that takes fully into account its analyticity properties. To be precise, we use the Omnès-Muskhelishvili method. We write

$$F_\pi (s)=G(s)J(s).$$

Here $`J`$ is expressed in terms of the P-wave $`\pi \pi `$ phase shift, $`\delta _1^1`$, as

$$J(s)=e^{1\delta _1^1(s_0)/\pi }\left(1\frac{s}{s_0}\right)^{[1\delta _1^1(s_0)/\pi ]s_0/s}\left(1\frac{s}{s_0}\right)^1\mathrm{exp}\left\{\frac{s}{\pi }_{4m_\pi ^2}^{s_0}𝑑s^{}\frac{\delta _1^1(s^{})}{s^{}(s^{}s)}\right\}.$$

$`s_0`$ is the energy at which inelasticity starts becoming important (in practice, above the percent level); we will take $`s_0=1.1GeV^2`$ in actual calculations.

Because of the equality of the phase of $`J(s)`$ and the phase of $`F_\pi (s)`$ below $`s=s_0`$, it follows that $`G(s)`$ will be an analytic function also for $`4m_\pi ^2ss_0`$, so in the whole $`s`$ plane except in a cut from $`s=s_0`$ to $`+\mathrm{}`$. If we now make the conformal transformation

$$z=\frac{\frac{1}{2}\sqrt{s_0}\sqrt{s_0s}}{\frac{1}{2}\sqrt{s_0}+\sqrt{s_0s}}$$

then, as a function of $`z`$, $`G`$ will be analytic in the unit disk and we can thus write a convergent Taylor series for it. Incorporating the condition $`G(0)=1`$, that follows from $`F_\pi (0)=1`$, and undoing the transformation, we have

$$G(s)=1+c_1\left[\frac{\frac{1}{2}\sqrt{s_0}\sqrt{s_0s}}{\frac{1}{2}\sqrt{s_0}+\sqrt{s_0s}}+\frac{1}{3}\right]+c_2\left[\left(\frac{\frac{1}{2}\sqrt{s_0}\sqrt{s_0s}}{\frac{1}{2}\sqrt{s_0}+\sqrt{s_0s}}\right)^2\frac{1}{9}\right]+\mathrm{},$$

$`c_1,c_2,\mathrm{}`$ free parameters. Actually, only two terms will be necessary to fit the data.

To obtain $`J`$, we use the effective range theory to write

$$\mathrm{cot}\delta _1^1(s)=\frac{s^{1/2}}{2k^3}(m_\rho ^2s)\widehat{\psi }(s),k=\frac{\sqrt{s4m_\pi ^2}}{2};$$

where we have extracted the zero corresponding to the rho resonance. The effective range function $`\widehat{\psi }(s)`$ is analytic in the full $`s`$ plane except for a cut for $`[\mathrm{},0]`$ and the inelastic cut $`[s_0,+\mathrm{}]`$. We can profit from this analyticity by making again a conformal transformation into the unit circle, which is now given by

$$w=\frac{\sqrt{s}\sqrt{s_0s}}{\sqrt{s}+\sqrt{s_0s}}.$$

We expand $`\widehat{\psi }`$ in a convergent series of powers of $`w`$; undoing the transformation we have

$$\delta _1^1(s)=\mathrm{Arc}\mathrm{cot}\left\{\frac{s^{1/2}}{2k^3}(m_\rho ^2s)\left[b_0+b_1\frac{\sqrt{s}\sqrt{s_0s}}{\sqrt{s}+\sqrt{s_0s}}+\mathrm{}\right]\right\}.$$

For the actual fits, only $`b_0,b_1`$, and $`m_\rho `$ are needed as parameters.

We give results both using only $`F_\pi `$ from $`e^+e^{}`$ annihilations, or involving also the decay $`\tau ^+\overline{\nu }_\tau \pi ^+\pi ^0`$, which last we consider to be our best estimates:

$$10^{11}\times a_\mu (2\pi ;t0.8GeV^2)=\{\begin{array}{cc}\mathrm{4\hspace{0.17em}774}\pm 31,(\mathrm{TY1}:e^+e^{}+\tau +\mathrm{spacelike}).\hfill & \\ \mathrm{4\hspace{0.17em}754}\pm 55,(\mathrm{TY2}:e^+e^{}+\mathrm{spacelike}).\hfill & \end{array}$$
(2)
As a byproduct of the calculation we can also give precise values for the $`\rho ^0`$, $`\rho ^+`$ masses and widths or the pion electromagnetic radius; see .

### 2.2 Hadronic light-by-light contributions

A contribution to $`a_\mu `$ in a class by itself is the hadronic light-by-light scattering. It can be evaluated only using models. One can make a chiral model calculation, in the Nambu–Jona-Lasinio version or the chiral perturbation theory variety, with a cut-off, or one can use a constituent quark model. The result depends on the cut-off (for the chiral calculation) or on the constituent mass chosen for the quarks. What is worse, for the chiral model, a recent calculation  has challenged even its sign. The variation from one model to the other is too large for comfort. We take as a typical value for the chiral calculation ,

$$10^{11}\times a_\mu (\text{Hadronic light-by-light})=86\pm 25\text{(Chiral calculation; BPP, HKS)}.$$

In the calculation by  one simply reverses the sign.

For the constituent quark model we use the results of Laporta and Remiddi . The contribution to $`a_\mu `$ of light-by-light scattering, with a loop with a fermion of charge $`Q_i`$, and mass $`m_i`$ larger than the muon mass, is

$$10^{11}\times a_\mu (\text{Hadronic light-by-light})=+92\pm 20\text{(Quark const. model)}$$

and the error is estimated by varying $`m_{u,d}`$ by 10%. This is essentially the same result following the chiral model calculation by .

Besides the sign problem, the variations are unfortunate; one expects the chiral calculation to be valid for small values of the virtual photon momenta, and the constituent model to hold for large values of the same. Thus, almost half of the contribution to $`a_\mu (\text{Hadronic light-by-light})`$ in the chiral calculation comes from a region of momenta above $`0.5`$ GeV, where the chiral perturbation theory starts to fail, while for this range of energies, and at least for the imaginary part of (diagonal) light-by-light scattering, the quark model reproduces reasonably well the experimental data. Because of the uncertainties, we will give results with the extreme choices.

### 2.3 Results

Our results for the hadronic part of the anomaly depend on which model one believes for the hadronic light-by-light contribution, discussed in the previous section. So we write,

$$10^{11}\times a_\mu (\text{Hadronic})=\{\begin{array}{cc}\mathrm{6\hspace{0.17em}993}\pm 69(\mathrm{TY1};\mathrm{Q}.\mathrm{c}.\mathrm{m}.)\hfill & \\ \mathrm{6\hspace{0.17em}815}\pm 71(\mathrm{TY1};\mathrm{Ch}.\mathrm{m}.(\mathrm{BPP},\mathrm{HKS}))\hfill & \\ & \\ \mathrm{6\hspace{0.17em}973}\pm 99(\mathrm{TY2};\mathrm{Q}.\mathrm{c}.\mathrm{m}.)\hfill & \\ \mathrm{6\hspace{0.17em}795}\pm 100(\mathrm{TY2};\mathrm{Ch}.\mathrm{m}.(\mathrm{BPP},\mathrm{HKS}))\hfill & \end{array}$$
(3)
Note that in $`a_\mu (\text{Hadronic})`$ we include all hadronic contributions, $`O(\alpha ^3)`$ as well as $`O(\alpha ^2)`$. In eqn. (3) “Q.c.m.” means that the light-by-light hadronic contribution was calculated with the quark constituent model, and “Ch.m.” that a chiral model was used. The errors include the statistical errors, as well as the estimated systematic and theoretical ones. This is to be compared with the value deduced from eqn. (1) and electroweak corrections

$$10^{11}\times a_\mu ^{\mathrm{exp}.}(\text{Hadronic})=\mathrm{7\hspace{0.17em}174}\pm 150,$$

from which eqn. (3) differs by something between $`1.1\sigma `$ and $`2.1\sigma `$. Our numbers are compared to other evaluations in Fig. 1. Whether these results should be interpreted as providing further experimental validation of the standard model, or one can consider them as “harbingers of new physics”, we will leave for the reader to decide.

## 3 Calculation of $`\overline{\alpha }_{\mathrm{Q}.\mathrm{E}.\mathrm{D}.}`$ on the $`Z`$

With a simple change of integration kernel the previous analysis can be extended to evaluate the hadronic contribution to the QED running coupling, $`\overline{\alpha }_{\mathrm{Q}.\mathrm{E}.\mathrm{D}.}(s)`$, in particular for $`s=M_Z^2`$; an important quantity that enters into precision evaluations of electroweak observables.

By using a dispersion relation one can write this hadronic contribution at energy squared $`s`$ as

$$\mathrm{\Delta }_{\mathrm{had}}\alpha (s)=\frac{s\alpha }{3\pi }_{4m_\pi ^2}^{\mathrm{}}𝑑t\frac{R(t)}{t(ts)}$$
(4)
with

$$R(t)=\frac{\sigma (e^+e^{}\mathrm{hadrons};t)}{\sigma ^{(0)}(e^+e^{}\mu ^+\mu ^{};t)},\sigma ^{(0)}(e^+e^{}\mu ^+\mu ^{};t)\frac{4\pi \alpha ^2}{3t}$$

and the integral in eqn. (4) has to be understood as a principal part integral. This is similar to the Brodsky–de Rafael expression for the hadronic contribution to the muon magnetic moment anomaly,

$$a_\mu (\text{h.v.p.})=_{4m_\pi ^2}^{\mathrm{}}𝑑tK(t)R(t).$$

Therefore, we can carry over all the work from the previous section with the simple replacement

$$K(t)\frac{s\alpha }{3\pi }\frac{1}{t(ts)}.$$

We find to next to leading order in $`\alpha `$,

$$10^5\times \mathrm{\Delta }_{\mathrm{had}}\alpha (M_Z^2)=2740\pm 12$$

or, excluding the top quark contribution,

$$10^5\times \mathrm{\Delta }_{\mathrm{had}}\alpha ^{(5)}(M_Z^2)=2747\pm 12.$$

Other recent determinations can be found in .

## 4 Conclusions

In conclusion, we have performed a detailed evaluation of the hadronic contributions to the muon $`g_\mu 2`$ , and to the running electromagnetic coupling .

In the case of the muon $`g_\mu 2`$ an important limitation of the attainable accuracy is the light-by-light scattering contribution; in particular, a third chiral model calculation is required to clear the sign controversy.


text:
# Size of the neutralino–nucleon cross–section in the light of a new determination of the pion–nucleon sigma term

## I Introduction

The evaluation of detection rates for relic neutralinos in WIMP direct experiments and in measurements of signals due to possible neutralino–neutralino annihilation in Earth and Sun  requires a reliable calculation of the neutralino–nucleon cross–section.

In Ref.  we stressed that sizeable uncertainties in the neutralino–nucleon cross–section are due to the estimate of the hadronic matrix elements, which describe the quark density contents in the nucleon; this point was also remarked in Refs. . In  we quantified the range of these uncertainties and showed that, within their ranges, the couplings of neutralinos to nucleon may be sizeably stronger than some of those currently employed in the literature. Subsequently, an analysis of this problem was also undertaken in Refs. .

One of the crucial ingredients in the calculation of the hadronic matrix elements is the so–called pion–nucleon sigma term. Its numerical derivation from experimental data of pion–nucleon scattering is rather involved, and thus the origin of considerable uncertainties. Recent experimental results in pion–nucleon scattering have now prompted a new determination of the pion–nucleon sigma term . In the present paper, we consider the implications of these new inputs for the neutralino–nucleon cross–section; in particular, we show that the new data reinforce our previous conclusions of Ref. . As compared to Ref. , the present paper, apart from the use of the new evaluation for the sigma term, contains also some other relevant updatings, mainly: i) radiatively corrected Higgs–quark couplings  and ii) recent experimental bounds on Higgs masses and supersymmetric parameters from LEP2  and CDF .

## II Neutralino–nucleon elastic cross–section

The neutralino–nucleon scalar cross–section may be written as

$$\sigma _{\mathrm{scalar}}^{(\mathrm{nucleon})}=\frac{8G_F^2}{\pi }M_Z^2m_{\mathrm{red}}^2\left[\frac{F_hI_h}{m_h^2}+\frac{F_HI_H}{m_H^2}+\frac{M_Z}{2}\underset{q}{}<N|\overline{q}q|N>\underset{i}{}P_{\stackrel{~}{q}_i}(A_{\stackrel{~}{q}_i}^2B_{\stackrel{~}{q}_i}^2)\right]^2.$$
(1)
The first two terms inside the brackets refer to the diagrams with exchanges of the two CP–even neutral Higgs bosons, $`h`$ and $`H`$, in the t–channel (the exchange diagram of the CP–odd one, $`A`$, is strongly kinematically suppressed and then omitted here)  and the third term refers to the graphs with squark–exchanges in the s– and u–channels . The mass $`m_{red}`$ is the neutralino–nucleon reduced mass. Since, for simplicity, in the present paper we explicitly discuss only the Higgs–mediated terms, we do not report the expressions for the squark propagator $`P_{\stackrel{~}{q}_i}`$ and for the couplings $`A_{\stackrel{~}{q}_i}`$, $`B_{\stackrel{~}{q}_i}`$, which may be found in Ref.. However, in the numerical results reported in this paper also the squark–exchange terms are included. The quark matrix elements $`<N|\overline{q}q|N>`$ are meant over the nucleonic state.

The quantities $`F_{h,H}`$ and $`I_{h,H}`$ are defined as follows

$`F_h`$ $`=`$ $`(a_1\mathrm{sin}\theta _W+a_2\mathrm{cos}\theta _W)(a_3\mathrm{sin}\alpha +a_4\mathrm{cos}\alpha )`$ (2)
$`F_H`$ $`=`$ $`(a_1\mathrm{sin}\theta _W+a_2\mathrm{cos}\theta _W)(a_3\mathrm{cos}\alpha a_4\mathrm{sin}\alpha )`$ (3)
$`I_{h,H}`$ $`=`$ $`{\displaystyle \underset{q}{}}k_q^{h,H}m_qN|\overline{q}q|N,`$ (4)
where the $`a_i`$’s are the coefficients in the definition of the neutralino as the lowest–mass linear superposition of photino ($`\stackrel{~}{\gamma }`$), zino ($`\stackrel{~}{Z}`$) and the two higgsino states ($`\stackrel{~}{H}_1^{}`$, $`\stackrel{~}{H}_2^{}`$)

$$\chi a_1\stackrel{~}{\gamma }+a_2\stackrel{~}{Z}+a_3\stackrel{~}{H}_1^{}+a_4\stackrel{~}{H}_2^{}.$$
(5)
The coefficients $`k_q^{h,H}`$ for up–type and down–type quarks are given in Table 1, in terms of the angle $`\beta `$, defined as $`\mathrm{tan}\beta =<H_2^0>/<H_1^0>`$ and the angle $`\alpha `$, which rotates $`H_1^{(0)}`$ and $`H_2^{(0)}`$ into $`h`$ and $`H`$. In Table I also included are the coefficients $`k_q^{h,H}`$ for the CP–odd Higgs boson $`A`$, which are important for the evaluation of the neutralino relic abundance. The entries include those radiative corrections which may be sizeable at large $`\mathrm{tan}\beta `$. These corrections affect the couplings to down–type quarks $`k_{dtype}`$, and are parametrized in terms of the quantity $`ϵ1/(1+\mathrm{\Delta })`$, where $`\mathrm{\Delta }`$ enters in the relationship between the fermion running masses $`m_d`$ and the corresponding Yukawa couplings $`h_d`$ :

$$m_d=h_d<H_0^1>(1+\mathrm{\Delta })$$
(6)
These corrections take contributions mainly from gluino–squark, chargino–squark and neutralino–stau loops .

## III Evaluation of the quantities $`m_qN|\overline{q}q|N`$

For the calculation of the quantities $`m_qN|\overline{q}q|N`$ it is useful to express them in terms of the pion–nucleon sigma term

$$\sigma _{\pi N}=\frac{1}{2}(m_u+m_d)<N|\overline{u}u+\overline{d}d|N>,$$
(7)
of the quantity $`\sigma _0`$, related to the size of the SU(3) symmetry breaking, and defined as

$$\sigma _0\frac{1}{2}(m_u+m_d)<N|\overline{u}u+\overline{d}d2\overline{s}s|N>,$$
(8)
and of the ratio $`r=2m_s/(m_u+m_d)`$.

Assuming isospin invariance for quarks $`u`$ and $`d`$, the quantities $`m_qN|\overline{q}q|N`$ for light quarks may be written as

$`m_u<N|\overline{u}u|N>`$ $``$ $`m_d<N|\overline{d}d|N>{\displaystyle \frac{1}{2}}\sigma _{\pi N}`$ (9)
$`m_s<N|\overline{s}s|N>`$ $``$ $`{\displaystyle \frac{1}{2}}r(\sigma _{\pi N}\sigma _0).`$ (10)
For the heavy quarks $`c`$, $`b`$, $`t`$, use of the heavy quark expansion provides

$`m_c<N|\overline{c}c|N>`$ $``$ $`m_b<N|\overline{b}b|N>m_t<N|\overline{t}t|N>`$ (11)
$``$ $`{\displaystyle \frac{2}{27}}\left[m_N\sigma _{\pi N}+{\displaystyle \frac{1}{2}}r(\sigma _{\pi N}\sigma _0)\right],`$ (12)
where $`m_N`$ is the nucleon mass. The quantities $`I_{h,H}`$ can then be rewritten as

$$I_{h,H}=k_{u\mathrm{type}}^{h,H}g_u+k_{d\mathrm{type}}^{h,H}g_d,$$
(13)
where

$`g_u`$ $``$ $`m_l<N|\overline{l}l|N>+\mathrm{\hspace{0.33em}2}m_h<N|\overline{h}h|N>`$ (14)
$``$ $`{\displaystyle \frac{4}{27}}(m_N+{\displaystyle \frac{19}{8}}\sigma _{\pi N}{\displaystyle \frac{1}{2}}r(\sigma _{\pi N}\sigma _0)),`$ (15)
$`g_d`$ $``$ $`m_l<N|\overline{l}l|N>+m_s<N|\overline{s}s|N>+m_h<N|\overline{h}h|N>`$ (16)
$``$ $`{\displaystyle \frac{2}{27}}(m_N+{\displaystyle \frac{23}{4}}\sigma _{\pi N}+{\displaystyle \frac{25}{4}}r(\sigma _{\pi N}\sigma _0));`$ (17)
$`l`$ stands for light quarks ($`l=u,d`$) and $`h`$ denotes the heavy ones ($`h=c,b,t`$).

We recall that from the values of $`\sigma _{\pi N}`$ and $`\sigma _0`$ one can derive the fractional strange–quark content of the nucleon $`y`$

$$y=2\frac{<N|\overline{s}s|N>}{<N|\overline{u}u+\overline{d}d|N>},$$
(18)
using the expression

$$y=1\frac{\sigma _0}{\sigma _{\pi N}}.$$
(19)
### A Values for $`\sigma _{\pi N}`$, $`\sigma _0`$ and $`r`$

The range of $`\sigma _{\pi N}`$ we used previously in Ref.  is

$$41\mathrm{MeV}<\sigma _{\pi N}<57\mathrm{MeV}.$$
(20)
This was derived from the pion–nucleon scattering amplitude, calculated at the so–called Cheng–Dashen point by Koch , and the evolution of the nucleon scalar form factor, as a function of the momentum transfer from $`t=2m_\pi ^2`$ to $`t=0`$, evaluated in Ref. .

We now consider the new determination of $`\sigma _{\pi N}`$ presented in Ref. . The George Washington University/TRIUMF group, using up–dated pion–nucleon scattering data  and a new partial–wave and dispersion relation analysis program, has derived a range for $`\sigma _{\pi N}`$

$$55\mathrm{MeV}<\sigma _{\pi N}<73\mathrm{MeV},$$
(21)
which turns out to be sizeably larger than the one of Eq. (20). Values of the nucleon scalar form factor at the Cheng–Dashen point higher than those of Ref.  were also reported in Ref. .

In the present paper we consider the effect of employing the new range for $`\sigma _{\pi N}`$, given in Eq. (21), in the evaluation of the quantities $`m_qN|\overline{q}q|N`$’s.

Here, $`\sigma _0`$ is taken in the range

$$\sigma _0=30÷40\mathrm{MeV}.$$
(22)
Also the mass ratio $`r=2m_s/(m_u+m_d)`$ may be affected by significant uncertainties . However, here, for consistency with some of the previous determinations, we use the default value $`r=25`$.

It has to be noted that combining together values of $`\sigma _{\pi N}`$ and $`\sigma _0`$ within their ranges in Eqs. (21) and (22) leads to rather large values for the fractional strange–quark content of the nucleon, as given by Eq. (19). This is a puzzle that urges further investigation in hadron physics. In what follows, in the variations of the quantities $`\sigma _{\pi N}`$ and $`\sigma _0`$ we impose the constraint that, anyway, $`y0.5`$.

### B Values for $`m_qN|\overline{q}q|N`$, $`g_u`$ and $`g_d`$

Inserting the numerical values of $`\sigma _{\pi N}`$, $`\sigma _0`$ and $`r`$ in the expressions given at the beginning of this section, we finally obtain estimates for the quantities $`m_qN|\overline{q}q|N`$, $`g_u`$ and $`g_d`$. Since $`m_s<N|\overline{s}s|N>`$ is the most important term among the $`m_q<N|\overline{q}q|N>`$’s , unless $`\mathrm{tan}\beta `$ is very small, the extremes of the range of the neutralino–nucleon cross–section are provided by the extremes of the range for $`m_s<N|\overline{s}s|N>`$. These, in turn, are given by: $`(m_s<N|\overline{s}s|N>)^{min}\frac{1}{2}r(\sigma _{\pi N}^{}{}_{}{}^{min}\sigma _{0}^{}{}_{}{}^{max})`$ and $`(m_s<N|\overline{s}s|N>)^{max}\frac{1}{4}r\sigma _{\pi N}^{}{}_{}{}^{max}`$ (to satisfy the constraint $`y0.5`$). We call Set b the one with $`m_s<N|\overline{s}s|N>=(m_s<N|\overline{s}s|N>)^{min}`$ and Set c the one with $`m_s<N|\overline{s}s|N>=(m_s<N|\overline{s}s|N>)^{max}`$; we denote by Set a the reference set of Ref. . The values of the quantities $`m_qN|\overline{q}q|N`$, $`g_u`$ and $`g_d`$ are given in Table II.

We note that, with the new values of $`\sigma _{\pi N}`$, the coefficient $`g_d`$, which usually dominates in the neutralino–nucleon cross–section, turns out to fall in the range

$$266\mathrm{MeV}<g_d<523\mathrm{MeV},$$
(23)
to be compared with the reference value $`g_d=241`$ MeV of Set a. Thus, for most supersymmetric configurations, one expects $`(\sigma _{scalar}^{(nucleon)})_{Setb}/(\sigma _{scalar}^{(nucleon)})_{Seta}1.5`$ and $`(\sigma _{scalar}^{(nucleon)})_{Setc}/(\sigma _{scalar}^{(nucleon)})_{Seta}6`$.

## IV Results and conclusions

The supersymmetric model employed in the present paper for the calculation of the neutralino–nucleon cross–section and for the neutralino relic abundance is the one defined in Ref.  and denoted there as effMSSM. We refer to  for details on the theoretical aspects and on the updated experimental bounds.

Inserting the values of Table II into Eq.(1), one obtains the results displayed in Fig. 1. The two ratios in the cross–sections are plotted as a function of $`\xi \sigma _{scalar}^{(nucleon)}`$, where $`\xi `$ is a rescaling factor for the neutralino local density. $`\xi `$ is taken to be $`\xi =\mathrm{min}\{1,\mathrm{\Omega }_\chi h^2/(\mathrm{\Omega }_mh^2)_{min}\}`$, in order to have rescaling, when $`\mathrm{\Omega }\chi h^2`$ turns out to be less than $`(\mathrm{\Omega }_mh^2)_{min}`$ (here $`(\mathrm{\Omega }_mh^2)_{min}`$ is set to the value 0.05).

The use of $`\xi \sigma _{scalar}^{(nucleon)}`$ instead of simply $`\sigma _{scalar}^{(nucleon)}`$ allows one to better identify the range of sensitivity in current WIMP direct experiments , in Fig. 1. Taking into account astrophysical uncertainties , this range may established to be

$$410^{10}\mathrm{nbarn}<\xi \sigma _{\mathrm{scalar}}^{(\mathrm{nucleon})}<210^8\mathrm{nbarn},$$
(24)
for WIMP masses in the interval $`40\mathrm{GeV}<m_W<200\mathrm{GeV}`$. In Fig. 1 we notice that, in the sensitivity range of Eq.(24), the cross–section ratios are actually of the sizes obtained by our previous estimate based on dominance of the term $`g_d`$, i. e. $`(\sigma _{scalar}^{(nucleon)})_{Setb}/(\sigma _{scalar}^{(nucleon)})_{Seta}1.5`$ and $`(\sigma _{scalar}^{(nucleon)})_{Setc}/(\sigma _{scalar}^{(nucleon)})_{Seta}6`$.

We turn now to the implications that our present analysis has for the cosmological properties of relic neutralinos, as explored by present WIMP direct detection experiments. Fig. 2 provides the essential information. The detection of relic neutralinos would be of great relevance even if these particles constitute only a subdominant dark matter population ; however, it is obvious that the most attractive case is when the neutralino relic abundance falls into the interval of cosmological interest: $`0.05<\mathrm{\Omega }_mh^2<0.3`$. Fig. 2 shows to what extent the most interesting region: $`0.05<\mathrm{\Omega }_mh^2<0.3`$ and $`410^{10}\mathrm{nbarn}<\xi \sigma _{\mathrm{scalar}}^{(\mathrm{nucleon})}<210^8\mathrm{nbarn}`$, is covered by neutralino configurations, depending on the values of the input parameters of Table II. The up–right frontier of the scatter plots moves progressively upward, as we move from Set a of inputs to Set b and to Set c.

The main results of the present paper may be summarized as follows:

* The new range of the pion–nucleon sigma term favours values of the neutralino–nucleon cross–section which are sizeably larger (a factor of 1.5 to 6) than some of the current estimates. However, a word of caution has to be said here: the new derivation of the pion–nucleon sigma term implies a rather high value for the fractional strange–quark content of the nucleon; this point requires further investigation in the framework of hadron physics. Thus, for instance, the consistency between the present indication for higher values of the pion–nucleon sigma term and the determination of the range of $`\sigma _0`$ has to be understood.
* Uncertainties implied by the hadronic quantities for the neutralino–nucleon cross–section are still quite sizeable; these have to be taken in due consideration in evaluations of the neutralino event rates.
* The previous remark applies not only to event rates for WIMP direct detection, but also to the evaluation of the neutrino fluxes expected at neutrino telescopes, as a consequence of possible neutralino–neutralino annihilations inside the Earth and the Sun. In fact, these fluxes depend on the capture rate of neutralinos by the macroscopic bodies; in turn, this rate depends on the neutralino–nucleon cross–section.
* A larger size of the neutralino–nucleon cross-section, as indicated by the new hadronic data, implies that present WIMP direct experiments (and some indirect measurements) explore larger sectors of the supersymmetric space where the neutralino may be of cosmological interest.


text:
# Radiation drag driven mass accretion in clumpy interstellar medium: implications for the supermassive black hole-to-bulge relation

## 1 Introduction

Recently, Kormendy $`\&`$ Richstone (1995) have pioneeringly suggested that the mass of a supermassive black hole (BH) does correlate linearly with the mass of the hosting bulge. (It is noted that the term of a bulge is used to mean a whole galaxy for an elliptical galaxy in this paper as is often so.) Further high-quality observations of the galactic center using stellar dynamics, gasdynamics, and maser dynamics (Miyoshi et al. 1995; Magorrian et al. 1998; Richstone et al. 1998; Ho 1999; Wandel 1999; Kormendy & Ho 2000; Ferrarese et al. 2000; Gebhardt et al. 2000a; Sarzi et al. 2001; Merritt & Ferrarese 2001a) allow us to make a detailed demography of supermassive BHs. The recent findings are the following: (1) The BH mass exhibits a linear relation to the bulge mass for a wide range of BH mass with a median BH mass fraction of $`f_{\mathrm{BH}}M_{\mathrm{BH}}/M_{\mathrm{bulge}}=0.0010.006`$ (Kormendy & Richstone 1995; Richstone et al. 1998; Magorrian et al. 1998; Gebhardt et al. 2000a; Ferrarese & Merritt 2000; Merritt & Ferrarese 2001a). (2) The BH mass correlates with the velocity dispersion of bulge stars with a power-law relation as $`M_{\mathrm{BH}}\sigma ^n`$, $`n=3.75`$ (Gebhardt et al. 2000a) or 4.72 (Ferrarese & Merritt 2000; Merritt & Ferrarese 2001a, 2001b). (3) $`f_{\mathrm{BH}}`$ tends to grow with the age of youngest stars in a bulge until $`10^9`$ yr (Merrifield, Forbes & Terlevich 2000). (4) In disc galaxies, the mass ratio is significantly smaller than 0.001 if the disc stars are included (Salucci et al. 2000; Sarzi et al. 2001). (5) For quasars, the $`f_{\mathrm{BH}}`$ is on a similar level to that for elliptical galaxies (Laor 1998; Mclure & Dunlop 2001a; Wandel 2001). (6) The $`f_{\mathrm{BH}}`$ in Seyfert 1 galaxies is under debate, which may be considerably smaller than 0.001 (Wandel 1999; Gebhardt et al. 2000a) or similar to that for ellipticals (McLure & Dunlop 2001a,2001b; Wandel 2001), while the BH mass-to-velocity dispersion relation in Seyfert 1 galaxies seems to hold good in a similar way to elliptical galaxies (Gebhardt et al. 2000b; Nelson 2000; Ferrarese et al. 2001). These BH-to-bulge correlations suggest that the formation of a supermassive BH is physically connected with the formation of a galactic bulge.

So far, very little is understood about the physical mechanism to produce such correlations, although some theoretical models has been proposed (Silk $`\&`$ Rees 1998; Ostriker 2000; Adams, Graff $`\&`$ Richstone 2001). Recently, as a possible mechanism to work in a spheroidal system, Umemura (2001) has considered the effects of radiation drag. The radiation drag is a relativistic effect, which may extract the angular momentum effectively in a spheroidal system like a bulge, so that plenty of interstellar medium (ISM) could accrete on to the galactic center. Obviously, the radiation drag is inefficient in present-day elliptical galaxies or galactic bulges, since they possess little ISM. If the contents of a supermassive BH are initially in the form of ISM, however, the bulge must have been optically thick in the early stage:

$$\tau \chi \rho r_\mathrm{b}=1.0\left(\frac{\chi }{100\mathrm{cm}^2\mathrm{g}^1}\right)\left(\frac{M_{\mathrm{gas}}}{10^9M_{}}\right)\left(\frac{r_\mathrm{b}}{3\mathrm{k}\mathrm{p}\mathrm{c}}\right)^2,$$

where $`\chi `$ is the mass extinction coefficient due to dust opacity, $`\rho `$ is the density of the ISM, $`M_{\mathrm{gas}}`$ is the mass of the ISM, and $`r_\mathrm{b}`$ is the bulge radius. If a considerable amount of gas is expelled by a galactic wind at some stage, the optical depth should be still larger before the wind. If the radiation drag works efficiently in an optically thick medium, the rate of mass accretion induced by the radiation drag is maximally $`L_{\mathrm{bol}}/c^2`$ (Umemura, Fukue & Mineshige 1997, 1998; Fukue, Umemura & Mineshige 1997), where $`L_{\mathrm{bol}}`$ is the bolometric luminosity. Umemura (2001) has found that, if the maximal drag efficiency is achieved, the resultant BH-to-bulge mass ratio is basically determined by the energy conversion efficiency of the nuclear fusion from hydrogen to helium, i.e., 0.007. However, it is not very clear whether this mechanism really works efficiently in realistic situations.

In this paper, we investigate in detail the efficiency of the radiation drag in an optically thick ISM to test whether the radiation drag model is promising to account for the putative BH-to-bulge correlations. In particular, we concentrate our attention on the effects of the inhomogeneity in the ISM. The model for the chemical evolution of elliptical galaxies suggests that an elliptical galaxy is initiated by a starburst in its early stage of ($`<10^7`$yr), and evolves passively after a galactic wind event at a few $`\times `$$`10^8`$ yr (Arimoto $`\&`$ Yoshii, 1986, 1987; Kodama $`\&`$ Arimoto, 1997; Mori et al. 1997). Also, in nearby starburst galaxies that have been studied, the ISM is highly clumpy (Sanders et al. 1988; Gordon et al. 1997). Thus, if we consider the radiation drag in the early phase of bulge evolution, we should consider an inhomogeneous optically thick ISM. In this paper, to elucidate the mutual effect between the clumpiness of ISM and the optical depth on the radiation drag efficiency, we build up a simple model of the bulge system and accurately solve the radiation transfer in a clumpy ISM.

The paper is organized as follows. In Section 2, we construct the model of a galactic bulge. In Section 3, the basic equations for the ISM are provided. In Section 4, the angular momentum transfer efficiency is assessed in a uniform ISM. In Section 5, we investigate the angular momentum transfer efficiency by solving the radiation transfer in a clumpy ISM, and elucidate the relationship between the clumpiness of the ISM and the angular momentum transfer efficiency. In Section 6, we give implications for the correlation between the supermassive BH mass and the bulge mass. In addition, we discuss further effects that would give significant influence on the BH mass. Section 7 is devoted to our conclusions.

## 2 Model

We assume that a spherical galactic bulge consists of three components, that is, dark matter, stars, and a dusty ISM. The bulge radius $`r_\mathrm{b}`$ is set to be 1-10kpc by taking account of the observed sizes of elliptical galaxies or bulges in spiral galaxies. Dark matter is distributed uniformly inside the bulge. The mass of dark matter component, $`M_{DM}`$, within the bulge size is equal to the stellar mass of the galactic bulge, $`M_{\mathrm{bulge}}`$. As for the stellar component, we assume star clusters with a specific stellar initial mass function (see below). The star clusters are distributed uniformly inside the galactic bulge. Hereafter, ‘a star’ in this paper means ‘a star cluster’. Here, $`N_{}(=100)`$ stars are distributed randomly. For a dusty ISM, we consider the two cases: one is a uniformly distributed ISM, and the other is a clumpy ISM. In the case of clumpy ISM, $`N_\mathrm{c}(=10^4)`$ identical clouds are distributed randomly. The density $`\rho _{\mathrm{gas}}`$ in a cloud is assumed to be uniform. The size of a gas cloud, $`r_\mathrm{c}`$, is a parameter. Then,the optical depth of a gas cloud is $`\overline{\tau }=\chi \rho _{\mathrm{gas}}r_\mathrm{c}`$, where $`\chi `$ is the mass extinction. We suppose that the stars and the ISM clouds corotate with the angular velocity corresponding to the angular momentum obtained by the tidal torque at the linear stage of density fluctuations. Quantitatively, the angular momentum is given by the spin parameter $`\lambda =(J_\mathrm{T}|E_\mathrm{T}|^{1/2})/(GM_\mathrm{T}^{5/2})=0.05`$, where $`J_\mathrm{T}`$, $`E_\mathrm{T}`$, and $`M_\mathrm{T}`$ are respectively the total angular momentum, energy, and mass (Barnes & Efstathiou 1987; Heavens & Peacock 1988). Here, the rigid rotation is assumed.

The mass range of galactic bulges is postulated to be $`10^{613}M_{}`$. (However, in the present analysis, it is not very important to specify $`M_{\mathrm{bulge}}`$, because the results are scaled with $`M_{\mathrm{bulge}}`$ as shown below.) The total mass of the ISM, $`M_{\mathrm{gas}}`$, and the mass of each gas cloud, $`m_\mathrm{c}`$, are parameters. If dust opacity as well as Thomson scattering is considered, the mass extinction is expressed by $`\chi =(n_\mathrm{e}\sigma _\mathrm{T}+n_\mathrm{d}\sigma _\mathrm{d})/(\rho _\mathrm{g}+\rho _\mathrm{d})`$, where $`\sigma _\mathrm{T}`$ is the Thomson scattering cross section, $`n_\mathrm{e}`$ is the electron number density, $`\rho _\mathrm{g}`$ is the gas density, and $`n_\mathrm{d}`$, $`\sigma _\mathrm{d}`$, and $`\rho _\mathrm{d}`$ are respectively the number density, cross-section, and mass density of dust grains. If we take the dust-to-gas mass ratio $`f_{\mathrm{dg}}`$ of $`10^2`$, then the opacity ratio is

$$\frac{n_\mathrm{d}\sigma _\mathrm{d}}{n_\mathrm{e}\sigma _\mathrm{T}}=1.9\times 10^2\left(\frac{a_\mathrm{d}}{1\mu \mathrm{m}}\right)^1\left(\frac{\rho _\mathrm{s}}{\mathrm{g}\mathrm{cm}^3}\right)^1\left(\frac{f_{\mathrm{dg}}}{10^2}\right),$$

where $`a_\mathrm{d}`$ is the grain radius and $`\rho _\mathrm{s}`$ is the density of solid material within the grain. Thus, we find $`n_\mathrm{d}\sigma _\mathrm{d}n_\mathrm{e}\sigma _\mathrm{T}`$ in the situations of interest. Hence, in this paper, we evaluate the mass extinction by $`\chi =n_\mathrm{d}\sigma _\mathrm{d}/\rho _\mathrm{g}`$.

Finally, as for the stellar evolution, we assume a Salpeter-type initial mass function (IMF) as $`\varphi =A(m_{}/M_{})^{1.35}`$ for a mass range of $`[m_l,m_u]`$. In the present analysis, we consider an initial starburst and the subsequent passive evolution of stars. The upper mass limit is inferred to be around $`40M_{}`$ in starburst regions (Doyon, Puxley & Joseph 1992). As for the lower mass limit, some authors suggest that the IMF in starburst regions is deficient in low-mass stars, with a cut-off of about $`23M_{}`$ (Doane & Mathews 1993; Charlot et al. 1993; Hill et al. 1994). In this paper, we assume $`m_l=2M_{}`$ and $`m_u=40M_{}`$ in the early stage of bulge formation. However, it should be also kept in mind that the lower mass limit is under debate; the claimed lower cut-off could be due to the magnitude limit effect (Selman et al. 1999), and also recently sub-solar-mass stars are found in a starburst region in our galaxy, NGC 3603 (Brandl et al. 1999). In order to incorporate the stellar evolution, we adopt the mass-luminosity relation $`(\mathrm{}_{}/L_{})=(m_{}/M_{})^{3.7}`$, and the mass-age relation $`\tau _{}=1.1\times 10^{10}(m_{}/M_{})^{2.7}\mathrm{yr}`$ (Lang 1974), where $`m_{}`$, $`\mathrm{}_{}`$, and $`\tau _{}`$ are respectively the stellar mass, luminosity, and age.

## 3 Basic equations

As a relativistic result of radiative absorption and subsequent re-emission, the radiation fields exert a drag force on moving material in resistance to its velocity. This radiation drag extracts angular momentum from the ISM, thereby allowing the gas to accrete on to the galactic center. The radiation drag is an effect of $`O(v/c)`$, but it could provide a key mechanism for angular momentum transfer in intense radiation fields. We put the origin at the center of the bulge, and adopt cylindrical coordinate $`r`$, $`\varphi `$, and $`z`$, where $`z`$-axis is the rotation axis of stars and gas. The components of the specific radiation force that is exerted on moving fluid elements with velocity $`𝐯`$ are given by

$$f_\mathrm{r}=\frac{\chi }{c}(F^\mathrm{r}v_\mathrm{r}Ev_\mathrm{r}P^{\mathrm{rr}}v_\varphi P^{\varphi \mathrm{r}}v_\mathrm{z}P^{\mathrm{zr}}),$$
(1)
$$f_\varphi =\frac{\chi }{c}(F^\varphi v_\varphi Ev_\varphi P^{\varphi \varphi }v_\mathrm{r}P^{\mathrm{r}\varphi }v_\mathrm{z}P^{\mathrm{z}\varphi }),$$
(2)
and

$$f_\mathrm{z}=\frac{\chi }{c}(F^\mathrm{z}v_\mathrm{z}Ev_\mathrm{z}P^{\mathrm{zz}}v_\mathrm{r}P^{\mathrm{rz}}v_\varphi P^{\varphi \mathrm{z}})$$
(3)
(Mihalas $`\&`$ Mihalas 1984) in the $`r`$-, $`\varphi `$-, and $`z`$\- directions respectively. Here, $`E`$ is the radiation energy density, $`F^\alpha `$ is the radiation flux, and $`P^{\alpha \beta }`$ is the radiation stress tensor where the non-diagonal components are null owing to the present symmetry.

Using equations (1)$``$(3), we have the radiation hydrodynamical equations to $`O(v/c)`$ as

$$\frac{dv_\mathrm{r}}{dt}=\frac{v_\varphi ^2}{r}f_\mathrm{g}^\mathrm{r}+\frac{\chi }{c}[F^\mathrm{r}(E+P^{\mathrm{rr}})v_\mathrm{r}]$$
(4)
$$\frac{1}{r}\frac{d(rv_\varphi )}{dt}=\frac{\chi }{c}[F^\varphi (E+P^{\varphi \varphi })v_\varphi ]$$
(5)
$$\frac{dv_\mathrm{z}}{dt}=f_\mathrm{g}^\mathrm{z}+\frac{\chi }{c}[F^\mathrm{z}(E+P^{\mathrm{zz}})v_\mathrm{z}]$$
(6)
where $`f_\mathrm{g}^\mathrm{r}`$ and $`f_\mathrm{g}^\mathrm{z}`$ are respectively $`r`$\- and $`z`$\- component of the gravitational force.

The azimuthal equation of motion (5) is the equation of angular momentum transfer. This equation implies that the radiation flux force (the first term on the right-hand side) makes fluid elements tend to corotate with stars, whereas the radiation drag (the second term on the right-side) works to extract the angular momentum from gas. Therefore, the gain and loss of total angular momentum is determined by equation (5).

## 4 Uniform ISM

In this section, we consider the angular momentum transfer by radiation drag in a uniform ISM. First, we analytically calculate the radiation fields produced by spherically distributed stars in an optically thin regime, and assess the angular momentum loss rate $`\dot{J}`$ and the mass accretion rate $`\dot{M}`$. Next, we extend the analysis to an optically thick regime. Then, the relationship between the optical depth of a dusty ISM and the angular momentum transfer efficiency is derived.

### 4.1 Optically thin regime

For uniform distributions of stars, the radiation fields inside the bulge are analytically integrated in an optically thin regime to be

$`cE`$ $`=`$ $`{\displaystyle \frac{3L_{\mathrm{bulge}}}{2\pi r_\mathrm{b}^2}}\left(1{\displaystyle \frac{\pi }{4}}{\displaystyle \frac{r}{r_\mathrm{b}}}\right),`$ (7)
$`F^\mathrm{r}`$ $`=`$ $`{\displaystyle \frac{L_{\mathrm{bulge}}}{4\pi r_\mathrm{b}^3}}r,`$ (8)
$`cF^\varphi `$ $`=`$ $`{\displaystyle \frac{3L_{\mathrm{bulge}}}{2\pi r_\mathrm{b}^2}}\omega r_\mathrm{b}\left({\displaystyle \frac{215}{216}}{\displaystyle \frac{49\pi }{192}}\right)\left({\displaystyle \frac{r}{r_b}}\right)^2`$ (9)
$`+{\displaystyle \frac{3L_{\mathrm{bulge}}}{2\pi r_\mathrm{b}^2}}\omega r_\mathrm{b}{\displaystyle \frac{8}{9}}{\displaystyle \frac{r}{r_\mathrm{b}}}\left(1{\displaystyle \frac{r}{r_\mathrm{b}}}\right),`$
$`F^\mathrm{z}`$ $`=`$ $`{\displaystyle \frac{L_{\mathrm{bulge}}}{4\pi r_\mathrm{b}^3}}z,`$ (10)
$`cP^{\mathrm{rr}}`$ $`=`$ $`{\displaystyle \frac{1}{3}}{\displaystyle \frac{3L_{\mathrm{bulge}}}{2\pi r_\mathrm{b}^2}}\left(1{\displaystyle \frac{3\pi }{16}}{\displaystyle \frac{r}{r_\mathrm{b}}}\right),`$ (11)
$`cP^{\varphi \varphi }`$ $`=`$ $`{\displaystyle \frac{1}{3}}{\displaystyle \frac{3L_{\mathrm{bulge}}}{2\pi r_\mathrm{b}^2}}\left(1{\displaystyle \frac{9\pi }{32}}{\displaystyle \frac{r}{r_\mathrm{b}}}\right),`$ (12)
$`cP^{\mathrm{zz}}`$ $`=`$ $`cP^{\varphi \varphi },`$ (13)
where $`L_{\mathrm{bulge}}`$ and $`r_\mathrm{b}`$ are the luminosity and radius of the bulge and $`\omega `$ is the angular velocity of stars. These quantities are equivalent to those obtained by Fukue, Umemura & Mineshige (1997), except that the radiation flux $`F^\varphi `$ in the azimuthal directions is generated by the rotation of the bulge. As seen in equation (5), the flux $`F^\varphi `$ works to corotate the ISM with stars in contrast to the drag force. In Fig. 1, we compare both forces exerted on gas per unit mass, $`f^{\mathrm{drag}}=\chi (E+P^{\varphi \varphi })v_\varphi /c`$ and $`f^\varphi =\chi F^\varphi /c`$, where $`v_\varphi =r\omega `$. It is found that $`f^{\mathrm{drag}}`$ overwhelms $`f^\varphi `$ everywhere. Thus, the optically thin ISM can always lose the angular momentum as a result of radiation drag. The angular momentum loss rate per unit volume per unit time $`\dot{j}`$ is evaluated by

$$\dot{j}=\frac{L_{\mathrm{bulge}}}{c^2}\frac{3\chi \omega }{2\pi }xf(x)$$
(14)
at the point of $`x=r/r_\mathrm{b}`$, where $`f(x)=(\frac{23}{216}+\frac{17\pi }{192})x^2+\frac{4}{9}x`$. Equation (14) is integrated over the volume of the bulge to give the total angular momentum loss rate $`\dot{J}`$ as

$$\dot{J}=\beta \frac{L_{\mathrm{bulge}}}{c^2}\frac{3\chi \omega }{2\pi }r_\mathrm{b}^3\rho _0,$$
(15)
where $`\rho _0`$ is the density of the ISM, and $`\beta =xf(x)𝑑V`$ with $`dV`$ being volume element in spherical coordinates. Noting that the initial angular momentum is

$$J_0=\frac{2}{5}M_{\mathrm{gas}}\omega r_\mathrm{b}^2,$$
(16)
the mass accretion rate is expressed as

$$\dot{M}_{\mathrm{gas}}=M_{\mathrm{gas}}\frac{\dot{J}}{J_0}=\eta \frac{L_{\mathrm{bulge}}}{c^2}\tau _\mathrm{T},$$
(17)
where $`\tau _\mathrm{T}=\chi \rho _0r_\mathrm{b}=(3\chi M_{\mathrm{gas}})/(4\pi r_\mathrm{b}^2)`$. The coefficient $`\eta `$ gives the radiation drag efficiency. By calculating the constant $`\beta `$ numerically, $`\eta `$ is found to be 0.32.

### 4.2 Optically thick regime

In optically thick ISM, the regions where the optical depth $`\tau _\mathrm{s}`$ from each star is less than unity are subject to the radiation drag (Tsuribe $`\&`$ Umemura 1997). In this situation, the flux for the region of $`\tau _\mathrm{s}1`$ is given by

$$F^\varphi =ER\omega 𝑑V,$$
(18)
and

$$F^{\mathrm{drag}}=(E+P^{\varphi \varphi })R_\mathrm{g}\omega 𝑑VER_\mathrm{g}\omega 𝑑V,$$
(19)
where $`R`$ are $`R_\mathrm{g}`$ are respectively the distance from the rotational axis to a star and to one volume element, $`dV=r^2\mathrm{sin}\theta drd\theta d\varphi `$ in spherical coordinates for the position of a star. $`E=\mathrm{}_{}/4\pi cr^2,P^{\varphi \varphi }0`$ with $`\mathrm{}_{}`$ being the luminosity of a star. Consequently, the local angular momentum transfer rate by equation (5) is

$$\dot{j}=\frac{\mathrm{}_{}}{c^2}\frac{\chi \omega }{6}r_\mathrm{s}^3\rho _0.$$
(20)
This is summed up over all stars to give

$$\dot{J}=\frac{1}{6}\frac{L_{\mathrm{bulge}}}{c^2}\chi \omega r_\mathrm{s}^3\rho _0,$$
(21)
where $`r_\mathrm{s}`$ is the size of the region of $`\tau _s=1`$, and $`L_{\mathrm{bulge}}=N_{}\mathrm{}_{}`$. Then, the mass accretion rate is given as

$$\dot{M}_{\mathrm{gas}}=M_{\mathrm{gas}}\frac{\dot{J}}{J_0}=\frac{5}{12}\frac{L_{\mathrm{bulge}}}{c^2}\tau _\mathrm{T}^2.$$
(22)
So far, we have considered a rigidly rotating system. We can also analytically estimate the accretion rate in a system with a different rotation law as $`v_\varphi r^n`$ (e.g. $`n=0.5`$ for the Keplarian rotation). Then, the mass accretion rate turns out to be

$$\dot{M}_{\mathrm{gas}}=M_{\mathrm{gas}}\frac{\dot{J}}{J_0}=\frac{n^2(n+4)}{9B}\left(\frac{\overline{R}_\mathrm{i}}{r_\mathrm{b}}\right)^{n1}\frac{L_{\mathrm{bulge}}}{c^2}\tau _\mathrm{T}^2,$$
(23)
where $`B=\sqrt{\pi }\mathrm{\Gamma }(\frac{n+3}{2})/\mathrm{\Gamma }(\frac{n+4}{2})`$, $`\overline{R}_\mathrm{i}`$ is the mean distance from the rotational axis to each star, and $`\mathrm{\Gamma }(n)`$ is the Gamma function. It is noted that the difference between (22) and (23) is merely a small factor for ordinary rotation laws. As a result, it is found that the mass accretion rate decreases according to $`\tau _\mathrm{T}^2`$ in an optically thick regime. The physical reason is that the larger is the optical depth of the bulge, the closer to zero is the difference between the velocity of stars and the gas velocity, so that the radiation drag efficiency falls.

Combined with the optically thin case, we have the mass accretion rate in the uniform case as

$$\dot{M}_{\mathrm{gas}}\{\begin{array}{cc}\frac{L_{\mathrm{bulge}}}{c^2}\tau _\mathrm{T}\hfill & (\tau _\mathrm{T}<1),\hfill \\ \frac{L_{\mathrm{bulge}}}{c^2}\tau _{\mathrm{T}}^{}{}_{}{}^{2}\hfill & (\tau _\mathrm{T}>1).\hfill \end{array}$$
(24)
As readily understood by these results, the angular momentum transfer efficiency by the radiation drag is maximum when the optical depth of the ISM is around unity.

## 5 Clumpy ISM

In the previous section, we considered the angular momentum transfer in a uniform ISM. In this section, we consider a clumpy ISM. First, we describe the treatment for the extinction by clumpy gas clouds. Next, by solving the radiative transfer numerically, we assess the total angular momentum loss rate $`\dot{J}`$ and the mass accretion rate $`\dot{M}_{\mathrm{gas}}`$. Then we elucidate the relation between the clumpiness of ISM and the angular momentum transfer efficiency quantitatively.

### 5.1 Extinction by a clumpy ISM

In the case of a clumpy ISM, $`N_\mathrm{c}(=10^4)`$ identical clouds with optical depth of $`\overline{\tau }=\chi \rho _{\mathrm{gas}}r_\mathrm{c}`$ are distributed randomly. In an optically thin regime, the radiation fields produced by a star are

$$dE_0=\frac{1}{c}\frac{\mathrm{}_{}}{4\pi r^2},\mathrm{𝐝𝐅}_0=\frac{\mathrm{}_{}}{4\pi r^3}𝐫,dP_0^{\mathrm{rr}}=dE_0,dP_0^{\varphi \varphi }=dP_0^{\mathrm{zz}}0,$$
(25)
where $`𝐫=𝐫_i𝐫_j`$ and $`r=|𝐫|`$ with $`𝐫_i`$ being the position of the $`i`$-th cloud and $`𝐫_j`$ being the position of $`j`$-th star. The physical quantities with suffix $`0`$ are the radiation fields without extinction. The corotating flux in the azimuthal direction and the flux contributing to the radiation drag are respectively

$$dF_0^{\mathrm{rot}}=(dE_0+dP_0^{\varphi \varphi })V_{},$$
(26)
$$dF_0^{\mathrm{drag}}=(dE_0+dP_0^{\varphi \varphi })v_{\mathrm{gas}},$$
(27)
where $`V_{}=r_{}\omega `$ and $`v_{\mathrm{gas}}=r_{\mathrm{gas}}\omega `$ are the rotational velocities of a star and a gas cloud, with $`r_{}`$ and $`r_{\mathrm{gas}}`$ being respectively the distances from the rotational axis to a star and to a gas cloud. Next, we consider the extinction by dust in clumpy gas clouds. We calculate the radiation fields by the direct integration of the radiation transfer (see Fig. 2), where the extinction by all intervening clouds from a star to a cloud is summed up. Then, the radiation flux $`𝐅`$, the radiation energy density $`E`$, and the radiation stress tensor $`P`$ are respectively given by

$`𝐅`$ $`=`$ $`{\displaystyle \underset{j=1}{\overset{N_{}}{}}}d𝐅_{0,j}\mathrm{exp}(\tau _j),`$ (28)
$`E`$ $`=`$ $`{\displaystyle \underset{j=1}{\overset{N_{}}{}}}dE_{0,j}\mathrm{exp}(\tau _j),`$ (29)
$`P`$ $`=`$ $`{\displaystyle \underset{j=1}{\overset{N_{}}{}}}dP_{0,j}\mathrm{exp}(\tau _j),`$ (30)
where $`\tau _j=2(1(b/r_\mathrm{c})^2)^{1/2}\overline{\tau }`$ is the optical depth for all intervening clouds along the light ray.

### 5.2 Angular momentum transfer in a clumpy ISM

By substituting (26)-(30) for (5), we can evaluate the total angular momentum loss rate as

$$\dot{J}=\frac{\chi }{c}\underset{\mathrm{i}=1}{\overset{N_\mathrm{c}}{}}r_\mathrm{i}(F_\mathrm{i}^{\mathrm{rot}}F_\mathrm{i}^{\mathrm{drag}}),$$
(31)
where $`F_i^\alpha =_{\mathrm{j}=1}^N_{}dF_{0,j}^\alpha \mathrm{exp}(\tau _j)`$ with $`\alpha `$ being ’rot’ or ‘drag’, and $`r_i`$ is the distance from the rotational axis to a gas cloud. If an optically thick cloud is irradiated by the radiation, only the optically thin surface layer is subject to radiation drag (Tsuribe & Umemura 1997). Then the surface layer is stripped by the drag and simultaneously loses angular momentum. We assume here that the stripped gas falls on to a central massive object. Then, with the above angular momentum loss rate, we estimate the total mass of the dusty ISM assembled on to a central massive object which may evolve into a supermassive BH. By the equation (31) and the relation $`\dot{M}_{\mathrm{gas}}/M_{\mathrm{gas}}=\dot{J}/J_0`$, the BH mass $`M_{\mathrm{BH}}`$ is assessed as

$$M_{\mathrm{BH}}=_0^{t(m_l)}\dot{M}_{\mathrm{gas}}𝑑t=_0^{t(m_l)}M_{\mathrm{gas}}\frac{\dot{J}}{J_0}𝑑t,$$
(32)
where $`t(m_l)`$ is the age of a star with mass $`m_l`$.

Here, we introduce $`\overline{N}_{\mathrm{int}}`$ as a measure of the clumpiness of ISM. The $`\overline{N}_{\mathrm{int}}`$ is defined by

$$\overline{N}_{\mathrm{int}}=n_\mathrm{c}\pi r_\mathrm{c}^2r_\mathrm{b}=\frac{3}{4}N_\mathrm{c}\left(\frac{r_\mathrm{c}}{r_\mathrm{b}}\right)^2,$$
(33)
where $`n_\mathrm{c}=N_\mathrm{c}/\frac{4}{3}\pi r_\mathrm{b}^3`$ is the number density of gas clouds. $`\overline{N}_{\mathrm{int}}`$ means the average number of gas clouds that are intersected by a light ray over a bulge radius. The degree of clumpiness is larger for smaller $`\overline{N}_{\mathrm{int}}`$, and $`\overline{N}_{\mathrm{int}}=\mathrm{}`$ corresponds to a uniform distribution. When $`\overline{N}_{\mathrm{int}}`$ and the optical depth $`\overline{\tau }`$ of each cloud are specified, the total optical depth $`\tau _\mathrm{T}`$ is given by $`\tau _\mathrm{T}=\overline{N}_{\mathrm{int}}\overline{\tau }`$. With changing $`\overline{N}_{\mathrm{int}}`$, we investigate the relationship between the clumpiness of ISM and the radiation drag efficiency.

In Fig. 3, for different $`\overline{N}_{\mathrm{int}}`$, we show the ratio of the resultant BH mass to the bulge mass ($`M_{\mathrm{BH}}/M_{\mathrm{bulge}}`$) against the total optical depth $`\tau _\mathrm{T}`$. In this figure, the dashed line is the analytic solution for a uniform ISM corresponding to $`\overline{N}_{\mathrm{int}}=\mathrm{}`$, i.e. equation (24). The solid lines show the numerical results for a clumpy ISM, where the arrows represent $`\overline{\tau }=1`$ for different $`\overline{N}_{\mathrm{int}}`$. The thick lines show the results for $`\overline{N}_{\mathrm{int}}1`$ and the thin lines are those for $`\overline{N}_{\mathrm{int}}<1`$. In Fig. 3, several points are clear regarding the effects of the clumpiness on the radiation drag efficiency. First, if $`\overline{N}_{\mathrm{int}}<1`$, the drag efficiency is saturated when $`\overline{\tau }>1`$. Since the ISM is highly clumpy in this case, the radiation from distant sources which have different velocities from the absorbing clouds can contribute to the radiation drag even if $`\tau _\mathrm{T}1`$. However, the covering factor of clouds is smaller than unity, and therefore a large fraction of photons escape from the bulge without intersecting clouds. The saturation is attributed to these two effects. Secondly, if $`\overline{N}_{\mathrm{int}}O(1)`$, $`M_{\mathrm{BH}}/M_{\mathrm{bulge}}`$ grows with $`\tau _\mathrm{T}`$ and maintains a high level in an optically thick ISM. Obviously, the behavior in the thick limit is different from the uniform case, although the situation is the same as the uniform case in that almost all of photons are consumed inside the bulge. This is again because the clumpiness increases the mean free paths of photons, so that the radiation from distant sources enhances the drag efficiency. Thirdly, when $`\overline{N}_{\mathrm{int}}`$ is much larger than unity, the distribution of ISM is closer to a uniform distributions and therefore the difference between the velocity of a star and the velocity of an absorbing cloud is closer to zero, so that the drag efficiency falls as the uniform case. To conclude, the combined effects of the high clumpiness of the ISM and the large covering factor of clouds are necessary for high drag efficiency.

## 6 Black hole-to-bulge mass ratio

The linear relation between the BH mass and the bulge mass is a direct consequence of the present radiation hydrodynamical mechanism. The possible mass accreted by the radiation drag is given by $`M_{\mathrm{max}}=_0^{t(m_l)}L_{\mathrm{bulge}}(t)/c^2𝑑t`$. Incorporating stellar evolution based on the present single starburst model, we find $`M_{\mathrm{max}}=1.2\times 10^3M_{\mathrm{bulge}}`$. As shown above, the drag efficiency is sensitively dependent on the density distributions of the ISM. In a highly clumpy ISM where $`\overline{N}_{\mathrm{int}}`$ is smaller than a few, the BH-to-bulge mass ratio $`M_{\mathrm{BH}}/M_{\mathrm{bulge}}`$ can be expressed by

$`{\displaystyle \frac{M_{\mathrm{BH}}}{M_{\mathrm{bulge}}}}`$ $`=`$ $`\eta {\displaystyle \frac{M_{\mathrm{max}}}{M_{\mathrm{bulge}}}}(1\mathrm{e}^{\tau _T})`$ (34)
$`=`$ $`1.2\times 10^3\eta (1\mathrm{e}^{\tau _T}),`$ (35)
where $`\eta `$ gives the radiation drag efficiency. From the above analysis, $`\eta `$ is found to be maximally 0.34 in the optically thick limit, and then $`M_{\mathrm{BH}}/M_{\mathrm{bulge}}=4.1\times 10^4`$.

The present model may be too simple in some respects to compare the results closely with the observations. For a more precise prediction of $`M_{\mathrm{BH}}/M_{\mathrm{bulge}}`$, it seems necessary to consider further effects which have not been incorporated in this simple model. In the present analysis, we assumed an initial coeval event of star formation, and the subsequent passive evolution without further star formation episodes. The radiation drag efficiency is basically determined by the total number of photons which are emitted from sources and absorbed by clouds during the whole history of the bulge. The recycling of the ISM for star formation generates more photons and therefore could enhance the mass ratio roughly by a factor of 2 (Umemura 2001). Also, in realistic situations, the radiation drag is not likely to remove completely the angular momentum of stripped gas, and also stripped gas may be mixed with ISM having appreciable angular momentum. Although the detailed processes are not very clear at present, a little leftover angular momentum may lead to the formation of a viscous accretion disk around a BH, which could ignite QSO activity with nearly Eddington luminosity (Umemura 2001). If QSO activity is triggered, the radiation from the QSO can also exert the drag force. This effect can enhance the mass ratio by maximally a factor of 1.7 (Umemura 2001). Furthermore, if a galactic wind occurs, it might cause a transition from an optically thick starburst phase to an optically thin QSO phase, and simultaneously reduce the drag efficiency.

If the above effects are taken into account, the mass ratio could be maximally $`M_{\mathrm{BH}}/M_{\mathrm{bulge}}=1.4\times 10^3`$. This is compared with the recent observational data in Fig. 4, combined with the result by a single-starburst model, $`M_{\mathrm{BH}}/M_{\mathrm{bulge}}=4.1\times 10^4`$. In Fig. 4, the hatched area is the prediction of the present analysis. As seen in this figure, the present radiation hydrodynamical model can roughly account for the observational data, although the prediction falls slightly short of the observed median mass ratio. The prediction is smaller by a factor of six than $`M_{\mathrm{BH}}/M_{\mathrm{bulge}}=0.006`$ by Magorrian et al. (1998), while it is comparable to $`M_{\mathrm{BH}}/M_{\mathrm{bulge}}=0.001`$ by Merritt & Ferrarese (2001a).

Another important effect is the geometrical dilution. In previous works, the radiation drag efficiency would be strongly subject to the effects of geometry (Umemura, Fukue, & Mineshige 1997, 1998; Ohsuga et al. 1999). A large fraction of emitted photons can escape from a disc-like system and thus the radiation drag efficiency is considerably reduced. Recent observation have revealed that the BH mass fraction is significantly smaller than 0.001 in disc galaxies (Salucci et al. 2000; Sarzi et al. 2001). Geometrical dilution may be a reason for the observational fact in disc galaxies, but the quantitative details are not clear before such an aspherical system is actually simulated.

## 7 Conclusions

By assuming a simple model of a bulge, we have investigated the mutual effect between the clumpiness of interstellar medium and the optical depth on the radiation drag efficiency for the angular momentum transfer. In a clumpy interstellar medium, we have accurately solved 3D radiation transfer to calculate the radiation drag force by the rotating bulge stars. We find that the radiation drag efficiency is sensitively dependent on the density distribution of the ISM in an optically thick regime. The efficiency drops by a factor of $`\tau _\mathrm{T}^2`$ in a uniform ISM, while the efficiency turns out to be almost constant at a high level in a clumpy ISM. Also, the radiation drag efficiency falls as the covering factor of clouds becomes smaller than unity. Hence, for the radiation drag to work effectively, it is necessary that the covering factor is larger than unity and that the distributions of ISM is highly clumpy. The present radiation hydrodynamical mechanism accounts for the linear relation between the mass of a supermassive black hole and the mass of a galactic bulge. The range of the predicted mass ratio is from $`M_{\mathrm{BH}}/M_{\mathrm{bulge}}=4.1\times 10^4`$ to $`M_{\mathrm{BH}}/M_{\mathrm{bulge}}=1.4\times 10^3`$. Although we have adopted a highly simplified model, the result roughly accounts for the recent observational trends. Thus, it seems that more elaborate calculations deserve to be performed in the present context. In the future, we should consider further several effects which have not been incorporated in the present analysis, and which can influence the mass ratio significantly: e.g. the realistic chemical evolution and geometrical dilution.

## Acknowledgments

We thank M. Mori, T. Nakamoto and H. Susa for helpful discussion. We are grateful to M. Hayashi and K. Ohsuga for many useful comments. We also thank the anonymous referee for valuable commments. Numerical simulations were performed with facilities at Center of Computational Physics, University of Tsukuba. This work was supported in part by the Grant-in-Aid of the JSPS, 11640225.


text:
# DESY 01-182 hep-ph/0111042 Collider versus Cosmic Ray Sensitivity to Black Hole Production

## Acknowledgements

We thank A. Brandenburg, S. Huber, M. Kowalski, K. Mannheim, F. Schrempp, E. Waxman, J. Wells, and E. Zas for fruitful discussions.