text:
# Fluctuations and Market Friction in Financial Trading

## Abstract

We study the relation between stock price changes and the difference in the number of sell and buy orders. Using a soft spin model, we describe the price impact of order imbalances and find an analogy to the fluctuation-dissipation theorem in physical systems. We empirically investigate fluctuations and market friction for a major US stock and find support for our model calculations.

The unpredictable up and down movements in the stock market have always captured the interest and imagination of investors. The scientific investigation of these phenomena started with Bachelier’s comparison of stock price dynamics with a random walk . This study has been refined in many respects . Interest has been devoted to the precise shape of the distribution of returns (difference of the logarithm of stock prices at different times), which is characterized by a high probability for large fluctuations .

After it was realized that large price fluctuations tend to cluster together in time, stock price returns were described by models of volatility (standard deviation of returns) changing in time . This effect is captured in time series models, in which the volatility at a given time depends on the magnitude of previous returns . It has been actively investigated how the stochastic properties of price dynamics can be related to the market microstructure, i.e. the rules and motivations according to which agents act in a financial market. Although the details of models differ  , they have been successful in reproducing the empirical observations.

Here, we follow a different approach motivated by the successful application of physics concepts to study the economy . We investigate the price dynamics on an intermediate or coarse grained time scale $`\mathrm{\Delta }t`$, which is long enough to average over the details of market microstructure, but short enough to resolve the trading dynamics. In physical systems, there exists a quantitative relation between the strength of fluctuations of an observable and the rate, at which energy is dissipated when the same observable is pulled by an external force. Einstein discovered this relation first, when he studied the fluctuations of small particles in a fluid . In this case, the variance of position fluctuations is proportional to temperature times the mobility of particles. This relation was later on generalized to a wide class of systems  and is known as fluctuation-dissipation theorem.

In this letter, we use an effective model for stock price dynamics to derive a similar relation between the time varying volatility of stock returns and the “friction constant for stock prices” or market liquidity. Market liquidity relates returns to the difference between buy and sell orders (order imbalance), which acts as an external force. We study the dependence of stock price changes on order imbalance empirically by using the method of data analysis and some of the results of . We find that the empirical results agree well with our model.

Model calculation: The observable quantity we are interested in is the logarithmic stock price changes within a time interval $`\mathrm{\Delta }t`$

$$G_{\mathrm{\Delta }t}(t)=\mathrm{ln}S(t)\mathrm{ln}S(t\mathrm{\Delta }t),$$
(1)
where $`S(t)`$ is the price of a given stock at time $`t`$. Transaction prices at a stock exchange lie usually in a finite interval between the bid price (the price a trader offers to pay for a stock) and the ask price (the price at which a dealer is willing to buy the stock). In addition, the historical prices studied take only discrete (tick) values. This motivates to model price changes by a spin model, and for the virtue of easier analytical calculations by a soft spin model .

The spin variable we use is the “instantaneous return” $`g(t)=\tau \frac{d}{dt}\mathrm{ln}S(t)`$, where the average time interval between trades $`\tau `$ sets the time scale of the problem. The observable returns on time scale $`\mathrm{\Delta }t`$ are related to them via

$$G(t)_{\mathrm{\Delta }t}=\frac{1}{\tau }_{t\mathrm{\Delta }t}^t𝑑t^{}g(t^{}).$$
(2)
In an analogous way, we define an instantaneous order imbalance $`q(t)=\tau \frac{d}{dt}Q(t)`$. Here, order imbalance is the number of buy orders minus the number of sell orders divided by the total number of outstanding shares of a given stock. We describe the dynamics of $`g(t)`$ by the following stochastic differential equation

$`\tau _tg(t)=rg(t)\kappa g^3(t)+\mu q(t)+\xi _i(t).`$ (3)
Mathematically, the $`r`$\- and $`\kappa `$\- term enforce the soft spin constraint of finite and discrete price changes. From a finance point of view, the linear term describes relaxation and controls the strength of fluctuations. The cubic term stabilizes the theory in a regime of strong fluctuations (small $`r`$). The order imbalance $`q(t)`$ is analogous to the magnetic field in the theory of spins, its coupling constant $`\mu `$ is in analogy to the magneton in magnetic systems. The Gaussian white noise $`\xi (t)`$ has the correlator $`\xi (t)\xi (t^{})=\tau \delta (tt^{})`$. It represents the influx of news not captured by the news content of the order imbalance. Stochastic differential equations for stock price changes with a quadratic instead of the cubic term were studied in . A generalization of Eq.(3) to many stocks was used to explain correlations between them .

The variance of price changes on a time scale $`\mathrm{\Delta }t\tau `$ is well approximated by

$`(G_i(\mathrm{\Delta }t))^2_0={\displaystyle \frac{2}{r^2}}{\displaystyle \frac{1}{\tau }}\mathrm{\Delta }t,`$ (4)
where the average is taken without the cubic term. This approximation is valid as long as the system is not in the vicinity of a critical point. Eq.(4) describes a diffusion process with diffusion constant $`D=\frac{2}{r^2\tau }`$. In order to describe the stochastic volatility and fat tails observed in empirical data, we must allow for a time dependence of both the frequency of trades $`\frac{1}{\tau }`$ and the standard deviation of the price of individual transactions. In Ref.  it was shown that price changes of financial transactions can be consistently interpreted in such a framework. When $`r`$ becomes a random variable, is causes multiplicative noise in Eq.(3), which can account for the power law tails of the return distribution .

The price impact function describes the functional relation between the expectation value $`G_{\mathrm{\Delta }t}_{Q_{\mathrm{\Delta }t}}`$ conditioned on the volume imbalance and the volume imbalance itself. We define the susceptibility $`\chi `$ as the slope of the price impact function close to zero volume imbalance

$`\chi =\underset{\mathrm{\Delta }Q0}{lim}{\displaystyle \frac{G_{\mathrm{\Delta }t}_{\mathrm{\Delta }Q}}{\mathrm{\Delta }Q}}.`$ (5)
As stocks are traded in discrete units, this limit can in an empirical study only be approximated. The inverse of this susceptibility can be interpreted as market liquidity. For an order imbalance which is constant in the time interval $`\mathrm{\Delta }t`$ the full price impact function can be calculated analytically for the model Eq.(3) and compares well to the empirical data in Fig. 1. Such a situation could be realized if an informed trader breaks up a large trade into many small ones and trades at a fixed frequency.

From Eq.(3) we find $`\chi =\mu /r`$. Combining this with Eq.(4) for the volatility, we obtain for the linear regime of small fluctuations

$`\sqrt{(G_{\mathrm{\Delta }t})^2}=\sqrt{2\mathrm{\Delta }t/\tau }\chi /\mu .`$ (6)
This relation and its empirical test are our central results. The square root on the r.h.s. is essentially the square root of the number of trades in the time interval $`\mathrm{\Delta }t`$ and in analogy to the temperature in physical systems. The relation Eq.(6) differs from the Einstein relation in that the standard deviation of fluctuations appears on the l.h.s. instead of the variance. The reason for this is the difference in “experimentally” relevant time scales: as $`\mathrm{\Delta }t`$ was assumed to be much larger than the microscopic time scale $`\tau `$ in our coarse grained description of stock price dynamics, $`G_{\mathrm{\Delta }t}(t)`$ is proportional to the zero frequency Fourier component of the instantaneous return, and not to the integral over all frequency components as in the derivation of the Einstein relation.

Empirical Results: The empirical test of Eq.(6) uses the method of data analysis and some of the results of . The present author got to know about the shape of the empirical price impact function and its similarity to the magnetization curve of spin systems from the authors of , and only later on learned about related work in the economics literature .-

In the framework of the efficient market hypothesis the price of a stock reflects present and discounted future earnings of the underlying company. The estimation of these earnings depends on information about that company. Information may be disseminated by a news release, a press article, or by rumors. In addition, traders extract information about a company (or at least information about the beliefs of other traders concerning that company) from the order imbalance. This is reasonable under the assumption that somebody selling the stock might do that because he has superior knowledge about the underlying company and expects its value to drop. Hence, one can assume that the order imbalance will have impact on the price of a stock.

We analyze stock prices of General Electric, a major American company for the year 1997. Price, volume, bid and ask price of all transactions of this stock are recorded in the Trades and Quotes (TAQ) data base published by the New York Stock exchange. We analyzed a total number of $`10^6`$ transactions contained in that data base. The question, whether a trade is buyer or seller initiated is difficult to answer, as there is a buyer and a seller for each trade. We try to answer it by finding out which one is more eager to do the trade, i.e. we compare the transaction price to the bid and ask price . First, we try to match the price of a trade with quotes at the same stock exchange with a time stamp at least two seconds prior to the transaction time. In this way we account for possible time lags in the reporting of prices. If the trade occurs at the bid price, it is classified as seller initiated, if it occurs between bid and ask price it is classified as neutral, and a trade at the ask price is classified as a buyer initiated. If the transaction price is lower than the bid price or higher than the ask price the trade is discarded as not matched.

Having classified all trades, we can give a sign to the volume of each trade. A buyer initiated trade has a positive sign, a seller initiated trade a negative sign. As our analysis does not focus on the market microstructure, we average both volume and price changes over a sampling time $`\mathrm{\Delta }t`$. This averaging is analogous to the concept of coarse graining in physics, which also averages over microscopic model details and allows the derivation of effective theories. In the finance context, the choice of the sampling time is influenced by the following considerations. If the sampling time is chosen too short, the effects of individual trades will not be smoothed out, if it is chosen too long it is not clear that there is still a causal relation between order imbalance and price change: price changes will influence traders and in this way have an effect on the order imbalance. For our analysis we choose $`\mathrm{\Delta }t`$ as five minutes. In this way we generate for each stock a time series of returns and order imbalances (measured relative to the total number of outstanding stocks of a company).

To visualize this information, we plot the average price change for a given order imbalance against the order imbalance $`\mathrm{\Delta }Q`$ (price impact function) in Fig.1. It agrees with the results of . Its characteristic form with a linear part for small $`\mathrm{\Delta }Q`$ and a less steep part for large $`\mathrm{\Delta }Q`$ is described by a hyperbolic tangent  and can be well approximated by the functional form resulting from Eq.(3).

Next, we sort the 5 minute intervals according to the average fluctuation ($`1/r`$ in terms of the Langevin model) per trade and calculate a price impact function for each fluctuation strength. We calculate the slope in its linear region and thus generate a series of susceptibilities $`\chi `$ depending on the local volatility. This data set is shown in Fig.2. The full line is a least square fit to the data points. The good quality of the fit is illustrated by the high regression coefficient $`R=0.975`$. Thus the empirical analysis supports well the linear relationship between fluctuation strength and susceptibility $`\chi `$ suggested by Eq.(6).

Discussion: The linear relation between fluctuation strength and liquidity implies that large price changes are not only due to a large volume of trades but also to a large price impact of a given volume. In the soft spin model Eq.(3), the time changing volatility is proportional to the square root of the number of trades times the liquidity. In Ref. it was shown that changes in the trading frequency cannot account for the appearance of large price fluctuations. Hence we conclude that large price fluctuations are caused by changes in the market liquidity.

The inverse susceptibility appears as the $`1/r`$-term in our model Eq.(3). For this reason, large price fluctuations associated with a large liquidity are described by the soft spin model in a critical state. In a critical spin model, even small changes in the external magnetic field cause huge fluctuations in magnetization. From this point of view, it is not large “external” influences which cause large price fluctuations, but the strong response of the system itself.

