Solar turbulence from sunspot records

ABSTRACT

It is generally assumed that coupling between the turbulent flow and the magnetic field at the top of the Sun's convection zone leads to a Kolmogorov cascade of kinetic to magnetic energy. An inertial range and a slope value close to −5/3 have been recognized in a log–log diagram of power spectral density versus frequency (or period). However, published values of the slope have large uncertainties and the inertial period range is limited to 0.1 s to 2 yr. We have applied an adapted version of the singular spectrum analysis (SSA) method to the series of (quasi-) daily sunspot numbers ISSN (an indirect way of monitoring solar activity) from 1868 to 2019. The log–log diagram of ISSN variance of SSA components versus frequency displays an inertial slope value of −1.66 ± 0.16 and an inertial range from about 4 to 100 yr. This is consistent with the existence of Kolmogorov turbulent behaviour in the Sun's convection zone.

turbulence, Sun: activity, Sun: evolution, Sun: fundamental parameters, sunspots

1 INTRODUCTION

Observation of sunspots has been an indirect way of monitoring solar activity for a long time. Besides well-known long periodicities (or rather pseudo-periodicities) at 11 yr, 22 yr (a result of the magnetic polarity swap of the 11 yr cycle), and longer cycles such as the ∼90 yr Gleissberg cycle (Le Mouël, Lopes & Courtillot 2017), and, on the shorter period side of the spectrum, 5.5, 3.6, 1 yr (Le Mouël, Lopes & Courtillot 2019a,c,d) and the solar rotation period at 27 d, the Sun is subject to complex, irregular variations that may reflect the non-linear coupling between the turbulent flow and the magnetic field in the solar convection zone (e.g. Browning et al. 2006). A relatively recent analysis of the international, daily and monthly sunspot numbers (SSN) by Plunian, Sarson & Stepanov (2009) has revealed an f^−2/3 scaling of their frequency spectrum, using wavelets for the spectral decomposition. These authors interpret that scaling as a buoyancy spectrum that they believe leads to a Kolmogorov magnetic energy spectrum of k^−5/3, k being the spatial wavenumber (Morfill et al. 1991; Lawrence, Cadavid & Ruzmaikin 1995; Plunian et al. 2009). Plunian et al. (2009) use wavelet analysis to determine the power spectral density (PSD) of ISSN.

Singular spectral analysis (SSA) has been applied with interesting results to the irregular El Niño–Southern Oscillation (ENSO) phenomenon, to geopotential height data, to a number of indicators of climate variability, to solar observations, and to cosmogenic isotopes. SSA provides at the same time a noise reduction technique, a detrending algorithm, and a way to identify oscillatory components. In a series of recent papers, we have used an adapted version of SSA to determine the trends and components of solar proxies (ISSN, sunspot group numbers, sunspot areas, number of polar faculae PF) and also a number of indices affected by solar variability, such as climate indices (Le Mouël et al. 2019a), length of day (Le Mouël et al. 2019b), the geomagnetic indices aa and Dst (Le Mouël et al. 2019c), and global surface temperatures (Le Mouël et al. 2019d).

In this paper, we try to go beyond the analysis of Plunian et al. (2009). We analyse ISSN using our adapted version of SSA, to see whether the value of the slope of the density spectrum (DS) can be confirmed and, if so, refined, and whether its spectral domain of validity can be extended.

2 AN SSA SPECTRUM OF ISSN

The series of (quasi-) daily sunspot number ISSN (one data point every 26.3 h) that we use starts on 1868 January 1 and ends on 2019 May 31 and is downloaded from http://www.sidc.be/silso/datafiles. We have started to analyse it using SSA in several earlier papers (Le Mouël et al. 2019a,b), and most recently in Le Mouël et al. (2019c,d). A short description of the method, as used in our previous papers, is given by Lopes, Le Mouël & Gibert (2017). We recall that SSA extracts from the data eigenvalues and eigenvectors that can be grouped to display a trend and components that are oscillatory but may be non-periodic. Eigenvectors of ISSN of particular interest are displayed in Le Mouël et al. (2019c) and the first eigenvalues/components are shown in this paper as Fig. 1 and listed in Table 1.

