text:
# Intersecting branes with an arbitrary excess angle

## I Introduction

In the seminal paper  gravity lives in five-dimensional anti de-Sitter ($`AdS_5`$) spacetime, infinite but strongly curved. Due to the presence of one graviton bound state, related to the finite proper size of “bulk” spacetime, gravity is effectively localized to a 4-dimensional “brane”. This mechanism however works only on codimension 1 branes but could be indirectly generalized  by placing $`D4`$ codimension 1 branes intersecting in one 3-brane (our world). The brane tension of each of these 4-branes localizes gravity in the respective dimension. Models have been built and their cosmological consequences analyzed in .

This geometric setup has also been used to address the cosmological constant. Static embeddings of 3-branes in $`AdS_5`$ spacetimes were found  where the induced metric is dS, Minkowski or AdS; the curvature of the brane could be made much smaller than the brane tension. The same solution, now for a 4-brane in 6 dimensions, was rediscovered in conformally flat form in . It was observed in this paper that it is possible to embed a static 4-brane, plane in conformal coordinates, at any angle $`\varphi `$ determined by the brane tension, $`\mathrm{sin}\varphi =\sigma /\sigma _c`$, where $`\sigma _c`$ is a critical value for the brane tension. A vanishing observed cosmological constant (i.e. flat induced metric) is produced when $`\sigma =\sigma _c`$. Cutting out wedges of the bulk $`AdS_5`$ with some angle $`\varphi `$ and pasting together $`n`$ copies can accommodate an excess angle $`n\varphi 2\pi `$. This must match the 3-brane vacuum energy in the intersection. The wedges are stabilized by imposing orbifold boundary conditions on each 4-brane. Whereas in the original Randall-Sundrum models  the 3-brane tension is fine tuned to the bulk cosmological constant, here it must be fine tuned to match the excess angle, a function of the 4-brane tensions.

A similar setup with different bulk cosmological constants in each wedge is considered in , although without the introduction of an excess angle and correspondingly, without a 3-brane tension. It has been speculated in this reference that an additional parameter of the solution could resolve the one remaining fine tuning. What we need is a flat direction; then a new dynamical mechanism could stabilize the cosmological constant at zero value. Any additional vacuum energy would result in the readjustment of the brane angles and the additional parameter.

In this paper we find such an additional parameter. Using Gaussian normal coordinates (GNC) of a 4-brane, we write down all vacuum solutions in the bulk (with a negative cosmological constant), satisfying the Ansatz with signature $`+++++`$***Allowing more than one timelike coordinate some restrictions on the brane tension disappear.

$$ds^2=dw^2+B(w,y)dy^2+A(w,y)\eta _{\mu \nu }dx^\mu dx^\nu .$$
(1)
The 4-brane is located at $`w=0`$. The stability of the brane with orbifold boundary conditions (i.e. $`ww`$ identified while $`y,x^\mu `$ are fixed) comes down to the Israel conditions  on the extrinsic curvature, $`K_{ab}=\frac{\sigma }{2}h_{ab}`$, where $`\sigma `$ is the 4-brane tension in appropriate units ($`h_{ab}`$ is the extrinsic curvature of the brane).

There is considerable freedom in transforming the coordinates, still within the Ansatz of Eq. (1). In addition, the same bulk metric is found starting with brane with different tensions, embedded at various angles. After eliminating all these ambiguities we are left with four different solutions (in each case the $`y`$ coordinate is chosen to render the induced metric conformally flat):

* A solution $`A`$ in which only a supercritical tension brane can be embedded, $`\sigma >\sigma _c`$, where the induced metric has two curvature singularities at finite $`y\pm y_{}`$ (these are also at finite proper distance)

and three solutions where only $`\sigma <\sigma _c`$ branes can be embedded, namely

* A solution $`B`$ where the induced metric has a curvature singularity at finite proper and coordinate) distance $`yy_{}`$ and an $`AdS`$ horizon at $`yy_+`$.
* A solution $`C`$ with $`AdS`$ horizon at both $`y\pm y_+`$. The $`w\mathrm{}`$ limit is also asymptotically AdS.
* A solution $`D`$ where both the bulk and the induced metric are AdS. This is the solution found previously .

(We do not consider the fine tuned case, $`\sigma =\sigma _c`$.)

All these solutions are vacuum solutions with a negative cosmological constant in the bulk, i.e. solutions with constant scalar curvature. With the exception of $`D`$, however, they are not maximally symmetric. The induced metric is not even a vacuum solution. It is worth noting that the corresponding embedding of 3-branes in 5 dimensional spacetimes can be similarly found, but the new solutions have little value because the induced 3-brane metric does not resemble our world at all.

All our solutions have four-dimensional Poincaré invariance and can be trivially generalized to any maximally symmetric 4-dimensional metric replacing $`\eta _{\mu \nu }`$ in Eq. (1). Eventually the function $`A(w,y)`$ is the same in that case, only $`B(w,y)`$ is changed. The solutions have a complicated mathematical form and are expressed in an implicit form as solutions of algebraic equations involving special functions.

The usual picture of localization of gravity to a codimension one brane acquires some peculiar features when the brane tension differs from the fine tuned value. The warp factor that effectively compactifies the bulk decreases up to a distance $`𝒪\left(\frac{\sigma }{\sigma _c\sigma }l_{AdS}\right)`$ in the extra dimension, but further away from the brane it starts to grow exponentially. The resulting “tail” of the graviton wave function eventually delocalizes brane gravity. This effect is absent for exact fine tuning and is also relevant for 3-branes embedded in a 5 dimensional bulk, e.g. for the models in . Models in which the brane tensions are far from their critical value do not seem to have brane gravity at all. We discuss in what extent this can be considered a coordinate artifact in Sec. II B.

At any point along the brane, another brane can be inserted at an angle determined by its tension. In coordinates where the first brane has a conformally flat induced metric, the second brane is curved. The (covariant) angle between the branes is $`\varphi =\pi \pm 2\mathrm{arcsin}(\sigma /V)`$. The quantity $`V`$ is independent of the tension of the second brane but changes as the intersection point is moved along the first brane ($`\pm `$ corresponds to two solutions and tells if $`A`$ increases or decreases as one moves away from the intersection along the brane). For example, close to the curvature singularity $`V\mathrm{}`$. Consequently, when several such wedges are sewn together in the same point, all branes will have the same parameter $`V`$. Now the induced metric on each brane is determined only by this $`V`$ (and the choice of the $`\pm `$ sign), so that the first and last branes can always be identified in order to build a manifold around the intersection. The orbifolding condition only requires, for global consistency, that every second brane around the intersection has the same tension (and $`\pm `$ sign).

The new parameter $`V`$ tells what part of the full 6-dimensional spacetime is cut away. Its emergence is due to the fact that the solution is not maximally symmetric. In the case of and $`AdS_6`$ solution, when the intersection point is shifted, the change in the metric can be compensated by a coordinate transformation, i.e. all points are physically equivalent. In our spacetime this is not so and $`V`$ parameterizes this difference. Because the angle of the wedges depends on $`V`$, so does the total excess angle. The topological requirement that a vacuum density produces an excess angle can be satisfied without fine tuning: the gravitational equations simply set the brane angles and the parameter $`V`$ to the consistent value.

When a brane tension is put on the 3-brane, another condition arises  which has not been discussed in . A positive 3-brane tension tries to minimize the 4-volume, which corresponds to minimizing the value of $`A(w,y)`$ at the intersection. There is such a minimum only in case $`C`$. Now at the minimum the value of $`V`$ is zero, so we loose the additional freedom we had when we only needed to stabilize the 4-branes. For purposes of illustration we consider the following unusual but consistent geometry. We take two branes (to be identified), of finite length. At each end the two are joined and the covariant intersection angle is taken to be more than $`2\pi `$. Upon identifying the two branes we find an asymptotically $`AdS_6`$ spacetime which satisfies the gravitational equations except the minimizing of $`A(w,y)`$. Introducing some additional repulsive interaction between the two 3-branes, one can keep them several (6-dimensional) Plank lengths apart to ensure the cancellation of the 4-dimensional cosmological constant. This requires a fine tuning of the new coupling constants.

The paper is organized as follows. In Sec. II we describe all possible embeddings of a 4-brane with 4-dimensional Poincaré invariance and detail our observations on the localization of gravity to the brane. Sec. III shows how intersecting branes can be embedded and explains the emerging additional parameter. In Sec. IV we show that these joining branes can be sewn together in a global manifold and illustrate on a simple model why the additional interaction strengths need to be fine tuned.

## II Static 4-branes in 6 dimensions

As explained in the introduction, we are looking for solutions to the six-dimensional Einstein equations:

$$\overline{G}_{ab}R_{ab}\frac{R}{2}g_{ab}g_{ab}=0.$$
(2)
These equations follow from the action with a number of branes included,

$$S[g]=_{bulk}(\mathrm{\Lambda }_62M_6^4R[g])+_{4brane(s)}\mathrm{\Lambda }_5,$$
(3)
while the units are so chosen that the bulk cosmological constant is set to $`1`$:

$$\frac{\mathrm{\Lambda }_6}{4M_6^4}1.$$
(4)
We are looking at metrics that satisfy the Ansatz of Eq. (1) and the Israel conditions  at $`w=0`$,

$$K_{ab}=\frac{\sigma }{2}h_{ab}.\text{ with }\sigma =\frac{\mathrm{\Lambda }_5}{8M_6^3}$$
(5)
These conditions ensure that keeping the $`w0`$ half plane and pasting back on another copy with $`ww`$ identified can accommodate a 4-brane with tension $`\sigma `$. Note that the critical (“fine tuned”) value of the 4-brane tension is at $`\sigma _c=\sqrt{\frac{2}{5}}`$ in our notation.

A direct calculation shows that the Einstein equations in the bulk become

$`0`$ $`=`$ $`\left(4{\displaystyle \frac{A_{yy}}{A}}+{\displaystyle \frac{A_y^2}{A^2}}2{\displaystyle \frac{A_y}{A}}{\displaystyle \frac{B_y}{B}}\right)+\left(3{\displaystyle \frac{A_w^2}{A^2}}+2{\displaystyle \frac{A_w}{A}}{\displaystyle \frac{B_w}{B}}2\right)B`$ (6)
$`0`$ $`=`$ $`\left(6{\displaystyle \frac{A_{yy}}{A}}3{\displaystyle \frac{A_y}{A}}{\displaystyle \frac{B_y}{B}}\right)+\left(6{\displaystyle \frac{A_{ww}}{A}}+3{\displaystyle \frac{A_w}{A}}{\displaystyle \frac{B_w}{B}}+2{\displaystyle \frac{B_{ww}}{B}}{\displaystyle \frac{B_w^2}{B^2}}4\right)B`$ (7)
$`0`$ $`=`$ $`3{\displaystyle \frac{A_y^2}{A^2}}\left(2{\displaystyle \frac{A_w^2}{A^2}}4{\displaystyle \frac{A_{ww}}{A}}\right)B`$ (8)
$`0`$ $`=`$ $`{\displaystyle \frac{A_yA_w}{A^2}}2{\displaystyle \frac{A_{yw}}{A}}+{\displaystyle \frac{A_w}{A}}{\displaystyle \frac{B_w}{B}}`$ (9)
while the Israel conditions are satisfied at $`w=0`$,

$$\frac{A_w}{A}=\sigma \text{ and }\frac{B_w}{B}=\sigma .$$
(10)
In the following we use the notation $`A(w,y)=e^{a(w,y)}`$. Observe from Eq. (8) that we either have $`a_y=0`$ and $`25a_w^24a_{ww}=0`$ simultaneously, or $`a_y0`$ and $`25a_w^24a_{ww}>0`$. By calculating the curvature tensor element of $`R_{wtt}^{}{}_{}{}^{w}`$ we observe that a necessary condition to have a maximally symmetric bulk solutionThat is, one satisfying $`R_{abcd}=\frac{R}{30}\left(g_{ac}g_{bd}g_{ad}g_{bc}\right)`$. is $`25a_w^210a_{ww}0`$.

In order to solve these equations in the $`a_y0`$ case we first express $`B`$ from Eq. (8),

$$B=\frac{3a_y^2}{25a_w^24a_{ww}},$$
(11)
and substitute this back into Eq. (9). The result contains only $`w`$-derivatives and can be solved as an ordinary differential equation for $`a_w`$ along the $`y=const`$ lines, which by construction are the brane orthogonal geodesics with affine parameter $`w`$:

$$_w\frac{\left(25a_w^210a_{ww}\right)^2}{\left(25a_w^24a_{ww}\right)^5}=0.$$
(12)
This is a first order differential equation on $`a_w`$ with boundary condition $`a_w(w=0)=\sigma `$, so that once the $`w`$-independent function $`c_1(y)=\left(25a_w^210a_{ww}\right)^2/\left(25a_w^24a_{ww}\right)^5`$ is specified it has one unique solution. It is not hard to see that Eq. (7) is a consequence of Eqs. (8,12). The boundary condition on $`B_w`$ trivially follows from Eq. (9). Now Eq. (6) relates the function $`c_1(y)`$ to $`a_0(y)a(0,y)`$. Substituting $`B`$ from Eq. (8) into Eq. (6) we find

$$0=_y\left(A\sqrt[5]{25a_w^210a_{ww}}\right)_y\left(\left[A\left(25a_w^24a_{ww}\right)\right]^5c_1(y)\right).$$
(13)
The term in the last parenthesis, $`e^{5c_0}\left[A\left(25a_w^24a_{ww}\right)\right]^5c_1(y)`$, is independent of $`w`$, as can be checked by showing that its $`w`$-derivative is proportional to Eq. (12). Then Eq. (6) is equivalent to $`c_0=const`$. We can calculate, with $`\eta =sign\left(25a_w^210a_{ww}\right)=\pm 1`$,

$$c_1(y)=\frac{\mathrm{exp}5[c_0a_0(y)]}{\left\{\frac{3}{5}(25\sigma ^2)+\frac{2}{5}\eta \mathrm{exp}\frac{5}{2}[c_0a_0(y)]\right\}^5}.$$
(14)
We have found all stable embeddings of a 4-brane with tension $`\sigma `$ satisfying the Ansatz in terms of the arbitrary constant $`c_0`$, the sign $`\eta `$ and an arbitrary function $`a_0(y)`$. The latter specifies the factor $`A`$ in the induced metric on the brane.

The unique solution of Eq. (12) can be written in an implicit form as follows. First introduce the quantity $`v(y,w)=\frac{3}{5}\frac{25a_w^2}{25a_w^24a_{ww}}`$, which can be found from solving the fifth order equation $`(25\tau ^2)^3\overline{c}_1(y)=v^3(v1)^2`$ in terms of $`\tau a_w`$ and $`\overline{c}_1(y)\frac{4\times 27}{5^5}c_1(y)`$; this is actually a rewriting of the definition of $`c_1(y)`$. Multiple solutions of this equation will lead to the different solutions mentioned in the introduction. Using the result we calculate the integral

$$w=_\sigma ^{a_w(w,y)}\frac{d\tau }{r(\tau ,y)}\text{ with }r(\tau ,y)a_{ww}=\frac{3}{20}(25\tau ^2)\left(\frac{5}{3}\frac{1}{v}\right),$$
(15)
which determines $`a_w`$ as a function of $`w,y`$. Integrating it, $`a(w,y)=a_0(y)+_0^wa_w(w,y)𝑑w`$, provides us with the function $`A(w,y)\mathrm{exp}a(w,y)`$. The function $`B(w,y)`$ is then found from Eq. (8).

The various solutions we are finding contain a great deal of redundancy that can be absorbed into coordinate transformations. A redefinition of the scale of the $`x^\mu `$ coordinates, $`x^\mu \lambda x^\mu `$ changes only $`A\lambda ^2A`$, $`a_0a_0+2\mathrm{log}\lambda `$, $`c_0c_0+2\mathrm{log}\lambda `$. This freedom can be used to set $`c_0=0`$.

The remaining freedom in arbitrarily choosing the function $`a_0(y)`$ locally corresponds to a coordinate transformation $`y\overline{y}(y)`$. This can be seen from the fact that as long as $`a_0^{}(y)0`$, the unique choice $`\overline{y}(y)=A(0,y)`$ provides $`\overline{A}(0,\overline{y})=\overline{y}`$ and the solution has no more continuous parameters left.

Instead of the above choice of the $`\overline{y}`$ coordinate, however, we find it more instructive to use coordinates in which the induced metric is conformally flat, i.e. $`A(0,y)=B(0,y)`$. Such coordinates always exist and can be found from any solution by the coordinate transformation to $`\overline{y}=^y\sqrt{B/A}𝑑y`$. Now expressing $`a_y`$ from Eq. (11) and using the definition of $`c_0`$ together with Eq. (14) we find (for w=0)

$$A_y=\pm \sqrt{\frac{25\sigma ^2}{5}A^3+\frac{2}{15}\eta e^{\frac{5}{2}c_0}\sqrt{A}}.$$
(16)
Using the variable $`u(y)=u_0/\sqrt{A(0,y)}`$ with $`u_0=\eta \frac{exp\frac{c_0}{2}}{\sqrt[5]{\frac{3}{2}(25\sigma ^2)}}`$ we write this in the form

$$\left(\frac{du}{dy}\right)^2=\left(\frac{25\sigma ^2}{20}u_0^2\right)\left(1u^5\right).$$
(17)
Note that, only on the brane,The general relationship is $`1\frac{1}{v}=\left(\frac{u_0^2}{A}\right)^{\frac{5}{2}}\frac{25\sigma ^2}{25a_w^2}`$. the variables $`u`$ and $`v`$ are related by $`1\frac{1}{v}=u^5`$.

### A Classification of the solutions

The various cases can be best understood by following the change of $`v`$ in Fig. 1 as we move along the brane. In the case $`A`$ when $`\sigma ^2>\frac{2}{5}`$, only $`u>1`$, i.e. $`v<0`$ can satisfy Eq. (17). It has one solution, unique up to a translation of the $`y`$-coordinate,

$$u(y)=F_s\left[y\times u_0\sqrt{\frac{5\sigma ^22}{20}}\right]\text{ i.e. }A(0,y)=\frac{u_0^2}{F_s^2\left[y\times u_0\sqrt{\frac{5\sigma ^22}{20}}\right]},$$
(18)
where the special function $`F_s(\overline{y})`$ is defined as the solution of $`|\overline{y}|=_0^{F_s}\frac{du}{\sqrt{u^51}}`$, see Fig. 2. The conformal factor of the induced metric, A(0,y) has a maximum at $`y=0`$ and vanishes at $`y=\pm y_{}=\pm \overline{y}_{}u_0\sqrt{\frac{5\sigma ^22}{20}}`$. By direct calculation of the curvature of the induced metric one can establish the asymptotics, $`R[g_{ab}^{(5)}]\left(y_{}y\right)^{\frac{10}{3}}`$. The induced metric has curvature singularities at these points which are easily seen to be at finite proper distance.

When $`\sigma ^2<\frac{2}{5}`$ we must have $`u<1`$, i.e. $`v>0`$ to satisfy Eq. (17). Again, we the solution can be written in terms of a special function as

$$u(y)=F_l\left[y\times u_0\sqrt{\frac{25\sigma ^2}{20}}\right]\text{ i.e. }A(0,y)=\frac{u_0^2}{F_l^2\left[y\times u_0\sqrt{\frac{25\sigma ^2}{20}}\right]},$$
(19)
where the special function $`F_l(\overline{y})`$ is now defined as the solution of $`|\overline{y}|=_{F_l}^1\frac{du}{\sqrt{1u^5}}`$, plotted in Fig. 3. In the two equivalent regions where $`F_l<0`$, case $`B`$, $`u`$ is negative and $`0<v<1`$. The conformal factor of the induced metric $`A(0,y)`$ tends to zero at a curvature singularity at $`yy_{}`$, located at finite proper distance, and the scalar curvature diverges as $`R[g_{ab}^{(5)}]\left(yy_{}\right)^5`$. The function $`A(0,y)`$ increases with $`y`$ and at $`yy_{}`$ diverges. The asymptotics is $`A(0,y)\frac{20}{25\sigma ^2}\frac{1}{(y+y_{})^2}`$, corresponding to an asymptotically $`AdS_5`$ induced metric. Consequently, the point $`y=y_{}`$ is at infinite proper distance away, and the regions $`B`$ and $`C`$ should be understood as different solutions. In the region $`y_{}<y<y_{}`$, case $`C`$, the conformal factor has a minimum and at the edges it diverges. The solution in this region sees $`AdS`$ horizons in both directions.

The definition of $`c_0`$ according to Eq. (13) supposes that the solution is not maximally symmetric, $`25a_w^210a_{ww}0`$. If it is, corresponding to $`c_0\mathrm{}`$, we have a fourth type of solution . This case $`D`$ corresponds to $`v(w,y)1`$ as one can see from the relation $`\frac{1}{v}1=\frac{2}{3}\frac{25a_w^210a_{ww}}{25a_w^2}0`$. In terms of Fig. 1, the solution is “sitting” in the $`AdS`$ point. In our coordinates its form is

$$ds^2=dw^2+\frac{25\sigma ^2}{2}\mathrm{cosh}^2\left[\frac{w}{\sqrt{10}}+\text{arctanh}\left(\sqrt{\frac{5}{2}}\sigma \right)\right]\times \left\{dy^2+\eta _{\mu \nu }dx^\mu dx^\nu e^{y\sqrt{\frac{25\sigma ^2}{5}}}\right\}.$$
(20)
The induced metric on the brane is obviously $`AdS_5`$, as it should be.

The solutions we were finding for fixed values of $`\sigma `$ represent 6-dimensional spacetimes with an embedded 4-brane. They have no parameters (other than $`\sigma `$) that could not be absorbed into coordinate redefinitions. We will see in the next section, that different values of $`\sigma `$ actually correspond to embedding branes with different tensions in the same bulk metric and there are only three solutions altogether, corresponding to the cases $`A,B,C`$ above. No 4-branes can be embedded (within the confines of the metric Ansatz of Eq. (1)) in any other 6-dimensional spacetime.§§§Excluding, of course, the case of a “fine tuned” brane, $`\sigma ^2=\frac{2}{5}`$, that we do not consider here.

In the exceptional case when $`A(0,y)=const.`$ along the brane it is not hard to see that there is no solution to our equations.

### B A note on localization of gravity to the brane

In the brane GNC’s we are using, the large-$`w`$ asymptotic behavior of the functions $`A`$ and $`B`$ can be better understood in the form

$$A(w,y)B(w,y)=\left(\frac{1\sqrt{\frac{5}{2}}\sigma }{2}\right)^2e^{\sqrt{\frac{2}{5}}w}+\frac{1}{2}\left(1\frac{5}{2}\sigma ^2\right)+\left(\frac{1+\sqrt{\frac{5}{2}}\sigma }{2}\right)^2e^{\sqrt{\frac{2}{5}}w}.$$
(21)
Far from the brane the first term dominates and evidently deconfines gravity. The form of the metric in the Randall-Sundrum model  is

$$ds^2=dw^2\eta _{\mu \nu }dx^\mu dx^\nu e^{2kw},$$
(22)
and the quick fall-off of the “warp factor” $`e^{2kw}`$ is essential in localizing gravity to the brane (at least in the limit where the Kaluza-Klein excited modes can be neglected). What is happening?

The resolution of the paradox hinges on two points. First, observe that a change in $`\sigma `$ corresponds to a shift in $`w`$, in addition to rescaling $`y`$ and $`x^\mu `$, so that a brane with a different tension could be placed parallel to the original one in the same bulk. The RS brane corresponds to $`\sigma ^2\frac{2}{5}`$, i.e. it is infinitely far away in $`w`$. The $`\sigma \sigma _c`$ limit is not smooth.That the limit is not smooth can also be seen in conformal coordinates, $`\rho >0`$, $`\mathrm{}<z<+\mathrm{}`$, $`ds^2=\frac{10}{\rho ^2}\left(\eta _{\mu \nu }dx^\mu dx^\nu +dz^2+d\rho ^2\right)`$. The original $`w`$ coordinate corresponds to the new angular coordinate. The branes are radial straight half lines in the $`(z,\rho )`$ plane while the fine tuned brane is along $`\rho =0,z>0`$. The limit is not smooth because the metric diverges as $`\rho 0`$. Now consider a brane with tension just smaller than the fine tuned value, $`\sigma _c\sigma \sigma _c`$. The behavior of the functions $`A,B`$ for increasing $`w`$ is first governed by the second term, the warp factor falls off and there is brane-localized gravity. By the time the exponential compensates for the smallness of its coefficient and the first two terms start to dominate, the warp factor is already $`𝒪\left(\frac{1\sigma /\sigma _c}{2}\right)1`$. It will grow again exponentially.

