Primary Parametric Amplification in a Weakly Forced Mathieu Equation

Since Eq. (1) is linear, the forced system is nonhomogeneous and its solution can be expressed as a sum of homogeneous and particular solutions, given as q = q_h + q_p, where q_h satisfies the damped Mathieu equation, $q¨h+2ϵμq˙h+(ω2+ϵγcosΩt)qh=0$⁠, and q_p satisfies Eq. (1). There is a large body of work aimed at understanding the unforced Mathieu equation, and hence q_h. It is well known to have instability regions in the parametric forcing amplitude–frequency parameter space, such that the fixed point at the origin can be stable or unstable, even if the damping and time-varying stiffness are positive. These stability wedges can be studied by applying the Floquet theory [12] with harmonic balance solutions of periodic transitional motions [13] and also by a higher-order perturbation expansion [14]. For lightly damped cases, these regions of instability are based near frequency ratios $Ω/ω=2/N$⁠, N being an integer. The dependence of the instability regions on the strength of damping has been further explored in detail in Ref. [15]. Free responses, whether stable or unstable, have also been analyzed based on the Floquet theory [16].

By linearity, if F(t) in Eq. (1) is scaled by a factor α, then the forced or particular solution q_p is also scaled by α. (However, scaling the level γ of parametric excitation does not simply scale the response.) The solution to the forced Mathieu equation can be approximated by perturbation techniques to reveal resonance conditions as well as instabilities. Our previous work [17] dealt with the superharmonic and subharmonic resonances for the forced Mathieu equation with cubic nonlinearity. A first-order perturbation analysis using the method of multiple scales revealed the existence of superharmonics at one-half and one-third the natural frequency and subharmonics at twice and thrice the natural frequency. The first-order analysis showed that the superharmonic at order 2 persists for the linear system. Numerical simulations of the linear system showed a noticeable stretch in the primary resonance (Fig. 1), suggesting a parametric amplification, as well as a superharmonic of order 3, which were not uncovered by the first-order perturbation multiple-scales analysis [17]. The aim of this study is to apply second-order multiple scales to unfold the parametric amplification at primary resonance.

Resonance conditions in parametrically excited systems with forcing and also nonlinearity have been examined [18], and there have been other studies on dynamics, stability, and bifurcations of nonlinear generalizations of forced Mathieu-type equations [19–24].

Although our interest began with wind turbines, such structures are unlikely to operate near primary resonance, unless there is an uncontrolled runaway incident. However, in this work, our interest is in the forced Mathieu equation as a fundamental dynamical system. An example area of study has been subharmonic parametric amplification (with primary direct excitation), which has undergone numerous analytical and experimental studies [8–10,25,26] in microelectromechanical systems (MEMS) and nanoelectromechanical systems with applications [27–30], including with arrays of resonators [31–33]. These examples involve a 2:1 ratio of parametric to direct excitation frequencies, and as such, the subharmonic parametric resonances are exploited. In contrast, Eq. (1) has a 1:1 ratio between the frequencies of the parametric and direct excitation terms. Conceivably, other configurations of MEMS resonators might involve 1:1 excitations and be conducive to primary parametric amplification.

We aim to apply higher-order perturbation techniques to unfold the existence and features of primary parametric amplification in the forced Mathieu equation. Higher-order perturbation techniques have been used extensively across various sciences to get more accurate expressions of the solutions to problems, to derive an understanding of the underlying physics [34,35]. Dynamical systems have been described by higher-order perturbation expansions [14,36,37] as well as other scaling techniques [38].

Hence, in the current work, we analyze Eq. (1) under hard forcing with two orders of multiple-scale expansion as formulated by Nayfeh [36] to capture the primary resonance and its features. The primary parametric amplification was introduced in Refs. [3,39] and is further developed here to describe the amplification and possible suppression, phase effects, and the effects of parameters.

We have three cases for Eq. (6), which can contribute secular terms:

1. No specific relationship between $Ω$ and ω

2. $Ω≈2ω$

3. $Ω≈ω$

