text:
# Self-organized critical network dynamics

## Abstract

We propose a simple model that aims at describing, in a stylized manner, how local breakdowns due unbalances or congestion propagate in real dynamical networks. The model converges to a self-organized critical stationary state in which the network shapes itself as a consequence of avalanches of rewiring processes. Depending on the model’s specification, we obtain either single scale or scale-free networks. We characterize in detail the relation between the statistical properties of the network and the nature of the critical state, by computing the critical exponents. The model also displays a non-trivial, sudden, collapse to a complete network.

Complex networks underlying many social and technological systems is a subject of booming recent interest AB ; Doro\_rev . On one side the structure of such networks has non-trivial properties AB , which influences dramatically the nature of processes taking place on them (see e.g. Optimal ; Epidemics ). On the other, the network’s structure constrains in a peculiar way the growth ABGrowth  and evolution Ebel ; Social  of the network itself. This calls for an extension of statistical physics, which traditionally studies collections of dynamical variables interacting through a fixed network, to systems where the network of interactions itself becomes a dynamical variable.

Here we focus on dynamical networks where links do not represent physical bonds but rather relationships or communication channels. Ref. Optimal  suggests that a structured network of communications aimed at solving problems or carrying out specific functions is a crucial feature of firms and organizations in general. Routers table in the Internet is also an example of a dynamic communication network. The network of financial institutions, linked by mutual contracts and loans, provides a further example of a dynamic network. Beyond “static” design problems, such as e.g. minimizing congestion Optimal  or redistributing optimally the loads Moreno , systems of this type also pose “dynamical” problems such as how and to what extent does congestion or breakdown events propagate through the system.

Here we address these problems in a dynamic network subject to two competing forces: on one side, there is a drive toward increasing complexity, by e.g. adding new links, because the system performs more efficiently its functions as it becomes more densely interconnected. On the other, the resulting increase in complexity may bring about conflicting constraints, imbalances or congestion problems which may cause a local breakdown of the network. A local breakdown may engender a re-adaptation in its neighborhood which may inadvertently cause the breakdown to propagate further on the network.

For example, a change in some router’s table in the Internet, which is meant to improve efficiency, may inadvertently cause congestion at some node downstream. This may trigger several other changes in that local neighborhood, as routers try to avoid the congested nodes. But these changes may, in their turn, cause further congestion elsewhere, and the problem may expand even further, as an avalanche, to a wider region. The stipulation of a contract between two institutions, which is in principle beneficial to both, may also increases their operative constraints, making them less adaptable to a changing environment and hence more exposed to the risk of bankruptcy. The failure of one institution, likely induces a rearrangement of the institutions linked with it and perhaps engenders effects which propagate further across the network Holyst . Similar phenomena may take place in social or trade networks.

Rather than trying to model in a realistic manner one of the problems just discussed, we focus on a simple model of network dynamics which captures the two main ingredients discussed above: a slow dynamics where links are added to the network and a fast relaxation dynamics of avalanche events. The motivation for this choice, is that, the detailed understanding of the behavior of a simple model with these features, may be the basis or at least a guide for addressing more complex and realistic situations, such as those discussed above. Our main finding is that such systems can self-organize close to a critical point where each modification of the network’s architecture can have unforeseable consequences which possibly affect a wide region of the system. This may have some bearing on the intermittence of Internet traffic or on the nature of financial crises and recessions.

We start from a empty network of $`N`$ nodes in which every node $`i`$ is assigned a fitness $`f_i`$ drawn from a probability distribution $`\rho (f)`$. Let $`F_i`$ be the set of neighbors of $`i`$ and $`k_i=|F_i|`$ be the number of neighbors of $`i`$. At every time step a link is added between two previously unconnected random nodes $`i`$ and $`jF_i`$. With probability $`f_j`$ nothing happens whereas with probability $`1f_j`$ the node $`j`$ becomes unstable or congested and it “topples”. As a result all its links (including that with $`i`$) are rewired to randomly chosen nodes, i.e. for any $`hF_j`$, a node $`lF_h`$ is chosen at random and the link $`jh`$ is rewired to $`hl`$. In its turn, with probability $`1f_l`$, also node $`l`$ may become unstable and topple. Hence toppling of node $`j`$ may start an avalanches of toppling events which propagates through the network rearranging it. Unstable nodes, after they topple, remain unconnected from the network notacon  and are assigned a new fitness value drawn from the distribution $`\rho (f)`$. Hence toppling is equivalent to replacing the unstable node with a new one.

In order to stabilize the network and reach a stationary state, we introduce dissipation of links: at each toppling event, with probability $`\lambda `$ all the links of the unstable node are removed from the network. Note that without dissipation the number of links would increase in time until the complete graph is reached. The complete graph ($`jF_i`$ $`i,j,ij`$) is an absorbing state of the dynamics, because no link can be added to it.

The distribution $`\rho (f)`$ is the only parameter of the model. An alternative class of models can be defined by specifying the probability $`u_k`$ that a node with $`k`$ neighbors becomes unstable upon addition of a further link. A relation between the two models is possible, along the lines of Ref. Polya , in the limit $`\lambda 0`$ on no dissipation. Then the probability to find a node with $`k`$ neighbors and $`f_i[f,f+df)`$ is $`\rho (f|k)dff^k\rho (f)df`$ where $`f^k`$ is the probability that a node with $`f_i=f`$ has $`k`$ neighbors. Then a model with

$$u_k=_0^1𝑑f(1f)\rho (f|k)=\frac{_0^1𝑑f(1f)f^k\rho (f)}{_0^1𝑑ff^k\rho (f)}$$
(1)
is completely equivalent to one specified in terms of $`\rho (f)`$, in the limit $`\lambda 0`$. The dependence on $`k`$ of $`u_k`$ reflects the fact that $`k_i`$ and $`f_i`$ are positively correlated because nodes with higher fitness have smaller chance of becoming unstable. For convenience we shall refer mostly to models specified in terms of $`u_k`$ using Eq. (1) to translate the results in the original model.

Our results can be summarized as follows: i) when $`u_k`$ decays faster than $`1/k`$ there is a critical $`\lambda _c`$ such that the network evolves toward a complete graph for $`\lambda <\lambda _c`$. The same happens for $`u_kb/k`$ ($`k1`$) and $`b<3/2`$ whereas when $`b>3/2`$ or when $`u_k`$ decays slower than $`1/k`$, the collapse takes place only in finite networks. Indeed we find complete graphs only for $`\lambda <\lambda _cN^\gamma `$, where $`\gamma =\frac{2b3}{b1}`$ for $`3/2<b<2`$ and $`\gamma =1`$ otherwise. ii) the non-collapsed phase $`\lambda >\lambda _c`$ is characterized by an uncorrelated random network NewmanRand  with finite average degree and a degree distribution $`p_k`$ which depends on $`\rho (f)`$ (or $`u_k`$) (see Fig. 1). In particular, if $`u_k`$ decays slower than $`1/k`$ then $`p_k`$ decays faster than any power, whereas if $`u_kb/k`$ then $`p_kk^b`$. iii) The dynamics converges to a stationary sequence of avalanches of rewiring processes with a power law distribution $`P(s)s^\tau `$ of sizes (see Fig. 3). As in Ref. Polya , we shall define the size $`s`$ of an avalanche as the number of toppling events that it causes. The exponent $`\tau `$ takes the mean field valueMF  $`\tau =3/2`$ when $`u_k`$ decays slower than $`2/k`$, whereas $`\tau =1+1/b`$ if $`u_kb/k`$ with $`3/2<b<2`$.

The collapse to a complete graph, where congestion is minimal, is reasonable, given that the network is trying to adapt by avoiding congestion. The non-trivial issue on which we shall concentrate mostly, is related to the self-organized critical state. In view of the special role played by the case $`u_k1/k`$, we focus on the following simple forms of $`\rho (f)`$ and to the corresponding $`u_k`$:

$`\rho (f)=b(1f)^{b1},`$ $`u_k={\displaystyle \frac{b}{b+k+1}}`$ (2)
In order to shed light on the model’s behavior, let us derive an equation for the avalanche distribution. Define $`s_k`$ as the avalanche size originating from an unstable node with $`k`$ neighbors. This is a random variable which can be decomposed as follows

$$s_k=1+d\underset{j=1}{\overset{k+1}{}}v_{k_j}s_{k_j}.$$
(3)
into the contribution of the unstable node and those of the avalanches $`s_{k_j}`$ ensuing from its neighbors, with $`k_j`$ being the number of neighbors of the $`j^{\mathrm{th}}`$ neighbor. Note that the sum runs over $`k+1`$ links as it includes the link which caused the instability and the $`k`$ pre-existing neighbors. In Eq. (3) $`d=0,1`$ describes the effect of dissipation with $`P(d=0)=\lambda =1P(d=1)`$ whereas $`v_{k_j}`$ takes value $`v_{k_j}=0`$ if rewiring of the link to the $`j^{\mathrm{th}}`$ neighbor causes no further toppling, and $`v_{k_j}=1`$ otherwise. Hence $`P(v_l=1)=1P(v_l=0)=u_l`$.

Now we can write the generating function $`\varphi _k(z)=z^{s_k}`$ of the probability $`P(s|k)`$ to have an avalanche of size $`s`$ given that the initiator node has $`k`$ neighbors. From this, it is easy to find the generating function $`\chi (z)`$ of the distribution $`P(s)`$ of avalanche sizes

$$\chi (z)\underset{s}{}P(s)z^s=\frac{1}{\overline{u}}\underset{k=0}{\overset{\mathrm{}}{}}p_ku_k\varphi _k(z).$$
(4)
with $`\overline{u}=_kp_ku_k`$. After some algebra, using the fact that the rewired nodes are chosen randomly in the network, Eq. (3) leads us to

