SEARCH FOR CHAOS IN NEUTRON STAR SYSTEMS: IS Cyg X-3 A BLACK HOLE?

X-ray binary systems vary on timescales ranging from months to milliseconds (see, e.g., Chen et al. 1997; Paul et al. 1997; Nowak et al. 1999; Cui 1999; Gleissner et al. 2004; Axelsson 2008). Detailed analysis of their temporal variability and fluctuation provides important insights into the geometry and physics of emitting regions and the accretion process. However, the origin of variability is still not clear. It could be due to varying external parameters, like the infalling mass accretion rate. It could also be due to possible instabilities in the inner regions of the accretion disk where the flow is expected to be nonlinear and turbulent. Uttley et al. (2005; see also Timmer et al. 2000; Thiel et al. 2001) argued that the nonlinear behavior of a system can be understood from the log-normal distribution of the fluxes and the rms–flux relation. This implies that the temporal behavior of the system may be driven by underlying stochastic variations. By studying the underlying nonlinear behavior, important constraints can be obtained on these various possibilities.

An elegant way of obtaining the constraint is to perform the nonlinear time series analysis of observed data and to compute the correlation dimension D₂ in a non-subjective manner. This technique has already been used in diverse situations (Grassberger & Procaccia 1983a, 1983b; Schreiber 1999; Aberbandel 1996; Serre et al. 1996; Misra et al. 2004, 2006; Harikrishnan et al. 2006, and references therein). By obtaining D₂ as a function of the embedding dimension M, one can infer the origin of the variability. For example, D₂ ≈ M for all M corresponds to the system having stochastic fluctuation which favors the idea that X-ray variations are driven by variations of some external parameters. On the other hand, a saturated D₂ to a finite (low) value, beyond a certain M, implies a deterministic chaos which argues in favor of inner disk instability. However, to implement the algorithm successfully, the system in question should provide enough data.

The technique was used earlier to understand the nonlinear nature of a black hole system Cyg X-1 (Unno et al. 1990) and an active galactic nucleus Ark 564 (Gliozzi et al. 2002), but due to insufficient data points the analyses were hampered and no concrete conclusions were made about D₂. Later on, another black hole system GRS 1915+105 was analyzed (Misra et al. 2004, 2006; Harikrishnan et al. 2006) which was shown to display low-dimensional chaos in certain temporal classes, while stochastic in other classes.

However, so far none of the neutron star systems have been analyzed in detail in order to understand the origin of nonlinearity. Decades back, Voges et al. (1987) attempted to understand the chaotic nature of Her X-1, but the analysis was hampered by low signal-to-noise ratio (Norris & Matilsky 1989). Since then the investigation of chaotic signature in neutron stars remains unattended. Can a neutron star system not be deterministic? Indeed, several features of X-ray binary systems consisting of a neutron star, such as their magnetosphere and twin kHz quasi-periodic oscillations (QPOs) and their possible relation to the spin frequency of the neutron star, favor the idea that they exhibit nonlinear resonance (e.g., Blaes et al. 2007; Mukhopadhyay 2009). While the QPO itself is a mysterious feature whose origin is still unclear, its possible link to the spin frequency of the neutron star¹ indicates the origin of QPO to be from nonlinear phenomena. Several low-mass X-ray binaries having a neutron star exhibit twin kHz QPOs (Mendez et al. 1998; van der Klis 2006). For some of them, e.g., 4U 1636-53 (Jonker et al. 2002), KS 1731-260 (Smith et al. 1997), 4U 1702-429 (Markwardt et al. 1999), 4U 1728-34 (van Straaten et al. 2002), the spin frequency of the neutron star has been predicted from observed data. However, for the source Sco X-1, which exhibits noticeable time variability (Mendez & van der Klis 2000), while we observe twin kHz QPOs, we do not know the spin frequency yet (but see Mukhopadhyay 2009). For another neutron star Cyg X-2, we do observe kHz QPOs (Wijnands et al. 1998) as well. Several black holes also exhibit QPOs, e.g., GRS 1915+105 (Belloni et al. 2001; McClintock & Remillard 2006) and Cyg X-1 (Angelini et al. 1994).

