text:
# 1 Introduction

## 1 Introduction

A maximally supersymmetric plane-wave background in ten dimensions is an exact solution of Type IIB superstring theory . Moreover, the string theory in this background is tractable, despite the fact that the background contains a nonvanishing RR field, provided that one uses the light-cone gauge Green–Schwarz (GS) formalism . In that approach the world-sheet theory consists of free massive bosons and fermions. Thus it is trivial to read off the complete spectrum of the noninteracting theory.

The string states and their interactions are holographically dual to certain operators and their correlation functions in $`𝒩=4`$ super Yang-Mills theory in a suitable limit . However, this paper will only consider the string side of the story and not discuss the duality.

The string interactions are encoded (using the formalism of light-cone-gauge string field theory) in a cubic interaction vertex. The three-string vertex has been formulated by Spradlin and Volovich  and explored further by other authors . These results generalize the flat-space light-cone gauge field theory results of  to the plane-wave geometry. For recent reviews see .

This paper reviews work that makes the formulas for the Neumann coefficients that enter in the interaction vertex more explicit. These coefficients are defined in the first instance in terms of the inverse of a certain infinite dimensional matrix. The first step is to express this inverse matrix in terms of a certain infinite component vector . The next step is to derive an expression for this vector in terms of a certain scalar quantity and then to derive an explicit formula for the scalar . Having obtained useful expressions for the Neumann coefficients (and hence the three superstring vertex), one can then explore the large $`\mu `$ (large curvature) limit, which is required for making contact with dual perturbative gauge theory computations.

## 2 Review of basic formulas

The type IIB superstring in the maximally supersymmetric plane-wave background is described in light-cone gauge by a free world-sheet theory. The eight bosonic and eight fermionic world-sheet fields each have mass $`\mu `$, a parameter that enters in the description of the plane-wave geometry and the RR five-form field strength. The mass term has two important consequences. One is that it leads to a mixing of left-movers and right-movers. The other is that the zero modes are also described by harmonic oscillators of finite frequency. Altogether, a convenient labelling of the bosonic lowering and raising operators arising from quantization of the free world-sheet theory is $`a_m^I`$ and $`a_m^I`$, where $`m`$ runs from minus infinity to plus infinity and $`I=1,\mathrm{},8`$. These satisfy ordinary oscillator commutation relations

$$[a_m^I,a_n^J]=\delta _{mn}\delta ^{IJ}.$$
(1)
There are also fermionic oscillators $`b_m^\alpha `$ and $`b_m^\alpha `$, which will not be discussed in this paper.

The spectrum of the free string theory is described by the light-cone Hamiltonian

$$H_2=\underset{m=\mathrm{}}{\overset{\mathrm{}}{}}\omega _mN_m,$$
(2)
where $`N_m`$ is the number of excitations of level $`m`$ oscillators

$$N_m=\underset{I=1}{\overset{8}{}}a_m^Ia_m^I+\mathrm{fermionic}\mathrm{terms},$$
(3)
and the frequencies are given by

$$\omega _m=\sqrt{m^2+\mu ^2\alpha ^2}.$$
(4)
The second term in the square root is actually $`(\alpha ^{}\mu p_{})^2`$, but we define $`\alpha =\alpha ^{}p_{}`$. (In flat space, we used the symbol $`p^+`$ rather than $`p_{}`$, but in curved space a subscript is more natural, since momenta are conjugate to coordinates that are defined with superscripts.) The physical spectrum is given by the product of all the oscillator spaces subject to one constraint

$$\underset{m=\mathrm{}}{\overset{\mathrm{}}{}}mN_m=0.$$
(5)
In flat space this constraint reduces to the usual level-matching condition for left-movers and right-movers.

The three-string interaction vertex for type IIB superstrings in flat space was worked out in  and generalized to the plane-wave geometry in . The formula can be written rather elegantly in terms of functionals, but to make its meaning precise and easily applicable to specific external states, it is desirable to expand it out in terms of oscillators. A convenient notation uses a tensor product of three string Fock spaces, labeled by an index $`r=1,2,3`$. Then the three string interaction vertex contains a factor

$$|V_B=\mathrm{exp}\left(\frac{1}{2}\underset{r,s=1}{\overset{3}{}}\underset{m,n=\mathrm{}}{\overset{\mathrm{}}{}}\underset{I=1}{\overset{8}{}}a_{mr}^I\overline{N}_{mn}^{rs}a_{ns}^I\right)|\mathrm{\hspace{0.17em}0}.$$
(6)
The quantities $`\overline{N}_{mn}^{rs}`$, called Neumann coefficients, are the main objects of concern in this paper.111Their definition here differs from that used in  by factors of $`\sqrt{mn}`$. The definition given here is more natural for the $`\mu 0`$ generalization. The three string vertex also contains a similar expression $`|V_F`$ made out of the fermionic oscillators and a “prefactor” that is polynomial in the various oscillators. We will not consider either of these in this paper. Suffice it to say that they are made out of the same basic objects, so that the results described here for the the bosonic factor are applicable to them as well.

In describing the Neumann matrices, it is convenient to consider separately the cases in which each of the indices $`m,n`$ are either positive, negative or zero. Henceforth, the symbols $`m,n`$ will always denote positive integers. Negative integers will be indicated by displaying an explicit minus sign. One result of , for example, using matrix notation for the blocks with positive indices, is

$$\overline{N}^{rs}=12(C_rC^1)^{1/2}A^{(r)T}\mathrm{\Gamma }_+^1A^{(s)}(C_sC^1)^{1/2}.$$
(7)
Here $`C_{mn}=m\delta _{mn}`$ and $`(C_r)_{mn}=\omega _{rm}\delta _{mn}`$, where

$$\omega _{rm}=\sqrt{m^2+(\mu \alpha _r)^2}.$$
(8)
The definitions of $`A^{(r)}`$ and $`\mathrm{\Gamma }_+`$, and other expressions that appear here, are collected in Appendix A.

The blocks with both indices negative are related in a simple way to the ones with both indices positive by

$$\overline{N}_{mn}^{rs}=\left(U_r\overline{N}^{rs}U_s\right)_{mn},$$
(9)
where

$$U_r=C^1(C_r\mu \alpha _r).$$
(10)
In the case of $`\overline{N}^{33}`$ these are the only nonvanishing terms. For the other Neumann coefficients the other nonvanishing terms are

$$\overline{N}_{m0}^{rs}=\overline{N}_{0m}^{sr}=\sqrt{2\mu \alpha _s}ϵ^{st}\alpha _t[(C_rC^1)^{1/2}A^{(r)T}Y]_mr=1,2,3s=1,2$$
(11)
$$\overline{N}_{00}^{rs}=(1+\mu \alpha k)ϵ^{rt}ϵ^{su}\sqrt{\alpha _t\alpha _u}r,s=1,2$$
(12)
$$\overline{N}_{00}^{r3}=\overline{N}_{00}^{3r}=\sqrt{\alpha _r}r=1,2.$$
(13)
Here we have introduced (see Appendix A)

$$Y=\mathrm{\Gamma }_+^1B,$$
(14)
$$k=B^T\mathrm{\Gamma }_+^1B,$$
(15)
$$y=\alpha _1/\alpha _3,$$
(16)
and (setting $`\alpha _3=1`$)

$$\alpha =\alpha _1\alpha _2\alpha _3=y(1y).$$
(17)
The asymmetry between string number three and the other two strings is a reflection of the fact that the $`\mu `$ dependence of the formula breaks the cyclic symmetry that is present in the flat space case.

To make the formulas useful for comparison with the dual gauge theory, it would be helpful to have explicit formulas for the various matrix multiplications and inversions that appear. The quantities that we especially would like to evaluate explicitly are the matrix $`\mathrm{\Gamma }_+^1`$, the vector $`Y=\mathrm{\Gamma }_+^1B`$, and the scalar $`k=B^T\mathrm{\Gamma }_+^1B`$. In the case of flat space ($`\mu =0`$) the results are known. Specifically

$$\overline{N}_{mn}^{rs}=\frac{mn\alpha }{m\alpha _s+n\alpha _r}\overline{N}_m^r\overline{N}_n^s\mathrm{for}\mu =0$$
(18)
where

$$\overline{N}_m^r=\frac{\sqrt{m}}{\alpha _r}f_m(\alpha _{r+1}/\alpha _r)e^{m\tau _0/\alpha _r}\mathrm{for}\mu =0,$$
(19)
$$f_m(\gamma )=\frac{\mathrm{\Gamma }(m\gamma )}{m!\mathrm{\Gamma }(m\gamma +1m)},$$
(20)
and

$$\tau _0=\underset{r=1}{\overset{3}{}}\alpha _r\mathrm{log}|\alpha _r|.$$
(21)
In particular, still for $`\mu =0`$, $`\mathrm{\Gamma }_+^1=\frac{1}{2}(1\overline{N}^{33})`$, $`Y_m=\overline{N}_m^3`$, and $`k=2\tau _0/\alpha `$.

## 3 Factorization theorem

In this section we will derive the generalization of eq. (18) that holds for the plane-wave geometry. The method of derivation is a fairly straightforward generalization of the one used for flat space in . We begin by defining

$$\stackrel{~}{\mathrm{\Gamma }}_+=\underset{r=1}{\overset{3}{}}A^{(r)}U_r^1A^{(r)T},$$
(22)
which differs from the $`\mathrm{\Gamma }_+`$ by the replacement of $`U_r`$ by $`U_r^1`$. Then we consider the product

$$\mathrm{\Gamma }_+C^1\stackrel{~}{\mathrm{\Gamma }}_+=(U_3+\underset{1}{\overset{2}{}}A^{(r)}U_rA^{(r)T})C^1(U_3^1+\underset{1}{\overset{2}{}}A^{(s)}U_s^1A^{(s)T})$$
(23)
Using various identities given in Appendix A, this simplifies to

$$\mathrm{\Gamma }_+C^1\stackrel{~}{\mathrm{\Gamma }}_+=U_3C^1\stackrel{~}{\mathrm{\Gamma }}_++\mathrm{\Gamma }_+C^1U_3^1\frac{1}{2}\alpha _1\alpha _2BB^T.$$
(24)
The next step is to use eqs. (67) and (63) to deduce that

$$\stackrel{~}{\mathrm{\Gamma }}_+=\mathrm{\Gamma }_++\mu \alpha BB^T$$
(25)
Substituting this into the previous equation and multiplying left and right by $`\mathrm{\Gamma }_+^1`$ gives

$$C^1U_3^1\mathrm{\Gamma }_+^1+\mathrm{\Gamma }_+^1U_3C^1=C^1+\frac{1}{2}\alpha _1\alpha _2YY^T+\mu \alpha ZY^T$$
(26)
where we have defined

$$Z=(1\mathrm{\Gamma }_+^1U_3)C^1B.$$
(27)
The next step is to eliminate $`Z`$ from eq. (26). This is achieved by multiplying the equation on the right with the vector $`B`$. This gives a linear equation for $`Z`$, whose solution is

$$Z=\frac{1}{1+\mu \alpha k}(C^1U_3^1\frac{1}{2}\alpha _1\alpha _2k)Y.$$
(28)
Substituting this back into eq. (26) and simplifying gives the formula

$$\{\mathrm{\Gamma }_+^1,C_3\}=C+\frac{1}{2}\frac{\alpha _1\alpha _2}{1+\mu \alpha k}CU_3^1YY^TCU_3^1.$$
(29)
In terms of components

$$(\mathrm{\Gamma }_+^1)_{mn}=\frac{m}{2\omega _m}\delta _{mn}+\frac{y(1y)(\omega _m\mu )(\omega _n\mu )Y_mY_n}{2[1\mu y(1y)k](\omega _m+\omega _n)},$$
(30)
where

$$\omega _m=\omega _{3m}=\sqrt{m^2+\mu ^2}.$$
(31)
This result be recast as a formula for the Neumann coefficient matrix $`\overline{N}_{mn}^{33}`$. The result is

$$\overline{N}_{mn}^{33}=\frac{mn\alpha _1\alpha _2}{1+\mu \alpha k}\frac{\overline{N}_m^3\overline{N}_n^3}{\omega _{3m}+\omega _{3n}}$$
(32)
where

