text:
# Untitled Document

Higgs induced light leptoquark-diquark mixing and proton decay

Uma Mahanta

Mehta Research Institute

Chhatnag Road, Jhusi

Allahabad-211019, India

email:mahanta@mri.ernet.in

Abstract

In low energy phenomenology to avoid the strong constraints of proton decay it is usually assumed that light ($``$ 250 Gev) leptoquarks couple only to quark-lepton pairs and light diquarks couple only to quark pairs. In this paper we present two specific examples where the higgs induced mixing between leptoquarks and diquarks through trilinear interaction terms reintroduces the troublesome couplings and gives rise to proton decay. The bound on the unknown parameters of this scenario that arise from proton life time has also been derived.

Leptoquarks (LQ)  and diquarks (DQ)  are colored scalar or vector particles that carry baryon numbers of $`\pm \frac{1}{3}`$ and $`\pm \frac{2}{3}`$ respectively. They occur naturally in many extensions of the standard model (SM) e.g. superstring inspired grand unified models based on $`E(6)`$ , technicolor models  and composite models  of quarks and leptons. If LQs and DQs couple both to quark-lepton pairs and quark pairs then they lead to too rapid proton decay. To be consistent with the proton life time such LQs and DQs must have a mass of the order of $`(10^{12}10^{15})`$ Gev. To avoid this strong constraint and to make them relevant for low energy phenomenology below 1 Tev it is usually assumed that LQs couple to quark-lepton pairs but not to quark pairs and DQs couple to quark pairs but not to quark-lepton pairs. However in the presence of the SM higgs doublet $`\varphi `$ the low energy effective Lagrangian will also contain a trilinear Higgs-LQ-DQ interaction term. After EW symmetry breaking (EWSB) this interaction term will induce a mixing between the LQ and the DQ which reintroduces the troublesome Yukawa couplings for LQ and DQ. If the physical LQ and DQ do not have the same mass then they lead to proton decay. There can be several different kinds \[1-3\] of LQ’s and DQ’s. However in this report we shall not go into a detailed discussion of every LQ-DQ mixing that can possibly occur. Rather we shall give two give specific examples of Higgs-LQ-DQ interaction term and the mixing between the LQ and DQ that takes place after EWSB. For one particular case we have also calculated the mixing angle in terms of the dimensional coupling that measures the strength of the Higgs-LQ-DQ interaction term. Finally we have estimated the proton decay rate in terms of the unknown parameters of our scenario and have derived the constraint that the parameters of the Lagrangian must satisfy in order to be consistent with proton lifetime.

Consider the SM to be extended by a light ($``$ 250 Gev ) chiral leptoquark D and a light chiral diquark S with the following assignments under $`SU(3)_c\times SU(2)_l\times U(1)_y`$: $`D(3^{},2,\frac{1}{6})`$ and $`S(3^{},1,\frac{1}{3})`$. The low energy effective Lagrangian of this extended scenario will contain besides the SM Lagranigian all possible renormalizable and gauge invariant interaction terms between D, S and the SM fields. Of particular importance for this work are the Yukawa like couplings of D and S to the SM fermions given by the following Lagrangian

$$\begin{array}{ccc}\hfill L_y& =g_R\overline{l}_Ld_RD+g^{}ϵ_{ijk}\overline{d}_{Ri}u_{Rj}^cS_k+h.c.\hfill & \\ & =g_R(\overline{\nu }_LD_1+\overline{e}_LD_2)d_R+g_R^{}ϵ_{ijk}\overline{d}_{Ri}u_{Rj}^cS_k+h.c.\hfill & (1)\hfill \end{array}$$

Here $`D_1`$ and $`D_2`$ are the isospin up and down components of D. i, j, k are the color indices. Note first that the leptoquark D couples only to quark-lepton pair and the diquark S couples only to quark pair which is required so that they do not lead to proton decay. Second both interaction terms are invariant under the SM gauge group $`SU(3)_c\times SU(2)_l\times U(1)_y`$. Finally the couplings of D and S are both chiral in nature i.e. they couple to quark fields of a particular chirality only. The Yukawa Lagrangian $`L_y`$ also implies that D and S should carry baryon numbers of $`\frac{1}{3}`$ and $`\frac{2}{3}`$ respectively so that $`L_y`$ conserves baryon number. This kind of structure of the Yukawa couplings could arise depending on the gauge symmetry of the high energy theory and its chiral matter representation. Besides the Yukawa couplings the low energy Lagrangian for describing physics much below the compositeness scale or grand unified scale $`\mathrm{\Lambda }`$ must also contain all possible renormalizable and gauge invariant interaction terms between $`\varphi `$, D and S. In principle such interaction terms cannot be neglected. They could arise from the scalar potential of the high energy theory after the gauge symmetry of the high energy theory breaks down into $`SU(3)_c\times SU(2)_l\times U(1)_y`$. The gauge invariant trilinear interaction term between $`\varphi `$, D and S is given by the Lagrangian

$$L_s=k_1(D^+\varphi _c)S+h.c.=k_1\frac{v+h}{\sqrt{2}}D_1^+S+h.c.$$
$`(2)`$
Here $`\varphi _c`$ is the charge conjugated Higgs doublet that gives mass to the up quarks in SM. The unknown mixing parameter $`k_1`$ carries the dimension of mass. After EW symmetry breaking this interaction term will lead to mixing between $`D_1`$ and S. Note that $`D_1`$ and S carries the same charge and color assignments which is necessary so that they could mix after EWSB. The above higgs-LQ-DQ interaction term violates both baryon number and lepton number by one unit ($`\delta B=\delta L=1`$). Here we are implicitly assuming that the breaking of baryon number symmetry is first communicated to the scalar sector of the high energy theory which in turn communicates it to other sectors e.g the Yukawa sector. Recall that in the SM baryon number conservation is realized as an accidental global symmetry . By that we mean that given the particle content and the gauge group $`SU(3)_c\times SU(2)_l\times U(1)_y`$ of the SM it is not possible to write a renormalizable and gauge invariant interaction term that violates baryon number. No adhoc global symmetry is required to explain the near abscence of proton decay. However as we have seen above that if we keep the gauge group the same but allow additional particles like color triplet scalars then it is possible to write down a renormalizable and gauge invariant interaction term that violates baryon number.

We shall now show that the mixing between $`D_1`$ and S induced by $`\varphi `$ can lead to proton decay in the “$`\nu `$ \+ any channel” where any refers to a positively charged non-strange meson. Consider the full scalar potential involving D and S

$$\begin{array}{ccc}\hfill V(D,S)& =\mu _1^2D^+D+\mu _2^2S^+S+\lambda _1(D^+D)^2+\lambda _2(S^+S)^2\hfill & \\ & +\lambda _1^{}(D^+D)(\varphi ^+\varphi )+\lambda _2^{}(S^+S)(\varphi ^+\varphi )+(k_1D^+\varphi _cS+h.c.)\hfill & (3)\hfill \end{array}$$

The quartic scalar interactions are not relevant to our present work. $`\mu _1^2`$ and $`\mu _2^2`$ are the mass parameters associated with the gauge eigenstates D and S. After EWSB the trilinear interaction term between $`\varphi `$, D and S induces a mixing between $`D_1`$ and S. It can be shown that the eigenvalues of the mass squared matrix $`M_{D_1S}^2`$ are given by

$$\lambda _+=M_{D_1}^2=\frac{1}{2}(\mu _1^2+\mu _2^2)+\frac{1}{2}\sqrt{ϵ^2+2k_1^2v^2}.$$
$`(4a)`$
and

$$\lambda _{}=M_S^2=\frac{1}{2}(\mu _1^2+\mu _2^2)\frac{1}{2}\sqrt{ϵ^2+2k_1^2v^2}.$$
$`(4b)`$
where $`ϵ=\mu _1^2\mu _2^2`$. The physical states correponding to $`\lambda _\pm `$ are given by $`D_1^{}=D_1\mathrm{cos}\theta S\mathrm{sin}\theta `$ and $`S^{}=D_1\mathrm{sin}\theta +S\mathrm{cos}\theta `$. Here $`\mathrm{cos}\theta =\frac{\sqrt{2}k_1v}{[2k_1^2v^2+(\sqrt{ϵ^2+2k_1^2v^2}ϵ)^2]^{\frac{1}{2}}}`$.

The Yukawa couplings of the LQ and the DQ written in terms of mass eigenstates are given by

