text:
# Relation between Magnitude Series Correlations and Multifractal Spectrum Width

\[

## Abstract

We study the correlation properties of long-range correlated time series, $`x_i`$, with tunable correlation exponent and built-in multifractal properties. For the cases we investigate we find that the correlation exponent of the magnitude series, $`|x_{i+1}x_i|`$, is a monotonically increasing function of the multifractal spectrum width of the original series.

\]

Many natural systems exhibit complex dynamics characterized by long-range power-law 2-point correlations . In some cases the output of such systems — a series of a fluctuating variable $`x_i`$ — may be characterized by a single 2-point correlation exponent, while in other cases the series may be quantified by a spectrum of exponents called the multifractal (MF) spectrum . The MF spectrum quantifies in details the long-range correlation properties of a series $`x_i`$. Nonetheless, additional way to measure the correlation properties of a complex series, $`x_i`$, is to analyze the long-range power-law correlations for the magnitudes of the increments $`|\mathrm{\Delta }x_i||x_{i+1}x_i|`$.

For example, complex fluctuations in financial indices are characterized by a broad MF spectrum . Moreover, recent reports showed that there are positive long-range power-law correlations in the volatility (magnitude of the increments) in such financial time series . Biomedical time series are also found to exhibit similar complex behavior — e.g., the human heartbeat dynamics generates fluctuations which are anti-correlated on a wide range of time scales  and are also characterized by a broad MF spectrum . Furthermore, it was recently reported that the magnitudes of these fluctuations are positively correlated and exhibit power-law scaling behavior over scales from seconds to hours . These and other examples  raise the possibility that there may be a relation between the scaling in the magnitude series and the MF spectrum of the original series which is independent of the type of long-range correlations present in the original series.

Here we study the relation between the scaling exponent which quantifies the power-law correlations in the magnitude series and the MF spectrum of artificially generated time series with analytically built-in MF properties. We find that the scaling exponent of the magnitude series is a monotonically increasing function of the MF spectrum width of the original series.

Long-range power-law correlated series are monofractal when they are characterized by exponents $`\tau (q)`$ for different moments $`q`$ being linearly depend on $`q`$,

$$\tau (q)=Hq+\tau (0),$$
(1)
with a single Hurst exponent ,

$$Hd\tau /dq=\mathrm{const}.$$
(2)
Other long-range correlated series are multifractal and are described by $`\tau (q)`$ that depend nonlinearly on $`q`$,

$$h(q)d\tau /dq\mathrm{const},$$
(3)
where $`h(q)`$ is the local Hurst exponent . The MF spectrum of the series is defined by

$$D(h)=h(q)q\tau (q),$$
(4)
where for monofractal series it collapses to a single point,

$$D(H)=\tau (0).$$
(5)
Here we show that the magnitude series scaling exponent is related to the width of the MF spectrum $`D(h)`$ of the original series,

$$\mathrm{\Delta }hh_{max}h_{min},$$
(6)
where $`h_{min}`$ and $`h_{max}`$ are the minimum and maximum of the local Hurst exponent.

We use the algorithm proposed in  to generate artificial MF series. Briefly, the random MF series is built by specifying discrete wavelet coefficients which are constructed recursively using a random variable $`W`$ (see  for details). Once the wavelet coefficients are obtained, we apply inverse wavelet transform to generate the MF random series; we use the 10-tap Daubechies discrete wavelet transform . The MF properties of the obtained artificial series depend on the form of distribution from which the random variable $`W`$ is chosen.

We consider two different types of probability distributions for the random variable $`W`$, the log-normal distribution and the log-Poisson distribution. For these two examples the MF properties are known analytically . The MF spectrum $`D(h)`$ is symmetric for the log-normal distribution and asymmetric for the log-Poisson distribution.

First we consider the case of a log-normal random variable $`W`$ (Fig. 1a) — i.e., $`\mathrm{ln}|W|`$ is normally distributed and $`\mu `$, $`\sigma ^2`$ are the mean and the variance of $`\mathrm{ln}|W|`$. We choose $`\mu =(\mathrm{ln}2)/4`$. In this case the scaling exponents, $`\tau (q)`$, are given by

$$\tau (q)=\frac{\sigma ^2}{2\mathrm{ln}2}q^2+\frac{q}{4}1,$$
(7)
and the MF spectrum $`D(h)`$ is

$$D(h)=\frac{(h1/4)^2\mathrm{ln}2}{2\sigma ^2}+1.$$
(8)
We calculate the MF spectrum width by solving $`D(h)=0`$ and find

$$\mathrm{\Delta }h=\left(\frac{\sqrt{2}\sigma }{\sqrt{\mathrm{ln}2}}+\frac{1}{4}\right)\left(\frac{\sqrt{2}\sigma }{\sqrt{\mathrm{ln}2}}+\frac{1}{4}\right)=\frac{2\sqrt{2}\sigma }{\sqrt{\mathrm{ln}2}}.$$
(9)
For $`\sigma =0`$, the series is monofractal, while for large value of $`\sigma `$ it has a broad MF spectrum. The MF spectrum is symmetric around $`D(1/4)=1`$ and collapses to a single point $`D(1/4)=1`$ when $`\sigma 0`$ \[Eq. (8) and Fig. 2a,b\].

We generate series with different $`\sigma `$ values ranging from 0 to 0.1 (20 realizations for each $`\sigma `$) and calculate the magnitude scaling exponent using the method of detrended fluctuation analysis (DFA)  (see Fig. 1c and Fig. 3). We find that the scaling exponent of the magnitude series increases monotonically with $`\sigma `$ from uncorrelated magnitude series with scaling exponent $`0.5`$ for $`\sigma =0`$ to strongly correlated magnitude series with exponent $`1`$ for $`\sigma =0.1`$. Further, the magnitude exponent converges to $`1`$.

Next we show that the scaling exponent of the magnitude does not only depend on the positive moments (and especially the fourth moment) but rather relates to the entire MF spectrum. For this purpose we study an example with an asymmetric MF spectrum where the exponents for negative moments, $`h(q<0)`$, are changed drastically when the MF spectrum width is changed, while the exponents for the positive moments $`h(q>0)`$ are less significantly changed. We tune the parameters in such a way that the fourth moment $`\tau (4)`$ is fixed and the second moment $`\tau (2)`$ is hardly changed.

To generate a series with an asymmetric MF spectrum we consider the case of a random variable $`W`$ from a log-Poisson distribution, where $`\lambda `$ is the mean and the variance of a Poisson distributed variable $`P`$, and $`\mathrm{ln}|W|`$ has the same distribution function as $`P\mathrm{ln}\delta +\gamma `$. We choose $`\lambda =\mathrm{ln}2`$ and $`\gamma =\delta ^4(\mathrm{ln}2)/4`$ where $`\delta `$ is an appropriately chosen positive parameter . Under this choice of parameters,

$$\tau (q)=\delta ^q+q\delta ^4/4,$$
(10)
and

$$D(h)=\frac{1}{\mathrm{ln}\delta }\left(h\frac{\delta ^4}{4}\right)\left[\mathrm{ln}\left(\frac{h\delta ^4/4}{\mathrm{ln}\delta }\right)1\right].$$
(11)
Thus the fourth moment $`\tau (4)`$ is independent of $`\delta `$

$$\tau (4)=0,$$
(12)
and the width of the MF spectrum depends on $`\delta `$

$$\mathrm{\Delta }h=\left(\delta ^4/4e\mathrm{ln}\delta \right)\left(\delta ^4/4\right)=e\mathrm{ln}\delta .$$
(13)
For $`\delta =1`$ the MF spectrum collapses to a single point (monofractal) with Hurst exponent $`H=1/4`$ and $`D(1/4)=1`$. For both cases, decreasing or increasing $`\delta `$, the MF spectrum becomes broader. Here we change the MF spectrum width by decreasing $`\delta `$. The major change in the MF spectrum occurs for negative moments $`q<0`$ (Fig. 2 c, and d). We find that the scaling exponent of the magnitude series increases monotonically with the MF spectrum width (Fig. 3) although the positive moments are hardly changed. Thus, the magnitude series scaling exponent is related to the entire MF spectrum of the original series and not just to the positive moments.

We summarize the results of the log-normal and the log-Poisson examples in Fig. 3 — the magnitude scaling exponent is plotted versus the width of the MF spectrum \[Eqs. (9),(13)\]. Surprisingly, these two examples collapse on the same curve. This collapse suggests a possible one-to-one relation between the magnitude scaling exponent and the MF spectrum width .

The relation between the scaling exponent of the magnitude series and the MF spectrum width of the original series is consistent with recent works on MF random walk models with built-in correlations . Recently  it was found that power-law correlations in the magnitude series are related to nonlinear  properties of the increment series $`\mathrm{\Delta }x_i`$ — nonlinear series have correlated magnitude series while linear series have uncorrelated magnitude series. On the other hand, monofractal series are linear series while MF series series are nonlinear . Thus the origin of the magnitude series correlations and the corresponding MF spectrum may be related through the nonlinearity of the original time series.

The results described here may be important from a practical point of view. The calculation of the MF spectrum from a time series involves advanced numerical techniques  and requires long time series. Our analysis of the magnitude series is less sophisticated than the MF analysis and is applicable to shorter time series than one need for the MF analysis. We note that a direct numerical calculation of the MF spectrum width from a time series may not follow the monotonic relation described in Fig. 3 due to overestimation of the MF spectrum width caused by (i) the numerical technique for calculating the MF spectrum \[and estimation of $`\tau (q\pm \mathrm{})`$\] and by (ii) finite series length. However, we expect the magnitude scaling exponent to increase monotonically with the width of the MF spectrum  since both are related to the nonlinearity of the original series.

Partial support was provided by the NIH/National Center for Research Resources (P41 RR13622), the Israel-USA Binational Science Foundation, and the German Academic Exchange Service (DAAD). We thank A.L. Goldberger and J.W. Kantelhardt for helpful discussions.


format: text