Second, the argument in  that leads to localization involves a Kaluza-Klein reduction, i.e. dropping the brane-orthogonal components of the metric (such as $`g_{ww},g_{w\mu }`$), in addition to the Kaluza-Klein expansion of the brane-tangential modes (i.e. $`g_{\mu \nu }`$). While close to the brane the separation of modes into brane-orthogonal and brane-tangential makes sense, far away it has no covariant meaning. The actual separation and large-$`w`$ behavior depends on the details of the model that explains the Kaluza-Klein reduction and cannot be expected to be faithfully represented by the general framework of . One may conjecture by the same token that the increase in the warp factor is unphysical and should be thrown away.

In any case the conclusion remains that in any scenario where the brane tensions are not close to their fine tuned values does not automatically provide localized gravity. In order to make a coordinate invariant statement one needs to restrict to such coordinates in which the thrown-away modes ($`g_{ww},g_{w\mu }`$) are massive. It cannot be decided without the explicit identification of these modes whether the warp factor along the resulting $`y=const`$ lines shows the same behavior as Eq. (21). This is equally relevant for codimension one branes in five dimensions.

In this respect our most interesting solution, $`C`$, which has no curvature singularities, can be seen to have the asymptotic behavior

$$A(w,y)function(y)\times e^{\sqrt{\frac{2}{5}}|w|}\left(1+2\frac{1+\sqrt{\frac{5}{2}}\sigma }{1\sqrt{\frac{5}{2}}\sigma }e^{\sqrt{\frac{2}{5}}|w|}+𝒪\left(e^{\sqrt{\frac{2}{5}}|w|}\right)^2\right),$$
(23)
the same as the asymptotics of the known $`AdS`$ solutions.

## III Intersecting branes

In the previous section we found all possible embeddings in GNC of the 4-brane. This still does not answer the question whether another brane can be inserted in the same bulk because for intersecting branes the two GNC’s certainly do not coïncide. In the CNC of one brane the other will not in general be straight. We were not able to solve the equations for the branes in closed form, but we could (almost) show their existence/unicity, and determine the relationship between the angle of intersection and the brane tensions.

We start with one solutions which permits a 4-brane with tension $`\sigma `$ located at $`w=0`$, keeping the $`W>0`$ half brane. If effect we are considering a half-brane at $`yy_0`$, and looking for another half-brane with tension $`\stackrel{~}{\sigma }`$, starting off at some angle $`0<\varphi <+\mathrm{}`$, and keeping a wedge between them as shown in Fig. 4 (the wedge metric will of course be $`2\pi `$-periodic if $`\varphi >2\pi `$.)

The condition for the stability of the $`\stackrel{~}{\sigma }`$ brane is again the Israel condition $`K_{ab}=+\frac{\stackrel{~}{\sigma }}{2}h_{ab}`$, where the opposite sign is due to the fact that this time we are keeping the bulk to the opposite side of the brane. We describe the brane with a pair of functions

$$\stackrel{~}{w}(\lambda ),\stackrel{~}{y}(\lambda )\text{ with }\stackrel{~}{w}(0)=0,\stackrel{~}{y}(0)=y_0\text{ and }0\lambda <\mathrm{}.$$
(24)
The Israel conditions become

$`2\left(\ddot{\stackrel{~}{y}}\dot{\stackrel{~}{w}}\ddot{\stackrel{~}{w}}\dot{\stackrel{~}{y}}\right)+b_y\dot{\stackrel{~}{w}}\dot{\stackrel{~}{y}}^2+b_w\dot{\stackrel{~}{w}}\dot{\stackrel{~}{y}}\left(2\dot{\stackrel{~}{w}}+B\dot{\stackrel{~}{y}}\right)`$ $`=`$ $`{\displaystyle \frac{\stackrel{~}{\sigma }}{\sqrt{B}}}`$ (25)
$`{\displaystyle \frac{a_y}{\sqrt{B}}}\dot{\stackrel{~}{w}}a_w\sqrt{B}\dot{\stackrel{~}{y}}`$ $`=`$ $`\stackrel{~}{\sigma }`$ (26)
with $`B=e^b`$, where we fixed the $`\lambda `$ variable by requiring

$$\dot{\stackrel{~}{w}}^2+B\dot{\stackrel{~}{y}}^2=1.$$
(27)
Now Eq. (26) can be viewed as a first order differential equation for the function $`\stackrel{~}{w}(\stackrel{~}{y})`$ with initial conditions set by $`\stackrel{~}{w}(y_0)=0`$ and by Eq. (27). Up to a discreet ambiguity, this system has a unique solution. We could not show in complete generality that Eq. (25) is a consequence of Eq. (26), we expanded the solution of Eq. (26) in a power series around $`\lambda =0`$ keeping terms up to $`𝒪(\lambda ^5)`$ and saw that they did indeed satisfy Eq. (25), up to $`𝒪(\lambda ^3)`$. (The loss of two powers of $`\lambda `$ is due to the three derivatives figuring in Eq. (25).) We take this as an indication that Eq. (25) is indeed not an independent equation.

We can define the angle between the in a covariant way by requiring $`\mathrm{cos}\varphi =n_ag^{ab}\stackrel{~}{n}_b`$, where $`n_a`$ and $`\stackrel{~}{n}_b`$ are the unit normals of the branes at the intersection. In our coordinates we find $`\left(\mathrm{sin}\varphi |\mathrm{cos}\varphi \right)=\left(\dot{\stackrel{~}{w}}|\sqrt{B}\dot{\stackrel{~}{y}}\right)_{\lambda =0}.`$ Then, substituting this definition into Eq. (26) we see that the brane tension and the position of the junction (i.e. $`y_0`$) determines the angle:

$$\frac{a_y}{\sqrt{B}}\mathrm{sin}\varphi +\sigma \mathrm{cos}\varphi +\stackrel{~}{\sigma }=0.$$
(28)
A trivial check of this calculation is that $`\varphi =\pi `$, $`\stackrel{~}{\sigma }=\sigma `$ is always a solution: the brane can always be continued across the point $`y_0=0`$. We also see that once Eq. (28) is satisfied, there is one unique solution to Eqs. (25-27).

Before discussing the solutions to Eq. 28 we look at the induced metric on the $`\stackrel{~}{\sigma }`$ brane. A simple substitution shows $`d\stackrel{~}{s}^2=d\lambda ^2+A\eta _{\mu \nu }dx^\mu dx^\nu `$, which can be written in a conformally flat from by introducing the new variable $`Y=_0^\lambda \frac{d\lambda }{\sqrt{A}}`$ along the new brane: $`d\stackrel{~}{s}^2=A\left(dY^2+\eta _{\mu \nu }dx^\mu dx^\nu \right)`$. Note that the $`x^\mu `$ coordinates did not have to be rescaled, and the conformal factor A is the same at the junction in the two induced metrics. The bulk metric can now be written in GNC’s of the $`\stackrel{~}{\sigma }`$ brane in the form $`ds^2=dW^2+\stackrel{~}{B}(W,Y)dY^2+\stackrel{~}{A}(W,Y)\eta _{\mu \nu }dx^\mu dx^\nu `$; we know $`\stackrel{~}{A}(0,Y)=\stackrel{~}{B}(0,Y)=A(\stackrel{~}{w},\stackrel{~}{y})`$. But this must be one of our previously found solutions which comes with its $`\stackrel{~}{v}(W,Y)`$ function. A direct calculation of $`\stackrel{~}{v}(0,0)`$ gives, using Eq. (28), and observing $`v=\frac{25\sigma ^2}{5}\frac{\stackrel{~}{B}}{\stackrel{~}{a}_Y^2}`$

$$\stackrel{~}{V}\stackrel{~}{\sigma }^2+\frac{25\stackrel{~}{\sigma }^2}{5\stackrel{~}{v}(0,0)}=V\sigma ^2+\frac{25\sigma ^2}{5v(0,0)}.$$
(29)
In other words, the quantity $`V`$ must be the same on the two sides of the wedge.

It will prove important to note now that once $`\stackrel{~}{\sigma },\stackrel{~}{v}(0,0)`$ and $`\stackrel{~}{A}(0,0)`$ are specified as they are, there are only two possible induced metrics on the new brane. The ambiguity results of the two ways to invert the functions $`F_{l,r}`$ and can be resolved by also specifying the sign $`\stackrel{~}{\kappa }=\text{sign}\stackrel{~}{a}_Y(0,0)=\pm 1`$. Then we will be able to argue that any two induced metrics with the same $`\sigma ,v(0,0),A(0,0)`$ and $`\kappa `$ coïncide. This will be the situation when a number of wedges is pasted together around one 3-brane.

The solutions of Eq. (28) can be visualized by the following auxiliary device. It is always possible to find a unique solution $`\varphi =\phi `$ to Eq. (28) with $`\stackrel{~}{\sigma }=0`$, also requiring $`\stackrel{~}{\kappa }=+1.`$ This corresponds to a “null brane” represented in Fig. 5. It is not hard to see that we could have started with the null brane and looked for the angles at which the two $`\sigma ,\stackrel{~}{\sigma }`$ branes to cut out the wedge. The equations for their angles become

$$\frac{a_y}{\sqrt{B}}\mathrm{sin}\stackrel{~}{\phi }+\stackrel{~}{\sigma }=0\text{ and }\frac{a_y}{\sqrt{B}}\mathrm{sin}\phi +\sigma =0,$$
(30)
where $`\frac{a_y}{\sqrt{B}}=\sqrt{V}`$ now refers to the null brane. The solutions that give the correct sign for the $`a_y`$’s are, with $`\frac{\pi }{2}\mathrm{arcsin}\frac{\pi }{2}`$,

$$\phi =\pi \frac{1\kappa }{2}\kappa \mathrm{arcsin}\frac{\sigma }{\sqrt{V}}+2\pi N\text{ and }\stackrel{~}{\phi }=\pi \frac{1\stackrel{~}{\kappa }}{2}\stackrel{~}{\kappa }\mathrm{arcsin}\frac{\stackrel{~}{\sigma }}{\sqrt{V}}+2\pi \stackrel{~}{N}.$$
(31)
Based on the “compatibility” of different branes, the possible setups can be classified as follows. It is not possible to move along any brane from one category to another.

* $`0<V<\frac{2}{5}`$, only branes with below critical tension, $`\sigma ^2<\frac{2}{5}`$ can be used, all type $`C`$ (i.e. $`v>1`$) with no curvature singularity. The angles are determined by the tensions and the parameter $`V`$ which is fixed in any one junction but varies from junction to junction
* $`V=\frac{2}{5}`$, the spacetime is $`AdS_6`$, only type $`D`$ branes with below critical tension, $`\sigma ^2<\frac{2}{5}`$ can be used. The angles are determined by the tensions and there is no additional parameter.
* $`V>\frac{2}{5}`$, one can accommodate both type $`B`$ branes with $`0<\sigma ^2<\frac{2}{5}`$ and type $`A`$ branes with $`\sigma ^2>\frac{2}{5}`$. There is again a parameter $`V`$ in the junctions. All induced metrics contain curvature singularities.

## IV Sawing together the wedges

Any combination of intersecting branes with angles between them set according to Eq. (31) represents a stationary point of the action in Eq. (3). The $`V`$ are parameters of the solution, one at each intersection. We must make sure, however, that a globally defined metric exists when the pieces we have discussed are sewn together.

Each brane carries orbifold boundary conditions, which are easiest to state in GNC’s: we identify points with the sign of the corresponding coordinate $`w`$ flipped. This allows at most two different sets of 4-brane tensions ($`\sigma `$) in each intersection: every second one most be equal. The corresponding signs $`\kappa `$ should also be chosen equal. In we join an odd number of branes in one intersection, then all must have the same $`\sigma `$ and $`\kappa `$.

When $`n=1,2,\mathrm{}`$ branes intersect along one 3-brane, the $`n+1^{th}`$ one must be identified with the first one. There are two conditions that must be met for consistency. (i) The induced metric on the 3-brane must be the same for all wedges. This condition is always satisfied because the functions $`A(0,0)`$ are always the same and no rescaling of the $`x^\mu `$ coordinates was necessary. (ii) The induced metric on the $`n+1^{th}`$ brane must be identical with the induced metric on the $`1^{st}`$. This condition is again automatically satisfied because, as we saw, the equality of $`V,A,\kappa `$ and $`\sigma `$ is sufficient for that.