$$\chi (z)=z\left[\lambda +\frac{1\lambda }{\overline{u}}\underset{k=0}{\overset{\mathrm{}}{}}p_ku_k\left[1\overline{u}+\overline{u}\chi (z)\right]^{k+1}\right]$$
(5)
which is a non-linear self-consistent equation for $`\chi (z)`$. Note that $`\chi (1)=1`$ as it should and that

$$\chi ^{}(1)=s=\frac{1}{1(1\lambda )(k+1)u}$$
(6)
Apart from the parameters $`\lambda `$ and $`u_k`$, Eq. (5) depends also on the stationary state degree distribution $`p_k`$. Hence, before discussing Eq. (5) further, we need to elaborate on the nature of the stationary state. A necessary condition in order to be at the stationary state is that the total number of links $`K(t)`$ present in the network is constant on average. $`K(t)`$ increases by one for the random addition of a new link and is reduced by the total number of links $`\mathrm{\Lambda }`$ dissipated during the avalanche that follows. Here

$$\mathrm{\Lambda }=v\underset{i=1}{\overset{s}{}}d_i(k_i+1)$$
(7)
where $`v=1`$ if the chosen site is unstable and $`v=0`$ otherwise, and $`d_i=1`$ if dissipation occurs at the toppling site $`i`$, otherwise $`d_i=0`$. Consequently $`v`$ and $`d_i`$ have average values $`\overline{u}`$ and $`\lambda `$ and thus $`\mathrm{\Lambda }=\lambda (k+1)u_ks`$. Then, stationarity $`\mathrm{\Lambda }=1`$ and Eq. (6) imply that

$$(k+1)u_k=1.$$
(8)
As a byproduct, we find $`s=1/\lambda `$, in perfect agreement with numerical simulations.

The stationary degree distribution $`p_k`$ can be found by quantifying the Markov chain of possible transitions $`k_ik_i^{}`$ during the dynamics. In the limit $`\lambda 0`$ and $`N\mathrm{}`$ where we can neglect dissipation and finite size effects, there are only two processes which take place on each node: $`k_ik_i+1`$, with probability $`1u_{k_i}`$ and $`k_i0`$ with probability $`u_{k_i}`$. Then, in the stationary state, $`p_k`$ satisfies

$$p_{k+1}=(1u_k)p_k.$$
(9)
Taking the sum of Eq. (9) on $`k`$ we find $`p_0=u_k`$, which means that the fraction of sites with no neighbor is equal to the probability that a node becomes unstable. Furthermore multiplying Eq. (9) by $`k+1`$ and taking the sum over $`k`$, we recover the stationary condition (8). In the simplest case $`u_k=u`$ for all $`k`$, we find $`p_k=u(1u)^k`$ whereas with Eq. (2) we find

$$p_k=(b1)\frac{\mathrm{\Gamma }(b)\mathrm{\Gamma }(k+1)}{\mathrm{\Gamma }(b+k+1)}k^b$$
(10)
where the asymptotic power law behavior holds for $`k1`$. Notice that $`k=1/(b2)`$ diverges when $`b2^{}`$ and that there is a finite fraction $`\overline{u}=11/b`$ of unconnected nodes. Still the network has a giant connected component for $`b<7`$ NewmanRand . Eq. (9) yields a $`p_k`$ which decays faster than any power if $`u_k`$ decays less slowly than $`1/k`$, or if it increases. Conversely, if $`u_k`$ decays faster than $`1/k`$, we find that $`p_k`$ is not normalizable for $`N\mathrm{}`$. Numerical simulations (see Fig. 1) fully supports this picture, even though the effects of dissipation and finite size are clearly evident for $`k1`$.

The neglect of dissipation, when $`N`$ is finite, is a reasonable approximation if nodes with maximal degree $`k_i=N1`$ are not stable. A node connected to all neighbors, cannot receive further links and hence cannot become unstable. Its degree decreases only if dissipation occurs at a node connected to them. The rate of this process, for a node with degree $`k`$, is $`\lambda k/k`$. Hence if $`\lambda N/k1`$, nodes with $`kN`$ decay very fast and the only effect of dissipation is to introduce a cutoff $`k_c1/\lambda `$ in the distribution $`p_k`$. When $`\lambda N/k1`$ we expect a transition (close) to the complete graph $`k_i=N1`$ for all $`i`$. When $`k`$ is finite, i.e. for $`b>2`$ or for $`u_k`$ which decays slower than $`1/k`$, the collapse to a complete graph takes place for $`\lambda <\lambda _cN^1`$. When $`3/2<b2`$ the average degree $`kN^{\frac{2b}{b1}}`$ diverges with the system size. Then the decay rate of totally connected nodes is $`\lambda N^\gamma `$ with $`\gamma =\frac{2b3}{b1}`$ and the collapse to a complete graph takes place for $`\lambda <\lambda _cN^\gamma `$. In both cases (i.e for $`b>3/2`$) for any $`\lambda <1`$, it is always possible to take $`N`$ large enough to make the decay rate $`\lambda N^\gamma `$ large enough so that finite size effects can be neglected. But when $`b<3/2`$ this is no longer true because $`kN`$ and the decay rate of completely connected nodes remains finite ($`\gamma =0`$) even when $`N\mathrm{}`$. Beyond a finite dissipation rate $`\lambda _c`$, the network collapses to the complete graph.

Fig. 2 fully confirms the theoretical insight discussed above. When $`\lambda >\lambda _c`$, where $`\lambda _cN^\gamma `$ (see inset), the dynamics reaches a network with finite average degree, whereas for $`\lambda <\lambda _c`$ a collapse to the complete graph is observed, with a transition which is sudden and discontinuous. In the case in which the $`u_k=u`$ we also observe a transition from a finite average connectivity to a average connectivity of order $`N`$. The transition occurs for values of $`\lambda _cN^1`$ but the transition is rather smooth.

Having discussed the stationary state, let us go back to the Eq. (5) for the avalanche size distribution. We focus on the region $`1\lambda N^\gamma `$ where the network is not densely connected and dissipation effects are weak. Anticipating that the avalanche size distribution acquires a cutoff $`s_c\lambda ^\sigma `$, for some exponent $`\sigma >0`$, we postulate the scaling form $`P(s)s^\tau \mathrm{\Phi }(s\lambda ^\sigma )`$ with finite $`\mathrm{\Phi }(0)`$ and $`\mathrm{\Phi }(x)0`$ faster than any power as $`x\mathrm{}`$. Such a scaling hypothesis is fully corroborated by numerical results, as shown in Fig. 3. This correspondsBMF  to an analogous scaling form

$$\chi (z)1(1z)^{\tau 1}h\left(\frac{1z}{\lambda ^\sigma }\right)$$
(11)
for the generating function for $`\lambda 1`$ and $`1z\lambda ^\sigma `$. Setting for convenience $`1z=x\lambda ^\sigma `$, asymptotic analysis for $`\lambda 1`$ shows that the leading orders of Eq. (5) are

$$\lambda ^{1+\sigma (\tau 1)}x^{\tau 1}h(x)\lambda ^\sigma x+c\lambda ^{\sigma \beta (\tau 1)}x^{\beta (\tau 1)}h^\beta (x)=0$$
(12)
where

$`\beta =2,`$ $`c={\displaystyle \frac{b1}{b(b2)}}`$  for $`b>2`$ (13)
$`\beta =b,`$ $`c={\displaystyle \frac{\pi b^{b1}(1b)^b}{\mathrm{sin}[\pi (2b)]}}`$ $`\text{ for }b<2.`$ (14)
Note that $`c1/|b2|`$ diverges when $`b2`$ is approached from both sides. Dividing Eq. (12) by $`\lambda ^\sigma `$ and taking the scaling limit $`\lambda 0`$ with $`x`$ finite, we find a non-trivial result with all three terms finite if we choose

$`\tau =3/2,`$ $`\sigma =2`$ for $`b>2`$ (15)
$`\tau =1+1/b,`$ $`\sigma =b/(b1)`$ $`\text{for }b<2.`$ (16)
and the scaling function is the inverse of $`x(h)=h^\beta /[1ch^\beta ]^\beta `$. In particular $`h(x)c^{1/\beta }`$ for $`x\mathrm{}`$ and $`h(x)x^{1/\beta }`$ for $`x1`$. For $`b>2`$ the solution coincides with that of other mean field models BMF  $`h(x)=[\sqrt{c/x+4}\sqrt{c/x}]/(2c)`$ and perfectly matches numerical simulations for a range of values of $`\lambda `$ (see the inset of Fig. 3). It is easy to check that the model with $`u_k`$ decaying slower than $`1/k`$ falls in the $`\beta =2,\tau =3/2`$ $`\sigma =2`$ universality class. In conclusion we have shown how a slow growth dynamics and a fast relaxation through avalanche events can generate a dynamical network with given degree distribution. The stationary state is critical in the sense that the avalanches of all sizes occur, and is reached spontaneously without fine tuning of external parameters as long as the dissipation rate is larger than a given threshold $`\lambda _cN^\gamma `$. For smaller dissipation rates, the network collapses to the complete graph. While the detailed solution depends on the particular simplicity of the model chosen, the generic picture may apply to a wider class of systems and capture some features of the non-linear and intermittent behavior of real systems. We acknowledge F. Vega-Redondo for useful comments and discussions. This work is supported in part by the European Community’s Human Potential Programme under contract HPRN-CT-2002-00319, STIPCO.


format: text

text:
# Comment on “Vacuum Photon Splitting in Lorentz-Violating Quantum Electrodynamics”

(August 2, 2003) preprint: KA–TP–06–2003

In a recent Letter KoPi , Kostelecký and Pickering argue that there may be observable effects from photon triple splitting in an extended version of quantum electrodynamics with Lorentz violation in the fermion sector. The argument is based on a supposed analogy with a known physical process. In this Comment, we point out that the analogy is misleading and that, at the order considered, the probability of on-shell photon triple splitting is strictly zero.

