text:
# Semiclassical quantization with bifurcating orbits

## Abstract

Bifurcations of classical orbits introduce divergences into semiclassical spectra which have to be smoothed with the help of uniform approximations. We develop a technique to extract individual energy levels from semiclassical spectra involving uniform approximations. As a prototype example, the method is shown to yield excellent results for photo-absorption spectra for the hydrogen atom in an electric field in a spectral range where the abundance of bifurcations would render the standard closed-orbit formula without uniform approximations useless. Our method immediately applies to semiclassical trace formulae as well as closed-orbit theory and offers a general technique for the semiclassical quantization of arbitrary systems.

The correspondence between atomic spectra and classical orbits has been of fundamental interest and importance since the early days of quantum mechanics. The “old” quantum theory suffered from the severe drawbacks that the Bohr-Sommerfeld quantization rules can only be applied to integrable systems, and, for atomic systems, the Heisenberg principle for matrix elements is silent about transition amplitudes between low-lying and highly excited states. An important step towards a deeper understanding of the intimate connection between classical orbits and the quantum spectra was achieved by the development of periodic orbit theory Gutzwiller (1990) and, as a variant for the photo-excitation of atomic systems, closed orbit theory Du and Delos (1988); Bogomolny (1989). In these theories, the density of states or photo-absorption spectra are given as the sum of two terms, one a smoothly varying part (as a function of energy) and the other a superposition of sinusoidal modulations. The frequencies, amplitudes, and phases of the modulations are directly given in terms of classical parameters of the orbits.

Closed-orbit theory has proven a powerful tool for the semiclassical interpretation of quantum spectra of, e.g., atoms in external magnetic and electric fields by explaining the peaks in the Fourier-transform recurrence spectra – qualitatively and even quantitatively – in terms of the closed orbits of the underlying classical system Main et al. (1994); Courtney et al. (1995); Kips et al. (1999). However, more than a decade after the development of closed-orbit theory the inverse procedure, i.e., the semiclassical calculation of the eigenenergies and transition amplitudes of individual states is still an unsolved problem. The reasons are twofold: Firstly, both closed-orbit and periodic-orbit theory suffer from fundamental convergence problems of the infinite orbit sums. Secondly, in generic systems the orbits undergo bifurcations when the energy is varied, and the semiclassical theories for isolated orbits exhibit unphysical divergences at the bifurcation points. Both problems have been addressed separately: Firstly, the harmonic inversion technique was introduced as a method for semiclassical quantization Main et al. (1997); Main (1999), which allows one to overcome the convergence problems of the closed orbit sum and to extract high-resolution spectra from a finite set of classical orbits. Secondly, in the vicinity of bifurcations the semiclassical approximation for isolated orbits was replaced with a uniform approximation describing all orbits involved in a bifurcation collectively Main and Wunner (1997); Gao and Delos (1997).

In this Letter both ideas are combined for the first time, i.e., we use both uniform approximations and harmonic inversion techniques for the semiclassical calculation of high-resolution spectra. In the presence of uniform approximations, the classical scaling laws which have been essential to all previous applications of the harmonic inversion technique Main (1999) are no longer valid. Therefore, the harmonic inversion method must be generalized to handle the non-scaling functional form of uniform semiclassical approximations. It then gains a degree of flexibility in the quantization of arbitrary systems no other semiclassical quantization scheme has been able to reach to date.

The novel method will be demonstrated by way of example of the hydrogen atom in an electric field. As is well known, the classical dynamics of this system is integrable, which means that semiclassical energy eigenvalues can be calculated with the help of the EBKM torus quantization rules Gao and Delos (1994); Kondratovich and Delos (1997). However, when closed-orbit theory is applied the hydrogen atom in an electric field exhibits properties typical of mixed regular-chaotic systems as, e.g., Rydberg atoms in a magnetic field or Hénon-Heiles type systems. In particular, the closed orbits starting at and returning to the nucleus undergo bifurcations as the energy is varied. Contrary to the torus quantization, the method introduced in this Letter is not restricted to the Stark effect but can be applied to a large variety of systems with chaotic or mixed regular-chaotic classical dynamics. Furthermore, it can be used in connection with periodic-orbit theory Gutzwiller (1990) as well as closed-orbit theory.

The classical dynamics of the Stark system has already been discussed in detail Gao and Delos (1994). For any energy, the electron can go “uphill” against the direction of the electric field until the external field and the Coulomb field make it return to the nucleus. Alternatively, the electron can leave the nucleus in the “downhill” direction of the external field. The downhill orbit is closed only for energies below the Stark saddle point energy, $`E_\mathrm{S}=2F^{1/2}`$, otherwise the electron will cross the Stark saddle and escape to infinity. In addition to these axial closed orbits, there are non-axial orbits returning to the nucleus after $`k`$ oscillations in the downhill direction and $`l>k`$ oscillations in the uphill direction. Each of these orbits is generated in a bifurcation off the downhill orbit at a critical energy $`E_{\mathrm{gen}}`$ and destroyed in a collision with the uphill orbit at $`E_{\mathrm{dest}}>E_{\mathrm{gen}}`$. Outside this energy range, they exist as complex “ghost” orbits.

Closed-orbit theory associates modulations observed in the quantum photo-absorption spectra of Rydberg atoms in external fields with the classical closed orbits. The quantum response function

$$g(E)=\underset{n}{}\frac{|i|D|n|^2}{EE_n+iϵ}=i|DG_E^+D|i,$$
(1)
where $`|i`$ is the initial state, $`D`$ the dipole operator and $`G_E^+`$ the retarded Green’s function, is given as a smooth background plus an oscillatory closed-orbit sum Du and Delos (1988); Bogomolny (1989); Gao et al. (1992)

$$g^{\mathrm{osc}}(E)=\underset{\mathrm{co}}{}𝒜_{\mathrm{co}}(E)e^{iS_{\mathrm{co}}(E)},$$
(2)
where $`S_{\mathrm{co}}`$ is the action of a closed orbit (co) and $`𝒜_{\mathrm{co}}`$ a recurrence amplitude calculated from the monodromy matrix of the orbit and its initial and final directions with respect to the electric field. It includes a complex phase given by the Maslov index. In the following we are using atomic units, with $`\mathrm{}=1`$ and $`F_0=5.14\times 10^9`$ V/cm the unit of the electric field strength.

The most convincing semiclassical interpretation of quantum spectra can be obtained by means of “scaled energy spectroscopy”: By rescaling the classical quantities with suitable powers of, e.g., the electric field strength $`F`$, the classical dynamics can be shown not to depend on the energy $`E`$ and the field strength $`F`$ separately, but only on the scaled energy $`\stackrel{~}{E}=EF^{1/2}`$. When recording quantum states at a fixed scaled energy $`\stackrel{~}{E}`$ as a function of the scaling parameter $`w=F^{1/4}`$, each isolated closed orbit contributes a sinusoidal modulation to the sum (2), which can be extracted by a Fourier transform of the quantum spectrum. Experimental scaled energy spectra of atoms in electric fields have been analyzed in this way Courtney et al. (1995); Kips et al. (1999). The analysis reveals strong evidence for closed orbit bifurcations.

The simple semiclassical approximation embodied in the closed-orbit formula (2) fails close to a bifurcation of closed orbits, resulting in the divergence of the recurrence amplitudes. To overcome this difficulty, the closed-orbit terms for isolated orbits in (2) must be replaced with a uniform approximation describing all orbits involved in a bifurcation collectively. A uniform approximation suitable for regularizing the bifurcation of a non-axial orbit off either the downhill or the uphill orbit was derived by Gao and Delos Gao and Delos (1997) as well as Shaw and Robicheaux Shaw and Robicheaux (1998). We will use a slightly modified version of their result that gives the collective contribution of the axial and non-axial orbits participating in a bifurcation in terms of their actions $`S_{\mathrm{ax}}`$ and $`S_{\mathrm{non}}`$ and recurrence amplitudes $`𝒜_{\mathrm{ax}}`$ and $`𝒜_{\mathrm{non}}`$ as

$$\mathrm{\Psi }(E)=\left[\frac{𝒜_{\mathrm{non}}}{(1+i)}I+\frac{1}{a}\left(a𝒜_{\mathrm{ax}}+\frac{1+i}{\sqrt{2\pi }}𝒜_{\mathrm{non}}\right)\right]e^{iS_{\mathrm{ax}}}$$
(3)
where $`I`$ is given in terms of the standard Fresnel integrals $`C(x)`$ and $`S(x)`$ Abramowitz and Stegun (1984),

$$I=e^{ia^2/4}\left[\frac{1+i}{2}C\left(\frac{a}{\sqrt{2\pi }}\right)iS\left(\frac{a}{\sqrt{2\pi }}\right)\right],$$
(4)
and

