- Email: [email protected]
Solution of the relativistic (an)harmonic oscillator using the harmonic balance method
ARTICLE IN PRESS JOURNAL OF SOUND AND VIBRATION Journal of Sound and Vibration 311 (2008) 1447–1456 www.elsevier.com/locate/jsvi
Short Communication
Solution of the relativistic (an)harmonic oscillator using the harmonic balance method A. Bele´ndez, C. Pascual, D.I. Me´ndez, C. Neipp Departamento de Fı´sica, Ingenierı´a de Sistemas y Teorı´a de la Sen˜al, Universidad de Alicante, Apartado 99, E-03080 Alicante, Spain Received 17 May 2007; received in revised form 3 October 2007; accepted 8 October 2007 Available online 7 November 2007
Abstract The harmonic balance method is used to construct approximate frequency–amplitude relations and periodic solutions to the relativistic oscillator. By combining linearization of the governing equation with the harmonic balance method, we construct analytical approximations to the oscillation frequencies and periodic solutions for the oscillator. To solve the nonlinear differential equation, firstly we make a change of variable and secondly the differential equation is rewritten in a form that does not contain the square-root expression. The approximate frequencies obtained are valid for the complete range of oscillation amplitudes, A, while the discrepancy between the second approximate frequency and the exact one never exceeds 0.82% and tends to 0.52% when A tends to infinity. Excellent agreement of the approximate frequencies and periodic solutions with the exact ones are demonstrated and discussed. r 2007 Elsevier Ltd. All rights reserved.
The simple harmonic oscillator and its extension to the relativistic case are important in physics because they are usually employed as the basis for analyzing more complicated motion. When the energy of a simple harmonic oscillator is such that the velocities become relativistic, the simple harmonic motion (linear oscillations) at low energy becomes anharmonic (nonlinear oscillations) at high energy [1]; hence the parentheses around the ‘‘an’’ in the title of this paper. Then, the strength of the nonlinearity increases as the total relativistic energy increases, and at the non-relativistic limit the oscillator becomes linear. Synge [2] gave an exact expression for the period in terms of an integral that Goldstein [3] identified as being expressible in terms of elliptical integrals. Mickens [4] showed that all the solutions to the relativistic (an)harmonic oscillator are periodic and determined a method for calculating analytical approximations to its solutions. Mickens considered the first-order harmonic balance method, but we think he did not apply the technique correctly and the analytical approximate frequency he obtained is not the correct one. The purpose of this paper is to determine the high-order periodic solutions to the relativistic oscillator by applying the harmonic balance method [5–14]. To do this, we use the analytical approach developed by Lim and Wu [6,7] in which linearization of the governing equation is combined with the harmonic balance method.
Corresponding author. Tel.: +34 96 5903651; fax: +34 96 5903464.
E-mail address: [email protected] (A. Bele´ndez). 0022-460X/$ - see front matter r 2007 Elsevier Ltd. All rights reserved. doi:10.1016/j.jsv.2007.10.010
ARTICLE IN PRESS 1448
A. Bele´ndez et al. / Journal of Sound and Vibration 311 (2008) 1447–1456
The governing non-dimensional nonlinear differential equation of motion for the relativistic oscillator is [4] " 2 #3=2 d2 x dx þ 1 x ¼ 0, (1) dt2 dt where x and t are dimensionless variables. The even power term in Eq. (1), (dx/dt)2, acts like the powers of coordinates in that it does not cause a damping of the amplitude of oscillations with time. Therefore, Eq. (1) is an example of a generalized conservative system [5]. At the limit when (dx/dt)251 Eq. (1) becomes (d2x/dt2)+xE0 the oscillator is linear and the proper time t becomes equivalent to the coordinate time t to this order. Introducing the phase space variable (x,y), Eq. (1) can be written in the system form dx dy ¼ y; ¼ ð1 y2 Þ3=2 x (2) dt dt and the trajectories in phase space are given by solutions to the first order, ordinary differential equation dy ð1 y2 Þ3=2 x ¼ . dx y
(3)
As Mickens pointed out, since the physical solutions of both Eqs. (1) and (3) are real, the phase space has a ‘‘strip’’ structure [4], i.e., 1oxo þ 1;
1oyo þ 1.
(4)
Then unlike the usual non-relativistic harmonic oscillator, the relativistic oscillator is bounded in the y variable. This is due to the fact that the dimensionless variable y is related with the relativistic parameter b ¼ v/c, where v is the velocity of the particle and c the velocity of light. In the relativistic case, the condition cov(t)o+c must be met, and so we obtain 1oy(t)o+1. Mickens proved that all the trajectories to Eq. (3) are closed in the open region of phase space given by Eq. (4) and then all the physical solutions to Eq. (1) are periodic. However, unlike the usual (non-relativistic) harmonic oscillator, the relativistic (an)harmonic oscillator contains higher-order multiples of the fundamental frequency [4]. The harmonic balance method can now be applied to obtain analytic approximations to the periodic solutions of Eq. (1). We make a change of variable, y-u, such that 1ouo þ 1.
(5)
u y ¼ pffiffiffiffiffiffiffiffiffiffiffiffiffi 1 þ u2
(6)
The required transformation is [4]
and the corresponding second-order nonlinear differential equation for u is d2 u u þ pffiffiffiffiffiffiffiffiffiffiffiffiffi ¼ 0. dt2 1 þ u2
(7)
We consider the following initial conditions in Eq. (7): du ð0Þ ¼ 0. (8) dt Eq. (7) is an example of a conservative nonlinear oscillatory system in which the restoring force has an irrational form. All the motions corresponding to Eq. (7) are periodic and the system will oscillate within symmetric bounds [B, B]. We can solve Eq. (7) approximately using the harmonic balance method. To do this, we first write this equation in a form that does not contain the square-root expression 2 2 d u ð1 þ u2 Þ ¼ u2 . (9) dt2 uð0Þ ¼ B
and
ARTICLE IN PRESS A. Bele´ndez et al. / Journal of Sound and Vibration 311 (2008) 1447–1456
1449
Since the restoring force is an odd function of u, the periodic solution u(t) has a Fourier series representation which contains only odd multiples of ot. The first-order harmonic balance solution takes the form [5] u1 ðtÞ ¼ B cos o t.
(10)
Observe that u1(t) satisfies the initial conditions, Eq. (8). The angular frequency of the oscillator is o, which is unknown and to be determined. Both the periodic solution u(t) and frequency o (thus period T ¼ 2p/o) depend on B. Substituting Eq. (10) into Eq. (9) gives ð1 þ B2 cos2 otÞo4 B2 cos2 ot ¼ B2 cos2 ot. Expanding and simplifying the above expression gives 1 4 2 1 4 1 1 4 3 4 2 1 1 4 2 2o þ 8o B 2 þ 2o B þ 2o 2 cos 2ot þ 8o B cos 4ot ¼ 0
(11)
(12)
and setting the coefficient of the resulting term cos(0ot) (the lowest harmonic) equal to zero gives the first analytical approximate frequency o1 as a function of B: 1=4 . (13) o1 ðBÞ ¼ 1 þ 34B2 In this limit, the present method gives exactly the same frequency [15] as the first-order homotopy perturbation method [16,17] applied to Eq. (9). The corresponding first analytical approximate periodic solution is given by u1 ðtÞ ¼ B cos½o1 ðBÞt.
(14)
The approximate frequency in Eq. (13) is the correct one when the harmonic balance method is applied to Eq. (11) and not the following frequency obtained by Mickens [4]: 1=4 oM ðBÞ ¼ 1 þ 12B2 . (15) This is due to the fact that Mickens divided Eq. (11) by cos2 ot and obtained ð1 þ B2 cos2 otÞo4 B2 ¼ B2 .
(16)
It is easy to verify that Eq. (15) can be obtained from Eq. (16). However, oM is not the first-order approximate frequency of Eq. (9). We can readily see that oM is the first-order approximate frequency of the following differential equation: 2 2 2 2 du ¼ u2 , (17) 1þ u 3 dt2 which can be obtained from d2 u u þ pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi ¼ 0. 2 dt 1 þ ð2=3Þu2
(18)
As we can see, Eq. (18) does not coincide with Eq. (7). The corresponding approximation to y is obtained from Eq. (6) u1 ðtÞ B cos o1 t y1 ðtÞ ¼ qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi ¼ pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi . 1 þ B2 cos2 o1 t 1 þ u21 ðtÞ
(19)
Likewise, x1(t) can be calculated by integrating equation y ¼ dx/dt subject to the restrictions x1 ð0Þ ¼ 0;
B y1 ð0Þ ¼ pffiffiffiffiffiffiffiffiffiffiffiffiffiffi , 1 þ B2
which can be easily obtained from Eqs. (6) and (7). This integration gives " # 3 2 1=4 1 B pffiffiffiffiffiffiffiffiffiffiffiffiffiffi sin½o1 ðBÞt . sin x1 ðtÞ ¼ 1 þ B 4 1 þ B2
(20)
(21)
ARTICLE IN PRESS 1450
A. Bele´ndez et al. / Journal of Sound and Vibration 311 (2008) 1447–1456
In order to apply the next level of harmonic balance method, firstly we express the periodic solution to Eq. (7) with the assigned conditions in Eq. (8) in the form of [6,7] u2 ðtÞ ¼ u1 ðtÞ þ Du1 ðtÞ,
(22)
where Du1(t) is the correction part. Linearizing the governing Equations (7) and (8) with respect to the correction Du1(t) at u(t) ¼ u1(t) leads to 2 2 2 2 2 2 2 2 d u1 d u1 d2 u1 d2 Du1 2 d u1 2 d u1 d Du1 u21 2u1 Du1 þ þ u þ 2u Du þ 2 þ 2u ¼0 (23) 1 1 1 1 dt2 dt2 dt2 dt2 dt2 dt2 dt2 and dDu1 ð0Þ ¼ 0. (24) dt The approximation Du1(t) in Eq. (22), which satisfies the initial conditions in Eq. (24), takes the form [6,7] Du1 ð0Þ ¼ 0;
Du1 ¼ c1 ðcos ot cos 3otÞ,
(25)
where c1 is a constant to be determined. Substituting Eqs. (10), (22) and (25) into Eq. (23), expanding the expression in a trigonometric series, and setting the coefficients of the resulting items cos(0ot) and cos(2ot) equal to zero, respectively, yield 12B þ 12Bo4 þ 38B3 o4 þ o4 c1 B2 o4 c1 c1 ¼ 0
(26)
2 4 12B þ 12Bo4 þ 12B3 o4 8o4 c1 11 2 B o c1 ¼ 0.
(27)
and
From Eq. (26) we can obtain c1 as follows: c1 ¼
4B þ 4Bo4 þ 3B3 o4 . 8ð1 o4 þ B2 o4 Þ
(28)
Substituting Eq. (28) into Eq. (27) and solving for the second analytical approximate frequency o2, we obtain pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi !1=4 40 þ 22B2 þ 2 256 þ 256B2 þ 71B4 o2 ðBÞ ¼ . (29) 72 þ 92B2 þ 25B4 Furthermore, c1 can be obtained by substituting Eq. (29) into Eq. (28) and the result is c1 ðBÞ ¼
B½64 þ 17B4 4D þ ð80 3DÞB2 , 4½32 þ 47B4 2D þ 2ð55 þ DÞB2
(30)
where DðBÞ ¼
pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 256 þ 256B2 þ 71B4 .
(31)
The corresponding second analytical approximate periodic solution is given by u2 ðtÞ ¼ ½B þ c1 ðBÞ cos½o2 ðBÞt c1 cos½3o2 ðBÞt.
(32)
The corresponding approximation to y is u2 ðtÞ B cos o2 t þ c1 ðcos o2 t cos 3o2 tÞ y2 ðtÞ ¼ qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi ¼ pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi . 2 1 þ ½B cos o2 t þ c1 ðcos o2 t cos 3o2 tÞ2 1 þ u ðtÞ
(33)
2
However, the analytical integration of Eq. (33) to obtain x2(t) is not possible. To obtain an analytical expression for x2(t), Eq. (25) is written as follows: B cos o2 t þ 4c1 ðcos o2 t cos3 o2 tÞ Bðcos o2 t þ ð4c1 =BÞ cos o2 t sin2 o2 tÞ y2 ðtÞ ¼ pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi ¼ qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi . 1 þ ½B cos o2 t þ 4c1 ðcos o2 t cos3 o2 tÞ2 1 þ B2 ðcos o2 t þ ð4c1 =BÞ cos o2 tsin2 o2 tÞ2
(34)
ARTICLE IN PRESS A. Bele´ndez et al. / Journal of Sound and Vibration 311 (2008) 1447–1456
1451
The above equation would be easily integrable if (4c1/B) cos o2t sin2o2t5cos o2t. Then Eq. (34) would have the same functional form as Eq. (19). From Eqs. (30) and (31) we can obtain the following limits: 4c1 ðBÞ ¼ 0, B!0 B pffiffiffiffiffi 4c1 ðBÞ 2337 þ 4283 71 pffiffiffiffiffi ¼ 0:12965. ¼ lim B!1 B 143867 þ 13822 71 lim
(35)
(36)
Then, 4c1(B)/B takes values between 0 and 0.12965 when B varies between 0 and N. For example, for B ¼ 1 and 10, 4c1(B)/B takes the values 0.0424 and 0.1271, respectively. We can write 4c1 2 2 cos o tsin o t (37) 2 2 p 0:12965 cos o2 tsin o2 t p0:05, B where we have taken into account that the maximum value of |cos o2t sin2 o2t| is 0.3849. In Fig. 1 we have plotted cos (2ph) and cos (2ph)+(4c1/B) cos(2ph) sin2(2ph) as a function of h ¼ o2t/2p ¼ t/T2 for B-N (4c1/B ¼ 0.12965). From Eqs. (34)–(37) and Fig. 1 we can conclude that Eq. (33) can be approximately written as follows: B cos o2 t y2 ðtÞ pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi . 1 þ B2 cos2 o2 t
(38)
Later we will verify that such a simple approximation gives very good results for x(t). Likewise, x2(t) can be calculated by integrating equation y ¼ dx/dt subject to the restrictions x2 ð0Þ ¼ 0;
B y2 ð0Þ ¼ pffiffiffiffiffiffiffiffiffiffiffiffiffiffi 1 þ B2
(39)
and this integration gives " # 3 2 1=4 1 B pffiffiffiffiffiffiffiffiffiffiffiffiffiffi sin½o2 ðBÞt . sin x2 ðtÞ ¼ 1 þ B 4 1 þ B2
(40)
We will show Eq. (40) gives good results for x2(t). However, we should not forget what we are really looking for is an approximate analytical solution to Eq. (1), that is, x(t). Moreover, it is convenient to express the approximate angular frequency and the solution in terms of oscillation amplitude A rather than as a function of B. It is now necessary to find a relation between oscillation amplitude A and parameter B used to solve Eq. (7) approximately.
Fig. 1. cos(2ph) (continuous line) and cos(2ph)+0.12965 cos(2ph)sin2(2ph) (dashed line) as a function of h ¼ o2t/2p ¼ t/T2.
ARTICLE IN PRESS 1452
A. Bele´ndez et al. / Journal of Sound and Vibration 311 (2008) 1447–1456
From Eq. (3) we get 1 ð1
y2 Þ1=2
1 þ x2 ¼ C, 2
(41)
where C is a constant to be determined as a function of initial conditions. From Eq. (20) we can easily obtain C ¼ (1+B2)1/2 and Eq. (41) can be written as follows: 1 ð1
y2 Þ1=2
1 þ x2 ¼ ð1 þ B2 Þ1=2 . 2
(42)
In addition, when x ¼ A, the velocity y ¼ dx/dt is zero. Taking this into account in Eq. (49), we obtain the following relation between amplitude A and parameter B: 1 þ 12A2 ¼ ð1 þ B2 Þ1=2 . From the above equation we can easily find that the solution for B is 1=2 B ¼ A 1 þ 14A2 .
(43)
(44)
Substituting Eq. (44) into Eqs. (13) and (21), the first analytical approximate periodic solution for the relativistic oscillator as a function of the oscillation amplitude A is given by 3 2 1=4 o1 ðAÞ ¼ 1 þ 34A2 þ 16 A , (45) 2sffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 3 2 4 1 4A þ A x1 ðtÞ ¼ sin1 4 sin½o1 ðAÞt5. o1 ðAÞ 4 þ 4A2 þ A4
(46)
Substituting Eq. (44) into Eqs. (29) and (40), the second analytical approximate periodic solution for the relativistic oscillator as a function of the oscillation amplitude A is given by pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi !1=4 640 þ 352A2 þ 88A4 þ 8 4096 þ 4096A2 þ 2160A4 þ 568A6 þ 71A8 o2 ðAÞ ¼ , (47) 1152 þ 1472A2 þ 768A4 þ 200A6 þ 25A8 2sffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 3 2 4 1 4A þ A sin1 4 sin½o2 ðAÞt5. x2 ðtÞ ¼ o2 ðAÞ 4 þ 4A2 þ A4
(48)
Now we illustrate the applicability, accuracy and effectiveness of the proposed approach by comparing the approximate analytical periodic solutions obtained in this paper with the exact ones. Calculation of the exact period frequency, Tex(A), proceeds as follows. Substituting Eq. (44) into Eq. (42) and integrating, we obtain Z A 1 þ 12ðA2 x2 Þ ffi dx, qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi T ex ðAÞ ¼ 4 (49) 0 A2 x2 þ 14ðA2 x2 Þ2 which allows us to obtain the exact angular frequency, oex(A), in terms of elliptical integrals as follows: " #1 pffiffiffiffiffiffiffiffiffiffiffiffiffiffi A2 2p 8 A2 2 ¼ 2p 4 4 þ A E , (50) oex ðAÞ ¼ pffiffiffiffiffiffiffiffiffiffiffiffiffiffi K T ex ðAÞ 4 þ A2 4 þ A2 4 þ A2 where K(q) and E(q) are the complete elliptic integrals of the first and second kind, respectively, defined as follows [18]: Z 1 Z 1 rffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi2 dz 1 qz pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi ; EðqÞ ¼ KðqÞ ¼ dz. (51) 2 2 1 z2 ð1 z Þð1 qz Þ 0 0
ARTICLE IN PRESS A. Bele´ndez et al. / Journal of Sound and Vibration 311 (2008) 1447–1456
1453
For small values of the amplitude A it is possible to take into account the following power series expansion of the angular frequencies oex ðAÞ 1
3 2 51 4 A þ A , 16 1024
(52)
o1 ðAÞ 1
3 2 42 4 A þ A , 16 1024
(53)
3 2 56 4 A þ A . (54) 16 1024 These series expansions were carried out using MATHEMATICA. As can be seen, in the expansions of the angular frequencies o1(A) and o2(A), the first two terms are the same as the first two terms of the equation obtained in the power-series expansion of the exact angular frequency, oex(A). If we compare the third term in Eqs. (53) and (54) with the third term in the series expansion of the exact frequency oex(A) (Eq. (52)), we can see that the third term in the series expansions of o2(A) is more accurate than the third term in the expansion of o1(A). For very large values of the amplitude A it is possible to take into account the following power series expansion of the exact angular frequency: o2 ðAÞ 1
p 1:5708 þ ¼ þ , 2A A
oex ðAÞ o1 ðAÞ
2 3
1=4
A
þ ¼
1:5197 þ , A
pffiffiffiffiffi1=4 88 þ 8 71 1 1:5790 o2 ðAÞ þ ¼ þ . 25 A A
(55) (56)
(57)
Once again we can see than o2(A) provides excellent approximations to the exact frequency oex(A) for very large values of oscillation amplitude. Furthermore, we have the following equations: lim oex ðAÞ ¼ lim o1 ðAÞ ¼ lim o2 ðAÞ ¼ 1,
(58)
lim oex ðAÞ ¼ lim o1 ðAÞ ¼ lim o2 ðAÞ ¼ 0,
(59)
A!0
A!1
A!0
A!1
lim
o1 ðAÞ
A!0
A!1
¼ 0:96745,
(60)
o2 ðAÞ ¼ 1:00523. A!1 oex ðAÞ
(61)
A!1 oex ðAÞ
lim
Eqs. (58)–(61) illustrate the very good agreement of the approximate frequency o2(A) with the exact frequency oex(A) for small as well as large values of oscillation amplitude. In Fig. 2 we plotted the percentage error of approximate frequencies o1 and o2, as a function of A. In this figure the percentage errors were calculated using the following equation: oj oex ; j ¼ 1; 2. (62) Relative error of oj ð%Þ ¼ 100 oex As we can see from Fig. 2, the relative errors for o1(A) are lower than 0.82% for all the range of values of amplitude of oscillation A (see Eq. (60)), and these relative errors tend to 0.52% when A tends to infinity (see Eq. (61)). The exact periodic solutions x(t) achieved by numerically integrating Eq. (1), and the proposed first-order approximate periodic solutions x1(t) in Eq. (46) and x2(t) in Eq. (48) are plotted in Figs. 3–6 for A ¼ 0.1, 1, 2 and 10 (B ¼ 0.1001, 1.12, 2.83 and 51), respectively. In these figures parameter h is
ARTICLE IN PRESS 1454
A. Bele´ndez et al. / Journal of Sound and Vibration 311 (2008) 1447–1456
Fig. 2. Relative error for approximate frequencies o1, (Eq. (45), dashed line) and o2 (Eq. (47), continuous line).
Fig. 3. Comparison of the normalized approximate analytical solutions x1/A (m) and x2/A (J) with the exact solution (continuous line) for A ¼ 0.1 (b0 ¼ v0/c ¼ 0.09963).
Fig. 4. Comparison of the normalized approximate analytical solutions x1/A(m) and x2/A (J) with the exact solution (continuous line) for A ¼ 1 (b0 ¼ v0/c ¼ 0.74536).
ARTICLE IN PRESS A. Bele´ndez et al. / Journal of Sound and Vibration 311 (2008) 1447–1456
1455
Fig. 5. Comparison of the normalized approximate analytical solutions x1/A (m) and x2/A (J) with the exact solution (continuous line) A ¼ 2 (b0 ¼ v0/c ¼ 0.94281).
Fig. 6. Comparison of the normalized approximate analytical solutions x1/A (m) and x2/A(J) with the exact solution (continuous line) A ¼ 10 (b0 ¼ v0/c ¼ 0.99981).
defined as follows: h¼
t ¼ 2poex ðAÞt. T ex ðAÞ
(63)
From Eqs. (20) and (44) we have v0 B b0 ¼ ¼ yð0Þ ¼ pffiffiffiffiffiffiffiffiffiffiffiffiffiffi ¼ c 1 þ B2
sffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 4A2 þ A4 , 4 þ 4A2 þ A4
(64)
where v0 is the maximum velocity of the particle and c is the velocity of light. Figs. 3–6 show that Eq. (48) provides a good approximation to the exact periodic solutions and that the approximation considered in Eq. (38) is adequate to obtain the approximate analytical expression of x2(t) given in Eq. (48). As we can see, for small values of A (Figs. 3 and 4) x(t) is very close to the sine function form of non-relativistic simple harmonic motion. For higher values of A the curvature becomes more concentrated at the turning points (x ¼ 7A). For these values of A, x(t) becomes markedly anharmonic and is almost straight between the turning points. Only in the vicinity of the turning points, where the magnitude of the
ARTICLE IN PRESS A. Bele´ndez et al. / Journal of Sound and Vibration 311 (2008) 1447–1456
1456
Hooke’s law force is maximum and the velocity becomes relativistic, is the force effective in changing the velocity [1]. Fig. 6 are a typical example of the motion in the ultra-relativistic region where b0-1. In summary, the harmonic balance method was used to obtain two approximate frequencies for the relativistic oscillator. To do this we rewrite the nonlinear differential equation in a form that does not contain an irrational expression. We can conclude that the approximate frequencies obtained are valid for the complete range of oscillation amplitude, including the limiting cases of amplitude approaching zero and infinity. Excellent agreement of the approximate frequencies with the exact one was demonstrated and discussed and the discrepancy between the second approximate frequency, o2, and the exact one never exceeds 0.82%. Some examples have been presented to illustrate the excellent accuracy of the approximate analytical solutions. Finally, we can see that the method considered here is very simple in its principle, and is very easy to apply. Acknowledgments This work was supported by the ‘‘Ministerio de Educacio´n y Ciencia’’, Spain, under Project FIS2005-05881C02-02, and by the ‘‘Generalitat Valenciana’’, Spain, under Project ACOMP-2007-020. We would like to thank the anonymous reviewers for their comments and suggestions, which have undoubtedly improved the original manuscript. References [1] [2] [3] [4] [5] [6] [7] [8] [9] [10] [11] [12] [13] [14] [15]
[16] [17] [18]
W. Moreau, R. Easther, R. Neutze, Relativistic (an)harmonic oscillator, American Journal of Physics 62 (1994) 531–535. J.L. Synge, in: S. Flu¨gge (Ed.), Classical Dynamics in Encyclopedia of Physics, Vol. III/1, Springer, Berlin, 1960, pp. 210–211. H. Goldstein, Classical Mechanics, second ed, Addison-Wesley, Reading, MA, 1981. R.E. Mickens, Periodic solutions of the relativistic harmonic oscillator, Journal of Sound and Vibration 212 (1998) 905–908. R.E. Mickens, Oscillations in Planar Dynamics Systems, World Scientific, Singapore, 1996. C.W. Lim, B.S. Wu, A new analytical approach to the Duffing-harmonic oscillator, Physics Letters A 311 (2003) 365–373. B.S. Wu, C.W. Lim, Large amplitude nonlinear oscillations of a general conservative system, International Journal of Non-linear Mechanics 39 (2004) 859–870. A. Bele´ndez, A. Herna´ndez, A. Ma´rquez, T. Bele´ndez, C. Neipp, Analytical approximations for the period of a simple pendulum, European Journal of Physics 27 (2006) 539–551. H. Hu, J.H. Tang, Solution of a Duffing-harmonic oscillator by the method of harmonic balance, Journal of Sound and Vibration 294 (2006) 637–639. H. Hu, Solution of a quadratic nonlinear oscillator by the method of harmonic balance, Journal of Sound and Vibration 293 (2006) 462–468. H.P.W. Gottlieb, Harmonic balance approach to limit cycles for nonlinear jerk equations, Journal of Sound and Vibration 297 (2006) 243–250. A. Bele´ndez, A. Herna´ndez, T. Bele´ndez, M.L. A´lvarez, S. Gallego, M. Ortun˜o, C. Neipp, Application of the harmonic balance method to a nonlinear oscillator typified by a mass attached to a stretched wire, Journal of Sound and Vibration 302 (2007) 1018–1029. G.R. Itovich, J.L. Moiola, On period doubling bifurcations of cycles and the harmonic balance method, Chaos, Solitons and Fractals 27 (2005) 647–665. Y.M. Chen, J.K. Lui, A new method based on the harmonic balance method for nonlinear oscillators, Physics Letters A 368 (2007) 371–378. A. Bele´ndez, C. Pascual, A. Ma´rquez, D.I. Me´ndez, Application of He’s homotopy perturbation method to the relativistic (an)harmonic oscillator. I: comparison between approximate and exact frequencies, International Journal of Non-linear Science and Numerical Simulations 8 (4) (2007) 483–491. J.H. He, Some asymptotic methods for strongly nonlinear equations, International Journal of Modern Physics B 20 (2006) 1141–1199. A. Bele´ndez, A. Herna´ndez, T. Bele´ndez, A. Ma´rquez, Application of the homotopy perturbation method to the nonlinear pendulum, European Journal of Physics 28 (2007) 93–104. L.M. Milne-Thomson, Elliptic integrals, in: M. Abramowitz, I.A. Stegun (Eds.), Handbook of Mathematical Functions, Dover Publications, Inc., New York, 1972.
Short Communication
Solution of the relativistic (an)harmonic oscillator using the harmonic balance method A. Bele´ndez, C. Pascual, D.I. Me´ndez, C. Neipp Departamento de Fı´sica, Ingenierı´a de Sistemas y Teorı´a de la Sen˜al, Universidad de Alicante, Apartado 99, E-03080 Alicante, Spain Received 17 May 2007; received in revised form 3 October 2007; accepted 8 October 2007 Available online 7 November 2007
Abstract The harmonic balance method is used to construct approximate frequency–amplitude relations and periodic solutions to the relativistic oscillator. By combining linearization of the governing equation with the harmonic balance method, we construct analytical approximations to the oscillation frequencies and periodic solutions for the oscillator. To solve the nonlinear differential equation, firstly we make a change of variable and secondly the differential equation is rewritten in a form that does not contain the square-root expression. The approximate frequencies obtained are valid for the complete range of oscillation amplitudes, A, while the discrepancy between the second approximate frequency and the exact one never exceeds 0.82% and tends to 0.52% when A tends to infinity. Excellent agreement of the approximate frequencies and periodic solutions with the exact ones are demonstrated and discussed. r 2007 Elsevier Ltd. All rights reserved.
The simple harmonic oscillator and its extension to the relativistic case are important in physics because they are usually employed as the basis for analyzing more complicated motion. When the energy of a simple harmonic oscillator is such that the velocities become relativistic, the simple harmonic motion (linear oscillations) at low energy becomes anharmonic (nonlinear oscillations) at high energy [1]; hence the parentheses around the ‘‘an’’ in the title of this paper. Then, the strength of the nonlinearity increases as the total relativistic energy increases, and at the non-relativistic limit the oscillator becomes linear. Synge [2] gave an exact expression for the period in terms of an integral that Goldstein [3] identified as being expressible in terms of elliptical integrals. Mickens [4] showed that all the solutions to the relativistic (an)harmonic oscillator are periodic and determined a method for calculating analytical approximations to its solutions. Mickens considered the first-order harmonic balance method, but we think he did not apply the technique correctly and the analytical approximate frequency he obtained is not the correct one. The purpose of this paper is to determine the high-order periodic solutions to the relativistic oscillator by applying the harmonic balance method [5–14]. To do this, we use the analytical approach developed by Lim and Wu [6,7] in which linearization of the governing equation is combined with the harmonic balance method.
Corresponding author. Tel.: +34 96 5903651; fax: +34 96 5903464.
E-mail address: [email protected] (A. Bele´ndez). 0022-460X/$ - see front matter r 2007 Elsevier Ltd. All rights reserved. doi:10.1016/j.jsv.2007.10.010
ARTICLE IN PRESS 1448
A. Bele´ndez et al. / Journal of Sound and Vibration 311 (2008) 1447–1456
The governing non-dimensional nonlinear differential equation of motion for the relativistic oscillator is [4] " 2 #3=2 d2 x dx þ 1 x ¼ 0, (1) dt2 dt where x and t are dimensionless variables. The even power term in Eq. (1), (dx/dt)2, acts like the powers of coordinates in that it does not cause a damping of the amplitude of oscillations with time. Therefore, Eq. (1) is an example of a generalized conservative system [5]. At the limit when (dx/dt)251 Eq. (1) becomes (d2x/dt2)+xE0 the oscillator is linear and the proper time t becomes equivalent to the coordinate time t to this order. Introducing the phase space variable (x,y), Eq. (1) can be written in the system form dx dy ¼ y; ¼ ð1 y2 Þ3=2 x (2) dt dt and the trajectories in phase space are given by solutions to the first order, ordinary differential equation dy ð1 y2 Þ3=2 x ¼ . dx y
(3)
As Mickens pointed out, since the physical solutions of both Eqs. (1) and (3) are real, the phase space has a ‘‘strip’’ structure [4], i.e., 1oxo þ 1;
1oyo þ 1.
(4)
Then unlike the usual non-relativistic harmonic oscillator, the relativistic oscillator is bounded in the y variable. This is due to the fact that the dimensionless variable y is related with the relativistic parameter b ¼ v/c, where v is the velocity of the particle and c the velocity of light. In the relativistic case, the condition cov(t)o+c must be met, and so we obtain 1oy(t)o+1. Mickens proved that all the trajectories to Eq. (3) are closed in the open region of phase space given by Eq. (4) and then all the physical solutions to Eq. (1) are periodic. However, unlike the usual (non-relativistic) harmonic oscillator, the relativistic (an)harmonic oscillator contains higher-order multiples of the fundamental frequency [4]. The harmonic balance method can now be applied to obtain analytic approximations to the periodic solutions of Eq. (1). We make a change of variable, y-u, such that 1ouo þ 1.
(5)
u y ¼ pffiffiffiffiffiffiffiffiffiffiffiffiffi 1 þ u2
(6)
The required transformation is [4]
and the corresponding second-order nonlinear differential equation for u is d2 u u þ pffiffiffiffiffiffiffiffiffiffiffiffiffi ¼ 0. dt2 1 þ u2
(7)
We consider the following initial conditions in Eq. (7): du ð0Þ ¼ 0. (8) dt Eq. (7) is an example of a conservative nonlinear oscillatory system in which the restoring force has an irrational form. All the motions corresponding to Eq. (7) are periodic and the system will oscillate within symmetric bounds [B, B]. We can solve Eq. (7) approximately using the harmonic balance method. To do this, we first write this equation in a form that does not contain the square-root expression 2 2 d u ð1 þ u2 Þ ¼ u2 . (9) dt2 uð0Þ ¼ B
and
ARTICLE IN PRESS A. Bele´ndez et al. / Journal of Sound and Vibration 311 (2008) 1447–1456
1449
Since the restoring force is an odd function of u, the periodic solution u(t) has a Fourier series representation which contains only odd multiples of ot. The first-order harmonic balance solution takes the form [5] u1 ðtÞ ¼ B cos o t.
(10)
Observe that u1(t) satisfies the initial conditions, Eq. (8). The angular frequency of the oscillator is o, which is unknown and to be determined. Both the periodic solution u(t) and frequency o (thus period T ¼ 2p/o) depend on B. Substituting Eq. (10) into Eq. (9) gives ð1 þ B2 cos2 otÞo4 B2 cos2 ot ¼ B2 cos2 ot. Expanding and simplifying the above expression gives 1 4 2 1 4 1 1 4 3 4 2 1 1 4 2 2o þ 8o B 2 þ 2o B þ 2o 2 cos 2ot þ 8o B cos 4ot ¼ 0
(11)
(12)
and setting the coefficient of the resulting term cos(0ot) (the lowest harmonic) equal to zero gives the first analytical approximate frequency o1 as a function of B: 1=4 . (13) o1 ðBÞ ¼ 1 þ 34B2 In this limit, the present method gives exactly the same frequency [15] as the first-order homotopy perturbation method [16,17] applied to Eq. (9). The corresponding first analytical approximate periodic solution is given by u1 ðtÞ ¼ B cos½o1 ðBÞt.
(14)
The approximate frequency in Eq. (13) is the correct one when the harmonic balance method is applied to Eq. (11) and not the following frequency obtained by Mickens [4]: 1=4 oM ðBÞ ¼ 1 þ 12B2 . (15) This is due to the fact that Mickens divided Eq. (11) by cos2 ot and obtained ð1 þ B2 cos2 otÞo4 B2 ¼ B2 .
(16)
It is easy to verify that Eq. (15) can be obtained from Eq. (16). However, oM is not the first-order approximate frequency of Eq. (9). We can readily see that oM is the first-order approximate frequency of the following differential equation: 2 2 2 2 du ¼ u2 , (17) 1þ u 3 dt2 which can be obtained from d2 u u þ pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi ¼ 0. 2 dt 1 þ ð2=3Þu2
(18)
As we can see, Eq. (18) does not coincide with Eq. (7). The corresponding approximation to y is obtained from Eq. (6) u1 ðtÞ B cos o1 t y1 ðtÞ ¼ qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi ¼ pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi . 1 þ B2 cos2 o1 t 1 þ u21 ðtÞ
(19)
Likewise, x1(t) can be calculated by integrating equation y ¼ dx/dt subject to the restrictions x1 ð0Þ ¼ 0;
B y1 ð0Þ ¼ pffiffiffiffiffiffiffiffiffiffiffiffiffiffi , 1 þ B2
which can be easily obtained from Eqs. (6) and (7). This integration gives " # 3 2 1=4 1 B pffiffiffiffiffiffiffiffiffiffiffiffiffiffi sin½o1 ðBÞt . sin x1 ðtÞ ¼ 1 þ B 4 1 þ B2
(20)
(21)
ARTICLE IN PRESS 1450
A. Bele´ndez et al. / Journal of Sound and Vibration 311 (2008) 1447–1456
In order to apply the next level of harmonic balance method, firstly we express the periodic solution to Eq. (7) with the assigned conditions in Eq. (8) in the form of [6,7] u2 ðtÞ ¼ u1 ðtÞ þ Du1 ðtÞ,
(22)
where Du1(t) is the correction part. Linearizing the governing Equations (7) and (8) with respect to the correction Du1(t) at u(t) ¼ u1(t) leads to 2 2 2 2 2 2 2 2 d u1 d u1 d2 u1 d2 Du1 2 d u1 2 d u1 d Du1 u21 2u1 Du1 þ þ u þ 2u Du þ 2 þ 2u ¼0 (23) 1 1 1 1 dt2 dt2 dt2 dt2 dt2 dt2 dt2 and dDu1 ð0Þ ¼ 0. (24) dt The approximation Du1(t) in Eq. (22), which satisfies the initial conditions in Eq. (24), takes the form [6,7] Du1 ð0Þ ¼ 0;
Du1 ¼ c1 ðcos ot cos 3otÞ,
(25)
where c1 is a constant to be determined. Substituting Eqs. (10), (22) and (25) into Eq. (23), expanding the expression in a trigonometric series, and setting the coefficients of the resulting items cos(0ot) and cos(2ot) equal to zero, respectively, yield 12B þ 12Bo4 þ 38B3 o4 þ o4 c1 B2 o4 c1 c1 ¼ 0
(26)
2 4 12B þ 12Bo4 þ 12B3 o4 8o4 c1 11 2 B o c1 ¼ 0.
(27)
and
From Eq. (26) we can obtain c1 as follows: c1 ¼
4B þ 4Bo4 þ 3B3 o4 . 8ð1 o4 þ B2 o4 Þ
(28)
Substituting Eq. (28) into Eq. (27) and solving for the second analytical approximate frequency o2, we obtain pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi !1=4 40 þ 22B2 þ 2 256 þ 256B2 þ 71B4 o2 ðBÞ ¼ . (29) 72 þ 92B2 þ 25B4 Furthermore, c1 can be obtained by substituting Eq. (29) into Eq. (28) and the result is c1 ðBÞ ¼
B½64 þ 17B4 4D þ ð80 3DÞB2 , 4½32 þ 47B4 2D þ 2ð55 þ DÞB2
(30)
where DðBÞ ¼
pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 256 þ 256B2 þ 71B4 .
(31)
The corresponding second analytical approximate periodic solution is given by u2 ðtÞ ¼ ½B þ c1 ðBÞ cos½o2 ðBÞt c1 cos½3o2 ðBÞt.
(32)
The corresponding approximation to y is u2 ðtÞ B cos o2 t þ c1 ðcos o2 t cos 3o2 tÞ y2 ðtÞ ¼ qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi ¼ pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi . 2 1 þ ½B cos o2 t þ c1 ðcos o2 t cos 3o2 tÞ2 1 þ u ðtÞ
(33)
2
However, the analytical integration of Eq. (33) to obtain x2(t) is not possible. To obtain an analytical expression for x2(t), Eq. (25) is written as follows: B cos o2 t þ 4c1 ðcos o2 t cos3 o2 tÞ Bðcos o2 t þ ð4c1 =BÞ cos o2 t sin2 o2 tÞ y2 ðtÞ ¼ pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi ¼ qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi . 1 þ ½B cos o2 t þ 4c1 ðcos o2 t cos3 o2 tÞ2 1 þ B2 ðcos o2 t þ ð4c1 =BÞ cos o2 tsin2 o2 tÞ2
(34)
ARTICLE IN PRESS A. Bele´ndez et al. / Journal of Sound and Vibration 311 (2008) 1447–1456
1451
The above equation would be easily integrable if (4c1/B) cos o2t sin2o2t5cos o2t. Then Eq. (34) would have the same functional form as Eq. (19). From Eqs. (30) and (31) we can obtain the following limits: 4c1 ðBÞ ¼ 0, B!0 B pffiffiffiffiffi 4c1 ðBÞ 2337 þ 4283 71 pffiffiffiffiffi ¼ 0:12965. ¼ lim B!1 B 143867 þ 13822 71 lim
(35)
(36)
Then, 4c1(B)/B takes values between 0 and 0.12965 when B varies between 0 and N. For example, for B ¼ 1 and 10, 4c1(B)/B takes the values 0.0424 and 0.1271, respectively. We can write 4c1 2 2 cos o tsin o t (37) 2 2 p 0:12965 cos o2 tsin o2 t p0:05, B where we have taken into account that the maximum value of |cos o2t sin2 o2t| is 0.3849. In Fig. 1 we have plotted cos (2ph) and cos (2ph)+(4c1/B) cos(2ph) sin2(2ph) as a function of h ¼ o2t/2p ¼ t/T2 for B-N (4c1/B ¼ 0.12965). From Eqs. (34)–(37) and Fig. 1 we can conclude that Eq. (33) can be approximately written as follows: B cos o2 t y2 ðtÞ pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi . 1 þ B2 cos2 o2 t
(38)
Later we will verify that such a simple approximation gives very good results for x(t). Likewise, x2(t) can be calculated by integrating equation y ¼ dx/dt subject to the restrictions x2 ð0Þ ¼ 0;
B y2 ð0Þ ¼ pffiffiffiffiffiffiffiffiffiffiffiffiffiffi 1 þ B2
(39)
and this integration gives " # 3 2 1=4 1 B pffiffiffiffiffiffiffiffiffiffiffiffiffiffi sin½o2 ðBÞt . sin x2 ðtÞ ¼ 1 þ B 4 1 þ B2
(40)
We will show Eq. (40) gives good results for x2(t). However, we should not forget what we are really looking for is an approximate analytical solution to Eq. (1), that is, x(t). Moreover, it is convenient to express the approximate angular frequency and the solution in terms of oscillation amplitude A rather than as a function of B. It is now necessary to find a relation between oscillation amplitude A and parameter B used to solve Eq. (7) approximately.
Fig. 1. cos(2ph) (continuous line) and cos(2ph)+0.12965 cos(2ph)sin2(2ph) (dashed line) as a function of h ¼ o2t/2p ¼ t/T2.
ARTICLE IN PRESS 1452
A. Bele´ndez et al. / Journal of Sound and Vibration 311 (2008) 1447–1456
From Eq. (3) we get 1 ð1
y2 Þ1=2
1 þ x2 ¼ C, 2
(41)
where C is a constant to be determined as a function of initial conditions. From Eq. (20) we can easily obtain C ¼ (1+B2)1/2 and Eq. (41) can be written as follows: 1 ð1
y2 Þ1=2
1 þ x2 ¼ ð1 þ B2 Þ1=2 . 2
(42)
In addition, when x ¼ A, the velocity y ¼ dx/dt is zero. Taking this into account in Eq. (49), we obtain the following relation between amplitude A and parameter B: 1 þ 12A2 ¼ ð1 þ B2 Þ1=2 . From the above equation we can easily find that the solution for B is 1=2 B ¼ A 1 þ 14A2 .
(43)
(44)
Substituting Eq. (44) into Eqs. (13) and (21), the first analytical approximate periodic solution for the relativistic oscillator as a function of the oscillation amplitude A is given by 3 2 1=4 o1 ðAÞ ¼ 1 þ 34A2 þ 16 A , (45) 2sffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 3 2 4 1 4A þ A x1 ðtÞ ¼ sin1 4 sin½o1 ðAÞt5. o1 ðAÞ 4 þ 4A2 þ A4
(46)
Substituting Eq. (44) into Eqs. (29) and (40), the second analytical approximate periodic solution for the relativistic oscillator as a function of the oscillation amplitude A is given by pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi !1=4 640 þ 352A2 þ 88A4 þ 8 4096 þ 4096A2 þ 2160A4 þ 568A6 þ 71A8 o2 ðAÞ ¼ , (47) 1152 þ 1472A2 þ 768A4 þ 200A6 þ 25A8 2sffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 3 2 4 1 4A þ A sin1 4 sin½o2 ðAÞt5. x2 ðtÞ ¼ o2 ðAÞ 4 þ 4A2 þ A4
(48)
Now we illustrate the applicability, accuracy and effectiveness of the proposed approach by comparing the approximate analytical periodic solutions obtained in this paper with the exact ones. Calculation of the exact period frequency, Tex(A), proceeds as follows. Substituting Eq. (44) into Eq. (42) and integrating, we obtain Z A 1 þ 12ðA2 x2 Þ ffi dx, qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi T ex ðAÞ ¼ 4 (49) 0 A2 x2 þ 14ðA2 x2 Þ2 which allows us to obtain the exact angular frequency, oex(A), in terms of elliptical integrals as follows: " #1 pffiffiffiffiffiffiffiffiffiffiffiffiffiffi A2 2p 8 A2 2 ¼ 2p 4 4 þ A E , (50) oex ðAÞ ¼ pffiffiffiffiffiffiffiffiffiffiffiffiffiffi K T ex ðAÞ 4 þ A2 4 þ A2 4 þ A2 where K(q) and E(q) are the complete elliptic integrals of the first and second kind, respectively, defined as follows [18]: Z 1 Z 1 rffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi2 dz 1 qz pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi ; EðqÞ ¼ KðqÞ ¼ dz. (51) 2 2 1 z2 ð1 z Þð1 qz Þ 0 0
ARTICLE IN PRESS A. Bele´ndez et al. / Journal of Sound and Vibration 311 (2008) 1447–1456
1453
For small values of the amplitude A it is possible to take into account the following power series expansion of the angular frequencies oex ðAÞ 1
3 2 51 4 A þ A , 16 1024
(52)
o1 ðAÞ 1
3 2 42 4 A þ A , 16 1024
(53)
3 2 56 4 A þ A . (54) 16 1024 These series expansions were carried out using MATHEMATICA. As can be seen, in the expansions of the angular frequencies o1(A) and o2(A), the first two terms are the same as the first two terms of the equation obtained in the power-series expansion of the exact angular frequency, oex(A). If we compare the third term in Eqs. (53) and (54) with the third term in the series expansion of the exact frequency oex(A) (Eq. (52)), we can see that the third term in the series expansions of o2(A) is more accurate than the third term in the expansion of o1(A). For very large values of the amplitude A it is possible to take into account the following power series expansion of the exact angular frequency: o2 ðAÞ 1
p 1:5708 þ ¼ þ , 2A A
oex ðAÞ o1 ðAÞ
2 3
1=4
A
þ ¼
1:5197 þ , A
pffiffiffiffiffi1=4 88 þ 8 71 1 1:5790 o2 ðAÞ þ ¼ þ . 25 A A
(55) (56)
(57)
Once again we can see than o2(A) provides excellent approximations to the exact frequency oex(A) for very large values of oscillation amplitude. Furthermore, we have the following equations: lim oex ðAÞ ¼ lim o1 ðAÞ ¼ lim o2 ðAÞ ¼ 1,
(58)
lim oex ðAÞ ¼ lim o1 ðAÞ ¼ lim o2 ðAÞ ¼ 0,
(59)
A!0
A!1
A!0
A!1
lim
o1 ðAÞ
A!0
A!1
¼ 0:96745,
(60)
o2 ðAÞ ¼ 1:00523. A!1 oex ðAÞ
(61)
A!1 oex ðAÞ
lim
Eqs. (58)–(61) illustrate the very good agreement of the approximate frequency o2(A) with the exact frequency oex(A) for small as well as large values of oscillation amplitude. In Fig. 2 we plotted the percentage error of approximate frequencies o1 and o2, as a function of A. In this figure the percentage errors were calculated using the following equation: oj oex ; j ¼ 1; 2. (62) Relative error of oj ð%Þ ¼ 100 oex As we can see from Fig. 2, the relative errors for o1(A) are lower than 0.82% for all the range of values of amplitude of oscillation A (see Eq. (60)), and these relative errors tend to 0.52% when A tends to infinity (see Eq. (61)). The exact periodic solutions x(t) achieved by numerically integrating Eq. (1), and the proposed first-order approximate periodic solutions x1(t) in Eq. (46) and x2(t) in Eq. (48) are plotted in Figs. 3–6 for A ¼ 0.1, 1, 2 and 10 (B ¼ 0.1001, 1.12, 2.83 and 51), respectively. In these figures parameter h is
ARTICLE IN PRESS 1454
A. Bele´ndez et al. / Journal of Sound and Vibration 311 (2008) 1447–1456
Fig. 2. Relative error for approximate frequencies o1, (Eq. (45), dashed line) and o2 (Eq. (47), continuous line).
Fig. 3. Comparison of the normalized approximate analytical solutions x1/A (m) and x2/A (J) with the exact solution (continuous line) for A ¼ 0.1 (b0 ¼ v0/c ¼ 0.09963).
Fig. 4. Comparison of the normalized approximate analytical solutions x1/A(m) and x2/A (J) with the exact solution (continuous line) for A ¼ 1 (b0 ¼ v0/c ¼ 0.74536).
ARTICLE IN PRESS A. Bele´ndez et al. / Journal of Sound and Vibration 311 (2008) 1447–1456
1455
Fig. 5. Comparison of the normalized approximate analytical solutions x1/A (m) and x2/A (J) with the exact solution (continuous line) A ¼ 2 (b0 ¼ v0/c ¼ 0.94281).
Fig. 6. Comparison of the normalized approximate analytical solutions x1/A (m) and x2/A(J) with the exact solution (continuous line) A ¼ 10 (b0 ¼ v0/c ¼ 0.99981).
defined as follows: h¼
t ¼ 2poex ðAÞt. T ex ðAÞ
(63)
From Eqs. (20) and (44) we have v0 B b0 ¼ ¼ yð0Þ ¼ pffiffiffiffiffiffiffiffiffiffiffiffiffiffi ¼ c 1 þ B2
sffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 4A2 þ A4 , 4 þ 4A2 þ A4
(64)
where v0 is the maximum velocity of the particle and c is the velocity of light. Figs. 3–6 show that Eq. (48) provides a good approximation to the exact periodic solutions and that the approximation considered in Eq. (38) is adequate to obtain the approximate analytical expression of x2(t) given in Eq. (48). As we can see, for small values of A (Figs. 3 and 4) x(t) is very close to the sine function form of non-relativistic simple harmonic motion. For higher values of A the curvature becomes more concentrated at the turning points (x ¼ 7A). For these values of A, x(t) becomes markedly anharmonic and is almost straight between the turning points. Only in the vicinity of the turning points, where the magnitude of the
ARTICLE IN PRESS A. Bele´ndez et al. / Journal of Sound and Vibration 311 (2008) 1447–1456
1456
Hooke’s law force is maximum and the velocity becomes relativistic, is the force effective in changing the velocity [1]. Fig. 6 are a typical example of the motion in the ultra-relativistic region where b0-1. In summary, the harmonic balance method was used to obtain two approximate frequencies for the relativistic oscillator. To do this we rewrite the nonlinear differential equation in a form that does not contain an irrational expression. We can conclude that the approximate frequencies obtained are valid for the complete range of oscillation amplitude, including the limiting cases of amplitude approaching zero and infinity. Excellent agreement of the approximate frequencies with the exact one was demonstrated and discussed and the discrepancy between the second approximate frequency, o2, and the exact one never exceeds 0.82%. Some examples have been presented to illustrate the excellent accuracy of the approximate analytical solutions. Finally, we can see that the method considered here is very simple in its principle, and is very easy to apply. Acknowledgments This work was supported by the ‘‘Ministerio de Educacio´n y Ciencia’’, Spain, under Project FIS2005-05881C02-02, and by the ‘‘Generalitat Valenciana’’, Spain, under Project ACOMP-2007-020. We would like to thank the anonymous reviewers for their comments and suggestions, which have undoubtedly improved the original manuscript. References [1] [2] [3] [4] [5] [6] [7] [8] [9] [10] [11] [12] [13] [14] [15]
[16] [17] [18]
W. Moreau, R. Easther, R. Neutze, Relativistic (an)harmonic oscillator, American Journal of Physics 62 (1994) 531–535. J.L. Synge, in: S. Flu¨gge (Ed.), Classical Dynamics in Encyclopedia of Physics, Vol. III/1, Springer, Berlin, 1960, pp. 210–211. H. Goldstein, Classical Mechanics, second ed, Addison-Wesley, Reading, MA, 1981. R.E. Mickens, Periodic solutions of the relativistic harmonic oscillator, Journal of Sound and Vibration 212 (1998) 905–908. R.E. Mickens, Oscillations in Planar Dynamics Systems, World Scientific, Singapore, 1996. C.W. Lim, B.S. Wu, A new analytical approach to the Duffing-harmonic oscillator, Physics Letters A 311 (2003) 365–373. B.S. Wu, C.W. Lim, Large amplitude nonlinear oscillations of a general conservative system, International Journal of Non-linear Mechanics 39 (2004) 859–870. A. Bele´ndez, A. Herna´ndez, A. Ma´rquez, T. Bele´ndez, C. Neipp, Analytical approximations for the period of a simple pendulum, European Journal of Physics 27 (2006) 539–551. H. Hu, J.H. Tang, Solution of a Duffing-harmonic oscillator by the method of harmonic balance, Journal of Sound and Vibration 294 (2006) 637–639. H. Hu, Solution of a quadratic nonlinear oscillator by the method of harmonic balance, Journal of Sound and Vibration 293 (2006) 462–468. H.P.W. Gottlieb, Harmonic balance approach to limit cycles for nonlinear jerk equations, Journal of Sound and Vibration 297 (2006) 243–250. A. Bele´ndez, A. Herna´ndez, T. Bele´ndez, M.L. A´lvarez, S. Gallego, M. Ortun˜o, C. Neipp, Application of the harmonic balance method to a nonlinear oscillator typified by a mass attached to a stretched wire, Journal of Sound and Vibration 302 (2007) 1018–1029. G.R. Itovich, J.L. Moiola, On period doubling bifurcations of cycles and the harmonic balance method, Chaos, Solitons and Fractals 27 (2005) 647–665. Y.M. Chen, J.K. Lui, A new method based on the harmonic balance method for nonlinear oscillators, Physics Letters A 368 (2007) 371–378. A. Bele´ndez, C. Pascual, A. Ma´rquez, D.I. Me´ndez, Application of He’s homotopy perturbation method to the relativistic (an)harmonic oscillator. I: comparison between approximate and exact frequencies, International Journal of Non-linear Science and Numerical Simulations 8 (4) (2007) 483–491. J.H. He, Some asymptotic methods for strongly nonlinear equations, International Journal of Modern Physics B 20 (2006) 1141–1199. A. Bele´ndez, A. Herna´ndez, T. Bele´ndez, A. Ma´rquez, Application of the homotopy perturbation method to the nonlinear pendulum, European Journal of Physics 28 (2007) 93–104. L.M. Milne-Thomson, Elliptic integrals, in: M. Abramowitz, I.A. Stegun (Eds.), Handbook of Mathematical Functions, Dover Publications, Inc., New York, 1972.