Linear Waves in Partially Ionized Plasmas in Ionization Non-equilibrium

We aim to investigate the properties of linear Alfv\'en and slow magnetoacoustic waves in a partially ionised plasma in ionisation non-equilibrium. The propagation characteristics of these waves are studied within the framework of a two-fluid plasma in terms of the collisional strength between heavy particles for different degrees of ionisation. In the ionisation non-equilibrium state the rates of ionisation and recombinations are not equal. For analytical progress we assume a background that is ionisation equilibrium, the non-equilibrium is driven by perturbations in the system, {\bf therefore, non-equilibrium effects are related to the perturbed state of the plasma}. Using simple analytical methods, we show that ionisation non-equilibrium can provide an additional coupling between ions and neutrals (implicitly a secondary damping mechanism in the collisionless limit) and this process is able to keep the neutrals in the system even in the collisionless limit. Due to the coupling between different species waves become dispersive. The present study improves our understanding of the complexity of dynamical processes partially ionised plasma in the lower solar atmosphere and solar prominences. Our results clearly show that the problem of partial ionisation and non-equilibrium ionisation introduce new aspects of plasma dynamics with consequences on the evolution waves and their dissipation.


INTRODUCTION
One of the most intriguing aspects that have been largely omitted so far in the field of dynamical process in solar plasmas is that the plasma in the lower atmosphere is partially ionized, with plasma made up from charged particles and neutrals that are interacting through collision. Although the exact degree of ionization is not fully known, the ratio of electron density to neutral hydrogen density covers a few orders of magnitude from the photosphere to the top of the chromosphere. The ionization state of the plasma is a very important factor, as the collision between various species will significantly enhance transport processes that control the appearance and evolution of instabilities in the presence of inhomogeneous flows. Depending on the range of frequencies we are interested in, the description of the dynamics in these regions can be described within the framework of single fluid magnetohydrodynamics (MHD) (for frequencies smaller than the ion-neutral frequency), two-fluid MHD (when frequencies of interest are larger than the ion-neutral collisional frequency), or multi-fluid description (in the high-frequency regime when waves' frequency is comparable to the electron-ion collisional frequency).
Research in the dynamical evolution of physical phenomena in the solar atmosphere is based on the assumption of ionization equilibrium and the equilibrium Maxwellian distribution of particles. However, this model is not accurate for rapidly changing phenomena (e.g., high frequency waves, shock waves), for rapid energy releases where high-energy tail of the electron distribution are observed. Non-equilibrium ionization can occur during heating or cooling events, significantly affecting line intensities and subsequently the plasma diagnostics (Bradshaw and Mason, 2003;Bradshaw et al., 2004). Departures from the equilibrium Maxwellian distribution have also been inferred from chromospheric and transition region line emission (Dzifčáková and Kulinová, 2011). In the chromosphere the ionization/recombination relaxation times scales are of the order of 10 3 − 10 5 s (Carlsson and Stein, 2002), meaning that dynamics occurring below this scales will be affected by non-equilibrium effects. Because the relaxation timescale is much longer than dynamic timescales, H ionization does not have time to reach its equilibrium value and its fluctuations are much smaller than the variation of its statistical equilibrium value appropriate for the instantaneous conditions. The problem of waves and oscillations in two fluid partially ionized plasmas is a relatively new area of solar physics, nevertheless some fundamental properties of such environments are already established (see, e.g., Soler et al., 2010Soler et al., , 2013Zaqarashvili et al., 2011;Ballester et al., 2018b).
The process of non-equilibrium ionization is very much related to the process of irreversible physics via inelastic collision between particles. In the present paper we will restrict our attention to collision impact ionization and radiative recombination. In addition, we assume that the collision between particles can lead to either ionization of recombination, while the processes of collisions leading to excited states of articles will be neglected. This approximation requires that the collision involving neutrals would imply an energy exchange that is at least equal to the first ionization potential, while excited states (if they appear) will have a lifetime that is much shorter than the dynamical scales involved in our problem.
The paper is structured as follows: in section 2 we introduce the necessary equations with their implications and limitations. Given the complexity of the problem, here we are employing a simplified model. In section 3 we study the propagation characteristics of Alfvén waves and investigate the effect of nonequilibrium ionization on the propagation speed of waves and their damping with respect to the collisional parameter between particles for different ionization degrees of the plasma. Using a simple configuration the wave characteristics of decoupled slow waves are studied in section 4. Finally our results are summarized in section 5.

BASIC EQUATIONS AND ASSUMPTIONS
The physical processes of ionization and recombination are far from trivial given the multitude of mechanisms that can result in one or more electrons being removed from a neutral atom and the reversed process of combination between positive ions and an energetic electrons to form neutral atoms. However, for simplicity here we are going to concentrate on collisions, as the main mechanism that can generate ions and neutrals. In order to describe the effects on non-equilibrium ionization and recombination processes let us introduce the quantities Ŵ r q , denoting the interaction rate of process r, affecting fluid q (Leake et al., 2012;Meier and Shumlak, 2012;Maneva et al., 2017). The radiative recombination rate of ions is proportional to the number of ions (with number density n i ) and electrons (with number density n e ) of the system, and it is given by, In the above equation the recombination frequency, R, is given by Cox and Tucker (1969) and Moore and Fung (1972) where X = Aǫ i /T(eV) is a quantity that is defined as the ratio between the ionization potential and thermal energy, and the constant A takes the value of 0.6. In can be shown that at temperatures we can find in the lower part of the chromosphere the electrons in hydrogen atoms are predominantly in the ground state (corresponding to n = 1 energy level), meaning that ǫ i is the first ionization potential, therefore ǫ i = 13.6 eV. Indeed, the relative population at the n = 2 energy level compared to the ground level in the hydrogen atom is given by the Boltzmann equation, where g n = 2n 2 is the statistical weight of each energy level, E 1 = ǫ i and E 2 = 3.4 eV are the energies of electrons on the two levels. In a H plasma at T = 10 4 K we can obtain that N 2 /N 1 = 2.5 × 10 −5 , while at T = 4 × 10 3 K, this ratio is of the order of 10 −13 . We should point out here that the radiative recombination should be treated with care. Under normal circumstances the recombination takes place as a result of the collision between a positive ion (proton in a H plasma) and an electron that has an energy, at least, equal to ǫ i . As a result of this interaction a photon is emitted that can further ionize the neutral H atoms. With a significant amount of neutral H, the emitted photon will have a higher probability of being absorbed by a neutral atom in the neighborhood of emission, with a creation of an ion. Therefore, recombination to the ground state has virtually no effect on the ionization state of the plasma. However, in our analysis we are going to consider that the photon emitted during recombination escapes, this corresponds to the response of optically thin plasma to ionizing radiation. Since we are dealing with spatial scales that are larger than the Debye radius in a H plasma, the plasma can be considered to be quasi-neutral, therefore, n e = n i . In addition, for simplicity we are going to deal with a uni-thermal plasma, where T e = T i = T n = T.
The impact ionization reaction takes place as a result of the collisional interaction between neutrals (with number density n n ) FIGURE 1 | The variation of the real and imaginary part of the dimensionless frequency for slow magnetoacoustic waves with respect to the dimensionless collisional frequency, a 1 , for two different values of the ionization degree of the plasma (color coded). Solid lines represent solutions obtained in the ionization non-equilibrium and dashed lines denote the slow waves in ionization equilibrium. and electrons with energies larger than the ionization potential. As a result, the ionization rate of neutrals is given by, Ŵ ion i = n n n e I = n n n i I, where I is the ionization rate and is given empirically by Cox and Tucker (1969) and Moore and Fung (1972) I = 2.34 × 10 −14 √ Xe X (m 3 s −1 ).
In order to determine the variation of the ionization/recombination rates in the solar atmosphere we use the VAL III C model (Vernazza et al., 1981) to plot the height dependence of these two quantities on logarithmic scale (see Figure 1). It is obvious that the variation with height of these rates follows closely the variation of the characteristic values typical for the atmospheric model employed here. Accordingly, the height-dependance of the recombination rate (R) resembles the height-dependance of the number density, while the ionization rate (I) shows similarity with the variation of the temperature, meaning that these physical quantities are the ones that determine the shape of variation with height. Both rates show an extreme value at the height of about 500 km (bottom of chromosphere, the point corresponding to the temperature minimum in the VALIII C model). The ionization rate overcomes the recombination rate at a height of 2 Mm, a height that corresponds to a temperature of approximately 7660 K.
The other effect we are going to consider is the collision between particles. Assuming a uni-thermal plasma, electrons will have their velocity (thermal velocity) √ m i /m e times larger than the thermal velocity of ions (and neutrals), where m e and m i are the masses of electrons and ions, respectively. During collisions between electrons and neutrals, electrons suffer a large change in momentum, but small change in their energy (given approximately by √ m e /m i ). In contrast, ion-neutral collisions are less frequent than electron-neutral collisions, however, the momentum exchange between these particles can influence the most propagation of waves, therefore, we are going to consider this effect as the dominant collisional mechanism. Collision between ions and neutrals will also ensure that neutrals (that are not influenced by the presence of magnetic field) are kept in the system.
In order to assess the importance of all physical effects included in our model, let us define the characteristic times involved connected to ionization, recombination, and collision as where ν in is the collisional frequency between ions and neutrals defined as, , with σ in = 1.16 × 10 −18 m 2 being the ion-neutral collisional cross section (Vranjes and Krstic, 2013), and k B is the Boltzmann constant. Obviously, the waves we are interested in must have periods that are larger than any of these times. It is clear that the smallest time scale is the collisional time scale (basically τ c gives the time between two consecutive collisions) and for most of lower solar atmosphere this characteristic time is several orders of magnitude smaller than ionization and recombination time. For temporal scales that are near or shorter than the ionneutral collisional time the plasma dynamics has to be described within the framework of two-fluid magnetohydrodynamics, where charged particles (here denoted by index i) and neutrals (denoted by an index n) can have separate behavior, depending on the relative strength of collision. Before embarking on finding the characteristics of waves in non-equilibrium plasmas we need to discuss one more limitation of the problem we are going to consider. In the presence of nonequilibrium ionization and recombination, the linearized mass conservation equations for the two species are written as, where ρ 0i,0n andρ i,n are the background and perturbation values of the densities for ions and neutrals, respectively and v i and v n are the velocity perturbations for the two species. In the above equations Ŵ ion i = −Ŵ ion n and Ŵ rec i = −Ŵ rec n . A direct consequence of the above two equation is that the processes of ionization and recombination are not balanced, particles are created and annihilated during the temporal length of the dynamics we are interested in. In addition, these two equations also imply that under normal circumstances (and in a static equilibrium) the equilibrium density of ions and neutrals should be time dependent quantities, i.e., the background state of the plasma is changing in time. In the absence of this temporal change the ionization and recombination rates must be equal and the system will be in equilibrium. Unfortunately, the description of a physical process when the background is a time-dependent is one of the most complex tasks, as it requires cumbersome mathematics.
We have three possibilities to deal with this problem. The first possibility is to leave the background as time dependent, and determine the evolution of perturbations. This task could be easily accomplished by numerical investigations. The second possibility to deal with the temporal variation of the background that makes the problem mathematically tractable would be to impose the condition that the ionization and recombination rates are much longer than the variation rate of the equilibrium density, that could be occurring over the scales presented by Carlsson and Stein (2002). In addition, temporal variation of perturbations (for instance their frequency) are of the same order as the ionization and recombination rates, so Finally, as a third possibility is to deal with this problem by assuming that non-equilibrium effects appear only in the perturbed state (and this is the possibility employed in the present paper). Similar to the assumption by, e.g., Brandenburg and Zweibel (1995), we can fix the rate of ionization and the recombination rate is chosen in such a way that in the unperturbed stage the right-hand sides of Equations (3, 4) are identically zero. The problem with this approach is that only one of the rates can be realistic, the chosen value of the other rate is artificially imposed. Assuming a VAL IIIC atmospheric model this assumption agrees with realistic values in the photosphere, i.e., in a weakly ionized region of the solar atmosphere. In our study we will employ the later assumption and choose, R = n 0n n 0i I.
With this assumption the characteristic times for ionization and recombination defined earlier become equal. This assumption also means that the empirical formula for R given earlier becomes redundant, as its value will be always given in terms of I, as above The dynamics of the coupled two fluids is given by a set of linearized equations that describe the conservation of mass (given earlier by Equations 3-4), together with a conservation of momentum, induction equation and energy equation (see, e.g., Zaqarashvili et al., 2011;Khomenko et al., 2014;Maneva et al., 2017;Martinez-Gomez et al. , 2017), is the magnetic field perturbation, and m in = m i m n /(m i + m n ) ≈ m i /2 is the reduced mass. Frictions between charged and neutral (close-range interaction) particles is ensured via collisional processes. The above equations must be supplemented by the equation of state for ions and neutrals p i,n = n i,n k B T i,n As explained earlier, the ionization non-equilibrium is present only in the perturbed state, the perturbations we are considering to take place will drive the system out of ionization equilibrium. Equations (5, 6) are the linearized momentum equations of the ion-electron fluid and neutrals, respectively. The terms on their right-hand side express the transfer of momentum between ions and neutrals through the diffusion of one species into the other. As a result of collisions, particles can loose energy and momentum.The same equations reveal an interesting aspect of non-equilibrium plasma. When the background is in equilibrium, in the limit of vanishing collisions (ν in = 0) the only dynamics that can be described is related to ions and the momentum equation reduces to the standard equation used in MHD. Of course, such a limit cannot exist as in that case there is nothing that can keep neutrals in the system and the twofluid description is meaningless. However, in the ionization nonequilibrium the two-fluid description can be applied even in the vanishing collision limit, as the creation and annihilation of ions generates a force per unit volume that acts to keep neutrals in the system.
The above system of equations will be used to study the properties of magnetic and slow magnetoacoustic waves propagating in non-equilibrium plasma and compare these with the values we obtain in the ionization equilibrium.

ALFVÉN WAVES
The simplest wave to study are Alfvén waves, for which the only restoring force is the Lorentz force and the driving force of Alfvén waves is the magnetic tension. We assume a homogeneous equilibrium magnetic field, B 0 , pointing in the z-direction. Alfvén waves will propagate plasma along the magnetic field and they will be polarized in the y-direction. Given the properties of Alfvén waves, the pressure terms in Equations (5-6) are neglected and, therefore, Equations (8, 9) are not needed. As a result, the dynamics of Alfvén waves is described by the system of Equations whereν in = ν in 2 + n 0n I, is the modified ion-neutral collisional frequency by the ionization non-equilibrium effects. Obviously, if I = 0, then we recover the case of Alfvén waves in ionization equilibrium, a case discussed earlier by Zaqarashvili et al. (2011). Let us Fourier-analyse the system of Equations (11-13) considering that all perturbations are proportional to exp[i(kz − ωt)], where k is the real wavenumber along the z-axis. The neutral fluid exerts a drag force against the motion of ions around the magnetic field and, therefore, Alfvén waves will decay. As a result, the frequency of waves, ω, will be a complex quantity, with negative imaginary part of it describing the damping of waves. After simple calculations the dispersion relation of Alfvén waves is obtained to be where χ = n 0i /n 0n ≈ ρ 0i /ρ 0n . In the absence of any partially ionized effects (e.g.,ν in = 0) the dispersion relation would reduce to the standard dispersion relation for Alfvén waves ω = ±kv A . The third order polynomial (14) describes the propagation of two Alfvén waves (propagating in opposite direction) and a third mode that is non-oscillatory (i.e., its frequency has a zero real part). Let us introduce the dimensionless quantities, where n 0 = n 0i + n 0n is the total number density of the plasma and ρ 0 is the total mass density. After some simple calculations the dispersion relation can be written in dimensionless form as, where N = a/2 + ξ nĨ . Luckily some analytical progress can be made by assuming that Alfvén waves will have a small damping. Accordingly, the dimensionless frequency of waves can be written as = r + i i (with both r , i real quantities) and we write | i | ≪ | r |. Focussing only on the forward propagating modes, the real and imaginary part of the frequencies can be easily obtained as, and i = −ξ 2 n a + 2ξ nĨ It is obvious that the imaginary part of the frequency is negative, meaning that waves will damp, regardless what is the strength of collisions between particles or degree of ionization. The corresponding values for ionization equilibrium can be found once the quantityĨ is set to zero. In ionization equilibrium the imaginary part of the frequency tends to zero when the collisional frequency, a, is set to zero, i.e., Alfvén waves in collisionless plasma do not damp. In this case the governing equations describe the dynamics of ions alone, while the dynamics of neutrals becomes undetermined. Since ions and neutrals do not interact through collisions, the possibility of having neutrals in such systems becomes a problem (technically speaking the plasma can be considered collisionless when the mean free path between collisions is much larger than the lengths over which the plasma quantities vary). In the limit of strong collisions (a ≫ 1) the imaginary part of the frequency becomes smaller and at a → ∞ the mixture of ions and neutrals behaves like a single fluid with no damping; this limit corresponding to the MHD limit. Before further discussion we need to clarify that the term "collisionless" used above refers to the case when the collisional frequency between ions and neutrals becomes zero. However, this does not exclude the possibility of having collisions between ions and electrons and neutrals and electrons, i.e., the collisions that affect the ionization state of the plasma. In contrast, in the ionization non-equilibrium, when we set the collisional frequency to zero we arrive to a critical damping rate, meaning that Alfvén waves will damp even in the absence of collisions between ions and neutrals and the ionization/recombination processes will ensure that the mixture of ions and neutrals stay coupled. The damping whose rate is given by Equation (19) constitute a new damping mechanism that roots itself in the additional drag force generated by different drift velocities of ions and neutrals during the process of ionization and recombination. For the sake of completeness, we should mention that Landau damping is also a mechanism that appears in collisionless plasmas, and this damping is not associated with an increase in entropy, and therefore is a thermodynamically reversible process. It remains to be seen whether the damping mechanism described in the present study has the same properties as Landau damping. The damping we discuss here can be due to the fact that with the number of ions changing in time, more and more ions become attached to FIGURE 2 | The variation of the recombination rate (R, shown in red line) and the ionization rate (I, blue line) with height (on logarithmic scale) in the lower solar atmosphere based on a VALIII C solar atmosphere model (Vernazza et al., 1981).
magnetic field lines, increasing their inertia, eventually causing the attenuation of Alfvén waves. Let us investigate the variation of the real and imaginary parts of the dimensionless frequency of forward propagating Alfvén waves in terms of the dimensionless collisional frequency, a for three different values of the relative neutral density (ξ n = 0.2, 0.5, 0.8). Clearly ξ n = 0 describes a fully ionized plasma, while the limit ξ n = 1 corresponds to a completely neutral fluid. The upper panel of Figure 2 shows the variation of the real part of the dimensionless frequency (solid lines) for the three values of ξ n and the corresponding dispersion curves corresponding to the ionization equilibrium (shown by dashed lines). Here the dimensionless quantity can be understood as the propagation speed of Alfvén waves in units of Alfvén speed, v A .
First of all it is clear that for very large collisional frequency the propagation speed of Alfvén waves and damping rate in the two regimes become identical and in a strongly collisional plasma the propagation speed of Alfvén waves and their damping rate are not influenced by collisions between ions and neutrals. Secondly, Alfvén waves in ionization non-equilibrium will propagate with a lower speed and the difference in propagation speed increases with the number of neutrals in the system. Furthermore, the more ionized the plasma is, the faster these Alfvén waves will propagate. Similarly, as long as the collisional frequency is larger than ( Ĩ2 + 4 −Ĩ)ξ n , Alfvén waves in ionization equilibrium will have a larger damping rate (smaller damping time) than those waves that propagate in a plasma in ionization non-equilibrium, however, these differences are not significant. For collisional frequency smaller than this threshold value, Alfvén waves damp much quicker in a non-equilibrium plasma. As pointed our earlier, the damping rate of Alfvén waves in equilibrium plasma tends to zero in the collisionless plasma. In contrast, in a nonequilibrium collisionless plasma, Alfvén waves damp with a rate given by cr i . In conclusion, in a strongly collisional plasma the propagation and attenuation of Alfvén waves is independent whether the plasma is in ionization equilibrium or not. Here the collisional time is at least one order of magnitude larger than the characteristic time for ionization and any non-uniformity that could potentially influence the characteristics of waves is smoothed out by collisions.
Finally we need to mention that Alfvén waves in the present configuration become dispersive, and dispersion is proportional toν 2 in . Figure 2 also shows that for a given collisional parameter, a, waves with shorter wavelength will propagate slower and the larger the amount of neutrals in the system, the more dispersive Alfvén waves are.

SLOW MAGNETOACOUSTIC MODES
To be able to decouple magnetoacoustic modes, we assume that slow magnetoacoustic waves propagate along the background magnetic field and the dominant dynamics occurs in the z direction. Since the species of the plasma have the same temperature, we can write Since slow waves propagate along the background magnetic field we assume that all perturbations are proportional to e i(kz−ωt) . This particular choice for the form of perturbations and background magnetic field can reduce the slow waves to acoustic modes. The dynamics of linear slow waves propagating along the magnetic field in a partially ionized plasma in the presence of non-equilibrium ionization is described by the system of Equations (3-10) where perturbations are considered to be proportional to the exponential factor introduced above. Accordingly the set of equations used are, These four equations can be reduced a set of coupled equations for the z-components of ion and neutral velocities. The compatibility condition of this system gives us the dispersion relation, where the coefficients A, B, C, and D depend on the parameters of the problem. This fourth-order polynomial in ω describes the propagation of two families of waves: one associated to ions, and the other one to neutrals. Due to collisions and nonequilibrium effects, the two kinds of modes are coupled. Again, let us introduce similar dimensionless quantities as in the case of Alfvén waves, but now we use the quantity kc Sn to write variables and constants in dimensionless form. As a result, the dispersion relation reduces to, where the constant coefficients are given by, where the parameters a 1 andĨ 1 are defined in the same way as parameters a andĨ given by Equation (15), but the quantity that is used to write them in dimensionless form is kc Sn . In the absence of collisions and ionization non-equilibrium the two pairs of sound waves are propagating with the sound speed of ions and neutrals, with the wave associated to ions propagating faster. Figure 3 shows the variation of the dimensionless frequency of slow waves (or the phase speed of slow waves in units of the neutral sound speed) with respect to the dimensionless collisional parameter, a 1 , for two distinct values of the ionization factor, ξ n (0.2 and 0.8, respectively). Due to the collision between particles and non-equilibrium effects, the quantity is complex, where the imaginary part describes the temporal modification of the amplitude of waves. The upper panel of Figure 3 shows the real part of , while the lower panel shows the imaginary part of this quantity. Similar to the plot obtained for Alfvén waves, the solid lines denote the dispersion curves for slow waves in ionization non-equilibrium, while dashed line correspond to the case when the plasma is in ionization equilibrium.
It is clear that for the whole spectrum of parameters, the imaginary part of the frequency will be negative, meaning that slow waves will damp; this result is in agreement with the conclusions of previous studies by, e.g., Braginskii (1965) and Zaqarashvili et al. (2011). The two sets of waves shown in Figure 3 have different behavior in terms of collisional frequency. Let us concentrate first on the dispersion curves obtained for the real part of the frequency. Depending on the amount of neutrals in the system, the ion-acoustic waves can have an enhanced propagation speed in the limit of weakly collisional plasma. For ξ n = 0.2 (black curves) the slow wave that corresponds to an ionization equilibrium travels slightly faster, however very quickly the propagation speed of these slow waves in the two regimes is equal. When the amount of neutrals is increased (red lines) we can see that the wave corresponding to the nonequilibrium state propagates faster but at the value of a 1 ≈ 2.5 the two speeds become identical. At this point the collisions between ions and neutrals become so frequent that their effect can overcome any modification due to the additional change in the ionization degree of the plasma.
On the other hand the slow waves associated to neutrals display a completely different behavior. In a plasma with ξ n = 0.2 these slow waves have a smaller phase speed, and the speed in a non-equilibrium plasma is slightly larger than the corresponding value obtained in an equilibrium state. Once the amount of neutrals is increased (the curves corresponding to ξ n = 0.8) the slow waves in an equilibrium plasma is larger, and this relation is maintained again, until the dimensionless collisional frequency, a, reaches the value of a 1 ≈ 2.5, after which, for a given value of a 1 , the neutral slow waves in equilibrium plasma propagates faster. It is also clear that at a 1 ≈ 2.5 the propagation speed of ion-slow waves and neutral-slow wave is equal. Similar to the results obtained by Zaqarashvili et al. (2011), the neutral-slow waves can exist up to a certain level of collisional rate and the value of collision where these waves cease to exist increases with the amount of neutrals in the system. At this point the collision between neutrals and ions is so frequent that the mixture of charged and neutral particles starts behaving like a single fluid. Increasing the frequency of collisions between ions and neutrals cause a slow down of neutral-slow waves. Now let us concentrate on the imaginary part of the frequency. It is clear that regardless what is the ionization degree of the plasma, the damping rate of ion-acoustic modes is very small, practically these waves can propagate with no attenuation for any value of collision. In contrast, the neutral-slow waves have a very different behavior. These waves show a very strong damping (with damping times of the order of the period of waves, or larger). For both ionization degrees chosen here, the damping time of neutral-slow waves is decreasing with the increase in the collisional frequency. Comparing the results we obtained in the two regimes, it is clear that the neutral-slow waves in an ionization equilibrium plasma have no damping in the absence of collision (as we would expect), however, in the presence of nonequilibrium these waves decay due to ionization/recombination processes and-in the absence of collisions-the damping time of these waves is independent on the amount of neutrals in the system. Once the collisional rate in increased, the mode that corresponds to smaller amount of neutrals will have a shorter damping time.
Similar to Alfvén waves, slow magnetoacoustic modes will also be dispersive, however the ion-acoustic modes will be practically non-dispersive (similar to the slow waves in the MHD description), while neutral-acoustic modes are strongly dispersive, again waves with larger wavelength traveling faster. Figure 3 also shows that in the case of neutral-acoustic modes there will be always a critical wavenumber above which these modes do not propagate.
It is very likely that the problem of propagation and decay of these waves will display a different behavior once the full polarization in the xz plane is considered. In that case slow waves will have a magnetic component and it remains to be seen how these waves will behave for different value of plasma-beta. However, this consideration requires numerical investigation.

CONCLUSIONS
Given the nature of partially ionized plasma in the lower solar plasma and prominences requires a different approach. For particular high frequencies range a two-fluid description is needed and neutrals can considerably influence the properties of waves. This framework was used to study the characteristics of Alfvén and slow waves propagation along a unidirectional homogeneous magnetic field. The novelty of our research resides in consideration of the effects of ionization non-equilibrium, i.e., the case when the rates at which neutrals are ionized through impact ionization and ions recombine with energetic electrons through radiative recombination are not equal. Our results show that this effect is more important when the collisional frequency is comparable with the frequency of waves. Any information about the existence of two fluids and the drag forces exerted by neutrals upon ions is lost in the case of strong collisions. Here the plasma behaves like a single fluid and the dynamics can be confidently described within the framework is MHD. We ought to mention that the nonequilibrium applies only to the perturbed state of the plasma; in the background state the plasma remains in ionization equilibrium. This assumption was needed to be able to make analytical progress.
Using a simple configuration we studied separately Alfvén and slow magnetoacoustic modes. The collision between heavy particles and non-equilibrium effects renders the frequency of waves to be complex, where the imaginary part of it describes damping. One of the main results of our investigation is that in a plasma that is in ionization non-equilibrium waves can damp even in the collisionless limit. In this case waves will damp due to the drag forces by neutrals that appear due to the ionization non-equilibrium. Given that waves are damped even in the collisionless limit makes us to think about the physical explanation of this effect to be similar as the theory of Landau damping. However, this statement needs to be investigated properly in the future. The processes of ionization and recombination will increase the degree of entropy in the system that could cause the additional damping.
Finally we should mention our results should treated with precaution when making far reaching conclusions valid for the whole photosphere and/or chromosphere. For simplicity the present study assumed collisional ionization and radiative recombination as the dominant mechanisms that determine the non-equilibrium, however, in reality this can change from region to region in the solar atmosphere. It is very likely that in the dense photosphere photoinduction and threebody recombination are more important effects. In addition, most of the partially ionized plasma is optically thick, in which case a similar treatment as presented by Ballester et al. (2018a) is needed. Moreover, the relatively low temperature plasma in these regions means that the gravitational scaleheight is short, therefore, in the case of waves traveling over long distances in the solar lower atmosphere, the gravitational stratification could influence the propagation of waves and work against the damping of waves. It remains to be seen how a realistic variation of ionization and recombination rate together with stratification will affect the characteristics of waves.