© ESO, 2013

1. Introduction

Observations show that magnetohydrodynamic (MHD) waves are ubiquitous in the solar atmosphere. Waves are routinely observed in the corona (e.g.., Tomczyk et al. 2007; Tomczyk & McIntosh 2009; McIntosh et al. 2011), the chromosphere (e.g.., Kukhianidze et al. 2006; De Pontieu et al. 2007; Zaqarashvili et al. 2007; He et al. 2009; Okamoto & De Pontieu 2011), and the photosphere (e.g., Jess et al. 2009; Fujimura & Tsuneta 2009). In addition, waves and oscillations in thin threads of solar prominences have been reported (e.g., Okamoto et al. 2007; Lin et al. 2007, 2009; Ning et al. 2009). These recent observations have motivated a number of theoretical works aiming at understanding the properties of the observed waves using the MHD theory (e.g., Van Doorsselaere et al. 2008; Pascoe et al. 2010, 2012; Terradas et al. 2010; Verth et al. 2010; Soler et al. 2011a,b, 2012a). It is believed that the observed waves may play an important role for the heating of the atmospheric plasma (see, e.g., Erdélyi & Fedun 2007; De Pontieu et al. 2007; McIntosh et al. 2011; Cargill & de Moortel 2011).

The study of waves in the solar corona is usually based on the linear ideal MHD theory for a fully ionized plasma (see Priest 1984; Goedbloed & Poedts 2004). In this context, the paper by Edwin & Roberts (1983) on wave propagation in a magnetic flux tube has been used as a basic reference for subsequent works in this field (see also the papers by, e.g., Wentzel 1979; Cally 1986; Goossens et al. 2009, 2012). The assumption of full ionization is adequate for the coronal plasma. However, the lower temperatures of the chromosphere and the photosphere cause this assumption to become invalid. Therefore, it is necessary to determine the impact of partial ionization on the wave modes discussed in Edwin & Roberts (1983) to perform a realistic theoretical modeling of wave propagation in partially ionized magnetic wave guides.

There are many papers that have studied the effect of neutral-ion collision on MHD waves (see, e.g., Braginskii 1965; De Pontieu et al. 2001; Khodachenko et al. 2004, 2006; Leake et al. 2005; Forteza et al. 2007). Since the present paper deals with wave propagation in a magnetic cylinder, we discuss here those works that studied waves in structured media. In solar plasmas the typical frequency of the observed waves is lower than the expected value of the neutral-ion collision frequency. For this reason the single-fluid theory is usually adopted as a first approximation. The single-fluid approximation assumes a strong coupling between ions and neutrals, so that both fluids behave in practice as a single fluid. Soler et al. (2009a) studied wave propagation in a partially ionized flux tube in the single-fluid approximation. The results of Soler et al. (2009a) show that the waves discussed in Edwin & Roberts (1983) are damped due to neutral-ion collisions. Although partial ionization may be relevant for plasma heating (see Khomenko & Collados 2012), the damping by neutral-ion collisions is in general too weak to explain the observed rapid attenuation. In addition, Soler et al. (2009a) found critical values for the longitudinal wavelength that constrains wave propagation. However, Zaqarashvili et al. (2012) have pointed out that these critical wavelengths may not be physical and may be an artifact of the single-fluid approximation. This is so because the single-fluid approximation misses important effects when short length scales and/or high frequencies are involved. In these cases, the general multifluid description is a more suitable approach (see, e.g., Zaqarashvili et al. 2011). In a multifluid description no restriction is imposed on the relative values of the wave frequency and the neutral-ion collision frequency. In this more general treatment the mathematical complexity of the equations is substantially increased compared to the single-fluid case.

A particular form of the multifluid theory is the two-fluid theory in which ions and electrons are considered together as an ion-electron fluid, i.e., the plasma, while neutrals form another fluid that interacts with the plasma by means of collisions. This approach was followed by Kumar & Roberts (2003), who studied surface wave propagation in a Cartesian magnetic interface. In the study by Kumar & Roberts (2003) the gas pressure of neutrals was neglected, meaning that the dynamics of neutrals was only governed by the friction force with ions. Recently, Soler et al. (2012a) used the two-fluid theory to study resonant Alfvén waves in flux tubes. Soler et al. (2012a) focused on transverse kink waves and, consequently, they also neglected gas pressure.

In the present paper we go beyond these previous studies. We derive the general dispersion relation for waves in cylindrical flux tubes in the two-fluid theory. To do so, we consider a consistent description of the neutral fluid dynamics that includes the effect of neutral gas pressure. This is necessary for a realistic description of magnetoacoustic waves. The use of cylindrical geometry and the consideration of neutral pressure are two significant improvements with respect to the Cartesian model of Kumar & Roberts (2003). The dispersion relation derived here is the two-fluid generalization of the well-known dispersion relation of Edwin & Roberts (1983) for a fully ionized single-fluid plasma. The collision frequency between ions and neutrals and the ionization degree are two important parameters of the model. The impact of these two parameters on the waves discussed by Edwin & Roberts (1983) in the fully ionized case is investigated.

This paper is organized as follows. Section 2 contains a description of the equilibrium configuration and the basic equations. In Sect. 3 we follow a normal mode analysis and derive the general dispersion relation for the wave modes. The case in which neutral gas pressure is neglected is explored in Sect. 4, both analytically and numerically. We compare our results with the previous findings of Kumar & Roberts (2003). Then, we incorporate the effect of neutral gas pressure in Sect. 5. Finally, we discuss the implications of our results and give our main conclusions in Sect. 6.

2. Model and two-fluid equations

We study waves in a partially ionized medium composed of ions, electrons, and neutrals. We use the two-fluid theory in which ions and electrons are considered together as an ion-electron fluid, while neutrals form another fluid that interacts with the ion-electron fluid by means of collisions (see, e.g., Zaqarashvili et al. 2011; Soler et al. 2012a,b; Díaz et al. 2012). Throughout, the subscripts “ie” and “n” refer to ion-electrons and neutrals.

We assume a cylindrically symmetric equilibrium and, for convenience, we use cylindrical coordinates, namely r, ϕ, and z for the radial, azimuthal, and longitudinal coordinates. The equilibrium is made of a straight magnetic flux tube of radius R embedded in an unbounded environment. The equilibrium magnetic field is straight and constant along the z-direction, namely B = B   ẑ. Gravity is ignored. As a consequence of the force balance condition, the gradients of the equilibrium gas pressures of the different species are zero. We also assume a static equilibrium so that there are no equilibrium flows. We restrict ourselves to the study of linear perturbations superimposed on the equilibrium state. Hence, the governing equations are linearized (see the general expressions in, e.g., Zaqarashvili et al. 2011). We adopt some additional simplifications that enable us to tackle the problem of wave propagation analytically. We neglect collisions of electrons with neutrals because of the low momentum of electrons. We also drop from the equations the nonadiabatic mechanisms and the magnetic diffusion terms since our purpose here is to determine the impact of neutral-ion collisions only. Thus, the set of coupled differential equations governing linear perturbations from the equilibrium state are $\begin{array}{ccc}{\mathit{\rho }}_{\mathrm{i}}\frac{\mathit{\partial }{{v}}_{\mathrm{i}}}{\mathit{\partial t}}& \mathrm{=}& \mathrm{-}\mathrm{\nabla }{\mathit{p}}_{\mathrm{ie}}\mathrm{+}\frac{\mathrm{1}}{\mathit{\mu }}\left(\mathrm{\nabla }\mathrm{\times }{b}\right)\mathrm{\times }{B}\mathrm{-}{\mathit{\rho }}_{\mathrm{n}}{\mathit{\nu }}_{\mathrm{ni}}\mathrm{(}{{v}}_{\mathrm{i}}\mathrm{-}{{v}}_{\mathrm{n}}\mathrm{)}\mathit{,}\\ {\mathit{\rho }}_{\mathrm{n}}\frac{\mathit{\partial }{{v}}_{\mathrm{n}}}{\mathit{\partial t}}& \mathrm{=}& \mathrm{-}\mathrm{\nabla }{\mathit{p}}_{\mathrm{n}}\mathrm{-}{\mathit{\rho }}_{\mathrm{n}}{\mathit{\nu }}_{\mathrm{ni}}\mathrm{(}{{v}}_{\mathrm{n}}\mathrm{-}{{v}}_{\mathrm{i}}\mathrm{)}\mathit{,}\\ \frac{\mathit{\partial }{b}}{\mathit{\partial t}}& \mathrm{=}& \mathrm{\nabla }\mathrm{\times }\mathrm{(}{{v}}_{\mathrm{i}}\mathrm{\times }{{B}}^{\mathrm{)}}\mathit{,}\\ \frac{\mathit{\partial }{\mathit{p}}_{\mathrm{ie}}}{\mathit{\partial t}}& \mathrm{=}& \mathrm{-}\mathit{\gamma }{\mathit{P}}_{\mathrm{ie}}\mathrm{\nabla }\mathrm{·}{{v}}_{\mathrm{i}}\mathit{,}\\ \frac{\mathit{\partial }{\mathit{p}}_{\mathrm{n}}}{\mathit{\partial t}}& \mathrm{=}& \mathrm{-}\mathit{\gamma }{\mathit{P}}_{\mathrm{n}}\mathrm{\nabla }\mathrm{·}{{v}}_{\mathrm{n}}\mathit{,}\end{array}$where v_i and v_n are the velocities of ions and neutrals, p_ie and p_n are the pressure perturbations of ion-electrons and neutrals, b is the magnetic field perturbation, ρ_i and ρ_n are the equilibrium densities of ions and neutrals, P_ie and P_n are the equilibrium gas pressures of ion-electrons and neutrals, μ is the magnetic permittivity, γ is the adiabatic index, and ν_ni is the neutral-ion collision frequency. For our following analysis it is useful to define the ionization degree as χ = ρ_n/ρ_i. This parameter ranges from χ = 0 for a fully ionized plasma to χ → ∞ for a neutral gas. Also for convenience, in the following expressions we replace the velocities of ions and neutrals by their corresponding Lagrangian displacements, ξ_i and ξ_n, given by ${{v}}_{\mathrm{i}}\mathrm{=}\frac{\mathit{\partial }{\xi }\mathrm{i}}{\mathit{\partial t}}\mathit{,} {{v}}_{\mathrm{n}}\mathrm{=}\frac{\mathit{\partial }{\xi }\mathrm{n}}{\mathit{\partial t}}\mathrm{·}$(6)

3. Normal modes