text:
# Bendings of radio jets in BL Lacertae objects

## 1 Introduction

While in most radio sources it is possible to see continuity in the direction of the jet, a number of blazars shows large misalignment angles between the jet direction on parsec and kiloparsec scale. The distribution of the misalignment angles ($`\mathrm{\Delta }`$PA), defined as the angle between the parsec and kiloparsec direction of the radio jet has been studied by several authors (i.e. Pearson & Readhead, 1988, Conway & Murphy, 1993, Hong et al., 1998). A comprehensive review on this argument can be found in Appl et al. (1996). They collect the data on 155 sources from six different blazar samples, obtaining a misalignment angle distribution ranging from 0$`\mathrm{°}`$ to 180$`\mathrm{°}`$. Restricting the study to the BL Lac objects the distribution show an expected peak around 0$`\mathrm{°}`$, but also a second peak around 90$`\mathrm{°}`$, harder to explain. Moreover this distribution is lacking in sources with $`20\mathrm{°}<\mathrm{\Delta }PA<40\mathrm{°}`$. Hong et al. (1998), have examined the misalignment angles for the blazars selected in the $`\gamma `$-rays by the EGRET experiment, and found that the BL Lac and quasars’ $`\mathrm{\Delta }PA`$ distributions are different, the quasars being more aligned than the BL Lac objects, with average $`\mathrm{\Delta }`$PA of 21$`\mathrm{°}`$ and 99$`\mathrm{°}`$, respectively.

Different explanations for the presence of the secondary peak have been proposed. Conway & Murphy (1993) claim that it would be impossible to explain the secondary peak through a simple bend or a double population of simply bent sources, proposing that the jet have a helical trajectory generated by Kelvin-Helmholtz instability, probably induced by precession of the rotation axis of the central black hole or of the accretion disk. Such a mechanism drives the jet in a fixed helical trajectory, a sort of tunnel in which the relativistic plasma is flowing. An alternative is the ballistic jet case, in which different components are emitted from the jet basis following straight paths whose direction depends on the instantaneous direction of the jet basis itself. In the case of a non-ballistic helical jet an interesting role is played by the resonance frequency $`\omega ^{}1.3c_{\mathrm{ext}}/R`$ (Conway & Wrobel, 1995), which depends on the external (to the jet) sound speed $`c_{\mathrm{ext}}`$ and on the transversal jet radius, $`R`$. This latter, indeed, grows until $`\omega ^{}`$ is greater than the precession frequency, and becomes constant (saturation) when the two frequencies become comparable. In this phase even the helix amplitude stops growing. In a FR II radio galaxy the intergalactic medium which the jet is going through has a sound speed larger than that required in a FR I, because of the different external density; moreover Bridle & Perley (1984) have determined that the average opening angle of the jets of the FR II is smaller than that of the FR I ones. These considerations lead them to postulate a greater resonance frequency of the gas in a FR II with respect to that of a FR I and hence the saturation of the helical jet would occur sooner in a FR I than in a FR II (Conway & Wrobel, 1995). In principle, this can be useful to distinguish the parent population of the blazar examined, but without adequate statistics it is impossible to use it with some success.

A different mechanism, based on a distortion of the accretion disk, has been proposed by Appl et al. (1996) to explain the misalignment. They have noted that some characteristics of the misalignments, like the lack of connection between the different scales and locally straight parts of the jet, can be compatible with a system of two independent jets, one of which originating in the inner part of the accretion disk, aligned with its rotation axis, and a second one departing from the outer part of the disk having a different rotation axis due to Lense-Thirring precession. Following Appl et al. (1996) the central jet would be the mas-scale jet while the arcsecond scale jet would be that originated from the outer region of the disk. This model could explain the independence of the Doppler factor from the $`\mathrm{\Delta }`$PA. If this is the general scenario, the sources with a small $`\mathrm{\Delta }`$PA would be due to the approximate coincidence of the inner and outer axes.

We observed seven BL Lac objects from the 1-Jy sample (Stickel et al. 1991), with the EVN and the MERLIN at 5 GHz to obtain information on the presence of helical jets on blazars and, if possible, to trace the path of these jets. Our aim was to image the jet from a few tens to few hundreds of mas and to study the evolution of the jet between the pc and the kpc-scale. Other VLBA data at 1.6 GHz have been obtained and are currently undergoing data processing.

## 2 Selection of the sample

We have chosen to restrict our research to the BL Lac objects of the 1-Jy sample because in a previous study on extended radio luminosity we have obtained images for most of the sources of this sample with arcsecond resolution (Cassaro et al. 1999). In order to select the objects to observe we derived the distribution of $`\mathrm{\Delta }`$PA of the sources of the 1-Jy sample, through the direct inspection of the available images on the arcsecond and the milliarcsecond-scales. We have been able to determine the $`\mathrm{\Delta }`$PA for 17 objects only, either because the lack of information on the pc-scale, or the difficulty in determining a clear jet or structure on the largest scale. We did not use the values of $`\mathrm{\Delta }`$PA reported in the literature because they are not obtained with a clearly defined and uniform method. We have compared our VLA images (Cassaro et al. 1999) with the images on the pc-scale available in literature and those obtained from the United States Naval Observatory (USNO) Radio Reference Frame Image Database (RRFID, Fey et al., 1996,1997,2000). We have determined the $`\mathrm{\Delta }PA`$ in the following way: when a clear jet was visible we took the jet axis as the jet direction in pc-scale images, otherwise we considered the straight line connecting the core with the furthest secondary component present in the VLBI image. We considered as kpc-scale direction of the jet the straight line connecting the core and the extended structure (component, elongated halo, jet) nearest to the core, which we believe should be the first region reached by the jet in its path. In the distribution reported in Fig. 1 two depletions are present: the first one, already known, between 45$`\mathrm{°}`$ and 60$`\mathrm{°}`$, and the second one around 120$`\mathrm{°}`$, less evident in analogous distributions by other authors (e.g. Appl et al. 1996). Obviously, this latter depletion in our distribution can be due to the poor statistics. We have chosen objects belonging to the three groups recognizable in the figure, that we call A,B, and C according to a $`0\mathrm{°}<\mathrm{\Delta }PA<45\mathrm{°}`$, $`60\mathrm{°}<\mathrm{\Delta }PA<90\mathrm{°}`$ and $`150\mathrm{°}<\mathrm{\Delta }PA<180\mathrm{°}`$ respectively, obtaining the list in Table 1.

## 3 Observations and data reduction

The radio observations were carried out on 24 November 1998 at 4.99 GHz dual polarization, with an effective bandwidth of 14 MHz, using the European VLBI Network (EVN) and the MERLIN, with the aim to obtain an angular resolution covering structures with angular size ranging from tens to hundreds of mas. The successful EVN antennas for this experiment were Effelsberg, Cambridge, Jodrell Bank (MK2), Medicina, Noto, Onsala, and Torun. The MERLIN observation made use of the antennas of Defford, Cambridge, Knockin, Darnhall, Jodrell Bank (MK2) and Tabley. Each source was observed for approximately one hour, with short scans distributed at different hour angle to improve the uv coverage. The calibration was carried out in a standard way, providing system temperature measurments and gain curves, at the Joint Institute for VLBI in Europe (JIVE) in Dwingeloo (The Netherlands), and the final self-calibration and imaging was performed at the Istituto di Radioastronomia in Noto. We have used the AIPS package for the whole data reduction and analysis.

First, we have obtained the images for the EVN and MERLIN separately and hence, if useful for our aim, we have combined the two datasets to obtain an intermediate resolution image. A small offset ($``$10%) was found between the amplitudes of the two datasets, as could be determined by comparing the data of the Cambridge-Jodrell baseline, common to both datasets. This offset was corrected accordingly.

## 4 Images and comments on individual sources

The VLBI images show jet components in good agreement with previous observations. The MERLIN images, instead, have detected emission at a few hundreds of mas from the core in five out of seven objects.

Secondary components are very weak if compared to those on the arcsecond structure, suggesting a high efficiency in carrying the energy at these intermediate scales. In Table 2 we have reported for each source some characteristics of the components on the EVN, MERLIN and (if any) the combined images: in particular the major and minor axis, the distance and the position angle with respect to the core and the flux density. The images are reported in figures 2 through 19. The reference for the arcsecond-scale images of each source is Cassaro et al. (1999), unless otherwise stated.

0814+425: it is one of the less distorted objects, with a $`\mathrm{\Delta }`$PA=28$`\mathrm{°}`$. The VLBI image (Fig. 2), shows a hint of a jet extended about 10 mas towards the SE direction, confirming previous observations (Polatidis et al., 1995). The MERLIN image (Fig. 3) shows that the jet continues its path for about 100 mas in a direction about 10$`\mathrm{°}`$ bent towards the East. This structure, despite its initial bending, is not a sufficient indication of the presence of a helical jet.

0954+658: this object has a relatively small size, even on the arcsecond-scale. With a $`\mathrm{\Delta }`$PA=74$`\mathrm{°}`$ it belongs to the B group of Fig. 1. The VLBI image (Fig. 4) confirms the presence of a jet in the W-NW direction, slightly resolved, but there are no evident components on the intermediate scale. The weak structure near the core in the MERLIN image (Fig. 5) appears not reliable.

The small size of the source on the arcsecond-scale and the absence of intermediate scale components lead to a jet slightly bent, with a high $`\mathrm{\Delta }`$PA due to projection effects. The lack of (visible) jet at the intermediate scale can be explained by beaming effects: we see the parsec and kpc-scale jet, while the intermediate scale jet is not detectable because it is highly Doppler boosted away from the line of sight.