The treatment of these cases will be laid out in the following sections. In each case, the second-order expansions will give us the slow time scale variations of the amplitude A of our fast time scale solution q₀. A second-order expansion will involve solvability conditions obtained by eliminating secular terms at both O(ε) and O(ε²).

But there are also specific combinations of $Ω$ and ω appearing at the second-order that could lead to resonance conditions in the system, and thereby add terms to the solvability condition at O(ε²). These can be identified by examining the remaining terms of the O(ε²) equation given by Eq. (10). We arrive at the following candidates:

• $Ω≈ω/2$

• $Ω≈2ω$

• $Ω≈ω$

The first case was treated with a first-order multiple-scales analysis in Ref. [17], in which hard excitation-induced secular terms at the O(ε) expansion. Similarly, in this case, the excitation, this time at O(ε), induced secular terms one order down at O(ε²). Hence, similar results are expected, and so this case is not emphasized in this article. The latter two cases actually produce secular terms at O(ε) and are therefore cannot be treated within case 1, since it is a nonresonant case at O(ε).

It is noticed that F does not appear in the solvability conditions and thus does not associate with a forced subharmonic response. There is a small linear forced-response term in Eq. (12), representing the forced component of the response. The solvability conditions thus provide the dynamic characteristics of the homogeneous solution and its dependence on the parametric excitation.

The expression for D is quadratic in γ², which when equated to zero gives the stability boundaries of the order 2 subharmonic under hard excitation. These agree with analyses in the existing literature [13,14], and we do not present further details of this analysis. A first-order analysis, in fact, is enough to approximate the stability boundaries (from Eq. (11)), and the second-order analysis refines the expressions of the stability boundaries.

We substitute the solution for q₀ and q₁ into the Eq. (4) at O(ε²). After expanding the right-hand side of the resulting q₂ differential equation, we seek to identify secular terms at this resonance condition.

The second-order analysis indicates that the parametric excitation affects the resonant amplitude at a very slow time scale, i.e., in the ε² terms in the reconstituted differential equation for A.

In this section, we study the amplitude and the phase of responses at primary resonance. Transforming Eq. (20) to polar coordinates via $A=12aeiβ$ and letting β = σT₁ − ϕ leads to slow-flow equations in a and ϕ. The steady-state equations include sin (ϕ), cos (ϕ) and sin (2ϕ), andcos (2ϕ) terms, which make them more difficult to solve than other familiar examples.

We have noted that γ affects the response amplitude at the second-order slow time scale. However, the steady-state response amplitude, via the expression for r², is at order O(ε⁰).

Equation (24) provides information related to the amplitude of the primary resonance as a function of parameters, including the detuning parameter σ and the phase angle θ between the parametric and direct excitation terms. Since α and δ depend on θ, the response amplitude is also dependent on θ.

Figure 2(a) plots the amplitude a = 2r, where r comes from Eq. (26) for θ = 0 (meaning the direct and parametric excitations are out of phase by π/2, since they are present as sine and cosine in Eq. (2)), for nondimensional values of parametric excitation level γ ranging from γ = 0 to γ = 4. The parameter values are ε = 0.1, ω = 1, μ = 0.25, and F = 5, and also nondimensional. The plot shows a clear resonant peak. This peak occurs at a value of σ_p (equivalently $Ωp=ω+ϵσp$⁠), which is very close to σ = 0 when γ is small, as γ increases, the peak moves to lower values of $Ω$⁠. When γ < 1.5, the resonance curve is indistinguishable from the curve of the regular forced linear oscillator. When γ = 1.5, the curves have a slight perceptible separation, which grows quadratically as γ increases further, indicating a parametric amplification of primary resonance. The case of γ = 4 is cropped for better visualization of the other curves, but remains bounded for all $Ω$⁠. For slightly larger values of γ, the primary resonance becomes unbounded for some values of $Ω$⁠. We will further discuss this shortly. For the cases of γ = 3 and γ = 3.5, simulation values are plotted. These were obtained by numerically solving the differential equation (using matlab ode45) and examining the fast Fourier transform (FFT) of 25 cycles of the response at the steady state. The sampling interval was 1/100 of the excitation period, such that leakage issues did not distort the FFT. The FFT indicated a strong peak at the primary resonance excitation frequency (close to the natural frequency) and a minor peak at its first harmonic. The amplitude of the primary term only is plotted in Fig. 2 and subsequent plots. The next harmonic (at twice the natural frequency) can be seen in the q₁ expression in Eq. (17). For example, when ε = 0.1, γ = 3, and $Ω=1$⁠, this term should be about 1/20 of the primary resonance amplitude; this is very consistent with the detected first and second harmonic amplitudes of 13.97 and 0.6988. The simulations shown have excellent agreement. Not shown are simulations samples for γ = 4. At this value of γ, some error can be detected. For example, when $Ω=ω$⁠, there is a 4.6% error between the simulated response amplitude and that predicted by the second-order perturbation expansion.

