Phase coherence among the Fourier modes and non-Gaussian characteristics in the Alfvén chaos system

Phase coherence among the Fourier modes and non-Gaussian characteristics in the Alfvén chaos system Yasuhiro Nariyuki1,∗, Makoto Sasaki, Naohiro Kasuya, Tohru Hada, and Masatoshi Yagi 1Faculty of Human Development, University of Toyama, 3190, Gofuku, Toyama City, Toyama 930-8555, Japan 2Research Institute for Applied Mechanics, Kyushu University, 6-1, Kasugakoen, Kasuga City, Fukuoka 816-8580, Japan 3Faculty of Engineering Sciences, Kyushu University, 6-1, Kasugakoen, Kasuga City, Fukuoka 816-8580, Japan 4National Institutes for Quantum and Radiological Science and Technology, 2-166, Omotedate, Obuchi, Rokkasho, Aomori 039-3212, Japan ∗E-mail: nariyuki@edu.u-toyama.ac.jp


Introduction
It is well known that finite amplitude Alfvénic fluctuations are frequently observed in the solar wind plasma (Refs. [1][2][3][4]). Since these fluctuations disappear with increasing heliocentric distances (Refs. [5,6]), their damping processes possibly play important roles in the generation of magnetic turbulence and/or in plasma heating. Damping of such finite amplitude Alfvén waves has been intensively discussed in terms of parametric decay (Refs. [7,8]) and modulational instabilities (Refs. [9,10]). When the system concerned contains the source of Alfvén waves and dissipation terms (sinks of energy), the statistical stationary Alfvénic states can be achieved approximately (Refs. [11,12]). Such a situation can be realized in space plasma, e.g., in the earth-foreshock region (Refs. [13,14]) and in the solar atmosphere (Ref. [15]), in which some sources of Alfvén waves may exist. It should also be noted that the magnetic fluctuations observed in the earth's foreshock region exhibit finite phase coherence among the Fourier modes (Refs. [16,17]). The result is in contrast to the often used ansatz of the "phase random approximation", which is assumed to be held in weak turbulence theories.
In the present study, we discuss the generation of non-Gaussian characteristics in nonlinear evolution of Alfvén waves by using the Alfvén chaos system proposed by Hada et al. (Ref. [11]). The PTEP 2017, 033J01 Y. Nariyuki et al.
Alfvén chaos system is derived from the driven-damped derivative nonlinear Schrodinger equation (DNLS) (Ref. [18]) with traveling wave solution. In Sect. 2, the phase coherence among the Fourier modes in the chaotic time series of the Alfvén chaos system is discussed. In Sect. 3, we discuss a linear Markovian equation by using the Mori projection operator method (Refs. [19,20]). The results are summarized and future issues are discussed in Sect. 4.

