Investigation of Parametric Instability of the Planetary Gear under Speed Fluctuations

Planetary gear is widely used in engineering and usually has symmetrical structure. As the number of teeth in contact changes during rotation, the time-varying mesh stiffness parametrically excites the planetary gear and may cause severe vibrations and instabilities. Taking speed fluctuations into account, the time-varyingmesh stiffness is frequencymodulated, and therefore sideband instabilities may arise and original instabilities are significantly affected. Considering two different speed fluctuations, original and sideband instabilities are numerically and analytically investigated. A rotational lumped-parameter model of the planetary gear is developed, in which the time-varying mesh stiffness, input speed fluctuations, and damping are considered. Closed-form approximations of instability boundaries for primary and combination instabilities are obtained by perturbation analysis and verified by numerical analysis.The effects of speed fluctuations and damping on parametric instability are systematically examined. Because of the frequency modulation, whether a parametric instability occurs cannot be simply predicted by the planet meshing phase which is applicable to constant speed. Besides adjusting the planet meshing phase, speed fluctuation supplies a new thought to minimize certain instability by adjusting the amplitude or frequency of the speed fluctuation. Both original and sideband instabilities are shrunken by damping, and speed fluctuation further shrinks the original instability.


Introduction
Planetary gears are widely used in power transmission because of their compact design, high efficiency, and reduced noise.As the number of teeth in contact changes during rotation, gear mesh stiffness varies periodically with time.This parametric excitation is a primary source of vibration and noise, causing severe vibrations and instabilities under certain operating conditions [1].When the parametric excitation interacts with clearance nonlinearity, complicated nonlinear behaviors such as jump phenomena and secondary resonances are observed in planetary gears [2,3].Moreover, when the planets are equally spaced or diametrically opposed, this structural symmetry will lead to structured vibration characteristics [4,5], which has a great effect on the parametric instability of the planetary gear caused by time-varying gear mesh stiffness [6].Therefore, it is of great significance to determine the operating conditions of parametric instability and identify parameters that minimize the occurrence.
Parametric instability in gear system with constant speed has been investigated extensively.Tordion and Gauvin [7] and Benton and Seireg [8] analyzed the instabilities of the same two-stage gear system but derived contradictory conclusions.Lin and Parker [9] clarified the conflict and derived simple design formulas to control particular instabilities.For planetary gears, the structural symmetry results in highly structured modal properties [4,5].Based on the unique properties and ignoring damping, Lin and Parker [6] obtained the expressions of instability boundaries with constant speed using the perturbation method.In practice, planet meshing phase rules are often applied to neutralize the resonant response where the mesh frequency is near a natural frequency [10][11][12][13].Considering the elastic deformation of the ring gear, Parker and Wu [14] investigated the parametric instability with an elastic-discrete model.
The investigations mentioned above assume that the input rotating speed is constant.However, speed fluctuation is unavoidable in practice such as the engine output speed and the wind speed and induces frequency modulation of the gear mesh stiffness in gear systems.Parametric instability of the single-mesh gear system under speed fluctuation has attracted more and more attention [15][16][17].The instability boundaries are numerically and analytically determined, and different speed fluctuation types are considered.However, investigations on the parametric instability of planetary gears under speed fluctuations are relatively scarce.Ignoring damping, Qiu et al. [18] numerically calculated a primary instability of the planetary gear under engine speed fluctuation.Analytical expressions of parametric instabilities were not obtained, and the influence rules of fluctuation parameters on the instabilities were not derived.
The objective of this investigation is to systematically analyze the parametric instabilities induced by two different speed fluctuations and present a new way to control the parametric instability by adjusting speed fluctuation parameters.Perturbation analysis is conducted to determine operating conditions leading to instabilities and the results are verified by numerical integration.A pure rotational model of the planetary gear considering speed fluctuation and damping is introduced first.The parametric instabilities are then numerically and analytically investigated.Finally, the influences of various parameters on instabilities are investigated, such as the frequency and amplitude of speed fluctuations and system damping.
The main differences between this paper and [18] are the following: (1) two speed fluctuation types are considered, and their different influences on parametric instability are investigated; (2) damping is included in this paper, and the combined effect of damping and speed fluctuation is indicated; (3) besides numerical investigation, analytical expressions of primary and combination parametric instabilities are derived; (4) the influence rules of fluctuation amplitude and frequency on the instabilities are concluded, and a new way to adjust instability is presented.