1147+245: this source has a $`\mathrm{\Delta }`$PA=80$`\mathrm{°}`$, and the jet at the VLBI scale (Fig. 6) is extended about 45 mas towards the West. In the MERLIN image (Fig. 7) a secondary component is present at about 150 mas from the core in the SW direction, at an angle of about 50$`\mathrm{°}`$ from the milliarcsecond-scale jet and 30$`\mathrm{°}`$ from the arcsecond-scale jet. There is no trace of the bending region between the pc-scale jet and the component on the MERLIN scale. The change of direction could be either due to a shock in the interstellar medium or to the passage through a region where the magnetic field changes direction abruptly (Koide, 1997), or, finally, to a curvature of the jet enhanced by projection effects. Even if we cannot exclude the first two hypotheses, we favour the last one, because of the absence, between the pc scale jet and the secondary component in the EVN+MERLIN image (Fig. 8) of a bright region near to the hypothetical shock in which the jet should dissipate energy and perhaps lose collimation. In Fig. 18 we present radio images of 1147+245 at different angular scales: from left to right our EVN+MERLIN image, the VLA A-array image at 8.4 GHz obtained by Dallacasa et al. (in preparation) and the VLA A+B array image at 1.36 GHz (Cassaro, 2000). The 8.4-GHz VLA image shows that the jet continues to the South then bends sharply toward the SE direction. In the 1.36-GHz VLA image we finally see a double structure with wide jets or maybe lobes where the southern emitting region appears to bend slightly towards the SE direction. The structure is compatible with a helical jet and the abrupt change of direction in the VLBI and VLA 8.4 GHz images would be due to the projection of the helicoid on the plane of the sky.

We note that our images do not cover the scale between 0.15″and 0.5″, that would be important to describe the jet behaviour in the first hundreds of parsec, where we expect to better see the jet behaviour. This angular scale should be adequately covered by the new VLBA data.
1308+326: on the arcsecond scale this source shows a diffuse halo around the core with a luminous region toward the East and a bright pointlike component at about 10 arcseconds toward the North (Murphy et al. 1993). The VLBI image (Fig. 9) shows a jet extended about 10 mas in the W direction, with a $`\mathrm{\Delta }`$PA of 180$`\mathrm{°}`$ or 90$`\mathrm{°}`$ depending on which arcsecond structure we consider as the end of the jet. We favour the first possibility, given that the MERLIN image (Fig. 10) reveals a well-defined jet which departs to the West and after about 100 mas bends sharply toward the SW direction, for a total extension of about 200 mas. The jet appears to follow a curved path that encounters the arcsecond-scale East component. In this case the hypothesis of a curvature due to a shock could be considered, but a helical structure can also easily explain the jet trajectory from mas to the arcsecond-scale.

The observed structures could be also compatible with a mechanism similar to that described for 0954+658. The jet bends smoothly, describing an arc from the mas-scale to the arcsecond-scale, which for projection effects would be visible only in the inner part, becoming invisible between 0.2″and 10″in which region the radio emission is boosted away from the line of sight; then it again becomes observable at $``$10 arcsecond of projected distance from the core at the opposite side of the core, when the velocity of the gas has slowed down and the Doppler factor decreased enough. To verify this scenario it would be important to image the uncovered range of distance from 0.2″to about 1″-2″, to search for possible components of the jet, which, if the latter description is correct, should not be found.
1418+546: in this case $`\mathrm{\Delta }`$PA=150$`\mathrm{°}`$. The VLBI image (Fig. 11) shows a straight jet directed toward the SE direction. On the intermediate scale (Fig. 12) this source shows a sharp bending perpendicular to the VLBI jet direction. For this source, as in the previous one, the misalignment and the structure are easily explained with a smoothly bending jet, whose curvature is enhanced by the projection.
1803+784: this is the largest object of the 1-Jy sample, with total projected linear size of 444 kpc. In the VLA image a low-brightness region is evident at about 9 arcseconds toward the SW direction with respect to the core; the total radio structure resembles a Wide Angle Tailed (WAT) source. Its $`\mathrm{\Delta }`$PA is 70$`\mathrm{°}`$ and belongs to the second group of Fig. 1. The VLBI image (Fig. 13) shows a collimated jet about 10 mas long which abruptly loses collimation and ends in a diffuse halo slightly bending toward the South, in the direction of the large scale structure. The MERLIN image (Fig. 14) does not reveal any further emission. The extended structure in the image of Fig. 13 could likely be originated by a shock with a denser medium, after which the jet loses collimation and changes its direction towards the arcsecond-scale structure. The enhancement of the apparent bending angle can be due to the projection onto the plane of the sky. The presence of a helical jet for this source is unlikely.
2131$``$021: we have determined for this source a $`\mathrm{\Delta }`$PA=36$`\mathrm{°}`$, although in the VLA image different structures can be considered the end of the VLBI jet. In order to determine the misalignment we have taken into account the arcsecond-scale component closest to the core, which presumably will be reached first. The VLBI image (Fig. 15) shows a component in the East direction extended for about 20 mas. In the MERLIN image (Fig. 16) a second component is present toward the South at about 280 mas from the core. The combined image (Fig. 17) shows a jet extended for about 30 mas in the SE direction, and a series of components extending up to 350 mas from the core. The jet appears well collimated, and seems to follow a swinging trajectory, indicating a helical motion. The jet seems collimated and the amplitude of the oscillation constant. If the new VLBA data confirm this behaviour, it would imply that the saturation takes place in the first tens of mas.

In Fig. 19 we present radio images of 2131$``$021 at different angular scales: from left to right our EVN+MERLIN image, the VLA A-array image at 8.4 GHz obtained by Dallacasa et al. (in preparation) and the VLA A+B array image at 1.36 GHz (Cassaro, 2000). The 8.4 GHz VLA image shows that the jet axis has a $``$45$`\mathrm{°}`$ change in its direction while the 1.36 GHz image shows a structure reminiscent of a Narrow Angle Tail (NAT) or a Wide Angle Tail (WAT) radio source. The NAT or WAT classification of this source is also consistent with the change in the direction of the jet axis, which likely occurs at the boundary of the host galaxy.

Even in this case the new VLBA data should cover the scale between 0.35″and about 1″missing in our images, and would be useful to properly trace the jet and possibly determine the saturation. Moreover, it could be possible to verify where the jet axis changes direction and test the NAT/WAT classification.

## 5 Discussion and conclusions

EVN+MERLIN observations at 5 GHz have been carried out to investigate the origin of the misalignment between the arcsecond and milliarcsecond-scale structures and to attempt to derive some information from the jet properties on intermediate scale.

A possible explanation for the misalignment, and in particular for the presence of the secondary peak in the $`\mathrm{\Delta }`$PA distribution has been provided by Conway & Wrobel (1995) who consider the presence of a helical jet. Unfortunately, the application of their conclusions has been only possible for Mrk 501, for which the images clearly showed the helical structure and the saturation region. Only two sources observed in this work have morphologies consistent with the presence of a helical jet (1147+245, 2131$``$021 Figs. 18 and 19). In no case it has been possible to clearly detect a helical jet and it is not possible to determine whether and where the saturation occurs.

Even if there is no sure identification of a helical jet, we can infer some qualitative conclusion from the inspection of our images. The likely helical jets belong to sources with different $`\mathrm{\Delta }`$PA: 80$`\mathrm{°}`$ for 1147+245, 36$`\mathrm{°}`$ for 2131$``$021. Moreover, in the radio source 1803+784, the strong loss of collimation favours an intrinsic bending due to a shock rather than a helical jet. In summary, one of the sources shows a helical jet structure and do not have a 90$`\mathrm{°}`$ $`\mathrm{\Delta }`$PA, and at least one other source with a misalignment of about 90$`\mathrm{°}`$ does not show a helical structure. This suggests that Conway & Murphy’s (1993) model, which postulates a population of helical jets to explain the 90$`\mathrm{°}`$ misalignment, may not work for all objects.

Appl et al. (1996) attribute the peak at 90$`\mathrm{°}`$ in the distribution of the misalignments to the presence of two jets departing from different regions of the accretion disk, without considering any intrinsic distortion of the jet. Several sources examined in this work show characteristics compatible with such a scenario. Some of them, however, show components on intermediate scale (i.e. 2131$``$021) which do not agree with the presence of only two jets. Moreover some sources have intermediate scale jets that begin to bend, as is the case of 1308+326 and 1418+546. This is in contrast with the conclusions of Appl et al. (1996) who consider straight jets. We note furthermore that such a mechanism can explain the sources with a misalignment around 0$`\mathrm{°}`$, but cannot do the same with sources belonging to the secondary peak in the misalignment distributions, or better, it cannot explain the lack of BL Lac objects with $`45`$°$`<\mathrm{\Delta }PA<60\mathrm{°}`$.

The misalignment distribution of Appl et al. (1996) restricted to the quarsars lacks the secondary peak around 90$`\mathrm{°}`$ which is present instead for the BL Lac objects. This suggests that the mechanism that produces the misalignment is likely to be different for these two classes of objects. In particular, for the BL Lacs we can conclude from the images obtained in the present work, that it is unlikely to invoke a single mechanism to explain the misalignment between the arcsecond and milliarcsecond scale structures and the different structures seen in the intermediate scale images.

While a helicoid is not surely proven to be present in all our misaligned sources, a less restrictive hypothesis of non-coplanar bending is required to explain a 90$`\mathrm{°}`$ peak in the distribution. By non-coplanar bending we mean a jet for which the directions in starting and ending regions of the path do not lie in the same plane; a helicoid itself or a branch of a helicoid resulting from a precession motion of only a fraction of a period, are examples of non-coplanar bending. It is also possible to explain the misalignment distribution of the BL Lacs considering that there may co-exist different kinds of objects. A first category with slightly distorted jets, and misalignments distributed from 0$`\mathrm{°}`$ to 180$`\mathrm{°}`$, with a higher probability of sources with $`\mathrm{\Delta }`$PA around 0$`\mathrm{°}`$. A second category has misalignments around 90$`\mathrm{°}`$ caused by non-coplanar jets, helical jets or strong ambient interactions (e.g. 1803+784), explaining the numerical eccess of BL Lac objects with such a misalignment.

###### Acknowledgements.

