Quantitative and qualitative analyses of the mKdV equation and modeling nonlinear waves in plasma

In this paper, nonlinear electrostatic structures on the ion time scale in plasma consisting of two populations of electrons (cold and hot), positrons, and warm adiabatic ions are investigated. The multiple scale method is used to derive the modified Korteweg–de Vries equation (mKdVE). The Jacobi elliptic function expansion method (JEFEM) is employed to find some exact analytical solutions such as periodic, solitonic, and shock solutions. It is shown that the variation in the plasma parameters of interest, for our model, allows the existence of solitary and periodic structures and no shocks. It is also shown that the most important plasma parameters for the plasma model under consideration are positron concentration, α, and the percentage of cold and hot electrons, represented by the parameters μ and ν, respectively. Additionally, the qualitative behavior of the mKdVE is studied using dynamical system theory. The topological structure of the solution is discussed in the phase plane. In this work, the phase plane analysis, which is restricted to the discrete values of the parameter, is extended to the continuous range of the parameter using a bifurcation diagram. Bifurcation diagrams are drawn to forecast the behavior of the solution for specifically chosen essential plasma parameters. The analytical solution and the qualitative behavior of the solution presented in this paper are shown to be compatible with each other. The results presented here are general and can be gainfully employed to study a variety of nonlinear waves in space, laboratory plasmas, and astrophysical plasmas.