Figure 2(b) shows the associated phase lag ψ behind the direct excitation for the same parameter values. The plot magnifies the phase transition region near the frequency ratio of $Ω/ω=1$⁠. For γ = 1, the phase lag overlaps with the nonparametric case (γ = 0). Increasing γ pushes the phase transition to a slightly lower frequency ratio, along with the trend of the response amplitude peak. Otherwise, the forced and parametrically excited system exhibits a phase lag, which approaches 0 for very low frequencies, approximately π/2 radians near resonance, and approaches π radians for very large frequencies.

Figure 3(a) plots the amplitude a = 2r, from Eq. (26) for θ = π/4 for nondimensional values of parametric excitation level γ ranging from γ = 0 to γ = 4, and for the same values of parameters as the previous figure. In this case, the plotted values of γ of 2 or less do not show a deviation from the regular forced linear oscillator, which is represented by the solid curve with a single local maximum near $Ω=1$⁠. As γ increases beyond a value of about 2, the resonance peak first decreases. Thus, this value of θ shows parametric deamplification. At the value of γ = 3.5, the resonance peak begins to split. The two peaks that form increase as γ continues to increase. Further slight increases in γ are accompanied by large growth in the amplitude before becoming unbounded. Simulations for the cases of γ = 3.5 and γ = 4 show strong agreement with the perturbation results. Figure 3(b) shows the associated phase lag ψ for θ = π/4 and the same parameter values. For γ = 1, the phase lag overlaps with the nonparametric case (γ = 0). Increasing γ causes the phase transition to become steeper, i.e., more abrupt with frequency changes near resonance. This trend, with increasing γ, is reminiscent of the trend seen in a regular (nonparametric) forced oscillator as damping decreases. Curiously, at this value of θ, the effect on the response amplitude (Fig. 3(a)) when increasing γ has a trend opposite that caused by decreasing damping. The phase steepening occurs with a crossover frequency, below which the phase is decreased and above which the phase is increased, with increasing γ. At the crossover point, the phase is independent of γ.

Insight can also be obtained from the plot of the complex amplitude in the complex plane. Figure 4 shows plots of X+iY, representing half complex amplitudes, phased relative to the direct excitation, for the cases of (a) θ = 0 and (b) θ = π/4, with various parametric levels γ. (In this case, the various levels of γ are plotted with solid curves to maintain a cleaner image.) The curves are parameterized by −1 < σ < 1 (equivalently $Ω$⁠). If the range σ were increased, the plots would approach closed curves. Figure 4(a) shows that the effect of γ is to extend the amplitudes from the γ = 0 reference, which is the smaller circle, to the oval curves that increase as γ increases. The plot shows that the maximum amplitude occurs with a phase lag ψ > π/2 radians (π/2 being the phase at resonance of the reference γ = 0 system). Figure 4(b) shows that as γ increases, the reference circle distorts into the oval curves that extend to both sides. From this, it can be seen that if γ is large enough, there can be two amplitude peaks associated with phase lags that are closer to 0 and π radians than to the reference phase of π/2. This feature can be cross checked in Fig. 3.

