Nonlinear Interaction of Low-frequency Alfvén Waves and Ions

The interaction between Alfvén waves and particles is a critical phenomenon in space and laboratory plasmas, and it has been observed that low-frequency Alfvén waves can accelerate and heat ions through subharmonic resonant interactions. In this study, we use test particle simulations to verify the nonlinear heating of parallelly propagating low-frequency Alfvén waves on ions and analyze the underlying process in terms of the Poincaré section. Our results demonstrate that low-frequency Alfvén waves can periodically pick up ions, leading to oscillations of average parallel velocity and temperature of the plasma by phase mixing, ultimately resulting in the stabilization of acceleration. Furthermore, we have developed an analysis that can estimate the time required for heating and accelerating. In the case of multiple waves, our findings indicate that the presence of more chaotic modes does not necessarily result in better wave heating. We have also discussed the effect of random phases on the heating process. Overall, this research sheds light on the crucial role played by the interaction between Alfvén waves and particles in astrophysics and provides new insights into the mechanisms underlying the heating and acceleration of ions through subharmonic resonant interactions with low-frequency Alfvén waves. These findings may have significant implications for the understanding of plasma dynamics in a range of astrophysical environments.


Introduction
The plasma heating by the Alfvén waves plays a vital role in both space and laboratory plasmas. The plasma heating by the Alfvén waves is generally believed to be one of the significant possible dissipative mechanisms for explaining the dynamics of coronal heating and solar wind acceleration (Nekrasov 1970;Lieberman & Lichtenberg 1973;Hollweg 1978;Abe et al. 1984;Markovskii & Vasquez Bernard 2011). The lowfrequency Alfvén waves are generally considered to have no strong interaction with ions to accelerate or heat plasma because their frequency, which is lower than the ion cyclotron frequency, does not meet the cyclotron resonant condition (ω − k ∥ v ∥ = nΩ p , where n is an integer; Wu et al. 1997;Chen et al. 2001;Guo et al. 2008). The theoretical analysis and numerical experiments done by Wu et al. (1997) show that ions can be picked up by the large-amplitude low-frequency Alfvén waves through nonlinear interactions even if the cyclotron resonance condition is not satisfied. Similarly, Chen et al. (2001) and White et al. (2002) studied heating by the obliquely propagating low-frequency Alfvén waves and found that significant perpendicular stochastic heating can be obtained with a sufficiently large wave amplitude. Kolesnychenko et al. (2005) demonstrated that the nonlinear subharmonic resonance exists for any k ⊥ and not for k ⊥ = 0 only, and it can decrease the stochasticity threshold when k ⊥ ≠ 0. Lu & Chen (2009) extended the work of Kolesnychenko et al. and showed that the increase of frequency, propagation angle, and number of modes can decrease the threshold when spectrum waves are taken into consideration. Physical processes in the space like the ion's velocity distribution and the temperature anisotropy in solar wind, heating of the corona, and acceleration of solar wind have been explained in many ways (Wu et al. 1998;Johnson & Cheng 2001;Yoon et al. 2009;Li et al. 2010;Wang et al. 2011;He et al. 2016;Liu et al. 2016). Karney was the first to use the Poincaré plot to study the stochastic state of particles during interactions with waves (Karney 1978). Chen et al. discussed the oblique propagating waves with the help of the Poincaré plot. However, no one has studied the parallel-propagating waves with the Poincaré plot. The motion of ions in the parallel-propagating wave field is simple and regular without carrying any chaos. For multiple modes, Lu & Chen (2009) discussed the influence of various parameters on the threshold of stochastic heating. However, the actual heating effect by multiple waves has not been discussed.
In the present study, we investigated the nonlinear interaction of ions with low-frequency Alfvén waves propagating along the background magnetic field by a test particle simulation. The Poincaré section plot demonstrates the trajectory in lower-dimensional phase space. For the multiple waves, we discuss the effect of the random phases considering the identical frequencies. These results will benefit the research and applications of the nonlinear interactions between the lowfrequency Alfvén waves and ions.
The structure of this paper is as follows. In Section 2, we briefly describe the simulation model used. The features of ions interacting with the single-and multimode waves are presented in Section 3. In Section 4, we summarize the results.

