Sensitivity of Helioseismic Measurements of Normal-mode Coupling to Flows and Sound-speed Perturbations

TL;DR

This paper derives the sensitivity kernels for helioseismic normal-mode coupling to flows and sound-speed perturbations, enabling more accurate inferences of solar interior dynamics using first-Born perturbation theory.

Contribution

It provides a detailed derivation and correction of the sensitivity kernels for mode coupling, including sound-speed effects, based on the formalism of Lavely and Ritzwoller (1992).

Findings

Derived and computed mode coupling sensitivity kernels to flows.

Reformulated the fluid-continuity condition for accuracy.

Established kernels for sound-speed perturbations.

Abstract

In this article, we derive and compute the sensitivity of measurements of coupling between normal modes of oscillation in the Sun to underlying flows. The theory is based on first-Born perturbation theory, and the analysis is carried out using the formalism described by \citet{lavely92}. Albeit tedious, we detail the derivation and compute the sensitivity of specific pairs of coupled normal modes to anomalies in the interior. Indeed, these kernels are critical for the accurate inference of convective flow amplitudes and large-scale circulations in the solar interior. We resolve some inconsistencies in the derivation of \citet{lavely92} and reformulate the fluid-continuity condition. We also derive and compute sound-speed kernels, paving the way for inverting for thermal anomalies alongside flows.

Equations

$$-\rho {\omega }_k^2{\mathbf{\xi }}_k-\nabla (\rho c^2\nabla \cdot {\mathbf{\xi }}_k-\rho g{\mathbf{\xi }}_k\cdot \hat{\mathbf{r}})-g\hat{\mathbf{r}}\nabla \cdot (\rho {\mathbf{\xi }}_k)=0,$$

$$(-\rho {\omega }_k^2+{\mathcal{L}}_0){\mathbf{\xi }}_k=0.$$

$${\int }_{\odot }d\mathbf{r}\rho {\mathbf{\xi }}_k^{\ast }\cdot {\mathbf{\xi }}_{k^{\prime }}={\delta }_{k{k}^{\prime }},$$

$$\text{mn@boldsymbol}\xi (\mathbf{r},\omega )=\sum _q{a}_q(\omega ){\mathbf{\xi }}_q(\mathbf{r})=\sum _q{a}_q^{\omega }{\mathbf{\xi }}_q(\mathbf{r}),$$

$$(-\rho {\omega }^2+{\mathcal{L}}_0)\text{mn@boldsymbol}\xi =\mathbf{J},$$

$$(-\rho {\omega }^2+{\mathcal{L}}_0+\delta \mathcal{L})(\text{mn@boldsymbol}\xi +\delta \text{mn@boldsymbol}\xi )=\mathbf{J},$$

$$(-\rho {\omega }^2+{\mathcal{L}}_0)\delta \text{mn@boldsymbol}\xi =-\int d{\omega }^{\prime }\delta {\mathcal{L}}_{\omega -{\omega }^{\prime }}\text{mn@boldsymbol}\xi (\mathbf{r},{\omega }^{\prime }),$$

$$\delta \text{mn@boldsymbol}\xi =\sum _j\delta a_j^{\omega }{\mathbf{\xi }}_j.$$

$$(-\rho {\omega }^2+{\mathcal{L}}_0)\delta \text{mn@boldsymbol}\xi =\rho \sum _j({\omega }_j^2-{\omega }^2)\delta a_j^{\omega }{\mathbf{\xi }}_j.$$

$$\sum _j\delta a_j^{\omega }({\bar{\omega }}_j^2-{\omega }^2){\int }_{\odot }d\mathbf{r}\rho {\mathbf{\xi }}_k^{\ast }\cdot {\mathbf{\xi }}_j=-\sum _{k^{\prime }}{\int }_{\odot }d\mathbf{r}\int d{\omega }^{\prime }{a}_{k^{\prime }}^{{\omega }^{\prime }}{\mathbf{\xi }}_k^{\ast }\cdot \delta {\mathcal{L}}_{\omega -{\omega }^{\prime }}{\mathbf{\xi }}_{k^{\prime }}.$$

$$\delta a_k^{\omega }=-\frac{1}{({\bar{\omega }}_k^2-{\omega }^2)}\sum _{k^{\prime }}\int d{\omega }^{\prime }{a}_{k^{\prime }}^{{\omega }^{\prime }}{\int }_{\odot }d\mathbf{r}{\mathbf{\xi }}_k^{\ast }\cdot \delta {\mathcal{L}}_{\omega -{\omega }^{\prime }}{\mathbf{\xi }}_{k^{\prime }}.$$

$$R_k^{\omega }=\frac{1}{{\bar{\omega }}_k^2-{\omega }^2}.$$

$${\Lambda }_{k^{\prime }}^k(\sigma )=-{\int }_{\odot }d\mathbf{r}{\mathbf{\xi }}_k^{\ast }\cdot \delta {\mathcal{L}}_{\sigma }{\mathbf{\xi }}_{k^{\prime }}.$$

