text:
# Dark Matter in Dwarf Galaxies: The First Dark Galaxy?

## 1. Introduction

The first hint of unusual goings-on in the neighborhood of LGS 3 came from Christian & Tully (1983), who mentioned that H i observations by Hulsbosch revealed that LGS 3 lies on the edge of a large cloud of gas (called HVC 127-41-330 by them) that contains a significant velocity gradient. However, with the 36 resolution of his Dwingeloo observations, Hulsbosch was apparently unable to draw any firm conclusions about the nature of this object. Over the next 17 years the situation remained murky as few observers paid attention to this part of the sky. The sole reference to the cloud during this period was in the HVC catalog of Hulsbosch & Wakker (1988), who assumed that it was part of LGS 3.

Three years ago the cloud was rediscovered by Blitz & Robishaw (2000) during their search for gas associated with Local Group dwarf spheroidals in the Leiden/Dwingeloo Survey (LDS) of Galactic Neutral Hydrogen (Hartmann & Burton 1997). Blitz & Robishaw (2000) pointed out the offsets in position and velocity between this cloud and the H i known to be directly associated with LGS 3 (Thuan & Martin 1979; Young & Lo 1997). They speculated that the cloud might have been ram-pressure stripped out of LGS 3 by hot gas in the halo of M31. They noted, however, that the velocity of the cloud should be less negative than the velocity of LGS 3 ($`v_{\text{HVC}}>287`$ km s-1), which it is not ($`v_{\text{HVC}}330`$ km s-1). Like Hulsbosch, their ability to consider other possible origins for the cloud was limited by the low angular resolution of the LDS data.

Robishaw, Simon, & Blitz (2002) shed new light on this object with a high-resolution, wide-field H i map of the region made at Arecibo. These observations showed that the cloud is completely separate from LGS 3, demonstrating that it is indeed an HVC.

## 2. Discussion

The Robishaw et al. (2002) observations raised several pressing questions. What is the nature of the LGS 3 HVC? Are LGS 3 and the HVC actually located at the same distance, or is their apparent proximity just a chance alignment? Even though there is no H i connection between LGS 3 and the HVC, could the high-velocity gas be associated with LGS 3 in any way?

### 2.1. Distance of the HVC

In Robishaw et al. (2002) we pointed out that the probability of a chance superposition between a compact HVC and a dwarf galaxy — if they are completely unrelated objects — is rather low. Using the statistics of HVCs in the LDS (de Heij et al. 2002) and known Local Group dwarfs, we showed that the likelihood of an HVC being located less than 30 from a dwarf is $`4\%`$ in this area of the sky. Given such an angular coincidence, the probability of the velocities also matching within 50 km s-1 is $`20\%`$, for a joint probability of less than $`1\%`$. This line of reasoning suggests, but does not prove, that the HVC and LGS 3 are physically associated.

A more powerful argument can be made based on the faint H i strips that our first Arecibo map revealed on either side of the HVC. These features appeared to be long, thin tails of gas that connect to the HVC. That appearance, along with their nearly symmetrical placement relative to the line between LGS 3 and the HVC, is a strong signal of a tidal interaction between the two objects. If this interpretation is correct, the HVC must be located physically close to LGS 3, at a distance of about 700 kpc from the Milky Way. This HVC would then be the first HVC known to lie more than $`50`$ kpc away.

### 2.2. HVC Rotation Curve

The HVC shows a surprisingly regular velocity gradient indicative of circular rotation. Robishaw et al. (2002) argued against the possibility of a chance alignment of two separate clouds moving at slightly different velocities, so if the gradient is *not* due to rotation, the only alternative is shear. A simple timing calculation shows that this explanation is unlikely. The cloud currently has a radius of about 3 kpc (D/1 Mpc), where D is the distance of the HVC. The velocity gradient is 15 km s-1, which would correspond to an expansion velocity of 7 km s-1. At that expansion rate, the time for the cloud to reach its current size is $`4\times 10^8`$ (D/1 Mpc) yr. In this scenario, distances of less than 100 kpc are highly implausible — the HVC would have to be extraordinarily young, and even then the nearest possible progenitors, the Magellanic Clouds, are too far away for the HVC to have reached its present position at any reasonable velocity. At 700 kpc the lifetime of the cloud is still uncomfortably shorter than a Hubble time, and the lack of a credible progenitor object is equally severe. We conclude that the velocity gradient is almost certainly due to rotation.

The rotation curve of the HVC, as measured by Robishaw et al. (2002), is symmetric and flattens out at a radius of about 10 (see Fig. 1). The total mass derived from the rotation curve is larger than the H i mass for any reasonable distance, implying that the HVC is dark matter-dominated. At 700 kpc, 82 % of the mass is in dark matter; at 50 kpc, 99 % of the mass is dark. The HVC can only be dynamically dominated by luminous material if its distance is $`>2`$ Mpc, which is extremely unlikely.

### 2.3. New Observations

In August 2002, we used another Arecibo observing run to make a new, more sensitive map of LGS 3 and the HVC. This map, displayed in Figure 2, confirms all of the features seen in our original data. The tidal tails are now visible in detail, and their connection to the main body of the HVC is evident. The tails are somewhat more extended than could be seen in the first map, and additional observations over a wider field have revealed that the southwest tail continues on beyond the edge of the map shown in Figure 2. The new data are fully consistent with the scenario proposed by Robishaw et al. (2002) of a tidal interaction between the HVC and LGS 3.

## 3. Conclusions

Given these findings, what can we conclude about the nature of this HVC? We have presented several arguments that it is probably located hundreds of kiloparsecs away, most likely at 700 kpc. The HVC appears to be undergoing a tidal interaction with LGS 3, which is stripping away a substantial portion of its neutral gas. If the HVC is in the Local Group it is dark matter-dominated, regardless of its exact distance. And, it does not contain any stars (Robishaw et al. 2002). This set of properties is exactly what is expected of the missing dark matter satellites that are predicted by simulations. We propose that this HVC is the first observed representative of this population of missing objects.

#### Acknowledgments.

This research was supported by NSF grant AST-9981308.

## References

Blitz, L., & Robishaw, T. 2000, ApJ, 541, 675

Christian, C. A., & Tully, R. B. 1983, AJ, 88, 934

de Heij, V., Braun, R., & Burton, W. B. 2002, A&A, 391, 159

Hartmann, D., & Burton, W. B. 1997, Atlas of Galactic Neutral Hydrogen (Cambridge: Cambridge Univ. Press)

Hulsbosch, A. N. M., & Wakker, B. P. 1988, A&AS, 75, 191

Robishaw, T., Simon, J. D., & Blitz, L. 2002, ApJ, 580, L129

Thuan, T. X., & Martin, G. E. 1979, ApJ, 232, L11

Young, L. M., & Lo, K. Y. 1997, ApJ, 490, 710


format: text

text:
# Logarithmic corrections to three-dimensional black holes and de Sitter spaces

## I Introduction

Recently there are several works which show that for a large class of AdS black holes, the Bekenstein-Hawking entropy receives logarithmic corrections due to thermodynamic fluctuations . The corrected formula takes the form

$$S=S_0\frac{1}{2}\mathrm{ln}\left(C_vT^2\right)+\mathrm{},$$
(1)
where $`C_v`$ is the specific heat of a given gravitational system at constant volume and $`S_0`$ denotes the Bekenstein-Hawking entropy. Here an important point is that for Eq.(1) to make sense, $`C_v`$ should be positive. However, the d-dimensional (dD) Schwarzschild black hole which is asymptotically flat has a negative specific heat of $`C_v^{Sch}=(d2)S_0`$ and thus its canonical ensemble is unstable . This means that the Schwarzschild black hole is never in thermal equilibrium with its environment. But introducing a negative cosmological constant $`\mathrm{\Lambda }=(d1)(d2)/2\mathrm{}^2`$ could make the specific heat positive and render Eq. (1) applicable. The 5D Schwarzschild-AdS black hole with the event horizon $`r_{EH}^2=\mathrm{}^2(1+\sqrt{1+4m/\mathrm{}^2})/2`$ belongs to this category . In the limit of large $`\mathrm{\Lambda }(\mathrm{}0,m>>\mathrm{}^2)`$, one has a large AdS black hole with $`C_v^{SAdS}3S_0`$, while in the limit of small $`\mathrm{\Lambda }(\mathrm{}\mathrm{},m<<\mathrm{}^2)`$, one finds a small AdS black hole which leads to the Schwarzschild black hole with $`C_v^{SAdS}3S_0`$. This implies that for a black hole with a large mass, the presence of a large negative cosmological constant is sufficient to make the specific heat positive.

In other words, a large black hole in asymptotically AdS space has a positive specific heat, while a small back hole in asymptotically AdS space has a negative specific heat. Roughly speaking, one considers the AdS space as a confining box. Then a large black hole in a small box (large AdS black hole) can be thermal equilibrium and gives a positive specific heat as a whole. On the contrary, if the box is large and unbounded as in the small AdS black hole ($``$ Schwarzschild black hole), its specific heat becomes negative and thus this system could not be in thermal equilibrium. Similarly one considers a cosmological horizon in 5D Schwarzschild-de Sitter space . In the limit of large $`\mathrm{\Lambda }(\mathrm{}0)`$ one finds a positive specific heat of the pure de Sitter space $`C_v^{SdS}3S_0`$, while in the limit of small $`\mathrm{\Lambda }(\mathrm{}\mathrm{})`$ one recovers a negative specific heat of the Schwarzschild black hole as $`C_v^{SdS}3S_0`$. The cosmological horizon of 5D hyperbolic topological de Sitter (HTdS) space  has the two limits . Also we expect a similar property in topological Reissner-Nordstrom de Sitter space .

However, the above two limits does not exist for a lower-dimensional gravity system. In other words, this system is simple in compared with 5D gravity systems. A class of 3D gravitational systems have no a singularity but a conical defect (excess) inside the single horizon. Also A(dS)/(E)CFT correspondences can be easily confirmed in connection with logarithmic corrections. In this work, we study 3D AdS black hole and de Sitter spaces which give positive specific heats and thus logarithmic corrections to the entropy could be achieved. These are the BTZ black hole with $`J=0`$ and Kerr-de Sitter (KdS) space with $`J=0`$ including the Schwarzschild-de Sitter (SdS) solution. All of them have single horizons.

## II 3D AdS black hole

It is believed that a higher-dimensional black hole in asymptotically flat spacetime should have a spherical horizon. When introducing a negative cosmological constant, a higher-dimensional black hole could have a non-spherical horizon. We call this the topological black hole. However, in three dimensions we could not construct various horizon geometries because one-dimensional horizon is locally flat and is equivalent to each other. Instead, we introduce an AdS black hole, BTZ black hole with $`J=0`$  as

$$ds_{AdS}^2=h_{AdS}(r)dt^2+\frac{1}{h_{AdS}(r)}dr^2+r^2d\theta ^2,$$
(2)
where a metric function $`h_{AdS}(r)`$ is given by

$$h_{AdS}(r)=\mu _{AdS}+\frac{r^2}{\mathrm{}^2}$$
(3)
with a reduced mass $`\mu _{AdS}=8G_3_{AdS}`$ for a massive point particle located at $`r=0`$ inside the event horizon. Here we choose a circle ($`S^1`$) for the horizon geometry. The event horizon is given by