We assume that the equilibrium densities ρ_i and ρ_n are functions of r alone, so that the equilibrium is uniform in both azimuthal and longitudinal directions. The temperatures of ion-electrons and neutrals vary with r accordingly to keep the equilibrium gas pressures of both fluids constant. Hence we can write the perturbed quantities proportional to exp(imϕ + ik_zz), where m and k_z are the azimuthal and longitudinal wavenumbers, respectively. We express the temporal dependence of the perturbations as exp(−iωt), with ω the frequency. Now we combine Eqs. (1)–(5) and arrive at four coupled equations for the radial components of the Lagrangian displacements of ions and neutrals, ξ_r,ie and ξ_r,n, the ion-electron total pressure Eulerian perturbation, P′ = p_ie + Bb_z/μ, and the neutral pressure Eulerian perturbation, p_n, namely $\begin{array}{ccc}\frac{\mathit{\partial }{\mathit{P}}^{\mathrm{\prime }}}{\mathit{\partial r}}& \mathrm{=}& \end{array}$(7)$\begin{array}{ccc}& & \frac{\mathit{\partial }{\mathit{p}}_{\mathrm{n}}}{\mathit{\partial r}}\mathrm{=}\mathrm{-}\mathrm{i}{\mathit{\rho }}_{\mathrm{n}}{\mathit{\nu }}_{\mathrm{ni}}\mathit{\omega }{\mathit{\xi }}_{\mathit{r,}\mathrm{ie}}\mathrm{+}{\mathit{\rho }}_{\mathrm{n}}\mathit{\omega }\left(\mathit{\omega }\mathrm{+}\mathrm{i}{\mathit{\nu }}_{\mathrm{ni}}\right){\mathit{\xi }}_{\mathit{r,}\mathrm{n}}\mathit{,}\\ & & \mathrm{𝒟}\frac{\mathrm{1}}{\mathit{r}}\frac{\mathit{\partial }\left(\mathit{r}{\mathit{\xi }}_{\mathit{r,}\mathrm{ie}}\right)}{\mathit{\partial r}}\mathrm{=}{\mathrm{𝒞}}_{\mathrm{1}}{\mathit{P}}^{\mathrm{\prime }}\mathrm{+}\mathrm{i}\frac{{\mathit{\nu }}_{\mathrm{ni}}}{\mathit{\omega }\mathrm{+}\mathrm{i}{\mathit{\nu }}_{\mathrm{ni}}}{\mathrm{𝒞}}_{\mathrm{2}}{\mathit{p}}_{\mathrm{n}}\mathit{,}\\ & & {\mathit{\rho }}_{\mathrm{n}}{\mathit{c}}_{\mathrm{n}}^{\mathrm{2}}\mathrm{𝒟}\frac{\mathrm{1}}{\mathit{r}}\frac{\mathit{\partial }\left(\mathit{r}{\mathit{\xi }}_{\mathit{r,}\mathrm{n}}\right)}{\mathit{\partial r}}\mathrm{=}\mathrm{i}\frac{{\mathit{\nu }}_{\mathrm{ni}}}{\mathit{\omega }\mathrm{+}\mathrm{i}{\mathit{\nu }}_{\mathrm{ni}}}{\mathit{\rho }}_{\mathrm{n}}{\mathit{c}}_{\mathrm{n}}^{\mathrm{2}}{\mathrm{𝒞}}_{\mathrm{2}}{\mathit{P}}^{\mathrm{\prime }}\mathrm{+}{\mathrm{𝒞}}_{\mathrm{3}}{\mathit{p}}_{\mathrm{n}}\mathit{,}\end{array}$with the coefficients , , , and defined as $\begin{array}{ccc}\mathrm{𝒟}& \mathrm{=}& {\mathit{\rho }}_{\mathrm{i}}\mathrm{(}{\mathit{c}}_{\mathrm{A}}^{\mathrm{2}}\mathrm{+}{{\mathit{c}}_{\mathrm{ie}}^{\mathrm{2}}}^{\mathrm{)}}\mathrm{(}\stackrel{}{{{\mathit{\omega }}^{\mathrm{2}}}_{\mathrm{˜}}}\mathrm{-}{{\mathit{\omega }}_{\mathrm{A}}^{\mathrm{2}}}^{\mathrm{)}}\mathrm{(}\stackrel{}{{{\mathit{\omega }}^{\mathrm{2}}}_{\mathrm{˜}}}\mathrm{-}{{\mathit{\omega }}_{\mathrm{c}}^{\mathrm{2}}}^{\mathrm{)}}\mathit{,}\\ {\mathrm{𝒞}}_{\mathrm{1}}& \mathrm{=}& \frac{{\mathit{m}}^{\mathrm{2}}}{{\mathit{r}}^{\mathrm{2}}}\mathrm{(}{\mathit{c}}_{\mathrm{A}}^{\mathrm{2}}\mathrm{+}{{\mathit{c}}_{\mathrm{ie}}^{\mathrm{2}}}^{\mathrm{)}}\mathrm{(}\stackrel{}{{{\mathit{\omega }}^{\mathrm{2}}}_{\mathrm{˜}}}\mathrm{-}{{\mathit{\omega }}_{\mathrm{c}}^{\mathrm{2}}}^{\mathrm{)}}\\ & & \mathrm{-}\mathrm{(}\stackrel{}{{{\mathit{\omega }}^{\mathrm{2}}}_{\mathrm{˜}}}\mathrm{-}{{\mathit{\omega }}_{\mathrm{A}}^{\mathrm{2}}}^{\mathrm{)}}\mathrm{(}\stackrel{}{{{\mathit{\omega }}^{\mathrm{2}}}_{\mathrm{˜}}}\mathrm{-}{\mathit{k}}_{\mathit{z}}^{\mathrm{2}}{{\mathit{c}}_{\mathrm{ie}}^{\mathrm{2}}}^{\mathrm{)}}\mathit{,}\\ {\mathrm{𝒞}}_{\mathrm{2}}& \mathrm{=}& {\mathit{k}}_{\mathit{z}}^{\mathrm{2}}{\mathit{c}}_{\mathrm{ie}}^{\mathrm{2}}\mathrm{(}\stackrel{}{{{\mathit{\omega }}^{\mathrm{2}}}_{\mathrm{˜}}}\mathrm{-}{{\mathit{\omega }}_{\mathrm{A}}^{\mathrm{2}}}^{\mathrm{)}}\mathrm{+}\frac{{\mathit{m}}^{\mathrm{2}}}{{\mathit{r}}^{\mathrm{2}}}\mathrm{(}{\mathit{c}}_{\mathrm{A}}^{\mathrm{2}}\mathrm{+}{{\mathit{c}}_{\mathrm{ie}}^{\mathrm{2}}}^{\mathrm{)}}\mathrm{(}\stackrel{}{{{\mathit{\omega }}^{\mathrm{2}}}_{\mathrm{˜}}}\mathrm{-}{{\mathit{\omega }}_{\mathrm{c}}^{\mathrm{2}}}^{\mathrm{)}}\mathit{,}\\ {\mathrm{𝒞}}_{\mathrm{3}}& \mathrm{=}& \mathrm{-}\{\frac{{\mathit{\nu }}_{\mathrm{ni}}^{\mathrm{2}}}{{\mathrm{(}\mathit{\omega }\mathrm{+}\mathit{i}{\mathit{\nu }}_{\mathrm{ni}}\mathrm{)}}^{\mathrm{2}}}{\mathit{\rho }}_{\mathrm{n}}{\mathit{c}}_{\mathrm{n}}^{\mathrm{2}}\mathrm{[}{\mathrm{𝒞}}_{\mathrm{2}}\mathrm{+}{\mathit{\omega }}_{\mathrm{A}}^{\mathrm{2}}{\mathrm{(}\stackrel{}{{{\mathit{\omega }}^{\mathrm{2}}}_{\mathrm{˜}}}\mathrm{-}{{\mathit{\omega }}_{\mathrm{A}}^{\mathrm{2}}}^{\mathrm{)}}}^{\mathrm{]}}\\ & & \mathrm{+}\frac{\mathrm{𝒟}}{\mathit{\omega }\mathrm{(}\mathit{\omega }\mathrm{+}\mathit{i}{\mathit{\nu }}_{\mathrm{ni}}\mathrm{)}}\left[\mathit{\omega }\mathrm{(}\mathit{\omega }\mathrm{+}\mathit{i}{\mathit{\nu }}_{\mathrm{ni}}\mathrm{)}\mathrm{-}{\mathit{c}}_{\mathrm{n}}^{\mathrm{2}}\left(\frac{{\mathit{m}}^{\mathrm{2}}}{{\mathit{r}}^{\mathrm{2}}}\mathrm{+}{\mathit{k}}_{\mathit{z}}^{\mathrm{2}}\right)\right]\}\mathit{,}\end{array}$where $\stackrel{}{{{\mathit{\omega }}^{\mathrm{2}}}_{\mathrm{˜}}}$, ${\mathit{\omega }}_{\mathrm{A}}^{\mathrm{2}}$, and ${\mathit{\omega }}_{\mathrm{c}}^{\mathrm{2}}$ are the squares of the modified frequency, the Alfvén frequency, and the cusp frequency given by $\begin{array}{ccc}\stackrel{}{{{\mathit{\omega }}^{\mathrm{2}}}_{\mathrm{˜}}}& \mathrm{=}& {\mathit{\omega }}^{\mathrm{2}}\left(\mathrm{1}\mathrm{+}\frac{\mathrm{i}\mathit{\chi }{\mathit{\nu }}_{\mathrm{ni}}}{\mathit{\omega }\mathrm{+}\mathit{i}{\mathit{\nu }}_{\mathrm{ni}}}\right)\mathit{,}\\ {\mathit{\omega }}_{\mathrm{A}}^{\mathrm{2}}& \mathrm{=}& {\mathit{k}}_{\mathit{z}}^{\mathrm{2}}{\mathit{c}}_{\mathrm{A}}^{\mathrm{2}}\mathit{,}\\ {\mathit{\omega }}_{\mathrm{c}}^{\mathrm{2}}& \mathrm{=}& {\mathit{k}}_{\mathit{z}}^{\mathrm{2}}\frac{{\mathit{c}}_{\mathrm{A}}^{\mathrm{2}}{\mathit{c}}_{\mathrm{ie}}^{\mathrm{2}}}{{\mathit{c}}_{\mathrm{A}}^{\mathrm{2}}\mathrm{+}{\mathit{c}}_{\mathrm{ie}}^{\mathrm{2}}}\mathrm{·}\end{array}$In addition, ${\mathit{c}}_{\mathrm{A}}^{\mathrm{2}}$ is the square of the Alfvén velocity and ${\mathit{c}}_{\mathrm{ie}}^{\mathrm{2}}$ and ${\mathit{c}}_{\mathrm{n}}^{\mathrm{2}}$ are the squares of the sound velocities of ion-electrons and neutrals, respectively, computed as ${\mathit{c}}_{\mathrm{A}}^{\mathrm{2}}\mathrm{=}\frac{{\mathit{B}}^{\mathrm{2}}}{\mathit{\mu }{\mathit{\rho }}_{\mathrm{i}}}\mathit{,} {\mathit{c}}_{\mathrm{ie}}^{\mathrm{2}}\mathrm{=}\frac{\mathit{\gamma }{\mathit{P}}_{\mathrm{ie}}}{{\mathit{\rho }}_{\mathrm{i}}}\mathit{,} {\mathit{c}}_{\mathrm{n}}^{\mathrm{2}}\mathrm{=}\frac{\mathit{\gamma }{\mathit{P}}_{\mathrm{n}}}{{\mathit{\rho }}_{\mathrm{n}}}\mathrm{·}$(18)Equations (7)–(10) are valid for all density profiles in the radial direction. For c_ie = c_n = 0 and p_n = 0 they revert to the equations discussed in Soler et al. (2012a). Equations (7)–(10) are singular when . The positions of the singularities are mobile and depend on the frequency, ω. The whole set of frequencies satisfying at some r form two continua of frequencies called the Alfvén and cusp (or slow) continua (see, e.g., Appert et al. 1974). A wave whose frequency is within any of the two continua is damped by resonant absorption (see, e.g., Goossens et al. 2011). Resonantly damped waves in partially ionized magnetic cylinders were studied by Soler et al. (2009b) in the single-fluid approximation and by Soler et al. (2012a) in the two-fluid formalism.

3.1. Piece-wise constant equilibrium