Kostelecký and Pickering KoPi  study the process of photon splitting in certain Lorentz- and CPT-violating extensions of quantum electrodynamics. For conventional quantum electrodynamics, the photon splitting amplitude is zero to all orders of perturbation theory, but in a Lorentz-violating extension of the theory the situation may be different. Concretely, the authors study a model where the photon sector is conventional and the Lorentz violation is due to additional terms in the fermion sector. The photon kinematics thus remains unchanged and the decay of an initial photon into any number of final photons is only possible if all final photons are collinear with the initial photon.

The main part of Ref. KoPi  is devoted to the calculation of the relevant fermionic one-loop diagrams, to first order in the Lorentz-violating terms. The authors still find a vanishing amplitude for the splitting of an initial photon into two final ones. But they do find a non-zero amplitude for splitting into *three* final photons, thereby establishing a difference with the case of conventional quantum electrodynamics.

The authors of Ref. KoPi  close with some remarks on the possible physical significance of their result, i.e., whether or not the finite nonzero amplitude for photon triple splitting could lead to observable effects. The obvious answer would appear to be negative, as the phase space volume is zero for decay into three or more photons; cf. Refs. Ha ; Ad ; WiSm . But the authors argue in favor of a possible physical effect by advocating an analogy to photon splitting via collinear parametric down-conversion in optically active crystals. Specifically, they quote an experiment reported in Ref. GiMi .

The analogy is, however, misleading, for the following two reasons. First, the experiment described in Ref. GiMi  is based on phenomena induced by *quadratically* nonlinear optics, as follows from their expressions for the nonlinear polarization. Microscopically, this corresponds to the decay of one initial photon with energy $`\mathrm{}\omega _\mathrm{p}`$ (“pump”) into *two* final photons with energies $`\mathrm{}\omega _\mathrm{s}`$ (“signal”) and $`\mathrm{}\omega _\mathrm{i}`$ (“idler”), where $`\omega _\mathrm{p}=\omega _\mathrm{s}+\omega _\mathrm{i}`$. The decay into two photons has, of course, a non-zero phase space volume Ha ; Ad ; WiSm .

Second, nonlinear optical phenomena like parametric down-conversion, parametric fluorescence, etc., almost always involve either the *spatially inhomogeneous structure* or the *nonstandard dispersion law* of the optical device used; cf. Ref. Yar . It may, for instance, happen that the photon momentum conservation condition $`\stackrel{}{p}_\mathrm{p}=\stackrel{}{p}_\mathrm{s}+\stackrel{}{p}_\mathrm{i}`$ does not hold exactly (momentum being absorbed by the crystal) or that the momenta of the initial and final photons are not perfectly collinear (consistent with the nonstandard dispersion law). These effects are, of course, absent for photons which propagate in a homogeneous vacuum and which have the standard vacuum dispersion law $`p_0=\omega \left(\stackrel{}{p}\right)=|\stackrel{}{p}|`$, as is the case for the model of Ref. KoPi .

A better physical analogy may be photon splitting in a constant external magnetic field Ad . As long as this background magnetic field is constant (i.e., time independent and homogeneous) and the photons obey the standard vacuum dispersion law, the only possibility is the splitting of a single photon into two. The reason is that energy and momentum cannot be extracted from the constant background field and that the phase space volume for higher splittings is zero Ha ; Ad ; WiSm . Since the vacuum of the model studied in Ref. KoPi  is perfectly homogeneous (although no longer isotropic) and the photons obey the standard dispersion law, the same reasoning applies also in this case. \[Possible quantum gravity effects which induce inhomogeneities (e.g., a spacetime foam) are not considered.\]

We conclude that the comparison to parametric down-conversion as advocated by Ref. KoPi  is inappropriate and that the probability for photon triple splitting is strictly zero, at least to quadratic order in the coefficients of the Lorentz-violating terms considered. Indeed, there appears to be a kind of “conspiracy” in quantum electrodynamics models with only the fermion sector modified, which keeps the photon stable by setting either the amplitude or the phase space volume equal to zero. For on-shell photon triple splitting, it may very well be necessary to have a modified photon sector; cf. Ref. AdKl .
*Acknowledgements:* CA acknowledges support from the Austrian START award project FWF-Y-137-TEC of N.J. Mauser.


text:
# Short-period AM CVn systems as optical, X-ray and gravitational wave sources

## 1 Introduction

AM CVn systems are binaries, with periods less than 1 hour, in which a white dwarf accretes matter from a low-mass, helium-rich, object. The mass transfer is driven by gravitational wave radiation (e.g. Paczyński, 1967). They are observed at optical wavelengths as faint blue and variable objects. For a reviews of observational work on AM CVn systems, see Warner (1995); Solheim (2003). A discussion of the formation and modelling of the population can be found in Tutukov & Yungelson (1996), Nelemans et al. (2001a), hereafter NPVY and Podsiadlowski, Han & Rappaport (2003). Hils & Bender (2000) and Nelemans, Yungelson & Portegies Zwart (2001b) discussed the contribution of the AM CVn population to the Galactic unresolved gravitational wave background.

The number of known or suspected systems is rapidly increasing: six new (candidate) AM CVn systems have been found in the last five years: CE 315, with a period of 65.1 minutes (Ruiz et al., 2001); KL Dra, with a period of 25 minutes (Jha et al., 1998; Wood et al., 2002); “2003aw”, with a possible period of 34 minutes Woudt & Warner (2003b); RX J1914.4$`+`$2456 (V407 Vul), with a period of 9.5 minutes (Cropper et al., 1998; Ramsay et al., 2000, 2002b; Marsh & Steeghs, 2002), KUV 01584$``$0939 (ES Cet) with a period of 10.3 minutes (Warner & Woudt, 2002) and RX J0806.3$`+`$1527 with a period of 5.3 minutes (Israel et al., 2002; Ramsay et al., 2002a). RX J1914.4$`+`$2456 and RX J0806.3$`+`$1527 were discovered as ROSAT X-ray sources (Motch et al., 1996; Israel et al., 1999). For these two systems the interpretation as AM CVn systems has been questioned (Wu et al., 2002; Norton et al., 2002), but not for the other four new systems.

The last three of the new systems mentioned above are particularly interesting as theoretical models of the evolution of AM CVn systems and their current Galactic population (e.g. Tutukov & Yungelson, 1996, NPVY) show that the expected number of systems at periods shorter that $``$10 minutes is small, and an optically selected sample will be dominated by longer period systems. However, at these short periods the mass transfer rate is as high as $`10^610^8`$ $`\mathrm{M}_{}`$ yr-1, and therefore systems might be discovered as X-ray, rather than optical sources, as is the case for the two ROSAT sources. At these short periods the systems are also promising sources for low-frequency gravitational wave (GW) detectors (e.g. Israel et al., 2002; Marsh & Steeghs, 2002).

In this paper we determine the range of orbital periods in which AM CVn systems could be discovered as X-ray sources, using simple, average temperature and blackbody emission models and taking into account interstellar absorption. We estimate the number of such systems in the Galaxy and their brightness in the optical band. We compute the number of short-period ($`P<`$ 1500 s) AM CVn systems that could be detected by the planned low-frequency gravitational wave radiation (GWR) detector in space, LISA, and focus on the subset of these systems that are expected to have optical and/or X-ray counterparts, as combined observations will lead to a better understanding of the physical state, formation and population of AM CVn systems. We will use the term “AM CVn system” throughout the paper for all double-degenerate interacting binaries with helium secondaries, irrespective of their optical properties.

In Sect. 2 we review our Galactic model of AM CVn systems, discuss the changes we made with respect to NPVY and the uncertainties in this model. The emission of optical, X-ray and gravitational wave signals from AM CVn systems is discussed in Sect. 3. The selection effects and the resulting expected population of sources detectable in the optical and X-ray bands, and as GWR sources are presented in Sect. 4. In Sect. 5 we consider the LISA sources that have optical and/or X-ray counterparts. Discussion (Sect. 6) and conclusions (Sect. 7) follow.

## 2 A model for the current Galactic population of AM CVn systems

### 2.1 Formation of AM CVn systems

For a detailed discussion of the formation of AM CVn systems we refer to NPVY, and references therein, so here we only give a brief summary. There are two main routes for the formation of AM CVn systems: (i) via a phase in which a double white dwarf loses angular momentum due to gravitational wave radiation and evolves to shorter and shorter periods to start mass transfer at periods of a few minutes and (ii) via a phase in which a non-degenerate helium star transfers matter to a white dwarf accretor.

Recently Podsiadlowski et al. (2003) proposed a third scenario to form AM CVn binaries, similar to the usual scenario for the formation of cataclysmic variables (CVs), via a phase of mass transfer in the course of which a low-mass star with initially almost exhausted hydrogen core or a small helium core becomes a helium white dwarf111Such a scenario has been proposed before (e.g. Tutukov et al., 1987) as a formation scenario for ultra-compact X-ray binaries. We found that with the assumptions used in this paper (see below) this scenario is relatively unimportant for the short-period systems discussed here (giving less than two per cent of the systems), so for the remainder of the paper we leave it out of our model, but briefly come back to it in the discussion (Sect. 6).

One of the main conclusions of NPVY was that the uncertainties in the models for the formation of AM CVn systems resulted in widely varying numbers of expected systems in the Galaxy. The most optimistic and most pessimistic cases of NPVY gave roughly 80 and 18 million systems in the Galaxy at present. The uncertainty is particularly severe for the shortest period AM CVn systems, which probably have to form through the double white dwarf scenario. In NPVY we showed that at the start of the mass transfer between double white dwarfs there is a phase of “direct impact accretion”, in which the gas stream flowing from the donor to the accretor hits the accretors surface directly without forming an accretion disc (see also Webbink, 1984). In this phase the stability of the mass transfer depends crucially on the tidal (or other) coupling between the spin of the accretor and its companion ( NPVY, see for a detailed discussion of the stability Marsh, Nelemans & Steeghs 2003).