$$a=\pm 2\sqrt{S_{\mathrm{ax}}S_{\mathrm{non}}}.$$
(5)
The negative sign for $`a`$ has to be chosen if the non-axial orbit is a complex ghost orbit.

The high-resolution quantization by harmonic inversion Main et al. (1997); Main (1999) is based on the observation that by equating the quantum recurrence function (1) to its semiclassical approximation (2) – the smooth part can be neglected – and taking the Fourier transform we obtain

$$i\underset{n}{}d_ne^{iE_nt}=C(t)$$
(6)
with $`d_n=|i|D|n|^2`$ and

$$C(t)=\frac{1}{2\pi }_{\mathrm{}}^{\mathrm{}}𝑑E\underset{\mathrm{co}}{}𝒜_{\mathrm{co}}(E)e^{iS_{\mathrm{co}}(E)}e^{iEt}.$$
(7)
The quantization problem has thus been recast as the problem of extracting the frequencies $`E_n`$ and amplitudes $`d_n`$ from a given time signal $`C(t)`$ of the form (6), provided the signal (7) can be calculated. In the case of a scaling system the signal is given as a sum of $`\delta `$ functions.

While the uniform approximation (3) successfully smoothes the divergences in (2), it spoils the classical scaling properties Shaw and Robicheaux (1998). Therefore, the Fourier transform of a spectrum including uniform approximations cannot be evaluated in terms of $`\delta `$ functions. In fact, for non-scaling systems, there seems to be no way at all to compute the integral (7) because, apart from the fact that the classical quantities can always be calculated in a finite energy interval only, the integral can in general not even be expected to converge. Therefore, neither by analytical nor by numerical means will one be able to compute a useful semiclassical signal from (7). The inclusion of uniform approximations in semiclassical quantization is thus a nontrivial and challenging task.

To solve the problem we resort to the observation made in Main et al. (2000) that a band-limited signal, which only contains the spectral information describing the quantum system in a finite energy interval $`[E_{\mathrm{min}},E_{\mathrm{max}}]`$, can be obtained by restricting the energy integral in (7) to this window. The resulting signal

$$C^{\mathrm{bl}}(t)=\frac{1}{2\pi }_{E_{\mathrm{min}}}^{E_{\mathrm{max}}}g^{\mathrm{osc}}(E)e^{iEt}𝑑E,$$
(8)
where $`g^{\mathrm{osc}}(E)`$ is in general a non-scaling function that includes both contributions of isolated closed orbits and uniform approximations, can be evaluated numerically by standard methods. As it contains only a small number of frequencies in the interval chosen, it can be processed in a numerically stable way by conventional high-resolution methods such as linear prediction or Padé approximants. Contrary to previous methods for scaling systems, all of which contained the analytic evaluation of an integral, the numerical integration imposes no restrictions on the semiclassical response function occurring in the integrand. Notice that our method cannot be applied in connection with the original filter diagonalization algorithm Main (1999). It is only the separation between a low-resolution frequency filtering stage and a high-resolution harmonic inversion stage introduced in Main et al. (2000) that allows for the present generalization to arbitrary non-scaling semiclassical signals.

To demonstrate our method we investigate Stark spectra of the hydrogen atom for transitions from the ground state $`|1s0`$ to highly excited Rydberg states with light polarized parallel to the electric field axis. The external field strength is $`F=10^8\text{a.u.}\widehat{=}51.4\text{V/cm}`$. The high-resolution Stark spectrum is obtained in two steps:

Firstly, a low-resolution semiclassical spectrum is obtained by truncating the closed-orbit sum (2) at a maximal period $`T_{\mathrm{max}}`$. The cut-off value for the periods should not be chosen smaller than the signal length of the band-limited signal (8) used for the harmonic inversion in the second step of the procedure. The low-resolution spectrum calculated with a scaled truncation time of $`T_{\mathrm{max}}=15\times 10^6`$ is presented in Fig. 1. In the lower energy range shown in Fig. 1(a), individual non-overlapping $`n`$-manifolds can be observed. (We have $`n=30`$ at $`E5.56\times 10^4`$.) In this region, the signal is sufficiently long to resolve individual spectral lines, although their precise determination from the plots remains difficult. In the higher energy range shown in Fig. 1(b), two, three, or even four different $`n`$-manifolds overlap, leading to a drastically increased spectral density. In this region the semiclassical signal is evidently too short to discriminate individual lines. It is important to note that the possibility of computing the low-resolution spectrum depends critically on the use of uniform approximations. If it was calculated from isolated-orbit contributions only, a dense sequence of bifurcation-induced divergences would cover even the large-scale structure of the spectra. With the low-resolution semiclassical spectrum (Fig. 1) at hand the band-limited time signal (8) is now obtained by a numerical Fourier transform of $`g^{\mathrm{osc}}(E)`$. In that calculation we used the signal length $`T_{\mathrm{max}}=40\times 10^6`$ in order to resolve individual levels in the region of overlapping $`n`$-manifolds.

In the second step, the high-resolution semiclassical spectrum is finally obtained by harmonic inversion of the band-limited time signal (8). In Fig. 2 the semiclassical spectrum is compared to the exact quantum spectrum. The overall agreement between the semiclassical and the quantum spectrum is excellent, although for a few levels the comparison reveals discrepancies between the semiclassical and the quantum matrix elements. Note, in particular, the region of high spectral density at $`E3.2\times 10^4`$. In this region, groups of 3 nearly degenerate levels exist, some of which are well resolved semiclassically. At $`E4.4\times 10^4`$, even closer lines exist – they can hardly be discerned in the quantum spectrum. These lines are not resolved semiclassically. Instead, the harmonic inversion yields single lines with amplitudes equal to the sum of the two quantum amplitudes. We are confident to fully resolve even these states in the semiclassical spectrum when applying the cross-correlation technique for harmonic inversion Main (1999).

In the future it should also be possible to extend the semiclassical spectrum to energies $`E>E_\mathrm{S}`$ where the classical motion is not completely bound, and to extract the semiclassical widths of the Stark resonances. A particular challenge is posed by the region around the Stark saddle point energy $`E_\mathrm{S}`$. Before the downhill orbit ceases to exist at $`E_\mathrm{S}`$, it undergoes an infinite sequence of bifurcations, giving birth to non-axial orbits with arbitrarily high uphill repetition numbers. If subsequent bifurcations of a single orbit are too close, the uniform approximation (3) is no longer appropriate. It must then be replaced with a uniform approximation describing several bifurcations collectively. The uniformization of an infinite bifurcation cascade, in particular, remains an open problem whose solution is required to semiclassically cross the saddle point energy.

In summary, we have extended the harmonic inversion approach to semiclassical quantization to the quantization of systems without a scaling property. The generalized method allows for the inclusion of uniform approximations into the quantization procedure. We have demonstrated the effectiveness of our method by calculating a high-quality semiclassical spectrum for the hydrogen atom in an electric field in a spectral region where the semiclassical approximation without uniform approximations would be completely useless due to the abundance of bifurcations. With the modifications presented here, the technique of quantization by harmonic inversion has reached a stage where it does not impose any conditions on the classical dynamics of the system under study except that a semiclassical approximation to the response function can be given. Besides uniform approximations, any other non-standard semiclassical contributions such as diffractive corrections can be included. Thus, the harmonic inversion can now be regarded as a truly universal tool for the semiclassical quantization of arbitrary systems.


format: text

text:
# Spontaneous heavy cluster emission rates using microscopic potentials

## Abstract

The nuclear cluster radioactivities have been studied theoretically in the framework of a microscopic superasymmetric fission model (MSAFM). The nuclear interaction potentials required for binary cold fission processes are calculated by folding in the density distribution functions of the two fragments with a realistic effective interaction. The microscopic nuclear potential thus obtained has been used to calculate the action integral within the WKB approximation. The calculated half lives of the present MSAFM calculations are found to be in good agreement over a wide range of observed experimental data.

Since the first experimental observation of cluster radioactivity , a lot of efforts, both experimental and theoretical, have gone into the understanding of the physics of cluster radioactivity. Lifetimes of the cluster radioactivities of radioactive nuclei have been predicted theoretically using various models and compared with existing experimental data from time to time. These models can be broadly classified as the superasymmetric fission model (SAFM)  and the preformed cluster model (PCM) . In the SAFM the barrier penetrabilities are calculated assuming two asymmetric clusters. In the PCM the cluster is assumed to be formed before it penetrates the barrier and its preformation probability is also included in the calculations. Though the physics of the two approaches is apparently different, but actually they are almost similar. Interpreting the cluster preformation probability within a fission model as the penetrability of the pre-scission part of the barrier, it was shown that the PCM is, in fact, equivalent to the fission model . However, the PCM has been found to be better applicable for lighter clusters while SAFM is more apt for all cluster decays .