In summary, we have studied both theoretically and experimentally the relation between the strength of stock price fluctuations and the friction constant, which relates order imbalance and price changes. We have modeled stock price dynamics on a coarse grained time scale by a soft spin model and found an analogy to the fluctuation-dissipation theorem in physics, which is supported by our empirical study.

Acknowledgements: We acknowledge several stimulating discussions with V. Plerou, P. Gopikrishnan, X. Gabaix, and H.E. Stanley and thank them for making available the results of  prior to publication. We acknowledge helpful conversations with A. Kempf, T. Nattermann, D. Stauffer, and L. Viceira.


format: text

text:
# Discovery of X-ray emission from the proto-stellar jet L1551 IRS5 (HH 154)

## 1. Introduction

During the final stages of the formation of low-mass stars (in the so-called classical T Tau phase) accretion of material from the proto-stellar nebula onto the Young Stellar Object (YSO) takes place through an accretion disk. Very often the presence of the accretion disk is correlated with the presence of energetic polar outflows, that is, collimated jets of material being ejected perpendicularly to the disk, along its axis. When these jets collide with the surrounding ambient medium – or with previously ejected material – they form a shock structure, which is directly observable in the form of so-called Herbig-Haro jets.

X-ray emission (and thus the presence of hot plasma, at temperatures in excess of several $`\times 10^5K`$, up to $`100`$ MK during energetic flares) has by now been observed in most stages of the formation of low-mass stars, ranging from the highly embedded, perhaps spherically accreting proto-stars (Class 0 objects) to the final stages of the pre-main sequence life of a star, the Weak-Line T Tau stage, during which the X-ray luminosity is thought to come from a “normal” (however very active) stellar corona (e.g. Feigelson & Montmerle 1999).

While accretion itself has been considered as a possible source of X-ray emission in classical T Tau stars, up to now no evidence of energetic phenomena associated with proto-stellar jets has been observed. Here we present the first observations of X-ray emission from a proto-stellar jet, obtained in a well-studied system in which the proto-star (and its immediate circumstellar environment) powering the outflow is so heavily obscured that the jet can be singled out as the source of emission of the X-rays without ambiguity. Our observations show that this jet is indeed an X-ray source with a luminosity equivalent to a significant fraction of the X-ray luminosity normally associated with YSOs. The observed X-ray spectrum is compatible with a thermal origin of the observed X-ray emission. The associated temperature is moderate, well matched to the shock velocities observed in this and other Herbig-Haro jets. This raises the question of whether the X-ray emission associated with jets could indeed be a common feature of stellar formation, so that in some cases a significant fraction of the X-ray luminosity associated with the star (YSO/accretion disk) is actually emanating from shocks in the jet.

## 2. The L1551 IRS5 outflow

The L1551 cloud is one of the nearest ($`d140`$ pc) sites of ongoing star formation, in which objects in several different stages of the process are clearly visible, from deeply embedded, actively accreting (proto-)stars to the final stages of star formation represented by the Weak-Line T Tau stars with no remaining circumstellar material. In this paper we are concerned with the jet associated with the IRS5 source embedded in the L1551 cloud and its associated outflow.

L1551 IRS5 is a deeply embedded proto-stellar binary system in the L1551 cloud, and it is effectively invisible at optical wavelengths as it is hidden behind some $`150`$ mag of visual extinction (White et al. 2000) which most likely originates in the circumstellar accretion disk. The two Class 0/1 stars have a total luminosity of $`30L_{}`$, and appear to be (jointly?) powering at least two observable outflows, a large (several arcmin) bipolar molecular outflow and a much smaller (with a length of $`10`$ arcsec) denser two-component jet (Fridlund & Liseau 1998), consisting of material at temperatures of $`T10^4`$ K, thus visible in the emission lines of e.g. H$`\alpha `$. The jet moves at transverse velocity of 200–400 km s-1 (Fridlund & Liseau 1994) and appears to end in a shock against the ambient medium (a “working surface”)

## 3. XMM-Newton observations

The X-ray observations presented here were obtained with the XMM-Newton observatory. A deep (50 ks exposure) of the star-forming region of the L1551 cloud was obtained starting on Sep. 9 2000 at 19:10 UTC. All three EPIC cameras were active at the time of the observation, in full-frame mode, with the “medium” filters. Here we present the data obtained with the EPIC-PN camera.

An X-ray source (the leftmost of the three sources in Fig. 1) is clearly visible at the position of the IRS5 jet, with two additional sources visible in the region of the molecular outflow. For each of these sources source photons have been extracted from a circular region of 45 arcsec diameter, while background photons have been extracted from a region on the same CCD chip and at the same off-axis angle as for the source region. The spectral analysis has been performed using the XSPEC package, after rebinning the source spectra to a minimum of 20 source counts per (variable width) bin. The background-subtracted count rate for the X-ray source associated with IRS5 is $`8.4\times 10^4`$ cts s-1 in the EPIC PN camera, so that 42 source cts are collected the $`50`$ ks exposure, allowing only a limited amount of spectral information to be derived for the source. The resulting spectrum (shown in Fig. 2) is soft, and can be reasonably described with a moderately absorbed thermal spectrum. The best-fit column density ($`1.4\pm 0.4\times 10^{22}`$ cm-2) corresponds to an extinction of $`A_V=7.3\pm 2.1`$ mag. The best-fit temperature is $`T=0.5\pm 0.3\times 10^6`$ K.

The full width at half energy (FWHE) of the XMM point-spread function (PSF) for EPIC PN camera is $`15`$ arcsec, significantly larger than the size of the jets (whose visible length is $`10`$ arcsec). Thus, it is not possible to locate the precise site of the X-ray emission within the jet structure.

The two additional X-ray sources visible in Fig. 1 have spectra and absorbing column densities compatible with their being background active galactic nuclei shining through the molecular outflow, and are thus most likely unrelated with IRS5 and its jet/outflow structure.

## 4. Discussion

The IRS5 jet has a number of shocks along it, and it is observed to end in a “working surface” against the ambient medium at $`10`$ arcsec from the presumed location of the source powering it. The absorbing column density toward the jet is estimated at 4–6 mag, a value compatible with the X-ray measured column density, making the association between the X-ray emission and the jet highly plausible. Since the IRS5 proto-stellar system is hidden behind a very thick layer of absorbing material, corresponding to $`A_V>150`$ mag, it can be excluded that the X-ray photons – given the small absorbing column density and the lack of high-energy photons in the spectrum – emanate from (or close to) the photosphere/chromosphere of the proto-stars powering the jet. We therefore draw the conclusion that this source is the result of thermal emission in the shocks whose recombination light is seen along the jet in the visual wavelength regime.

As discussed in detail by Favata et al. (2001), the fluid velocity in the shock is $`270`$ km s-1, with an inclination, $`\mathrm{i}`$, of $`45`$ deg. The immediate post-shock temperature can then be estimated at $`T_{\mathrm{p}s}0.67`$ MK, fully compatible with the observed X-ray temperature of $`0.5\pm 0.3`$ MK. Thus, the observed X-ray emission is likely to be due to material heated at the interface shock (the working surface) between the jet and the ambient medium medium, or possibly in shocks along the cavity excavated by the jet. The X-ray luminosity of the emission associated with the jets is $`L_\mathrm{X}3\times 10^{29}`$ erg s-1 (assuming a distance of 140 pc for the L1551 complex).

While the jet’s X-ray luminosity is at the low end of the X-ray luminosity distribution of typical PMS stellar X-ray sources, the location of the emitting source some 10 AU *above* the star and thus the disk may result in a relevant effect of the jet’s X-ray emission on the conditions of the disk even in the presence of a higher stellar X-ray luminosity. The X-ray emission coming from the jet’s shock will illuminate the disk from above (and below), and can thus change the ionization of the disk’s material at large distances from the proto-star, in regions which are normally (unless the disk is very strongly flared) shielded from the influence of the stellar X-ray emission.

### 4.1. Energetics

The mass of the jet (see Favata et al. 2001 for the details) can be estimated at $`1`$$`2\times 10^6`$ M, which, with a shock velocity of $`200`$ km s-1, results in a mechanical luminosity of the jet of $`10^{41}`$$`10^{42}`$ erg s-1, so that a very low conversion efficiency between mechanical and radiant luminosity is sufficient to justify the observed X-ray luminosity from the shock. The H$`\alpha `$ luminosity is $`4.7\times 10^{28}`$ erg s-1, comparable with the X-ray luminosity derived here.

## 5. Conclusions

While the energetic nature of the collimated jets observed to be originating from proto-stellar sources has been evident for some time, no high-energy photons have up to now been observed from these phenomena. Here we report the first convincing evidence of X-ray emission from the proto-stellar jet associated with the IRS5 proto-star(s) in the L1551 cloud. The X-ray source and the proto-star and related jets are positionally coincident, and the small absorbing column density observed for the X-ray spectrum ($`7`$ mag, fully compatible with the absorbing column density observed towards the jet) allow us to exclude that the X-ray emission is associated with the proto-stellar sources (which are hidden behind $`150`$ mag of obscuration). The size of the jets ($`10`$ arcsec) originating at L1551 IRS5 is smaller in angular extent than the XMM EPIC PSF ($`15`$ arcsec), so that no inference is possible on spatial grounds about the possible detailed location of the origin of the X-ray emission.

The emission from the IRS5 jet is compatible with being caused by thermal emission from a plasma heated to a moderate temperature ($`T0.5`$ MK). This is equivalent to the shock temperature that is expected at the interface (“working surface”) between the jet and the surrounding circumstellar medium, on the basis of the observed jet velocity.

## References

Favata, F., Fridlund, C. V. M., Micela, G., Sciortino, S., Kaas, A. A. 2001, A&A, in press (astro-ph/0110112)

Feigelson, E. D. & Montmerle, T. 1999, ARA&A, 37, 363

Fridlund, C. V. M. & Liseau, R. 1994, A&A, 292, 631

Fridlund, C. V. M. & Liseau, R. 1998, ApJ, 499, L75

White, G. J., Liseau, R., Men’shchikov, A. B., et al. 2000, A&A, 364, 741


text:
# Constancy of the bilayer splitting as a function of doping in 𝐵𝑖₂𝑆𝑟₂𝐶𝑎𝐶𝑢₂𝑂_{8+𝛿}

## Abstract

Using high energy resolution angle resolved photoemission spectroscopy, we have resolved the bilayer splitting effect in a wide range of dopings of the bilayer cuprate $`Bi_2Sr_2CaCu_2O_{8+\delta }`$. This bilayer splitting is due to a nonvanishing intracell coupling $`t_{}`$, and contrary to expectations, it is not reduced in the underdoped materials. This has implications for understanding the increased c-axis confinement in underdoped materials.