electrons and positrons having a small fraction of ions, whereas the converse is true for the laboratory plasmas where the positron concentration is very low.
Many approaches have been proposed to create laboratory EP plasmas over the past few decades [5][6][7]. In initial experiments, relativistic positrons from radioactive neon gas were trapped directly into a magnetic mirror (MM) [5,6]. Tsytovich and Wharton [7] suggested trapping positrons in an MM from a LINAC source. Later on, Boehmer [8] employed cyclotron heating to trap moderated positrons from a radioactive source and to heat the trapped positrons to relativistic energies. The most successful experimental approach, to obtaining positron plasmas, however, has turned out to be the scattering from a buffer gas into a Penning trap. In this case, the positron density is large enough to study the collective modes in plasmas [9]. It is important to note that pair annihilation can take place in EP plasmas. It is, therefore, a prerequisite that, for the study of collective modes in such plasma, the annihilation time scale ought to be much longer than the time scale for plasma effects, which is typically the inverse of the plasma frequency [9]. It is quite common in astrophysical and laboratory plasmas to have ions in addition to electrons and positrons, and, therefore, it is important to understand the linear and nonlinear behavior in electron-positron-ion (EPI) plasmas. Many studies have been carried out in the past few decades to investigate the linear and nonlinear wave propagation in EPI plasmas that highlight the importance of the inclusion of positrons in the system [10][11][12][13][14][15][16][17][18].
The observations of geophysical and laboratory plasmas have indicated the presence of two populations of electrons, namely, cold and hot electrons. Examples include hot cathode discharge plasmas, solar wind at about 1 AU, turbulent plasma for thermonuclear interest, and strongly interacting beam plasma systems [19]. In collisionless plasmas of space, Earth's bow shock [20], and in the up/ downstream of interplanetary shocks [21], two-temperature plasmas have been observed. Moreover, this type of plasmas has been observed in heliospheric termination shock [22], planetary magnetospheres [23], and space plasmas [24]. Additionally, different missions such as Polar [25], Geotail [26], FAST [27], and S3 − 3 [28] observed and discovered the two populations of electrons in the auroral zone and magnetosphere. The investigation of the linear ion acoustic wave (IAW) in plasma of two population electrons has a significant impact on its characteristics [29]. Many authors studied the ion acoustic solitons in plasma of two population electrons and reported many interesting facts about them [30][31][32].
No well-established methods have been developed, so far, that can solve all nonlinear differential equations arising in mathematics, physics, engineering, and biological sciences. Concerted efforts have thus been made to develop various analytical techniques to study nonlinear phenomena. Some of these methods include the homogeneous balance method [33][34][35], hyperbolic tangent expansion method [36][37][38], trial function method [39,40], nonlinear transformation method [41,42], and sine-cosine method [43]. Although, soliton and shock solutions have been obtained, periodic solutions could not be found by these methods. Porubov et al. [44,45] obtained some exact periodic solutions for some nonlinear wave equations using the Weierstrass elliptic function which involved complicated deductions. Fu et al. [46] proposed the Jacobi elliptic function expansion method (JEFEM) and applied it to some nonlinear wave problems. In addition to shock and solitary wave solutions, various periodic solutions, based on the Jacobi elliptic sine function, were also found [46]. To discuss all three possible waveforms (for traveling wave solution), we use JEFEM to present all possible analytic solutions of the mKdVE. These include both bounded and unbounded periodic solutions. In recent years, theory of the dynamical system has generated a lot of interest in the study of nonlinear equations arising in plasma physics and fluid mechanics [47][48][49][50][51][52][53][54][55][56].
Saha and Chatterjee [57] used bifurcation of the phase portrait to analyze the qualitative behavior of the dust ion acoustic traveling wave solution of the modified Kadomtsev-Petviashvili equation (mKPE). The authors showed the existence of periodic, solitary, and kink wave solutions using the theory of the dynamical system by investigating the nature of the critical points. Analytical traveling wave solution of the mKPE is also obtained, and the two approaches are shown to be consistent. This work [57] is a generalization of the work presented by Samanta et al [52] who identified two solutions instead of three, using bifurcation analysis. Tamang et al; [53] used the planar dynamical system approach to present nonlinear homoclinic and nonlinear periodic trajectories of the modified mKdVE and generalized Gardner equation (GGE). A few articles, in which Asit Saha [52,53,57] is one of the authors, bifurcation has been used to identify the wave profiles and their properties depending on the parameter values. It may be noted that the bifurcation of the wave structure, as mentioned in these articles, should not be confused with the bifurcation diagram. We know that the qualitative behavior of the phase portraits such as periodic, solitary, or kink depends on the characterizing parameters of the problem. The phase portraits are drawn for some discrete values of the parameters without any prior knowledge of their behaviors. The qualitative change occurring in the phase diagrams is observed visually depending on the closed, homoclinic, or heteroclinic orbits.
Let us first recall the procedure adopted in these studies. A model nonlinear partial differential equation (PDE) is initially transformed into a nonlinear ordinary differential equation (ODE) using a traveling wave solution. The ODE is then formulated in terms of the corresponding dynamical system, and the phase plane diagram is drawn near the critical points using linear stability analysis. The solution behavior is based on the nature of the critical points, and the wave behavior is predicted from the phase diagram. For example, the center in the phase space corresponds to the closed orbit representing the periodic wave, the homoclinic orbit representing the solitary wave, and heteroclinic orbits corresponding to the kink (or shock) wave front. The closeddomain traveling wave solution generally corresponds to these waveforms. The stability of the wave solution is determined from the critical point and the Jacobian matrix. The phase diagram and the corresponding wave nature are determined for specifically choosing the parameters.
The aforementioned analysis raises a few questions for further thinking: a) Can we choose a set of premeditated parameters which gives rise to the wave structure of our choice? b) Can we divide the whole range of a parameter into sub-intervals with different topological behaviors (if so)? Can we know the behavior of the solution for some choice of the parameters without actually experimenting it through the phase plane diagram? Can we guess Frontiers in Physics frontiersin.org the behavior for a continuous spectrum of the parameter of interest instead of its discrete values? The answer to these questions lies in stability and bifurcation theory, which, to our understanding, has not been given due importance. The issues raised previously are addressed in this work besides finding the exact analytical solution.
After taking up the analytical solution and the solution behavior of mKdVE, we show that the two approaches lead to same conclusions and corroborate each other. We are lucky to find the exact solution here; however, in more complex problems involving higher-order non-linearities, where analytical solution is not possible, the qualitative solution gives us all the physics we need to know. In the present paper, no heteroclinic orbit appears, and, therefore, no shock solution exists for all possible parameters of the mKdVE. The existence of periodic (bounded and unbounded) and solitary wave solutions is found both mathematically and qualitatively. The bifurcation points are identified at which the topological behavior of the wave form changes, that is, the switch of the dynamical behavior from periodic to solitary waves or the other way round. A complete range of the important parameters like α, β, δ, γ, μ, and ] have been specified from the bifurcation diagram predicting the behavior in advance without any fear of wrong conclusion or incomplete conclusion (leaving some values of the parameters unattended). In passing, we note that the homoclinic behavior mentioned in [57] is actually from 0.8 < q < 1.0, instead of, q > − 1. The important conclusions and observations highlighting important physical features are discussed in the concluding section. The analysis presented here will be useful to explore the physical behaviors in the problems with greater complexities in some future works where the analytical solution is not a possibility.
The paper is organized as follows: in Section 2 and Section 3, we present a set of model equations and derive the mKdVE using the reductive perturbation technique (RPT) [58-60]. In Section 4, we present the analytical solutions of the mKdVE using the JEFEM. All possible analytic solutions including bounded periodic solutions, solitary solutions, and unbounded periodic solutions are found. Section 5 explores the qualitative behavior of the mKdVE using the planar theory of the dynamical system. In Section 5, we present the results and discussion of both the quantitative and qualitative analyses using the plasma parameters that are consistent with the satellite observations of space plasmas. The main findings of the paper are recapitulated in Section 7.