The formation of systems formed via the helium star scenario is also uncertain because of the possibility of so called edge-lit-detonation (ELD) of the accreting white dwarf before the system reaches the AM CVn phase (see NPVY, ). In the current study we use the most optimistic model (Model II of NPVY, ), in which we assume efficient coupling in the direct impact phase, and that accretion of at least 0.3 $`\mathrm{M}_{}`$ of He is necessary for ELD. But we note that recent calculations, including rotation, of helium accreting white dwarfs suggest that ELD might not happen at all (Yoon & Langer, 2003). This would increase the formation rate of AM CVn systems via the helium star scenario in our most optimistic model by about 15 per cent. For the short period systems discussed here the total uncertainty is probably of the order of a factor 10 at least.

We note that we do not take into account the unstable helium burning on the surface of carbon-oxygen white dwarfs accreting helium at high rates and the possible associated mass loss (see, e.g. Kato & Hachisu, 1999).

### 2.2 Galactic model, star formation history and stellar evolution

Instead of the simple Galactic model we used before, we implemented a star formation rate (SFR) as function of time and position in the Galaxy based on the model of Boissier & Prantzos (1999). We obtain the SFR as function of $`R`$ and $`t`$ by linear interpolation between values in a table (S. Boissier, private communication), and assume that the probability of a star being born with a radius $`R`$ or smaller goes with the integrated SFR, i.e.

$$P(R,t)=\frac{_0^R\text{SFR}(R^{},t)\mathrm{\hspace{0.17em}2}\pi R^{}dR^{}}{_0^{R_{\mathrm{max}}}\text{SFR}(R^{},t)\mathrm{\hspace{0.17em}2}\pi R^{}dR^{}},$$
(1)
where $`R_{\mathrm{max}}`$ is the maximum extent of the disc (19 kpc) and we assume that the stars do not migrate radially in time. The age of the Galaxy in this model is 13.5 Gyr. As the Boissier & Prantzos (1999) model gives the SFR projected onto the plane of the Galaxy, we assume a $`z`$-distribution,

$$P(z)\text{sech}(z/z_h)^2,$$
(2)
as in Nelemans et al. (2001c), where $`z_h`$ = 200 pc, neglecting its age and mass dependence. We have added a bulge to the model, simply by doubling the SFR in the inner 3 kpc of the Galaxy compared to Boissier & Prantzos (1999). As a result the total mass in this region is 2.6 1010 $`\mathrm{M}_{}`$ at the current age of the Galaxy, consistent with kinematic and micro-lensing results (e.g. Klypin et al., 2002). We distribute the stars in the bulge according to

$$\rho _{\mathrm{bulge}}=\mathrm{exp}((r/0.5\text{ kpc})^2),$$
(3)
where $`r=\sqrt{x^2+y^2+z^2}`$.

In Fig. 1 we plot the resulting distribution of stars in the $`xz`$ plane of the Galaxy, showing the $`z`$-distribution of stars and the change in the radial extent of the star formation as function of time. For the SFR integrated over $`R`$ as function of time, see Fig. 2 (left panel).

The important effect resulting from the different choice of Galactic model is that for most of the time in the Galactic history, the SFR is lower in the Boissier & Prantzos (1999) model than in the ‘exponential’ model used in NPVY (see Fig. 2); in particular, the current SFR is lower than the one we obtained previously (1.9 vs 3.6 $`\mathrm{M}_{}`$ yr-1). This has a strong effect on the current birth rates of AM CVn systems, since the latter follow star formation history with a typical delay of several Gyr (Tutukov & Yungelson, 1996; Nelemans et al., 2001a). However, the integrated SFR remains quite similar ($`8.2\times 10^{10}`$ vs. $`8\times 10^{10}`$ $`\mathrm{M}_{}`$ ), partly due to the different Galactic ages (13.5 vs. 10 Gyr). Another characteristic of our new model is that the Galactic distribution is more centrally concentrated than with an exponential radial density distribution , mainly due to the inclusion of the bulge. In Fig. 3 we show the distribution of distances from the Earth of the population of short-period AM CVn systems in the Galaxy for the two models. The concentration in the Galactic centre with the new model is clearly visible. As the longer period AM CVn systems are older, their distribution is even more centrally concentrated

Another new element in our current study is the use of the IMF as proposed by Kroupa, Tout & Gilmore (1993), rather than the Miller & Scalo (1979) IMF as used in NPVY. Because the Kroupa et al. (1993) IMF is steeper for masses above 1 $`\mathrm{M}_{}`$ , the same star formation history leads to fewer AM CVn systems, especially for the helium star formation scenario, the efficiency of which is reduced by about 25 per cent (see Sect. 2.3).

We use the population synthesis code SeBa, developed by Portegies Zwart & Verbunt (1996, see also ). We assume a 50 per cent binary fraction, an initial distribution of separations in binary systems ($`a`$) flat in $`\mathrm{log}a`$ between contact and 106 R, a flat mass ratio distribution and a thermal eccentricity distribution.222 All binaries discussed here will be circularized to somewhat smaller separation before the first phase of mass transfer and thus the shape of the initial eccentricity distribution hardly influences the number of AM CVn systems. The current study is based on a model with 250,000 ZAMS binaries. For the treatment of the stellar evolution and the mass transfer we refer to our earlier papers, but recall that we use a common envelope formalism based on the angular momentum balance for the first episode of unstable mass transfer in the system, as detailed in Nelemans et al. (2000, 2001c). Two small changes that turn out to be important are in the formation and evolution of white dwarfs. We updated our treatment of the growth of the core in giants as described in Nelemans et al. (2001c), by including a simple model of the second dredge-up to be more in line with detailed stellar evolution models (e.g. Hurley et al., 2000). The second change is in the treatment of white dwarf cooling (see Appendix A) and white dwarf mass – radius relation, for which we now use the fit by Eggleton, as given in Verbunt & Rappaport (1988), rather than the Zapolsky & Salpeter (1969) relation, which is not appropriate for massive white dwarfs.

The last change we made compared to NPVY is in the maximum rate at which a main sequence star can accrete before it expands and mass will be lost from the system. We used to limit this to accretion at the thermal time scale of the accretor. However, as is shown e.g. by Kippenhahn & Meyer-Hofmeister (1977), higher accretion rates do not necessary lead to expansion of the accretor. We therefore use the treatment of Pols & Marinus (1994) which takes the reaction of the radius of the accreting star upon mass gain into account.

### 2.3 Total population of AM CVn systems in the Galaxy

With the changes described in the previous section, we obtain, for the two formation channels discussed here, the current birth rates and total number of AM CVn systems in the Galaxy as shown in Table 1 and Fig. 2 (compare Table 1 in NPVY). The difference between the numbers presented here and in NPVY are due to the use of the new Galactic mode, different IMF and the changes in the treatment of stellar evolution described in Sect. 2.2. All these change the birth rate of AM CVn systems at any time, e.g. at the current Galactic age the birth rate of the white dwarf family goes down by $``$60 per cent due to the Galactic model and is reduced further by another $``$40 per cent due to the changes in treatment of white dwarf formation and mass – radius relation. The birth rate of the helium star family goes down by $``$60 per cent due to the Galactic model, $``$25 per cent due to the IMF and $``$40 per cent due to the higher accretion rate limit for main-sequence stars (because less systems pass through common envelopes and become close enough to produce AM CVn systems). Partly compensated by the change in the age of the Galaxy this results in 55 and 65 per cent fewer systems present in the Galaxy for the two families, respectively.

In Fig. 4 we show a histogram of the shortest period systems as function of orbital period, indicating the systems in the direct impact phase and the periods of the observed (candidate) AM CVn systems.

## 3 Modelling the optical, X-ray and GW emission of AM CVn systems

In order to be able to compare our results with (future) observations we have to simulate the emission produced by our model systems. We determine the temperature and luminosity of the different emission sources and determine how much energy is emitted in a certain band, assuming the emitted spectrum is a blackbody. For the present study we only discuss emission in $`V`$-band and the ROSAT 0.1–2.4 keV X-ray band plus gravitational wave radiation. We chose the ROSAT band, because the two candidate systems were found as ROSAT sources and because we can compare our results to the ROSAT all-sky survey. However, in Sect. 6 we briefly discuss the expected results for the Chandra and XMM.

There are four main emission sites: the accretion disc and boundary layer between disc and accreting white dwarf, the impact spot in the case of direct impact accretion, the accreting star and the donor star.

### 3.1 Optical emission

We assume that half of the accretion luminosity is radiated by the disc (the other half in the boundary layer), i.e.

$$L_{\mathrm{disc}}=\frac{1}{2}GM\dot{m}\left(\frac{1}{R}\frac{1}{R_{\mathrm{L1}}}\right)\text{erg s}^1,$$
(4)
where $`M`$ and $`R`$ are the mass and radius of the accretor and $`R_{\mathrm{L1}}`$ is the distance of the first Lagrangian point to the centre of the accretor and $`\dot{m}`$ is the mass transfer rate. Because the two white dwarfs are so close together the usual assumption that the matter falls in from infinity is not valid (Han & Webbink, 1999). We model the optical emission from the disc like in NPVY, according to Wade (1984), as a single temperature disc that extends to 70 per cent of the Roche lobe of the accretor, radiates as a blackbody and use the bolometric correction appropriate for its temperature.

We model the emission from the donor as the cooling luminosity of the white dwarf, using approximations to the models of Hansen (1999), see the Appendix A.

We model the emission from the accretor as being the unperturbed cooling luminosity of the white dwarf.

In our simple model we do not take into account such effects as irradiation of the donor by the impact spot or disc, heating of the accretor by radiation from the disc/boundary layer (e.g. Pringle, 1988), compressional heating (e.g. Townsley & Bildsten, 2002) or tidal heating of the two white dwarfs prior to contact (Iben et al., 1998).