It has long been thought that the high temperature superconducting cuprates can to first order be treated as two dimensional systems since the transport properties within the ab-plane ($`CuO_2`$ plane) and along the c-axis are very different. Such anisotropy shows up for example in the resistivity ratio $`\rho _c/\rho _{a,b}`$ and their temperature dependence $`\rho _{(a,b),c}(T)`$. Other probes such as the optical conductivity and interlayer tunneling measurements also indicate that the charge transport along the c-axis is severely reduced, especially in underdoped samples. This leads to the general notion that the physics of these materials is almost entirely two dimensional, and can be well described by a single $`CuO_2`$ plane. However, this picture seems to be incompatible with the fact that the transition temperature $`T_c`$ is closely related to the number of $`CuO_2`$ planes per unit cell, with single layer compounds of a family generically having lower $`T_c`$’s than bilayer or trilayer compounds. Additionally, recent experiments on the organic-chain intercalated $`Bi_2Sr_2CaCu_2O_{8+\delta }`$ (Bi2212) showed that isolated $`CuO_2`$ planes do not possess superconductivity. These results highlight the importance of the c-axis couplings for the superconductivity. Therefore, although they have been largely neglected so far, it appears clear that studies of the c-axis coupling effects on the electronic structure will be critical for a clear understanding of the high $`T_c`$ phenomenon.

A large advance was made recently by angle-resolved photoemission (ARPES) experiments on the overdoped bilayer cuprates $`Bi_2Sr_2CaCu_2O_{8+\delta }`$ (Bi2212) which showed for the first time the existence of the normal state c-axis intracell coupling. The coupling between adjacent $`CuO_2`$ planes within the bilayer, $`t_{}`$, was observed as a splitting of the near Fermi energy ($`E_F`$) bands into a Bonding (B) and Antibonding (A) set. This effect was found to be maximal at the $`(\pi ,0)`$ point of the zone and to go to zero along the $`(0,0)(\pi ,\pi )`$ nodal line, consistent with theoretical predictions. The extracted coupling value $`t_{}`$ of $`55meV`$ is reduced from band theory predictions by about a factor of 3, indicating sizable correlation effects.

Considering the highly unusual physics in the underdoped regime including in particular the pseudogap and the non-Fermi liquid behavior due to some sort of increased correlations, it is natural to wonder whether the correlations which reduce $`t_{}`$ are enhanced in the underdoped regime as well. General expectations are that they would, since (1) the c-axis conductivity (DC and optical conductivity) is greatly reduced in the underdoped samples, indicating greater two dimensional confinement, (2) correlations are expected to be stronger as the parent Mott insulator is approached, and (3) the bilayer splitting has not yet been directly observed in the underdoped regime. In this letter, we present the first direct evidence that bilayer splitting exists generically in the Bi2212 family irrespective of doping, and in fact, the energy splitting is essentially constant as a function of doping.

We performed high resolution ARPES on $`Bi_2Sr_2CaCu_2O_{8+\delta }`$ over a wide range of doping from heavily overdoped (OD, $`T_c=55K`$), through optimal (OpD, $`T_c=91K`$) and underdoped (UD, $`T_c=78K`$), where the doping levels corresponding to different oxygen concentrations were achieved by high temperature anneals . The experiments were carried out at the Advanced Light Source (ALS), Berkeley using the HERS Scienta endstation at the undulator BL10.0.1. Experimental details can be found elsewhere. An important aspect of the experiments was the choice of photon energy, which we found adjusts the relative intensities of the bonding and antibonding bands, allowing an easier experimental deconvolution of the various features.

Figure 1 shows experimental spectra along the $`(\pi ,0)(\pi ,\pi )`$ line (green line of the Brillouin zone in panel (g), which also shows a schematic of the bonding FS (red) and antibonding FS (black)). False color images of the energy and momentum spectra along this cut are shown in panels (a)-(f) for the three dopings and for two photon energies ($`22eV`$, $`47eV`$). The bilayer-split bands labeled A and B are seen very clearly in panel (a) (OD, $`22eV`$), in the (red) energy distribution curve (EDC) in panel (h) and (red) momentum distribution curve (MDC) in panel (i). In the OpD and UD samples however, the intrinsic linewidth of the individual features is increased to a value equivalent to or even greater than the bilayer splitting energy, effectively masking the splitting.

Additional insight can be obtained by using the complementary photon energy of $`47eV`$, which we empirically found enhances the antibonding band A relative to the bonding band B (see panels (b), (h) and (i)). A recent first-principles calculation of the photoemission matrix elements has independently shown this same effect. Aside from the matrix element effect changing the intensity ratio of the two bands, we find that the lineshape and energy/momentum position of the A band is to within experimental uncertainty unchanged. This large matrix element modulation allows us to deconvolve the contribution from each of the bilayer-split bands in the optimal and underdoped samples, where the broadening of the individual features is relatively large.

Similar to what we have observed in panels (a) and (b), the $`47eV`$ images of (d) and (f) show a larger concentration of spectral weight near $`E_F`$ than do the $`22eV`$ images of (c) and (e). This is seen more clearly in the EDCs at $`(\pi ,0)`$ shown in (j) and (l) which show a sharper feature at $`47eV`$(blue) than at $`22eV`$(red). We argue that this is because the $`47eV`$ data shows significant weight from the A band, while the 22 eV data additionally shows weight from the B band. The MDCs at $`E_F`$ (panel (i), (k) and (m)) tell a similar story - the $`47eV`$ data is more concentrated at $`k_y=0.0`$ corresponding to the A band, while the $`22eV`$ data extends farther out in $`\stackrel{}{k}`$ space because it contains weight from the B bands as well.

Figure 2 shows the Fermi Surfaces plots extracted from the data, including the effects of superstructure bands and some photon energy dependences. Panels (a) -(f) show the the locus of normalized near-$`E_F`$ ( in $`\pm 10meV`$ window) spectral weight throughout much of the Brillouin zone for the same OD and OpD samples and photon energies as in figure 1. The data were taken in the normal state ($`T=80K`$ for OD and $`100K`$ for OpD and UD) and the $`\stackrel{}{k}`$ space location is indicated by the green box in panel (g). A few of the additional superstructure bands possible are also indicated in panel (g). Modulo matrix element effects which may be significant, the high intensity locus of this data should represent the Fermi Surface topology.

As in figure 1, the bilayer splitting is most clear in panel (a) of figure 2 (OD, $`22eV`$) where we see the A (black) and B bands (red) plus some of their superstructure replicas (light black and red lines). Also similar to figure 1, the $`47eV`$ data emphasizes the A band and so panel (b) (OD, $`47eV`$) almost exclusively shows the A band with a great amount of weight near $`(\pi ,0)`$ (in fact the matrix elements at $`47eV`$ favor the second zone somewhat, causing the high weight portion to be slightly off of the $`(\pi ,0)`$ point). Panels (d), (f) show the FS topology very similar to that in (b), i.e. they all show the A FS topology. Panels (c) and (e) at $`22eV`$ are dominated by the standard hole-like FS corresponding to the B band.

The above data is relevant for the issue of the Fermi Surface topology of Bi2212, which has been an intensely debated issue recently. In particular, there has been a discussion about whether the Fermi Surface may sometimes (depending upon photon energy) exhibit a more electron-like character , or whether there is a single unchanging hole-like piece plus superstructures, with no effects other than matrix elements which vary as a function of photon energy. Figures 1 and 2 show that both FS topologies can be seen on any one sample depending upon the photon energy used, and this can be clearly understood to be simply due to the bilayer splitting effect. We believe that this should close the debate of the possibility of two types of Fermi Surfaces in Bi2212- they both exist. Whether one, the other, or both are observed depends upon the matrix elements which naturally vary with photon energies.

We note that the current data can not explicitly determine which side of $`(\pi ,0)`$ the antibonding band crosses, as it approaches exceedingly close to this point (this is the reason for the shading in panel (g) of figures 1 and 2). In other words, the van Hove singularity of the antibonding band is very near the Fermi level at this point. Typical notion is that if the band crosses $`E_F`$ along $`(0,0)(\pi ,0)`$ line, the van Hove singularity is above the $`E_F`$ whereas if it crosses along $`(\pi ,0)(\pi ,\pi )`$ then it is below $`E_F`$. In general it had been known that the van Hove singularity was near $`E_F`$, with a well known parameterization being $`34meV`$ below $`E_F`$. While this is quite near $`E_F`$, it is still much larger than the temperature scale of $`T_c`$. On the other hand, the current data indicates that for the overdoped sample the antibonding band is at binding energy $`9\pm 2meV`$ at the $`(\pi ,0)`$ point, which is comparable to the temperature scale of the problem. This implies it should be quite active electronically. The van Hove singularity of the standard “hole-like” Fermi Surface (bonding component) is on the other hand about $`100meV`$ below $`E_F`$, making it much less relevant. This is important, as there have been many theoretical suggestions that a van Hove singularity near $`E_F`$ may be critical to understand the physics of these compounds .

The existence of the bilayer splitting and its effect on the Fermi Surface topology is qualitatively new and was obtained directly from the raw data with no fitting or modeling required. This is the ideal most trustworthy case. There is much new quantitative information as well, and this can be extracted through a simple fitting of the data. Near $`E_F`$ where the dispersion is steep we fit each set of MDCs with two Lorentzian peaks and a linear background (see for example figures 3(a)-(d)). Near the bottom of the bands where the dispersion is flat we fit the EDCs with two Lorentzian peaks plus a linear background multiplied by a Fermi cutoff (figure 3(e)-(g)). Usually we can only fit a portion of the spectra (positive $`k_y`$ portion in (a), negative $`k_y`$ in (b),(c),(d)) because of contamination by superstructure bands. The easiest data to fit is $`22eV`$ OD, which shows both A and B bands clearly and well separated (figure 3(a) and 3(e)). However, data from the optimal and underdoped samples can also be well fit if we tune to $`47eV`$ and use two bands A and B. It is much more difficult or impossible to fit these spectra with one peak.

The overall fitting results are compiled in figures 3(h) - 3(k), each set of which corresponds to about 20 fits of the type in panels (a)-(g). These were each taken along a cut about $`27\%`$ of the zone away from $`(\pi ,0)`$ where the superstructure contamination was minimal ($`k_x=0.73\pi `$ for $`22eV`$ data and $`1.27\pi `$ for $`47eV`$ data). The open circles are the EDC peak centroids and the closed circles are the MDC peak centroids. The dashed lines are quadratic fits to the centroids, implying a simple parabolic band dispersion. The fact that both the MDC and EDC centroids are very well fit by the same parabolic bands gives high confidence for the quality of the data and fits for all sample types.

Panels (i) - (k) show the parabolic fits from all three sample types measured under identical conditions ($`47eV`$, $`k_x=1.27\pi `$ cut, normal state). While the different doping levels give an energy offset of the parabolas from one sample to another ($`30meV`$ from OD to OpD and another $`10meV`$ from OpD to UD) due to the band filling effect, the shape and energy splitting of the parabolic curves are essentially unchanged with doping. In practice, we found the shape and splitting ($`70meV`$) of the parabolic curves from the OD sample and then overlaid this parabola (shifted in energy but with the same splitting) with the data for the other samples. The good agreement is clear from the fits.

The bilayer splitting as a function of k has been theoretically parameterized as $`\mathrm{\Delta }(\stackrel{}{k})=0.5t_{}(cos(k_xa)cos(k_ya))^2`$, which is zero along the node lines and maximal at the $`(\pi ,0)`$ points. Using this formula to extrapolate the current data to $`(\pi ,0)`$ we get a maximum splitting of $`114meV\pm 8meV`$ or a $`t_{}57meV\pm 4meV`$, obtained from the $`47eV`$ data. The $`22eV`$ data shown here (panel (h)) gives an estimate of $`54meV`$ also. This is in excellent agreement with our previous estimate of $`55meV`$ obtained on $`25eV`$ data. Thus we conclude that $`t_{}`$ is independent of the incident photon energies and has at the most a very weak dependence on the doping level.

