© ESO 2021

1. Introduction

Solar prominences are fascinating inhabitants of the solar corona that are made of relatively cool and dense plasma whose physical properties are akin to those in the chromosphere (see, e.g., Vial & Engvold 2015). High-resolution observations show that prominences are formed by a myriad of thin and long sub-structures, usually called threads, which are believed to outline particular field lines of the prominence global magnetic structure that is anchored at the photosphere (see, e.g., Lin 2011; Martin 2015). While it has been well established that the magnetic field provides the necessary force to support the prominence material against gravity, the energy balance in the prominence plasma is less understood. This balance requires an understanding of the processes involved with heating and cooling (see, e.g., Parenti 2014; Gilbert 2015; Heinzel 2015). According to Gilbert (2015), radiative heating due to the illumination from the surrounding solar atmosphere may be the dominant heating mechanism, however, determining the radiation field within prominences is not an easy task (Heinzel 2015). In addition, the prominence models that are obtained when balancing incident radiation and cooling (e.g., Heasley & Mihalas 1976; Anzer & Heinzel 1999; Heinzel et al. 2010; Heinzel & Anzer 2012) suggest that an additional source of heating must be present to explain the temperatures at prominence cores, which are estimated to be around 6000–8000 K (see Parenti 2014). Models of prominence threads under mechanical and thermal equilibrium typically require ad hoc heating to reproduce the expected density and temperature profiles (see the recent work by Terradas et al. 2020). So, the problem of energy balance remains open.

Observations show the ubiquitous presence of waves and oscillations in the fine structure of prominences (see the review by Arregui et al. 2018), which are interpreted as magnetohydrodynamic (MHD) waves (Ballester 2015). Hillier et al. (2013) suggested that transverse MHD waves in prominence threads are driven by horizontal motions in the photosphere. In this scenario, MHD waves can transport energy from the photosphere up to prominences suspended in the corona (see Soler et al. 2019). If an efficient enough dissipation mechanism is at work in the prominence plasma, MHD wave energy could be thermalized. Therefore, wave heating arises as another viable heating mechanism in prominences whose actual relevance needs to be determined.

The prospects of wave heating in prominences were previously discussed by Pécseli & Engvold (2000) and Parenti & Vial (2007), but their conclusions were not in agreement. Pécseli & Engvold (2000) estimated that the heating associated with the dissipation of Alfvén waves due to ion-neutral collisions would only contribute a tiny fraction to the bulk radiative energy. However, Pécseli & Engvold (2000) only considered Alfvén waves with a fixed frequency, while observational evidence has pointed to the presence of a broadband spectrum (Hillier et al. 2013). Conversely, Parenti & Vial (2007) estimated the Alfvén wave energy flux from measured nonthermal velocities and suggested that a large fraction of radiative losses could be compensated for. In their estimation, Parenti & Vial (2007) assumed that the entire energy flux would be dissipated, while the actual dissipated energy should depend on the efficiency of the damping mechanism at work.

To our knowledge, the first detailed investigation of wave heating in prominences was performed in Soler et al. (2016) using a very idealized configuration. The prominence was represented by a homogenous slab with a constant transverse magnetic field embedded in a homogeneous corona. No fine structure was considered. The prominence plasma was assumed to be partially ionized with an arbitrarily prescribed ionization degree, while the action of ion-neutral collisions was adopted as a dissipation mechanism. Alfvén waves were launched towards the prominence slab and the plasma heating rate due to wave dissipation was computed. Soler et al. (2016) estimated the volumetric heating integrated over the range of periods between 0.1 s and 100 s to be as large as 10% of the bulk radiative energy of the cool prominence plasma. Despite the important idealizations of the model of Soler et al. (2016), their results provide evidence supporting the potential of Alfvén wave heating in prominences.

The purpose here is to go further and deeper into the problem of prominence wave heating by improving the treatment of Soler et al. (2016). Although it is still simple, we shall consider a more elaborated model for a 1D thin thread with a longitudinal magnetic field. After assuming a constant pressure along the thread and imposing an ad hoc density profile, the temperature, and ionization profiles will be self-consistently computed considering both LTE and non-LTE approximations for the hydrogen ionization degree. The considered dissipation mechanisms will be Ohm’s diffusion due to magnetic resistivity and ambipolar diffusion due to ion-neutral collisions. A broadband driver for Alfvén waves will be considered at one end of the magnetic field line, aiming to mimic photospheric excitation. The plasma heating rate due to the dissipation of the waves will be computed and their relevance for the energy balance will be discussed.

This paper is organized as follows. Section 2 contains a description of the thin thread model. Section 3 includes the basic equations for the discussion of Alfvén waves and a description of the numerical method. The results are presented and discussed in Sect. 4, while some concluding remarks and prospects for future works are given in Sect. 5.

2. Background configuration

2.1. Thin thread model

We consider a simple 1D model for a prominence thin thread made of a single magnetic field line of length L = 10⁸ m. We neglect curvature and the effect of gravity, so that the magnetic field line is straight and aligned with the z-axis for convenience. The constant magnetic field strength is denoted by B. A typical value in quiescent prominences is B = 10 G. The center of the thread is located at z = 0 and its ends are at z = ±L/2. The pressure and density variation along prominence threads is linked to the thermal instability process that leads to the condensation of cool and dense plasma along the prominence magnetic field lines (see, e.g., Luna et al. 2012). Recently, Terradas et al. (2020) explored models of 1D prominence threads under hydrostatic and thermal equilibrium. The consideration of such detailed models is beyond the aim and purpose of the present work. Instead, we assume for simplicity that the thread has a constant pressure of p = 5  ×  10⁻³ Pa and that the density variation along the thread follows a Lorentzian dependence that aims to mimic the expected density profile along prominence threads (see Soler et al. 2015b) and is a plausible choice for the density profile according to Bayesian model comparison (see Arregui & Soler 2015), namely:

$$\begin{array}{c}\hfill \rho (z)=\frac{{\rho }_0}{1+4(\chi -1)z^2/L^2},\end{array}$$(1)

where ρ₀ = 10⁻¹⁰ kg m⁻³ is the density at the center of the thread and χ = 100 is the ratio between the density at z = 0 and that at z = ±L/2. The density profile is represented in Fig. 1a.

Fig. 1. Background magnitudes along the thread; (a) adopted density profile; (b) computed temperature; (c) computed ionization fraction; (d) Ohm’s diffusion coefficient; (e) ambipolar diffusion coefficient; (f) Cowling’s diffusion coefficient. Results for the LTE and non-LTE methods in panels b–f.

Once we have defined the pressure and density of the thread, the next step is obtaining the profiles along z for the temperature, T, and the ionization fraction, ξ_i, which is defined as the ratio of the ion density to the total density of the plasma. Due to the relatively low temperatures of prominences, the plasma is only partially ionized (see, e.g., Labrosse et al. 2010). Hence, the ionization fraction can vary from fully ionized plasma, ξ_i = 1, to almost neutral gas, ξ_i ≈ 0. Both ξ_i and T are related through the ideal equation of state, which for a plasma only composed by hydrogen reads

$$\begin{array}{c}\hfill p=(1+{\xi }_i)\rho RT,\end{array}$$(2)

where R is the ideal gas constant. However, we need another independent relation between ξ_i and T to close the system. This additional relation is the one that determines the ionization state of the plasma. Computing the ionization degree of a prominence thread in a consistent way requires the full solution of the radiative transfer equations taking into account the incident radiation (see, e.g., Gouttebroze & Labrosse 2009). To avoid that complicated calculations, here we use two different approximate approaches. The first one uses the Saha equation, which assumes local thermodynamic equilibrium (LTE) and the second one uses the tabulated values of the ionization fraction as a function of the temperature and pressure given by Heinzel et al. (2015), which are based on non-LTE radiative transfer computations in 1D prominence slabs¹. One of the purposes of this work is to compare the results obtained for both approaches. Throughout the paper, we use “LTE case” to refer to the method that involves the Saha equation and “non-LTE case” to refer to the one based on the computations of Heinzel et al. (2015).