Simulation Model
The left-hand circularly polarized waves that propagate parallel to the ambient magnetic field are used in this work.
Moreover, the heating and acceleration of plasma caused by waves are studied. We assume that the background magnetic field B 0 = B 0 i z , and the wave is nondispersive with ω is the local Alfvén wave speed. Thus, the electromagnetic field of the wave can be expressed as follows: here i x , i y , and i z are the unit vectors of the three dimensions. ψ j = k jz z − ω j t + f j is the phase of the wave, and f j denotes the random phase of the mode. The particle equation of motion in the wave field can be summed up as follows: where the subscript i indicates ions labeled as i. The coupling between the variables of the entire system shows a strong nonlinear effect, which cannot be solved by the general analytical methods. The acceleration and heating mechanism of ions in the parallel-propagating Alfvén wave field is explained by the mapping on the Poincaré section of the system mentioned above.
To obtain the Poincaré section of the ions, we used a test particle simulation to study the differential equation system by ignoring the influence of particle motion on waves. The Boris algorithm has been used to solve equations. The time step was t 0.025 P 1 D = W -. Initially, all particles were placed in the same position, which was the center of the simulation space. We also have constrained the initial position distribution of ions to a fixed range centered around the simulation center, to investigate the effect of the initial position distribution of ions on the subsequent perpendicular heating. The length of the space in each direction was v 40,000 p A 1 W -. For the sake of simplicity, though, the simulation is one-dimensional consider-ingonly changes in positionalong the background magnetic field, which is alsothe direction of wave propagation.