In the present paper we are not interested in studying the resonant behavior of the waves. We refer to Soler et al. (2009b, 2012a) for studies of resonant waves in partially ionized plasmas. Instead, here we focus on the modification of the wave frequencies due to neutral-ion collisions. For this reason, we avoid the presence of the resonances by choosing both ρ_i and ρ_n to be piece-wise constant functions, namely ${\mathit{\rho }}_{\mathrm{i}}\mathrm{\left(}\mathit{r}\mathrm{\right)}\mathrm{=}\{\begin{array}{c}\\ {\mathit{\rho }}_{\mathrm{i}\mathit{,}\mathrm{0}}& \mathrm{if}& \mathit{r}\mathrm{\le }\mathit{R,}\\ {\mathit{\rho }}_{\mathrm{i}\mathit{,}\mathrm{ex}}& \mathrm{if}& \mathit{r}\mathit{>}\mathit{R,}\end{array} {\mathit{\rho }}_{\mathrm{n}}\mathrm{\left(}\mathit{r}\mathrm{\right)}\mathrm{=}\{\begin{array}{c}\\ {\mathit{\rho }}_{\mathrm{n}\mathit{,}\mathrm{0}}& \mathrm{if}& \mathit{r}\mathrm{\le }\mathit{R,}\\ {\mathit{\rho }}_{\mathrm{n}\mathit{,}\mathrm{ex}}& \mathrm{if}& \mathit{r}\mathit{>}\mathit{R,}\end{array}$(19)where R is the radius of the cylinder with constant densities ρ_i,0 and ρ_n,0 embedded in an external environment with constant densities ρ_i,ex and ρ_n,ex. This is the density profile adopted by Edwin & Roberts (1983) in the fully ionized case. Hereafter the indices “0” and “ex” represent the internal and external media. Because the magnetic field and equilibrium pressure of both fluids are constant, the rest of equilibrium quantities are piece-wise constants as well. For simplicity we also take ν_ni as a constant free parameter.

For a piece-wise constant equilibrium we can combine Eqs. (7) and (10) to eliminate ξ_r,ie and ξ_r,n and obtain two coupled equations involving P′ and p_n only, namely $\begin{array}{ccc}\frac{{\mathit{\partial }}^{\mathrm{2}}{\mathit{P}}^{\mathrm{\prime }}}{\mathit{\partial }{\mathit{r}}^{\mathrm{2}}}\mathrm{+}\frac{\mathrm{1}}{\mathit{r}}\frac{\mathit{\partial }{\mathit{P}}^{\mathrm{\prime }}}{\mathit{\partial r}}\mathrm{+}\left({\mathit{k}}_{\mathrm{ie}}^{\mathrm{2}}\mathrm{-}\frac{{\mathit{m}}^{\mathrm{2}}}{{\mathit{r}}^{\mathrm{2}}}\right){\mathit{P}}^{\mathrm{\prime }}& \mathrm{=}& {\mathit{q}}_{\mathrm{1}}^{\mathrm{2}}{\mathit{p}}_{\mathrm{n}}\mathit{,}\\ \frac{{\mathit{\partial }}^{\mathrm{2}}{\mathit{p}}_{\mathrm{n}}}{\mathit{\partial }{\mathit{r}}^{\mathrm{2}}}\mathrm{+}\frac{\mathrm{1}}{\mathit{r}}\frac{\mathit{\partial }{\mathit{p}}_{\mathrm{n}}}{\mathit{\partial r}}\mathrm{+}\left({\mathit{k}}_{\mathrm{n}}^{\mathrm{2}}\mathrm{-}\frac{{\mathit{m}}^{\mathrm{2}}}{{\mathit{r}}^{\mathrm{2}}}\right){\mathit{p}}_{\mathrm{n}}& \mathrm{=}& {\mathit{q}}_{\mathrm{2}}^{\mathrm{2}}{\mathit{P}}^{\mathrm{\prime }}\mathit{,}\end{array}$with $\begin{array}{ccc}{\mathit{k}}_{\mathrm{ie}}^{\mathrm{2}}& \mathrm{=}& \frac{\mathrm{(}\mathit{\omega }\left(\mathit{\omega }\mathrm{+}\mathrm{i}\mathit{\chi }{\mathit{\nu }}_{\mathrm{ni}}\right)\mathrm{-}{{\mathit{\omega }}_{\mathrm{A}}^{\mathrm{2}}}^{\mathrm{)}}}{\mathrm{(}{\mathit{c}}_{\mathrm{A}}^{\mathrm{2}}\mathrm{+}{{\mathit{c}}_{\mathrm{ie}}^{\mathrm{2}}}^{\mathrm{)}}\mathrm{(}\stackrel{}{{{\mathit{\omega }}^{\mathrm{2}}}_{\mathrm{˜}}}\mathrm{-}{{\mathit{\omega }}_{\mathrm{c}}^{\mathrm{2}}}^{\mathrm{)}}}\\ & & \mathrm{\times }\left(\stackrel{}{{{\mathit{\omega }}^{\mathrm{2}}}_{\mathrm{˜}}}\mathrm{-}{\mathit{k}}_{\mathit{z}}^{\mathrm{2}}{\mathit{c}}_{\mathrm{ie}}^{\mathrm{2}}\frac{\stackrel{}{{{\mathit{\omega }}^{\mathrm{2}}}_{\mathrm{˜}}}\mathrm{-}{\mathit{\omega }}_{\mathrm{A}}^{\mathrm{2}}}{\mathit{\omega }\mathrm{(}\mathit{\omega }\mathrm{+}\mathrm{i}\mathit{\chi }{\mathit{\nu }}_{\mathrm{ni}}\mathrm{)}\mathrm{-}{\mathit{\omega }}_{\mathrm{A}}^{\mathrm{2}}}\right)\mathit{,}\\ {\mathit{k}}_{\mathrm{n}}^{\mathrm{2}}& \mathrm{=}& \frac{\mathit{\omega }\mathrm{(}\mathit{\omega }\mathrm{+}\mathit{i}{\mathit{\nu }}_{\mathrm{ni}}\mathrm{)}\mathrm{-}{\mathit{k}}_{\mathit{z}}^{\mathrm{2}}{\mathit{c}}_{\mathrm{n}}^{\mathrm{2}}}{{\mathit{c}}_{\mathrm{n}}^{\mathrm{2}}}\\ & & \mathrm{+}\frac{\mathit{\chi }{\mathit{\nu }}_{\mathrm{ni}}^{\mathrm{2}}}{\mathit{\omega }\mathrm{+}\mathrm{i}{\mathit{\nu }}_{\mathrm{ni}}}\frac{\mathit{\omega }{\mathit{\omega }}_{\mathrm{A}}^{\mathrm{2}}}{\mathrm{(}{\mathit{c}}_{\mathrm{A}}^{\mathrm{2}}\mathrm{+}{{\mathit{c}}_{\mathrm{ie}}^{\mathrm{2}}}^{\mathrm{)}}\mathrm{(}\stackrel{}{{{\mathit{\omega }}^{\mathrm{2}}}_{\mathrm{˜}}}\mathrm{-}{{\mathit{\omega }}_{\mathrm{c}}^{\mathrm{2}}}^{\mathrm{)}}}\mathit{,}\\ {\mathit{q}}_{\mathrm{1}}^{\mathrm{2}}& \mathrm{=}& \mathrm{i}\frac{{\mathit{\nu }}_{\mathrm{ni}}}{\mathit{\omega }\mathrm{+}\mathrm{i}{\mathit{\nu }}_{\mathrm{ni}}}[\frac{\mathit{\omega }\left(\mathit{\omega }\mathrm{+}\mathrm{i}{\mathit{\nu }}_{\mathrm{ni}}\right)\mathrm{-}{\mathit{k}}_{\mathit{z}}^{\mathrm{2}}{\mathit{c}}_{\mathrm{n}}^{\mathrm{2}}}{{\mathit{c}}_{\mathrm{n}}^{\mathrm{2}}}\mathrm{+}\frac{{\mathit{k}}_{\mathit{z}}^{\mathrm{2}}{\mathit{c}}_{\mathrm{ie}}^{\mathrm{2}}\mathrm{(}\stackrel{}{{{\mathit{\omega }}^{\mathrm{2}}}_{\mathrm{˜}}}\mathrm{-}{{\mathit{\omega }}_{\mathrm{A}}^{\mathrm{2}}}^{\mathrm{)}}}{\mathrm{(}{\mathit{c}}_{\mathrm{A}}^{\mathrm{2}}\mathrm{+}{{\mathit{c}}_{\mathrm{ie}}^{\mathrm{2}}}^{\mathrm{)}}\mathrm{(}\stackrel{}{{{\mathit{\omega }}^{\mathrm{2}}}_{\mathrm{˜}}}\mathrm{-}{{\mathit{\omega }}_{\mathrm{c}}^{\mathrm{2}}}^{\mathrm{)}}}\\ & & \mathrm{+}\frac{\mathit{\chi }{\mathit{\nu }}_{\mathrm{ni}}^{\mathrm{2}}}{\mathit{\omega }\mathrm{+}\mathrm{i}{\mathit{\nu }}_{\mathrm{ni}}}\frac{{\mathit{\omega }}_{\mathrm{A}}^{\mathrm{2}}\mathit{\omega }}{\mathrm{(}{\mathit{c}}_{\mathrm{A}}^{\mathrm{2}}\mathrm{+}{{\mathit{c}}_{\mathrm{ie}}^{\mathrm{2}}}^{\mathrm{)}}\mathrm{(}\stackrel{}{{{\mathit{\omega }}^{\mathrm{2}}}_{\mathrm{˜}}}\mathrm{-}{{\mathit{\omega }}_{\mathrm{c}}^{\mathrm{2}}}^{\mathrm{)}}}]\mathit{,}\\ {\mathit{q}}_{\mathrm{2}}^{\mathrm{2}}& \mathrm{=}& \mathrm{i}\frac{\mathit{\chi }{\mathit{\nu }}_{\mathrm{ni}}\mathit{\omega }\stackrel{}{{{\mathit{\omega }}^{\mathrm{2}}}_{\mathrm{˜}}}}{\mathrm{(}{\mathit{c}}_{\mathrm{A}}^{\mathrm{2}}\mathrm{+}{{\mathit{c}}_{\mathrm{ie}}^{\mathrm{2}}}^{\mathrm{)}}\mathrm{(}\stackrel{}{{{\mathit{\omega }}^{\mathrm{2}}}_{\mathrm{˜}}}\mathrm{-}{{\mathit{\omega }}_{\mathrm{c}}^{\mathrm{2}}}^{\mathrm{)}}}\mathrm{·}\end{array}$Equations (20) and (21) are the basic equations of this investigation. They are two coupled Bessel-type differential equations representing the coupled behavior of ion-electrons and neutrals due to neutral-ion collisions. Note that the terms with q₁ and q₂ are the ones that couple the equations. These terms vanish in fully ionized and fully neutral cases. We need to solve Eqs. (20) and (21) to find the general dispersion relation.

For later use we also give the expressions of ξ_r,ie and ξ_r,n in terms of P′ and p_n, namely $\begin{array}{ccc}{\mathit{\xi }}_{\mathit{r,}\mathrm{ie}}& \mathrm{=}& \frac{\mathrm{1}}{{\mathit{\rho }}_{\mathrm{i}}\mathrm{(}\stackrel{}{{{\mathit{\omega }}^{\mathrm{2}}}_{\mathrm{˜}}}\mathrm{-}{{\mathit{\omega }}_{\mathrm{A}}^{\mathrm{2}}}^{\mathrm{)}}}\left(\frac{\mathit{\partial }{\mathit{P}}^{\mathrm{\prime }}}{\mathit{\partial r}}\mathrm{+}\mathrm{i}\frac{{\mathit{\nu }}_{\mathrm{ni}}}{\mathit{\omega }\mathrm{+}\mathrm{i}{\mathit{\nu }}_{\mathrm{ni}}}\frac{\mathit{\partial }{\mathit{p}}_{\mathrm{n}}}{\mathit{\partial r}}\right)\mathit{,}\\ {\mathit{\xi }}_{\mathit{r,}\mathrm{n}}& \mathrm{=}& \left(\frac{\mathrm{1}}{{\mathit{\rho }}_{\mathrm{n}}\mathit{\omega }\left(\mathit{\omega }\mathrm{+}\mathrm{i}{\mathit{\nu }}_{\mathrm{ni}}\right)}\mathrm{-}\frac{{\mathit{\nu }}_{\mathrm{ni}}^{\mathrm{2}}}{{\left(\mathit{\omega }\mathrm{+}\mathrm{i}{\mathit{\nu }}_{\mathrm{ni}}\right)}^{\mathrm{2}}}\frac{\mathrm{1}}{{\mathit{\rho }}_{\mathrm{i}}\mathrm{(}\stackrel{}{{{\mathit{\omega }}^{\mathrm{2}}}_{\mathrm{˜}}}\mathrm{-}{{\mathit{\omega }}_{\mathrm{A}}^{\mathrm{2}}}^{\mathrm{)}}}\right)\frac{\mathit{\partial }{\mathit{p}}_{\mathrm{n}}}{\mathit{\partial r}}\\ & & \mathrm{+}\mathrm{i}\frac{{\mathit{\nu }}_{\mathrm{ni}}}{\mathit{\omega }\mathrm{+}\mathrm{i}{\mathit{\nu }}_{\mathrm{ni}}}\frac{\mathrm{1}}{{\mathit{\rho }}_{\mathrm{i}}\mathrm{(}\stackrel{}{{{\mathit{\omega }}^{\mathrm{2}}}_{\mathrm{˜}}}\mathrm{-}{{\mathit{\omega }}_{\mathrm{A}}^{\mathrm{2}}}^{\mathrm{)}}}\frac{\mathit{\partial }{\mathit{P}}^{\mathrm{\prime }}}{\mathit{\partial r}}\mathrm{·}\end{array}$Note that for ν_ni = 0 ions and neutrals are decoupled and the magnetic field has no influence on the dynamics of neutrals.

3.2. Strong thermal coupling and characteristic velocities

Up to here, no restriction on the values of the characteristic Alfvén and sound velocities has been made. As a result of the pressure balance condition at r = R, the six characteristic velocities in the equilibrium are related by ${\mathit{\rho }}_{\mathrm{i}\mathit{,}\mathrm{0}}\left(\frac{\mathit{\gamma }}{\mathrm{2}}{\mathit{c}}_{\mathrm{A}\mathit{,}\mathrm{0}}^{\mathrm{2}}\mathrm{+}{\mathit{c}}_{\mathrm{ie}\mathit{,}\mathrm{0}}^{\mathrm{2}}\mathrm{+}\mathit{\chi }{\mathit{c}}_{\mathrm{n}\mathit{,}\mathrm{0}}^{\mathrm{2}}\right)\mathrm{=}{\mathit{\rho }}_{\mathrm{i}\mathit{,}\mathrm{ex}}\left(\frac{\mathit{\gamma }}{\mathrm{2}}{\mathit{c}}_{\mathrm{A}\mathit{,}\mathrm{ex}}^{\mathrm{2}}\mathrm{+}{\mathit{c}}_{\mathrm{ie}\mathit{,}\mathrm{ex}}^{\mathrm{2}}\mathrm{+}\mathit{\chi }{\mathit{c}}_{\mathrm{n}\mathit{,}\mathrm{ex}}^{\mathrm{2}}\right)\mathit{.}$(28)Now for simplicity we assume a strong thermal coupling between the ionized and neutral fluids, so that both fluids have the same temperature locally. However, the temperature can still be different in the internal and external media. Using the ideal gas law for both species, this requirement has the consequence that the local sound velocity of the ionized and neutral fluids are related by ${\mathit{c}}_{\mathrm{ie}}^{\mathrm{2}}\mathrm{=}\mathrm{2}{\mathit{c}}_{\mathrm{n}}^{\mathrm{2}}$. Hence we can define an effective sound velocity of the whole plasma, c_s,eff, as ${\mathit{c}}_{\mathrm{s}\mathit{,}\mathrm{eff}}^{\mathrm{2}}\mathrm{=}\frac{{\mathit{c}}_{\mathrm{ie}}^{\mathrm{2}}\mathrm{+}\mathit{\chi }{\mathit{c}}_{\mathrm{n}}^{\mathrm{2}}}{\mathrm{1}\mathrm{+}\mathit{\chi }}\mathrm{=}\frac{\mathrm{2}\mathrm{+}\mathit{\chi }}{\mathrm{2}\mathrm{(}\mathrm{1}\mathrm{+}\mathit{\chi }\mathrm{)}}{\mathit{c}}_{\mathrm{ie}}^{\mathrm{2}}\mathit{,}$(29)so that the parameter space of characteristic velocities can be reduced from six to four different velocities, namely c_A,0, c_s,eff,0, c_A,ex, and c_s,eff,ex.

3.3. Coupled solutions

We continue the mathematical study of the normal modes. We look for solutions to the coupled Eqs. (20) and (21). In the internal medium, i.e., r ≤ R, physical solutions imply that both P′ and p_n are regular at r = 0. The general solutions of P′ and p_n satisfying this condition in the internal medium are $\begin{array}{ccc}& & {\mathit{P}}_{\mathrm{0}}^{\mathrm{\prime }}\mathrm{=}{\mathit{A}}_{\mathrm{1}}{\mathit{J}}_{\mathit{m}}\left({\mathit{k}}_{\mathrm{1}\mathit{,}\mathrm{0}}\mathit{r}\right)\mathrm{+}{\mathit{A}}_{\mathrm{2}}{\mathit{J}}_{\mathit{m}}\left({\mathit{k}}_{\mathrm{2}\mathit{,}\mathrm{0}}\mathit{r}\right)\mathit{,}\\ & & {\mathit{p}}_{\mathrm{n}\mathit{,}\mathrm{0}}\mathrm{=}\mathrm{-}{\mathit{A}}_{\mathrm{1}}\frac{{\mathit{q}}_{\mathrm{2}\mathit{,}\mathrm{0}}^{\mathrm{2}}}{{\mathit{k}}_{\mathrm{1}\mathit{,}\mathrm{0}}^{\mathrm{2}}\mathrm{-}{\mathit{k}}_{\mathrm{n}\mathit{,}\mathrm{0}}^{\mathrm{2}}}{\mathit{J}}_{\mathit{m}}\left({\mathit{k}}_{\mathrm{1}\mathit{,}\mathrm{0}}\mathit{r}\right)\mathrm{-}{\mathit{A}}_{\mathrm{2}}\frac{{\mathit{q}}_{\mathrm{2}\mathit{,}\mathrm{0}}^{\mathrm{2}}}{{\mathit{k}}_{\mathrm{2}\mathit{,}\mathrm{0}}^{\mathrm{2}}\mathrm{-}{\mathit{k}}_{\mathrm{n}\mathit{,}\mathrm{0}}^{\mathrm{2}}}{\mathit{J}}_{\mathit{m}}\left({\mathit{k}}_{\mathrm{2}\mathit{,}\mathrm{0}}\mathit{r}\right)\mathit{,} \end{array}$where J_m is the Bessel function of the first kind of order m, A₁ and A₂ are constants, and k_1,0 and k_2,0 are the two possible values of the radial wavenumber, k, in the internal medium given by the solution of the following equation $\mathrm{(}{\mathit{k}}^{\mathrm{2}}\mathrm{-}{{\mathit{k}}_{\mathrm{ie}\mathit{,}\mathrm{0}}^{\mathrm{2}}}^{\mathrm{)}}\mathrm{(}{\mathit{k}}^{\mathrm{2}}\mathrm{-}{{\mathit{k}}_{\mathrm{n}\mathit{,}\mathrm{0}}^{\mathrm{2}}}^{\mathrm{)}}\mathrm{-}{\mathit{q}}_{\mathrm{1}\mathit{,}\mathrm{0}}^{\mathrm{2}}{\mathit{q}}_{\mathrm{2}\mathit{,}\mathrm{0}}^{\mathrm{2}}\mathrm{=}\mathrm{0.}$(32)Equivalently, in the external medium, i.e., r > R, the solutions of P′ and p_n are of the form $\begin{array}{ccc}{\mathit{P}}_{\mathrm{ex}}^{\mathrm{\prime }}& \mathrm{=}& {\mathit{A}}_{\mathrm{3}}{\mathit{H}}_{\mathit{m}}^{\mathrm{\left(}\mathrm{1}\mathrm{\right)}}\left({\mathit{k}}_{\mathrm{1}\mathit{,}\mathrm{ex}}\mathit{r}\right)\mathrm{+}{\mathit{A}}_{\mathrm{4}}{\mathit{H}}_{\mathit{m}}^{\mathrm{\left(}\mathrm{1}\mathrm{\right)}}\left({\mathit{k}}_{\mathrm{2}\mathit{,}\mathrm{ex}}\mathit{r}\right)\mathit{,}\\ {\mathit{p}}_{\mathrm{n}\mathit{,}\mathrm{ex}}& \mathrm{=}& \mathrm{-}{\mathit{A}}_{\mathrm{3}}\frac{{\mathit{q}}_{\mathrm{2}\mathit{,}\mathrm{ex}}^{\mathrm{2}}}{{\mathit{k}}_{\mathrm{1}\mathit{,}\mathrm{ex}}^{\mathrm{2}}\mathrm{-}{\mathit{k}}_{\mathrm{n}\mathit{,}\mathrm{ex}}^{\mathrm{2}}}{\mathit{H}}_{\mathit{m}}^{\mathrm{\left(}\mathrm{1}\mathrm{\right)}}\left({\mathit{k}}_{\mathrm{1}\mathit{,}\mathrm{ex}}\mathit{r}\right)\\ & & \mathrm{-}{\mathit{A}}_{\mathrm{4}}\frac{{\mathit{q}}_{\mathrm{2}\mathit{,}\mathrm{ex}}^{\mathrm{2}}}{{\mathit{k}}_{\mathrm{2}\mathit{,}\mathrm{ex}}^{\mathrm{2}}\mathrm{-}{\mathit{k}}_{\mathrm{n}\mathit{,}\mathrm{ex}}^{\mathrm{2}}}{\mathit{H}}_{\mathit{m}}^{\mathrm{\left(}\mathrm{1}\mathrm{\right)}}\left({\mathit{k}}_{\mathrm{2}\mathit{,}\mathrm{ex}}\mathit{r}\right)\mathit{,}\end{array}$where ${\mathit{H}}_{\mathit{m}}^{\mathrm{\left(}\mathrm{1}\mathrm{\right)}}$ is the Hankel function of the first kind of order m, A₃ and A₄ are constants, and k_1,ex and k_2,ex are the two possible values of the radial wavenumber, k, in the external medium also given by the solution of Eq. (32) but now using external values for the parameters. The expressions of P′ and p_n for r > R are written generally in terms of Hankel functions instead of the usual modified Bessel functions, K_m, used to represent trapped waves. Thus, our formalism takes into account the possibility of wave leakage, i.e., wave radiation in the external medium. Although we do not explore leaky waves in the present work, the obtained dispersion relation will also be valid for leaky waves. For ideal trapped waves, the external radial wavenumber is purely imaginary and the function ${\mathit{H}}_{\mathit{m}}^{\mathrm{\left(}\mathrm{1}\mathrm{\right)}}$ consistently reverts to the modified Bessel function of the second kind, K_m, used by Edwin & Roberts (1983). In the partially ionized case the radial wavenumber is a complex quantity and the use of the function ${\mathit{H}}_{\mathit{m}}^{\mathrm{\left(}\mathrm{1}\mathrm{\right)}}$ causes us to be very careful when choosing the branch of the external radial wavenumber. To avoid nonphysical energy propagation coming from infinity, we need to select the appropriate branch of the external radial wavenumber so that the condition for outgoing waves is satisfied (see details in Cally 1986; Stenuit et al. 1999).

To understand why there are two different values of the radial wavenumber, k, it is instructive to assume a weak coupling between the species, i.e., ν_ni ≪ |ω |, so that the quadratic terms in ν_ni can be neglected in Eq. (32). The neglected terms involve the coupling coefficients q₁ and q₂. Then, the two independent wavenumbers in Eq. (32) simplify to ${\mathit{k}}_{\mathrm{1}}\mathrm{\approx }{\mathit{k}}_{\mathrm{ie}}\mathit{,} {\mathit{k}}_{\mathrm{2}}\mathrm{\approx }{\mathit{k}}_{\mathrm{n}}\mathit{,}$(35)where we have dropped the indices “0” or “ex” because the same result is valid in both internal and external media. The two values of k for weak coupling reduce to the wavenumbers of the ionized and neutral fluids, respectively. For high collision frequencies, these two wavenumbers contain terms with q₁ and q₂ that couple the two fluids.

3.4. Dispersion relation

To find the dispersion relation we match the internal solutions to the external solutions by means of appropriate boundary conditions at r = R. The boundary conditions at the interface between two media in the two-fluid formalism are discussed in Díaz et al. (2012). In the absence of gravity, these boundary conditions reduce to the continuity of P′, p_n, ξ_r,ie, and ξ_r,n at r = R. After applying the boundary conditions at r = R we find a system of four algebraic equations for the constants A₁, A₂, A₃, and A₄. The condition that the system has a non-trivial solution provides us with the dispersion relation. The dispersion relation is ${\mathrm{𝒟}}_{\mathit{m}}\left(\mathit{\omega ,}{\mathit{k}}_{\mathit{z}}\right)\mathrm{=}\mathrm{0}\mathit{,}$(36)with the full expression of given in Appendix A. This dispersion relation is the two-fluid generalization of the dispersion relation of Edwin & Roberts (1983) for a fully ionized single-fluid cylindrical flux tube. In the following sections we discuss the corrections due to neutral-ion collisions on the wave modes described by Edwin & Roberts (1983) in the fully ionized case.

4. Pressureless neutral fluid

Before exploring the solutions of the general dispersion relation we first consider the case that neutral pressure is neglected and only the gas pressure of the ionized fluid is taken into account. In the absence of neutral pressure the only force acting on the neutral fluid is the friction force due to neutral-ion collisions. This is the situation studied by Kumar & Roberts (2003) in planar geometry. We set P_n = 0 so that c_n = 0. After some algebraic manipulations the general dispersion relation simplifies to $\frac{{\mathit{k}}_{\mathrm{0}}}{{\mathit{\rho }}_{\mathrm{i}\mathit{,}\mathrm{0}}\mathrm{(}\stackrel{}{{{\mathit{\omega }}^{\mathrm{2}}}_{\mathrm{˜}}}\mathrm{-}{{\mathit{\omega }}_{\mathrm{A}\mathit{,}\mathrm{0}}^{\mathrm{2}}}^{\mathrm{)}}}\frac{{\mathit{J}}_{\mathit{m}}^{\mathrm{\prime }}\mathrm{\left(}{\mathit{k}}_{\mathrm{0}}\mathit{R}\mathrm{\right)}}{{\mathit{J}}_{\mathit{m}}\mathrm{\left(}{\mathit{k}}_{\mathrm{0}}\mathit{R}\mathrm{\right)}}\mathrm{=}\frac{{\mathit{k}}_{\mathrm{ex}}}{{\mathit{\rho }}_{\mathrm{i}\mathit{,}\mathrm{ex}}\mathrm{(}\stackrel{}{{{\mathit{\omega }}^{\mathrm{2}}}_{\mathrm{˜}}}\mathrm{-}{{\mathit{\omega }}_{\mathrm{A}\mathit{,}\mathrm{ex}}^{\mathrm{2}}}^{\mathrm{)}}}\frac{{\mathit{K}}_{\mathit{m}}^{\mathrm{\prime }}\mathrm{\left(}{\mathit{k}}_{\mathrm{ex}}\mathit{R}\mathrm{\right)}}{{\mathit{K}}_{\mathit{m}}\mathrm{\left(}{\mathit{k}}_{\mathrm{ex}}\mathit{R}\mathrm{\right)}}\mathit{,}$(37)

with k₀ and k_ex given by $\begin{array}{ccc}{\mathit{k}}_{\mathrm{0}}^{\mathrm{2}}& \mathrm{=}& \frac{\mathrm{(}\stackrel{}{{{\mathit{\omega }}^{\mathrm{2}}}_{\mathrm{˜}}}\mathrm{-}{{\mathit{\omega }}_{\mathrm{A}\mathit{,}\mathrm{0}}^{\mathrm{2}}}^{\mathrm{)}}\mathrm{(}\stackrel{}{{{\mathit{\omega }}^{\mathrm{2}}}_{\mathrm{˜}}}\mathrm{-}{\mathit{k}}_{\mathit{z}}^{\mathrm{2}}{{\mathit{c}}_{\mathrm{ie}\mathit{,}\mathrm{0}}^{\mathrm{2}}}^{\mathrm{)}}}{\mathrm{(}{\mathit{c}}_{\mathrm{A}\mathit{,}\mathrm{0}}^{\mathrm{2}}\mathrm{+}{{\mathit{c}}_{\mathrm{ie}\mathit{,}\mathrm{0}}^{\mathrm{2}}}^{\mathrm{)}}\mathrm{(}\stackrel{}{{{\mathit{\omega }}^{\mathrm{2}}}_{\mathrm{˜}}}\mathrm{-}{{\mathit{\omega }}_{\mathrm{c}\mathit{,}\mathrm{0}}^{\mathrm{2}}}^{\mathrm{)}}}\mathrm{·}\\ {\mathit{k}}_{\mathrm{ex}}^{\mathrm{2}}& \mathrm{=}& \mathrm{-}\frac{\mathrm{(}\stackrel{}{{{\mathit{\omega }}^{\mathrm{2}}}_{\mathrm{˜}}}\mathrm{-}{{\mathit{\omega }}_{\mathrm{A}\mathit{,}\mathrm{ex}}^{\mathrm{2}}}^{\mathrm{)}}\mathrm{(}\stackrel{}{{{\mathit{\omega }}^{\mathrm{2}}}_{\mathrm{˜}}}\mathrm{-}{\mathit{k}}_{\mathit{z}}^{\mathrm{2}}{{\mathit{c}}_{\mathrm{ie}\mathit{,}\mathrm{ex}}^{\mathrm{2}}}^{\mathrm{)}}}{\mathrm{(}{\mathit{c}}_{\mathrm{A}\mathit{,}\mathrm{ex}}^{\mathrm{2}}\mathrm{+}{{\mathit{c}}_{\mathrm{ie}\mathit{,}\mathrm{ex}}^{\mathrm{2}}}^{\mathrm{)}}\mathrm{(}\stackrel{}{{{\mathit{\omega }}^{\mathrm{2}}}_{\mathrm{˜}}}\mathrm{-}{{\mathit{\omega }}_{\mathrm{c}\mathit{,}\mathrm{ex}}^{\mathrm{2}}}^{\mathrm{)}}}\mathrm{·}\\ & & \end{array}$In Eq. (37) we have used the modified Bessel function of the second kind, K_m, instead of the Hankel function of the first kind, ${\mathit{H}}_{\mathit{m}}^{\mathrm{\left(}\mathrm{1}\mathrm{\right)}}$, in order to write Eq. (37) in the same form as the dispersion relation of Edwin & Roberts (1983). Indeed, the comparison between Eq. (37) and the dispersion relation of Edwin & Roberts (1983) reveals that both equations are exactly the same if we replace $\stackrel{}{{{\mathit{\omega }}^{\mathrm{2}}}_{\mathrm{˜}}}$ by ω² in Eq. (37). We can take advantage of this result to easily study the modifications due to neutral-ion collisions of the solutions of the fully ionized case.

4.1. Analytic approximate study

4.1.1. Standing waves

We focus on standing waves so that we fix the longitudinal wavenumber, k_z, to a real value. Due to temporal damping by neutral-ion collisions ω is complex. i.e., ω = ω_R + iω_I, where ω_R and ω_I are the real and imaginary parts of the frequency. Therefore the wave amplitude is damped in time due through the exponential factor exp(−|ω_I|t).

We assume that ω = ω₀ is a trapped solution of the dispersion relation of Edwin & Roberts (1983) in the fully ionized case. Here we only consider trapped waves so that ω₀ is real. For leaky waves see, e.g., Cally (1986). Then $\stackrel{}{{{\mathit{\omega }}^{\mathrm{2}}}_{\mathrm{˜}}}\mathrm{=}{\mathit{\omega }}_{\mathrm{0}}^{\mathrm{2}}$ is a solution of Eq. (37). Using the expression of $\stackrel{}{{{\mathit{\omega }}^{\mathrm{2}}}_{\mathrm{˜}}}$ (Eq. (15)) we expand the equation $\stackrel{}{{{\mathit{\omega }}^{\mathrm{2}}}_{\mathrm{˜}}}\mathrm{=}{\mathit{\omega }}_{\mathrm{0}}^{\mathrm{2}}$ as ${\mathit{\omega }}^{\mathrm{3}}\mathrm{+}\mathrm{i}\left(\mathrm{1}\mathrm{+}\mathit{\chi }\right){\mathit{\nu }}_{\mathrm{ni}}{\mathit{\omega }}^{\mathrm{2}}\mathrm{-}{\mathit{\omega }}_{\mathrm{0}}^{\mathrm{2}}\mathit{\omega }\mathrm{-}\mathrm{i}{\mathit{\nu }}_{\mathrm{ni}}{\mathit{\omega }}_{\mathrm{0}}^{\mathrm{2}}\mathrm{=}\mathrm{0.}$(40)The study of the modifications due to neutral-ion collisions of the modes of Edwin & Roberts (1983) reduces to the study of the solutions of Eq. (40). This study is general and is independent of the particular mode considered since we have not specified what mode ω₀ corresponds to. All modes are affected in the same way by neutral-ion collisions when neutral pressure is absent.

Equation (40) is a cubic equation, hence it has three solutions. To determine the nature of the solutions we perform the change of variable ω = −is, so that Eq. (40) becomes a cubic equation in s with all its coefficients real. We compute the discriminant, Δ, of the resulting equation as $\mathrm{\Delta }\mathrm{=}\mathrm{-}{\mathit{\omega }}_{\mathrm{0}}^{\mathrm{2}}\mathrm{[}\mathrm{4}{\left(\mathrm{1}\mathrm{+}\mathit{\chi }\right)}^{\mathrm{3}}{\mathit{\nu }}_{\mathrm{ni}}^{\mathrm{4}}\mathrm{-}\mathrm{(}{\mathit{\chi }}^{\mathrm{2}}\mathrm{+}\mathrm{20}\mathit{\chi }\mathrm{-}{\mathrm{8}}^{\mathrm{)}}{\mathit{\nu }}_{\mathrm{ni}}^{\mathrm{2}}{\mathit{\omega }}_{\mathrm{0}}^{\mathrm{2}}\mathrm{+}\mathrm{4}{{\mathit{\omega }}_{\mathrm{0}}^{\mathrm{4}}}^{\mathrm{]}}\mathit{.}$(41)The discriminant, Δ, is defined so that (i) Eq. (40) has one purely imaginary zero and two complex zeros when Δ < 0; (ii) Eq. (40) has a multiple zero and all its zeros are purely imaginary when Δ = 0; and (iii) Eq. (40) has three distinct purely imaginary zeros when Δ > 0. This criterion points out that oscillatory solutions of Eq. (40) are only possible when Δ < 0.

To determine when oscillatory solutions are not possible we set Δ = 0 and find the corresponding relation between ω₀ and ν_ni in terms of χ as $\frac{{\mathit{\nu }}_{\mathrm{ni}}}{{\mathit{\omega }}_{\mathrm{0}}}\mathrm{=}{\left[\frac{{\mathit{\chi }}^{\mathrm{2}}\mathrm{+}\mathrm{20}\mathit{\chi }\mathrm{-}\mathrm{8}}{\mathrm{8}{\left(\mathrm{1}\mathrm{+}\mathit{\chi }\right)}^{\mathrm{3}}}\mathrm{\pm }\frac{{\mathit{\chi }}^{\mathrm{1}\mathit{/}\mathrm{2}}{\left(\mathit{\chi }\mathrm{-}\mathrm{8}\right)}^{\mathrm{3}\mathit{/}\mathrm{2}}}{\mathrm{8}{\left(\mathrm{1}\mathrm{+}\mathit{\chi }\right)}^{\mathrm{3}}}\right]}^{\mathrm{1}\mathit{/}\mathrm{2}}\mathrm{·}$(42)\\Since ν_ni/ω₀ must be real, Eq. (42) imposes a condition on the minimum value of χ which allows Δ = 0. This minimum value is χ = 8 and the corresponding critical ν_ni/ω₀ is ${\mathit{\nu }}_{\mathrm{ni}}\mathit{/}{\mathit{\omega }}_{\mathrm{0}}\mathrm{=}\mathrm{1}\mathit{/}\mathrm{3}\sqrt{\mathrm{3}}\mathrm{\approx }\mathrm{0.19}$. When χ < 8 oscillatory solutions are always possible. When χ > 8 we have to take into account the + and the − signs in Eq. (42), so that Eq. (42) defines a range of values of ν_ni/ω₀ in which oscillatory solutions are not possible. We call this interval the cut-off region. To our knowledge, Kulsrud & Pearce (1969) were the first to report on the existence of a cut-off region of MHD waves in a partially ionized plasma, although Kulsrud & Pearce (1969) only investigated Alfvén waves in a homogeneous medium. The existence of cut-off regions for MHD waves in structured media is also deduced from the plots of Kumar & Roberts (2003), but these authors did not explore this phenomenon in detail.

We look for approximate analytic expressions of the frequency of the oscillatory solutions. To do this, we assume that ν_ni/ω₀ is outside the cut-off interval defined by Eq. (42), meaning that Eq. (40) has two complex (oscillatory) solutions and one purely imaginary (evanescent) solution. First we focus on the oscillatory solutions. We write ω = ω_R + iω_I and insert this expression into Eq. (40). We assume weak damping and take | ω_I| ≪ |ω_R |. Hence we neglect terms with ${\mathit{\omega }}_{\mathrm{I}}^{\mathrm{2}}$ and higher powers. It is crucial for the validity of this approximation that ν_ni/ω₀ is not within or close to the cut-off region in which ω_R = 0. After some algebraic manipulations we derive approximate expressions for ω_R and ω_I. For simplicity we omit the intermediate steps and give the final expressions, namely $\begin{array}{ccc}{\mathit{\omega }}_{\mathrm{R}}& \mathrm{\approx }& {\mathit{\omega }}_{\mathrm{0}}\sqrt{\frac{{\mathit{\omega }}_{\mathrm{0}}^{\mathrm{2}}\mathrm{+}\left(\mathrm{1}\mathrm{+}\mathit{\chi }\right){\mathit{\nu }}_{\mathrm{ni}}^{\mathrm{2}}}{{\mathit{\omega }}_{\mathrm{0}}^{\mathrm{2}}\mathrm{+}{\left(\mathrm{1}\mathrm{+}\mathit{\chi }\right)}^{\mathrm{2}}{\mathit{\nu }}_{\mathrm{ni}}^{\mathrm{2}}}}\mathit{,}\\ {\mathit{\omega }}_{\mathrm{I}}& \mathrm{\approx }& \mathrm{-}\frac{\mathit{\chi }{\mathit{\nu }}_{\mathrm{ni}}}{\mathrm{2}\mathrm{[}{\mathit{\omega }}_{\mathrm{0}}^{\mathrm{2}}\mathrm{+}{\left(\mathrm{1}\mathrm{+}\mathit{\chi }\right)}^{\mathrm{2}}{{\mathit{\nu }}_{\mathrm{ni}}^{\mathrm{2}}}^{\mathrm{]}}}{\mathit{\omega }}_{\mathrm{0}}^{\mathrm{2}}\mathit{.}\end{array}$A mode with the same ω_I but with ω_R of opposite sign is also a solution. When ν_ni = 0, so that neutral-ion collisions are absent, we find ω_R = ω₀ and ω_I = 0. Hence we recover the ideal undamped modes of Edwin & Roberts (1983). Now these modes are damped due to neutral-ion collisions and the frequency has an imaginary part. In addition, the real part of the frequency is also modified.

It is instructive to study the limit values of ν_ni. For a low collision frequency, ${\mathit{\nu }}_{\mathrm{ni}}^{\mathrm{2}}\mathrm{\ll }{\mathit{\omega }}_{\mathrm{0}}^{\mathrm{2}}$. Equations (43) and (44) become $\begin{array}{ccc}{\mathit{\omega }}_{\mathrm{R}}& \mathrm{\approx }& {\mathit{\omega }}_{\mathrm{0}}\mathit{,}\\ {\mathit{\omega }}_{\mathrm{I}}& \mathrm{\approx }& \mathrm{-}\frac{\mathit{\chi }{\mathit{\nu }}_{\mathrm{ni}}}{\mathrm{2}}\mathrm{·}\end{array}$Thus, the real part of the frequency when ${\mathit{\nu }}_{\mathrm{ni}}^{\mathrm{2}}\mathrm{\ll }{\mathit{\omega }}_{\mathrm{0}}^{\mathrm{2}}$ is just the frequency found in the fully ionized case (Edwin & Roberts 1983) and it does not depend on the amount of neutrals. The imaginary part of the frequency is independent of ω₀, meaning that it is the same for all wave modes. On the other hand, for a high collision frequency, ${\mathit{\nu }}_{\mathrm{ni}}^{\mathrm{2}}\mathrm{\gg }{\mathit{\omega }}_{\mathrm{0}}^{\mathrm{2}}$ and Eqs. (43) and (44) simplify to $\begin{array}{ccc}{\mathit{\omega }}_{\mathrm{R}}& \mathrm{\approx }& \frac{{\mathit{\omega }}_{\mathrm{0}}}{\sqrt{\mathrm{1}\mathrm{+}\mathit{\chi }}}\mathit{,}\\ {\mathit{\omega }}_{\mathrm{I}}& \mathrm{\approx }& \mathrm{-}\frac{\mathit{\chi }{\mathit{\omega }}_{\mathrm{0}}^{\mathrm{2}}}{\mathrm{2}{\left(\mathrm{1}\mathrm{+}\mathit{\chi }\right)}^{\mathrm{2}}{\mathit{\nu }}_{\mathrm{ni}}}\mathrm{·}\end{array}$Now the expression of ω_R involves χ in the denominator, so that the larger the amount of neutrals, the lower ω_R compared to the frequency in the fully ionized case. The frequency is reduced by the factor $\mathrm{1}\mathit{/}\sqrt{\mathrm{1}\mathrm{+}\mathit{\chi }}$ compared to the fully ionized value. This is the same factor as found by Kumar & Roberts (2003) and is equivalent to replace the ion density by the sum of densities of ions and neutrals or, equivalently, to replace ρ_i by (1 + χ)ρ_i. The imaginary part of the frequency depends now on ω₀, meaning that the various modes have different damping rates.

In addition to the oscillatory solutions, Eq. (40) has one purely imaginary solution whose approximation is given by $\mathit{\omega }\mathrm{\approx }\mathrm{-}\mathrm{i}{\mathit{\nu }}_{\mathrm{ni}}\frac{{\mathit{\omega }}_{\mathrm{0}}^{\mathrm{2}}\mathrm{+}{\left(\mathrm{1}\mathrm{+}\mathit{\chi }\right)}^{\mathrm{2}}{\mathit{\nu }}_{\mathrm{ni}}^{\mathrm{2}}}{{\mathit{\omega }}_{\mathrm{0}}^{\mathrm{2}}\mathrm{+}\left(\mathrm{1}\mathrm{+}\mathit{\chi }\right){\mathit{\nu }}_{\mathrm{ni}}^{\mathrm{2}}}\mathrm{·}$(49)\\The perturbations related to these modes are evanescent in time. There exists a different evanescent solution related to each mode, ω₀, of the fully ionized case. Purely imaginary solutions were also found by Zaqarashvili et al. (2011), who studied waves in a partially ionized homogeneous medium. These modes are not present in the fully ionized case. These evanescent perturbations may be relevant during the excitation of the waves, since part of the energy used to excite the waves may go to these modes instead of being used to excite the oscillatory solutions. This could be investigated by going beyond the present normal mode analysis and solving the initial value problem.

4.1.2. Propagating waves

Now we turn to propagating waves. In the fully ionized case the study of propagating waves is equivalent to that of standing waves. Here we shall see that neutral-ion collisions break this equivalence and propagating waves are worth being studied separately. Propagating waves were not investigated by Kumar & Roberts (2003). For propagating waves, we fix the frequency, ω, to a real value and solve the dispersion relation for the complex k_z, i.e., k_z = k_z,R + ik_z,I, where k_R and k_I are the real and imaginary parts of k_z. Therefore the wave amplitude is damped in space due to the exponential factor exp(−k_z,Iz).

Equation (37) also holds for propagating waves. We assume that k_z = k_z,0, with k_z,0 real, is a solution of the ideal dispersion relation of Edwin & Roberts (1983). Then, ${\mathit{k}}_{\mathit{z}}^{\mathrm{2}}\mathrm{=}{\mathit{k}}_{\mathit{z,}\mathrm{0}}^{\mathrm{2}}\frac{\mathit{\omega }\mathrm{+}\mathrm{i}\left(\mathrm{1}\mathrm{+}\mathit{\chi }\right){\mathit{\nu }}_{\mathrm{ni}}}{\mathit{\omega }\mathrm{+}\mathit{i}{\mathit{\nu }}_{\mathrm{ni}}}\mathit{,}$(50)\\is a solution of Eq. (37) following the same argument as explained for standing waves. We write k_z = k_z,R + ik_z,I and insert this expression in Eq. (50). We obtain exact expressions for k_z,R and k_z,I, namely $\begin{array}{ccc}{\mathit{k}}_{\mathit{z,}\mathrm{R}}& \mathrm{=}& {\mathit{k}}_{\mathit{z,}\mathrm{0}}\sqrt{\frac{{\mathit{\omega }}^{\mathrm{2}}\mathrm{+}\mathrm{(}\mathrm{1}\mathrm{+}\mathit{\chi }\mathrm{)}{\mathit{\nu }}_{\mathrm{ni}}^{\mathrm{2}}}{\mathrm{2}\mathrm{(}{\mathit{\omega }}^{\mathrm{2}}\mathrm{+}{{\mathit{\nu }}_{\mathrm{ni}}^{\mathrm{2}}}^{\mathrm{)}}}}\\ & & \\ {\mathit{k}}_{\mathit{z,}\mathrm{I}}& \mathrm{=}& \frac{{\mathit{k}}_{\mathit{z,}\mathrm{0}}^{\mathrm{2}}}{\mathrm{2}{\mathit{k}}_{\mathit{z,}\mathrm{R}}}\frac{\mathit{\chi }{\mathit{\nu }}_{\mathrm{ni}}\mathit{\omega }}{{\mathit{\omega }}^{\mathrm{2}}\mathrm{+}{\mathit{\nu }}_{\mathrm{ni}}^{\mathrm{2}}}\mathrm{·}\end{array}$The same values of k_z,R and k_z,I but with opposite signs are also solutions. In principle, Eq. (51) allows two different values of k_z,R because of the ± sign. However, the solution with the − sign is not physical since it corresponds to k_z,R imaginary, which is an obvious contradiction. For this reason we discard the solution with the − sign and take the + sign in Eq. (51). Now we realize that $\frac{{\mathit{\chi }}^{\mathrm{2}}{\mathit{\nu }}_{\mathrm{ni}}^{\mathrm{2}}{\mathit{\omega }}^{\mathrm{2}}}{{\mathrm{(}{\mathit{\omega }}^{\mathrm{2}}\mathrm{+}\mathrm{(}\mathrm{1}\mathrm{+}\mathit{\chi }\mathrm{)}{{\mathit{\nu }}_{\mathrm{ni}}^{\mathrm{2}}}^{\mathrm{)}}}^{\mathrm{2}}}\mathrm{\ll }\mathrm{1.}$(53)Thus we can approximate k_z,R and k_z,I as $\begin{array}{ccc}{\mathit{k}}_{\mathit{z,}\mathrm{R}}& \mathrm{\approx }& \\ {\mathit{k}}_{\mathit{z,}\mathrm{I}}& \mathrm{\approx }& \end{array}$These expressions are equivalent to Eqs. (43) and (44) obtained for standing waves. However, contrary to Eqs. (43) and (44), Eqs. (54) and (55) are valid for all values of the ratio ν_ni/ω. There is no cut-off region for propagating waves because k_z,R never vanishes. There are no purely imaginary solutions of k_z, i.e., all the solutions always have an oscillatory behavior in z. This is an important difference to standing waves.

4.2. Comparison with numerical results

Here we numerically solve the dispersion relation (Eq. (37)). For simplicity we focus on standing waves and show the results for transverse kink (m = 1) waves only, although we have checked that equivalent results are found for other modes as, e.g., slow magnetoacoustic modes. In the fully ionized case, the frequency of the transverse kink waves in the thin tube (TT) limit, i.e., k_zR ≪ 1, is ω = ω_k, with ω_k the kink frequency given by ${\mathit{\omega }}_{\mathrm{k}}\mathrm{=}{\mathit{k}}_{\mathit{z}}\sqrt{\frac{{\mathit{\rho }}_{\mathrm{i}\mathit{,}\mathrm{0}}{\mathit{c}}_{\mathrm{A}\mathit{,}\mathrm{0}}^{\mathrm{2}}\mathrm{+}{\mathit{\rho }}_{\mathrm{i}\mathit{,}\mathrm{ex}}{\mathit{c}}_{\mathrm{A}\mathit{,}\mathrm{ex}}^{\mathrm{2}}}{{\mathit{\rho }}_{\mathrm{i}\mathit{,}\mathrm{0}}\mathrm{+}{\mathit{\rho }}_{\mathrm{i}\mathit{,}\mathrm{ex}}}}\mathrm{·}$(56)Figure 1 displays the real and imaginary parts of the frequency of the kink mode as functions of the ratio ν_ni/ω_k for a particular choice of parameters given in the caption of the figure. We compare the numerical results with the analytic approximations. In Fig. 1 we have considered χ < 8 so that no cut-off region is present. Regarding the real part of the frequency (Fig. 1a), we obtain that for ν_ni/ω_k ≪ 1 the results are independent of the ionization degree, with ω_R ≈ ω_k as in the fully ionized case. When the ratio ν_ni/ω_k increases, plasma and neutrals are more coupled and the ionization degree becomes a relevant parameter, so that ω_R decreases until de value ${\mathit{\omega }}_{\mathrm{R}}\mathrm{\approx }{\mathit{\omega }}_{\mathrm{k}}\mathit{/}\sqrt{\mathrm{1}\mathrm{+}\mathit{\chi }}$ is reached. This behavior is consistent with the analytic Eq. (43) although the details of the transition are not fully captured by the approximation. Regarding the imaginary part of the frequency (Fig. 1b), we find that ω_I tends to zero in the limits ν_ni/ω_k ≪ 1 and ν_ni/ω_k ≫ 1. The damping is most efficient when ν_ni and ω_k are, approximately, on the same order of magnitude. This result is also consistent with the analytic Eq. (44), although the approximation underestimates the actual damping rate when ω_I is minimal. This is so because Eq. (44) was derived in the weak damping approximation, i.e., | ω_I| ≪ |ω_R |, while | ω_I | and | ω_R | are on the same order when ω_I is minimal, meaning that the damping is strong.

Fig. 1 a)ωR/ωk and b) ωI/ωk versus νni/ωk for the transverse kink (m = 1) mode with kzR = 0.1, ρi,0/ρi,ex = 3, cie,0/cA,0 = 0.2, and χ = 4. The solid line is the result obtained by numerically solving the dispersion relation in the absence of neutral pressure (Eq. (37)). The symbols are the approximate analytic results in the weak damping approximation (Eqs. (43) and (44)).