The European VLBI Network is a joint facility of European and Chinese radio astronomy institutes funded by their national research councils. MERLIN is operated as a National Facility by the University of Manchester at Jodrell Bank Observatory on behalf of the UK Particle Physics & Astronomy Research Council. The National Radio Astronomy Observatory is a facility of the National Science Foundation operated under cooperative agreement by Associated Universities, Inc. PC acknowledge the Joint Institute for VLBI in Europe (JIVE) for the support and hospitality during the Summer Research Program 1999. PC is also very grateful to Denise Gabuzda and Michael Garrett for the very useful help in data reduction.


text:
# Theoretical Limb Darkening for Pulsating Cepheids.

## 1 Introduction

Classical Cepheids play a central role in setting the Extragalactic Distance Scale. In fact, Cepheid high intrinsic luminosity, and the existence of a well studied relation between their pulsation period and luminosity (first discovered by Leavitt 1906), makes them crucial primary distance indicators. These properties have been exploited extensively, from early surveys to the recently completed HST Key Project to measure the Hubble constant (Freedman et al., 2001).

The reliability of Cepheids as distance indicators depends on the accurate calibration of the period-luminosity relation zero-point. This is difficult to accomplish, due to the low density distribution of Cepheids in our galaxy. Direct determination of the trigonometric parallaxes of the nearest Cepheids has been attempted using the database of the HIPPARCOS mission (see e.g. Feast & Catchpole 1997), but large errors do not allow an unambiguous interpretation of the results (Madore & Freedman, 1998). An accurate measurement of the parallaxes of even the closest galactic Cepheids does not appear feasible, even with the space based surveys proposed for the near future (Horner et al., 2000). Most promising, among the techniques which can provide an independent calibration of the period luminosity zero-point, is the Baade-Wesselink (BW) method (Baade, 1926; Wesselink, 1946).

Recent progress in ground-based interferometry finally allows direct detection of pulsations in nearby Cepheids (Lane et al., 2000). The measurement of such radial angular variations, combined with spectroscopic observations yielding the radial velocity of the pulsating photosphere, can in principle allow a precise determination of nearby Cepheid distances (Sasselov & Karovska, 1994, hereafter SK94). As described in SK94, the feasibility of such a measurement is conditioned on the availability of accurate hydrodynamic models for the Cepheid atmospheres. These are required (1) to derive the projection factor necessary to convert the photospheric line velocities into radial motion, and (2) to compute reliable time-dependent limb darkening profiles. The first of these requirements is described in SK94 and addressed in detail by Sabbey et al. (1995). In this paper we describe a new approach to solve the second of these requirements.

Currently, limb darkening tables specifically computed for Cepheid atmospheres are not available, so coefficients derived for static atmospheres of non-variable yellow supergiants are generally used (Parsons, 1971; Manduca at al., 1977; Manduca, 1979; Kurucz, 1993a; Claret et al., 1995). This follows the assumption that a Cepheid atmosphere at a given phase can be approximated with a static plane-parallel atmosphere having the same $`T_{\mathrm{eff}}`$ and $`\mathrm{log}g`$, as estimated by fitting the depth ratio of weak photospheric lines. This is generally true as long as the observational constraints do not demand a level of precision that is affected by the time-dependent hydrodynamics and by the spherical geometry of the real pulsating star.

The assumptions of hydrostatic equilibrium and plane parallel geometry lead to an underestimation of the temperature gradients in yellow supergiants, and thus fail to reproduce the correct IR and UV fluxes (Fieldus & Lester, 1990). They also cannot reproduce the relative separation between the continuum photosphere and the line-forming region, and thus underestimate the difference in the radial-velocity amplitudes (Karp, 1975; Sasselov & Raga, 1992). All these and other limitations of static plane-parallel atmospheres seriously affect the precision of Cepheid distance determination with the BW method. Current demands on a Cepheid distance scale accuracy of 10% and better require us to abandon the approximation of hydrostatic equilibrium and plane-parallel geometry (see SK94).

In Sasselov & Lester (1994) new time-dependent, non-LTE dynamic models of Classical Cepheids were presented, and Sabbey et al. (1995) used them to study spectral line formation. In our paper, updated models are used to further explore the consequences of the dynamic treatment of Cepheids pulsation in the resulting spectral brightness distribution and limb profiles. This last step is complex, since a consistent model, able to treat at the same time non-LTE radiative transfer and spherical geometry with all the necessary opacities and line sources in a time-dependent hydrodynamic framework, is still not feasible. Section 2 illustrates our solution to this problem.

In section 3 the typical results for a classical Cepheid are analyzed, and the differences between our computed limb darkened profile and the one from an equivalent static atmosphere are discussed. The wavelength dependence of our limb profiles is also analyzed. In section 4 we finally discuss the applications of our models, and their implication on the analysis of presently available and future interferometric measurements.

## 2 Spectral Intensities for Dynamic Atmospheres

The computation of a realistic spectral intensity distribution for a stellar model requires a complete set of opacities for the atmospheric atomic species. If $`T_{\mathrm{eff}}`$$`7000`$ K (which is generally true for classical Cepheids), molecules should also be included. Ideally, this requires solving the full radiative transfer problem and the hydrodynamic equations in a self consistent way. However, the needed computing power is excessive, and therefore the hydrodynamic structure of the atmosphere is computed by considering only a limited set of elements (see discussion in section 2.1). The full set of elements and their opacities should somehow be reintroduced when deriving the emergent spectra.

A possible approach consists in using the Rosseland opacities derived in hydrostatic equilibrium. The hydrostatic equivalent of the dynamic model is chosen to match as close as possible the dynamic atmosphere, as it changes during the Cepheid pulsation. For each snapshot of the hydrodynamic simulation, one has to search for a static LTE model producing the same total flux of the star, and having a similar atmospheric structure. A “hybrid” model is thus artificially constructed, having the depth-dependent opacities of the static atmosphere, and the temperature, electron density and pressure of the dynamic model atmosphere.

The main limitation of this procedure is that the opacities are still computed in LTE conditions. However, the intensity distributions thus derived are much improved, since the thermodynamical state is set by the hydrodynamical calculations. The accuracy and caveats of this technique are analyzed and discussed in the following sections.

### 2.1 Dynamic Models

Hydrodynamic simulations show fundamental departures in the structure of pulsating atmospheres from hydrostatic equilibrium and LTE (Willson & Bowen, 1985; Cuntz, 1987). Our hydrodynamical calculations only include model ions of H, He I and II, Ca II and Mg II (see Sasselov & Lester, 1994, for details). Additional non-LTE calculations are then performed with model ions of C I and II, O I, Mg I, Ca I and Fe I. These additional species do not provide any feedback for the radiative hydrodynamic simulations and are thus considered only for the starting models, and after the dynamic calculations have converged. For galactic Cepheids, Solar composition is assumed for all elements.

The pulsation is introduced by perturbing the atmosphere with a “piston” having suitable period, amplitude and shape. These parameters are specific for each Cepheid under study, and are derived from observational constraints. The period of the piston is the period of the Cepheid. The shape of the piston is taken from linear non-adiabatic model computations of Buchler et al. (1990). The piston amplitude is the main free parameter, and is derived by matching the radial velocity curve, the phase lags, the estimated extension of the atmosphere, the observed color and light variations and other basic observable features. When the oscillations introduced by the piston stabilize, and no transients remain, a set of time-dependent pulsating photospheric models are computed for the entire pulsational cycle. The details of this procedure are explained in Sasselov & Lester (1994).

Given a good model for the Cepheid piston, it is thus possible to fit the main stellar observables with a set of dynamic models, each one representing a snapshot of the atmosphere at any pulsational phase. We have applied this procedure for $`\zeta `$ Gem, which is a typical classical Cepheid.

We note that in Sasselov & Lester (1994) a further step is performed, which consists of building a semi-empirical chromosphere to derive the phase-dependent profile of individual chromospheric lines. Since we are mostly interested in deriving low spectral resolution intensities and visibilities, the computed chromosphere is ignored in this work, due to its negligible contribution to the overall limb darkening profiles.

### 2.2 Dynamic Models and Stellar Observables

The Cepheid pulsation leads to the extension and contraction of the atmosphere, responsible for the changes of the stellar observables ($`T_{\mathrm{eff}}`$, photospheric lines velocity, radius, luminosity, etc). These quantities, which change with the phase of a pulsating star, parameterize the computation of its phase-dependent spectrum. Let us define a consistent reference framework for the phases of a Cepheid cycle.

Traditionally, the zero-point phase of variable stars is set to coincide with maximum luminosity. This choice, however, is not convenient when dealing with time-dependent hydrodynamic models, since the emergent flux in the visible is a derived quantity, not always available in the numerical computation. A better approach is to use hydrodynamic quantities directly related to observables. The common choice is the pulsational velocity of the Cepheid photosphere, which can be observed by measuring the Doppler shift of suitable photospheric lines.

Figure Theoretical Limb Darkening for Pulsating Cepheids. (left) shows the radial velocity and lightcurve of $`\zeta `$ Gem. These were derived by tables published by Bersier et al. (1994a, b), where these parameters are measured for a large set of Cepheids, and expressed in terms of their Fourier coefficients. The conversion of radial velocities into pulsational motion requires us to perform the following geometric transformation:

$$v_pp\left(v_r\gamma \right)$$
(1)
where $`\gamma `$ is the systemic velocity taking into account the radial velocity of the star relative to the solar system, the minus sign inverts Doppler shifts into velocities relative to the Cepheid atmosphere, and $`p`$ is the projection factor. The pulsational velocity of $`\zeta `$ Gem derived using Equation 1 is plotted in Figure Theoretical Limb Darkening for Pulsating Cepheids. (right-top panel), using a constant projection factor of 1.36 (Burki et al., 1982). The systemic velocity $`\gamma =5.91`$ km/s was determined by requiring conservation of radius after the completion of one pulsational cycle, and is consistent with the value of $`6.5(\pm 1.0)`$ km/s derived by Wilson (1953) for the radial velocity of the star.

