text:
# Intense field stabilization in circular polarization: 3D time-dependent dynamics

\[

## Abstract

We investigate the stabilization of a hydrogen atom in circularly polarized laser fields. We use a time-dependent, fully three dimensional approach to study the quantum dynamics of the hydrogen atom subject to high intensity, short wavelength laser pulses. We find enhanced survival probability as the field is increased under fixed envelope conditions. We also confirm wavepacket dynamics seen in prior time-dependent computations restricted to two dimensions.

\]

The advent of lasers producing electric fields at or above inter-atomic electric fields has led to the discovery of many new and highly nonlinear phenomena. As nonlinear laser-atom physics has matured, computational approaches employed by researchers in the field have naturally been applied to scenarios not easily realized in experiments. One such scenario is the interaction of an intense, *high frequency* laser with the single-electron hydrogen atom. Commercially available pulsed laser systems readily produce intensities at or above the atomic unit $`I=3.51\times 10^{16}W/cm^2`$, but available photon energies $`\mathrm{}\omega `$ are well below the atomic unit frequency, which is given by the ground state ionization energy $`|E_o|`$. In the case where the driving frequency is near or above the ground state ionization energy, $`\mathrm{}\omega |E_o|`$, a remarkable phenomena known as “stabilization” may occur. Stabilization is characterized by a decrease in the ionization probability as the laser intensity increases. For constant driving field, stabilization is manifest as an increase in dressed state *lifetime*. However, due to the initial increase in ionization rate for increased fields as the interaction turns on, increased steady-state lifetimes do not necessarily translate into increased end-of-pulse survival probabilities. This poses a significant challenge to experimental observation of stabilization and is commonly known as the “death-valley” problem.

Nevertheless, indications of stabilization have been seen by tailoring the initial state to act as an effective ground state. In particular, the high frequency condition has been realized using dipole forbidden circular Rydberg states ($`n=4,5`$) of neon ($`|E_o|<1eV`$), and a driving laser of $`\mathrm{}\omega =2eV`$ . In this experiment, a decrease in total ionization yield was observed with increasing peak laser intensity, while fluence was held constant. Fermi’s golden rule predicts the total ionization depends only on the fluence, so this result clearly indicates non-perturbative stabilization . With this notable exception, however, stabilization remains experimentally unconfirmed, and detailed studies have been confined to the realm of simulation.

The vast majority of these studies have concentrated on the case of linearly polarized (LP) fields, where much has been learned about the dynamics of stabilization. Studies considering the case of circularly polarized (CP) fields have also noted comparable or enhanced stabilization . However, the physical mechanism for stabilization in the CP case remains open, as the electron dynamics are dramatically different from the dynamics in LP fields. Recent time-dependent investigations of the CP Rydberg system, restricted to the two-dimensional plane of polarization, have found “ring-like” probability distributions. This ring structure is generally found to rotate in phase with the laser field, *about a point displaced from core*. Thus the motion resembles the motion of a hula-hoop. A probability maximum may also occur at the center of the rotating ring. All prior fully three dimensional studies have been limited to the case of LP systems .

Computational Model – In this Letter, we investigate the quantum dynamics of stabilization in the CP hydrogen system using a three dimensional time-dependent approach. In this system, the electron wavefunction moves under the influence of the Coulomb force and circularly polarized radiation field. For the fields and frequencies considered herein, the interaction is appropriately described in the dipole approximation with classical radiation field. Adopting atomic units ($`\mathrm{}=1=e=m`$), the Hamiltonian for the electron may be written

$$H=\frac{1}{2}^2\frac{1}{r}+F(t)(\widehat{x}\mathrm{cos}(\omega t)+\widehat{y}\mathrm{sin}(\omega t))$$

where $`\widehat{x}`$ and $`\widehat{y}`$ are position operators, and $`F(t)`$ describes the laser field envelope. The Schrödinger equation is integrated using Crank-Nicholson finite differencing and absorbing boundary conditions. In order to avoid the singularity at the origin, we offset the position of points on the grid by half the distance between adjacent gridpoints. This places gridpoints symmetrically about the origin, and poses no problem provided the spatial discretization level is sufficient to resolve the maximal components of momentum $`k_{max}`$ attained by the wavefunction. This can be estimated with confidence from the high harmonic generation (HHG) cutoff law since harmonics produced under stabilization conditions are known to be suppressed . The HHG cutoff energy is given by $`\mathrm{}\omega |E_o|+3U_p`$, where $`U_p=F^2/4\omega ^2`$ is the pondermotive energy of an electron in a LP field. Thus $`k_{max}`$ is essentially proportional to $`F/\omega `$. Setting the laser frequency $`\omega =1.2`$, and considering fields up to $`F4.0`$, we find a spatial discretization of $`\delta x_i0.1667`$ is sufficient to resolve all components of the wavefunction. To our knowledge this is the largest field considered to date in studies of stabilization for circular polarization.

Survival Probability – To determine the survival probability, we start with a stationary gaussian wavepacket centered at the origin with width $`\sigma =1.5`$, and introduce the laser interaction with a pulse envelope $`F(t)=F\mathrm{sin}^2(\frac{\pi t}{6\tau })`$, over the time interval $`0<t<6\tau `$, where $`\tau `$ is the optical cycle. A portion of the wavefunction ionizes off the grid during the laser turn-on, and at the end of the pulse, remaining probability essentially occupies the atomic ground state. As the field strength is increased, we find the total probability remaining on the grid at the end of the pulse increases, characteristic of stabilization. This is demonstrated in Fig. 1, where we show the time dependent survival probability on the grid over the 6-cycle pulse, for the field strengths $`F=2.5`$, $`F=3.0`$, $`F=3.5`$, and $`F=4.0`$. The probability on the grid remains close to unity during the first half of the pulse for all field strengths. This corresponds to the time it takes components of the wavefunction to traverse the grid and ionize. At a point near two laser cycles the probability begins to drop significantly, with lower total probability occuring for larger field strengths. Therefore, an increase in total ionization occurs during the laser ramp-up for higher fields. However, during the following two cycles when the laser field is closest to its maximum, ionization occurs more slowly, or at a reduced rate, for more intense fields. The total ionization for the lower fields then catches up to and surpasses that of the higher fields. This occurs just beyond four laser cycles for each of the fields we consider, and this trend continues to the end of the laser pulse, resulting in enhanced survival probability as a function of field strength.

Wavepacket Dynamics – While Fig. 1 confirms the presence of the stabilization for the fields and frequency we consider, the fundamental stabilization mechanism remains obscure. In order to address this issue we must examine the time dependent wavepacket dynamics. Now, when a CP laser field of sufficient intensity to induce direct ionization is turned on, the electron probability will generally flow away from the core along a spiral path to eventual ionization. However, for the fields we consider, an interesting twist on this effect is seen. The spiral path begins to close back upon itself. This produces a probability trapping effect, as the ionization channel is closed. This is evidenced by the formation a ring-like probability distribution. In Fig. 2 we show a series of snapshots near the peak of the pulse exhibiting this behavior. The field strength is $`F=3.5`$ and the time is $`t=15.544`$, approximately $`3`$ laser cycles. We find a ring-like structure of enhanced probability. Probability on this ring is asymmetrically distributed with maximum occuring along the direction of the rotating laser field. This peak follows a circular orbit about the core, with approximate radius $`F/\omega ^2`$. The entire ring structure also rotates about the core, with one edge intersecting the core. Thus, the ring is displaced along the field, and rotates, locked in phase with the field, about the origin. This behavior cannot be predicted or observed in time-averaged (Kramers-Henneberger) approaches. Specifically, if one takes a time average of the probability, it may result in an apparent “ring-like” distribution, but this will necessarily be rotationally invariant about the origin, and cannot exhibit the displacement along the laser field evident in time dependent-studies. Thus the actual “hula-hoop” motion will be concealed by the time averaging. The last frame of Fig. 2 also shows the end-of-pulse probability distribution, which has overlap with the hydrogen ground state to within $`10^5`$.

Shifts in Ring Structure – In Fig. 3, we illustrate the displacement and closure of the spiral arm explicitly with a cross sectional slice of the probability distribution. We take a cross section for each field strength, near the peak of the pulse. We choose $`t3\tau `$, so that the laser field is aligned along the $`x`$-axis. Each profile exhibits three maxima: one near the core, one along the laser field at positive $`x`$, and one on the “back-side” at negative $`x`$. These maxima correspond to probability remaining near the core, and the two locations at which the spiral tail intersects the $`x`$-axis. As the field increases, the maxima at positive $`x`$ shifts outward along the field direction, and the “back-side” maxima shifts toward the core. This corresponds to the displacement of the ring along the field, and closing of the tail, which forms the ring-like structure. These displaced ring structures have been noted in time dependent two dimensional studies of stabilization . Another interesting feature of Fig. 3 is the decrease in probability amplitude seen at the spiral tail crossing points, relative to the probability amplitude at the core, as the field increased. This a laser ramping effect, *i.e.* the limited time available for the probability to spread over the ring leads to unequal population distribution on the ring. In this case, as we increase the field, the maximum probability shifts along the ring, from a point located along the laser field, backwards along the spiral arm, toward the core. This is indicative of the fact that the probability distribution on the stabilized ring structure may be adjusted by the pulse ramping.

In conclusion, for a hydrogen atom subject to an intense, high frequency laser pulse of circular polarization, we have found a decrease in total ionization yield as the field strength is increased under fixed pulse envelope conditions. For this stabilization phenomena in circularly polarized fields, the characteristic waveform is a ring-like structure, displaced along the laser electric field. This entire displaced ring structure is locked in phase with the rotating field, executing a “hula hoop” type motion. This behavior can only be observed in time-dependent approaches. Our fully three dimensional investigation confirms prior results obtained using a two-dimensional time dependent approach, and indicates the two dimensional time-dependent approach is likely to provide the necessary physical insight into stabilization while maintaining reasonable computational effort. These ring-like waveforms are typically asymmetrically populated, and population of these stabilized structures may be achieved from the ground state with sufficiently short pulses.

The authors would like to thank Linda Reichl, Joe Eberly, and Peter Knight for useful conversations. We also thank The University of Texas at Austin High Performance Computing Center for use of their facilities.


format: text

text:
# Multiscaling in Ising quantum chains with random Hilhorst–van Leeuwen perturbations

## 1 Introduction

In semi-infinite two-dimensional (2D) Ising models, extended perturbations of the coupling constants of the form $`\mathrm{\Delta }K_l=Al^\kappa `$, decaying as a power of the distance $`l`$ from a free surface, were first studied by Hilhorst and van Leeuwen (HvL) at the begining of the 1980s . Since then they have been the subject of continuous interest . Similar extended perturbations associated with line defects in the bulk  or point defects  have also been considered (see Ref.  for a review).