Both the theoretical approaches described above use either phenomenological potentials or the proximity type potentials to calculate nuclear interaction between the two fragments. The SAFM calculations using proximity type potentials or semiempirical heavy ion potentials obtained by fitting the elastic scattering data or other phenomenological nuclear potentials for interaction between the fragments do not reproduce the observed cluster radioactivity lifetimes successfully. The SAFM using a parabolic potential approximation for the nuclear interaction potential, which is a rather unusual fragment interaction potential, however, has been found to provide reasonable estimates for the lifetimes of cluster radioactivity . The PCM with various nuclear potentials have also been tried with some succes for the alpha radioactivity but was not much succesful even for a very limited number of heavier cluster decays. In the present work microscopically calculated nuclear interaction potentials have been used in the SAFM approach with reasonable success for calculating the lifetimes of cluster radioactive decays over a wide range of emitted heavy clusters from a large number of parent nuclei. The microscopic nuclear potentential obtained by double folding the cluster density distributions with realistic effective interaction is also very fundamental in nature. Moreover, the use of a single microscopic nuclear potential over a wide range of daughter and emitted cluster interaction is also aesthetically appealing.

In the SAFM the half life of the parent nucleus against the split into a cluster and a daughter is calculated using the WKB barrier penetration probability. The assault frequecy $`\nu `$ is obtained from the zero point vibration energy $`E_v=(1/2)\mathrm{}\omega =(1/2)h\nu `$. The half life $`T`$ of the parent nucleus $`(A,Z)`$ against its split into a cluster $`(A_e,Z_e)`$ and a daughter $`(A_d,Z_d)`$ is given by

$$T=[(h\mathrm{ln}2)/(2E_v)][1+\mathrm{exp}(K)]$$
(1)
where the action integral $`K`$ within the WKB approximation is given by

$$K=(2/\mathrm{})_{R_a}^{R_b}[2\mu (E(R)E_vQ)]^{1/2}𝑑R$$
(2)
Here $`\mu =mA_eA_d/A`$ is the reduced mass, m is the nucleon mass, and $`E(R)`$ is the total interaction energy of the two fragments separated by the distance R between the centres, which is equal to the sum of nuclear interaction energy, Coulomb interaction energy and the centrifugal barrier. The amount of energy released in the process is $`Q`$ and $`R_a`$ and $`R_b`$ are the two turning points of the WKB action integral determined from the equations

$$E(R_a)=E(R_b)=Q+E_v$$
(3)
Energetics allow spontaneous emission of cluster only if the released energy

$$Q=M(M_e+M_d)$$
(4)
is a positive quantity, where $`M`$, $`M_e`$ and $`M_d`$ are the atomic masses of the parent, the emitted cluster and the daughter nuclei, respectively, expressed in the units of energy. Correctness of predictions for possible decay modes therefore rests on the accuracy of the ground state masses of nuclei while the reliability of the half life calculations requires proper zero point vibration energies and nuclear interaction energies.

In the present work the total interaction energy $`E(R)`$ has been evaluated using microscopic nuclear potential along with the Coulomb potential over the entire domain of interaction. The microscopic nuclear potentials have been obtained by double folding in the densities of the fragments with the finite range realistic M3Y effective interacion as

$$V(R)=\rho _1(\stackrel{}{r_1})\rho _2(\stackrel{}{r_2})v[|\stackrel{}{r_2}\stackrel{}{r_1}+\stackrel{}{R}|]d^3r_1d^3r_2$$
(5)
The density distribution used for the clusters has been chosen to be of the spherically symmetric form given by

$$\rho =\rho _0/[1+exp((rc)/a)]$$
(6)
where

$$c=R(1\pi ^2a^2/3R^2),R=1.13A^{1/3}anda=0.54fm$$
(7)
and the value of $`\rho _0`$ is fixed by equating the volume integral of the density distribution function to the mass number of the cluster. The finite range M3Y effective interaction $`v(s)`$ appearing in the eqn. (5) is given by

$$v(s)=7999.\mathrm{exp}(4s)/(4s)2134.\mathrm{exp}(2.5s)/(2.5s)$$
(8)
For the direct part of the M3Y effective interaction the long range one-pion exchange potential is exactly equal to zero. As the cluster decays involve only very low energies, the finite range exchange interaction has not been considered because it is important only at higher energies . This microscopic nuclear potential energy is then used to calculate the total interaction energy $`E(R)`$ for use inside the WKB action integral. The two turning points of the action integral have been obtained by solving eqns.(3) using microscopic double folding potential given by eqn.(5) along with the Coulomb potential. Then the WKB action integral between the two turning points has been evaluated numerically for calculating the half lives of the cluster decays. The zero point vibration energies used in the present calculations are same as that described by eqns.(5) in reference . The shell effects for every cluster radioactivity are implicitly contained in the zero point vibration energy due to its proportionality with the Q value, which is maximum when the daughter nucleus has a magic number of neutrons and protons. A normalisation factor of 0.9 for the microscopic nuclear potential has been used to obtain the optimum fit. The present calculation uses the experimental ground state masses for calculating the released energy Q. Whenever the experimental ground state masses are not available, it uses the theoretically calculated ground state masses from the latest mass table .

It is important to mention here that in the analytical superasymmetric fission model (ASAFM)  calculations, the entire interaction region is divided into two distinct zones. In the overlapping zone, where the distances of separation between the centres of the two fragments are below the touching radius, a parabolic form for the nuclear interaction potential has been used. And for distances beyond the touching radius only the Coulomb potential plus the centrifugal barrier for the separated fragments have been considered within a framework of a liquid drop model \[LDM\] two centre spherical parametrization. Treating the region beyond the touching radius as a nuclear force free zone and approximating the nuclear interaction potential to a parabolic form in the overlapping region yield analytical expression for the WKB action integral . Although the overall uncertainty of this analytical superasymmetric fission model (ASAFM) was found to be small, neither the division of the interaction region into two distinct domains is justifiable nor the use of parabolic nuclear potential has much physical basis.

In Fig. 1, Fig. 2, Fig. 3, and Fig. 4 the experimental data for logarithmic half lives  have been plotted against the mass numbers of parent nuclei along with the results of the present calculations for zero angular momentum of the fragments. In all the figures, the open circles depict the experimental data while the continuous line with solid circle represents the present calculations (MSAFM). The upward arrows to some experimental data points indicate that those are only the lower limits of the decay half lives determined experimentally. Fig. 1 contains the results of the present (MSAFM) theoretical calculations and the data points for carbon-14, oxygen-20 and fluorine-23 cluster emissions. Fig. 2, Fig. 3, and Fig. 4 represent the data and theoretical results of MSAFM calculations for cluster emissions of neon, magnesium and silicon isotopes, respectively. The decay modes and the experimental values for their half lives have been presented in Table 1. Those data that represent only the lower limits for the decay half lives have been placed at the bottom. The corresponding results of the present calculations of superasymmetric fission model with microscopic potentials (MSAFM) are also presented along with the results of ASAFM calculations of 1986  and 1991  so as to facilitate the comparison of the results of older calculations  with the present one.

The results of the present calculations of the MSAFM have been found to predict the general trend very well for a wide range of experimental data. The quantitave agreement with experimental data for lighter cluster emissions is excellent while that for heavier clusters is reasonable. The degree of reliability of the MSAFM predictions for cluster decay lifetimes are comparable to that of ASAFM , although they are not exactly the same. It is worthwhile to mention that all the ASAFM results of 1986 and of 1991 listed in the Table 1 have been recalculated using zero-point vibration energies given by eqn. of reference  and eqns.(5) of reference  respectively.

The half lives for cluster-radioactivity have been analyzed with microscopic nuclear potentials which are based on profound theoretical basis. The results of the present calculations with MSAFM are in good agreement over a wide range of experimental data and are comparable to the best available theoretical calculations  of ASAFM which used parabolic interaction potentials that did not have any microscopic basis. Present calculations certainly put the SAFM on a firm theoretical basis. Refinements such as introduction of dissipation while tunneling through the barrier or incorporating the dynamic shape deformations in the density distributions of the clusters may further improve results. It may, however, be realised that as the first illustrative calculations using realistic microscopic cluster interaction potentials, the results of cluster radioactive decay lifetimes obtained without adjusting parameters are remarkable. In future, such calculations may therefore be extended to provide reasonable estimates of the lifetimes of nuclear decays by cluster emissions for the entire domain of exotic nuclei.

