text:
# Lyman Continuum Escape from Inhomogeneous ISM

## 1 Introduction

Recent studies suggest that the ionising radiation escaping from galaxies could give a substantial (possibly dominant) contribution to the ultraviolet background radiation (UVB) during a large cosmological time span (Giallongo, Fontana & Madau 1997; Giroux & Shull 1997; Bianchi, Cristiani & Kim 2001). Also, at large redshifts, radiation from the first stellar objects has very likely driven the process of cosmic reionisation (Gnedin & Ostriker 1998; Ciardi et al. 2000; Miralda-Escudé, Haehnelt & Rees 2000; Gnedin 2000; Benson et al. 2000; Ciardi et al. 2001). In spite of this extensive body of work on both aspects, their predictive power is jeopardized by the persisting (theoretical and experimental) ignorance on the value of $`f_{esc}`$, the fraction of hydrogen-ionising photons that escapes from the parent galaxy into the intergalactic medium (IGM). Obviously, this quantity enters the modeling of the UVB and reionisation process and the results depend quite sensibly on the assumptions made about this poorly constrained parameter.

A promising way to make theoretical progresses on this issue is to improve the degree of realism of the modeling and the treatment of physical processes to make predictions that can be directly compared with available local and intermediate redshift data.

Dove & Shull (1994b) initially tackled the problem of determining $`f_{esc}`$ by assuming smoothly varying H I distributions in the Galactic disk. They concluded that about 14% of Lyman continuum (Lyc) photons are able to escape the disk. Dove, Shull & Ferrara (2000) later improved the calculation by solving the time-dependent radiation transfer problem of stellar radiation through evolving superbubbles. Their main result is that the shells of the expanding superbubbles quickly trap ionising photons, so that most of the radiation escapes shortly after the formation of the superbubble. This results in a value of $`f_{esc}`$ roughly a factor of 2 lower than obtained by Dove & Shull (1994b). Additional theoretical works (Ricotti & Shull 2000; Wood & Loeb 2000) have extended the analysis to include high-redshift galaxies, for which the escape of Lyc photons can be modified by geometrical effects (i.e. disk vs. spheroidal systems) and by the higher mean galactic interstellar medium (ISM) density. Generally speaking, the escape fraction decreases as the object virialization redshift or mass become larger; for a typical 2-$`\sigma `$ fluctuation in a CDM model at redshift $`10`$ these studies derive even lower values ($`f_{esc}\stackrel{<}{}1`$%).

Observationally, a wide range of values for the escape fraction has been deduced, but it appears that some data favor larger values than expected from theory. For example, most of the detections of starburst galaxies with the Hopkins Ultraviolet Telescope and the Far Ultraviolet Spectroscopic Explorer (Leitherer et al. 1995; Hurwitz, Jelinsky & Dixon 1997; Heckman et al. 2001) are consistent with $`f_{esc}<10\%`$, although objects have been observed with $`f_{esc}<57\%`$ (Hurwitz, Jelinsky & Dixon 1997). However, absorption from undetected interstellar components could allow the true escape fractions to exceed these upper limits. Bland-Hawthorn & Maloney (1999; see Erratum 2001) used optical line emission data for the Magellanic Stream to derive $`f_{esc}45`$%. Even more puzzling are the recent results (Steidel, Pettini & Adelberger 2001; Haenhelt et al. 2001) who detected flux beyond the Lyman limit (with significant residual flux at $`\lambda <912\AA `$) in a composite spectrum of 29 LBGs at $`z=3.4`$. Also, this implies that at these early epochs galaxies were much more transparent to ionising radiation than at present time, contrary to the trend found by the theoretical works mentioned before.

Can the low values of $`f_{esc}`$ predicted by theory and the high values suggested by recent observations be brought into agreement? In this paper we try to ascertain if the dichotomy can be removed by considering the effects of inhomogeneities in the ISM structure and a more realistic distribution of massive stars, the primary sources of ionising photons. As already pointed out before, essentially all studies (with the partial exception of Wood & Loeb 2000) to date have derived the escape fraction assuming smoothly varying gas density distributions. However, this position is clearly untenable in the light of the large number of observations showing that the ISM in the Galaxy and nearby ones has a hierarchical, very likely fractal, structure (Elmegreen & Falgarone 1996 and references therein). This view is also supported by theoretical works (see for example Norman & Ferrara 1996), which predict a fractal behavior over a large dynamical range. This type of organization naturally arises from turbulent pressure, which plays a crucial role in the various phases of the ISM (see for example Kulkarni & Fich  for the cold neutral HI component and Reynolds  for the warm ionised medium): Norman & Ferrara (1996) showed that the ISM turbulent pressure is roughly 30 times higher than the thermal one. Turbulence is pumped into the ISM primarily by multi-SN explosions; however, reaching the high energy density of present day galaxies requires a considerable fraction of the Hubble time at high redshift. Hence the ISM of these primordial systems, differently from their local counterparts, is likely to be more quiescent and smooth.

For these reasons, it seems worthwhile to examine a model for the escape fraction in which a fractal, turbulence dominated ISM is taken into account. This is what we explore in the rest of the paper.

## 2 Model

In this section we describe the model adopted for the Milky Way (MW) gas density and the stellar distribution, as well as the source emission properties. We use the MW as a reference template, but our results should hold with good approximation for the entire class of MW-type galaxies, as long as the properties of their ISM are similar.

To study the escape of (ionising) Lyman continuum (Lyc) photons from the MW, we have chosen a representative cubic volume of 1 kpc3 located above the disk midplane and centered at the Solar position. This allows us to exploit the high quality data available for the spatial distribution of the ionising OB stars (e.g. Garmany, Conti & Chiosi 1982). Also, the same argument applies to the exquisite detail with which the gas density distribution has been derived (Dickey & Lockman 1990; Reynolds & Haffner 2001 and references therein).

Comparison between the soft X-ray background and absorption line experiments toward stars near the Sun, have suggested that we live in a rarefied, hot ISM cavity, usually known as the Local Bubble, whose average radial extent is $``$ 100 pc (Sanders et al. 1977; Snowden et al. 1990; Sfeir et al. 1999). As the volume of the Local Bubble is only 0.2% of the simulation volume, we neglect for simplicity the effects of such structure on our results.

### 2.1 ISM Density Distribution

In order to compare the effects of inhomogeneities on the escape of Lyc photons with previous results, we have first considered the commonly adopted Gaussian HI density distribution (GDD):

$$n_{HI}(z)=n_0\mathrm{exp}(z^2/2H^2),$$
(1)
where $`z`$ is the height above the midplane of the disk. Such distribution is a solution of hydrostatic equation for the gas in the gravitational field of a disk galaxy, including a dark matter component. The Gaussian distribution that most closely resembles the three-component Dickey-Lockman (1990) vertical distribution has parameters $`n_0=0.367`$ cm-3 and $`H=0.184`$ kpc (Dove, Shull & Ferrara 2000). This gives a total neutral hydrogen column density perpendicular to the disk of $`N_{HI}=5.2\times 10^{20}`$ cm-2.

Observations of molecular clouds have shown that with increasing spatial resolution the gas distribution breaks up into substructures of yet smaller scales. The self-similarity of the relations linking some of the molecular cloud properties (such as the mass-size relation and the distribution functions of size and mass) suggests that the density distribution very likely has a fractal structure induced by turbulence (Falgarone, Phillips & Walker 1991; Elmegreen & Falgarone 1996; Elmegreen 1999). As a more realistic description, we have thus alternatively adopted a fractal model.