On the one hand, the expression of ξ_i according to the Saha equation in the case of a hydrogen plasma is

$$\begin{array}{c}\hfill {\xi }_i=\frac{1}{2}\mathcal{M}(\sqrt{1+\frac{4}{\mathcal{M}}}-1),\end{array}$$(3)

where, in MKS units,

$$\begin{array}{c}\hfill \mathcal{M}\approx 4\times {10}^{-6}{\rho }^{-1}{T}^{3/2}exp(-T^{\ast }/T),\end{array}$$(4)

with T^* = 1.578  ×  10⁵ K. On the other hand, in the case of the numerical results by Heinzel et al. (2015) we use the values given in Table 1, which relate the temperature and the ionization fraction for the constant pressure considered here. In Fig. 2, we plot the ionization fraction as a function of temperature for both LTE and non-LTE approaches. In the case of the non-LTE method, we interpolated the values of Table 1 using a third-order polynomial and have imposed ξ_i = 1, that is, full ionization, for T ≥ 20 000 K. Extrapolation has been used when T < 6000 K. When the temperature increases, the ionization fraction in the LTE case reaches the completely ionized state much faster than in the case of the non-LTE case.

Fig. 2. Ionization fraction, ξi, as a function of temperature, T, for the LTE and non-LTE methods.

The resulting profiles of temperature and ionization fraction along the thread computed for both approaches can be seen in Figs. 1b and c, respectively. For the temperature, both methods give a similar result: a slightly modified parabolic-like profile with the minimum at the center of the thread and the maximum at the ends. The central temperature is ∼5000 K and the temperature at the ends is ∼3  ×  10⁵ K. Such values are consistent with what is expected in prominences. For the ionization fraction, both approaches provide a partially ionized zone around the thread center, where the temperature is lowest and the density is largest. The LTE method gives a narrower partially ionized region and the minimum value of ξ_i is less than half that obtained in the non-LTE method, for which there is a wider partially ionized region but the minimum ξ_i is larger. The ionization fraction at the center is ∼0.1 for the LTE case and ∼0.37 for the non-LTE case.

2.2. Diffusion coefficients

In partially ionized plasma, Ohm’s and ambipolar diffusion are two mechanisms that can dissipate Alfvén waves (see, e.g., Ballester et al. 2018). Ohm’s diffusion (also called electrical resistivity) is caused by the collisions of electrons with ions and neutrals, while ambipolar diffusion is caused by the collisions of neutrals with ions and electrons. The Ohm’s, η, and ambipolar, η_A, coefficients are calculated as:

$$\begin{array}{cc}& \eta =\frac{{\alpha }_e}{{\mathrm{\mu }}_0{e}^2{n}_e^2},\hfill \end{array}$$(5)

$$\begin{array}{cc}\hfill & {\eta }_{\mathrm{A}}=\frac{{\xi }_n^2}{{\mathrm{\mu }}_0{\alpha }_n},\hfill \end{array}$$(6)

where μ₀ is vacuum magnetic permeability, n_e and e are the electron number density and charge, respectively, ξ_n = 1 − ξ_i is the neutral atoms fraction, and α_e and α_n are the electron and neutral friction coefficients, respectively. The friction coefficients are computed as:

$$\begin{array}{c}\hfill {\alpha }_{\beta }={\displaystyle \sum _{{\beta }^{\prime }\ne \beta }{\alpha }_{\beta {\beta }^{\prime }},}\end{array}$$(7)

where β and β′ denote either electrons (e), ions (i), or neutrals (n), and α_ββ′ = α_β′β is the symmetric friction coefficient for collisions between two particular species (see, e.g., Braginskii 1965). For collisions between two charged particles, the friction coefficient is:

$$\begin{array}{c}\hfill {\alpha }_{\beta {\beta }^{\prime }}=\frac{n_{\beta }{n}_{{\beta }^{\prime }}{e}^{4}ln{\mathrm{\Lambda }}_{\beta {\beta }^{\prime }}}{6\pi \sqrt{2\pi }{ϵ}_0^2{m}_{\beta {\beta }^{\prime }}{(k_{B}T/m_{\beta {\beta }^{\prime }})}^{3/2}},\end{array}$$(8)

where $m_{\beta {\beta }^{\prime }}=m_{\beta }{m}_{\beta }^{\prime }/\left(m_{\beta }+m_{\beta }^{\prime }\right)$ is the reduced mass, e is the electron charge, k_B is Boltzmann’s constant, ϵ₀ is the permittivity of free space, n_β is the number density, and lnΛ_ββ′ is Coulomb’s logarithm, given by (see, e.g., Vranjes & Krstic 2013):

$$\begin{array}{c}\hfill ln{\mathrm{\Lambda }}_{\beta {\beta }^{\prime }}=ln\left(\frac{24\pi {ϵ}_0^{3/2}{k}_B^{3/2}{T}^{3/2}}{e^3\sqrt{n_{\beta }+n_{{\beta }^{\prime }}}}\right).\end{array}$$(9)

If at least one of the colliding particles is a neutral, so that β′ = n, the friction coefficient is

$$\begin{array}{c}\hfill {\alpha }_{\beta n}=n_{\beta }{n}_n{m}_{\beta n}{\left[\frac{8k_{B}T}{\pi m_{\beta n}}\right]}^{1/2}{\sigma }_{\beta n},\end{array}$$(10)

where σ_βn is the collision cross-section.

The profiles of the diffusion coefficients are presented in Fig. 1d for Ohm’s diffusion and in Fig. 1e for ambipolar diffusion. The profiles of the Ohm’s coefficient are similar in the two cases. The Ohm’s coefficient for the LTE case has a lower maximum but its profile is slightly wider than in the non-LTE case. In both cases, the maximum value is situated at the center of the thread and is ∼3  ×  10³ m² s⁻¹ for the LTE case and ∼4  ×  10³ m s² for the non-LTE case. Conversely, the profiles of the ambipolar coefficient are significantly different in the LTE and non-LTE cases. The ambipolar coefficient for the LTE case takes the maximum value, ∼1.6  ×  10¹⁴ m² s⁻¹ T⁻², at the center of the thread, whereas for the non-LTE case the maximum, ∼5.4  ×  10¹³ m² s⁻¹ T⁻², is displaced from the center and instead of a maximum there is a relative minimum at the center of the thread. Comparing the relative numerical values of both diffusion coefficients and assuming typical values of the magnetic field strength in prominences, it turns out that the effect of ambipolar diffusion should be much larger than that of Ohm’s diffusion, so that the ambipolar term is expected to be dominant in the calculation of wave dissipation. This simple estimation agrees with previous results from, for instance, Barceló et al. (2011) regarding the efficiency of Ohm’s and ambipolar diffusion in prominence conditions.

3. Method

3.1. Basic equations

To study the propagation of Alfvén waves along the thread, we used the MHD equations for a partially ionized plasma under the single-fluid approximation (see details in, e.g., Ballester et al. 2018). Alfvén waves are incompressible transverse waves driven by the magnetic tension force. We neglected all non-ideal terms in the basic single-fluid MHD equations, with the exception of Ohm’s and ambipolar diffusion. In order to study linear Alfvén waves, we linearize the equations by writing each variable as the sum of the background value and a small perturbation. Then, we only retain linear terms in the perturbations. After some calculations, the linearized equations that govern Alfvén waves are the transverse components of the momentum and induction equations, namely:

$$\begin{array}{cc}& \rho \frac{\partial v_{\perp }}{\partial t}=\frac{B}{{\mathrm{\mu }}_0}\frac{\partial B_{\perp }}{\partial z},\hfill \end{array}$$(11)

$$\begin{array}{cc}\hfill & \frac{\partial B_{\perp }}{\partial t}=B\frac{\partial v_{\perp }}{\partial z}+\frac{\partial }{\partial z}\left({\eta }_{\mathrm{C}}\frac{\partial B_{\perp }}{\partial z}\right),\hfill \end{array}$$(12)

where v_⊥ and B_⊥ are the velocity and magnetic field perturbations perpendicular to the background magnetic field, respectively, and η_C = η + B²η_A is Cowling’s (or total) diffusion coefficient. We plot in Fig. 1f Cowling’s coefficient, which varies many orders of magnitude along the thread. In the cool partially ionized region, η_C is dominated by the ambipolar coefficient multiplied by the magnetic field strength squared, whereas in the hot fully ionized region η_C = η.

We study Eqs. (11) and (12) in considering the stationary state of wave propagation. In this situation, it is assumed that the driver that excites the waves has been working for a sufficiently long time that is, at least, twice the Alfvén travel time along the thread. So, we can express the temporal dependence of the perturbations as exp(−iωt), where ω is the angular wave frequency. The linear wave frequency is simply f = ω/2π. Hence, we can combine Eqs. (11) and (12) to obtain a single partial differential equation for B_⊥ alone, namely

$$\begin{array}{c}\hfill \frac{{\partial }^2{B}_{\perp }}{\partial z^2}+\frac{\frac{\partial }{\partial z}(v_{\mathrm{A}}^2-\mathrm{i}\omega {\eta }_{\mathrm{C}})}{(v_{\mathrm{A}}^2-\mathrm{i}\omega {\eta }_{\mathrm{C}})}\frac{\partial B_{\perp }}{\partial z}+\frac{{\omega }^2}{(v_{\mathrm{A}}^2-i\omega {\eta }_{\mathrm{C}})}{B}_{\perp }=0,\end{array}$$(13)

where ${\mathit{v}}_{\text{A}}=B/\sqrt{{\mathrm{\mu }}_0\rho }$ is the Alfvén speed. In turn, v_⊥ is related to B_⊥ as

$$\begin{array}{c}\hfill v_{\perp }=\frac{\mathrm{i}}{\omega }\frac{v_{\mathrm{A}}^2}{B}\frac{\partial B_{\perp }}{\partial z}.\end{array}$$(14)

Equation (13) is our main equation and will be integrated using standard numerical methods after assuming boundary conditions at z = ±L/2.

3.2. Boundary conditions

To numerically solve Eq. (13), we assume that waves are driven at z = −L/2, with a prescribed amplitude and a given frequency, f. The driver imposes a magnetic field perturbation amplitude of the form:

$$\begin{array}{c}\hfill B_{\perp }=A\left(f\right)exp\left(\mathrm{i}\mathrm{\Phi }\left(f\right)\right)\mathrm{at}z=-\frac{L}{2},\end{array}$$(15)

where A(f) is the spectral weighting function and 0 < Φ(f) < 2π is a random phase different for each frequency. The result of waves excited by a broadband driver is constructed by varying the driver frequency in a wide range and superposing the solutions obtained for individual, discrete frequencies. We consider 1001 discrete driver frequencies between 0.1 and 1000 mHz in logarithmic spacing. The superposition is done after assuming a particular form for A(f). The spectral weighting function for Alfvén waves in solar prominences is largely unknown. Apart from Hillier et al. (2013), who performed a statistical study of transverse oscillations in a quiescent prominence, no other observational work has addressed this issue. Here, the spectral weighting function is assumed to be of a double power-law form, namely:

$$\begin{array}{c}\hfill A\left(f\right)=A_0\{\begin{array}{ccc}{\left(\frac{f}{f_b}\right)}^{1/4},\hfill & \mathrm{if}\hfill & f\le f_b,\hfill \\ {\left(\frac{f}{f_b}\right)}^{-5/6},\hfill & \mathrm{if}\hfill & f>f_b,\hfill \end{array}\end{array}$$(16)

where f_b is the break frequency of the spectrum and A₀ is a constant. The low-frequency exponent of 1/4 is considered after the results of Hillier et al. (2013). They showed that for periods, P, larger than 50 s, the velocity amplitude of the transverse oscillations approximately scales as ∼P^{−0.25 ± 0.04}. We set the break frequency to f_b = 20 mHz, which corresponds to the smallest period reported by Hillier et al. (2013) owing to instrumental constrains. Conversely, since we have no information about the scaling for higher frequencies, we simply assume a Kolmogorov-like exponent of −5/6. Finally, the value of the constant A₀ depends on the prescribed value of the driven energy flux.

For an Alfvén wave, the energy flux averaged over one full period is given by (see Walker 2004):

$$\begin{array}{c}\hfill ⟨\mathit{\pi }⟩=-\frac{1}{2{\mathrm{\mu }}_0}\mathrm{Re}\left[v_{\perp }{B}_{\perp }^{\ast }\right]\mathbf{B},\end{array}$$(17)

where * denotes the complex conjugate. The energy flux can be separated into the parallel and anti-parallel contributions with respect to the magnetic field direction, which are denoted by ⟨π⟩^↑ and ⟨π⟩^↓, respectively. In order to write these contributions, we use the Elsässer variables (Elsasser 1950), which are defined as

$$\begin{array}{cc}& Z^{\uparrow }=v_{\perp }-\frac{B_{\perp }}{\sqrt{{\mathrm{\mu }}_0\rho }},\hfill \end{array}$$(18)

$$\begin{array}{cc}\hfill & Z^{\downarrow }=v_{\perp }+\frac{B_{\perp }}{\sqrt{{\mathrm{\mu }}_0\rho }},\hfill \end{array}$$(19)

where Z^↑ and Z^↓ represent parallel-propagating and anti-parallel-propagating Alfvénic disturbances, respectively. With these new variables, we can write the energy flux as ⟨π⟩=⟨π⟩^↑ − ⟨π⟩^↓, where

$$\begin{array}{cc}& {⟨\mathit{\pi }⟩}^{\uparrow }=\frac{1}{8}\sqrt{\frac{\rho }{{\mathrm{\mu }}_0}}{Z}^{\uparrow }{Z}^{\uparrow \ast }\mathit{B},\hfill \end{array}$$(20)

$$\begin{array}{cc}\hfill & {⟨\mathit{\pi }⟩}^{\downarrow }=\frac{1}{8}\sqrt{\frac{\rho }{{\mathrm{\mu }}_0}}{Z}^{\downarrow }{Z}^{\downarrow \ast }\mathit{B}.\hfill \end{array}$$(21)

The fluxes ⟨π⟩^↑ and ⟨π⟩^↓ represent propagation of energy in the directions parallel and anti-parallel to the magnetic field, respectively.

At the driver location, that is, at z = −L/2, ⟨π⟩^↑ corresponds to the energy flux injected by the driver, while ⟨π⟩^↓ is the reflected energy flux. The energy flux of Alfvén waves driven at the photosphere is believed to be of the order of 10⁴ W m⁻² (see, e.g., De Pontieu et al. 2001; Goodman 2011; Tu & Song 2013; Arber et al. 2016). However, due to the filtering effect of the chromosphere, only about 1% of the photospheric flux is capable of reaching coronal heights (see a detailed study in Soler et al. 2019). Therefore, we set ∑_f⟨π⟩^↑ = 10² W m⁻² at z = −L/2, which allows us to determine the value of A₀.

On the other hand, the boundary condition at z = L/2 is set according to the imposed relation between ⟨π⟩^↑ and ⟨π⟩^↓. The paradigmatic cases are the perfectly reflecting boundary and the perfectly transparent boundary. In the perfectly reflecting boundary ⟨π⟩^↑ = ⟨π⟩^↓, then ⟨π⟩ = 0. This condition requires that ${\mathit{v}}_{\perp }{B}_{\perp }^{\ast }=0$. If B_⊥ is arbitrary, then v_⊥ = 0. Using Eq. (14), we obtain

$$\begin{array}{c}\hfill \frac{\partial B_{\perp }}{\partial z}=0\mathrm{at}z=\frac{L}{2},\end{array}$$(22)

for the perfectly reflecting boundary. In the case of a perfectly transparent boundary, we have ⟨π⟩^↑ ≠ 0 and ⟨π⟩^↓ = 0, which is equivalent to Z^↓ = 0. Using Eqs. (14) and (19) we obtain:

$$\begin{array}{c}\hfill \frac{\partial B_{\perp }}{\partial z}=\frac{i\omega }{v_A}{B}_{\perp }\mathrm{at}z=\frac{L}{2},\end{array}$$(23)

for the perfectly transparent boundary. In general, a realistic boundary condition would be neither perfectly reflective nor perfectly transparent. Our model does not include the chromosphere and photosphere where the magnetic field of prominences is actually anchored. It is the properties of the plasma in those layers that determine the reflection or transmission of the waves (see, e.g., Leroy 1980; Hollweg 1981; Similon & Zargham 1992; Soler et al. 2017, 2019). Based on the above results, we can define a general boundary condition as

$$\begin{array}{c}\hfill \frac{\partial B_{\perp }}{\partial z}=\epsilon \frac{i\omega }{v_A}{B}_{\perp }\mathrm{at}z=\frac{L}{2},\end{array}$$(24)

where ε ∈ [0,1] is a parameter. Checking the extreme values of ε we recover the perfectly reflecting boundary for ε = 0 and the perfectly transparent boundary for ε = 1. Values of ε in between those limit cases correspond to intermediate situations.

3.3. Wave absorption and heating rate

We define the reflectivity, ℛ, and transmissivity, 𝒯, coefficients, which physically describe the fraction of the driven energy that is reflected back to the driver location or transmitted all the way through the opposite boundary of the model, respectively. We compute the incoming, reflected, and transmitted fluxes as (see Soler et al. 2019):

$$\begin{array}{cc}& {⟨\pi ⟩}_{\mathrm{inc}.}=\left|{⟨\mathit{\pi }⟩}^{\uparrow }\right|\mathrm{at}z=-\frac{L}{2},\hfill \end{array}$$(25)

$$\begin{array}{cc}\hfill & {⟨\pi ⟩}_{\mathrm{ref}.}=\left|{⟨\mathit{\pi }⟩}^{\downarrow }\right|\mathrm{at}z=-\frac{L}{2},\hfill \end{array}$$(26)

$$\begin{array}{cc}\hfill & {⟨\pi ⟩}_{\mathrm{tra}.}=|⟨\mathit{\pi }⟩|\mathrm{at}z=\frac{L}{2}.\hfill \end{array}$$(27)

With the help of these fluxes, the coefficients are readily computed as

$$\begin{array}{cc}& \mathcal{R}=-\frac{{⟨\pi ⟩}_{\mathrm{ref}.}}{{⟨\pi ⟩}_{\mathrm{inc}.}},\hfill \end{array}$$(28)

$$\begin{array}{cc}\hfill & \mathcal{T}=\frac{{⟨\pi ⟩}_{\mathrm{tra}.}}{{⟨\pi ⟩}_{\mathrm{inc}.}}.\hfill \end{array}$$(29)

In a stationary system with no dissipation, ℛ + 𝒯 = 1 because of energy conservation. However, the presence of dissipation causes ℛ + 𝒯 < 1. Then, we can define the absorption coefficient, 𝒜, which represents the fraction of the driven wave energy that is absorbed or deposited in the plasma because of dissipation, namely:

$$\begin{array}{c}\hfill \mathcal{A}=1-(\mathcal{R}+\mathcal{T}).\end{array}$$(30)

These three coefficients are functions of the driver frequency but are independent of the driver amplitude. They are intrinsic properties of the model.

The consequence of the wave energy absorption is the heating of the plasma. The heating is calculated using the Joule heating function, which represents the dissipation of electric currents associated with the Alfvén waves. The heating function, Q, can be written as

$$\begin{array}{c}\hfill Q={\mathrm{\mu }}_0\eta {\left|{\mathit{j}}_{\mathbf{\Vert }}\right|}^2+{\mathrm{\mu }}_0{\eta }_{\mathrm{C}}{\left|{\mathit{j}}_{\mathbf{\perp }}\right|}^2,\end{array}$$(31)

where j_∥ and j_⊥ are the perturbations of the parallel and perpendicular components of the electric current with respect to the background magnetic field. We see that Cowling’s diffusion, that is, the joint effect of Ohm’s and ambipolar diffusion, is responsible for the dissipation of perpendicular currents, whereas parallel currents are dissipated by Ohm’s diffusion alone (see, e.g., Khomenko & Collados 2012). In our case, we have that

$$\begin{array}{cc}& {\mathit{j}}_{\Vert }=0,\hfill \end{array}$$(32)

$$\begin{array}{cc}\hfill & {\mathit{j}}_{\perp }=\frac{1}{{\mathrm{\mu }}_0}\frac{\partial B_{\perp }}{\partial z}{\widehat{e}}_{\perp },\hfill \end{array}$$(33)

where ${\widehat{e}}_{\perp }$ denotes the unit vector in the perpendicular direction. There are no parallel currents in our model, so that all the heating is produced by the dissipation of perpendicular currents. Therefore, the heating function averaged over one full period of the wave is:

$$\begin{array}{c}\hfill ⟨Q⟩=\frac{{\eta }_{\mathrm{C}}}{2{\mathrm{\mu }}_0}{\left|\frac{\partial B_{\perp }}{\partial z}\right|}^2.\end{array}$$(34)

The computation of the heating rate in the case of a broadband spectrum is done in the same fashion as for the energy flux. Namely, the heating function is calculated for all the individual frequencies in the spectrum using Eq. (34) and then the results are added together according to the spectral weighting function.

4. Results

4.1. Magnetic field and velocity perturbations

We start by studying the magnetic field and velocity perturbations along the thread because of the presence of the broadband spectrum of Alfvén waves. These results are displayed in Fig. 3. The perturbations are complex quantities and we plot their real part only.

Fig. 3. Real part of the (a) magnetic field and (b) velocity perturbations along the thread for the non-LTE case with ε = 0, 0.5, and 1; (c) and (d) are same as (a) and (b) but with ε = 0 and comparing the results for LTE and non-LTE cases.

Figure 3a and b compare the perturbations for ε = 0, 0.5, and 1 in the non-LTE case. Before discussing these results, we note that the random phase assumed during superposition of the different frequencies in the spectrum is different for the three values of ε. This implies that the fine details of the obtained perturbations are also intrinsically different.

The driven broandband spectrum naturally results in a mixture of many spatial scales (wavelengths) along the thread. However, it is evident from Figs. 3a and b that smaller wavelengths are present in the central, denser part of the thread than in the surrounding evacuated zones. The reason for this result is that the density enhancement at the thread center causes the local decrease of the Alfvén speed. The values of the Alfvén speed at the tube center and at its ends are ∼89 km s⁻¹ and ∼892 km s⁻¹, respectively. For an Alfvén wave, the wavelength is proportional to the Alfvén speed. The consequence of this is the decrease of the wavelength in the central part (for a similar result but in the case of chromospheric waves see, e.g., De Pontieu et al. 2001; Zaqarashvili et al. 2013; Soler et al. 2017).

The amplitude of the magnetic field perturbation is larger at the thread center than at the ends, while the opposite behavior is found in the case of the velocity perturbation. Again, this can be explained by the dependence of the density along the thread and the fact that in an Alfvén wave there is equipartition between kinetic and magnetic energies (Ferraro & Plumpton 1958). The maximum value of the magnetic field perturbation found near the thread center is ∼1 G, which is one order of magnitude below the equilibrium magnetic field strength. In turn, the maximum values of the velocity found near the tube ends are ∼25 km s⁻¹, which is again much smaller than the local Alfvén speed. Therefore, the obtained amplitudes justify the use of linearized theory. In addition, the obtained velocity amplitudes are of the same order as the observationally reported Doppler velocities in small-amplitude prominence oscillations (see the review by Arregui et al. 2018), which indicates that the considered influx of energy by the driver is consistent with the observed amplitudes.

The results obtained for the three different values of ε are qualitatively similar. The case with ε = 0, that is, total reflection at the right boundary, is the one that differs most from the others. For this case, larger magnetic field perturbations are found at the right end. Clearly, the condition that the waves are forced to reflect at the right boundary has an impact on the magnetic field perturbation. Also, slightly larger velocity perturbations are obtained in the evaluated part of the tube when ε = 0. There are no other outstanding differences.

On the other hand, Figs. 3c and d compare the perturbations obtained for the LTE and non-LTE cases for the particular value of ε = 0. Both magnetic field and velocity perturbations have similar behaviors in the two cases, as already discussed above. The main difference between the two cases resides in the velocity perturbation, which is larger for the non-LTE case near the tube ends. It is difficult to determine whether this disparity is actually related to the different profiles used in the LTE and non-LTE cases or, on the contrary, it is simply a consequence of the different random phases considered in the superposition of frequencies.

4.2. Reflectivity, transmissivity, and absorption

Here, we show the computations of the coefficients of reflectivity, ℛ, transmissivity, 𝒯, and absorption, 𝒜. These coefficients are calculated using Eqs. (28)–(30).

In Fig. 4 we plot the three coefficients as functions of the driver frequency, f, for the non-LTE case and compare the results of ε = 0, 0.5, and 1. The reflection and transmission coefficients, shown in Figs. 4a and b respectively, are remarkably different for ε = 0 from those obtained for ε = 0.5 and 1, so we start by discussing the results with ε = 0. In this case, 𝒯 = 0 because of the imposed boundary condition: the waves are not allowed to be transmitted through the right end of the tube. As a consequence, the reflectivity is found to be maximal, that is, ℛ ≈ 1, for f ≲ 10 mHz, approximately. There is no energy deposition into the plasma for those low frequencies since all the energy returns back to the driver location. For higher frequencies the reflectivity starts to decrease as a result of the increase of the absorption (see Fig. 4c). For f ≳ 300 mHz, approximately, 𝒜 ≈ 1, so that all the wave energy is absorbed into the thread plasma. These findings are consistent with the well-known result that the efficiency of the damping due to Cowling’s diffusion increases with the wave frequency (see the review by Ballester et al. 2018, and references therein).

Fig. 4. Reflectivity (a); transmissivity (b); and absorption as functions of the frequency (c) in log scale for ε = 0, 0.5, and 1 in the non-LTE case. The purple dotted line in panel a corresponds to the sum of reflectivity and transmissivity for ε = 0.5.

The behavior of the reflection and transmission coefficients when ε = 0.5 and 1 is very different from that when ε = 0. The most obvious disparity is the nonzero value of 𝒯. Indeed, the transmission coefficient is found to be significantly larger than the reflection coefficient in most parts of the frequency range, which indicates that the condition of no transmission imposed when ε = 0 affects dramatically wave behavior. As expected, the larger ε, the larger 𝒯 for a given frequency. Another remarkable result is the nonmonotonic behavior of ℛ and 𝒯 for low and intermediate frequencies. Both coefficients display oscillations, with the first one being of very large amplitude. In the appendix, we show that the presence of these oscillations is caused by resonances associated with the standing modes of the thread. The oscillations in 𝒯 are highlighted in Fig. 5 for a particular case. The oscillations of the coefficient have their minimums at the standing mode frequencies, approximately, which are indicated by vertical lines. This clearly confirms that the nature of the oscillations is related to the existence of resonances with the standing modes. Physically, the result that 𝒯 is minimum and ℛ is maximum at the standing mode frequencies can be explained by the fact that a standing wave, represented as the superposition of two waves propagating in opposite directions, has a null net energy flux. We refer readers to the appendix for a study of the standing modes.

Fig. 5. Transmissivity coefficient in the non-LTE case for ε = 1 in the range of frequencies corresponding to the first ten standing modes of the thread (vertical black lines).

Putting the oscillations aside, the global trend of the coefficients shows the decrease (increase) of the reflectivity (transmissivity) with the frequency, until the absorption starts to play a role for sufficiently high frequencies. No oscillations are present once the absorption begins to rise. This is so because dissipation starts to dominate and diminishes the role of the resonances (see also Soler et al. 2016). Thus, regarding the oscillations in the reflectivity or transmissivity, they are not present for frequencies lower than those of the standing mode frequencies because of the absence of resonances. They are neither present for very high frequencies because dissipation is strong enough to suppress the oscillations.

Concerning the absoption coefficient, there are no important differences between the results for various values of ε. The case with ε = 0 reaches total absorption, that is, 𝒜 ≈ 1, faster than the other cases, which is consistent with the fact that, in this situation, no fraction of wave energy can escape through the right end of the thread, and so there is more energy available in the system to be dissipated.

To further understand why for ε = 0 there are no oscillations related to resonances with standing modes in the reflectively, we have overplotted in Fig. 4a the sum of reflectivity and transmissivity, namely ℛ + 𝒯, for the ε = 0.5 case. Surprisingly, ℛ + 𝒯 does not display oscillations and the result matches quite well with the reflectively of the ε = 0 case. The reason is that when ε ≠ 0 the oscillations in ℛ are out of phase with respect to those of 𝒯 in a fashion so that they cancel out when ℛ and 𝒯 are added. This implies that resonances must surely be present in the ε = 0 case too, but the reflectively is not the appropriate parameter to show them. As shown later in Sect. 4.3, resonances clearly appear in the integrated heating rate when ε = 0.

The results discussed above correspond to the non-LTE case. In Fig. 6, we compare the three coefficients obtained in the LTE and non-LTE cases for ε = 0.5. The reflectivity is almost identical for both cases in the whole frequency range. The transmissivity is also almost the same for low and intermediate frequencies until the absorption begins the grow. In this model, the wave reflection and transmission is essentially governed by the longitudinal variation of density, which is the same in both cases. This explains why the reflectivity and transmissivity are identical in the LTE and non-LTE cases. Conversely, the wave damping depends on the values of Cowling’s diffusion, which are different for the two cases. In the LTE case total absorption is reached faster than in the non-LTE case, which indicates that wave damping is more efficient in the LTE case. This is consistent with the larger values of Cowling’s diffusion obtained in the LTE case (see again Fig. 1f).

Fig. 6. Same as Fig. 4 but for a comparison between LTE and non-LTE cases with ε = 0.5.

A remarkable difference with the results of Soler et al. (2016) is that here the absorption coefficient is a smooth, monotonic function of the frequency, while in Soler et al. (2016) the absorption was characterized by the presence of narrow peaks associated with resonances with standing modes. We do find resonances in the reflection and transmission coefficient, but not in the absorption coefficient. Conversely, the shape of our absorption coefficient is very similar to that computed in the chromosphere (see Soler et al. 2017). In addition, the absorption values obtained by Soler et al. (2016) for nonresonant frequencies were significantly smaller than the absorption computed here, so that globally the absorption of wave energy is larger in the present case. We believe that the distinct feature that causes our results to be different from those of Soler et al. (2016) is the nonuniform thread model considered here, whereas Soler et al. (2016) considered a uniform prominence slab embedded in a uniform corona. The piecewise model of Soler et al. (2016) favors the trapping of waves within the prominence slab, and this may artificially enhance the role of resonances. A similar argument was used by Soler et al. (2017) to explain the absence of resonances in their nonuniform chromospheric model as compared with the resonances obtained in piecewise models (see Hollweg 1981). This highlights the importance of using nonuniform models to correctly describe the reflection, transmission, and absorption properties of the waves.

4.3. Heating rates

The total volumetric heating rates associated with wave dissipation are plotted in Fig. 7. The heating profiles for the LTE and non-LTE profiles have almost the same global shape, which is consistent with the previous results. Heating is maximum at the thread center and decreases by many orders of magnitude in the evacuated, coronal part. This result is related with the relative values of Cowling’s diffusion coefficient in these regions. In particular, in the LTE case the central region where heating is significant is narrower than in the non-LTE case, and the maximum value of the heating rate is slightly larger for the LTE case than for the non-LTE case. In both cases, the order of magnitude of the volumetric heating rate at the center of the thread is ∼10⁻⁷ W m⁻³, while at the thread ends it is ∼10⁻¹⁶ W m⁻³.

Fig. 7. Wave heating function (for ε = 0, 0.5, and 1), cooling function, and thermal conduction term along the thread for the (a) LTE and (b) non-LTE cases.

The heating profiles have similar forms for the three different boundary conditions considered. In the cases of total and partial transmission, the results almost overlap and are barely distinguishable in Fig. 7. Somewhat smaller heating rates are obtained in the case of total reflection (ε = 0). This is a counter-intuitive result, since we showed before in Sect. 4.2 that, for the same wave frequency, the absorption coefficient is larger in the ε = 0 case. In order to understand this result, we study the contribution to the heating rate of the different frequencies present in the broadband spectrum. To do so, we compute the heating rate integrated along the whole domain, namely:

$$\begin{array}{c}\hfill {⟨Q⟩}_{\mathrm{int}}={\displaystyle {\int }_{-L/2}^{L/2}⟨Q⟩\mathrm{d}z,}\end{array}$$(35)

where the integration is performed for each individual frequency in the spectrum. Figure 8 shows the integrated heating as a function of the frequency for the three different boundary conditions in the non-LTE case. The results for the three boundary conditions display a similar trend. The heating rate increases with the frequency in the low frequency region until reaching the spectrum peak frequency at f = 20 mHz. Then, the heating rate saturates to a roughly constant value until it starts to decrease for f ≳ 200 mHz. This behavior can be related to the form of the spectral weighting function, which increases in the low frequency region, has its maximum at f = 20 mHz, and decreases afterwards. The result that for an intermediate region of frequencies, corresponding to 20 mHz ≲ f ≲ 200 mHz, the heating rate is roughly constant despite of the decrease of the weighting function can be understood by the fact that the efficiency of wave damping due to Cowling’s diffusion increases with the frequency. These two opposite effects, namely the decrease of the spectral weighting function and the increase of the damping efficiency, seem to balance each other almost perfectly.

Fig. 8. Integrated heating rate along the thread as a function of the wave frequency in the non-LTE case with ε = 0, 0.5, and 1.

In addition to the global trend discussed above, the curves in Fig. 8 have also some oscillations that can again be associated with resonances with the standing modes of the system (see the appendix for more). The presence of resonances in the heating rate is a result already present in the previous paper by Soler et al. (2016). Importantly, the amplitude of these oscillations depends on the value of ε. When ε = 0.5 and 1, the oscillations have a small amplitude and they barely affect the global trend. However, in the case of total reflection (ε = 0) the oscillations have large amplitudes that heavily distort the shape of the curve. The fact that the presence of the resonances have a more significant effect in the ε = 0 case is related to the condition that for standing modes the boundaries of the domain must behave as perfectly reflecting boundaries. The integrated heating rate in the ε = 0 case is larger than that for ε = 0.5 and 1 between the resonant frequencies, but this only happens in very narrow intervals. Generally, the heating rate for ε = 0 is lower than that for ε = 0.5 and 1, and this happens in wider intervals. The net result is that the total heating for ε = 0 is smaller, as shown in Fig. 7.

To shed more light on the oscillatory behavior of the integrated heating rate with the frequency, Fig. 9 shows the integrated heating in the range of frequencies corresponding to the first ten standing modes. The integrated heating has relative minimums at the standing modes frequencies, which are marked with vertical lines in Fig. 9. This result reinforces the idea of the existence of resonances between the propagating waves and the standing modes of the thread, as already shown in Fig. 5. The fact that the heating is minimum at the standing modes frequencies appears to stand in contradiction to the general idea of a resonance. This result is also opposite to that of Soler et al. (2016), where the maximums in the heating rate were found at the prominence slab resonance frequencies. The reason why the heating rate is minimum and not maximum at the resonance frequencies is explored in the appendix. In short, the conclusion is that the excitation of standing modes via resonances is disadvantageous with regard to the heating rate. Instead of being directly dissipated, part of the energy of the propagating waves is transferred to the standing modes at the resonance frequencies, but standing modes are very weakly damped by Cowling’s diffusion and their contribution to the heating is negligible. The net result is that the heating rate is smaller than the would-be heating in the absence of resonances. In turn, the discrepancy with the results of Soler et al. (2016) can be explained by the fact that in their piecewise slab model, the transmission of waves to the prominence was only efficient at the resonance frequencies because of the so-called cavity resonances (see, e.g., Hollweg 1981). So, in the case of Soler et al. (2016) only the wave energy transmitted into the slab at the resonance frequencies was available to be dissipated efficiently. Conversely, in the present nonuniform model, wave transmission to the dense part of the thread is much more efficient than in Soler et al. (2016), as evidenced by the larger absorption obtained in the present nonuniform model. Another consequence is that the heating rates obtained here are larger than those obtained in Soler et al. (2016).

Fig. 9. Integrated heating rate along the thread as a function of the wave frequency in the non-LTE case with ε = 0, 0.5, and 1. Same as Fig. 8 but in the range of frequencies corresponding to the first ten standing modes of the thread (vertical black lines).

4.4. Discussing the energy balance

Finally, we discuss the role that wave heating may have in the plasma energy balance. To do so, we evaluate the heat-loss function, ℒ, which is given by:

$$\begin{array}{c}\hfill \mathcal{L}=-\mathrm{\nabla }·\mathit{q}-L+⟨Q⟩+C,\end{array}$$(36)

where q = −κ∇T is the heat flux due to thermal conduction with κ thermal conductivity, L is the radiative cooling function, ⟨Q⟩ is the wave heating function calculated before, and C is an additional unspecified source of heating. The term C may represent the heating function due to incident coronal radiation. To appropriately treat this term, we should solve the radiative transfer problem (see Heinzel 2015). Since this is beyond the scope of the present paper, for simplicity we set C = 0 and omit the role of incident coronal radiation in these estimations.

In a partially ionized plasma, the parallel thermal conductivity to the magnetic field can be approximated by κ ≈ κ_e + κ_n, where κ_e and κ_n are the electron and neutral thermal conductivities, respectively, given by (see, e.g., Soler et al. 2015a)

$$\begin{array}{cc}& {\kappa }_{\mathrm{e}}=3.2\frac{n_e^2{k}_B^{2}T}{{\alpha }_{\mathrm{e}}+{\alpha }_{\mathrm{ee}}},\hfill \end{array}$$(37)

$$\begin{array}{cc}\hfill & {\kappa }_{\mathrm{n}}=\frac{5}{3}\frac{n_n^2{k}_B^{2}T}{{\alpha }_{\mathrm{n}}+{\alpha }_{\mathrm{nn}}}.\hfill \end{array}$$(38)

The cooling function is adopted after Athay (1986). The Athay function is frequently used in the prominence literature and may provide a more realistic representation of the cooling in the densest parts of the thread than other optically-thin cooling functions, which may overestimate radiative losses (see, e.g., Soler et al. 2012). Hence, we write the cooling function as:

$$\begin{array}{c}\hfill L(\rho ,T)=f_{\mathrm{p}}(T)\frac{{\rho }^2{T}^2}{m_{\mathrm{p}}^2},\end{array}$$(39)

where m_p is the proton mass and f_p(T) is an analytic function of the temperature given in Athay (1986), which in MKS units is:

$$\begin{array}{cc}& f_{\mathrm{p}}(T)={10}^{-35}{T}^{-2}\{0.4exp[-30{({log}_{10}T-4.6)}^2]\hfill \\ \hfill & +4exp[-20{({log}_{10}T-4.9)}^2]\hfill \\ \hfill & +4.5exp[-16{({log}_{10}T-5.35)}^2]\hfill \\ \hfill & +2exp[-4{({log}_{10}T-6.1)}^2]\}.\hfill \end{array}$$(40)

To compare with the wave heating term, we overplot in Fig. 7 the thermal conduction term and the cooling term in absolute values. We find that radiative cooling is the largest term in the hot, evacuated part of the thread, whereas wave heating dominates in the central, cool region. Near the ends of the thread, although radiative cooling is still dominant, thermal conduction reaches its largest contribution. In the central part, thermal conduction can become more important than radiation because of the enhancing effect of the neutral conductivity. The radiative cooling and the thermal conduction terms depend on the temperature profile. When comparing the LTE and non-LTE cases we do not see important differences in the totally ionized region, where the temperature profiles are almost identical in the LTE and non-LTE cases. We note that owing to the influence of the density profile, the maximum of the cooling function is not located at the ends of the thread, but at both sides of the densest region where the conditions are equivalent to those of prominence-corona transition region (PCTR). Because of thermal conduction, there is a flow of energy from the hot part of the thread towards the cool central part. We find that the thermal conduction term changes sign within the cool region at some specific locations. There, this term actually behaves as a heating term.

In a situation of thermal equilibrium, ℒ = 0. To determine how far we are from thermal equilibrium, we plot the complete heat-loss function in Fig. 10. In the hot part of the thread, the heat-loss function has almost the same form for both LTE and non-LTE cases and for the three boundary conditions due to the dominance of radiative cooling. The heat-loss function is negative along the whole thread (again due to the radiative cooling) except at the central, densest zone. The minimum value of the heat-loss function is ℒ ≈ −3  ×  10⁻⁴ W m⁻³ and is located at z ≈ ±20 Mm, which corresponds to the location where radiative cooling is most efficient. This is located within the PCTR of the thread. Since we have considered no heating due to incident coronal radiation in this simple estimation, cooling cannot be compensated for in the PCTR and coronal part of the thread. This is a logical and expected result.

Fig. 10. Heat-loss function along the thread for the LTE method: (a) complete domain and (b) close-up view of the central part of the thread; (c) and (d) same as (a) and (b) but for the non-LTE method.

Conversely, in the central, densest part of the thread, corresponding to conditions of prominence core, the wave heating function reaches its highest values, while there is much less cooling than in the evacuated parts. As a consequence, the cooling can be compensated for. There is even a slight excess of heating around the thread center, so that ℒ becomes positive in that region. In the LTE case, we find ℒ ∼ 10⁻⁶ W m⁻³, while in the non-LTE case, an even smaller value of ℒ ∼ 10⁻⁷ W m⁻³ is present. Therefore, in the thread core the energy balance situation is very close to that of thermal equilibrium, especially in the non-LTE case. The effect of the different boundary conditions is minor in this regard. We finally note that the heating rates calculated here depend on the injected energy flux by the driver. As explained in Sect. 3.2, we considered an injected flux of 10⁻² W m⁻² that is consistent with previous results on the transmission of photospherically-driven Alfvén through the chromosphere (Soler et al. 2019). If the injected flux is increased or decreased, the resulting heating rates would be larger or smaller. This would consequently affect the energy balance estimation.

5. Concluding remarks

In this exploratory paper, we study the prospects of Alfvén wave dissipation as a heating mechanism in thin threads of solar prominences. The present work builds upon the previous paper by Soler et al. (2016), where a more rudimentary prominence model was used. Here, we consider a 1D model for a thin prominence thread, which is still simple but incorporates several important ingredients that were missing in Soler et al. (2016). For instance, the model is fully non-uniform so that the density, temperature, ionization degree, etc. vary along the thread. In addition, a more consistent representation of the Alfvén wave driver is considered here, including an observationally-motivated spectral weighting function.

The magnetic field and velocity perturbations associated with the Alfvén waves have been obtained by solving the linearized MHD equations for a partially ionized plasma in the single-fluid approximation. Assuming an injected energy flux consistent with the studies of the photospherically-driven Alfvén waves that are transmitted through the chromosphere (see Soler et al. 2019), the obtained amplitudes for the Alfvénic perturbations in the thread are similar to those reported in the observation of prominence oscillations (Arregui et al. 2018). Moreover, the perturbations amplitudes are much smaller than the background values, which justifies the linearization adopted at the beginning of the study. During the solution of the linearized MHD equations, we considered three different boundary conditions, namely total reflection, partial reflection, and total transmission, at the right end of the thread. However, the results are very similar in the three cases, which highlights the robustness of the conclusions.

Ambipolar diffusion due to partial ionization of the plasma is a dissipation mechanism for Alfvén waves. The efficiency of this process depends on the ionization degree. An important part of the paper has been devoted to compare the results obtained when the temperature and ionization profiles are computed considering either LTE or non-LTE approximations for the hydrogen ionization. In general, the results are very similar in the two cases, although a somewhat larger efficiency of wave dissipation is found in the LTE case. The reason for this result is that, for core prominence temperatures, the ionization degree computed in the LTE case is much lower than that computed in the non-LTE approximation. The consequence of this is an enhanced efficiency of ambipolar diffusion.

Alfvén wave heating in the model is caused by Cowling’s diffusion, that is, the joint effect of Ohm’s and ambipolar diffusion. The computed Alfvén wave heating is largest at the thread center, where the density is highest and the temperature and ionization degree are lowest. At this spot, ambipolar diffusion is dominant. Conversely, in the hot coronal part of the model, the plasma is fully ionized and the Alfvén wave heating is negligible. In the fully ionized part, ambipolar diffusion is not present and the dissipation of wave energy is owing to Ohm’s diffusion alone. For coronal conditions, the coefficient of Ohm’s diffusion is so small that no appreciable wave heating is obtained.

In order to determine the relevance of wave heating, we estimated the energy balance by taking into account radiative cooling and thermal conduction in addition to wave heating. We assumed the cooling function after Athay (1986) and the thermal conduction coefficient, including the contributions of electrons and neutrals. We found that in the coronal and PCTR parts, the model is far from being in thermal balance because of the predominance of cooling. We have not included a coronal heating function in this simple estimation, so that this result was expected because there is no mechanism in the corona that could counterbalance the lost of energy due to radiation. On the contrary, a situation of almost equilibrium is achieved in the central part of the thread, especially in the non-LTE case. There is even a slight excess of heating around the thread center. So, in the present case, the efficiency of wave heating in prominence conditions is found to be larger than that anticipated by Soler et al. (2016) in an over-simplified model. This last result encourages us to further investigate the role of wave heating in the prominence energy balance using more realistic models.

Another relevant conclusion that can be drawn is that the ad hoc density profile assumed in this work should evolve towards a new profile that would satisfy the energy balance if the background model was allowed to vary with time. Then, a new wave heating function could be computed for the new density profile, which would in turn drive the system towards a different equilibrium. This process could iteratively be repeated all over again until convergence towards a self-consistent model would be achieved. The process would be very similar to that of Ofman et al. (1998) in the case of heating by resonant absorption. This interesting study would be a natural follow-up to the present paper and may be tackled in a forthcoming work.

Acknowledgments

We acknowledge the support from grant AYA2017-85465-P (MINECO/AEI/FEDER, UE).

References

Appendix A: Standing modes

Here we briefly study the standing modes of the thread, which are used to explain the existence of the resonances with the propagating waves. We use Eqs. (11) and (12) with the temporal dependence specified as exp(−iωt). The equations become

$$\begin{array}{cc}& \omega v_{\perp }=\frac{\mathit{iB}}{\rho {\mathrm{\mu }}_0}\frac{\partial B_{\perp }}{\partial z},\hfill \end{array}$$(A.1)

$$\begin{array}{cc}\hfill & \omega B_{\perp }=iB\frac{\partial v_{\perp }}{\partial z}+i\frac{\partial }{\partial z}\left({\eta }_C\frac{\partial B_{\perp }}{\partial z}\right).\hfill \end{array}$$(A.2)

Equations (A.1) and (A.2) form an eigenvalue problem, where ω is the eigenvalue and B_⊥ and v_⊥ are the eigenfunctions. The eigenvalue problem is numerically solved after imposing the boundary conditions

$$\begin{array}{cc}& v_{\perp }=0\mathrm{at}z=\pm \frac{L}{2},\hfill \end{array}$$(A.3)

$$\begin{array}{cc}\hfill & \frac{\partial B_{\perp }}{\partial z}=0\mathrm{at}z=\pm \frac{L}{2}.\hfill \end{array}$$(A.4)

The eigenfunctions and eigenfrequencies are obtained for both LTE and non-LTE cases.

Figure A.1 shows the eigenfunctions of the fundamental mode and the first four harmonics in the non-LTE case. Both perturbations have a clear harmonic-like behavior with a spatially-dependent amplitude that is determined by the Alfvén speed profile along the thread. The amplitude of the magnetic field perturbations are larger in the densest part of the thread. On the contrary, except for the fundamental mode the velocity perturbations tend to have their maximums located further away from the center, especially as the mode number increases. The shapes of the eigenfunctions in the LTE case are almost identical to those of the non-LTE case and are not shown here for simplicity.

Fig. A.1. Absolute value of the eigenfunctions of the fundamental mode and the first 4 harmonics of standing modes in the non-LTE case: (a) magnetic field perturbation and (b) velocity perturbation. The eigenfunctions are normalized to the maximum of the velocity perturbation in all cases. The eigenfunctions in the LTE case are indistinguishable from those displayed here.

Following Joarder & Roberts (1992) and Oliver et al. (1993), the standing modes can be classified according to the location of the maximum of the velocity perturbation. If the maximum is in the dense prominence region, the mode is called an internal mode. If the maximum is in the evacuated coronal part, the mode is called an external mode. If the velocity is large in both regions, the mode is called a hybrid or string mode. According to Fig. A.1b, the only mode that could be considered as an internal or hybrid mode is the fundamental one. All overtones should be classified as external modes. However, we note that the classification of Joarder & Roberts (1992) and Oliver et al. (1993) applies to a slab model with a piecewise constant density, while here we use a continuous density profile.

The eigenfrequencies for the first 10 modes are shown in Fig. A.2a, where f = ω_R/2π with ω_R the real part of the eigenfrequency. The eigenfrequencies show an almost linear dependence with the mode number and almost the same result is obtained for both LTE and non-LTE cases. This can be explained by the lack of influence of Cowling’s diffusion on the real part of the eigenfrequencies. The different profiles of Cowling’s coefficient is the distict feature between both cases (Fig. 1f).

Fig. A.2. (a) Real part of the eigenfrequency as function of the mode number for the first 10 modes. (b) Damping per period as a function of the eigenfrequency for the same modes as in (a).

The damping per period measures the efficiency of the damping and is defined as D_p = ω_R/2πω_I, where ω_I is the imaginary part of the eigenfrequency. Figure A.2b shows the damping per period of the first 10 modes as a function of their corresponding f. In both LTE and non-LTE cases, the damping per period has extremely large values. The damping per period decreases as the eigenfrequency increases, but we should consider very high harmonics for the damping to be significant. The damping per period is larger in the non-LTE case than in the LTE case, which is again related to the relative values of the Cowling’s coefficient in the two cases.

The result that the damping of the standing modes is very inefficient implies that the heating associated with the dissipation of their energy is unimportant for all practical purposes. If standing modes are resonantly excited by the propagating waves, as the results of the paper evidence, part of the energy of the propagating waves should be transferred to the standing modes. However, since the standing modes are practically undamped, their contribution to the plasma heating would be negligible. This explains why the heating rates computed in the paper have relative minimums at the resonance frequencies.

All Tables

All Figures

	Fig. 1. Background magnitudes along the thread; (a) adopted density profile; (b) computed temperature; (c) computed ionization fraction; (d) Ohm’s diffusion coefficient; (e) ambipolar diffusion coefficient; (f) Cowling’s diffusion coefficient. Results for the LTE and non-LTE methods in panels b–f.
In the text

	Fig. 2. Ionization fraction, ξi, as a function of temperature, T, for the LTE and non-LTE methods.
In the text

	Fig. 3. Real part of the (a) magnetic field and (b) velocity perturbations along the thread for the non-LTE case with ε = 0, 0.5, and 1; (c) and (d) are same as (a) and (b) but with ε = 0 and comparing the results for LTE and non-LTE cases.
In the text

	Fig. 4. Reflectivity (a); transmissivity (b); and absorption as functions of the frequency (c) in log scale for ε = 0, 0.5, and 1 in the non-LTE case. The purple dotted line in panel a corresponds to the sum of reflectivity and transmissivity for ε = 0.5.
In the text

	Fig. 5. Transmissivity coefficient in the non-LTE case for ε = 1 in the range of frequencies corresponding to the first ten standing modes of the thread (vertical black lines).
In the text

	Fig. 6. Same as Fig. 4 but for a comparison between LTE and non-LTE cases with ε = 0.5.
In the text

	Fig. 7. Wave heating function (for ε = 0, 0.5, and 1), cooling function, and thermal conduction term along the thread for the (a) LTE and (b) non-LTE cases.
In the text

	Fig. 8. Integrated heating rate along the thread as a function of the wave frequency in the non-LTE case with ε = 0, 0.5, and 1.
In the text

	Fig. 9. Integrated heating rate along the thread as a function of the wave frequency in the non-LTE case with ε = 0, 0.5, and 1. Same as Fig. 8 but in the range of frequencies corresponding to the first ten standing modes of the thread (vertical black lines).
In the text

	Fig. 10. Heat-loss function along the thread for the LTE method: (a) complete domain and (b) close-up view of the central part of the thread; (c) and (d) same as (a) and (b) but for the non-LTE method.
In the text

	Fig. A.1. Absolute value of the eigenfunctions of the fundamental mode and the first 4 harmonics of standing modes in the non-LTE case: (a) magnetic field perturbation and (b) velocity perturbation. The eigenfunctions are normalized to the maximum of the velocity perturbation in all cases. The eigenfunctions in the LTE case are indistinguishable from those displayed here.
In the text

	Fig. A.2. (a) Real part of the eigenfrequency as function of the mode number for the first 10 modes. (b) Damping per period as a function of the eigenfrequency for the same modes as in (a).
In the text