Such extended perturbations are known to modify the local critical behaviour, for an arbitrarily small value of the perturbation amplitude $`A`$, as soon as the decay exponent $`\kappa `$ is sufficiently small. Analytical results for a surface defect in the 2D Ising model  show that with a temperature-like perturbation (the case considered later) marginal behaviour is obtained when $`\kappa =1`$. Then the surface exponents vary continuously with the perturbation amplitude $`A`$ and the surface transition becomes first order when $`A`$ is greater than a critical value $`A_\mathrm{c}`$. When $`\kappa <1`$, the surface properties display essential singularities when $`A<0`$ and the transition is first order when $`A>0`$.

This rich critical behaviour has been explained independently by Cordery and Burkhardt  through scaling considerations. Under a change of the length scale $`l^{}=l/b`$, the perturbation transforms as

$$[\mathrm{\Delta }K_l]^{}=\frac{A^{}}{l_{}^{}{}_{}{}^{\kappa }}=b^{1/\nu }\frac{A}{l^\kappa },$$
(1)
where $`\nu `$ is the bulk correlation length exponent. Comparing both sides of the last equation, one obtains the following scaling behaviour for the perturbation amplitude:

$$A^{}=b^{\kappa +1/\nu }A.$$
(2)
Thus, for the 2D Ising model, $`\kappa _\mathrm{c}=1/\nu =1`$ corresponds to the critical value of the decay exponent below which the extended perturbation becomes relevant.

An interesting aspect of these smoothly inhomogeneous models is that in 2D their correlation functions are covariant under conformal transformations when the perturbation is marginal. Thus gap-exponent relations are satisfied  and conformal profiles are obtained in the strip geometry , provided the form of the perturbation is also properly transformed as first noticed in Ref. .

Recently, aperiodic  and random  versions of the HvL model have been studied for the quantum Ising chain in a transverse field, which corresponds to the extreme anisotropic limit  of the 2D classical Ising model and belongs to the same universality class.

In Ref.  random extended perturbations decaying towards a pure bulk system were considered with two types of quenched randomness.

In model #1, the random couplings take the following form:

$$J_l=J\left[1+a_1\frac{(1)^{f_l}}{l^\kappa }\right],f_l=\{\begin{array}{c}1\hfill \\ 0\hfill \end{array}\mathrm{with}\mathrm{probability}\begin{array}{c}1/2\hfill \\ 1/2\hfill \end{array}.$$
(3)
Thus the average perturbation $`J_lJ`$ vanishes.

In model #2

$$J_l=J\left[1+f_la_2\right],f_l=\{\begin{array}{c}1\hfill \\ 0\hfill \end{array}\mathrm{with}\mathrm{probability}\begin{array}{c}l^\kappa \hfill \\ 1l^\kappa \hfill \end{array}.$$
(4)
Then the average coupling $`J_l=J(1+a_2l^\kappa )`$ has the same form as in the HvL model.

In both models the transverse fields are kept constant, $`h_l=h`$.

A perturbation expansion of the average free energy and scaling considerations allowed us to obtain a relevance criterion for both types of perturbations . With model #1, the first-order correction to the average free energy vanishes and the second-order correction leeds to the following scaling behaviour for the amplitude $`a_1`$ of a thermal perturbation:

$$a_1^{}=b^{1/\nu \kappa 1/2}a_1.$$
(5)
With model #2, the scaling behaviour of the perturbation amplitude $`a_2`$ is governed by the non-vanishing first-order correction to the average free energy so that

$$a_2^{}=b^{1/\nu \kappa }a_2,$$
(6)
as for the non-random HvL model in (2). According to (5) and (6) marginal behaviour is obtained for $`\kappa =1/2`$ with model #1 and $`\kappa =1`$ with model #2.

Analytical expressions for the finite-size scaling behaviour at criticality of the average and typical surface magnetization were obtained in Ref. , in the case of marginal and relevant perturbations, in agreement with the results of numerical simulations. We gave also some conjectured expressions for the surface energy exponents in the marginal case.

In the present work, we extend our previous study by considering the finite-size scaling behaviour of the average of the moments of the surface magnetization $`m_\mathrm{s}^q(L)`$ and of the singular part of the surface energy $`e_\mathrm{s}^q(L)`$ in the marginal situation. In Section 2, we briefly recall the techniques used to diagonalize the corresponding quantum Hamiltonians and to obtain the scaling behaviour of the surface magnetization and surface energy. Analytical expressions for the average of the moments of the surface magnetization and the surface energy are deduced from scaling considerations in Sections 3. The results are compared to numerical finite-size scaling data in Section 4 and their range of validity discussed in Section 5.

## 2 Free fermions techniques

### 2.1 Diagonal Hamiltonian

We consider the semi-infinite transverse-field Ising model (TIM) with Hamiltonian

$$=\frac{1}{2}\underset{l=1}{\overset{\mathrm{}}{}}(J_l\sigma _l^z\sigma _{l+1}^z+h_l\sigma _l^x),$$
(7)
where $`\sigma _l^x`$ and $`\sigma _l^z`$ are Pauli matrices, $`J_l`$ is the random coupling and $`h_l`$ the transverse field.

The TIM Hamiltonian can be diagonalized using standard methods . A Jordan–Wigner transformation  changes the original Hamiltonian (7) into a quadratic form in fermions which is diagonalized through a canonical transformation leading to

$$=\underset{q=1}{\overset{L}{}}ϵ_q\left(\eta _q^{}\eta _q\frac{1}{2}\right)$$
(8)
for a chain with length $`L`$. The $`\eta _q^{}`$ ($`\eta _q`$) are fermion creation (annihilation) operators and the non-negative excitation energies $`ϵ_q`$ satisfy the set of equations

$`ϵ_q\psi _q(l)`$ $`=`$ $`h_l\varphi _q(l)J_l\varphi _q(l+1),`$
$`ϵ_q\varphi _q(l)`$ $`=`$ $`J_{l1}\psi _q(l1)h_l\psi _q(l)`$ (9)
with the boundary conditions $`J_0=J_L=0`$. Introducing $`2L`$-dimensional vectors $`𝑽_q`$ with components

$$V_q(2l1)=\varphi _q(l),V_q(2l)=\psi _q(l)$$
(10)
the relations (9) lead to an eigenvalue problem for the following tridiagonal matrix:

$$𝖳=\left(\begin{array}{cccccc}0& h_1& & & & \\ h_1& 0& J_1& & & \\ & J_1& 0& h_2& & \\ & & \mathrm{}& \mathrm{}& \mathrm{}& \\ & & & J_{l1}& 0& h_L\\ & & & & h_L& 0\end{array}\right).$$
(11)
When $`𝖳`$ is squared, odd and even components of $`𝑽_q`$ decouple and one obtains two separate eigenvalue problems for $`\mathit{\varphi }`$ and $`𝝍`$. The relations (9) are invariant under the transformation $`\mathit{\varphi }_q\mathit{\varphi }_q`$, $`ϵ_qϵ_q`$ so that, changing $`\mathit{\varphi }_q`$ into $`\mathit{\varphi }_q`$ in $`𝑽_q`$, one obtains the eigenvector of $`𝖳`$ with eigenvalue $`ϵ_q`$. Thus the physically relevant information is contained in that part of the spectrum of $`𝖳`$ with $`ϵ_q0`$.

### 2.2 Surface magnetization and surface energy

The surface critical properties can be obtained through finite-size scaling at criticality, i.e., by working on a finite system with length $`L`$ and $`J=h=1`$.

The imaginary time spin-spin autocorrelation function is given by

$$𝒢_l^{\sigma \sigma }(\tau )=0|\sigma _l^z(\tau )\sigma _l^z(0)|0=\underset{n}{}|n|\sigma _l^z|0|^2\mathrm{exp}[\tau (E_nE_0)],$$
(12)
where $`|0`$ and $`|n`$ are the ground state and the $`n`$th excited state of $``$, $`E_0`$ and $`E_n`$ the corresponding eigenvalues. With fixed boundary conditions, $`\sigma _L^z=\pm 1`$, at $`l=L`$, the ground state is degenerate and asymptotically, one obtains $`lim_\tau \mathrm{}𝒢_l^{\sigma \sigma }(\tau )=m_l^2`$ where

$$m_l=\sigma |\sigma _l^z|0$$
(13)
is the local magnetization. It is given by an off-diagonal matrix element involving the first excited state with one fermionic excitation, $`|\sigma =\eta _1^{}|0`$, which is degenerate with the ground state. For the surface spin at $`l=1`$ a simple expression is obtained :

$$m_\mathrm{s}(L)=m_1=\varphi _1(1)=\left[1+\underset{l=1}{\overset{L1}{}}\underset{k=1}{\overset{l}{}}\left(\frac{h_k}{J_k}\right)^2\right]^{1/2}.$$
(14)
This expression can be rewritten under the following form :

$$m_\mathrm{s}(L)=\overline{m}_\mathrm{s}^\mathrm{d}(L)\underset{l=1}{\overset{L1}{}}\left(\frac{J_l}{h_l}\right),$$
(15)
where

$$\overline{m}_\mathrm{s}^\mathrm{d}(L)=\left[1+\underset{j=1}{\overset{L1}{}}\underset{k=Lj}{\overset{L1}{}}\left(\frac{J_k}{h_k}\right)^2\right]^{1/2}$$
(16)
is the surface magnetization at $`l=L1`$ on the dual chain with transverse fields $`J_l`$, couplings $`h_l`$ $`(l=0,L1)`$ and fixed boundary conditions $`J_0=0`$ at $`l=0`$.

The scaling dimension of the surface energy density can be obtained by considering the finite-size behaviour with free boundary conditions of the off-diagonal matrix element

$$e_\mathrm{s}(L)=\epsilon |\sigma _1^x|0=(ϵ_2ϵ_1)\varphi _1(1)\varphi _2(1),$$
(17)
where $`|\epsilon =\eta _1^{}\eta _2^{}|0`$ is the lowest eigenstate leading to a non-vanishing matrix element. This matrix element enters into the expression of the connected energy-energy surface autocorrelation function $`𝒢_1^{\epsilon \epsilon }(\tau )`$ and scales like the singular part of the surface energy density.

For the critical HvL model with a marginal decay of the perturbation, $`h_l=1`$, $`J_l=1+a/l`$, the product in (15) behaves asymptotically as $`L^a`$. When $`L1`$, the couplings on the right side of the chain are asymptotically unperturbed, thus one expects that $`\overline{m}_\mathrm{s}^\mathrm{d}(L)L^{1/2}`$ as for the homogeneous Ising model at the ordinary transition. This leads to the scaling behaviour:

$$m_\mathrm{s}(L)L^{x_\mathrm{m}^\mathrm{s}},x_\mathrm{m}^\mathrm{s}=\frac{1}{2}a,a\frac{1}{2}.$$
(18)
This expression cannot be valid beyond $`a=1/2`$ where the surface remains ordered at the bulk critical point . The surface magnetization then displays a first order transition and $`x_\mathrm{m}^\mathrm{s}=0`$ when $`a>a_c=1/2`$. Thus our assumption concerning the behaviour of $`\overline{m}_\mathrm{s}^\mathrm{d}(L)`$ is wrong when $`a>1/2`$ and then one must have $`\overline{m}_\mathrm{s}^\mathrm{d}(L)L^a`$.

For the surface energy in Eq. (17), when $`a1/2`$ both $`\varphi _1(1)`$ and $`\varphi _2(1)`$ scale as $`m_\mathrm{s}(L)L^{x_\mathrm{m}^\mathrm{s}}`$. The excitations $`ϵ_1`$ and $`ϵ_2`$ scale as $`L^1`$ leading to

$$e_\mathrm{s}(L)L^1m_\mathrm{s}^2(L)L^{x_\mathrm{e}^\mathrm{s}},x_\mathrm{e}^\mathrm{s}=2(1a),a\frac{1}{2},$$
(19)
in the second-order regime. When $`a>1/2`$, the onset of surface order leads to a more complicated scaling behaviour. The first gap $`ϵ_1`$ is associated with a localized state which reflects the surface behaviour. It scales as $`L^{2a}`$ while $`ϵ_2L^1`$ as usual. We have $`\varphi _1(1)=m_\mathrm{s}(L)L^0`$ but $`\varphi _2(1)`$ behaves anomalously as $`L^{a+1/2}`$. Finally, $`x_\mathrm{e}^\mathrm{s}=a+1/2`$ in the first order regime, $`a>1/2`$.

## 3 Calculation of the moments

In this section we evaluate the moments of the surface magnetization and the surface energy density at criticality on a chain with length $`L`$ in the regime of second-order surface transition. We restrict our study to the case of marginal perturbations, i.e., with $`\kappa =1/2`$ for model #1 and $`\kappa =1`$ for model #2.

### 3.1 Surface magnetization

We assume that the surface magnetization at $`l=L1`$ on the dual chain, $`\overline{m}_\mathrm{s}^\mathrm{d}(L)`$ in (15), is still scaling as $`L^{1/2}`$ for the random models so that

$$m_\mathrm{s}^q(L)=L^{q/2}\underset{l=1}{\overset{L1}{}}J_l^q=L^{q/2}\underset{l=1}{\overset{L1}{}}J_l^q,$$
(20)
where we used the statistical independence of the couplings on different bonds. Using the distribution of the couplings given in (3) for model #1, with $`J=1`$ we obtain

$$J_l^q=\frac{1}{2}\left(1+\frac{a_1}{l^{1/2}}\right)^q+\frac{1}{2}\left(1\frac{a_1}{l^{1/2}}\right)^q1+\frac{q(q1)}{2}\frac{a_1^2}{l}.$$
(21)
The last expression gives the asymptotic behaviour when $`l1`$ which governs the critical behaviour. Taking the logarithm of the product in (20) and a continuum approximation leads to

$$\mathrm{ln}\underset{l=1}{\overset{L1}{}}J_l^q_1^L\mathrm{ln}\left(1+\frac{q(q1)}{2}\frac{a_1^2}{l}\right)l\frac{q(q1)a_1^2}{2}\mathrm{ln}L.$$
(22)
Thus the $`q`$th moment of the surface magnetization scale as

$$m_\mathrm{s}^q(L)L^{qx_\mathrm{m}^\mathrm{s}(q)}L^{(q/2)[1(q1)a_1^2]}$$
(23)
and the surface magnetization displays multiscaling with

$$x_\mathrm{m}^\mathrm{s}(q)=\frac{1}{2}[1(q1)a_1^2](\mathrm{model}\mathrm{\#}1).$$
(24)
Similarly for model #2, using (4) with $`J=1`$ and $`\kappa =1`$, we obtain

$$J_l^q=\frac{1}{l}(1+a_2)^q+1\frac{1}{l}=1+\frac{(1+a_2)^q1}{l}.$$
(25)
Since

$$\mathrm{ln}\underset{l=1}{\overset{L1}{}}J_l^q_1^L\mathrm{ln}\left[1+\frac{(1+a_2)^q1}{l}\right]l[(1+a_2)^q1]\mathrm{ln}L,$$
(26)
the $`q`$th moment of the surface magnetization scales as

$$m_\mathrm{s}^q(L)L^{(q/2)+(1+a_2)^q1},$$
(27)
so that

$$x_\mathrm{m}^\mathrm{s}(q)=\frac{1}{2}\frac{(1+a_2)^q1}{q}(\mathrm{model}\mathrm{\#}2).$$
(28)
The exponents in Eqs. (24) and (28) are in agreement with the average and typical exponents obtained previously in Ref.  when $`q=1`$ and $`q0`$, respectively:

$$\begin{array}{ccc}[x_\mathrm{m}^\mathrm{s}]_{\mathrm{av}}=\frac{1}{2},\hfill & [x_\mathrm{m}^\mathrm{s}]_{\mathrm{typ}}=\frac{1}{2}(1+a_1^2)\hfill & (\mathrm{model}\mathrm{\#}1)\hfill \\ [x_\mathrm{m}^\mathrm{s}]_{\mathrm{av}}=\frac{1}{2}a_2,\hfill & [x_\mathrm{m}^\mathrm{s}]_{\mathrm{typ}}=\frac{1}{2}\mathrm{ln}(1+a_2)\hfill & (\mathrm{model}\mathrm{\#}2)\hfill \end{array}$$
(29)
### 3.2 Surface energy density

For the scaling behaviour of the moments of the surface energy density, one can start from the expression giving $`e_\mathrm{s}(L)`$ in Eq. (19). The factor $`L^1`$ is associated with the scaling of the excitations and reflects the isotropy of the bulk critical behaviour. It should not be modified with a random surface perturbation, hence we have

$$e_\mathrm{s}^q(L)L^qm_\mathrm{s}^{2q}(L).$$
(30)
The $`q`$th moment of the surface energy density depends on the moment of order $`2q`$ of the surface magnetization.

For model #1, making use of Eq. (23), we obtain

$$e_\mathrm{s}^q(L)L^{qx_\mathrm{e}^\mathrm{s}(q)}L^{q[2(2q1)a_1^2]},$$
(31)
so that

$$x_\mathrm{e}^\mathrm{s}(q)=2(2q1)a_1^2(\mathrm{model}\mathrm{\#}1).$$
(32)
In the same way Eq. (27) leads to

$$e_\mathrm{s}^q(L)L^{2q+(1+a_2)^{2q}1}$$
(33)
for model #2 so that

$$x_\mathrm{e}^\mathrm{s}(q)=2\frac{(1+a_2)^{2q}1}{q}(\mathrm{model}\mathrm{\#}2).$$
(34)
The average and typical exponents are then given by:

$$\begin{array}{ccc}[x_\mathrm{e}^\mathrm{s}]_{\mathrm{av}}=2a_1^2,\hfill & [x_\mathrm{e}^\mathrm{s}]_{\mathrm{typ}}=2+a_1^2\hfill & (\mathrm{model}\mathrm{\#}1)\hfill \\ [x_\mathrm{e}^\mathrm{s}]_{\mathrm{av}}=2a_2(2+a_2),\hfill & [x_\mathrm{e}^\mathrm{s}]_{\mathrm{typ}}=2[1\mathrm{ln}(1+a_2)]\hfill & (\mathrm{model}\mathrm{\#}2)\hfill \end{array}$$
(35)
The typical exponents in Eqs. (29) and (35) are in agreement with the values conjectured in Ref. . They keep the same form as for the pure HvL model in Eqs. (18) and (19) with $`a`$ replaced by an effective amplitude, $`a_1^2/2`$ for model #1 and $`\mathrm{ln}(1+a_2)`$ for model #2. These typical effective amplitudes are deduced from the relation $`\mathrm{ln}J_l(a_{\mathrm{eff}})\mathrm{ln}J_l`$ for $`l1`$. A similar conjecture proposed in Ref.  for the average exponents is actually verified only up to the first order in the perturbation amplitude.

## 4 Numerical results

### 4.1 Surface magnetization

The surface magnetization for a given configuration of the couplings $`\{J_k\}`$ was obtained using Eq. (14) at criticality, i.e., with $`h_k=1`$ and $`J=1`$ in (3) and (4).

The moments were averaged either by performing exact enumerations of the random configurations in the case of chains with lengths $`L=4`$$`24`$ or through Monte Carlo samplings (over $`10^6`$ samples) for longer chains with lengths $`L=2^5`$$`2^{14}`$. The exponents $`x_\mathrm{m}^\mathrm{s}(q)`$ were deduced from an extrapolation of two-point approximants using the BST algorithm . The marginal exponents are shown in Fig. 1 for model #1 and in Fig. 2 for model #2.

### 4.2 Surface energy density

The surface energy density was obtained using the expression given in (17) where $`ϵ_1`$, $`ϵ_2`$, $`\varphi _1(1)`$ and $`\varphi _2(1)`$ are obtained through the diagonalisation of a $`L\times L`$ matrix constructed from the elements of $`𝖳^2`$ with odd indices.

The average moments were deduced from exact enumerations for chains with lengths $`L=4`$$`24`$ and from Monte Carlo samplings (over $`10^4`$ samples) for chains with lengths $`L=2^2`$$`2^9`$. The exponents $`x_\mathrm{e}^\mathrm{s}(q)`$ were deduced from an extrapolation of two-point approximants using the BST algorithm in the case of exact enumerations. The approximants deduced from the Monte Carlo data were too noisy to use the BST extrapolation process. In this case the exponents were deduced from a non-linear fit of the following expression:

$$\mathrm{ln}e_\mathrm{s}^q(L)\mathrm{ln}Ax_\mathrm{e}^\mathrm{s}(q)\mathrm{ln}L+BL^\omega .$$
(36)
The continuously varying exponents are compared to the analytical expressions in Fig. 3 for model #1 and in Fig. 4 for model #2.

## 5 Discussion

A good agreement between numerical and conjectured analytical results have been obtained for the surface magnetization and the surface energy density, except in the vicinity of first-order surface transitions.

Already for the pure HvL model a slow convergence to the exactly known behaviour was observed close to the first-order regime . One may verify on Figs. 2 and 4 for model #2 that the Monte Carlo data, obtained on larger systems, are closer to the analytical results when $`a_2>0`$. Nevertheless, like for the the pure system, some of the scaling assumptions leading to the analytical expressions for the surface exponents cannot be valid beyond some critical value of the perturbation amplitude since $`x_\mathrm{m}^\mathrm{s}(q)`$ and $`x_\mathrm{e}^\mathrm{s}(q)`$ have to remain non-negative.

The dual magnetization in (15) was assumed to scale as $`L^{1/2}`$, like for the unperturbed system at the ordinary transition. We performed a numerical study of the moments of $`\overline{m}_\mathrm{s}^\mathrm{d}(L)`$ for the two models using Eq. (16), exact enumerations for chains with lengths $`L=4`$$`24`$ and BST extrapolations. The results are shown in Fig. 5.

For model #1, the ordinary surface behaviour is obtained for any value of the perturbation amplitude which gives some confidence in the conjectured analytical expressions in the regime of second-order transition. For model #2, the ordinary surface behaviour is only obtained when $`a_2.2`$. Beyond this value deviations from the ordinary surface behaviour are observed and the analytical expressions for the exponents are probably no longer valid.

Useful discussions with Dragi Karevski and Ferenc Iglói are gratefully acknowledged.


text:
# Relic Black Holes and Type 2 Quasars

## 1. Introduction

It has long been suspected that the most luminous Active Galactic Nuclei (AGN) are powered by accretion of matter onto massive black holes (BH). This belief, combined with the observed evolution of the space-density of AGN suggests that massive BHs should be present in local galaxies as relics of past activity. Indeed, the current observational evidence indicates that most, possibly all, luminous galaxies host a massive BH at their centres and its mass correlates with the luminosity and stellar velocity dispersion of the bulge in which it resides (e.g. Richstone et al. 1998, Ferrarese & Merritt 2000). Soltan (1982) and Chokshi & Turner (1992) derived the local mass density in relic BHs,  $`\rho _{\mathrm{BH}}`$ , expected from past activity using optical number counts of AGNs. More recently, Fabian & Iwasawa (1999) estimated  $`\rho _{\mathrm{BH}}`$  from the energy density of the hard X-Ray Background (XRB) and this value is in agreement with the observed local BH density derived by Merritt & Ferrarese (2001) using the BH-mass vs Bulge-mass relation. However it seems that the AGNs found in optical and soft X-ray surveys are not enough to account for the local  $`\rho _{\mathrm{BH}}`$  and this fact prompted Fabian & Iwasawa to talk about ”obscured accretion”, i.e. obscured AGN activity. The existence of obscured AGNs (the so-called ”type 2” objects) is required to account for the XRB spectrum with integrated emission from AGNs (e.g. Gilli, Salvati & Hasinger 2001 and references therein): the ratio between type 2 and unobscured type 1 objects (hereafter $`\mathrm{R}_{21}`$) is usually in the range 4-10, implying that both optical and Xray surveys are missing a large fraction of the AGN population. With a completely independent argument, the missing Lyman edge in blue quasars implies that $`\mathrm{R}_{21}`$ is at least $`0.5`$ (Maiolino et al. 2001). Locally, low luminosity obscured AGNs are well known (Seyfert 2 galaxies) and in the ratio of about 4:1 with respect to type 1 objects (Maiolino & Rieke 1995). However, at higher redshifts few type 2 quasars (high luminosity objects) are known, a number which is largely insufficient to account for the spectral shape of the XRB. What is the relationship between type 2 AGNs and massive black holes? The arguments based on  $`\rho _{\mathrm{BH}}`$  indicates that some obscured AGNs must exist in order to account for the fossil mass but nothing can be told about the luminosity of these obscured objects: are many low luminosity AGNs needed or could fewer high luminosity quasars suffice? Is it really true that AGN activity can explain the local BHs, not only globally, but also in terms of their mass distribution? This paper will try to answer these questions by comparing the local BH mass function (BHMF) expected from AGN activity with the local BHMF inferred from galaxy bulges. Throughout this paper we assume $`H_{}=70`$km/s/Mpc-1, $`\mathrm{\Omega }_\mathrm{M}=1`$ and $`\mathrm{\Omega }_\mathrm{\Lambda }=0`$.

## 2. The continuity equation

We consider the continuity equation which regulates the time evolution of the BH mass function with the Small & Blandford (1992) formalism. If $`N(M,t)`$ is the BH mass function ($`N(M,t)\mathrm{d}M`$ is the number of BHs per unit comoving volume (Mpc-3) at cosmic time $`t`$), the continuity equation can be written as:

$$\frac{N(M,t)}{t}+\frac{}{M}\left(N(M,t)<\dot{M}(M,t)>\right)=S(M,t)$$
(1)
where $`<\dot{M}(M,t)>`$ is the ”average” accretion rate and $`S(M,t)`$ is the source function. As in Small & Blandford (1992) we assume that a massive BH is ”born” when its mass is $`MM_{\mathrm{min}}`$ and that no BH is created or destroyed (i.e. we neglect any merging process). With these assumptions the source function is simply $`S(M,t)=S_{}\delta (MM_{\mathrm{min}})`$ where $`\delta `$ is the Dirac delta function. If AGNs are powered by accretion onto BHs, we can relate the AGN Luminosity Function $`\varphi (L,t)`$ ($`\varphi (L,t)\mathrm{dln}L`$ is the number of AGNs per unit comoving volume (Mpc-3) at cosmic time $`t`$) to the BH mass function

$$\varphi (L,t)\mathrm{d}\mathrm{ln}L=\delta (M,t)N(M,t)\mathrm{d}M$$
(2)
where $`\delta (M,t)`$ is the fraction of BHs with mass M active at time $`t`$. If a BH is accreting at a fraction $`\lambda `$ of the Eddington rate its emitted luminosity is then $`L=\lambda Mc^2/t_E=\epsilon \dot{M}c^2`$ where $`t_E`$ is the Eddington time, $`\epsilon `$ is the accretion efficiency and $`\dot{M}`$ is the accretion rate. Since $`<\dot{M}(M,t)>=\delta (M,t)\dot{M}(M,t)`$, we can write:

$$N(M,t)<\dot{M}(M,t)>=\frac{\lambda }{\epsilon t_E}\left[\varphi (L,t)\right]_{L=\lambda \frac{Mc^2}{t_E}}$$
(3)
$$\frac{N(M,t)}{t}=\frac{c^2\lambda ^2}{\epsilon t_E^2}\left[\frac{\varphi (L,t)}{L}\right]_{L=\lambda \frac{Mc^2}{t_E}}$$
(4)
which can be easily integrated given the AGN luminosity function AND the initial conditions. When $`\lambda `$ is dependent on $`M`$ the equation is more complex but still easily solvable. As initial conditions, we assume that at a given redshift $`z_{}`$ ($`t_{}=t(z_{})`$) ALL Black Holes are active, i.e. that $`\delta (M,t(z_{}))=1`$. The only constraint is that $`\delta (M,t)=\varphi (L,t)/MN(M,t)1`$.

A fundamental ingredient is the AGN luminosity function (LF) expressed in terms of the total AGN luminosity, $`L`$. We first consider the LF of optically selected quasars (OLF in the B band; Boyle et al. 2000) and the LF of soft X-ray selected AGNs (XLF in 0.5-2 keV band; Gilli et al. 2001, it represents only unabsorbed AGNs and is derived using the XLF by Miyaji, Hasinger & Schmidt 2000). Both LFs consider only type 1 (unobscured) AGNs. In principle the LF in one spectral band should be enough to obtain $`\varphi (L,t)`$ but in this paper we consider more instructive to combine the above LFs using Boyle’s LF to describe high $`L`$ AGNs and Gilli’s one to describe low $`L`$ AGNs. A more refined treatment will be presented elsewhere. The AGN luminosity function for type 1 objects is plotted in Fig. 1a at two redshift values. We have assumed $`L/L_\mathrm{X}=30`$, for low L objects and $`L/\nu _\mathrm{B}L_{\nu _\mathrm{B}}=15`$ as typical of quasars (e.g. Elvis et al. 1994).

## 3. The local BH mass function

To estimate the local BH mass function we follow the method by Salucci et al. (1999). The BHMF is derived from the Galaxy LF through the following steps: (1) we consider the LF of galaxies in the B band divided by morphological type (Folkes et al. 1999); (2) for each morphological type we derive the LF of bulges applying the bulge/disk ratio by Simien & de Vaucouleurs (1985); (3) we obtain the mass function of bulges using the mass-to-light ratio by Fukugita, Hogan, & Peebles (1998); (4) we convolve the bulge mass function with the probability distribution of the ratio $`M_{\mathrm{BH}}/M_{\mathrm{Bulge}}`$ to obtain the BH mass function. An important assumption in the last point is that the observed scatter in the $`M_{\mathrm{BH}}`$\- $`M_{\mathrm{Bulge}}`$ relation is real and not due to observational uncertainties. We then assume a gaussian probability distribution $`P(\mathrm{log}x)`$ ($`x=M_{\mathrm{BH}}/M_{\mathrm{Bulge}}`$) with $`<\mathrm{log}x>=2.93`$ and $`\sigma _x=0.45`$ as derived by Merritt & Ferrarese (2001). The BHMF is presented in Fig. 1b and the ”grey” region represent an error estimate on the BH mass function obtained by varying the parameters used in each of the above steps. The BHMF has $`\text{ }\rho _{\mathrm{BH}}\text{ }=(5\pm 2)\times 10^5\mathrm{M}_{}\mathrm{Mpc}\text{-3}`$ and this value is in agreement with $`\text{ }\rho _{\mathrm{BH}}\text{ }=(3÷6)\times 10^5\mathrm{M}_{}\mathrm{Mpc}\text{-3}`$ obtained by Fabian & Iwasawa (1999) and Salucci et al. (1999) from the hard XRB.

## 4. Model Results

Model results are shown in Fig. 2a. The solid line in the grey area is the BH mass function derived from the galaxy mass function as indicated in the previous section. The solid line is the BH mass function expected from type 1 AGNs only with a starting redshift of $`z_{}=3`$. The starting redshift at which we make the assumption of $`\delta (M,t_{})=1`$ (or alternatively the redshift at which we select the active BHs whose evolution we will follow) is not important for the final results. Indeed the dotted line representing the BHMF obtained assuming $`z_{}`$=5 and $`\delta =0.1`$, is not significantly different from the one with $`z_{}`$=3 and $`\delta =1`$. The only constraint is that $`z_{}>z_{\mathrm{max}}`$, the redshift of maximum AGN activity. From Fig. 2a it is clear that Optical QSOs (high $`L`$ AGNs) produce black holes with masses $`M_{\mathrm{BH}}>10^8\mathrm{M}_{}`$ while X-ray AGNs (low $`L`$ AGNs) produce the lower mass ones. Moreover, comparing the BHMF obtained from AGN activity at $`z=0`$ with the one at $`z_{}`$ it can be clearly seen that BH masses are produced during the AGN activity: this is the main reason why the results are little sensitive to the initial conditions. The model presented in the figure has $`\mathrm{R}_{21}`$=0, i.e. it considers only unobscured AGNs. Thus it is not unexpected to derive a low  $`\rho _{\mathrm{BH}}`$ =1.8$`\times 10^5`$$`\mathrm{M}_{}`$$`\mathrm{Mpc}`$-3. Let’s first consider the high $`L`$ AGNs. It is not yet clear if type 2 quasars exist or not or what is their ratio $`\mathrm{R}_{21}`$ relative to type 1 objects. The comparison between the ”observed” local BHMF with the one deduced from AGN activity can yield a constraint on $`\mathrm{R}_{21}`$. From Fig. 2b, $`1<\mathrm{R}_{21}<2`$ is required to match the local BHMF in absolute level (the slope is nicely reproduced). Correspondingly $`\text{ }\rho _{\mathrm{BH}}\text{ }=3.2÷4.6\times 10^5M_{\mathrm{BH}}\mathrm{Mpc}\text{-3}`$. Let’s now focus on the low $`L`$ AGNs. The existence of low $`L`$ type 2 AGNs is well know and the local ratio is $`\mathrm{R}_{21}`$=4 (e.g. Maiolino & Rieke 1995). As shown in Fig. 3, with that ratio the local BHMF at low masses ($`M_{\mathrm{BH}}<10^8\mathrm{M}_{}`$) is overproduced. A change in the bolometric correction does not improve the situation because the BHMF is just shifted along the M axis. A possible solution to this problem is that low luminosity AGNs are not emitting at the Eddington Luminosity (Fig. 3a), a possibility also considered in other works (e.g. Salucci et al. 1999 and references therein), and indeed, apart for some ”wiggles”, the BHMF is now well reproduced. One must keep in mind that this might also mean that there are some problems with the AGN luminosity functions at low luminosities or that the $`M_{\mathrm{BH}}`$-$`M_{\mathrm{Bulge}}`$ relation is different at low masses.

## 5. Summary and Conclusions

We have used the continuity equation to follow the evolution of the BH mass function in the assumption that mass accretion onto a BH is the powering mechanism of AGN activity. AGNs accrete at a fraction of their Eddington luminosity and at a starting redshift $`z_{}`$ all BHs were active. By comparing the expected BHMF with the one derived from local galaxy luminosity functions we can make the following conclusions. (1) The majority of the mass of a local BH was produced during AGN activity. (2) Quasar activity was responsible for the growth of BHs with $`M>10^8\mathrm{M}_{}`$. (3) The final results are not sensitive to the value of the redshift $`z_{}`$ at which all BH were active. (4) The ratio between type 2 and type 1 quasars is $`\mathrm{R}_{21}12`$, i.e. $`30\%50\%`$ of  $`\rho _{\mathrm{BH}}`$  is made by type 2 QSOs. (5) The required $`\mathrm{R}_{21,B}`$ can be reduced if we assume that $`L<L_{Eddington}`$ or $`\epsilon <0.1`$. In any case it is clear that (a) $`\mathrm{R}_{21,B}`$ is much lower than the $`410`$ value required to fit the XRB and (b) no small efficiencies are required to explain the local BHMF. (6) $`\mathrm{R}_{21,X}>0`$ gives a problem at low masses overpredicting the number of relic BHs. This could indicate that either $`L<L_{Eddington}`$ at low luminosities or that there might be uncertainties on the BHMF and/or the AGN LF at low M and L. (7) No merging of BHs is required by our model.

The main conclusion of this work is that with ”reasonable” assumptions, compatible with current knowledge, mass accretion onto massive BH during AGN activity can explain the local BH mass function both in shape and normalization. There is not much room for a large ($`R>>1`$) population of high luminosity obscured AGNs (i.e. type 2 quasars).

## References

Boyle, B.J., et al. 2000, MNRAS, 317, 1014

Chokshi, A., & Turner, E.L. 1992, MNRAS, 259, 421

Elvis, M., et al. 1994, ApJS, 95, 1

Fabian, A.C., & Iwasawa, K. 1999, MNRAS, 303, L34

Ferrarese, L., & Merritt, D. 2000, ApJ, 539, L9

Folkes, S., et al. 1999, MNRAS, 308, 459

Fukugita, M., Hogan, C.J., & Peebles, P.J.E. 1998, ApJ, 503, 518

Gilli, R., Salvati, M., & Hasinger, G. 2001, A&A, 366, 407

Macchetto, F., et al. 1997, ApJ, 489, 579

Maiolino, R., & Rieke, G.H. 1995, ApJ, 454, 95

Maiolino, R., Salvati, M., Marconi, A., & Antonucci, R. 2001, A&A, 375, 25

Marconi, A., et al. 2001, ApJ, 549, 915

Merritt, D., & Ferrarese, L. 2001, MNRAS, 320, L30

Miyaji, T., Hasinger, G., & Schmidt, M. 2000, 353, 25

Richstone, D. et al. 1998, Nature, 395, 14

Salucci, P., et al. 1999, MNRAS, 307, 637

Simien, F., & de Vaucouleurs, G. 1985, ApJ, 302, 564

Small, T.A., & Blandford, R.D. 1982, MNRAS, 259, 725

Soltan, A. 1982, MNRAS, 200, 115


text:
# Energy Release and Transport Processes in the Centres of Galaxies

## 1. Introduction

This is an exciting time in the study of active galactic nuclei. Radio, optical and X-ray observations are not only persuading the last doubters that most normal galaxies contain massive black holes, but are measuring of their masses \[Dressel, Baan, Greenhill\] and probing their spins. Furthermore, the high angular resolution accessible to ground-based and space-based observation of relativistic jets is allowing us to trace these outflows down to radii less than $`100m`$ vindicating the view that they are at least powered by physical processes located fairly close to the central black hole \[Falcke, Hirabayashi\]. Finally, radio polarimetric and $`\gamma `$-ray \[Edwards\] observations of these relativistic jets are being interpreted in a manner that allows us to make arguments concerning their composition, which, in turn, impact our understanding of how they are formed. Many of these topics are described at greater length elsewhere in these proceedings and designated with square brackets. I shall just emphasize some more recent theoretical ideas.

## 2. Accretion Power vs Black Hole Power

As is often the case, the most important factor controlling the evolution of accretion disks is the flow of energy. There are two power sources. Firstly, accreting gas will release its binding energy before it makes its final adiabatic plunge towards the event horizon. However, the energy will generally not be liberated all at the point of release because the magnetic torques which transport angular momentum, also transport energy outward radially through the disk. (There can be additional, convective and radiative energy transport which is ineffective in transporting angular momentum.)

There are several possible channels through which this energy can escape and their relative importance surely depends upon the accretion rate. It can be radiated from a disk photosphere. It can emerge from the disk surface as a hydromagnetic energy flux and heat an active corona, which can in turn radiate non-thermally (for example by Comptonising disk photons) or drive a wind. The coronal heating will probably occur via magnetic reconnection or a turbulent cascade. Energy and angular momentum can also be extracted from the disk if there is a large scale field that launches a hydromagnetic wind. There will then be an external torque per unit radius $`dG/dr`$ which will extract energy at a rate $`𝑑G\mathrm{\Omega }`$ where $`\mathrm{\Omega }`$ is the disk angular velocity without dissipation.

The second source of energy is the spin of the hole. As is well known, up to 29 percent of the mass of a spinning black hole is, in principle, extractable (eg Lee, Wijers & Brown 1999). This is $`6\times 10^{62}(M/10^9\mathrm{M}_{})`$ erg. In practice no more than, say, a quarter of this energy is likely to be available, but still it is ample to account for the non-thermal emission associated with the most profligate quasars. It is important to realize that, although energy is not strictly localisable in general relativity, the spin energy is just as much associated with the so-called ergosphere outside the event horizon as with the interior of the hole. (The ergosphere can be defined as the region within which a distant observer must rotate with respect to the non-rotating frame at infinity.) This is a crucial and important difference from a spinning, newtonian star.

Again there are many possibilities \[Meier\]. It is possible to extract work from orbiting matter which can be dragged by magnetic stress on to lower or even negative energy (including rest mass) orbits within the ergosphere (e. g. Ruffini & Wilson 1975; Agol & Krolik 2000). When this matter enters the horizon, the net gravitational mass of the hole will increase by less than the accreted rest mass and will actually decrease if the orbit has negative energy. If the matter is attached to open magnetic field lines, then work is done on the base of the jet and drives an outflow, either continuously or intermittently. Alternatively, if the infalling matter is magnetically attached to the disk and still in causal contact, that is to say outside the relevant, ingoing critical points, then orbital energy and angular momentum can be transferred to the disk from within its innermost stable circular orbit, thereby enhancing the disk power. In either case, it is an arbitrary choice whether to charge the energy extracted to the hole or the accreting matter. Magnetic flux that threads the horizon can transfer energy directly from the hole to the jet or the disk at a rate $`a^2B^2c`$, where $`a`$ is the specific angular momentum of the hole and $`B`$ is the field strength. To an observer hovering just above the horizon and, necessarily, orbiting with the hole angular velocity, the electromagnetic energy flux will appear to enter the hole. However energy flux is conserved in the non-rotating frame and transforming into this frame produces an outward energy flux at the horizon. In other words, the black hole, a spinning conductor, acts like a giant battery. Under the conditions invoked for Cygnus A, say, an EMF of $`10^{20}`$ V is induced, with a field $`B10^4`$ G, and the current is $`10^{18}`$ A. (The resistance $`100\mathrm{\Omega }`$ turns out to be similar to the impedance of free space. The merit of this process as a candidate for powering relativistic jets is that the electromagnetic power is likely to be relatively uncontaminated with baryons, the disadvantage is that if the field strength above the hole is comparable with that above the inner disk then, the relatively small area of the horizon will ensure that the hole power be less than the disk power, which, although adequate for most jets may not be sufficient for the most luminous examples. Note, though, that this conclusion does not follow if, as commonly envisaged, the disks associated with powerful radio galaxies are thick.

There is actually a second way to extract energy from the hole. If the disk plane wobbles about the equatorial plane of the hole, then there will be a “Bardeen-Petterson” torque that acts between the disk and the hole. The nature of this torque has been quite controversial over the years, but the net result is that the disk will be brought into alignment with the hole within some radius, $`r_{\mathrm{BP}}`$. There will be a reaction torque operating on the hole and an associated dissipation in the disk near $`r_{\mathrm{BP}}`$. This could provide a significant local heating of the disk under some models (eg Natarajan & Pringle 1998)

## 3. Jets

So far, I have emphasized the (electro)magnetic extraction of energy from the disk and the hole because most models are now of this type. However, this does not imply that observed jets mostly carry Poynting flux. Just as is the case with the Crab Nebula, there is believed to be a transformation of the energy flux from mechanical to electromagnetic to particle. The issue now is where and how do these transitions occur. One extreme possibility is that a hot plasma forms through magnetic dissipation in the corona itself and that this flows hydrodynamically into the jet. (This would probably have to be an electron-ion plasmas as it is very hard to see how a pair plasma could be accelerated to ultrarelativistic bulk speed close to the hole without suffering catastrophic radiative drag.) Another possibility is that ultrarelativistic ions are formed close to the hole and somehow directed along a jet and that these cascade down to lower energy particles that will produce the radiation we observe. In many regards the most natural scheme is that the energy is carried away electromagnetically from the vicinity of the hole and that the pair fraction gradually builds up at a rate that may be limited by annihilation. Provided that this happens after the radiative drag is ignorable, an ultrarelativistic jet will ensue.

Extracting the energy is only half of the challenge in forming a jet. The outflow must also be tightly collimated. For a long while it has been thought that jets required disks. Howevever, with the discovery of jets associated with the Crab and Vela pulsars, perhaps this is not even necessary. Magnetic collimation is most commonly invoked simply because it is hard to sustain the necessary gas pressure against large radiative loss at the location where the collimation must be occuring. (Note that this does not imply that gas pressure is irrelevant. Magnetic stress must be ultimately confined transversely at large cylindrical radius and static gas pressure or ram pressure is surely responsible.)

There are still choices that have to be made in magnetic collimation. The confining stress can be attributed to the pressure of poloidal magnetic field or the tension associated with the toroidal field. Alternatively, the field may be quite disorganised but show an anisotropic stress tensor which somehow achieves the necessary collimation.

One argument that has often been made against magnetic collimation, particularly by toroidal field is that such configurations are probably unstable. However, this may not be such an unwelcome conclusion because observed jets often do appear to corkscrew in real space and this is a natural instability mode for a pinched jet. At an even more extreme level, many jets appear to be episodic and to exhibit strong internal shocks \[Klare\]. This, too, can be regarded as the nonlinear development of a jet instability. One possibility is that the natural increase in the ratio of the toroidal to the poloidal field with distance along the jet is balanced by steady magnetic reconnection, which limits the growth of toroidal field and lead to the steady acceleration of electrons and positrons.

In order to construct more quantitative models of jets, it is necessary to make better estimates of the power that they carry (as distinct from the power that they radiate). There are several approaches that have been followed. One is to assume that the pressure is well-estimated by the equipartition pressure and that the jet expands on the Mach cone. In this case the power will be $`L_\mathrm{j}P_{\mathrm{eq}}r^2V_\mathrm{j}`$, a formula that is roughly valid relativistically. Alternatively, if the external pressure is known and the jet is thought to be in pressure equilibriium with its surroundings, then we can replace the equipartition pressure with the external pressure. A third approach is to estimate the power supply to large double lobes. Using arguments like this it has been possible to argue that the jet power in some sources, including M87 \[Junor\] and some blazars, exceeds the bolometric power. If this were generally true, then it would effectively rule out many jet formation models. The mid infrared observations of SIRTF will be crucial for measuring the total nuclear luminosity of radio galaxies.

## 4. Brightness Temperature and Polarisation

One of the more interesting and potentially significant developments in recent years has been the accumulation of evidence for high brightness temperatures in compact extragalactic radio sources. Intraday variable sources exhibit measured brightness temperatures that are at least as high as $`3\times 10^{13}`$ K and could be significantly larger \[Jauncey, Dennett-Thorpe, Kedziora-Chudczer\]. This is conventionally interpreted in terms of the traditional inverse Compton limit. For the observed ratios of Compton to synchrotron power in jets, the comoving brightness temperature has an upper limit of $`10^{12}`$ K. One power of the Doppler factor, which is twice the relevant Lorentz factor converts to the comoving frame suggesting that $`\mathrm{\Gamma }>20`$. (Large Lorentz factors are not so intimidating now as values of $`\mathrm{\Gamma }300`$ are derived for $`\gamma `$-ray bursts although the radiative efficiency may be low at radio wavelengths.) If it turns out that comoving brightness temperatures well in excess of the Compton limit are required then it may be necessary to invoke a coherent emission mechanism. This is not, by itself, unreasonable; after all, the Sun and Jupiter support coherent emission under relatively benign conditions. A cyclotron maser is a particularly enticing possibility. What may be a bit more problematic is the escape of very high brightness radiation from the source when it is subject to nonlinear induced Compton and Raman scattering (Levinson & Blandford 1995).

Now, in a self-absorbed synchrotron source, close to the inverse Compton limit, the energies of the emitting electrons at the radio photosphere must be $`300`$ MeV. We do not see the emission from lower energy electrons and it is interesting to ask if they are present and if the electrons that we do see have proton or positron partners. As argued by Wardle et al. (1998 and references cited therein), in the former case, the lower cut off in the electron distribution function must be $`50`$ MeV, otherwise there will be too much jet kinetic energy, carried mainly by the protons. In the latter case, the pair distribution function will extend down to mildly relativistic energy Similar conclusions are drawn from the observation of strong linear polarization which limits the internal Faraday rotation within the source. Finally the detection of persistent circular polarisation \[Macquart\] has been used to rule in favor of a pair plasma, on the grounds that the frequency dependence is better fit by Faraday conversion in a pair plasma than intrinsic circular polarisation in a proton plasma. Note that if the radiation mechanism is coherent cyclotron radiation, then the fact that the polarisation is not even stronger may be an indication that there is a pair plasma present.

## 5. Galactic Centre

The 2.6 million solar mass black hole at our Galactic center is in many respects the best example to study because it is so inactive \[Eckart\]. It has recently been detected in X-rays with a very soft 2–10 keV luminosity (Baganoff 2000) of $`2\times 10^{33}`$ erg s-1, which may be variable. This is a strikingly low luminosity. The independently estimated mass supply rate is very roughly $`10^{22}`$ g s-1 and yet the bolometric luminosity is $`10^{36}`$ erg s-1 corresponding to an apparent radiative efficiency $`10^7c^2`$ \[Coker\]. As has been discussed elsewhere by Blandford & Begelman (1999), this is consistent with the theoretical prediction that adiabatic (in the sense that there is no loss of heat) accretion flows lose a large fraction of their mass before crossing the horizon. This is in contrast to the prediction of the widely-discussed and conservative ADAF models. (A naive interpretation of the X-ray observations suggests a maximum mass accretion rate of $`10^{20}(T_e/10^{10}\mathrm{K})^{1/4}`$ g s-1. It could plausibly be much lower.) The reported spectrum suggests that the X-ray emission is nonthermal, which is really expected, because it would be remarkable if trans-sonic, trans-Alfvénic, trans-relativistic flows avoided accelerating relativistic particles.

If the mass accretion rate really is very low, then this opens up the possibility that we may be seeing radio or, more likely, mm emission from very close to the hole and that polarization observations might be quite diagnostic of the physical conditions. Now, in the case of Sgr A, there is no linear polarization reported at radio wavelengths (although it has been reported, quite surprisingly, at submm wavelengths). The lack of linear polarization is quite consistent with the high Faraday rotation expected from an accretion flow onto a black hole. However, variable circular polarization has been reported with a persistent sign. There is the intriguing possibility that this reflects the orbital motion of the plasma around the central black hole. To see how this might come about, consider a magnetoactive plasma with $`X=\omega _P^2/\omega ^2,Y=\omega _G/\omega <1`$. Under quasi-longitudinal conditions, essentially $`\mathrm{cos}\theta <Y`$, where $`\theta `$ is the angle between $`\stackrel{}{k}`$ and $`\stackrel{}{B}`$, the electromagnetic eigenmodes are elliptically polarized with axis ratio $`r=1\pm Y\mathrm{sin}\theta \mathrm{tan}\theta `$ and phase velocity difference $`\mathrm{\Delta }V=cXY\mathrm{cos}\theta `$.

Now consider synchrotron or cyclotron radiation emitted from within an accretion disk of thickness $`H`$. The major axis of the polarization ellipse will be Faraday rotated at a rate $`\mathrm{\Delta }\mathrm{\Phi }/ds=\mathrm{\Delta }V\omega /2c^2`$. Now, we expect that the magnetic field will reverse often along a ray and as $`X|Y|\omega H/c`$ is anticipated to be much greater than unity, then we expect that the emergent linear polarisation will be vanishingly small. (The limiting polarisation along an individual ray is likely to be determined by the decrease in the density rather than the magnetic field strength.)

The circular polarisation is a bit more problematic. If we suppose that there is a net magnetic field normal to the disk, as is true of some models, then there should be a preferred sense of the circular polarization that is emitted in cyclotron or low energy synchrotron radiation. This will be largely preserved in propagating out of the disk. However, eventually, the sign of the net vertical field must change and we might even be able to see this change propagate from low to high frequency if, as seems reasonable, there is a radius to frequency mapping in the disk. However, it is also possible that there is essentially no large scale field. In this case, there is a second mechanism which may create circular polarisation.

Consider a single eigenmode propagating out of the disk through a spatially varying magnetic field. We suppose that the variation happens relatively slowly on the scale of the wavelength so that the polarization adjusts adiabatically (with no mode crossings). We ignore refraction and recall that the electric vector is normal to $`\stackrel{}{k}`$. Next, suppose that the two eigenmodes are launched with equal amplitude and that the field is uniform. The beating between the two eigenmodes will result in a circular polarisation of amplitude $`Y\mathrm{sin}\theta \mathrm{tan}\theta `$ that changes sign as the plane of linear polarisation rotates. The limiting circular polarisation that emerges from the disk is equally likely to have either sign and so there will be no net circular polarisation. Now let the magnetic field direction vary along a ray. If the angle $`\theta `$ varies, then there will be a corresponding change in the axis ratio of the polarization ellipse, but still no preference for one sign over the other. However, if the azimuthal angle $`\varphi `$ relating $`\stackrel{}{B}`$ to $`\stackrel{}{k}`$ changes in a systematic fashion along all rays then it appears that the circular polarisation will have one sign as opposed to the other for slightly more than half the time. The average circular polarization along a ray will be $`<Y\mathrm{sin}\theta \mathrm{tan}\theta d\varphi /d\mathrm{\Phi }>|Y\mathrm{ln}Y|d\varphi /d\mathrm{\Phi }`$.

Now this is just what we might anticipate in a magnetised accretion disk. In the disk interior, the typical field direction will trail to reflect the differential rotation in the disk. However, the field will be swept back by progressively smaller angles as the ray approaches the disk surface, corresponding to a net rotation of the average azimuthal angle $`\varphi `$. There will only be a preferred sense of limiting circular polarization, if the magnetostatic field is still changing in this systematic manner over the last radian of Faraday rotation. The net circular polarization will be $`c\mathrm{ln}Y/\omega HX`$, where $`H`$ is the disk thickness. However, if $`X>10^6`$, as seems assured in Sgr A, then this mechanism is unlikely to create the observed circular polarisation at radio wavelengths, although it may be relevant in the sub mm. I conclude that the measurement of persistent circular polarization is more naturally interpreted in terms of a long term, ordered magnetic field, as discussed above. There may be a related effect that can occur in the ergosphere, where the dragging of inertial frames is reponsible for a rotation of the polarisation ellipse. These ideas are currently under study.

## 6. Black Hole Formation

The foregoing considerations apply to sub-critical accretion. Similar principles are relevant to super-critical accretion. Here the problem is not creating the photons, but allowing them to escape. Electron scattering opacity almost certainly dominates the radiative transfer and the photons will be trapped out to a radius $`r_{\mathrm{tr}}=\dot{M}\kappa _T/4\pi c`$. Within this radius, we can regard the inflow as adiabatic, with an outward transport of energy and angular momentum. I have conjectured elsewhere that supercritical accretion is “demand-limited” rather than “supply-driven”. In other words black holes only grow at roughly the critical rate, $`\dot{M}=4\pi GM/ϵ\kappa _Tc`$, where $`ϵ0.1`$ is the overall efficiency of energy release. If a black hole is supplied at a greater rate, it will radiate at a few times the Eddington luminosity and the excess mass and energy will be carried off in an outflow with a similar power. The mass may be liberated at all radii so that the local net mass accretion rate grows roughly linearly with radius. Alternatively, there could be strong, internal circulation which carries the energy released near the black hole to some outer radius from where the wind is launched. This behaviour is, arguably, observed in Broad Absorption Line Quasars and some of the most luminous, Galactic, transient sources \[Spencer\].

If this principle were correct, it would raise an important question. How do black holes grow to a billion solar masses, when the universe was less than $`1`$ Gyr old (as observations of high redshift quasars now suggest) if there were only time for a few e-foldings in mass? One possible answer is that large holes might have formed by dynamical collapse. However, this seems a bit unlikely as they should have made a sort of mega gamma ray burst, and there is no evidence for such prolonged, luminous explosions on our past light cone. Alternatively, what may be different about accretion during the early epochs of galaxy formation is that the black hole should be surrounded by a chaotic, self-gravitating reservoir of gas, extending from the Bondi radius, $`r_\mathrm{B}GM/\sigma ^2`$ where $`\sigma `$ is the velocity dispersion determined by the dark matter, stars and gas to some outer radius, $`r_{\mathrm{gal}}`$. We expect the gas initially to dominate the mass in the inner parts of the galaxy.

Now, modifying a suggestion of Silk & Rees (1998), suppose that a small black hole forms. The maximum, net mass supply rate is the Bondi (or Jeans’) rate, $`\dot{M}\sigma ^3/G`$ under these conditions. If this were to accrete onto the black hole and release energy efficiently, the resulting wind would quickly blow away all the surrounding gas and shut off the supply. What seems most likely is that the feedback limits the accretion rate to a value such that the wind power is just sufficient to drive away the accreting gas, ie to $`\dot{M}\sigma ^5/ϵGc^22\times 10^3(\sigma /100\mathrm{k}\mathrm{m}\mathrm{s}^1)^5`$ M yr-1. The hole grows until it has a mass $`M2\times 10^5(\sigma /100\mathrm{k}\mathrm{m}\mathrm{s}^1)^5`$ M after which point the accretion rate will probably be Eddington-limited and the hole mass will grow exponentially until the wind associated with a sub-critical accretion (formed with an efficiency $`ϵ_W`$) can limit the accretion again. If we allow this phase to continue for a dynamical time, $`3r_{\mathrm{gal}}/\sigma `$ then the final black hole mass will be $`M3\sigma ^4r_{\mathrm{gal}}/ϵ_WGc^26\times 10^6(\sigma /100\mathrm{k}\mathrm{m}\mathrm{s}^1)^4(r_{\mathrm{gal}}/10\mathrm{k}\mathrm{p}\mathrm{c})(ϵ_W/10^3)^1`$ M. Choosing plausible values $`ϵ_W10^3;r_{\mathrm{gal}}10\mathrm{k}\mathrm{p}\mathrm{c}`$ give a relation similar to that found by Ferrarese et al. (2000) and Gebhardt et al. (2000), \[Marconi, van der Marel\].

Obviously, this is a gross over-simplification of the many complex processes at work, simultaneously, in a nascent galaxy. Radiative cooling and star formation will be stimulated in the surrounding gas and this may cause additional heating. The direct and indirect influence of the hole may soften the potential. Rotation and mergers are surely relevant. Most importantly, I have not specified whether or not the gas falls in towards the nucleus or moves out towards $`r_{\mathrm{gal}}`$ – whether galaxies form from the outside in or the inside out. This raises a radical possibility. Perhaps early-type galaxies are so called because their central black holes grow quickly and expel enough gas to prevent the subsequent formation of disks. This scheme, which is subject to hydrodynamic simulation, has clear, observable implications.

### Acknowledgments.

I am indebted to many collaborators, including Avery Broderick, Marco Salvati, Todd Small and especially Mitch Begelman for discussion of the above. Support under NSF grant AST 99-00866 and NASA grant 5-2837 is gratefully acknowledged.

## References

Agol. E. & Krolik, J. 2000, ApJ, 528, 161

Baganoff, F. et al. 2000 Paper presented at November 2000 HEAD meeting, Hawaii

Blandford, R. D. & Begelman, M. C. 1999, MNRAS, 303, L1

Ferrarese, L. & Merritt, D. 2000, ApJ, 539, L9

Gebhardt, K. et al. 2000 ApJ, 539, L13

Lee, H-K, Wijers, R. A. M. J. & Brown, G. E. 1999, Phys. Rep., 325, 83

Levinson, A. & Blandford, R. D. 1995, MNRAS, 274, 717

Nandra, P. 2000 Proc. X-ray Astronomy ’99 - Stellar Endpoints, AGN and the Diffuse Background, in press

Natarajan, P. & Pringle, J. 1998, ApJ, 506, L97

Ruffini, R. & Wilson, J. R. 1975, Phys.Rev.D, 12, 2959

Silk, J. & Rees, M. 1998, A&A, 331, L1

Wardle, J. F. C. et al. 1998, Nature, 395, 457

Wehrle, A. et al. 1998, ApJ, 497, 178

Wilson, A., Shopbell, P. & Young, A. 2000, ApJ, in press


text:
# Untitled Document

A Simple Solution to the Strong CP Problem

Sheldon Lee Glashow

Department of Physics

Boston University

Boston, MA 02215

We propose a minimal modification of the standard model, remarkable in its simplicity, which may solve the strong CP problem. It employs three Higgs doublets with interactions taken to be invariant under a flavor symmetry. Both CP and the flavor symmetry are softly broken by Higgs boson mass terms. In tree approximation, quark mass matrices are triangular and $`\mathrm{arg}detM`$ vanishes. Radiative corrections lead to a tiny and tolerable value of $`\overline{\theta }`$.

10/01

The problem is simply stated. The minimal three-family standard model involves two distinct CP-violating phases: the $`\delta `$ parameter in the Kobayashi-Maskawa matrix and $`\overline{\theta }`$, the measure of strong CP violation, both of which can be generated by complex terms in the quark mass matrices. The former angle is of order unity, whereas the latter is known to be less than $`10^9`$. The strong CP problem arises if this disparity is regarded as other than coincidental.

There have been many proposed resolutions to this problem of which three classics are: a massless up quark , an invisibilized version of the Peccei-Quinn axion , and the Barr-Nelson implementation of spontaneous CP violation . Although none of these $`\overline{\theta }`$–avoidance mechanisms are decisively excluded, neither is convincingly true. Consequently there have been many other suggested remedies which would be too tedious to enumerate. Some models appeal to supersymmetry , others appeal to technicolor , and still others invoke extra dimensions . In addition, many solutions have been put forward in which CP is softly broken . The last-cited models (of which I am a coauthor) are decisively ruled out by currently available experimental data. In this note I present yet another model making use of soft CP violation, one which is both simpler than its predecessors and seemingly compatible with experiment.

We begin by assigning a flavor quantum number $`F`$ to each quark family. We assign $`F=+1`$ to the right- and left-handed quarks of a nominal first family, $`F=0`$ to those of a second family, and $`F=1`$ to those of the third family. (Similar $`F`$ assignments could be made to the leptons.) In addition, we introduce three doublets of Higgs bosons, $`h_0`$, $`h_1`$ and $`h_2`$, where the subscripts indicate their $`F`$ quantum numbers. Each of these doublets contributes both to the up and down quark mass matrices. Although the Glashow-Weinberg constraint  is not respected, unacceptable flavor-changing effects should be avoidable with a modest degree of fine tuning.

All dimension-4 terms in the Lagrangian are assumed to be both CP and flavor invariant. As a result of this hypothesis, the Yukawa couplings of the $`h_i`$ to quarks (as well as their quartic self-couplings) must be real and must conserve the flavor quantum number $`F`$. Furthermore, the CP-violating Chern-Simons term cannot be present in the Lagrangian. Quadratic expressions in the Higgs fields are not constrained by the above hypothesis. These dimension-2 mass terms are permitted to be complex and $`F`$-violating. This is our proposed mechanism for the explicit but soft violation of both CP invariance and the flavor symmetry.

The parameters of the Higgs portion of the Lagrangian are chosen so that each of the three electrically neutral Higgs bosons develops a complex vacuum expectation value (vev). An appropriate weak-hypercharge rotation lets us choose $`h_0`$ to be real with no loss of generality. The sum of the squared magnitudes of the vevs is constrained to take its conventional electroweak value. We assume that all three vevs are similar in magnitude to avoid the appearance of light surviving bosons (i.e., with masses less than $`100`$ GeV) such as are known not to exist. Observed violations of CP invariance are caused by the complexity of the Higgs vevs, which lead to the irremoveable complexity of the quark mass matrices $`M_U`$ and $`M_D`$. Furthermore, the $`F`$ invariance of the Yukawa couplings ensures that in tree approximation $`M_D`$ is an upper triangular matrix whilst $`M_U`$ is lower triangular. All diagonal entries of these matrices are real and the phases of the off-diagonal entries are constrained in a manner to be discussed.

The triangular mass matrices we obtain are sufficiently general to produce any spectrum of quark masses and can result in any desired Kobayashi-Maskawa matrix. Furthermore, the determinants of these mass matrices are evidently real, so that there is no strong CP problem in tree approximation.

With one possible identification of the observed quarks with the flavor quantum number, we obtain:

$$M_D=\left(\begin{array}{ccc}m_d& ϵ_{12}& ϵ_{13}\\ 0& m_s& ϵ_{23}\\ 0& 0& m_b\end{array}\right)\mathrm{and}M_U=\left(\begin{array}{ccc}m_u& 0& 0\\ ϵ_{21}^{}& m_c& 0\\ ϵ_{31}^{}& ϵ_{32}^{}& m_t\end{array}\right),$$

where the $`m_i`$ are the quark masses in the absence of mixing. Quark mixing is induced by smaller off-diagonal entries, whose complex phases are constrained as follows:

$$\mathrm{arg}ϵ_{12}=\mathrm{arg}ϵ_{21}=\mathrm{arg}ϵ_{23}=\mathrm{arg}ϵ_{32},\mathrm{and}\mathrm{arg}ϵ_{13}=\mathrm{arg}ϵ_{31}.$$

Note that the complexity of the mass matrices is irremoveable except under the coincidental circumstance that

$$\chi \mathrm{arg}(ϵ_{12}^2ϵ_{13}^{})=0.$$

It is illustrative to estimate the Kobayashi-Maskawa parameters in the small-mixing approximation. We obtain:

$$\theta _{12}|ϵ_{12}/m_sϵ_{21}^{}m_u/m_c^2|,\theta _{23}|ϵ_{23}/m_bϵ_{32}^{}m_c/m_t^2|,\theta _{31}|ϵ_{13}/m_bϵ_{31}^{}m_u/m_t^2|,$$

with the $`\delta `$ parameter assuming any desired value. Note that the KM parameters are essentially determined by the off-diagonal terms of $`M_D`$, provided that the $`ϵ_{ij}`$ with $`j>i`$ are comparable to those with $`i>j`$. Putting in numbers for quark masses and mixings, we obtain the rough estimates:

$$ϵ_{12}25\mathrm{MeV},ϵ_{13}13\mathrm{MeV},ϵ_{23}150\mathrm{MeV}.$$

While our model does not suffer a strong CP problem in tree approximation, radiative corrections to quark masses will modify the quark mass matrices and can thereby lead to a non-vanishing value for $`\overline{\theta }`$. Are these corrections small enough to avert a problem? An examination of one-loop corrections to the quark mass matrices reveals that the only such terms are those explicitly involving $`\chi `$, as defined in Eq. (1). The most threatening term by far is a complex contribution to the up-quark mass mediated by charged Higgs bosons, which may be estimated to yield:

$$\mathrm{\Delta }\overline{\theta }=\left(\frac{1}{4\pi }\right)^2\frac{ϵ_{13}ϵ_{23}^{}ϵ_{21}^{}}{(\mathrm{vev})^2m_u}K,$$

where $`K`$ is a dimensionless integral which could be small for an appropriate choice of the spectrum and mixing pattern of the Higgs bosons. In any case, its prefactor is of order $`10^9`$.

Several distinctive features of this model may be testable. In particular, it requires the existence of two surviving charged Higgs bosons, which should be readily detectable at a future linear collider. Furthermore, should one or both of these particles lie much below the top quark mass, we should expect a significant (although not readily calculable) branching ratio for the decay $`tb+h^+`$. The predicted existence of five neutral bosons should offer interesting challenges to experimenters. Although flavor violation via their exchange can be made small, it could be large enough to yield measureable departures from the standard-model description of CP violation.

This work was supported in part by the National Science Foundation, under grant NSF-PHY-0099529.

References

relax For a recent analysis, see: A.G. Cohen, D.B. Kaplan and A.E. Nelson, JHEP 9911 (1999) 27. relax R.D. Peccei and H.R. Quinn, Phys. Rev. Lett. 38 (1977) 387. relax S.M. Barr, Phys. Rev. Lett. 53 (1984) 329; A.E. Nelson, Phys. Lett. B136 (1984) 387. relax E.g., G. Hiller and M. Schmaltz, hep-ph/0105254. relax E.g., K.D. Lane, hep-ph/0106328. relax E.g., M. Chaichan and A.B. Kobakhidze, hep-ph/0011376. relax E.g., P. H. Frampton and S.L. Glashow, Phys. Rev. Lett., 87 (2001) 011801; H. Georgi and S.L Glashow, Phys. Lett., B451 (1999) 372. relax S.L. Glashow and S. Weinberg, Phys. Rev. D15 (1977) 1958.


text:
# BI-TP 2001/24 DESY 01-153 October 2001 Meson Spectral Functions at finite Temperature

## 1 INTRODUCTION

The maximum entropy method (MEM)  has been applied successfully to reconstruct hadronic spectral functions from correlation functions in Euclidean time calculated at zero temperature in quenched QCD . It was demonstrated that the spectral shape reflects ground and excited state contributions, which yields reliable results for the respective masses and decay constants.

At finite temperature, where very little is known about the spectral shape, MEM would provide new insight in the thermal changes of hadronic properties. Moreover, this approach requires no a priori assumptions or specific ansaetze and thus allows a reconstruction of the spectral functions from first principles. A first application of MEM at finite temperature could show that a reconstruction of the continuum part of the free thermal meson spectral function is possible . In the following a similar test is reported for free thermal meson correlators calculated on the lattice. On the basis of such knowledge the MEM analysis is extended to meson correlation functions at several temperatures below and above the deconfinement transition.

## 2 APPLICABILITY OF MEM AT <br>FINITE TEMPERATURE

Instead of sharp ground state peaks at $`T=0`$ broad bumps and continuum-like structures are expected to dominate the spectral shape at sufficiently high temperature. In order to test the applicability of MEM at finite temperature we have calculated free thermal meson correlation functions $`G(\tau )=d^3xJ(\tau ,\stackrel{}{x})J^{}(0,\stackrel{}{0})`$ in the (pseudo-)scalar channel for different lattice sizes projected to zero momentum. A Gaussian noise with the variance $`\mathrm{\Delta }(\tau )\tau G(\tau )`$ was added to the exact values to simulate the statistical error for correlators at finite temperature . The correlator is then related to the respective spectral function $`\sigma (\omega )`$ through the integral equation

$`G(\tau ,T)={\displaystyle _0^{\mathrm{}}}𝑑\omega \sigma (\omega )K(\tau ,\omega ),`$
where the continuum integral kernel is given by

$`K(\tau ,\omega )={\displaystyle \frac{\mathrm{cosh}(\omega (\tau 1/2T))}{\mathrm{sinh}(\omega /2T)}}.`$
The thermal correlators and reconstructed spectral functions are illustrated in figure 1.

Strong finite size effects are visible for the correlators at small time separations compared to the continuum curve (solid line). In the reconstructed spectral function this results in unphysical bumps (dotted curve).

The application of a lattice adapted kernel

$`K_L(\tau ,\omega )={\displaystyle \frac{2\omega }{T}}{\displaystyle \underset{n=0}{\overset{N_\tau 1}{}}}{\displaystyle \frac{\mathrm{exp}(i2\pi n\tau T)}{(2N_\tau \mathrm{sin}(n\pi /N_\tau ))^2+(\omega /T)^2}}`$
is absolutely mandatory to obtain the correct spectral shape. It can be observed in figure 1 that this choice of the kernel combined with a restriction of the $`\omega /T`$-range up to about $`4N_\tau `$ absorbs the cut-off effects and yields an almost perfect reconstruction of the free spectral function $`\sigma (\omega )=\frac{3}{8\pi ^2}\omega ^2\mathrm{tanh}(\omega /4T)`$ for $`N_\tau =8,12`$ and 16.

### 2.1 Smearing at finite Temperature

At zero temperature the application of fuzzing  and/or smearing techniques to optimize the projection onto the ground state leads to a reduction of the peak width in the spectral function and an almost complete elimination of excited states. Therefore one might be tempted to use this method also at finite temperature, where the peaks are less pronounced compared to the continuum. However, the entire concept of smearing has to fail when a single state is not well separated from higher excited states. The situation will naturally be even worse when only a continuum exists. In the case $`T>T_c`$ one can no longer be sure to project on an actual ground state, since even the fuzzing of the free meson correlation function leads to sharp peaks instead of the broad continuum. This is illustrated in figure 2 for different fuzzing radii R. In order not to be biased in the MEM analysis we therefore use unmodified point-point correlators only, which preserve the full information about ground and excited states as well as the continuum contribution.

## 3 MESON SPECTRAL FUNCTIONS

Having tested the applicability of MEM at finite temperature, we use the approach to analyze meson spectral functions in the temperature range 0.4 to 3 $`T_c`$. Gauge field configurations on lattice sizes $`(2464)^3\times 16`$ were generated with the plaquette gauge action, while the Clover action with non-perturbatively improved coefficients  was used in the fermion sector. Temporal and spatial meson correlators were calculated for four quark mass values below $`T_c`$ and at almost zero quark mass (in the vicinity of $`\kappa _c`$) above $`T_c`$.

In general a broadening of the peaks and a reduction in height can be observed at finite temperature due to the short extent in the temporal direction. Nevertheless, a good agreement of the meson ground state masses is obtained comparing the MEM results at 0.4 $`T_c`$ with zero temperature data ($`24^3\times 32`$ lattice ) and conventional exponential fits results of the spatial correlator. The advantage of MEM is most obvious for the temporal vector meson correlator, where the exponential fits lead to overestimated mass values. MEM allows a distinction between the pole and continuum contributions (see fig. 3) and thus detects precisely the actual ground state mass. Error bars are indicated for the average spectral functions in the given $`\omega /T`$-range.

While the situation at0.4 $`T_c`$ and 0.6 $`T_c`$ is similar and almost no thermal effects are visible, the situation is clearly different at 0.9 $`T_c`$. Here we observe a broadening of the spectral functions as well as a shift in the location of the peaks. To what extent this is a physical or statistical effect has to be clarified with larger statistics.

Increasing the temperature above $`T_c`$ a gradual approach towards the free curve is found in the (pseudo-)scalar channel (see fig. 4a). The vector meson correlator is instead much closer to free quark behavior. In this channel we have summed up the contributions with all four $`\gamma `$-matrices, which yields twice the (pseudo-)scalar correlator in the free case. These properties are similarly evident from the reconstructed spectral shape in figure 4b. However, the vector meson spectral function shows an enhancement around $`\omega /T6`$, while a suppression compared to the free curve is found at smaller energies. This observation is in contrast to HTL-resummed perturbation theory , where the spectral function is IR-divergent. The deviations in the (pseudo-)scalar channel cannot be sufficiently explained by HTL medium effects like thermal quark masses and Landau damping. This leads to the conclusion that the behavior in the deconfined plasma phase is characterized by strongly correlated quarks and gluons in the temperature range up to 3 $`T_c`$.