A fractal ISM distribution is obtained from hierarchically clustered points, along the prescriptions given by Elmegreen (1997). In brief, starting from a point with coordinates $`(i_0,j_0)=(0.5,0.5)`$, fractals on the $`xy`$ plane were made with $`K=6`$ hierarchical levels by choosing continuous coordinates such that for random numbers $`R_i`$ in the interval , $`i_1=i_0+2(R_10.5)/L^1`$ at the first level, $`i_2=i_1+2(R_20.5)/L^2`$ at the second level, and so on up to level $`K`$; $`L=2`$ is a geometric factor for subdivision of one level into the next. The same procedure was used for $`j`$. At each level, $`N=5`$ new positions were substituted for each position in the previous level, leaving $`N^K`$ positions after $`K`$ levels. At each height $`z`$, a $`L^{2K}=64^2`$ square grid was placed around all the points, and the number of points inside each grid cell was counted. We thus recovered the gas density in each grid cell constraining its average value to be the same as for the Gaussian distribution of eq. 1. The resulting gas distribution has thus a fractal dimension $`D=\mathrm{log}N/\mathrm{log}L=2.32`$, which is consistent with the experimentally derived value (Elmegreen & Falgarone 1996). This procedure is repeated at each of the 64 heights that have discrete values in the range $`0<z<1`$ kpc; this range represents a good compromise between the need to encompass most of the vertical extent of the HI disk and to achieve a sufficient spatial sampling. The clumping factor of this distribution, $`Cn^2/n^2`$, calculated as a function of distance from the midplane, ranges between 4 and 8.

In principle, it would have been possible to use the HI distribution directly obtained from 21-cm line observations. However, resolving power limitations make such maps of little use to study the effects of inhomogeneities. For example, the spatial sampling of the Dwingeloo survey (Burton & Hartmann 1994) is $`0.5^{}`$ (corresponding to $`0.03`$ kpc) but the spectral sampling is 1.03 km s-1, or 0.1 kpc. This is about a factor of 6 lower than the resolution required by our simulation. Surveys with higher spatial resolution (1 arcmin) do exist, but they do not cover a sufficiently large latitude range (the best is probably the Canadian Galactic Plane Survey which covers $`3.5^{}`$ to $`+5.5^{}`$, equivalent to half the simulation box) at roughly the same spectral resolution. Thus, we have restricted our analysis to the previous two model cases. A summary of the adopted density distributions in the various simulation runs is shown in Table 1. For illustration, in Fig. 1 we show a slice of the GDD (left panel) and the FDD (right panel) cut through a plane perpendicular to the Galactic disk.

Measurements of faint optical interstellar emission lines and pulsar dispersion measures, show that in the MW approximately 1/3 of the HI mass is contained in the Warm Ionized Medium (Reynolds 1991a; Reynolds 1993). To account for the presence of such component, we have multiplied the GDD and FDD described above by a factor of 1.3. We assume that initially the gas has a temperature of 100 K, as derived for the Cold Neutral Medium.

### 2.2 Stellar Distribution

Garmany, Conti & Chiosi (1982) compiled a catalog of O-type stars in the Solar neighborhood. The catalog is supposed to be complete to a distance of 2.5 kpc from the Sun. Within this radius, the surface density of such component is about $`\mathrm{\Sigma }_{}=24`$ stars kpc-2. Therefore, we have reproduced the local stellar population by randomly extracting positions and photon ionisation rates for 24 stars in our 1 kpc3 computational volume as follows.

We have adopted for the stars a vertical Gaussian distribution with a scaleheight $`H_{}=63`$ pc, as recently derived by Maíz-Apellániz (2001) from a sample of O-B5 stars from the Hipparcos catalog. By randomly sampling this distribution, we have derived the $`z`$-coordinate of each star; the $`x`$ and $`y`$ coordinates have been selected randomly. All stars above the galactic plane have been included in the simulations.

As the star formation process should preferentially take place in the densest regions of the ISM, we have also produced stellar distributions assuming that the probability of finding a star in a certain region is proportional to a power of the gas density, $`P\rho ^{1.5}`$, i.e. a Schmidt-type law (Kennicut 1998).

Finally, we have also produced stellar distributions with a surface density two times larger than the local one, $`\mathrm{\Sigma }_{}=48`$ kpc-2, to study the dependence of $`f_{esc}`$ on a larger range of ionisation rates. A summary of the adopted stellar distributions is shown in Table 1.

### 2.3 Source Emission Properties

A spectral type and luminosity class has been randomly assigned to each star, using the number frequencies of objects of a given type in the catalog of Garmany, Conti & Chiosi (1982). Within 2.5 kpc from the Sun, about 50% of the stars are of O9.0 - O9.5 types, half of which of luminosity class V. We have then used the stellar models of Schaerer & de Koter (1997) to derive the rate of ionising photons from the given spectral type and luminosity class.

Following this recipe, we have produced 20 different stellar distribution realizations for each case in Table 1, to get some semblance of the scatter among different realizations. The scatter is produced by the random sampling of the number, position and ionising photon rate of individual stars in the simulation volume.

For simplicity, we have adopted a single black-body spectrum for each star. We have assumed a temperature of 40000 K, a mean value for the effective temperatures of O-type stars ($`T_{eff}=3000050000`$ K; Schaerer & de Koter, 1997). The pattern of ionisation does not heavily depend on the detailed shape of the spectrum, as long as most of the ionising photons are emitted close to the HI ionisation limit.

Early B-type stars also produce ionising radiation. Assuming a Salpeter Initial Mass Function (IMF) and the ionising photon rate of Schaerer & de Koter (1997), we have estimated that B0-B0.5 stars contribute on average to less than 10% of the total number of ionising photons. We have run a few simulations to check the influence of B stars. The difference between models with and without B stars are masked by the random scatter in the results. Therefore, we did not include B-type stars in the stellar distributions.

## 3 Radiative Transfer Simulations

To study the propagation of ionising radiation produced by the stars through the given ISM density distribution we use the Monte Carlo (MC) radiative transfer code CRASH (Cosmological RAdiative transfer Scheme for Hydrodynamics), described in Ciardi et al. (2001). The code has been originally developed for cosmological applications, but below we show that it can correctly handle also the rather different conditions prevailing in the interstellar medium of galaxies. For clarity, we briefly summarize the main features of the numerical scheme relevant to the present study.

In the application of a MC scheme to radiative transfer problems, the radiation intensity is discretized into a representative number of monochromatic photon packets. The processes involved (e.g. packet emission and absorption) are then treated statistically by randomly sampling the appropriate distribution function. The 1 kpc3 simulation volume has been discretized in $`64^3`$ grid cells; we consider absorption by HI and dust; we neglect the contribution of He because of the paucity of the He-ionising photons produced. We have estimated the dust optical depth in the optical using the canonical relation between the HI column density, $`N_{HI}`$, and the color excess E(B-V) (Bohlin, Savage & Drake 1978). This has been converted into an optical depth at the ionisation limit by adopting the far-UV parametrization of the mean Galactic extinction law (Fitzpatrick &Massa 1988). We obtain:

$$\tau _d(912\AA )=1.9\times 10^{21}(N_H/\mathrm{cm}^2),$$
(2)
where we have excluded the contribution of dust scattering to the opacity, by multiplying for $`(1\omega )`$; $`\omega =0.4`$ is the average albedo of dust in the UV close to the ionisation limit (Witt et al. 1993). Given the above equation, as the final mean ionisation fraction for run A (run B) is $`x=0.69`$ (0.40) (see Section 4), for the frequencies of interest here, we find that $`\tau _d/\tau _{HI}<2.46\times 10^2`$ ($`1.28\times 10^2`$). Thus, dust contribution to absorption is negligible.

For the range of densities considered here, a number of photon packets equal to $`𝒩_p=5\times 10^7`$ has to be emitted in order to reach numerical convergence.