In both examples of Fig. 4, there are two crossover points, which indicate, for any γ, there are two frequencies that produce a common amplitude and phase. We have observed that these events occur, however, at different frequencies $Ω$ for different values of γ (not shown here).

Figure 5 delves more into the effect of θ on the resonance behavior, for a fixed value of γ = 4. The top left shows how the half complex amplitude (X+iY) varies with 0 < θ < π at increments of π/8 radians. θ = 0 is represented by the dashed curve and also shows up in Fig. 4(a). The oval rotates clockwise with increase in θ. At θ = π/4, the oval is small and horizontal, represented by the tightly hatched curve. Continued increase in θ rotates the oval and increases it size. The dotted curve represents θ = 3π/4 and leads to the maximum resonance amplification of the cases depicted. Continued increase to θ = π brings the ovals back to the dashed case, and the process repeats. The solid black circular curve shows the case of γ = 0 for reference.

The right column of Fig. 5 separates the amplitude and the phase of the response. The dotted, hatched, and dashed curves are consistent with those in the top left image, and the other curves show the transitions with varying θ. The solid black curve (lowest single maximum) shows γ = 0 for reference. The top right curves show that the resonance occurs at a slightly lower frequency than the regular (γ = 0) case. θ has a strong effect on the character of the peak. The case of θ = π/4 is the hatched curve with twin peaks. It can be cross checked in Fig. 3(a). The lower right image shows the associated phase lags. The case of θ = π/4 is the hatched (looks dotted in this image, however) curve that is steepest. It can be cross checked in Fig. 3(b). The γ = 0 reference is the solid black curve that transits π/2 radians at $Ω=1$⁠. The general, but not strict trend, is that for most θ, the phase transition occurs at slightly lower frequencies than $Ω=1$⁠, consistent with the amplitude graph. Also, θ can be associated with steeper transitions (e.g., θ = π/4) or softer transitions (e.g., the dotted curve of θ = 3π/4).

Steepening the phase curve, and amplifying the amplitude peak, might be suspected to be “equivalent damping reduction” effects. When θ = π/4, the amplitude peak is suppressed (like increased damping), while the phase curve is more steep (like decreased damping). When θ = 3π/4, both effects are reversed. Hence, if one tries to look for “equivalent damping effects,” the parametric resonance has conflicting patterns.

Thus, θ can have a significant effect on the parametric amplification [8,10,26]. In this case, depending on the value of θ, parametric-plus-direct excitation can lead to parametric amplification or suppression. We can define a gain as the ratio of the response amplitude a = 2r to that of the system without parametric excitation (γ = 0) [8]. For the latter, we use the resonant (⁠$Ω=ω$⁠) value for the simple harmonic oscillator, namely, $ah=F/2μω2$⁠. Figure 6 shows the dependence of the maximum gain on the phase θ, for a family of values of γ ranging from γ = 0 at the horizontal line (unit gain) to γ = 4 at the highest curve, with the same parameters as mentioned earlier. This figure was produced by choosing a value of θ and then sampling Eq. (24) to determine the maximum amplitude. The maximum amplitude is then normalized by a_h to produce the gain and plotted versus θ. The figure shows that the maximum amplification occurs near θ ≈ 3π/4, and the minimum amplification, indeed suppression, occurs near θ ≈ π/4. However, for most values of θ, the parametric effect amplifies the resonance. The curve is periodic in θ, with one period depicted. As γ increases, the amplification becomes larger (and in cases of deamplification, the suppression becomes less pronounced until it eventually shows slight amplification). At some point, γ becomes large enough that the maximum response becomes uniformly unbounded. (For these parameters, this occurs around γ =4.5.)

Figure 6(b) shows that as damping decreases, the maximum system gain uniformly increases. Much like when γ achieved a high enough value, it turns out that when μ reaches a small enough value, the system becomes unbounded.