$$\begin{array}{ccc}\hfill L_y& =g_R[\overline{\nu }_L(D_1^{}\mathrm{cos}\theta +S^{}\mathrm{sin}\theta +\overline{e}_LD_2]d_R\hfill & \\ & +ϵ_{ijk}g_R^{}\overline{d}_{Ri}u_{Rj}^c[S_k^{}\mathrm{cos}\theta D_{1k}^{}\mathrm{sin}\theta +h.c.]\hfill & (5)\hfill \end{array}$$

. The above Yukawa Lagrangian does lead to proton decay in the “$`\nu `$ \+ any” channel via the exchange of D and S particles. The effective four fermion Lagrangian for this decay is given by

$$L_{eff}=g_Rg_R^{}\mathrm{sin}\theta \mathrm{cos}\theta ϵ_{ijk}(\overline{\nu }_Ld_{Rk})(\overline{u}_{Rj}^cd_{Ri})(\frac{1}{M_S^2}\frac{1}{M_D^2})+h.c.$$
$`(6)`$
This is the effective Lagrangian at a scale $`M^2=M_D^2M_S^2(250Gev)^2`$. In order to use it for proton decay it has to be renormalized down to a scale $`\mu ^2=(1Gev)^2`$ under the unbroken QCD and EM interactions. The EM corrections are small because the coupling itself is small. The QCD corrections are not that large either because $`\mathrm{ln}\frac{M^2}{\mu ^2}`$ is not large in our case. It can be shown that the proton decay rate arising from the above effective Lagrangian is given by

$$\mathrm{\Gamma }(p\nu _e+any)\frac{1}{17\pi ^2R^3}G_{eff}^2m_q^2A^2(\frac{M}{\mu })$$
$`(7)`$
.

Here $`A(\frac{M}{\mu })`$ includes the effects of renormalization group evolution from M to $`\mu `$. We shall assume it to be of order one. $`R=\frac{3}{4}`$ fm, $`m_q\frac{1}{3}`$ Gev and $`G_{eff}=\frac{1}{2\sqrt{2}}g_Rg_R^{}\mathrm{sin}\theta \mathrm{cos}\theta (\frac{1}{M_S^2}\frac{1}{M_D^2})`$. The present lower bound on $`\tau (p\nu +any)`$ is 25$`\times 10^{30}`$ yr . It is clear from the expression of $`G_{eff}`$ that the tight constraint of proton lifetime can be satisfied either with $`\mathrm{sin}\theta 1`$ or $`\mathrm{cos}\theta 1`$. It can be shown that for $`\mathrm{sin}\theta 1`$ we require $`k_1v\sqrt{2}ϵ`$. On the other hand for $`\mathrm{cos}\theta 1`$ we need $`\sqrt{2}k_1v(\sqrt{ϵ^2+2k_1^2v^2}ϵ)`$. For $`M_S=200`$ Gev and $`M_D=300`$ Gev  the product combination $`g_Rg_R^{}\mathrm{sin}\theta \mathrm{cos}\theta `$ must be less than $`3\times 10^{26}`$ in order that the decay rate in eqn(7) is consistent with the proton lifetime. If we assume that the Yukawa couplings $`g_R`$ and $`g_R^{}`$ are of the order of .1, then for $`\mathrm{sin}\theta 1`$ we must adjust the parameters of the Lagrangian so as to satisfy $`k_1v<4\times 10^{24}ϵ`$. On the other hand for $`\mathrm{cos}\theta 1`$ the parameters must be adjusted so as to satisfy $`\sqrt{2}k_1v<3\times 10^{24}(\sqrt{ϵ^2+2k_1^2v^2}ϵ)`$. These constraints may not be stable against radiative corrections requiring extreme fine tuning to each order in loop expansion. We would like to reiterate that there can be several different species of LQ’s as well as DQ’s. However the mixing between $`D_1`$ and S and hence proton decay occurs only for specific $`SU(3)_c\times SU(2)_l\times U(1)_y`$ assignments of the LQ and the DQ. Unless these assignments are achieved there is no mixing between the LQ and DQ and hence no contribution to proton decay.

Having shown that higgs induced mixing between the doublet leptoquark D and the singlet diquark S does lead to proton decay we shall now present another concrete example where the same phenomenon also takes place. An EW triplet diquark $`T_{ak}`$ (a is the $`SU(2)_l`$ index and k the color index) can also mix with the doublet leptoquark D and contribute to both proton decay and neutron decay. Let the $`SU(3)_c\times SU(2)_l\times U(1)_y`$ assignments of $`T_{ak}`$ be given by (3, 3, $`\frac{1}{3}`$). The triplet diquark $`T_{ak}`$ can couple to LH quark pair according to the following Lagrangian:

$$\begin{array}{ccc}\hfill L_y^{}& =g_R^{\prime \prime }ϵ_{ijk}\overline{q}_{Li}^c\tau _ai\tau _2T_{ak}+h.c.\hfill & \\ & =g_R^{\prime \prime }ϵ_{ijk}[\sqrt{2}\overline{u}_{Li}^cu_{Lj}T_k+\sqrt{2}\overline{d}_{Li}^cd_{Lj}T_{+k}+(\overline{u}_{Li}^cd_{Lj}+\overline{d}_{Li}^cu_{Lj})T_{0k}]\hfill & (6)\hfill \end{array}$$

where $`T_{+k}=\frac{T_{1k}+iT_{2k}}{\sqrt{2}}`$, $`T_k=\frac{T_{1k}iT_{2k}}{\sqrt{2}}`$ and $`T_{0k}=T_{3k}`$. $`T_{+k}`$ and $`T_k`$ carry electromagnetic charges of $`\frac{2}{3}`$ and $`\frac{4}{3}`$ units respectively. Note that $`T_k`$ is not the antiparticle of $`T_{+k}`$ since they carry different electromagnetic charges. The gauge invariant Higgs-LQ-DQ interaction term in this case will be given by

$$\begin{array}{ccc}\hfill L_s^{}& =k_2(D^+\tau _a\varphi _c)T_a^{}+h.c.\hfill & \\ & =k_2\frac{v+h}{\sqrt{2}}[D_2^+T_+^{}+D_1^+T_0^{}+h.c.]\hfill & (7)\hfill \end{array}$$

The mixing between $`D_1`$ and $`T_0^{}`$leads to $`p\pi ^+\nu `$ and the mixing between $`D_2`$ and $`T_+^{}`$ leads to $`n\pi ^+e^{}`$ both of which violate baryon number.

The low energy effective Lagrangian will also contain trilinear Higgs-LQ-LQ and Higgs-DQ-DQ interaction terms. Afetr EWSB these terms give rise to mixing between different multiplets of LQ and DQ. Such interaction terms do not violate baryon number but they could violate lepton number. Of particular importance is the Higgs-LQ-LQ ineraction for a pair of chiral leptoquarks belonging to different weak SU(2) multiplets. The mixing between two different LQ multiplets will give rise to helicity unsuppressed contribution to $`\pi ^{}e^{}\overline{\nu }_e`$ The helicity suppressed SM contribution to $`\pi ^{}e^{}\overline{\nu }_e`$ is in excellent agreement with the experimental data. Therefore any helicity unsuppressed contribution arising from new physics must be strongly constrained. This leads to stringent bounds on the mixing parameter between the two LQ multiplets. The mixing between different LQ multiplets can also lead to majorana mass matrix for neutrinos. A complete discussion of these topics can be found in Ref.  and will not be taken up here.

To conclude in this report we have shown that the usual assumption of low energy phenomenology that LQ’s couple only to q-l pairs and DQ’s couple only to quark pairs is not sufficient to stabilize the proton. We have given two specific examples where the Higgs-LQ-DQ interaction term induces a mixing between the LQ and DQ after EWSB. This mixing reintroduces the troublesome couplings for the LQ and the DQ and leads to proton decay. We find that in order to be consistent with the bound on proton lifetime the parameters of the Lagrangian must either satisfy $`k_1v<4\times 10^{24}ϵ`$ or $`\sqrt{2}k_1v<3\times 10^{24}(\sqrt{ϵ^2+2k_1^2v^2}ϵ)`$. In deriving these constraints we have assumed that $`g_Rg_R^{}.1`$, $`M_{D_1}=300`$ Gev and $`M_S`$=200 Gev. These constraints may not be stable against radiative corrections requiring extreme fine tuning of the parameters to each order in loop expansion. This naturalness problem can perhaps be resloved by means of some new symmetry that remains unbroken in the low energy theory and prevents the Higgs-LQ-DQ interaction term from occuring in the effective Lagrangian.

References

1. W. Buchmuller Acta Phys. Austriaca, Suppl 27, 517 (1985); W. Buchmuller and D. Wyler, Phys. Lett. B 177, 377 (1986).

2. E. Farhi and L. Susskind, Phys. Rep. 74, 277 (1981); K. Lane and M. V. Ramanna, Phys. Rev. D 44, 2678 (1991).

3. J. Hewett and T. Rizzo, Phys. Rep. 183, 193 (1989).

4. A. Cohen, A. Nelson and D. Kaplan, Phys. Lett. B 388, 588 (1996).

5. Review of Particle Properties Euro. Phy. Jour. C3, 50 (1998).

6. The values of $`M_S`$ and $`M_D`$ assumed by us are consistent with the latest bounds on LQ and DQ masses: Euro. Phys. Jour. C3, 260 (1998); B. Abbot et al. (DO Collaboration), Phys. Rev. Lett. 80, 2051 (1998); F. Abe et al. (CDF Collaboration), Phys. Rev. D 55, R5263 (1997).

7. M. Hirsch et al., Phys. Rev. D 54, 4207(1996); Phys. Lett. B 378, 17 (1996); U. Mahanta, MRI-PHY/990930, hep-ph/9909518.


format: text

text:
# Angular Momentum Transport in Particle and Fluid Disks

## 1 Introduction

When turbulence in unmagnetized Keplerian disks is driven by thermal convection, angular momentum appears to be transported inwards (Ryu & Goodman 1992; Kley, Papaloizou, & Lin 1993; Stone & Balbus 1996; Cabot 1996). If generically true, this implies that convection cannot provide the “enhanced” angular momentum transport (relative to collisional viscosity) needed to account for observations of accreting systems in astrophysics.

Despite appearances, inward transport does not contradict the fact that an isolated Keplerian disk relaxes to its minimum energy state by transporting angular momentum outwards (Lynden-Bell & Pringle 1974). The reason is that simulated convective disks are not completely isolated. They tap a source of free energy other than their intrinsic differential rotation: an externally imposed, superadiabatic entropy gradient.

A more mechanistic understanding of the behavior of turbulent fluid disks can be obtained by considering disks composed of discrete macroscopic particles. Transport properties and stability criteria obtained under one picture often carry over to the other, to within dimensionless factors of order unity. The fluid/particle dualism has been employed with considerable success in studies of planetary rings (see, e.g., Goldreich & Tremaine 1982) and galactic structure (see, e.g., Chapter 6 of Binney & Tremaine 1987).

Indeed, one-to-one correspondences can be established between terms in the Reynolds stress tensor equation describing turbulent fluid disks, and terms in the Boltzmann stress tensor equation describing particle disks. This paper explores the transport properties of the latter system, with the aim of shedding light on the former. We demonstrate how simple particle disk models can reproduce many of the transport phenomena found in numerical simulations of vertically or radially convective fluid disks. The particle viewpoint also provides insight into the stability properties of isolated, unmagnetized, fluid disks. An analysis similar to ours was carried out by Kumar, Narayan, & Loeb (1995); this is discussed in §5.

In §2, we write down the Boltzmann stress tensor equation governing particle disks and discuss its intimate relationship with the Reynolds stress tensor equation governing hydrodynamic turbulence in a fluid disk. In §3, we prove a general theorem that relates the direction of angular momentum transport to the behavior of the velocity dispersion tensor in the presence of external driving. Section 4 is the most important of this paper; there we consider a simple model for the collision term in the Boltzmann equation that allows one to easily solve for the stress tensor of the disk and, we believe, capture much of the relevant physics regarding both stability and the direction of transport. In §5, we summarize and discuss our results.

## 2 Stress Tensor Equations for Particle Disks

The second moments of the Boltzmann equation describe the evolution of the viscous stress tensor. For an axisymmetric disk, in cylindrical coordinates, these read (see, e.g., Borderies, Goldreich, & Tremaine 1983):

$$_tp_{rr}4\mathrm{\Omega }p_{r\varphi }=\left(_tp_{rr}\right)_c+\dot{H}_{rr}$$
(1)
$$_tp_{\varphi \varphi }+\frac{\kappa ^2}{\mathrm{\Omega }}p_{r\varphi }=\left(_tp_{\varphi \varphi }\right)_c+\dot{H}_{\varphi \varphi }$$
(2)
$$_tp_{zz}=\left(_tp_{zz}\right)_c+\dot{H}_{zz}$$
(3)
$$_tp_{r\varphi }+\frac{\kappa ^2}{2\mathrm{\Omega }}p_{rr}2\mathrm{\Omega }p_{\varphi \varphi }=\left(_tp_{r\varphi }\right)_c,$$
(4)
where we have neglected third moments. In equations (1)–(4), $`\mathrm{\Omega }`$ is the angular frequency, $`\kappa ^2=r^3d(\mathrm{\Omega }r^2)^2/dr`$ is the square of the epicyclic frequency, and $`p_{ij}`$ are the components of the stress tensor. The terms $`\left(_tp_{ij}\right)_c`$ describe changes in velocity dispersions due to interparticle collisions. Collisions perform three tasks: they (1) extract energy from the mean shear into random motions, (2) transfer energy in random motions associated with one direction to those associated with another, and (3) convert kinetic energy into heat. The terms $`\dot{H}_{ij}`$ are inserted to simulate external stirring of the disk, with $`\dot{H}_{r\varphi }=0`$ so as not to explicitly induce angular momentum transport.

### 2.1 Boltzmann vs. Reynolds

Since the fluid equations are derived from the Boltzmann equation it is not surprising that hydrodynamic turbulence in fluid disks satisfies energy conservation equations very similar to equations (1)–(4). The Reynolds stress equation can be obtained by taking the Navier-Stokes equation, decomposing the velocity and pressure into mean and fluctuating components, and multiplying the entire equation by the velocity. This yields (e.g., Speziale 1991; Kato & Yoshizawa 1997):

$$_tT_{ij}=T_{ik}_kU_jT_{jk}_kU_i+\pi _{ij}ϵ_{ij}_kC_{ijk},$$
(5)
where $`U`$ is the mean velocity of the fluid, $`u`$ is the fluctuating velocity, $`T_{ij}=u_iu_j`$ is the stress tensor, and $``$ denotes an ensemble average. The quantity $`C_{ijk}`$ is a sum of turbulent and viscous fluxes,

$$\pi _{ij}=u_i_jp+u_j_ip$$
(6)
and

$$ϵ_{ij}=2\nu _ku_i_ku_j,$$
(7)
where $`p`$ is the fluctuating gas pressure and $`\nu `$ is the kinematic viscosity. The first two terms on the right hand side of equation (5) describe how turbulence extracts energy from the shear in the mean flow. The terms $`\pi _{ij}`$ are known as the “pressure-strain” correlation in the fluid dynamics literature;444To be precise, our $`\pi _{ij}`$ is related to the usual pressure-strain correlation tensor by a total divergence, which we have absorbed into $`C_{ijk}`$. the diagonal terms of $`\pi _{ij}`$ are simply the work done by fluctuating pressure forces associated with the turbulence; among their other tasks, these terms attempt to enforce isotropy of the turbulence by transferring energy from velocity fluctuations associated with one direction to those associated with another. The terms $`ϵ_{ij}`$ represent the viscous dissipation of energy on small scales (note that $`ϵ_{ii}`$ is positive definite).

It is straightforward to evaluate the above expression for $`_tT_{ij}`$ in cylindrical coordinates. We neglect the flux divergence555This assumes “homogeneous” turbulence, i.e., that the flux varies on a length scale $`r`$, much greater than the vertical scale height of the disk. to obtain

$$_tT_{rr}4\mathrm{\Omega }T_{r\varphi }=\pi _{rr}ϵ_{rr}$$
(8)
$$_tT_{\varphi \varphi }+\frac{\kappa ^2}{\mathrm{\Omega }}T_{r\varphi }=\pi _{\varphi \varphi }ϵ_{\varphi \varphi }$$
(9)
$$_tT_{zz}=\pi _{zz}ϵ_{zz}$$
(10)
$$_tT_{r\varphi }+\frac{\kappa ^2}{2\mathrm{\Omega }}T_{rr}2\mathrm{\Omega }T_{\varphi \varphi }=\pi _{r\varphi }ϵ_{r\varphi }.$$
(11)
Equations (8)–(10) are identical to the energy equations given by, e.g., Balbus & Hawley (1998; see their equations 58–60).

Equations (1)–(4) for the evolution of the stress tensor in a particle disk are very similar to equations (8)–(11) for the evolution of the stress tensor in a fluid disk. In particular, the rotational profile of the disk enters identically in the two cases; both $`\mathrm{\Omega }`$ and $`\kappa ^2`$ couple to the stress tensor components in exactly the same way for particle and fluid disks. Moreover, the “source” and “sink” terms in the energy equations are analogous; the collision term, $`\left(_tp_{ij}\right)_c`$, in the Boltzmann equation plays the same role as $`\pi _{ij}ϵ_{ij}`$ in the Reynolds stress equations.

This comparison indicates that the dynamical role of epicycles is the same in both fluid and particle disks, and that any difference in the angular momentum transport properties of the two disks can only be due to differences in the sources and sinks in the energy equation. We demonstrate this formally in the next section by proving a theorem that relates the direction of angular momentum transport to the properties of the collisional and stirring terms. This theorem is very useful for interpreting the model problem of §4.

## 3 The “Long Axis Must Get Longer” Theorem

Theorem: For thin axisymmetric disks in which $`d\mathrm{\Omega }/dr<0`$, steady inward transport of angular momentum requires that the net tendency of collisions and external stirring be to further lengthen the principal equatorial long axis of the velocity ellipsoid and to further shorten the principal equatorial short axis. This theorem applies to fluid disks as well as to particle disks so long as the interparticle collision term, $`(_tp_{ij})_c`$, is replaced by $`\pi _{ij}ϵ_{ij}`$.

Proof: We transform equations (1)–(4) to the principal axis basis of the velocity ellipsoid. Let $`p_{11}`$ and $`p_{22}`$ be the greater and lesser of the two components of the stress tensor, respectively, in the equatorial plane:

$`p_{11}`$ $`=`$ $`[p_{rr}+p_{\varphi \varphi }+\sqrt{(p_{rr}p_{\varphi \varphi })^2+4p_{r\varphi }^2}]/\mathrm{\hspace{0.17em}2}`$ (12)
$`p_{22}`$ $`=`$ $`[p_{rr}+p_{\varphi \varphi }\sqrt{(p_{rr}p_{\varphi \varphi })^2+4p_{r\varphi }^2}]/\mathrm{\hspace{0.17em}2},`$ (13)
where we shall always take the positive square root. By symmetry about the midplane, $`p_{33}=p_{zz}`$. The orientation of the velocity ellipsoid is given by $`\delta `$, the angle between the longer equatorial axis and the radial direction:

$`\mathrm{tan}\delta `$ $``$ $`{\displaystyle \frac{p_{r\varphi }}{p_{11}p_{\varphi \varphi }}}`$ (14)
$`=`$ $`{\displaystyle \frac{2p_{r\varphi }}{(p_{rr}p_{\varphi \varphi })+\sqrt{(p_{rr}p_{\varphi \varphi })^2+4p_{r\varphi }^2}}},`$
valid for $`|\delta |<\pi /2`$ and nonisotropic $`p_{ij}`$.666If $`\delta =\pm \pi /2`$ or $`p_{ij}`$ is isotropic, $`p_{r\varphi }=0`$ and there is no transport. Since the denominator in equation (14) is positive,

$$\text{sgn}p_{r\varphi }=\text{sgn}\delta .$$
(15)
Outward (inward) transport corresponds to positive (negative) $`p_{r\varphi }`$, which in turn corresponds to positive (negative) $`\delta `$.

In the principal axis basis, the viscous stress equations (1)–(4) read (cf. Goldreich & Tremaine 1978):

$$_tp_{11}S\mathrm{\Omega }(\mathrm{sin}2\delta )p_{11}=\left(_tp_{11}\right)_c+\dot{H}_{11}$$
(16)
$$_tp_{22}+S\mathrm{\Omega }(\mathrm{sin}2\delta )p_{22}=\left(_tp_{22}\right)_c+\dot{H}_{22}$$
(17)
$$_tp_{33}=\left(_tp_{33}\right)_c+\dot{H}_{33}$$
(18)
$$(p_{11}p_{22})_t\delta (2S\mathrm{sin}^2\delta )\mathrm{\Omega }p_{22}+(2S\mathrm{cos}^2\delta )\mathrm{\Omega }p_{11}=(p_{11}p_{22})(_t\delta )_c+\dot{H}_{12},$$
(19)
where $`S=d\mathrm{ln}\mathrm{\Omega }/d\mathrm{ln}r`$.

Equations (16) and (17) state that for $`S>0`$, the steady ($`_t=0`$) inward transport ($`\delta <0`$) of angular momentum demands that the net tendency of collisions and external stirring be to further lengthen the principal long axis of the velocity ellipsoid and to further shorten the principal short axis.

Collisions between hard, imperfectly elastic spheres tend to isotropize the velocity ellipsoid: $`(_tp_{11})_c<0<(_tp_{22})_c`$. Thus, in the absence of stirring, transport is normally outwards in axisymmetric particle disks.

## 4 Angular Momentum Transport in Externally Stirred Disks

In this section we solve for the steady state structure of a particle disk using a simple model for the effects of interparticle collisions and external stirring. We also discuss the fluid analogue of our particle disk model in some detail.

### 4.1 The Particle Scattering Function

Following, e.g., Cook & Franklin (1964), Shu & Stewart (1985), and Narayan, Loeb, & Kumar (1994) we take the scattering function in the Boltzmann equation to be given by the Krook approximation:

$$(_tf)_c=\frac{f_0f}{\tau },$$
(20)
where $`\tau `$ is the collisional mean free time, $`f`$ is the particle distribution function, and $`f_0`$ is the “equilibrium” or “post-scattering” distribution function. The corresponding collision term in the viscous stress equations is $`\left(_tp_{ij}\right)_c=(p_{ij}^0p_{ij})/\tau `$. The physical idea behind the Krook model is that when particles collide they are removed from phase space at a rate $`f/\tau `$ and are returned to phase space with the new distribution function $`f_0`$.

We assume that collisions do not generate off diagonal stresses so that $`p_{ij}^0=0`$ if $`ij`$. We account for the fact that particle collisions are inelastic, converting some fraction of the kinetic energy of the colliding particles into heat. With this model, the relevant second moments of $`f_0`$ are given by (see Kumar et al. 1995)

$$p_{rr}^0=(1\xi _r)\frac{\sigma ^2}{3},p_{\varphi \varphi }^0=(1\xi _\varphi )\frac{\sigma ^2}{3},\mathrm{and}p_{zz}^0=(1\xi _z)\frac{\sigma ^2}{3},$$
(21)
where $`\sigma `$ is the velocity dispersion ($`\sigma ^2=p_{rr}+p_{\varphi \varphi }+p_{zz}`$) and the parameters $`\xi _iϵ[0,1]`$ measure the inelasticity of the collisions; $`\xi 0`$ is elastic, while $`\xi 1`$ is highly inelastic.

Equation (21) allows the inelasticity of colliding particles to be anisotropic, i.e., to be different for the $`r`$, $`\varphi `$, and $`z`$ directions. Such an anisotropic scattering function is unrealistic for real particle disks in nature (e.g., planetary rings); collisions between hard spheres at the microscopic level cannot be expected to respect the cylindrical coordinate system. Our strategy is to first solve for the structure of “real” externally stirred particle disks, in which the $`\xi _i`$’s are equal. This captures much of the physics in which we are interested. We will show, however, that the assumption of equal $`\xi _i`$’s does not furnish a fully adequate model of fluid disks and that the introduction of an “anisotropic inelasticity” can remedy this defect.

### 4.2 The Fluid Analogue of the Scattering Function

Our particle scattering function has a mathematically rigorous analogue in fluid turbulence. Namely, it is equivalent to a particular closure scheme for the a priori unknown turbulent quantities $`ϵ_{ij}`$ and $`\pi _{ij}`$:

$$ϵ_{ij}=\frac{T_{ij}}{\tau }\mathrm{and}\pi _{ij}=\frac{(1\xi _i)}{3\tau }\sigma ^2\delta _{ij}.$$
(22)
Interpreted as such, $`\tau `$ is roughly the timescale on which free turbulence would decay (the eddy turnover time) and the inelasticity parameter, $`\xi _i`$, measures the efficiency with which fluctuating pressure forces perform work in the $`i`$th direction. The closure model proposed in equation (22) differs from conventional second-order closure models of the Reynolds stress equation (e.g., Speziale 1991; see Kato & Yoshizawa 1997 for an application to Keplerian disks). We will show, however, that this simple model provides considerable insight into the angular momentum transport properties of accretion disks.

In contrast to $`\tau `$, it is not obvious what values to assign $`\xi _i`$ in hydrodynamic turbulence. Although in one sense turbulence is “inelastic,” with energy cascading to small scales in the absence of source terms, the physical meaning of $`\xi _i`$ defined by equation (22) is more precise: fluid turbulence is “inelastic” if and only if fluctuating pressure forces are inefficient at doing work. Throughout the remainder of this paper it is worth keeping in mind that numerical simulations suggest that the highly “inelastic” limit is appropriate for Keplerian disks. Note also that, unlike in the particle disk case, we would not expect the $`\xi _i`$’s to be identical in the fluid case. The fluctuating pressure forces are macroscopic in origin and thus can plausibly be different for the $`r`$, $`\varphi `$, and $`z`$ directions.

### 4.3 The Stress Tensor Equations

With our model for the interparticle collision term, the steady state equations for the stress tensor become

$$p_{rr}(1\xi _r)\frac{\sigma ^2}{3}4\mathrm{\Omega }\tau p_{r\varphi }=H_{rr},$$
(23)
$$p_{\varphi \varphi }(1\xi _\varphi )\frac{\sigma ^2}{3}+\frac{\kappa ^2\tau }{\mathrm{\Omega }}p_{r\varphi }=H_{\varphi \varphi },$$
(24)
$$p_{zz}(1\xi _z)\frac{\sigma ^2}{3}=H_{zz},$$
(25)
and

$$\frac{\kappa ^2\tau }{2\mathrm{\Omega }}p_{rr}2\mathrm{\Omega }\tau p_{\varphi \varphi }+p_{r\varphi }=0,$$
(26)
where $`H_{ij}\tau \dot{H}_{ij}`$, i.e., the net stirring over the mean free time.

An “energy equation” can be derived by summing equations (23)–(25):

$$HH_{rr}+H_{\varphi \varphi }+H_{zz}=\xi \sigma ^2+2p_{r\varphi }\mathrm{\Omega }^{}\tau ,$$
(27)
where $`\mathrm{\Omega }^{}=d\mathrm{\Omega }/d\mathrm{ln}r`$ and $`3\xi =\xi _r+\xi _\varphi +\xi _z`$. For $`\mathrm{\Omega }^{}<0`$, equation (27) shows that, if there is no external stirring, transport must be outwards or zero in steady state. This is because any energy lost due to inelasticity must be extracted from the shear in the mean flow, which requires $`p_{r\varphi }>0`$. No such thermodynamic constraint on the direction of angular momentum transport exists if there is external stirring.

### 4.4 Transport in a Keplerian Disk

The general solution to equations (23)–(26) is given in the Appendix. Here we focus on the case of a Keplerian disk, for which $`\kappa =\mathrm{\Omega }`$; in the next subsection we briefly contrast the Keplerian disk with a constant angular momentum disk. We define $`x=\mathrm{\Omega }\tau `$ and solve for the $`r\varphi `$ component of the stress tensor:

$$p_{r\varphi }=\frac{x\sigma ^2}{2}\left(\frac{\frac{4}{3}\left[1\xi _\varphi \right]\frac{1}{3}\left[1\xi _r\right]\stackrel{~}{H}_{rr}\xi +4\stackrel{~}{H}_{\varphi \varphi }\xi }{1+x^2\left[2.5\stackrel{~}{H}_{rr}+10\stackrel{~}{H}_{\varphi \varphi }+4\stackrel{~}{H}_{zz}\right]}\right),$$
(28)
where $`\stackrel{~}{H}_{ii}=H_{ii}/H`$. The velocity dispersion of the disk is not an independent parameter but is determined by solving the equation

$$\sigma ^2=\frac{H+x^2\left(2.5H_{rr}+10H_{\varphi \varphi }+4H_{zz}\right)}{(1+4x^2)\xi 0.5\xi _rx^2+2\xi _\varphi x^21.5x^2}.$$
(29)
#### 4.4.1 Real Particle Disks: Isotropic Inelasticity

For a real particle disk with isotropic inelasticity ($`\xi _r=\xi _\varphi =\xi _z=\xi `$) equations (28) and (29) become

$$p_{r\varphi }=\frac{x\sigma ^2}{2}\left(\frac{1\xi \stackrel{~}{H}_{rr}\xi +4\stackrel{~}{H}_{\varphi \varphi }\xi }{1+x^2\left[2.5\stackrel{~}{H}_{rr}+10\stackrel{~}{H}_{\varphi \varphi }+4\stackrel{~}{H}_{zz}\right]}\right)$$
(30)
and

$$\sigma ^2=\frac{H+x^2\left(2.5H_{rr}+10H_{\varphi \varphi }+4H_{zz}\right)}{(1+5.5x^2)\xi 1.5x^2}.$$
(31)
In the limit of no external stirring, we recover the usual results of particle disk theory. In such an isolated disk, the inelasticity adjusts to the optical depth of the disk,

$$\xi =\frac{1.5x^2}{1+5.5x^2}.$$
(32)
This follows from equation (31) by setting $`H=0`$ and requiring finite $`\sigma ^2`$. Since the collisional inelasticity depends on the relative impact velocity, equation (32) determines a velocity dispersion (Goldreich & Tremaine 1978); call this $`\sigma _0`$. Inserting equation (32) for $`\xi `$ in equation (30), we confirm that the transport is always outwards, with777To derive eq. , write eq.  as $`p_{r\varphi }=0.5x\sigma ^2N/D`$. Rewritting the denominator as $`D=1+4x^2+x^2(6\stackrel{~}{H}_{\varphi \varphi }1.5\stackrel{~}{H}_{rr})`$ and substituting eq.  for $`\xi `$ into the numerator, $`N`$, all terms $`H_{ij}`$ cancel, yielding eq. .

$$p_{r\varphi }=\frac{\xi \sigma ^2}{3x}=\frac{x\sigma ^2}{2+11x^2}.$$
(33)
The fact that isolated fluid disks appear to be hydrodynamically stable can be explained mechanistically in the sense that they do not satisfy equation (32). For any $`x`$, equation (32) has a solution only if $`\xi ϵ[0,0.27]`$. Simulations suggest that the work done by fluctuating pressure forces is far less efficient than this. Therefore, for $`H=0`$, the only solution to equation (31) is $`\sigma _0^2=0`$, i.e., no turbulence.

What happens if the disk is externally stirred? For $`H\xi \sigma _0^2`$ the disk structure is only weakly affected by the stirring; corrections to $`p_{r\varphi }`$, e.g., are $`O[H/\xi \sigma _0^2]1`$. Physically, this is because, for $`\sigma _00`$, heating must become significant relative to the “natural” velocity dispersion in order to modify the disk structure. For a hydrodynamically stable fluid disk with $`\sigma _0=0`$, however, the “natural” velocity dispersion vanishes; thus, for any $`H`$, the properties of the external stirring (together with the “collisions”) determine the angular momentum transport properties of the disk.

In the presence of significant stirring, i.e., H
>σ02
>𝐻subscriptsuperscript𝜎20H\mathrel{\mathchoice{\vbox{\offinterlineskip\halign{$\m@th\displaystyle\hfil#\hfil$\cr>\crcr\sim\crcr} }}{\vbox{\offinterlineskip\halign{$\m@th\textstyle\hfil#\hfil$\cr>\crcr\sim\crcr} }}{\vbox{\offinterlineskip\halign{$\m@th\scriptstyle\hfil#\hfil$\cr>\crcr\sim\crcr} }}{\vbox{\offinterlineskip\halign{$\m@th\scriptscriptstyle\hfil#\hfil$\cr>\crcr\sim\crcr} }}}\sigma^{2}_{0}, the direction of angular momentum transport can reverse; in particular, the sign of the angular momentum transport is then determined by the sign of the numerator in equation (30). Note that the sign of $`p_{r\varphi }`$ does not explicitly depend on $`x=\tau \mathrm{\Omega }`$. This is consistent with the arguments given in §2 and §3; only the properties of the collisions and the stirring determine the direction of angular momentum transport, not the degree to which the disk is collisional.

From the numerator in equation (30) it follows that (1) energy input along $`r`$ tends to lead to inward transport ($`H_{rr}`$ enters with a negative sign), (2) energy input along $`\varphi `$ tends to lead to outward transport ($`H_{\varphi \varphi }`$ enters with a positive sign), and (3) elastic scattering tends to lead to outward transport ($`p_{r\varphi }>0`$ if $`\xi 1`$). All of these results are readily understood using the theorem of §3. The principal long axis of the velocity ellipsoid generally has a much larger projection along the radial direction than the azimuthal direction.888In the collisionless unstirred limit, equation (19) states that $`p_{rr}/p_{\varphi \varphi }=2/(2S)=4`$ for steady-state Keplerian disks. Thus $`r`$ stirring corresponds to elongating the long axis of the velocity ellipsoid and so should (and does) lead to inward transport (and vice-versa for $`\varphi `$ stirring). The transport only respects the direction of external stirring, however, if the collisions are fairly inelastic. In the elastic limit, collisions efficiently transfer energy from one direction to another, isotropizing the distribution function; transport is then always outwards.

Two special cases of equation (30) are of particular interest:

1. $`R`$ stirring only: This is relevant for interpreting the work of Stone, Pringle, & Begelman (1999) and Igumenshchev & Abramowicz (1999) who carried out simulations of advection-dominated accretion flows, which are radially convective. In these flows, one can treat the convection as “externally” driven by the viscous dissipation associated with magnetic fields. As discussed by Narayan, Igumenshchev, & Abramowicz (2000) and Quataert & Gruzinov (2000), the convective angular momentum transport is inwards in these simulations, with $`\alpha `$ as large as $`0.1`$ ($`\alpha p_{r\varphi }/\sigma ^2`$ is the dimensionless viscosity). As noted above, it is particularly natural for radial convection to produce inward transport since radial convection tends to directly expand the long axis of the velocity ellipsoid. More concretely, equation (30) predicts that the transport should be inwards in the presence of $`r`$ stirring so long as $`\xi >1/2`$. Moreover, for $`\xi 1`$ and $`x1`$,999This means that the eddy turnover time is comparable to the rotational period of the disk, as is indeed the case. $`\alpha `$ should be $`0.1`$ in the presence of $`r`$ stirring, consistent with simulations.
2. $`Z`$ stirring only: This is relevant for interpreting the simulations of Stone & Balbus (1996) in which vertical convection in an otherwise stable Keplerian disk generated low levels of inward angular momentum transport ($`\alpha 10^4`$). For a particle disk, equation (30) shows that transport is always outwards in the presence of $`z`$ stirring: $`p_{r\varphi }(1\xi )`$. This can be understood by noting that the $`z`$ dynamics in thin reflection-symmetric disks is largely decoupled from the $`r`$ and $`\varphi `$ dynamics which determine the direction of transport. In particular, energy input in the $`z`$ direction is transferred to the $`r`$ and $`\varphi `$ directions by particle collisions, which attempt to isotropize the distribution function. By the theorem of §3 the transport must then be outwards.

#### 4.4.2 “Anisotropic Inelasticity”: Inward Transport with Vertical Stirring

The analysis of the previous section applies to real particle disks in which the inelasticity is isotropic. Interestingly, it shows that vertically stirred particle disks would transport angular momentum outwards, opposite to what is found in simulations of hydrodynamic turbulence. As noted in §4.2, the fluid analogue of inelasticity is inefficient fluctuating pressure forces; such forces are macroscopic in origin and therefore need not be the same in the $`r`$, $`\varphi `$, and $`z`$ directions. We now discuss how such “anisotropic inelasticity” can produce inward transport in the presence of vertical stirring.

In the presence of vertical stirring only, equation (28) reduces to

$$p_{r\varphi }=\frac{x\sigma ^2}{2}\left(\frac{\frac{4}{3}\left[1\xi _\varphi \right]\frac{1}{3}\left[1\xi _r\right]}{1+4x^2}\right).$$
(34)
From the numerator in equation (34) it follows that elasticity in the $`r`$ direction tends to lead to inward transport ($`1\xi _r`$ has a negative prefactor), while elasticity in the $`\varphi `$ direction tends to lead to outward transport ($`1\xi _\varphi `$ has a positive prefactor).

Inward transport in the presence of $`z`$ stirring requires $`\xi _\varphi >0.75+0.25\xi _r`$. Physically, this means that fluctuating pressure forces are particularly inefficient at doing work in the $`\varphi `$ direction (cf. Stone & Balbus 1996). In fact, in the absence of explicit $`\varphi `$ stirring, equation (28) shows that the transport is always inwards if $`\xi _\varphi 1`$.101010This is analogous to Stone & Balbus’s argument that axisymmetric simulations of hydrodynamic turbulence will always find inward transport. In order to reproduce $`\alpha 10^4`$ in the presence of $`z`$ stirring, as found in the simulations, equation (34) requires (taking $`x1`$) $`1\xi _r10^3`$, and $`1\xi _\varphi 10^3`$.

### 4.5 Transport in a Constant Angular Momentum Disk

Balbus and collaborators have emphasized the interesting point that a constant angular momentum, $`\kappa =0`$, disk is analogous to a linear shear flow (e.g., Balbus, Hawley, & Stone 1996). Indeed, numerical simulations show that the $`\kappa =0`$ disk is unstable to nonlinear perturbations and produces outwards angular momentum transport. This feature is also readily captured by our model. By equation (A1) the $`r\varphi `$ stress in a $`\kappa =0`$ disk is given by

$$p_{r\varphi }=2x\sigma ^2\left(\frac{\frac{1\xi _\varphi }{3}+\xi \stackrel{~}{H}_{\varphi \varphi }}{1+8x^2\stackrel{~}{H}_{\varphi \varphi }}\right),$$
(35)
which is positive definite. For unstirred $`\kappa =0`$ disks in which the particle inelasticity is isotropic,

$$\xi =\frac{8x^2}{3+8x^2}.$$
(36)
In contrast to equation (32) for a Keplerian disk, equation (36) can be satisfied for some $`x`$ for any $`\xi ϵ[0,1]`$. For a fluid disk this means that, regardless of the efficiency of fluctuating pressure forces, outward transport can be sustained in an unstirred disk. This is consistent with the fact that shear flows are nonlinearly unstable.

## 5 Discussion

This paper explores the angular momentum transport properties of axisymmetric, externally stirred particle disks as a means of better understanding hydrodynamic turbulence in fluid disks. We are not, however, advocating that turbulence is a sort of macroscopic kinetic theory, with colliding turbulent “blobs” literally replacing colliding particles. Rather, we argue that it is physically revealing to exploit the similar mathematics governing turbulent fluid disks and particle disks.

There are one-to-one correspondences between terms in the Boltzmann stress tensor equation governing particle disks and terms in the Reynolds stress tensor equation governing fluid disks (§2). This analogy reveals that the role of differential rotation and epicyclic motion is the same in both fluid and particle disks; both the rotation rate and the epicyclic frequency couple to the components of the stress tensor identically in the two cases. Consequently, the transport properties of fluid disks can differ from those of particle disks only to the extent that the “scattering function” for fluid elements differs from that of particles. The formal, and very useful, analogy is that interparticle collisions in particle disks play the same role as viscous dissipation and fluctuating pressure forces in turbulent disks. Both transfer energy in random motions associated with one direction to those associated with another, and convert kinetic energy into heat.

Motivated by the above considerations, we solve for the structure of a particle disk using a simple model for the collision term in the Boltzmann equation (§4). We allow different collisional inelasticities along the $`r`$, $`\varphi `$, and $`z`$ directions. Although an anisotropic inelasticity is not realistic for actual particle collisions, its introduction serves to make the analogy with fluid turbulence complete. The fluid analogue of particle inelasticity is the inefficacy with which fluctuating pressure forces do work, and there is no a priori reason for the latter to be isotropic. A final ingredient of our model is the insertion of an external stirring mechanism, modeled as constant source terms in the $`r`$, $`\varphi `$, and $`z`$ energy equations. External stirring is relevant for interpreting numerical simulations in which either vertical or radial convection is present (e.g., Stone & Balbus 1996; Stone, Pringle, & Begelman 1999).

We believe that this particle model captures the transport properties of fluid disks. With a single equation for the stress tensor (eq. ) we can reproduce results from standard particle disk theory and from a number of numerical simulations of hydrodynamic turbulence.

In the limit of no external stirring, differential rotation is the only source of free energy in our disks. The angular momentum transport induced by turbulence or particle collisions must then be outwards or zero in steady state. This is because both particle collisions and turbulence generate heat, which can only be compensated for by a positive $`r\varphi `$ stress (see eq. ). To have outward transport in a Keplerian disk requires that “collisions” be relatively efficient at isotropizing the distribution function (see eq.  and the associated discussion). This tendency towards isotropization is normally satisfied in axisymmetric particle disks such as planetary rings; these systems therefore transport momentum outwards while maintaining some degree of random motion (turbulence). By contrast, the hydrodynamic stability of isolated fluid disks can be understood by noting that, according to simulations, fluctuating pressure forces are extremely inefficient at isotropization. Thus there can be no self-sustaining turbulence in these systems.111111Kato & Yoshizawa (1997) have argued that if pressure fluctuations in turbulent Keplerian disks behave like they do in shear flows and other turbulent laboratory fluids, the isotropization is sufficient to yield self-sustaining turbulence. This assumes, however, that Reynold’s stress closure models developed for laboratory problems are applicable to Keplerian disks.

In the presence of significant external stirring, there are no thermodynamic constraints on the direction of angular momentum transport. Instead, the sign of the transport is determined by the nature of the stirring and the properties of the collision term in the Boltzmann equation (or its suitably defined analogue in the fluid problem); formally, a necessary and sufficient condition for inward transport is that the combined tendency of collisions and stirring be to make the long equatorial axis of the velocity ellipsoid yet longer and the short equatorial axis yet shorter (§3).

More physically, for any mechanism of external stirring (e.g., convection, spiral shocks, etc.), our model predicts the direction of transport as a function of the “inelasticity” of collisions (eq. ). The latter is, of course, not known a priori for the fluid problem, which underscores the parametrized approach of our model. Nonetheless, it is gratifying that, so long as collisions are suitably inelastic, our model finds inward transport for both vertical and radial stirring (e.g., convection), with the transport naturally much more efficient for the latter than for the former. Inward transport in the presence of purely vertical stirring requires, however, an anistropic “inelasticity,” for which “collisions” are more inelastic in the azimuthal direction than in the radial direction (§4.4.2). By equating inelasticity with inefficient fluctuating pressure forces, as advocated above, we find results that are consistent with analytic arguments and simulations of hydrodynamic turbulence (e.g., Balbus & Hawley 1998).

To conclude, we note that our analysis expands on several aspects of the Boltzmann equation model for radial convection in rotating media considered by Kumar, Narayan, & Loeb (1995). These authors treated “collisions” between convective eddies using the same scattering function we have considered (eq. ) and found that for sufficiently inelastic collisions, radial convection would produce inward angular momentum transport. Interestingly, we now see that this inward transport has little to do with the details of convectively unstable fluids (e.g., the superadiabatic entropy gradient, correlated temperature and velocity fluctuations, etc. analyzed by Kumar et al.). Rather, it is a generic property of disks in the presence of external “stirring” of the radial velocity dispersion. Moreover, we have emphasized the analogy between inelastic collisions in particle disks and inefficient fluctuating pressure forces in turbulent disks, better establishing the relevance of kinetic theory models to the fluid problem.

###### Acknowledgements.

We thank Steve Balbus, Roger Blandford, Peter Goldreich, Jeremy Goodman, Pawan Kumar, Ramesh Narayan, Scott Tremaine, and an anonymous referee for useful discussions and comments. EQ is supported by NASA through Chandra Fellowship PF9-10008, awarded by the Chandra X–ray Center, which is operated by the Smithsonian Astrophysical Observatory for NASA under contract NAS 8-39073. EC acknowledges partial support from a Caltech Kingsley Foundation Fellowship.

## Appendix A Solution for General $`\kappa `$

Defining $`x=\mathrm{\Omega }\tau `$, $`y=\kappa \tau `$, and $`x^{}=dx/d\mathrm{ln}r=(y^24x^2)/2x`$, one can solve equations (23)–(26) to find:

$$p_{r\varphi }=x\sigma ^2\left(\frac{\frac{2}{3}\left[1\xi _\varphi \right]\frac{y^2}{6x^2}\left[1\xi _r\right]\stackrel{~}{H}_{rr}\frac{\xi y^2}{2x^2}+2\stackrel{~}{H}_{\varphi \varphi }\xi }{1+4y^2+\frac{y^2x^{}}{x}\stackrel{~}{H}_{rr}4xx^{}\stackrel{~}{H}_{\varphi \varphi }}\right)$$
(A1)
and

$$\sigma ^2=\frac{H+4Hy^2+\frac{y^2x^{}}{x}H_{rr}4xx^{}H_{\varphi \varphi }}{(1+4y^2)\xi \frac{x^{}}{3}\left[\frac{y^2}{x}\left(1\xi _r\right)4x\left(1\xi _\varphi \right)\right]}.$$
(A2)


text:
# On the formation and evolution of the globular cluster 𝜔 Centauri

## 1 Introduction

$`\omega `$ Centauri is one of the most interesting globular cluster in the Milky Way, and maybe the most studied one (Majewski et al. 1999, Lee et al. 1999, Pancino 1998). A compilation of its main properties is given in Table 1.
The most striking feature of this cluster is the measured metallicity spread, which has been interpreted as the evidence of a multiple stellar population inside it (Norris et. al 1996, Suntzeff & Kraft 1996). The precise metallicity distribution function (MDF) is nonetheless still disputed. Norris et al (1996) suggest the presence of a secondary peak at $`[Ca/H]0.9`$ roughly 5 times smaller that the main peak at $`[Ca/H]1.4`$. This trend is not confirmed by Suntzeff & Kraft (1996), who claim for a more regular MDF.
On the base of a large photometric survey Lee et al. (1999) (but see also Ortolani et al 1999) analyzed the color distribution of a sample of bright stars, showing that on the average it has an e-folding trend, with the presence of several significant metallicity peaks. However Majewski et al. (1999) arrive at a somewhat different conclusion, showing that the MDF has a gaussian shape with the maximum at $`[Fe/H]1.7`$, and with some evidences of a secondary peak. Although different, all these analyses point to the common picture of an object which experienced an irregular self-enrichment over its evolution.
Putting together the chemical and kinematical properties Majewski et al. (1999) claim that $`\omega `$ Centauri might be a possible dwarf galaxy relict.
In this paper we propose a N-body/gasdynamical model for the formation and evolution of $`\omega `$ Centauri, suggesting that this globular cluster can actually be the remnant of a dwarf elliptical galaxy, formed and evolved avoiding strong mergers, and eventually captured by the Milky Way.
To this aim, the plan of this paper is as follows. In Sect. 2 we describe our model and the initial conditions setup; Sect. 3 to 5 are dedicated to the analysis of the structure, chemistry and internal kinematics, respectively, of our model, and the comparison with $`\omega `$ Centauri; in Sect. 6 we investigate about the possible presence of an extended DM halo around the cluster. Finally Sect. 7 summarizes our results.

## 2 Technique and Initial Conditions

The simulation presented here has been performed using the Tree–SPH code developed by Carraro at al. (1998) and Buonomo et al. (2000). In this code, the properties of the gas component are described by means of the Smoothed Particle Hydrodynamics (SPH) technique, whereas the gravitational forces are computed by means of a hierarchical tree algorithm using a tolerance parameter $`\theta =0.8`$ and expanding tree nodes to quadrupole order. We adopt a Plummer softening parameter. In SPH each particle represents a fluid element whose position, velocity, energy, density etc. are followed in time and space. The properties of the fluid are locally estimated by an interpolation which involves the smoothing length $`h_i`$. In our code each particle possesses its own time and space variable smoothing length $`h_i`$, and evolves with its own time-step. This renders the code highly adaptive and flexible, and suited to speed-up the numerical calculations. Radiative cooling is described by means of numerical tabulations as a function of temperature and metallicity taken from Sutherland & Dopita (1993) and Hollenbach & McKee (1979). This allows us to account for the effects of variations in the metallicity among the fluid elements and for each of these as a function of time and position. The chemical enrichment of the gas-particles caused by Star Formation (SF) and stellar ejecta is described by means of the closed-box model applied to each gas-particle (cf. Carraro et al. 1998 for more details). SF and Feed-back algorithm are described in great detail by Buonomo et al. (2000). In brief, SF is let depend on the total mass density - baryonic (gas and stars) and dark matter \- of the system and on the metal–dependent cooling efficiency. Moreover we consider the effects of energy (and mass) feed-back from SNe of type II and Ia, and stellar winds from massive stars.

Instead of selecting halos from large cosmological N-body simulations, we opted to set up a protogalaxy as an isolated virialized DM halo, whose density profile is:

$$\rho (r)\frac{1}{r}.$$

DM particles are distributed inside the sphere according to an acceptance–rejection procedure. Along with positions, we assign also velocities according to:

$$v(r)\sqrt{(}r\times ln(\frac{1}{r})).$$

which is the solution of the Jeans equation for spherical isotropic collisionless system for the adopted density profile.

The system is let evolve until virial equilibrium is reached. Gas particles are then distributed homogeneously inside the DM halo with zero velocity field, thus mimicking the cosmological gas infall inside the DM potential well (White & Rees 1978).
For the simulation presented here we consider a $`9\times 10^7M_{}`$ DM halo, sampled by 10,000 particles, homogeneously filled with $`1\times 10^7M_{}`$ baryonic mass in form of metal poor ($`Z=10^4`$) gas, sampled by 10,000 particles. This way each gas particle has a mass of $`10^3M_{}`$.

## 3  The structure

The evolution of our model is shown in Fig. 1, where we plot the SF as a function of time. Starting from a $`10^4`$ oK metal poor gas, cooling drives baryons towards the center of the DM potential well. Part of this gas (about 20%) is transformed into stars with a very irregular SF history, whereas the remaining 80% is thrown away via a galactic wind mechanism, for which SNe and stellar winds are responsible.
The final stars distribution resembles a triaxial object, with $`a/b=0.859`$ and $`a/c=0.967`$. Correspondingly the ellipticity $`(ϵ=1b/a)`$ turns out to be 0.141, not far from the mean value reported by Geyer et al. (1983). The stars distribution is shown in Fig. 2 (onto the X-Z plane) and 3 (onto the X-Y plane). For the sake of an easier comparison with $`\omega `$ Centauri the model has been placed at the cluster actual distance (5.3 kpc).
The total mass in stars amount to about $`2\times 10^6M_{}`$, sampled by about 2,000 star particles. The spherically averaged mass density profile derived from our simulation is shown in Fig. 4 . Each filled circle represents the stars density $`\rho _i`$ computed as the mass inside a spherical shell $`i`$ with size the softening $`ϵ`$, divided by the shell volume. For each density estimate $`\rho _i`$ we plot as error bar the uncertainty of the density estimate evaluated as the poissonian error related to the shell population.
We perform a fit with a King (1962) profile to check whether our model recovers the structural properties of $`\omega `$ Centauri. The cluster indeed (see Table 1) has a core and tidal radius $`r_c`$ and $`r_t`$ equal to 2.6 and 45 arcmin, respectively.
By using a King (1962) model (solid overimposed line), we found the best fit for $`r_c`$= 3.2 and $`r_t`$= 65 arcmin, which implies a concentration $`c=log(\frac{r_t}{r_c})`$ = 1.30, somewhat higher than the value reported in Table 1. While the value of $`r_c`$ we obtain agrees with the observed one, the value of $`r_t`$ turns out to be much larger than the measured one for $`\omega `$ Centauri. In other words we over-estimate the tidal radius. In our opinion this is due to the fact that in our simplified model we neglect the effect of the tidal stripping that the cluster probably did experience when it was accreted by the Milky Way. This way the outermost stars were probably removed reducing the tidal radius.
Finally, the central density, computed assuming the $`r_c`$ derived above, is about 2660 $`M_{}/pc^3`$.

## 4 Internal kinematics

The internal kinematics of $`\omega `$ Centauri has been analyzed in great detail by Merritt et al. (1997) using the radial velocity catalogue obtained with CORAVEL by Mayor et al. (1997). From this study it emerges that the cluster has a peak rotational speed of about 8 km/sec at 11 pc from the center. Moreover the cluster has a central velocity dispersion $`\sigma `$ = 17 km/sec, the higher one among globular clusters (see Fig. 8 in Merritt et al 1997).
From the analysis of our data we find evidences of smaller rotational speed (about 4 km/sec at 15 pc), whereas the velocity dispersion profile $`\sigma (r)`$ (see Fig. 5) shows a central value of about 13 km/sec, denoting that our model is somewhat colder than $`\omega `$ Centauri.

## 5 Metallicity distribution

The stars MDF at the end of the simulation is presented in Fig. 6. The bulk of stars has a metallicity in the range between $`\text{[Fe/H]}=1.9`$ and $`\text{[Fe/H]}1.6`$. A secondary peak is observed at $`\text{[Fe/H]}=1.2`$, about 4 times smaller than the major peak.
Finally a significant number of stars has a relatively high metallicity around $`\text{[Fe/H]}=0.50`$. The global shape resembles the e-folding SF history shown in Fig. 1, with secondary peaks which mark successive bursts of star formation.
The good agreement we find with the data presented by Majewski et al. (1999) confirms their suggestion that $`\omega `$ Centauri experienced an irregular self-enrichment over its evolution and may actually be the core of a larger dwarf elliptical galaxy.

## 6 Dark Matter around $`\omega `$ Centauri?

Our simulation starts with a virialised DM halo, in whose center gas concentrates due to cooling instability, forming stars. The final stars and DM distribution is shown in Fig. 7, where open squares refer to DM, and filled circles represent stars. The inner region of the cluster is dominated by stars up to 20 arcmin, which we call the transition radius $`r_{tr}`$. In the outer region DM dominates over a distance 6 times larger than the transition radius.
We expect that the internal kinematics of stars in the inner region is unlikely to be influenced by DM. On the contrary, stars outside the transition radius $`r_{tr}`$ (20 arcmin) should show a kinematics strongly influenced by DM. A detail spectroscopic study of the stars kinematics in the cluster envelope should be able to confirm or deny the presence of DM. There is indeed in our model a trend of the velocity dispersion to weakly increase out of the transition radius.

## 7 Conclusions

We have presented a N-body/gasdynamical simulation of the formation and evolution of the globular cluster $`\omega `$ Centauri. We are able to reproduce the bulk properties of the cluster, namely structure, kinematics and chemistry assuming that it formed and evolved in isolation, and eventually fell inside the Milky Way potential well.
According to our results and to the dwarf galaxies taxonomy proposed by Roukema (1999), $`\omega `$ Centauri can actually be a cosmological dwarf by mass, formed in a high redshift low mass halo, which escaped significant merging up to the present time. Finally we stress that in order to obtain these results, an extended DM halo should surround the present day $`\omega `$ Centauri.

### Acknowledgments

G.C. expresses his gratitude to Prof. S. Ortolani for fruitful discussions. The referee, George Meylan, is acknowledged for important suggestions which led to the improvement of the paper. This work has been financed by italian MURST and ASI.

## References

Buonomo F., Carraro G., Chiosi C., Lia C., 2000, MNRAS 312, 371

Carraro G., Lia C., Chiosi C., 1998, MNRAS 291, 1021

Geyer E.H. Hopp U., Nelles B., 1983, A&A 125, 359

Harris W. E. 1996, AJ 112, 1487

Hollenbach D., McKee C.F., 1979, ApJS 41, 555

King I., 1962, AJ 67, 471

Lee Y.-W., Joo J.-M., Sohn Y.-J., Lee H.-C., Walker A.R., 1999, Nat., in press (astro-ph/9911137)

Majewski S.R., Patterson R.J., Dinescu D.I. et al., 1999, in ”The Galactic Halo: from Globular Clusters to Field Stars”, 35th Liege Astrophysical Colloquium (astro-ph/9910278)

Mayor M., Meylan G., Udry S. et al., 1997, AJ 114, 1087

Merritt D., Meylan G., Mayor M., 1997, AJ 114, 1074

Ortolani S., Covino S., Carraro G., 1999, in ”The Third Stromlo Simposium: The Galactic Halo”, Eds Gibson B.K., Axelrod T.S., Putman M.E., ASP Conference Series Vol. 165, 316

Norris J.E., Freeman K.C., Mighell K.J., 1996, ApJ 462, 251

Pancino E., 1998, Master Thesis, Padova University

Roukema B.F., 1999, in ”XVI Rencontre de Moriond: Dwarf Galaxies and Cosmology”, in press (astro-ph/9806111)

Searle L., Zinn R., 1978, ApJ 225, 357

Suntzeff N.B., Kraft R.P., 1996, AJ 111, 1913

Sutherland R.S., Dopita M.A., 1993, ApJS 88, 253

White S.D.M., Rees M., 1978, MNRAS 183, 341


text:
# Spin-Gap States of a Periodic Mixed Spin Chain

## 1 Introduction

Since the discovery of cuprate superconductors, quantum spin systems have been studied to understand the magnetic properties of their undoped mother materials. In particular quantum spin systems with spin gap are interesting since the spin gap is possibly related to the superconducting gap in the doped case. Although the materials are two-dimensional, it is basically important to study magnetic properties of spin gap systems in any dimensions.

In one dimension, Haldane  predicted that the spin excitation spectrum of a uniform spin chain is gapless if the magnitude of a spin is a half-odd-integer, and is gapful if it is an integer. The prediction has been confirmed theoretically and experimentally in many standpoints. Since Haldane’s prediction is based on a mapping of a spin Hamiltonian to the nonlinear $`\sigma `$ model (NLSM), it is recognized that the NLSM is a useful tool to investigate spin systems.

The NLSM has been developed for periodic inhomogeneous spin chains. Affleck  derived and examined an NLSM for a spin chain with bond alternation. The application of the NLSM to spin chains with more than one spin species in the period of more than two lattice spacings are also considered . In particular we have succeeded to derive an NLSM, with accurately keeping the degrees of freedom of the spin variables, for a general mixed spin chain with arbitrary finite period .

We applied the NLSM method to systems with period 4 and spin species 2 (we defines the magnitudes as $`s_a`$ and $`s_b`$) . Imposing the condition that the ground state is singlet, only the possible order of spin magnitudes in a unit cell is of $`s_a`$-$`s_a`$-$`s_b`$-$`s_b`$. We obtained phase diagrams in the parameter space of the exchange couplings. In a few special cases, numerical calculations have been performed . The phase diagrams obtained by the NLSM are qualitatively agree with the numerical results.

In this paper we examine another mixed spin chain with period 4 which includes three species of spins. The sequential order of magnitudes of spins in a unit cell is 1/2, 1, 3/2 and 1. There are four kinds of nearest-neighbour exchange couplings among them. If we consider a simplified case, the exchange parameters on the both sides of a 1/2 spin are the same value $`J`$, and those on the both sides of a 3/2 spin are of the same value $`J^{}`$. We in fact consider the distortion of $`J^{}`$ on the both sides of a 3/2 spin. We map the spin Hamiltonian describing the spin chain to an NLSM following Ref. . From the value of the topological angle in the NLSM, we determine the phase diagram of this model in the space of the exchange parameters.

## 2 Mapping to the nonlinear $`\sigma `$ model

The Hamiltonian for the present spin chain is written as

$`H`$ $`=`$ $`{\displaystyle \underset{j}{}}(J_1𝐒_{4j+1}𝐒_{4j+2}+J_2𝐒_{4j+2}𝐒_{4j+3}`$ (1)
$`+`$ $`J_3𝐒_{4j+3}𝐒_{4j+4}+J_4𝐒_{4j+4}𝐒_{4j+5}),`$
where the magnitudes of $`𝐒_{4j+1}`$, $`𝐒_{4j+2}`$, $`𝐒_{4j+3}`$ and $`𝐒_{4j+4}`$ are

$`(s_1,s_2,s_3,s_4)=({\displaystyle \frac{1}{2}},1,{\displaystyle \frac{3}{2}},1).`$ (2)
A unit cell is illustrated in Fig. 1. To reduce the number of parameters, we restrict and reparametrize the antiferromagnetic exchange parameters as

$`J_1=J_4=J,J_2={\displaystyle \frac{J^{}}{1\gamma }},J_3={\displaystyle \frac{J^{}}{1+\gamma }},`$ (3)
where $`\gamma `$ ($`1<\gamma <1`$) is the distortion parameter describing the asymmetry between the couplings of the both sides of a spin with magnitude 3/2.

Following Ref. , we can map a general periodic mixed spin chain satisfying a restriction to the NLSM. The mapped NLSM action is given by

$`S_{\mathrm{eff}}=`$ $`{\displaystyle }d\tau {\displaystyle }dx\{i{\displaystyle \frac{J^{(0)}}{J^{(1)}}}𝐦({\displaystyle \frac{𝐦}{\tau }}\times {\displaystyle \frac{𝐦}{x}})`$ (4)
$`+{\displaystyle \frac{1}{2aJ^{(1)}}}\left({\displaystyle \frac{J^{(1)}}{J^{(2)}}}{\displaystyle \frac{J^{(0)}}{J^{(1)}}}\right)\left({\displaystyle \frac{𝐦}{\tau }}\right)^2`$
$`+{\displaystyle \frac{a}{2}}J^{(0)}\left({\displaystyle \frac{𝐦}{x}}\right)^2\},`$
where $`𝐦`$ is an O(3) unit-vector field and $`a`$ is the lattice constant. For a general spin chain with period $`2b`$, the constants in the action are given as follows:

$`{\displaystyle \frac{1}{J^{(n)}}}`$ $`=`$ $`{\displaystyle \frac{1}{2b}}{\displaystyle \underset{q=1}{\overset{2b}{}}}{\displaystyle \frac{(\stackrel{~}{s}_q)^n}{J_qs_qs_{q+1}}}(n=0,1,2)`$ (5)
with accumulated spins

$`\stackrel{~}{s}_q`$ $`=`$ $`{\displaystyle \underset{k=1}{\overset{q}{}}}(1)^{k+1}s_k.(q=1,2,\mathrm{},2b)`$ (6)
The action (4) is of the standard form of the NLSM:

$`S_{\mathrm{st}}`$ $`=`$ $`{\displaystyle }d\tau {\displaystyle }dx\{i{\displaystyle \frac{\theta }{4\pi }}𝐦({\displaystyle \frac{𝐦}{\tau }}\times {\displaystyle \frac{𝐦}{x}})`$ (7)
$`+{\displaystyle \frac{1}{2gv}}\left({\displaystyle \frac{𝐦}{\tau }}\right)^2+{\displaystyle \frac{v}{2g}}\left({\displaystyle \frac{𝐦}{x}}\right)^2\},`$
where $`\theta `$ is the topological angle, $`g`$ is the coupling constant and $`v`$ is the spin wave velocity.

In the present model, the accumulated spins (6) are

$`(\stackrel{~}{s}_1,\stackrel{~}{s}_2,\stackrel{~}{s}_3,\stackrel{~}{s}_4)=({\displaystyle \frac{1}{2}},{\displaystyle \frac{1}{2}},1,0)`$ (8)
and then the constants (5) are calculated as

$`{\displaystyle \frac{1}{J^{(0)}}}={\displaystyle \frac{1}{J}}+{\displaystyle \frac{1}{3J^{}}},`$
$`{\displaystyle \frac{1}{J^{(1)}}}={\displaystyle \frac{1}{4}}\left({\displaystyle \frac{1}{J}}+{\displaystyle \frac{1}{3J^{}}}+{\displaystyle \frac{\gamma }{J^{}}}\right),`$
$`{\displaystyle \frac{1}{J^{(2)}}}={\displaystyle \frac{1}{8}}\left({\displaystyle \frac{1}{J}}+{\displaystyle \frac{5}{3J^{}}}+{\displaystyle \frac{\gamma }{J^{}}}\right).`$ (9)
Comparing Eq. (4) with Eq. (7), we have the topological angle

$`\theta =4\pi {\displaystyle \frac{J^{(0)}}{J^{(1)}}}=\pi +{\displaystyle \frac{\pi \gamma }{{\displaystyle \frac{J^{}}{J}}+{\displaystyle \frac{1}{3}}}},`$ (10)
the spin wave velocity

$`v={\displaystyle \frac{4aJ^{}}{\sqrt{\left({\displaystyle \frac{J^{}}{J}}+{\displaystyle \frac{1}{3}}\right)\left({\displaystyle \frac{J^{}}{J}}+3\right)\gamma ^2}}}`$ (11)
and the coupling constant

$`g={\displaystyle \frac{4\left({\displaystyle \frac{J^{}}{J}}+{\displaystyle \frac{1}{3}}\right)}{\sqrt{\left({\displaystyle \frac{J^{}}{J}}+{\displaystyle \frac{1}{3}}\right)\left({\displaystyle \frac{J^{}}{J}}+3\right)\gamma ^2}}}.`$ (12)
These quantities are shown in Figs. 2, 3 and 4, as functions of $`J^{}/J`$ for several values of $`\gamma `$.

## 3 Phase diagram

It is well known that an NLSM has gapless excitations if the topological angle $`\theta `$ is an odd-integer multiple of $`\pi `$. In the general periodic spin model, this condition is turned into the gapless equation

$$\frac{2J^{(0)}}{J^{(1)}}=\frac{2l1}{2}$$
(13)
with arbitrary integer $`l`$. In the present case of Eq. (10), it becomes of a simple form

$`\gamma =2l({\displaystyle \frac{J^{}}{J}}+{\displaystyle \frac{1}{3}}).(l=1,0,1)`$ (14)
The equation for each $`l`$ determines a phase boundary between gapful disordered phases. The restricted values of integer $`l`$ is due to the condition $`1<\gamma <1`$. In the symmetric case of $`\gamma =0`$, the system is gapless irrespective of the ratio $`J^{}/J`$, because it satisfies Eq. (14) with $`l=0`$.

The phase diagram of the ground state is obtained by identifying the gapless lines (14) as the phase boundaries. We show the phase diagram in Fig. 5. There exist four phases A+, A-, B+ and B- in the ($`J^{}/J`$, $`\gamma `$) parameter space. The phase structure is symmetric for $`\gamma \gamma `$ according to the symmetry of the Hamiltonian (1) with Eqs. (2) and (3).

In phase B+, $`\gamma `$ is close to 1 and $`J_2`$ is selectively large. The ground state of B+ may be explained by a VBS picture . That is, if we decompose each spin into spins with magnitude 1/2, the ground state is approximately a direct product of dimers. Two dimers in a unit cell are formed between 1/2-spins at the both sides of a coupling with exchange parameter $`J_2`$. Phase A+ is of an intermediate character and spreads over almost all area in the parameter space. The ground state of A+ does not seem to be explained by the conventional VBS picture.

## 4 Summary

We have examined the 1/2-1-3/2-1 spin chain with a distortion of the exchange parameters by using the nonlinear $`\sigma `$ model method. We obtained the phase diagram of the ground state in the parameter space. There exist four gapful phases. In the case of no distortion, the system is always gapless and is on a phase boundary. When the distortion is strong, there are two phases to be explained by a VBS picture. It is a future work to explain the other phases by the singlet-cluster-solid (SCS) picture developed in Ref. .


text:
# Quantum Dot as Spin Filter and Spin Memory

\[

## Abstract

We consider a quantum dot in the Coulomb blockade regime weakly coupled to current leads and show that in the presence of a magnetic field the dot acts as an efficient spin-filter (at the single-spin level) which produces a spin-polarized current. Conversely, if the leads are fully spin-polarized the up or down state of the spin on the dot results in a large sequential or small cotunneling current, and thus, together with ESR techniques, the setup can be operated as a single-spin memory.

\]

An increasing number of spin-related experiments  show that the spin of the electron offers unique possibilities for finding novel mechanisms for information processing—most notably in quantum-confined semiconductors with unusually long spin dephasing times approaching microseconds , and where spins can be transported coherently over distances of up to 100 micrometers . Besides the intrinsic interest in spin-related phenomena, spin-based devices hold promises for future applications in conventional  as well as in quantum computer hardware . One of the challenging problems for such applications is to obtain sufficient control over the spin dynamics in nanostructures. In the following we address this issue and propose a quantum-dot setup which can be either operated as a spin-filter (spin diode) to produce spin-polarized currents or as a device to detect and manipulate single-spin states (single-spin memory). Both effects occur at the single-spin level and thus represent the ultimate quantum limit of a spin-filter and spin-memory. In both cases, we will work in the Coulomb blockade regime and consider sequential and cotunneling processes. A new feature of our proposal is that the spin-degeneracy is lifted with different Zeeman splittings in the dot and in the leads which then results in Coulomb blockade peaks which are uniquely associated with a definite spin state on the dot.

Formalism.– Our system consists of a quantum dot (QD) connected to two Fermi-liquid leads which are in equilibrium with reservoirs kept at the chemical potentials $`\mu _l`$, $`l=1,2`$, where outgoing currents can be measured, see Fig. 1. Using a standard tunneling Hamiltonian approach, we write for the full Hamiltonian $`H_0+H_T`$, where $`H_0=H_L+H_D`$ describes the leads and the dot, with $`H_D`$ including the charging and interaction energies of the dot electrons as well as the Zeeman energy $`g\mu _BB`$ of their spins in the presence of an external magnetic field $`𝐁=(0,0,B)`$, where $`g`$ is the effective g-factor. We concentrate first on unpolarized lead currents and assume that the Zeeman splitting in the leads is negligibly small compared to the one in the QD. This can be achieved e.g. by using InAs for the dot ($`g=15`$) attached to GaAs 2DEG leads ($`g=0.44`$), or by implanting a magnetic impurity (say Mn) inside a GaAs dot (again attached to GaAs 2DEG leads) with a strongly enhanced electron g-factor due to exchange splitting with the magnetic impurity . \[Below we will also consider the opposite situation with a fully spin polarized lead current, and a much smaller Zeeman splitting on the dot.\] The tunneling between leads and the QD is described by the perturbation $`H_T=_{l,k,p,s}t_{lp}c_{lks}^{}d_{ps}+\mathrm{h}.\mathrm{c}.`$, where $`d_{ps}`$ and $`c_{lks}`$ annihilate electrons with spin $`s`$ in the dot and in the $`l`$th lead, resp. While the orbital $`k`$-dependence of the tunneling amplitude $`t_{lp}`$ can be safely neglected, this is not the case in general for the QD orbital states $`p`$. From now on we concentrate on the Coulomb blockade (CB) regime, where the charge in the QD, $`\widehat{N}=_{p,s}d_{ps}^{}d_{ps}`$, is quantized, $`\mathrm{Tr}\rho \widehat{N}=N`$.

Next, turning to the dynamics induced by $`H_T`$, we introduce the reduced density matrix for the dot, $`\rho _D=\mathrm{Tr}_L\rho `$, where $`\rho `$ is the full stationary density matrix, and $`\mathrm{Tr}_L`$ is the trace over the leads. To describe the stationary limit, we use a standard master equation approach  formulated in terms of the dot eigenstates and eigenenergies, $`H_D|n=E_n|n`$, where $`n=(𝐧,N)`$. Denoting with $`\rho (n)=n|\rho _D|n`$ the stationary probability for the dot to be in the state $`|n`$, and with $`W(n^{},n)`$ the transition rates between $`n`$ and $`n^{}`$, the stationary master equation to be solved reads $`_n\left[W(n^{},n)\rho (n)W(n,n^{})\rho (n^{})\right]=0`$.

The rates $`W`$ can be calculated in a standard “golden rule” approach  where we go up to 2nd order in $`H_T`$, i.e. $`W=_lW_l+_{l^{},l}W_{l^{}l}`$, where $`W_lt^2`$ is the rate for a tunneling process of an electron from the $`l`$th lead to the dot and back, while $`W_{l^{}l}t^4`$ describes the simultaneous tunneling of two electrons from the lead $`l`$ to the dot and from the dot to the lead $`l^{}`$. Thus, two regimes of transport through the QD can be distinguished: Sequential tunneling (ST) and cotunneling (CT). The ST regime is at the degeneracy point, where $`\widehat{N}`$ fluctuates between $`N`$ and $`N^{}=N\pm 1`$, and 1st order transitions are allowed by energy conservation with the explicit rates

$`W_l(n^{},n)`$ $`=`$ $`2\pi \nu [f_l(\mathrm{\Delta }_{n^{}n})|A_{lnn^{}}|^2\delta _{N^{},N+1}`$ (1)
$`+`$ $`[1f_l(\mathrm{\Delta }_{nn^{}})]|A_{ln^{}n}|^2\delta _{N^{},N1}],`$ (2)
where $`\nu =_k\delta (\epsilon _F\epsilon _k)`$ is the lead density-of-states per spin at the Fermi energy $`\epsilon _F`$, $`f_l(\epsilon )=[1+\mathrm{exp}((\epsilon \mu _l)/k_BT)]^1`$ the Fermi function at temperature $`T`$, $`\mathrm{\Delta }_{n^{}n}=E_n^{}E_n`$ is the level distance, and we have introduced the matrix elements $`A_{ln^{}n}=_{ps}t_{lps}n^{}|d_{ps}|n`$. In the ST regime the current through the QD can be written as $`I_s=\pm e_{n,n^{}}W_2(n^{},n)\rho (n)`$, where $`\pm `$ stands for $`N^{}=N1`$. We emphasize that the rates $`W(n,n^{})`$ and thus the current depend on the spin state of the dot electrons via $`n,n^{}`$. The ST current takes a particularly simple form if bias $`\mathrm{\Delta }\mu =\mu _1\mu _2>0`$ and temperature are small compared to the level distance on the dot (the case of interest here), $`\mathrm{\Delta }\mu ,k_BT<|\mathrm{\Delta }_{mn}|,m,n`$, and, thus only the lowest energy levels participate in the transport . The solution of the master equation gives in this case for the ST current

$$I_s=\frac{e\gamma _1\gamma _2}{\gamma _1+\gamma _2}\left[f_1(\mathrm{\Delta }_{n^{}n})f_2(\mathrm{\Delta }_{n^{}n})\right],$$
(3)
where $`\gamma _l=2\pi \nu |A_{lnn^{}}|^2`$ is the tunneling rate through the $`l`$th barrier.

In the CT regime the only allowed processes are 2nd order transitions with initial and final electron number on the QD being equal, i.e. $`N=N^{}`$, and with the rate

$`W_{l^{}l}(n^{},n)=2\pi \nu ^2{\displaystyle 𝑑\epsilon f_l(\epsilon )[1f_l^{}(\epsilon \mathrm{\Delta }_{n^{}n})]}`$ (4)
$`\times \left|{\displaystyle \underset{n_1}{}}{\displaystyle \frac{A_{l^{}n^{}n_1}A_{lnn_1}^{}}{\mathrm{\Delta }_{nn_1}+\epsilon }}+{\displaystyle \underset{n_2}{}}{\displaystyle \frac{A_{l^{}n_2n}A_{ln_2n^{}}^{}}{\mathrm{\Delta }_{n^{}n_2}\epsilon }}\right|^2,`$ (5)
where, $`N_1=N+1`$, and $`N_2=N1`$, and thus the two terms in Eq. (5) differ by the sequence of tunneling. Our regime of interest here is elastic CT where $`E_n^{}=E_n`$, which holds for $`|\mathrm{\Delta }_{mn}|>\mathrm{\Delta }\mu ,k_BT,mn`$. That means the system is always in the ground state with $`\rho (n)=1`$, and thus the CT current is given by $`I_c=eW_{21}(n,n)eW_{12}(n,n)`$. In particular, close to a ST resonance (but still in the CT regime) Eq. (5) considerably simplifies—only one term contributes—and for $`\mathrm{\Delta }\mu ,k_BT<|\mu \pm \mathrm{\Delta }_{nn_i}|`$, we obtain

$$I_c=\frac{e}{2\pi }\frac{\gamma _1\gamma _2\mathrm{\Delta }\mu }{(\mu \pm \mathrm{\Delta }_{nn_i})^2},$$
(6)
where $`+`$ stands for $`i=1`$, and $``$ for $`i=2`$, and $`\mu =(\mu _1+\mu _2)/2`$. From Eqs. (3) and (6) it follows that $`I_s\gamma _i`$, while $`I_c\gamma _i^2`$, and therefore $`I_cI_s`$. Thus, in the CB regime the current as a function of $`\mu `$ (or gate voltage) consists of resonant ST peaks, where $`\widehat{N}`$ on the QD fluctuates between $`N`$ and $`N\pm 1`$. The peaks are separated by plateaus, where $`N`$ is fixed, and where the (residual) current is due to CT.

We note that the tunneling rates $`\gamma _l`$ depend on the tunneling path through the matrix elements $`A_{lmn}`$. In general, this can lead to a spin-dependence of the current, which, however, is difficult to measure. In contrast to this, we will show now that a much stronger spin dependence can come from the resonance character of the currents $`I_s`$ and $`I_c`$, when the position of a resonance (as function of gate voltage) depends on the spin orientation of the tunneling electron. To proceed we first specify the energy spectrum of the QD more precisely. In general, the determination of the spectrum of a QD is a complicated many-electron problem . However, it is known from experiment that especially away from orbital degeneracy points (which can be easily achieved by applying magnetic fields) the spectrum is formed mainly by single-particle levels, possibly slightly renormalized by exchange interaction .

For a QD with $`N`$ odd there is one unpaired electron in one of the two lowest energy states, $`|`$ and $`|`$, with energies $`E_{}`$ and $`E_{}`$, which become Zeeman split due to a magnetic field $`B`$, $`\mathrm{\Delta }_z=|E_{}E_{}|=\mu _B|gB|`$. Let us assume that $`|`$ is the ground state, and set $`E_{}=0`$ for convenience. For $`N`$ even, the two topmost electrons (with the same orbital wave function), form a spin singlet, $`(||)/\sqrt{2}`$, with energy $`E_S`$. This is the ground state, while the other states, such as three triplets $`|T_+=|`$, $`|T_{}=|`$, and $`|T_0=(|+|)/\sqrt{2}`$ with energies $`E_{T_\pm }`$ and $`E_{T_0}`$ are excited states, because of higher (mostly) single-particle orbital energy.

Sequential tunneling. First we consider the ST peak, which separates two plateaus with $`N`$ and $`N+1`$ electrons on the dot, where $`N`$ is odd (odd-to-even transition). In the regime $`E_{T_+}E_S,\mathrm{\Delta }_z>\mathrm{\Delta }\mu `$, $`k_BT`$, only ground-state transitions are allowed by energy conservation. The tunneling of spin-up electrons is blocked by energy conservation, i.e. $`I_s()=0`$, because it involves excited states $`|T_+`$ and $`|`$. The only possible process is tunneling of spin-down electrons as shown in Fig. 2, which leads to a spin-polarized current, $`I_s()`$, given by Eq. (3), where $`\mathrm{\Delta }_{n^{}n}=E_S>0`$ (since $`E_{}=0`$). Thus, we have

$`I_s()/I_0=\theta (\mu _1E_S)\theta (\mu _2E_S),k_BT<\mathrm{\Delta }\mu ,`$ (7)
$`I_s()/I_0={\displaystyle \frac{\mathrm{\Delta }\mu }{4k_BT}}\mathrm{cosh}^2\left[{\displaystyle \frac{E_S\mu }{2k_BT}}\right],k_BT>\mathrm{\Delta }\mu ,`$ (8)
where $`I_0=e\gamma _1\gamma _2/(\gamma _1+\gamma _2)`$. Hence, in the specified regime the dot acts as spin filter through which only spin-down electrons can pass .

Next, we consider the ST peak at the transition from even to odd, i.e. when $`N`$ is even. Then the current is given by Eq. (3) with $`\mathrm{\Delta }_{n^{}n}=E_S>0`$. The spin-down current is now blocked, $`I_s()=0`$, while spin-up electrons can pass through the QD, with the current $`I_s()`$ given by (7) and (8), where $`E_S`$ has to be replaced by $`E_S`$. Because this case is very similar to the previous one with $``$ replaced by $``$, we shall concentrate on the odd-to-even transition only. Next, we will demonstrate that although CT processes can in general lead to a leakage of unwanted current, this turns out to be a minor effect, and spin filtering works also in the CT regime.

Cotunneling. Above or below a ST resonance, i.e. when $`E_S>\mu _{1,2}`$ or $`E_S<\mu _{1,2}`$, the system is in the CT regime. Close to this peak the main contribution to the transport is due to two processes (a) and (c) see Fig. 3, where the energy deficit of the virtual states, $`|\mu E_S|`$, is minimal. According to Eq. (6) we have

$$I_c()=\frac{e}{2\pi }\frac{\gamma _1\gamma _2\mathrm{\Delta }\mu }{(\mu E_S)^2}.$$
(9)
Thus, we expect the spin filtering of down electrons to work even in the CT regime close to the resonance. However, there are additional CT processes, (b) and (d), which involve tunneling of spin-up electrons and lead to a leakage of up-spin. If $`N`$ is odd (below the resonance), the dot is initially in its ground state $`()`$, and an incoming spin-up electron can only form a virtual triplet state $`|T_+`$ (process (b) in Fig. 3). This process contributes to the rate (3) with an energy deficit $`E_{T_+}\mu `$, so that for the efficiency of spin filtering \[defined as $`I()/I()`$\] we obtain in this regime,

$$I_c()/I_c()\left(1+\frac{E_{T_+}E_S}{E_S\mu }\right)^2,N\text{ odd.}$$
(10)
Above the resonance, i.e. when $`N`$ is even and the ground state is the spin-singlet $`|S`$, the tunneling of spin-up electrons occurs via the virtual spin-down state (process (d) in Fig. 3) with an energy deficit $`\mathrm{\Delta }_z+\mu E_S`$, which has to be compared to the energy deficit $`\mu E_S`$ of the main process (c). Thus, we obtain for the efficiency of the spin filtering in the CT regime

$$I_c()/I_c()\left(1+\frac{\mathrm{\Delta }_z}{\mu E_S}\right)^2,N\text{ even.}$$
(11)
We see that in both cases, above and below the resonance, the efficiency can be made large by tuning the gate voltage to the resonance, i.e. $`|\mu E_S|0`$. Eventually, the system goes to the ST regime, $`|\mu E_S|k_BT,\mathrm{\Delta }\mu `$. Now, using Eqs. (6), (7), and (8) we can estimate the efficiency of spin filtering in the ST regime,

$$I_s()/I_c()\frac{\mathrm{\Delta }_z^2}{(\gamma _1+\gamma _2)\mathrm{max}\{k_BT,\mathrm{\Delta }\mu \}},$$
(12)
where we have assumed that $`\mathrm{\Delta }_z<|E_{T_+}E_S|`$. In the ST regime considered here we have $`\gamma _i<k_BT,\mathrm{\Delta }\mu `$ . Therefore, if the requirement $`k_BT,\mathrm{\Delta }\mu <\mathrm{\Delta }_z`$ is satisfied, filtering of spin-down electrons in the ST regime is very effecient. In the quantum regime, $`\gamma _i>k_BT,\mathrm{\Delta }\mu `$, tunneling occurs as Breit-Wigner resonance, and $`\mathrm{max}\{k_BT,\mathrm{\Delta }\mu \}`$ in Eq. (12) has to be replaced by $`\gamma _i`$. Finally, we note that the spin polarization of the transmitted current oscillates between up and down as we change the number of dot electrons $`N`$ one by one.

The functionality of the spin filter can be tested e.g. with the use of a p-i-n diode  which is placed in the outgoing lead 2. Via excitonic photoluminescence, the diode transforms the spin polarized electrons (entering lead 2) into correspondingly circularly polarized photons which can then be detected.

Spin memory. We consider now the opposite case where the current in the leads is fully spin polarized. Recent experiments have demonstrated that fully spin-polarized carriers can be tunnel-injected from a spin-polarized magnetic semiconductor (III-V or II-VI materials) with large effective g-factor into an unpolarized GaAs system . Another possibility would be to work in the quantum Hall regime where spin-polarized edge states are coupled to a quantum dot. To be specific, we consider the case where $`E_{T_+}E_s+\mathrm{\Delta }_z>\mathrm{\Delta }\mu ,k_BT`$ with $`E_{T_+}>E_s`$ ($`\mathrm{\Delta }_z>k_BT`$ is not necessary as long as the spin relaxation time is longer than the measurement time, see below). We assume that the spin polarization of both leads is, say, up and $`N`$ is odd. There are now two cases for the current, $`I^{}`$ or $`I^{}`$, corresponding to a spin up or down on the QD. First, we assume the QD to be in the ground state with its topmost electron-spin pointing up. According to previous analysis \[see paragraph before Eqs. (7,8)\], the ST current vanishes, i.e. $`I_s^{}=0`$, since the tunneling into the level $`E_{T_+}`$ (and higher levels) is blocked by energy conservation, while the tunneling into $`E_s`$ is blocked by spin conservation (the leads can provide and take up only electrons with spin up). However, there is again a small CT current, $`I_c^{}`$, which is given by Eq. (9). Now we compare this to the second case where the topmost dot-spin is down with additional Zeeman energy $`\mathrm{\Delta }_z>0`$. Here, the ST current $`I_s^{}`$ is finite, and given by Eqs. (7, 8) with $`E_s`$ replaced by $`E_s\mathrm{\Delta }_z`$. In the limit $`E_s>\mathrm{\Delta }_z`$, we get $`I_s^{}=I_s()`$, and thus the ratio $`I_s^{}/I_c^{}`$ is again given by Eq. (12). Hence, we see that the dot with its spin up transmits a much smaller current than the dot with spin down. This fact allows the read-out of the spin-state of the (topmost) dot-electron by comparing the measured currents. Furthermore, the spin state of the QD can be changed (“read-in”) by ESR techniques, i.e. by applying a pulse of an ac magnetic field (perpendicular to B) with resonance frequency $`\omega =\mathrm{\Delta }_z`$. Thus the proposed setup functions as a single-spin memory with read-in and read-out capabilities, the relaxation time of the memory given by the spin relaxation time $`\tau _S`$ on the QD (which can be expected to exceed 100’s of nanoseconds). We note that this $`\tau _S`$ can be easily measured since it is the time during which $`I_s^{}`$ is finite before it strongly reduces to $`I_c^{}`$. Finally, this scheme can be upscaled: In an array of such QDs where each dot separately is attached to in- and outgoing leads (for read-out) we can switch the spin-state of each dot individually by locally controlling the Zeeman splitting $`\mathrm{\Delta }_z`$. This can be done  e.g. by applying a gate voltage on a particular dot that pushes the wave function of the dot-electrons into a region of, say, higher effective g-factor (the induced level shift in the QD can be compensated for by the chemical potentials).

In conclusion, we have shown that quantum dots in the Coulomb blockade regime and attached to leads can be operated as efficient spin filters at the single-spin level. Conversely, if the leads are spin-polarized the spin state of the quantum dot can be read-out by a traversing current which is (nearly) blocked for one spin state while it is unblocked for the opposite spin state.

Acknowledgements. This work has been supported by the Swiss NSF.