In Fig. 5 we show the range of absolute $`V`$-magnitudes (MV) and temperatures of the different optical components described above, for systems forming from double white dwarfs. The top panel shows the range of absolute magnitudes of the discs as function of period (in grey). The two symbols show the observed absolute magnitudes for AM CVn and GP Com \[using distances of 235 pc for AM CVn (C. Dahn333On behalf of the USNO CCD Parallax Team., private communication) and 70 pc for GP Com (Thorstensen, 2003)\], while the bars in the bottom panel show the range of effective temperatures for these systems. The upper and lower limit of the shaded regions are for two combinations of masses that bracket the range of masses of the majority of the systems: 0.3 + 1.0 and 0.15 + 0.25 $`\mathrm{M}_{}`$ respectively.

Because for the accretor and donor we assume that their luminosity is the cooling luminosity, their magnitudes depend on their ages, rather than on the orbital period. However, during the AM CVn phase there is a direct relation between time and period (e.g. figure 4 in NPVY) so for given cooling luminosity at the start of the mass transfer we can plot the magnitude of the cooling white dwarf as function of orbital period. In Fig. 5 we plot that curve for the extreme case in which the cooling luminosity at the starts of the mass transfer is that of a newly born, very hot and luminous white dwarf. In reality the donor and accretor are typically of the order of a Gyr old, so they evolve along a track that lies below these curves. Towards the longer periods, the age of the system is a few Gyr anyway, so the curves for all initial luminosities converge.

The few observed systems fall in the grey area in Fig. 5 indicating that our simple model is at least able to understand these systems satisfactorily except the temperature of GP Com444Indeed the continuum emission of GP Com originates at the reheated accreting white dwarf, rather than in the disc (see Bildsten et al., in preparation). The same model can be adopted for the systems that form via the helium-star donors. In these cases, however, the curves shift to longer periods (see NPVY).

### 3.2 X-ray emission

In modelling the X-ray emission of AM CVn systems we distinguish two cases: with and without disc.

If there is no disc, the gas stream hits the companion directly and a small area of the surface will be heated. We assume that the full accretion luminosity is radiated as a blackbody with a temperature given by

$$\left(\frac{T_{\mathrm{BB}}}{T_{}}\right)^4=\frac{1}{f}R^2L_{\mathrm{acc}},$$
(5)
where $`L_{\mathrm{acc}}`$ and $`R`$ are in solar units and $`L_{\mathrm{acc}}`$ is defined by Eq. (4). The fraction $`f`$ of the surface that is radiating depends on the details of the accretion, but will in general be small. We used $`f=0.001`$, which is consistent with the observed X-ray emission of V407 Vul and simple arguments regarding the impact stream (Marsh & Steeghs, 2002).

If there is an accretion disc, the X-ray emission can come from the boundary layer between the inner disc and the accretor, where we assume half of the accretion luminosity is radiated, as suggested e.g. by Pringle (1977), or it might have a non-thermal origin. The average temperature of the disc (Wade, 1984), which we use to estimate the optical emission is never high enough to produce X-rays (Bath et al., 1974). Even though from observations of cataclysmic variables it is clear that the real X-ray emission mechanism is much more complicated (e.g. Hoare & Drew, 1991; Mauche et al., 1991; van Teeseling & Verbunt, 1994; van Teeseling et al., 1996), we model the X-ray emission of systems with a disc for the moment as coming from the boundary layer, with temperature (Pringle, 1977)

$$T_{\mathrm{BL}}=5\times 10^5\left(\frac{\dot{m}}{10^{18}\text{ g/s}}\right)^{\frac{2}{9}}\left(\frac{M}{M_{}}\right)^{\frac{1}{3}}\left(\frac{R}{5\times 10^8\text{ cm}}\right)^{\frac{7}{9}}\mathrm{K}.$$
(6)
In Fig. 6, top panel, we show how the luminosity in X-rays compares to the total accretion luminosity. For the X-ray luminosity we consider here the emission in the ROSAT band, from 0.1 – 2.4 keV. The temperature of the impact spot or boundary layer is shown in the bottom panel. The sharp drop in the X-ray luminosity and temperature at the bottom of the darker region at $`\mathrm{log}P2.9`$ is caused by the change from direct impact to disc accretion in the low-mass system. We assume that the transition from the impact temperature to the boundary layer temperature is instantaneous. This transition for the more massive system occurs just left off the plot, where the evolution is so fast that hardly any systems are expected to be observed.

At periods larger than 1000 – 2000 s the mass transfer rates drop to below 10-10 $`\mathrm{M}_{}`$ yr-1, the temperature drops below 105 K and hardly any X-ray emission is produced in the boundary layer. The inferred temperatures and X-ray luminosities of RX J0806 and V407 Vul are plotted as the bars in the figures. For RX J0806, we use the X-ray flux ($`4\times 10^{10}`$ erg cm-2 s-1), distance (0.1 – 1 kpc) and blackbody temperature ($`6_4^{+13}\times 10^5`$ K) as given by Israel et al. (2002). For V407 Vul blackbody temperatures of 43 and 56 eV are derived by Motch et al. (1996) and Ramsay et al. (2002b) respectively, while the latter authors derive an X-ray luminosity of $`1.2_{0.3}^{+4.1}\times 10^{33}`$ erg s-1, for a distance of 100 pc and we plot these limits for a distance range 100 – 400 pc (Ramsay et al., 2000). The X-ray luminosity of AM CVn is taken from Ulla (1995). Except for the luminosity of RX J0806, our very simple model agrees well with the observations.

### 3.3 Gravitational wave emission

We assume that AM CVn systems emit GWR as two point sources, i.e. neglecting any effect of the deformation of the donor, which fills its Roche lobe and the mass in the accretion disc, if present (see Rezzolla et al., 2001, for a discussion of the small effect of deformation for cataclysmic variables). The amplitude of the GWR signal, as can be detected by LISA, is given by (e.g. Evans et al., 1987)

$$h=0.5\times 10^{21}\left(\frac{}{\mathrm{M}_{}\text{ }}\right)^{5/3}\left(\frac{P_{\mathrm{orb}}}{1\mathrm{hr}}\right)^{2/3}\left(\frac{d}{1\mathrm{k}\mathrm{p}\mathrm{c}}\right)^1,$$
(7)
where $`=(Mm)^{3/5}/(M+m)^{1/5}`$ is the chirp mass, $`P_{\mathrm{orb}}`$ is the orbital period and $`d`$ is the distance to the system.

## 4 Results I: the detectable population of short-period AM CVn systems

### 4.1 Selection effects

To compare our models with the number of observed systems, we have to decide which of the model systems would be detectable.

For the optical emission we focus on model systems that are brighter than $`V=20`$ mag, which is similar to the optical magnitude of the observed short period AM CVn candidates. To take the interstellar reddening into account we use the Sandage (1972) model

$`A_\mathrm{V}(\mathrm{})`$ $`=`$ $`0.165{\displaystyle \frac{(\mathrm{tan}(50\mathrm{°})\mathrm{tan}(b))}{\mathrm{sin}(b)}}\text{ for }b<50\mathrm{°}`$ (8)
$`=`$ $`0\text{ for }b50\mathrm{°},`$
and $`A_\mathrm{V}(d)`$ follows from our Galactic model \[Eq. (2)\]

$$A_\mathrm{V}(d)=A_\mathrm{V}(\mathrm{})\mathrm{tanh}\left(\frac{d\mathrm{sin}(b)}{z_h}\right).$$
(9)
For the sources detectable in X-rays we select the systems with an X-ray flux in the 0.1 – 2.4 keV (ROSAT) band higher than $`10^{13}`$ erg s-1 cm-2. To estimate this flux we use the intrinsic flux in this band (Sect. 3.2), the distance and an estimate of the Galactic hydrogen absorption (Morrison & McCammon, 1983). We modelled the hydrogen column density after Predehl & Schmitt (1995) as

$$N_\mathrm{H}=0.179A_\mathrm{V}(d)\mathrm{\hspace{0.17em}10}^{22}\mathrm{cm}^2.$$
(10)
### 4.2 Systems detectable in the optical and X-ray bands

In Fig. 7 we show the sample of the systems that would be detectable in the optical and/or X-ray bands, taking into account the selection effects described above. The relevant numbers are given below in the text, between square brackets (see also Table 3). (We should recall that the results we discuss represent one possible random realization of the model for current Galactic population of AM CVn systems, so all numbers given are subject to Poisson noise). We limited ourselves to the short period systems, with $`P<1500`$ s. At longer periods the optical emission from the disc will completely dominate and the detectable sample has been discussed in NPVY. Fig. 7, top panel, shows the systems detectable in X-rays. Pluses are  systems only detectable in X-rays, while the filled symbols mark systems also detectable in the $`V`$-band.

At the shortest periods a significant fraction of the small number of systems that exist (see Fig. 4) produces enough X-rays that they are detectable even though the sources are close to the Galactic centre (Fig. 7). At longer periods the X-rays get weaker (and softer) so only the systems close to the Earth can be detected. Several of these systems are also detectable in the optical band (filled symbols). The  filled circles show the systems in the direct impact phase, where the optical emission comes solely from the donor as the impact spot is too hot to emit in the $`V`$-band. Systems with optical emission from a disc are shown as triangles . There are some  systems that are close enough to the Earth that the donor stars can be seen as well as the discs (filled squares). Above 600 s the systems descending from helium star donors appear and have high enough mass transfer rate to be X-ray sources, the closer ones also being visible in the $`V`$-band (filled triangles/squares, as these systems always have a disc).

The X-ray sources at larger distances (pluses) are almost all in the bulge and thus are concentrated in a small area on the sky. The closer systems, which often are also detectable in the $`V`$-band are more evenly distributed over the sky.

The bottom panel of Fig. 7 shows the 1230 ‘classical’ AM CVn systems, detectable only by the optical emission, for most systems only from their accretion disc (open triangles ). Of this population the ones closest to the Earth in some cases  also have a visible donor (open squares). A few  direct impact systems only have a visible donor, but are not X-ray sources (open circles). The majority of the optically detectable systems represented in the lower panel of Fig. 7 (with orbital periods between 1000 and 1500 s) are expected to show outbursts due to the viscous-thermal disc instability (Tsugawa & Osaki, 1997) which could enhance the chance of their discovery.