In the following we consider a case with one positive tension 4-brane with two junctions around which the brane is identified with its own reverse side, shown in Fig. 6. We require, due to the effectively two dimensional gravity on the $`w,y`$ plane, an excess angle $`0<\lambda <\pi `$. This is possible only with $`\kappa =1`$; and in order to have no curvature singularities, we will use type $`C`$ branes. Because then $`a_0(y)`$ must increase towards both “junctions”, it must have a minimum (at $`y=0`$) on the brane between then. Using the excess angles $`\lambda _{1,2}`$ as parameters, Eq. (31) tells the parameter (on each end separately) $`V=\frac{\sigma ^2}{\mathrm{sin}^2\frac{\lambda }{2}}`$. Because for the $`C`$ type branes $`V<\frac{2}{5}`$, we can accommodate only an excess angle at least as large as $`\mathrm{sin}\frac{\lambda }{2}>\sqrt{\frac{5}{2}}\sigma `$. The 3-brane tension cannot be much smaller than the 4-brane tension. In brane GNC’s the conformal factor is $`A(0,y)=F_l^2\left[y\sqrt{\frac{25\sigma ^2}{20}}\right]`$, shown in Fig. 7. We calculated the position of the two ends of the brane, with the result

$$y_j=\frac{(1)^j}{\sqrt{\frac{25\sigma ^2}{20}}}_{u_j}^1\frac{du}{\sqrt{1u^5}}\text{ with }j=1,2.$$
(32)
The proper distance of the endpoints from the origin is

$$l_j=\frac{(1)^j}{\sqrt{\frac{25\sigma ^2}{20}}}_{u_j}^1\frac{du}{u\sqrt{1u^5}}\text{ with }u_j=\left[\frac{2\mathrm{sin}^2\frac{\lambda _j}{2}5\sigma ^2}{(25\sigma ^2)\mathrm{sin}^2\frac{\lambda _j}{2}}\right]^{\frac{1}{5}}.$$
(33)
Consistency requires that this distance should be $`l1`$ (otherwise the underlying field theory description of gravity is not justified). This requires that the 4-brane tension not be much smaller than $`\lambda `$ and also either $`25\sigma ^21`$ or $`0<\mathrm{sin}\frac{\lambda }{2}\sqrt{\frac{5}{2}}\sigma 1`$.

The presence of a vacuum energy density on the 3-brane requires the proper amount of excess angle. When this relationship is not satisfied, our Ansatz does not give a solution: the induced metric on the 3-brane is not flat any more. It is a trivial exercise to derive the equations of motion when a maximally symmetric but not flat 4-dimensional metric is allowed, i.e. $`dt^2+d𝐱^2`$ is replaced by $`dt^2+d𝐱^2e^{2ht}`$ in Eq. (1). We find that Eqs. (9,12) are not modified, so that the function $`A(w,y)`$ is not modified. The modification occurs because Eq. (11) gets replaced by

$$B=\frac{3a_y^2}{25a_w^24a_{ww}+12\frac{h^2}{A}},$$
(34)
so that $`B(w,y)`$ changes. This change also affects the intersection angles,

$$\frac{\mathrm{tan}^2\phi _h}{\mathrm{tan}^2\phi }=\frac{V}{V+4\frac{h^2}{A}}.$$
(35)
If the “mistuning” is little, $`\delta \lambda =\lambda _h\lambda =2(\varphi _h\varphi )1`$, the result is a small $`h^2=\frac{AV}{2\mathrm{sin}\lambda }\delta \lambda `$. This translates to a small observable Hubble constant (the induced metric must be brought to the form $`ds_4^2=dt^2+d𝐱^2e^{2Ht^{}}`$, $`H=h/\sqrt{A}`$),

$$H^2=\frac{\sigma ^2}{\mathrm{sin}^2\frac{\lambda }{2}}\times \frac{\delta \lambda }{2\mathrm{sin}\lambda }.$$
(36)
Given that $`\lambda `$ cannot be made much smaller than $`\sigma `$, any sizable “mistuning” of $`\lambda `$ would result in a Hubble constant on the 6-dimensional Planck scale.

Up to this point however, we have not required any fine tuning. The 4-dimensional vacuum density, through the topologically required conical singularity of the metric enforces the adjustment of the parameter $`V`$ to its required value. However, when a 3-brane tension is added to the action,

$$S[g]=_{bulk}(\mathrm{\Lambda }_62M_6^4R[g])+_{4brane(s)}\mathrm{\Lambda }_5+_{3brane}\mathrm{\Lambda }_4,$$
(37)
the added term can be written as

$$_{3brane}\mathrm{\Lambda }_4=\mathrm{\Lambda }_4A^2(0,0)_{3brane}d^4x.$$
(38)
In order to avoid that the spacetime collapses onto the point where $`A(0,0)`$ is minimum, i.e. $`y=0`$, additional repulsive interactions have to be introduced. Their strength has to be fine tuned to ensure that the equilibrium length of the 4-brane is that of Eq. (33).

## V Conclusions

The main result of this paper is that there exists embeddings of 4-branes in less than maximally symmetric 6-dimensional spacetimes that preserve 4-dimensional Poincaré invariance. Then different points on the 4-brane are not physically equivalent as they were in the $`AdS`$ case. The angle between intersecting branes, in any fixed intersection point, is always determined by the nondynamical parameters through the Israel conditions $`K_{ab}=\frac{\sigma }{2}h_{ab}`$. We saw that there is no obstacle to building a global metric around the intersection and the additional parameter can be translated to an arbitrary excess angle. We also saw that there is no solution (without fine tuning the 4-brane tensions) with constant warp factor along the 4-brane. As a consequence, as we illustrated on a simple example, a 3-brane with tension tries to “move” to the mininum of that warp factor and additional interactions are needed to stabilize the 3-brane. We note that in all suggested solutions to the cosmological constant problem  in brane contexts some sort of tuning of the bulk interactions has been required.

## VI Acknowledgments

The author thanks Bob Holdom for many useful discussions and stressing the importance of finding flat directions, as well as Hael Collins for discussions.


format: text

text:
# Hall state quantization in a rotating frame

## Abstract

We derive electromagnetomotive force fields for charged particles moving in a rotating Hall sample, satisfying a twofold U(1) gauge invariance principle. It is then argued that the phase coherence property of quantization of the line integral of total collective particle momentum into multiples of Planck’s quantum of action is solely responsible for quantization in the Hall state. As a consequence, the height of the Hall quantization steps should remain invariant in a rapidly rotating Hall probe. Quantum Hall particle conductivities do not depend on charge and mass of the electron, and are quantized in units of the inverse of Planck’s action quantum.

Modern molecular beam epitaxy enables the preparation of modulation-doped semiconductor heterostructures in which, at low enough temperatures, a high mobility two-dimensional electron gas is formed. This system is characterized by a long Thouless dephasing length $`l_\varphi `$, the distance within which phase coherence of mobile electrons is maintained. In good samples, and at low enough temperatures, the length $`l_\varphi `$ reaches several micrometers, exceeding the magnetic length $`l_B=\sqrt{\mathrm{}/eB}`$ for applied magnetic fields of order one Tesla. Under these conditions, it should be possible to detect noninertial effects due to rotation or acceleration of the sample as a result of the change of quantum interference conditions, since the gauge potentials of the electromagnetic and noninertial fields experienced by the electrons both appear in their collective phase. In what follows, we shall argue that quantum coherence under the influence of noninertial force fields is directly observable in the quantum Hall effect . The quantum of Hall conductivity for particle transport is given by the inverse of Planck’s quantum of action alone, and involves no properties specific to the electron. The arguments used to prove this result crucially rely on the existence of a collective particle momentum expressing quantum coherence. By including gauge fields other than the electromagnetic one, we promote the idea that Hall quantization is of necessity derivable from a topological quantum number related to this coherence. It will be shown that this prediction about the nature of the quantum Hall effect is verifiable within current technological means.

The main quantity of interest to us is the total particle momentum

$$𝒑=m𝒗+m𝛀\times 𝒓+q𝑨.$$
(1)
Here, $`𝑨`$ is the electromagnetic vector potential, $`q`$ the charge of the electron in vacuo, $`𝒗`$ the particle velocity, and $`m`$ the inertial mass. The body is rigidly rotating with respect to the laboratory frame at a (time dependent) angular velocity $`𝛀`$. The noninertial force on a massive test particle inside a rotating Hall probe, as measured in the rotating sample frame, is then given by the standard expression

$$𝑭=2m𝛀\times 𝒗m𝛀\times 𝛀\times 𝒓m_t𝛀\times 𝒓.$$
(2)
The first term on the right hand side represents (minus) the Coriolis force, the second one the centripetal force, and the last term is due to temporal changes of the angular velocity. The presence of this last term will prove to be crucial for our argument on proper Hall quantization in a rotating frame presented below. We omit possible additional terms on the right hand side of equation (2) due to potential forces (e.g. gravity) or externally imposed linear acceleration. Vector and scalar potentials associated to rotation are defined as follows

$$𝒂=𝛀\times 𝒓,a_0=\frac{1}{2}\mathrm{\Omega }^2𝒓_{}^2,$$
(3)
where $`𝒓_{}`$ is the distance vector perpendicular to the axis of rotation. For a charged massive particle like the electron, we merge these potentials and the electromagnetic potentials into a generalized vector potential, incorporating the coupling constants charge $`q`$ and mass $`m`$,

$$𝒜=q𝑨+m𝒂,$$
(4)
and a generalized scalar potential

$$\chi =qA_0ma_0.$$
(5)
The sum of the generalized electromotive and magnetomotive forces, acting on an electron , consisting of noninertial plus proper Lorentz and electric forces, then takes on the form

$`𝑭_{}=+𝒗\times ,`$ (6)
where the generalized electric and magnetic fields are

$``$ $`=`$ $`\chi _t𝒜,`$
$``$ $`=`$ $`\times 𝒜.`$ (7)
As a consequence of this relation for the total force, the usual expression for the drift velocity of the charge carriers, resulting from zero total force in perpendicular electric and magnetic fields, experiences the obvious modification that $`𝑬`$ and $`𝑩`$, so that $`𝒗_D=\times /^2`$.

The force fields displayed in equations (6) and (7) give a theory possessing in effect two U(1) gauge symmetries. The standard U(1) from electromagnetism, with coupling constant $`q`$ (charge), and another U(1) gauge symmetry, with coupling constant $`m`$ (inertial rest mass). The gauge potential of this second U(1) has a scalar part $`a_0`$ and a vectorial part $`𝒂`$. The homogeneous Maxwell equations $`\mathrm{rot}=_t`$ and $`\mathrm{div}=0`$ then follow from the existence of the potentials $`𝒜`$ and $`\chi `$ in (4) and (5). That the Faraday law holds is due to our admitting a variation of the angular velocity with time and the resulting last force term in (2). This important term, leading to gauge invariance in explicitly time dependent situations, has not been considered in , where the Hall effect under rotation, without the simultaneous existence of a magnetic field, was investigated.

The gauge invariant particle current induced by the electromotive force field $``$ is in linear response ($`a,b\{x,y\}`$):

$`𝑱_a^{\mathrm{ind}}`$ $`=`$ $`\stackrel{\mathbf{~}}{𝝈}_{ab}_b`$ (8)
$`=`$ $`\stackrel{\mathbf{~}}{𝝈}_{ab}\left(q𝑬+m𝒈\right)_b,`$
where $`𝒈=a_0_t𝒂`$ represents the ‘electric’ part of the mechanical acceleration experienced by the electron. In the present case this acceleration is purely caused by rotation, and takes the form $`𝒈=𝛀\times 𝛀\times 𝒓_t𝛀\times 𝒓`$.

Observe that the left hand side of equation (8) contains the number current density rather than the electric current density and that, dimensionally, the particle conductivity $`[\stackrel{~}{\sigma }]=[\sigma _{el}/q^2]`$. In the case of two coupling constants, $`m`$ and $`q`$, it is the number of particles crossing (in two spatial dimensions) a line of unit length per unit time, which is the relevant observable. This quantity is proportional to the electromotive force field $``$, which causes these particles to move. If we were to use the transport coefficient $`\sigma _{\mathrm{el}}`$ and the transport equation $`𝑱_{\mathrm{el}}^{\mathrm{ind}}=𝝈_{\mathrm{el}}𝑬`$, a rotating Hall sample does not yield the sharp conductivity quantization steps observed in nonrotating samples. We will show below that Hall state quantization is, according to equation (16), to be expressed in the particle Hall conductivity $`\stackrel{~}{\sigma }_{xy}`$ occuring in relation (8).

