Yair Zarmi
Jacob Blaustein Institutes for Desert Research, Be Gurion University of the Negev, Midreshet Ben Gurion, Israel.
DOI: 10.4236/am.2026.176023   PDF    HTML   XML   9 Downloads   75 Views   Citations

Abstract

The dissipative term in the Rayleigh oscillator contains a nonlinear (cubic) contribution. Usually, the coefficient, ε, of the dissipative term is assigned a small value. The solution evolves into a limit cycle with O(1) maximal amplitude, maximal velocity and period. In this paper, the case, in which this coefficient is extremely large (ε ≥ 1), is studied. As expected, for long times, the solution tends to a limit cycle. The characteristics of the solution are different from those encountered in energy conserving oscillatory systems in which the magnitude of the nonlinear term is increased indefinitely. Based on numerical solutions of the equation, it is found that the amplitude and period grow linearly with ε, whereas the maximal velocity remains O(1). The consequences of the linear growth with ε of amplitude and period are studied through the definition a scaled solution and a scaled time variable, an angle θ, in which the scaled form of the limit cycle is periodic.

Keywords

Rayleigh Oscillator, Large Dissipative Term

Share and Cite:

1. Introduction

The equation for Rayleigh oscillator,

$\ddot{x}\left(t\right)+x\left(t\right)=\epsilon \dot{x}\left(t\right)\left(1-\dot{x}{\left(t\right)}^2\right)$ (1)