Figure 1.

The first 65 SSA eigenvalues of ISSN as a function of rank.

Table 1.

List of the 65 first eigenvalues of the SSA analysis of ISSN ranked as a function of decreasing period (see Figs 1 and 2).

Let us replace for brevity ISSN by N(t), available as we have seen from 1868 January 1 to 2019 May 31. The eigenvalues n₁,n₂,…,n_L shown in Fig. 1 are ranked as a function of their magnitude. To each n_j is associated an eigenfunction (of time),

$$\begin{eqnarray*} {n_j}\left( t \right) = {{\it{n}}_{{j}}}\exp{\rm{ }}(2\mathrm{ i}\pi t/{T_j}), \end{eqnarray*}$$

(1)

where n_j is the amplitude, and T_j the period of the oscillation/component of rank j (see remarks below). We have

$$\begin{eqnarray*} N\left( t \right) = {\sum _j}{{\it{n}}_{{j}}} \exp{\rm{ }}(2\mathrm{ i}\pi t/{T_j}). \end{eqnarray*}$$

(2)

We display the spectral density (SD) obtained from the n_j as a function of f_j = 1/T_j in log–log coordinates in Fig. 2 (see also Table 1).

Figure 2.

Same as Fig. 1 but with the log of the SD as a function of the log of their period. Inset: histogram of values of the slope of the central inertial segment, calculated using a bootstrap technique, with one million estimates for samples with 50 per cent of the data included. The slope is found to be −1.66 with a standard deviation of 0.50 and a standard error on the mean of 0.08.

Note: PSD is defined as the square of the modulus of a signal's Fourier transform. As is the case for most matrix decompositions, SSA uses singular value decomposition (SVD; Golub & Reinsch 1971), ensuring orthonormality on the basis of eigenvectors and unicity of the decomposition. A parallel can be made between the two methods: the only thing that differs is the nature of the matrix that allows to pass from data to ‘parameters’ (Vandermonde for Fourier – with an imposed base − and Hankel for SSA – where we construct the base). In both cases the bases are orthogonal. There is a link between Fourier related PSD and what we call SD in the case of SSA. Thus we use a similar terminology, as is done by most authors (e.g. Plunian et al. 2009).

The SSA spectrum of Fig. 2 displays two successive period domains: a linear portion from ∼100 to ∼4 yr, and a noisy range with (possibly) steeper slope for periods shorter than a few years. There is not enough data on the longer period side to indicate another spectral domain for periods longer than 100 yr (the trend, a 105-yr and the 89-yr Gleissberg cycle, corresponding to the first component in Fig. 1).

We determine the slope of the linear segment in Fig. 2 using a bootstrap technique (Efron 1981; Efron & Tibshirani 1994), involving one million estimates for samples with 50 per cent of the data included (Fig. 2, inset). The slope is found to be −1.66 with a standard deviation of 0.50 and a standard error on the mean of 0.08. −1.66 ± 0.16 (with a confidence of 95 per cent) is the Kolmogorov value of −5/3 to within 10 per cent.

By analogy with the vocabulary used in turbulence theory, we see that the sunspot number, therefore solar activity, displays an inertial zone and a dissipative zone typical of turbulence. The large period domain (∼4–100 yr) and the accuracy of the −5/3 exponent (to 10 per cent) are the central results of this paper (see discussion below).

We note that three points in the period domain of the inertial range lie well above the linear spectral segment: they correspond to quasi-periodic components of solar activity at 11, 5.5, and 4.8 yr.

3 THE KOLMOGOROV EXPONENT