System Model
The analysis deals with the parametric instability of the planetary gear subjected to input speed fluctuations.A rotational lumped-parameter model of the planetary gear is shown in Figure 1.All components are modeled as rigid bodies with moments of inertia   ,   ,   , and   ( = 1, 2, . . ., ,  is the number of planets).The subscripts , , , and  denote the carrier, ring, sun, and the th planet, respectively.The circumferential angle of the th planet is represented by   .Only gear rotational displacements  ℎ =  ℎ  ℎ (ℎ = , , , 1, . . ., ) are considered, where  ℎ is the base circle radius and  ℎ are the rotations in radian.Sun-planet and ring-planet gear meshes are modeled as linear springs acting along the line of action and are denoted by   () and   (), respectively.
For spur gears, gear mesh stiffness is usually approximated as rectangular wave and expressed in Fourier series as [6] ( ()  sin  +  ()  cos ) .
( The calculation of  under input speed variation is stated in detail as follows.In general, the input speed for rigid body conditions can be introduced via the Fourier series as [18] where Ω 0 is the nominal input speed and a small parameter  is defined to indicate the amplitude of the speed fluctuation: In (3),  = Ω 0 is the nominal mesh frequency without considering the speed fluctuation, and  = /  depends on the amplitude of the speed fluctuation and the ratio between the nominal input speed and the fluctuation frequency.In planetary gear,  is determined by the configuration and the tooth number of the central gears.With a fixed ring,  =   when the carrier is the input element, and  =     /(  +  ) when the sun is the input element.Substituting (3) into (1), sun-planet and ring-planet mesh stiffness can be rewritten as where the two small parameters  1 =  V /  and  2 =  V /  are defined to indicate the amplitudes of the gear mesh stiffness variations.It can be clearly seen that the gear mesh stiffness is frequency modulated by the speed fluctuation.The time history and the spectrum of an example gear mesh stiffness under speed fluctuation are shown in Figure 2. Because of the frequency modulation induced by speed fluctuations, sideband frequencies are introduced and symmetrically distributed on both sides of the harmonics of the nominal Frequency (Hz) gear mesh frequency , and the amplitude at nominal gear mesh frequency  changes as well.In contrast to the single frequency excitation without considering speed fluctuations, more instability regions will arise due to sideband frequencies, and the original instability (which is similar to that of constant speed) will be affected.
Based on the relationship between the fluctuation frequency   and the nominal gear mesh frequency , the input speed fluctuations can be classified into two categorizations: in the first type, the fluctuation frequency   is constant; in the second type, fluctuation frequency   is proportional to the nominal mesh frequency ; that is,  =   , such as the engine speed fluctuation [15,16].
Through force analysis, the equations of motion of the planetary gear can be derived, and the system stability is governed by the free vibration equation [19].The equation is applicable for general 2K-H planetary gear, and input and output component is not restricted: x = [  ,   ,   ,  1 , . . .,   ]  , Table 1: Simulation parameters of an example planetary gear with three equally spaced and in-phase planets.

Number of planets 𝑁 = 3
Inertia (kg) Substituting (4) into (9) and letting  =  1 =  2   /  (  and   are contact ratios), the time-varying stiffness matrix can be rewritten as where K 0 is the time-invariant stiffness matrix with the inclusion of the average gear mesh stiffness and K () V1 and K ()

V2
are the Fourier coefficient matrices.Damping is introduced via the modal damping ratio and expressed as where   ( = 1, 2, . . .,  + 3) are the modal damping ratios.The modal matrix V and the natural frequencies   are calculated by solving the eigenvalue problem The vibration matrix V is normalized as

Numerical Analysis
In the case of steady speed, parametric instability occurs in the vicinity of the critical frequencies defined as /(  +   ) = 2/ ( = 1, 2, 3, . ..),where   and   are the natural frequencies of the planetary gear [6].If  = , the situation is defined as primary instability.If  ̸ = , the situation is defined as combination instability.As high-order instabilities have much smaller instability regions and are unlikely to occur in practice [6], the following analysis focuses on  = 1.
An example planetary gear with equally spaced and inphase planets is used for numerical simulation.The ring is fixed to the gearbox housing, and its vibration is constrained to be zero.Damping is ignored for the example numerical analysis.Simulation parameters are listed in Table 1.The operating conditions leading to parametric instabilities can be derived by calculating the free vibrations under nontrivial initial conditions using numerical integration.As shown in Figure 3, if the amplitude of the response diverges, the response is unstable.Otherwise, the response is stable.Based on this, the instability regions of the planetary gear under two different speed fluctuations can be numerically obtained.
Because of the unique cyclic symmetry, the vibration modes of the planetary gear with three equally spaced planets can be classified into 1 rigid body mode ( 1 ), 2 distinct modes ( 2 ,  5 ), and 2 degenerate modes ( 3 =  4 ).As shown in Figure 4, the planets have identical motions in the distinct modes, and the central components have no motions and the motions of the planets differ in the degenerate modes.
Because of this unique modal property, whether certain instability occurs under constant speed can be directly estimated by planet meshing phase [6]: when the planet meshes are sequentially phased, primary instability and combination instability of the distinct modes are suppressed; when the planet meshes are in-phase, combination instability of distinct and degenerate mode is suppressed.Therefore, for the example planetary gear with in-phase planets, parametric instability does not occur in the vicinity of  2 + 3 and  3 + 5 without considering speed fluctuation, as shown in Figure 5.
Responses of the planetary gear at point A ( = 3787 Hz which is in the vicinity of  3 +  5 ) in the presence of speed fluctuation are compared with those for constant speed.For constant speed, the vibration displacement and spectral amplitudes converge, as shown in Figure 6.That is, parametric instability does not occur at point A for constant speed.Because of the influence of the frequency modulation induced by speed fluctuations, instability may occur in the frequency range where no instability occurs under constant speed.As shown in Figure 7, parametric instability occurs at point A with the inclusion of speed fluctuation.Compared with the constant speed case shown in Figure 6(b), more frequency components occur in the power spectrum, and the spectral amplitude increases with time.Therefore, whether instability occurs at certain nominal gear mesh frequency cannot be simply predicted by the planet meshing phase.
The stability diagrams in the vicinity of 2 5 for two different speed fluctuations are shown in Figure 8.The width of original instability which is similar to that of the constant  speed case is changed, and sideband instabilities arise as well.These changes are closely related to the frequency modulation of the gear mesh stiffness [17]: (1) the amplitude at nominal mesh frequency changes; (2) sideband frequencies have influence on original instability and may be in the frequency range of parametric resonance, generating sideband instabilities.Comparing Figure 8(a) with Figure 8(b), it can be clearly seen that the original instability is approximately symmetrical, while the sideband instability at higher frequencies leans to high frequency direction and that at lower frequencies leans to low frequency direction.Moreover, the distributions of instabilities for the two fluctuations are different.For constant speed fluctuation frequency   , the sideband instabilities are equally spaced around the original instabilities and the interval is   , while for  =   , sideband instabilities distribute more densely in the low frequency region.
The characteristics of the speed fluctuation have a great effect on the widths of the instability regions when planet meshing phases are determined, and the influence rule is closely related to the amplitude and frequency of the speed fluctuation.As shown in Figure 9(a), the sideband instabilities increase with the speed fluctuation amplitude , while the original instability decreases with the speed fluctuation amplitude .As shown in Figure 9(b), the original instability is relatively small and first decreases and then increases with the speed fluctuation amplitude .The phenomena shown in Figure 9 inspire a new thought to decrease the instabilities of the planetary gear with the inclusion of the speed fluctuation.In order to systematically investigate the parametric instability of the planetary gear under speed fluctuation and present a new way to control instabilities, analytical investigation is conducted using the method of multiple scales (MMS).

Perturbation Analysis
Substituting the modal transformation x = Vz into (5), the free vibrations of the planetary gear with a fixed ring are transformed into modal response and are expressed as where the matrices D () = V  K () V1 V and E () = V  K () V2 V and   = 2    .Because  = 1 responds to the rigid body mode which does not affect the parametric instability of the planetary gear, it is not considered in the following analysis.Using the method of multiple scales, the solutions of ( 12) can be expressed as [19,20] where  0 =  and  1 = .Substituting (13) into ( 12) and making the coefficients of the same power in  of both sides equal yield

Shock and Vibration
The general solutions of ( 14) are where c.c. represents the complex conjugate of the preceding terms.Because of speed fluctuations, sin() and cos() are no longer the standard Fourier series.With the aid of the Bessel function of the first kind [21], they can be expanded in generalized Fourier series as Substitution of ( 16) and ( 18) into (15) yields In the following perturbation analysis, primary instabilities of distinct modes, primary instabilities of degenerate modes, and combination instabilities of distinct modes are considered.Other instabilities can be obtained following the similar procedure listed below.
(1) Primary Instability of Distinct Modes.Let  +   = 2  +  ( = 0, ±1, ±2, . ..),where   is distinct and  is the detuning parameter to be determined.Elimination of secular terms in (19) requires Substitution of   ( 1 ) = (1/2)( 1 ) ( 1 ) +c.c.into (20) yields where  =  1 − .To obtain the steady-state motion, let  1  =  1  = 0, and then  can be determined: If the speed fluctuation frequency   is proportional to the nominal gear mesh frequency ,  =  does not change with , and thus the instability boundaries can be easily expressed as If   is constant,  = /  varies with , and the instability boundaries can be obtained by solving the following equation: The case  = 0 corresponds to the original instability, and the cases  ̸ = 0 correspond to sideband instabilities.From ( 23) and ( 24), it can be seen that the intersections of instability boundaries and the abscissa axis are 2 0 /(1 + /) and 2 0 −   , respectively, which explains the different distributions of instabilities for the two fluctuation types.Moreover, with the increase of   , the interval of original and sideband instability becomes larger.With the increase of , the original and sideband instability region becomes more intense.
Once the detuning parameter  is determined, the calculations of the instability boundaries are similar for all instabilities.Therefore, only the calculation of  is emphasized in the following analysis.
(2) Combination Instability of Distinct Modes.Let  +   =   +   +  ( = 0, ±1, ±2, . ..),where   and   are both distinct.Elimination of secular terms in (19) requires Substituting   ( 1 ) = (1/2)  ( 1 )   ( 1 ) + c.c. and   ( 1 ) = (1/2)  ( 1 )   ( 1 ) + c.c. into (25), the detuning parameter  is derived as (3) Primary Instability of Degenerate Modes.Let  +   =   +  + ( = 0, ±1, ±2, . ..),where   =   are degenerate.Elimination of secular terms in (19) requires Because modal coupling between the degenerate modes is stronger than that between distinct modes, the determination of the instability boundaries is more complex, as shown in (27).Following the standard procedure, the detuning parameter is expressed as To validate the analytical results, comparisons of the numerical and analytical results using the same parameters are shown in Figure 10.It can be clearly seen that the perturbation solution for the original instabilities matches well with the numerical results for the two speed fluctuation types.22), (26), and (29), planet meshing phase changes the instability regions by influencing the elements of matrices D and E.Besides planet meshing phase, it can be clearly seen that speed fluctuations have a significant influence on the width of instability region by changing the value of   ().That is, speed fluctuation amplitude  and frequency   have a combined effect on the parametric instability width, and the influence rule is closely related to the properties of the Bessel function of the first kind.