Differently from Ciardi et al. (2001), here we deal with more than one ionising source. We have found that, if the number of emitted packets is sufficiently high, a more reliable solution of the discretized time-dependent ionisation equation is obtained if the time step $`\mathrm{\Delta }t`$ used in their eq. 10 is not treated statistically. Thus, here $`\mathrm{\Delta }t`$ is the time elapsed since a photon packet has gone through the cell for which the equation is being solved. Only photons escaping from the top side of the cubic region contribute to the escape fraction, while those escaping from side faces or the bottom one are subject to a reflecting boundary condition, thus simulating photons coming from adjacent regions.

The same numerical tests of the code described in Ciardi et al. (2001) have been performed in an interstellar (rather than intergalactic) environment. We do not report them here; instead, in Fig. 2 we show the evolution of ionised regions (in black) produced in the GDD described in the previous Section, by sources with different ionisation rate, located on the disk midplane. The slices show a plane perpendicular to the Galactic disk through each source location. Different columns refer to different source ionisation rates ($`\dot{𝒩}_\gamma =1.8,16,34\times 10^{48}`$ s-1 respectively, from left to right), while from the top to the bottom we show the ISM ionisation structure at times $`t=0.05,1,1.5,5`$ Myr after the source has been turned on. The solid lines correspond to the stationary analytical solution given in Dove & Shull (1994b). The ionisation rates were chosen in order to reproduce the various possible shapes of the ionised regions: spherical, elongated and funnel-like.

To prevent misunderstandings, it is useful here to give our working definition of the escape fraction. We define the global escape fraction, $`f_{esc,g}`$, as:

$$f_{esc,g}(t)\frac{_0^t𝑑t\dot{𝒩}_\gamma ^e(t)}{_0^t𝑑t\dot{𝒩}_\gamma (t)},$$
(3)
where $`\dot{𝒩}_\gamma ^e`$ is the escaping photon rate (i.e. the number of photons reaching $`z1`$ kpc per unit time) and $`\dot{𝒩}_\gamma `$ is the total photon rate production at (simulation) time $`t`$. The instantaneous escape fraction, $`f_{esc,i}`$, is instead given by:

$$f_{esc,i}(t)\frac{\dot{𝒩}_\gamma ^e(t)}{\dot{𝒩}_\gamma (t)}.$$
(4)
In order to derive $`f_{esc}`$, we assume that the stars have a constant ionisation rate and we run the simulation until convergence (i.e. the escape fraction and the ionisation structure do not change) is reached. As we find that, after a time of $`5\times 10^6`$ yr the value of $`f_{esc}`$ remains roughly constant, we stop the simulations after $`10^7`$ yr. This is expected as the recombination time at the midplane is 0.3 Myr.

Following the method described above, we have run 20 simulations for each of the cases indicated in Table 1.

## 4 Results

In Fig. 3 we show the evolution of the global escape fraction, $`f_{esc,g}`$, as a function of the total ionisation rate, $`\dot{𝒩}_\gamma `$, for runs A (filled triangles in Fig. 3a) and B (filled circles in Fig. 3b). In both cases $`f_{esc,g}`$ increases with increasing $`\dot{𝒩}_\gamma `$, although for the FDD with a flatter slope. The reason for this different behavior is that photons in a GDD can escape only if an ionised channel can be produced by the stars themselves. Consequently, a low ionisation rate results in a low escape fraction. In a FDD instead, photons can travel along clear sigthlines through low density ionised channels and escape more easily. Thus, for $`\dot{𝒩}_\gamma \stackrel{<}{}2\times 10^{50}`$ s-1 $`f_{esc,g}`$ is always higher for a FDD, while for larger ionisation rates a GDD becomes more transparent. In both cases, the scatter in $`f_{esc,g}`$ is mainly due to differences in the star position and ionisation rate.

To study the effect of a larger total ionisation rate, we have increased the stellar surface density, as explained in Section 2, obtaining a maximum value $`\dot{𝒩}_\gamma 5\times 10^{50}`$ s-1. These additional runs are represented in Fig. 3 as open triangles (run C) and circles (run D). The scatter is reduced with respect to runs A and B as, increasing the number of stars, their position affects less sensibly the value of $`f_{esc}`$, which instead depends more strongly on the total ionisation rate. In order to check if the curve for the GDD would reach a plateau or rather keep growing (as the analytical derivation by Dove & Shull 1994a suggests), we have increased the total ionisation rate up to a value of $`8\times 10^{50}`$ s-1 (corresponding to a stellar surface density of $`\mathrm{\Sigma }_{}120`$ kpc-2). As we expect the scatter to be reduced, in order to limit the computational time, we have only run 6 different realizations.

The evolution of $`f_{esc,g}`$ is linear over the considered range of $`\dot{𝒩}_\gamma `$ and it can be fitted by the following function:

$$f_{esc,g}=\alpha (\dot{𝒩}_\gamma /10^{50}\mathrm{s}^1)+\beta ,$$
(5)
where for a GDD (FDD) $`\alpha =0.089\pm 0.003`$ ($`0.023\pm 0.003`$) and $`\beta =0.07\pm 0.01`$ ($`0.07\pm 0.01`$). In both cases the fit is represented as a dotted line in Fig. 3. We find that the GDD curve eventually flattens and departs from the fit in eq. 5; on the contrary the results for the FDD are well described by the same fit also for high values of $`\dot{𝒩}_\gamma `$.

Dove & Shull (1994a) derived analytically the value of $`f_{esc}`$ expected in a GDD from a point source of given ionisation rate, located on the disk of the Galaxy. The value they obtained (represented by a dashed line in Fig. 3a) should be compared with the asymptotic value of $`f_{esc,i}`$. However, as for late evolutionary times $`f_{esc,g}`$ approaches $`f_{esc,i}`$ (see below), we can compare their results with our runs A and C directly from Fig. 3a. Clearly, the evolution of $`f_{esc,g}`$ differs from the analytic calculation and higher values for the escape fraction are obtained, especially for large ionisation rates. This is due to the fact that in our case the ionisation rate is distributed among various sources located at different heights above the midplane.

In Fig. 4 the time evolution of the escape fraction is shown, for a Gaussian (lower curves) and a fractal (upper curves) density distribution. As a reference, we have chosen a run that produces, for both distributions and given the same set of stars (with $`\dot{𝒩}_\gamma =2\times 10^{50}`$ s-1), a comparable value of $`f_{esc,g}`$. The solid (dotted) lines indicate the evolution of $`f_{esc,g}`$ ($`f_{esc,i}`$). As in a GDD photons can escape only through ionised channels produced by the stars themselves, the photon escape is retarded and the evolution of $`f_{esc}`$ is steeper compared to the case of a FDD. Initially, $`f_{esc,i}`$, which becomes constant after $`10^6`$ yr, is higher than $`f_{esc,g}`$. This is due to the fact that, when the first photons escape (at a time $`t10^3[10^5]`$ yr for a FDD \[GDD\]), a large number of photons has already been emitted and this results in a very low value of $`f_{esc,g}`$. Then $`f_{esc,g}`$ increases, approaching $`f_{esc,i}`$.

In Figs. 5 and 6 illustrative slices extracted from the simulation box for the same run of Fig. 4, show the evolution of the ionised gas (black regions) produced by the sources in the given density field. The top six panels are for a GDD and show the ionised regions at a time $`t=0.1`$ Myr from the source turn on, at different vertical distances from the midplane. The bottom six panels are the same for a FDD. As already pointed out above, although in the GDD the upper gas layers at this time are still neutral, some ionising photons have already escaped the FDD. In Fig. 6 the top (bottom) six panels refer to the GDD (FDD). Here we show the time evolution (columns, from left to right, refer to $`t=0.005,0.5,5`$ Myr respectively) of slices at two heights above the midplane. Again, the ionised regions are more extended in a fractal than in a Gaussian medium.