$$\overline{N}_m^3=\left[(C^1C_3)^{1/2}U_3^1Y\right]_m.$$
(33)
Some further simple manipulations give the generalization

$$\overline{N}_{mn}^{rs}=\frac{mn\alpha }{1+\mu \alpha k}\frac{\overline{N}_m^r\overline{N}_n^s}{\alpha _s\omega _{rm}+\alpha _r\omega _{sn}}$$
(34)
where

$$\overline{N}_m^r=\left[(C^1C_r)^{1/2}U_r^1A^{(r)T}Y\right]_m.$$
(35)
This is the desired generalization of the flat-space formula eq. (18). However, we still require a generalization of eq. (19) as well as an explicit formula for $`k`$. Note that combining eq. (35) with eq. (11) gives

$$\overline{N}_{m0}^{rs}=\overline{N}_{0m}^{sr}=\sqrt{2\mu \alpha _s}ϵ^{st}\alpha _tU_r\overline{N}_m^rr=1,2,3s=1,2.$$
(36)
## 4 Determination of $`Y`$ and $`k`$

To complete the explicit determination of the Neumann matrices, and thus the three string vertex, we need useful formulas for

$$(\mathrm{\Gamma }_+^1)_{mn},Y_m=(\mathrm{\Gamma }_+^1B)_m,k=B^T\mathrm{\Gamma }_+^1B$$

as functions of $`y`$ and $`\mu `$. In view of of eq. (30), if we knew $`Y_m`$ and $`k`$ we would know $`(\mathrm{\Gamma }_+^1)_{mn}`$. The strategy that we will use for obtaining them is to derive first-order differential equations (in $`\mu `$) and input the known values at $`\mu =0`$ as initial conditions. Using the various definitions and identities in Appendix A, ref.  derived the differential equation

$$\frac{Y_m}{\mu }=\left[\frac{1}{2}\frac{F}{\mu }\left(1\frac{\mu }{\omega _m}\right)\frac{\mu }{\omega _m^2}\right]Y_m,$$
(37)
where

$$F(\mu ,y)=\mathrm{log}[1\mu y(1y)k(\mu ,y)].$$
(38)
The derivation of eq. (37) is sketched in Appendix B.

Eq. (37) has the solution

$$Y_m(\mu ,y)=\frac{m}{\omega _m}\mathrm{exp}\left[\frac{1}{2}_0^\mu \frac{F}{\mu }\left(1\frac{\mu }{\omega _m}\right)𝑑\mu \right]Y_m(0,y).$$
(39)
Thus, since $`Y_m(0,y)`$ is known, if we knew $`F(\mu ,y)`$, we would know $`k(\mu ,y)`$ and $`Y_m(\mu ,y)`$ and hence all the Neumann coefficients. Since, we have one fewer equations than unknowns, we need to input one additional piece of information. One that is easy to obtain and does the job is the asymptotic formula

$$Y_m\frac{m}{2\mu }B_m+O(\mu ^2)$$
(40)
for large $`\mu `$. Combining this condition with eq. (39) implies that

$$\mathrm{exp}\left\{\frac{1}{2}_0^{\mathrm{}}\frac{F}{\mu }\left(1\frac{\mu }{\omega _m}\right)𝑑\mu \right\}=\frac{B_m}{2Y_m(0)}.$$
(41)
Taking the logarithm, integrating by parts, and substituting the known value of the right-hand side gives

$$_0^{\mathrm{}}(m^2+\mu ^2)^{3/2}F(\mu ,y)𝑑\mu =G(m,y),$$
(42)
where

$$G(z,y)=\frac{2\tau _0}{z}+\frac{2}{z^2}\mathrm{log}\left(\frac{\mathrm{\Gamma }(1+z)}{\mathrm{\Gamma }(1+zy)\mathrm{\Gamma }(1+z(1y))}\right).$$
(43)
This formula must hold for $`m=1,2,\mathrm{}`$ and $`0y1`$. The inverse integral transform, which does not seem to exist in the mathematical literature, determines $`F(\mu ,y)`$. Ref.  proves that for a function $`G(m,y)`$ that is holomorphic in the right half $`m`$ plane and vanishes at infinity in that half plane, the inverse integral transform that solves the integral equation eq. (42) is

$$F(\mu ,y)=\frac{i\mu ^2}{\pi }_0^\pi \mathrm{cos}\theta G(i\mu \mathrm{cos}\theta ,y)𝑑\theta .$$
(44)
The proof is sketched in Appendix C. Using this, one can show that for our specific choice of $`G(m,y)`$ in eq. (43)

$$F(\mu ,y)=2\mu \tau _0+2\underset{r=1}{\overset{3}{}}\varphi (\mu \alpha _r),$$
(45)
where

$$\varphi (x)=\underset{n=1}{\overset{\mathrm{}}{}}\left[\mathrm{log}\left(\frac{\sqrt{n^2+x^2}+x}{n}\right)\frac{x}{n}\right].$$
(46)
We now have sufficiently explicit formulas to construct large $`\mu `$ asymptotic expansions.222This analysis was initiated in  and completed in . Large $`\mu `$ corresponds to small $`\lambda ^{}=1/\mu ^2`$, which is the effective coupling constant in the dual gauge theory . Differentiating eq. (46) twice, introducing a contour integral (Sommerfeld-Watson) representation of the series, and integrating by parts one can show that

$$\varphi ^{\prime \prime }(x)=x\underset{n=1}{\overset{\mathrm{}}{}}\frac{1}{(x^2+n^2)^{3/2}}=\frac{1}{x}+\frac{1}{2x^2}\pi _1^{\mathrm{}}\frac{zdz}{\sqrt{z^21}}\frac{1}{\mathrm{sinh}^2(\pi xz)}.$$
(47)
Integrating back, this allows us to deduce that

$$F(\mu ,y)=\mathrm{ln}[4\pi \mu y(1y)]+J(\mu y)+J(\mu (1y))J(\mu ),$$
(48)
where

$$J(x)=\frac{2}{\pi }_1^{\mathrm{}}\frac{\mathrm{ln}(1e^{2\pi xz})}{z\sqrt{z^21}}𝑑z.$$
(49)
Note that $`J\mathrm{exp}(2\pi x)`$ for large $`x`$.

## 5 Asymptotic behavior

We can now derive asymptotic results that include all inverse powers of $`\mu `$ and have leading corrections of order $`\mathrm{exp}[2\pi y\mu ]`$ (if $`y1/2`$). For example,

$$F(\mu ,y)\mathrm{ln}[4\pi \mu y(1y)].$$
(50)
The amazing thing is that a remarkably simply expression captures the result to all finite orders in $`1/\mu `$ and therefore one can easily read off predictions for the dual gauge theory that are valid (at this order of nonplanarity) to all orders in perturbation theory! The exponentially suppressed corrections encoded by the $`J`$ functions correspond to nonperturbative effects in the dual gauge theory. What these effects are, however, is mysterious. They don’t seem to be related to instantons or any other familiar nonperturbative phenomena.

Combining eqs. (50) and(38) one finds that

$$k(\mu ,y)\frac{1}{\mu y(1y)}\frac{1}{4\pi \mu ^2y^2(1y)^2}.$$
(51)
Similarly, substituting eq. (50) in eq. (39) gives

$$Y_m(\mu ,y)\sqrt{\frac{\mu +\omega _m}{2\mu }}\frac{m}{2\omega _m}B_m,$$
(52)
or in matrix notation

$$Y\frac{1}{2\sqrt{2\mu }}U_3^{1/2}C^{3/2}C_3^1B.$$
(53)
From this it follows that

$$\overline{N}_n^r\frac{(1)^{r(n+1)}}{2\pi y(1y)}\sqrt{\frac{|\alpha _r|}{2\mu n\omega _{rn}U_{rn}}}s_{rn},r\{1,2,3\}$$
(54)
where

$$s_{1m}=s_{2m}=1,s_{3m}=2\mathrm{sin}(\pi my).$$
(55)
The asymptotic expansions of the Neumann matrices are then given by substitution in eqs. (34) and (36).

## Acknowledgments

I am grateful to M. Spradlin and A. Volovich for helpful discussions. This work was supported in part by the U.S. Dept. of Energy under Grant No. DE-FG03-92-ER40701.

## Appendix A. Useful definitions and identities

The light-cone momenta in the minus direction that appear in the three-string vertex are defined to be $`\alpha _r/\alpha ^{}`$, $`r=1,2,3`$. For the process in which string #3 splits into strings #1 and #2, we take $`\alpha _1,\alpha _2>0`$, $`\alpha _3<0`$. Momentum conservation implies that $`\alpha _1+\alpha _2+\alpha _3=0`$. We also define the momentum fraction carried by string #1 to be $`y=\alpha _1/\alpha _3`$, which satisfies $`0<y<1`$. It follows that $`1y=\alpha _2/\alpha _3`$ is the momentum fraction carried by string #2. It is sometimes convenient to set $`\alpha _3=1`$, which can always be achieved by a suitable Lorentz boost. Only then does the mass parameter $`\mu `$ have an invariant meaning.

The matrices $`A_{mn}^{(r)}`$, which appear in the Neumann coefficients, are given by

$$A_{mn}^{(1)}=\frac{2}{\pi }(1)^{m+n+1}\sqrt{mn}\frac{y\mathrm{sin}(m\pi y)}{n^2m^2y^2},$$
(56)
$$A_{mn}^{(2)}=\frac{2}{\pi }(1)^m\sqrt{mn}\frac{(1y)\mathrm{sin}(m\pi y)}{n^2m^2(1y)^2},$$
(57)
and $`A_{mn}^{(3)}=\delta _{mn}`$. The indices $`m,n`$ range from 1 to infinity. Additional quantities that we need are

$$B_m=\frac{2\alpha _3}{\pi \alpha _1\alpha _2}(1)^m\frac{\mathrm{sin}(m\pi y)}{m^{3/2}}$$
(58)
and

$$C_{mn}=m\delta _{mn}.$$
(59)
These quantities are all $`\mu `$ independent and were defined already in the flat space analysis.

The infinite matrices $`A_{mn}^{(r)}`$ and the infinite vector $`B_m`$ satisfy a number of useful relations, which we record here

$$A^{(r)T}CA^{(s)}=\frac{\alpha _3}{\alpha _r}C\delta ^{rs}r,s=1,2$$
(60)
$$A^{(r)T}C^1A^{(s)}=\frac{\alpha _r}{\alpha _3}C^1\delta ^{rs}r,s=1,2.$$
(61)
The symbol $`T`$ means matrix transpose. Some additional useful identities are

$$\underset{r=1}{\overset{3}{}}\frac{1}{\alpha _r}A^{(r)}CA^{(r)T}=0,$$
(62)
$$\underset{r=1}{\overset{3}{}}\alpha _rA^{(r)}C^1A^{(r)T}=\frac{\alpha }{2}BB^T,$$
(63)
where we have introduced

$$\alpha =\alpha _1\alpha _2\alpha _3.$$
(64)
Additional matrices that involve the mass parameter $`\mu `$ of the plane-wave geometry, introduced in , are

$$(C_r)_{mn}=\omega _{rm}\delta _{mn}=\sqrt{m^2+\mu ^2\alpha _r^2}\delta _{mn}$$
(65)
and

$$U_r=C^1(C_r\mu \alpha _r).$$
(66)
Note that

$$(U_r)^1=C^1(C_r+\mu \alpha _r)=U_r+2\mu \alpha _rC^1.$$
(67)
A crucial construct is the infinite matrix

$$\mathrm{\Gamma }_+=\underset{r=1}{\overset{3}{}}A^{(r)}U_rA^{(r)T}.$$
(68)
Explicit formulas for the inverse of $`\mathrm{\Gamma }_+`$ are a main goal of our work. Related quantities that also are needed are the infinite vector

$$Y=\mathrm{\Gamma }_+^1B$$
(69)
and the scalar

$$k=B^T\mathrm{\Gamma }_+^1B.$$
(70)
## Appendix B. Derivation of the differential equation

This appendix sketches the derivation of eq. (37), which we copy here

$$\frac{Y_m}{\mu }=\left[\frac{1}{2}\frac{F}{\mu }\left(1\frac{\mu }{\omega _m}\right)\frac{\mu }{\omega _m^2}\right]Y_m.$$
(71)
The matrix $`\mathrm{\Gamma }_+=_rA^{(r)}U_rA^{(r)\mathrm{T}}`$, introduced in Appendix A, only depends on $`\mu `$ through the dependence of $`U_r`$ on $`\mu `$. Its derivative can be written in the form