The zero phase in the velocity curve, and in the subsequent discussion, is set at the time in which $`v_p`$ passes from negative to positive, e.g. when the star begins to expand after reaching the minimum radius. This point is marked in Figure Theoretical Limb Darkening for Pulsating Cepheids. with a solid vertical line. In the case of $`\zeta `$ Gem, the maximum luminosity occurs at phase $`\varphi 0.28`$, which is 3 days after the zero-point phase (dashed line), given its 10.15 days period. Minimum luminosity is at $`\varphi 0.76`$ (dotted line), while maximum radius is at $`\varphi 0.59`$ (dot-dashed line).

The radius R is obtained by integrating the velocity $`v_p`$ over time:

$$R(t)=R_0+_{t_0}^tv_p(t^{})dt^{}$$
(2)
where $`R_0`$ is the mean radius. An example is shown in Figure Theoretical Limb Darkening for Pulsating Cepheids. (bottom-right) for $`\zeta `$ Gem, which has $`R_069`$ R (Krockenberger, 1996).

Note, however, that the above expression is strictly correct only if $`v_p`$ is associated with an individual mass element in the stellar atmosphere. This is not the case for the photospheric line velocities (as the lines are produced by different photospheric layers during the stellar pulsation), nor for the velocity field of the hydrodynamic model (as the radial grid is not comoving with the atmosphere). This inconsistency, which translates in the non-conservation of the path integral in equation 2 over one pulsational period, should be taken into account in order to have a consistent reference frame for the observed and model velocities.

The condition of constant average radius for the pulsating atmosphere can however be achieved by adding an extra factor to the $`\gamma `$ constant in Equation 1. This was done implicitly when deriving the systemic velocity for $`\zeta `$ Gem. We did the same with the model pulsational velocity, deriving a correction factor of $`0.475`$ km/s. With this condition, the model and observed $`v_p`$ are directly comparable, allowing to set the zero point phase in the sequence of models.

The pulsational velocity of the hydrodynamic model was defined as the average velocity of the photospheric layers, weighted over the region where the line formation occurs. This quantity, shown in Figure Theoretical Limb Darkening for Pulsating Cepheids., is very similar to the observational $`v_p`$. Note the perturbation starting at $`\varphi 0.9`$, which is associated with a shock-wave crossing the photosphere.

The bottom panel of Figure Theoretical Limb Darkening for Pulsating Cepheids. shows the effective temperature of $`\zeta `$ Gem as a function of phase, derived by fitting the depth ratio of weak metallic lines as Fe I, Fe II, Si I and a few other elements insensitive to the changes in metallicities (Krockenberger et al., 1997). The effective temperature here defined is the $`T_{\mathrm{eff}}`$ of the plane-parallel LTE atmosphere used to fit the observed line ratios. We will maintain this convention through the following discussion, since this quantity allows direct comparison between a model parameter and the observations.

The effective temperature plotted in Figure Theoretical Limb Darkening for Pulsating Cepheids. is related to the luminosity curve of the Cepheid, and is thus an important observable that we use to constrain our model intensity spectra, as explained in the following section.

### 2.3 Rosseland Opacities for Dynamic Atmospheres

The largest collection of stellar atomic and molecular opacities is available as part of the ATLAS9$``$12 numeric code (Kurucz, 1970, 1979, 1993b), which incorporates no less than 58 million lines for different stellar metallicities and values of microturbulent velocity.

The Kurucz code initially computes static LTE model atmospheres in plane-parallel geometry and radiative equilibrium. The input parameters which define each model are the effective temperature $`T_{\mathrm{eff}}`$, the gravity $`\mathrm{log}g`$, the convective mixing length $`L/H`$, the turbulent velocity $`V_{\mathrm{turb}}`$ and the metallicity. The structure of a model atmosphere is evaluated on a radial grid in the mass coordinate $`\mathrm{d}M=\rho \mathrm{d}r`$. For each depth the model computes the local temperature, total pressure, electron density, Rosseland opacity, radiation pressure and turbulent velocity (which is not allowed to change along the radial grid). The SYNTHE program, included in the ATLAS package, then allows the computation of intensity spectra $`I_\nu (\lambda ,\mu )`$, from which limb darkened profiles are derived.

In this work we adopted the ATLAS9 version of the Kurucz code, which incorporates a full treatment for line-blanketing; we do not expect that adopting ATLAS12 model atmospheres with opacity sampling would have produced significant differences in the atmospheric structure and resulting intensity profiles, at least in our wavelength range of interest (Kurucz, private communication).

Figure Theoretical Limb Darkening for Pulsating Cepheids. shows the main differences between the snapshot of a dynamic atmosphere and a family of ATLAS9 models having $`\mathrm{log}g`$ from 0.5 to 2.0 and $`T_{\mathrm{eff}}`$ $`5600`$ K. The effective temperature for the static model has been chosen to match the observational estimate for the Cepheid $`T_{\mathrm{eff}}`$ at the same phase. Contrary to the hydrostatic case, the dynamic model shows a temperature inversion for $`\mathrm{log}\mathrm{d}M1.5`$. This is due to the chromosphere attached to the Cepheid atmosphere, which is not reproduced by the static radiative equilibrium models. The deeper layers of the atmosphere, on the contrary, can be better approximated by hydrostatic equilibrium, with an appropriate choice in the $`\mathrm{log}g`$ parameter.

The agreement between the temperature distribution in the static and dynamic atmospheres, below the temperature inversion point, depends on how significant are the departures from hydrostatic and radiative equilibrium. A minor role is also played by departures from plane-parallel geometry. As shown in Sasselov & Lester (1994), however, spherical geometry becomes an important factor only for long period Cepheids, and can be neglected for our model of $`\zeta `$ Gem and the other short period Cepheids that can be observed with the available interferometers.

Consider phase $`\varphi 0.40`$ (figure Theoretical Limb Darkening for Pulsating Cepheids., top panel), which is close to average radius, when the star is freely expanding. This is well after the effects of the shock-wave (occurring at $`\varphi 0.9`$) dissipated. The dynamic and static thermal structures are very similar, at least where the convective overshooting incorporated in ATLAS9 (absent in the dynamic models) is not important. The overshooting only affects the deepest layers in the stellar atmosphere, and is responsible for the lower temperature gradient in the hydrostatic temperature profile.

At other phases (e.g. for $`\varphi 0.94`$), the hydrodynamic effects (the shock-wave in particular) are stronger, causing a larger departure from hydrostatic equilibrium. (Figure Theoretical Limb Darkening for Pulsating Cepheids., bottom). The spread in the temperature distribution in this case can be large, and a single hydrostatic model which is a good approximation of the hydrodynamic atmosphere cannot be identified.

In a typical Cepheid atmosphere, the contribution function for the photospheric lines is negligible both in the chromosphere, and in the deep layers subjected to convective overshooting, as shown in Sabbey et al. (1995). In fact, the region where the hydrodynamic thermal structure is best approximated by an equivalent hydrostatic atmosphere, is precisely where the photospheric contribution functions are at maximum. Our assumption is that the opacities computed in hydrostatic equilibrium in this region can be used to compute emergent spectral intensities with the dynamic model.

On these basis, we fit our dynamic model for $`\zeta `$ Gem with a grid of ATLAS atmospheres having $`V_{\mathrm{turb}}`$ $`=4`$ km/s and $`\mathrm{log}g`$ from 0.5 to 2.0 in steps of 0.1. Since the effective temperature for our hydrodynamic models is not well defined, we have set the ATLAS fitting grid to the $`T_{\mathrm{eff}}`$ value determined by Krockenberger et al. (1997) for each pulsational phase. The fit was done by minimizing the normalized $`\chi ^2`$, computed on the “fitting region” in which the contribution functions are at maximum, between each dynamic model and all static atmospheres in the ATLAS grid. This technique allows to find the best match between the dynamic models and the static atmospheres, even in the cases in which the dynamic model temperature structure cuts across several static models (as shown in Figure Theoretical Limb Darkening for Pulsating Cepheids. for $`\varphi =0.94`$). The result is a family of static atmospheres, each having a Rosseland opacity computed with the full set of relevant atomic species, to be used instead of the unknown opacities for the dynamic models. The error associated with this procedure can be visualized as the separation between the opacities of two consecutive models in the fitting grid, which is $`\pm 0.3`$ in the logarithmic scale (see figure Theoretical Limb Darkening for Pulsating Cepheids.).

Note that, once the best fit static model is determined for each phase, it is the whole atmospheric thermodynamic status that is constrained, e.g. not only the temperature distribution $`T`$, but also the electron density $`x_{\mathrm{ne}}`$ and the pressure $`P`$. Differences in the $`x_{\mathrm{ne}}`$ and $`P`$ distributions between static and dynamic models are to be expected, given the different equation of state used in the two cases. These inconsistencies reflect in the accuracy of the Rosseland opacity adopted for the dynamic atmosphere, which is computed with the pressure and electron density relative to the hydrostatic case.

### 2.4 Intensities and Limb Profiles

The spectral intensity distribution for each phase of our dynamic atmospheres are finally computed by adapting SYNTHE/ATLAS to process our “hybrid” models which combine the hydrodynamic atmospheres with hydrostatic Rosseland opacities.

Each “hybrid” model is characterized by the following parameters: (1) $`T_{\mathrm{eff}}`$ is the effective temperature of the best fit ATLAS model, which was set to the value derived by Krockenberger et al. (1997); (2) $`\mathrm{log}g`$ determined by the best fit of the dynamic model with the static grid; (3) $`V_{\mathrm{turb}}`$, $`L/H`$, metallicity and Rosseland opacities taken from the best fit hydrostatic model; (4) $`T`$, $`x_{\mathrm{ne}}`$ and $`P`$ obtained by splicing the hydrodynamic model (below the temperature inversion point) with the static model (above $`T`$ inversion, e.g. effectively removing the chromosphere). The removal of the chromosphere is necessary, because SYNTHE is in essence a radiative equilibrium code and cannot treat temperature inversions and discontinuities in the radial grid. This step is justified, as noted before, by the consideration that the only part of the atmosphere relevant in the photospheric emission, is where the contribution function is maximum. This is the same region used in section 2.3 to perform the fitting between static and dynamic atmospheres. Eventual discontinuities around the temperature inversion point in the two segments of the $`T`$, $`x_{\mathrm{ne}}`$ and $`P`$ profiles have been eliminated by rescaling the static part with a suitable normalization factor.