The author is grateful to Dr. A.K. Chaudhury, Dr. K. Krishan, Dr. S. Bhattacharya and Dr. J.N. De for many helpful discussions and suggestions.


text:
# Beam Models for Gamma-Ray Bursts Sources: Outflow Structure, Kinematics and Emission Mechanisms

## 1 Introduction

The discovery of afterglows in recent years (Costa et al. 1997; van Paradijs et al. 1997; Frail et al. 1997) has moved the study of gamma-ray bursts (GRB) to a new plane. It not only has extended observations to longer timescales and other wavelength bands, enabling the measurement of redshift distances, but has also provided a direct determination of the source size and confirmation of its relativistic expansion (see Piran 1999 and Mészáros 2001 for recent reviews). Detailed predictions of the afterglow properties, made in advance of the observations (Mészáros & Rees 1997; Vietri 1997), agreed well with subsequent detections at these wavelengths, followed over timescales of months. While progress has been made in understanding how the GRB and afterglow radiation arises in terms of an unsteady relativistic outflow followed by the development of a blast wave moving into the external medium, interest continues to grow as new observations provided new challenges of interpretation.

Many GRBs are now found to be associated with star forming regions (Kulkarni et al. 1998; Fructher et al. 1999; Berger, Kulkarni & Frail 2001), the remnants of a massive stellar progenitor system (Piro et al. 2000; Amati et al. 2000) and possibly supernovae (Galama et al. 1998; Bloom et al. 1999; Reichart 1999; Björnsson et al. 2001; Lazzati et al. 2001b). These observations give support to the idea that the most common GRBs could be linked to the collapse of massive stars (Woosley 1993; Paczyński 1998; MacFadyen & Woosley 1999). This speculation is made more intriguing by the recent report that GRBs could be highly collimated, so that there may be many more events Doppler boosted in directions other than ours, bringing the GRB rate to within a factor of 100 or so of the supernova rate and their derived explosive bulk energies down to $`10^{52.5}`$ erg (Frail et al. 2001, Panaitescu & Kumar 2001; Piran et al. 2001; assuming a conversion efficiency of 1%). The apparent dispersion in isotropic energies is then caused by a distribution of jet opening angles – bursts with higher energy per solid angle are more collimated (or maybe viewed closer to the jet axis; see Rossi, Lazzati & Rees 2002; Salmonson & Galama 2002; Zhang & Mészáros 2002). New evidence and new puzzles were added when the total isotropic luminosity was found to be also related with a variety of temporal and spectral behaviors exhibited by $`\gamma `$-ray light curves - namely the degree of variability or “spikiness” (Fenimore & Ramirez-Ruiz 2002a), the differential time lags for the arrival of burst pulses at different energies (Norris, Marani & Bonnell 2000, Norris 2002), and the rest frame GRB peak energy (Lloyd-Ronning & Ramirez-Ruiz 2002).

These correlations are still tentative, but if confirmed they could be used to relate the properties of material moving almost directly towards us (e.g. the variable activity of the $`\gamma `$-ray light curve) with the ejection in directions away from our line of sight (e.g. degree of collimation; see Salmonson 2000; Ioka & Nakamura 2001; Kobayashi et al. 2002; Plaga 2001). This could turn out to be very important, since at present we can only infer the energy per solid angle; as yet the constraints on the angle-integrated $`\gamma `$-ray energy are not strong. The gamma-rays we receive come only from material whose motion is directed within $`1/\mathrm{\Gamma }`$ of our line of sight. They therefore provide no information about the ejecta in other directions. At observer times of more than a week, the beaming and aberration effects are less extreme so we observe afterglow emission from a wider range of angles. At these later times, however, all memory of the initial time-structure of the outflow responsible for the $`\gamma `$-ray emission would be lost.

Information regarding how much energy and momentum has been injected, its distribution in angle, the mass fractions in shells with different Lorentz factors and its efficiency in radiated gamma-rays are at the forefront of attention. We address some aspects of these issues here by studying the physical conditions under which simulated GRBs arising from internal shocks in relativistic winds can account for this burst phenomenology - degree of collimation, variable $`\gamma `$-ray activity, the rest frame peak energy and the total burst luminosity. To that effect, we consider anisotropic relativistic outflows where changes in the energy per solid angle are caused either by an increase in the bulk Lorentz factor or in the baryon loading. We apply these models to interpret the large variety of behaviors exhibited by both the afterglow and the prompt emission, and discuss their possible use for predicting the internal structure of the collimated outflow. Some possible scenarios in which large variations of the energy per solid angle are naturally expected are highlighted, along with the types of observation that would discriminate among the various models. We assume $`H_0=65\mathrm{km}\mathrm{s}^1\mathrm{Mpc}^1`$, $`\mathrm{\Omega }_{\mathrm{matter}}=0.3`$, and $`\mathrm{\Omega }_\mathrm{\Lambda }=0.7`$.

## 2 Anisotropic or Highly Beamed Outflows?

The high Lorentz factors and energies seen in GRBs are consistent with the catastrophic formation of a stellar black hole of a few $`M_{}`$, with $``$ 1 % going to a relativistic outflow. This could be the extreme example of the asymmetric explosion produced by supernova (Khokhlov et al. 1999), in which instead of halting at the neutron star stage, the collapse continues to the black-hole stage, producing an even faster jet in the process (MacFadyen & Woosley 1999). GRBs arising from a very small fraction of stars that undergo this type of catastrophic energy release are likely to produce collimated outflows.

Even if the outflow is not highly collimated, some beaming is expected because energy would channeled preferentially along the rotation axis. Also, one would expect baryon contamination to be lowest near the axis, because angular momentum flings material away from the axis and material with low-angular momentum falls into the black hole. The dynamics, however, are complex. While numerical simulations of collapse scenarios can address the fate of the bulk of the matter (MacFadyen & Woosley 1999; MacFadyen, Woosley & Heger 2001; Aloy et al. 2000; Zhang et al. 2002, in preparation), higher resolution simulations of the outer layers of the stellar mantle seem to be required since even a very small amount of baryons polluting the outflow could severely limit the attainable Lorentz factor. It is quite possible, for instance, that the stellar pressure will tend to collimate the fireball into a jet. A broad spread of Lorentz factors is thus expected – close to the rotation axis $`\mathrm{\Gamma }`$ may be high. At larger angles away from the jet axis, there may be an increasing degree of entrainment, with a corresponding decrease in $`\mathrm{\Gamma }`$. So in this model, the measured flux is more intense when observed closer to the jet axis (with an angle $`\theta _v0`$ from the jet axis). A large spread in the inferred luminosities thus originates from the differences in viewing angles (see Figure 1a). This interpretation may be appealing because a large spread of both Lorentz factors and luminosities could be attributable to a standard type of event viewed from different orientations (Rees 1999). On the other hand, the large variations in the fluences of GRB could also originate from differences in the opening angle $`\theta _j`$ of an isotropic outflow (see Figure 1b; Frail et al. 2001; Panaitescu & Kumar 2001; Piran et al. 2001). The dynamics of the blast-wave in the off-axis anisotropic jet are almost indistinguishable from the uniform jet model, provided that the energy per unit solid angle varies as $`\theta _{v,j}^2`$ (Kobayashi et al. 2002; Rossi et al. 2002; Salmonson & Galama 2002; Zhang & Mészáros 2002).

The correlations between variability, the luminosity per unit solid angle and the characteristic photon energy in the cosmological rest frame (Lloyd-Ronning & Ramirez-Ruiz 2002) impose severe restrictions on any emission model. The question at hand is whether this burst phenomenology, combined with the intrinsic variations in energy per unit solid angle discussed above, can be understood in terms of internal shocks and if so, under what conditions. To this end, the results presented in §3 are independent of whether the variations in energy per unit solid angle are caused either by the viewing angle of a standard event (provided that $`d\mathrm{\Gamma }/d\theta _v`$ is small for $`d\theta _v\mathrm{\Gamma }^1`$) or by collimated jets, in which the typical Lorentz factor (or energy per unit solid angle) is higher when the jet is highly collimated.

### 2.1 Highly collimated outflows

