Dust-acoustic solitary and periodic waves in magnetized self-gravito-electrostatic opposite polarity dusty plasmas

Dust-acoustic (DA) solitary and periodic waves investigations were performed in a magnetized self-gravitating dusty plasma consisting of negatively and positively charged dust grains in the presence of inertialess ions and electrons. The Korteweg–de Vries–Burger (KdVB) equation has been derived. The numerical investigations revealed the compressive or rarefactive DA solitons depending on the plasma parameters. The nonlinear homoclinic and periodic trajectories from the KdVB equation were obtained for the phase portrait profiles when employing the phase plane theory of dynamical systems. The periodic wave solution depends also on the system parameters. The present results are considered to be beneficial in understanding the nonlinear structures in experimental devices and different astrophysical environments such as the Earth’s mesosphere, cometary tails, and Jupiter’s magnetosphere.


Introduction
When positive and negative dust exists in plasmas, a new dusty plasma system termed "opposite polarity dusty plasma" (OPDP) arises [1][2][3][4][5][6][7][8][9]. Generally, the collection of electrons from the surrounding medium and the thermal speeds of electrons that far exceed that of ions cause dust grains to be negatively charged. However, there are some observations of positively charged dust grains in different parts of space (Jupiter's magnetosphere [10], upper mesosphere [11], cometary tails [12], thick molecular clouds [13], and in the Martian atmosphere [3]). The charging processes of the positive dust are achieved by different mechanisms, such as the secondary emission of electrons from the surface of the dust grains [3], the thermionic emission from particles heated by radiative sources [3], and the photoemission under ultraviolet radiation [14].
Usually, the size of positively charged particles compared to negatively charged dust might be smaller or larger than unity in the dusty plasma, and it can even be equal to unity. This important feature of dust particles makes the OPDP distinguished from electron-ion and electron-positron plasmas. Depending on the dust size, the mass and the magnitude of the charge of the dust are determined. Chow et al. [15] studied theoretically the effect of dust species size (the size of positively charged dust grains is smaller than that of the negatively charged dust grains) on the charging of insulated dust grains by the secondary electron emission. Also, Mendis and Rosenberg [16] and Zhao et al. [17] studied the opposite situation (the size of positively charged dust grains is larger than that of the negatively charged dust grains).
Dust-acoustic waves (DAWs) propagating in a multicomponent, collisional, magnetized dusty plasma have been investigated, recently by Tolba [18], and found that relative densities and temperatures parameters affect the behavior of the obtained waves. The extremely massive charged dust particles may be affected by both the electromagnetic force [19] and the gravitational field [20] in astrophysical environments when compared to charged electrons and ions, which are controlled by electromagnetic fields only. The gravitational and the electric forces for a massive dust grain in this case become nearly comparable in their values. Therefore, they play an important role in dust particles dynamics [20,21] in space plasmas containing dust grains with a radius larger than 1 µm [22]. Accordingly, for the dust species with radius r d > 1µm the gravitational force is important for the dynamical properties. On the other hand, for r d < 1µm, the electrostatic (repulsion) force must be considered [23,24]. Both electrostatic and gravitational forces govern the massive and highly charged dust grain dynamics since the dust sizes vary from 1 to 15 µm. In this situation,

The model
We consider the propagation of the DA waves in a self-gravitating OPDP consisting of negatively and positively charged dust grains and inertialess ions (following the Maxwellian distribution). We neglect the contribution of the free electrons since the latter were almost completely depleted to the surface of the negative dust species; also, we neglect the effect of the self-gravitational field on the ions since the number density of ions is neglected in comparison with that of positive or negative dust species. Since the charged dust grains are highly massive particles, the self-gravitational, electrical forces should be taken into account and in the presence of an external static magnetic field B 0ẑ which is aligned in the z-direction, where B 0 is its magnitude andẑ is the unit vector along the z-direction. The nonlinear gravito-electrodynamics [19] of the DA waves are described by the normalized equations in the form [24,[43][44][45] where n 1 (n 2 ) represents the positive (negative) dust number density normalized by its equilibrium value n 10 (n 20 ); u 1 (u 2 ) is the positive (negative) dust fluid speed normalized by the negative DA speed C d2 = √ Z 2 k B T i /m 2 (T i being the ion temperature, and k B being the Boltzmann constant); α = Z 1 /Z 2 with Z 1 being the number of electrons depleted from the positive dust and Z 2 being the number of electrons residing in the negative dust. β = m 1 /m 2 with m 1 (m 2 ) being the mass of positive (negative) dust; η 1,2 being the normalized longitudinal viscosity coefficient which is normalized by m 1,2 n 10,20 ω p2 λ 2 D2 ; φ is the electrostatic wave potential normalized by k B T i /e (with e being the charge of an electron); ψ is the self-gravitational potential normalized by C 2 d2 ; the space variables (x, y, z) in ∇ =x ∂ ∂ x +ŷ ∂ ∂ y +ẑ ∂ ∂z are normalized by Debye length λ D2 = K B T i /4π Z 2 n 20 e 2 ; and t is the time variable normalized by ω −1 p2 = m 2 /4π Z 2 2 n 20 e 2 ; ω c1,2 = Z 1,2 eB 0 /m 1,2 ω p2 is the cyclotron frequency for the positive and negative dust which are normalized by negative dust plasma frequency ω p2 ; μ = n 10 /n 20 ; and γ = ω 2 J 2 /ω 2 p2 (with ω J 2 = √ 4π Gn 20 m 2 being the Jeans frequency for the negative dust component and G being the universal gravitational constant).

Derivation and solution of the KdVB equation
To study the characteristics of the nonlinear DAWs propagating in this magnetized plasma system, we derive the nonlinear KdVB equation by using the reductive perturbation method [46]. So, we first introduce the stretched coordinates as The viscosity coefficient is rescaled as η i = 1 2 η i0 , i.e., low viscosity is needed, where is the smallness parameter measuring the perturbation amplitude. V p is the phase speed of the perturbation wave. The wave vector k has directional cosines given by L x , L y , and L z along x, y, and z axes, respectively, satisfying the relation L 2 x + L 2 y + L 2 z = 1. We expand the dependent parameters n 1,2 , u 1,2x,1,2y,1,2z , ψ and φ in power series of as [47] Now, substituting Eqs. (2) and (3) into the set of basic Eqs. (1) and collecting equations having the same power, the lowest order coefficients of lead to (1) .
Also, we obtain the higher-order coefficients of for the momentum equation components along x-and y-directions as Eliminating the second-order perturbed quantities n (2) 1,2 , u 1,2y , u 1,2z , φ (2) , and ψ (2) by solving Eqs. (7) and (8) with the aid of Eqs. (4) and (6), finally, we get the following KdVB equation: where A, B, and C represent the nonlinearity, dispersion, and dissipation coefficients, respectively, which are given by Both the coefficients of nonlinear and dispersive terms are affected also by the same parameters as the phase speed. Figure 2 shows that the nonlinear term can be positive or negative depending on the system parameters. For positive A, A decreases (increases) as θ, μ, and β (α) increase (s), while when A becomes negative, its absolute value increases (decreases) as β and μ (α and θ ) increase. Figure 3 shows that the dispersive term is affected by the same parameters in an opposite way, where B increases (decreases) as θ, μ, and β (α) increase(s).
These are physical for the DA waves to exist in this magnetized plasma system due to the reasons of both the restoring force that comes from the inertialess electrons pressure, while the inertia is provided by the dust mass. This leads to an increase in the phase velocity as shown in Fig. 1, which in situ leads to affect both the nonlinear and dispersive terms as presented in Figs. 2 and 3, respectively.
To find the solution of KdVB Eq. (9), we use the variable χ = η − U T , where χ is the transformed coordinate relative to a reference frame moving with the shock speed U . Integrating Eq. (9) with respect to χ leads to We use the tanh method to determine the analytical solution of this equation. According to this method, the solution of the KdVB equation takes the form [48] φ (1) Moreover, the associated electric field (E = −∇φ (1) ) [49]with the solution (12) is given by Figures 4 and 5 present the potential and the associated electric field against χ, respectively, for different system parameters. The potential includes solitary potential due to the sech 2 (χ) term and shock potential due to the tanh (χ) term. The shock term is dominant in the presence of viscosity effect and large obliquity angle as shown in Fig. 4a, while the solitary term is dominant as θ decreases. Also, as α increases the potential changes from rarefactive solitary to compressive shock, and finally, it becomes compressive solitary with double layer part as shown in Fig. 4b. Figure 5a, b depicts the associated electric field of the potentials corresponding to Fig. 4a, b, respectively. The electric field spreads as θ and α increase.
The results obtained from Figs. 4 and 5 are physical because the obliquity angle increase results in the decrease in the nonlinearity with the increase in dispersion which leads to decreasing the soliton amplitude. Small values of protons reside in the positive dust, while larger values of electrons reside in the negative dust; the phase velocity decreases and the rarefactive soliton amplitude increases while the situation is reversed for larger values of protons since the dispersion decreases and the compressive soliton amplitude increases. The associated electric field is affected by the angle, electrons, and protons in the same manner.
These results coincide with those obtained in Refs. [24,30] where α and β have the same effect on the shock wave phase speed, and in the same manner larger β changes the nonlinear term to be negative as is obtained in Fig. 2. Similarly, α affects the nonlinear term in the same way and consequently affects the perturbed potential. Also, the results obtained by Sumi et al. [32] represent the same behavior of the same parameter on phase velocity and the perturbed potential in the case of studying DA shock waves in a dusty self-gravitating opposite polarity plasma in the presence of trapped ions.
Also, we obtain the next higher-velocity components along x-and y-directions on employing the next higher-order coefficients of as ∂η + L y α β ∂φ (3) ∂η + L y (3) ∂η + L y Now, combining Eqs. (21)-(22) after simplifying by using Eqs. (19)- (20), we obtain the well-known mKdV equation as follows: where Fig. 6 The evolution of a the potential φ (1) represented by Eq. Now, taking the same frame χ = (η − U T ) as that employed in the KdV solution, the stationary solution of Eq. (24) can be directly given as [51] φ (1) where the amplitude φ (1) m = 6U F and the width δ = B U and the electric field associated with this solution is given by (1) , (26) Figure 6 represents the changes of the potential (φ (1) ) and the associated electric field (E (1) ) against χ for various values of α (= Z 1 /Z 2 ). It is depicted that, the solution, Eq. (25) gives only compressive waves, and the amplitude increases as α increases. This is physical as demonstrated before due to the effect of the electrons on the restoring force and in situ affects the phase velocity which in turn results in the change of the mKdV nonlinear term. (1) , (27) where V φ (1) is the Sagdeev potential, which is given by

From Eq. (14), the total energy equation becomes
The variation of V φ (1) with different values of the positive-to-negative dust charges ratio for compressive and rarefactive solitons, depending on β values, is plotted in Fig. 7, where the other variables (μ, θ, ω c1 , ω c2 , η 10 , η 20 ) are kept fixed. It is shown by increasing α , the width of the potential well V φ (1) is enhanced for the rarefactive soliton and shrinks for the compressive soliton. This demonstrates that the amplitude of the DA compressive solitary wave decreases by increasing α as shown in Fig. 7a. It is also seen that the depth of the potential well decreases with increasing values of α. In fact, steeper solitary waves are completely associated with a deeper potential well. The opposite effects appear in the case of rarefactive solitons associated with Fig. 7b. This effect takes place since the nonlinear term behaves in opposite way as it becomes negative (in case of rarefactive solitons) against the increase in α.
The phase portraits associated with the system of Eq. (11) can be obtained by employing the following transformations: There are two critical points E 1 (0, 0) and E 2 ( 2U A , 0) for this phase plane. For a chosen set of physical parameters, viz. α = 0.2, μ = 0.6, θ = 0, U = 0.1, γ = 0.3, and ω c1 = ω c2 = 0.5, the coefficient matrix of the linearized system (Jacobian matrix) at the critical point, E 1 (0, 0), satisfies J E 1 (0,0) = −U/B < 0 which gives a saddle node. The coefficient matrix at at the critical point, , so it represents a center point. There are two qualitatively different phase trajectories, namely nonlinear homoclinic trajectory (NHT) and nonlinear periodic (NPT) trajectory when we take η 10 = η 20 = 0. These phase trajectories NHT and NPT represent solitary and periodic wave solutions of Eq. (27), respectively. The parameters α, β, μ, and other system parameters affect the formed solitary waves that can be rarefactive or compressive as shown in Fig. 8a, b, respectively. This is physical due to the effect of β on the phase velocity as shown in Fig. 1 and consequently affects both the magnitude and sign of the nonlinear term. The positive region corresponds to a compressive solitary solution, and the negative region corresponds to a rarefactive solitary solution of the KdV Eq. (14). Therefore, we obtain both compressive and rarefactive DAW solutions of the KdV Eq. (14) for the plasma system.
One can note from Fig. 8 that the phase curve (NPT) is repeated on the same path and one complete cycle corresponds to a wavelength in the physical space. Physically, the pseudoparticle's speed becomes zero (dφ (1) /dχ = 0); it is reflected as a result of the potential force since(−dV φ (1) /dφ (1) = 0) and it oscillates between the two points. Therefore, the NPT in the phase plane is a periodic orbit. On the other side, the NHT curve in Fig. 8a emerges from the saddle point, circling anticlockwise around φ (1) axis, and again stops at the saddle point, entering from the upper side. As the analogy, the pseudoparticle starts with zero speed and attains some speed along the φ (1) axis. After maximum speed, its speed decreases and reaches zero at φ (1) = φ (1) max . Furthermore, the pseudoparticle bounces back toward the saddle point as a result of the potential force, where −dV φ (1) /dφ (1) = 0, and gains speed in the opposite direction and returns to rest at the saddle point. In physical space, the electric potential increases from zero (atχ = −∞), attains a maximum value at χ = 0, and then decreases until it becomes zero (at χ = ∞). This potential structure represents the DA compressive soliton and does not repeat [23], while Fig. 8b corresponds to rarefactive DA solitons in which the nonlinear term becomes negative.

Nonlinear periodic wave solutions
We follow the same procedure in Refs. [50,53] and [54] to find the cnoidal wave solution of the obtained KdV Eq. (14) for DA waves. The solution is assumed to be φ (1) (χ), in which we introduce the same variable transformation χ = η − U T . Accordingly, Eq. (14) can be written as where = φ (1) (χ). We get some algebraic calculations for evaluating the integration of Eq. (30), which takes the energy balance equation as where the Sagdeev potential, W ( ), in this case can be defined as When the potential vanishes, both ρ 0 and E 0 are the density of charge and the electric field, respectively. The potential function W ( ) has two points of extreme for , which are defined at ∂ W/∂ = 0; these are 1,2 and can be given as Moreover, for real values, U 2 /A 2 > 2ρ 0 B/A must hold. The zero's of the potential energy (32), i.e., = z 1 , z 2 , z 3 , are given by Now, U 2 /4 − 2 Aρ 0 B/3 > 0 must hold to get the shape of the potential well. The positive cnoidal waves are associated with the case of > 0, while the case of negative cnoidal waves is associated with < 0. Equation (31) is associated with the first integral of the energy that can be inspired by Using the initial conditions (0) = φ 0 and d (0) /dχ = 0, we can find Substituting Eqs. (32) and (36) into Eq. (35) and after the factorization process, we have where Further, the following inequalities 2 ≤ φ 0 ≤ 1 or 1 ≤ φ 0 ≤ 2 in the last relations should also be kept. Moreover, we can get from Eqs. (35) and (37) the following relation: The periodic wave solution of Eq. (35) is given by [50,55] (χ) = φ 1 + ψ cn cn 2 (Dχ, m) , where cn is the Jacobian elliptic function and the amplitude ψ cn = (φ 0 − φ 1 ) and the modulus m (0 < m < 1) and D are defined as and Physically, the parameter m is the nonlinearity indicator, with the linear limit being m → 0 and the extreme nonlinear limit being m → 1. The cnoidal solution of Eq. (39) requires the conditions φ 1 ≤ ≤ φ 0 and φ 0 > φ 1 ≥ φ 2 . Furthermore, the cnoidal wave frequency ν and the cnoidal wave wavelength λ are defined as where K (m) is is the first kind of complete elliptic integrals. Figure 9 shows the influence of parameters α, β, μ, and θ on basic characteristics of periodic wave solution corresponding to NPT enduring critical point C 0 . The amplitude of DA periodic waves grows with the increase in α and μ but depletes as θ enhances. The periodic wave widths increase also as α, μ, and θ increase. This means physically that for lower values of depleted electrons from the positive dust (decrease in α) or lower values of positive dust density (decrease in μ); the restoring force of nonlinear periodic traveling wave decreases which is provided by the electrons pressure. The variation of the profile of the periodic traveling waves with θ is examined. It is obvious that the amplitude decreases and the width increases as θ increases. Physically, one can predict that as the nonlinear periodic traveling wave approaches the direction perpendicular to the magnetic field, the nonlinear periodic traveling wave disappears where the increase in θ causes the increase in the dispersion term coefficient in the nonlinear KdV equation and the periodic traveling wave becomes fatter. Also, the obtained periodic waves are affected by the obliqueness and species densities similarly as presented by Selim et al. [56]. It is noticed that the used plasma system is entirely contained by the magnetic force, an arrangement referred to as magnetic confinement (i.e., the magnetic field makes the nonlinear periodic traveling wave structures in the field direction spikier as demonstrated in Ref. [57]. Now, the soliton solution case [58] can be achieved at φ 1 = φ 2 = 0, in which m → 1 and ρ 0 = E 0 = 0. In this case, ψ cn = φ 0 = 3u 1 /A = ψ m , D = √ Aφ 0 /12B = √ u 1 /4B = 1/ , Jacobi elliptic function, cn, degenerates to the sech function [59,60], and then the cnoidal solution, Eq. (24), is reduced to the DA solitary wave solution, Eq. (15).

Discussion and conclusions
The reductive perturbation method has been employed to scrutinize small-amplitude DA waves in self-gravitating strongly coupled opposite polarity dust plasma system consisting of strongly coupled inertial positive and negative dust species and inertialess weakly coupled Maxwellian distributed ions, while the electrons are completely depleted to the surface of dust species to make them negatively charged. The KdVB equation has been derived, and its analytical solution is determined using a tanh method. The associated electric field has been also derived and illustrated. In this paper, phase plane analysis of the planar dynamical systems is applied to deal with the KdVB equation. The KdVB equation is reduced to a dynamical system with perturbation on a twodimensional invariant manifold so that some relations and dynamical system bifurcation methods can be applied to investigate its phase space geometry in detail, which can be used to determine the wave speed for various traveling waves such as periodic and solitary waves in the system. These results will promote the study of numerical solutions and help to seek new analytic solutions of the KdV equation. We point out that the KdV equation may have other oscillatory traveling waves. When a closed trajectory (or a homoclinic trajectory) comes from the Poincaré bifurcation (or a homoclinic bifurcation), the region enclosed by the homoclinic trajectory is compact. It implies that all trajectories and the corresponding traveling waves are bounded. Since these solutions depend on the system parameters, spiral orbits are found to approach the equilibrium within the compact region and its boundary (the periodic orbit), which exactly correspond to oscillatory traveling waves.
Our study may apply to understanding the finite-amplitude localized electrostatic disturbances and nonlinear periodic wave features in different astrophysical environments such as the Earth's mesosphere, cometary tails, and Jupiter's magnetosphere [11,61,62] where the self-gravitational field effect is presented.
Funding Open access funding provided by The Science, Technology & Innovation Funding Authority (STDF) in cooperation with The Egyptian Knowledge Bank (EKB).

Data availability
The data used to support the findings of this study are included within the article and available in Ref. [24].
Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons.org/licenses/by/4.0/.