Now we take χ > 8 and repeat the previous computations. The results are shown in Fig. 2. We notice the presence of a cut-off region for a certain range of ν_ni/ω_k. In this cut-off region ω_R = 0. The location of the cut-off region agrees very well with the range given by Eq. (42). This zone exists for relatively low ν_ni/ω_k. In the solar atmosphere the expected value of the collision frequency (see, e.g., De Pontieu et al. 2001, Figs. 1 and 2) is much higher than the frequency of the observed waves (e.g., De Pontieu et al. 2007; Okamoto & De Pontieu 2011). The minimum value of the neutral-ion collision frequency in the solar atmosphere is on the order of 10 Hz, while the dominant frequency in the observations of chromospheric kink waves by Okamoto & De Pontieu (2011) is around 22 mHz. Therefore, the cut-off region found here is not in the range of ν_ni/ω_k consistent with solar atmospheric parameters and observed wave frequencies. However, in other astrophysical situations the presence of the cut-off region may be relevant (see Kulsrud & Pearce 1969). As expected, the analytic approximations (Eqs. (43) and (44)) completely miss the presence of the cut-off region, although they are reasonably good far from the location of the cut-off.

Fig. 2 Same as Fig. 1 but with χ = 20. The shaded area denotes the cut-off region given by Eq. (42).