In the present paper, we first aim at analyzing the time series of two neutron star sources Sco X-1 and Cyg X-2 to understand if a neutron star is a deterministic nonlinear (chaotic) system. Then we try to manifest the knowledge of nonlinear (chaotic/random) property of compact sources to distinguish a black hole from a neutron star. Subsequently, knowing their difference based on the said property, we try to identify the nature of a unknown source (whether it is a black hole or a neutron star). While the nature of some sources, as mentioned above, has already been predicted based on alternate method, for some others, e.g., Cyg X-3, SS433, it has not yet been confirmed.

For the present purpose, we therefore concentrate on three additional sources Cyg X-1, Cyg X-2, and Cyg X-3. While Cyg X-1 has been predicted to be a black hole and Cyg X-2 be a neutron star, the nature of Cyg X-3 is not confirmed yet. Some authors (Ergma & Yungelson 1998; Schmutz et al. 1996; Szostek & Zdziarski 2008) argued for a black hole nature of Cyg X-3, on the basis of its jet, the time variations in the infrared emission lines, the BeppoSAX X-ray spectra, and so on. However, earlier it was argued for a neutron star (Chadwick et al. 1985) by measuring its 1000 GeV γ-rays which suggests a pulsar period of 12.5908 ± 0.0003 ms. By analyzing the time series and computing the correlation dimension D₂, here we aim at pinpointing the nature of Cyg X-3: whether a black hole or a neutron star.

In the next section, we briefly outline the procedure to be followed in understanding the nonlinear nature of a compact object from observed data and to implement it to analyze the neutron star source Sco X-1. In Section 3, we then describe nonlinear behaviors of Cyg X-1, Cyg X-2, and Cyg X-3. Subsequently, in Section 4, we compare all the results and argue for a black hole nature of Cyg X-3. Finally, we summarize in Section 5.

The method to obtain D₂ is already established (see, e.g., Grassberger & Procaccia 1983b; Misra et al. 2004; Harikrishnan et al. 2006). Therefore, here we discuss it briefly. We consider PCA data of the RXTE satellite (see Table 1 for the observations IDs) from the archive for our analysis. We process the data using the FTOOLS software. We extract a few continuous data streams of 2500–3500 s long. The time resolution used to generate light curves is ∼0.1–1 s. This is the range of optimum resolution, at least for the sources we consider, to minimize noise without losing physical information of the sources. A finer time resolution would be Poisson noise dominated and a larger binning might give too few data points to derive physical parameters from it (see Misra et al. 2004, 2006, for details).

Then we calculate the correlation dimension according to the Grassberger & Procaccia (1983a, 1983b) algorithm. From the time series s(t_i) (i = 1, 2, ..., N), we construct an M-dimensional space (called embedding space), in which any vector has the following form:

$$\begin{equation} x(t_{i})=[s(t_{i}),s(t_{i}+\tau),\ldots, s(t_{i}+(M-1)\tau)], \end{equation} \tag{ 1 }$$

where τ is the time delay chosen in such a way that each component of the vector x(t_i) is independent of each other. For a particular choice of embedding dimension M, we compute the correlation function:

$$\begin{equation} C_M(r)=\frac{1}{N(N_c-1)}\sum _{i=1}^{N} \sum _{j=1,j\ne i}^{N_c}\Theta (r-|x_{i}-x_{j}|), \end{equation} \tag{ 2 }$$

which is basically the average number of points within a hypersphere of diameter r, where Θ is a Heaviside step function, N is the total number of points, and N_c is the number of centers. If the system has a strange attractor, then one can show that for a small value of r

$$\begin{equation} {D_{2}(M)}= \frac{d \,{\rm log} C_M(r)}{d\, {\rm log}\, r}. \end{equation} \tag{ 3 }$$

In this numerical calculation, we divide the whole phase space into M cubes of length r around a point and we count the average number of data points in these cubes to calculate C_M(r). The edge effects, which come due to the finite number of data points, have been avoided by calculating C_M(r) in the range r_min < r < r_max, where r_min is the value of r for C_M(r) just greater than 1 and r_max can be found by restricting the M cubes to be within the embedding space. In Figure 1, we show the variation of log(C_M(r)) with log(r) for different values of M for Sco X-1 data.