As discussed above, the gamma-rays we receive come only from material whose motion is directed within an angle $`\mathrm{\Gamma }^1`$ of our line of sight. At observer times of more than a week, the shocked outflow starts to slow down. As $`\mathrm{\Gamma }`$ drops after the deceleration shock, the causal angle includes an increasing amount of the solid angle toward the jet as well as towards the equator, so we observe emission from a wider range of angles. The afterglow is thus a direct probe of the geometry of the ejecta - if the outflow is beamed, we expect a downturn in the light curve as the edge of the jet becomes visible (Rhoads 1997). Collimation factors of $`\mathrm{\Omega }_j/4\pi <0.01`$ have been derived from the observation of such steepenings111The above argument assumes that the breaks observed in many GRB afterglow light-curves are due to a geometrical beam effect and not to either a transition to non-relativistic expansion (Huang, Dai & Lu 2000) or an environmental effect such as a sharp density gradient (Chevalier & Li 2000; Ramirez-Ruiz et al. 2001). An important feature produce by this jet collimation is the achromaticity of the afterglow break, which will clearly distinguishes it from the steepening that may be produced by the passage of a spectral break through the observing band. (Kulkarni et al. 1999; Castro-Tirado et al. 1999; Harrison et al. 1999; Frail et al. 2001; Panaitescu & Kumar 2001). Although several methods of producing and collimating jets have been proposed, these highly collimated outflows of plasma with velocities close to the speed of light may be a unique feature of rapidly rotating, gravitationally confined plasmas that are threaded by a strong magnetic field (e.g. Uchida & Shibata 1985). An important result from present numerical simulations of Poynting flux-dominated outflows, which was not predicted by semi-analytical steady state studies, is that the jets arising from these outflows are likely to be highly collimated (see Meier, Koide & Uchida 2001 for a recent review). Because of the strong pinch that develops, a narrow jet that delivers its thrust in a narrow solid angle may be a common ingredient of strong rotating magnetic fields, not only for accretion disks but also perhaps for many jet-producing objects.

The range of parameter space where high jet speeds are expected (i.e. jets propelled at speeds much larger than the local escape velocity) is determined by the condition $`L_{\mathrm{MHD}}L_\mathrm{C}`$. The MHD luminosity

$$L_{\mathrm{MHD}}=\frac{B_p^2R_0^4\mathrm{\Omega }_0^2}{32c}$$
(1)
depends only on the poloidal magnetic field $`B_p`$ protruding from the jet-production region at radius $`R_0`$ and on the angular velocity of the field $`\mathrm{\Omega }_0`$ (Blandford & Znajek 1977), while the critical luminosity $`L_\mathrm{C}`$ is defined to be that given by the ratio between the energy needed for the plasma to reach the escape velocity and the free-fall time,

$$L_\mathrm{C}=4\pi R_0^2\rho _0\left(\frac{GM}{R_0}\right)^{3/2},$$
(2)
where $`\rho _0`$ is the density of the plasma. The condition $`L_{\mathrm{MHD}}L_\mathrm{C}`$ is a necessary requirement to explain the $`\mathrm{\Gamma }>100`$ outflow seen in GRBs, even from the vicinity of black holes. This criterion is similar to that of LeBlanc and Wilson 1970, in which the magnetic field grows to a dynamically important value in less than a dynamical time, which seems to be a condition of forming MHD jets in collapsing supernova cores (Meier et al. 1976; Wheeler et al. 1999). Although carrying more power, these highly collimated outflows will be much less efficient in imparting energy and momentum to the outer stellar layers and may not explode (Khokhlov et al. 1999). These failed explosions may be similar to the failed supernova model for GRBs (MacFadyen & Woosley 1999), continuing to accrete much of the stellar mantle and thus collapsing into a black hole.
While the very narrow $`2^o`$ and fast jets $`\mathrm{\Gamma }>100`$ that are speculated to exist in $`\gamma `$-ray bursts (Panaitescu & Kumar 2001) may be a general feature of strong rotating magnetic fields, some observations may be difficult to reconcile with this interpretation. For example, several afterglows have shown evidence of large amount of X-ray line emitting material, possibly arising from ionized iron (Piro et al. 1999; Yoshida et al. 1999; Piro et al. 2000; Antonelli et al. 2000). The large Fe K$`\alpha `$ equivalent widths inferred from the X-ray observations favor models in which the line is produced when the primary X-ray emission from the source strikes Thomson-thick material and Compton scatters into our line of sight (Lazzati, Campana & Ghisellini 1999; Böttcher & Fryer 2000; Vietri et al. 2000; Rees & Mészáros 2000; Ballantyne & Ramirez-Ruiz 2001, Lazzati et al. 2001). These conditions impose strong constraints on the location and geometry of the optically thick reflecting material, as well as on the structure of the outflow producing the primary X-rays - illumination of the stellar remnant may be complicated if the GRB emission is highly beamed. This problem may be overcome in particular source geometries for which lower luminosities and softer spectra are expected at the edges of the relativistic outflow (Ballantyne & Ramirez-Ruiz 2001). More importantly, highly collimated outflows bring the GRB rate to within a factor of 100 or so of the supernova rate (Frail et al. 2001). Those stars that produce GRBs are likely to have core masses $`>6M_{}`$ at the time they collapse (Woosley 1993). Stars satisfying this criterion are $`30-50`$ times rarer than those producing SN Ib/c (Izzard, Ramirez-Ruiz & Tout 2002, in preparation). Because burst formation is also likely to be favored by rapid rotation and low metallicity (MacFadyen & Woosley 1999), highly collimated outflows could face some difficulties in explaining large event rates (i.e. $`1/100`$ of the supernova rate). On the other hand, one appealing aspect of massive star progenitors is that the large variety of stellar parameters could explain the diversity of jet opening angles. If however the GRB 980425/1998bw association is real, then we may have a new class of GRB with lower energy $`E_\gamma 10^{48}(\mathrm{\Omega }_j/4\pi )`$ erg, which is only rarely observable even though its comoving density could be substantial.

### 2.2 Anisotropic relativistic outflows

Broader jets can also be produced by a magneto rotational mechanism, however, some concern has been raised that these outflows may not be easy accelerated to $`\mathrm{\Gamma }>10`$ (Meier et al. 2001 and references therein). Detailed predictions are nonetheless still in need of three dimensional, high-resolution core-collapse calculations. Other mechanisms of jet production such as neutrino radiation in the context of the collapsar model (MacFadyen & Woosley 1999) or intense radiation of a newly born pulsar (Blackman & Yi 1998) could be less constraining in their beaming requirements.

If the large variety of behaviors exhibited from burst to burst is produced by viewing an anisotropic universal beam configuration from different directions, the opening angle of the jet $`\theta _j>\theta _v`$ should still affect the light curve at the time when the edge of the jet becomes visible (this break would be additional to the one predicted from the viewing angle dependence; Rossi et al. 2002, Zhang & Mészáros 2002). If the critical Lorentz factor is less than 3 or so (the opening angle exceeds 20o) such a transition might be masked by the transition from ultra relativistic to mildly relativistic flow, so quite generally it would be difficult to limit the opening angle of this anisotropic outflow if it exceeds 20o. Although, the jet sideways expansion in this model is likely to be very important and more detailed calculations are required in order to understand the exact shape of the afterglow light curve, it is nonetheless important to mention that some afterglows are unbroken power laws for over 100 days (e.g. GRB 970228), implying that the opening angle of the late-time afterglow at long wavelengths is probably $`>20^o`$.

Also, under this interpretation, the $`\gamma `$-rays are likely to be more narrowly beamed than the optical afterglow and therefore there should be many “orphan” afterglows. The transient sky at faint magnitudes is poorly understood, but there are two major searches to find supernova down to R=23 (Garnavich et al. 1998; Perlmutter et al. 1998). These searches are sensitive to afterglows of the brightness levels observed to date and have covered a few tens of square degree years of exposure. It then follows that the afterglow rate is not more than a few times 0.1/sq deg/yr. Considering that the magnitude limit of these searches allows optical counterparts brighter than 1 ph cm-2 s-1 to be detected, the ratio of orphan afterglows to GRB is unlikely to exceed $`10`$ (Mészáros, Rees & Wijers 1999). This condition is not in disagreement with the few number of orphan afterglows that are predicted from an anisotropic universal beam configuration (Rossi et al. 2002). Studies of the transient sky at magnitudes down to R $``$ 25 should be able to disprove or consolidate this anisotropic model. The location of GRBs within their parental galaxy could help settle this question, since recently it was suggested that bursts located closer to the center of their parent galaxies have smaller isotropic equivalent energies (Ramirez-Ruiz, Lazzati & Blain 2002). If confirmed in further host observations, this correlation will strongly complicate this universal beam interpretation. A direct test of this model remains, however, the time dependent measurements of polarization.

## 3 Internal Shocks: A Phenomenological Study

### 3.1 Model Outline