Basic set of equations
We consider multicomponent plasma comprising hot and cold electrons, warm adiabatic ions, and isothermal positrons. In order to study the nonlinear behavior of the IAWs, we employ the following set of normalized fluid equations [32]: where n e = n h + n c .
The cold and hot populations of electrons and positrons are assumed to be inertialess on the ion time scale and assumed to follow the Maxwell-Boltzmann distribution function and are given as [32] and In the aforementioned equations, n i is the density of ions normalized by n i0 ; u i is the fluid velocity of the ion species normalized by ion acoustic speed C i ; n c0 , n h0 , and n i0 are the equilibrium densities of two electron populations and the ion component, respectively; ϕ is the electrostatic potential normalized by thermal potential T eff /e; α is the equilibrium density ratio of positron to electron species; ] is the equilibrium density ratio of hot electrons to the total electron number density; μ is the equilibrium density ratio of cold electrons to the total electron number density; β is the ratio of cold to hot electron temperatures; and δ is the ratio of ion to electron effective temperature. In the aforementioned equations, time and space coordinates have been normalized by the inverse of the ion-plasma frequency ω −1 pi and Debye length λ D , respectively, whereas electron densities n h and n c are normalized by n e0 . Here, (4πn eo e 2 /m i ) −1/2 , and λ D (T eff /(4πn eo e 2 )) 1/2 .

Derivation of the mKdVE
Before deriving the mKdV equation, we remark that in some cases, it so happens that the nonlinear term in the KdV equation becomes identically zero because of the value of plasma parameters. As a consequence, no solitary wave is possible. This necessitates to consider the next-order (higher-order) non-linearity to balance dispersion in the governing equation leading to the mKdV equation. Moreover, it is worthwhile to explore the effect of cubic non-linearity on the formation of nonlinear structures, which becomes possible with the mKdV equation. For that, the following stretched coordinates ξ and τ are introduced [32] where λ indicates the phase velocity and ε is a small and real perturbation parameter. We next expand the physical variables in Eqs. (4)-(7) as follows: Frontiers in Physics frontiersin.org Substituting Eqs (8)-(9) into Eqs (1)-(3) and equating the terms with same power of ϵ will yield the following equations in the lowest order, that is, O (ϵ 2 ) For the next higher-order equations, that is, O (ϵ 3 ), and after substituting the values of n i1 , u i1 , and λ, we obtain whereλ (3λ 2 − 5δ).
Substituting Eq. (23) into Eq. (21) and equating the terms of the same orders to zero, we obtain Note that the velocity u is a free parameter. Thus, from Eq. 24, we can obtain the value of k as follows: leading to the following new form of a 1 Thus, the periodic solution to the mKdVE can be expressed as However, knowing that B > 0 and u > 0, this solution can only exist for A < 0. Eq. (28), therefore, takes the following form: where m 1 = 1 − m. However, this is the unbounded periodic solution which is of little physical significance. Thus, we can consider a new ansatz Applying Eq. (29) to Eq. (21) and by applying the balancing principle, we obtain n = 1, which leads to Accordingly, the periodic Cn solution reads The solution given by Eq. (33) demands u > 0, A > 0, and 2m 2 − 1 > 0 or u > 0, A < 0, and 2m 2 − 1 < 0. The first choice gives the bounded Cn periodic solution, whereas the second choice leads to the solution where nC = 1/Cn. Equation 34 reveals unbounded periodic solutions for 0.29 < m < 0.7071 and bounded periodic solutions for 0 < m < 0.29. Remember that the wave profile depends on the phase and the modulus of the JEF. Note that hitherto invisible behavior for the second choice emerges for some range of the parameter values of m, which is also of physical significance (bounded periodic solution).