Although in this work we focus on the variation of $`f_{esc}`$ with the different density distribution, it is interesting to check the values we obtain for the column density of the neutral gas. For runs A and B, we should be able to recover the observed Galactic HI column density ($`N_{HI}=5.2\times 10^{20}`$ cm-2), i.e. about 30% of the total gas should be ionised (Sect. 2.1). For a GDD (run A) we obtain a final mean HI column density of $`N_{HI}=2.1\times 10^{20}`$ cm-2, i.e. about 40% lower; for a FDD (run B), $`N_{HI}=4.1\times 10^{20}`$ cm-2 is closer to the observed value. In both cases, however, some of the simulations have column density closer to the observed. Despite the larger escape fraction for a FDD, the resulting mean HI column density is larger for a fractal medium, where the highest density regions are more difficult to ionise and recombine faster. Also with an increased stellar surface density (runs C and D) the mean $`N_{HI}`$ in a FDD remains as high as $`2.7\times 10^{20}`$ cm-2, while in a GDD it is $`5.9\times 10^{19}`$ cm-2, with values as low as $`10^{17}`$ cm-2. The reduced ionisation fraction in the case of an inhomogeneous medium was also noted in the model for the diffused ionised medium of Miller & Cox (1993), which took into account the absorption of UV radiation by clouds with density larger than the surrounding medium in a statistical way.

In Fig. 7 the vertical distribution of the mean HI (solid lines) and HII (dotted lines) number densities is shown, for runs A and B with the above $`N_{HI}`$, in the case of a Gaussian (top panel) and fractal (bottom panel) density field at different times after the source turn on. For both distributions, the HI profile does not change until $`t10^5`$ yr, when the inner regions, where the stars reside, have an HII density of $`n_{HII}0.1`$ cm-3. The final HI density profile has a similar shape for both distributions, with a slightly lower value in the inner regions for a GDD case. The HII number density increases gradually in a FDD, while in a GDD, as previously mentioned, the ionisation of the upper gas layers is postponed until the stars have been able to produce an ionised channel for the escaping of photons. At the final stages, the HII distribution appears more extended than the HI one. This is consistent with observations of the Reynolds layer, which seems to have a vertical scale height of $`1`$ kpc (Mezger 1978; Reynolds 1991a, 1991b). However, a proper modelling of the diffuse ionised gas should take into account the dynamical effects on gas, which would make the ionised gas to expand, as well as a larger vertical size for the simulation box.

In Fig. 8 we show the evolution of $`f_{esc,g}`$ as a function of the total ionisation rate for run B (filled circles) and E (stars), the only difference being the stellar distribution. The distribution of the escape fraction looks similar in the two cases and the obtained values are comparable: in fact, in run E the stars, following the Gaussian distribution of the gas, have a larger vertical extension, but, as they are located in denser regions, more photons are needed to ionised them and the two effects roughly balance.

Finally, we present in Fig. 9 the spectrum of the emerging ionising radiation, for five realizations of runs A and B, respectively. In both plots, the uppermost curve is the intrinsic spectrum of stars, i.e. a black body with T=40000 K. As already discussed, each simulation is characterised by the rate of emitted ionising photons, $`\dot{𝒩}_\gamma `$. Since we are mainly interested in the shape of the spectrum, we normalise all simulations to the value $`\dot{𝒩}_\gamma =2\times 10^{50}\text{s}^1`$. As expected from the $`\nu ^3`$ frequency dependence of the HI photoionisation cross section, the emerging filtered spectrum is harder than the intrinsic one. For run A, the shape of the spectrum depends on the amount of escaping photons (i.e. on the escape fraction) in each realization: the lower and upper emerging spectra in the left panel of Fig. 9 being those for the lower and higher values of $`\dot{𝒩}_\gamma `$. Different realizations of run B yield much more similar emerging spectra, this reflecting the nearly constant value of the escape fraction. As a reference, we plot in the right panel of Fig. 9 a blackbody with T=45000 K. We remind that here we do not take into account the absorption of photons with $`h\nu >24.6`$ eV by HeI, which may be the cause for the low HeI ionisation observed in the diffuse ionised medium (Reynolds & Tufte 1995).

## 5 Summary and Conclusions

We have studied the effects of gas density inhomogeneities on the escape of ionizing Lyman continuum photons from Milky Way-type galaxies via radiative transfer numerical simulations. To this aim a comparison between a smoothly stratified and an inhomogeneous, fractal distribution have been made with realistic assumptions for the ionizing stellar sources based on available data in the solar neighborhood. The main results obtained can be summarized as follows.

(1) The global escape fraction, $`f_{esc,g}`$, in the case of a Gaussian Density Distribution (GDD) rapidly increases with increasing total ionization rate, $`\dot{𝒩}_\gamma `$, and it flattens for high values of $`\dot{𝒩}_\gamma `$; for a Fractal Density Distribution (FDD) the dependence on $`\dot{𝒩}_\gamma `$ is milder.

(2) For $`\dot{𝒩}_\gamma \stackrel{<}{}2\times 10^{50}`$ s-1 $`f_{esc,g}`$ is always higher for a FDD, while for larger ionization rate a GDD becomes more transparent.

(3) For low stellar surface densities, $`\mathrm{\Sigma }_{}`$, there is a large scatter in the results, due to differences in the star positions and ionization rates. The scatter is reduced for larger $`\mathrm{\Sigma }_{}`$, as the star position affects less sensibly the value of $`f_{esc}`$, which instead depends more strongly on $`\dot{𝒩}_\gamma `$.

(4) In a GDD the photon escape is retarded and the time evolution of the escape fraction is steeper compared to the FDD case.

(5) For a GDD we obtain a final mean HI column density of $`N_{HI}=2.1\times 10^{20}`$ cm-2, lower than the observed one; for a FDD, $`N_{HI}=4.1\times 10^{20}`$ cm-2 closely matches the experimental value. In all cases, the HII distribution appears to be more extended than the HI one.

(6) The value of the escape fraction is more sensitive to the gas density rather than the stellar distribution.

In order to draw a more general conclusion from the results summarized above, it is useful to translate the total ionization rate into a star formation rate for a given galaxy. To do so we use the results of the stellar population synthesis code Starburst99 (Leitherer et al. 1999). Assuming a continuous star formation mode, a Salpeter IMF with $`M_{down}=1M_{}`$ and $`M_{up}=100M_{}`$ and a metallicity equal to 1/20 of solar, we find that the expected ionization rate inside our simulated kpc3 produced by the corresponding stellar population is $`\dot{𝒩}_\gamma =4.5\times 10^{50}(\dot{M}_{}/M_{}\mathrm{yr}^1)`$. From Fig. 3 we then note that for $`\dot{M}_{}2M_{}\mathrm{yr}^1`$ the escape fraction is already equal to roughly 0.6 (0.28) for the GDD (FDD) case.

Star formation rates of a few tens of solar masses are commonly derived from observations of LBGs (Pettini et al. 2001). This result may suggest a galaxy-dominated UV background. In fact, Bianchi, Cristiani & Kim (2001) have recently shown that estimates of the local and high-$`z`$ meta-galactic ionizing flux are consistent with a galaxy-dominated background if $`f_{esc}10\%`$. Hence it appears that these type of sources, as long as they could be modelled as relatively normal disk galaxies, can strongly influence the ionization history of the intergalactic medium and, possibly, the galaxy formation process.