The conclusion that $`t_{}`$ depends weakly on the hole doping levels is quite surprising since the normal state resistivities show very different behavior, both in the magnitude and the temperature dependence. For example, in the UD regime, $`\rho _c`$ shows insulating behavior while $`\rho _{ab}`$ is still metallic. In contrast, both $`\rho _{ab}`$ and $`\rho _c`$ are metallic in the overdoped regime. This type of data has been widely taken to indicate an increased two-dimensional confinement in the underdoped samples, a phenomenology that would be consistent with other exotic behaviors such as spin-charge separation or marginal Fermi Liquid behavior. While coupling between the unit cells is also relevant for the confinement, this coupling should be much weaker than the bilayer coupling discussed here, and in fact this coupling has not yet been directly observed. Our data directly implies a lack of increase in two dimensional confinement from the strong intracell $`t_{}`$ and hints at a lack of change from the intercell $`t_{}`$ effect. We suggest that other factors such as the increase in the scattering rate near $`(\pi ,0)`$, the onset of the pseudogap, or carrier fractionization should be more responsible for these anomalies in the transport behavior.

We acknowledge sample preparation help from J. Koralek and M. Varney and beamline support from S. Kellar, X.J. Zhou, P. Bogdanov, and Z. Hussain. This work was supported by the NSF Career-DMR-9985492 and the DOE DE-FG03-00ER45809. ALS is operated by the DOE, Office of Basic Energy Sciences.


text:
# Singular Density of States of Disordered Dirac Fermions in the Chiral Models

## Abstract

The Dirac fermion in the random chiral models is studied which includes the random gauge field model and the random hopping model. We focus on a connection between continuum and lattice models to give a clear perspective for the random chiral models. Two distinct structures of density of states (DoS) around zero energy, one is a power-law dependence on energy in the intermediate energy range and the other is a diverging one at zero energy, are revealed by an extensive numerical study for large systems up to $`250\times 250`$. For the random hopping model, our finding of the diverging DoS within very narrow energy range reconciles previous inconsistencies between the lattice and the continuum models.

72.15.Rn,71.23.-k,05.30.Fk

Dimensionalities and symmetries play central roles for universalities in the Anderson localization problem. More than two decades ago, Abrahams, Anderson, Licciardello and Ramakrishnan  presented the well-known scaling theory, predicting that electron wave functions always localize in one and two dimensions and metal-insulator transition occurs in three dimensions. However, even in two dimensions, delocalized states are also marginally allowed to appear when systems possesses some special symmetries. For symmetries, it is often convenient to borrow terminologies from the random matrix theory which was first introduced by Wigner and Dyson. Recently, Altland and Zirnbauer  have reported seven new symmetry classes in conjunction with mathematical classification scheme of the Riemannian symmetric spaces.

Amongst these new symmetry classes, chiral models have attracted much attention as a novel exception for the scaling theory. Compared to the conventional models where randomness enters as an on-site potential, it resides on links  (i.e., as a gauge field) in chiral models. This type of randomness may play important roles for the composite fermion theory of fractional quantum Hall effects and vortex states of dirty superconductors. In view of a localization problem, models with this randomness have a special symmetry, referred as chiral symmetry, and thereby belong to a new universality class . This symmetry is expected to affects localization properties of the systems drastically.

Several chiral models on a lattice, which are convenient for numerics, have been studied and interesting physics around zero-energy have been revealed. DoS of these models show singularities at zero-energy and the corresponding wave functions exhibit a delocalized multifractal behavior. Examples of these include the Gade’s model , the random flux model , and the $`\pi `$-flux model with link disorders . These models, defined on a two dimensional square lattice, have the chiral symmetry which is conveniently stated as $`\{H,\gamma \}=0`$ where $`\gamma `$ is a matrix which changes the sign of wave functions on a sub lattice. Consequently, for any realization of the disorder, the energy spectrum is invariant under the transformation $`EE`$. Therefore, given an eigenstate $`\psi `$ with energy $`E`$, $`\gamma \psi `$ is also an eigenstate with energy $`E`$. This symmetry is responsible for the existence of delocalized states at zero-energy.

Chiral models in a continuum space have also been investigated extensively. Especially, models which include the Dirac fermions have attracted an interest. The Dirac fermion is a quasi particle which appears in several interesting situations in condensed matter physics such as $`d`$-wave superconductivity, graphite sheets, the gap-closing transition in quantum Hall effects, the Chalker-Coddington network model, and the mean field theory of the $`tJ`$ model. Effects of randomness are also of fundamental interest in these contexts. Since many analytic approaches such as a field theoretic one are applicable for the Dirac fermion, several interesting results have been obtained so far. Amongst them, there exist simple models which allow us to construct an explicit zero-mode wave function for any realization of disorder . Due to this advantage, it has been revealed that the zero-mode wave functions are not localized and exhibit a multifractal behavior like ones in the lattice models. Moreover, there exists a transition in the multifractal spectrum  and the density-density correlation  as the disorder strength varied.

Although the delocalized multifractal nature of the exact zero-energy states has been well established now, our knowledge for finite energy states, especially for DoS around zero-energy, is still in confusion. For example, a continuum model where the Dirac fermions feel random gauge fields was studied in Ref.  and it was found that DoS exhibits a power-law dependence on energy. On the other hand, a similar model with species doubling was studied in Ref. and diverging DoS was found. From the lattice point of view, the corresponding model including Dirac fermions has not been studied so far. However, the random flux model, a lattice model where flux is randomly distributed for each plaquette, shows diverging DoS . A clear link between them is missing.

Another example is the case where the Dirac fermion feels imaginary vector potentials. For this case, there even exist inconsistent results between lattice and continuum models. Field theoretic studies  predicts diverging DoS at zero-energy for any randomness strength. However, a numerical result for the random $`\pi `$-flux model, where a random hopping amplitude act as a imaginary vector potential, shows that DoS behaves as a power-law with its exponent dependent on the disorder strength for weak randomness .

In this paper, we make an attempt to clarify the relationship between lattice models and continuum Dirac fermions, reconcile some inconsistent results in the previous studies, and thereby give a whole perspective for the localization problem of the chiral models. For these purposes, we consider lattice models which recover the random Dirac fermions in the continuum limit. We focus on DoS, especially around zero energy where quantum interference play an important role and fully quantum mechanical treatments are necessary. We use the transfer-matrix method developed in Ref.. It allows us to handle large enough systems up to $`250\times 250`$ , which is indispensable to give reliable results for localization problems.

We realize Dirac fermions on a 2D square lattice via the $`\pi `$-flux model:

$`_{pure}={\displaystyle \underset{<ij>}{}}c_i^{}t_{ij}^{pure}c_j+h.c.`$ (1)
where $`t_{j+\widehat{x},j}^{pure}=()^{j_x}`$ , $`t_{j+\widehat{y},j}^{pure}=1`$ and flux piercing a plaquette is $`\pi `$. The energy spectrum is given by $`E=\pm 2\sqrt{\mathrm{cos}^2k_xa+\mathrm{cos}^2k_ya}`$ where $`𝐤`$ belongs to the magnetic Brillouin zone $`(\pi /a,\pi /a]\times (0,\pi /a]`$ and $`a`$ is a lattice constant. In the continuum limit $`a0`$, this includes the doubled massless Dirac fermions around $`𝐤=(\pi /2a,\pi /2a)`$ and $`(\pi /2a,\pi /2a)`$. This realization is a minimum model for our purposes.

In the following, we will consider two types of disorder which live on a link: the random gauge field and the random hopping. We implement the random gauge field as

$`t_{ij}^{rg}=t_{ij}^{pure}\mathrm{exp}[iaA_{ij}].`$ (2)
Taking a Coulomb gauge, we determine the random gauge field $`A_{ij}`$ via a scalar potential $`\mathrm{\Phi }`$ on a dual lattice as $`A_{j+\widehat{x},j}=\left(\mathrm{\Phi }_{j+(\widehat{x}+\widehat{y})/2}\mathrm{\Phi }_{j+(\widehat{x}\widehat{y})/2}\right)/a,A_{j+\widehat{y},j}=\left(\mathrm{\Phi }_{j+(\widehat{x}+\widehat{y})/2}\mathrm{\Phi }_{j(\widehat{x}\widehat{y})/2}\right)/a`$ where $`\mathrm{\Phi }`$ is randomly chosen from a simple Gaussian distribution $`P[\mathrm{\Phi }]\mathrm{exp}[\frac{a^2}{2g}_{<ij>}\left(\mathrm{\Phi }_i\mathrm{\Phi }_j\right)^2/a^2]`$. In the continuum limit, this reduces to $`P[\mathrm{\Phi }]\mathrm{exp}[\frac{1}{2g}d^2x\left(\mathrm{\Phi }\right)^2]`$, a natural choice for effective field theoretic treatments. As stated above, this Hamiltonian possesses the chiral symmetry $`\{,\gamma \}=0`$ due to the special nature of randomness. Note also that time-reversal invariance is broken by the random magnetic flux.

In the continuum limit, the Hamiltonian is expressed as

$`_{cont}={\displaystyle d^2x\mathrm{\Psi }^{}\left(\begin{array}{cc}0& D\\ D^{}& 0\end{array}\right)\mathrm{\Psi }}`$ (3)
where $`D:=i_x\sigma _xi_y\sigma _y+A_x(𝐱)\sigma _x+A_y(𝐱)\sigma _y+M(𝐱)\sigma _z+V(𝐱)\mathrm{𝟏}`$ , $`\sigma _{i=x,y,z}`$ are the $`2\times 2`$ Pauli matrices, and $`\mathrm{\Psi }`$ is a four-component spinor. The coefficients $`A_x,A_y,M`$ and $`V`$ are arbitrary complex fields. This is the most general chiral symmetric form of a Dirac Hamiltonian in a continuum space. The chiral symmetry is, in the present basis, expressed as $`\{,\sigma _z\mathrm{𝟏}\}=0`$. For the random gauge field model, $`(A_x,A_y)`$ is real and serves as the random gauge field. The other coefficients $`M`$ and $`V`$ do not appear at the first order in the lattice constant $`a`$, but they are non-zero for the lattice models in general.

For strong enough randomness, this model is naturally related to another lattice model, the random flux model. (By the random flux model, we mean the model discussed in Ref. which does not include Dirac fermions, while we simply call the model with Dirac fermions defined above as the random gauge field model.) To see this, note that the spatial correlation of the flux per plaquette $`\varphi _i:=a^2\mathrm{\Delta }_i\mathrm{\Phi }_i`$ is expressed as $`\varphi _i\varphi _j=a^4\mathrm{\Delta }_i\mathrm{\Delta }_jG(|ij|)`$ where $`G(|ij|):=\mathrm{\Phi }_i\mathrm{\Phi }_j`$ is a lattice Green’s function for $`\mathrm{\Phi }`$, and $`\mathrm{\Delta }_i`$ is lattice Laplacian at a site $`i`$. Taking a continuum limit $`a0`$, we find $`\varphi _i\varphi _j=0`$ for finite $`a|ij|`$, i.e., the spatial correlation of flux is short ranged. Moreover, for strong enough randomness, lattice fermions do not remember original flux $`\pi `$ anymore, and we expect this model to exhibit a similar behavior to that of the random flux model where the flux for each plaquette is an independent random variable.