Shock solution to the mKdVE
First, we note that the elliptic Jacobi functions become hyperbolic functions at m = 1. When m = 1, the solution Eq. (27) reduces to the shock wave solution For the real solution to this equation, we have two choices i) B < 0 and A > 0 and ii) B > 0 and A < 0. The first choice is not possible for the values of the plasma parameters of this problem, and the second choice gives the unbounded periodic solution. Thus, it eliminates the possibility of the shock wave structure.

Solitary wave solution to the mKdVE
For m = 1, the solution Eq. 33 reduces to the solitary wave solution We note that based on the solution Eq. (33), the solitary wave solution is only possible for A > 0 and 2m 2 − 1 > 0. Thus, the analytical analysis using the JEFEM shows periodic and solitary solutions and no shock wave solution. The method also reveals unbounded periodic solutions which are presented for mathematical interest.

Qualitative analysis of the mKdVE (traveling wave solution)
For brevity, we take BC = A here and rewrite Eq. (20) in the following form: According to the Lie symmetry method, Eq. (37) where k = 1. Differential Eq. (38) is now expressed as an autonomous system of first-order equations (dynamical system) The Hamiltonian function (total energy) of the aforementioned system reads as follows: The dynamical system Eq. 39 reveals three equilibrium points (stationary solutions): The nature of the equilibrium can be determined from the eigenvalues of the Jacobian matrix at the location of the critical points. The eigenvalues of the Jacobian matrix of stationary solutions turn out to be The first stationary solution is a saddle, while the second and third are centers. The local behavior of the solution at the saddle point (hyperbolic point) matches with the global behavior of the Hartman-Grobman theorem. For the centers (non-hyperbolic point), the Hartman-Grobman theorem does not apply [67], and we use the energy method. Since the Hamiltonian function Eq. (40) is conserved dH dψ 0, closed orbits are represented by the center. We make a remark that for the traveling wave solution, we expect only three types of the waves, namely, periodic, solitary, or the shock  6 Results and discussion