$$r_{EH}=\mathrm{}\sqrt{\mu _{AdS}}.$$
(4)
Relevant thermodynamic quantities: free energy ($`F`$), Bekenstein-Hawking entropy ($`S_0`$), Hawking temperature ($`T_H`$), specific heat ($`C_v=(dE/dT)_V`$) and energy (ADM mass: $`E=M`$) are given by

$`F={\displaystyle \frac{\mu _{AdS}}{8G_3}},S_0={\displaystyle \frac{\pi \mathrm{}\sqrt{\mu _{AdS}}}{2G_3}}=C_v,`$ (5)
$`T_H={\displaystyle \frac{\sqrt{\mu _{AdS}}}{2\pi \mathrm{}}}=\left[{\displaystyle \frac{G_3}{\pi ^2\mathrm{}^2}}\right]S_0,E=F+T_HS_0={\displaystyle \frac{\mu _{AdS}}{8G_3}}=_{AdS}=M_{AdS},`$ (6)
where $`G_3`$ is the 3D Newton constant. Here we obtain a positive specific heat. In this case we have the ADM mass $`M_{AdS}=_{AdS}`$ which is consistent with that from the Brown-York approach in Ref..

## III 3D Schwarzschild-de Sitter spacetime

The 3D Schwarzschild-de Sitter spacetime is given by

$$ds_{SdS}^2=h_{SdS}(r)dt^2+\frac{1}{h_{SdS}(r)}dr^2+r^2d\theta ^2$$
(7)
where $`h_{SdS}(r)`$ is given by

$$h_{SdS}(r)=1\mu _{SdS}\frac{r^2}{\mathrm{}^2}$$
(8)
with a reduced mass $`\mu _{SdS}=8G_3_{SdS}`$. The cosmological horizon is located at

$$r_{CH}=\mathrm{}\sqrt{1\mu _{SdS}}.$$
(9)
In this case we have to put a particle with a small mass at $`r=0`$ ($`r_{EH}<\mathrm{},0<\mu _{SdS}<1`$). This describes a conical defect spacetime with a deficit angle $`2\pi (1\sqrt{1\mu _{SdS}})`$, indicating a world with a positive cosmological constant and a pointlike object with mass $`_{SdS}`$ . In the case of $`\mu _{SdS}>1`$, any cosmological horizon does not exist.

Interesting thermodynamic quantities for the cosmological horizon are given by

$`F=\left[{\displaystyle \frac{1\mu _{SdS}}{8G_3}}\right],S_0={\displaystyle \frac{\pi \mathrm{}\sqrt{1\mu _{SdS}}}{2G_3}}=C_v,`$ (10)
$`T_H={\displaystyle \frac{\sqrt{1\mu _{SdS}}}{2\pi \mathrm{}}}=\left[{\displaystyle \frac{G_3}{\pi ^2\mathrm{}^2}}\right]S_0,E=F+T_HS_0={\displaystyle \frac{1\mu _{SdS}}{8G_3}}=M_{SdS}.`$ (11)
We obtain a positive specific heat only for a small massive particle ($`\mu <1`$) inside the cosmological horizon. In other words, this gravitating system composed of the cosmological horizon $`r_{CH}<\mathrm{}`$ and a massive particle with mass $`_{SdS}`$ is thermodynamically stable because $`F<0`$ when $`\mu <1`$. If $`_{SdS}=0`$, we recover an ADM mass of $`M_{dS}=1/8G_3`$ for the pure de Sitter space with a cosmological horizon at $`r_{CH}=\mathrm{}`$. In the case of $`_{SdS}0`$, the ADM mass ($`M_{SdS}=1/8G_3_{SdS}`$) of this system is less than that of the pure de Sitter space. Here we address why a particle with $`_{SdS}>0`$ inside the cosmological constant contributes to $`M_{SdS}`$ as a negative one. The reason is that even if the matter of a pointlike object has a positive energy, the binding energy to the gravitational de Sitter background can make it negative. That is, a conical defect at $`r=0`$ swallows up a part of the spacetime and thus appears to reduce the net amount of the positive energy stored in the cosmological constant. In this sense we no longer take $`_{SdS}`$ as a true mass (energy) in SdS space. The true mass is given by the ADM mass $`M_{SdS}`$. As a result, the presence of a pointlike object inside the cosmological horizon decreases the size of the cosmological horizon in compared with the pure de Sitter case.

## IV 3D Kerr-de Sitter space

When introducing a negative mass in higher-dimensional de Sitter space one could find a non-spherical horizon. We call this the topological de Sitter space. However, in three dimensions we could not construct various horizon geometries because one-dimensional horizon is locally flat and is equivalent to each other. Here we introduce the Kerr-de Sitter space with $`J=0`$

$$ds_{KdS}^2=h_{KdS}(r)dt^2+\frac{1}{h_{KdS}(r)}dr^2+r^2d\theta ^2,$$
(12)
where $`h_{KdS}(r)`$ is given by

$$h_{KdS}(r)=\stackrel{~}{\mu }_{KdS}\frac{r^2}{\mathrm{}^2}=\mu _{KdS}\frac{r^2}{\mathrm{}^2}$$
(13)
with a negative reduced mass $`\stackrel{~}{\mu }_{KdS}\mu _{KdS}=8G_3_{KdS}`$ in compared with the SdS case in Eq. (8). The cosmological horizon is given by

$$r_{CH}=\mathrm{}\sqrt{\mu _{KdS}}.$$
(14)
Thermodynamic quantities are calculated as

$`F={\displaystyle \frac{\mu _{KdS}}{8G_3}},S_0={\displaystyle \frac{\pi \mathrm{}\sqrt{\mu _{KdS}}}{2G_3}}=C_v,`$ (15)
$`T_H={\displaystyle \frac{\sqrt{\mu _{KdS}}}{2\pi \mathrm{}}}=\left[{\displaystyle \frac{G_3}{\pi ^2\mathrm{}^2}}\right]S_0,E=F+T_HS_0={\displaystyle \frac{\mu _{KdS}}{8G_3}}=_{KdS}=M_{KdS}.`$ (16)
Here we obtain a positive specific heat, which means that the KdS space is thermodynamically stable. We observe that there is no difference between Eqs. (5) and (15).

Let us ask of how to understand the positivity of ADM mass ($`M_{KdS}>0`$) in the presence of a pointlike object with a negative energy. The topological de Sitter space in five-dimensions has a positive mass and a timelike singularity which remains within a single cosmological horizon. Hence it satisfies the mass bound conjecture: any asymptotically de Sitter space whose mass exceeds that of pure de Sitter space contains a cosmological singularity . In the case of 3D KdS space, there is no cosmological singularity. However, the negative energy of a pointlike object at $`r=0`$ gives the ADM mass of $`M_{KdS}=_{KdS}`$ in asymptotically de Sitter space. As a result, one finds that $`r_{CH}^{KdS}=\mathrm{}\sqrt{\mu _{KdS}}>r_{CH}^{dS}=\mathrm{}`$ if $`\mu _{KdS}>1`$. This corresponds to a conical excess. If $`\mu _{KdS}<1`$, one finds a conical defect.

## V correction to Bekenstein-Hawking entropy

We summarize our result in TABLE I. In this section we make corrections to the Bekenstein-Hawking entropy according to the formula of Eq.(1). For all cases, one finds $`C_v=S_0,T_H=\left[\frac{G_3}{\pi ^2\mathrm{}^2}\right]S_0`$ without any approximation. Here we require a small SdS solution for obtaining a positive specific heat. All logarithmic corrections to the Bekenstein-Hawking entropy are given by a single relation as

$$S=S_0\frac{3}{2}\mathrm{ln}S_0+\mathrm{}.$$
(17)
This means that any thermodynamically stable system gets a logarithmic correction to the Bekenstein-Hawking entropy due to thermal fluctuation around the equilibrium state. Here we confirm that this correction is universal for all 3D gravitational systems with $`J=0`$. Further we wish to mention that all thermal quantities of the KdS spaces take the same form as those from the BTZ black holes. This implies that there is no distinction between the event horizon of BTZ black hole and the cosmological horizon of Kerr-de Sitter space. That is, a pointlike object with a positive energy inside the event horizon is equivalent to a pointlike object with a negative energy inside the cosmological horizon. We note that even for $`\mu _{SdS}=0`$, the SdS space has a cosmological horizon. If $`\mu _{SdS}0`$, the size of the cosmological horizon is decreased.

## VI AdS/CFT and dS/ECFT correspondence

The holographic principle means that the number of degrees of freedom associated with the bulk gravitational dynamics is determined by its boundary spacetime without gravity. Let us first discuss the AdS black hole. The AdS/CFT correspondence represents a realization of this principle . Brown and Henneaux have shown that asymptotic symmetries in 3D AdS spacetime are described by a pair of Virasoro algebras with central charges $`c=\overline{c}=3\mathrm{}/2G_3`$ . The generators of Brown-Henneaux Virasoro algebras are simply Hamiltonian and momentum constraints of 3D gravity smeared against appropriate vector fields. For BTZ black hole with $`J=0`$, these are $`L_0=M\mathrm{}/2=\overline{L}_0`$ with the ADM mass $`M=_{AdS}`$ . In order to calculate the entropy on the boundary spacetime, we use the Cardy formula for a pair to yield

$$S_{CFT}=2\pi \sqrt{\frac{cL_0}{6}}+2\pi \sqrt{\frac{\overline{c}\overline{L}_0}{6}}=\frac{\pi r_{EH}}{2G_3}=S_0$$
(18)
which establishes the AdS/CFT correspondence for the entropy. According to the Carlip’s work , the density of states can be approximated as

$$\rho (L_0,\overline{L}_0)=\frac{8G_3\mathrm{}}{r_{EH}^3}\mathrm{exp}\left[\frac{\pi r_{EH}}{2G_3}\right].$$
(19)
Thus the entropy correction due to thermal fluctuation of the CFT system leads to

$$S=S_03\mathrm{ln}S_0+\mathrm{}.$$
(20)
Comparing the above with the bulk correction in Eq. (17), one finds logarithmic term that differs from it by a factor of 2. However, this difference does not seem to be serious because other approaches suggest a correct factor of 3/2 .

For SdS and KdS spaces, according to Strominger , the dS/ECFT proposal is useful for calculating the logarithmic corrections to the boundary entropy. We note that 3D de Sitter space can be represented by the group manifold of SL(2,C)/SL(2,R). The asymptotic symmetry group of de Sitter space, subject to the boundary conditions, is the two-dimensional Euclidean conformal group SO(3,1), which contains the isometry group SL(2,C) as a subgroup. Indeed, the Brown-Henneaux analysis can be formally continued to arrive at Virasoro algebras with central charges $`c=\overline{c}=3\mathrm{}/2G_3`$ . Its dual theory might be a Lorentzian CFT with Euclidean signature Virasoro algebras. Actually the ADM masses $`M_{SdS},M_{KdS}`$ in Eqs.(10) and (15) are consistent with those computed with the Euclidean signature on the spacelike boundary. Thus we have generators as $`L_0=M\mathrm{}/2=\overline{L}_0`$ where the ADM mass is given by $`M_{SdS}=1/8G_3_{SdS}`$ for the SdS case and $`M_{KdS}=_{KdS}`$ for the KdS space. At this stage, we remind the reader that for an observer near infinity, there is no actual difference between the event horizon of AdS black hole and the cosmological horizon of SdS and KdS space. The difference results in the precise value of a conserved quantity, the ADM mass. Hence we expect that their dual ECFTs exist on the boundary. The Cardy formula provides