The Heating and Acceleration of Ions by Monochromatic Alfvén Wave
From the work of Wu et al. (1997), different initial velocities of ions lead to different acceleration effects. In their work, only the average parallel velocity versus time was shown. It is still an open question what the nature of ion motion at the microscale level is. We used the Poincaré plot to see the motion of ions at the microscale level by taking the points v y = 0 and  v 0 y ¢ . The Poincaré plot can reveal many properties of orbits of the system in low-dimensional phase space (Karney 1978;Karimabadi et al. 1990). The visible trajectory of each ion means that the motion of ions is regular. If the motion of ions is irregular, the Poincaré plot is in chaos. Figure 1 is a typical Poincaré plot of ions in parallelpropagating Alfvén waves with different amplitudes. The k z t z y w ¢ = denotes the particle's phase in the wave field. A region of acceleration, but with no net energy gain for each individual ion, exists on the Poincaré plot between phases and 1 2 1 2 ( ) ( ) p p -, irrespective of the amplitudes. Considering that the initial velocities v x0 and v y0 of ions are zero, the ions maintain their initial velocity when they are in the phase interval 1 2 ( ) p p --and 1 2 ( ) p p -, but this rule will be breached while the initial perpendicular velocity is not zero. When ions move into the phase interval 0 acceleration in both the x-and z-directions can be observed. Additionally, the value of acceleration is proportional to the difference of local Alfvén wave speed and initial parallel velocity, which can be represented mathematically as Δv ∝ v A − v z0 , and then the velocity obtained decelerates to the initial state in the phase interval 0 The trajectories in the Poincaré plot reflect the changes of ion velocities in the phase of the wave by selecting the point v y = 0, which are projections of phase-space trajectories of particles on the Poincaré section we selected. Therefore, they can be discrete curves, such as the top one in panel (a1) of Figure 1. This indicates a deeper connection between the variations of perpendicular and parallel velocities.
The motion of a single ion in the lab frame is more complex as shown in Figure 2. Two typical characteristic frequencies appear with the velocity along the x-and y-direction owing to the simultaneous modulation of the wave field and background magnetic field. Furthermore, the higher characteristic frequency is five times greater than the lower one. The trajectory within the plane perpendicular to the background magnetic field comprises two distinct cyclotron motion paths with differing radii. These circular motion paths correspond to two frequencies of perpendicular velocity mentioned above. The largeradius one is controlled by the wave field, and the small-radius one is controlled by cyclotron motion corrected by the wave field. Consistent with the Poincaré plot for the parallel motion, the velocity in the z-direction periodically changes with only one frequency in response to the cyclotron motion modified by the wave field, not controlled by the wave field directly. It leads an ion's pickup motion by waves in the direction along background magnetic field, which refers to the periodic acceleration of ions in the longitudinal direction due to the Lorentz force. Regions corresponding to acceleration, deceleration, and constant speed on the Poincaré plot can be found separately. Further analysis tells us that the changes of zdirection velocity stop during the period of small radius cyclotron motion. The points, which we selected for printing on the Poincaré plot, are the end of cyclotron motion with the initial state when they are located at the upper semicircle, while they are undergoing the cyclotron motion with an accelerated state at the lower semicircle. That would explain the shape of the trajectories in the Poincaré plot. It should be noted that the Poincaré plot represents a low-dimensional projection of a higher-dimensional phase space, capturing the state of particles as they traverse a selected section along a given direction. In Figure 2, these states ( ) are denoted by blue circles. These circles are distributed homogeneously along a large circle in panel (a), corresponding to the phases of the wave field. As more points are recorded over a sufficient period, the Poincaré plot appears as a solid line owing to the homogeneously distributed circles. However, the density of the circles in panel (a) may sometimes appear sparse, resulting in a dotted line in the Poincaré plot, as seen in Figure 1. This effect can be mitigated by increasing the duration of the simulation.
For the situation when the initial perpendicular velocity is not equal to zero, a more exciting result in the Poincaré plot is demonstrated in Figure 3. The most distinct feature is that some ion trajectories in the v ⊥ /v A versus y¢ plot are broken at the parts v x0 → 0 − and v y0 > 0 for v x0 ≠ 0 and v y0 ≠ 0 cases, respectively. Particles with initial velocities v x0 > 0 and v y0 < 0 have unbroken lines in the Poincaré plot. The phase island arises for v y0 > 0, and chaos arises at the edge for v x0 < 0. The  relationship between the states of v x0 and v y0 can be likened to conjugation, as their respective Poincaré plot patterns appear to be upside down. If v x0 ≠ 0, the ions with initial velocity above zero will obtain stronger acceleration, but the contrary situation occurs when v y0 ≠ 0. With the further analysis, we discovered that ions in the wave phase interval tend to keep the initial velocity state. The wave field affects the motion of ions in phase intervals similar to the abovementioned nonzero parallel velocity case. However, ions with a nonzero perpendicular velocity can be trapped in the phase space to trace a circle at 2π phase in the Poincaré plot. This suggests that the subharmonic resonance appears in the direction perpendicular to the background magnetic field. This gives a likely reason why stochasticity always occurs with the presence of vertical wavevectors in Kolesnychenko et al. (2005). Combining Figure 1 with Figure 3, we can conclude that initial velocity along the background magnetic field controls the energy obtained from waves. The waves parallelly propagating affect the ion's perpendicular velocity first. Then, the change in vertical velocity is transferred to parallel velocity owing to the conservation of energy in the wave frame. In real space, parallel-propagating waves have a perpendicular electric field that induces a perpendicular velocity change proportional to |v A − v ∥ | because of the Doppler shift and Faraday's law (Equation (2)). The magnetic field of waves and the changed perpendicular velocity modify the cyclotron motion, which can accelerate parallel motion periodically for ions in the system. In the absence of a wave field, the ion's cyclotron motion would complete a full circle. However, in the presence of a wave field, the ion just undergoes a semicircular trajectory upon returning to its initial state. This results in a breaking of the conservation of magnetic moment and leads to parallel acceleration.
All the above analyses are carried out on a single particle at the microscale, and these changes have a significant impact on the average velocity and temperature of plasma. Figure 4 shows that the periodic variation of the velocity of individual ions during wave−particle interactions drives the transition of plasma macroscopic parameters from oscillatory behavior to equilibration. The average energy is utilized in this study to demonstrate plasma heating by wave, which is as follows: where T k is the "equivalent temperature," k B is Boltzmann's constant, and the · á ñ denotes an average over all particles. E k was defined as the "ion kinetic temperature T kin " by Wang et al.  Asymptotic values of the E ⊥/∥ and U ∥ can be found in Figure 4 for all the cases. This is consistent with the previous work. Because of the different initial velocities of ions, each ion's phase in the wave field is much more distinguishing after a while. Hence, a phase mixing emerges (Kanekar 2014), and the oscillations of average parallel velocity and temperature decay with the phase mixing. We can infer that the acceleration of plasma by a wave is proportional to the value of v v A  -á ñ based on the single ion behavior in the wave field. The oscillation periods of plasma macroscopic parameters in Figures 4(b) and (d) are obviously shorter than those in Figures 4(a) and (c), while panels (b) and (d) have a higher wave frequency. Additionally, it can be concluded that the increase in wave amplitude leads to an increase in the oscillation periods in comparison with the panels of Figure 4 within each column. After the oscillation, a significant increase in plasma temperature can be found, and we can define this oscillating process of macroscopic quantity as the period of "nonresonant heating." During nonresonant heating, the heating is mainly caused by the phase mixing process and periodic pickup of ions with the existence of the wave field.
When most ions have experienced several periods, the parameters of ions in the parallel direction get stable. However, the time to get the asymptotic value of E ⊥ is longer than that of U ∥ and E ∥ . This time difference can be observed in Figure 4. Furthermore, a wave of lower frequency and smaller amplitude leads to a longer lagging time. Since the lower-frequency and small-amplitude wave indicates slower and weaker movement of ions, the lagging time will be longer naturally. Besides that, we also guess that lagging time is relevant to the initial position distribution, which also can be the initial phase distribution in the wave frame. A modified version simulation has been carried out to verify this idea, where we distribute the ions in a range of position (phases). Results are demonstrated in Figure 5. Plasma energy in the orthogonal direction has different behaviors with the initial position range of ions distributed. The overlapped parallel energy curves further prove that the velocity variation in the parallel direction is not affected by the wave field mentioned above, while a broader position region that ions distributed leads to a better heating in the perpendicular direction during the nonresonant heating process. There is a threshold of initial position range that leads to a sufficient nonresonant heating energy to make the lagging time vanish. In this example, a wavelength is v v 2 120 p A 1 A 1 p w » W --, which equals the initial position range of the best perpendicular heating case. Consequently, we can infer that one wavelength is the threshold value of the initial position range to produce the distinct behaviors of the two orthogonal directions. Furthermore, regardless of any lagging time, all cases shown in Figure 5 converge to the same asymptotic value. Combined with the result of Figure 2, it appears that the motion of each ion is determined as it enters the wave field, leading to predetermined plasma energy acquired through wave−particle interaction. If the initial position range is sufficiently small, the heating process from the end of the oscillation period to the equilibrium of the predetermined energy in the perpendicular direction can be referred to as "stochastic heating," which is dominated by stochastic scattering of ions in velocity space. It is evident in Figures 6 and 7 that particles are scattered by the wave in v ∥ -v ⊥ space during the period of stochastic heating. FromFigure 2,we observe that the two vertical velocities of a single ion appear to be conjugate, with their time-averaged value being zero. In contrast, the time-averaged value of the We put all the particles at the same position (center of the simulation space) initially, and the wave's random phase is zero. parallel velocity is finite. As a result, the velocity distribution function eventually becomes Maxwellian, with a zero average velocity in the perpendicular direction and a finite average velocity in the parallel direction. This can be proved by panel (b2) in Figures 6 and 7. Regardless of the presence of stochastic heating, the shapes of distribution functions with the same parameters of waves are consistent during the heating process, while a shift arises in the distribution of v ⊥ that makes the averaged velocity nonzero. This suggests that nonresonant heating with a small range of ions' initial positions produces a nonzero v ⊥ , and then this unstable distribution causing the scattering of ions in velocity space leads to "stochastic heating" and the average velocity of zero in the perpendicular direction.
The initial position range is sufficiently small, which may be rare in actual phenomena such as the solar wind. Instead, what can predominantly be observedin most cases is the process of "nonresonant heating".However, here we discuss the "stochastic heating," which can help us understand how the heating in the two different directions takes place. These results suggest that the interactions that take place at the two directions are  different processes, which verified our previous conclusion that the heating and acceleration of the parallel direction are caused by the periodic pickup during the modified cyclotron motion, but perpendicular heating is not only influenced by the modified cyclotron motion but also affected by the wave field directly. The parallel direction experiences a higher-order effect compared to the perpendicular direction. This can beexplained in physical terms as to why nonresonant heating results in a temperature anisotropy (Liu et al. 2016). Nonresonant heating is a vital process during the wave −particle interaction that can be inferred from the above analyses. We study the relationship between the nonresonant heating ending time and amplitude and frequency of waves. The nonresonant heating ending time is inversely proportional to frequency and directly proportional to amplitude (from Figure 8(c)). An empirical formula (Equation (7)) has been carried out by PySR (Cranmer 2023), a high-performance symbolic regression in Python, as follows: The simple formula alignswell with the real dataand the mean absolute percentage error (MAPE) is 10.6706%. This study can help us understand the process of wave−particle interaction more deeply.