An effective sigma model via supersymmetry (SUSY) technique was studied in Ref. , predicting diverging DoS for all range of randomness strength. A similar model without species doubling $`=d^2x\mathrm{\Psi }^{}\left[\sigma 𝐩+\sigma 𝐀\right]\mathrm{\Psi }`$ has also been studied  and power-law DoS with a disorder-dependent exponent was proposed. On the other hand, the random flux model, which is expected to be the strong randomness limit of the present case, shows a diverging DoS .

We also consider the case where fermions on a 2D square lattice with $`\pi `$-flux per plaquette feel random hopping amplitudes:

$`t_{ij}^{rh}=t_{ij}^{pure}+\delta t_{ij}`$ (4)
where $`\delta t_{ij}`$ is a real random variable. As for a probability distribution of $`\delta t`$, we consider a Gaussian distribution $`P[\delta t]\mathrm{exp}[(\delta t)^2/2g^2]`$.

In this case, $`(A_x,A_y)`$ in Eq.(3) is purely imaginary, and Dirac fermions feel an imaginary random gauge potential. Note that the total four-component Hamiltonian is hermitian. Unlike the random gauge field model, time-reversal symmetry is not broken for this model since we can take all matrix elements to be real.

As is the case of the random gauge field model, this model is connected to another lattice model for sufficiently strong randomness. For strong enough randomness, again, one expects that fermions do not remember the original flux anymore and shows a similar behavior to that of the Gade’s model.

Fukui  studied the continuum model by a replica nonlinear sigma model with a large number of fermion flavors and proposed that DoS at zero-energy diverges and delocalized states exist at the band center for any randomness strength. One-loop renormalization group (RG) study with SUSY method was also applied for this model in Ref. , predicting diverging DoS at zero energy as $`\rho (E)\frac{E_R}{E}e^{c\sqrt{\mathrm{ln}(E_R/E)}}`$ where $`E_R`$ is a constant and $`c`$ depends on randomness. A conjectural RG flow beyond one-loop order was presented in Ref.  and similar results were obtained. The Gade’s model, which we expect to be the strong randomness limit of the present model, also shows diverging DoS . On the other hand, according to a numerical study for the corresponding lattice model, vanishing DoS $`\rho (E)|E|^\alpha `$ was observed with a disorder-dependent exponent $`\alpha >0`$ for weak randomness . We re-examined this problem by the transfer-matrix method.

A quantity of interest we investigate numerically is random averaged DoS: $`\rho (E):=\frac{1}{L^2}_i\delta (E_iE)`$. For the pure $`\pi `$-flux model, DoS vanishes linearly around zero-energy due to the relativistic dispersion of Dirac fermions. Since for the strong enough disorder, the present models are expected to exhibit similar properties to those of the random flux model or the Gade’s model, both of which exhibit divergent DoS at zero-energy, a natural question is how the vanishing DoS at zero-energy becomes divergent as we increase the disorder strength $`g`$. It should be contrasted to the usual case where DoS of pure systems are also divergent at zero-energy due to the van-Hove singularities of the tight binding model on a two-dimensional square lattice.

First, we discuss the random gauge field model, i.e., $`\pi `$-flux model with a random gauge field as in Eq. (2). In Fig. 1 and 2, we present numerically calculated DoS for several $`g`$. In the transfer matrix method, the Green’s function for a given energy $`E`$ is calculated at $`E+i\delta `$. We chose $`\delta 0.010.0005`$ which gives us enough resolution of energy. We used a square geometry rather than a quasi one-dimensional one since we are mainly interested in pure two-dimensional properties.

As shown, there exist two structures. For over the wide range of the energy scale, the global feature of DoS changes from V-shaped to flat one, the latter is characteristic to the random flux model as expected. In addition, if we look more precisely, for very tiny region around zero-energy, there exists another structure, i.e., diverging DoS for sufficiently large $`g`$. Note that since we used a lattice with even number of sites and adopted periodic boundary condition, there is no exact zero energy state, i.e., the diverging DoS found here is not an artifact by the special choice of boundary conditions. The numerical results found here are consistent with the sigma model study via SUSY technique in Ref. , where the divergence of DoS is predicted.

Now let us go on to the $`\pi `$-flux model with random hopping amplitudes ( Eq.(4) ). Calculated DoS for systems up to $`250\times 250`$ are shown in Fig.3 and 4. Again, one recognizes two different structures as in the random gauge field model. Away from zero energy, DoS exhibits a power law behavior with a disorder-dependent exponent as in Ref. . The whole DoS profile becomes divergent as the disorder strength increased, which is natural since the present model is expected to show similar behaviors to that of the Gade’s model. On the other hand, as is the case with the random gauge field model studied above, within a very narrow range near zero-energy, there exists another structure. As the disorder strength increased, DoS divergence at zero-energy occurs when $`g`$ reaches about the order of the bandwidth. This behavior is independent of the global power-law profile away from the zero energy and before it turns to be divergent. These findings are consistent with the field theoretic analyses, especially with the SUSY approach where the divergence of DoS is predicted and a power-law behavior appears in the intermediate region of the RG flow. We also investigated the case where $`\delta t`$ is uniformly distributed between $`[w/2,w/2]`$. The results for this case are qualitatively similar to those of a Gaussian distribution.

We should compare the present case to dirty $`d`$-wave superconductors which also include Dirac fermions, however, the underlying symmetry class is different in the random matrix theory . For dirty $`d`$-wave superconductors, all quasiparticle states are found to be localized, and a small energy scale around zero energy appearers where quantum interference effects produce several critical DoS profiles depending on the details of randomness . This energy scale is determined by the diffusion constant in the diffusive regime and the localization length of quasiparticles . For the present case, on the other hand, zero-energy states are delocalized due to the chiral symmetry, and quantum interference effects give rise to diverging DoS at zero energy. The emergence of a power-law behavior away from zero energy is, from the effective field theory point of view, well described by the intermediate regime of the RG flow .

In conclusion, we have investigated DoS of the chiral models where the chiral symmetry plays crucial roles for the existent of the delocalized states at zero energy. We considered the lattice counter parts of continuum models including the Dirac fermions, the random gauge field model and the random hopping model, which have been well studied by field theoretic methods. Large-scale calculations by the transfer matrix method up to $`250\times 250`$ lattices revealed the existence of two distinct structures in DoS, one is a power-law behavior in the intermediate energy range and the other is a diverging one around zero energy. Quantum interference plays crucial roles for these fine structures, which naive semiclassical treatments may miss. Our finding of the diverging DoS at zero-energy for the random hopping model reconciles inconsistencies between lattice and continuum models and thus establishes a clear connection between them.

We thank Y. Morita for fruitful discussions. The computation in this work has been partly done at the YITP Computing Facility and at the Supercomputing Center, ISSP, University of Tokyo.


text:
# References

## Abstract

Black hole solutions with nonspherical event horizon topology are shown to exist in an Einstein-Yang-Mills theory with negative cosmological constant. The main characteristics of the solutions are presented and differences with respect to the spherically symmetric case are studied. The stability of these configurations is also addressed.

NEW HAIRY BLACK HOLES WITH NEGATIVE COSMOLOGICAL CONSTANT

J.J. van der Bij and Eugen Radu111corresponding author:
Albert-Ludwigs-Universität Freiburg, Fakultät für Physik,
Hermann-Herder-Straße 3, D-79104 Freiburg, Germany
email: radu@newton.physik.uni-freiburg.de
Telephone: +49-761/203-7630
Fax: +49-761/203-5967

Albert-Ludwigs-Universität
Fakultät für Physik
Freiburg Germany

PACS: 04.40; 04.70.-s; 98.80.E
Although anti-de Sitter (AdS) spacetime does not seem to correspond to the world in which we live , its importance has been noticed in many occasions. As shown by Hawking and Page, the presence of a negative comological constant $`\mathrm{\Lambda }`$ makes it possible for a black hole to reach stable thermal equilibrium with a heat bath .

The black holes discussed by Hawking and Page have a spherically symmetric event horizon. The topological structure of the event horizon of a black hole is an intriguing subject in black hole physics. When asymptotic flatness and the energy conditions are given up, there are no fundamental reasons to forbid the existence of black holes with nontrivial topologies. In particular, for a negative cosmological constant, black holes for which the topology of the horizon is an arbitrary genus Riemann surface have been considered by many authors (see $`e.g.`$ -). The thermodynamics of these solutions has been discussed in - while higher dimensional generalizations were obtained in -. The theorems about spherical horizon topology do not apply here because the negative cosmological constant can be interpreted as a negative vacuum energy density. They generalize the known solutions replacing the round two-sphere by a two-dimensional space $`\mathrm{\Sigma }`$ of negative or vanishing curvature.

These black holes embedded in ’locally AdS’ background spacetimes (background locally isometric to spacetimes of constant negative curvature) have been seminal to recent developments in black hole physics.

All these investigations are mainly based on Einstein(-Maxwell) theory (note however the inclusion of a dilaton field in ). It may be of interest to generalize these solutions for a nonabelian matter content. Here we restrict attention to $`SU(2)`$ Einstein-Yang-Mills (EYM) theory and investigate black hole solutions with locally flat or hyperbolic horizons, which asymptotically approach a locally AdS spacetime.

The case of a spherically symmetric horizon has been discussed in  with some surprising results. For example, there are finite energy solutions for several intervals of the shooting parameter (the value of gauge function at event horizon), rather than discrete values. Solutions exist for all values of $`\mathrm{\Lambda }<0`$. There are also stable solutions in which the gauge field has no zeros, and the asymptotic value of the gauge potential is not fixed.

This behavior is drastically different from those observed for EYM black holes in asymptotically flat  or de Sitter spacetime . Their regular conterparts are discussed in  and present also interesting properties.

In this letter we show numerically that, in the EYM system, the topology of the event horizon is not restricted to the spherically symmetric case. The family of black hole solutions considered in  is extended to geometries with locally flat and hyperbolic horizons.

We consider spacetimes whose metric can be written locally in the form

$$ds^2=\frac{dr^2}{H(r)}+r^2d\mathrm{\Omega }_k^2\frac{H(r)}{p^2(r)}dt^2$$
(1)
where

$$H(r)=k\frac{2m(r)}{r}\frac{\mathrm{\Lambda }r^2}{3}.$$
(2)
Here $`d\mathrm{\Omega }_k^2=d\theta ^2+f^2(\theta )d\phi ^2`$ is the metric on a two-dimensional surface $`\mathrm{\Sigma }`$ of constant curvature $`2k`$. $`r`$ is the radial coordinate for which $`r\mathrm{}`$ defines the asymptotic region. The discrete parameter $`k`$ takes the values $`1,0`$ and $`1`$ and implies the form of the function $`f(\theta )`$