An example of final temperature, electron density and pressure radial profiles of our $`\zeta `$ Gem model is shown in figure Theoretical Limb Darkening for Pulsating Cepheids. for two representative phases. At $`\varphi 0.40`$, the $`T`$ and $`x_{\mathrm{ne}}`$ dynamic model distributions closely match their hydrostatic approximation, at least in the region used for the fit (marked in the figure by the vertical dotted lines). The pressure of the dynamic atmosphere, on the contrary, is lower than the static equivalent, reflecting the change in the equation of state. At phase 0.94, the extra energy dissipated by the traveling shock wave is responsible for an increase in the dynamic model electron density and pressure, which becomes higher than in the hydrostatic atmosphere with the same temperature profile. The increased electron density of the “hybrid” model, reflecting the structure of the hydrodynamic atmosphere, is in part responsible for the different limb profile of the Cepheid at this phase.

The spectral intensity distributions (and the associated limb profiles) are then computed (together with a consistent set of molecular opacities) by running SYNTHE on the “hybrid” models. The results are discussed in the following section.

## 3 Discussion

Limb profiles for our model Cepheid are obtained by plotting the intensity $`I_\nu (\lambda ,\mu )`$ as a function of the projected radial coordinate $`\mathrm{sin}\alpha `$ (where $`\alpha `$ is the angle between the line of sight and the direction of the emergent flux, and $`\mu =\mathrm{cos}\alpha `$). Given the detailed treatment of stellar opacities in our computation, the profiles thus obtained can be trusted analyzing the wavelength dependence of limb darkening. Also, the explicit time dependence of our dynamic models allows us to follow the variations in the limb profile as the Cepheid star pulsates.

### 3.1 Differences between static and dynamic models

Figure Theoretical Limb Darkening for Pulsating Cepheids. compares static and dynamic models at three representative phases. At $`\varphi 0.4`$ the limb darkening profiles computed with the hydrodynamic models are very similar to the hydrostatic equivalent. This is not surprising, since at this phase the atmosphere is quasi-static, having long dissipated the transient hydrodynamic effects caused by the passing shock-wave. Even though the pressure is very different in the static and dynamic case, the spectral emission (and thus the limb darkening) is similar, confirming that $`T`$ and $`x_{\mathrm{ne}}`$ are the main parameters constraining the photospheric emission. This is an important validation of the procedure described in section 2.

Dramatic evidence that the limb darkening is substantially affected by the hydrodynamic effects is given in the two lower panels of Figure Theoretical Limb Darkening for Pulsating Cepheids.. The two panels compare the static and dynamic limb profiles at $`\varphi 0.92`$ (shock-wave crossing the photosphere) and at $`\varphi 1`$ (minimum radius). In the first case the newborn shock-wave reduces $`T`$ with respect to the hydrostatic atmosphere, and the limb darkening is lower. At minimum radius the effects of the shock-wave have dissipated, and the compression of the atmosphere gives rise to a larger temperature gradient. This is responsible for a larger limb darkening. In both cases the changes in the amount of limb darkening with respect to the static model is mostly visible at shorter wavelengths, where the opacity is lower and the limb darkening larger.

### 3.2 Changes of limb darkening with phase and wavelength

Figure Theoretical Limb Darkening for Pulsating Cepheids. puts in evidence the full amount of monochromatic limb profile variations as a function of phase, for $`\lambda =500`$ nm. As explained before, the limb darkening for our model of $`\zeta `$ Gem is greater at minimum radius, and lower when the shock-wave is crossing the photosphere.

The effect of phase on the limb profiles is not only restricted to the amount of limb darkening, but is also reflected in the shape of the profile itself. Consider the two profiles for $`\varphi 0.76`$ and $`\varphi 0.91`$. Despite the small difference in phase, they cross each other at $`\mathrm{sin}\alpha 0.75`$, which indicates that the gradient of $`I_\nu `$ with respect of $`\mathrm{sin}\alpha `$ is also a function of the pulsational phase. This example highlights the intrinsic limits, in the parameterization of the limb profiles, that result from a single limb darkening parameter.

A closer examination of figure Theoretical Limb Darkening for Pulsating Cepheids. shows that, even when the dynamic atmosphere is close to hydrostatic equilibrium ($`\varphi 0.4`$), static models are slightly less darkened for $`\mathrm{sin}\alpha 0.85`$ and more darkened on the limb. Even in the most favorable phases, the dynamic limb profiles cannot be obtained by scaling suitable profiles computed with a static atmosphere. This is a direct consequence of the similar difficulties, mentioned in section 2.3, of fitting the hydrodynamic thermal structure with a single hydrostatic model. Again, the different functional shape between the static and dynamic profiles depends strongly on the phase and the wavelength, and requires detailed modeling.

As expected, our limb profiles have larger darkening at shorter wavelengths. This is mainly a consequence of the source function and the optical depth dependence on wavelength. At 900 nm and in the near-IR, different lines of sight along the stellar disk probe regions of the stellar atmosphere where the source function is similar. This results in a flat limb intensity profile, and thus low limb darkening. At 500 nm the source function in the atmospheric layers probed by the different directions change more rapidly towards the limb, producing a larger limb darkening. This effect is even larger when going to shorter wavelengths, in the UV (300–400 nm).

This wavelength dependence of limb darkening is also affected by phase dependent hydrodynamical effects. Figure Theoretical Limb Darkening for Pulsating Cepheids. presents limb profiles for $`\lambda =500`$, 700 and 900 nm at two different phases. As shown in the figure, the limb profiles derived for dynamic atmospheres are in most cases very different from those computed using static models, or from linear limb darkening curves (Milne, 1921).

## 4 Conclusions

The increasing amount of available interferometric data requires accurate stellar models to convert the raw visibilities into meaningful angular diameters. This is particularly true in the case of Cepheids, whose radial pulsations add another degree of complexity in modeling their spectral emission. The Cepheid variability is being used as a means to accurately determine their distance, and thus the scale of the entire universe; at the same time there is a tendency to neglect their pulsations when computing the properties of their atmospheres, as they change during the Cepheid luminosity cycle.

We presented here a new approach that is used to derive time dependent center-to-limb brightness distributions for classical Cepheids. The adoption of second-order accurate hydrodynamical calculations in spherical geometry, and the inclusions of a complete set of atomic and molecular opacities, improve our limb darkening profiles over previous computations. Conventional limb darkened interferometric measurements are done by adopting a suitable stellar model having the same spectral type and average $`T_{\mathrm{eff}}`$ of the Cepheid (see e.g. Nordgren et al., 1999). In the past, however, none of the available models, has been specifically computed for a Cepheid, and the variations of $`T_{\mathrm{eff}}`$ and $`\mathrm{log}g`$ as the star pulsates are usually not taken into account (the case of $`\mathrm{log}g`$ is particularly dire, given the difficulty of a precise measurement of this parameter, and $`T_{\mathrm{eff}}`$ may suffer from the ambiguities in its definition).

These limitations are absent in our models, which can be computed for individual Cepheids; the only requirements are an accurate model for the pulsational piston (see Buchler et al., 1990) and the temporal variation of $`T_{\mathrm{eff}}`$. Both conditions can be satisfied for most of the nearby Cepheids, and constrained with available high resolution spectral observations.

The main differences between our dynamic models and conventional static limb darkening profiles are analyzed in detail in the previous section. However, it is important to remark once more that the amount of limb darkening between static and dynamic models is significantly different, and can be either larger or smaller according to the pulsational phase. The actual shape of the limb profile is also a function of phase, and cannot be reproduced by a single family of static models, or a single limb darkening parameter.

These results highlight the importance of hydrodynamic calculations in the determination of center-to-limb brightness distribution for pulsating Cepheids. In a future paper we will demonstrate that the effects here described are already measurable with the current interferometric technology, and affect the determination of Cepheid angular diameters. This issue is becoming very important when measuring the variations of such diameters as the star pulsates, which is the first step to provide a reliable determination of its distance through the geometric Baade-Wesselink method.

We wish to thank Bob Kurucz for valuable discussions and the anonymous referee for comments which helped us to improve the paper. This work was partially supported by NSF grant AST9876734. M.K. is a member of the Chandra Science Center, which is operated under contract NAS8-39073, and is partially supported by NASA.


text:
# Searching and studying clusters with the SZ effect

## 1. Introduction

Galaxy cluster catalogs are useful for cosmological studies for a number of reasons: cluster abundance, distribution and evolution reflect the structure formation process, which depends sensitively on certain cosmological parameters, such as the matter density $`\mathrm{\Omega }_\mathrm{M}`$ (Oukbir & Blanchard 1992, 1997; Bahcall 1998; Bartlett 1997; Borgani et al. 1999 & these proceedings; Eke et al. 1998; Henry 1997; Viana & Liddle 1999); and cluster internal properties provide clues to the nature of dark matter and on certain mechanisms associated with galaxy formation, such as heating and feedback.

The construction of well–defined catalogs appropriate for statistical studies demands good understanding of the cluster selection criteria (see Postman, these proceedings). These criteria in turn depend on the type of observations – optical, X–ray, radio – as well as on the data analysis proceedures. Each type of observation has its own advantages and disadvantages that we may refer to as its intrinsic biases. For example, X–ray cluster selection eliminates the problem of projection effects suffered by optical observations, but on the other hand is significantly influenced by the rather poorly understood physics taking place in the core region.