$$S_{ECFT}=2\pi \sqrt{\frac{cL_0}{6}}+2\pi \sqrt{\frac{\overline{c}\overline{L}_0}{6}}=\frac{\pi r_{CH}}{2G_3}=S_0$$
(21)
which establishes the dS/ECFT correspondence for the entropy. Following the dS/ECFT correspondence, their logarithmic corrections will be given by the same form as in Eq. (20).

On the other hand, a complete Cardy formula for the asymptotic density of states of a unitary (E)CFT is given by

$$S_{CCFT}=2\pi \sqrt{\frac{c}{6}\left(L_0\frac{c}{24}\right)}+2\pi \sqrt{\frac{\overline{c}}{6}\left(\overline{L}_0\frac{\overline{c}}{24}\right)}.$$
(22)
A naive application of this formula with the same central charge and generators fails to find the Bekenstein-Hawking entropy $`S_0`$, in contrast with the use of the Cardy-formula in Eq. (18). However, in the case of $`M_{AdS/KdS}\mathrm{}/2>>c/24=\mathrm{}/16G_3`$, the complete Cardy formula is given approximately

$$S_{CCFT}\frac{\pi \mathrm{}\sqrt{\mu _{AdS/KdS}}}{2G_3}=S_0.$$
(23)
Hence we find that the Cardy formulae in Eqs. (18) and (21) are considered as approximate forms to the complete Cardy formula in Eq. (22).

## VII Discussion

First of all, we summarize our main result. We find that all thermal quantities of KdS space take the same form as those from BTZ black hole. This means that there is no difference between the event horizon and cosmological horizon in 3D gravity systems with the single horizon.

We have a few of comments in order. Any ADM mass for a gravitational system with a single horizon including the BTZ black hole with $`J=0`$, the Kerr-de Sitter space with $`J=0`$, and SdS solution is always positive. In three dimensions, any system with a single horizon has a positive specific heat and its canonical ensemble is thermodynamically stable. Furthermore, this system has an equality of $`C_v=S_0`$, which dictates that it has a single horizon in three dimensions. Thus we find the same form for logarithmic corrections to the Bekenstein-Hawking entropy and the Cardy formula.

Concerning A(dS)/(E)CFT correspondences, we note that these are applicable for the event horizon for AdS black hole and the cosmological horizons for KdS and SdS spaces. Boundary systems have the same central charge of $`c=\overline{c}=3\mathrm{}/2G`$ with generators $`L_0=\overline{L}_0=M\mathrm{}/2`$ where $`M`$ is one ADM mass for AdS black holes, KdS and SdS spaces. In general the boundary entropy is equal to the Bekenstein-Hawking entropy and its logarithmic correction is given by Eq. (20).

Finally we wish to remark the difference between 3D and 5D gravity systems. In 5D gravity, the parameter $`k`$ is introduced to classify the 3D horizon geometry: $`k=0,1,1`$ denote the flat, elliptic, and hyperbolic spaces with constant curvature. That is, it is necessary to classify topological AdS black holes and topological de Sitter spaces. In 3D gravity, however, the horizon is one-dimensional space and it is locally flat. Hence the parameter like $`k`$ is irrelevant to 3D gravity systems. Even if $`k`$ is introduced, it could be absorbed into $`\mu `$ by the redefinition. Therefore, in this work we investigated thermal properties of three known examples: BTZ black hole with $`J=0`$, Kerr-de Sitter space with $`J=0`$, and Schwarzschild-de Sitter space. In other words, the BTZ black hole is representative of 3D topological AdS black holes, while the Kerr-de Sitter is representative of 3D topological Sitter spaces.

## Acknowledgments

The author thanks Rong-Gen Cai and Mu-In Park for helpful discussions. This work was supported in part by KOSEF, Project No. R02-2002-000-00028-0.


text:
# The long reach of DNA sequence heterogeneity in diffusive processes

## I Introduction

Diffusion appears in most basic processes in the living matter and therefore has been studied extensively by theoretical and experimental biophysicists for many decades. At the macroscopic scale, the phenomena are adequately described by continuum models that form a well established methodology finding many applications in science and technology Murray (2002). Advanced experimental methods, such as nanoprobing and single-molecule techniques, provide us with a wealth of data at the microscopic level. Theoretical description of the observed phenomena at such scales is often a considerable challenge, since many irregular features that average out on the macroscopic scale cannot be ignored anymore. Sometimes, however, rather simple characteristics emerge, allowing for exact analytic treatment.

One-dimensional (1D) transport is rarely found on the macroscopic scale; at the molecular level though, one can find several examples, e.g. kinesin motion along microtubules Bustamante et al. (2001); R. D. Vale and R. A. Milligan (2000); R. D. Astumian (1997) or DNA translocation through a nanopore D. K. Lubensky and D. R. Nelson (1999); P. G. de Gennes (1999); Salman et al. (2001); Meller et al. (2001). Usually, in such problems, the underlying potential profile is considered to be constant or at least regular. However, as we show in this paper, DNA sequence heterogeneity and the resulting random energy landscape can have a considerable influence on the diffusion up to biologically relevant length scales at room temperatures.

### I.1 Protein-DNA interaction

The first example we study here arises in the context of protein-DNA interaction. As proposed by von Hippel and Berg O. G. Berg et al. (1981); P. H. von Hippel and O. G. Berg (1989), and recently observed in many systems Shimamoto (1999), 1D “sliding” of proteins along the DNA molecule is an important component of protein specific site location; at least in prokaryotes. The “sliding” is viewed as an unbiased, thermally activated process. The actual rules of motion for sliding depend on the details of interaction between the protein and the DNA. The general belief is that there are two protein-DNA binding modes: a strong “specific” mode that characterizes binding of operator sites, and a much weaker “non-specific” mode in which binding of non-cognate DNA occurs P. H. von Hippel and O. G. Berg (1989); O. G. Berg and P. H. von Hippel (1987); Gerland et al. (2002); R. F. Bruinsma (2002). In the “non-specific” or “search” mode, the interaction energy is usually assumed to be independent of the DNA sequence that the protein is bound to, though not much experimental evidence beside relatively fast observed search times favors this strictly “equipotential” picture. On the other hand, scanning force microscopy experiments by Erie et. al. D. A. Erie et al. (1994) clearly demonstrate DNA bending by Cro repressor protein, both at operator and at non-operator sequences 111DNA bending by transcription factors is a well-known phenomenon, though practically all the available experimental data focus on proteins bound to operator sequences.. Since local DNA elasticity is known to be highly sequence-dependent M. G. Munteanua et al. (1998), the energy of protein bound at random locations should have a random component, correlated at length scales of the order of the protein binding domain size; see Fig. 1a. This sequence-dependent interaction energy component appears in addition to possible local uncorrelated sequence-dependent contributions from amino acid - base pair contacts.

To estimate the significance of the random component of the elastic energy, we use DNA elasticity data supplied by the BEND.IT server Vlahoviček et al. (2003), that incorporates DNase I based bendability parameters Brukner et al. (1995) and the consensus bendability scale Gabrielian and Pongor (1996). We assume that the protein-DNA complex in Fig. 1a has a fixed geometry, i.e. the protein is “hard.” Then, the random component of the binding energy $`\delta U`$ is proportional to the random component of the Young’s modulus $`\delta E`$

$$\delta U=\frac{\delta E}{\overline{E}}\left(\frac{\mathrm{}_p\theta ^2}{2L}\right)k_BT,$$
(1)
where $`\mathrm{}_p50\mathrm{nm}`$ is the DNA persistence length, $`\theta 60^{}`$ is the curvature angle D. A. Erie et al. (1994), $`L=1020\mathrm{bp}`$ is the bent sequence length and $`\overline{E}3.4\times 10^8\mathrm{N}/\mathrm{m}`$ is the average Young’s modulus. The resulting potential profile is plotted in Fig. 2a. The standard deviation of the random component is $`(\delta U)^2^{1/2}0.51.5k_BT`$, so that disorder appears to be relevant for this problem.

Another interesting example, also from the field of protein-DNA interaction, was considered recently by Schiessel et. al. Schiessel et al. (2001), and deals with nucleosome repositioning by DNA reptation. It was argued that chromatin remodeling Widom (1998, 2001) can be readily understood in terms of intranucleosomal loop diffusion, the size of the loop resulting mainly from a compromise between elastic energy and nucleosome-DNA binding energy. Here again, for a given size of the loop, the elastic energy is sequence-dependent Widom (2001), and therefore has a random component with finite correlation length; see Fig. 1b. For nucleosome repositioning, this effect may be even more pronounced than for prokaryotic protein-DNA interaction; the bending angles $`\theta `$ and the sequence lengths $`L`$ are 2-3 times larger so that the net effect may be twice as strong as for the $`Cro`$ repressor Schiessel et al. (2001).

It is known that DNA can have an $`intrinsic`$ curvature arising from the stacking interactions between base pairs. Such sequence-dependent curvature can play a role similar to sequence-dependent DNA bendability in providing a correlated landscape. The bending energy of an intrinsically curved region is easier, requiring a smaller angular deformation $`\theta =\theta _{\mathrm{complex}}\theta _{\mathrm{intrinsic}}`$ by the DNA-protein complex. Such sequence-dependent intrinsic curvature was suggested to be involved in positioning nucleosomes Kiyama and E. N. Trifonov (2002).

Aside from DNA bendability and curvature, local correlations in nucelotide composition, known to be present in eukaryotic genomes, (AT/GC-rich isochores) can result in a correlated landscape of the protein-DNA binding energy. This effect becomes especially pronounced when a DNA-binding protein has a strong preference toward a particular AT/GC composition of its site. However, in this case, variations take place over much longer scales, and are not quantitatively relevant in the specific contexts addressed in this paper.

Both above examples can be viewed as specific cases of DNA reptation by means of a propagating defect (or “slack”) of a fixed size. Elastic energy associated with the slack creation is sequence-dependent and correlated on the scale of the slack size. The propagating defect is well localized and samples the energies of well-defined subsequent DNA segments. As was pointed out by Cule and Hwa Cule and Hwa (1998), short-range correlated randomness of this kind has no effect on the scaling of the reptation time. However, as we show below, the defect motion itself is strongly influenced by the disorder and has nontrivial behavior at different length scales.

### I.2 DNA translocation through a nanopore

Consider a piece of single-stranded DNA (ssDNA) passing through a large membrane channel. If the potential difference across the membrane is zero, the motion of the ssDNA is governed by thermal fluctuations. Since the channel width differs from the ssDNA external diameter only by few Ångstroms 222For $`\alpha `$-haemolysin, the diameter of the limiting aperture is about 15 Å., it is reasonable that local interactions between the nucleotides and the amino acids of the channel take place. These interactions may have a local base-dependent component. In addition, longer-range terms are likely to appear in the presence of a voltage difference. In the cytoplasm, the DNA negative charge is almost completely screened out at distances of few nanometers by the counterion cloud. When the DNA molecule enters the pore, most of the counterions are likely to be “shaven off,” though some of them may remain stuck to the DNA; see Fig. 3. Thus, the linear charge density inside the pore acquires a random and basically uncorrelated component:

$$q(x)=\overline{q}(x)+\delta q(x),\delta q(x)\delta q(y)=\rho ^2d\delta (xy).$$
(2)
The potential energy of the DNA segment inside the pore in the presence of a voltage difference of $`V_0`$ is

$$U(x)=\frac{V_0}{h}_x^{x+h}x^{}q(x^{})𝑑x^{}.$$
(3)
Since the average charge density $`\overline{q}(x)`$ is nonzero, DNA transport is driven by the average force $`V_0\overline{q}(x)/h`$. The correlation function of the random component of $`U(x)`$ is readily calculated to be

$$\delta U(x)\delta U(x+y)=\frac{V_0^2\rho ^2d}{3h^2}(h|y|)^2\left(h+\frac{|y|}{2}\right)H(h|y|),$$
(4)
where $`H(x)`$ is the Heaviside function. Thus, the potential profile for DNA motion has a random component with correlation length of $`h`$. Taking $`V_0100\mathrm{mV}`$, $`\rho e/h`$ ($`e`$ is the elementary charge), $`h10\mathrm{nm}`$, we obtain $`\delta Uk_BT`$.

Although this example differs from the above ones in that a nonzero average driving force is present, large random fluctuations of the energy landscape may have significant effect on the distribution of translocation times - a problem that has attracted much interest lately Meller (2003).

## II Diffusion in a random potential

### II.1 The model

The problems described above map onto a one-dimensional random walk with position-dependent hopping probabilities $`p_i`$, $`q_i=1p_i`$ to the right and to the left, respectively; it is most natural to assume the regular activated transport form

$$p_ie^{\beta \left(U_{i+1}U_i\right)},q_ie^{\beta \left(U_{i1}U_i\right)},$$
(5)
where $`\beta \left(k_BT\right)^1`$ and $`U_i`$ is the sequence-dependent component of the potential energy. The latter is basically a sum of many random contributions and can therefore be considered to be normally distributed Gerland et al. (2002). Thus, in the absence of correlations, the probability for realization of a certain profile $`U(x)`$ of length $`L`$ is (in the continuum limit)

$$P[U(x)]\mathrm{exp}\left[\alpha _0^L𝑑xU^2(x)\right].$$
(6)
This is the well-known Random-Energy Model Derrida (1981) that was applied successfully to various biophysical problems, from protein folding Bryngelson and Wolynes (1987) to protein-DNA interaction Gerland et al. (2002). It assumes no correlations between energies of different sites. One can think of a more general form of potential profile

$$P[U(x)]\mathrm{exp}\left[_0^L_0^L𝑑y𝑑xU(x)G(xy)U(y)\right].$$
(7)
Taking for example, $`G(xy)_{xy}^2\delta (xy)`$, we obtain the Random-Force Model J. P. Bouchaud et al. (1990) that describes an energy landscape appearing as a random walk with linearly growing correlations. This model was studied during the last decades in the context of heteropolymer dynamics De Gennes (1975); Cule and Hwa (1998), glassy systems Le Doussal and Vinokur (1995); D. S. Fisher et al. (1998) and quite recently - to describe DNA denaturation dynamics Hwa et al. (2003). Characteristic features of the Random-Force Model are logarithmically slow (“Sinai’s”) diffusion Sinai (1982); S. H. Noskowicz and I. Goldhirsh (1988) and aging D. S. Fisher et al. (1998); Hwa et al. (2003). More generally, $`G`$ is related to the correllator of $`U`$ by $`U(x)U(y)=G^1(xy)`$.

To include finite-range correlations into Eq. (6), we must incorporate a limitation on the acceptable forces. The ensemble of energy profiles is therefore naturally described by the following probability density

$$P[U(x)]e^{[U]},$$
(8a)
with pseudoenergy
$$[U]=_0^L𝑑x\left[\alpha U^2(x)+\gamma \left(\frac{dU}{dx}\right)^2\right].$$
(8b)
Energy level statistics for this kind of potential profile is also Gaussian, as can be seen from the average

$$e^{ikU}=\frac{𝒟[U]e^{ikU}e^{[U]}}{𝒟[U]e^{[U]}}=\mathrm{exp}\left(\frac{k^2}{8\sqrt{\alpha \gamma }}\right),$$
(9)
which is the characteristic function for Gaussian distribution with zero mean and variance

$$\sigma ^2=\frac{1}{4\sqrt{\alpha \gamma }}.$$
(10)
The correlator of the potential profile is readily calculated as

$$g(r)\frac{1}{2}[U(x)U(x+r)]^2=\sigma ^2\left(1e^{|r|/\xi _c}\right),$$
(11)
where $`\xi _c=\sqrt{\gamma /\alpha }`$ is the correlation length.

### II.2 Mean First Passage Time

A convenient formalism for analyzing diffusion in a random one-dimensional potential profile is that of mean first-passage time S. H. Noskowicz and I. Goldhirsh (1988); K. P. N. Murthy and K. W. Kehr (1989). For a given set of probabilities $`\left\{p_i\right\}`$, the mean first-passage time (MFPT) from $`i=0`$ to $`i=N`$ (in terms of number of steps) is

$$\overline{t}_{0,N}=N+\underset{k=0}{\overset{N1}{}}\omega _k+\underset{k=0}{\overset{N2}{}}\underset{i=k+1}{\overset{N1}{}}\left(1+\omega _k\right)\underset{j=k+1}{\overset{i}{}}\omega _j,$$
(12)
where $`\omega _iq_i/p_i`$ (see Appendix A for derivation). The MFPT given by this expression is for a fixed realization of probabilities, i.e. for a given potential energy profile; as such, it is itself a random variable. The disorder-averaged version of the MFPT is readily obtained after we note that the sequential products in Eq. (12) reduce to

$$\underset{j=k}{\overset{i}{}}\omega _j=\mathrm{exp}\left[\beta (U_{i+1}+U_iU_kU_{k1})\right].$$
(13)
For an uncorrelated potential profile, this exponential factorizes into independent exponentials; after the ensemble averaging and the summations are carried out, we obtain for $`N1`$

$$\overline{t}_{0,N}=N^2e^{2\beta ^2\sigma ^2},$$
(14)
where, for the uncorrelated potential ($`\gamma =0`$)

$$\sigma ^2=\frac{1}{2\alpha d},$$
(15)
where $`d`$ is the lattice spacing. Note that this expression cannot be obtained by simply putting $`\gamma =0`$ in Eq. (10). The reason is that when $`\gamma `$ becomes small, the discrete nature of the underlying lattice starts to matter. The integration in the momentum space extends only up to $`|q_{\mathrm{max}}|=\pi /d`$, and thus,

$$\sigma ^2|_{\gamma 0}=_{\pi /d}^{\pi /d}\frac{dq}{4\pi \alpha }=\frac{1}{2\alpha d}.$$
(16)
Returning to the case of a finite correlation length, we note that in the limit of $`\xi _cd`$, variations of the potential between neighboring sites can be neglected compared to variations between sites separated by distances of order $`\xi _c`$ or larger. Since the main contribution to the MFPT comes from the double sum in Eq. (12), we can write the continuum version as

$$\overline{t}_{0,N}2_0^N𝑑x_x^N𝑑ye^{2\beta (U(x)U(y))}.$$
(17)
To average over all possible realizations of $`\left\{U(x)\right\}`$, we calculate

$`e^{2\beta (U(x)U(y))}={\displaystyle \frac{𝒟[U]e^{2\beta (U(y)U(x))}e^{[U]}}{𝒟[U]e^{[U]}}}`$
$`=\mathrm{exp}\left[{\displaystyle \frac{\beta ^2\xi _c}{\gamma }}(1e^{|xy|/\xi _c})\right].`$ (18)
For $`|xy|\xi _c`$, Eq. (II.2) reduces to $`\mathrm{exp}(\beta ^2|xy|/\gamma )`$, so that for $`N\xi _c`$ we have

$$\overline{t}_{0,N}N^2\mathrm{exp}(4\beta ^2\sigma ^2N/\xi _c).$$
(19)
(Here and in what follows, we measure distances in units of $`d`$, unless specified otherwise.) This kind of exponential creep is quite expected, since for $`\alpha 0`$, $`\xi _c\mathrm{}`$ our model (8a) reduces to the Random-Force Model.

In the opposite limit $`|xy|\xi _c`$, we can neglect the exponent $`e^{|xy|/\xi _c}`$ , so that Eq. (17) produces an ordinary diffusion law, with a disorder-renormalized diffusion coefficient:

$$\overline{t}_{0,N}=N^2e^{4\beta ^2\sigma ^2}.$$
(20)
Comparing Eqs. (20) and (14), we see that diffusion in a correlated potential profile proceeds more slowly than in an uncorrelated profile. It is straightforward to obtain an expression for the disorder-averaged MFPT for arbitrary correlation length. If we keep all four terms in the exponential in Eq. (13) while going to the continuum limit, we obtain

$$\overline{t}_{0,N}2_0^N𝑑x_x^N𝑑ye^{\beta (U(x+d)+U(x)U(y)U(yd))}.$$
(21)
Averaging this expression over the disorder as in Eq. (II.2) yields for $`N\xi _c`$

$$\overline{t}_{0,N}=N^2\mathrm{exp}[2\beta ^2\sigma ^2(1+e^{d/\xi _c})],$$
(22)
which has the obvious limits of Eqs. (14) and (20) for $`\xi _c0`$ and $`\xi _cd`$, respectively.

## III Typical vs average

Large deviations from the average are characteristic to many disordered systems. In this section, we therefore explore the typical properties of random walks as compared to the disorder-averaged ones.

### III.1 Quantifying fluctuations

After the potential profile is generated (see Appendix B), we calculate the MFPT using Eq. (12). Fig. 4a presents the mean first passage times calculated for various realizations of $`U(x)`$ at biologically relevant temperature ($`\sigma k_BT`$). It is clear that although the ensemble-averaged MFPT does behave as prescribed by Eq. (22), typical MFPT exhibits high variability from one profile to another. The stepwise shape of typical curves suggests that a random walk in such a profile consists of regions characterized by subdiffusion (vertical “steps”) and superdiffusion (plateaus), appearing intermittently.

Uncorrelated potential profiles, as Fig. 4b shows, also lead to a certain disorder-induced variability, though of a considerably smaller magnitude. To quantify the sample dependence of the MFPT, we calculate its variance over the ensemble of potential profiles. Fig. 5 presents the standard deviation in $`\overline{t}_{0,N}`$ as a function of $`N`$ for correlated as well as uncorrelated potential profiles. We observe that the variance scales as $`N^3`$ for all profiles. This dependence can be obtained analytically in a quite straightforward fashion. Consider the average of the square of MFPT in a potential profile with correlation length $`\xi _c`$. The leading term is obviously $`N^4\mathrm{exp}(8\beta ^2\sigma ^2)`$ and it comes from independent $`\{i,k,l,m\}`$. The next largest contribution comes from terms with $`i=m`$ or $`k=l`$. There are $`N^3`$ of such terms, each contributing $`\mathrm{exp}(12\beta ^2\sigma ^2)`$. Next, we note that in order to make a contribution of the same order of magnitude, the two indexes ($`i,m`$ or $`k,l`$) should not necessarily coincide exactly; it is sufficient that they are less than one correlation length apart. Hence, after the leading $`O(N^4)`$ term is cancelled by $`\overline{t}_{0,N}^2`$, the variance is

