Solar convective velocities: Updated helioseismic constraints

Solar subsurface convection is thought to be an important participant in global-scale dynamics as it determines the Reynolds stresses that play a role in driving differential rotation (e.g., Hanasoge , 2016; Hotta , 2023). The strong stratification of the solar interior together with solar rotation causes the convection to be neither spatially homogeneous nor isotropic. Due to the inherent non-linearity of solar convection, numerical simulations have been essential for our understanding of this topic (see, e.g., Nordlund , 2009 for a review). A long-standing puzzle is that simulations of the convection in the Sun tend to produce much larger amplitude convective flows than what is seen in most observational studies (e.g., Hanasoge , 2012—hereafter HDS2012; Lord , 2014; and O'Mara , 2016). This puzzle is one aspect of what has become known as the “convective conundrum” (O'Mara , 2016).

At small spatial scales (⁠ $≈$ one Mm), solar convection manifests as the well-observed granulation pattern. There is good agreement between observations and numerical simulations of convection at this length scale (e.g., Nordlund , 2009). At intermediate scales (⁠ $≈$ 30 Mm), the dominant feature is termed supergranules. Though these features are well observed near the solar surface, there is not yet general agreement on their physical origin (for a review, see Rieutord and Rincon, 2010).

At scales larger than the supergranulation scale, there is not yet, however, a consensus in the literature on the strength and spatial spectrum of the convective flows in the solar interior. HDS2012 used local helioseismology (see Gizon , 2010 for a review) applied to Helioseismic and Magnetic Imager on the Solar Dynamics Observatory (SDO/HMI) Doppler observations (Schou , 2012) to obtain an upper limit on the root mean square horizontal convective flows in the depth range of 20–30 Mm below the photosphere. This result suggested that the simulations of solar convection of Miesch (2008) (hereafter M2008) overestimate the strength of convective flows in this depth range by roughly two orders of magnitude at large horizontal length scales (see also Gizon and Birch, 2012, hereafter GB2012). Using three-dimensional inversions of helioseismic measurements, Greer (2015) (hereafter GHFT2015), however, found that the observed amplitudes of the flows in this depth range are (order of magnitude) compatible with the simulations. The stark differences between the observational inferences from HDS2012 and GHFT2015 and also between the result of HDS2012 and simulations motivate a study of the current state of understanding of solar subsurface convection.

The overall goal of the current work is to present inferences from the literature of the strength of east-west convective flows in the upper layers of the convection zone. In particular, we aim to provide a consistent presentation of the results from GHFT2015, HDS2012, GB2012, and the measurements of Proxauf (2021) (hereafter P2021). Our specific goal is to produce a consistent plot of the existing results in a way that facilitates comparisons. Consistent comparisons will hopefully elucidate some paths toward reducing any remaining discrepancies. We are not aiming to reproduce the data analysis from the original publications.

GHFT2015 used a 3D inversion (three spatial dimensions) of densely packed ring-diagram measurements of horizontal flow velocity parameters. This approach is different than the one used in the HMI/SDO Ring-Diagram Pipeline, which employs independent one-dimensional (depth only) inversions for each ring tile (for an introduction to the SDO/HMI pipeline, see Bogart , 2011a; 2011b).

HDS2012 used deep-focusing time-distance measurements (Duvall , 1993), together with the assumption that convective flows around the depth of 0.96  $R⊙$ cover a vertical extent of about 1.8 local density scale heights [this corresponds to about 17 Mm at this depth, see Sec. 2.4.2 of Böning (2021) for a detailed description] to obtain an upper limit on east-west convective flows.

Böning (2021) applied the analysis of HDS2012 to synthetic measurements for a simulation of solar convection from Lord (2014). This work concluded that, for the particular simulation of convection that was studied, the method of HDS2012 produced a correct order of magnitude upper limit on the strength of the subsurface convection as a function of horizontal scale. In addition, this work showed that the method employed by HDS2012 could be made more accurate by accounting for a noise-to-signal ratio that depends on horizontal scale.

P2021 applied correlation tracking and ring-diagram analysis to measure the large-scale horizontal flows near the solar surface. These flow maps can be used to obtain constraints on the spectrum of near-surface convection, as described in the next two subsections.