Comparing these results with the tentative classification of V407 Vul and RX J0806 as direct impact systems, detectable in X-rays and (just) in the optical band, shows that a few of such systems (filled circles) could be expected in the Galaxy, given all the uncertainties in population synthesis model, the emission model and the selection effects. Furthermore it shows that many more short period AM CVn systems might be detectable only as X-ray sources. Some have periods as short as three minutes. About one third of the X-ray sources is also detectable in the optical band.

### 4.3 Gravitational waves from short period AM CVn systems

In Nelemans et al. (2001b) we presented a study of the gravitational wave signal from Galactic binaries in which both stars are compact objects. In that study we estimated the contribution of AM CVn systems to the double white dwarf noise background and found that AM CVn systems hardly contribute, even though at frequencies between 0.3 and 1 mHz (periods between 33 and 105 minutes) they outnumber the detached double white dwarfs. This is because at these frequencies their chirp mass has fallen far below that of a typical detached system (see also Hils & Bender, 2000). We therefore did not consider them anymore. However at shorter periods, where the number of AM CVn systems is much smaller their chirp mass is similar to that of the detached systems from which they originate.

In Fig. 8 we show the Galactic population of detached double white dwarfs and AM CVn systems that can be resolved by LISA. *There are $``$11 000 resolved double white dwarfs as well as $``$11 000 resolved AM CVn systems.* Their distribution is shown by the grey shading. The 200 strongest sources of each type are shown as dots to enhance their visibility. Figure 8 clearly shows that AM CVn systems may be important sources for LISA. Compared to our previous study the ‘confusion limit’, i.e. the frequency at which the unresolved background disappears (see Nelemans et al. (2001b) for details) moves to slightly higher frequency because in our current model there are 213 million detached close double white dwarfs, instead of 110 million before, due to the use of different IMF, Galactic model and the changes in stellar evolution.

The crucial difference between AM CVn systems and detached double white dwarfs in the short period range is the sign of their period derivatives. The periods of AM CVn systems increase due to mass transfer, while the periods of detached systems decrease due to GW losses555Note that systems formed from helium stars and evolved secondaries in CVs have decreasing periods before they reach their period minimum (e.g. Tutukov & Fedorova, 1989; Podsiadlowski et al., 2003). However, in this phase they are (much) brighter than after the period minimum and are not consistent with being classified as “classical” AM CVn systems.. Assuming conservative mass transfer, the change of the GW frequency ($`f=2/P`$) of AM CVn systems formed from double white dwarfs can be derived from

$$\frac{\dot{f}}{f}=\frac{\dot{P}}{P}=\frac{3}{2}\frac{\dot{a}}{a}=3\frac{\dot{J}_{\mathrm{GWR}}}{J_{\mathrm{orb}}}3(1q)\frac{\dot{m}}{m},$$
(11)
where $`m`$ is the mass of the donor star and $`q=m/M`$. We can use the equilibrium mass transfer rate as derived in Marsh et al. (2003) to get

$$\dot{f}=5.8\times 10^7f^{11/3}\left(\frac{}{\mathrm{M}_{}\text{ }}\right)^{5/3}\left(\frac{\zeta _{r_L}\zeta _2}{2+\zeta _2\zeta _{r_L}2q}\right)\text{Hz s}\text{-1},$$
(12)
where $`\zeta _2`$ is the logarithmic derivative of the secondary radius with respect to mass and $`\zeta _{r_L}=\mathrm{d}\mathrm{ln}(R_L/a)/\mathrm{d}\mathrm{ln}m1/3`$ (see Marsh et al., 2003). We wrote the equation in the same form as the corresponding equation for frequency change in detached systems (e.g. Eq. (6) in Nelemans et al., 2001b), with the (positive) term in parentheses giving the change with respect to the detached case. This term is of order unity at the onset of the mass transfer (i.e. at highest frequency, above $``$5 mHz), decreasing to about 0.1 at 2 mHz.

LISA will be able to measure the frequency change of a number of the resolved binaries. In Nelemans et al. (2001b) we estimated the number of systems for which this can be done by requiring that in the course of the experiment the frequency would change by one frequency bin. However, it really is the phase shift resulting from a frequency change that will be measured (e.g. Webbink & Han, 1998), so the requirement has to be that over the course of the experiment the phase shifts by a certain amount, compared to a system without frequency change. Writing the phase as an expansion in $`f,\dot{f}`$ and $`\ddot{f}`$, this means $`\dot{f}`$ can be measured if

$$\mathrm{\Delta }\mathrm{\Phi }=\dot{f}\frac{T^2}{2}>\mathrm{\Delta }\mathrm{\Phi }_{\mathrm{min}},$$
(13)
where $`\mathrm{\Delta }\mathrm{\Phi }_{\mathrm{min}}`$ is the minimum phase shift that can be measured and $`T`$ is the duration of the experiment. Analogously $`\ddot{f}`$ can be measured if

$$\mathrm{\Delta }\mathrm{\Phi }=\ddot{f}\frac{T^3}{3}>\mathrm{\Delta }\mathrm{\Phi }_{\mathrm{min}}.$$
(14)
Measuring the frequency change provides important information, as the chirp mass can be measured directly for the detached systems. Combined with the wave amplitude $`h`$ \[Eq. (7)\], this yields the distance. For AM CVn systems this argument does not hold, as $`\dot{f}`$ is determined by the mass transfer \[Eq. (12)\]. However, for an assumed mass–radius relation of the donor star, both the model mass transfer rate and the mass of the donor star (and thus the chirp mass of the system to within a factor of $``$2) as function of period are fixed (e.g. Vila, 1971; Faulkner, 1971). Thus a measurement of $`\dot{f}`$ can test mass–radius models, which is especially interesting because more realistic mass–radius relations for low-mass white dwarfs become available (e.g. Bildsten, 2002). In Table 2 we show the number of systems for which we can measure the frequency change as function of the minimum phase shift.

Finally it would be very interesting to measure the second derivative of the frequency, as that gives a direct check whether the change in frequency is due to GWR only, or whether other effects (like tidal torques etc.) are important (Webbink & Han, 1998, Phinney, private communication). For pure GW losses we have666Note that Webbink & Han (1998) give 5/3 as numerical factor, which is appropriate for the corresponding equation for the second derivative of the orbital period.

$$\ddot{f}=\frac{11}{3}\frac{\dot{f}^2}{f}.$$
(15)
As can be seen from Table 2, the number of systems for which the second derivative can be measured is small, even for a more optimistic assumption of a total LISA operation time of 5 yr. For none of the AM CVn systems we expect its second derivatives to be measurable. This is due to the fact that the AM CVn systems form from pairs of low-mass ($`M0.5\mathrm{M}_{}\text{ }`$) white dwarfs (e.g. Nelemans et al., 2001a) which have larger radii and thus do not penetrate into the highest frequency range where $`\ddot{f}`$ is large (see Fig. 8). In this range there are only massive ($`M>1\mathrm{M}_{}\text{ }`$) detached double white dwarf pairs that merge rather than become AM CVn systems.

## 5 Results II: Short-period AM CVn systems detectable by LISA with optical and/or X-ray counterparts

In Figure 9 we show the number distributions vs. orbital periods for LISA sources that have optical and/or X-ray counterparts, direct impact systems that are GWR and X-ray sources, and the total number of systems that are optical and/or X-ray sources (see also Table 3). For comparison, we indicate also the distribution of the total number of AM CVn systems with $`P<1500`$ s as the white histogram in each panel. The systems with different properties are marked by different hatching styles. We will discuss these groups now in more detail.

Resolved LISA sources (GWR)
The top left panel shows the 11,000 systems that are expected to be resolved by LISA (grey histogram, also shown in Fig. 8), most of which do not have optical or X-ray counterparts. There are 151 short-period systems that have X-ray counterparts (GWR+X) and 143 systems, mainly with periods between 500 and 1000 s, that have optical counterparts (GWR + Opt). Only 31 systems have both optical and X-ray counterparts (GWR + Opt + X).

Direct impact accretors
The top right panel shows the 2680 systems that are in the direct impact phase of accretion (grey histogram, as in Fig. 4). Most of these systems (2057) are expected to be resolved by LISA \[DI (GWR)\], while at the shortest periods there are 105 X-ray systems that are in the direct impact phase \[DI (X)\].

Optical sources
The bottom left panel shows (in grey) that for periods below 25 minutes there are 1336 systems that have $`V<`$ 20 mag (grey histogram). These are the systems shown in the bottom panel of Fig. 7, and are the local population of AM CVn systems. Some (106) of these systems are also X-ray systems (Opt + X) and are shown as the filled symbols in the top panel of Fig. 7. Many of the short period systems will show eclipses (even more than the detached systems which are discussed in Cooray et al., 2003) which would compliment the GWR data and provide information on the radii of the components and the inclinations of the systems. E.g. a crude estimate for a typical system with initial parameters $`0.25+0.6\mathrm{M}_{}\text{ }`$ gives a probability for eclipsing the accreting star of about 30 per cent at $`P=1000s`$, and even higher for eclipsing (part of) the accretion disc.

X-ray sources
Finally, the bottom right panel shows the 326 X-ray sources that are also shown in the top panel of Fig. 7. The majority is only detectable in X-rays and many are in the direct impact phase. Only a few direct impact systems are detectable in the optical band, by the radiation from the cooling white dwarf (the filled circles in Fig. 7). At longer periods we expect detectable X-ray emission from the boundary layers of systems formed from helium stars, some of which are also detectable in the $`V`$-band (Opt + X). These X-ray sources in the ROSAT band in our model are concentrated in the Galactic plane, especially towards the Galactic centre. There are 193 systems (60 per cent of the total) in an area between Galactic longitudes $``$50 and +50 degrees and Galactic latitude between $``$10 and +10 degrees. The ROSAT all-sky survey777http://www.xray.mpe.mpg.de/rosat/survey/rass-fsc/ found $``$4500 sources in this region, so the short-period AM CVn systems can only account for a small fraction of these sources.