We simulate GRB light curves by adding pulses radiated in a series of internal shocks that occur in a transient, unsteady relativistic wind (Rees & Mészáros 1994). Several authors have modeled this process by randomly selecting the initial conditions at the central site (Kobayashi et al. 1997; Daigne & Mochkovitch 1998; Pilla & Loeb 1998; Panaitescu, Spada, Mészáros 1999, Spada, Panaitescu & Mészáros 2000; Ramirez-Ruiz, Merloni & Rees 2001). The model used here is similar to that described by Ramirez-Ruiz et al. (2001) but differs in the following aspects: (i) Our treatment of the radiation emission takes into account synchrotron emission and Compton scattering; (ii) The effect of photon diffusion through the optically thick wind is included, which may be very important in obscuring pulses that occur at small radii; (iii) We use the shock jump equations to determine the physical conditions in the shocked fluid, and (iv) The duration of the pulse width $`\mathrm{\Delta }t_0`$ is calculated by adding in quadrature the angular spreading time $`\mathrm{\Delta }T_\mathrm{a}R/(2\mathrm{\Gamma }^2c)`$, the radiative cooling time $`\mathrm{\Delta }T_\mathrm{r}t_\gamma /\mathrm{\Gamma }`$ (where $`t_\gamma `$ is the lab-frame radiative timescale), the shell crossing time $`\mathrm{\Delta }T_\mathrm{\Delta }t_\mathrm{\Delta }/\mathrm{\Gamma }^2`$ (where $`t_\mathrm{\Delta }\mathrm{\Delta }/|v_{\mathrm{sh}}-v_0|`$ is the lab-frame shock’s crossing time, $`v_{\mathrm{sh}}`$ is the shell preshock flow velocity and $`\mathrm{\Delta }`$ is the shell thickness) and the diffusion time $`\mathrm{\Delta }T_\mathrm{d}`$ through the optically thick wind as done by Panaitescu et al. 1999 and Spada et al. 2000.

The wind is discretized into a sequence of $`N=T_{\mathrm{dur}}/\overline{\delta t_i}`$ shells with a range of initial thickness $`\mathrm{\Delta }_i`$, where $`T_{\mathrm{dur}}`$ is the duration time of the wind ejection from the central source and $`\overline{\delta t_i}T_{\mathrm{dur}}`$ is the average interval between consecutive ejections. The time interval $`\delta t_i`$ between two consecutive ejections $`i`$ and $`i+1`$ is assumed to be proportional to the $`i+1`$th shell energy. This implies that longer periods of quiescence, during which the engine accumulates fuel, are followed by more energetic shells (this is motivated by the observations of Ramirez-Ruiz & Merloni 2001 and Nakar & Piran 2002). The different values of $`\delta t_i`$ are selected from the log normal distribution of pulse intervals found by Norris et al. (1996).

We calculate the radii where shells collide and determine the emission features for each pulse. If some inner shell moves faster than an outer one ($`\mathrm{\Gamma }_i>\mathrm{\Gamma }_j`$), it will overtake the slow one at a radius $`R_i(t_{ij})=R_j(t_{ij})=R_c\mathrm{\Gamma }_i^2`$. For each collision between two shells there is a reverse and a forward shock. The shock jump equations determine the physical parameters of the shocked fluids: the velocities of the shock fronts, the internal energy in the shocked fluid frame $`u^{}`$, and both the thickness $`\mathrm{\Delta }_{ij}`$ and Lorentz factor $`\mathrm{\Gamma }_{ij}`$ of the merged shell at the end of the collision. We assume that in between two consecutive collisions the thickness of the shell increases proportionally to the fractional increase of its radius $`\frac{d\mathrm{\Delta }_i}{\mathrm{\Delta }_i}\frac{dR}{R}`$ (e.g. Spada et al. 2000). The ejection parameters determine the dynamics of the wind and the dynamical efficiency, $`ϵ_{ij}`$. This efficiency reflects the differences between the Lorentz factors of a pair of colliding shells. The efficiency for an individual collision can be calculated from the initial and final bulk energies,

$$ϵ_{ij}=1-\frac{M_{ij}\mathrm{\Gamma }_{ij}}{M_i\mathrm{\Gamma }_i+M_j\mathrm{\Gamma }_j}$$
(3)
where

$$\mathrm{\Gamma }_{ij}^2=\mathrm{\Gamma }_i\mathrm{\Gamma }_j\frac{M_i\mathrm{\Gamma }_i+M_j\mathrm{\Gamma }_j}{M_i\mathrm{\Gamma }_j+M_j\mathrm{\Gamma }_i},$$
(4)
and the resulting mass is $`M_{ij}=M_i+M_j`$. The first collisions remove the initial random differences between the Lorentz factors of successive shells. If the mean Lorentz factor $`\overline{\mathrm{\Gamma }}`$ remains steady for the entire burst duration, then the efficiency steadily decreases during the wind expansion. If $`\overline{\mathrm{\Gamma }}`$ is modulated on a timescale much smaller than the overall duration of the wind, dynamically efficient collisions at large radii are still possible.

A fraction $`\zeta 1`$ of the shock-accelerated electrons is assumed to have a power law distribution $`p`$ in the electron Lorentz factor $`\gamma _e`$, starting from a low random Lorentz factor $`\gamma _m\frac{1}{3}\frac{m_p}{m_e}\frac{ϵ_e}{\zeta }(\mathrm{\Gamma }_{ij}^{}-1)`$, where the energy stored in electrons is a fraction $`ϵ_e`$, $`\mathrm{\Gamma }_{ij}^{}`$ is the Lorentz factor for internal motions in the shocked frame and we assume $`p=2.5`$. The magnetic field is assumed to be turbulent and parameterized through the fraction $`ϵ_B`$ of the internal energy it contains: $`B^28\pi ϵ_Bu^{}`$ (primed quantities are measured in the comoving frame). The synchrotron spectrum of each pulse is approximated as three power-law segments, with the high energy slope(s) depending on $`p`$ and the relative values of the peak and cooling frequencies. The typical energy of synchrotron photons as well as the synchrotron cooling time depend on the Lorentz factor $`\gamma _e`$ of the accelerated electrons and on the strength of the magnetic field. The two energy release parameters $`\zeta `$ and $`ϵ_B`$ alter the peak energies of the synchrotron and inverse Compton spectra (they are higher for lower $`\zeta `$) and the Comptonization parameter (which is larger for smaller $`ϵ_B`$).
The shock-accelerated electrons radiate, and the emitted photons can be upscattered on the hot electrons ($`\gamma _e>>1`$) or downscattered by the cold ones ($`\gamma _e1`$). The optical depth to upscattering, far from the Klein-Nishina regime, is $`\tau _{\mathrm{ic}}\sigma _{\mathrm{Th}}\zeta n_e^{}\mathrm{min}(ct_\gamma ^{},\mathrm{\Delta }^{})`$ (Spada et al. 2000), where $`n_e^{}`$ is the electron density and $`t_\gamma ^{}=t_{\mathrm{sy}}^{}/(1+y)`$ is the radiative timescale with $`t_{\mathrm{sy}}^{}`$ the synchrotron cooling time and $`y`$ the Comptonization parameter ($`y=\gamma _m^2\tau _{\mathrm{ic}}`$ for $`\tau _{\mathrm{ic}}<1`$). A fraction min(1,$`\tau _{\mathrm{ic}}`$) of the synchrotron photons is inverse Compton scattered $`n_{\mathrm{ic}}`$=max(1, $`\tau _{\mathrm{ic}}^2`$) times unless the Klein-Nishina regime is reached (Panaitescu et al. 1999). The energy of the upscaterred photons is $`h\nu _{\mathrm{ic}}=\mathrm{min}[\gamma _mm_ec^2,h\nu _{\mathrm{syn}}(4\gamma _m^2/3)^{n_{\mathrm{ic}}}]`$. The optical thickness $`\tau _c`$ for the cold electrons within the emitting shell is calculated by taking into account the cold electrons within the hot fluid - those that were accelerated but have cooled radiatively while the shocked crossed the shell- and those within the yet unshocked part of the shell. When $`\tau _c>1`$, photons are downscattered by the cold electrons before they escape the emitting shell, leading to a decrease in photon energy and increase in photon duration. We approximate the increase in pulse duration due to the diffusion through optically thick shells by the time $`\mathrm{\Delta }T_\mathrm{d}`$ it takes for a photon to diffuse through them (Spada et al. 2000). The addition of all pulses gives both the burst spectra and the gamma-ray light curve. The latter is binned on a timescale of 64 ms and is used to computed the burst variability, that is, the average mean-square of the count variations relative to a smoothed time profile (calculated using a boxcar function with a timescale equal to 30% of the duration that contains 90% of the total counts; see Fenimore & Ramirez-Ruiz 2002a). A simulated profile depends on the specific values of the relevant model parameters describing the wind ejection ($`N50`$, $`\delta t_i`$, $`\mathrm{\Gamma }_{\mathrm{max}}`$, and $`\mathrm{\Gamma }_{\mathrm{min}}`$) and the radiative efficiency of the pulses ($`ϵ_e`$, $`\zeta `$ and $`ϵ_B`$).