It is worth mentioning again that the effect of resonant absorption on the damping of transverse kink waves is absent because we used a piece-wise constant density profile (Ruderman & Roberts 2002; Goossens et al. 2002). Soler et al. (2012a) showed that damping by resonant absorption is more efficient than damping by neutral-ion collisions when ν_ni/ω_k ≫ 1. Hence, the actual damping of transverse kink waves would be stronger than the damping shown here if resonant absorption were taken into account (see details in Soler et al. 2012a).

Fig. 3 a)ωR/ωc,0 and b) ωI/ωc,0 versus νni/ωc,0 for the slow magnetoacoustic mode with χ = 4 and the rest of parameters the same as in Fig. 1. The solid line is the full result taking neutral pressure into account, while the dashed line is the result when neutral pressure is absent.

5. Role of neutral pressure

We incorporate the effect of neutral pressure. We anticipate that neutral pressure is relevant for slow magnetoacoustic waves, while its effect on transverse kink waves is minor. Goossens et al. (2009, 2012) showed that transverse kink waves in flux tubes have the typical properties of surface Alfvén (or Alfvénic) waves. In thin tubes, i.e., k_zR ≪ 1, the main restoring force of these waves is magnetic tension and the gas pressure force is negligible. We computed the kink mode frequency in the presence of neutral pressure and found no significant differences with the result in the absence of neutral pressure. Hence, here we focus on the study of slow magnetoacoustic waves and explore the modification of the slow mode frequency due to neutral-ion collisions when neutral pressure is taken into account. In the fully ionized case (Edwin & Roberts 1983), slow magnetoacoustic waves in the TT limit have frequencies ω ≈ ω_c,0, with ω_c,0 the internal cusp frequency given by ${\mathit{\omega }}_{\mathrm{c}\mathit{,}\mathrm{0}}\mathrm{=}{\mathit{k}}_{\mathit{z}}\sqrt{\frac{{\mathit{c}}_{\mathrm{A}\mathit{,}\mathrm{0}}^{\mathrm{2}}{\mathit{c}}_{\mathrm{ie}\mathit{,}\mathrm{0}}^{\mathrm{2}}}{{\mathit{c}}_{\mathrm{A}\mathit{,}\mathrm{0}}^{\mathrm{2}}\mathrm{+}{\mathit{c}}_{\mathrm{ie}\mathit{,}\mathrm{0}}^{\mathrm{2}}}}\mathrm{·}$(57)This result is independent of the azimuthal wavenumber, m. In a low-β plasma, i.e., when the magnetic pressure is much more important than the gas pressure, ω_c,0 ≈ k_zc_ie,0. Due to the complexity of the general dispersion relation (Eq. (36)) it is not possible to obtain simple analytic approximations of the slow mode frequency. For this reason we perform this investigation in a numerical way.