It is clear from Fig. 9 that LISA is expected to be very effective in detecting the short-period AM CVn systems. A large fraction of these will be in the direct impact phase. At short periods, a significant fraction of the LISA sources might be detectable in X-rays, while at longer periods a small fraction might be optical sources.

LISA sources with optical and/or X-ray counterparts are important, as LISA will measure a combination of all the parameters that determine the GWR signal (frequency, chirp mass, distance, position in the sky and inclination angle, e.g. Hellings, 2003), so if some of these parameters (period, position) can be obtained from optical or X-ray observations the other parameters can be determined with higher accuracy. This is particularly interesting for the distances, inclinations and masses of the systems, which are very difficult to measure with other methods.

## 6 Discussion

We stress here that AM CVn systems form from binaries with a rather small range of initial parameters and thus that their birth rate is very sensitive to small changes in assumptions in binary population synthesis calculations. One has to keep in mind that a lot of these assumptions are badly constrained, which means that the number of progenitor systems is rather uncertain (e.g. Nelemans et al., 2001c). On top of that, there are uncertainties connected with the formation of AM CVn systems from their immediate progenitors as discussed in Sect. 2.1 and NPVY, making the AM CVn birth rates uncertain by at least an order of magnitude. The models we use for the X-ray and optical emission are very simplified as well. The interstellar absorption model might very well underestimate the X-ray absorption for the systems in the Galactic centre, which comprise a large fraction of the systems that are detectable only in X-rays (pluses in Fig. 7). We furthermore recall that a number of (possibly important) effects like irradiation and compressional heating are not taken into account in our simple emission models.

For comparison we also computed the results for X-ray detectors sensitive in the range 0.1 – 15 keV (similar to the Chandra and XMM bands). Because most of our spectra are rather soft, the flux ($``$) in this band is generally not much larger than in the ROSAT band: 80 per cent has $`_{\mathrm{XMM}}/_{\mathrm{ROSAT}}<`$ 1.5; 96 per cent has $`_{\mathrm{XMM}}/_{\mathrm{ROSAT}}<3`$. However, the improved sensitivity of Chandra and XMM obviously will have an impact. For instance the Muno et al. (2003) Chandra observations of the Galactic Centre have a completeness limit of 3 10-15 erg cm-2 s-1, almost two orders of magnitude deeper than our assumed ROSAT limit. The number of X-ray sources in the 0.1 – 15 keV band detectable down to $`10^{14}`$ erg cm-2 s-1 is 644 and 1085 down to $`10^{15}`$ erg cm-2 s-1. In the Chandra mosaic image of the Galactic centre (Wang et al., 2002), roughly down to $`10^{14}`$ erg cm-2 s-1 there are $``$1000 point sources, while our model gives 16 X-ray systems in this region. In the tiny Muno et al. (2003) region, where they identify more than 2000 sources, we predict 3 AM CVn systems with X-ray flux above $`10^{15}`$ erg cm-2 s-1.

As briefly discussed in Sect. 2, we did not include AM CVn systems formed like CVs in our models. The reason is that for the production of the short period systems, which we consider in this paper, this channel is not effective. An evolutionary track for a donor with an initial mass of 0.96 $`\mathrm{M}_{}`$ and a hydrogen-exhausted core $`0.05M_{}`$ (sequence nr. 4 in Fedorova & Ergma, 1989) reaches a minimum period of 8.4 minutes and spends 3.8 107 yr in the 25 – 8.4 minute period range (A.V. Fedorova, private communication). Taking this time as typical, and assuming a birth rate of 6.4 10-5 yr-1 for systems evolving to $`P<25`$ minutes, as in the standard model of Podsiadlowski et al. (2003), we expect 2400 post-CV ultra-compact systems in the Galaxy, compared to the 140,000 from the channels discussed here. However, if these channels are ineffective, e.g., due to the instability of the mass transfer and/or destruction of accretors by edge-lit detonations (NPVY), the AM CVn systems from CVs may become more important even at short periods.

## 7 Conclusion

We discussed a model for the Galactic population of short-period AM CVn systems and the uncertainties in this model. We find a total population of AM CVn systems in the Galaxy with periods shorter than 1500 s of $``$140,000 systems. We estimate that the final uncertainty on this number is at least a factor of 10. The model presented here gives an optimistic estimate of the total number of systems, i.e. all quoted numbers could be (substantially) smaller. We calculated the optical, X-ray and gravitational wave signals that could be detectable, using simple emission models and selection effects, and discussed the overlap between the subpopulations which may be discovered by different detection methods.

In this optimistic model we find that there are a few systems expected that are in the direct impact phase and are detectable as optical and X-ray sources, where the optical emission comes from the cooling of the (possibly heated) white dwarfs. Thus if RX J0806 and V407 Vul turn out to be short period AM CVn systems, that would argue for an optimistic model as in a pessimistic model there would be no such systems that are observable. Strohmayer (2002) and Hakala et al. (2003) have reported a shortening of the periods in V407 Vul and RX J0826, arguing against this interpretation despite the fact that these stars might be helium-donor or CV systems before the period minimum (see also Woudt & Warner (2003a) who argue that there is not sufficient evidence for a shortening of the period in RX J0806).

At the observed periods of RX J0806, V407 Vul and ES Cet, the fraction of systems that is in the direct impact phase is 80, 40 and 30 per cent respectively. The total number of direct impact systems in the Galaxy is 2680. We further find that there are 220 (very) short period AM CVn systems that are only detectable in the X-ray band, some with orbital periods as short as 3 minutes.

We conclude that the future gravitational wave detector LISA is expected to discover thousands of short-period AM CVn systems, about ten per cent of the total number of short-period systems in the Galaxy. The number of systems with periods below 1500 s that is detectable in the optical band is 1336 (about one per cent of the total). Only 326 systems (0.25 per cent) are detectable in the ROSAT X-ray band. For the typical sensitivities of Chandra and XMM, the number of detectable sources over the whole sky is $``$ 1000, still only a tiny fraction of all X-ray sources.

For many of the LISA sources the change in frequency can be measured, which would provide a test of the mass – radius relation for the donor stars in these systems. LISA is expected to detect a similar number of detached double white dwarfs, for many of which the frequency change can be determined, yielding a direct measurement of the distances to these sources. Maybe for a few the second derivative of the frequency can be measured, which would give a direct test whether the frequency change is purely due to GWR, or whether tidal effects play a role.

Our results indicate that many short-period AM CVn systems may be detectable in more than one way. A total of 106 systems is detectable in the optical and X-ray bands. For the systems that are detectable with LISA, about 300 (3 per cent) have optical and/or X-ray counterparts, giving the opportunity to independently determine the position and orbital period. This leads to a better determination of the other parameters that can be derived from the GWR signal (i.e. inclination and masses) and would lead to unprecedented accuracy in the study of these extreme binaries.

## Acknowledgments

We thank the referee for comments that improved the paper We thank Samuel Boissier for providing the table of star formation rates and helpful discussion about the Galactic model. GN is supported by PPARC. LRY is supported by RFBR grant no. 03-02-16254 and Federal Science and Technology Program “Astronomy” (contract no. 40.022.11.1102). SPZ is supported by the KNAW and NWO.

## Appendix A White dwarf cooling

We model the cooling of white dwarf based on the curves of Hansen (1999). We abandoned the previous use of the Driebe et al. (1998) models, since their models with thick hydrogen envelopes cool slower than observations of detached double white dwarfs suggest and their thick hydrogen envelopes actually may be lost in hydrogen shell flashes (e.g. Nelemans et al., 2001c).

Our simple models have carbon-oxygen white dwarfs that cool faster than helium white dwarfs, but within the groups the most massive objects are brighter at given age. We model the curves as (see Fig. 10):

$$\mathrm{log}L=L_{\mathrm{max}}1.33\mathrm{log}(t/10^6\mathrm{yr}),$$
(16)
where $`L_{\mathrm{max}}`$ is a linear fit given by

$$L_{\mathrm{max}}=1(0.9M_{\mathrm{WD}})\text{ for }M_{\mathrm{WD}}>0.5\mathrm{M}_{}\text{ },$$
(17)
and

$$L_{\mathrm{max}}=1.41.33(0.45M_{\mathrm{WD}})\text{ for }M_{\mathrm{WD}}<0.5\mathrm{M}_{}\text{ }$$
(18)
(mass and luminosity in solar units). For white dwarf masses below 0.5 $`\mathrm{M}_{}`$ the luminosity is constrained to be below $`\mathrm{log}L/L_{\mathrm{}}=0.5`$, for more massive white dwarfs below $`\mathrm{log}L/L_{\mathrm{}}=2`$.


text:
# Long and short gamma-ray bursts, and the pulsar kicks

## Abstract

One of the mysteries that surround gamma-ray bursts (GRB) is the origin of two classes of events: long and short GRB. The short GRB are similar to the first second of a long GRB. We suggest that the short bursts are interrupted long bursts, we point out a plausible mechanism for the interruption, and we explain the observed time scales. The supernova-like central engine may contain a neutron star or a black hole surrounded by an accretion disk and jets. In the case of a neutron star, the same mechanism that is responsible for the pulsar kicks can disrupt the central engine, thus producing an interrupted, short GRB. If a black hole is produced, the absence of the kick ensures the long duration of the GRB. The time delay of the kick and the absence of the kick for a black hole are natural consequences of a model based on neutrino oscillations.