was proposed by Rayleigh as the dynamical description of stable sound oscillations in musical instruments [1]. Traditionally, the coefficient, ε, is assigned small values (ε  1). Approximations to the limit-cycle solution are obtained by perturbation methods, such as the Method of Multiple Time Scales (see, e.g., [2]-[4] or the Method of Normal Forms (see. e.g., [5]-[8]]. Some works study the limit cycle solution (or of the equivalent Van der Pol oscillator) for ε = O (1) (see. e.g., [9] [10]. However, the effect on the solution in the limit (ε  1) has not been studied. It is analyzed in this paper.

As ε is increased, the amplitude, x_max, and the period of oscillations, T, grow linearly with ε, whereas the maximal velocity v_max, remains O (1). The variation of the amplitude near the turning point is very fast and the velocity varies very rapidly between its positive and negative extreme values (both O (1)).

The motivation for studying the ε  1 limit for Equation (1) is that this limit is the parallel of “nonlinear violence” in conservative oscillatory systems. In the latter case, the nonlinear term is increased indefinitely. The effect on the solutions in the two cases is qualitatively different. The case of conservative oscillatory systems has been studied in [11]. For the sake of showing the qualitative difference of the effect of nonlinear violence on conservative oscillatory systems and on the Rayleigh oscillator, the Duffing Equation:

$x\left(t\right)+\ddot{x}\left(t\right)\text{}+\epsilon x{\left(t\right)}^3=0$ (2)

is discussed in Appendix I.

2. Numerical Analysis

All numerical solutions of Equation (1) are obtained through Wolfram Mathematica NDSolve, which uses the Explicit Runge-Kutta method. The boundary conditions are x (0) = 1 and v (0) = 0. The program is run up to long times, where the limit-cycle pattern emerges and prevails.

The numerical solutions of Equation (1) have been obtained for 0 ≤ ε ≤ 2000. Figure 1 and Figure 2 present, respectively, the long-time dependence of x (t) and dx/dt (t) for ε = 2000.

The numerical solutions for a range of values ε provide the ε - dependence of the period of oscillations, T (Figure 3), maximal oscillation amplitude, x_max, (Figure 4) and the maximal velocity, v_max (Figure 5). Both are obtained through

Figure 1. Long time dependence of x (t), solution of Equation (1); ε = 2000. Dashed lines: ±x_max.

Figure 2. Long time dependence of dx/dt for solution of Equation (1); ε = 2000. Dashed lines: extreme velocity, ±v_max; dotted lines – velocity (±1) at x (t) = 0.

Figure 3. ε - dependence of period of oscillations, T.

Figure 4. ε - dependence of maximal oscillation amplitude, x_max.

Figure 5. ε - dependence of maximal velocity, v_max.

Wolfram Mathematica instructions: T by finding the zeroes of the solution for x (t); x_max - by numerically finding maxima of the solution for x (t); v_max - by numerically finding maxima of the solution for dx/dt.

A study of numerical fits for the dependence of v_max and T on ε yields that, except for small values of ε, both T and x_max grow linearly with ε. Writing

$T=T_0\epsilon ,x_{\max }=X_0\epsilon$ , (3)

For 10 ≤ ε ≤ 1000, the results of the numerical solutions presented in Figures 3-5 yield that

$T_0=1.61492,X_0=0.38513,v_{\max }=\mathrm{1.15532.}$ (4)

The standard deviations the linear fit T from the numerical data vary from 0.0006 at ε =10, to 0.00002 at ε =1000. The standard deviations the linear fit from X₀ from the numerical data vary data vary respectively, from.0001 to 0.00005. The deviations of the fit for a constant v_max are always of the order of than 0.05.

3. Scaled Variables

It pays to use scaled variables as follows:

$x\left(t\right)=x_{\max }\xi \left(\theta \right),\theta =\left(\left(2\pi t\right)/T\right)\Rightarrow \dot{x}\left(t\right)=\left(\left(2\pi x_{\max }\right)/T\right){\xi }^{\prime }\left(\theta \right)$. (5)

In addition, as ε is large, it pays to define a small parameter:

$\mu =2\pi /\left(T_0\epsilon \right)$ . (6)

Using Equations (5) and (6), Equation (1) leads to

$\xi \left(\theta \right)-\left(\left(2\pi \right)/T_0\right){\xi }^{\prime }\left(\theta \right)+\left(\left(8{\pi }^3{X}_0^2\right)/T_0^3\right){\xi }^{\prime }{\left(\theta \right)}^3+{\mu }^2{\xi }^{″}\left(\theta \right)=0$ . (7)

For t  1, x (t) tends to a limit cycle; ξ becomes periodic in θ with period 2π. Solutions of Equation (1) for different values of ε differ. Figure 6 shows phase-space plots of solutions with ε = 2000 (full line) and 850 (dashed line).

Once the transformation, Equation (5), is employed, the phase-space plots of the two scaled solutions fall on top on another, as shown in Figure 7.

The plots in Figure 6 are for very long times, where the limit cycle has been reached.

Figure 6. Phase-space plots of solutions of Equation (1); ε = 2000 (full line) and 850 (dashed line).

Figure 7. Asymptotic phase-space plot of scaled solution (Equations (5) and (6)).

Plots of ξ and dξ/dθ are shown in Figure 8 and Figure 9, respectively.

Figure 8. Plot of ξ. Full line- solution of Equation (1) scaled by Equations (3), (5) and (6); Dashed line- approximate solution (Equation (12)) discussed in Section 4.A.

Figure 9. Plot of dξ/dθ vs. θ: solution of Equation (1) scaled by Equations (3), (5) and (6).

While ξ is periodic and vanishes at θ = 0, π and 2π, the extreme values (±1) occur at points that are slightly shifted from π/2 and 3π/2 (see Figure 9). As the shifts from π/2 and 3π/2 are the same, let us focus on one of the points, θ₀ in Figure 8. The numerical solution for ε = 2000 yields

${\theta }_0=\text{1}\text{.75685}\cdots$ . (8)

The dashed line in Figure 9 is a plot of an approximate solution for ξ (θ) discussed in Section 4.A (Equation (12)).

4. Matching Approximate Solutions

A. Approximation starts around θ = 0

Consider the solution of Equation (7) near θ = 0, with ξ (0) =0. As ${\xi }^{″}$ is bounded there (it is extremely large only for θ near θ₀), in the limit μ → 0, the O (μ²) term in Equation (7) can be neglected as long as (θ₀ − θ) exceeds O (μ) (see Section 4.B), yielding

$\xi \left(\theta \right)-\left(\left(2\pi \right)/T_0\right){\xi }^{\prime }\left(\theta \right)+\left(\left(8{\pi }^3{X}_0^2\right)/T_0^3\right){\xi }^{\prime }{\left(\theta \right)}^3=0$ . (9)

We look for an approximate solution for Equation (9) in powers of θ:

$\xi \left(\theta \right)\approx {\displaystyle \sum _{k=1}^n{c}_k{\theta }^k}\equiv {\theta }_n\left(\theta \right)$ (10)

Substituting Equation (10) in Equation (9), one solves for c_k power by power in θ. Up to n = 5, one finds:

$\begin{array}{l}{c}_1=T_0/(2\pi X_0)\\ c_2=-T_0^2/\left(16{\pi }^2{X}_0\right)\\ c_3=-T_0^3/\left(96{\pi }^3{X}_0\right)\\ c_4=-13T_0^4/\left(3072{\pi }^4{X}_0\right)\\ c_5=-143T_0^5/\left(61440{\pi }^5{X}_0\right)\end{array}$ (11)

Using the values of T₀ and, X₀ given in Equation (4), through n = 5, Equation (10) becomes:

$\begin{array}{c}\xi \left(\theta \right)=0.667221\theta -0.0428289{\theta }^2-0.00366558{\theta }^3\\ -0.000764703{\theta }^4-0.000215979{\theta }^5.\end{array}$ (12)

The matching procedure is based on finding an approximation for θ₀ by requiring that, for a polynomial with a given n in Equation (10),

$\xi \left({\theta }_0\right)=1$ . (13)

The value of θ₀, computed from the numerical solution of Equation (1), for ε = 2000 (μ = 0.001947), is 1.75685 (see Equation (8)). The values computed from Equation (12) depend on the highest power, n. The deviation from 1.75685 is about 15% for polynomials ending with n =1 and 2. It becomes smaller as the number of terms in ξ (θ) (Equation (9)) is increased and is reduced to about 1% for n = 5. A plot of the n = 5 polynomial (Equation (12)) is included in Figure 9 (dashed line).

B. Approximation around θ = θ₀

In Figure 8, the approximate solution developed around θ = 0 (Equation (12), dashed line) agrees very well with the full numerical solution (full line). They grossly disagree beyond θ = θ₀ because ${\xi }^{″}$ for the polynomial approximation misses completely the violent behavior of ${\xi }^{″}\left(\theta \right)$ for the full numerical solution. This can be understood by taking into account that ${\xi }^{\prime }\left(\theta \right)\approx 0$ in the vicinity of θ = θ₀ ${\xi }^{\prime }\left(\theta \right)\approx 0$ . Neglecting ${\xi }^{\prime }\left(\theta \right)$ in Equation (7), the latter is reduced to

$\xi \left(\theta \right)+{\mu }^2{\xi }^{″}\left(\theta \right)=0$ , (14)

For

$\xi \left({\theta }_0\right)=1,{\xi }^{\prime }\left({\theta }_0\right)=0$ , (15)

Equation (13) is solved by

$\xi \left(\theta \right)=\cos \left(\left({\theta }_0-\theta \right)/\mu \right)\Rightarrow {\xi }^{″}\left(\theta \right)=-\frac{1}{{\mu }^2}\cos \left(\left({\theta }_0-\theta \right)/\mu \right)$ . (16)

${\xi }^{″}\left(\theta \right)$ is extremely large, of O (1/μ²) for (θ₀ − θ) = O (μ). It goes down rapidly once (θ₀ − θ) exceeds O (μ). The polynomial approximation of Equation (10) cannot reproduce this violent behavior. This is demonstrated in Figure 10, where the solution of Equation (1), scaled by Equations (3)-(5), the approximate solution developed around θ = 0 (Section 4.A, Equation (12)), and the approximate solution around θ = θ₀ (Equation (16)), are plotted in the vicinity of θ = θ₀.

Figure 10. Comparison of solution of Equation (1), scaled by Equations (3)-(5) (full line) and approximate solutions generated in vicinity of θ = 0 (Equation (12), dotted line)) and in vicinity of θ = θ₀ (Equation (16), dashed line).

Whereas the solution generated in the vicinity of θ = 0 (Equation (12)) is almost constant over the narrow range of θ around θ₀, the solution generated in vicinity of θ = θ₀ (Equation (16) varies dramatically in this narrow range of θ.

Appendix I Nonlinear Violence in Duffing Equation

Nonlinear Violence occurs when the nonlinear term in the equation of motion is much greater than the linear part of the driving force (x (t)) [12]. In the case of the Duffing Equation (Equation (2)), this may occur by increasing either ε or x_max. It can be analyzed through studying, E, the total oscillation energy:

$\begin{array}{l}E=1/2x^2+1/2{\dot{x}}^2+1/4\epsilon x^4\\ =1/2x_{\max }^2+1/4\epsilon x_{\max }^4=1/2{\dot{x}}^2\\ \Rightarrow x_{\max }=\sqrt{\frac{\sqrt{1+4E\epsilon }-1}{\epsilon }}\\ v_{\max }=\max \left(\dot{x}\right)=\sqrt{2E}\end{array}$ (17)

As the total energy may be chosen to have any value, increasing the total energy goes hand in hand with increasing the amplitudes of oscillation and velocity. The increase in velocity leads to a reduction in the period of oscillations. This pattern consequence is independent of the magnitude of ε.

For fixed E, as ε → ∞, x_max → 0. For fixed, ε and x_max → ∞, E → ∞. When E → ∞, v_max → ∞.

A direct calculation yields for the period:

$T=\frac{8K\left(-\frac{\epsilon x_{\max }^2}{2+\epsilon x_{\max }^2}\right)}{\sqrt{4+2\epsilon x_{\max }^2}}$ , (18)

where K (x) is the complete elliptic integral of the first kind.

T is reduced as x_max or ε are increased and tends to zero as either x_max or ε → ∞.

In summary, the effect of nonlinear violence in conservative oscillatory systems is qualitatively quite different from the effect on the Rayleigh oscillator.

Author Responsibility

1) The author made all contributions to the conception or design of the work. No data acquisition was required. No new software was created;

2) The author drafted the work or revised it critically for important intellectual content;