Figure 3 shows the numerically obtained real and imaginary parts of the slow mode frequency as functions of the ratio ν_ni/ω_c,0 for χ < 8. The rest of parameters are given in the caption of the figure. For comparison, we overplot the results in the absence of neutral pressure, which are obtained by solving Eq. (37) with the same parameters. First we discuss the behavior of the ω_R (Fig. 3a). When ν_ni/ω_c,0 ≪ 1 we recover the result in the fully ionized case, i.e., ω_R ≈ ω_c,0. When ν_ni/ω_c,0 increases, the real part of the slow mode frequency decreases until a plateau is reached at the value ${\mathit{\omega }}_{\mathrm{R}}\mathrm{\approx }{\mathit{\omega }}_{\mathrm{c}\mathit{,}\mathrm{0}}\mathit{/}\sqrt{\mathrm{1}\mathrm{+}\mathit{\chi }}$. This is equivalent to performing the replacements ${\mathit{c}}_{\mathrm{A}}^{\mathrm{2}}\mathrm{\to }{\mathit{c}}_{\mathrm{A}}^{\mathrm{2}}\mathit{/}\mathrm{(}\mathrm{1}\mathrm{+}\mathit{\chi }\mathrm{)}$ and ${\mathit{c}}_{\mathrm{ie}}^{\mathrm{2}}\mathrm{\to }{\mathit{c}}_{\mathrm{ie}}^{\mathrm{2}}\mathit{/}\mathrm{(}\mathrm{1}\mathrm{+}\mathit{\chi }\mathrm{)}$ in Eq. (57) and is the same behavior as that found in the absence of neutral pressure, so that both results are superimposed in Fig. 3a. However, as the ratio ν_ni/ω_c,0 continues to increase, the real part of the frequency rises again until a second plateau is finally reached. The presence of this second plateau of ω_R for high collision frequencies owes its existence to the effect of neutral pressure and is absent when neutral pressure is neglected.

We now analyze in detail the second plateau of ω_R at high collision frequencies. Since we have no simple analytic expressions for the slow-mode frequency we performed a parameter study. We considered various values of c_ie and numerically computed the slow mode ω_R at the second plateau. The result of this parameter study (not shown here for simplicity) points out that the value of ω_R at the second plateau is, approximately, ${\mathit{\omega }}_{\mathrm{R}}\mathrm{\approx }{\mathit{\omega }}_{\mathrm{c}\mathit{,}\mathrm{0}}\sqrt{\frac{\mathrm{2}\mathrm{+}\mathit{\chi }}{\mathrm{1}\mathrm{+}\mathit{\chi }}}{\left(\mathrm{2}\mathrm{+}\frac{\mathit{\chi }{\mathit{c}}_{\mathrm{ie}\mathit{,}\mathrm{0}}^{\mathrm{2}}}{{\mathit{c}}_{\mathrm{A}\mathit{,}\mathrm{0}}^{\mathrm{2}}\mathrm{+}{\mathit{c}}_{\mathrm{ie}\mathit{,}\mathrm{0}}^{\mathrm{2}}}\right)}^{\mathrm{-}\mathrm{1}\mathit{/}\mathrm{2}}\mathrm{·}$(58)This expression can be obtained from Eq. (57) by performing the replacements ${\mathit{c}}_{\mathrm{A}}^{\mathrm{2}}\mathrm{\to }{\mathit{c}}_{\mathrm{A}}^{\mathrm{2}}\mathit{/}\mathrm{(}\mathrm{1}\mathrm{+}\mathit{\chi }\mathrm{)}$ and ${\mathit{c}}_{\mathrm{ie}}^{\mathrm{2}}\mathrm{\to }{\mathit{c}}_{\mathrm{s}\mathit{,}\mathrm{eff}}^{\mathrm{2}}$, with the expression of the effective sound velocity, c_s,eff, given in Eq. (29). In a low-β plasma, ${\mathit{c}}_{\mathrm{A}}^{\mathrm{2}}\mathrm{\gg }{\mathit{c}}_{\mathrm{ie}}^{\mathrm{2}}$ and ω_c,0 ≈ k_zc_ie,0 so that Eq. (58) simplifies to ${\mathit{\omega }}_{\mathrm{R}}\mathrm{\approx }{\mathit{k}}_{\mathit{z}}{\mathit{c}}_{\mathrm{ie}\mathit{,}\mathrm{0}}\sqrt{\frac{\mathrm{2}\mathrm{+}\mathit{\chi }}{\mathrm{2}\mathrm{(}\mathrm{1}\mathrm{+}\mathit{\chi }\mathrm{)}}}\mathrm{=}{\mathit{k}}_{\mathit{z}}{\mathit{c}}_{\mathrm{s}\mathit{,}\mathrm{eff}\mathit{,}\mathrm{0}}\mathit{.}$(59)Hence, for practical purposes we obtain that at high collision frequencies the sound velocity of the ionized fluid has to be replaced by the effective sound velocity in the expressions of the frequency. The effective sound velocity takes into account the sound velocity of neutrals in addition to the sound velocity of ions. In contrast, for intermediate collision frequencies it is enough to replace c_ie by ${\mathit{c}}_{\mathrm{ie}}\mathit{/}\sqrt{\mathrm{1}\mathrm{+}\mathit{\chi }}$, so that the sound velocity of neutrals can be ignored.

We turn to the imaginary part of the frequency (Fig. 3b). As for the real part of the frequency, there are striking differences between the results with and without neutral pressure. The full result shows the presence of two different minima of ω_I, whereas only the first minimum is present in the absence of neutral pressure. The second minimum takes place in the range of high collision frequencies compared to the ideal cusp frequency. As a consequence of the presence of this second minimum for high collision frequencies, the efficiency of neutral-ion collisions for the damping of the slow mode is higher than for transverse kink modes when ν_ni approaches realistic values. The stronger damping of the slow modes due to neutral-ion collisions compared to the weak damping of the transverse kink modes was already obtained by Soler et al. (2009a) in the single-fluid approximation.

Fig. 4 Same as Fig. 3 but with χ = 20. The shaded area denotes the cut-off region given by Eq. (42).

Now we consider χ > 8 and repeat the previous computations (see Fig. 4 for χ = 20). As in the case without neutral pressure, we notice the presence of a cut-off region that agrees well with the range given by Eq. (42). The existence of the cut-off region is not affected by neutral pressure, i.e., the location of the cut-off region is the same as when neutral pressure is neglected. For the set of parameters considered in Fig. 4 the cut-off region appears in the vicinity of ν_ni/ω_c,0 ≈ 10^-1. We increase the ionization fraction to χ = 100 (Fig. 5). In addition to the cut-off region described above, now we find a second cut-off region. The second cut-off region appears for relatively high collision frequencies, ν_ni/ω_c,0 ≈ 10¹ for the set of parameters used in Fig. 5. This new cut-off region is only present when neutral pressure is taken into account. The result in the absence of neutral pressure completely misses this second cut-off region. The existence of a second cut-off region is a result of slow modes only, since it is not found for transverse kink modes. Unfortunately, unlike for the first cut-off region, we have no analytic expression that gives us the location of the second cut-off region. Instead, out numerical study informs us that the second cut-off appears when χ ≳ 24, although this value may be affected by the particular choice of parameters used in the model.