$$(\mathrm{\Delta }\overline{t}_{0,N})^2\xi _cN^3\mathrm{exp}(12\beta ^2\sigma ^2).$$
(23)
Similar reasoning yields for the uncorrelated case

$$(\mathrm{\Delta }\overline{t}_{0,N})^2N^3\mathrm{exp}(6\beta ^2\sigma ^2).$$
(24)
We see that for given $`\sigma `$ and $`\beta `$, the correlated energy landscape produces stronger fluctuations in MFPT than uncorrelated ones, in agreement with Fig. 4.

Comparing the expressions for the variance with the corresponding expressions for disorder-averaged MFPT, we see that for any temperature, there is a characteristic distance $`N_c`$, below which there is no self-averaging and the typical MFPT is determined by fluctuations. This length is

$$N_c\xi _ce^{4\beta ^2\sigma ^2}$$
(25)
for correlated profiles, and

$$N_ce^{2\beta ^2\sigma ^2}$$
(26)
for uncorrelated ones. This effect is akin to “freezing” in the Random-Energy Model Derrida (1981): for low enough temperatures, typical passage times for distances below $`N_c`$ are dominated by high barriers. This is more pronounced for correlated profiles since in addition to stronger temperature dependence, there is amplification by a factor of $`\xi _c`$, as sites within a correlation length give similar contributions. Figure 6 demonstrates the lack of self-averaging for uncorrelated potential profiles at short distances and low temperatures: the median MFPT (defined as the 50th percentile of a sample) shows large deviations from the average at distances shorter than $`N_c`$ and coincides with it at distances larger than $`N_c`$.

### III.2 Anomalous diffusion

The lack of self-averaging in the region $`\xi _cNN_c`$ can be quantified by estimating the typical MFPT. Consider Eq. (13) for an uncorrelated potential and define the following coarsening procedure: $`\stackrel{~}{U}_i=U_{2i}+U_{2i+1}`$. Then, in the “freezing regime,” the double sum

$$\underset{k}{}\underset{i}{}\mathrm{exp}[\beta (\stackrel{~}{U}_i\stackrel{~}{U}_k)],$$
(27)
is dominated by $`(i,k)`$ producing the largest exponent. For a finite sample $`\{U_i\}`$ of size $`N`$ and variance $`\sigma ^2`$, the corresponding sample $`\{\stackrel{~}{U}_i\}`$ contains $`N/2`$ values distributed with a variance $`2\sigma ^2`$. The minimum and the maximum of $`\{\stackrel{~}{U}_i\}`$ have therefore characteristic values of $`\pm 2\sigma \sqrt{\mathrm{ln}[N/(2\sqrt{2\pi })]}`$, respectively. Thus, a typical MFPT for an uncorrelated potential reads

$$\overline{t}_{0,N}\mathrm{exp}\left[4\beta \sigma \sqrt{\mathrm{ln}\frac{N}{2\sqrt{2\pi }}}\right].$$
(28)
For the purposes of estimating the extreme values of a correlated energy landscape, the sample size is effectively reduced by a factor of $`\xi _c`$, therefore, the extrema of $`\{U_i\}`$ are approximately $`\pm \sigma \sqrt{2\mathrm{ln}[N/(\xi _c\sqrt{2\pi })]}`$. Noting that sites within a correlation length around the extrema contribute similarly to the MFPT, for a correlated potential we write

$$\overline{t}_{0,N}\xi _c^2\mathrm{exp}\left[4\beta \sigma \sqrt{2\mathrm{ln}\frac{N}{\xi _c\sqrt{2\pi }}}\right].$$
(29)
Figure 7 compares typical values of $`\overline{t}_{0,N}`$ calculated from Eqs. (28) and (29) with numerically calculated median values of MFPT. We see that our analytical estimates produce a correct order of magnitude for $`\overline{t}_{0,N}`$. As expected, for uncorrelated profiles, the agreement is better at lower temperatures; for higher temperatures, Eq. (28) is an underestimation since we do not include contributions from second-lowest, second-highest, etc., energy levels. Eq. (29), on the other hand, turns out to be a slight overestimation, since we have replaced the average of $`\xi _c^2`$ terms by their maximum value.

Large difference between the median and the average values is a signature of a broad probability distribution. The insets of Fig. 8 present two probability density functions for MFPT, at $`NN_c`$ and $`NN_c`$. For the short distance, the distribution is very broad and spans several orders of magnitude. For $`NN_c`$, the system is self-averaging, in the sense that the MFPT distribution is much narrower with almost coinciding median and average values.

### III.3 Characteristics of random walk

To complete the picture, we perform direct simulations of random walks in correlated and uncorrelated potential profiles; typical results are depicted in Fig. 9. One can see a clear qualitative difference between the two cases: random walks in the uncorrelated profile look very much like standard walks with $`p_i=q_i=1/2`$, whereas motion of a particle in a correlated profile has a somewhat different nature.

As above, we see that macroscopic motion of a particle in a correlated potential consists of subdiffusive as well as superdiffusive segments. It also appears that the particle tends to be localized near the bottom of “valleys” of few $`\xi _c`$ in extent, whereas in an uncorrelated profile, there are no preferable sites for localization. Obviously, when the time is measured in real-time units, rather than in number of steps, the particle is more likely to be found at the minima of the energy landscape in both cases. In terms of the number of steps though, all sites of the uncorrelated landscape are revisited more or less uniformly.

## IV Biological implications

### IV.1 Transcription Factors

Consider a DNA-binding protein searching for its target site on the genome. As explained in the introduction, a correlated random energy landscape can arise from the interplay of sequence-dependent flexibility, and the bending contribution to the total DNA-binding energy. Diffusion on such a landscape may then lead to localization in the energy ‘valleys,’ i.e. the protein will reside preferentially in specific (favorable) areas of the genome. Such nonuniform sampling has important implications for biological strategies of transcription factor bindings: First, if a valley contains several binding sites, the rapid (superdiffusive) scanning of the valley leads to quick equilibration between these sites (while equilibration for similarly spaced sites outside a valley will take much longer). This is important when the protein binds nearby sites with distinct binding energies, and the strongest one has to be occupied first to provide correct regulation (as in the case of the cro repressor). Second, several proteins bind their specific sites only when activated by ligands (e.g. $`PurR`$, $`GalS`$ etc), spending the rest of the time in an inactive form “waiting” for the ligand. These proteins can benefit from staying close the site in the waiting mode, since they can then quickly find their target upon activation.

One of the results of this study was that inhomogeneities significantly reduce the overall diffusion rate, as in Eq. (22). While this may be beneficial in confining a protein to favorable regions, it severely restricts the ability to search large portions of the genome by one dimensional diffusion. Since we argue that a portion of the inhomogeneity originates from variations in the bending energy of the DNA, a potential strategy is for the binding protein to switch between two states which bend the DNA weakly or strongly. The weak bending state is subject to reduced variations in the energy landscape and can diffusive more freely (search mode), compared to the strongly bending state which is more likely to be confined in the vicinity of favorable energy valleys (waiting mode). One potential candidate for exploiting this strategy is the tertarmeric $`LacI`$ protein that consists of two DNA-binding dimeric subunits. Each subunit binds DNA and bends it slightly; when both subunits are bound, DNA is deformed into an extended loop. Several experimental results suggest that $`LacI`$ binds DNA with only one subunit while searching for its target site (“holding DNA with one arm”). Only when both both subunits find their site, the DNA is bent into a loop. Very few structural data are available for proteins bind to DNA non-specifically (search mode). The above strategy suggests that DNA is less deformed in such complexes.

Another potential source for a correlated inhomogeneous energy landscape is an extended protein–DNA interface with net interactions that are the sum of several local contributions. (The addition of such correlated contributions leads to a much larger variance of energy than if they were uncorrelated.) This can be a significant effect for large multi-protein complexes (such as polymerases, TFIID, TFIIB complexes in yeast, etc.). To avoid slow-down by such inhomogeneities, protein complexes can avoid scanning DNA in the fully assembled state when the protein–DNA interface is extensive. Individual components of the complex can search for their sites independently, assembling the whole complex only on the right site. In fact, most of large protein–DNA complexes follow this strategy of assembly on the site, while many dimers and tetramers are assembled in the solution.

### IV.2 Nucleosomes

Other implications concern nucleosome positioning and dynamics. Wrapping of the DNA around these large multi-protein complexes is essential for packing DNA in the small volume of the cell nucleus. Nucleosomes, however, prevent transcription factors and other proteins from accessing DNA. To allow a transcription factor to access its target, nucleosomes close to that site have to be removed from the DNA or re-positioned. While removal of nucleosomes is made by specific enzymes that chemically modify them (e.g. by histon methilation), re-positioning relies in part on nucleosome mobility. In general, nucleosomes have to be (i) positioned at specified locations, and (ii) be able to move along the DNA in the vicinity to the initial placement site allowing access to this region of the DNA.

Nucleosome positioning is determined by specific sequences on the DNA. Such sequences are also known to provide DNA flexibility and/or internal curvature Thastrom et al. (1999); Kiyama and E. N. Trifonov (2002). As discussed above, local DNA flexibility and curvature create a correlated energy landscape for binding. We suggest that inhomogeneous diffusion on such landscapes is an important element that provides both (i) preferential positioning of the nucleosomes due to DNA flexibility and curvature, and (ii) relatively rapid diffusion within the confines of the energy valley. Conversely, uncorrelated landscapes cannot achieve both objectives, since strong nucleosome binding sites prevent local diffusion along the DNA, while weak sites are not able to localize these proteins, leading to their random placement. In fact, experiments Thastrom et al. (1999) have shown that nucleosome positioning sites are extended and are fairly weak. Such structure of positioning sites creates an extended valley on the correlated binding landscape, supporting our hypothesis.

This mechanism can also explain how certain proteins (such as HMGB) can reposition nucleosomes by binding to the DNA in their proximity. It has been suggested that such proteins alter the local mechanical properties of the DNA (such as its flexibility, curvature, or super-coiling) leading to repositioning of the nucleosome A. A. Travers (2003). If the nucleosome is indeed preferentially localized by being trapped in a valley of the binding landscape, HMGB proteins may well alter the shape of the valley (e.g. by shrinking it on one side). Mobile nucleosomes, rapidly diffusing within the boundaries of the valley, will then reposition themselves in the new landscape.

### IV.3 Translocation

In Sec. I.2 we described how slow (activated) passage of ssDNA through a nanopore can be modeled by diffusion over a correlated landscape. In particular, we demonstrated that if there are inhomogeneities in the charge of the DNA inside the channel, there will be variations in the potential energy landscape that are proportional to the applied voltage difference $`V`$. There is in fact scant structural information about the reconfigurations of charges (both free and bound) as DNA passes through a channel. Examining the variations in the MFPT of DNA as a function of the applied voltage Meller et al. (2001), may provide an indirect probe of any inhomogeneities in the charge passing through a channel.

## V Conclusions