Analytical approach
In this section, we will discuss the graphs of the numerically investigated parametric regimes for which we obtain solitonic, shock, and periodic solutions of the mKdVE admitted by the JEFEM. We mention at the outset that for the present model, the mKdVE does not admit any shock solutions; however, it allows the formation of solitonic and periodic solutions. It is also apt to mention here that the most important plasma parameters of the model are the positron concentration, α, and the percentage of cold and hot electrons present in the system represented by the parameters μ and ], respectively.
The JEF, Sn, gives unbounded periodic solutions for A < 0 for all values of 0 < m < 1 (an outcome of Sn property). When m = 1, the corresponding tanh solution (shock wave) does not exist for the parametric values of the present model. The function, Cn, clearly gives bounded periodic solutions for A > 0 and 2m 2 − 1 > 1, that is, m > 0.707. Finally, when A < 0 and 2m 2 − 1 < 0, Cn gives the unbounded periodic solution for 0.29 < m < 0.707 and bounded periodic solution for 0 < m < 0.29. The unexpected periodic solution suddenly emerges depending on the value of the modulus of Jacobi elliptic function (eccentricity of the ellipse) that strongly influences the behavior of the wave. Once again, this conclusion comes from the property of the Cn function with imaginary argument. The corresponding solitonic solution, sec h (limiting the solution of Cn when m reaches 1), exists for A > 0. From the aforementioned discussion, we conclude that for A > 0, bounded periodic Cn and solitary (sec h) wave solutions exist. For A < 0, there are unbounded periodic solutions because Sn, Cn, and tanh have imaginary arguments. An exception in the case of A < 0 comes from nC for which the bounded periodic solution exists for a range of small values of m. As m reaches closer to zero, elliptic Jacobi functions tend to yield trigonometric functions and the unbounded periodic solutions become bounded periodic functions. The mathematical reason is that in nC, it is the complimentary of m, that is, m 1 that matters instead of m. The fate of the waveforms depends on the sign of the non-linearity coefficient A and the value of m. For all the cases, the amplitude needs to be real for physical solutions. We once again emphasize that unbounded solutions are not physical. The shock solution does not appear for any choice of parameters of the problem. We now proceed to the discussion of the figures and emphasize that the reader should look at them in the light of aforementioned explanation. Here, we plot only the positive solutions.
The Sn periodic solution for μ = 0.1 and α = 0.1 is shown in Figure 1. These values yield A < 0, which corresponds to the unbounded periodic solution. If we keep α the same and choose μ > 0.17, the non-linearity   Frontiers in Physics frontiersin.org coefficient A becomes greater than zero and Sn function offers no solution. Physically, this means that the change in values of the positron concentration will change the range of A for which it goes from positive to negative, turning the unbounded periodic solution into no solution. The aforementioned analysis is true for all values of m. Before discussing Figure 2, we recall that the periodic solution, Cn, exists, in principle, for both A ≶ 0 subject to (2m 2 − 1) ≶ 0. For A > 0 (μ > 0.17), the Cn solution gives the bounded periodic solution-the physically plausible solution (Figure 3). For A < 0 (μ = 0.1) and 0.29 < m < 0.707, we obtain the unbounded periodic solution (Figure 2), whereas for 0 < m < 0.29, the solution turns into the bounded periodic solution (Figure 4). The former solution is not physical, whereas the latter is physical. The interesting physical solution emerges suddenly because of small values of m, and the Jacobi function approaches the trigonometric function which is certainly bounded. The Jacobi functions turn into hyperbolic functions for m = 1. In this case, the bounded Cn function becomes sech function representing the solitary wave solution. This is substantiated in Figure 5  All solutions that have been discussed so far theoretically are supported by the dynamical system via the phase plane analysis and the bifurcation diagrams. The topological behavior (local behavior) near the equilibrium points in the phase plane provides the wave patterns representing the solution to the mKdVE, and the bifurcation diagrams predict the behavior for the continuous spectrum of the physical parameters of the problem. The dynamical system approach is presented in the following section.
6.2 Dynamical system and bifurcation diagram 6.2

.1 Variation with respect to α
The phase portrait in the (ψ, z) plane is given by Figure 8A. In this figure, the center corresponds to the closed orbits having a finite   Frontiers in Physics frontiersin.org observation made previously. Figure 8B is the solitary wave representing the homoclinic orbit that satisfies the boundary conditions at infinity for mKdVE. Figure 8F gives the bifurcation diagram for α in the α − ψ plane, keeping other parameters fixed. The bifurcation diagram defines the solution behavior over a continuous range 0 < α < 1. For 0 < α < 0.3, there is only one saddle point at (0,0), whereas for α ≥ 0.3, there is a saddle at (0,0) and centers at ± 3u BC ( Figures 8E, F). The phase portraits give the solution behavior for the specific choice of α parameter, whereas the bifurcation diagrams predict the solution behavior for all possible choices of this parameter. The qualitative behavior for other values of α (including α = 0.3) is determined from the bifurcation diagram.
For instance, at α = 0.5, the bifurcation diagram shows the existence of two centers (red branches) and one saddle (blue branch). Similarly, for α = 0.1, there is only one singular point, saddle, which is unstable equilibrium. The value of α = 0.3 corresponds to two bifurcation points (−0.64, 0) and (0.64,0). The behavior for this value of α is also supported by the bifurcation diagram, that is, Figure 8F. A single saddle point (which does exist in the phase plane) has two branches. One branch moves toward the critical point (stable branch), and the other (unstable branch) moves away from the critical point. The one which moves away from the saddle point gives the periodic wave of infinite amplitude (the singular solution). Moreover, the inset branches of the saddle show the stable manifold, and outset branches are the unstable manifold (see Figure 8E). From the foregoing analysis, the sensitivity of the positron concentration is evident. The switching of the dynamical behavior takes place at the parametric value α = 0.3. For the range 0 < α < 0.3, there is only saddle (showing unbounded periodic solutions) and two centers (closed orbits) for α ≥ 0.3. It may further be noted that for 0 < α < 0.3 and A < 0, unbounded periodic solutions are obtained as shown in Figure 8E, whereas for α ≥ 0.3, the bounded periodic solutions are obtained as shown in Figures 8C, D. These observations are in conformity with the analytical results given previously in Section 4.