Fig. 5 Same as Fig. 4 but with χ = 100.

As for the transverse kink modes, it is worth mentioning that some mechanisms not considered in the present analysis may produce damping of the slow mode comparable to or even stronger than that due to neutral-ion collisions. This is the case for non-adiabatic mechanisms, e.g., thermal conduction and radiative losses (see, e.g., De Moortel & Hood 2003; Soler et al. 2008). Partial ionization increases the efficiency of thermal conduction due to the additional contribution of neutral conduction (see Forteza et al. 2008; Soler et al. 2010), so that the combined effect of thermal conduction and partial ionization produces a very efficient damping of the slow mode (see, e.g., Khodachenko et al. 2006; Soler 2010).

6. Discussion

In this paper we have explored the impact of partial ionization on the properties of the MHD waves in a cylindrical magnetic flux tube using the two-fluid formalism. Unlike previous works (e.g., Kumar & Roberts 2003; Soler et al. 2012a), we considered a consistent description of the dynamics of the neutral fluid that takes gas pressure into account. We derived the dispersion relation for the wave modes, that is, the two-fluid generalization of the dispersion relation of Edwin & Roberts (1983). Instead of performing particular applications, we considered the neutral-ion collision frequency as an arbitrary parameter. The modifications of the wave frequencies due to neutral-ion collisions were explored.

First, we investigated the case that neglects neutral pressure. This case was previously investigated for surface waves in a Cartesian interface by Kumar & Roberts (2003). For ν_ni/|ω| ≪ 1 we recovered the frequencies in fully ionized plasma (Edwin & Roberts 1983). As the ratio ν_ni/|ω | increases, both ion-electrons and neutrals become more and more coupled. For ν_ni/|ω| ≫ 1 ion-electrons and neutrals behave as a single fluid. Compared to the fully ionized case, when ν_ni/|ω| ≫ 1 the wave frequencies are lower and depend on the ionization ratio, χ. The frequencies are reduced by the factor (χ + 1)^−1/2 compared to their values in the fully ionized case. This result is the same for all wave modes and is equivalent to replacing the density of ions by the sum of densities of ions and neutrals. The same conclusion was reached in the previous work by Kumar & Roberts (2003). Concerning the damping by neutral-ion collisions, we found that damping is weak when ν_ni/|ω| ≪ 1 unless the plasma is very weakly ionized, i.e., χ ≫ 1. Damping is most efficient when the wave frequency is on the same order of magnitude as the collision frequency. Damping is again weak when ν_ni/|ω| ≫ 1. This last result is always true regardless of the value of χ. Again, our results agree well with the previous findings of Kumar & Roberts (2003).

Then we incorporated the effect of neutral gas pressure. Since transverse kink waves are largely insensitive to the sound velocity, we obtained the same results as in the absence of gas pressure. However, the consideration of neutral pressure has dramatic consequences for slow magnetoacoustic modes. For high collision frequencies the slow-mode frequency is higher than the value obtained in the absence of neutral pressure. This is so because the slow-mode frequency depends on an effective sound velocity that corresponds to the weighted average of the sound velocities of ions-electrons and neutrals. When χ ≫ 1, the effective sound velocity tends to be like the sound velocity of neutrals. Concerning the imaginary part of the frequency, the results with and without neutral pressure also show significant differences for high collision frequencies. The imaginary part of the slow-mode frequency has an additional minimum at high collision frequencies when neutral pressure is included. This causes the slow mode to be more efficiently damped by collisions when neutral pressure is taken into account.

Here we go back to the discussion given in the introduction about the necessity of using a multi-fluid theory instead of the simpler single-fluid MHD approximation. In the solar atmosphere the expected value of the neutral-ion collision frequency (see De Pontieu et al. 2001, Figs. 1 and 2) is much higher than the frequency of the observed waves (e.g., De Pontieu et al. 2007; Okamoto & De Pontieu 2011). In this case, the results of this paper point out that for practical purposes, it is enough to use the single-fluid MHD approximation, but taking into account the following simple recipe to adapt the results of fully ionized models to the partially ionized case. First of all, the ion density has to be replaced by the total density, i.e., the sum of densities of ions and neutrals. Second, the effective sound velocity should be used instead of the sound velocity of the ionized fluid. This effective sound velocity is the weighted average of the sound velocities of ion-electrons and neutrals. This approach is appropriate for works that are not interested in the details of the interaction between ions and neutrals, but only in the result of this interaction on the wave frequencies. Hence, this recipe may be useful to perform magneto-seismology of chromospheric spicules and other partially ionized structures as prominence threads. For example, in their seismological analysis of transverse kink waves in spicules, Verth et al. (2011) considered a definition of the kink velocity that incorporates the density of neutrals. Although our recipe provides a good approximation to obtain the wave frequencies, it misses the effect of damping due to neutral-ion collisions, which might have an impact on the seismological estimates using the amplitude of the waves. Taking into account the range of observed wave frequencies, the damping of kink waves due to neutral-ion collisions may be negligible compared to other effects such as resonant absorption (Soler et al. 2012a). However, the damping of slow magnetoacoustic waves may be important even for high collision frequencies. Other damping effects such as resonant absorption and non-adiabatic mechanisms should be considered in addition to neutral-ion collisions to explain the damping of MHD waves in the solar atmosphere (see, e.g., Khodachenko et al. 2006; Soler et al. 2012a).

On the other hand, a relevant result obtained in this paper is the presence of cut-offs for certain combinations of parameters (see Kulsrud & Pearce 1969). This result cannot be obtained in the single-fluid MHD approximation. At the cut-offs the real part of the frequency vanishes, meaning that perturbations are evanescent in time instead of oscillatory. For transverse kink waves we found a cut-off region for low values of ν_ni/ω_k when χ > 8. For slow magnetoacoustic waves there are two cut-off regions. The first one is the same as that found for transverse kink waves, whereas the second one is only present in the case of slow modes when neutral pressure is considered. The second cut-off region takes place when χ ≳ 24 and relatively high collision frequencies. The relevance of these cut-offs for wave propagation in the solar atmosphere should be investigated in forthcoming works.

Finally, we are aware that the model used in this paper is simple and misses effects that may be of importance in reality. The model should be improved in the future by considering additional ingredients not included in the present analysis. Among these effects, the variation of physical parameters, e.g., density, ionization degree, etc., along the magnetic wave guide, ionization and recombination processes, and the influence of a dynamic background are worth being explored in the near future.

Acknowledgments

We acknowledge the anonymous referee for his/her constructive comments. We thank T. V. Zaqarashvili for reading an early draft of this paper and for giving helpful comments. RS, JLB, and MG acknowledge support from MINECO and FEDER funds through project AYA2011-22846. RS and JLB acknowledge support from CAIB through the “Grups Competitius” scheme and FEDER funds. MG acknowledges support from KU Leuven via GOA/2009-009. AJD acknowledges support from MINECO through project AYA2010-1802.

References

Appendix A: expression of the dispersion relation