Observational evidence suggests that gamma-ray bursts are supernova-like phenomena. This evidence is stronger for the long GRB, which last longer than a few seconds and for which the optical afterglows have been detected. It is now established that the long GRB occur in star forming regions. The central engine is probably a sort of a supernova explosion, which is different from ordinary supernovae in that a large fraction of energy is concentrated in highly relativistic jets. On the other hand, for short GRB, the lack of optical afterglows leaves open several other possibilities. For example, neutron star and black hole mergers are probably not ruled out.

In this Letter we will argue that both long and short GRB can originate from the same supernova type central engine that produces either a black hole or a neutron star, respectively. The main reason why the short bursts are shorter is that the jet, pointing originally along the line of sight, is disrupted by the same dynamical mechanism that is responsible for the “pulsar kicks”. A model for such a mechanism based on neutrino oscillations predicts a delay on the time scale of order of a few seconds, which can naturally explain the durations of the interrupted bursts. The same kick mechanism does not act on a black hole. Thus, if the central engine is a black hole, the jet exists and maintains its direction along the line of sight for (much) longer than a second, which results in a long GRB. We suggest, therefore, that the long GRB come from central engines containing black holes, while the short ones originate from those with neutron stars.

The pulsar kick may also be the origin of unusually high rotational velocities of pulsars at birth. The kick is not expected to be a central force. Hence, it simultaneously kicks the pulsar and spins it up. At the time of the kick both the position and the angular momentum of the neutron star change dramatically.

Since the accretion disk and the jet are gravitationally bound to the neutron star, the pulsar kick can disrupt the jet entirely. An average pulsar velocity is 200-500 km/s, although about 15% of pulsars have velocities in excess of 1000 km/s. The neutron star is strongly gravitationally bound to a part of the accretion disk of the size, roughly, 10 neutron star radii, $`D10R_{_{\mathrm{NS}}}150`$ km. Therefore, it takes from 0.1 to 1 second after the kick for the neutron star to leave the central region of the accretion disk. This can either terminate the jet entirely, or can change its direction. The kick changes the angular momentum of neutron star, which affects the total angular momentum of the system and, thus, the direction of the jet. In either case, whether the jet is disrupted entirely or it is kicked off the line of sight, the kick terminates the GRB.

The two essential elements of this picture, illustrated in Fig.1, are the presence of a compact star in the central engine and the existence of a pulsar kick mechanism which acts on a neutron star within one second of the onset, and which does not act on a black hole. We will first discuss the phenomenological consequences of such a mechanism, and then we will argue that a microscopic model based on neutrino oscillations is suitable.

This simple explanation can account for the following observational features of GRB:

* short-time variability on the millisecond time scale for all bursts
* temporal and spectral similarities between the short bursts and the first seconds of the long bursts
* lack of the optical afterglows following the short GRB
* approximate fluence-duration proportionality
* characteristic time scale of about $`1`$s separating the two classes of bursts
* supershort bursts

We will discuss each of these items below.

Short-time variability is a natural consequence of a small size of the compact star at the center of the engine. The observed gamma rays are generated in the jet by some dynamics that may involve shocks and flow irregularities. However, the millisecond time scale in the temporal profile of bursts is suggestive of a small size dynamics that modulates the jet at its origin. Causality forces one to associate the millisecond time intervals with length scales of the order of a kilometer. Hence, a compact central object is called for. In addition, the smaller the central engine, the stronger is the field of gravity in the region where the jet originates. This can help accelerate the accreting material to higher speeds before it is ejected out into the jet. Given the tight energetic constraints on the central engine, it is desirable to extract as much energy as possible in the form of a narrow jet. Therefore, the presence of a compact star in the central engine is definitely a desirable feature of the model. An improvement in temporal sensitivity of observations in the millisecond range, in combination with improved theoretical modeling, may help distinguish black holes from neutron stars by measuring short-time variability, which may be different for those two cases.

Temporal similarities between the short and long bursts are to be expected if the short bursts are merely the interrupted long bursts. This is strongly supported by recent observations, which demonstrate that the first seconds of the long bursts are almost identical to the short bursts. Improved future analyses may reveal a softening of the short GRB toward the end, when the central engine is disrupted.

The lack of optical afterglows for the short bursts may be the consequence of a disruption of the central engine and the jet by the kick force. The isotropic emission of light in a supernova is insufficient to produce an observable afterglow at large red shift. Thus, after the jet is disrupted, one does not expect an afterglow of observable luminosity. In the black hole case, the jet remains stable for the entire duration of the (long) GRB, and the same jet can be detected at a later time as an afterglow. However, an early optical afterglow for a short GRB is also possible if it results from an interaction of the interrupted jet with ambient matter.

The approximate fluence-duration proportionality is also explained in our model because the short bursts are supposed to be the interrupted long ones. In the absence of the interruption by the kick force, the short burst would have continued and would have produced the total luminosity typical of a long burst. The disruption of the jet ends the emission of energy in the direction of the observer.

The fluence – duration power law is softer for the short GRB then it is for the long ones. This observation fits in our model because the black holes can make more powerful engines than the neutron stars. This is because the black holes can have a larger mass and a smaller size simultaneously.

The characteristic time scale of order 1 s arises naturally in a model that explains the pulsar kicks by neutrino oscillations. There are actually two viable mechanisms, based on the resonant and non-resonant oscillations, respectively. The unknown neutrino parameters determine which of the two (mutually exclusive) mechanisms is a work. In the case of resonant oscillations, the force is applied from the onset of the neutron star cooling, with no delay. However, for different neutrino parameters, the off-resonant mechanism is operative, and the kick starts only after some time $`\tau 1s`$, which is a natural time scale to explain the difference between the short and the long GRB. This pulsar kick mechanism requires that a singlet fermion, having a keV mass and mixed with the electron neutrino, be present in the neutrino spectrum. Serendipitously, the same singlet neutrino can be the dark matter in the universe. The kick begins after the effective matter potential for neutrino oscillations is reduced by neutrino conversions. The time it takes before the start of the kick is estimated to be

$`\tau `$ $``$ $`{\displaystyle \frac{4\sqrt{2}\pi ^2m_n}{G__F^3\rho }}{\displaystyle \frac{(V_m^{(0)})^3}{(\mathrm{\Delta }m^2)^2\mathrm{sin}^22\theta }}{\displaystyle \frac{1}{\mu ^3}}`$
$``$ $`1\mathrm{s}\left({\displaystyle \frac{10^8}{\mathrm{sin}^22\theta }}\right)\left({\displaystyle \frac{50\mathrm{M}\mathrm{e}\mathrm{V}}{\mu }}\right)^3\left({\displaystyle \frac{10\mathrm{keV}^2}{\mathrm{\Delta }m^2}}\right)^2,`$
where $`\rho `$ and $`\mu `$ are the density and the electron chemical potential inside the neutron star, $`\mathrm{\Delta }m^2`$ is the difference in neutrino masses squared, $`\theta `$ is the singlet-active neutrino mixing angle, and $`V_m^{(0)}0.1\mathrm{eV}`$ is the initial value of the matter potential. The mixing $`\mathrm{sin}^22\theta 10^{11}10^7`$ for a $`110`$keV singlet neutrino is consistent with dark matter. Clearly, the one-second time scale is a natural prediction of this mechanism.

Even if there is no delay, and the pulsar kick occurs at the time of the formation of the jet, the observed pulsar velocities of the order of a few hundred km/s can explain the the $`1`$s time scale. Indeed, it may be necessary to remove the pulsar at least a distance $`D10R_{_{NS}}150`$km away from its original position to disrupt the jet. Thus, the jet may remain intact for the time $`D/v`$, where $`v`$ is the pulsar velocity. The range of time delays that arise this way is determined by the distribution of pulsar velocities. Since most pulsars move faster than 100 km/s, we do not expect a burst that lasts much longer than a second as long as the engine is powered by a neutron star.

Supershort bursts. Gamma-ray bursts with a very short duration $`\tau <0.1`$ second are sometimes considered as a separate class of GRB. It is possible that such short bursts appear because a jet, which was not pointing at Earth initially, changes its direction and quickly passes through the line of sight when the neutron star is kicked. Alternatively, the supershort bursts may correspond to the neutron stars with the highest velocities, which move out of the central region in much less than a second.

Our model can, in principle, relate the numbers of the long and the short gamma-ray bursts to the birth rates of black holes and neutron stars, respectively. Simple estimates of the integrated galactic birth rates for compact stars suggest that neutron stars are born at two or three times the rate of black holes. However, neutron stars and black holes are born from progenitors of different masses. This may affect the angular size of the jet, the brightness, and, ultimately, the detectability of the GRB. Therefore, we cannot relate the numbers of the short and long GRB directly to the birth rates of compact objects.

Our model agrees with non-observation of late optical afterglows for short GRB. If such afterglows will be observed in the future, interpretation of short GRB as interrupted jets will become less convincing.

To summarize, we have proposed that the long GRB originate from supernova-like explosions that result in a formation of a black hole, while the short GRB originate in the same way, except a neutron star forms instead of a black hole. The neutron stars undergo the same dynamical impacts that cause the pulsar kicks. As a result, the momentum and the angular momentum of the neutron star and of the gravitationally bound accretion disk change dramatically. This impact either disrupts the jet altogether or changes its direction away from the line of sight. Central engines that contain neutron stars can, therefore, produce short GRB. At the same time, a central engine with a black hole, which is not a subject to the kick force (at least as long as the pulsar kick is due to neutrino oscillations), would produce a long GRB. At least one plausible model for the pulsar kick mechanism predicts that the kick force would be delayed by a period of time of order one second. However, regardless of the microscopic mechanism, since most pulsars have velocities in excess of 100 km/s, the pulsar kick should remove the neutron star from the central region within about a second, which should destabilize the jet. Our model can be generalized to include any kick mechanism that acts on a neutron star but does not act on a black hole. This can be the explanation for the time separation between the two populations of bursts.

Acknowledgments This work was supported in part by the U.S. Department of Energy and NASA grants.

Competing interests statement The authors declare that they have no competing financial interests.

Correspondence to D.S. (email: semikoz@physics.ucla.edu).