The Heating and Acceleration of Ions by Multiple Alfvén Waves
Based on Kolesnychenko's work, particle motion is stochastic at lower wave amplitudes with a nonmonochromatic wave spectrum, and it is true even when the waves propagate along the ambient magnetic field. Because of the different ψ k of each wave, the original method cannot be used to compute the phase in waves: (8) has been used to compute the phases of ions in the wave field according to the superposition principle of waves, where N is the number of wave modes. We take the v ∥ /v A versus y¢ and v ⊥ /v A versus y¢ plots to construct the Poincaré plots depicted in Figure 6. The frequencies of the wave modes follow the uniform distribution in ω 1 ∼ ω N . The amplitude satisfies the relation B B =is the generally accepted value for the power spectrum of magnetic fluctuations in the solar wind. Figures 9(a2), (b2), and (c2) show the Poincaré plot of the multimode wave case with the same frequencies. The pattern is similar to that of the single-wave case. The superposition of multiple waves with the same frequency is a single wave with the original frequency but with an increased amplitude. This similarity indicates that our computation of y¢ is practical. From Figures 9(a2)-(c2), the maximum parallel velocity that a W -, which indicates no "stochastic heating" process. The other parameters were set to the same values as those used in Figure 6.  single ion can achieve is increased compared to the N = 2 case, but it is almost the same between the N = 20 and N = 50 cases. This means that the heating effect has a saturation for the mode number. When a broadband frequency spectrum is taken under consideration, as shown in Figures 9(a1), (b1), and (c1), chaos occurs in the Poincaré plot. In addition, the region where chaos exists expands with the increment of the mode number N. The maximum parallel velocity that a single ion can reach is close for different wave mode numbers. However, the increased randomization with increasing modes leads to better heating than the two modes based on previous work. It is not consistent with the temperature evolution presented in Figure 10.
Distinguishingthe time evolution in the single-wave case presented in Figure 9, the oscillations turn to random amplitudes similar to turbulent wave fluctuations. A notable frequency can be found in the evolution process, which corresponds to the lowest frequency of the superposed wave. The stochastic heating period in the single-wave case has been covered by intensive turbulence. Plasma energy and acceleration obtained by the interaction when N = 5 are more than when N = 2 presented in Figure 10(a). This is consistent with the work by Lu & Chen (2009), except the heating with N = 10 and N = 20, which is inconsistent with the intuition that the energy obtained increases with the number of modes. Lu's work reported in 2009 showed that the heating effect with N = 5 is much better than with N = 2 presented in Figure 11 of Lu & Chen (2009). However, the random phases of the waves are not shown there. Only the phases of the two modes 0°and 30°are given in that report. When we adjust the f N from 2π to π /6, we find that the energy obtained from the waves by heating decreases with the increased number of modes. The phenomenon of a greater number of wave modes lowering the heating energy can be interpreted as follows. The increase in wave modes with a characteristic power spectrum index leads to more intensive turbulence. Thus, the motion of ions in the wave field is more stochastic. However, the stochastic region in the Poincaré plot has no significant enhancement in heating. This means that the plasma heating is less effective in more chaotic cases than in regular cases. Figures 9 and 10 show that the heating depends on the achievable maximum parallel velocity of the ions, which is controlled by the intensity of the wave field. Moreover, the magnetic field with more modes is weaker than the two modes at most positions along the background magnetic field, as shown in Figure 11. Figure 11 illustrates that, in most positions along the background magnetic field, magnetic fields with more modes are weaker than those with two modes. This phenomenon is due to the fact that as the number of wave modes increases, the distance between adjacent wave packets becomes greater. Therefore, although the magnetic field at the wave packet location may be strong, the majority of the magnetic field intensity in space is very weak when the number of wave modes is sufficiently large. This trend of the magnetic field spatial distribution can lead to chaos in the Poincaré plot, and it is one of the reasons for the inefficient heating.
The contents of Figure 10 can be partly explained by the different field strength mentioned above. However, the random phases also influence the heating process and can be observed from the comparison of Figures 10(a) and (b). The effect of random phases on heating and acceleration should be discussed in detail. It is assumed that the effects of the two wave modes are entirely independent. This means that there is no wave −wave interaction. The phases of wave modes determine the position of particles in each different wave field. Moreover, the single-wave Poincaré plot indicates that the different position in the field leads to different velocity states. Waves with the same frequency have been studied because the heating of waves with a spectrum has many influential factors to distinguish the effects of random phases.
A linear relationship between parallel plasma velocity and the number of modes can be found in Figure 12. This shows that more waves mean better heating and acceleration when f N π . A higher amplitude of an equivalent single wave generated by superposed fields of multiple waves makes the relationship valid. This means that the opposite phenomenon observed in the case of N 6 f = p with a frequency spectrum presented in Figure 10(b) is mainly due to the nonlinear wave −wave interaction. This interaction is coherent and linear when the two identical frequencies waves superpose. Therefore, the curve for N = 2 is symmetrical and the axis of symmetry is f N = π . Nonlinearity appears for the N > 2 case. The curve is no longer symmetrical, and the inflection point has an apparent rightward shift with the increase of wave modes. The nonlinearity of multiple-wave interaction is an essential topic in plasma physics. The significance of nonlinearity shown in our work means that the further study of nonlinear wave−wave interaction should be taken into the particle−wave interactions.