Zoom In Zoom Out Reset image size

D₂(M) can be calculated from the linear part of the log(C_M(r)) versus log(r) curve and its value depends on the value of M. For a stochastic system, D₂ ≈ M for all M. On the other hand, for a chaotic or deterministic system, initially D₂(M) increases linearly with the increase of M, then it reaches a certain value and saturates. This saturated value of D₂ is taken to be the correlation dimension of the system which is a non-integer. The standard deviation gives the error in D₂.

We first concentrate upon the neutron star source Sco X-1. In Figures 2(a) and (b), we show respectively the light curve and the variation of D₂ as a function of M and find that D₂ saturates to a value 2.6(±0.8). As this is a non-integer, the system might be chaotic. On the other hand, we know that the Lorenz attractor is an example of an ideal chaos with D₂ = 2.05. Therefore, Sco X-1 may be like a Lorenz system. But due to noise its D₂ seems appearing slightly higher (Misra et al. 2004, 2006) than the actual value. However, one should be cautious about the fact that Sco X-1 is a bright source (much brighter than other sources considered later). Hence, the dead time effect on the detector might affect the actual value of saturated D₂ and the computed value might be slightly different than the actual one. However, this cannot rule out the signature of chaos in Sco X-1, particularly because the corresponding count rates are confined in the same order of magnitude and hence the dead time effect, if any, is expected to affect all the count rates in a similar way.

Zoom In Zoom Out Reset image size

However, a saturated D₂ is necessary but not a sufficient evidence for chaos. Existence of color noise (for which the power spectrum P(ν) ∝ ν^−α, where the power spectral indices α = 0, 1, and 2 correspond to "white," "pink," and "red" noise, respectively) into a stochastic system might lead to a saturated D₂ of low value as well (e.g., Osborne & Provenzale 1989; Theiler et al. 1992; Misra et al. 2006; Harikrishnan et al. 2006). Therefore, it is customary to analyze data by alternate approach(s) to distinguish it from a pure noisy time series (Kugiumtzis 1999). One of the techniques is the surrogate data analysis (e.g., Schreiber & Schmitz 1996), which has been described earlier in detail and implemented for a black hole (Misra et al. 2006; Harikrishnan et al. 2006). In brief, surrogate data are random data generated by taking the original signal and reprocessing it so that data have the same/similar Fourier power spectrum and autocorrelation along with the same distribution, mean, and variance as of the original data, but has lost all deterministic characters. Then the same analysis is carried out to the original data and the surrogate data to identify any distinguishable feature(s) between them. The scheme proposed by Schreiber & Schmitz (1996), known as Iterative Amplitude-Adjusted Fourier Transform (IAAFT), is more consistent to generate surrogate data.

Figures 2(c)–(f) compare results for the original data with the surrogate data. It is clear that while distributions and power spectra are same/similar for both the data sets, D₂ is much higher for the surrogate data which suggest existence of low-dimensional chaos in Sco X-1 with D₂ ∼ 2.6. This confirms, for the first time to the best of our knowledge, a neutron star source to display chaotic behavior. As the existence of chaos is a plausible signature of instability in the inner region of accretion flows which is nonlinear and turbulent, as mentioned in Section 1, the corresponding QPO, which is presumably an inner disk phenomenon as well, is expected to be governed by nonlinear resonance mechanisms (e.g., Mukhopadhyay 2009).

In Table 1, we enlist the average counts 〈S〉, its rms variation $\sqrt{\langle S^2\rangle -\langle S\rangle ^2}/\langle S\rangle$, the expected Poisson noise $\langle PN\rangle \equiv \sqrt{\langle S\rangle }$, and the ratio of the expected Poisson noise to the rms value for all sources. It clearly shows a strong correlation between the inferred behavior of the systems and the ratio of the expected Poisson noise to the rms value.