$$\delta a_k^{\omega }=R_k^{\omega }\sum _{k^{\prime }}\int d{\omega }^{\prime }{a}_{k^{\prime }}^{{\omega }^{\prime }}{\Lambda }_{k^{\prime }}^k(\omega -{\omega }^{\prime })=R_k^{\omega }\sum _{k^{\prime }}\int d{\omega }^{\prime }{a}_{k^{\prime }}^{{\omega }^{\prime }}{\Lambda }_k^{k^{\prime }\ast }({\omega }^{\prime }-\omega ),$$

$$\delta a_k^{{\omega }^{\ast }}=R_k^{\omega \ast }\sum _{k^{\prime }}\int d{\omega }^{\prime }{a}_{k^{\prime }}^{{\omega }^{\prime }\ast }{\Lambda }_k^{k^{\prime }}({\omega }^{\prime }-\omega ).$$

$$⟨a_q^{\omega +\sigma }\delta a_k^{\omega \ast }⟩=R_k^{\omega \ast }\sum _{k^{\prime }}\int d{\omega }^{\prime }⟨a_q^{\omega +\sigma }{a}_{k^{\prime }}^{{\omega }^{\prime }\ast }⟩{\Lambda }_k^{k^{\prime }}({\omega }^{\prime }-\omega )=N_q{R}_k^{\omega \ast }|R_q^{\omega +\sigma }{|}^2{\Lambda }_k^q(\sigma ).$$

$$⟨\delta a_q^{\omega +\sigma }{a}_k^{\omega \ast }⟩=N_k{R}_q^{(\omega +\sigma )}|R_k^{\omega }{|}^2{\Lambda }_q^{k\ast }(-\sigma ),$$

$$⟨a_q^{\omega +\sigma }\delta a_k^{\omega \ast }+\delta a_q^{\omega +\sigma }{a}_k^{\omega \ast }⟩=(N_q{R}_k^{\omega \ast }|R_q^{\omega +\sigma }{|}^2+N_k{R}_q^{(\omega +\sigma )}|R_k^{\omega }{|}^2){\Lambda }_k^q(\sigma )=H{\Lambda }_k^q(\sigma ),$$

$$H(k,q,\sigma )=N_q{R}_k^{\omega \ast }|R_q^{\omega +\sigma }{|}^2+N_k{R}_q^{(\omega +\sigma )}|R_k^{\omega }{|}^2.$$

$${\xi }_{k,r}$$

$${\xi }_{k,\theta }={\gamma }_{\ell }V(r){\partial }_{\theta }{Y}_{\ell }^{0m}$$

$${\xi }_{k,\varphi }={\gamma }_{\ell }V(r)\frac{1}{\sin \theta }{\partial }_{\varphi }{Y}_{\ell }^{0m}$$

$${\mathbf{u}}_0$$

$$u_s^{-t}(r,\omega )=(-1{)}^t[u_s^t(r,-\omega ){]}^{\ast },v_s^{-t}(r,\omega )=(-1{)}^t[v_s^t(r,-\omega ){]}^{\ast },w_s^{-t}(r,\omega )=(-1{)}^t[w_s^t(r,-\omega ){]}^{\ast }.$$

$$u_r$$

$$u_{\theta }$$

$$u_{\varphi }$$