The –5/3 exponent was obtained by Kolmogorov (1941) in the frame of an analysis of (homogeneous and isotropic) turbulence theory. According to Chorin (1994), a major goal was to provide an intuitive explanation of what happens in a turbulent flow (see also Yaglom 2001). The basic starting physical idea is that the small-scale turbulent flow is produced by a cascade process where energy is transferred from the mean flow and the largest structures of the turbulent flow to smaller and smaller structures. Let the turbulent flow be expressed as

$$\begin{eqnarray*} {\it{u}}\left( {\boldsymbol{x}} \right) = {\sum _k}{{\it{u}}_k}\,\exp\left( {\mathrm{ i}{\boldsymbol{k}}\cdot {\boldsymbol{x}}} \right), \end{eqnarray*}$$

(3)

where exp(i|${\boldsymbol{k}} \cdot {\boldsymbol{x}}$|⁠) can be associated with an ‘eddy’ of size k⁻¹, where k is the wavenumber (Chorin 1994, p. 51; Frisch 1995). Let < u_k² > be the energy spectrum, k₁ be a characteristic value of the domain where energy is given to the (large scale) flow, and k₂ be a value of k that characterizes the domain where k²E(k) becomes very large (dissipation range). Between k₁ and k₂ (the inertial domain), energy passes from the region around k₁ to the region around k₂. What could E(k) in this intermediate region depend on? Kolmogorov makes the hypothesis that E(k) depends on only two parameters: k and the rate of dissipation of energy Ɛ = (d/dt) < u² >. k has the dimension of L⁻¹ (L being a length) and Ɛ of L²T⁻³ (T being a duration). Kolmogorov shows that the only dimensionally correct combination is

$$\begin{eqnarray*} E\left( k \right) = C{\varepsilon ^{2/3}}{\rm{ }}{k^{ - 5/3}}, \end{eqnarray*}$$

(4)

where C is a dimensionless constant. This is where the (−5/3) exponent of the wavenumber k comes from.

4 THE CASE OF THE SUNSPOT (ISSN) SERIES

The ISSN index is a priori a pure number and cannot be given a dimension. It could however be given a quasi-dimension, for instance that of an energy. This is rather natural when one takes into consideration the origin of sunspots, see e.g. Plunian et al. (2009), who put forward the coupling between the turbulent flow and the magnetic field at the top of the Sun's convection zone (Browning et al. 2006): kinetic energy is exchanged with magnetic energy. One can observe, comparing equations (1) and (3), that their dimensions are not the same: the wavenumber k in equation (3) has the dimension of an inverse length, whereas the T in equation (1) has the dimension of a time. Recall that the −5/3 law comes from a dimensional argument. In general, without delving too deeply in the question, one calls for the Taylor hypothesis (e.g. Frisch 1995): the time variation of u (in equation 3) at a fixed point is reinterpreted as a space variation of the small-scale flow v' (measured in the frame of reference of the mean flow). The correspondence between space increments l for v' and time increments t for u is then simply l = Ut, where U is the velocity of the mean flow. In the many experiments in wind tunnels that provided much of the information on developed turbulence, results obtained in the time domain were easily recast in the space domain (and vice versa) via the Taylor hypothesis (Frisch 1995, p. 59; Plunian et al. 2009). In fact, the time (or frequency) is often relabelled as a position (wavenumber) without mentioning the Taylor hypothesis. The key element one looks for in the analysis of a potentially turbulent flow/time series is simply the −5/3 exponent of log(SD) versus frequency.

5 THE KOLMOGOROV EXPONENT AND ISSN PRIOR TO THIS STUDY

Of course, the importance of a power estimate from an observational series depends on closeness of a log–log slope value to −5/3 and to the range of the straight line segment (inertial zone) with that slope that is actually observed. We have indeed shown in Section 2 that the SSA analysis of ISSN leads to a density spectrum of eigenvalues that follows the inertial Kolmogorov slope of −5/3 in a period domain extending from roughly 4 to 100 yr. A number of authors have undertaken this kind of work in different period ranges and with different uncertainties. We now summarize these results.