### 3.2 On the internal structure of the relativistic wind

We study the effects of varying the the energy per unit solid angle of the expelled outflow by an increase either in the bulk Lorentz factor or in the baryon loading. One simple possibility is that the central engine ejects consecutive shells in which the Lorentz factors $`\mathrm{\Gamma }_i`$ are randomly selected between $`\mathrm{\Gamma }_{\mathrm{min}}`$ and $`\mathrm{\Gamma }_{\mathrm{max}}`$. The shell mass $`M`$ is then assumed constant and determined by the wind isotropic equivalent energy by requiring that $`Mc^2=E_{4\pi }/_{i=1}^N\mathrm{\Gamma }_i`$. In this scenario, narrow jets (or jets viewed close to the jet axis) should have higher Lorentz factors. Assuming $`\mathrm{\Gamma }_{\mathrm{max}}\theta _{j,v}^2`$, then $`E_{4\pi }\mathrm{\Gamma }_{\mathrm{max}}`$ (Kobayashi et al. 2002). Nonetheless, it could be also possible that changes in the energy per unit solid angle are mainly caused by variations in the baryon loading of the wind. In this case, the average bulk Lorentz factor $`\overline{\mathrm{\Gamma }}=\mathrm{\Gamma }_{\mathrm{min}}/2+\mathrm{\Gamma }_{\mathrm{max}}/2`$ remains unchanged throughout the entire wind and equal to $`\overline{\mathrm{\Gamma }}=E_{4\pi }/_{i=1}^NM_ic^2`$ with $`M_i`$ randomly drawn from the interval $`[M_{\mathrm{min}},M_{\mathrm{max}}]`$ (the shell Lorentz factors $`\mathrm{\Gamma }_i`$ are random between $`\overline{\mathrm{\Gamma }}`$-$`\mathrm{\Delta }\overline{\mathrm{\Gamma }}`$ and $`\overline{\mathrm{\Gamma }}`$+$`\mathrm{\Delta }\overline{\mathrm{\Gamma }}`$). Assuming a correlation between $`M_{\mathrm{max}}`$ and $`\theta _{j,v}`$, then $`E_{4\pi }M_{\mathrm{max}}`$. These two extreme scenarios will be referred to in the following as impulsive and lethargic  outflows, respectively.

### 3.3 Impulsive outflows

To simulate the effects of an impulsive outflow, we calculate the temporal and spectral profiles of 102 realizations for a given $`\mathrm{\Gamma }_{\mathrm{max}}`$, which we used to evaluate the peak of the spectrum in $`\nu F_\nu `$ and the mean variability. The wind (isotropic equivalent) energy, in all cases, is drawn from a log normal distribution with an average value $`\overline{E_{4\pi }}=10^{54}(\mathrm{\Gamma }_{\mathrm{max}}/10^3)`$ erg and a dispersion $`\sigma _E=10^{1/3}`$ erg. Numerous collisions happen during the evolution of the wind. Each collision produces a pulse, whose strength strongly depends on $`\mathrm{\Gamma }_{\mathrm{max}}`$, since the wind luminosity increases with $`\mathrm{\Gamma }_{\mathrm{max}}`$.

The optical thickness $`\tau _c`$ is determined mainly by the wind luminosity, $`L_w=E_{4\pi }/T_{\mathrm{dur}}`$, and the range of Lorentz factors in the wind. For small collision radii $`R\mathrm{\Gamma }^2`$, the pulse duration is dominated by the time it takes the photon to propagate through all the shells of optical thickness above unity. Thus, photons are downscaterred by the cold electrons before they escape the wind, leading to an increase of the pulse width and a decrease in the photon energy - photon diffusion through the wind widens pulses that would otherwise would appear narrow since they occur at small radii (Figure 2). This would be the case when the typical Lorentz factor $`\mathrm{\Gamma }_i`$ is relatively small. The smoothing effect is expected to be stronger when the the density of scattering is amplified by new $`e^\pm `$ pairs are formed in the originally thick scattering medium (Thompson & Madau 2000; Guetta, Spada & Waxman 2001; Mészáros, Ramirez-Ruiz & Rees 2001; Beloborodov 2002). For larger radii, the pulse duration is dominated by both shell crossing and radiative timescales (Figure 2), which decrease as $`\mathrm{\Gamma }_i`$ increases ( $`\mathrm{\Delta }T_\mathrm{\Delta }`$ increases as a result of the shell widening, while $`\mathrm{\Delta }T_\mathrm{r}`$ is larger for later collisions because $`B`$ and $`\gamma _m`$ are lower). The latter timescale is found to be comparable to the angular spreading one during the whole wind expansion, for the assumed linear shell-broadening between collisions. This implies that for large Lorentz factors, internal collisions will produce pulses with small $`\mathrm{\Delta }T_0`$ and high $`E_{4\pi }`$ – giving rise to a very spiky and luminous light curve.

By running models with $`10^2`$ different values of $`\mathrm{\Gamma }_{\mathrm{max}}`$ ranging from 100 to $`3\times 10^3`$ (assuming $`\mathrm{\Gamma }_{\mathrm{max}}/\mathrm{\Gamma }_{\mathrm{min}}`$=20), we have found that more variable profiles ($`>\mathrm{\Gamma }_{\mathrm{max}}`$) have higher luminosities and can reasonably fit the observational data (see Figure 3, where the shaded region depict the 1$`\sigma `$ region). The variability reaches an asymptotic value around $`V0.1`$ because the light curve is binned on a timescale of 64 ms to compare with BATSE temporal resolution. On the other hand, for very low radii (small $`\mathrm{\Gamma }_{\mathrm{max}}`$), the first collisions remove the initial random differences between the Lorentz factors of successive shells before reaching the radius of transparency. As a result, the efficiency steadily decreases during the wind expansion giving rise to very under luminous burst in which the variability measurement saturates due to count statistics. In this particular example, the values of the ejecta Lorentz factors are such (i.e. $`\mathrm{\Gamma }_{\mathrm{max}}/\mathrm{\Gamma }_{\mathrm{min}}`$ is set to be constant) that the average efficiency of individual collisions $`\overline{ϵ_{ij}}`$ is independent of the intrinsic luminosity of the outflow (or $`\mathrm{\Gamma }_{\mathrm{max}}`$). Within our model more efficient collisions can alternatively be produced by increasing $`\mathrm{\Gamma }_{\mathrm{max}}`$, while leaving $`\mathrm{\Gamma }_{\mathrm{min}}`$ unchanged. The difference in the shell Lorentz factors (i.e. $`\mathrm{\Gamma }_{ij}`$) in this case will be much larger, giving rise to more luminous pulses than in the previous case. A smaller range of $`\mathrm{\Gamma }_{\mathrm{max}}`$ values is thus required to produce similar changes in the burst isotropic luminosity (see the banded region in Figure 3). The numerical results shown in Figure 3 assume that the fraction $`ϵ_B`$ is distributed uniformly in logarithmic space between $`10^2`$ and $`10^1`$, $`n_{\mathrm{ISM}}=1\mathrm{c}\mathrm{m}^3`$ (in this case most collisions take place before the deceleration radius of the wind; see Fenimore & Ramirez-Ruiz 2002b), $`ϵ_e=0.25`$ and that $`\zeta 0.1`$.