$$\frac{\mathrm{\Gamma }_+}{\mu }=\frac{1}{2}\alpha BB^\mathrm{T}+\mu N,$$
(72)
where

$$N=\underset{r=1}{\overset{3}{}}\alpha _r^2A^{(r)}C^1C_r^1A^{(r)\mathrm{T}}.$$
(73)
It follows that

$$\frac{Y}{\mu }=\frac{1}{2}k\alpha Y\mu \mathrm{\Gamma }_+^1NY.$$
(74)
The product $`NY`$ can be recast in the form

$$NY=g_1C_3^2B+g_2B,$$
(75)
where we define the coefficients $`g_1`$ and $`g_2`$ to be the scalar quantities

$$g_1=\frac{2(1+\mu \alpha k)}{2+\mu \alpha k+\mu ^2\alpha k_1},$$
(76)
$$g_2=\left(\frac{\alpha }{2}\right)\frac{\alpha k^2+\mu \alpha kk_1+2k_1}{2+\mu \alpha k+\mu ^2\alpha k_1},$$
(77)
and

$$k_i=B^\mathrm{T}C_3^iY.$$
(78)
The above equations imply that

$$\frac{k}{\mu }=B^\mathrm{T}\frac{Y}{\mu }=\frac{1}{2}\alpha k^2\mu g_2k\mu g_1k_2.$$
(79)
This is not very useful as it stands, since there is no other apparent way to determine $`k_2`$. ($`k_1`$ could be determined, but that will turn out not to be necessary.) Substituting the equation for $`NY`$ and an identity for $`[C_3^2,\mathrm{\Gamma }_+^1]`$ deduced from eq. (29), one can recast eq. (74) in the form

$$\frac{Y}{\mu }=(F_0+F_1C_3^1+F_2C_3^2)Y,$$
(80)
where the scalar functions $`F_i`$ are given by

$$F_0=\frac{1}{2}\alpha k\mu g_2+\frac{1}{2}\mu g_1\frac{\alpha }{1+\mu \alpha k}(k_1\mu k_2),$$
(81)
$$F_1=\frac{1}{2}\mu g_1\frac{\alpha }{1+\mu \alpha k}(k\mu ^2k_2),$$
(82)
$$F_2=\mu g_1+\frac{1}{2}\mu ^2g_1\frac{\alpha }{1+\mu \alpha k}(k\mu k_1).$$
(83)
Using the equations above to eliminate $`g_1`$, $`g_2`$, $`k_1`$, and $`k_2`$, we find

$$F_2=\mu ,F_1=\mu F_0,F_0=\frac{\alpha }{2}\frac{1}{1+\mu \alpha k}(k+\mu k^{}).$$
(84)
This allows us to rewrite eq. (80) in the desired form eq. (71).

## Appendix C. The integral transform

The analysis in Sec. 4 of the text required solving the following integral equation for $`f(x)`$

$$g(z)=_0^{\mathrm{}}\frac{f(x)}{(x^2+z^2)^{3/2}}𝑑x,$$
(85)
where $`g(z)`$ is a given function that is holomorphic in the right-half plane and vanishes at infinity in that half plane. We claim  that the unique solution is

$$f(x)=i\frac{x^2}{\pi }_0^\pi g(ix\mathrm{cos}\theta )\mathrm{cos}\theta d\theta .$$
(86)
The proof that elimination of $`f`$ from eqs. (85) and (86) gives $`g=g`$ is a consequence of elementary integration. This proves the existence of a solution. The proof of uniqueness requires that elimination of $`g`$ should give $`f=f`$. After making changes of variables and deforming integration contours, one can argue that this requires the identity

$$\delta (yy^{})=\frac{\sqrt{y}}{4\pi i}_𝒞\frac{dw}{\sqrt{1+w}}\frac{1}{(wy+y^{})^{3/2}},$$
(87)
where $`y,y^{}>0`$. The contour can be taken to be the unit circle $`|z|=1`$, in the counterclockwise sense, starting and ending at the point $`z=1`$. It is an elementary application of Cauchy’s theorem to show that this integral vanishes for $`y<y^{}`$ and for $`y^{}<y`$. That it is has the right singularity at $`y=y^{}`$ can be verified by showing that the Laplace transform of both sides give an identity.


format: text

text:
# Guiding optical flows by photonic crystal slabs made of dielectric cylinders

## I Introduction

When propagating through periodically structured media such as photonic crystals (PCs), optical waves will be modulated with the periodicity. As a result, the dispersion of waves will no longer behave as in a free space, and so called frequency band structures appear. Under certain conditions, waves may be prohibited from propagation in certain or all directions, corresponding to partial and complete bandgaps respectively. The photonic crystals revealing bandgaps are called bandgap materials.

Photonic crystals and band gap materials have a broad spectrum of applications, ranging from computing to digital communication and from laser cavities to optical transistorsBook . The possibilities are unlimited. In fact, applications have well gone beyond expectation, and are so far reaching that a fruitful new field called photonic crystals has come into existence. Most updated information about the research of photonic crystals and related materials can be found in the comprehensive webpage web .

So far, most applications are associated with the properties of the complete bandgaps of PCs. On one hand, the bandgaps confine optical propagation within certain frequency regimes. On the other, when encountering the complete bands, optical waves can be guided into desired directions. For example, one of the main applications of PCs is to control optical flows, so that they can be used for such as telecommunications. A comprehensive survey of phonic crystal research can be referred to Refs. Book ; web ; Mono1 ; Mono2 ; Mono3 . To our knowledge, however, there have been very few attempts in the literature to explore possible usage of partial bandgaps. In this paper, we wish to discuss a previously undiscussed phenomenon associated with partial bandgaps, that is, deflection of optical waves. That is, the partial bandgap can collimate wave propagation into certain directions. This property may allow for novel applications in manipulating optical flows.

The paper is organized as follows. The systems and the theory will be outlined in the next section. The results and discussion will be presented in Section III, followed by a short summary.

## II The systems and formulation

The systems considered here are two dimensional photonic crystals made of arrays of parallel dielectric cylinders placed in a uniform medium, which we assume to be air. Such systems are common in both theoretical simulations or experimental measurements of two dimensional PCsBook ; web . For brevity, we only consider the E-polarized waves (TM mode), that is, the electric field is kept parallel to the cylinders. The following parameters are used in the simulation. (1) The dielectric constant of the cylinders is 14, and the cylinders are arranged to form a square lattice. (2) The lattice constant is $`a`$ and the radius of the cylinders is 0.3$`a`$; in the computation, all lengths are scaled by the lattice constant. (3) The unit for the angular frequency is $`2\pi c/a`$. After scaling, the systems become dimensionless; thus the features discussed here would be applicable to a wider range of situations.

While the frequency band structure in the systems can be calculated by the plane-wave expansion methodBook , the propagation and scattering of electromagnetic (EM) waves in such systems can be studied by the standard multiple scattering theory. The theory originated from the self-consistent idea first discussed by FoldyFoldy , and then made maturity through the significant efforts by LaxLax , Waterman et al.Waterman , and particularly by TwerskyTwersky .

The essence of the theory is summarized as follows. In response to the incident wave from the source and the scattered waves from other scatterers, each scatter will scatter waves repeatedly, and the scatterered waves can be expressed in terms of a modal series of partial waves. When this scattered wave serves as an incident wave to other scatterers, a set of coupled equations can be formulated and computed rigorously. The total wave at any spatial point is the summation of the direct wave from the source and the scattered waves from all scatterers. The intensity of the wave is represented by the square of the wave field.

For the reader’s convenience we present briefly the general multiple scattering theory. Consider that $`N`$ straight cylinders of radius $`a^i`$ located at $`\stackrel{}{r_i}`$ with $`i=1,2,\mathrm{},N`$ to form an array. A line source transmitting monochromatic waves is placed at $`\stackrel{}{r}_s`$. Here we take the standard approach with regard to the source. That is, the transmission from the source is calculated from the multiple scattering theory, and assume that the source is not affected by the surroundings. If some other sources such as a line of atoms are used, the reaction between the source and the backscattered waves should be taken into account.

The scattered wave from each cylinder is a response to the total incident wave composed of the direct wave from the source and the multiply scattered waves from other cylinders. The final wave reaching a receiver located at $`\stackrel{}{r}_r`$ is the sum of direct wave from the source and the scattered waves from all the cylinders.

The scattered wave from the $`j`$-th cylinder can be written as

$$p_s(\stackrel{}{r},\stackrel{}{r_j})=\underset{n=\mathrm{}}{\overset{\mathrm{}}{}}i\pi A_n^jH_n^{\left(1\right)}\left(k\left|\stackrel{}{r}\stackrel{}{r_j}\right|\right)e^{in\varphi _{\stackrel{}{r}\stackrel{}{r_j}}},$$
(1)
where $`k`$ is the wavenumber in the medium, $`H_n^{\left(1\right)}`$ is the $`n`$-th order Hankel function of first kind, and $`\varphi _{\stackrel{}{r}\stackrel{}{r_j}}`$ is the azimuthal angle of the vector $`\stackrel{}{r}\stackrel{}{r_j}`$ relative to the positive $`x`$ axis. The total incident wave around the $`i`$-th cylinder $`\left(i=1,2,\mathrm{},N;ij\right)`$ is the summation of the direct incident wave from the source and the scattered waves from all other scatterers, can be expressed as

$$p_{in}^i\left(\stackrel{}{r}\right)=\underset{n=\mathrm{}}{\overset{\mathrm{}}{}}B_n^iJ_n\left(k\left|\stackrel{}{r}\stackrel{}{r_i}\right|\right)e^{in\varphi _{\stackrel{}{r}\stackrel{}{r_i}}}.$$
(2)
In this paper, $`p`$ stands for the electrical field in the TM mode and the magnetic field in the TE mode.

The coefficients $`A_n^i`$ and $`B_n^i`$ can be solved by expressing the scattered wave $`p_s(\stackrel{}{r},\stackrel{}{r_j})`$, for each $`ji`$, in terms of the modes with respect to the $`i`$-th scatterer by the addition theorem for Bessel function. Then the usual boundary conditions are matched at the surface of each scattering cylinder. This leads to

$$B_n^i=S_n^i+\underset{j=1,ji}{\overset{N}{}}C_n^{j,i},$$
(3)
with

$$S_n^i=i\pi H_n^{\left(1\right)}\left(k\left|\stackrel{}{r_i}\right|\right)e^{in\varphi _{\stackrel{}{r}_i}},$$
(4)
and

$$C_n^{j,i}=\underset{l=\mathrm{}}{\overset{\mathrm{}}{}}i\pi A_l^jH_{ln}^{\left(1\right)}\left(k\left|\stackrel{}{r_i}\stackrel{}{r_j}\right|\right)e^{i\left(ln\right)\varphi _{\stackrel{}{r_i}\stackrel{}{r_j}}},$$
(5)
and

$$B_n^i=i\pi \tau _n^iA_n^i,$$
(6)
where $`\tau _n^i`$ are the transfer matrices relating the properties of the scatterers and the surrounding medium and are given as

$$\tau _n^i=\frac{H_n^{\left(1\right)}\left(ka^i\right)J_n^{}\left(ka^i/h^i\right)g^ih^iH_{n}^{\left(1\right)}{}_{}{}^{}\left(ka^i\right)J_n\left(ka^i/h^i\right)}{g^ih^iJ_n^{}\left(ka^i\right)J_n\left(ka^i/h^i\right)J_n\left(ka^i\right)J_n^{}\left(ka^i/h^i\right)},$$
(7)
where