Evidence for the necessity of using the particle transport equation (8) comes from the existence of the London field in superconductors. Complete expulsion of the field $``$ deep inside in a superconductor requires the particle conductivity $`\stackrel{\mathbf{~}}{𝝈}`$ to have a contribution proportional to $`1/i\omega `$, which yields a term on the right hand side of (8), proportional to the generalized vector potential $`𝒜`$. Corresponding to complete Meissner type screening, $`=\mathrm{rot}𝒜=q𝑩+2m𝛀=0`$, the London spontaneous field $`𝑩_L`$ then takes the value

$$𝑩_L=2\frac{m}{q}𝛀.$$
(9)
This relation corresponds to zero winding number of the phase $`\theta `$, cf. equations (10)–(12) below. Equation (9) has been verified experimentally already 35 years ago , in an experiment in which it was used to infer the Compton wavelength of superconducting electrons. For the linear in velocity (nonrelativistic) limit and in a superconductor, the Cooper pair mass $`m`$ equals twice the electron inertial rest mass in vacuo, $`2m_e`$ (the outcome of a more recent high precision experiment using a rotating superconducting niobium ring  has been $`m/2m_e=1.000084(21)`$). If we insert on the left hand side of the equation (9) the bare electron values $`m=2m_e`$ and $`q=2e`$ ($`e>0`$), we have $`𝑩_L=(1.1410^{11}\mathrm{Tesla}\mathrm{sec})𝛀`$. Only the ratio of $`m`$ and $`q`$ enters the London induced magnetic flux strength. In a quantum Hall liquid, where the “elementary” quanta are $`q=e`$, $`m=m_e`$ and $`\varphi _0=2\pi \mathrm{}/e`$, instead of $`q=2e`$, $`m=2m_e`$ and $`\varphi _0=2\pi \mathrm{}/2e`$ in a superconductor, the London flux strength thus takes for a given $`𝛀`$ the same value (9) like in the superconductor. That the mass is exactly the bare inertial mass to linear order in the particle velocity is independently substantiated by a recent discussion of the London equation in , by using (thermodynamic) arguments different from our gauge invariance argument. If the quantum Hall fluid is described in a (relativistic) theory with interactions, the particle mass is to be replaced by the chemical potential, cf. , but the identity of $`m`$ with the bare mass to linear order in the velocity still persists.

Quantum coherence properties enter the stage if we require for the line integral of collective particle momentum along a closed path

$$𝒑𝑑𝒓=N_v2\pi \mathrm{},$$
(10)
where $`N_v`$ is the winding number of phase $`\theta `$, such that the total canonical momentum

$`𝒑`$ $``$ $`\mathrm{}\theta `$ (11)
$`=`$ $`m𝒗+m𝛀\times 𝒓+q𝑨`$
$`=`$ $`m𝒗+𝒜.`$
The uniqueness condition of the collective phase represented in (10) then leads, if we take a path in the bulk of the electron liquid, for which the integral of $`m𝒗`$ may be neglected, to the quantization of the sum of a Sagnac flux  and the magnetic flux:

$`\mathrm{\Phi }`$ $`=`$ $`q{\displaystyle 𝑨𝑑𝒓}+m{\displaystyle 𝛀\times 𝒓𝑑𝒓}`$ (12)
$`=`$ $`{\displaystyle 𝑑𝑺}=N_v2\pi \mathrm{}.`$
This flux quantization rule associated with the field $``$ corresponds to the fact that a vortex is fundamentally characterized by the winding number $`N_v`$ alone . No properties of the medium in which it lives, in particular the mass and charge of the medium’s constituents, enter the quantum of generalized flux.

Consider now the (purely magnetic) filling factor of a nonrotating two-dimensional electronic system of areal density $`n_{2D}`$ in a large magnetic field $`B`$ at low temperatures. The ratio $`\nu =n_{2D}/(B/\varphi _0)`$, where $`\varphi _0=2\pi \mathrm{}/e`$, gives the inverse of the number of singly quantized ($`N_v=1`$) magnetic flux quanta available per electron. The Hall resistance is quantized into $`R_H=(2\pi \mathrm{}/e^2)\nu ^1=R_K/\nu `$ (the von Klitzing constant $`R_K=25812.807`$ $`\mathrm{\Omega }`$), with integer or fractional $`\nu `$ . We have seen that the dynamical, generalized magnetic force field occuring in the Hamiltonian in the noninertial rotating state is $``$. The solution of the Landau problem for the electronic energy levels in a magnetic field thus refers to a Hamiltonian $`H=(𝒑𝒜)^2/2m+\chi `$, depending on the generalized vector potential $`𝒜`$ and magnetic field $``$. Hence, in the noninertial case, the Landau level degeneracy per unit area is $`/(2\pi \mathrm{})`$, and contains the generalized magnetomotive field $``$ instead of the magnetic flux strength $`𝑩`$. By assigning this value to the degeneracy, use is made of the fact that the “magnetic length” associated with rotation, $`l_\mathrm{\Omega }\sqrt{\mathrm{}/(2m\mathrm{\Omega })}`$, is much larger than the Thouless length as well as the proper magnetic length, $`l_\mathrm{\Omega }l_\varphi >l_B`$ . As a consequence, the linear dependence on position of the part of $``$ associated with rotation, given by $`m𝛀\times 𝛀\times 𝒓m_t𝛀\times 𝒓`$, does not lift the Landau level degeneracy beyond the broadening of the levels already taking place due to scattering. Under the condition $`l_\mathrm{\Omega }l_\varphi >l_B`$, the filling factors assigned to the Landau levels are

$$\nu _{}=\frac{n_{2D}}{/(2\pi \mathrm{})}.$$
(13)
The Faraday law $`\mathrm{rot}=_t`$, telling us how the flux strength corresponding to the vector potential $`𝒜`$ changes in time, is in its integrated form

$$_\mathrm{\Gamma }𝑑𝒓=\frac{d\mathrm{\Phi }}{dt}.$$
(14)
Consider the adiabatic process of slowly turning on a flux quantum $`2\pi \mathrm{}`$ in a nondissipative Hall state, which has $`\stackrel{~}{\sigma }_{xx}=\stackrel{~}{\sigma }_{yy}=0=\stackrel{~}{\rho }_{xx}=\stackrel{~}{\rho }_{yy}`$, and the antisymmetry property $`\stackrel{~}{\sigma }_{xy}=\stackrel{~}{\sigma }_{yx}`$ . The path $`\mathrm{\Gamma }`$ is led around the flux tube. The induced current then obeys $`\widehat{z}\times 𝑱^{\mathrm{ind}}=\stackrel{~}{\sigma }_{xy}`$, and the number of particles $`N`$ inside the area enclosed by $`\mathrm{\Gamma }`$ changes according to

$$\frac{dN}{dt}=\stackrel{~}{\sigma }_{xy}\frac{d\mathrm{\Phi }}{dt}.$$
(15)
After a single quantum of generalized flux $`2\pi \mathrm{}`$ has been added, and because the number of particles is integral, this implies that the off-diagonal part of the particle transport conductivity defined in (8) is quantized according to

$$\stackrel{~}{\sigma }_{xy}=\nu _{}/2\pi \mathrm{},$$
(16)
with $`\nu _{}`$ an integer . Hence, to summarize this derivation, the generalized Faraday law in a rotating frame gives in a nondissipative Hall state the quantization of the Hall conductivity if and only if the moving generalized flux in (12) is quantized in units of the action quantum. This argument is similar to the one given by Laughlin , extended to a rotating frame.

The quantization of the Hall resistance into $`R_K/\nu `$ for nonrotating samples has been measured to an absolute accuracy of a few parts in $`10^8`$ for an individual, specific sample, whereas in a comparative study of different materials, a relative accuracy of about $`10^{10}`$ of the ratio of Hall resistances has been achieved . Considering that magnetic fields in quantum Hall experiments cover a range $`B1\mathrm{}30`$ Tesla, this implies that with a rotation rate of the Hall sample of, say, $`\mathrm{\Omega }=10^3\mathrm{sec}^1\mathrm{}10^4\mathrm{sec}^1`$, the Sagnac contribution in (12) is large enough to verify if Hall quantization experiences a change if the sample rotates. For the comparative measurement, the rotation rates required are lower by about two orders of magnitude. At the very low temperatures of order $`10^3`$ K needed for superfluid $`^3`$He, rotation rates of the cryostat of order $`\mathrm{\Omega }10`$ sec<sup>-1</sup> already have been realized . In relation to superfluid $`^3`$He, it is also worthwhile to mention here that for thin $`^3`$He-A films, an electrically neutral system, a half integer quantum Hall effect, owing to a topological invariant of the $`p`$-wave order parameter in this system, has been suggested .

We point out that with respect to the interpretation of experimental results, it should be borne in mind that the force fields in (6) and (7) refer to the rotating sample frame, and not to the laboratory frame.

The experiment proposed, then, consists in a comparison of the quantum Hall resistances in the reference frame of a rotating sample as well as in a nonrotating sample. If identical quantization results are obtained, this yields direct proof that what is actually observed in the Hall experiment is the quantization according to (16), rather than quantization into $`e^2/2\pi \mathrm{}`$. The electric charge of the electron, the U(1) coupling constant of electromagnetism, enters if we count the number of particles arriving at the Hall contacts, by ascribing to their transport properties an electric conductance in an electric circuit. What is invariantly measured, though, is the induced number current in (8). Particle number currents can be distinguished from conventional electric currents as follows. Whereas electrical currents are measured by (classical) means of an impedance and the ensuing voltage drop, intrinsic particle currents can be measured if we use tunneling contacts, for which the wave-particle duality embodied in the tunneling process ensures that particles are counted.

We stress that if we were to use a Chern-Simons effective field theory for a description of the quantum Hall effect , such an effective description is expressible within our approach. For that purpose, one uses a term in the action density of the form

$$(J^{\mathrm{ind}})^\mu 𝒜_\mu +\frac{1}{4}\stackrel{~}{\sigma }_{xy}ϵ^{\alpha \beta \gamma }𝒜_\alpha _{\beta \gamma },$$
(17)
where $`ϵ^{\alpha \beta \gamma }=\pm 1`$ (with the sign convention $`ϵ^{0xy}=+1`$) is the unit antisymmetric symbol in three space-time dimensions and $`_{\mu \nu }=_\mu 𝒜_\nu _\nu 𝒜_\mu `$ is the field tensor constructed from $``$ and $``$. Such a term reproduces the induced current (8) in the nondissipative Hall state as a result of varying the action with respect to $`𝒜_i`$. Integration of the zeroth component of the complete field equations $`(J^{\mathrm{ind}})^\mu =\frac{1}{2}\stackrel{~}{\sigma }_{xy}ϵ^{\mu \alpha \beta }_{\alpha \beta }`$, giving the induced density

$$\rho ^{\mathrm{ind}}=\stackrel{~}{\sigma }_{xy}_z,$$
(18)
tells us that each particle associated with $`\rho ^{\mathrm{ind}}`$ carries generalized flux $`1/\stackrel{~}{\sigma }_{xy}=2\pi \mathrm{}/\nu _{}`$.

In conclusion, the quantization of the Hall particle conductivity under rotation has been derived by invoking the following basic arguments. (i) The generalized Lorentz force equation (6), containing the invariant field strengths $``$ and $``$, and describing the motion of the charge carriers, is valid. (ii) The flux conservation law for the magnetomotive field $``$ in (14) holds true. (iii) We consider a nondissipative Hall state, which has $`\stackrel{~}{\sigma }_{xx}=\stackrel{~}{\sigma }_{yy}=0=\stackrel{~}{\rho }_{xx}=\stackrel{~}{\rho }_{yy}`$. (iv) A necessary condition for Hall quantization to hold is that the total collective canonical particle momentum is derivable from a collective quantum phase, $`𝒑\mathrm{}\theta `$, such that the Bohr-Sommerfeld type integral of this momentum is quantized into units of Planck’s action quantum. The translation of the conventional quantum Hall problem into a rotating frame elucidates and strengthens the point of view that it is essentially a basic Hamiltonian quantity in phase space, the closed action integral of the collective particle momentum $`𝒑`$, which yields Hall state quantization.

An experimental proof or disproof of the quantization rule (16) in a rotating system will show if our assertion about the nature of the quantum Hall phenomenon is true and that it indeed relies upon the existence of collective particle momentum and the motion of the associated generalized flux quanta.

We thank Grisha Volovik for a discussion of the contents of this work. U. R. Fischer acknowledges financial support by the DFG (FI 690/1-1).


text:
# A Modified Parallel Tree Code for N-body Simulation of the Large Scale Structure of the Universe