In principle, then, photons can escape rather efficiently (particularly if their ISM is smooth) from high redshift galaxies, provided their vertical gas distribution is not too different from the Milky-Way type galaxies investigated here. Although galactic disk models based on simple semi-analytical prescriptions (Fall & Efstathiou 1980; Mo, Mao & White 1998; Weil, Eke & Efstathiou 1998; Ferrara, Pettini & Shchekinov 2000; Wood & Loeb 2000) predict that the disk column density should increase with redshift by a factor $`10`$ between the present value and redshift 5, this estimate is still subject to several uncertainties. For example it is not clear how both the gas temperature (governing the disk scaleheight) and the fraction of baryons that are able to cool and settle in the disk (determining the midplane gas density) evolve with cosmic time. Even if the column density is larger, we do expect that the galactic ISM is smoother than the present day one, due to the short time interval available for the build up of a fully developed turbulent spectrum and hence a fractal gas distribution. This process might require several (say, 10) eddy turnover times, $`t_e`$, on the scale at which turbulence is injected by superbubbles, i.e. $`\mathrm{}1`$ kpc (Norman & Ferrara 1996). If the characteristic eddy velocity is $`v_t1`$ km s-1, this implies that $`t_e10`$ Gyr is much longer than the Hubble time at redshift 3. Hence, a final conclusion from models like the one presented here extended to encompass high redshift galaxies have to await additional data (both observational and from galaxy formation models) on the distribution and dynamical state of the gas in primordial galaxies.

## Acknowledgments

We would like to thank A. Lorenzani and M. Normandeau for help with the HI density distributions and the referee, K. Wood, for useful comments. This work was partially supported by the Research and Training Network ”The Physics of the Intergalactic Medium” set up by the European Community under the contract HPRN-CT2000-00126 RG29185.


format: text

text:
# 1 Threshold entropy for cooling within the age of the universe. The solid lines show the specific entropy 𝑆=𝑇𝑛_𝑒^{-2/3}, where 𝑛_𝑒 is the electron density, at which gas of temperature 𝑇 will cool within 15 Gyr (red) and 10 Gyr (blue) for elemental abundances of 30% solar. (The usual thermodynamic entropy for an ideal monatomic gas is proportional to ln𝑆.) Note that the core entropy levels in clusters measured at the core radius 0.1𝑟_𝑣 (ref. 5) track the cooling threshold for typical cluster heavy-element abundances of 30% solar. Within that core radius many clusters contain gas whose entropy lies below the cooling threshold, but the luminosity of that relatively small amount of low-entropy gas usually contributes only modestly to the overall luminosity of the cluster.18 The dashed green line indicates the entropy threshold 𝑆_𝑐 used to determine the 𝐿:𝑇 relation in Figure 2.

slugcomment: Published in Nature 22 November 2001

Regulation of the X-ray luminosity of clusters of galaxies by cooling
and supernova feedback

Clusters of galaxies are thought to contain about ten times as much dark matter as baryonic matter1. The dark component therefore dominates the gravitational potential of a cluster, and the baryons confined by this potential radiate X-rays with a luminosity that depends mainly on the gas density in the cluster’s core.2 Predictions of the X-rays’ properties based on models of cluster formation do not, however, agree with the observations. If the models ignore the condensation of cooling gas into stars and feedback from the associated supernovae, they overestimate the X-ray luminosity because the core gas is too high. An early episode of uniformly distributed supernova feedback could rectify this by heating the uncondensed gas and therefore making it harder to compress into the core3-11, but such a process seems to require an implausibly large number of supernovae6,8,12-14. Here we show how radiative cooling of intergalactic gas and subsequent supernova heating conspire to eliminate highly compressible low-entropy gas from the intracluster medium. This brings the core entropy and X-ray luminosities of galaxy clusters into agreement with the observations, in a way that depends little on the efficiency of supernova heating in the early Universe.

Numerical simulations of cosmological structure formation show that the dark-matter potential wells of clusters are quite similar in shape across a wide mass range15. A cluster’s mean mass density within the virial radius $`r_v`$ is roughly 200 times $`\rho _{cr}`$, the critical density for closing the universe, with a dark-matter density distribution $`\rho _{DM}[r(1+cr)^2]^1`$, where $`r`$ is the radius in units of $`r_v`$ and $`c410`$ is a parameter governing the concentration of dark matter toward the centre of the cluster. In the absence of radiative cooling and supernova heating, the gas density in such simulated clusters follows nearly the same profile16, with some additional flattening for $`r1/c`$ .

If all clusters had similar gas density distributions, then cluster temperature would scale as $`T\rho _{cr}r_v^2`$ and X-ray luminosity would scale as $`L\mathrm{\Lambda }\rho _{cr}^2r_v^3`$, where $`\mathrm{\Lambda }`$ is the cooling rate per unit density of the hot gas. At typical cluster temperatures ($`T>2\mathrm{keV}`$), free-free emission dominates the cooling, so $`\mathrm{\Lambda }T^{1/2}`$ and we would expect $`LT^2`$. Instead, the observed power-law index17-19 in the relation $`LT^b`$ is $`b2.62.9`$, meaning that low-temperature clusters and groups of galaxies are far less luminous than expected. This luminosity deficit is also evident in the low level of the unresolved $`1\mathrm{keV}`$ X-ray background, which is an order of magnitude smaller than predicted by cosmological simulations without radiative cooling or supernova heating20-23.

The X-ray luminosities of low-temperature clusters are unexpectedly small because their gas is less centrally concentrated than in hotter clusters, an effect that has been attributed to a universal minimum entropy level in intracluster gas resulting either from supernova heating5,9-11, from heating by active galactic nuclei8, or from radiative cooling8,12,24,25. A lower limit to intracluster entropy steepens the relationship between luminosity and temperature because the shallower potential wells of low-temperature clusters are less able to overcome the resistance to compression owing to that minimum entropy. Both observations5 and theoretical models6,10,11,20,21 have established that a core entropy level corresponding to $`S100200\mathrm{keV}\mathrm{cm}^2`$ can account for the observed slope of the $`L`$:$`T`$ relation and the low level of the $`1\mathrm{keV}`$ background.

This entropy scale emerges most naturally from considerations involving radiative cooling. For any cluster of temperature $`T`$, one can compute the specific entropy level at which a gas parcel would radiate away its thermal energy in a time equivalent to the age of the Universe. Figure 1 shows that locus in the entropy-temperature plane for a typical cluster heavy-element abundance of 30% of the solar value. The locus for this typical abundance lies in the $`100150\mathrm{keV}\mathrm{cm}^2`$ range for $`T<2\mathrm{keV}`$ and rises in proportion to $`T^{2/3}`$ at higher temperatures. Note that the measured core entropies of clusters and groups5 closely track the cooling locus.

The remarkable correspondence between the threshold entropy for cooling ($`S_c`$) and observations of core entropy suggests the following picture for the evolution of the intracluster medium. We can consider the medium to consist of independent gas parcels, each with its own entropy history. The luminosity of a sufficiently relaxed cluster, in which entropy increases monotonically with radius, is then completely determined by the shape of its potential well $`\rho _{DM}(r)`$ and the distribution of specific entropy among its gas parcels. In the context of hierarchical merging, a parcel’s entropy history proceeds like this: each merger event affecting the parcel raises its entropy to some new value $`S_i`$, and as the gas relaxes after the merger, the parcel’s temperature will approach $`T_i`$, the characteristic temperature of its new dark-matter halo. A parcel’s trajectory in the $`ST`$ plane is thus defined by a series of points $`(S_i,T_i)`$ that proceed upward and to the right in a diagram like Figure 1.