$$h^i=\frac{1}{\sqrt{ϵ^i}},\text{and}g^i=\{\begin{array}{cc}ϵ^i\hfill & \text{for TE waves}\hfill \\ 1\hfill & \text{for TM waves}\hfill \end{array},$$

in which $`ϵ^i`$ is the dielectric constant ratio between the $`i`$-th scatterer and the surrounding medium.

The coefficients $`A_n^i`$ and $`B_n^j`$ can then be inverted from Eq. (3). Once the coefficients $`A_n^i`$ are determined, the transmitted wave at any spatial point is given by

$$p\left(\stackrel{}{r}\right)=p_0\left(\stackrel{}{r}\right)+\underset{i=1}{\overset{N}{}}\underset{n=\mathrm{}}{\overset{\mathrm{}}{}}i\pi A_n^iH_n^{\left(1\right)}\left(k\left|\stackrel{}{r}\stackrel{}{r_i}\right|\right)e^{in\varphi _{\stackrel{}{r}\stackrel{}{r_i}}},$$
(8)
where $`p_0`$ is the field when no scatterers are present. The transmitted intensity field is defined as $`\left|p\right|^2`$.

## III Results and discussion

The frequency band structure is plotted in Fig. 1, and the qualitative features are similar to that obtained for a square array of alumina rods in air. A complete band gap is shown between frequencies of 0.22 and 0.28. Just below the complete gap, there is a regime, sandwiched by two horizontal lines, of partial band gap in which waves are not allowed to travel along the $`\mathrm{\Gamma }X`$ or  direction. We will consider waves whose frequency is within this partial bandgap. In particular, we choose the frequency to be 0.192.

First we consider the propagation of EM waves through two rectangular slabs of arrays of dielectric cylinders. Figure 2 shows the images of the fields. The left panel shows the real parts of the fields $`E_z`$, while the right panel presents the images of the intensity fields $`\left|E_z\right|^2`$. In (a1) and (a2), the slab measures 14x45, and the slab is oriented such that the  direction, i. e. the $`\mathrm{\Gamma }M`$ direction, is along the horizontal level. The size of the slab in (b1) and (b2) is 10x45, and the  direction is titled upwards and makes an angle of 22.5 degree with respect to the horizontal direction. The frequency is chosen at 0.192. A transmitting point source is placed at 2 lattice constant away from the left side of the slabs. The detailed geometrical information can be referred to in Fig. 2.

A few observations can be drawn from Fig. 2. First, there is a focused image across the slab in (a1) and (a2). Earlier, this focused image was attributed to the effect of negative refractionYe1 , inferred from the group velocity calculation. If this conjecture were valid, another focused image would be expected inside the slab as well. Our result does not support this conjecture. As seen from (a1) and (a2), there is no focused image inside the slab. Rather, the waves are mostly confined in a tunnel and travel to the other side of the slab, then release to the free space. This is understandable, because the forbidden direction in (a1) and (a2) is along $`\mathrm{\Gamma }X`$, which makes an angle of 45 degree from the $`\mathrm{\Gamma }M`$ direction that lies horizontally. The passing band in the $`\mathrm{\Gamma }M`$ direction thus acts as a transportation carrier that moves the source to the other side of the slab. The waves on the right hand side of the slab look as if they were radiated by an image that has been transported across the slab within a narrow guide. Second, the waves tend to bend towards the $`\mathrm{\Gamma }M`$ direction, as evidenced by Fig. 2 (b1) and (b2). Third, the decay of the transported intensity along the travelling path is not obvious, an indication of efficient guided propagation.

The results in Fig. 2 are promising. They show that in the presence of partial bandgaps and when incident upon a slab of photonic crystals, waves tend to bend toward directions which are mostly away from forbidden directions. This would indicate that partial bandgaps may be considered as a candidate for guiding wave flows. To verify this conjecture, we have further explored the guiding phenomenon associated with partial bandgaps.

In Fig. 3, we show the EM wave propagation through stacks of photonic crystal slabs. Two setups are considered. In (a1) and (a2), two slabs of dielectric cylinders are stacked together. The first (left) slab is oriented such that the  direction is horizontal, while the second (right) slab is arranged to make the  direction tilted upward, making 22.5 degree with respect to the horizontal direction. The two slabs measure as 9x44 and 14x44 respectively. In (b1) and (b2), two slabs are adjacently attached. The  direction is tilted upward by 10 degree for the first (left) slab, while it is along the horizontal direction in the second (right) slab. The sizes of the two slabs are 8x40 and 10x40 respectively. In both situations (a) and (b), the point source is placed at a distance of 1.5 away the left side of the stacks. The purpose here is to show how the light would travel when two adjacent slabs have different orientations.

Here it is clearly shown that the EM waves indeed always tend to travel along the  direction. In the case of (b1) and (b2), an image has been formed on the right hand side. Compared to the source, the image is uplifted by a distance about $`8\mathrm{tan}\left(\pi /18\right)`$. In the cases considered here, the first slabs (left) serve as a collimating device, and then the collimated waves will be guided by subsequent photonic crystal slabs. This consideration can be extended to multiple consecutive slabs so that the wave flows can be guided into desired orientations, making possible alternative ways of controlling optical flows.

In Fig. 4, we consider two other situations of stacked photonic crystal slabs. The geometrical parameters are indicated in the figure. Again, the waves tend to move along the  direction. Here the amphoteric diffraction is observed. It draws analogy with the amphoteric refraction observed when waves propagate from an isotropic to an anisotropic mediumYau .

The results from Figs. 2, 3 and 4 clearly indicate that the partial bandgaps can be indeed used as a guiding channel for optical flows. It can be also inferred that the guided transport is efficient. We have carried out further simulations against variations of frequencies, filling factors, and dielectric constants, the results are quantitatively the same for waves within partial bandgaps. The observation presented here has also been confirmed by FDTD simulations. The controlled wave transport due to partial bandgaps of PCs should be interpreted in terms of diffraction or scattering rather than refraction; in fact, no refraction index can be determined for the phenomenon.

An immediate question may thus arise: Why the waves of frequencies within the partial bandgap tend to bend to particular directions? To answer this question, we have examined the properties of the energy flow of the eigenmodes which correspond to the frequency bands. While details will be published elsewhere, here we only outline our thoughts. The usual approach mainly relies on the curvatures of frequency bands to infer the energy flow. As documented in Ref. Yeh , an energy velocity is defined as $`\stackrel{}{v}_e=\frac{\frac{1}{V}{\scriptscriptstyle \stackrel{}{J}_\stackrel{}{K}d^3r}}{\frac{1}{V}{\scriptscriptstyle U_\stackrel{}{K}d^3r}},`$ where $`\stackrel{}{J}_\stackrel{}{K}`$ and $`U_\stackrel{}{K}`$ are the energy flux and energy density of the eigenmodes, and the integration is performed in a unit cell. It can be shown that thus defined energy velocity equals the group velocity obtained as $`\stackrel{}{v}_g=_\stackrel{}{K}\omega \left(\stackrel{}{K}\right).`$ Therefore it is common to calculate the group velocity to infer the energy velocity and subsequently the energy flows or refraction of waves. A few questions, however, may arise with regard to this approach. First, when the variation in the Bloch vector, i. e. $`\delta \stackrel{}{K}`$, is small, the changes in $`\omega `$, $`\stackrel{}{E}_\stackrel{}{K}`$ and $`\stackrel{}{H}_\stackrel{}{K}`$ should also be small. Second, the variation operation should be exchangeable with the partial differential operations. When these two conditions fail, the energy velocity will become illy defined. Third, even if the two conditions hold, whether the net current flow through a unit cell really follows the direction of $`\stackrel{}{v}_e`$ remains unclear. We note here that the average flux through a surface may be defined as $`\stackrel{}{J}=\frac{\widehat{n}}{S}𝑑\stackrel{}{S}\stackrel{}{J}`$, where $`\widehat{n}`$ is the unit normal vector of the surface $`S`$. Clearly, the volume averaged current within a unit cell does not necessarily correspond to the actual current flow. We will publish verifications elsewhere.

To avoid possible ambiguities, here we consider the energy flow based upon its genuine definition. One advantage of this approach is that we are also able to examine the local properties of energy flows. By Bloch’s theorem, the eigenmodes corresponding to the frequency bands of PCs can be expressed as $`E_\stackrel{}{K}\left(\stackrel{}{r}\right)=e^{i\stackrel{}{K}\stackrel{}{r}}u_\stackrel{}{K}\left(\stackrel{}{r}\right)`$, where $`\stackrel{}{K}`$ is the Bloch vector, as the wave vector, and $`u_\stackrel{}{K}\left(\stackrel{}{r}\right)`$ is a periodic function with the periodicity of the lattice. When expressing $`E_\stackrel{}{K}\left(\stackrel{}{r}\right)`$ as $`\left|E_\stackrel{}{K}\left(\stackrel{}{r}\right)\right|e^{i\theta _\stackrel{}{K}\left(\stackrel{}{r}\right)}`$, the corresponding energy flow is derived as $`\stackrel{}{J}_\stackrel{}{K}\left(\stackrel{}{r}\right)\left|E_\stackrel{}{K}\left(\stackrel{}{r}\right)\right|^2\theta _\stackrel{}{K}\left(\stackrel{}{r}\right)`$; clearly $`\theta _\stackrel{}{K}`$ combines the phase from the term $`e^{i\stackrel{}{K}\stackrel{}{r}}`$ and the phase from the function $`u_\stackrel{}{K}\left(\stackrel{}{r}\right)`$. To explore the characteristics of the partial bandgap, we have computed the eigen-field $`E_\stackrel{}{K}\left(\stackrel{}{r}\right)`$ and also the energy flow of the eigenmodes. The results are shown in Fig. 5. Fig. 5(a) shows that the energy eventually tends to flow into the direction of $`\mathrm{\Gamma }M`$, i. e. the  direction, while the Bloch vector points to an angle of 22.5o that lies exactly between $`\mathrm{\Gamma }X`$ and $`\mathrm{\Gamma }M`$.

Although the above features are only investigated for the first partial bandgap in this paper, we have found that they are also valid for other partial bandgaps. For example, we have also considered the second partial bandgap which is located between 0.283 and 0.325. All above features remain quantitatively valid. Within this second gap, however, the waves are collimated to travel along the  direction rather than the  direction. In addition, we have also carried out simulations for various slab sizes, all the features are the same, thus excluding the boundaries as the possible cause.

## IV Summary

We have considered EM wave propagation through slabs of photonic crystals which are made of arrays of dielectric cylinders. Properties of partial bandgaps are investigated. It was shown that the partial bandgaps may act as a guiding channel for wave propagation inside the photonic crystals. Such a feature may lead to novel applications in manipulating optical flows.

## Acknowledgments

The work received support from National Science Council, National Central University, and Fudan University.

## Figure Captions

The band structure of a square lattice of dielectric cylinders. The lattice constant is $`a`$ and the radius of the cylinders is $`0.3a`$. $`\mathrm{\Gamma }M`$ and $`\mathrm{\Gamma }X`$ denote the  and  directions respectively. A partial gap is between the two horizontal lines.

Imaging of the transmitted fields across two slabs of dielectric cylinders. The black circles denote the cylinders (for clarity, not all cylinders are plotted).

Imaging of the intensity fields across two-consecutive slabs of arrays of dielectric cylinders in two arrangements. The left and right panels respectively plot the real part of the field and the intensity.

Imaging of the intensity fields across two-consecutive slabs of arrays of dielectric cylinders in two arrangements. The left and right panels respectively plot the real part of the field and the intensity. Here is show the amphoteric diffraction at the interfaces between two adjacent slabs: (a) positive and (b) negative. In both (a) and (b), the adjacent slabs measure as 10x50 and 12x50.

Left panel: the field pattern of eigenmodes. Right panel: the energy flow of the eigenmodes. The eigenmodes along two directions are considered: (a) $`\stackrel{}{K}=(0.9\pi /a,0.37\pi /a)`$, i. e. along an angle of 22.5o exactly between $`\mathrm{\Gamma }X`$ and $`\mathrm{\Gamma }M`$ directions; the corresponding frequency is 0.185; (b) $`\stackrel{}{K}=(0.7\pi /a,0.7\pi /a)`$, i. e. along $`\mathrm{\Gamma }M`$; the corresponding frequency is 0.192. The direction of the Bloch vectors are denoted by the blue arrows, while the red arrows denote the local energy flow including the direction and the magnitude. The circles refer to the cylinders. Both frequencies in (a) and (b) lie within the partial gap. Due to the periodicity, we only plot the energy flow within one unit cell. Note that although the features shown by (a) also hold for other Bloch vectors for which the corresponding frequencies lie within the partial gap regime, we only plot here for the case of 22.5o.


text:
# Scaling of the anomalous Hall effect in Sr1-xCaxRuO3

## Abstract

The anomalous Hall effect (AHE) of ferromagnetic thin films of Sr1-xCaxRuO3 (0 $`x`$ 0.4) is studied as a function of $`x`$ and temperature $`T`$. As $`x`$ increases, both the transition temperature $`T_c`$ and the magnetization $`M`$ are reduced and vanish near $`x`$ 0.7. For all compositions, the transverse resistivity $`\rho _H`$ varies non-monotonously with $`T`$, and even changes sign, thus violating the conventional expression $`\rho _H=R_oB+4\pi R_sM(T)`$ ($`B`$ is the magnetic induction, while $`R_o`$ and $`R_s`$ are the ordinary and anomalous Hall coefficients). From the rather complicated data of $`\rho _H`$, we find a scaling behavior of the transverse conductivity $`\sigma _{xy}`$ with $`M(T)`$, which is well reproduced by the first-principles band calculation assuming the intrinsic origin of the AHE.

It has been known that the Hall resistivity $`\rho _H`$ in a ferromagnet Hall  has some extra contribution originated from the spontaneous magnetization, which is assumed and often observed experimentally to be fitted by $`\rho _H=R_0B+4\pi R_sM`$, where $`B`$ is the magnetic induction and $`M`$ is the magnetization of the material. $`R_0B`$ represents the ordinary Hall contribution, which is related to the nature and amount of charge carriers. It is a linear function of the applied magnetic field $`H`$ as in the Hall measurement geometry, $`B`$=$`H`$. $`R_sM`$ is referred to as the anomalous Hall term, and is usually associated with the spin polarization of the conduction carriers and the relativistic spin-orbit interaction. According to the above definition of $`\rho _H`$, the anomalous Hall term is proportional to the magnetization of the material. However the quantitative analysis of the AHE has seldom been completed since its mechanism has not yet been established, and theories give much smaller values compared with the experiments. Furthermore, most of the theories regard the AHE as from extrinsic origins, involving processes such as skew scatterings skew  and side-jump mechanismsside . Therefore the magnitude of the AHE depends on the concentration and scattering strength of impurities, thermal spin-agitation, etc. In contrast to these extrinsic mechanisms, several workskl ; ye ; nagaosa ; salamon ; jungwirth  regard AHE of the intrinsic origin. Namely the phase of the Bloch wave-function in the momentum space determines the Hall conductivity $`\sigma _{xy}`$, which is largely determined by the band crossing points acting as “magnetic monopoles”Fang . As an explicit example, the AHE and magneto-optical effect of SrRuO3 have been studied, and a good agreement was obtained between theory and experimentsFang . In the present paper, we report the extensive study of the AHE in Sr1-xCaxRuO3 (0 $`x`$ 0.4) as a function of $`T`$ and $`x`$ to reveal its systematics. We have found the scaling behavior of the transverse conductivity $`\sigma _{xy}`$ in terms of the $`T`$\- and $`x`$-dependence of the magnetization $`M(x,T)`$, which is in fairly good agreement with the first-principles band calculation. This gives a firm evidence for the intrinsic origin of the AHE.

Several members of the Ruddelson-Popper-type Srn+1RunO3n+1 series show metallic properties, as well as superconductivity and magnetic orderSROxyz . SrRuO3 ($`n=\mathrm{}`$, perovskite) is ferromagneticSRO1 ; SRO2  with a Curie temperature ($`T_c`$) around 160 K, and a fairly large spin orbit coupling. The $`4d`$ orbitals of Ru4+ are rather extended, and the Coulomb repulsion is small compared to the band width. The ferromagnetic properties of SRO are usually associated with a narrow itinerant band resulting from the hybridization between the Ru($`t_{2g}`$) and O($`2p`$) orbitals.

While similar structurally, and as well metallic, CaRuO3 (CRO) does not exhibit ferromagnetism and its magnetic state is still under discussionCRO . In Sr1-xCaxRuO3, the ferromagnetic interaction becomes weaker with increasing $`x`$SCRO . For the compounds with larger Ca-doping ($`x`$ $``$ 0.7), no clear phase transition is discerned, and only some irreversibility is observed in the magnetization curves of these materials. The disappearance of the long range magnetic order is commonly related to the distortion of the RuO6 octahedra associated with the partial or total replacement of Sr by Ca, and the corresponding narrowing of the $`4d`$ bandwidthSCRO2 .

It is difficult to grow clean single crystals of Ca-doped SRO, although it is possible to prepare high quality single crystals of the end compounds, SRO and CRO. On the other hand, it is nowadays possible to grow high quality epitaxial films of SrRuO3 and Sr1-xCaxRuO3eom . In the present article, we study the anomalous Hall resistivity of epitaxial films of Sr1-xCaxRuO3 (0 $`x`$ 0.4), and its evolution as the ferromagnetic interaction decreases with $`x`$. The Hall resistivity contains an anomalous component, associated with the ferromagnetic ordering of the samples at low temperatures. This component is not simply proportional to the $`M`$, as usually assumed. The results reveal a close relation between $`\sigma _{xy}`$ and the spin polarization of the system. Such a correlation was predicted by first-principles calculations taking account of the spin-orbit interaction in terms of the Berry phase connectionFang . The calculations successfully reproduce the non-monotonous variation of $`\sigma _{xy}`$ with temperature (via its magnetization), as well as its sign change. The effects of disorder and structural changes on the anomalous conductivity are discussed.

Thin ($``$ 500 Å) films of Sr1-xCaxRuO3 ($`x`$ = 0, 0.1, 0.2, 0.3, and 0.4) were epitaxially grown on the (001) surfaces of high quality SrTiO3 single-crystal substrateskawasaki  by pulsed laser deposition (PLD). Bulk single crystal SrRuO3 and CaRuO3 were prepared using a flux method for comparison. The quality of the films and phase-purity of the single crystals were confirmed by x-ray diffraction. Both single crystals have orthorhombic structure; SrRuO3 has a relatively small orthorhombicity ($`c/a`$=1.003, $`c/b`$=0.996), while it is larger for CaRuO3 ($`c/a`$=1.01, $`c/b`$=0.980). The epitaxial thin films are coherently strained by the SrTiO3 substrate, yielding a tetragonal distortion in the  direction ($`c/a`$=$`c/b`$=1.01 for $`x`$=0). As a result, the out-of-plane lattice constants of the films are elongated (c.f. inset of Fig. 1), and due to the spin-orbit interaction, the easy axis of magnetization is perpendicular to the film planestrain2 . Magnetic and transport measurements were performed on the Sr1-xCaxRuO3 thin films and the single crystals of SRO and CRO. The magnetization data was recorded on a MPMS5S SQUID magnetometer using a magnetic field applied normal to the plane of the films. The films were then patterned in a six-lead Hall bar geometry using conventional photo-lithography and Ar ion etching for transport measurements. The Hall resistivity $`\rho _H`$ was measured with a PPMS6000 system together with the longitudinal resistivity $`\rho _{xx}=\rho `$ as a function of $`H`$ and $`T`$. The anomalous resistivity $`\rho _{xy}`$ was extrapolated to $`H`$ = 0 from $`\rho _H`$ vs $`H`$ measurements up to $`H`$ = $`\pm `$ 9 T at constant temperatures (from 2 K to 200 K) after subtraction of the ordinary Hall contribution, and the transverse conductivity $`\sigma _{xy}`$ was estimated as -$`\rho _{xy}/\rho _{xx}^2`$. A small (as the patterned leads are nearly symmetric) magnetoresistance was removed by subtracting $`\rho _H`$($`H`$) to $`\rho _H`$($`H`$). First-principles calculations of $`\sigma _{xy}`$ were performed assuming orthorhombic and cubic crystal structures. The plane-wave pseudo-potential calculations were performed based on the local spin density approximation (LSDA), and the spin-orbit coupling was treated self-consistently by using the relativistic fully separable pseudo-potentials in the framework of non-collinear magnetism formalism. The finite life-time broadening was estimated from the experimental residual resistivity and the extended Drude analysis of the longitudinal conductivityFang .

All the investigated samples are metallic as seen in Fig. 1 and inset. The residual resistivity at low temperatures is very low ($``$ 2.5 $`\mu \mathrm{\Omega }`$cm) for the single crystals, and increases for the films from 26 ($`x`$ = 0) to 232 $`\mu \mathrm{\Omega }`$cm ($`x`$ = 0.4). A kink is observed in the resistivity curves, around the paramagnetic-to-ferromagnetic transition temperatures ($`T_c`$) of the films. No long range magnetic order is observed in CRO, and as seen in inset, the resistivity curve has no anomaly in the measured range of temperature. Due to the above mentioned strain effects, the epitaxial film of SRO has a slightly lower $`T_c`$ than the single crystal, near 150 Kstrain2 ; strain . The substitution of Sr by Ca in Sr1-xCaxRuO3 weakens the ferromagnetic interaction, and $`T_c`$ is greatly reduced. It is reduced to $``$ 110 K for $`x`$ = 0.2, and $``$ 70 K for $`x`$ = 0.4.

The top panel of Fig. 2 shows the temperature dependence of the 7 T-cooled magnetization of the films probed in a small field of $`H`$ = 0.05 T applied normal to the plane, i.e. in the same direction as in the Hall measurements, along the easy axis of magnetization. As seen in inset, the coercivity $`H_c`$ of the films monotonously increases with increasing Ca doping. In addition to possible disorder, the Ca substitution yields structural changes, affecting the magnetic properties and, as we will discuss below, the AHE. At low temperatures, $`M`$($`H`$ = 7 T) amounts to $``$ 1.5 $`\mu _B`$/Ru for the undoped SRO film. Due to the itinerant character of the magnetism, the obtained moment is smaller than expected according to Hund’s rule (2 $`\mu _B`$/Ru; $`S`$=1). $`M`$($`H`$ = 7T) decreases with increasing Ca doping, and amounts to $``$ 0.75 $`\mu _B`$/Ru for $`x`$ = 0.4. The middle panel of Fig. 2 shows the temperature dependence of the anomalous Hall resistivity $`\rho _{xy}`$ for all the films and the SRO single crystal. At a constant temperature, $`\rho _H`$ is a linear function of $`H`$ at high fields, with a negative proportionality constant ($`R_o<`$ 0), indicating charge carriers of electron-like nature.

As seen in the figure, $`\rho _{xy}`$ of the films varies non-monotonously with $`T`$. No AHE is observed at high temperatures. The AHE appears just above $`T_c`$; $`\rho _{xy}`$ is negative, and show a maximum near $`T_c`$. At lower temperatures, $`\rho _{xy}`$ changes sign (near 120 K for $`x`$ = 0, and 60 K for $`x`$ = 0.2), and remain positive down to the lowest temperature. For $`x`$ = 0.4 with a lower $`T_c`$, however, the sign change is not observed. Similar features are observed for the single crystal of SRO; the resistivity of the single crystals was too low below 40 K to estimate $`\rho _{xy}`$ (and thus $`\sigma _{xy}`$). As seen in the bottom panel of Fig. 2, the anomalous conductivity remains fairly large at low temperatures, amounting to $``$ -100 Scm-1 at 2 K for the undoped and $`x`$ = 0.1 films. $`\sigma _{xy}`$($`T`$) shows a similar high-temperature peak, which, as $`T_c`$, is shifted to lower temperatures as the Ca doping increases. If the anomalous conductivity data in Fig. 2 is plotted against the magnetization instead of the temperature, as in the top of Fig. 3, one observes, within the measurement uncertainties, a similar or universal behavior of the Hall conductivity for all the samples.

Now we discuss the physical origin of these behaviors. The scaling $`\rho _{xy}\rho _{xx}^2`$ itself is often observed experimentally Hall  and is even derived considering extrinsic mechanismsside . In these conventional theories, a simple proportionality relation $`\rho _{xy}M`$ is derived in terms of the perturbative expansion in the spin-orbit coupling $`\lambda `$ and $`M`$. However, this simple relation is violated because the band crossing occurs in the band structure, and $`\lambda M`$ lifts this degeneracynagaosa ; Fang . This degeneracy point is known to act as a monopole for the gauge field representing the Berry phase curvature, producing its singular distribution. This non-perturbative feature causes the rapid change of $`\sigma _{xy}`$ including sign reversal as a function of $`M`$ because the Fermi energy crosses this monopole energy as $`M`$ changes. First-principles calculations confirm this scenario as shown in the inset of Fig. 3 for SRO (undoped case), using the orthorhombic structure obtained for the bulk SRO crystal. As seen in this inset, the calculations reproduce closely the non-monotonous variation, as well as the sign change of $`\sigma _{xy}`$ with $`M`$ (or $`T`$). The first-principles calculations (bottom panel of Fig. 3) show how the anomalous Hall conductivity depends on the crystal structure and lifetime of the electrons. Results obtained by considering a cubic structure are indeed quite different from those obtained in the orthorhombic case, even though they qualitatively show a similar non-monotonous behavior. It is also shown in the figure how $`\sigma _{xy}`$ is reduced upon increasing the scattering rate in the calculations; the broadening parameters were chosen so as to reflect the increase of longitudinal resistivity shown in Fig. 1. The Ca doping of the Sr1-xCaxRuO3 films also induces slight structural changes, such as a smaller orthorhombicity (or more correctly tetragonality, c.f. inset of Fig. 1). It is thus expected that the temperature dependence of $`\sigma _{xy}`$ of the doped films should differ, more or less, from that of the undoped SRO film, reflecting the changes in the local lattice-structure of the system and their effect on the band structure. Nevertheless, if the additional electron scattering effects arising from the Ca doping are taken into account as the broadening-induced reduction of the magnitude of $`\sigma _{xy}`$, the anomalous Hall conductivity of the Sr1-xCaxRuO3 films shows a good scaling to $`M`$, while changing $`T`$ and $`x`$. This indicates that the AHE is mainly of intrinsic origin, as described in terms of the Berry phase connection.Fang

In summary, the anomalous Hall effect was investigated for thin films of Sr1-xCaxRuO3, in which the ferromagnetic interaction is weakened with increasing Ca content $`x`$. The Hall resistivity of the films vary in a similar fashion with the temperature $`T`$. The obtained anomalous Hall conductivity varies non-monotonously and non-trivially with $`T`$, and even changes sign. The results however, show a good scaling solely to the $`T`$\- and $`x`$-dependent magnetization $`M`$, which can be reproduced by first-principles calculations. The anomalous Hall effect appears, as asserted by Fang et alFang , as a hallmark of the presence of magnetic monopoles in the momentum space of the crystal.

Z. F. acknowledges the support from the NSF of China (N. 10334090 and 90303022).


text:
# The Boson peak and the phonons in glasses

## 1 Introduction

$`X`$-ray and neutron scattering techniques allow to obtain very detailed physical insight into the high-frequency (0.1–10 THz) vibrational dynamics of supercooled liquids and glasses. Within this range of frequencies their spectra reveal several universal properties exp , related with the presence of sound-like excitations even for momenta $`p`$ of the same order of magnitude of $`p_0`$, the first maximum of the static structure factor (typically corresponding to wave numbers of a few $`nm^1`$). This *high-frequency sound* is revealed as Brillouin-like peaks in the Thz region of the dynamic structure factor. An accessible quantity to experiments is the vibrational density of states (VDOS), $`g(\omega )`$, whose most striking feature is the presence of an excess of states over the Debye $`\omega ^2`$ law in the “low” frequency region, (i.e.where the dispersion relation is linear, but still in the Thz region)debate . This excess of states is seen as a peak when plotting $`g(\omega )/\omega ^2`$ and has been named Boson peak (BP) 111There exist alternative ways of defining the boson peak from experiments, for example as a peak in Raman scattering data or as a peak in the difference between the observed VDOS of the glass and that of the corresponding crystal. The peak position $`\omega _{BP}`$ usually shifts to lower frequency on heating bpvsheat , except for the case of silica bpsilica1 . In this material the shift is seen on lowering the density bpsilica2 .

Due to its universality, the relevant physics underlying the Boson peak can be hopefully captured by some simple model. Furthermore, several recent numerical simulations have shown that a model of harmonic vibrations is wholly adequate to describe this frequency range harmonic  and that anharmonicity need not be invoked. Given the presence of well formed local structures (SiO2 tetrahedra, for instance) a natural approximation is to consider that the oscillation centers form a crystalline structure, the disorder in the atomic positions being mimicked by randomness in their interaction potential schirmacher98 ; taraskin  (disordered lattice models elliott ). The main drawback of such a models is that they dramatically underestimate the scattering of sound waves mayor00 . A different approach studies vibrations around a topologically disordered elliott  (liquid like) structure. It is followed by two different theories: modified mode-coupling theory mct  (which is not limited to harmonic excitations), and euclidean random matrix theory (ERMT) ermZee ; erm2 ; erm2bis . ERMT owes its name to the fact that it formulates the vibrational problem as random matrix problem metha . The matrices involved are called Euclidean random matrices erm , and their study has required the development of new analytical tools. Both MCT and ERMT predict an enhanced scattering of sound waves as compared to disordered crystals.

On the other hand, even within the harmonic framework the nature of the extra low frequency modes giving rise to the BP is still an open point. At a qualitative level, the frequency $`\omega _{BP}`$ is close to the Ioffe-Regel iofferegel  frequency $`\omega _{IR}`$, suggesting the possibility that the excess BP modes are localized alexander89 . However, numerical simulations have shown that the localization edge is at frequencies greater than $`\omega _{BP}`$ and $`\omega _{IR}`$ diffusoni . The Ioffe-Regel criterion signals rather a crossover to a region where the harmonic excitations are not longer propagating, due to the strong interaction with the disorder. We call these modes glassons (since they do not propagate but “diffuse”, they have also been called *diffusons* diffusoni ). A large bump of glassons is generally found around the Ioffe-Regel frequency, due to the flattening of the dispersion relation. This can be considered as the glass counterpart of the van Hove singularity of crystals taraskin ; erm3 . All the recently proposed theoretical frameworks predict that this peak of glassons should move to lower frequencies when approaching an instability transition, where negative eigenvalues (imaginary frequencies) appear. The aim of the paper is showing that the ERM theory makes sharp predictions about the values of universal critical exponents describing the approach to this singularity and comparing them with numeric results. The emerging scenario describes the BP modes as given by the hybridization between the phonons and the low-energy tail of the glasson peak which softens when the system approaches the instability transition nature ; vetto .

## 2 The Euclidean Random Matrix Theory

The starting approximation is that particles can only oscillate around fixed random positions, so that the position of particle $`i`$ at time $`t`$ is $`𝐱_i(t)=𝐱_i^{eq}+\phi _i(t)`$; the $`𝐱_i^{eq}`$ are quenched equilibrium positions (whose distribution must be specified) and $`\phi _i(t)`$ are the displacements. Hence the Hamiltonian is

$$H\left[𝐱\right]=\underset{i,j}{\overset{1,N}{}}V(𝐱_i𝐱_j)\frac{1}{2}\underset{i,j}{\overset{1,N}{}}\underset{\mu ,\nu }{\overset{1,3}{}}M_{i\mu ,j\nu }[𝐱^{eq}]\phi _i^\mu \phi _j^\nu $$
(1)
where the dynamical matrix M is an Euclidean Random Matrix:

$$M_{i\mu ,j\nu }[𝐱^{eq}]f_{\mu \nu }(𝐱_i^{eq}𝐱_j^{eq})+\delta _{ij}\underset{k=1}{\overset{N}{}}f_{\mu \nu }(𝐱_i^{eq}𝐱_k^{eq}),$$
(2)
with $`f_{\mu \nu }(𝐱)_{\mu \nu }V(𝐱)`$.

In the one-excitation approximation the dynamic structure factor is

$$S^{(1)}(𝐩,\omega )=\frac{k_\mathrm{B}T}{m\omega ^2}\overline{\underset{n}{}\left|\underset{i}{}𝐩𝐞_{n,i}e^{i𝐩𝐱_i^{eq}}\right|^2\delta (\omega \omega _n)},$$
(3)
where $`𝐞_n`$ are the eigenvectors of the dynamical matrix and $`\omega _n`$ its eigenfrequencies (= square root of eigenvalues). The overline means average over the disordered quenched positions, whose distribution $`P[𝐱^{eq}]`$ has to be specified. The density of states (VDOS) is obtained in the limit of large momenta:

$$g(\omega )=\underset{p\mathrm{}}{lim}\frac{m\omega ^2}{k_BTp^2}S^{(1)}(p,\omega ).$$
(4)
We can obtain $`S^{(1)}(𝐩,\omega )`$ through the resolvent $`G(𝐩,z)`$:

$`G_{\mu \nu }(𝐩,z){\displaystyle \frac{1}{N}}{\displaystyle \underset{jk}{}}\overline{e^{i𝐩(𝐱_j^{\mathrm{e}q}𝐱_k^{eq})}\left[{\displaystyle \frac{1}{zM}}\right]}_{j\mu ,k\nu }`$
$`G_L(p,z){\displaystyle \frac{p_\mu p_\nu }{p^2}}+G_T(p,z)\left(\delta _{\mu \nu }{\displaystyle \frac{p_\mu p_\nu }{p^2}}\right)`$ (5)
separating the axial tensor in a longitudinal term and a transversal one. The dynamic structure factor is then:

$$S^{(1)}(𝐩,\omega )=\frac{2k_\mathrm{B}Tp^2}{\omega \pi }\mathrm{Im}G_L(𝐩,\omega ^2+\mathrm{i0}^+).$$
(6)
A transverse dynamic structure factor (not measurable in experiments) can be defined in an analogous way. However, a most important and general result is that for $`p\mathrm{}`$ the resolvent becomes isotropic:

$$G_{\mu \nu }^{\mathrm{}}(z)=\frac{1}{N}\underset{j}{}\overline{\left[\frac{1}{zM}\right]}_{j\mu ,j\nu }=\delta _{\mu \nu }\frac{1}{N}\overline{\text{Tr}[zM]^1}.$$
(7)
So both longitudinal and transverse structure factors tend to a common limit (the VDOS, see eq. 4) at infinite momentum. 222Conseguently, both the dispersion relations saturate at the same value. However due to the broadening of the line, they are rather ill-defined when $`\omega \omega _{IR}`$ Leaving the potential $`V(r)`$ unspecified and taking the simplified case $`P[𝐱^{eq}]=1/V^N`$ ($`V`$ being the volume), one finds that:

$$G_{\mu \nu }(𝐩,z)=\left[\frac{1}{z\rho \widehat{f}(0)+\rho \widehat{f}(𝐩)\mathrm{\Sigma }(𝐩,z)}\right]_{\mu \nu },$$
(8)
The self-energy $`\mathrm{\Sigma }(𝐩,z)`$ is a matrix with the standard form

$$\mathrm{\Sigma }_{\mu \nu }(𝐩,z)=\mathrm{\Sigma }_L(p)\frac{p_\mu p_\nu }{p^2}+\mathrm{\Sigma }_T(p)\left(\delta _{\mu \nu }\frac{p_\mu p_\nu }{p^2}\right).$$
(9)
which vanishes at $`\rho =\mathrm{}`$ and that can be computed in a series expansion in $`1/\rho `$. The main point is that the sum of all the infinite diagrams obtained recursively starting from this next-to-leading order result gives a self-consistent integral equationvetto :

$$\mathrm{\Sigma }_{\mu \nu }(𝐩,z)=\frac{1}{\rho }\frac{d^3q}{(2\pi )^3}V_{\mu \lambda }(𝐪,𝐩)G_{\lambda \sigma }(𝐪,z)V_{\sigma \nu }(𝐪,𝐩).$$
(10)
where the vertices have the form $`V_{\mu \nu }(𝐪,𝐩)=\rho (\widehat{f}_{\mu \nu }(𝐪)\widehat{f}_{\mu \nu }(𝐩𝐪))`$. Let us remark that the self energy renormalizes the dispersion relations and gives a finite width to the Brillouin peaks:

$`\omega _{L,T}^2(p)`$ $`=`$ $`\left(\omega _{L,T}^0\right)^2(p)+\text{Re}\mathrm{\Sigma }_{L,T}(p,\omega _{L,T}(p)),`$
$`\mathrm{\Gamma }_{L,T}(p)`$ $`=`$ $`\text{Im}\mathrm{\Sigma }_{L,T}(p,\omega _{L,T}(p))/\omega _{L,T}(p).`$ (11)
The correlations between the equilibrium positions of the particles can be taken into account quite easily at the level of the superposition approximation in the above approach. The results derived above for the case without correlations are translated to the correlated case by replacing the functions $`f(𝐱)`$ by $`g^{(2)}(𝐱)f(𝐱)`$. In this way the usual power law divergence of the pair potential for $`|𝐱|0`$ is balanced by the exponential behaviour of the pair distribution function, and this ensures the existence of the Fourier transform of the product $`f(𝐱)g^{(2)}(𝐱)`$.

### 2.1 The phase transition

From equation (10) it is possible to derive a few analytic model-independent results about the arising of the Boson Peak. These results are expressed in form of scaling laws, whose exponents are predicted in this approximation. Simulations (see below) and experiments will allow to clarify the dependence of the exponents on the approximation. The VDOS can be obtained from

$$g(\omega )=\frac{2\omega }{\pi }\text{Im}G^{\mathrm{}}(\omega ^2+i0^+),$$
(12)
where $`z=\omega ^2+i0^+`$ and $`G^{\mathrm{}}(z)lim_p\mathrm{}G(𝐩,z)`$.

Hence one have to solve the integral equation (10) in the $`p\mathrm{}`$ limit:

$$\frac{1}{G^{\mathrm{}}(z)}=z\rho \widehat{f}(0)\rho AG^{\mathrm{}}(z)\rho \frac{d^3q}{(2\pi )^3}\widehat{f}^2(𝐪)G(𝐪,z)$$
(13)
where $`G^{\mathrm{}}(z)`$, $`A(2\pi )^3d^3q\widehat{f}^2(𝐪)`$ and the last term are matrices proportional to the identity.

The solution of the above integral equation yields a VDOS which contains both the phonons, since $`g(\omega )\omega ^2`$ at $`\omega 0`$, and the extended but not propagating glassons, described by a semicircle with center at $`\omega =\rho \widehat{f}(0)`$ and radius $`2\sqrt{\rho A}`$ erm3 . If we limit (for pedagogical purposes) to the case where the VDOS changes because of changes in the density, the key point is the existence of a phase transition in the space of the eigenvalues of the Hessian matrix. In fact $`G^{\mathrm{}}(0)`$ develops an imaginary part when $`\rho <\rho _c`$ ($`\rho _c`$ being a critical density), then the transition separates the stable phase (all positive eigenvalues) and the unstable phase (negative and positive eigenvalues). The glassons are the modes which move towards the negative zone of the spectrum (hybridizing the phonons) when approaching the transition. The order parameter is $`\phi =\text{Im}G^{\mathrm{}}(i0^+)`$ which vanishes as $`\phi |\mathrm{\Delta }|^\beta `$, with $`\beta =1/2`$ and $`\mathrm{\Delta }\rho \rho _c`$

The relation with the Boson Peak becomes quite clear when one writes down the VDOS in the stable phase arising from the theory without any reference to the control parameter, which then does not need to be $`\rho `$. In fact, one has

$$g(\omega ,\mathrm{\Delta })=\omega ^\gamma h(\omega \mathrm{\Delta }^\rho ),h(x)\{\begin{array}{ccc}x^{2\gamma }\hfill & & \hfill x1\\ \text{const.}\hfill & & \hfill x1\end{array},$$
(14)
with $`\mathrm{\Delta }`$ defined in terms of an arbitrary control parameter.The ERMT (in the cactus approximation) predicts $`\rho =1`$, $`\gamma =3/2`$. Hence it exist a crossover frequency (in the region where the dispertion relation is still linear) between a $`\omega ^2`$ and a $`\omega ^\gamma `$ region. We shall identify that with the BP frequency $`\omega _{BP}`$. This implies that $`\omega _{BP}\mathrm{\Delta }^\rho `$ and $`g(\omega _{BP},\mathrm{\Delta })/\omega _{BP}^2\mathrm{\Delta }^\eta `$, with

$$\eta =\rho (2\gamma ).$$
(15)
Let us note that the BP is indicated from a peak in the function $`g(\omega )/\omega ^2`$m not in $`g(\omega )`$. Summarizing, according to ERMT the BP frequency moves linearly toward $`0`$ when approaching the transition (from the stable side) and its height diverges as a power law whose exponent is $`\eta =1/2`$. Eq. (4) shows that at the level of the one-phonon approximation it can also be detected in the large $`p`$ limit of the dynamic structure factor $`S(p,\omega )`$.

## 3 Boson Peak in a gaussian model

In order to confirm that the saddle-phonon transition described by the Euclidean Random Matrix theory is not an artifact of the approximation involved (cactus resummation), we solved numerically the cactus equation for the case where $`f(p)`$ has a Gaussian form and compare with direct numerical results for the same modelvetto . The model is described by

$`\widehat{f}_{\mu \nu }(𝐩)`$ $`=`$ $`\widehat{f}_L(p){\displaystyle \frac{p_\mu p_\nu }{p^2}}+\widehat{f}_T(p)\left(\delta _{\mu \nu }{\displaystyle \frac{p_\mu p_\nu }{p^2}}\right),`$
$`\widehat{f}_{L,T}(p)`$ $`=`$ $`\left({\displaystyle \frac{2\pi }{\sigma _0^2}}\right)^{3/2}\mathrm{exp}(p^2/2\sigma _{L,T}^2).`$ (16)
This choice for $`\widehat{f}(p)`$ is mainly due to its simplicity. However the behaviour of the Boson peak close to the saddle-phonon transition have to be independent of the details of the model. Moreover the superposition approximation takes $`\widehat{f}_{\mu \nu }(p)=[g(r)v_{\mu \nu }(r)]`$ yielding a $`\widehat{f}_{L,T}(p)`$ finite at $`p=0`$, like in the Gaussian model. We shall consider various values of the density, which is here the control parameter, comparing the analytical (cactus) results with the numerical spectra and dynamic structure factor obtained from the method of moments mom . In the high density regime the approximations used in deriving the integral equation (13) are quite good since the analytic solution reproduces the numerical spectrum (and in particular the Debye behaviour) rather accurately (fig.1).

However, the crucial check regards the exponents of the transition. Figs. 2a and 2c show that the position of the BP is linear with respect to $`\mathrm{\Delta }(\rho \rho _c)`$ and that the height of BP diverges as $`\mathrm{\Delta }^{1/2}`$. This confirms the theoretical predictions $`\nu =1`$ and $`\eta =1/2`$. In Fig. 2b we determine the value of $`\gamma `$ by studying the fraction of unstable modes. In fact, in the region of parameters where $`\rho <\rho _c`$ the fraction of unstable modes, defined as $`f_u=_{\mathrm{}}^0g_\lambda (\lambda )𝑑\lambda `$, is given by

$$f_u(\mathrm{\Delta })=_0^{\mathrm{}}𝑑\omega \omega ^\gamma \stackrel{~}{g}(\omega /|\mathrm{\Delta }|)|\mathrm{\Delta }|^{1+\gamma }.$$
(17)
We find numerically (Figs. 3b) that $`f_u(\rho _c\rho )^{5/2}`$, i.e. $`\gamma =3/2`$. Finally, the order parameter $`\phi `$ vanishes as $`(\rho _c\rho )^\beta `$ with $`\beta =1/2`$ (Fig. 2d).

Hence our analytic treatment based on euclidean random matrix theory describe quite well the vibrational features of simple topologically disordered systemsvetto . The following step is understanding to what degree of accurateness ERMT could describe the high frequency properties of more realistic systemsnature .

## 4 Boson Peak in a fragile glass

Starting from the hypothesis that the Thz region of supercooled liquids and glasses can be described in terms of purely harmonic excitations, the origin of the Boson peak in glasses can be understood if we consider the ensemble of generalized inherent structures (GIS). For each equilibrium configuration the associated GIS is the nearest stationary point of the Hamiltonian. If we start from an equilibrium configuration at low temperature, the GIS is a local minimum, and it coincides with the more frequently used inherent structures (IS) INM (i.e. the nearest minimum of the Hamiltonian). On the contrary, if we start from high temperature, the GISs are saddle points. In the GIS ensemble there is a sharp phase transition separating these two regimes. It takes place in glass-forming liquids saddles  at the Mode Coupling temperature tmc  ($`T_{MC}`$), above which liquid diffusion is no longer ruled by rare “activated” jumps between ISs but by the motion along the unstable directions of saddles. Phonons are present in the spectrum of the VDOS in the low temperature phase (IS dominated) but are absent in the saddle phase. The key point is that the minima obtained starting from configurations below $`T_{MC}`$ and the saddles obtained starting above $`T_{MC}`$ join smoothly at $`T_{MC}`$. Thus we can study GIS as a single ensemble parametrized by their energy saddles .

Since this transition separates a phase where all the eigenvalues are positive from another one where even negative eigenvalues exist, we expect that ERMT is able to describe correctly this phenomenum. Hence we measured numerically the values of some exponents predicted by the theory in a simple model of a fragile glass nature . We simulated a soft-spheres binary mixture Bernu87  in the stable (phonon) phase with the Swap Monte Carlo algorithm GrigeraParisi01 , and computed the VDOS of the ISs obtained starting from equilibrium configurations at temperatures below $`T_{MC}`$333At very low $`T`$, where equilibrium is not achieved, runs were followed until $`e_{IS}`$ got very close to its asymptotic value.

In Fig.3 we show that the theoretical predictions agree with the numerical data. Taking the IS’s energy as the relevant parameter for describing the spectral properties, one has $`\mathrm{\Delta }=e_ce_{IS}`$, $`e_{IS}`$ being the energy of the ISs and $`e_c`$ the critical value. In fact, plotting $`g(\omega )/\omega ^2`$ a peak is clearly identified, which is seen to grow in height and shift to lower frequency on rising the IS’s energy. Using all the spectra for which the peak position can be clearly identified, we find that the relationship between $`\omega _{BP}`$ and the energy of the IS is linear (Fig. 3a). The energy at which $`\omega _{BP}`$ becomes zero, $`e_c`$, is found from a linear fit as $`e_c=1.74ϵ`$ ($`ϵ`$ is the energy scale), quite close to the value where the GIS stop to be minima (IS) and become saddles saddles . As for the height of the peak (Fig. 3b), the results are compatible with a power-law divergence. Fixing $`e_c`$ at the value $`1.74ϵ`$ arising from the linear fit of $`\omega _{BP}`$ vs. $`e_{IS}`$, a power-law fit yields an exponent $`\beta =0.40(15)`$, while if one fixes the exponent at $`\beta =1/2`$, then the critical value is $`e_c=1.752(2)ϵ`$. Thus the numerical data are compatible with the theoretically predicted scaling, although we have not been able to work close to the critical point, and thus cannot get a great accuracy on the critical exponents or the critical point.

## 5 Conclusions

In summary, we have shown that the saddle (negative eigenvalues)-phonon (no negative eigenvalues) transition, a common feauture of vibrating topologically disordered systems, is well described by the euclidean random matrix theory. It provides a coherent scenario for the arising of a Boson Peak, which results from the hybridization of acoustic modes with high-energy modes that soften upon approaching the transition. Hence we applied the theory to describe the saddle-phonon transition and the BP in supecooled liquids, comparing the predicted scaling laws with the numeric results obtained for a simple fragile glass former. The agreement found is quite encouraging. The present discussion applies to experiments as long as one is in the regime where the inverse frequency is much larger than the structural relaxation time, when the harmonic approximation makes sense. We expect that the saddle-phonon transition point of view will be able to bridge the realms of experiment and numerical studies of the energy landscape. As a matter of fact, a recent experiment on the poly(methyl methacrylate) (PMMA) glass gave the first experimental confirmation of the ERMT predictionsDuval03 .

V.M.-M. is a Ramón y Cajal research fellow (MCyT, Spain). P.V. was supported through the European Community’s Human Potential Programme under contract HPRN-CT-2002-00307, DYGLAGEMEM.


text:
# Jet-Induced Star Formation

## 1. Introduction

Observations suggest that large scale, shock-induced star formation is an integral part of the galaxy evolution process. It is important in colliding and merging galaxies such as the Ultra-Luminous Infrared Galaxies (ULIRGs) and in forming galaxies. It may also occur to a lesser extent near the central dominant galaxies in cooling flow clusters and in galaxies moving through dense cluster atmospheres. Perhaps the most spectacular shock-induced star formation occurs as the result of the collision of extragalactic radio jets with over-dense ambient material. Here we review existing observational evidence in support of jet-induced star formation in both nearby and distant radio galaxies. We will then discuss numerical simulations performed at the University of California Lawrence Livermore National Laboratory which are aimed to explore the conditions under which such star formation may occur, and apply this to the nearby ’Minkowski’s Object’ (M.O.).

## 2. Environments of radio galaxies

Radio galaxies do not live in a ’vacuum’ but are surrounded by gaseous halos and/or debris from recent merger events with may have triggered the radio galaxy activity in the first place. Broadly speaking one can discriminate between the low-luminosity Fanaroff and Riley FR-I types (such as 3C 31), which have expanding, turbulent jets without terminal shocks or hotspots, and the high luminosity FR-II types (such as Cygnus A) which have very well collimated jets ending at high surface brightness hotspots and large diffuse lobes (Fanaroff and Riley 1974). FR-I’s may be found in clusters of galaxies where they can interact with cooling clouds in the Inter Galactic Medium (Ferland et al. 2002). At high redshift ($`z>2`$), because of sensitivity limits, only the FR-II types are easy to detect. Since they inhabit young, forming galaxies their jets propagate through relatively dense and clumpy media. Observational evidence for this dense and cool gas in high redshift radio galaxies has been plentiful and includes the detection of dust, HI, extended line emission, associated absorption line systems and molecular gas (e.g. De Breuck et al. 2003).

### 2.1. Jet-induced star formation in nearby FR-I type radio galaxies

One of the first radio galaxies where evidence was found for jet-induced star formation was the nearest radio galaxy, Centaurus A (Blanco et al 1975). Further examples have been found as the sensitivity and spatial resolution of radio and optical telescopes has improved. In the case of Centaurus A, recent observations with the Hubble Space Telescope, when compared to radio images obtained with the Very Large Array, have confirmed that there are about half a dozen young ($`<15`$ Myr) OB associations near filaments of ionized gas located between the radio jet and a large HI cloud (Mould et al. 2000). Another nearby example is ’Minkowski’s Object’, a peculiar small starburst system at the end of a radio jet emanating from the elliptical galaxy NGC 541, located near the center of the cluster of galaxies Abell 194 (van Breugel et al. 1985). Star forming regions associated with radio sources have also been found in cooling flow clusters, with the best example being in Abell 1795 (McNamara 2002).

### 2.2. Jet-induced star formation in distant FR-II type radio galaxies

Radio galaxies have now been identified up to $`z5.2`$. They are the most massive galaxies at any redshift and, at high redshift, are usually associated with extended, clumpy systems. Many have giant Ly-a halos, dust and CO molecular gas. One of the most striking correlations in $`z>0.6`$ radio galaxies is that their optical line and continuum emission is aligned with the radio sources (McCarthy et al. 1987; Chambers et al. 1987). The best studied example here is 4C41.17 at $`z=3.8`$, where deep spectroscopic observations have shown that the bright, spatially extended, rest-frame UV continuum emission aligned with the radio axis of this galaxy is unpolarized and shows P Cygni-like features similar to those seen in star-forming galaxies (Dey et al. 1977).

Collectively, these observations are best explained by models in which shocks generated by the radio jet propagate through an inhomogeneous medium and trigger gravitational collapse in relatively overdense regions (Begelman et al. 1989; De Young 1989; Rees 1989). A detailed analysis of the jet-induced star formation in 4C41.17 has been presented by Bicknell et al. (2000). In that object, Hubble Space Telescope images showed a bimodal optical continuum structure parallel to the radio jet, strongly supporting the idea that the star formation was triggered by sideways shocks.

## 3. Minkowski’s Object

One of the most spectacular observed jet-induced starbursts is ’Minkowski’s Object’ associated with the elliptical galaxy NGC 541 in the cluster of galaxies Abell 194 (van Breugel et al. 1985). Its morphology is strongly suggestive of a collision between the FR-I type jet from NGC 541 and a dense cloud: M.O. has the same overall diameter as the jet, appears wrapped around the end of the jet, has bright emission in the upstream direction and filamentary structure down-stream where the jet appears disrupted. Spectroscopically M.O. looks like an HII region, resembling starburst galaxies. The H-a luminosity suggests a modest star formation rate of 0.3 M/yr. VLA observations (van Breugel and van Gorkom, unpublished) show two detections down-stream from the jet-cloud collision site, indicating a total HI mass of $`3\times 10^8`$ M.

## 4. Numerical simulations

The basic idea of shock-induced star-formation is that when a strong shock passes through a clumpy medium it may trigger many smaller-scale compressive shocks in overdense clumps. These shocks increase the density inside the clumps and make it possible for them to radiate more efficiently. If the radiative efficiency of the gas has a sufficiently shallow dependence upon the temperature, then radiative emissions are able to cool the gas rapidly, in a runaway process, producing even higher densities as the cooling gas attempts to re-attain pressure balance with the surrounding medium. Both the reduction in temperature and the increase in density act to reduce the Jeans mass, above which gravitational forces become important. Any clump that was initially close to this instability limit will be pushed over the edge and forced into gravitational collapse. Thus, the passing of a shock through a clumpy medium may trigger a burst of star formation. Several processes, acting on different timescales, govern whether or not cooling and star formation are able to proceed.

To investigate under which conditions shock-induced star formation can occur requires numerical simulations. This also provides information about environments that may be observationally out of reach, such as in forming galaxies at high redshift, and which could be important for understanding the role of jets in the feedback of active galactic nuclei (AGN) on their environment. Numerical simulations can also give us insight into the role of shock-induced star formation in other environments such as supernova shocks and cloud-cloud collisions.

The interaction of a strong shock with non-radiative clouds has been the subject of many numerical studies (e.g. Klein et al 1994; Poludnenko et al 2002). For a non-radiative cloud the passing shock ultimately destroys the clouds within a few dynamical timescales. Destruction results primarily from hydrodynamic instabilities at the interface between the cloud and the post-shock background gas.

The effects of strong shocks interacting with radiative clouds are very different and studies of this have only recently begun (Mellema et al. 2002; Fragile et al. 2003). Instead of re-expanding and quickly diffusing into the background gas, the compressed cloud instead breaks up into numerous dense, cold fragments. These fragments survive for many dynamical timescales and are presumably the precursors to star formation. Previous work has focused primarily on FR-II type radio galaxies even though, because of their proximity, much better observational data can be obtained for FR-I type jet-induced starbursts such as M.O. Fragile et al (2003) used the LLNL developed, multi-dimensional, multi-physics, massively parallel ’COSMOS’ numerical simulations package (Anninos et al. 2003; Anninos and Fragile 2003) to investigate both FR-I and FR-II type jet-induced star formation systems.

### 4.1. Application to Minkowski’s Object

To simulate the jet induced star formation in M.O. Fragile et al (2003) assumed that NGC 541 is surrounded by a multi-phase medium resembling cluster atmospheres (e.g. Ferland et al. 2002). Specifically, it was assumed that the FR-I jet interacts with a moderate density, hot ‘mother’-cloud with a semimajor axis of 10 kpc, a semiminor axis of 5 kpc, a density $`n_{mcl}=0.1`$ cm-3, and temperature $`T=10^6`$ K. This corresponds to an initial total cloud mass of $`10^9M_{}`$. Within this mother-cloud denser, warm clouds were assumed to be embedded with typical sizes of 100 pc, $`n_{cl}=10`$ cm-3, and temperature $`T=10^4`$ K.

The detailed radio and X-ray study of the proto-typical FRI-type radio galaxy 3C31 by Laing and Bridle (2002) was used to estimate plausible jet parameters near M.O., at $`15`$ kpc from the NGC 541 AGN. We assumed that the long axis of the mother-cloud is aligned with a jet of high velocity ($`0.3c`$), low density ($`10^4`$ cm-3) gas flowing onto the grid. The diameter of the jet nozzle is equal to half the diameter of the cloud along the semiminor axis.

Our numerical simulations consisted of two parts. To investigate the overall structure of M.O. we followed the evolution of the collision of the jet with the low density mother-cloud. To determine whether star formation can be triggered we followed the evolution of shocks interacting with one or more of the dense, warm clouds embedded in the mother-cloud.

The collision with the mother-cloud triggers a nearly planar shock down the long axis of the cloud. As the bow shock from the jet wraps around, it also triggers shocks along the sides of the cloud. A similar shock structure may explain the filamentary nature of the star-forming region in M.O. (Fig. 4). To determine whether the clouds embedded within the mother-cloud would indeed collapse and form stars we ran a number of simulations to explore the density / velocity parameter space. These results are summarized in Figure 5.

## 5. Conclusions

From our numerical simulations, when applied to M.O., a number of important conclusions can be drawn. First, its peculiar morphology \- bright star forming region orthogonal to the jet, and fainter filamentary features downstream from there - can be easily reproduced by our simulations. Second, the modest amount star formation required - $`0.3M_{}`$ yr-1 for the entire object, is also easily achieved for the plausible parameter space explored by our simulations. Third, and most interestingly, we conclude that the star formation in M.O. could be induced by a moderate velocity jet ($`9\times 10^4`$ km s-1) interacting with a collection of slightly overdense ($`10`$ cm-3), warm ($`10^4`$ K) clouds, i.e. it is NOT necessary to assume that this was an accidental collision between a jet and a preexisting gas-rich galaxy. This also suggests that the neutral hydrogen associated with M.O. ($`3\times 10^8M_{}`$; W. van Breugel & J. van Gorkom 2003, private communication) may have cooled from the warm gas phase as a result of the radiative cooling triggered by the radio jet. An update on the observations of M.O. will be presented in a forthcoming paper (S. D. Croft et al., in preparation).

Our models can be used to infer the importance of shock-induced star formation in other regimes by simply adjusting the criteria used to develop Figure 5. They have confirmed the reality of the division between regions II and V (cooling and non-cooling) derived using approximate analytical methods, while the other region boundaries are set by physical limits for the clouds. Our key conclusion is that shocks associated with jets may indeed trigger the collapse of clouds to form stars. Whether this occurs at the impact area, or along the sides of expanding lobes depends on jet power and the ambient gas density distribution. In forming galaxies, where both dense gas and AGN activity are likely, jet-induced star formation might help trigger even more star formation, at an earlier stage, then otherwise might have occured.

In future work, we shall expand our numerical simulations package to include Adaptive Mesh Refinement and magnetic fields and will investigate regimes appropriate for other possible shock-cloud interaction scenarios.

#### Acknowledgments.

This work was performed under the auspices of the U.S. Department of Energy by University of California, Lawrence Livermore National Laboratory under Contract W-7405-Eng-48. W.v.B. also acknowledges NASA grants GO 9779 and GO3-4150X in support of high-redshift radio galaxy research with HST and Chandra.

## References

Anninos, P., & Fragile, P. C. 2003, ApJS, 144, 243

Anninos, P., et al. 2003, ApJS, 147, 177

Bicknell, G. V. et al. 2000, ApJ, 540, 678

Begelman, M. C. & Cioffi, D. F. 1989, ApJ, 345, L21

Blanco, V. M. et al. 1975, ApJ, 198, L63

Chambers, K. C. et al. 1987, Nature, 329, 604

De Breuck, C. et al. 2003, å, 401, 911

De Young, D. S. 1989, ApJ, 342, L59

Dey, A. et al. 1997, ApJ, 490, 698

Fanaroff, B. L. & Riley, J. M. 1974, MNRAS, 167, 31P

Ferland, G. J. et al. 2002, MNRAS, 333, 876

Fragile, P. C. et al. 2003, ApJ, in press

McKee, C. F. et al 1994, ApJ, 420, 213

Laing, R. A. & Bridle, A. H. 2002, MNRAS, 336, 1161

Mould, J. R. et al. 2000, ApJ, 536, 266

McCarthy, P. J. et al. 1987, ApJ, 321, L29

McNamara, B. R. 2002, New Astronomy Review, 46, 141

Mellema, G. et al. 2002, A&A, 395, L13

Poludnenko, A. Y. et al. 2002, ApJ, 576, 832

Rees, M. J. 1989, MNRAS, 239, 1P

van Breugel, W. et al. 1985, ApJ293, 83