## 1 INTRODUCTION

N-body codes are one of the most important tools of theoretical cosmology  because they offer the possibility of simulating most of the gravitational processes driving the formation of the large scale structure of the universe (hereafter LSS) . These simulations are often used to check cosmological models and to constrain the free parameters of these models which cannot be fixed either theoretically or observationally.
The typical mass scale for gravitational instability, the Jeans mass, has a value of $`10^{6.5}`$ solar masses (1 solar mass $`1.9\times 10^{33}g`$) at the recombination epoch, and it gives the approximate size of the first objects forming by gravitational collapse at that epoch. On the other hand, the largest structure we see in our Universe today, the “Supercluster” of galaxies, has a mass in excess of $`10^{18}`$ solar masses. Moreover, the gravitational force has a truly long-range character, which makes it impossible to introduce reasonable upper cutoffs in the mass range. For all these reasons, one would like to be able to perform simulations spanning more than 12 orders of magnitude in mass, but present-day state-of-the-art software and hardware technology does not allow simulations with more than $`10^9`$ bodies. For these reasons, the quest for increasingly efficient algorithms is still in progress. However, the importance of making N-body simulations is clear to several authors . During the past years N-body codes have been much improved and applied successfully to various problems in galaxy dynamics, galaxy formation, and cosmological large structure formation. Nevertheless, the computational expense has remained prohibitive for $`N>10^9`$, even using tree-based algorithms on the most powerful computers.
The situation is even worse for other N-body algorithms. The N-body direct evolution method scales as $`O(N^2)`$, which makes it impossible to run simulations with more than $`10^4`$ particles. To overcome this difficulty, and when high accuracy is required, alternative numerical methods based on hierarchical force-computation algorithms are widely used. The recent effort has addressed the production of new software and algorithms for the new generation of high-performance computer systems. The ultimate target is an implementation of the tree N-body algorithm to run simulations with higher accuracy and particle number, decreasing the cost of the simulation in terms of CPU time and increasing performance in terms of number of particles/second elaborated when running on MPP systems.
Among the tree algorithms designed to compute the gravitational force in N-body systems, one of the most used and powerful in modern cosmology is that by Barnes and Hut (BH) . The BH octal-tree recursive method is inherently adaptive and allows one to achieve a higher mass resolution even if parallel implementation of this algorithm  suffers from a serious drawback: it can easily run into imbalance as soon as the configuration evolves, causing performance degradation. In this paper we present a modified version of the BH algorithm in which we have introduced an enhanced grouping strategy. We will show how this feature allows an increase in performance when we consider N-body simulation with a large number of particles ($`N10^6`$). The code we present incorporates fully periodic boundary conditions using the Ewald method, without the use of fast Fourier transform techniques .
In Section 2 we give a brief description of our N-body parallel code, based on the BH tree algorithm, and the dynamic load balance (DLB) policy adopted. In Section 3 we describe our enhanced grouping strategy. In Section 4 we show the results of our tests and in Section 5 we report our conclusions.

## 2 THE PARALLEL CODE

Since the publication of the monograph by Hockney and Eastwood in 1981 , a new class of particle simulation methods  has emerged as an alternative to particle-particle (PP)  particle-mesh (PM) (for a review of this method see ) and particle-particle-particle-mesh ($`P^3M`$)  methods. These new methods are characterized by the particles being arranged into a hierarchy of clusters, which span the full range of length scales from the minimum interparticle spacing up to the diameter of the entire system. These methods are usually known as tree methods or tree codes because of the data structures used. With these new methods, the short-range force on a particle $`p`$ is calculated as a direct sum over nearby particles. Remote bodies are organized into groups which become progressively larger with the distance from the particle; then a multipole expansion of the potential of each cluster about its center of mass is performed. The long-range contribution to the acceleration is given by the sum of the particle-cluster interactions.

### 2.1 The Barnes-Hut Tree Algorithm

The BH algorithm works using a hierarchy of cubes arranged in an octal-tree structure; that is, each node in the tree has eight siblings and each node represents a physical volume of the space. The total mass of all particles within a given volume and their centers of mass are stored at the corresponding node. Thus, the system is first surrounded by a single cell (cube) encompassing all the particles. This main cell (called root) is subdivided into eight subcells of equal volume, each containing its own subset of particles. Each subcell in turn is subdivided into eight new subcells and so on. This procedure is repeated until each cell at the lowest level contains only one particle.
. The force on any given $`p`$ is then the sum of the forces by the nearby particles plus the force by the distant cells whose mass distributions are approximated by multipole series truncated typically at the quadrupole order . The criterion for determining whether a cell is sufficiently distant for a multipole force evaluation (that is, for approximating the cell as a multipole) is based on an opening angle parameter $`\theta `$ given by

$$\frac{C_l}{d}\theta ,$$
(1)
where $`C_l`$ is the size of the cell and $`d`$ is the distance of $`p`$ from the center of mass of the cell. Smaller values of $`\theta `$ lead to more cell opening and more accurate forces (for $`\theta =1`$ we have an error lower than 1% on the accelerations ). The equations of the dynamics are solved using the Leapfrog integrator.

### 2.2 Data Distribution and DLB

In our parallel implementation of the BH tree algorithm, using the PGHPF/CRAFT (an implementation of High Performance Fortran by the Portland Group) programming environment for the Cray T3E system, we have exploited both the Data Sharing and the Work Sharing programming models. The flexibility of the PGHPF/CRAFT environment allows one to mix these two modes in order to gain the maximum efficiency and speed-up. We can distinguish two main phases in our code structure: the Tree\_formation (TF) and the Force\_compute (FC). A data distribution in contiguous blocks

!HPF$ DISTRIBUTE PARTICLE\_ATTRIBUTE(BLOCK,\*)

and alternatively, a fine grain distribution

!HPF$ DISTRIBUTE TREE\_ATTRIBUTE(CYCLIC,\*)

were adopted to distribute the particle data properties and the tree data properties.
The !HPF$ DISTRIBUTE directive of the PGHPF/CRAFT compiler allows us to consider an array like PARTICLE\_ATTRIBUTE (or TREE\_ATTRIBUTE) as a unique large array, accessible from all the processors, the array being physically distributed in the local memory of all the processors. We used two different sets of initial conditions, namely uniform and clustered distributions having 2 million particles each, and they were carried out using from 16 to 128 PEs. Our results show that the higher code performances are obtained using a fine grain tree data distribution and a coarse grain bodies data distribution. A detailed description can be found in . The static array distribution, fixed as described above, allows each PE to cooperate during the TF phase by using principally the DO INDEPENDENT structure of PGHPF that is a synchronous mechanism, and then to execute the FC phase in asynchronous mode. To minimize the communication overhead, each PE executes the FC phase mainly on the local residing bodies. The BLOCK distribution arranges bodies with the nearest logical number (near in the space) in the same PE local memory, or in the nearest PEs. Using the above mentioned data distribution, each PE has a block of closed bodies in the local memory ($`Np=N_{bod}/N\$PEs`$, where N$PEs is the number of processors used for the simulation);in an initial condition with a uniform distribution, the PEs having extreme numeration in the pool of available PEs have a lower load at each time-step. The load imbalance is enhanced when a clustered situation occurs during the system evolution. The PEs having bodies in clustered regions have a greater workload since the load of the FC phase increases as the mass density grows. The technique we follow to perform a load redistribution among the PEs is to assign each PE to execute this phase only for a fixed portion of the bodies residing in the local memory $`NB_{lp}`$ given by

$$NB_{lp}=(N_{bod}/N\$PEs)P_{lp},$$
(2)
where $`P_{lp}=`$ const. ($`0P_{lp}1`$).

The FC phase for all the remaining bodies

$$N_f=N\$PEs(N_{bod}/N\$PEs)(1P_{lp})$$
(3)
is executed by all the PEs that have concluded the FC phase for the assigned $`NB_{lp}`$ bodies. No correlation is considered between the PE memory location of the body belonging to the $`N_f`$ set and the PE that computes the FC phase on it. The results imply that it is possible to fix a $`P_{lp}`$ value that allows the best code performances. Data already presented in  show that it is convenient to fix the $`P_{lp}`$ value near $`0.25`$, which is the value that maximizes the load balance for N-body simulations of the LSS both in uniform and clustered situation.

## 3 GROUPING

Our work and data sharing-parallel tree (WDSH-PT) code is principally aimed at running LSS cosmological simulations with a number of particles as high as possible using supercomputers such as Cray T3E systems. In order to increase the code efficiency, we adopt initially the grouping method proposed by Barnes  and introduce a modified implementation of his grouping policy yielding very high gains in the code performances with the same accuracy.

### 3.1 B90 Grouping

To compute the force on a body, the BH algorithm needs to build an interaction list ($`IL`$) for each particle $`p`$. Starting from the root cell, a tree inspection is done and the opening angle parameter $`\theta `$ is used to evaluate whether a cell must be opened or closed as mentioned above. If a cell has dimension $`C_l`$ and distance $`d`$ from the particle $`p`$ so that Eq. (1) is verified, the cell is closed, it is added to the $`IL`$, and its subcells are not investigated further. Otherwise the cell is opened and its subcells are investigated in the same way. Bodies belonging to an opened cell are added to the $`IL`$.
Next, the force on the body is computed using the monopole and quadrupole momenta for all the cells in the list.

##### BH timing

The tree inspection phase represents a sizeable task to compute the force because the cell opening criterion is applied many times for each particle. The CPU time $`T_o`$ to compute the force in a time-step for all the $`N`$ particles is

$$T_o=NT_l+NT_f,$$
(4)
where $`T_l`$ is the average time to build an $`IL`$ and $`T_f`$ is the average time to compute the force on each particle using the interaction list.

##### B90 timing

The basic idea of B90 was to build a unique interaction list that allows the force for a group of particles inside a region; i.e., a cell $`C_{group}`$ of the tree (grouping cell), to be computed reducing the number of tree inspections to build the $`ILs`$. B90 builds an IL that applies everywhere within $`C_{group}`$ and reuses this IL for each particle $`pC_{group}`$ in turn. In this way it is possible to reduce the tree inspection phase. The CPU time $`T_g`$ for B90 may be written as

$$T_g=N_{gc}T_{gl}+NT_{gf},$$
(5)
where

* $`N_{gc}`$ is the number of grouping cells (assuming that each body is inside a group region);
* $`T_{gl}`$ is the average time to build an interaction list for a group;
* $`T_{gf}`$ is the average time to compute the force on a particle using the list formed for the group.

In the following paragraphs we will compare the $`T_g`$ time with the $`T_o`$ time considering the generic case $`\theta =0.8`$ . We notice that different values of $`\theta `$ give similar results, as shown by the accompanying figures.

##### B90 opening criterion

The original BH algorithm adopts an opening criterion $`\theta _{BH}`$, based on the distance between the position of the $`p`$ particle and the center of mass of the remote cells, the $`IL`$ length ($`L_{IL}`$) being proportional to $`(\theta _{BH}^3)^1logN`$. In order to have the same accuracy as the original algorithm, the interaction list of the grouping cell is formed using Eq. (1) but now the $`d`$ term is computed in terms of the distance from the center of mass of an inspected cell and the edge of the grouping cell, as shown in Fig. 1 ($`d_{B90}`$ is used instead of $`d_{BH}`$).

This implies that the $`IL`$ formed using the grouping cell contains more elements than the $`IL`$ formed by applying the original BH algorithm.

##### B90 interaction list increment

The B90 adopts an opening criterion ($`\theta _{B90}`$) based on the distance between the edge of $`C_{group}`$ and the center of mass of the remote cells. In this case the $`IL_g`$ length will be proportional to $`\theta _{B90}^{3}{}_{}{}^{1}logN`$. Moreover, the B90 criterion uses $`\theta _{B90}`$ numerically equal to $`\theta _{BH}`$ when using the original BH algorithm. We found experimentally the relation between $`\theta _{B90}`$ and $`\theta _{BH}`$, using 2 million particles in a uniform distribution (see Fig. 2): this relation agrees with data in Salmon . A typical value used as opening criterion to run simulations for the LSS is $`\theta _{BH}=0.8`$. Consequently we consider $`\theta _{B90}=0.8`$ , which, in terms of $`IL`$ length increment, corresponds to running a simulation with $`\theta _{BH}=0.6`$.

##### B90 timing vs BH timing