We studied one-dimensional diffusion in a random energy landscape with short-range correlations. We found that disorder with short correlation length $`\xi _c`$ leads to a strong sample dependence of diffusion characteristics. The diffusive transport is influenced up to length scales exceeding $`\xi _c`$ by orders of magnitude. Three diffusion regimes can be identified:

1. For distances smaller than the correlation length ($`N\xi _c`$), the disorder-averaged Mean First Passage Time (MFPT) is

$$\overline{t}_{0,N}N^2\mathrm{exp}(4\sigma ^2\beta ^2N/\xi _c).$$

At biologically relevant temperatures, the $`N^2`$ factor prevails; however, at low temperatures ($`k_BT2\sigma /\sqrt{\xi _c}`$), we obtain exponential creep (Sinai’s diffusion).
2. For distances $`N`$ much larger than the characteristic value $`N_c`$, MFPT exhibits some variability from sample to sample. However, the typical value of the MFPT is given by the disorder-averaged MFPT

$$\overline{t}_{0,N}=N^2\mathrm{exp}[2\beta ^2\sigma ^2(1+e^{d/\xi _c})].$$

The variance of MFPT over the ensemble of potential profile realizations scales as $`N^3`$ with distance above $`N_c`$. The characteristic distance $`N_c`$ equals $`\xi _ce^{4\beta ^2\sigma ^2}`$ for correlated profiles and $`e^{2\beta ^2\sigma ^2}`$ for uncorrelated ones.
3. In the intermediate case $`\xi _cNN_c`$, the disorder-averaged MFPT behaves as described by Eq. (22). However, the MFPT distribution over the ensemble of profile realizations is much broader below $`N_c`$ than above it, as Fig. 8 demonstrates. As a result, a typical sample yields diffusion times orders of magnitude shorter than the average. This effect can be qualitatively understood in terms of the Random Energy Model. Below $`N_c`$, diffusion times are mostly influenced by high barriers and deep valleys that are at the extrema of energy landscape histogram. The typical diffusion times are given by

$$\overline{t}_{0,N}\mathrm{exp}\left[4\beta \sigma \sqrt{\mathrm{ln}\frac{N}{2\sqrt{2\pi }}}\right]$$

for an uncorrelated profile, and

$$\overline{t}_{0,N}\xi _c^2\mathrm{exp}\left[4\beta \sigma \sqrt{2\mathrm{ln}\frac{N}{\xi _c\sqrt{2\pi }}}\right],$$

for a correlated one. Above $`N_c`$, most obstacles to the particle motion lie in the central region, so that Eq. (22) produces a valid estimation for a typical diffusion time: the system becomes self-averaging.

These regimes appear to be relevant for biological systems and provide qualitative insight into the kinetics of protein-DNA interaction.

###### Acknowledgements.

This work was supported by the National Science Foundation through grant No. DMR-01-18213 (M.K.). L.M. is an Alfred P. Sloan Research Fellow and is also supported by John F. and Virginia B. Taplin Award. M.S. is supported by NEC research fund.

## Appendix A Mean First-Passage Time Derivation.

The mean first passage time (MFPT) from site $`\mathrm{\#}0`$ to site $`\mathrm{\#}N`$ is defined as the mean number of steps the particle has to make in order to reach site $`\mathrm{\#}N`$ for the first time. The derivation here follows the one in Ref. K. P. N. Murthy and K. W. Kehr (1989).

Let $`P_{i,j}\left(n\right)`$ denote the probability to start at site $`\mathrm{\#}i`$ and to reach the site $`\mathrm{\#}j`$ in exactly $`n`$ steps. Then, for example,

$$P_{i,i+1}\left(n\right)=p_iT_i\left(n1\right),$$
(30)
where $`T_i\left(n\right)`$ is defined as the probability of returning to the $`i`$-th site after $`n`$ steps without stepping to the right of it. Now, all the paths contributing to $`T_i\left(n1\right)`$ should start with the step to the left and then reach the site $`\mathrm{\#}i`$ in $`n2`$ steps, not necessarily for the first time. Thus, the probability $`T_i\left(n1\right)`$ can be written as

$$T_i\left(n1\right)=q_i\underset{m,l}{}P_{i1,i}\left(m\right)T_i\left(l\right)\delta _{m+l,n2}.$$
(31)
We now introduce generating functions

$$\stackrel{~}{P}_{i,j}\left(z\right)=\underset{n=0}{\overset{\mathrm{}}{}}z^nP_{i,j}\left(n\right),\stackrel{~}{T}_i\left(z\right)=\underset{n=0}{\overset{\mathrm{}}{}}z^nT_i\left(n\right).$$
(32)
One can easily show (see e.g. Ref. Goldhirsh and Gefen (1986)) that

$$\stackrel{~}{P}_{0,N}\left(z\right)=\underset{i=0}{\overset{N1}{}}\stackrel{~}{P}_{i,i+1}\left(z\right).$$
(33)
Knowing $`\stackrel{~}{P}_{i,i+1}\left(z\right)`$, one calculates the MFPT straightforwardly as

$`\overline{t}_{0,N}={\displaystyle \frac{_nnP_{0,N}\left(n\right)}{_nP_{0,N}\left(n\right)}}=\left[{\displaystyle \frac{d}{dz}}\mathrm{ln}\stackrel{~}{P}_{0,N}\left(z\right)\right]_{z=1}`$
$`={\displaystyle \underset{i=0}{\overset{N1}{}}}\left[{\displaystyle \frac{d}{dz}}\mathrm{ln}\stackrel{~}{P}_{i,i+1}\left(z\right)\right]_{z=1}.`$ (34)
Using Eqs. (30) and (31), we obtain the following recursion relation for $`\stackrel{~}{P}_{i,i+1}\left(z\right)`$:

$$\stackrel{~}{P}_{i,i+1}\left(z\right)=\frac{zp_i}{1zq_i\stackrel{~}{P}_{i1,i}\left(z\right)}.$$
(35)
To solve for $`\overline{t}_{0,N}`$, we must introduce boundary conditions. Let $`p_0=1,q_0=0`$, which is equivalent to introducing a reflecting wall at $`i=0`$. This boundary condition clearly influences the solution for short times and distances. However, as numerical simulations suggest, its influence relaxes quite fast, so that for longer times, the result is clearly independent of the boundary. The benefit of setting $`p_0=1`$ becomes clear when we observe that

$$\stackrel{~}{P}_{0,1}\left(1\right)=1,i\stackrel{~}{P}_{i,i+1}\left(1\right)=1.$$
(36)
Hence,

$$\overline{t}_{0,N}=\underset{i=0}{\overset{N1}{}}\stackrel{~}{P}_{i,i+1}^{}\left(1\right).$$
(37)
The recursion relation for $`P_{i,i+1}^{}\left(1\right)`$ is readily obtained from Eq. (35):

$$\stackrel{~}{P}_{i,i+1}^{}\left(1\right)=\frac{1}{p_i}+\frac{q_i}{p_i}\stackrel{~}{P}_{i1,i}^{}\left(1\right)=1+\omega _i\left[1+\stackrel{~}{P}_{i1,i}^{}\left(1\right)\right],$$
(38)
with $`\omega _iq_i/p_i`$. Thus, the expression for $`\overline{t}_{0,N}`$ is obtained in the closed form as

$$\overline{t}_{0,N}=N+\underset{k=0}{\overset{N1}{}}\omega _k+\underset{k=0}{\overset{N2}{}}\underset{i=k+1}{\overset{N1}{}}\left(1+\omega _k\right)\underset{j=k+1}{\overset{i}{}}\omega _j.$$
(39)
## Appendix B Potential Profile Generation

Given the pseudoenergy partition function

$$Z(\lambda )=𝒟[U]e^{\lambda [U]},$$
(40)
the average pseudoenergy is

$$=\frac{}{\lambda }\mathrm{ln}Z(\lambda )|_{\lambda =1},$$
(41)
and the variance is

$$(\mathrm{\Delta })^2=^2^2=\frac{^2}{\lambda ^2}\mathrm{ln}Z(\lambda )|_{\lambda =1}.$$
(42)
Straightforward calculation for the pseudoenergy given by Eq. (8a) yields

$$=L/2,(\mathrm{\Delta })^2=L/2.$$
(43)
Hence, typical potential profiles have pseudoenergies in the range $`L/2\pm \sqrt{L/2}`$. This result together with Gaussian statistics of energy levels of Eq. (9) forms the basis of the algorithm we employ for building the energy profiles. First, a random and uncorrelated potential profile obeying Gaussian statistics with the required variance $`\sigma ^2`$ is generated on a one-dimensional lattice. Next, we look for a permutation of lattice sites that produces a typical pseudoenergy $`[U]`$ for a given correlation length $`\xi _c`$ (or, equivalently, for given values of $`\alpha `$ and $`\gamma `$). This is accomplished by a Metropolis-type algorithm that converges to a prescribed value of pseudoenergy picked at random from Gaussian distribution around $``$; see Fig. 10.


text:
# Magnetic Spin-Up of Line-Driven Stellar Winds

## 1. Background and Context

Within the broad theme of this conference on “Magnetic Fields in O, B, and A Stars”, a key issue is the role of magnetic fields in channeling the line-driven stellar wind outflows from such hot, luminous stars. Initial results of full MHD simulations for models without stellar rotation (ud-Doula and Owocki 2002) are reviewed in these proceedings through the write-up for the talk presented by A. ud-Doula. Here we summarize further simulation results for rotating hot-stars with a magnetic dipole aligned to the stellar rotation axis. These simulations provide a timely dynamical test of the “Magnetically Torqued Disk” (MTD) paradigm, which, as reviewed in the talk by J. Brown in these proceedings, has recently been promoted as a central mechanism for producing the disks inferred in Be stars (Cassinelli et al. 2002).

To allow greater emphasis on our MHD results, the broader discussion of alternative Be-disk models given in the oral version of the present paper will not be reproduced here, as most of this has already been included in the write-up for a related talk presented at the recent IAU Symposium 215 on Stellar Rotation (Owocki 2003). A key general result of these MHD simulations regards the tendency for magnetic torquing of wind material to result in centrifugal mass ejection rather a circumstellar disk.

## 2.  Alfven Radius vs. Keplerian Rotation Radius

Recent MHD simulations (ud-Doula 2002; ud-Doula and Owocki 2002) indicate that the effectiveness of magnetic fields in channeling a stellar wind outflow can be characterized by the ratio of the magnetic to wind energy densities

$$\eta (r)\frac{B^2/8\pi }{\rho v^2/2}=\eta _{}\frac{(r/R_{})^{22q}}{(1R_{}/r)^\beta }.$$
(1)
Here $`\eta _{}B_{}^2R_{}^2/(\dot{M}v_{\mathrm{}})`$ defines an overall “magnetic confinement parameter” in terms of the strength of the equatorial field $`B_{}`$ at the stellar surface radius $`R_{}`$, and the wind terminal momentum $`\dot{M}v_{\mathrm{}}`$. The latter equality thus isolates the radial variation in terms of a magnetic power-law index $`q`$ ($`=3`$ for a dipole) and a velocity index $`\beta `$ ($`1`$ for a standard CAK wind). If, for simplicity, we ignore the wind velocity variation (i.e. by taking $`\beta =0`$), we can easily solve for an “Alfven radius” $`\eta (R_A)1`$ at which the magnetic and wind energy densities are equal