P2021 applied the correlation tracking algorithm of Löptien (2017) to measure the horizontal flows at the surface from the motion of the granulation pattern seen in the SDO/HMI intensity observations. Maps of the east-west component of the flow were constructed by averaging over 27.5 h. Carrington maps covering the latitude range $±59.8°$ were made by assembling strips of $14.4°$ in longitude centered on the central meridian. The quantity $Eℓ$ was computed for each Carrington map. A final average $Eℓ$ was obtained by averaging over the 25% of Carrington rotations with the lowest sunspot numbers.

In contrast with the approach of HDS2012, these spectra were not corrected for the contributions of realization noise. In this sense, they can been seen as upper limits.

P2021 began from the 1D (depth only) ring-diagram inversions supplied by the SDO/HMI ring-diagram pipeline. Each individual ring-diagram measurement covers a time period of 27.2753 h. The construction of the Carrington maps was made using strips of width $15°$ in longitude centered on the central meridian. From these Carrington maps, P2021 computed power spectra of the flows. The averaging of the power spectra over the Carrington rotations was the same as in Sec. II C 1. As for the granulation tracking described earlier, these flow maps were used to estimate $Eℓ$⁠. Again, there was no effort made to remove the contribution of realization noise to the power spectra. Figure 2.8 of P2021 shows the ranges in depth to which the inversions are sensitive. For the case of a target depth near $r=0.98R⊙$⁠, the vertical resolution is roughly 6 Mm.

Figure 1 summarizes the results of the analyses described in the previous sections, together with the unmodified curves for the Anelastic Spherical Harmonic (ASH) and Stagger simulations described in GB2012. This figure is a consolidated baseline for our current understanding of the spectrum of surface and subsurface east-west velocities from observations and simulations.

The granulation tracking and ring-diagram measurements by P2021 provide some important points of comparison. The granulation tracking measurements (these are sensitive to near-surface flows) are qualitatively similar to the result of HDS2012, in the sense that $Eℓ$ increases with angular degree below the supergranulation peak around $ℓ≈120$⁠. As discussed in Sec. II C 1, the granulation tracking measurements include the contribution of noise and are, thus, upper limits. The ring-diagram estimates from P2021 are roughly an order of magnitude (in $Eℓ$⁠) above the upper limit from HDS2012 at the lowest angular degrees. As for the granulation tracking measurements, the ring-diagram estimates do not remove the contribution of noise to the measured $Eℓ$⁠. In this sense, they are also upper limits, and in that generous sense, they are compatible with HDS2012.

The ring-diagram measurements are sensitive to flows at roughly similar depths as the GHFT2015 result but are incompatible in the region around $ℓ≈10$⁠. The source of this discrepancy is not known. The GHFT2015 result depends on two main steps: (1) estimation of the $ux$ and $uy$ ring-diagram parameters using a multi-ridge fitting method (Greer , 2014; Greer, 2019a) and (2) three-dimensional inversion of the measured $ux$ and $uy$ flow velocity parameters (Greer, 2019b). Greer (2015) and Nagashima (2020) both concluded that the multi-ridge fitting method returns $ux$ and $uy$ parameters that are qualitatively compatible with those from traditional single-ridge fitting, as used in the SDO/HMI ring-diagram pipeline providing the input data for the analysis of P2021. This suggests that the significant discrepancy between the results of GHFT2015 and P2021 is due to the use of a three-dimensional (GHFT2015) vs a one-dimensional inversion (P2021).

Figure 1, following GB2012, shows the simulations from M2008. Recently Hotta and Kusano (2021) reported simulations with a much higher spatial resolution than those of M2008. In the highest-resolution simulation, the quantity $Eℓ$ was reduced by about an order of magnitude compared with the result from M2008 (see left panel of Fig. 5 from Hotta , 2023). This change leads to the simulation results being closer to GHFT2015, but the results remain inconsistent with HDS2012 and P2021. It remains to future work to determine if improvements in the numerical simulations (e.g., increased spatial resolution, see Hotta and Kusano, 2021) will bring them into agreement with the observational inferences of HDS2012 and P2021.

The consistent presentation of the convective velocity power spectra shown in Fig. 1 does not substantially resolve any of the main inconsistencies in the literature. To be specific:

• The observational results of HDS2012 and GHFT2015 remain incompatible.