Variation with respect to β
Similar analysis is repeated for β = 0.1 in Figure 9. The analysis shows two centers and one saddle in the phase plane. According to the physics of these critical points, they correspond to close orbits and homoclinic orbit, respectively ( Figure 9A). The homoclinic orbit represents the solitary wave structure of the wave form shown in Figure 9B. The value β = 0.01 reveals only one critical point-saddle point (unstable equilibrium). These observations are authenticated from the bifurcation diagram ( Figure 9F), showing a saddle at (0,0) for 0 < β < 0.02 and centers at ± 3u BC for β ≥ 0.02. The value of β = 0.026 corresponds to two bifurcation points (−0.64, 0) and (0.64,0) ( Figure 9F). We further observe that the non-linearity coefficient A < 0 for 0 < β < 0.02 gives rise to unbounded periodic solutions ( Figure 9E), and for A > 0, the bounded periodic solution is obtained, as shown in Figures 9C, D. We conclude that the nature of these critical points, for the continuous range of their values, can be forecasted without actually experimenting for the specific values of these parameters from the phase plane analysis or analytic considerations.

Variation with respect to μ
First, we recall that variation of positron concentration (μ) in the system changes the range of percentage of cold electrons for which the periodic or solitary waves exist. The bifurcation diagram between μ and ψ is shown in Figures 10A, B) for different values of α (positron concentration). It is observed that the region for the existence of solitary wave increases with α. From the bifurcation diagram, we infer that an increase in positron concentration gives solitary behavior for lesser percentage of cold electrons μ. The branch points with respect to μ for α = 0.0 and α = 0.1 are (±0.97, 0) and (±0.94, 0), respectively ( Figures 10A, B). This observation is in agreement with the analytic results obtained previously. The wave forms for saddles in Figures 8E, 9E are not presented since they correspond to unbounded periodic solutions and are of little physical significance.

Some significant observations of bifurcation analysis
The bifurcation diagrams provide the solution behavior for a continuous range of physically important parameters, namely, α, β, and μ. Another important aspect of even greater significance is that the bifurcation theory predicts and identifies the switching of the solution behavior from one topological space to another when the value of a parameter crosses some threshold value. For example, the bifurcation diagram predicts the change in behavior from periodic to solitary at α = 0.3 and β = 0.026. Similarly, a cut off from the periodic to solitary waves with respect to μ for α = 0.0 and α = 0.1 comes at (±0.97, 0) and (±0.94, 0), respectively. Moreover, the bifurcation theory plays an important role in the analysis of highly nonlinear models whose analytical solution is not possible. It is effective to probe the stability of the solution by investigating the equilibrium curves drawn versus characterized parameters in advance without plotting the topological structures.

Summary and conclusion
The nonlinear ion acoustic waves in unmagnetized collisionless plasma consisting of two populations of electrons (cold and hot), positrons, and warm adiabatic ions are investigated. The mKdVE is derived using multi-scale analysis and solved by the Jacobi elliptic function expansion method. The bounded periodic and solitary solutions are obtained. In addition to that, unbounded periodic solutions (of little physical significance) also appear in the solution. No shock wave is found for any value of plasma parameters.
The dynamical system theory and the bifurcation analysis are applied to determine the qualitative behavior of the mKdVE. The behaviors are determined from the nature of the critical points and the corresponding phase plane diagrams. The centers and homoclinic orbits ensure periodic and solitary wave behaviors, whereas the absence of heteroclinic orbit confirms the nonexistence of the shock wave and corroborates the analytical results. A complete convergence between the analytical results and the qualitative behavior is established. In the end, we conclude with the remark that the qualitative behavior is useful when the closed-form solutions of the nonlinear differential equations are known and crucial and when analytic solutions are not possible. The results presented here provide a general framework and can be used to study a variety of plasma waves in laboratory, space, and astrophysical plasmas.
Frontiers in Physics frontiersin.org