The process by which the dissipated energy is finally radiated depends on the energy distribution of protons and electrons in the shocked material and on the values of the comoving density and magnetic field. The internal shocks heat the expanding ejecta, amplifying the preexisting magnetic field or generating a turbulent one, and accelerate electrons, leading to synchrotron emission and inverse Compton scattering. By varying the injection fraction $`\zeta `$ one can study the relative intensity of the synchrotron and inverse Compton components and determine those values that maximize the received flux in the BATSE range. For $`\zeta >10^2`$ the inverse Compton emission occurs in the Thomson regime and carries most of the burst energy in the BATSE window, while for $`\zeta <10^3`$ synchrotron dominates over inverse Compton scattering as the latter takes place in the Klein-Nishina regime (Panaitescu et al. 1999). Figure 4 shows the evolution of the synchrotron (which for $`\zeta `$=1 lies mainly below the BATSE window) and inverse Compton peaks before and after downscattering for both thick and a thin winds. For smaller collision radii (i.e. small $`\mathrm{\Gamma }_{\mathrm{max}}`$), the emission takes place where the wind is optically thick, leading to a decrease in the photon energy - diffusion through the wind decreases the photon energy that would otherwise would appear large since they are produced at small radii where $`B`$ is larger. For the set of parameters considered here, the Thompson limit is usually a good approximation to treat the downscattering of the synchrotron photons. The inverse Compton emission, however, peaks at large comoving frame energies, and the general cross section has to be considered. For these energetic photons, the cross section strongly depends on the photon energy and varies after each photon interaction. Figure 5 shows the shifting of the burst emission towards higher energies with increasing variability, due to the increase of $`\mathrm{\Gamma }_{\mathrm{max}}`$ from 100 to $`3\times 10^3`$ and assuming $`\mathrm{\Gamma }_{\mathrm{max}}/\mathrm{\Gamma }_{\mathrm{min}}`$=20 (corresponding to the internal shock parameters of the simulations shown in Figure 3). This increase in energy stops once the average bulk Lorentz factor is large enough that collisions take place where the wind is optically thin (see the shaded region in Figure 5). The typical radiation energy of each pulse also strongly depends on the resulting Lorentz factor $`\mathrm{\Gamma }_{ij}`$. Thus it is possible to produce a larger dynamic range of spectral peak energies by increasing $`\mathrm{\Gamma }_{\mathrm{max}}`$, while leaving $`\mathrm{\Gamma }_{\mathrm{min}}`$ unchanged (see banded region in Figure 5). We found that more variable profiles have higher spectral peaks and can reasonably fit the observational data (Lloyd-Ronning & Ramirez-Ruiz 2002) provided that $`\zeta >10^2`$. In this case, the synchrotron emission lies mainly below the BATSE window and the inverse Compton component carries most of the energy of burst. When synchrotron emission dominates over inverse Compton scattering, as the latter takes place in the Klein Nishina regime (i.e. $`\zeta <\mathrm{\hspace{0.17em}10}^2`$), we find that simply varying $`\mathrm{\Gamma }_{\mathrm{max}}`$ of the shells does not produced the dynamical range of peak energies that is observed. This is mainly because the typical observed synchrotron frequency, which at $`R_c`$ gives

$$h\nu _{\mathrm{syn}}\frac{ϵ_e^2ϵ_B^{1/2}(\mathrm{\Gamma }_{ij}^{}-1)^2L_w^{1/2}}{R_c(1+z)},$$
(5)
is unlikely to increase by simply incrementing $`\mathrm{\Gamma }`$ as $`R_c\mathrm{\Gamma }^2`$ gives $`h\nu _{\mathrm{syn}}\mathrm{\Gamma }^2`$ for a constant $`L_w`$, $`ϵ_B`$, $`ϵ_e`$ and $`\mathrm{\Gamma }_{ij}^{}`$. Assuming the correlation between the Lorentz factor and $`\theta _{j,v}`$, the isotropic energy itself is related to the source expansion velocity, $`L_wE_{4\pi }\mathrm{\Gamma }_{\mathrm{max}}`$, and the typical observed synchrotron frequency is then222This scaling is valid if $`ϵ_B`$, $`ϵ_e`$ and $`\mathrm{\Gamma }_{ij}^{}`$ are constant. However, it is possible to produce a larger dynamical range of spectral energies by increasing $`\mathrm{\Gamma }_{\mathrm{max}}`$ while leaving $`\mathrm{\Gamma }_{\mathrm{min}}`$ unchanged. In this case, $`\mathrm{\Gamma }_{ij}^{}\mathrm{\Gamma }`$ and thus $`h\nu _{\mathrm{syn}}\mathrm{\Gamma }^{1/2}`$ $`\mathrm{\Gamma }^{3/2}`$. This may be overcome by implicitly assuming that the electrons are continually accelerated in the shocked region (see Lloyd & Petrosian 2000) or by postulating that the equipartition energy acquired by the electrons and the magnetic field ($`ϵ_e`$ and $`ϵ_B`$) are dependent on $`\mathrm{\Gamma }`$. While these are strong constraints in the emission properties, it is worthy of notice that the correlation between variability and spectral energy could be significantly affected by selection effects and more GRBs with independent spectroscopic redshifts are need it in order to test its validity (Lloyd-Ronning & Ramirez-Ruiz 2002). This relationship -if true- can help shed light on the relevant emission mechanisms in unsteady outflows. We would also like to emphasize that the wind optical thickness is one of the factors to which the variability and spectral peak energy is most sensitive. Indeed, it is reassuring that the inclusion of this effect is of the utmost importance to reproduce the power density spectrum of gamma-ray bursts arising from internal shocks in unsteady outflows (Panaitescu et al. 1999, Spada et al. 2000).

### 3.4 Lethargic outflows

Until this point we have assumed that the changes in the energy per unit solid angle are caused solely by an increase in the bulk Lorentz factor ($`M`$=constant). Nonetheless, it could be possible that changes in $`E_{4\pi }`$ may be mainly caused by variations in the baryon loading of the wind. We have run such lethargic models and find that simply varying the baryon loading of the wind of the ejected shells does not reproduce the dynamic range of variability and spectral energies that is observed. We also study intermediate behaviors by assuming that the shell baryon loading and its Lorentz factor are related: $`\mathrm{\Gamma }M^\alpha `$. We conclude that, within such models, the burst phenomenology can be accommodated provided that the changes in the energy per unit solid angle within different beams is dominated by $`\mathrm{\Gamma }`$ with $`\alpha <0.2`$. This is not surprising since the timescales that determine the pulse structure and thus the variability are highly dependent on $`R_c\mathrm{\Gamma }_i^2`$. An interesting consequence is that large variations in the energy per unit solid angle within different beam models should be strongly dominated by changes in $`\mathrm{\Gamma }`$.

## 4 Discussion

The main conclusion we can draw from the simulations presented above is that, in the framework of the internal shocks model, there are realistic assumptions in which the variable activity of GRB is found to be correlated with absolute peak luminosity. We attribute this correlation to a variation of the energy per unit solid angle caused mainly by an increase in the emitting regions’ Lorentz factors in the context of impulsive relativistic outflows. This change of the Doppler factor could be attributable either to a standard type of event viewed from different orientations or by changes in the collimation of GRBs jets. Within such models, collisions between shells with higher Lorentz factors, which are intrinsically more luminous, produce more variable profiles. The existence of a correlation between the characteristic photon energy in the cosmological rest frame and the GRB variability, on the other hand, imposes more constraining requirements. The main characteristics of the modeled bursts with the above mentioned features are: a high electron injection fraction – which implies that the inverse Compton scattering dominates over synchrotron emission – required to increase the spectral peak energy of the larger collision radii; and a wind optical thickness to scattering on cold electrons above unity is required to increase the duration of the pulses and decrease the photon energy of the collisions that take place at smaller radii.

We must also remain aware of other possibilities. For instance, we may be wrong in supposing that the main radiation mechanism assumed to be responsible for the burst event is synchrotron or its inverse Compton component. It could be that the GRB emission originates from a fireball moving out through an extended stellar envelope, along a funnel that is empty of matter but pervaded by thermal radiation from the funnel walls, leading to Compton drag (Lazzati et al. 1999; Ghisellini et al. 2001). The fireball itself remains optically thick until it expands beyond the stellar surface (see also Ramirez-Ruiz, MacFadyen & Lazzati 2002). A burst with complex time-structure could then be produced by a series of expanding shells. In this scenario, faster moving shells will produce greater peak energies, as the burst emission simply reflects the temperature of the funnel photons up-scattered by the square of the bulk Lorentz factor. Furthermore, these fast shells will produce both smaller and more luminous pulses. This is because the observed variability timescale, which is related to the typical size of the region containing the dense seed photon divided by the time compression factor, increases with $`\mathrm{\Gamma }`$ and the Compton drag is more efficient for faster shells. To this end, ejection of very rapid shells will lead to bursts produced by the ejection of very rapid shells are both more variable and have larger peak spectral energies.

Much progress has been made in understanding how afterglows can arise from a forward shock or blast-wave moving into the external medium ahead of the ejecta, and in deriving the generic properties of the long wavelength afterglows that follow from this (Mészáros 2001). There still remain a number of mysteries, especially regarding the prompt emission, in particular, the formation of the ultra relativistic outflow, its structure and the radiation process. As we have shown, the correlations between the characteristic photon energy, the gamma-ray burst variability and its luminosity, as well as the jet collimation, can help shed light on all three of these issues.

We are indebted to M. J. Rees and D. Lazzati for very useful discussions. We thank D. Lazzati for highlighting the relevance of the Compton drag scenario in reproducing the burst phenomenology. We would also like to thank E. Fenimore, S. Kobayashi, P. Kumar, A. Panaitescu, E. Rossi and the referee P. Mészáros for helpful insight regarding internal-shock and non-thermal radiation calculations. ERR acknowledges support from CONACYT, SEP and the ORS foundation.