Figs. 3 and 4 show the relationship $`T_{gl}T_l`$ and $`T_{gf}T_f`$ with $`\theta `$ ranging between $`0.4`$ and $`1.2`$. Considering $`T_{gl}=1.3T_l`$ and $`T_{gf}=2.2T_f`$ with $`\theta =0.8`$, Eq. (5) may be rewritten using the above relations as follows:

$$T_g=1.3\frac{N}{N_{gp}}T_l+2.2NT_f,$$
(6)
where $`\frac{N}{N_{gp}}=N_{gc}`$, and $`N_{gp}`$ is the average number of particles in a grouping cell.
For a large number of systems and in particular for our WDSH-PT code running on the T3E system, $`T_f`$ ranges from $`1.2T_l`$ to $`1.5T_l`$ at $`\theta =0.8`$. Considering these figures, we obtain from Eqs. (4) and (6), respectively,

$$T_o=2.2NT_l$$
(7)
and

$$T_g=2.2NT_l(\frac{1.3}{2.2N_{gp}}+1.2)$$
(8)
and then $`T_g>T_o`$.

A real gain could be obtained using B90 if the CPU time spent to form the interaction list were longer than the phase to compute the force on the particle. The results reported in  demonstrate that if $`T_lT_f`$ the code performance is between two and three times faster of the BH algorithm, and a good choice for the $`N_{crit}`$ value is about 32 .

### 3.2 The Modified Grouping Strategy for LSS Simulations

We will now describe the modification we introduce in the 1999 version of our WDSHPT code (WD99) to increase performance even if $`T_lT_f`$. The basic idea is to assign the same $`IL_g`$ to each particle within a cell $`C_{group}`$, containing a maximum of $`N_{crit}`$ particles.
We will not use the B90 criterion to build the interaction list. Instead, we will use the same $`\theta _{BH}`$ criterion used in the original BH algorithm. This criterion is applied to a hypothetical particle placed in the center of mass of the $`C_{group}`$, hereafter VB (Virtual Body) (Fig. 5). Moreover, we consider the $`IL_g`$ as formed by two parts given by

$$IL_g=IL_{far}+IL_{near}$$
(9)
$`IL_{far}`$ and $`IL_{near}`$ being two subsets of the interaction list. An element is included in one of the two subsets, using the following Sphere criterion for all the elements that satisfy Eq. (1).

Define $`Sphere_{radius}=3\frac{Cellsize(C_{group})\sqrt{3}}{2}`$
If $`Distance(IL_g(element),VB)>Sphere_{radius}`$
Add element to $`IL_{far}`$
Else
Add element to $`IL_{near}`$
Endif

Moreover all $`pC_{group}`$ are included in $`IL_{near}`$.

Using the two subsets it is possible to compute the force $`𝐅_p`$ on a particle $`pC_{group}`$ as the sum of two components,

$$𝐅_p=𝐅_{far}+𝐅_{near},$$
(10)
where $`𝐅_{far}`$ is a force component due to the elements listed in $`IL_{far}`$ and $`𝐅_{near}`$ is the force component due to the elements in $`IL_{near}`$. We assume the component $`𝐅_{far}`$ to be the same for each particle $`pC_{group}`$ and compute it considering the gravitational interaction between the VB and only the elements listed in $`IL_{far}`$, while the $`𝐅_{near}`$ component is computed separately for each $`p`$ particle by the direct interaction with the elements listed in $`IL_{near}`$.

Moreover, $`𝐅_{near}`$ contains a restricted number of elements in comparison with the $`𝐅_{far}`$ list, so we expect a net gain in performance even if $`T_lT_f`$. The gain that is possible depends on several parameters ($`N_{crit}`$, the size of the $`C_{group}`$ and the Sphere radius), whose ranges of variation are constrained by the maximum allowed value of the overall error of the method, as we will describe in the following sections.

### 3.3 Errors Analysis and Performance considerations

Before showing the performance of our WD99 procedure in an N-body simulation of the large scale structure of the universe, it is important that we perform an error analysis of the procedure itself. Considering that the cumulative error, when simulations for the LSS studies are run using the original BH algorithm, is lower than $`1\%`$ , fixing the opening criterion $`\theta =1`$, we will give some constraint concerning the size of $`C_{group}`$ , the $`N_{crit}`$ value, and the Sphere radius needed to have negligible cumulative error. The following sub-sections discuss the two main sources of error.

#### 3.3.1 The Differences in the Interaction List

The first error source is that WD99 uses the Sphere criterion and the VB to create an interaction list $`IL_gIL_{VB}`$ and WD99 applies the $`IL_{VB}`$ to all bodies $`pC_{group}`$. This approximation could introduce an error in the force value on the $`p`$ particle if the $`IL_p`$, created using the original BH algorithm, and the $`IL_{VB}`$ have a difference in the elements greater then $`1\%`$. As we found with our tests, in order to decrease this difference it is necessary to limit the size of $`C_{group}`$ that is equivalent to fixing a critical level of the tree structure: cells above the critical level cannot form a grouping cell. The user has to fix the critical level considering the density of local bodies in the box where bodies are arranged: the critical level must be chosen in such a way as to make the difference between $`IL_p`$ and $`IL_{VB}`$ negligible (no more than $`1\%`$ of the elements). It seems reasonable for a LSS simulation in a 50 Mpc box with more than 2 million particles to fix the critical level as the sixth level of the tree. The next section shows in detail the obtained results. In any case, the cumulative error is very small considering the increase in accuracy compared with the original BH algorithm, due to the inclusion of all $`pC_{group}`$ in the interaction list.

#### 3.3.2 Approximation of the Force Component

The second error source is due to the assignment of $`𝐅_{far}`$, computed for the VB, to each $`pC_{group}`$. We found that the Sphere criterion allows us to reduce this error to values much lower than 0.01% for $`N10^6`$ if the dimension of the $`C_{group}`$ cell is fixed with the critical level as mentioned above, and the Sphere radius is three times the radius of the sphere enclosing the $`C_{group}`$ cell (Fig. 5).
Another important constraint to be fixed is the value of $`N_{crit}`$. All the elements $`pC_{group}`$ are listed in the $`IL_{near}`$ list and there is a direct body-body interaction among the $`N_{gp}`$ ($`N_{gp}N_{crit}`$) elements forming the group. This introduces a term $`O(N_{gp}N_p)`$ in the algorithm complexity. In order to avoid a decrease of the code efficiency and to maintain a good code accuracy, as with the original BH algorithm, it seems reasonable, running LSS simulation with more than 1 million particles, to maintain $`N_{crit}32`$. We adopted, in our runs, a safe value $`N_{crit}=16`$ even if we obtained good results with $`N_{crit}=32`$.
In the next section we show the errors obtained using the above-mentioned constraints when applying WD99 to LSS simulation with both uniform and clustered distributions.

## 4 TESTS AND RESULTS

We carried out many tests to estimate the error introduced in the WD99 and obtained increased performances using several values of $`N_{crit}`$. Therefore, this section is subdivided as follows: first we test whether our algorithm increases the average length of the interaction list, then we measure the resulting percentage error, and we conclude with an overall performance analysis. As a test case we ran a simulation using 2 million particles for LSS in a cubic region of 50 Mpc, starting from a homogeneous initial condition (redshift $`Z=50`$) and reaching a clustered configuration (redshift $`Z=0`$). We used an opening angle parameter $`\theta `$ ranging from $`0.8`$ to $`1.2`$. Our tests were executed on a Cray T3E system and the results will be shown in the following sections.

### 4.1 Measuring the Interaction List Length

The aim of this first test is to verify that the WD99 algorithm does not introduce a significant computational cost when the force for a generic particle is computed. This measurement is substantially performed on the average length of the $`IL`$ we form adopting our code. Fig. 6 reports the result we obtain when the simulation evolves at redshift $`Z=50`$.

Tests were executed for several values of redshift, but the differences between BH and our algorithm was computed only at the end of the run. The curves were obtained by fixing $`N_{crit}=32`$ and varying the critical level from 5 to 8. In all cases the differences we obtained are negligible, which means that the computed $`IL`$ for the VB (with the adopted Sphere criterion) is about equal to the $`IL`$ we obtain for a generic particle with the original BH algorithm. The first important result is that WD99 does not produce any increment in the IL length and consequently $`T_l=T_{gl}`$.

### 4.2 Error Measurement

We carry out this measurement in two phases. First we run a single time-step of the 2-million-particle simulation at redshift $`Z=50`$. We compare the values we obtain running the BH original algorithm and the WD99. As a reference case, we adopt the critical level equal to 6.

A similar comparison is made at $`Z=0`$ and the BH and WD99 histograms of the forces of each component are compared. The comparison shows a negligible difference in the force distribution in a single time-step, at least an order of magnitude less than the error of the original BH algorithm.
The second measurement is made analysing an entire system evolution. We start with the initial condition of 2 million particles with redshift $`Z=50`$, $`\mathrm{\Delta }t=0.001`$, $`\theta =0.8`$, and particle mass about $`1.65510^{10}`$ solar masses. The system evolution is carried out up to redshift $`Z=0`$. The evolution is executed with the original BH algorithm and with the WD99 code. As reference case, we adopt the critical level equal to 6. We measure the absolute error $`ϵ`$ in the position of particles and in the velocities of particles in the mean square sense,

$$ϵ_{pos}=\frac{1}{N}\sqrt{\underset{i=1}{\overset{N}{}}(X_{BH}X_{WD99})^2+(Y_{BH}Y_{WD99})^2+(Z_{BH}Z_{WD99})^2}$$
(11)
and

$$ϵ_{vel}=\frac{1}{N}\sqrt{\underset{i=1}{\overset{N}{}}(Vx_{BH}Vx_{WD99})^2+(Vy_{BH}Vy_{WD99})^2+(Vz_{BH}Vz_{WD99})^2}.$$
(12)
The study of the final evolution is described in the next sub-section. Here we give the measured value at the end of the simulation: $`ϵ_{pos}=0.003`$ and $`ϵ_{vel}=0.01`$.

Similar values are measured when running simulations with more than 2 million particles and with a critical level equal to 6. The obtained $`ϵ`$ values lead us to conclude that the WD99 procedure does not introduce significant errors in comparison with the BH algorithm.

#### 4.2.1 Simulation Analysis

The final stages obtained running simulations with the BH algorithm and the WD99 code are very similar. The two-point correlation function is defined as

$$\xi (r)=\frac{N_q}{\rho N_cV}1,$$
(13)
$`N_q`$ being the number of pairs of particles with separations between $`r`$ and $`r+\mathrm{\Delta }r`$, $`V`$ the volume considered, $`N_c`$ the particle number taken as centres, and $`\rho `$ the mean particle density. We calculate this function at redshift $`Z=0`$ for WD99 and BH algorithms. The values we obtain are perfectly equal, and the substructures we form (number and size) are identical.

### 4.3 WD99 Performances

To conclude our WD99 description, we report the performances measured for the WD99 code (Fig. 7) (including the boundary periodic conditions using the Ewald method ) and the performances of the original BH algorithm.

The measured performances lead us to the conclusion that when the system evolution is clustered ($`Z=0`$) the WD99 does not decrease the performance as the BH algorithm. This effect is due to the nature of the WD99 algorithm, which has a structure that increases the efficiency when clusters of particles are well closed. This important effect allows us to run simulations with very clustered systems, obtaining very good performance and negligible errors. Moreover, the efficiency of the WD99 increases by a factor of up to five at the redshift $`Z=0`$. The gain is enhanced when bigger simulations are run: a recent simulation with 16 million particles performed on the Cray T3E system using WD99 showed an increase in performance by a factor of 7 at the redshift $`Z=0`$. We note that the gain obtained, in comparison with that obtained by the original BH algorithm, is greater using a lower critical level (5 or 6). The gain is incremented using a Sphere criterion with Sphere radius lower than the value we consider, having only a small increment in the global error.

## 5 CONCLUSIONS AND FUTURE

The code WD99 is mainly used for LSS studies, but it could be tested and used for other applications where accuracy not higher than 1% is necessary. Considering the high performances we obtained, the WD99 method may be very successfully applied when clustered configurations such as galaxies or clusters of galaxies have to be studied. The new approach could be applied also to other fields of physics where collisionless systems are to be simulated, as in plasma and hydrodynamic studies.
The code is written in Fortran 90 with PGHPF/CRAFT, but the latest version (written in F90 and C languages) uses the one-side communication library SHMEM, allowing it to run on the ORIGIN 2000 systems. A new version will be implemented using dynamical array allocation, and we are studying the implementation of the parallel out-of-core , moving data in the disk. This version will be developed for a CC-NUMA machine with MPI-2. We plan to have a freely available version of WD99 in October 2000.