Phase coherence of Alfvén chaos
Let us first define the Alfvén chaos system (Refs. [11,21,22]). By incorporating the monochromatic growth term and the dissipation term into the DNLS with traveling wave ansatz, Hada et al. (Ref. [11]) obtained the ordinary differential equation seṫ The ordinary differential equation set (1), (2) is derived by using the coordinate transformation ξ = X − VT , where X and T are spatial coordinate and time respectively, and V is a constant. If there is no driving force, Eqs. (1), (2) describe the stationary solutions of the DNLS system (Ref. [11]). In accordance with the past study in Ref. [21], here we use the variable t = ξ as the conventional time variable. The terms b y and b z are the transverse magnetic field, ν corresponds to the coefficient of the Burgers-type dissipation term, λ indicates the velocity V normalized to the ambient magnetic field, and a and are the intensity and frequency of the driving term, respectively. The ambient magnetic field is assumed to be in the x-y plane.
Past studies (Refs. [11,21]) carried out direct numerical simulations of the system (1)-(3) and discussed the bifurcation diagram of the Alfvén chaos. In the present study, we also carry out direct numerical simulation of the Alfvén chaos system (1)-(3) with the parameters that lead to the chaotic time series observed by Chian et al. (Ref. [21]). The runs are presented in Table 1. In a similar way to Chian et al. (Ref. [21]), we vary only the dissipative parameter ν and fix the other parameters as λ = 0.25, = −1, and a = 0.3. The hodogram (b y − b z ) is shown in Fig. 1. The reason for the gap between runs 1-6 and runs 7-10 is a window in the bifurcation diagram (Ref. [21]), in which period-doubling bifurcation of the stable periodic orbit with decreasing ν is observed. At the end of the bifurcation (ν = 0.06212), an interior crisis is observed, which occurs due to the collision between the unstable orbit from the saddle-node bifurcation and the bounded chaotic attractor (Ref. [21]). Detailed discussion on the bifurcation was given in the past study in Ref. [21]. As for numerical integration, a fourth-order Runge-Kutta scheme is used for the time integration, with t = 5 × 10 −4 .
To check the statistical stability, runs are repeated fifty times using different initial conditions each time, i.e., by giving b y (t = 0) = b z (t = 0) = c 0 , where c 0 is a uniform random number within the Table 1. The parameter (ν) used in simulation runs. The other parameters are fixed as λ = 0.25, = −1, and a = 0.3 so that the system shows the chaotic behavior (Ref. [21]  To examine the non-Gaussian characteristics in Alfvén chaos, we use the phase coherence index (Refs. [16,26]) as a measure for evaluating the phase coherence among the Fourier modes contained in a given time series (Ref. [23]). Here we use the modified version of the phase coherence index (Ref. [26]) where L i (α, τ ) = t |f (t + τ ) − f (t)| α is the structure function, f = b y , b z , i = or (original data), rs (randomized surrogate data), and cs (coherent surrogate data). The value α = 1 is used in the present study. Figure 3 shows the time series of (a) the original data of b y , (b) the randomized surrogate   Table 1 data, and (c) the coherent surrogate data in run 3. As seen in Fig. 3(c), the wave form of the coherent surrogate (randomized surrogate) becomes localized and the variance measured by the structure function becomes smaller (larger) than the one for the original data. Here we use one time series for each run, since the dependence of the statistical properties of the runs on the initial conditions is almost negligible, as mentioned above.
In the present study, the phase shuffle surrogate data (Refs. [16,17,23,25,26]), in which the correlation function (power spectrum) of the original data is conserved, is used as the randomized surrogate data. This can be done by first discrete Fourier transforming (DFT) the target time series to obtain the complex DFT spectra, and then by rewriting its phase part by random numbers within the range [0, 2π ]. To evaluate the significance, multiple surrogate data are generated by using different sets of random numbers. In the present study, the number of random surrogate data for each window (one DFT window) is 100. The number of windows is 5, so we have 5 windows to check the variance of C φ . As shown in Fig. 5, the error bar is negligibly small. The coherent surrogate data are obtained by substituting a uniform number for the original phase.
When we apply the phase shuffle surrogate to the Monte Carlo significance testing, the null hypothesis is the linear Gaussian process (color noise) (Refs. [23,25]). While the histogram of the original data of b y ( Fig. 4(a)) is far from the Gaussian distribution, one of the surrogate data sets ( Fig. 4(b)) becomes close to the Gaussian distribution (Ref. [25]). Here we use 16 384 data points (data from t = 100.1 to ∼ 1738.4 with interval of 0.1) for the DFT with box windows similar to the past studies in Refs. [17,27]. The histograms are also calculated by using the same time series (16 834 data points). As mentioned above, the number of random surrogate data for each window is 100 and the number of windows is 5. The level of significance for the one-sided rank-order test is defined as (Ref.    surrogate data. In the present study, the significance of C φ is defined by the structure function L i . If L rs > L or for all the random surrogate data, the level of the significance is larger than about 99 percent (= (1 − β) × 100). If L rs < L or is observed, we treat C φ (τ ) as insignificant. The time difference between neighboring windows is 5. The validity of the present calculation is checked by using the relative error of the energy between the original data ( b 2 y + b 2 z ) and the surrogate data. In runs except for run 7, the maximum energy difference in the surrogate data is 6.6 percent of that of the original data, except for run 7, where a value of 16.6 percent is obtained. The effect of the windows is also discussed at the end of this section.
Figures 5(a)-(f) shows the dependence of C φ on τ in runs 1-6. Error bars are added to the plot but are not visible due to small standard deviation. We find that C φ decreases monotonically with increasing τ , and also with increasing ν except for runs 3 and 4. In run 4 (ν = 0.06212), an interior crisis occurs (Ref. [21]). The weakly chaotic time series is observed at the parameter (Ref. [21]). On the other hand, the chaotic time series is observed in run 3 (ν = 0.0616) (Ref. [21]).
Except for run 1 (Fig. 5(a)), C φ (τ ) at the maximum τ values is insignificant, i.e., L rs > L or . For instance, in Fig. 5(e), C φ at τ = 0.1, 0.5 is significant but C φ at τ = 1.0 is insignificant. In runs 7-10, C φ is insignificant for all τ . The dependence on ν and τ indicates that the non-Gaussian characteristics measured by C φ are made by the fluctuations within the short period that is damped with increasing ν. Figure 6 shows the amplitude of the DFT spectra in (a) run 3 and (b) run 8. Although the side-band mode with a higher frequency than the driven mode (ω = 1) can clearly be observed in run 3, the spectrum in run 8 does not show the side-band mode at higher frequencies. It is noteworthy that as shown in Fig. 7, the time series in the runs with insignificant C φ is also intermittent and the histogram 6   is not Gaussian. This clearly indicates that there are some non-Gaussian characteristics that are not evaluated by the procedure using the structure function. In other words, the type II error (Ref. [24]), which is the error that the null hypothesis is adopted although not valid, occurs in the present Monte Carlo significant testing. Improvement of the present method is a future subject. Instead, here we discuss a different aspect of the non-Gaussian characteristics in the next section. At the end of this section, we comment on the effect of the window used in the DFT. A simple box window, in which a rectangular shape is multiplied by the time series without a change of the energy, is used in the present study. As demonstrated by Suzuki et al. (Ref. [37]), when the number of data points is sufficiently large (larger than about 2000), the false rejection rate of the box window becomes smaller than those of the data windows such as the Welch window. On the other hand, when we apply the cosine cube-tapered rectangle function used in the past study in Ref. [26], a significant finite phase coherence appears even in runs 7-10. However, in our calculation using the cosine cube-tapered rectangle function, the normalized maximum energy difference of the surrogate data is more than 10 percent of that of the original data in all the runs except for run 3 (9.6 percent). When we limit our discussion to runs 7-10, the difference is more than 14 percent. We chose the box window for the better performance in energy conservation. Since such a characteristic can depend on other factors such as the choice of the structure function and the original time series, further experiments are needed to evaluate detailed characteristics of the phase coherence index.

Evaluation of the chaos-induced friction and the model of time correlation
In the previous section, we find that the runs with large ν have insignificant phase coherence (i.e., L rs > L or ). Although statistical testing provides us with an evaluation criteria to regard the fluctuation as weak turbulence, the consequences do not ensure that the fluctuation is produced by a linear Gauss process. This point is mentioned as a type II error in the previous section. In this section, we discuss the physics-based modeling of the stochastic process in the Alfvén chaos system by using the Mori projection operator method (Refs. [19,20,[28][29][30][31][32]), which enables us to decompose the macro behavior from the chaotic system. The effects of the random motion appear as the difference between the nonlinear terms and the projected term.

P(H (A(t))) = H (A(t))A
where · is the long-time average, and the symbol † denotes the Hermitian conjugate. In the same manner as the past studies in Refs. [20,28,30], where is the evolution operator (Refs. [20,28,30]), Q = 1 − P, and is the fluctuating force, and (t) = r(t)r † bb † −1 is the memory function that is assumed to be a diagonal matrix in the present study. Equation (5) can be rewritten as From Eq. (11), we can obtain an evolution equation for the time correlation function as (Refs. [20,28,30]) C zy (t) − νĊ yy (t) = L tzy C yy + L tzz C zy + C Fzy , where C F = F(t)b † , L t = L+L P . By performing the Laplace transformation, the relation between the time correlation spectrum and the memory spectrum can be obtained (Refs. [29,30]). Here we consider the linear Markovian stochastic equation as (Refs. [20,28,29]) γ zz ) is the chaos-induced friction coefficient (Refs. [20,29]). From Eqs. (12) and (13), we obtain where κ = (γ yy − ν(L tyz + L tzy ) − (L tyy + L tzz ))/(1 − ν 2 ), 2 0 = (|L t | − γ yy L tzz )/(1 − ν 2 ), and Direct numerical simulations indicate that F tyy becomes the sinusoidal function  Numerically evaluated coefficients are presented in Table 2. Thus, Eq. (14) models a simple forcedamped oscillator with a special solution C yy,s = I 1 cos(t + ) + I 2 sin(t + ), where or where From Eq. (21), we can have the quadratic equation for γ yy , where Past studies evaluated γ by using mode-coupling theory (Ref. [20]) and the numerical results of the nonlinear force (Ref. [29]). Here we evaluate the unique γ yy using the multiple root of Eq. (23). From the condition that the discriminant of Eq. (23)) is zero, we obtain the equation In Table 2, γ yy and the coefficients of Eq. (14) for each run are tabulated. Figure 8 shows the dependence of C yy on τ . Due to the existence of the side-band mode as shown in Fig. 6, the modulation is observed in C yy calculated by direct numerical simulation. We remark that the chaos-induced friction γ yy does not exist in the original system (1)-(3) and appears as a consequence of the Markov approximation. The chaos-induced friction is derived by using the Markovian approximation, which corresponds to the null hypothesis of the Monte Carlo testing in the previous section. As shown in Fig. 9, the phases of C yy,s are decorrelated to C yy with increasing ν. It is noteworthy that even if Eq. (14) assumes a linear Markovian process, C yy,s in runs 7-10, which has no significant phase coherence, does not agree with direct numerical simulation. This suggests that the surrogate data method can overlook non-Gaussian characteristics, as shown in the previous section. In contrast, even if finite phase coherence is observed, C yy,s in the runs with smaller ν agrees well with direct numerical simulation. In contrast to the phase, the amplitude of C yy,s is not consistent with that of direct numerical simulation. We also note that the cases with κ < 0 are not realistic, since their general solutions diverge in time.

Summary and Discussion
In the present paper, we discuss the non-Gaussian characteristics of the time series data produced by the reduced nonlinearAlfvén system (1)-(3) by using the phase coherence index (Sect. 2) and the Mori projection operator method (Sect. 3). Through Monte Carlo significance testing, phase coherence among Fourier modes appears when the dissipative coefficient ν is small, while it is not significant when the dissipative coefficient is large. On the other hand, the probability distributions (histogram) of the time series are not similar to a Gaussian distribution even if there is no significant phase coherence. This indicates that evaluation using the phase coherence index can overlook evidence of non-Gaussian characteristics, which are also suggested through the evaluation of chaos-induced friction. As shown in Fig. 5, C φ decays within a shorter time scale than the wave period of the drive mode (= 2π ). This is in contrast to the past observational study in Ref. [17], which suggested that the peak of C φ corresponds to the frequency of the waves excited by the ion beam instability. To resolve the contradiction, more detailed observational study of the turbulence in the earth-foreshock region is necessary.
To discuss non-Gaussian characteristics from the aspect of a physics-based stochastic model, we apply the Mori projection operator method to the Alfvén chaos system to derive the force-damped oscillator model, and by using this, we determine the chaos-induced friction coefficient. When the dissipation is relatively small, the model equation of the time correlation function (the special solution of the force-damped oscillator) indicates that the phase evolution of the model equation agrees well with the phase evolution computed by direct numerical simulations. This means that the derived Markovian equation is not consistent with consequences of the present Monte Carlo significance testing (linear Gauss processes), although the Markovian approximation applied to the model agrees with the null hypothesis of the Monte Carlo testing. Notice that even if the histogram is clearly non-Gaussian as shown in Fig. 7, it is not trivial whether the processes concerned can be regarded PTEP 2017, 033J01 Y. Nariyuki et al. as a linear Gaussian process. The phase coherence index and the chaos-induced friction react to the different aspects of non-Gaussian characteristics. On the other hand, the theoretical model of the anisotropic spectra of magnetohydrodynamic (MHD) turbulence (Ref. [33]) is based on the closure derived by using a small deviation from the Gaussian distribution (quasi-normal approximation).  Practically, the applicability condition of the Gaussian distribution is significant in understanding MHD turbulence in the solar wind. Such an applicability condition should be justified from a plurality of viewpoints. As a result, the present result suggests that even if the phase coherence is significantly low, the physics-based model derived by using the Mori-projection operator method declines to apply the Gaussian distribution to the Alfvén chaos system. When we discuss turbulence observed by the spacecraft, a lesson from the present analysis should be absorbed.
Some recent studies also discussed random (thermal) fluctuations in the solar wind plasma by using the fluctuation-dissipation theorem (Refs. [34][35][36]), while the physical origin of the fluctuations is unclear at the present time, since the collisional effects of the solar wind plasma are usually believed to be negligible. As discussed in the present study, even if the time series is non-Gaussian, the Markov approximation can partially come into effect in one aspect. From this point of view, the randomization due to non-Gaussian processes can become one origin of the random fluctuation in the solar wind plasma.