Coleman (1968) is among the first to have found in the power spectra of the magnetic field of the solar wind evidence of Kolmogorov turbulence at 1 au in the period range from ∼1 month to 27 d, the solar rotation period. An incomplete list of references includes Burlaga & Goldstein (1984), who studied the evolution of plasma flow at 4–5 au and found intermittency in the 5/3 law; they covered a range of 6 months. Marsh & Tu (1990) analysed magnetic field and plasma data from the Helios 1 and Helios 2 spacecrafts between 0.3 and 1 au at the activity minimum of solar cycle 21; they found a −5/3 spectral law, with turbulence better developed in low-speed flows. Espagnet et al. (1993) computed the power spectra of solar granulation; their energy decreases as a Kolmogorov −5/3 power law in the granulation range (1000–2000 km). Lawrence et al. (1995) studied ISSN and the 10.7 cm microwave flux; using Fourier and wavelet spectral analysis, they uncovered two regimes of magnetic variability. Cascade (−5/3) behaviour was observed from 2 yr down to 2 d, whereas low dimensional chaotic behaviour was observed at time-scales longer than 8 yr. Abramenko (2005) worked on images of the photosphere, and extracted the magnetic field of the Sun at times of low flare activity; the study extends over only 60 h. Alexandrova et al. (2008) computed the power density spectrum using wavelet transforms over a range from 0.1 s to 3 months. McAteer, Gallagher & Conlon (2010) studied turbulence and complexity in order to predict solar flares; using wavelets, these authors found that the PSD of flares is indeed characterized by a −5/3 law. Ashwanden (2019) recently reviewed models of the solar wind. He noted that, by introducing the reflection of Alfvénic waves, Cranmer & van Ballegooijen (2005) produced ‘a Kolmogorov-like spectrum that does not change dramatically from the photosphere to the solar wind’.

The analysis by Plunian et al. (2009) is the closest to ours and among the most recent ones. As recalled in the Introduction, these authors analysed several kinds of time series of sunspot numbers, concentrations of ¹⁴C in tree rings and ¹⁰Be in polar ice and found two frequency scalings, f^−2/3 when the Sun was active and f⁻¹ when it was quiescent. These values of exponents are not that expected for Kolmogorov turbulence that should be exactly (within uncertainty) −5/3. Plunian et al. (2009) suggest that these exponents reflect buoyancy spectra implying values of k^−5/3 and k⁻² for the magnetic energy spectra of the Sun. Plunian et al's (2009) fig. 3, which displays the PSD of daily and monthly SSN and 10 yr averages of ¹⁴C, ranges from 0.1 to over 10⁴ months, shows significant dispersion of the curves, slope changes, and energy peaks at 1 month (27 d) and 10 yr (Schwabe). The −2/3 slope is based on a small range of frequencies (2–20 months). A value of −5/3 (between 1 and 104 months) is not such a bad fit given these uncertainties. But our own results with SSA are much more compelling.

6 CONCLUSION

The diverse proxies of solar activity, and in particular solar wind velocity, that have been studied to conclude that Kolmogorov turbulence is present in their power density spectra lead to slopes close to −5/3, but with large uncertainties (well in excess of our 10 per cent standard error on the value −1.66) and over a period range that is limited to the range from 0.1 s to 2 yr (values from the references in Section 5). The improvement of accuracy (10 per cent) and the extended (longest) period over which we identify the inertial domain (∼4–100 yr) in the sunspot time series ISSN have been made possible by use of our new SSA analysis. It strengthens the hypothesis that a Kolmogorov cascade is present and the evidence of inertial behaviour is very significantly extended towards longer, multidecadal to secular periods.

ACKNOWLEDGEMENTS

The ISSN data used in the paper are accessible at the following site: http://www.sidc.be/silso/datafiles. IPGP contribution 4095.