We now look into three additional compact sources: Cyg X-1 (black hole), Cyg X-2 (neutron star), and Cyg X-3 (nature is not confirmed yet), and apply the same analysis as in the case of Sco X-1. Figure 3(b) shows that D₂ for Cyg X-1 seems not to saturate and appears very high. However, there is no surprise in it because its variability is similar to the temporal class χ of the black hole GRS 1915+105 which was shown to be Poisson noise dominated and stochastic in nature (Misra et al. 2004). Indeed, earlier analysis of Cyg X-1 data, while it could not conclusively quantify the underlying chaotic behavior due to insufficient data, revealed very high dimensional chaos. Moreover, a large 〈PN〉 (as well 〈PN〉/rms) for Cyg X-1, compared to that for Sco X-1 given in Table 1, reveals the system to be noise dominated. This ensures Cyg X-1 to be different from Sco X-1. However, the variation of D₂ as a function of M for the original data does not deviate noticeably from that of corresponding surrogate data, as shown in Figure 3(c), which argues that Cyg X-1 is not a chaotic system.

Zoom In Zoom Out Reset image size

Figures 4(b) and (c) show that D₂ for Cyg X-2 saturates to a low value ∼4, which is significantly different than that of corresponding surrogate data. The power spectra and distributions, on the other hand, for original and surrogate data are same/similar (as shown in Sco X-1, not repeated further). The saturated D₂ for Cyg X-2 is almost double than that of a Lorenz system, possibly due to high Poisson noise to rms ratio (see Table 1). This suggests the corresponding system to be a low-dimensional chaos.

Zoom In Zoom Out Reset image size

From Figures 5(b) and (c), we see that for Cyg X-3 the variations of D₂ as a function of M for original and surrogates data are similar to that of Cyg X-1. This confirms that the behavior of the unknown source Cyg X-3 is similar to that of the black hole source Cyg X-1 (see, however, Axelsson et al. 2008). Note from Table 1 that Cyg X-1, X-2, and X-3 are significantly noise dominated compared to Sco X-1. Although noise could not suppress the chaotic signature in the neutron star Cyg X-2, its saturated D₂ is higher than that of Sco X-1. On the other hand, even though the Poisson noise to rms ratio in Cyg X-1 is lower than that in Cyg X-2 (but Poisson noise itself is higher in Cyg X-1), its D₂ never saturates, which confirms the source to be nonchaotic; the apparent stochastic signature is not due to Poisson noise present into the system.

Zoom In Zoom Out Reset image size

Finally, we compare the variations of D₂ for all three cases of Cygnus in Figure 6. Remarkably we find that D₂ values for Cyg X-1 and Cyg X-3 practically overlap, appearing much larger compared to that for Cyg X-2 which is shown to be a low-dimensional chaotic source. On the other hand, Cyg X-2 is a confirmed neutron star and Cyg X-1 a black hole. Therefore, Cyg X-3 may be a black hole.

Zoom In Zoom Out Reset image size

The source Cyg X-3, whose nature is not confirmed yet, seems to be a black hole based on the analysis of its nonlinear behavior. On the other hand, we have shown, for the first time to the best of our knowledge, that neutron star systems could be chaotic in nature. The signature of deterministic chaos, which argues in favor of inner disk instability, into an accreting system has implications in understanding its transport properties particularly in Keplerian accretion disk (Winters et al. 2003). Note that in Keplerian accretion disks transport is necessarily due to turbulence in the absence of significant molecular viscosity. The signature of chaos confirms instability and then plausible turbulence. On the other hand, for a rotating neutron star having a magnetosphere, the signature of chaos suggests their QPOs to be nonlinear resonance phenomena (Mukhopadhyay 2009). The absence of chaos and related/plausible signature of instability in Cyg X-1 and Cyg X-3 suggests the underlying accretion disk to be sub-Keplerian (Narayan & Yi 1995; Chakrabarti 1996) in nature which is dominated significantly by gravitational force.

This work is partly supported by a project (grant no. SR/S2/HEP12/2007) funded by Department of Science and Technology (DST), India. B.B.K. thanks CSIR, India for financial support. Also the financial support to one of the authors (J.D.) has been acknowledged. The authors thank Arnab Rai Choudhuri and Debbijoy Bhattacharya of IISc and the anonymous referee for carefully reading the manuscript, constructive comments, and suggestions.