Parametric excitation and chaos through dust-charge fluctuation in a dusty plasma
aa r X i v : . [ n li n . C D ] A ug Parametric excitation and chaos through dust-charge fluctuation in a dusty plasma
Madhurjya P Bora and Dipak Sarmah
Physics Department, Gauhati University, Guwahati, India. ∗ We consider a van der Pol-Mathieu (vdPM) equation with parametric forcing, which arises in asimplified model of duty plasma with dust-charge fluctuation [1]. We make a detailed numericalinvestigation and show that the system can be driven to chaos either through a period doublingcascade or though a subcritical pitchfork bifurcation over an wide range of parameter space. We alsodiscuss the frequency entrainment or frequency-locked phase of the dust-charge fluctuation dynamicsand show that the system exhibits 2:1 parametric resonance away from the chaotic regime.
I. INTRODUCTION
The subject of parametric excitation can be traced back to Faraday in 1831 [2], when he observed that surfacewaves in a fluid-filled cylinder under vertical excitation exhibited twice the period of the excitation itself. The mostsimplified version of parametric excitation was given by Mathieu in 1868 [3] related to the vibrations of an ellipticalmembrane, which now has become a model equation for response of many systems to sinusoidal parametric excitation.The simplest Mathieu equation can be stated as [4],¨ x + ( δ + ǫ cos t ) x = 0 , (1)with δ and ǫ as constants. There have been extensive investigations of parametric excitation and resonance relatedMathieu equation by several authors [5, 6, 7]. One of the very common examples of parametric forcing modeled bythe nonlinear Mathieu equation is the forced and unforced inverted pendulum [8].In this work, we have studied the parametric excitation and resonance of the van der Pol-Mathieu (vdPM) equation,which arises in a simplified model of dusty plasma with dust-charge fluctuation. Dusty plasmas are characterized bypresence of massive dust (impurities) particles embedded in an electron-ion plasma [9]. Immersed in the plasma,the dust particles acquire charges by collecting the electrons and ions on their surfaces, which are mostly negative.However, the charge on a dust particle is never a constant and varies temporally. Thus, along with other theirdynamical properties, the charge of the dust particles becomes a dynamic variable, which can severely modify theplasma properties. The presence of dust particles in a plasma may modify the dynamics of the plasma in many ways,of which the most prominent is the appearance of low frequency dust-acoustic waves [9]. The dust-charge fluctuationis known to damp the acoustic waves in a dusty plasma [9, 10]. Recently, Momeni, Kourakis, and Shukla [11] hasstudied a simplified model of nonlinear dust-charge fluctuation based on a vdPM equation, where they have discussedthe stability regions of the vdPM equation and have shown the existence of stable and unstable periodic orbits indifferent parameter space. This vdPM equation is originally derived by Saitou and Honzawa [1], where they haveshown that in a very restricted region of parameter space, the system exhibit chaotic behavior. They argue that theoscillations leading to chaos basically stems out of the balancing between the van der Pol (vdP) and Mathieu-liketerms.We report, in this paper, a detailed investigation of the vdPM equation for dust-charge fluctuation and show thatthe vdPM equation proposed by Saitou [1] can exhibit chaotic behavior over an wide range of relevant parameterspace, which is, in many instances, preceded by period doubling cascades. We have found that both period doublingand pitchfork bifurcations take place as one varies the bifurcation parameters, both leading to chaos. The range ofparameters for a chaotic regime comes out to be not as restrictive as pointed out by Saitou [1]. We have discussedthe stability and bifurcation of the nonlinear system and found that the system can be completely deterministic inbetween chaotic regimes. The paper is organized as follows. In Section II, we formulate the nonlinear dust-chargefluctuation model yielding the vdPM equation and discuss about the parametric forcing. In Section III, we havediscussed about frequency entrainment (frequency-locked phase) of the nonlinear oscillator, driven by the parametricforcing term. We have shown that far away from the chaotic regime, the system can be driven by a 2:1 parametricresonance leading to a stable limit cycle. Away from the resonance, it displays quasi-periodic behavior. We discussthe stability and bifurcation of the periodic orbit of the nonlinear system in Section IV with the help of Floquetstability theory and show that the instability of a limit cycle may manifest through a period doubling bifurcation. ∗ Electronic address: [email protected]
In Section V, we explore the chaotic regime of the vdPM equation, where we show that the system can be driven tochaos over an wide range of parameters. We have found that the route to chaos may be through a period doublingcascade. We draw the conclusions in Section VI.
II. DUST-CHARGE FLUCTUATION MODEL
We consider an unmagnetized collisionless plasma consisting of electrons, ions, and massive dust particles, whichbecome charged by acquiring charged particles (ion or electrons) on their surfaces. The subject of dust-chargefluctuation in a dusty plasma is a well studied process which, primarily, has a damping effect on the acoustic waves[10]. Here, we assume that the number densities of the charged particles are considerably larger than that of thedust particles, so that the effect of dust-charge fluctuation on the dynamics of the electrons and ions is negligible andcharge neutrality is always satisfied [1]. So, the charge on the dust grains q ( t ) becomes a time dependent function.Assuming the equilibrium (unperturbed) state to be static, we can write the nonlinear continuity, momentumbalance, and the Poisson equation for dust-charge fluctuation as [1, 11], ∂n∂t + n ∇ · v = αn − βn , (2) m d d v dt = q E , (3) ε ∇ · E = qn, (4)where m d , n , ε are dust mass, equilibrium dust number density, and permittivity of free space. The variables n, v , E are perturbed dust density, velocity, and electric field. In the first equation, Eq.(2), the terms on the right hand sidedenote the rate of production and loss of charged dust grains, where α and β are constants of proportionality. Inwriting these terms, we have assumed that the production rate of charged dust particles is proportional to the dustdensity. The cubic loss term appears mainly due to the loss of dust grains through a three-body recombination process[1]. In the momentum equation, Eq.(3), we have assumed the dust particles to be cold which basically eliminatesany variant of dust-acoustic waves. In writing these equations, we have assumed that the average dust velocity, v , isfairly uniform in space and its spatial gradient is considerably smaller so that the convective derivative term in themomentum equation, namely, the term ( v · ∇ ) v can be neglected. Further, we have approximated the term ∇ · ( n v )with n ∇ · v , assuming a uniform distribution of the charged dust particles in space ( ∇ n ≈
0) [11]. We assume thedust-charge q ( t ) to be changing harmonically with time with a frequency ν and use the ansatz [1, 11], q ( t ) = q (1 − ǫλ cos νt ) / , (5)where the term ( ǫλ ) denotes the strength of charge fluctuation. Note that, in principle, there is no need to restrictthe value of the term ( ǫλ ) to a smaller value, which gives us freedom to explore an wider parameter space of α - ǫ .Without loss of generality, we consider only one dimension, z and write Eqs.(2-4) as, ∂n∂t + n ∂v z ∂z = αn − βn , (6) m d ∂v z ∂t = qE z , (7) ε ∂E z ∂z = qn. (8)By taking a z -derivative of Eq.(7), we can eliminate the terms involving perturbed velocity and electric field usingEqs.(8) and (6) to get a coupled differential equation in perturbed density, d ndt − ( α − βn ) dndt + nω d (1 − ǫλ cos νt ) = 0 , (9)where ω d = ( n q /m d ε ) / is the plasma frequency corresponding to the dust particles. The above equation, Eq.(9)can be classified as van der Pol-Mathieu (vdPM) equation [1], owing to the nonlinear term ( α − βn ) which is like avan der Pol (vdP) term [4] and parametric forcing term (1 − ǫλ cos νt ) which like the parametric term of a classicalMathieu equation [4]. A. Parametric forcing
As is well known from the theory of classical Mathieu equation, the parametric forcing term in Eq.(9) makes thedynamics of the dust-charge fluctuation prone to chaos [8, 12, 13, 14]. As the vdP equation has a stable limit cycle,we can see that Eq.(9) should show vdP-type behavior for large α . However the parametric forcing term may stilldrive the system unstable. In all probability, we expect the onset of chaotic behavior as ǫ increases, which should bemore pronounced when α ≪ ν is close to twice the frequency of the unforced oscillator. So, when both the vdP and Mathieuterms are present, as in Eq.(9), we expect a frequency entrainment at 2:1 [15] and the system represented by Eq.(9)must exhibit some sort of quasi-periodic and frequency-locked (entrainment) behaviors in the parameter space of ǫ - α before it can be driven to chaos. The route to chaos should be through a series of quasi-periodic regime or through aperiod-doubling cascade rather than the other universal route i.e. through intermittancy [16]. III. FREQUENCY ENTRAINMENT
Entrainment dynamics plays an important part in design engineering and many other dynamical systems [17].Recently, frequency entrainment is shown to exist in nonautonoums chaotic oscillators [18]. In this work, we considerthe possible entrainment by the parametric term, which can lead to quasi-periodicity and finally to chaos. We considerthe following dynamical equation for this dust-charge dynamics,¨ x − ( α − βx ) ˙ x + ω d x (1 − ǫλ cos νt ) = 0 , (10)where we have replaced the dust density n by the variable x . In order to facilitate the multiple scales in the problemand entrainment, we assume that α ∼ β = δ ≪
1, a small number. As the entrainment is possible only when thesystem as far away from the chaotic regime, when ǫ is small, we assume that ǫ <
1. We further re-scale the time by t → ω d t and write Eq.(10) as ¨ x − µǫ (1 − x ) ˙ x + x (1 − ǫλ cos 2 ωt ) = 0 , (11)where µǫ = δ/ω d . The strength of the parametric forcing term is given by ǫλ with ǫ < ν = 2 ω . We expect thatthe parametric forcing should result in a 2:1 subharmonic resonance, when the parametric frequency ν is close to 2 or ω ∼
1. Note that in absence of the parametric forcing term ( ǫλ = 0), the natural frequency of the oscillator is unity.We now introduce two different time-scales [15, 16], the stretched time ξ = ωt and the slow time η = ǫt and expandthe forcing frequency ω about the natural frequency of the oscillator i.e. 1 with ǫ as the expansion parameter, ω = 1 + kǫ + O ( ǫ ) , (12)where k is a detuning parameter at order ǫ . The variable x now is expanded in a power series x = x ( ξ, η ) + ǫx ( ξ, η ) + O ( ǫ ) . (13)Substituting Eqs.(12) and (13) in Eq.(11) and collecting terms at the order ǫ = 0 and 1, we have, x ξξ + x = 0 , (14) , (15)where the subscripts refer to derivatives with respect to ξ and η . The solution to Eq.(14) can be taken as x ( ξ, η ) = A ( η ) cos ξ + B ( η ) sin ξ, (16)where the coefficients A and B are functions of the slow time-scale. Substituting Eq.(16) into Eq.(15) and removingthe secular terms [16], we have the following coupled differential equations for the slow time-scale, A ′ = − kB + 12 µA − µA ( A + B ) + 14 λB, (17) B ′ = kA + 12 µB − µB ( A + B ) + 14 λA. (18) (cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1) −1 −0.5 PSfrag replacements λ = 4 | k | λ = 4 √ k k λ
00 0 . − . Figure 1: Bifurcation diagram of Eqs.(17) and (18) in the k - λ plane. The line λ = 4 | k | along which a Hopf bifurcation occursat the origin. Below this line, in the shaded region, there exists a limit cycle. The dashed line denote a saddle-node bifurcationat the origin. PSfrag replacements
QPQP EN R p forward sweepbackward sweep0 0 . Figure 2: Entrainment by the parametric dust-charge fluctuation for α = 0 . , β = 0 . , λ = 1 . , ω d = 1 . , ǫ = 0 .
5, where theamplitude R of the oscillation is plotted against the period p of the parametric forcing term. The two vertical dashed linesat p = 2 . . We note that the hyperbolic fixed points of the slow flow correspond to the periodic motion of the original equationEq.(11) i.e. an entrainment and limit cycles of the slow flow correspond to quasi-periodic motions of Eq.(11) [19].From simple resonance dynamics it is evident that the entrainment region of Eq.(11) by the parametric forcing termshould increase as the parametric forcing amplitude λ increases, allowing an wider range in the detuning parameter k during which the entrainment is observed. Therefore it is worthwhile to study the flow of the slow variables throughEqs.(17) and (18) as we vary the parameters k and λ .In Fig.1, the bifurcation diagram of the slow flow in the k - λ plane is shown. Along the line λ = 4 | k | , a Hopfbifurcation occurs at the origin, below which, in the shaded region, a limit cycle appears. As λ falls below 4 | k | , thefixed point at the origin becomes an unstable spiral [complex eigenvalues of the linearized Jacobian of Eqs.(17) and(18)] from an unstable node. Along the line λ = 4 √ k , saddle node bifurcation occurs for sufficiently high λ whenthe origin becomes a saddle point from an unstable node after two other saddle points coalesce at the origin. Fromwhat we have observed, one can conclude that entrainment by parametric forcing occurs above the shaded region ofFig.1. In the shaded region, we expect a quasi-periodic behavior. These results are confirmed from the numericalsolutions of Eq.(11), which are shown in Fig.2. The entrainment region for λ = 1 can be calculated from Fig.1 as − . < k < +0 .
25, which in terms of period of Eq.(11) for the parameters of Fig.1, is given by 2 . < p < .
59 whichvery well agrees with the numerical results (see Fig.2).As we increase ǫ , the expansion parameter, the prediction from Fig.1 however agrees less and less as the systemmakes a transition to chaos. IV. STABILITY AND BIFURCATIONS
The stability of a periodic orbit can be effectively studied using Floquet theory [19, 20]. In this section, we brieflyreview the Floquet theory with reference to Eq.(10). We first write the second order dynamical equation Eq.(10) astwo first order equations, ˙ x = y, (19)˙ y = ( α − βx ) y − ω d x (1 − ǫλ cos νt ) . (20)We have already shown in the previous section that a frequency-locked (entrained) phase with a single periodic limitcycle can exist for Eqs.(19,20). In this section, we are going to study the stability of this periodic orbit and bifurcationsleading to the onset of chaos.The Poincar´e map of an initial point z = ( x , y ) on the periodic orbit (limit cycle) can be obtained by samplingthe orbit points z n at discrete time interval t = t n , n = 1 , , , . . . . So the transformation which successively mapsthe Poincar´e section P ( z ) is z n +1 = P ( z n ). The linear stability of a q -periodic orbit with P q ( z ) = z can now bedetermined from the linearized map given by the matrix DP q of P q at an orbit point z , where P q is the q -timesiterated Poincar´e map. The linearized matrix M = DP q can be obtained by integrating the linearized equationscorresponding to Eqs.(19,20) for small perturbations along the q -periodic orbit [19].Assume that z ⋆ ( t ) = z ⋆ ( t + q ) is a point lying on the q -periodic limit cycle of Eqs.(19,20). We perturb the orbitwith a small perturbation δz = ( δx, δy ) and linearize Eqs.(19,20) about the closed orbit, (cid:18) ˙ δx ˙ δy (cid:19) = J ( t ) (cid:18) δxδy (cid:19) , J ( t ) = (cid:18) f x f y (cid:19) ( x,y )=( x ⋆ ,y ⋆ ) , (21)where J ( t ) is the q -periodic linearized Jacobian. The partial derivatives f x,y are given by, f x = − βxy − ω d (1 − ǫλ cos νt ) , f y = ( α − βx ) . (22)We now assume that W ( t ) = [ w ( t ) , w ( t )] is a fundamental solution matrix with W (0) = I [19]. The general solutionof the q -periodic system, Eq.(21), is then given by, (cid:18) δx ( t ) δy ( t ) (cid:19) = W ( t ) (cid:18) δx (0) δy (0) (cid:19) . (23)We then substitute Eq.(23) into Eq.(21) to obtain the initial value problem,˙ W ( t ) = J ( t ) W ( t ) , W (0) = I, (24)where W ( q ) is the linearized map DP q ( z ). So, the matrix DP q can, in principle, be obtained from numericalintegration of Eq.(24) over period q . However, the numerical procedure is not very straight forward and requiressophisticated techniques to determine the exact form of the matrix DP q which is very sensitive to initial conditions( x ⋆ , y ⋆ ). The matrix DP q is known as the monodromy matrix for the q -periodic orbit and the eigenvalues of thismonodromy matrix, popularly known as the Floquet multipliers [19, 20], indicate the stability of the q -periodic orbit.Therefore, the values of the Floquet multipliers have to be determined with considerable precision for understandingthe true nature of the stability of the nonlinear system. The characteristic equation of the linearized map M = DP q is given by, ζ − τ ζ + ∆ = 0 , (25)where the eigenvalues ζ , are the Floquet multipliers and τ = tr( M ) , ∆ = det( M ). The determinant ∆ is given by[19, 20, 21] ∆ = e R q tr( J ) dt = e ( α − βx ⋆ ) q , (26) τ = 1 q Z q tr( J ) dt (cid:18) mod 2 πiq (cid:19) . (27)We know from Floquet theory that the periodic orbit is stable only if the pair of Floquet multipliers lie inside theunit circle. The bifurcations of the periodic orbit occur on the unit circle. From the expression for ∆, Eq.(26), it canbe seen that we do have bifurcations depending on a balancing of the parameters α, β , and the periodic orbit, which is αǫ Figure 3: Stability diagram of Eqs.(19,20) in the α - ǫ plane. The period doubling bifurcations occur along the lines denoted bythe points ‘ × ’. The points denoted by a ‘ ◦ ’ along the dashed line indicate pitchfork bifurcations. The chaotic regime lies abovethe two lines. As we can see that some pitchfork bifurcations are followed by a period doubling cascade. The other parametersare β = 0 . ν = ω d = λ = 1 .
0. Table I: Scaling of period doubling cascades. k α k δ k ǫ k δ k determined by the parameter ǫ , the magnitude of the parametric driving force. In all probability, the unstable regionshould lie in the region of large α . We numerically determine the Floquet multipliers for a range of periodic orbit inthe parameter space ( α, ǫ ) with β = 0 . ω d = λ = 1 . × ’ sign and subcritical pitchfork bifurcationsoccur along the dashed line denoted by a ‘ ◦ ’. The line joining the ‘ × ’ points in the figure denotes the accumulationpoints or the limiting points of the period doubling bifurcations before transition to chaos. The chaotic regime liesabove these lines. As the two lines seem to intersect, when extended, we see that some of the period doubling cascadesare preceded by pitchfork bifurcations.As the pair of Floquet multipliers decreases through − q -periodic orbit loses itsstability to jump to a 2 q -periodic stable limit cycle. So, in the observed parameter regime of Fig.3, all period doublingbifurcations are supercritical [16]. In case of the pitchfork bifurcations, the Floquet multipliers increases through+1 and are subcritical as there are no stable limit cycles after the bifurcations and the system becomes aperiodic[16]. In Figs.3(a) and (b), two successive period doubling orbits are shown. In Figs.5(a) and (b), we have shown thebifurcation diagrams with the bifurcation parameters as α and ǫ (for details, please see the captions in the figure).These bifurcation diagrams are obtained with the help of AUTO [22] as a part of the XPPAUT package [23]. A. Scaling of the period doubling cascades
It is interesting to investigate the scaling behaviour of the period doubling cascades in light of the scaling ofthe period doubling sequences in 1-D maps. As usually observed in any period doubling cascades, the bifurcationparameter, which in our case are α and ǫ , converges geometrically to a limiting value with the convergence ratioapproaching a unique value, analogous to Feigenbaum number in case of 1-D maps [16, 24]. In Table.I, we havelisted the values of the bifurcation parameters as these converge to their respective limiting values. The ratio of thisconvergence δ k , expressed as, δ k = A k − A k − A k +1 − A k , lim k →∞ δ k = δ, (28) (a) xy (b) xy Figure 4: Phase portrait of two successive period doubling orbits for α = 0 . ǫ = 1 . ǫ = 1 . PSfrag replacements αy (a)00 . . . . . . PSfrag replacements ǫy (b)1 . . . . . . . . . . . Figure 5: Bifurcation diagrams for the period doubling sequences. In (a), the bifurcation parameter is α with ǫ = 1 . ǫ with α = 0 .
15. All other parameters are same as in Fig.3. The open boxes in the diagramsare unstable orbits and the solid lines indicate stable periodic orbits (limit cycles) with period doubling bifurcations. Blow-upregions of the second period double cascade is shown in the insets. approaches a unique limiting value δ as the bifurcation parameter A k → const . , with k denoting each successiveperiod doubling point. The value of δ characterizes the scaling property of the bifurcation, which agrees well with theFeigenbaum number 4 . . . . for 1-D maps [8, 14, 24]. V. TRANSITION TO CHAOS
In this section, we study the chaotic behavior of Eq.(10). The single most prominent feature of chaos is its sensitivedependence on initial conditions, which is measured by the Lyapunov characteristic exponents (LCEs) [25]. Theseexponents are invariant global indicators of the non-linear system. For a continuous dynamical system described bya set of autonomous ordinary differential equations, the number of LCEs is equal to the dimension of the system. Bydiscretizing the temporal dimension as ∆ t , the LCEs can be defined as, σ i = lim N →∞ lim ∆ t → N ∆ t ln[ S i ( M N )] , (29)where M N = N Y i =0 e J ( i ∆ t )∆ t . (30)In the above relations, S i are the singular values of the matrix M N and J is the Jacobian. Numerically, the number N denotes the number of integration steps of length ∆ t . Here, we employ a numerical algorithm, based on the Wolf’s (a) αǫ (b) ǫx max Figure 6: (a) Plot of the maximal Lyapunov exponent, σ , for Eq.(10) in the parameter space of α - ǫ . The Lyapunov exponentis calculated in the entire parameter space indicated in the figure. The blackened portions indicate a positive and the blankportions of the figure indicate a negative Lyapunov exponent. As can be seen from the figure, the fractal nature of the chaosis indicated by repetitive appearance of the whole figure in the smaller regions close to the horizontal axis. It also seems thatthe entire figure is only a part of a large figure covering the entire α - ǫ plane. (b) The orbit diagram for Eq.(31), where thesuccessive local maxima are plotted against the bifurcation parameter, for α = 0 . . Other parameters are same as in (a). well known method to calculate the LCEs [26]. Note that, Wolf’s original method is not the best approach to calculatethe LCEs and it can be modified to contain the non-uniformity-factors of the LCEs [27]. An important necessary stepin Wolf’s algorithm is to re-orthogonalize of the set of vectors, which is carried out, usually, through Gram-Schmidorthogonalization procedure. In our numerical routine too, we have used the Gram-Schmid orthogonalization tocalculate the LCEs. In particular, we have calculated the LCEs for the following autonomous system,˙ x = y, ˙ y = ( α − βx ) y − ω d x (1 − ǫλ u ) , ˙ u = u (1 − u − v ) − πv/T, ˙ v = v (1 − u − v ) + 2 πu/T, (31)where ν = 2 π/T . Note that, we have made the original non-autonomous equation, Eq.(10), a system of autonomousequations by introducing the last two equations of Eqs.(31), which have stable and unique solutions, u ( t ) = cos(2 πt/T ) , v ( t ) = sin(2 πt/T ) , (32)for the initial values [ u (0) , v (0)] = [1 , α - ǫ for the range shown.The parameters are β = 0 . , ω d = 1, and λ = 1. The period of the orbit is chosen as T = 6 . ν = 0 . blackened portions of the figure indicate a positive Lyapunov exponent indicating chaos. Inall other places the maximal Lyapunov exponent is negative, signifying stable oscillations of the system. We can seethe fractal behavior of the chaos from the figure. It also seems that the whole figure is only a part of a larger figurecovering the entire domain of the α - ǫ plane. In Fig.7, the maximal LCE, σ , is plotted along with a blow up of theregion of period doubling cascades. The period doubling points are marked with ‘arrows’ in Fig.7(b).We have constructed an orbit diagram [16] [Fig.6(b)] for the system, Eq.(31), where we have plotted the successivelocal maxima for the variable x of the oscillation against the bifurcation parameter ǫ . The period doubling bifurcationsoccur near 1 . . ǫ = 2 . . VI. CONCLUSION
In this work, we have carried out a detailed investigation of the stability, bifurcation leading to chaos for the vdPMsystem arising out dust-charge fluctuation in a a dusty plasma. We have shown that the system can be highly chaoticdepending the chosen parameters and not as restrictive as has been pointed out by Saitou and Honzawa [1]. In fact,an wide range of chaotic region exists for the parametric driving strength ( ǫλ ) as low as 0.05, as shown by the plot of −0.7−0.6−0.5−0.4−0.3−0.2−0.1 0 0.1 0.2 0 1 2 3 4 5 6 7 8 PSfrag replacements ǫσ (a) −1 PSfrag replacements ǫσ (b)01 . .
95 2 2 .
05 2 . .
15 2 . .
25 2 . − − . − . − . − . . −0.3−0.25−0.2−0.15−0.1−0.05 0 0.05 0.1 1.16 1.18 1.2 1.22 1.24 1.26 1.28 1.3(c) PSfrag replacements ασ Figure 7: (a) Plot of the maximal Lyapunov exponent σ against the bifurcation parameter ǫ for α = 0 .
15. Other parametersare same is in Fig.6. (b) A blow up of the previous figure is shown where the period doubling points are marked with arrows.Note that after each period doubling bifurcation, the system transits to a stable limit cycle (characterized by large negative σ before getting unstable for the next period doubling point when σ becomes close to zero. (c) The Maximal LCE in a regionof pitchfork bifurcation. The first arrow (up) indicates the pitchfork bifurcation after which the system transits to a chaoticregime before becoming periodic (denoted by the negative LCE) and becomes chaotic again as α increases, following a perioddoubling cascade, denoted by the second arrow (down). the maximal Lyapunov exponent σ . It has also been found that the system can be completely deterministic in themiddle of two chaotic regions and exhibit quasi-periodicity, as shown by the orbit diagram where a 5-period windowappears in the middle of a chaotic region. We have further shown that when the parametric forcing term related todust-charge fluctuation is small, away from the chaotic regime, the system can be driven in a frequency-locked statewhen a harmonic resonance of 2:1 takes place between the driving frequency and the fundamental frequency of thesystem. We have found that in most of the cases, the system transits to chaos through a cascade of period doublingbifurcations and the scaling of the period doubling cascades closely agrees to that of 1-D maps [24]. [1] Y. Saitou Y and T. Honzawa, in Proceedings of the Int. Cong. Plasma Phys. and 25th EPS Conf. Fusion Plasma Phys. (Prague, Czech Republic, 1998), Vol 22, pp. 2521-2524.[2] M. Faraday, Phil. Trans. Royal Soc. , 299 (1831).[3] E. Mathieu, J. Math. , 137 (1868).[4] A. H. Nayfeh and D. T. Mook, Nonlinear Oscillations (Wiley Interscience, NY, 1979).[5] R. H. Rand, A. Barcilon, and T. M. Morrison, Nonlin. Dyn. , 411 (2005).[6] R. H. Rand and T. M. Morrison, Nonlin. Dyn. , 195 (2005).[7] R. H. Rand, A. Barcilon, and T. M. Morrison, Nonlin. Dyn. , 411 (2005).[8] S.-Y. Kim and B. Hu, Phys. Rev. E , 3028 (1998).[9] P. K. Shukla and A. A. Mamun, Introduction to Dusty Plasma Physics (IOP Publishing Ltd., 2002).[10] M. R. Jana, A. Sen, and P. K. Kaw, Phys. Rev. E , 3930 (1993). (a) xy x = − . , y = − . x = − . , y = − . tx (c) t mod Tx (d) xy Figure 8: (a) Projection of the phase portrait in the x - y plane in the chaotic regime of Eq.(10). (b) The corresponding timeevolution of the variable x . Note the divergence of the two curves for change in the initial condition ( x , y ) at t = 0 at thefifth decimal place, showing the sensitivity of the system on initial conditions in the chaotic region. (c) The phase portrait inthe cylindrical space x × y × ( t mod T ), where T is the fundamental period of the system, which is 6 . α = 0 . , ǫ = 3 . , λ = 1 , ω d = 1, and ν = 0 . T = 6 . , F473 (2007).[12] D. Armbruster, M. George, and I. Oprea, Chaos , 52 (2001).[13] R. Bl¨umel, E. Bonneville, and A. Carmichael, Phys. Rev. E , 1511 (1998).[14] J. Jeong and S.-Y. Kim, J. Korean Phys. Soc. , 393 (1999).[15] M. Pandey, R. Rand, and A. Zehnder, in Proceedings of Proceedings of ASME 2005 International Design EngineeringTechnical Conferences (Long Beach, California, 2005), DETC2005-84018.[16] S. H. Strogatz,
Nonlinear Dynamics and Chaos (Addison-Wesley, 1994).[17] M. Zalalutdinov et al. , App. Phys. Lett. , 3281 (2003).[18] I. Bove, S. Boccaletti, J. Bragard, J. Kurths, and H. Mancini, Phys. Rev. E , 016208 (2004).[19] J. Guckenheimer and P. Holmes, Nonlinear Oscillations, Dynamcal Systems, and Bfurcations of Vector Fields (Springer-Verlag, NY, 1983).[20] J. Hale and H. Ko¸cak,
Dynamics and Bifurcations (Springer-Verlag, NY, 1991).[21] V. I. Arnold,
Ordinary Differential Equations (MIT Press, Cambridge, 1973), pp. 114.[22] E. J. Doedel et al. , AUTO 2000 : Continuation and Bifurcation Software for Ordinary Differential Equations . available online at : http://indy.cs.concordia.ca/auto/ [23] B. Ermentrout,
Simulating, Analyzing, and Animating Dynamical Systems : A Guide to XPPAUT for Researchers andStudents (Cambrdige, 1987).[24] M. J. Feigenbaum, J. Stat. Phys. , 25 (1978); , 669 (1979). [25] I. Shimada and T. Nagashima, Prog. Theor. Phys. , 1605 (1979).[26] A. Wolf, J. B. Swift, H. L. Swinney, and J. A. Vastano, Physica D , 285 (1985).[27] F. Grond et al. , Chaos. Sol. Fract.16