$$R_A=\eta _{}^{1/4}R_{}.$$
(2)
As shown by simulation results summarized below, this Alfven radius provides a reasonable estimate for the maximum radius of closed loops in a wind outflow. Moreover, since in rotating models such closed loops tend to keep the outflow in rigid-body rotation with the underlying star, it also defines the radius of maximum rotational spin-up of the wind azimuthal speed.

To characterize such rotational effects, let us next define a Keplerian corotation radius $`R_K`$ at which rigid-body rotation would yield an equatorial centrifugal acceleration that just balances the local gravitational acceleration from the underlying star,

$$\frac{GM}{R_K^2}=\frac{v_\varphi ^2}{R_K}=\frac{V_{eq}^2R_K}{R_{}^2},$$
(3)
where $`V_{eq}`$ is the stellar surface rotation speed at the equator. This can be solved to yield

$$R_K=W^{2/3}R_{},$$
(4)
where $`WV_{eq}/V_{crit}`$, with $`V_{crit}\sqrt{GM/R_{}}`$ the critical rotation speed.

Finally, it is also worth noting here that corotation out to an only slightly higher “escape radius”,

$$R_E=2^{1/3}R_K=2^{1/3}W^{2/3}R_{},$$
(5)
would imply an azimuthal speed that equals the local escape speed from the star’s gravitational field.

## 3. MHD Models with Parameters Optimized for Keplerian Spin-up

The above scalings suggest that a likely necessary condition for propelling outflowing material into a Keplerian disk is to choose a combination of parameters for magnetic confinement vs. stellar rotation such that $`R_K<R_A<R_E`$. In the parameter plane of $`\sqrt{\eta _{}}B_{}`$ vs. $`W`$ defined in figure 1, the required combination is represented by the relatively narrow gray band. The dark region below this represents parameter combinations for which the azimuthal speed at the Alfven radius is expected to be sub-Keplerian, while the white region above represents cases for which the rotation speed at the Alfven radius should exceed the local gravitational escape speed.

As a sample test case, we focus here on the specific combination $`\eta _{}=10`$ and $`W=1/2`$, which as shown by the dot in figure 1, lies in the middle of the gray domain, and thus should represent an optimal case for magnetic spin-up into Keplerian orbit.

Figures 2 and 3 illustrate results of 2D MHD simulations for this case, using the approach and general model assumptions described in ud-Doula and Owocki (2002), but now extended to include field-aligned rotation. Figure 2 shows that conditions at a time 90 ksec after introduction of the field do superficially resemble a magnetically torqued disk. In the upper panel the darkest regions – which represent a near rigid-body rotation with $`V_\varphi /V_{rigid}1`$ – include essentially the entire closed magnetic loop. The log-density grayscale in the lower panel shows moreover that magnetic channeling of the outflowing stellar wind leads to formation of a dense equatorial compression that might be interpreted as the expected circumstellar disk.

Closer examination shows, however, that most of this equatorial compression does not have the appropriate velocity for a stable, stationary, Keplerian orbit. Thus in just a few ksec of subsequent evolution, this putative “disk” becomes completely disrupted, characterized generally by infall of the material in the inner region, i.e. below the Keplerian radius $`R_K`$, and by outflow in the outer region above this Keplerian radius. Figure 3 illustrates the irregular form of the dense compression at an arbitrarily chosen later time (390 ksec from the initial start). The arrows emphasize the flow divergence of the dense material both downward and upward from the Keplerian radius. This evolution is most vividly illustrated through animations, which can be viewed on the web at:

www.bartol.udel.edu/$``$owocki/animations/den4wp5eta10.avi

We have carried out similar MHD simulations for a moderately extensive set of combinations for the rotation and magnetic confinement parameters. In all cases we find that any equatorial compressions are dominated by radial inflows and/or outflows, with no apparent tendency to form a steady, Keplerian disk. In weak magnetic spin-up cases with $`R_A<R_K`$, material trapped on closed loops generally falls back on the star, much as in the non-rotating models discussed in the talk by A. ud-Doula (see also ud-Doula 2002 and ud-Doula and Owocki 2002). The intermediate magnetic spin-up cases with $`R_K<R_A<R_E`$ (lying in the gray area of figure 1) generally show a combination of infall and outflow, much as described above. Finally, for the strong magnetic spin-up models with $`R_A>R_E`$ the material density tends to build up at the tops the highest loops before breaking out in semi-regular epsiodes of discrete “mass ejections”.

Unfortunately technical issues make it difficult to carry out simulations with very high magnetic fields and/or very rapid stellar rotation. For rotations near the critical speed, the star becomes distorted into an oblate spheroid, thus making it difficult to formulate a static atmosphere lower boundary condition within our spherically symmetric version of the ZEUS MHD code. For very strong magnetic fields (i.e. $`\eta _{}>100`$), the large associated Alfven speed implies a small Courant time-step within the explicit time-stepping of the Zeus MHD code, and so such cases become very computationally expensive to evolve through the required few characteristic flow evolution times. Further work is thus needed to explore such cases of near-critical rotation and/or very strong magnetic confinement.

## 4. Magnetically Confined “Rigid-Body” Disks

In lieu of detailed simulations, it is however possible to infer some likely attributes of magnetically torqued outflows in the strong field limit, for which the magnetic lines-of-force can be idealized as “rigid pipes” that channel the wind outflow toward a collision near fixed loop tops. For the case of a rotation aligned dipole, this effectively concentrates material at the combined magnetic/rotational equator. Below the Keplerian radius, such material will again tend to fall back to the star, but above this radius, it will be centrifugally pushed against the loop top, gradually building in mass until it finally has a sufficient rotational energy to break open the field. By comparing the gravitational-rotational energy of this building disk to the magnetic field energy, one can derive an estimate for this “break-out” time. For typical stellar parameters, the ratio of this to a characteristic wind flow time $`\tau _{flow}=R_{}/V_{\mathrm{}}`$ scales roughly as

$$\frac{\tau _{break}}{\tau _{flow}}\frac{0.01\eta _{}W^{4/3}}{(r/R_K)^4(r/R_K)}.$$
(6)
Such a disk would thus have a roughly constant inner edge at the Keplerian radius $`R_K`$, but a highly variable outer edge. Indeed, having any sort of disk of accumulated material that extends much beyond $`R_K`$ requires that the numerator in eqn. (6) exceed unity. In that case, we can estimate a maximum outer cutoff radius to be where $`\tau _{break}\tau _{flow}`$, which solves to

$$R_{out,max}\frac{R_A}{\sqrt{10}W}.$$
(7)
A key point here, however, is that such disks would be quite distinct from a Keplerian disk, both physically and in terms of likely obervational signatures. In particular, they seem likely to be much more variable, characterized by intervals with substantial radial flow speeds in the outer disk. Perhaps even more notably, their rotation would follow a rigid-body law, $`v_\varphi r`$, instead of the Keplerian form, $`v_\varphi 1/\sqrt{r}`$.

In terms of a proposed application to Be-star disks, both properties seem to run counter to traditional interpretations of Be-disk emission line-profiles (Hanuschik 1995), which in many Be stars exhibit features (e.g. “central-quasi emissions”; Rivinius et al. 1999) that suggest a quite low upper limit (ca. 10-20 km/s) to radial outflow speeds, and moreover seem generally consistent with an azimuthal speed that follows a Keplerian law. The problems becomes particularly acute for extended rigid-body disks, since these imply more rapid rotation speeds than are typically inferred from Be emission-line-widths. For disks with a limited radial extent, the differences from Keplerian rotation are less dramatic, and as explored in the poster contribution by D. Telfer et al., the resulting profiles could indeed be consistent with existing observational analyses.

However, the long-term V/R variations seen in a substantial fraction of Be stars do still seem best explained by a Keplerian disk undergoing a precession of elliptical orbits that characterize an one-arm disk oscillation (Savonije and Heemskerk 1993; Telting et al. 1994; Savonije 1998). In a rigid-body disk in which the individual fluid elements of the disk are tied to the rotation period of the star, it is difficult to see how any processes could reproduce the year-long timescales of V/R asymmetries that require an associated long-term distinction in the physical properties of emitting material in a specific fixed-frame direction.

On the other hand, such rigid-body magnetic disks may represent a good paradigm for the observed circumstellar emission from Bp stars (see, e.g. talk write-up by Groote and references therein). The observed ca. 10 kG fields of these stars combined with their generally modest mass loss rates $`\dot{M}10^9M_{}`$/yr imply a huge magnetic confinement parameter, of order $`\eta _{}10^6`$. With the associated Alfven radii of $`R_A30R_{}`$, even the generally modest rotations of stars with, say, $`W0.1`$, would allow centrifugal disk support above $`R_K10^{2/3}R_{}4.6R_{}`$, extending perhaps out to $`R_{out,max}100R_{}`$. Such low-density winds are moreover subject to ion runaway instability (Krticka and Kubat 2000, 2001; Owocki and Puls 2002), and this likely plays an important role in the abundance variations often associated with the circumstellar emission from Bp stars. Future efforts should thus explore further such a wind-fed “magnetically rigid disk” (MRD) paradigm for Bp stars.

### Acknowledgments.

The research herein was supported in part by NASA grant NAG5-3530 and NSF grant AST-0097983 to the University of Delaware, and by an UK-PPARC Visitor Fellowship to SPO. SPO is grateful to J. Brown at University of Glasgow and A. Willis at University College London for their hopitality during his sabbatical visit to these institutions.

## References

Cassinelli, J. P., Brown, J. C., Maheswaran, M., Miller, N. A., and Telfer, D. C. 2002, ApJ, 578, 951.

Castor, J. I., Abbott, D. C., & Klein, R. I. 1975, ApJ, 195, 157 (CAK).

Hanushik, R. W. 1995, A&A, 295, 423.

Krtička, J. & Kubát, J. 2000, A&A, 359, 983.

Krtička, J. & Kubát, J. 2001, A&A, 369, 222.

Owocki, S. P. and Puls, J. 2002, ApJ, 568, 965.

Owocki, S. P. 2003, in “Stellar Rotation”, IAU Symposium 215, A. Maeder and P. Eenens, eds., in press.

Rivinius, T., Stefl, S., and Baade, D. 1999, A&A, 348, 831.

Savonije, G. J. 1998, in Cyclical Variability in Stellar Winds, L. Kapers and A. Fullerton, eds., Springer: Berlin, p. 337.

Savonije, G. J., and Heemskerk, M. H. M. 1993, A&A, 276, 409.

Telting, J. H., Heemskerk, M. H. M., Henrichs, H. F., and Savonije, G. J. 1994, A&A, 288, 558.

ud-Doula, A. 2002, Ph. D. Thesis, University of Delaware.

ud-Doula, A. and Owocki, S. P. 2002, ApJ, 576, 413.


text:
# Enhanced transmission and beaming of light via photonic crystal surface modes

## Abstract

Surface modes are generally believed to be an undesirable feature of finite photonic crystals (PC), unlike point or line defect modes. However, it is possible to make the surface mode radiate by appropriate corrugation of the PC interface. In this paper we show theoretically that the coherent action of these surface indentations can be engineered to collimate within a few degrees the light exiting a PC waveguide, or to funnel light coming from free space into the waveguide.