The dispersion relation is , with given by the solution of the following determinant, ${\mathrm{𝒟}}_{\mathit{m}}\left(\mathit{\omega ,}{\mathit{k}}_{\mathit{z}}\right)\mathrm{=}\left|\begin{array}{c}\\ {\mathit{a}}_{\mathrm{11}}& {\mathit{a}}_{\mathrm{12}}& {\mathit{a}}_{\mathrm{13}}& {\mathit{a}}_{\mathrm{14}}\\ {\mathit{a}}_{\mathrm{21}}& {\mathit{a}}_{\mathrm{22}}& {\mathit{a}}_{\mathrm{23}}& {\mathit{a}}_{\mathrm{24}}\\ {\mathit{a}}_{\mathrm{31}}& {\mathit{a}}_{\mathrm{32}}& {\mathit{a}}_{\mathrm{33}}& {\mathit{a}}_{\mathrm{34}}\\ {\mathit{a}}_{\mathrm{41}}& {\mathit{a}}_{\mathrm{42}}& {\mathit{a}}_{\mathrm{43}}& {\mathit{a}}_{\mathrm{44}}\end{array}\right|\mathit{,}$(A.1)with $\begin{array}{ccc}{\mathit{a}}_{\mathrm{11}}& \mathrm{=}& {\mathit{a}}_{\mathrm{12}}\mathrm{=}{\mathit{a}}_{\mathrm{13}}\mathrm{=}{\mathit{a}}_{\mathrm{14}}\mathrm{=}\mathrm{1}\mathit{,}\\ {\mathit{a}}_{\mathrm{21}}& \mathrm{=}& \frac{{\mathit{q}}_{\mathrm{2}\mathit{,}\mathrm{0}}^{\mathrm{2}}}{{\mathit{k}}_{\mathrm{1}\mathit{,}\mathrm{0}}^{\mathrm{2}}\mathrm{-}{\mathit{k}}_{\mathrm{n}\mathit{,}\mathrm{0}}^{\mathrm{2}}}\mathit{,} {\mathit{a}}_{\mathrm{22}}\mathrm{=}\frac{{\mathit{q}}_{\mathrm{2}\mathit{,}\mathrm{0}}^{\mathrm{2}}}{{\mathit{k}}_{\mathrm{2}\mathit{,}\mathrm{0}}^{\mathrm{2}}\mathrm{-}{\mathit{k}}_{\mathrm{n}\mathit{,}\mathrm{0}}^{\mathrm{2}}}\mathit{,}\\ {\mathit{a}}_{\mathrm{23}}& \mathrm{=}& \frac{{\mathit{q}}_{\mathrm{2}\mathit{,}\mathrm{ex}}^{\mathrm{2}}}{{\mathit{k}}_{\mathrm{1}\mathit{,}\mathrm{ex}}^{\mathrm{2}}\mathrm{-}{\mathit{k}}_{\mathrm{n}\mathit{,}\mathrm{ex}}^{\mathrm{2}}}\mathit{,} {\mathit{a}}_{\mathrm{24}}\mathrm{=}\frac{{\mathit{q}}_{\mathrm{2}\mathit{,}\mathrm{ex}}^{\mathrm{2}}}{{\mathit{k}}_{\mathrm{2}\mathit{,}\mathrm{ex}}^{\mathrm{2}}\mathrm{-}{\mathit{k}}_{\mathrm{n}\mathit{,}\mathrm{ex}}^{\mathrm{2}}}\mathit{,}\end{array}$$\begin{array}{ccc}{\mathit{a}}_{\mathrm{31}}& \mathrm{=}& \frac{{\mathit{k}}_{\mathrm{1}\mathit{,}\mathrm{0}}}{{\mathit{\rho }}_{\mathrm{i}\mathit{,}\mathrm{0}}\mathrm{(}\stackrel{}{{{\mathit{\omega }}^{\mathrm{2}}}_{\mathrm{˜}}}\mathrm{-}{{\mathit{\omega }}_{\mathrm{A}\mathit{,}\mathrm{0}}^{\mathrm{2}}}^{\mathrm{)}}}\left(\mathrm{1}\mathrm{-}\frac{\mathrm{i}{\mathit{\nu }}_{\mathrm{ni}}}{\mathit{\omega }\mathrm{+}\mathrm{i}{\mathit{\nu }}_{\mathrm{ni}}}\frac{{\mathit{q}}_{\mathrm{2}\mathit{,}\mathrm{0}}^{\mathrm{2}}}{{\mathit{k}}_{\mathrm{1}\mathit{,}\mathrm{0}}^{\mathrm{2}}\mathrm{-}{\mathit{k}}_{\mathrm{n}\mathit{,}\mathrm{0}}^{\mathrm{2}}}\right)\\ & & \mathrm{\times }\frac{{\mathit{J}}_{\mathit{m}}^{\mathrm{\prime }}\mathrm{\left(}{\mathit{k}}_{\mathrm{1}\mathit{,}\mathrm{0}}\mathit{R}\mathrm{\right)}}{{\mathit{J}}_{\mathit{m}}\mathrm{\left(}{\mathit{k}}_{\mathrm{1}\mathit{,}\mathrm{0}}\mathit{R}\mathrm{\right)}}\mathit{,}\end{array}$(A.5)$\begin{array}{ccc}{\mathit{a}}_{\mathrm{32}}& \mathrm{=}& \frac{{\mathit{k}}_{\mathrm{2}\mathit{,}\mathrm{0}}}{{\mathit{\rho }}_{\mathrm{i}\mathit{,}\mathrm{0}}\mathrm{(}\stackrel{}{{{\mathit{\omega }}^{\mathrm{2}}}_{\mathrm{˜}}}\mathrm{-}{{\mathit{\omega }}_{\mathrm{A}\mathit{,}\mathrm{0}}^{\mathrm{2}}}^{\mathrm{)}}}\left(\mathrm{1}\mathrm{-}\frac{\mathrm{i}{\mathit{\nu }}_{\mathrm{ni}}}{\mathit{\omega }\mathrm{+}\mathrm{i}{\mathit{\nu }}_{\mathrm{ni}}}\frac{{\mathit{q}}_{\mathrm{2}\mathit{,}\mathrm{0}}^{\mathrm{2}}}{{\mathit{k}}_{\mathrm{2}\mathit{,}\mathrm{0}}^{\mathrm{2}}\mathrm{-}{\mathit{k}}_{\mathrm{n}\mathit{,}\mathrm{0}}^{\mathrm{2}}}\right)\\ & & \mathrm{\times }\frac{{\mathit{J}}_{\mathit{m}}^{\mathrm{\prime }}\mathrm{\left(}{\mathit{k}}_{\mathrm{2}\mathit{,}\mathrm{0}}\mathit{R}\mathrm{\right)}}{{\mathit{J}}_{\mathit{m}}\mathrm{\left(}{\mathit{k}}_{\mathrm{2}\mathit{,}\mathrm{0}}\mathit{R}\mathrm{\right)}}\mathit{,}\end{array}$(A.6)$\begin{array}{ccc}{\mathit{a}}_{\mathrm{33}}& \mathrm{=}& \frac{{\mathit{k}}_{\mathrm{1}\mathit{,}\mathrm{ex}}}{{\mathit{\rho }}_{\mathrm{i}\mathit{,}\mathrm{ex}}\mathrm{(}\stackrel{}{{{\mathit{\omega }}^{\mathrm{2}}}_{\mathrm{˜}}}\mathrm{-}{{\mathit{\omega }}_{\mathrm{A}\mathit{,}\mathrm{ex}}^{\mathrm{2}}}^{\mathrm{)}}}\left(\mathrm{1}\mathrm{-}\frac{\mathrm{i}{\mathit{\nu }}_{\mathrm{ni}}}{\mathit{\omega }\mathrm{+}\mathrm{i}{\mathit{\nu }}_{\mathrm{ni}}}\frac{{\mathit{q}}_{\mathrm{2}\mathit{,}\mathrm{ex}}^{\mathrm{2}}}{{\mathit{k}}_{\mathrm{1}\mathit{,}\mathrm{ex}}^{\mathrm{2}}\mathrm{-}{\mathit{k}}_{\mathrm{n}\mathit{,}\mathrm{ex}}^{\mathrm{2}}}\right)\\ & & \mathrm{\times }\frac{{\mathit{H}}_{\mathit{m}}^{\mathrm{\prime }\mathrm{\left(}\mathrm{1}\mathrm{\right)}}\mathrm{\left(}{\mathit{k}}_{\mathrm{1}\mathit{,}\mathrm{ex}}\mathit{R}\mathrm{\right)}}{{\mathit{H}}_{\mathit{m}}^{\mathrm{\left(}\mathrm{1}\mathrm{\right)}}\mathrm{\left(}{\mathit{k}}_{\mathrm{1}\mathit{,}\mathrm{ex}}\mathit{R}\mathrm{\right)}}\mathit{,}\end{array}$(A.7)$\begin{array}{ccc}{\mathit{a}}_{\mathrm{34}}& \mathrm{=}& \frac{{\mathit{k}}_{\mathrm{2}\mathit{,}\mathrm{ex}}}{{\mathit{\rho }}_{\mathrm{i}\mathit{,}\mathrm{ex}}\mathrm{(}\stackrel{}{{{\mathit{\omega }}^{\mathrm{2}}}_{\mathrm{˜}}}\mathrm{-}{{\mathit{\omega }}_{\mathrm{A}\mathit{,}\mathrm{ex}}^{\mathrm{2}}}^{\mathrm{)}}}\left(\mathrm{1}\mathrm{-}\frac{\mathrm{i}{\mathit{\nu }}_{\mathrm{ni}}}{\mathit{\omega }\mathrm{+}\mathrm{i}{\mathit{\nu }}_{\mathrm{ni}}}\frac{{\mathit{q}}_{\mathrm{2}\mathit{,}\mathrm{ex}}^{\mathrm{2}}}{{\mathit{k}}_{\mathrm{2}\mathit{,}\mathrm{ex}}^{\mathrm{2}}\mathrm{-}{\mathit{k}}_{\mathrm{n}\mathit{,}\mathrm{ex}}^{\mathrm{2}}}\right)\\ & & \mathrm{\times }\frac{{\mathit{H}}_{\mathit{m}}^{\mathrm{\prime }\mathrm{\left(}\mathrm{1}\mathrm{\right)}}\mathrm{\left(}{\mathit{k}}_{\mathrm{2}\mathit{,}\mathrm{ex}}\mathit{R}\mathrm{\right)}}{{\mathit{H}}_{\mathit{m}}^{\mathrm{\left(}\mathrm{1}\mathrm{\right)}}\mathrm{\left(}{\mathit{k}}_{\mathrm{2}\mathit{,}\mathrm{ex}}\mathit{R}\mathrm{\right)}}\mathit{,}\end{array}$(A.8)$\begin{array}{ccc}{\mathit{a}}_{\mathrm{41}}& \mathrm{=}& {\mathit{k}}_{\mathrm{1}\mathit{,}\mathrm{0}}[\frac{\mathrm{i}{\mathit{\nu }}_{\mathrm{ni}}}{\mathit{\omega }\mathrm{+}\mathrm{i}{\mathit{\nu }}_{\mathrm{ni}}}\frac{\mathrm{1}}{{\mathit{\rho }}_{\mathrm{i}\mathit{,}\mathrm{0}}\mathrm{(}\stackrel{}{{{\mathit{\omega }}^{\mathrm{2}}}_{\mathrm{˜}}}\mathrm{-}{{\mathit{\omega }}_{\mathrm{A}\mathit{,}\mathrm{0}}^{\mathrm{2}}}^{\mathrm{)}}}\\ & & \mathrm{\times }\left(\mathrm{1}\mathrm{-}\frac{\mathrm{i}{\mathit{\nu }}_{\mathrm{ni}}}{\mathit{\omega }\mathrm{+}\mathrm{i}{\mathit{\nu }}_{\mathrm{ni}}}\frac{{\mathit{q}}_{\mathrm{2}\mathit{,}\mathrm{0}}^{\mathrm{2}}}{{\mathit{k}}_{\mathrm{1}\mathit{,}\mathrm{0}}^{\mathrm{2}}\mathrm{-}{\mathit{k}}_{\mathrm{n}\mathit{,}\mathrm{0}}^{\mathrm{2}}}\right)\\ & & \mathrm{-}\frac{\mathrm{1}}{{\mathit{\rho }}_{\mathrm{n}\mathit{,}\mathrm{0}}\mathit{\omega }\left(\mathit{\omega }\mathrm{+}\mathrm{i}{\mathit{\nu }}_{\mathrm{ni}}\right)}\frac{{\mathit{q}}_{\mathrm{2}\mathit{,}\mathrm{0}}^{\mathrm{2}}}{{\mathit{k}}_{\mathrm{1}\mathit{,}\mathrm{0}}^{\mathrm{2}}\mathrm{-}{\mathit{k}}_{\mathrm{n}\mathit{,}\mathrm{0}}^{\mathrm{2}}}]\frac{{\mathit{J}}_{\mathit{m}}^{\mathrm{\prime }}\mathrm{\left(}{\mathit{k}}_{\mathrm{1}\mathit{,}\mathrm{0}}\mathit{R}\mathrm{\right)}}{{\mathit{J}}_{\mathit{m}}\mathrm{\left(}{\mathit{k}}_{\mathrm{1}\mathit{,}\mathrm{0}}\mathit{R}\mathrm{\right)}}\mathit{,}\end{array}$(A.9)$\begin{array}{ccc}{\mathit{a}}_{\mathrm{42}}& \mathrm{=}& {\mathit{k}}_{\mathrm{2}\mathit{,}\mathrm{0}}[\frac{\mathrm{i}{\mathit{\nu }}_{\mathrm{ni}}}{\mathit{\omega }\mathrm{+}\mathrm{i}{\mathit{\nu }}_{\mathrm{ni}}}\frac{\mathrm{1}}{{\mathit{\rho }}_{\mathrm{i}\mathit{,}\mathrm{0}}\mathrm{(}\stackrel{}{{{\mathit{\omega }}^{\mathrm{2}}}_{\mathrm{˜}}}\mathrm{-}{{\mathit{\omega }}_{\mathrm{A}\mathit{,}\mathrm{0}}^{\mathrm{2}}}^{\mathrm{)}}}\\ & & \mathrm{\times }\left(\mathrm{1}\mathrm{-}\frac{\mathrm{i}{\mathit{\nu }}_{\mathrm{ni}}}{\mathit{\omega }\mathrm{+}\mathrm{i}{\mathit{\nu }}_{\mathrm{ni}}}\frac{{\mathit{q}}_{\mathrm{2}\mathit{,}\mathrm{0}}^{\mathrm{2}}}{{\mathit{k}}_{\mathrm{2}\mathit{,}\mathrm{0}}^{\mathrm{2}}\mathrm{-}{\mathit{k}}_{\mathrm{n}\mathit{,}\mathrm{0}}^{\mathrm{2}}}\right)\\ & & \mathrm{-}\frac{\mathrm{1}}{{\mathit{\rho }}_{\mathrm{n}\mathit{,}\mathrm{0}}\mathit{\omega }\left(\mathit{\omega }\mathrm{+}\mathrm{i}{\mathit{\nu }}_{\mathrm{ni}}\right)}\frac{{\mathit{q}}_{\mathrm{2}\mathit{,}\mathrm{0}}^{\mathrm{2}}}{{\mathit{k}}_{\mathrm{2}\mathit{,}\mathrm{0}}^{\mathrm{2}}\mathrm{-}{\mathit{k}}_{\mathrm{n}\mathit{,}\mathrm{0}}^{\mathrm{2}}}]\frac{{\mathit{J}}_{\mathit{m}}^{\mathrm{\prime }}\mathrm{\left(}{\mathit{k}}_{\mathrm{2}\mathit{,}\mathrm{0}}\mathit{R}\mathrm{\right)}}{{\mathit{J}}_{\mathit{m}}\mathrm{\left(}{\mathit{k}}_{\mathrm{2}\mathit{,}\mathrm{0}}\mathit{R}\mathrm{\right)}}\mathit{,}\end{array}$(A.10)$\begin{array}{ccc}{\mathit{a}}_{\mathrm{43}}& \mathrm{=}& {\mathit{k}}_{\mathrm{1}\mathit{,}\mathrm{ex}}[\frac{\mathrm{i}{\mathit{\nu }}_{\mathrm{ni}}}{\mathit{\omega }\mathrm{+}\mathrm{i}{\mathit{\nu }}_{\mathrm{ni}}}\frac{\mathrm{1}}{{\mathit{\rho }}_{\mathrm{i}\mathit{,}\mathrm{ex}}\mathrm{(}\stackrel{}{{{\mathit{\omega }}^{\mathrm{2}}}_{\mathrm{˜}}}\mathrm{-}{{\mathit{\omega }}_{\mathrm{A}\mathit{,}\mathrm{ex}}^{\mathrm{2}}}^{\mathrm{)}}}\\ & & \mathrm{\times }\left(\mathrm{1}\mathrm{-}\frac{\mathrm{i}{\mathit{\nu }}_{\mathrm{ni}}}{\mathit{\omega }\mathrm{+}\mathrm{i}{\mathit{\nu }}_{\mathrm{ni}}}\frac{{\mathit{q}}_{\mathrm{2}\mathit{,}\mathrm{ex}}^{\mathrm{2}}}{{\mathit{k}}_{\mathrm{1}\mathit{,}\mathrm{ex}}^{\mathrm{2}}\mathrm{-}{\mathit{k}}_{\mathrm{n}\mathit{,}\mathrm{ex}}^{\mathrm{2}}}\right)\\ & & \mathrm{-}\frac{\mathrm{1}}{{\mathit{\rho }}_{\mathrm{n}\mathit{,}\mathrm{ex}}\mathit{\omega }\left(\mathit{\omega }\mathrm{+}\mathrm{i}{\mathit{\nu }}_{\mathrm{ni}}\right)}\frac{{\mathit{q}}_{\mathrm{2}\mathit{,}\mathrm{ex}}^{\mathrm{2}}}{{\mathit{k}}_{\mathrm{1}\mathit{,}\mathrm{ex}}^{\mathrm{2}}\mathrm{-}{\mathit{k}}_{\mathrm{n}\mathit{,}\mathrm{ex}}^{\mathrm{2}}}]\frac{{\mathit{H}}_{\mathit{m}}^{\mathrm{\prime }\mathrm{\left(}\mathrm{1}\mathrm{\right)}}\mathrm{\left(}{\mathit{k}}_{\mathrm{1}\mathit{,}\mathrm{ex}}\mathit{R}\mathrm{\right)}}{{\mathit{H}}_{\mathit{m}}^{\mathrm{\left(}\mathrm{1}\mathrm{\right)}}\mathrm{\left(}{\mathit{k}}_{\mathrm{1}\mathit{,}\mathrm{ex}}\mathit{R}\mathrm{\right)}}\mathit{,}\end{array}$(A.11)$\begin{array}{ccc}{\mathit{a}}_{\mathrm{44}}& \mathrm{=}& {\mathit{k}}_{\mathrm{2}\mathit{,}\mathrm{ex}}[\frac{\mathrm{i}{\mathit{\nu }}_{\mathrm{ni}}}{\mathit{\omega }\mathrm{+}\mathrm{i}{\mathit{\nu }}_{\mathrm{ni}}}\frac{\mathrm{1}}{{\mathit{\rho }}_{\mathrm{i}\mathit{,}\mathrm{ex}}\mathrm{(}\stackrel{}{{{\mathit{\omega }}^{\mathrm{2}}}_{\mathrm{˜}}}\mathrm{-}{{\mathit{\omega }}_{\mathrm{A}\mathit{,}\mathrm{ex}}^{\mathrm{2}}}^{\mathrm{)}}}\\ & & \mathrm{\times }\left(\mathrm{1}\mathrm{-}\frac{\mathrm{i}{\mathit{\nu }}_{\mathrm{ni}}}{\mathit{\omega }\mathrm{+}\mathrm{i}{\mathit{\nu }}_{\mathrm{ni}}}\frac{{\mathit{q}}_{\mathrm{2}\mathit{,}\mathrm{ex}}^{\mathrm{2}}}{{\mathit{k}}_{\mathrm{2}\mathit{,}\mathrm{ex}}^{\mathrm{2}}\mathrm{-}{\mathit{k}}_{\mathrm{n}\mathit{,}\mathrm{ex}}^{\mathrm{2}}}\right)\\ & & \mathrm{-}\frac{\mathrm{1}}{{\mathit{\rho }}_{\mathrm{n}\mathit{,}\mathrm{ex}}\mathit{\omega }\left(\mathit{\omega }\mathrm{+}\mathrm{i}{\mathit{\nu }}_{\mathrm{ni}}\right)}\frac{{\mathit{q}}_{\mathrm{2}\mathit{,}\mathrm{ex}}^{\mathrm{2}}}{{\mathit{k}}_{\mathrm{2}\mathit{,}\mathrm{ex}}^{\mathrm{2}}\mathrm{-}{\mathit{k}}_{\mathrm{n}\mathit{,}\mathrm{ex}}^{\mathrm{2}}}]\frac{{\mathit{H}}_{\mathit{m}}^{\mathrm{\prime }\mathrm{\left(}\mathrm{1}\mathrm{\right)}}\mathrm{\left(}{\mathit{k}}_{\mathrm{2}\mathit{,}\mathrm{ex}}\mathit{R}\mathrm{\right)}}{{\mathit{H}}_{\mathit{m}}^{\mathrm{\left(}\mathrm{1}\mathrm{\right)}}\mathrm{\left(}{\mathit{k}}_{\mathrm{2}\mathit{,}\mathrm{ex}}\mathit{R}\mathrm{\right)}}\mathrm{·}\end{array}$(A.12)