• The upper limits on the flow strengths from HDS2012 and P2021 are both well below those seen in the simulation of M2008.

The last point is highlighted in Fig. 2, which is a simplified version of Fig. 1. This simplified figure shows the upper limit from HDS2012, the granulation tracking measurement of P2021, and the ASH simulation of M2008. Taken together, these curves are a minimum set that highlights the differences between observations and simulations at the large spatial scales.

To summarize the situation with regard to observations: it appears that the result of GHFT2015 is substantially different than the other observational results. As discussed earlier in this section, we speculate that this is due to the 3D inversion used in that work to infer flows beneath the surface. To summarize the situation for simulations: simulations of solar convection produce east-west flows at the largest spatial scales that are stronger than what is observed on the Sun.

The results shown here reveal some important directions for future work. Additional study is needed to determine the causes of the discrepancy in the helioseismic measurements of subsurface flows between P2021 and GHFT2015. The way forward for bringing the simulations into agreement with the observations is not clear. The current situation is that models need relatively high velocities to transport sufficient heat and also to transport enough angular momentum to explain solar differential rotation (Hotta , 2023). We emphasize that the heat and angular momentum do not need to be carried by the same convective spatial scales. For the case of heat transport, one possibility is that scales that are too small to be resolved in the global simulations in reality carry the necessary heat flux. We also note that the very low Prandtl number of the convective fluid in the sun makes a huge difference to its transport properties (e.g., Schumacher and Sreenivasan, 2020; Animasaun , 2022). The resolution of the “convective conundrum” remains an important topic for simulation and theoretical work.

The SDO/HMI observations used here are courtesy of NASA/SDO and the HMI science team. L.G. acknowledges partial support from ERC Synergy under Grant No. WHOLE SUN 810218.

The authors have no conflicts to disclose.

A. C. Birch: Conceptualization (equal); Data curation (equal); Formal analysis (equal); Methodology (equal); Writing – original draft (equal). B. Proxauf: Data curation (equal); Formal analysis (equal); Validation (equal); Writing – original draft (supporting); Writing – review & editing (equal). T. L. Duvall, Jr.: Methodology (equal); Writing – review & editing (equal). L. Gizon: Conceptualization (equal); Methodology (equal); Writing – review & editing (equal). S. Hanasoge: Writing – review & editing (equal). B. W. Hindman: Writing – review & editing (equal). K. R. Sreenivasan: Conceptualization (equal); Writing – review & editing (equal).

The curves shown in Figs. 1 and 2 are available online at the Edmond data repository of the Max Planck Society: https://doi.org/10.17617/3.DFU3SQ.

In this section, we describe how we compute $Eℓ$ from the upper limit on the flow strength obtained in HDS2012. In the course of this description, we elaborate on a few steps that were not clearly described in the original publication. In addition, we describe and correct two minor issues: an error in the area normalization and an error in the estimation of the S/N ratio in the flow maps.

HDS2012 used time-distance helioseismology (Duvall , 1993) to measure east-west travel-time differences (i.e., differences in the travel time for west-going and east-going waves). The output of the deep-focusing time-distance analysis was 17 daily travel-time maps. These 17 travel-time maps were then used to construct a Carrington map covering the latitude range $±58.3°$⁠.

HDS2012 used time-distance helioseismology to compute maps of east-west travel-time differences for consecutive time intervals of 0.1 day covering a total time interval of 27 day. Travel-time maps corresponding to other time intervals T longer than 0.1 day were constructed by averaging in time. This procedure was followed for T running from 0.1 to 1.0 in steps of 0.1 day. The variance $σ2(T)$ of the resulting maps was then computed.

We used the travel-time maps from HDS2012 to carry out the fitting procedure described in the previous paragraph. Figure 3 shows the standard deviation $σ$ in the travel-time maps as a function of the averaging time T. As expected for realization noise, the variance falls essentially as $1/T$⁠. HDS2012 used the small deviation from a $1/T$ dependence to place an upper limit on the contribution of a time-independent signal to the variance. A least squares fit using the model $σ2(T)=S2+N2/T$ yields $N=11.77±0.02$ and $S=2.5±0.3$⁠. This leads to $N/S=4.7±0.8$⁠. This is in contrast with the value of $N/S=2$ used in HDS2012.