3) The author approved the version to be published;

4) The author agrees to be accountable for all aspects of the work in ensuring that questions related to the accuracy or integrity of any part of the work are appropriately investigated and resolved.

Data availability statement:
data sharing is not applicable to this article, as no datasets were generated or analyzed during the current study.

No data associated with this study have been deposited in a publicly available repository.

No data were used for the research described in the article.

Conflicts of Interest

The authors declare no conflicts of interest regarding the publication of this paper.

References

[1]	Rayleigh, J.W.S. (1945) The Theory of Sound. Dover Publications.
[2]	Kevorkian, J. and Cole, J.D. (1996) Multiple Scale and Singular Perturbation Methods. Springer.[CrossRef]
[3]	Bender, C.M. and Orszag, S.A. (1999) Perturbation Series. In: Bender, C.M. and Orszag, S.A. Advanced Mathematical Methods for Scientists and Engineers I, Springer 319-367.[CrossRef]
[4]	Nayfeh, A.H. (2004) Perturbation Methods. Wiley-VCH Verlag.
[5]	Bruno, A. (1971) Analytical Form of Differential Equations. Trudy Moskovskogo Matematicheskogo Obshchestva, 25, 132-288.
[6]	Arnold, V.I. (1983) Geometrical Methods in the Theory of Ordinary Differential Equations. Springer.[CrossRef]
[7]	Guckenheimer, J. and Holmes, P. (1983) Nonlinear Oscillations, Dynamical Systems, and Bifurcations of Vector Fields. Springer.[CrossRef]
[8]	Kodama, Y. and Mikhailov, A.V. (1997) Obstacles to Asymptotic Integrability. In: Fokas, A.S. and Gelfand, I.M., Eds., Algebraic Aspects of Integrable Systems, Birkhäuser, 173-204.[CrossRef]
[9]	Georgescu, A., Bazavan, P. and Sterpu, M. (2008) Approximate Limit Cycles for the Rayleigh Model. ROMAI J, 4, 73-780.
[10]	del Pino, J., Košata, J. and Zilberberg, O. (2024) Limit Cycles as Stationary States of an Extended Harmonic Balance Ansatz. Physical Review Research, 6, Article 033180.[CrossRef]
[11]	Zarmi, Y. (2023) On the Relativistic Harmonic Oscillator. Applied Mathematics, 14, 1-20.[CrossRef]
[12]	Zarmi, Y. (2024) Nonlinear Violence in Nonlinear Oscillations at High Energy. Applied Mathematics, 15, 65-95.[CrossRef]