## 6 ACKNOWLEDGMENTS

All the tests were carried out using the CRAY T3E 1200/256 machine at the CINECA (Casalecchio di Reno (BO), Italy), a 256 PE system, using the financial support of the Italian Consortium CNAA (Consorzio Nazionale per l’Astronomia e l’Astrofisica). We thank Dr. G. Erbacci of CINECA and Dr. A. F. Lanza of Catania Astrophysical Observatory for their useful help.


text:
# Relaxation properties of small-world networks

## 1 Introduction

Networks of the real world often seem to combine aspects from regular and from completely random lattices. Social networks, neural networks, electrical powergrids, and traffic networks  are all examples of patterns not described satisfactorly by conventional regular lattices, nor by completely random lattices. Social structures, for instance, do not behave as regular lattices, since (as is well known) randomly chosen people are connected in general by a small number of intermediary bilateral ties. Here, as in random graphs, the minimal (chemical) distance between any two points in the system scales logarithmically with the system size .

To combine these two properties, Watts and Strogatz recently introduced the idea of small-world networks . This construction is a superposition of a regular lattice with a random lattice, and includes simultaneously well defined local clusters and short global connections. As we will demonstrate, these systems also display properties intermediate between those of regular and tree-like (loop-less) lattices, already under a small number of global connections, provided the system size is large enough.

Much work has already been done on the properties of small-world networks  but most of it has focused on static (geometric) properties. We shall not address these issues, but rather concentrate on a dynamical model defined on the structure. Treatments of the dynamics of small-world networks include for instance the study of an Ising model defined on the lattice , spectral properties of the small-world Laplacian, , percolation , spreading of diseases  and neural networks . In the following we will examine the properties of random walks on small-world networks, in particular the relaxation, exemplified by the probability for a random walker of being at the original site at a later time. This is a simple quantity to extract numerically, and very relevant for various physical properties: It is sensitive to the topology of the network, and is related to its vibrational modes.

## 2 Definition of the model and presentation of the results.

The small-world networks we consider are built as follows: We start from a regular lattice with $`L`$ vertices in $`1`$ dimension under periodic boundary conditions, each site being connected symmetrically to its $`2k`$ nearest neighbours, i.e. having as coordination number $`z=2k`$. Then we add to each of the sites with probability $`p`$ a new bond. The other end gets attached with equal probability to any of the lattice sites; this allows also the possibility of vertices to become connected to themselves. In this way we add, independent of $`k`$, on the average $`pL`$ new bonds to the underlying regular lattice.

This construction follows  for $`k=1`$ and is simpler than the original procedure , by which one rewires with probability $`p`$ each of the original $`kL`$ bonds randomly.

A step-wise diffusion process is now defined by specifying all the transition probabilities $`W_{i,j}`$ entering the master-equation:

$$P(i,n+1)P(i,n)=\underset{j}{}W_{i,j}P(j,n)P(i,n)\underset{j}{}W_{j,i}$$
(1)
The $`W_{f,i}`$ is the probability to go from site $`i`$ to site $`f`$ during one time step, and the probability $`P(i,n)`$, $`i=1\mathrm{}L`$ is just the probability of being at site $`i`$ after the $`n`$th step. The process defined in Eq. (1) is the discrete variant of diffusion on an arbitrary lattice, a topic interesting in its own right. Diffusion on regular lattices is ubiquitous, and diffusion on random graphs has (among other things) also been studied in the context of glassy relaxation . We are therefore inspired to investigate what happens on the small-world model, which interpolates between these two extremes. Previously a lot of interest has also been seen in the related problem of diffusion on fractals (see for example  and references therein). As we proceed to show, diffusion on Cayley trees  shows also features closely related to the present problem. Furthermore, the motion of charge carriers or of excitons over polymer chains, where steps between spatially close sites can connect regions far apart along the chemical backbone, also involves global shortcuts .

The transition probabilities $`W_{i,j}`$ in Eq. (1) are as follows: First $`W_{i,j}=0`$ if there are no bonds between $`i`$ and $`j`$. For $`i`$ connected to $`j`$ by one or more direct bonds, $`W_{i,j}`$ is proportional to the number of such bonds. The same holds for the probability of remaining at the same site after one time unit, i.e. we allow “sticking”. Formally

$$W_{i,j}=\frac{z_{i,j}+\delta _{i,j}}{z_j+1}$$
(2)
In this equation, $`z_{i,j}`$ is the number of bonds between the two sites $`i`$ and $`j`$, and $`z_i`$ is the total number of bonds emanating from vertex $`i`$, i.e. the coordination number of the site. Hence $`z_i=_jz_{i,j}`$. Note that the $`z_{i,j}`$-values are determined both by the additional wiring as well as by the underlying lattice. The $`\delta _{i,j}`$ and the $`1`$ in the denominator appear because we allow for the possibility of the walker to remain at site $`i`$ during a time step. This procedure renders the numerically determined $`P(i,n)`$ smoother in $`n`$. We remark that the rates defined according to Eq. (2) are not symmetrical in $`i`$ and $`j`$, i.e. in general $`W_{i,j}W_{j,i}`$.

The algorithm we have used is the exact (cellular automaton) enumeration of random walks , corresponding to the implementation of Eq. (1). All the results plotted are averaged over $`200`$ disorder configurations. We have worked mostly with the value $`k=1`$. This is also the value implied if we do not state otherwise.

We focus on the probability $`P(i,n|i,0)`$ that a particle initially at site $`i`$ is found at the same site just after the $`n`$th step. In the figures below we plot $`P(i,n|i,0)`$, i.e. the average of $`P(i,n|i,0)`$ over the different realisations of the small-world lattice. Since all sites are equivalent in an ensemble of small-world networks, this quantity does not depend the particular site $`i`$ chosen, and we hereafter denote it by $`P_n(0)`$. In Fig. 1 we have chosen $`p=0.05`$ and plotted $`P_n(0)`$ in double-logarithmic scales for system sizes ranging from $`L=1000`$ to $`L=10000`$. This allows us to examine the dependence of $`P_n(0)`$ on the size of the system.

From Fig. 1 we infer that initially all the curves fall on one curve and that for large $`n`$ they saturate at their respective equilibrium values, $`1/L`$. However, Eq. (2) implies an inhomogenous equilibrium distribution

$$P^{eq}(i)(z_i+1).$$
(3)
Therefore $`P(i,\mathrm{}|i,0)`$ depends on the specific small-world realisation, and will fluctuate from realisation to realisation around its average value $`1/L`$.

To find out how much of the behavior is due to finite size effects, we subtract from each average curve in Fig. 1 its corresponding average equilibrium value $`P_{\mathrm{}}(0)1/L`$, and replot $`P_n(0)P_{\mathrm{}}(0)`$ in Fig. 2.

From Fig. 2 we see that all curves collapse nicely onto what we view as representing $`P_n(0)`$ on small-world networks in the limit $`L\mathrm{}`$. Both Figs. 1 and 2 display initially a quasi-linear decay in the chosen double logarithmic scales, and this may be viewed as being an approximate power law decay. Depending on $`p`$, the exponents range from around $`0.5`$ for the smallest $`p`$ to around $`0.6`$ for the largest. This regime is followed by a steeper decay at larger $`n`$. To highlight the power law character we have plotted in Fig.3 $`P_n(0)`$ for $`p=0.01`$. As is evident from the figure, the power-law domain extends well over two orders of magnitude in $`n`$.

The results can be understood qualitatively in the following way: For a fractal one has

$$P_n(0)n^{d_s/2},$$
(4)
where $`d_s`$ is the spectral dimension. Thus the initial decay in Figs. 1 to 3 follows that of a fractal with a $`d_s`$ close to $`1`$, i.e. that of a quasi $`1`$-d system. This is reasonable given our construction: for sufficiently small $`p`$ and small $`n`$, only relatively few random walkers encounter any long-range connections (shortcuts). Therefore in the beginning the behavior of $`P_n(0)`$ closely reflects the character of the underlying $`1`$-d lattice. However for larger $`n`$, the random walkers probe larger and larger portions of the graph, and thus follow more and more shortcuts. This speeds up progressively the decay of $`P_n(0)`$ as more regions at larger and larger length scales are visited, and the fractal picture is lost. One would thus expect that the concept of a $`d_s`$ begins to be invalid when the random walkers visit enough short cuts, i.e. when the $`1`$-d diffusion extends longer than the typical distance between shortcuts. This is the fundamental length scale $`\xi `$ of small-world networks, besides the lattice constant, which is less important here. In our case we have

$$\xi =p^1,$$
(5)
measuring $`\xi `$ in units of the lattice constant. For diffusion on scales smaller than $`\xi `$ one furthermore has in terms of the diffusion constant $`D`$ of the regular lattice $`\xi ^22Dn`$, so that

$$n\frac{1}{2Dp^2}.$$
(6)
Given that we allow random walkers to stay at a site during a time step, $`D=1/3`$ and thus $`n=2/3p^2`$. However some walkers do encounter shortcuts at length-scales below $`p^1`$, and numerically the cross-over to a region that does not have approximate power-law character is seen to take place earlier than $`np^2`$.

We turn now to the analysis of this region. To be able to follow the it more closely, we replot the results of Fig. 2 for $`L=10000`$ on semi-logarithmic scales in Fig. 4. Evidently, the decay for larger $`n`$ is slower than exponential.

The decay of $`P_n(0)`$ is hence quicker than a power-law, but slower than the one for Cayley-trees, for which one has (for coordination numbers greater than $`2`$)

$$P_n(0)n^{3/2}\mathrm{exp}\left[Cn\right],$$
(7)
where $`C`$ is a constant. Comparing this behavior to the one displayed in Figs. 1 and 4 we remark that in our case, for a relatively small number of steps, the decay goes approximately as $`n^\alpha `$ with $`\alpha 1/2`$, whereas at larger $`n`$, a more adequate description would be a stretched exponential, $`\mathrm{exp}[Cn^\beta ]`$. One may even suspect that the decay in Fig. 4 obeys $`P_n(0)n^\alpha \mathrm{exp}[Cn^\beta ]`$. A fit of the data in Fig. 4 to this functional form (keeping $`\alpha 0.5`$ fixed) is also shown; the fit turns out to be indistinguishable from the numerical data, when we choose $`\beta =0.56`$ and $`C=0.04`$, as plotted in Fig. 4.

We now consider the dependence of the decay on the value of $`p`$. For this we plot in Fig. 5 the decay law $`P_n(0)`$ for $`L=2000`$ and $`p`$ ranging from $`p=0.01`$ to $`p=0.8`$.

We note that the initial power-law like region diminishes with increasing $`p`$. Furthermore, the plateau-region $`P_n(0)1/L`$ is reached earlier for larger $`p`$. This is in accordance with our argument above, that the long-range connections (short-cuts) interrupt the simple diffusion on the underlying lattice, such that the crossover length decreases with increasing $`p`$ (c.f. also Eq. (6)). As $`p`$ becomes large enough, the power law regime practically disappears. This is so because the random walker rapidly meets a short-cut. As before, the influence of finite size effects can be reduced by plotting, as in Fig. 2, $`P_n(0)P_{\mathrm{}}(0)`$.

We have also performed simulations of the random walk on small-world networks where the underlying lattice has a $`k`$ value larger than $`1`$. In Fig. 6 we plot the results for $`p=0.1`$ and $`L=2000`$ in the cases of $`k=1`$, $`k=2`$, $`k=3`$ and $`k=4`$. The findings reproduce the general picture: $`P_n(0)`$ behaves like a power law for small $`n`$, while decaying more rapidly as $`n`$ gets larger. The curves for different $`k`$ are mainly shifted with respect to each other, and the network with the largest coordination number (largest $`k`$) also displays the quickest relaxation. To be noted, however, is that the case $`k=1`$ has the largest dynamical range and thus shows best the decay forms, while also being the one simplest to implement; hence $`k=1`$ may be the ideal small-world model.

## 3 Conclusions

In this work we have studied numerically the behavior of random walks on small-world lattices. Our work has focused on the probability of being at the initial site $`P_n(0)`$ as a function of the number of steps $`n`$. This quantity is found to show a complex, very interesting pattern: initially $`P_n(0)`$ displays a power-law, “quasi-fractal” regime. At larger $`n`$ a quicker decay takes over, reminiscent of stretched exponentials. In this respect the $`P_n(0)`$ decay is intermediate between the decays found for fractal structures and the ones found for tree-like (loop-less) structures, exemplified here by Cayley trees.