$$f(\theta )=\{\begin{array}{cc}\mathrm{sin}\theta ,\hfill & \mathrm{for}k=1\hfill \\ \theta ,\hfill & \mathrm{for}k=0\hfill \\ \mathrm{sinh}\theta ,\hfill & \mathrm{for}k=1.\hfill \end{array}$$
(3)
In case of vanishing $`\mathrm{\Lambda }`$, these metrics describe black holes only for spherically symmetric solutions. The topology of spacetime is $`R^2\times H_g^2`$, where $`H_g^2`$ is the topology of the surface $`\mathrm{\Sigma }`$. When $`k=+1`$, the universe takes on the familiar spherically symmetric form, and the $`(\theta ,\phi )`$ sector has constant positive curvature. The cosmological constant may assume either sign. When $`k=0`$, the $`\mathrm{\Sigma }`$ is a flat surface, with $`\mathrm{\Lambda }`$ less than zero. When $`k=1`$, $`\mathrm{\Lambda }<0`$ and the $`(\theta ,\phi )`$ sector is a space with constant negative curvature, also known as a hyperbolic plane. A discussion of these cases is available in .

When $`\mathrm{\Sigma }`$ is closed, we denote its area by $`V`$.

The EYM coupled system is described by the action

$$S=d^4x\sqrt{g}[\frac{1}{16\pi G}(2\mathrm{\Lambda })\frac{1}{4}F_{\mu \nu }^aF^{a\mu \nu }],$$
(4)
The most general expression for the appropriate $`SU(2)`$ connection is obtained by using the standard rule for calculating the gauge potentials for any spacetime group . Taking into account the symmetries of the line element (1) we find

$$A=\frac{1}{2e}\{u(r,t)\tau _3dt+\nu (r,t)\tau _3dr+\left(\omega (r,t)\tau _1+\stackrel{~}{\omega }(r,t)\tau _2\right)d\theta +\left(\frac{d\mathrm{ln}f}{d\theta }\tau _3+\omega (r,t)\tau _2\stackrel{~}{\omega }(r,t)\tau _1\right)fd\phi \},$$
(5)
where $`\tau _a`$ are the Pauli spin matrices. To simplify the general picture we set $`u=0`$ $`i.e.`$ no dyons. However, nontrivial solutions exist also for a nonvanishing electric part (see  for a discussion of black hole dyon solutions in the $`k=1`$ case). The existence of dyon solutions without a Higgs field is a new feature for AdS spacetime; if $`\mathrm{\Lambda }0`$ the electric part of the gauge fields is forbidden .

For static configurations, it is convenient to take the $`\nu =0`$ gauge and eliminate $`\stackrel{~}{\omega }`$ by using a residual gauge freedom. The remaining function $`\omega `$ depends only on the coordinate $`r`$.

As a result, we obtain a simplified YM curvature $`F=dAieAA`$

$$F=\frac{1}{2e}\left(\omega ^{}\tau _1drd\theta +f\omega ^{}\tau _2drd\phi +(w^2k)f\tau _3d\theta d\phi \right)$$
(6)
where a prime denotes a derivative with respect to $`r`$. For this purely magnetic YM ansatz, we obtain the equations

$`m^{}`$ $`=`$ $`\omega ^2H+{\displaystyle \frac{(\omega ^2k)^2}{2r^2}},`$ (7)
$`\omega ^{\prime \prime }`$ $`=`$ $`{\displaystyle \frac{\omega ^{}}{H}}\left({\displaystyle \frac{2\mathrm{\Lambda }r}{3}}{\displaystyle \frac{2m}{r^2}}+{\displaystyle \frac{(\omega ^2k)^2}{r^3}}\right)+{\displaystyle \frac{\omega (\omega ^2k)}{r^2H}},,`$ (8)
where we have used dimensionless variables rescaled as $`r(\sqrt{4\pi G}/e)r,\mathrm{\Lambda }(e^2/4\pi G)\mathrm{\Lambda }`$ and $`m(eG/\sqrt{4\pi G})m`$.
The equation for $`p`$ decouples from the rest. Given a solution to (8), $`p`$ is obtained from

$`p=\mathrm{exp}\left(2{\displaystyle _r^{\mathrm{}}}𝑑r{\displaystyle \frac{\omega ^2}{r}}\right).`$ (9)
For black hole solutions having a regular event horizon at $`r=r_h>0`$,

$`m(r_h)={\displaystyle \frac{r_h}{2}}\left(k{\displaystyle \frac{\mathrm{\Lambda }r_h^2}{3}}\right).`$ (10)
We find the following expansion near the event horizon

$`m(r)`$ $`=`$ $`m(r_h)+m^{}(r_h)(rr_h),`$ (11)
$`\omega (r)`$ $`=`$ $`\omega _h+\omega ^{}(r_h)(rr_h),`$ (12)
where

$`m^{}(r_h)`$ $`=`$ $`{\displaystyle \frac{(\omega _h^2k)^2}{2r_h^2}},`$ (13)
$`\omega ^{}(r_h)`$ $`=`$ $`{\displaystyle \frac{r_h\omega _h(\omega _h^2k)}{(k\mathrm{\Lambda }r_h^2)r_h^2(\omega _h^2k)^2}}.`$ (14)
The condition for a regular event horizon is

$`{\displaystyle \frac{d}{dr}}\left(k{\displaystyle \frac{2m}{r}}{\displaystyle \frac{\mathrm{\Lambda }r^2}{3}}\right)|_{r_h}>0.`$ (15)
For given ($`r_h,\mathrm{\Lambda }`$), this places a bound on $`\omega _h`$

$`2m^{}(r_h)={\displaystyle \frac{\left(\omega _h^2k\right)^2}{r_h^2}}<k\mathrm{\Lambda }r_h^2,`$ (16)
and implies the positiveness of the quantity $`\omega ^{}(r_h)`$. In the $`k=1`$ case, (16) implies the existence of a minimal value of $`|\mathrm{\Lambda }|`$, $`i.e.`$ for a given $`r_h`$

$`|\mathrm{\Lambda }|>{\displaystyle \frac{1}{r_h^2}}(1+{\displaystyle \frac{1}{r_h^2}}).`$ (17)
At $`r\mathrm{}`$, the spacetime is locally isometric to AdS spacetime. In this limit, the field variables have the following asymptotic forms

$`m(r)`$ $`=`$ $`M+{\displaystyle \frac{M_1}{r}}+O\left({\displaystyle \frac{1}{r^2}}\right),`$ (18)
$`\omega (r)`$ $`=`$ $`\omega _{\mathrm{}}+{\displaystyle \frac{C_1}{r}}+O\left({\displaystyle \frac{1}{r^2}}\right),`$ (19)
where $`\omega _{\mathrm{}}`$, $`M`$ and $`C_1`$ are parameters and $`M_1=\frac{\mathrm{\Lambda }C_1^2}{3}\frac{1}{2}(\omega _{\mathrm{}}^2k)^2.`$ For $`k=1`$, $`m(r)`$ is the mass contained inside the radius $`r`$, while $`M`$ is the Arnowitt-Deser-Misner (ADM) mass. In the general case, $`M`$ is related to mass of the hole $`𝐌_{ADM}`$ (appropriately generalized for a negative cosmological constant ) by $`𝐌_{ADM}=MV/4\pi `$. However, for $`k=1`$ this identification is contentious. At issue is the zero of mass-energy, since even in pure Einstein gravity (with negative $`\mathrm{\Lambda }`$) one has a black hole solution with $`M<0`$. Surprisingly, this negative mass black hole can also be formed by gravitational collapse . Vanzo  gives explicit mass formulæ both for the solutions in an EM theory and for more general stationary black holes with the same asymptotic behaviour. In Ref. , the shifting of the zero of energy is advocated so that the total energy is proportional to $`MM_0`$, assuming there is suitable shift $`M_0`$ that works not only for the special solutions considered by Vanzo, but for all reasonable solutions with the same asymptotic behaviour.

The black hole solutions discussed in this paper contain a YM magnetic charge

$`𝐐_𝐌={\displaystyle \frac{1}{4\pi }}{\displaystyle 𝑑S_k\sqrt{g}\stackrel{~}{F}^{kt}}=Q_M{\displaystyle \frac{V}{4\pi }}{\displaystyle \frac{\tau _3}{2}},`$ (20)
where $`Q_M=(k\omega _{\mathrm{}}^2)`$. We remark that all $`k=1`$ black holes posses a nonzero magnetic charge. $`𝐌_{ADM}`$ and $`𝐐_𝐌`$ turn out to be finite if $`\mathrm{\Sigma }`$ is closed.

We numerically solve Eqs. (7) and(8) with boundary conditions (18), (19) using a standard shooting method. As mentioned above, the case $`k=1`$ is discussed in Ref. . We shall therefore focus on the cases $`k=0`$ and $`k=1`$.

Using the initial conditions on the event horizon (11), (12), the equations were integrated for a range of values of $`(r_h,\mathrm{\Lambda })`$ and varying $`\omega _h`$. Since the equations (7), (8) are invariant under the transformation $`\omega \omega `$, only values of $`\omega _h>0`$ are considered. In the $`k=1`$ case, for $`\omega _h=0`$, $`\omega ^{}(r_h)=0`$ and we have the corresponding EM solutions discussed in . For sufficiently small $`\omega _h`$, all field variables remain close to their values for the abelian configuration with the same $`r_h`$. Significant differences occur for large enough values of $`\omega _h`$ and the effect of the nonabelian field on the geometry becomes more and more pronounced. A similar behavior is noticed for $`k=0`$; here $`\omega _h=0`$ corresponds to a vacuum solution.

The solutions we obtain have many common properties with their spherically symmetric $`\mathrm{\Lambda }<0`$ counterparts. Also, most of the analytical results found in  for $`k=1`$ can easily be generalized for the cases $`k=0,1`$. As expected, the black holes properties are rather different from the corresponding ones for asymptotically flat and de Sitter spacetimes.

Solutions are classified by the parameters $`M`$ and $`Q_M`$. The behavior of metric functions $`m`$ and $`p`$ is similar to the spherically symmetric black hole solutions. For every considered value of $`(\mathrm{\Lambda },r_h)`$, we find regular black hole solutions for only one interval $`0<\omega _h<\omega _h^c`$. For $`\omega _h>\omega _h^c`$ solutions blow up or the function $`H(r)`$ becomes negative. The value of $`\omega _h^c`$ increases as $`|\mathrm{\Lambda }|`$ increases. There are also solutions for which $`\omega _{\mathrm{}}>1`$ although $`\omega _h<1`$.

In contrast to the spherically symmetric case, we find only nodeless solutions. This can be analytically proven by integrating the equation for $`\omega `$, $`(H\omega ^{}/p)^{}=\omega (\omega ^2k)/pr^2`$ between $`r_h`$ and $`r`$; thus we obtain $`\omega ^{}>0`$ for every $`r>r_h`$. For $`k=1`$, both nodeless solutions and solutions where $`\omega `$ crosses the axis can exist.

For $`k=0`$ we have always $`m>0`$ and the black holes therefore only occur with positive values of mass. Typical solutions in this case are presented in Figure 1. In Figure 2, $`M`$ is plotted as a function of $`Q_M`$ for $`k=0`$ black holes and several values of $`\mathrm{\Lambda }`$.

Not unexpectedly, when $`\mathrm{\Sigma }`$ has a negative curvature, we have found for suitable values of ($`\mathrm{\Lambda },r_h,\omega _h`$), a class of solutions with $`M<0`$ or even $`M=0`$, apart from configurations with positive $`M`$ (see Figure 3). Typical monopole spectra for $`k=1`$ case are presented in Figure 4.

Black hole solutions in asymptotically Minkowski and de Sitter space, which necessarily have at least one node in $`w(r)`$, are unstable . In contrast, the spherically symmetric monopole black hole solutions in the asymptotically AdS spacetime with no node in $`\omega (r)`$ are stable ,,.

We shall show that there are also black hole solutions with nonspherical event horizon stable against linear fluctuations (this is likely since we have found only nodeless solutions).

In examining time-dependent fluctuations around black hole solutions it is convenient to work in the $`u(r,t)=0`$ gauge. All field variables are written as the sum of the static equilibrium solution whose stability we are investigating and a time dependent perturbation. By following the standard methods well-known from the spherically symmetric case, we derive linearized equations for $`\delta w(r,t)`$, $`\delta \stackrel{~}{w}(r,t)`$, $`\delta \nu (r,t)`$, $`\delta p(r,t)`$ and $`\delta m(r,t)`$. The fluctuations decouple in two groups: $`\delta w(r,t)`$, $`\delta p(r,t)`$ and $`\delta m(r,t)`$ form the gravitational sector, whereas $`\delta \stackrel{~}{w}(r,t)`$ and $`\delta \nu (r,t)`$ form matter (sphaleronic) sector. The linearized equations imply that $`\delta p(r,t)`$ and $`\delta m(r,t)`$ are determined by $`\delta w(r,t)`$, and $`\delta \stackrel{~}{w}(r,t)`$ by $`\delta \nu (r,t)`$. For a harmonic time dependence $`e^{i\mathrm{\Omega }t}`$, the linearized system of the matter sector implies a standard Schrödinger equation

$`\left\{{\displaystyle \frac{d^2}{d\rho ^2}}+U_\beta (\rho )\right\}\beta `$ $`=`$ $`\mathrm{\Omega }^2\beta ,`$ (21)
where $`\beta (r,t)=r^2p\delta \nu /\omega =e^{i\mathrm{\Omega }t}\beta (r)`$ and a new tortoise radial coordinate is introduced $`d\rho /dr=p/H`$.
Here

$`U_\beta `$ $`=`$ $`{\displaystyle \frac{H}{r^2p^2}}(\omega ^2+k)+{\displaystyle \frac{2}{\omega ^2}}({\displaystyle \frac{d\omega }{d\rho }})^2,`$ (22)
with $`U_\beta (r_h)=0`$. Since $`\omega `$ has no nodes, this potential is regular everywhere. For $`k=0`$, $`U_\beta >0`$ and the potential is greater than the lower of its two asymptotic values. Standard results from quantum mechanics  imply that there are no negative eigenvalues for $`\mathrm{\Omega }^2`$ and no unstable modes. The same condition is satisfied by those $`k=1`$ solutions with $`\omega _{\mathrm{}}>1`$.
Similarly, the gravitational sector perturbation satisfies the equation

$`\left\{{\displaystyle \frac{d^2}{d\rho ^2}}+U_\omega (\rho )\right\}\delta \omega `$ $`=`$ $`\mathrm{\Omega }^2\delta \omega ,`$ (23)
with the potential

$`U_\omega `$ $`=`$ $`{\displaystyle \frac{H}{r^2p^2}}(3\omega ^2k)+4{\displaystyle \frac{d}{d\rho }}\left[{\displaystyle \frac{p}{rH}}({\displaystyle \frac{d\omega }{d\rho }})\right]`$ (24)
or, by using the equations (7)-(8)

$`U_\omega `$ $`=`$ $`{\displaystyle \frac{H}{r^2p^2}}\left[3\omega ^2k+{\displaystyle \frac{8\omega \omega ^{}}{r}}(\omega ^2k)+4\mathrm{\Lambda }r^2\omega ^24k\omega ^2+{\displaystyle \frac{4\omega ^2(\omega ^2k)^2}{r^2}}\right].`$ (25)
This potential is not positive definite, but is regular in the entire range $`\mathrm{}<\rho <\mathrm{}`$. Near the event horizon, $`U_\omega 0`$; at infinity $`U_\omega `$ takes always a positive value. It follows that eq. (23) will have no bound states if the potential $`U_\omega `$ is everywhere greater than the lower of its two asymptotic values $`i.e.U_\omega >0`$. We have found numerically $`U_\omega >0`$ for all $`k=0`$ considered configurations. This condition holds also for all $`k=1`$ solutions stable in the sphaleronic sectors (with $`\omega _{\mathrm{}}>1`$). Thus, these solutions are linearly stable in both sectors of the theory.

This fact can also be proven more rigorously. The analytical arguments presented in  for the existence of linearly stable black hole solutions can be directly generalized for every $`k`$. The only negative term in (25) is $`4\mathrm{\Lambda }r^2\omega ^2`$. Following  we can show that this term can be made as small as we like by taking $`|\mathrm{\Lambda }|`$ sufficiently large. Thus we see that a nonspherical topology of the event horizon has a tendency to stabilize the unperturbed solution.

The thermodynamic stability of these black holes will not be considered here and we will restrict ourselves on deriving the Hawking temperature. This can easily be done by using the standard Euclidean method. For the line element (1), if we treat $`t`$ as complex, then its imaginary part is a coordinate for a nonsingular Euclidean submanifold iff it is periodic with period

$`\beta ={\displaystyle \frac{4\pi r_hp(r_h)}{k2m^{}(r_h)\mathrm{\Lambda }r_h^2}}.`$ (26)
Then continuous Euclidean Green functions must have this period, so by standard arguments the Hawking temperature is (with $`k_B=\mathrm{}=1`$)

$`T_H={\displaystyle \frac{\left(k2m^{}(r_h)\mathrm{\Lambda }r_h^2\right)}{4\pi r_hp(r_h)}}{\displaystyle \frac{k\mathrm{\Lambda }r_h^2}{4\pi r_h}}.`$ (27)
This is a general result, for every matter content satisfying the weak energy condition. Thus the Hawking temperature of such systems appears to be suppressed relative to that of a vacuum black hole of equal horizon area. This relation generalizes for a negative cosmological constant the result found in  for a vanishing $`\mathrm{\Lambda }`$. Note also that, for given $`k`$, $`T_H`$ is independent of the topology of two-space $`\mathrm{\Sigma }`$.

In Figure 5 we plot the temperature $`T_H`$ $`vs.`$ the parameter $`M`$ for black holes with fixed ($`\mathrm{\Lambda },r_h`$) and varying the parameter $`\omega _h`$. On these curves, the maximal values of $`T_H`$ correspond to $`\omega _h=0`$.

The configurations we have found broaden the family of black holes with nonabelian matter content and manifest some new features. Also, all known hairy EYM black hole solutions have particle-like, regular counterparts. However, there are no regular counterparts of the solutions discussed in this work. This fact may be of some relevance when generalizing the isolated horizon formalism for a negative cosmological constant . We conjecture the existence of nontrivial black hole solutions with nonspherical event horizon for a more general nonabelian matter content $`e.g.`$ when including a Higgs field.
Acknowledgement
We gratefully thank Dr. B. Tausk for a careful reading of this manuscript.
This work was performed in the context of the Graduiertenkolleg of the Deutsche Forschungsgemeinschaft (DFG): Nichtlineare Differentialgleichungen: Modellierung,Theorie, Numerik, Visualisierung.

Figure Captions

Figure 1: The gauge function $`\omega `$ and the mass $`m`$ are shown as a function of the radial coordinate $`r`$ for three different solutions; here $`k=0`$, $`r_h=1`$.

Figure 2: Black hole monopole spectrum for $`k=0`$ black hole solutions with $`r_h=1`$ and different values of $`\mathrm{\Lambda }`$.

Figure 3: The gauge function $`\omega `$ and the mass $`m`$ are shown as a function of the radial coordinate $`r`$ for three different solutions; here $`k=1`$, $`r_h=5`$.

Figure 4: Black hole monopole spectrum for $`k=1`$ black hole solutions with $`r_h=5`$, and different values of $`\mathrm{\Lambda }`$. The existence of solutions with negative values of $`M`$ for small values of parameter $`\omega _h`$ is explicit.

Figure 5: The variation of the Hawking temperature $`T_H`$ with $`M`$ for topological black holes as a function of the parameter $`\omega _h`$. The radius of event horizon radius $`r_h`$ is fixed and several values of cosmological constant are considered. On these curves, the maximal values of $`T_H`$ correspond to $`\omega _h=0`$.

Figure 1

Figure 2.

Figure 3.

Figure 4.

Figure 5a.

Figure 5b.


text:
# Large scale star formation in galaxies. II. The spirals NGC 3377A, NGC 3507 and NGC 4394.

## 1 Introduction

Global star formation in a spiral galaxy is influenced by several factors. Different authors have investigated the dependence of the star formation rate and efficiency on environment. Star formation in interacting galaxies is generally enhanced, while it is not well assessed its dependence upon galaxy morphological type (Young et. al 1996 and references therein). Some recent papers have shown the importance of the determination of empirical correlations between large scale star forming regions and parent galaxy properties in understanding the star formation process. In this framework Elmegreen et al. (1994, 1996) and Elmegreen & Salzer (1999), analyzing different samples of spirals, have found a clear correlation between the size, the total B magnitude of a galaxy and the size of its largest stellar complex. These authors suggest, indeed, that these scaling laws can be related to the properties of the interstellar medium and of the instability process itself. However, it is worth to note that Selman & Melnick (2000) have recently explained that the correlation between the size of the largest superassociation and the luminosity of the parent galaxy, found by Elmegreen et al. (1994), is most likely due to a size-of-sample effect.

The problem of the identifications of regions of star formation in distant, unresolved, galaxies is quite difficult, because the identification suffers from several biases, like those introduced by differences in observational data and/or in the identification criterion adopted. For instance, Hodge (1986) clearly discussed these difficulties when trying to derive general properties of star forming regions.

To overcome these complications we developed an automatic method for the identification of the star forming regions (Young Star Groupings, or YSGs hereafter) in unresolved galaxies (Adanti et al. 1994). The availability of various sets of colours and fluxes allows every pixel of the galaxy image to be represented as a point in the space of the variables (i.e. of fluxes and colours). Principal Component Analysis (PCA) and Cluster Analysis (CA) result in an artificial image of the galaxy, where pixels are grouped into classes according to their relative distances in the space of the variables. As in Battinelli et al. (2000, hereafter Paper I), we performed our classification using U, U-B, B-V and B-R as variables. A first step to a deeper understanding of the link between the YSG’s and parent galaxy’s properties is to build a homogeneous and ‘objective’ database. Target galaxies were chosen according to the criteria discussed in Paper I.

In this paper we present data and results of YSGs in three spiral galaxies: NGC 3337A, NGC 3507 and NGC 4394. In Sect. 2 observational data acquisition and reduction are described; in Sect. 3 we present and discuss our identifications in each galaxy. Finally, in Sect. 4, we discuss some properties of the sample of YSGs in all the six galaxies we have studied so far.

## 2 Observations and data reduction

All data were obtained with the SPIcam camera at the 3.5 m telescope of the Apache Point Observatory (New Mexico, USA) 111APO is privately owned and operated by the Astrophysical Research Consortium (ARC), consisting of the University of Chicago, the Institute for Advanced Study, Johns Hopkins University, New Mexico State University, Princeton University, the University of Washington and Washington State University., in four different broadbands: U, B, V and R. SPIcam is a backside illuminated thinned SITe CCD, with a size of $`1024\times 1024`$ pixels, with an angular scale of $`0.28\mathrm{}`$/pixel and with a field of view of about $`4.8\mathrm{}\times 4.8\mathrm{}`$.

Three exposures per filter were taken for each galaxy. The CCDPROC routine in IRAF222IRAF is a free software distributed by National Optical Astronomy Observatories, which are operated by the Association of Universities for Research in Astronomy Inc., under contract to the National Science Foundation. was used to subtract the bias, as measured in the overscan region of each image, as well as to divide by the appropriate flat field for that night and filter. The spatial offsets among the images of a particular galaxy were obtained using the positions of stars within each field. The images were shifted to a common reference frame using the IMALIGN routine in IRAF. Finally, the IRAF task IMCOMBINE was used to average together the three exposures in each filter while at the same time rejecting cosmic rays.

NGC 3377A was observed on the night of 1998 November 21, a photometric night, while the NGC 3507 and NGC 4394 images were taken on 1998 December 21 and on 17-18 January 1999 respectively. These last three nights were not photometric, so we obtained shallower exposures of NGC 3507 and of NGC 4394 on a clear night, 16 April 1999, which we used to calibrate these images. The ratio in each filter between the averaged images of the non-photometric nights and the shallower images taken on the last photometric night was determined using stars within the field, if present, or using the radial profile of the galaxy itself. Total exposure times are always 1200 s in the U images and 600 s in the B, V and R bands, with the exception of the U image of NGC 4394, whose exposure time is 2000 s. The seeing of each image is reported in Table 1.

As mentioned in Paper I, the flat fields (one for each night and for each band) may be contaminated by scattered light from the lack of baffling of the 3.5m APO telescope, leading to an artificial spatial gradient in the sky background. This is most evident in the B and R band images than in the U and B ones. The one exception is NGC 4394 where some of the gradient in the background is due to the presence of the galaxy M85, just off the field of view.

Calibration was computed using observation of several standard stars from Landolt (1992). An aperture magnitude for each observation of each standard star was obtained with an aperture radius of $`7\mathrm{}`$, chosen to match the aperture used by Landolt. Standards included 3 stars in the field of PG 1323-086 and 5 stars in the field of PG 1633+099, with color index in the range from -0.9 to 1.1 for U-B, from -0.2 to 1.1 for B-V and from -0.1 to 0.6 for V-R. The standards were mostly observed at low airmass, less than 1.2, with the entire range being from 1.1 to 2.2.

## 3 Results

### 3.1 NGC 3377A

NGC 3377A is an almost face-on dwarf spiral galaxy, about $`7\mathrm{}`$ north-west from the E6 giant elliptical NGC 3377, in the Leo I (M96) group. The radial velocity difference of 120 km/s (de Vaucouleurs et al. 1991, hereafter RC3) between these two galaxies suggests that NGC 3377A may be a companion of the giant elliptical. According to Tonry et al. (1997), NGC 3377A is about 10.7 Mpc distant, implying a linear separation from NGC 3377 of only 20 kpc. Sandage et al. (1991) and Knezek et al. (1999) showed that this very low surface brightness galaxy has various peculiarities. In particular, in spite of its high MHII/LB ratio ($`0.30`$ in solar units, which is a typical value for a Sc type spiral), it has a very low star formation rate, about 0.003 M/yr, so that star formation can continue for more than a Hubble time (Knezek et al. 1999).

The application of our algorithm to this galaxy led us to the identification of 83 YSGs, whose positions in the galaxy are shown in Fig. 1 and whose main characteristic parameters are given in Table 2. We remind that, in order to reduce the contamination by random groups in this list of YSGs candidates, we adopted (as in Paper I) the procedure described by Battinelli & Demers (1992). This procedure results into the introduction of a threshold, Nlim, in the number of pixels in each groups, such that the contamination by random groups in the sample of all the groups composed by at least Nlim pixels is less than 10$`\%`$. The computed Nlim are 9, 4, 7 for NGC 3377A, NGC 3507 and NGC 4394, respectively.  Blue magnitudes and surface brightnesses have been corrected for a total extinction of 0.06 mag (0.04 mag galactic extinction and 0.02 mag internal extinction) as given in RC3. As explained in Paper I (Sect. 3), a precise determination of integrated colours of YSGs results hard since they often lie in areas with a lot of ”structure”.

The size distribution of the YSGs is shown in Fig. 4 for the three galaxies studied in this paper. For NGC 3377A the average size is 87 pc with standard deviation 59 pc. Fitting the high-size tail of the size distribution with the power law $`dND_{xy}^\alpha dD_{xy}`$ we find $`\alpha =2.3\pm 0.3`$ (correlation coefficient $`r=0.93`$. The shaded strips shown in Fig. 4 are the Dxy intervals directly effected by the introduction of the Nlim threshold for each galaxy. Such shaded areas are therefore certainly uncomplete.

Knezek et al. (1999) give a multi-color map of NGC 3377A where several H$`\alpha `$ ”knots” are evident over the whole optical disk. A comparison of our YSG map with the H$`\alpha `$ knots shows that a large percentage of the latter (about $`80\%`$) correspond to YSGs. Two out of the six H$`\alpha `$ sources that have no counterpart in our YSG sample are very close to a bright field star, three are in the very outskirts of the optical galaxy and just one lies along a spiral branch. The latter is a very tiny H$`\alpha `$ region perhaps related to a YSG too small to be resolved by us. The most significant difference between our YSGs and Knezek ’s sources is in two big clumps near to the center that don’t show well localized H$`\alpha `$ emission, but are well characterized in our work.

### 3.2 NGC 3507

The spiral NGC 3507 is a typical barred galaxy (SBb type), with two long regular symmetric arms circling around the nucleus and the bar. The inner parts of each arm are more regular than the outer parts. The isophotal major axis is $`D_{25}=3.4\mathrm{}`$ (RC3), the inclination is 26  (Schulman et al. 1997) while the position angle is about 90  (de Jong & van der Kruit 1994). According to Garcia et al. (1993), NGC 3507 is the brightest member of a small group of galaxies. From the radial velocity of 906 km/s (RC3) and $`H_0`$=75 km/s/Mpc, a distance of 12.1 Mpc is obtained, that corresponds to a scale of 16.4 pc/pixel. Schulman et al. (1997) studied the HI distribution and dynamics in a small sample of galaxies, including NGC 3507. This paper concludes that in NGC 3507 there are no High Velocity Clouds (HVCs, i.e. HI regions with peculiar motion). Since HVCs are thought to be star formation indicators (being probably originated by Supernova explosions) the lack of HVCs in NGC 3507 suggests a low star formation rate in this galaxy. Fig. 2 shows the map of the 90 YSGs detected: most of the objects are located along the inner part of the arms. In Table 3 we give the main properties of the YSGs. Galactic and internal B-absorptions are taken 0.00 mag and 0.10 mag, respectively, from RC3. The YSGs average size is 121 pc, with a standard deviation of 147 pc and the slope of the fitted power law to the size distribution is the smallest one with $`\alpha =1.6\pm 0.3`$ ($`r=0.80`$).

### 3.3 NGC 4394

NGC 4394 is a barred and ringed spiral galaxy, with integrated B magnitude 11.53 mag and with $`D_{25}=3.6\mathrm{}`$ (RC3). Its proximity to M85, just $`6\mathrm{}`$ East, and their similar radial velocities suggest that NGC 4394 may belong to the M85 subgroup of the Virgo Cluster. Thus we assume as distance of this galaxy the same (about 16 Mpc) deduced by Ferrarese et al. (2000) for the M85 subcluster basing on Cepheid’s indicators. The inclination and the position angle are taken from Chapelon et al. (1999) and are 25 and 103, respectively.

The distribution of the YSGs, shown in Fig. 3, is clearly concentrated in the ring region and just few YSGs are detected along the spiral arms. This is supported by the radial distribution of the surface luminosity density, $`\sigma _B`$, defined as the total B luminosity of the YSGs found in an annulus divided by its area (see Fig. 5). This quantity is clearly related to the high–mass star formation rate. The concentrated distribution of the YSGs probably reflects the high abundance of cold gas detected in the rings of galaxies (Wong et al. 2000). Moreover, HII regions were studied by Hodge (1974) in NGC 4394 with a 2.1 m telescope. He also found a clear concentration of the HII regions in the ring, few of them being distributed in the outer arms.

Table 4 gives the YSG parameters; the average YSG diameter results to be 114 pc with a standard deviation of 77 pc, while the exponent $`\alpha `$ of the power law is $`2.7\pm 0.3`$ ($`r=0.92`$).

## 4 Discussion and conclusions

With our semi-automatic identification algorithm we have determined the population and main characteristics of Young Star Groupings in the spiral galaxies NGC 3377A, NGC 3507 and NGC 4394. The data collected constitute not only an enlargement of the YSG data-base but also a direct source of astrophysical information about the properties of the regions of star formation in spirals.

The analysis of the YSG differential B luminosity functions (Fig. 6) shows that their high luminosity tails can be well fitted with power laws $`dNL_B^\beta dL_B`$. We find $`\beta =1.79\pm 0.10`$, $`1.31\pm 0.08`$ and $`1.90\pm 0.18`$, for NGC 3377A, NGC 3507 and NGC 4394, respectively (the correlation coefficients r are equal to 0.94, 0.99 and 0.99). These values are in good agreement to those found by Elmegreen & Salzer (1999) for the star forming complexes in a sample of 11 galaxies. Again, the introduction of the Nlim threshold described in Sect. 3.1 implies a correspondent detection limit in the YSG luminosity. We evaluated such detection limits for each galaxy, by plotting the luminosity as a function of number of pixels of the YSGs and then determining the value of the luminosity corresponding to Nlim. The detection limits given in Fig. 6 are close to the peaks of the luminosity functions, implying statistical uncompleteness in those regions.

The sample of six spirals studied so far (three in Paper I and three in this paper) allows us a preliminary investigation of the existence of correlations among the YSG and parent galaxy properties. In Fig. 7 we show some evidence of positive correlations among the number, dimension, integrated B magnitude and the total area of the YSG populations with global characteristics of the galaxies (integrated B magnitude and dimension).

Incidentally, we note that in Fig. 7 the two ringed galaxies, NGC 7217 and NGC 4394, are close each other and, at least in panel b) and d), they seem to behave differently from the rest of the sample. Of course, a larger number of ringed galaxies is necessary to check how real is such a behaviour. In particular, the presence of these two galaxies in our small sample results into a shallower slope of the correlation between the size of the largest YSG in a galaxy and the galaxy $`M_B`$ respect to the one obtained by Elmegreen et al. (1994), who claimed their correlation closely matching the expected variation in the characteristic length of the gaseous gravitational instability with $`M_B`$. Such correlation has been suggested by Selman & Melnick (2000) as a mere result of size-of-sample effect acting on a universal size distribution ($`dND_{xy}^{4.2}dD_{xy}`$). In order to check the Selman & Melnick (2000) suggestion we can compare the $`4.2`$ slope with those of the fitted power laws to the size distributions of the YSGs in NGC 3377A, NGC 3507 and NGC 4394, limited to the statistically reliable size ranges. The slopes found ($`2.3`$,$`1.6`$,$`2.7`$ for NGC 3377A, NGC 3507, and NGC 4394, respectively) are quite close to the value, $`2`$, for star forming molecolar cloud size distribution (see Solomon et al. 1987 and Larson 1981). They are always too shallow to explain the Elmegreen et al. (1994) relation by mean of the size-of-sample effect. We stress, however, that our six points in Fig. 7a suggest a slope not as steep as Elmegreen et al’s relation. Anyway, we remark that our sample is yet too small to draw firm conclusion about this question.