The proposal of photonic crystals (PC) Yablonovitch87 ; John87  in dielectric materials has stimulated the search of photonic band gap effects in periodic metallic structures McGurn93 ; Smith94 ; Kitson96 ; FJ96 ; Fan96 ; Sievenpiper98 ; Bozhe01 ; Fleming02 . These latter structures have shown remarkable properties, such as extraordinary optical transmission both through subwavelength hole arrays Ebbesen98  and single apertures Thio01 , and beaming of light from single apertures flanked by periodic surface corrugations Lezec03 . The root of all these phenomena lies in the existence of surface electromagnetic (EM) modes at metal-dielectric interfaces LMM01 ; LMM03 ; FJ03 . PCs can also support surface EM modes for appropriate crystal termination surfaces Meade91 ; winn94 ; Ramos99 . Additionally, for frequencies within the band gap, the PC behaves as a mirror for incoming EM radiation. These two facts suggest that the cited optical phenomena found in metals could have counterparts in PC systems. However, this analogy is not straightforward as the field penetration inside PCs (of the order of one wavelength) is much larger than the skin depth in good metals in the optical regime, which could strongly modify the interference phenomena we are referring to. In this paper we explore the possibility of exporting some of the results found in metallic structures to dielectric PCs. We will show that, through appropriate corrugation of the PC surface, it is possible to obtain both enhanced transmission through slits in PC slabs, and strong beaming of light coming out of a PC waveguide.

Let us first discuss the basics of the beaming phenomena in metallic systems, considering the simpler one-dimensional case of a guided mode coming out of one single slit flanked by a finite array of grooves LMM03 . The grooves act as centers that re-scatter the EM radiation diffracted by the slit exit. For certain resonant conditions a leaky surface EM wave develops, which emits into vacuum as it propagates along the surface. The phase accumulated during this surface propagation can be compensated, after emission, by optical paths in free space. Therefore, for certain angles, radiation coming from all indentations interferes constructively, leading to beaming. This mechanism suggests the strategy shown in this paper in order to pursue the analogy: first a waveguide in a PC is constructed and then, by modifying the PC-vacuum interface, a surface EM mode is created. Finally, in order to achieve the desired beaming effect, this mode is connected to the radiative modes by adding a periodic modulation to the PC surface. Additionally, these surface leaky modes could also be used to funnel light from free space into the waveguide, therefore enhancing the transmission, as has been previously demonstrated in single apertures surrounded by a periodic corrugation in a metallic film FJ03 .

In this paper, for proof of principle purposes, we concentrate in two-dimensional PCs. As photonic band gap (PBG) system, we consider a square lattice of cylinders in vacuum. We choose cylinder radius $`r=0.18a`$ ($`a`$ being the lattice parameter) and dielectric constant of the cylinders $`ϵ_\text{r}=11.56`$ (we take $`ϵ_\text{r}`$ for GaAs at a wavelength of $`1.5\mu m`$ as a representative value). For this set of parameters the system presents, for E-polarization (electric field pointing along the cylinders axis), a frequency gap for reduced frequencies $`\omega a/2\pi c`$ in the range $`(0.30,0.44)`$, where $`c`$ is the speed of light. All results presented are for E-polarization and have been obtained with the multiple multipole (MMP) method Hafner99 . The computational technique to treat PC waveguide discontinuities within the MMP method is thoroughly described in Refs. Moreno02a, ; Moreno02b, . A waveguide is constructed by removing one row of cylinders in the $`z`$-direction (see inset of Fig. 2). For most frequencies in the PBG this waveguide supports one propagating mode. The PC is infinitely extended towards $`z\mathrm{}`$ but, in order to study the coupling of the waveguide mode into free space, we terminate the PC in a set of cylinders parallel to the $`x`$-direction (in the coordinate system we have chosen, the center of these last cylinders is at $`z=7.5a`$). This regular PC termination does not support surface states winn94 . Fig. 1(a) shows the real space electric field amplitude when a guided mode with frequency $`\omega _0=0.408(2\pi c/a)`$, launched at $`z\mathrm{}`$, is diffracted by the waveguide exit. The radiated electric field shows a broad angular distribution, as corresponds to an aperture of dimensions comparable to the wavelength (the wavelength corresponding to $`\omega _0`$ is $`\lambda _0=2.45a`$). This is more clearly seen in Fig. 2 (black line), rendering the far-field radial component of the Poynting vector, $`S_r`$, as a function of azimuthal angle, $`\theta `$. In order to create a surface state, one possibility would be to cut the outer part of the interface cylinders, thus leaving one monolayer of hemicylinders winn94 . The origin of this surface state is that a mode that would be in the dielectric band for complete cylinders shifts its energy up when part of the dielectric is removed Meade91 . In this work the surface state is induced in a similar manner by decreasing the radius $`r_\text{i}`$ of the interface monolayer of cylinders, which is a more convenient system to treat by the MMP method. For the ulterior illustration of beaming in the forward direction, it is convenient to choose $`r_\text{i}`$ in such a way that a surface mode appears with a modal wavelength $`\lambda _{\text{sm}}=2a`$. For the chosen $`\omega _0`$ this occurs approximately for $`r_\text{i}=0.09a`$. Fig. 1(b) renders a similar plot to Fig. 1(a) but, this time, all cylinders in the interface monolayer have that particular $`r_\text{i}`$. The presence of a surface mode being fed by the waveguide is clear in Fig. 1(b), but as this mode does not couple to radiative modes, the radiated power again spreads in all directions as corresponds to an aperture size of the order of $`\lambda _0`$. This point is corroborated by Fig. 2 (red curve) showing an almost uniform $`S_r(\theta )`$ for the case considered. In order to alter the radiation pattern, the surface mode should couple to the continuum of radiative modes. This can be achieved by introducing a periodic modulation to the PC surface. In this way the surface state can “pick up” a reciprocal lattice vector of the superimposed lattice and connect to waves inside the light cone. Out of the different possibilities for this modulation, we present here results for a simple choice with a corrugation period $`\mathrm{\Lambda }=2a`$: at either side of the waveguide exit, the cylinders at the PC-vacuum interface in the first $`N`$ even positions (we order the cylinders consecutively according to their distance to the waveguide) are shifted a distance $`\mathrm{\Delta }z`$ along the $`z`$-axis. Fig. 1(c) renders the real space electric field amplitude for $`\mathrm{\Delta }z=0.3a`$ and $`N=9`$. This plot is radically different from the one obtained for systems without surface mode or with a non-radiative surface mode, distinctly showing the re-emission process. The beaming characteristics are clearly seen in the corresponding far-field transmission angular distribution shown in Fig. 2 for $`N=1,3,6,9`$. These curves present a strong redistribution of the transmitted radiation into the forward direction, with a beam angular width of a few degrees.

In the corrugated metal case, a focus in the Fresnel regime was found associated to the beaming in the far-field region FJAPL . This is also the situation in the present case. Fig. 3(a) shows the real space dependence of the Poynting vector modulus, $`S`$, for a modulated surface with $`N=3`$, where already the presence of a collimated beam and a focus is apparent. In similar way to the metallic counterpart, the focus characteristics depend on the number of periods of modulation in the surface. As illustration, Fig. 3(b) shows $`S_z(z,x=0)`$ for $`N=0,1,3,6,9`$. The maxima of these curves provide the focal distance of the, very elongated, focus. In order to further assign these effects (beaming and lensing) to the existence of a leaky surface mode, and to rule out that they are only due to the modification of the waveguide opening induced by the corrugation, the $`x`$-component of the Poynting vector, $`S_x`$, can be scrutinized along the interface. Fig. 4 shows $`S_x(x,z=7.5a)`$ for the cases of no surface mode, surface mode ($`N=0`$), and surface leaky mode for $`N=3,9`$. For $`N=0`$, after an initial dependence on $`x`$ for small distances to the waveguide (due to the coupling between the waveguide mode and the surface mode), the modal carried power remains constant, as expected for the non-radiative surface mode present in the system. On the contrary, cases with $`N0`$ show that the field decreases as the mode travels along the surface, due to the radiation losses responsible for the beaming and lensing phenomena. Notice that the case with a surface termination that does not support surface modes or resonances, presents an almost negligible intensity along the surface. It is worth stressing here that, although we have chosen a favorable range of parameters in order to show that the cited effects exist, no attempt has been made in order to optimize the structures. This is already apparent in Figs. 1(b, c), where the strong electric field modulation inside the waveguide is due to a high reflection \[$`R=52.3\%`$ and $`R=24.5\%`$ in (b) and (c), respectively\] caused by the large impedance mismatch between waveguide and free space regions. This reflection can be largely modified by very local changes of the dielectric structure at the opening (such as tapering Moreno02b ). Similarly, the coupling between waveguide mode and surface resonance surely has a strong influence on the transmission angular distributions. Which geometry is optimal is an interesting question that, however, falls outside the scope of this paper.

Beaming at other angles is also possible. For instance, for the chosen frequency $`\omega _0`$, shifting normally to the surface one every three cylinders ($`\mathrm{\Lambda }=3a`$), results in two collimated beams at $`\theta =\pm 24^{}`$. For the system presented in this paper, coupling to a surface state with an in-plane wave vector different from $`2\pi /\mathrm{\Lambda }`$ (which would occur at certain frequencies different from $`\omega _0`$) would also result in beaming at $`\theta 0^{}`$.

As previously said, another interesting effect found in metallic structures is that a structure presenting beaming when corrugated in the exit surface, is capable of presenting enhanced transmission when that corrugation is induced in the input surface. Fig. 5 shows the transmittance through a waveguide in a PC slab when illuminated by a normally incident plane wave. The black curve corresponds to a regular unstructured interface \[same as in Fig. 1(a)\], whereas the red curve corresponds to the corrugated interface shown in Fig. 1(c). The transmittance has been normalized to the power that impinges onto a length $`a`$ (i.e., roughly the waveguide cross section). In both curves rendered in Fig. 5, transmission outside the PBG is due to the EM field permeating through the crystal, while the region of transmission of order unity inside the gap is due to propagation through the waveguide. Interestingly, Fig. 5 shows a large resonant transmission inside the gap for the PC with a modulated surface. The resonant frequency is the same that was previously obtained for beaming of light out of a waveguide. Notice that, even in this case without full optimization of the geometry, the energy flux collected into the waveguide corresponds to the one impinging in a cross section around $`14a`$. This effective cross section could be increased, for instance, by considering a shallower corrugation with a larger lateral size. Unlike the analog situation in metallic structures, here the system is lossless, and therefore dissipation does not set a limit to the possible device lateral size.

To summarize, in this paper we have demonstrated the use of surface leaky modes, appearing in PCs under appropriate conditions, to collimate EM radiation out of PC waveguides. Very narrow beams with angular distributions of a few degrees can be achieved by modulating the PC exit surface. Furthermore, these surface resonances may also be used for injecting light efficiently into these waveguides. This enhanced transmission property could have interesting applications in the efficient coupling between conventional optical fibers and PC-based waveguides. (Funding by the Swiss National Science Foundation/NCCR Quantum Photonics and financial support by the Spanish MCyT under contracts MAT2002-01534 and MAT2002-00139 are gratefully acknowledged).