Now consider the effect of the cooling locus, which creeps vertically up the diagram as time progresses. If a parcel’s entropy trajectory always remains above this threshold, it will never cool. However, parcels of gas that find themselves below the cooling threshold as structure develops begin to condense and form stars. Shortly thereafter, supernovae heat the neighbouring gas with an efficiency that remains unknown, adding entropy that tends to suppress further cooling and condensation. Stars will continue to form and supernovae will continue to explode until there is no more gas below the cooling threshold; then both processes must cease. Thus cooling and whatever feedback accompanies it act in tandem to eliminate gas parcels with $`S<S_c`$, inevitably creating an entropy floor at the cooling threshold.

This model leads to an $`L`$:$`T`$ relation that closely matches cluster observations. Suppose that both the gas and dark matter density are proportional to $`[r(1+cr)^2]^1`$ in the absence of cooling. Then we can compute the unmodified entropy distribution $`M_g(S)`$, the mass of gas with specific entropy less than $`S`$, and the inverse relation $`S(M_g)`$. When cooling is allowed to operate, gas parcels with $`S<S_c`$ are subject to condensation, and any subsequent feedback is targeted directly at those low-entropy gas parcels. If feedback is inefficient, condensation will remove this gas from the intracluster medium. If feedback is highly efficient, supernova heating will raise the entropy of this gas to $`SS_c`$, and it will convect to the outer regions of the cluster. In either limit, the modified entropy distribution approaches $`S^{}(M_g)=S(M_g+M_c)`$, where $`M_cM_g(S_c)`$. Figure 2 illustrates the $`L`$:$`T`$ relation that we compute when we allow gas with entropy distribution $`S^{}(M_g)`$ to be in hydrostatic equilibrium with the dark matter potential, assuming a Hubble constant of $`100h\mathrm{km}\mathrm{s}^1\mathrm{Mpc}^1`$ with $`h=0.65`$ and a baryon fraction $`0.02h^2\mathrm{\Omega }_M^1`$ with a total matter density $`\mathrm{\Omega }_M=0.33`$. (Using an alternative modified entropy distribution, $`S^{}(M_g)=S(M_g)+S_c`$, leads to similar results.)

The $`L`$:$`T`$ relation we derive from the cooling threshold agrees not only with observations but also with cluster simulations26 that include radiative cooling but no supernova heating. That agreement underlines an important feature of this model: the $`L`$:$`T`$ relation it predicts is insensitive to the efficiency of supernova heating. Any combination of cooling and feedback that truncates the unmodifed entropy distribution at the entropy scale $`S_c`$ will lead to a similar $`L`$:$`T`$ relation. In that sense, the relation determined by the cooling threshold represents the upper envelope to all $`L`$:$`T`$ relations consistent with the age of the universe. Extreme amounts of feedback could further lower the luminosity, but apart from relatively small amounts of low-entropy gas that have not yet completely cooled27, the luminosity of a cluster cannot substantially exceed this relation.

Despite the spherical symmetry of the end state, the cooling, condensation, and feedback processes that regulate the $`L`$:$`T`$ relation are not restricted to the cluster’s core. In hierarchical models of structure formation, the spatial distribution of specific entropy is initially quite complex. If cooling is not allowed to operate, the highest-density, lowest-entropy gas parcels will eventually find their way into the cluster core. However, in the real Universe, this low-entropy gas condenses and forms galaxies long before it reaches the centre of the cluster. Much of the condensation and feedback that ultimately determines a cluster’s luminosity therefore pre-dates the epoch of cluster formation. Only modest quantities of gas remain below the cooling threshold in most present-day clusters, where they collect into centrally-focused cooling flows.

If this interpretation of the $`L`$:$`T`$ relation is correct then the level of the entropy floor in clusters ($`100200\mathrm{keV}\mathrm{cm}^2`$) conveys little information about the efficiency of supernova heating in the early Universe. That information is more directly reflected in the fraction of baryons that have formed stars. A large proportion of the baryons that are now in clusters were associated with objects of less than $`10^{12}`$ solar masses at $`z2`$. At that time, the cooling threshold in a dark-matter halo of $`10^{12}`$ solar masses would have been $`10\mathrm{keV}\mathrm{cm}^2`$, and virtually all the baryons in those haloes were below that threshold28. Some form of feedback, probably supernovae, prevented most of these baryons from condensing. Once feedback drives uncondensed gas to high entropy and low density, merger shocks can more easily maintain it above the cooling threshold. Detailed computations tracing the entropy history of each gas parcel through the series of cooling, feedback, and merger events leading to the final cluster will therefore be needed to relate the supernova heating efficiency to the fraction of baryons now in stars.

Received 2 July; accepted 26 September 2001.

1. Evrard, A. E. The intracluster gas fraction in X-ray clusters: constraints on the clustered mass density. Mont. Not. R. Astron. Soc., 292, 289-297 (1997).

2. Sarazin, C. L. X-ray emissions from clusters of galaxies (Cambridge University Press, Cambridge, 1988)

3. Evrard, A. E. & Henry, J. P. Expectations for X-ray cluster observations by the ROSAT satellite. Astrophys. J., 383, 95-103 (1991).

4. Kaiser, N. Evolution of clusters of galaxies. Astrophys. J., 383, 104-111 (1991).

5. Ponman, T. J., Cannon, D. B., & Navarro, J. F. The thermal imprint of galaxy formation on X-ray clusters. Nature, 397, 135-137 (1999).

6. Balogh, M. L., Babul, A., & Patton, D. R. Pre-heated isentropic gas in groups of galaxies. Mont. Not. R. Astron. Soc., 307, 463-479 (1999).

7. Wu, K. K. S., Fabian, A. C., & Nulsen, P. E. J. The effect of supernova heating on cluster properties and constraints on galaxy formation models. Mont. Not. R. Astron. Soc., 301, L20-L24 (1998).

8. Wu, K. K. S., Fabian, A. C., & Nulsen, P. E. J. Non-gravitational heating in the hierarchical formation of X-ray clusters. Mont. Not. R. Astron. Soc., 318, 889-912 (2000).

9. Cavaliere, A., Menci, N., & Tozzi, P. Hot gas in clusters of galaxies: the punctuated equilibria model. Mont. Not. R. Astron. Soc., 308, 599-608 (1999).

10. Tozzi, P. & Norman, C. The evolution of X-ray clusters and the entropy of the intracluster medium. Astrophys. J., 546, 63-84 (2001).

11. Bialek, J. J., Evrard, A. E., & Mohr, J. J. Effects of preheating on X-ray scaling relations in galaxy clusters. Astrophys. J., 555, 597-612 (2001).

12. Bryan, G. L. Explaining the entropy excess in clusters and groups of galaxies without additional heating. Astrophys. J., 544, L1-L5 (2001).

13. Kravtsov, A. V. & Yepes, G. On the supernova heating of the intergalactic medium. Mont. Not. R. Astron. Soc., 318, 227-238 (2000).

14. Bower, R. et al. The impact of galaxy formation on the X-ray evolution of clusters. Mont. Not. R. Astron. Soc., 325, 497-508 (2001).

15. Navarro, J. F., Frenk, C. S., & White, S. D. M. A universal density profile from hierarchical clustering. Astrophys. J., 490, 493-508 (1997).

16. Frenk, C. S. et al. The Santa Barbara cluster comparison project: a comparison of cosmological hydrodynamics simulations. Astrophys. J., 525, 554-582 (1999).

17. Edge, A. C. & Stewart, G. C. EXOSAT observations of clusters of galaxies. I - The X-ray data. Mont. Not. R. Astron. Soc., 252, 414-427 (1991).

18. Markevitch, M. The $`L_X`$-$`T`$ relation and temperature function for nearby clusters revisited. Astrophys. J., 504, 27-34 (1998).