Summary
Extending the work of Wang et al. (2011), we studied the motion of ions interacting with the parallel-propagating Alfvén wave. Using the Poincaré plots to examine the ion motions in low-dimensional phase space, the calculated regular trajectories agree with the corresponding trajectories reported by Kolesnychenko et al. (2005). Most previous works use the Poincaré plot to study the chaos in interaction of waves and particles, but if we treat the Poincaré plot as a projection of ions' motion in phase space on a selected section to study the motion of ions, it also can be helpful for the parallel-propagating wave research. There is no obvious chaos structure in the Poincaré plot of parallel-propagating waves; however, the process of plasma acceleration and heating is indicated. Some significant features of the macroscopic physical quantity in this process are inferred from the Poincaré plot and the motion of a single ion on amicroscopic scale.
The waves with a spectrum have been considered in this study as well. The stochastic regions appear with the spectrum, which does not effectively enhance the heating in this work. The stochastic region or "chaos structure" in the Poincaré plot of oblique propagating waves can obtain a significant perpendicular stochastic heating, which is different from our results. The stochasticity in our multiwave study is caused by the turbulent fluctuation of the magnetic field. However, the chaos in the oblique propagating wave case can undergo a complete transition from regular to chaos with the increasing of wave amplitude. Although the chaos existed in the motion of ions, a sufficiently weak magnetic field generated by multiple waves superposed still can make the heating effectless. The effect of the random phases of each wave should be discussed as well. However, the heating of plasma by waves with a spectrum is so complicated that we discussed the identical frequency wave modes instead. An apparent regularity and a nonlinearity were found in this case.
Distinguishing from the works of Chen et al. and Lu & Chen (2009), we just discussed the parallel-propagating wave, whereno global chaos will appear on the Poincaré plot. In our situation, there is no threshold amplitude for the "stochastic heating" caused by the stochastic ions' trajectory; the relevant theory has been shown in the work of Kolesnychenko et al. (2005). Indeed, the inclusion of "stochastic heating" in their work only accounts for the heating process, while significant acceleration in the parallel direction is not considered. Through a comparative analysis of Figures 3 and 4 from the work of Kolesnychenko et al. with the Poincaré plot presented in our study, one can also observe the presence of the "nonresonant heating" process when the waves propagate obliquely, and it significantly contributes to the energy gained by the plasma from the waves.
As a potential mechanism for corona heating, it has been considered that Alfvén waves can only heat low-energy plasma through the cyclotron resonance when the frequency satisfies certain conditions. But the widespread low-frequency waves on the Sun surface suggest that interaction of ions and lowfrequency Alfvén waves is vital for the corona heating. Our results will be meaningful in helping to understand the process of plasma heating by low-frequency Alfvén waves through the interaction of ions and waves.