$$\delta ⟨a_{k^{\prime }}^{\omega +\sigma }{a}_k^{\omega \ast }⟩=H(k,k^{\prime },\sigma )\sum _{s=0}^{s_{\max }}\sum _{t=-s}^s{\int }_{\odot }dr[i{u}_s^t(r,\sigma ){\mathcal{K}}_u^{st}(r;k,k^{\prime })$$

$${\Lambda }_k^{k^{\prime }}(\sigma )=-\int d\mathbf{r}{\mathbf{\xi }}_{k^{\prime }}^{\ast }\cdot \delta {\mathcal{L}}_{\sigma }{\mathbf{\xi }}_k=2i\omega {\int }_{\odot }d\mathbf{r}\rho {\mathbf{\xi }}_{k^{\prime }}^{\ast }\cdot [{\mathbf{u}}_0(\mathbf{r},\sigma )\cdot \nabla ]{\mathbf{\xi }}_k,$$

$${\Lambda }_k^{k^{\prime }}={\Lambda }^u+{\Lambda }^v+{\Lambda }^w,$$

$$T_{rr}$$

$$T_{r\theta }$$

$$T_{r\varphi }$$

Full text

††pubyear: 2017††pagerange: Sensitivity of Helioseismic Measurements of Normal-mode Coupling to Flows and Sound-speed Perturbations–LABEL:lastpage

1 Introduction

The use of wavefield correlations to infer internal structure has a long history in seismology (e.g. Dahlen, 1968; Woodhouse, 1980; Lavely & Ritzwoller, 1992; Dahlen & Tromp, 1998; Woodard, 2006, 2014, 2016). The normal modes of oscillation are typically calculated with respect to a given model of the Sun (generally non-magnetic, non-rotating, spherically symmetric but not necessarily so), and in this state the modes are deemed “uncoupled”. However, given that our knowledge of the Sun’s interior is incomplete, there will be deviations in our model from reality, causing the modes to become coupled. For a temporally static perturbation to a given model, scattering can only occur across wavenumber, redistributing modal power in a manner that depends on the type of perturbation. The redistributed power is quantified by correlating the modes (i.e. spherical-harmonic time-series filtered around the corresponding resonant frequency) and relating the coupling to internal structure.

The study of convective flows in the solar interior is a topic of major interest in contemporary helioseismology. Recent measurements using time-distance helioseismogy (Duvall et al., 1993) have suggested that convective-flow amplitudes are significantly weaker than theory and numerical simulations suggest (Hanasoge et al., 2012). However, in contradiction to these results, Greer et al. (2015) measured large convective velocities using ring-diagram analysis (Hill, 1988), in agreement with convection simulations (ASH, Miesch et al., 2008). In particular, Woodard (2016) used mode-coupling theory to interpret correlations of Michelson Doppler Imager data (MDI; Scherrer et al., 1995) of spherical harmonic time series across different azimuthal orders $m$ at the same harmonic degrees $\ell$. The conclusion of the Woodard (2016) study was that the amplitudes are larger by an order in magnitude than those Hanasoge et al. (2012) obtained. In addition to understanding the physics of solar thermal transport, the resolution of this issue of convective-flow amplitudes is a critical to progress in related areas such as angular momentum transport and the solar dynamo (e.g. Hanasoge et al., 2016).

The simplicity of the data-handling procedure and elegance of the underlying interpretive theory make mode-coupling analysis a very attractive technique to study this flow inference problem. In the Woodard (2016) study, the inference of flow velocity magnitude was empirical and inversions to obtain flow magnitudes were not carried out. To that end, this manuscript aims to build the theory for performing inversions using these correlations.

Although set out in Lavely & Ritzwoller (1992), typographical errors appear in the expressions for kernels, specifically the lack of consistency between the toroidal-flow kernel in their equations (D20) and (C34) and to a lesser extent, a missing factor $i$ in their equation (C27). Moreover, the condition on flow conservation should limit the number of free scalar functions to two (as opposed to three), but Lavely & Ritzwoller (1992) write the inverse problem in terms of three scalar functions. To test the expressions for all kernels and address the issue of flow conservation, we set out to derive them here in the context of normal-mode measurements that are currently available from the Helioseismic and Magnetic Imager (HMI; Schou et al., 2012) and MDI missions.

2 Definitions and helioseismic measurements

The Lagrangian-displacement eigenfunction ${\boldsymbol{\xi}}_{k}$ of a non-rotating, non-attenuating, non-magnetic, spherically symmetric solar model, obeys (e.g. Christensen–Dalsgaard, 2003)

[TABLE]

where $\omega_{k}$ is the real resonant frequency, $\rho$ the density, $c$ the sound speed, ${\bf\hat{r}}$ the radially outward unit vector and $g$ the gravity. In the following analysis, we assume that the density is only a function of radius and is entirely determined (therefore not a function of horizontal coordinate or time). We have invoked Cowling’s approximation (e.g. Christensen–Dalsgaard, 2003) whereby perturbations to the gravitational potential of the Sun induced by waves are ignored. This is justified because we are primarily interested in waves that propagate in the convection zone which are well modelled in the Cowling limit (only very low-degree modes that traverse the radiative interior and core significantly perturb the gravitational potential). However, the present treatment may also be extended to the full oscillation equations, e.g. Lavely & Ritzwoller (1992). The operator (1) is Hermitian (e.g. Lynden-Bell & Ostriker, 1967) and therefore the eigenfrequencies and eigenfunctions for that operator are real. Note that for the sake of notational convenience, $k=(\ell,m,n)$ is used to denote a specific mode, typically categorised by three quantum numbers, spherical harmonic degree $\ell$, azimuthal order $m$ and radial order $n$. Although for non-rotating, spherically symmetric models, the modes are degenerate in $m$, i.e. $\omega_{\ell,m,n}=\omega_{\ell,0,n}$, we retain the dependence on $m$ for the sake of completeness here. We rewrite equation (1) for the sake of brevity as

[TABLE]

The eigenfunctions ${\boldsymbol{\xi}}_{k}$ satisfy the following orthonormality relationship

[TABLE]

where $d{\bf r}=r^{2}\,\sin\theta\,dr\,d\theta\,d\phi$ is an infinitesimal volume element, $r$ is radius, $\theta$ is co-latitude and $\phi$ is longitude. This is an important property for the analysis to follow.

Because the eigenfunctions are a complete set (with zero Lagrangian-pressure perturbation on the outer boundary, e.g. Christensen–Dalsgaard, 2003), a general wavefield $\mn@boldsymbol{\xi}({\bf r},\omega)$ may be written as

[TABLE]

where $\omega$ is temporal frequency, $a_{q}$ are coefficients representing the phase and amplitude of the $q$th mode and we write $a^{\omega}_{q}=a_{q}(\omega)$ for the sake of notational convenience. The governing equation for driven, undamped oscillations is given by

[TABLE]

where ${\bf J}$ is the unperturbed source, $\mn@boldsymbol{\xi}$ the zeroth-order wavefield and ${{\mathcal{L}}}_{0}$ is the unperturbed operator (whose properties are derived from a basic solar model such as model S; Christensen-Dalsgaard et al., 1996) as defined in equation (2). Note that the MDI and HMI missions provide time series of spherical-harmonic coefficients of oscillations, i.e. $a_{k}(\omega)$, measured at the solar photosphere. In practice, because we only observe the side of the Sun that faces us, it is not possible to precisely spatially decompose the measured oscillation signals into individual spherical harmonics. This partial visibility leads to leakage between modes across $\ell$, introducing systematical errors in eventual inferences (e.g. Larson & Schou, 2015; Woodard, 2016). Further, at a given frequency (between ridges for instance), many different radial orders contribute, leading to diminished ability in distinguishing between them. Temporal gaps in observations, due to instrumental or other issues, lead to mixing across temporal frequencies (and hence radial orders), thereby exacerbating this issue. To mitigate the problem of radial-order leakage, we use measurements at frequencies within a linewidth of resonance, ensuring that the associated radial order contributes a dominant fraction of observed mode power. Nevertheless, leakage across radial orders and spherical-harmonic degrees will systematically influence the results of mode-coupling analyses and need to be carefully considered.

2.1 Perturbation theory

The effects of rotation, flows, magnetic fields or thermal asphericities are modelled by adding a perturbation to the operator $\delta{{\mathcal{L}}}$ to ${{\mathcal{L}}}_{0}$ in equation (5), in turn leading to a perturbation to the wavefield, $\delta\mn@boldsymbol{\xi}$. Ignoring perturbations to the source ${\bf J}$, we obtain

[TABLE]

for the wavefield in a perturbed Sun. For time-varying perturbations, $\delta{{\mathcal{L}}}$ behaves as a convolution operator. To first order, we have therefore

[TABLE]

where the convolution over frequency captures the time variability of the perturbation operator $\delta{{\mathcal{L}}}$ (Eq. [85] of Appendix D). Perturbing equation (4) while assuming that the eigenfunctions are unmodified by the perturbation, we may express the change to the wavefield as a linear combination of the original set of eigenfunctions,

[TABLE]

Using equation (2), the left hand side of equation (6) may be simplified thus

[TABLE]

Solar modes are subject to a small amount of attenuation $\gamma_{k}\ll\omega_{k}$ that perturb the eigenfunctions and eigenfrequencies (obtained from Eq. [1]) by introducing small imaginary components. However, because the attenuation term is generally much smaller than the resonant frequency, we only consider it as a perturbation to the resonant frequency and ignore changes to the eigenfunctions. The inclusion of damping is driven to a large extent by the desire to apply the theoretical model to the analysis of oscillation data. Therefore, we substitute $\omega_{k}$ by $\bar{\omega}_{k}=\omega_{k}-i\gamma_{k}/2$ and continue to treat ${\boldsymbol{\xi}}_{k}$ as a real quantity. The sign of the attenuation term is negative owing to the Fourier convention we adopt (Appendix D). Dotting both sides of equation (7) by ${\boldsymbol{\xi}}^{*}_{k}$, using equation (9) and integrating over the solar volume, we obtain

[TABLE]

This gives us the coupling between modes $k$ and all of $k^{\prime}$,

[TABLE]

We quantify the spectral profile of the mode thus,

[TABLE]

The term in the denominator of equation (12) suggests that the coupling is strongest when the temporal frequency $\omega$ is very close to the resonant mode frequency and falls rapidly as the frequency difference grows. The modification of the frequency through the introduction of a small imaginary component (as described above) removes the singularity associated with the point $\omega=\omega_{k}$ in equation (12). The coupling integral between the two eigenstates $k$ and $k^{\prime}$ is defined thus,

[TABLE]

Equation (11) may thus be written as

[TABLE]

where for the perturbation operators of interest here, i.e. for sound-speed and flows, it is shown in Appendix C that $\Lambda^{k^{\prime}*}_{k}(\omega^{\prime}-\omega)=\Lambda^{k}_{k^{\prime}}(\omega-\omega^{\prime})$. We therefore consider the complex conjugate of equation (14)

[TABLE]

We now multiply both sides of equation (15) by $a^{{\omega+\sigma}}_{q}$, ensemble average and note that the excitation of waves in the Sun is well approximated as uncorrelated across modes (e.g. Woodard, 2007). We have $\langle a^{\omega+\sigma}_{q}\,a^{\omega^{\prime}*}_{k}\rangle=N_{q}\,|R^{\omega^{\prime}}_{q}|^{2}\delta(\omega+\sigma-\omega^{\prime})\,\delta_{qk}$, where the angular brackets indicate ensemble averaging and $N_{j}$ is the mode-amplitude normalization. Equation (15) is simplified thus

[TABLE]

Similarly,

[TABLE]

and using equations (16) and (17), we model the expectation value of the cross-spectral correlation signal,

[TABLE]

where

[TABLE]

Equations (18) and (19) define the essence of normal-mode-coupling problem in helioseismology. In the remaining part of this article, we focus on calculating the coupling term $\Lambda^{q}_{k}$ and expressing it in terms of the underlying perturbations.

2.2 Basis expansion of perturbations and eigenfunctions

Making use of terminology originally set out in Lavely & Ritzwoller (1992), we first introduce spherical harmonics $Y_{\ell}^{m}(\theta,\phi)$, the basis set on which eigenfunctions are projected. The generalised spherical harmonic $Y_{\ell}^{Nm}=D^{\ell}_{Nm}(\phi,\theta,0)$, where $D^{\ell}_{Nm}$ is the Wigner rotation matrix that relates spherical harmonics in rotated frames, will also play an important role in the analysis.

The mode eigenfunction describing oscillations in a non-rotating, non-magnetized, spherically symmetric model of the Sun may be written using spheroidal harmonics thus ${\boldsymbol{\xi}}_{k}=U(r)\,Y^{m}_{\ell}\,{\bf\hat{r}}+V(r)\,{\boldsymbol{\nabla}}_{h}Y^{m}_{\ell}$, where $k=(\ell,m,n)$ denotes a specific mode, and ${\boldsymbol{\nabla}}_{h}$ is the horizontal covariant derivative (e.g. Lavely & Ritzwoller, 1992; Christensen–Dalsgaard, 2003). Rewriting the eigenfunction in terms of generalised spherical harmonics and simplifying (see appendix A, Eq. [76] and Eq. [LABEL:specuse]), we obtain expressions for the three components $(r,\theta,\phi)$ of the eigenfunction,

[TABLE]

where $U(r)$ and $V(r)$ are functions of radius alone. Using vector spherical harmonics, a general time-varying (and therefore frequency dependent) flow ${\bf u}_{0}({\bf r},\omega)$ in a sphere may be written as

[TABLE]

where $s$ is the degree and $t$ is the azimuthal order. The first two terms capture the poloidal flows and the remaining term represents toroidal flow. In the inverse problem of determining flows, the goal is to infer the best-fit coefficients $u^{t}_{s}(r,\omega),v^{t}_{s}(r,\omega),w^{t}_{s}(r,\omega)$. The symmetries associated with ${\bf u}_{0}$ being real in the spatio-temporal domain give us the following relationships for the coefficients,

[TABLE]

The three components of the flow are (applying Eq. [76] and [LABEL:specuse]),

[TABLE]

Recalling equation (18), we seek to write an inverse problem (also the forward problem) of the form

[TABLE]

where the kernels ${{\mathcal{K}}}$ relate the sensitivity of the flow coefficients $u,v,$ and $w$ to the cross-spectral measurements and $H(k,k^{\prime},\sigma)$ was defined in equation (19). Equation (2.2) states the inverse problem for cross-spectral measurements (see also Eq. 35 of Woodard, 2014).

3 Mode coupling due to flows

We now analyse the coupling integral for flows that occurs in equation (18), i.e. where $\delta{{\mathcal{L}}}=-2i\omega\rho{\bf u}_{0}\cdot{\boldsymbol{\nabla}}$ (appropriate in the limit of anelastic flows, Gough, 1969), which must be evaluated in order to compute the kernels ${{\mathcal{K}}}$,

[TABLE]

where we had defined $\Lambda^{k^{\prime}}_{k}$ in equation (13). The perturbation operator for flows $\delta{{\mathcal{L}}}=-2i\omega\rho{\bf u}_{0}({\bf r},\sigma)\cdot{\boldsymbol{\nabla}}$ is seen to depend on two frequencies, $\omega$ and $\sigma$. Since modes couple only when their frequencies are very similar (see Eq. [12] and subsequent discussion), we invoke the approximation $\omega\approx\omega_{\rm ref}$, where $\omega_{\rm ref}$ is a reference frequency close to the frequencies of the two coupled modes. Note that this approximation applies because the degree of mode coupling weakens rapidly as the difference between the two resonant frequencies increases. The latter, $\sigma$, is the temporal frequency at which the flow evolves, i.e. ${\bf u}_{0}={\bf u}_{0}({\bf r},\sigma)$. Anelastic flows, i.e. ${\boldsymbol{\nabla}}\cdot(\rho\,{\bf u})=0$ vary on long time scales (Gough, 1969), implying that $\sigma\ll\omega_{\rm ref}$.The coupling matrix $\Lambda^{k^{\prime}}_{k}$ experiences contributions from poloidal $(u,v)$ and toroidal flows $w$. We therefore split the matrix into three components to make the computation of kernels easier,

[TABLE]

where dependencies on $k$ and $k^{\prime}$ are not explicitly noted on the right hand side.

Defining the tensor ${\bf T}={\boldsymbol{\nabla}}{\boldsymbol{\xi}}_{k}$ and using the expressions from (78) along with the definitions of mode eigenfunctions (20), we obtain equation (LABEL:Tsimple). Using recursion relations (76), we further simplify equation (LABEL:Tsimple),

[TABLE]

where the dot over a symbol indicates its radial derivative, i.e. ${\dot{V}}=\partial_{r}V(r)$.

4 Kernel for poloidal flow ($u^{t}_{s}$ coefficients)

The algebra to compute each of these coefficients is lengthy so we limit the discussion in the main body of the text to two coefficients, $u^{t}_{s}$ and $w^{t}_{s}$. The derivation of the kernel for $v^{t}_{s}$ follows the same procedure as that for $w^{t}_{s}$. To compute ${{\mathcal{K}}}_{u}$, we consider radial flows in equation (27),

[TABLE]

where the components $\xi_{k^{\prime},}$ are described in equation (20). Since the eigenfunctions are real (attenuation is assumed to only affect the resonant frequency), we cease to explicitly note the conjugate symbols on the $U$ and $V$ functions (Eq. [20]). Substituting expressions for the eigenfunction (Eq. [20]) and tensor $T$ (Eq. [29]) and expanding,

[TABLE]

where $U^{\prime},V^{\prime}$ belong to the eigenfunction ${\boldsymbol{\xi}}_{k^{\prime}}$. The Wigner 3-$j$ coefficient with $N=0$ is zero unless $\ell^{\prime}+s+\ell$ is even, justifying the appearance of the $1+(-1)^{\ell^{\prime}+s+\ell}$ term in the definition of the $B$ coefficient (see Appendix A).

[TABLE]

Next we simplify the following terms thus

[TABLE]

and

[TABLE]

where in equation (LABEL:maniple), we have use a standard manipulation of the Wigner symbol (see, for instance appendix C, Section d of Lavely & Ritzwoller, 1992). The sum of these terms gives the expression for the kernel for $u^{t}_{s}$,

[TABLE]

and the inverse problem for the $u$ coefficients based on equations (18) and (28) is

[TABLE]

where ${{\Lambda}}^{u}$ is the contribution by radial flows to the cross-spectral measurement between modes $k,k^{\prime}$. This is the first term of the expression (2.2) that we were seeking.

5 Kernel for poloidal flow ($v^{t}_{s}$ coefficients)

We do not derive the kernel here but state it for the sake of completeness. Note that it is identical to that derived by Lavely & Ritzwoller (1992) in their equation (C33),

[TABLE]

The inverse problem is similar to that stated in equations (18), (28) and (36)

[TABLE]

This is the second term of the target expression (2.2).

6 Kernel for toroidal flow ($w^{t}_{s}$ coefficients)

We collect contributions from two sets of terms

[TABLE]

We study some of the terms here.

6.1 Analyzing $u_{\theta}T_{\theta r}\xi^{*}_{k^{\prime},r}$

We study the first term in the integral in equation (39),

[TABLE]

[TABLE]

[TABLE]

[TABLE]

6.2 Analyzing $u_{\theta}T_{\theta\theta}\xi^{*}_{k^{\prime},\theta}$

Now we focus on the second term in the integral in equation (39),

[TABLE]

6.3 Overall sum

Albeit tedious, this analysis may be applied to the remaining terms in equation (39). Summing all terms together, we obtain

[TABLE]

The kernel for the $w$ coefficients is given by

[TABLE]

and the inverse problem is similar to that stated in equations (18), (28), (36) and (38)

[TABLE]

Note the difference in sign in the $UU^{\prime}-VU^{\prime}$ term for the $w^{t}_{s}$ kernel given in equation (44) and equation (C34) of Lavely & Ritzwoller (1992) - the latter appear to have made typographical errors. Equation (45) is the third term of the inverse problem that was posed in (2.2).

7 Kernels for anelastic flows

The anelastic condition on flows (Gough, 1969), ${\boldsymbol{\nabla}}\cdot[\rho{\bf u}_{0}({\bf r},\sigma)]=0$ applies to small Mach-number flows in stratified environments, relevant to convective flow in the solar interior (away from the surface layers). In principle, a time varying flow will cause the continuity equation to possess a time derivative of density as well. We ignore this on the basis that seismic signals are minimally influenced by density changes. The solenoidal condition on $\rho{\bf u}_{0}$ has the following consequences

[TABLE]

where ${\bf T}^{t}_{s}$ and ${\bf P}^{t}_{s}$ are the toroidal and poloidal components of the flow.

Consider the kernel for the $u^{t}_{s}$ coefficients, ${{\mathcal{K}}}_{u}$ (Eq. [35]). Introduce the symbol ${{\mathcal{F}}}$ that includes all factors not dependent on $r$ such that

[TABLE]

Now the inverse problem, given in equation (36) is

[TABLE]

The term (48) is symmetric between the primed and non-primed eigenfunctions whereas term (49) is anti-symmetric. However, the presence of the $i$ pre-multiplying both terms indicates that the terms are Hermitian (49) and anti-Hermitian (48) respectively. The modified anelatic-flow kernel for $u$, denoted by $\tilde{{\mathcal{K}}}^{st}_{u}$, will only contain the Hermitian term, and is given by

[TABLE]

Equation (50), the anelastic-flow kernel for $u$, is identical to the expression (D18) in Lavely & Ritzwoller (1992). The anti-Hermitian term is removed through integration by parts,

[TABLE]

where we have assumed that $u^{t}_{s}$ is zero on the boundaries (causing the total derivative term to drop out), and where we invoke the divergence-free condition (46). Therefore the symmetric term, i.e. equation (48), contributes to the kernel for the $v$-coefficients. Now to correct ${{\mathcal{K}}}^{st}_{v}$ to account for term (52), we first introduce the following relationships (see Appendix D, Eq. D7 of Lavely & Ritzwoller, 1992),

[TABLE]

Consider the definition of $B^{(1)+}_{s\ell\ell^{\prime}}$ (Appendix A),

[TABLE]

In the equation for ${{\mathcal{K}}}_{v}$ in equation (37), we manipulate the pair of terms

[TABLE]

Substituting expression (54) into equation (37), and adding term (52), we obtain an expression for the anelastic form of the kernel for $v$, denoted by $\tilde{{\mathcal{K}}}_{v}$,

[TABLE]

Applying the relationship $B^{(1)+}_{s\ell^{\prime}\ell}=\frac{1}{2}[s(s+1)+\ell(\ell+1)-\ell^{\prime}(\ell^{\prime}+1)]B^{(0)+}_{s\ell^{\prime}\ell}$ and $B^{(0)+}_{s\ell^{\prime}\ell}=B^{(0)+}_{\ell^{\prime}s\ell}$ (Eq. C45 and C46 of Lavely & Ritzwoller, 1992), we simplify equation (55) further to obtain the anelastic kernel for $v$,

[TABLE]

Equation (56), the anelastic-flow kernel for $v$, is identical to the expression (D19) in Lavely & Ritzwoller (1992). Note that $B^{(1)+}_{s\ell\ell^{\prime}}=B^{(1)+}_{\ell^{\prime}\ell s}$ and $B^{(1)+}_{s\ell^{\prime}\ell}=B^{(1)+}_{\ell\ell^{\prime}s}$ always hold.

Next, we manipulate the terms involving $VV^{\prime}$ in the expression for ${{\mathcal{K}}}_{w}$ in equation (44),

[TABLE]

We use the identities $B^{(1)-}_{\ell^{\prime}s\ell}=B^{(1)-}_{\ell^{\prime}\ell s}=B^{(1)-}_{s\ell\ell^{\prime}}$ (e.g., Eq. 46 of Lavely & Ritzwoller, 1992, and references therein) and

[TABLE]

we obtain the following relationship

[TABLE]

Rearranging equating (44), we obtain a simplified expression for the kernel for the $w$-coefficients, denoted as $\tilde{{\mathcal{K}}}^{st}_{w}$,

[TABLE]

Equations (60) and (44) are necessarily identical since toroidal flows automatically satisfy continuity. The anelastic kernel (60) is identical to expression (D20) of Lavely & Ritzwoller (1992) in two places; indeed the typographical error appears to not have carried forward from their equation (C34). Finally, noting again from the continuity equation that

[TABLE]

we substitute expression (61) for $v^{t}_{s}$ to obtain the following constrained inverse problem for anelastic flows that has only two scalar parameters,

[TABLE]

where we assume that the function $u^{t}_{s}\,\tilde{{\mathcal{K}}}_{v}$ vanishes at the boundary and

[TABLE]

Note that our formulation of the inverse problem differs from that of Lavely & Ritzwoller (1992) in that we only invert for two scalar functions, one a poloidal component $u^{t}_{s}$ and the other the toroidal, $w^{t}_{s}$. From $u^{t}_{s}$, we obtain the other poloidal component $v^{t}_{s}$ using the continuity equation, i.e. specifically equation (61). Computing Wigner-3$j$ symbols for high degrees and orders (i.e. $\ell>80$ for instance) represents a non-trivial challenge. We use an implementation of an algorithm described in Schulten & Gordon (1975) downloaded from the SLATEC library for these calculations. In Figures 1 and 3, we show plots of continuity-constrained kernels for $u^{t}_{s}$ and $w^{t}_{s}$ respectively.

Figure 4 shows kernels for anelastic poloidal flows based on equation (63), termed “fully anelastic” since it reduces the inference to a one-term inversion, compared with expressions (50) and (56), which result in a (redundant) two-term inversion problem for $u^{t}_{s}$ and $v^{t}_{s}$ (also described in Lavely & Ritzwoller, 1992).

In Figures 5 through  7, we compare the kernels for poloidal flow, computed using equation (56) and for toroidal flows using equation (60), with corresponding asymptotic expressions (see equations 43 through 47 of Woodard, 2014),

[TABLE]

and

[TABLE]

where equations (64) and (65) are asymptotic kernels sensitive to poloidal and toroidal flow components $v$ and $w$ respectively, $f$ and $g$ are asymptotic expressions for Wigner-3$j$ symbols in the limit $s\ll\ell,\ell^{\prime}$ (Vorontsov, 2011; Dahlen & Tromp, 1998) and $|\ell^{\prime}-\ell|\ll s$. These expressions also assume $|\ell-\ell^{\prime}|\ll\ell$ and $n=n^{\prime}$, which allow for the approximation $U\approx U^{\prime}$ and $V\approx V^{\prime}$. Since only eigenfunctions for a single $(n,\ell)$ mode appear in the expression, these asymptotic kernels cannot account for cross $\ell,n$ coupling. The asymptotics match well with exact expressions even in cases where $s$ is comparable $\ell,\ell^{\prime}$ and $|\ell^{\prime}-\ell|$.

8 Kernels for Sound-speed perturbations

Recalling the wave equation (1), sound-speed anomalies affect the wavefield through the perturbation operator $\delta{{\mathcal{L}}}_{\omega}=-{\boldsymbol{\nabla}}[\delta(\rho c^{2})\,{\boldsymbol{\nabla}}\cdot]$, where $c=c({\bf r},\omega)$ may be time varying. The coupling integral due to this perturbation, given by equation (13) is,

[TABLE]

where sound speed is positive-definite and time varying and where we have ignored the boundary terms (i.e. assuming that $\delta c=0$ on the outer boundary). We now expand $\delta\ln c$ in spherical harmonics,

[TABLE]

where because $c$ is real in the spatio-temporal domain, $c^{-t}_{s}(r,\sigma)=(-1)^{t}[c^{t}_{s}(r,-\sigma)]^{*}$, and

[TABLE]

where the tensor ${\bf T}$ is defined in equation (29) and the trace may be calculated as $Tr({\bf T})=T_{rr}+T_{\theta\theta}+T_{\phi\phi}$. Using equation (29) to compute the trace, substituting it into equation (66) and simplifying, we obtain the following expression for the kernel for the logarithm of sound speed,

[TABLE]

Note that we may take into account time variations of sound-speed anomalies in exactly the same as we did for flows (see Eq. [18]). Thus the inverse problem may be written as

[TABLE]

where $c^{t}_{s}(r,\sigma)$ represents the coefficient associated with the sound-speed anomaly temporally filtered at frequency $\sigma$. We reconstruct the sound speed thus

[TABLE]

where $c_{0}(r)$ is the fiducial sound speed of the spherically symmetric solar model for which the eigenfunctions have been computed. Examples of kernels for sound-speed perturbations are shown in Figure 8.

9 Discussion

The use of cross-spectral normal-mode measurements in helioseismology represents an exciting avenue to investigate the structure of convective flows, large-scale circulations and thermal anomalies in the solar interior. A description of the measurement process applied to obtain the cross-spectra is found in e.g. Woodard (2016). In this article, we have derived and computed kernels for flows and sound-speed perturbations. We have also resolved minor inconsistencies in the mode-coupling theory for helioseismology set out by Lavely & Ritzwoller (1992). Taking into account the constraint that the continuity equation imposes on flows, we also reformulate the inverse problem for two scalar functions. To summarise, the inverse problem for continuity-preserving flows is given in equation (62), where the kernels and related quantities are defined in equations (63), (50), (56), and (60). Other relevant symbols and a means to recompute the flow field in the spatial domain are described in equations (61), (19), and (21). The inverse problem for sound-speed perturbations is described in equation (70), and related expressions for kernels and reconstruction of the sound speed in spatial domain are described in equations (69), (19) and (71).

The kernels computed here can be used to study a variety of time-varying flow structures such as supergranulation and large-scale convection, meridional flow and associated sound-speed perturbations (should they be measurable). The formalism here can also be adopted for a starting model that is differentially rotating using perturbation theory outlined in, e.g. Gough (1990). Indeed, the theory may be extended to infer interior magnetism although the theory is complicated by the fact that these measurements are sensitive not to the magnetic field directly but rather the Lorentz stress tensor. One way to proceed is to reformulate linear magnetohydrodynamics as an analogous elastodynamics problem and use geophysical literature to obtain expressions for sensitivity kernels for the fourth-rank wavespeed tensor (e.g. Woodhouse, 1980; Dahlen & Tromp, 1998).

Acknowledgments

SMH, LG and KRS thank NYUAD’s Center for Space Science. SMH acknowledges support from the Max-Planck partner group program and Ramanujan fellowship SB/S2/RJN-73. MFW acknowledges support from NASA grant NNX14AH84G.

Appendix A List of Symbols

[TABLE]

We note that $\Omega^{\ell}_{0}=\Omega^{\ell}_{1}$, and $\Omega^{\ell}_{-1}=\Omega^{\ell}_{2}.$

A useful identity involving Wigner-3$j$ symbols is

[TABLE]

The Wigner-3$j$ symbol has a number of symmetry properties and the constituent parameters ($N,\ell$ etc.) must satisfy a variety of conditions for it to be non-zero. We do not describe these here but refer the reader to e.g. appendix C, Section d of Lavely & Ritzwoller (1992) for a detailed list of properties. The $B$-coefficient symbol that Woodhouse (1980) and Lavely & Ritzwoller (1992) significantly reduces the notational burden so we define and apply it in our analysis,

[TABLE]

Appendix B Useful Relations

[TABLE]

The following relations are particularly useful

[TABLE]

We note the following relationships,

[TABLE]

This allows us to write tensor ${\bf T}={\boldsymbol{\nabla}}\mn@boldsymbol{\xi}_{k}$ so

[TABLE]

Appendix C Properties of the Adjoints

Consider the flow operator $\delta{{\mathcal{L}}}=-2i\omega{\bf u}\cdot{\boldsymbol{\nabla}}$, and keeping in mind that because ${\bf u}$ is real in the spatio-temporal domain, ${\bf u}^{*}(\omega)={\bf u}(-\omega)$. Anelasticity (Gough, 1969) allows us to write ${\boldsymbol{\nabla}}\cdot(\rho{\bf u})=0$. Recalling the definition from equation (13), we have for $\Lambda^{k^{\prime}}_{k}(\omega)$,

[TABLE]

Similarly, for the scattering due to thermal anomalies, where $\delta{{\mathcal{L}}}=-{\boldsymbol{\nabla}}[\delta(\rho c^{2})\,{\boldsymbol{\nabla}}\cdot]$, and keeping in mind that $c({\bf r},-\omega)=c^{*}({\bf r},\omega)$, we have

[TABLE]

Appendix D Fourier convention

The following forward and inverse transforms connect the function $f(t)$ and its Fourier transform $\hat{f}(\omega)$,

[TABLE]

and

[TABLE]

The $\delta$ function is defined thus,

[TABLE]

A product in time domain $h(t)=f(t)g(t)$ transforms to the following convolution in Fourier domain,

[TABLE]