19. Arnaud, M. & Evrard, A. E. The $`L_X`$-$`T`$ relation and intracluster gas fractions of X-ray clusters. Mont. Not. R. Astron. Soc., 305, 631-640 (1999).

20. Pen, U. Heating of the intergalactic medium. Astrophys. J., 510, L1-L5 (1999).

21. Voit, G. M. & Bryan, G. L. On the distribution of X-ray surface brightness from diffuse gas. Astrophys. J., 551, L139-L142 (2001).

22. Bryan, G. L. & Voit, G. M. The X-ray surface brightness distribution from diffuse gas. Astrophys. J., 556, 590-600 (2001).

23. Wu, K. K. S., Fabian, A. C., & Nulsen, P. E. J. The soft X-ray background: evidence for widespread disruption of the gas halos of galaxy groups. Mont. Not. R. Astron. Soc., 324, 95-107 (2001).

24. Knight, P. A. & Ponman, T. J. The properties of the hot gas in galaxy groups and clusters from 1D hydrodynamical simulations - I. Cosmological infall models. Mont. Not. R. Astron. Soc., 289, 955-972 (1997).

25. Pearce, F. R., Thomas, P. A., Couchman, H. M. P., & Edge, A. C. The effect of radiative cooling on the X-ray properties of galaxy clusters. Mont. Not. R. Astron. Soc., 317, 1029-1049 (2000).

26. Muanwong, O., Thomas, P. A., Kay, S. T., Pearce, F. R., & Couchman, H. M. P. The effect of radiative cooling on scaling laws of clusters and groups. Astrophys. J., 552, L27-L30 (2001).

27. Fabian, A. C., Crawford, C. S., Edge, A. C., Mushotsky, R. F. Cooling flows and the X-ray luminosity-temperature relation for clusters. Mont. Not. R. Astron. Soc., 267, 779-784 (1994).

28. Balogh, M., Pearce, F. R., Bower, R. G., & Kay, S. T. Revisiting the cosmic cooling crisis. Mont. Not. R. Astron. Soc. 326, 1228-1234 (2001)

29. Eke, V.,Navarro, J. F., & Steinmetz, M. The power spectrum dependence of dark matter halo concentrations. Astrophys. J., 554, 114-125 (2001).

30. Helsdon, S. F. & Ponman, T. J. The intragroup medium in loose groups of galaxies. Mont. Not. R. Astron. Soc., 315, 356-370 (2000).

Acknowledgements. We acknowledge R. Bower and M. Balogh for discussions. G.M.V. receives partial support from NASA, and G.L.B. was supported by NASA through a Hubble Fellowship.


text:
# Theoretical developments on Quartic Gauge Boson Couplings at LEPBased on work done in collaboration with G. Montagna, M. Moretti and O. Nicrosini.

## 1 Introduction

Despite the striking success of the Standard Model (SM) in accommodating the precision data collected at high-energy colliders, important tests of the theory, such as the non-abelian nature of the gauge symmetry and the mechanism of electroweak symmetry breaking, are still at a beginning stage. To this end, gauge-boson self interactions play a key role. At present, triple gauge couplings are being probed at LEP  and the Tevatron , while direct measurements of quartic couplings only very recently became available through the study of radiative events at LEP . Actually, events with one or two isolated, hard photons are analysed at LEP, to search for anomalies in the sector of quartic gauge-boson couplings. Only vertices involving at least one photon can be constrained, since quadrilinear interactions containing four massive gauge bosons give rise to a final state of three massive gauge bosons and are therefore beyond the potential of LEP, because of the lack of phase space. The processes considered in the experimental analyses are $`e^+e^{}W^+W^{}\gamma `$, $`e^+e^{}Z\gamma \gamma `$ and $`e^+e^{}\nu \overline{\nu }\gamma \gamma `$ . The $`W^+W^{}\gamma `$ signature, which yields a four-fermion plus gamma final state, is interesting in order to test $`WW\gamma \gamma `$ and $`WWZ\gamma `$ vertices. For the $`Z\gamma \gamma `$ events, the final states due to the hadronic decays of the $`Z`$ boson with two jets and two visible photons are selected to probe the purely anomalous vertex $`ZZ\gamma \gamma `$, which is of particular interest, being absent from the SM at tree level. The final state with two neutrinos and two acoplanar photons allows a study of the quartic $`WW\gamma \gamma `$ and $`ZZ\gamma \gamma `$ interactions.

In the light of these experimental analyses, the aim of the present contribution is to review the state of the art of exact SM calculations for the processes $`e^+e^{}4`$ fermions ($`4f`$) + $`\gamma `$, $`e^+e^{}q\overline{q}\gamma \gamma `$ and $`e^+e^{}\nu \overline{\nu }\gamma \gamma `$, including the effects of quartic anomalous gauge couplings (QAGC).

## 2 Theoretical approach

The theoretical framework of interest is the formalism of electroweak chiral Lagrangians. In such a scenario, QAGC involving four massive gauge bosons emerge as operators of dimension four at next-to-leading order, while QAGC with at least one photon originate from six (or higher) -dimensional operators at next-to-next-to-leading order. They are said genuinely anomalous if they do not induce new trilinear gauge interactions.

Anomalous $`WW\gamma \gamma `$ and $`ZZ\gamma \gamma `$ vertices were originally introduced in ref. . In this paper, the authors show that, by assuming $`C`$ and $`P`$ conservation and further imposing $`U(1)_{em}`$ gauge invariance and $`SU(2)_c`$ custodial symmetry, two independent Lorentz structures contribute to $`WW\gamma \gamma `$ and $`ZZ\gamma \gamma `$ interactions according to the following Lagrangians

$`_0`$ $`=`$ $`{\displaystyle \frac{e^2}{16}}{\displaystyle \frac{a_0}{\mathrm{\Lambda }^2}}F_{\mu \nu }F^{\mu \nu }\stackrel{}{W}^\alpha \stackrel{}{W}_\alpha `$
$`_c`$ $`=`$ $`{\displaystyle \frac{e^2}{16}}{\displaystyle \frac{a_c}{\mathrm{\Lambda }^2}}F_{\mu \alpha }F^{\mu \beta }\stackrel{}{W}^\alpha \stackrel{}{W}_\beta ,`$ (1)
where $`F_{\mu \nu }`$ is the electromagnetic field tensor, and $`\stackrel{}{W}`$ is a $`SU(2)`$ triplet describing the $`W`$ and $`Z`$ physical fields, $`\mathrm{cos}\theta _w`$ being the cosine of the weak mixing angle. In eq. (1) $`a_0`$ and $`a_c`$ are (dimensionless) anomalous couplings, divided by an energy scale $`\mathrm{\Lambda }`$, which has the meaning of scale of new physics. Generally speaking, $`\mathrm{\Lambda }`$ is in principle unknown and model-dependent. However, the ratios $`a_i/\mathrm{\Lambda }^2`$ entering the phenomenological Lagrangians can be meaningfully extracted from the data in a model-independent way.

The anomalous $`WWZ\gamma `$ vertex has been analysed in ref. , showing that, under the assumption of $`C`$, $`P`$ and $`U(1)_{em}`$ invariance, five additional Lorentz structures with respect to eq. (1) contribute. It is further demonstrated that, by embedding all the structures related to $`WW\gamma \gamma `$, $`ZZ\gamma \gamma `$ and $`WWZ\gamma `$ vertices in $`SU(2)\times U(1)`$ gauge-invariant and $`SU(2)_c`$ symmetric combinations, fourteen $`C`$\- and $`P`$\- conserving operators are allowed, with $`k_j^i`$ parameters, which parametrize the strength of anomalous couplings.