Influence of Speed Fluctuations. As shown in (
As shown in Figure 11,  0 () shows cosine oscillation attenuation, while the rest shows sine oscillation attenuation.This property explains the variation tendency of the original and sideband instabilities in the numerical analysis.For the first speed fluctuation type, where the speed fluctuation frequency   is constant, the influences of  and   on original instabilities in the vicinity of 2 3 ,  3 +  5 , and 2 5 are shown in Figure 12.All original instabilities decrease with the increase of  and increase with the increase of   under the simulation conditions.Because   is constant, the ratio of the nominal mesh frequency  and   varies with .The change of instability width differs in different frequency range.With the definition of instability ratio  = width ( ̸ = 0)/width ( = 0), the instability ratios of 2 3 ,  3 +  5 , and 2 5 for   = 200 Hz and  = 0.05 are 0.852, 0.826, and 0.715, respectively.
For the second speed fluctuation type, where  =   , the influences of fluctuation amplitude  and frequency ratio  on original instabilities of 2 3 ,  3 + 5 , and 2 5 are shown in Figure 13.Similarly, all original instabilities decrease with the increase of  and  under the simulation conditions.Because the frequency ratio  does not change with the nominal gear mesh frequency, instability ratios of 2 3 ,  3 +  5 , and 2 5 are identical once the fluctuation amplitude is determined, which is different from the first fluctuation type.
Through the above parameter analysis, it can be found that the properties of the Bessel function of the first kind with different order supply a new possibility to minimize certain parametric instability by adjusting the parameters of the speed fluctuation.The adjusting rule is to minimize the value of   ():  = 0 corresponds to the original instability, and  = ±1, ±2, . . .correspond to the sideband instabilities.It should be noted that the influence trends of speed fluctuation on original and sideband instabilities may be different.Based on the first-order perturbation solution, the instability width with different speed fluctuation parameters can be easily obtained.For the second speed fluctuation type, variations of original and sideband instability width of 2 5 with speed fluctuation parameters are shown in Figure 14, for example.Based on these figures, speed fluctuation parameters can be easily determined to minimize certain instability.Comparing point B ( = 24.05 and  = 0.1) with point C ( = 24 and  = 0.05), although speed fluctuation amplitude of point B is larger than point C, the corresponding instability region of point B is much smaller.The numerical results in Figure 15 show the similar decrease when the speed fluctuation is changed from point C to point B.

Combined Effect of Damping and Speed
Fluctuations.In all of the above calculations, damping is constrained to be zero.Taking damping into account, the combined effect of damping and speed fluctuations on original instabilities is shown in Figure 16.It is well known that damping shrinks the instability region for the constant speed.With the inclusion of speed fluctuations, original instability regions are further shrunken.This is because  0 () obtains its maximum value when the speed is constant ( = 0).The minimum value  min of the mesh stiffness variance corresponding to the occurrence of original instability is affected by damping and speed fluctuation.As shown in Figure 16, the minimum value for fluctuating speed is larger than that for constant speed.This phenomenon can be easily explained by analytical solutions (see ( 22), (26), and (29)), and speed fluctuation offers a new way to suppress the occurrence of parametric instability under determined mesh stiffness variation besides increasing the system damping.Taking point D ( = 3393 Hz and  = 0.13) in Figure 16 as an example, parametric instability occurs at point D for constant speed ( = 0), while the instability is suppressed when speed fluctuation is added ( = 0.05).Comparing the minimum stiffness variation  min in three instability regions, it can be found that speed fluctuation has greater influence on the larger instability region.
The influence of speed fluctuation amplitude  on the original and sideband instabilities in the presence of damping is shown in Figure 17.It can be clearly seen that the original instability shrinks, while the sideband instabilities expand with the increase of  for the two speed fluctuations.Because  − () = (−1)    () and  =   , the minimum values of  for the corresponding sideband instabilities for the second fluctuation type are equal, while for the first fluctuation type, the ratio of the nominal gear mesh frequency  and fluctuation frequency   changes with , and therefore the minimum values of  for the corresponding sideband instabilities are not the same.

Conclusions
In this study, parametric instabilities of the planetary gear under two different speed fluctuations are systematically investigated.A rotational model of the planetary gear and the mesh stiffness modeling are introduced first.Perturbation analysis is then conducted to determine operating conditions leading to instabilities and the results are verified by numerical integration.Finally, the influences of speed fluctuation parameters and damping on instabilities are investigated, and a new way to control instabilities by adjusting speed fluctuation is proposed.Main conclusions are summarized as follows: (1) Speed fluctuations induce frequency modulation of the gear mesh stiffness and then cause sideband instabilities on both sides of original instabilities and greatly influence the widths of original instabilities.
(2) Because of the influence of sideband frequencies, whether parametric instability occurs at certain nominal gear mesh frequency cannot be simply predicted by the planet meshing phases that are applicable to constant speed.
(3) The influence of speed fluctuations on instabilities is determined by the Bessel function of the first kind with different order.Original and sideband instabilities can be controlled by changing the value of the corresponding Bessel function, which is realized by adjusting the amplitude and frequency of speed fluctuations.
(4) Damping and speed fluctuations have a combined effect on the occurrence of original instabilities.When the mesh stiffness variation is determined, original instabilities can be suppressed by introducing proper speed fluctuation besides increasing damping.

Figure 1 :
Figure 1: Rotational lumped-parameter model of the planetary gear.

Figure 2 :
Figure 2: (a) Time history and (b) spectrum of an example sun-planet mesh stiffness under input speed fluctuation with  = 2000 Hz,   = 200 Hz, and  = 0.5.

Figure 4 :Figure 5 :
Figure 4: Modal properties of the planetary gear with three equally spaced planets: (a) distinct mode; (b) degenerate mode.Dashed lines are the equilibrium positions and solid lines are the deflected positions.

Figure 6 :
Figure 6: Dynamic response of the planetary gear at point A: (a) time history; (b) power spectrum.

Figure 7 :Figure 8 :Figure 9 :
Figure 7: Dynamic response of the planetary gear at point A for  = 0.05 and   = 200 Hz: (a) time history; (b) power spectrum.