Figure 7 shows the amplification as a function of parametric excitation level γ, for various values of θ. The graph indicates that for γ near zero, the gain has a zero slope, and toward a value of γ slightly larger than 4, the gain becomes unbounded (as it enters a parametric instability). Similar γ dependence has been observed for 2:1 parametric pumps [8,26,28]. For values of θ for which the system is parametrically de-amplified, with large enough γ, there is a parametric amplification up to a point at which the system becomes unbounded. In the figure, the lowest curve corresponds to θ = π/4, and the highest curve corresponds to θ = 3π/4, with equal increments between. Similar phase dependence has been observed for 2:1 parametric pumps [8,26,28]. For small γ, there is little amplification effect. As γ increases, the amplification increases sharply (or suppression develops gently) before becoming unbounded. There may be an opportunity to introduce nonlinearity to broaden the range of γ for significant amplification and provide robustness in applications [10].

The features of Figs. 2 and 3 show that the resonance may have a single peak or two peaks or may be unbounded. With damping present, we might in principle seek the presence of one or more peaks by finding zeros of d(r²)/dσ. However, this is a complicated calculation, and thus, we do not pursue it here. Furthermore, since the peaks can stray from $Ω=ω$ (σ = 0), it may not be useful to approximate the peak value with an $Ω=ω$ (or σ = 0), as is equivalently done to estimate the resonance peak of a lightly damped regular linear forced vibration system. The upshot is that it is not easy to find a simple expression for approximating the peak response under parametric amplification.

The equilibrium x₀ = A⁻¹F was found previously, and its amplitude r was studied. Since Eq. (29) is linear, dynamics of X and Y, and hence the leading-order slowly varying response amplitude $A=12a(T1,T2)eiβ(T1,T2)$⁠, about the equilibrium, and hence the stability, are governed by matrix A. In addition, the homogeneous solution has its own equilibrium at the origin. The homogeneous solution is that for which F = 0 in Eq. (2) and consequently F = 0 in Eq. (29). Thus, matrix A simultaneously governs the stability of the homogeneous solution and the characteristics of the particular solutions of Eq. (2).

The calculated value of γ_cr is the value of γ in Fig. 7 for which the maximum amplitude becomes unbounded for specific values of detuning parameter. The increase in amplitude until that point is a result of parametric amplification, analytically approximated by Eq. (26). The analytically predicted values of γ_cr are plotted versus σ in Fig. 8. Numerical simulations were also checked for stability in this parameter space, and stable (dots) and unstable (×’s) events are plotted along with the predicted stability boundary from Eq. (31) with good agreement. For the parameters of the plot, the minimum γ for which the system destabilizes is near γ = 4.5, which is consistent with the implied asymptote in Fig. 7.

Looking at Eq. (26), for the lightly damped case with sufficiently large γ, there are two values of σ (equivalently $Ω$⁠) for which the denominator is zero and the amplitude is unbounded. Thus, there are two frequencies, $Ω1$ and $Ω2$⁠, which meet the stability-wedge boundary. For $Ω1<Ω<Ω2$⁠, the response is unstable and unbounded. Fixing γ and increasing $Ω→Ω1$⁠, the response amplitude plot approaches a vertical asymptote, while decreasing $Ω→Ω2$⁠, the response amplitude plot approaches another vertical asymptote. The sides of the unbounded resonance peak surround an interval of values of $Ω$⁠, in contrast to a single frequency of unbounded responses in the regular (γ = 0) oscillator.

We have seen that the linear forced Mathieu equation has a homogeneous and particular solution. The homogeneous solution has the known behavior of the Mathieu equation. The particular solution is a forced response with the potential for a primary resonance, which was analyzed in this section by multiple scales. A second-order analysis shows that the primary resonance is stretched, or amplified (for most values of θ), by the parametric excitation. Given the parameter set, as the level of parametric excitation γ grows, the resonance amplitude is at first hardly affected by small γ, but then begins to grow sharply for larger γ, and approaches infinity at a critical γ_cr. The bounded resonance response is stable. As the response becomes infinite at γ_cr, the response destabilizes. At the same time, the homogeneous part of the response also becomes unstable as the parameters cross into the unstable part of the Arnold tongue.