By allowing the violation of at least one discrete symmetry, but retaining $`U(1)_{em}`$ invariance and global custodial $`SU(2)_c`$ symmetry, additional terms can be introduced in the Lagrangian. In particular, in the literature three different contributions have been considered: $`_n`$ , which violates $`C`$ and $`CP`$, $`\stackrel{~}{_0}`$ , which violates $`P`$ and $`CP`$, and $`\stackrel{~}{_n}`$ , which violates both $`C`$ and $`P`$, thus conserving $`CP`$.

In the approach of ref. , all the different operatorial structures contributing to the vertices analysed at LEP have been implemented directly at the Lagrangian level in the Monte Carlo codes NUNUGPV  and WRAP  with arbitrary $`a_i`$ coefficients. By means of appropriate relations between the $`a_i`$ parameters, both parametrizations available in the literature for QAGC, namely the parametrization in terms of $`a_0,a_c,a_n`$ couplings and the one in terms of $`k_j^i`$ coefficients, can be obtained, as shown explicitly in ref. .

By means of the above event generators, theoretical predictions for the processes $`e^+e^{}\nu \overline{\nu }\gamma \gamma `$, $`e^+e^{}q\overline{q}\gamma \gamma `$ and $`e^+e^{}4f+\gamma `$ can be obtained. As an option of WRAP, predictions for the inclusive final states $`WW\gamma `$ and $`Z\gamma \gamma `$ can also be obtained, especially in order to compare them with results existing in the literature treating the $`W`$ and $`Z`$ bosons in the on-shell approximation, as those of refs. . Recently the code RacoonWW has been upgraded to include QAGC for the set of $`4f+\gamma `$ processes .

## 3 Discussion

A detailed numerical investigation of the potentialities of the above-mentioned processes in the search for QAGC at LEP can be found in refs. . Here only the main results are summarized.

The options of $`WW\gamma `$ and $`Z\gamma \gamma `$ real production allowed a comparison of the predictions given by the Monte Carlo WRAP with the results already present in the literature. In particular, such a comparison showed a very satisfactory agreement with the numerical results of ref. , while discrepancies were found with the results published in ref. , for the dependence of the $`WW\gamma `$ and $`Z\gamma \gamma `$ cross sections on $`a_0`$ and $`a_c`$ parameters. More precisely, an opposite sign is present for the relative effects of $`a_0`$ and $`a_c`$ on the $`WW\gamma `$ cross section, as noted in ref. , whereas this is not the case for the $`Z\gamma \gamma `$ final state.

An important issue, which can be addressed with complete matrix element calculations, is the reliability of the calculations performed in the narrow-width approximation, used so far in the experimental search for QAGC.

As far as the neutral QAGC sector is concerned, the cross sections obtained with complete calculations for the final states $`q\overline{q}\gamma \gamma `$ and $`\nu \overline{\nu }\gamma \gamma `$ have been compared in ref.  with the $`Z\gamma \gamma `$ approximation as functions of the parameters $`a_0`$ and $`a_c`$ at the c.m. energy of 200 GeV. In the case of the $`q\overline{q}\gamma \gamma `$ channel, a cut on the invariant mass of the jet–jet system around the $`Z`$ mass (80 GeV $`M_{q\overline{q}}100`$ GeV) has been imposed, while for the $`\nu \overline{\nu }\gamma \gamma `$ final state a cut on the recoil mass again around the $`Z`$ mass (80 GeV $`M_{recoil}120`$ GeV) has been adopted, in order to perform a consistent comparison with the $`Z\gamma \gamma `$ approximation. It is worth noting that these selection criteria are also adopted by the real experimental analysis. This comparison allows the effects due to the $`\gamma `$$`Z`$ interference in the $`q\overline{q}\gamma \gamma `$ channel and to $`W`$$`Z`$ interference in the $`\nu \overline{\nu }\gamma \gamma `$ one to be quantified, as well as the effects of the off-shellness of the $`Z`$ boson, in the extraction of limits on the QAGC. The numerical investigation of ref.  shows an agreement between the integrated cross sections at the per cent level as a function of $`a_0`$ and $`a_c`$ variations inside the currently allowed experimental constraints, thus illustrating the reliability of the $`Z\gamma \gamma `$ approximation in view of the expected experimental precision. This conclusion also holds true for the differential distributions mostly sensitive to QAGC and considered in the experimental studies.

A similar analysis has been performed in ref.  for a $`4f+\gamma `$ final state ($`e^+e^{}\mu \overline{\nu }_\mu u\overline{d}\gamma `$), as computed by means of the exact calculation of WRAP, in comparison with the $`WW\gamma `$ approximation considered in the literature . The sensitivity to the anomalous couplings $`a_0,a_c,a_n`$ has been studied for the $`4f+\gamma `$ final state according to two different event selections: no cuts on the invariant masses of the decay products, and cuts on the invariant masses of the decay products around the $`W`$ mass, i.e. 75 GeV $`M_{u\overline{d},\mu \overline{\nu }_\mu }`$ 85 GeV, in order to disentangle, as much as possible, the contributing Feynman graphs with two final-state resonant $`W`$ bosons. The $`WW\gamma `$ approximation predicts a quite different sensitivity with respect to the complete $`4f+\gamma `$ calculation.

By considering variations of the anomalous couplings within the allowed experimental bounds, differences at the ten per cent level are registered between the $`WW\gamma `$ approximation and the $`4f+\gamma `$ prediction, when invariant mass cuts are imposed in the $`4f+\gamma `$ calculation. Notice that, if invariant mass cuts are not considered, the $`4f+\gamma `$ cross section grows up by a factor of 2 with respect to the cross section in the presence of cuts. Therefore, the $`WW\gamma `$ approximation should be employed with due caution in QAGC studies, especially if we take into account that exact $`4f+\gamma `$ generators, such as WRAP  and RacoonWW , incorporate the effects of QAGC and are at our disposal for such experimental studies.

As a last remark, it is worth mentioning that, in realistic experimental analyses, the effects of ISR (properly simulated by means of the available event generators) should be considered, since they tend to diminish the sensitivity on the QAGC at the level of some per cent.

## 4 Conclusions

The search for QAGC in radiative events at LEP demands precise predictions for the processes $`e^+e^{}\nu \overline{\nu }\gamma \gamma `$, $`e^+e^{}q\overline{q}\gamma \gamma `$ and $`e^+e^{}4`$ fermions+$`\gamma `$. To this end, the available exact calculations for such processes, including the contribution of QAGC (and related event generators), have been reviewed.

After a brief description of the possible Lagrangians giving rise to QAGC in vertices involving at least one photon, and of their implementation in existing event generators, comparisons between exact calculations and approximate results existing in the literature have been discussed. It turns out that, for the $`\nu \overline{\nu }\gamma \gamma `$ and $`q\overline{q}\gamma \gamma `$ final states, the $`Z\gamma \gamma `$ approximation works well, the interference effects present in the complete calculations being confined at the per cent level, if appropriate cuts around the $`Z`$-boson mass are required. The $`\nu \overline{\nu }\gamma \gamma `$ final state with appropriate cuts on the recoil mass can be successfully exploited to extract limits on neutral QAGC, as a complementary channel to the $`q\overline{q}\gamma \gamma `$ one. As far as the 4 fermions+$`\gamma `$ final states are concerned, significant differences are seen between the exact calculation and the $`WW\gamma `$ approximation, even in the presence of invariant mass cuts around the $`W`$-boson mass.

The authors are grateful to G. Montagna for his careful reading of the manuscript. F. P. thanks the Conveners of the Parallel Session “Test of the Standard Model” for their kind invitation.