The question that arises is, what kind of selection is best? The answer clearly depends on the utlimate use of the catalog, and for cosmological studies we may argue that the most ideal selection would be based on cluster mass. This is because cluster mass and redshift are the two primary theoretical cluster descriptors. From this line of argument, understanding cluster selection criteria is reduced to relating the observational selection criteria to cluster mass and redshift. This is a useful viewpoint, for it provides a framework for quantifying the biases associated with each type of observation for the purposes of cosmological studies. We do not argue that one particular observational method is the best; rather, that the important point is to use several methods with different biases.

The Sunyaev–Zel’dovich (SZ) effect (Sunyaev & Zel’dovich 1972) will soon open the way to a new method of surveying for clusters. Several ground–based instruments optimized for SZ surveying are currently under construction, and the Planck satellite will supply an almost full–sky catalog of several $`10^4`$ clusters at the turn of the decade. In anticipation and in the spirit just outlined, we shall now examine in more detail this new kind of survey (the ‘SZ band’), and in particular the nature of its cluster section.

This potential of surveying for clusters with the SZ effect has been appreciated for some time. To my knowledge, the first to give source count predictions were Korolyov, Sunyaev, & Yakubstev (1986), and many authors have since considered the issue. An indicative, useful, though by no means exhaustive, list of references focused on SZ surveying (i.e, counts and redshift distributions of SZ detected clusters) would include: Bond & Meyers (1991), Bartlett & Silk (1994), Markevitch et al. (1994), Barbosa et al. (1996), Eke, Cole, & Frenk (1996), Aghanim et al. (1997), Colafrancesco et al. (1997), da Silva et al. (2000), Holder et al. (2000), Lo et al. (2000), Benson, Reichardt, & Kamionkowski (2001), Fan & Chuieh (2001), Kay, Liddle, & Thomas (2001), Kneissl et al. (2001), Springel, White, & Hernquist (2001), Xue & Wu (2001)111The reader will appreciate the precarious position of the author in giving such a list. In offering my apologies to those that I have unjustly, but without malice, forgotten, I ask no more than to be corrected by a kind email in reminder..

## 2. SZ Cluster Selection

In the non–relativistic limit (which applies in practice to all but the hottest clusters), the SZ effect may be written as

$`\delta i_\nu (\widehat{n})`$ $`=`$ $`y(\widehat{n})j_\nu `$
$`y`$ $`=`$ $`{\displaystyle _{los}}{\displaystyle \frac{kT}{mc^2}}n\sigma _T𝑑l`$
where $`\delta i_\nu `$ is the brightness change induced in the cosmic microwave background (CMB) by compton scattering off electrons (mass $`m`$) in the hot intracluster medium (ICM) at temperature $`T`$ and electron density $`n`$. The Compton $`y`$–parameter measures the amplitude of the energy transfer to the photons ($`\sigma _T`$ is the Thompson cross section); and the function $`j_\nu `$ describes the spectral form of the distortion. It is convenient to consider the effect integrated over the entire cluster face:

$$Y=𝑑\mathrm{\Omega }y(\widehat{n})f_{\mathrm{gas}}\frac{MT}{D_{\mathrm{ang}}^2}$$
(1)
where now $`M`$ is the cluster (virial) mass, $`f_{\mathrm{gas}}`$ is the ICM gas mass fraction, $`T`$ is the (true) mean electron temperature (i.e., particle weighted), and $`D_{\mathrm{ang}}(z)`$ is the angular–size distance. Notice that $`Y`$ is measured in, for instance, square arcmins.

It is from these equations that we begin to understand the nature of SZ cluster selection (or detection). Let us emphasize four aspects:

1. The spectral signature $`j_\nu `$ is unique – negative at low frequency, positive at high frequency, with a fixed zero point at $`\lambda 1.4`$ mm – and the same for all clusters ($`j_\nu `$ depends only on the observation frequency). There is therefore no k–correction;
2. The integrated effect $`D_{\mathrm{ang}}^2`$, which is the statement that the surface brightness is independent of redshift $`z`$. This is in contrast to other emission mechanisms, that vary as the luminosity distance $`D_{\mathrm{lum}}=(1+z)^2D_{\mathrm{ang}}`$, and which suffer from ‘cosmic dimming’ (the redshift);
3. The integrated signal is directly proportional to the total ICM thermal energy $`f_{\mathrm{gas}}MT`$ (this is essentially the definition of the mean temperature $`T`$). We may therefore expect a tight correlation between the measurable, $`Y`$ and the cluster mass, e.g., $`Yf_{\mathrm{gas}}M^{5/3}`$, especially in the more massive systems where non–graviational effects appear less important.
4. Cluster detection and subsequent modeling may be based on the same measurable quantity $`Y`$.

A particularily useful quantity for characterizing the selection properties of a survey is the minimum detectable mass as a function of redshift, $`M_{\mathrm{det}}(z)`$. For example, it determines the source counts as

$$\frac{dN}{d\mathrm{\Omega }}=_0^{\mathrm{}}𝑑z\frac{dV}{dzd\mathrm{\Omega }}_{M_{\mathrm{det}}(z)}^{\mathrm{}}𝑑M\frac{dn}{dM}(M,z)$$
(2)
where $`dn/dM`$ is the mass function. A comparison of $`M_{\mathrm{det}}`$ from all–sky X–ray and SZ (the Planck curve) surveys is made in Figure 1. The X–ray survey flux limit corresponds to that of the RASS. An interferometer–like similar to those proposed (see below) and capable of surveying several square degrees from the ground is also shown. To make the figure, a non–evolving $`L_xT`$ relation of the form given by Arnaud & Evrard (1999) and a $`TM`$ relation from Evrard, Metzler, & Navarro (1996) have been adopted. The cosmological model is an open, low–density model, and Eq. (1) was used with a constant gas mass fraction $`f_{\mathrm{gas}}=0.06h^{1.5}`$ (Evrard 1997)

The figure is illustrative, especially because we have no idea of the evolution of the $`L_xT`$ relation to such redshifts. It is nevertheless instructive to note the different behavior of the X–ray and SZ curves, and in particular the almost flat asymptote of the SZ curves; in fact, they demonstrate a slow decrease towards high redshift (this is independent, of course, of the $`L_xT`$ relation). In other words, SZ selection corresponds essentially to a mass selection. The conclusion depends in detail on the constancy of $`f_{\mathrm{gas}}`$ with redshift, but the potential consequences are interesting: for example, this implies that there is no $`V_{\mathrm{max}}`$ for a SZ detected cluster. Even more significant is that SZ selection therefore finds the same kinds of objects (mass–wise) over a large range of redshifts. For evolutionary studies of clusters, this is an important advantage over other selection techniques that force one to compare ever more massive objects detected at high redshift to less massive local samples.

## 3. Aspects of Surveying

The discussion so far as assumed that clusters are unresolved – that the integrated SZ flux of Eq. (1) is directly measured. This situation applies in practice to the Planck survey, with angular resolution of at best 5 armins. Ground–based instruments, on the other hand, will operate at higher resolution and clusters will appear extended, complicating the selection criteria: some clusters may be resolved out, relative to a hypothetical point–source approximation. Detection depends now not only on $`Y`$, but also on the detailed SZ cluster profile; and photometry becomes an issue, as with any extended source.

Many ground–based instruments will operate, at least at first, at a single frequency. In this case, primary CMB fluctuations (in the background) become an important source of noise, and one must extract extended sources over a fluctuating background. Another potentially important effect is confusion by point sources, either randomly distributed on the sky, or preferentially in cluster cores. Without spectral information, their removal requires higher resolution imaging; for this reason, some ground–based interferometer projects plan to control this confusion with additional telescopes placed on long baselines.

There are certainly other issues, for example other foregrounds, that we could list as problematic. However, I would instead like to address a specific question that often arises when preparing a survey: In a given amount of observing time, should one go (integrate) ‘deep’, or instead ‘wide’, covering more sky to less depth? Actually, the real approach is to try to do several surveys, i.e., in some sense to do both. All the same, it is useful to try to answer the question in terms of maximizing the number of detected objects. This depends on the slope of the source counts as a function of survey sensititivity222For resolved sources, these counts are not necessarily the same as the counts as a function of integrated cluster signal $`Y`$..

### 3.1. Deep or Wide?

We calculate the counts using Eq. (2). For unresolved sources, the detection mass is only a function of sensitivity, $`Y_{\mathrm{noise}}`$, and redshift – $`M_{\mathrm{det}}(Y_{\mathrm{noise}},z)`$. Additional factors are important in the case of resolved source detection, such as angular resolution ($`\mathrm{\Omega }_{\mathrm{pix}}`$, the solid angle of a ‘pixel’) and source profile: $`M_{\mathrm{det}}(Y_{\mathrm{noise}},z,\mathrm{\Omega }_{\mathrm{pix}},profile)`$. A simple analytic model of extended source detection in a map may be employed to examine cluster selection (Bartlett 2000). It supposes that one can identify and then integrate over all the pixels covering a cluster (out to the virial radius) to access the total cluster signal $`Y`$ for detection; it is therefore the best one can hope to do.

Using this approach, the effects of resolution on $`M_{\mathrm{det}}`$ are presented in Figure 1, where the dashed lines (for the SZ curves) give $`M_{\mathrm{det}}(z)`$ once cluster extension has been taken into account. Clearly, clusters at low/intermediate redshifts are resolved out relative to the point–source calculation, because they are covered by several noisy pixels; this is most evident for the higher resolution interferometer–like curve. The difference for a low resolution survey, such as the Planck survey, is relatively minor, as can been seen.

This point is important when it comes to predicting the number of objects to be expected. Clearly the predicted number will be lower when correctly accounting for resolution effects. Even more relevant to the question at hand (deep or wide), is the fact that, as seen in Figure 2 (and in Figure 4 below), the actual counts as a function of sensitivity are steeper than the point–source calculation would suggest. In fact, at the bright end they are steeper than $`Y_{\mathrm{noise}}^2`$, the critical slope for answering our question. Counts steeper than $`2`$ imply that more objects will be found by integrating deep on a small patch of sky (as long as $`Y_{\mathrm{noise}}t^{1/2}`$, which usually applies). At the faint end, the counts flatten out and rejoin the point–source calculation, which is always flatter than $`Y_{\mathrm{noise}}^2`$. Two important conclusions follow:

1. the point–source calculation gives the wrong answer in suggesting that one should go wide;
2. the correct answer is found by accounting for resolution effects and indicates that, for each experimental set–up, there is an optimal sensitivity corresponding to the point where the counts reach a slope of $`2`$.

## 4. Expectations

A number of dedicated instruments optimized for SZ surveying have recently been proposed. At the time of writing, many of these projects are funded and under construction, including:

* Interferometer arrays

+ Arcminute MicroKelvin Imager (AMI)
    (http://www.mrao.cam.ac.uk/telescopes/ami/index.html)
    + Array for Microwave Background Anisotropy (AMiBA)
    (http://www.asiaa.sinica.edu.tw/amiba/)
    + SZ Effect Imaging Array
    (http://www.astro.uiuc.edu/ jmohr/SZE/index.html)
* Bolometer cameras

+ ACBAR
    (http://cfpa.berkeley.edu/ swlh/research/acbar.html)
    + BOLOCAM
    (http://www.astro.caltech.edu/ lgg/bolocam/bolocam.html)
* Planck satellite
(http://astro.estec.esa.nl/Planck/)

Several authors have recently made predictions for these various instruments (Aghanim et al. 1997; Holder et al. 2000; Lo et al. 2000; Benson et al. 2001; Kay et al. 2001; Kneissl et al. 2001). In Figures 3 & 4 I show some expectations for the Planck Survey and for the AMI interferometer, both for two extreme cosmologies: a critical model and a purely open model with $`\mathrm{\Omega }_\mathrm{M}=0.3`$. A flat low–density model with non–zero cosmological constant will fall between these two cases, although somewhat closer to the purely open model.

In the following, it should be kept in mind that the predicted total number of detected clusters is subject to large modeling uncertainties and can vary by as much as a factor of $`5`$. A particularly important source of uncertainty comes from the uncertainty in the amplitude of the density perturbations $`\sigma _8`$. To be conservative, the predictions are given for ‘low $`\sigma _8`$’: $`(\sigma _8=0.8,\mathrm{\Omega }_\mathrm{M}=0.3)`$ and $`(\sigma _8=0.5,\mathrm{\Omega }_\mathrm{M}=1)`$ (e.g., Eke et al. 1998). Use of higher values quoted in the literature (e.g., Blanchard et al. 2000; Pierpaoli, Scott, & White 2001) would substantially increase the predicted total number of clusters. Another source of uncertainty lies in the potential evolution of the cluster gas mass fraction $`f_{\mathrm{gas}}`$. The lack of evidence for strong evolution of $`f_{\mathrm{gas}}`$ with redshift or for an important dependence on cluster mass may simply reflect a lack of sufficient data. In this light, it is important to note that the SZ signal for unresolved clusters (e.g., Planck detections) depends only on the global quantity $`f_{\mathrm{gas}}`$, and not on the details of the gas distribution; X–ray observations, as a contrast, are heavily influenced by the physics of the core, which is very poorly understood and difficult to model.

The Planck predictions are given in terms of redshift distributions: shown as a function of redshift $`z`$ are the total number of clusters with $`Y>10^4`$ armin2 at redshifts larger than $`z`$333i.e., the integrated counts from $`z=\mathrm{}`$ down to $`z`$, as a function of $`z`$ and at the fixed sensitivity of $`10^4`$ arcmin2.. The numbers at $`z=0`$ therefore give the total number of detectable clusters, and we see that on the order of 1 cluster/sq.deg. is expected for the survey, depending on the cosmology. At $`z1`$ the two models differ wildly. There are essentially no clusters expected in the critical model. On the other hand, about 10% of clusters are expected to lie at $`z>1`$ in the open scenario. These conclusions assume that $`f_{\mathrm{gas}}`$ remains relatively constant with redshift and mass. The presence of an important population of clusters at high redshift would be a strong indication of a low–density universe; but the lack of such a population could be accommodated in the open model by decreasing $`f_{\mathrm{gas}}`$ with $`z`$. The result in any case would be an important piece of cosmological information.

The predictions for an AMI–like instrument (resolution of a few arcminutes or better) are given in terms of integrated source counts as a function beam sensitivity. Both cosmological models are once again shown, each one for two different low–mass cutoffs ($`10^{13}`$ and $`10^{14}`$ solar masses) in order to separate the contribution of low and high mass objects. The dashed curves show the equivalent counts in a point–source approximation, and we observe that the latter fails rather badly because most clusters are resolved at this resolution. The dotted lines in the background follow the critical slope of $`Y^2`$. This kind of instrument reaches its optimal sensitivity at $`few\times 10^5`$ arcmin2/beam, where the predicted counts attain a slope of $`2`$. At this sensitivity, one may expect, as a conservative round number and within the modeling uncertainties mentioned, on the order of several clusters/sq.deg.

## 5. Conclusion

It will shortly be possible to perform purely SZ–based surveys, opening up what may be refered to as the ‘SZ band’ for cluster studies. A number of ground–based instruments will realize such surveys within a few years over patches of sky covering several square degrees; and by the end of the decade, the Planck mission will produce a catalog of several $`10^4`$ SZ clusters over the majority of the sky. Although predictions for the number of clusters expected for these surveys suffer from modeling uncertainties444Which is not surprising: if there were no uncertainties the surveys would be of little value!, various studies indicate roughly on the order of $`1`$ cluster/sq.deg. for Planck and $`few`$ clusters/sq.deg. for the deeper ground–based efforts under way.

These surveys represent a new method for studying clusters that will find its place alongside X–ray and optical methods. As a population, clusters provide a unique statistical sample for probing cosmology and the large–scale structure formation process. Each observational approach to surveying for clusters has its proper advantages and disadvantages (collectively, its biases). In light of cosmological studies, it is most instructive to quantify observational biases in terms of cluster mass and redshift, the two primary theoretical cluster descriptors. The specific advantages inherent to SZ surveying are listed above under Section 2 and include the fact that the integrated SZ signal is proportional to a physical quantity – the ICM thermal energy.

For a more complete understanding of the cluster population in view of its application as a cosmological probe it is essential to combine a variety of surveying (catalog selection) methods with different biases. SZ surveying will soon add an important new, independent and complementary approach in this endeavor.

I would like to thank the organizers for their invitation to a most interesting and enjoyable meeting.

## References

Aghanim, N., de Luca, A., Bouchet, F.R., Gispert, R., & Puget, J. L. 1997, A&A, 325, 9

Arnaud, M., & Evrard, A.E. 1999, MNRAS, 305, 631

Bahcall, N., & Fan, X. 1998, ApJ, 504, 1 &

Barbosa, D., Bartlett, J.G., Blanchard, A., & Oukbir, J. 1996, A&A, 314, 13

Bartlett, J.G. 1997, in ASP Conf. Ser. Vol. 126, From Quantum Fluctuations to Cosmological Structures, ed. D. Valls–Gabaud et al. (San Francisco: ASP) 365

Bartlett, J.G. 2000, astro–ph/0001267

Bartlett, J.G., & Silk, J. 1994, ApJ, 423, 12

Benson, A.J., Reichardt, C., & Kamionkowski, M. 2001, astro–ph/0110299

Blanchard, A., Sadat, R., Bartlett, J.G., & Le Dour, M. 2000, A&A, 362, 809

Bond, J.R. & Meyers, S.T. 1991, in Trends in Astroparticle Physics, ed. D. Cline (Singapore: World Scientific)

Borgani, S., Rosati, P., Tozzi, P., & Colin, N. 1999, ApJ, 517, 40

Colafrancesco, S., Mazzotta, P., Rephaeli, Y., & Vittorio, N. 1997, ApJ, 433, 454

da Silva, A.C., Barbosa, D., Liddle, A.R., & Thomas, P.A. 2000, MNRAS, 317, 37

Eke, V.R., Cole, S., & Frenk, C.S. 1996, MNRAS, 282, 263

Eke, V.R., Cole, S., Frenk, C.S., & Henry, P.J. 1998, MNRAS298, 1145

Evrard, A.E. 1997, MNRAS, 292, 289

Evrard, A.E., Metzler, C.A., & Navarro, J.F. 1996, ApJ, 469, 494

Fan, Z. & Chiueh, T. 2001, ApJ, 550, 547

Henry, J.P. 1997, ApJ, 489, L1

Holder, G.P., Mohr, J.J., Carlstrom, J.E., Evrard, A.E., & Leitch, E.M. 2000, ApJ, 544, 629

Kay, S.T., Liddle, A.R., & Thomas, P.A. 2001, MNRAS, 325, 835

Kneissl, R., Jones, M.E., Saunders, R., Eke, V.R., Lasenby, A.N., Grainge, K., & Cotter, G. 2001, astro–ph/0103042

Korolyov, V.A., Sunyaev, R.A., & Yakubtsev, L.A. 1986, Sov. Astron. Lett., 12, L141

Lo, K.Y., Chiueh, T.H., Martin, R.N., Kin–Wang Ng, Liang, H., Ue-li Pen, & Chung-Pei Ma 2000, astro–ph/0012282

Markevitch, M., Blumenthal, G.R., Forman, W., Jones, C., & Sunyaev, R.A. 1994, ApJ, 426, 1

Oukbir, J., & Blanchard, A. 1992, A&A, 262, L21

Oukbir, J., & Blanchard, A. 1997, A&A, 317, 1

Pierpaoli, E., Scott, D., & White M. 2001, MNRAS,325, 77

Springel, V., White, M., & Hernquist, L. 2001, ApJ, 549, 681

Sunyaev, R.A. & Zel’dovich, Ya.B. 1972, Comm. Astrophys. Space Phys., 4, 173

Viana, P.T.P., & Liddle, A.R. 1999, MNRAS, 303, 535

Xue, Y.–J., & Wu, X.–P. 2001, ApJ, 552, 452