We have done an expansion based on powers of ε and derived a predictive amplitude equation (24), which is a function of parameters. Equation (24) has a denominator which is a function of ε. There may be a question as to whether this expression for the amplitude should be expanded into powers of ε, retaining up to ε², for consistency with the perturbation expansion. It turns out, in our example, that such a treatment reduces the accuracy of the amplitude expression.

To check this, we have expanded Eq. (24) and omitted terms of O(ε³) and higher (not shown). For example, for μ = 0.25, ε = 0.1, F = 5, γ = 3 and σ = 0, the amplitude obtained from Fig. 1 at primary resonance (⁠$Ω≈ω$⁠) is a = 13.95. The numerical simulation of Eq. (2), as described earlier, led to the resonant harmonic amplitude of a_s = 13.97. However, the quadratic expansion of Eq. (24) produced an estimation of a_x = 13.3.

Thus, it seems that the reconstituted equation (20) has information of a second-order expansion, and the approximate steady-state solution (26) uses all of this information. A truncated expansion of (26) loses information and deviates from the true solution.

We have looked at the linear Mathieu equation with combined parametric and direct excitation at the same frequency. The forced linear Mathieu equation has homogeneous and particular solutions relative to consideration of the direct excitation term. The homogeneous solution has all of the characteristics shown in classical studies of the Mathieu equation. Thus, our attention was on the forced responses. In an earlier study, numerical simulations revealed resonance phenomena that were not uncovered by a first-order multiple-scales analysis, including a parametric amplification at primary resonance, noticeable to a trained eye.

In this study, we performed a second-order multiple-scales analysis with weak excitation to uncover the parametric amplification effects on primary resonance. The work complements parametric amplification studies exploiting subharmonic resonance. With a parametric excitation of $γcosΩt$ and a direct excitation of $Fsin(Ωt+θ)$⁠, the parametric amplification shows a stretching of the primary resonance peak with the increasing strength γ of parametric excitation, until a critical value γ_cr at which the resonance becomes unbounded. Amplitudes of the numerical solutions were compared with the analytical expression to very high agreement. The response phase lag relative to the direct excitation varies with the excitation frequency from zero to π radians as with a standard linear oscillator. The parametric term tends to push the resonance peak and phase-lag transition to a frequency slightly less than the natural frequency.

We looked at the role of a phase θ in the direct excitation term. Depending on parameters, it may be possible to choose a phase such that the parametric term suppresses the resonance. Under most phases, the parametric term has an amplifying effect. Variation in the response phase lag with respect to the direct-parametric phase difference was also examined. It is possible to have an excitation phase such that there is a frequency for which the response phase lag is independent of the parametric excitation level γ.

The amplified primary resonance goes unstable coincidentally as it becomes unbounded. This instability coincides with that of the 1:1 wedge in the unforced Mathieu equation. Thus, at primary resonance, the forced response and the homogeneous solution destabilize together.

Parametric amplification could be applied by design to accentuate (or suppress) resonance quality. However, care should be taken regarding the possibility of instability associated with the slender instability wedge associated with the unforced Mathieu equation near the natural frequency. The next study could address effects of stiffness nonlinearity on bounding the unstable behavior or bending the resonance, as done for the case when the parametric excitation has twice the frequency of the direct excitation [10].

The problem studied here shows an example where the second-order perturbation analysis uncovers behavior of nominally O(ε⁰) (or greater based on numerical results), but is not uncovered in the first-order perturbation analysis. The parametric effect at primary resonance occurs with a second-order slow time scale.

This material is relevant to projects supported by the National Science Foundation (CBET-0933292 and CMMI-1335177). Any opinions, findings and conclusions or recommendations expressed are those of the authors and do not necessarily reflect the views of the NSF.

There are no conflicts of interest.

The datasets generated and supporting the findings of this article are obtainable from the corresponding author upon reasonable request. The data and information that support the findings of this article are freely available. The authors attest that all data for this study are included in the paper.