- Email: [email protected]
A method of calculation of steady state oscillations in autonomous non-linear systems
Journal ofSound
and Vibration (1990) 136(l), 35-44
A METHOD OSCILLATIONS
OF
CALCULATION
IN AUTONOMOUS
OF STEADY
STATE
NON-LINEAR
SYSTEMS
S. S. QIU South China Universit_vof Technology, Guangzhou, People’s Republic of China AND 1.
Department of Electrical Engineering,
M.
FILANOVSKY
Universit_yof Alberta, Edmonton, Alberta, Canada T6G 2Gl
(Received 20 December 1988, and in revised form 13 April 1989)
An extension of the harmonic balance method is presented for calculation of the steady state (periodic) solutions in autonomous non-linear systems. The non-linear differential equation is expanded into a system of linear differential equations which are solved consecutively. An approximate solution which includes a few harmonics (the main oscillation) is found from the first equation of this system by the ordinary harmonic balance method. Other harmonics and the corrections of the main oscillation are found as the forced solutions of the second, third and the following equations of the system. To improve the efficiency of correction calculations a systematic method for derivation of the forcing term in these equations is given. As an example, a non-linear conservative system with a non-symmetric non-linearity of the polynomial type is considered.
1. INTRODUCTION method of harmonic balance is one of the most frequently used methods available for the analysis of periodic oscillations in non-linear systems. The method is a relatively simple one since it leads to a system of algebraic equations, even though non-linear ones. The method is especially effective for oscillating systems for which the first harmonic of the oscillation, even if only approximately calculated, provides adequate information on the system [ 1,2]. However, this method becomes inefficient very quickly when one tries to apply it for calculation of additional harmonics, even with a computer at hand. The system of non-linear algebraic equations becomes too complicated [3] for the analysis and the solutions can be overlooked. Besides, it has been observed that the method may also lead to inaccurate and inconsistent results in some applications. The method was criticized [4] but the criticism, in our opinion, should be directed not on the fact that “one needs to know a great deal about the solution a priori or to carry enough terms in the solution and check the order of the coefficients of all neglected harmonics” [2] but on the inability of the method, in its original form, to provide sufficient precision for calculation of a relatively small number of some initial harmonics (if the solution is represented as a series) and, thus, to handle effectively these neglected harmonics. In many cases the terms which should be included in the solution are known. In view of the relative conceptual simplicity of this classical method it was recently modified [5] by addition of the perturbation equations from the system differential equation. This allows the consideration of additional harmonics in a uniform and a simple 35 The
0022-460X/90/010035
+ 10 !SO3.00/0
0
1990 Academic
Press Limited
S. S.
36
QIU
AND
1. M.
FILANOVSKY
way. Here we describe another modification of the method for the same purpose, which proceeds as follows. The original non-linear differential equation is expanded into a system of linear differential equations which are solved consecutively. An approximate solution which includes a few harmonics (the main oscillation) is found from the first equation of this system by the ordinary harmonic balance method. This equation gives only the non-linear algebraic equations for the amplitudes in the main oscillation. Corrections for the main oscillation and other harmonics are found as the forced solutions of the second, third and the following equations of the system: i.e., the amplitudes of the corrections are found from linear algebraic equations. The harmonic content of each consecutive correction is fully determined by the terms neglected in the solution of the previous equation. To improve the efficiency of correction calculations a systematic method for derivation of the forcing terms in the differential equations is given. As a result a previous correction is calculated with maximum precision and the subsequent correction in many cases does not even include the harmonics of the previous correction but involves only new harmonics. At first glance, the proposed method requires the notorious a priori knowledge of the harmonic content of the main oscillation. The approach, however, has been designed to increase the efficiency of calculations when the main oscillation includes the first harmonic only. Thus, the proposed method, which is a further refinement (a better precision of calculation is achieved) of a previous approach [6,7], is in the same category as many other methods [8,9] in which only one harmonic is assumed in the initial approximation. Yet the efficiency of correction calculations in the proposed method is such that the first two corrections are, in many cases, sufficient to obtain a solution.
2. EXPANSION
We consider
the autonomous
OF THE
systems
DIFFERENTIAL
which are described
a,w”x’“‘+a,,_,w”-lx’“-I’+. Here x(I) = d’x/dr’,
where
EQUATION
by the differential
. ~+a”x+K(w)f=O.
T = or is the normalized
time. We assume
f=f[x, x(I), XC>), . . x(P’] .)
equation (1)
that
(2)
is a function which can be represented by a Taylor’s type series expansion when the variable x is represented as a sum of successively diminishing terms (functions of a polynomial type, for example, satisfy this condition) and K(w) is a polynomial. The periodic solution of equation (1) is obtained as a sum x=x,+x,+x,+.
. .)
(3)
is called the main oscillation and where the term x0, which includes a few harmonics, the terms x, , x2,. . . are corrections. Each subsequent correction may include the same harmonics as the previous ones and additional harmonics of other frequencies. The frequency w is obtained as a sum w=6Jo+w,+w~+...,
(4)
of the approximate value w,,. where w,, w2,... are small corrections The approximate periodic solution of equation (1) is expressed as x,, = A,,, cos 7 + C (Ak,, cos k~+ &, sin kr), ktEi,
(5)
AUTONOMOUS
NON-LINEAR
SYSTEMS
37
where the set E0 identifies the harmonics included in x0, which are different from the first one (thus, k # 1). If one substitutes equation (5) into equation (2) one obtains the result
.hl=f[.%, 4’),XC’.,. . , xl?1 =hm + uil -.I& 1,
(6)
where f;,m = C,, cos T+ D,,, sin 7 + C (C,, cos kT+ Dko sin kT) kEE,,
(7)
includes the same harmonics as the main oscillation and the function (fo-A,,,,) includes the harmonics of other frequencies. Now one can rewrite equation (1) as a,o”x’“‘+
a,_,W*-‘X(n-‘)+.
. ~+~ox+K(~l.hm
=-K(w)[(f,-f,,)+(f-f,,)l.
(8)
Assuming that the higher harmonics with large amplitudes are included in the second term of equation (5) one can consider that the terms in the brackets, on the right side of equation (8), represent small quantities. Following reference [lo], one can consider, instead of equation (8), a more general equation,
. *+a,x+~(wIh,
a,w”x’“‘+a,_,w”~‘x’“-“+.
= -K(W)[&(fO-fOm)+(f-f0)1, (9)
with O< F G 1. When E = 1 equation (9) becomes the same as equation (8) (and, thus the same as the initial equation (1)). It is natural in this case to seek a solution of equation (9) as a power series x=x,f&X,+E*Xz+~**.
(10)
One can substitute this series into equation (2) and obtain the expansion f=fo+&f+&Ej2+...,
(11)
where f0 depends on x0 and its derivatives only, f, depends on x,, x, and their derivatives only, f2depends on x0, x, , x2 and their derivatives only, etc. To obtain precision of calculation better than in references [6,7] the procedure which was done for x0 is repeated for x, . As was mentioned before, x, is a sum of harmonics as well: i.e., x, = C (Ak, cos kT+ Bk, sin k7), keE,
(12)
where the set E, identifies the harmonics included in x,. If this sum is substituted into
f,one obtains (assuming that x0 is known at this stage) that f, =f,(xo,
x,) =f,m -c(f,
-f,m)t
(13)
where f,,,, =
c
htE,
(CA, cos kT+Dk,
sin kT),
(14)
with C,, and Ok, being linear functions of Ak, and BI, (they can be, of course, non-linear functions of the amplitudes of harmonics which are included in x0). The residual term f, -f,m includes the harmonics with numbers k F?E, or the harmonics with k E E, but with amplitudes which are non-linear functions of AL, and Bk,. The same calculation is carried out for all the f;: i.e., for each f; one finds the form
f;(X”,
XI
7..
‘,x,)=.L+(f;-.L)
(15)
38
S. S. QIU
as soon as the identifying
AND
I. M. FILANOVSKY
set E, of x, = 1 ( Ak, cos kT + B,, sin k7) A\~ F
is shown. Equation (13) and all the equations (15) are substituted in their direct form but in the form of
./-I=f,,n +&(.A-f,mL
(16)
into equation
f2=fim+4fi-.Ll),
(1 l), except not
...> (17)
f;=.L+e(_L-.L), where
.s is the previously
introduced
....
parameter.
One obtains
that
One can now substitute o=w,+&W,+& into K (w ) and obtain
wz+.
..
(19)
the expansion
K(w) = K(w,)+FK(w~,)w,+F~K~(w~,
w,)w2+.
9..
(20)
Finally, one can substitute equations (10) and (18)-(20) into equation (9) and equate the terms which have the same exponents of E. Then the following system of equations is obtained u,w~x~~‘+a,
,w;;-‘xl;‘~ ‘I+. . ~+a,,x,,+K(w,,)fc
a,w~x(ln’+a,,~,w;l~‘.~~*~‘)
+. .
,,,, =o,
.+ua,,x,+K,(w,,)w,f,,,,+K(wol.fi,,,
=w,P,(x,,w,,)-K(w,,)(.f,,-f;,,), a,WgnX:n’+a,_,Wgn-‘X~,l--“+~. . + a,+,+K,(w,,, w,h.f;~m+ K,(wohfim + K(wol_L
= wJ’z(xo, x,, w,) - K,(w,b,(h,-An)
- K(wo)(f, -h,n),
. ..>
(21)
where P, and PI are polynomials of their arguments (they are obtained from the terms additional to the terms of the left parts in equation (21) when the expansions (10) and (19) are substituted in the linear part of equation (3)). Equations (21) should be solved consecutively to obtain x,,, x,, x2,. . . and w,,, w,, w2,. . . . . The second equation includes the additional term K ( wo).f,mr the third the additional term K,(w,,)w,f,, + K (w,)frm, etc. These terms did not appear in the approach used earlier [6,7], but they essentially increase the precision of calculations of the first few harmonics. It is not difficult to see that to return from equation (9) to equation (8) one has to put E = 1. The same can be done in equations (17), but this means that one can solve equations (21) and find the required solution, using equation (3) (and the required frequency from equation (4)), by summation of the solutions obtained. The first equation in system (21) is the usual equation of the harmonic balance method, the second and other consecutive equations are linear differential equations. The solution of the second and other consecutive equations does not present any difficulties because one has to find only their steady state solutions: i.e., one uses the harmonic balance method again and solves the linear algebraic equations to find the amplitudes in x,, w,, w?, . . . x2,. . . and the frequency corrections
AUTONOMOUS
3. SIMPLIFlCATION
NON-LINEAR
OF THE LINEARIZED
39
SYSTEMS
EQUATIONS
The following steps for simplification of the system of calculations. (1) If the final form of the solution is set as
AND CALCULATION
OF f;,,
(21) can be used from the onset
x=A,cos7+A2cos2T+Brsin27+*..,
(22)
then A, =A,o+A,,+A,,7+.
*. ,
n=l,2,3,...
i.e., B,i=O (i=1,2,3 ,... ). (2) In conservative systems the harmonic in this case, can have the form
(23)
. . 1.’
n = 2,3,4,.
1 B,=B,o+B,,+B,z+...,
phases are independent
[5,8] and the solution,
~=A,,+A,cos7+A~cos2~+~~~.
(24)
(3) In conservative systems and systems with odd non-linearities w, = 0. (4) In systems with odd non-linearities A0 = 0 and A,, = 0. Additionally, the calculation of the terms l;,n in equations (21) can be carried out in a systematic way. If one considers that in equation (2) x,, .x2,. . . are small perturbations of x0 one can find that the perturbation Sf; off caused by the perturbation x, is
- af 4’l \ :\,\*, a/;=$.=x,+fd[X’
14)
\
(25)
9
q=l
.?I
where summation embraces all the derivatives included in equation (2). It is not difficult to see that Sf; can be obtained by the substitution of x, instead of x, in f,.Then if x, = 1 (Ak, cos kT+ BL, sin kT), ACE, the substitution
of equation
(26) into equation
(26)
(25) will give (27)
where 8kf; is the total contribution of the kth harmonic content of x,, is determined one can find that
fikf; = 1 ( Cki COSk7 i- Dki Combining harmonics
the parts of contributions included in x, one obtains
Sill
k7) + c
included
( ck,
COS
in Xi. When the harmonic
kT + &
from all harmonics having the required expression
Sill k7).
(28)
the frequencies
J,,, = 1 1 (C,; cos k7+ AL, sin kT) hk t, 1.. I;, Thus,
all the functions
f;“? are calculated
uniformly,
of
(29)
from the same function
f,.
4. EXAMPLE Consider equation
the oscillations
in a non-linear
conservative
“~x’~‘+x+cu~x’+a,x~=o.
system described
by the differential
(30)
40
S. S. QIU
AND
I. M.
FILANOVSKY
This is the same example which forced some authors [4] to abandon the harmonic balance techniques and was used for verification of the modification of the harmonic balance method in reference [5]. The initial conditions for the system considered are taken as the following: i-1,3,....
x,(O)= A,,xi(O)+x,+,(O)=O,
(31)
For conservative systems w, = 0. As a result the system (21) is simplified and becomes o;x:“+x,,+f”,
=o,
w;x:“+xz+fim
wW+x,+“L
= -(fo-fom), (32)
= -2x)i’w,,wz - (f, --j-i.,),
. . *.
For the example considered f= (u2x2+ (Y3x3
(33)
and to stress the possibilities of the proposed method we will carry out the calculations assuming only one harmonic in the main oscillation: i.e., we set x,, = A,” cos T.
(34)
If one substitutes equation (34) into equation (33) one obtains f0 = cxZA&cos* 7 + a3A;,, cos3 T = &,A:,
+ acy3A:, cos T+ fa2Af0 cos 2r + acu3Aiocos 37,
(35) so that j& = &,A;,
cos T
(36)
and f0 -f&, = $x,A$ + fa,A:,
cos 27 + &,A;,
cos 3 7.
(37)
Substituting equation (36) in the first equation of system (32) one obtains -w;A~,,+A,,+&~A~~=O.
Following the initial conditions
(38)
(31) A,, = A, and in this case one finds that w0= 41 +&A:.
(39)
The substitution of x = x0+ FX,+ e2x7. . . into f given by equation (33) gives f, = 2azx,x, + 3a,x;x,
)
(40)
and the substitution of x, instead of x, in this result gives sf; = (2~ax,~+3cz3x;)x,.
(41)
This is a conservative system and for it one can choose x,=
1 Ak, cos
hi E,
kr.
(42)
The substitution of the obtained x0 = A, cos T and xi given by equation (42) into equation (41) results in S,f; = (2azA, cos 7+3a,Af
cos’ T)A,, cos kT
=~(YJA~&COS(~-~)T+(Y~A,A~;COS(~-~)T
+;cu,AfAk,
cos kT+ + cuIA,Ak, cos (k+ 1)T
+&AfAki
cos (k +
2)T,
(43)
AUTONOMOUS
NON-LINEAR
SYSTEMS
41
and this simple trigonometric expansion is the key formula of the proposed method. Equation (37) gives the right side for the second equation in the system (32). As a result one has to set x, =A,,
+Az,
cos 27+A3,
cos 37.
(44)
E, = {0,2,3}
and using equation (43) (i.e., one substitutes i = 1 and then, for each the terms which give the contribution into the zeroth, or into the first, or into the third harmonic and then combines these terms together) one obtains
Thus,
k = 0, 1, 3 one selects
fin, =$,A:(A,,
+:A~,)+[&A:(A~,
+(cx2A,A2,+$r3A:A3,)
+A~,)+(Y~A,A~,]
cos 27
cos 37.
(45)
Further calculations will show that AZ, << AO, and for a more simple representation final results we will use _f ,m =&A:A~,
+ (&A:A~,
Substituting equation one can find that
+ a2~,~3,)
(46) and equation
cos 2T+(a,A,~,,
(37) into the second
+&A;)
equation
cos 3T.
of the
(46)
of the system (32)
where A = 24 +,A;
- &A;.
(48)
To solve the third equation of the system (32) we, first of all, will find f, . The substitution of x, given by equation (44) into equation (40) gives f, =;a3Af(Ao,
+A2,)
+ (2azA,A,,,
+ a>A,A,,
+&r3A:A3,)
cos T
+[~~,A:(A~,,+A~,)+(Y~A,A~,]COS~T+(((Y~A,A~,+~~~A:A~,)COS~T + (ct2A,A3, +:(Y~A;A~,)
and the subtraction
of equation
f, -f,,,,
cos
4~+~cqA~x‘t3,
cos
(45) from equation
= (2qA,A,,, +((Y,A,A,,
(49)
(49) gives
+ (Y:A,A~, +:cx3A:A3,) +:qA,AZ,)
5T,
cos
cos 7
~T+&,A;A~,
cos
57.
(50)
Ideally, the calculation (49) should not be done at all because the terms of the right side of equation (50) represent the terms discarded in the calculation ofJim (just as the terms discarded in the calculation off”,,, in equations (35)-(37) werefo-f&,,); we have to retain them only. Now one can conclude that x,=A,~cosT+A~~cos~T+A~~cos~T.
(51)
We will not carry out the calculation of the fourth and fifth harmonics (even though it is easy to do here and not so easy a task when using other methods) and instead will restrict ourselves to the calculation of the second harmonic and the frequency correction. Thus, we set E2 = (2) and, using equation (43) with k = i = 2, one can find that f?,, = &A:A,Z
COS
T+ .
+ . .
(52)
S.S.
42 Substitution equation
of equations
QIU
AND
I.
M.
FILANOVSKY
(50) and (52) in the third equation
&,A:A,~= fw,,w,A, - (2a:A,A,,, From the initial
+ a?A,A,,
of the system
+&,A;A,,).
(53)
(31) for i = 1 one obtains
conditions
A ,2=-(Ao,+Ax+A3,), and substitution
of equation
(54) into equation
w2 = (1/2w,)[(2a,-~cu,~,)A,, A brief analysis
of the results
obtained
(54)
(53) gives
+~cQ-$c+~,~A~,
-&+,A~,].
can be done by assuming
allows
(55)
that
0<(~~A:
(32) gives the
(56)
one to use the approximation
&A;
oo=l+
(57)
and to compare the results with results obtained by other methods. are of small amplitude in a Case (a): ff2” (Yeand A,<<1. In this case the oscillations system with equal influence of the non-linear terms, and one is close to the left side of the inequality (56). Thus, = -;a,A:+A,
x=x,+x,+x,
COS
d&+4;
COS27f’
. .
(58)
and w=oo+oz
~(Y~A~-~~A~.
=l+
(59)
Case (b): CY~<< 1 but (Y~> 1 and a,Af= 1. Here the oscillations are of a larger amplitude in a system with asymmetric influence of the non-linearity terms, and one is close to the right side of the inequality (56). In this case X=
-&,A:+A,
cos
7+&~3~:COS37
(60)
and w-l+
&qAf-~+x~A;.
(61)
The results (58) and (59) were obtained by other methods as well [4,5] but their modification into (60) and (61) caused by the change of system parameters and the oscillation amplitude is not so easily traced in other methods. To conclude this section we would like to stress that the corrections of the harmonics which have just been found will be corrected again not in the next equation but in the equation which follows it. This means that the calculated results have good precision.
5. DISCUSSION
The proposed method should be applied, as many other methods, in conditions when it is known that a periodic solution exists. The usual approach is to solve the first equation of the system (21) and to find the approximate solution x,). If this solution exists (i.e., one is able to find real values for A,,,", l& and wO) one assumes that the exact periodic solution exists as well. Even though this assumption requires an additional mathematical proof, it is widely used in many other methods of calculations; for example, in the describing function method [ 1,2]. As in many other methods, the calculation of the solution should be followed by verification of its stability [ 111.
AUTONOMOUS
NON-LINEAR
SYSTEMS
43
The question of the spectral content of the main oscillation should not be considered as a major defect of the harmonic balance method. Of course, the modifications providing efficient calculations of harmonics are desirable. All other methods should be considered from this point of view as well, and we think that there is no universal method which can predict all possible situations in non-linear systems. The efficient calculations of harmonics can help to discover a situation that was overlooked earlier. The correct content of the main oscillations is a desirable feature improving the convergence of the calculations. The answer is sometimes determined by the problem itself (for example, in the calculation of subharmonic oscillations the main oscillation will include a subharmonic and the oscillation with the frequency of the external force), some help can be obtained from the investigation of resonance properties of the linearized second equation in the system (21). But in any case the convergence of equation (9) and (in the limit) of its initial form (8) was proved [ 121 for the mentioned limitations of xi. In calculation practice it is reasonable to check that Si + 0, where S;= C
(A:,+@,)
(62)
k
for the correction xi. Besides, it is convenient to compare S, with Si_, to reach a conclusion about the smallness of a correction. 6. CONCLUSIONS The features of the approach presented can be summarized the following way. 1. The proposed modification of the harmonic balance method allows the efficient calculation of a small amount (up to ten) harmonics in the oscillation. This is usually enough for many engineering problems. 2. The method can be extended for systems described by differential equations of order higher than the second. In such cases, one obtains only more extended systems of linear algebraic equations. 3. The method requires that the non-linearity function be an analytical function of the variable and its derivatives. Even though many systems can be reduced to this case, this is an essential limitation of the proposed approach. 4. The method is especially efficient when the main oscillation x0 includes one harmonic only and the corrections x, and x2 are sufficient to obtain the solution.
ACKNOWLEDGMENT This work was supported in part by the National Science Foundation under Grant 6872011.
of China (PRC)
REFERENCES 1. D. P. ATHERTON 1982 Nonlinear Control Engineering. New York: Van Nostrand. 2. D. P. ATHERTON and H. T. DORRAH 1980 International Journal of Control 31(6), 1041-1105. A survey on nonlinear oscillations. 3. C. HAYASHI 1964 Nonlinear Oscillations in Physical Systems. New York: McGraw-Hill. 4. A. H. NAYFEH and D. T. MOOK 1979 Nonlinear Oscillations. New York: John Wiley. Vibrations 95,525-530. An intrinsic 5. A. S. ATADAN and K. HUSEYIN 1984 JournalofSoundand method of harmonic analysis for nonlinear oscillations (a perturbation technique). 6. SHUI-SHENG QUI and I. M. FILANOVSKY 1987 IEEE Transactions on Circuits and Systems, CAS-34, 913-918. Periodic solutions of Van der Pol equation with moderate values of damping coefficient.
44
S. S. QIU AND
I.M. FILANOVSKY
7. SHUI-SHENG QIU, I.M. FILANOVSKY and K. A. STROMSMOE 1987 ControL77wory and Advanced Technology 3(3), 189-195. On one combination of the harmonic balance method and perturbation techniques. 8. N. MINORSKY 1962 Nonlinear Oscillafions. Princeton, N.J.: .Van Nostrand. 1952 Self-oscillating Systems. Moscow: Gostechizdat. (In Russian). 9. K. F. THEODORCHIK 1954 Oscillations. Moscow: Gostechizdat. (In Russian). 10. 8. V. BULGAKOV 11. D. P. ATHERTON 1981 Stability of Nonlinear Systems. New York: John Wiley. 1942 Applied Mathematics and Mechanics (Prikladnaia Matematika i 12. B. V. BULGAKOV Mechanika, in Russian) 6(4), 263-280. On the application of Poincare’s method to free pseudolinear oscillatory systems.
and Vibration (1990) 136(l), 35-44
A METHOD OSCILLATIONS
OF
CALCULATION
IN AUTONOMOUS
OF STEADY
STATE
NON-LINEAR
SYSTEMS
S. S. QIU South China Universit_vof Technology, Guangzhou, People’s Republic of China AND 1.
Department of Electrical Engineering,
M.
FILANOVSKY
Universit_yof Alberta, Edmonton, Alberta, Canada T6G 2Gl
(Received 20 December 1988, and in revised form 13 April 1989)
An extension of the harmonic balance method is presented for calculation of the steady state (periodic) solutions in autonomous non-linear systems. The non-linear differential equation is expanded into a system of linear differential equations which are solved consecutively. An approximate solution which includes a few harmonics (the main oscillation) is found from the first equation of this system by the ordinary harmonic balance method. Other harmonics and the corrections of the main oscillation are found as the forced solutions of the second, third and the following equations of the system. To improve the efficiency of correction calculations a systematic method for derivation of the forcing term in these equations is given. As an example, a non-linear conservative system with a non-symmetric non-linearity of the polynomial type is considered.
1. INTRODUCTION method of harmonic balance is one of the most frequently used methods available for the analysis of periodic oscillations in non-linear systems. The method is a relatively simple one since it leads to a system of algebraic equations, even though non-linear ones. The method is especially effective for oscillating systems for which the first harmonic of the oscillation, even if only approximately calculated, provides adequate information on the system [ 1,2]. However, this method becomes inefficient very quickly when one tries to apply it for calculation of additional harmonics, even with a computer at hand. The system of non-linear algebraic equations becomes too complicated [3] for the analysis and the solutions can be overlooked. Besides, it has been observed that the method may also lead to inaccurate and inconsistent results in some applications. The method was criticized [4] but the criticism, in our opinion, should be directed not on the fact that “one needs to know a great deal about the solution a priori or to carry enough terms in the solution and check the order of the coefficients of all neglected harmonics” [2] but on the inability of the method, in its original form, to provide sufficient precision for calculation of a relatively small number of some initial harmonics (if the solution is represented as a series) and, thus, to handle effectively these neglected harmonics. In many cases the terms which should be included in the solution are known. In view of the relative conceptual simplicity of this classical method it was recently modified [5] by addition of the perturbation equations from the system differential equation. This allows the consideration of additional harmonics in a uniform and a simple 35 The
0022-460X/90/010035
+ 10 !SO3.00/0
0
1990 Academic
Press Limited
S. S.
36
QIU
AND
1. M.
FILANOVSKY
way. Here we describe another modification of the method for the same purpose, which proceeds as follows. The original non-linear differential equation is expanded into a system of linear differential equations which are solved consecutively. An approximate solution which includes a few harmonics (the main oscillation) is found from the first equation of this system by the ordinary harmonic balance method. This equation gives only the non-linear algebraic equations for the amplitudes in the main oscillation. Corrections for the main oscillation and other harmonics are found as the forced solutions of the second, third and the following equations of the system: i.e., the amplitudes of the corrections are found from linear algebraic equations. The harmonic content of each consecutive correction is fully determined by the terms neglected in the solution of the previous equation. To improve the efficiency of correction calculations a systematic method for derivation of the forcing terms in the differential equations is given. As a result a previous correction is calculated with maximum precision and the subsequent correction in many cases does not even include the harmonics of the previous correction but involves only new harmonics. At first glance, the proposed method requires the notorious a priori knowledge of the harmonic content of the main oscillation. The approach, however, has been designed to increase the efficiency of calculations when the main oscillation includes the first harmonic only. Thus, the proposed method, which is a further refinement (a better precision of calculation is achieved) of a previous approach [6,7], is in the same category as many other methods [8,9] in which only one harmonic is assumed in the initial approximation. Yet the efficiency of correction calculations in the proposed method is such that the first two corrections are, in many cases, sufficient to obtain a solution.
2. EXPANSION
We consider
the autonomous
OF THE
systems
DIFFERENTIAL
which are described
a,w”x’“‘+a,,_,w”-lx’“-I’+. Here x(I) = d’x/dr’,
where
EQUATION
by the differential
. ~+a”x+K(w)f=O.
T = or is the normalized
time. We assume
f=f[x, x(I), XC>), . . x(P’] .)
equation (1)
that
(2)
is a function which can be represented by a Taylor’s type series expansion when the variable x is represented as a sum of successively diminishing terms (functions of a polynomial type, for example, satisfy this condition) and K(w) is a polynomial. The periodic solution of equation (1) is obtained as a sum x=x,+x,+x,+.
. .)
(3)
is called the main oscillation and where the term x0, which includes a few harmonics, the terms x, , x2,. . . are corrections. Each subsequent correction may include the same harmonics as the previous ones and additional harmonics of other frequencies. The frequency w is obtained as a sum w=6Jo+w,+w~+...,
(4)
of the approximate value w,,. where w,, w2,... are small corrections The approximate periodic solution of equation (1) is expressed as x,, = A,,, cos 7 + C (Ak,, cos k~+ &, sin kr), ktEi,
(5)
AUTONOMOUS
NON-LINEAR
SYSTEMS
37
where the set E0 identifies the harmonics included in x0, which are different from the first one (thus, k # 1). If one substitutes equation (5) into equation (2) one obtains the result
.hl=f[.%, 4’),XC’.,. . , xl?1 =hm + uil -.I& 1,
(6)
where f;,m = C,, cos T+ D,,, sin 7 + C (C,, cos kT+ Dko sin kT) kEE,,
(7)
includes the same harmonics as the main oscillation and the function (fo-A,,,,) includes the harmonics of other frequencies. Now one can rewrite equation (1) as a,o”x’“‘+
a,_,W*-‘X(n-‘)+.
. ~+~ox+K(~l.hm
=-K(w)[(f,-f,,)+(f-f,,)l.
(8)
Assuming that the higher harmonics with large amplitudes are included in the second term of equation (5) one can consider that the terms in the brackets, on the right side of equation (8), represent small quantities. Following reference [lo], one can consider, instead of equation (8), a more general equation,
. *+a,x+~(wIh,
a,w”x’“‘+a,_,w”~‘x’“-“+.
= -K(W)[&(fO-fOm)+(f-f0)1, (9)
with O< F G 1. When E = 1 equation (9) becomes the same as equation (8) (and, thus the same as the initial equation (1)). It is natural in this case to seek a solution of equation (9) as a power series x=x,f&X,+E*Xz+~**.
(10)
One can substitute this series into equation (2) and obtain the expansion f=fo+&f+&Ej2+...,
(11)
where f0 depends on x0 and its derivatives only, f, depends on x,, x, and their derivatives only, f2depends on x0, x, , x2 and their derivatives only, etc. To obtain precision of calculation better than in references [6,7] the procedure which was done for x0 is repeated for x, . As was mentioned before, x, is a sum of harmonics as well: i.e., x, = C (Ak, cos kT+ Bk, sin k7), keE,
(12)
where the set E, identifies the harmonics included in x,. If this sum is substituted into
f,one obtains (assuming that x0 is known at this stage) that f, =f,(xo,
x,) =f,m -c(f,
-f,m)t
(13)
where f,,,, =
c
htE,
(CA, cos kT+Dk,
sin kT),
(14)
with C,, and Ok, being linear functions of Ak, and BI, (they can be, of course, non-linear functions of the amplitudes of harmonics which are included in x0). The residual term f, -f,m includes the harmonics with numbers k F?E, or the harmonics with k E E, but with amplitudes which are non-linear functions of AL, and Bk,. The same calculation is carried out for all the f;: i.e., for each f; one finds the form
f;(X”,
XI
7..
‘,x,)=.L+(f;-.L)
(15)
38
S. S. QIU
as soon as the identifying
AND
I. M. FILANOVSKY
set E, of x, = 1 ( Ak, cos kT + B,, sin k7) A\~ F
is shown. Equation (13) and all the equations (15) are substituted in their direct form but in the form of
./-I=f,,n +&(.A-f,mL
(16)
into equation
f2=fim+4fi-.Ll),
(1 l), except not
...> (17)
f;=.L+e(_L-.L), where
.s is the previously
introduced
....
parameter.
One obtains
that
One can now substitute o=w,+&W,+& into K (w ) and obtain
wz+.
..
(19)
the expansion
K(w) = K(w,)+FK(w~,)w,+F~K~(w~,
w,)w2+.
9..
(20)
Finally, one can substitute equations (10) and (18)-(20) into equation (9) and equate the terms which have the same exponents of E. Then the following system of equations is obtained u,w~x~~‘+a,
,w;;-‘xl;‘~ ‘I+. . ~+a,,x,,+K(w,,)fc
a,w~x(ln’+a,,~,w;l~‘.~~*~‘)
+. .
,,,, =o,
.+ua,,x,+K,(w,,)w,f,,,,+K(wol.fi,,,
=w,P,(x,,w,,)-K(w,,)(.f,,-f;,,), a,WgnX:n’+a,_,Wgn-‘X~,l--“+~. . + a,+,+K,(w,,, w,h.f;~m+ K,(wohfim + K(wol_L
= wJ’z(xo, x,, w,) - K,(w,b,(h,-An)
- K(wo)(f, -h,n),
. ..>
(21)
where P, and PI are polynomials of their arguments (they are obtained from the terms additional to the terms of the left parts in equation (21) when the expansions (10) and (19) are substituted in the linear part of equation (3)). Equations (21) should be solved consecutively to obtain x,,, x,, x2,. . . and w,,, w,, w2,. . . . . The second equation includes the additional term K ( wo).f,mr the third the additional term K,(w,,)w,f,, + K (w,)frm, etc. These terms did not appear in the approach used earlier [6,7], but they essentially increase the precision of calculations of the first few harmonics. It is not difficult to see that to return from equation (9) to equation (8) one has to put E = 1. The same can be done in equations (17), but this means that one can solve equations (21) and find the required solution, using equation (3) (and the required frequency from equation (4)), by summation of the solutions obtained. The first equation in system (21) is the usual equation of the harmonic balance method, the second and other consecutive equations are linear differential equations. The solution of the second and other consecutive equations does not present any difficulties because one has to find only their steady state solutions: i.e., one uses the harmonic balance method again and solves the linear algebraic equations to find the amplitudes in x,, w,, w?, . . . x2,. . . and the frequency corrections
AUTONOMOUS
3. SIMPLIFlCATION
NON-LINEAR
OF THE LINEARIZED
39
SYSTEMS
EQUATIONS
The following steps for simplification of the system of calculations. (1) If the final form of the solution is set as
AND CALCULATION
OF f;,,
(21) can be used from the onset
x=A,cos7+A2cos2T+Brsin27+*..,
(22)
then A, =A,o+A,,+A,,7+.
*. ,
n=l,2,3,...
i.e., B,i=O (i=1,2,3 ,... ). (2) In conservative systems the harmonic in this case, can have the form
(23)
. . 1.’
n = 2,3,4,.
1 B,=B,o+B,,+B,z+...,
phases are independent
[5,8] and the solution,
~=A,,+A,cos7+A~cos2~+~~~.
(24)
(3) In conservative systems and systems with odd non-linearities w, = 0. (4) In systems with odd non-linearities A0 = 0 and A,, = 0. Additionally, the calculation of the terms l;,n in equations (21) can be carried out in a systematic way. If one considers that in equation (2) x,, .x2,. . . are small perturbations of x0 one can find that the perturbation Sf; off caused by the perturbation x, is
- af 4’l \ :\,\*, a/;=$.=x,+fd[X’
14)
\
(25)
9
q=l
.?I
where summation embraces all the derivatives included in equation (2). It is not difficult to see that Sf; can be obtained by the substitution of x, instead of x, in f,.Then if x, = 1 (Ak, cos kT+ BL, sin kT), ACE, the substitution
of equation
(26) into equation
(26)
(25) will give (27)
where 8kf; is the total contribution of the kth harmonic content of x,, is determined one can find that
fikf; = 1 ( Cki COSk7 i- Dki Combining harmonics
the parts of contributions included in x, one obtains
Sill
k7) + c
included
( ck,
COS
in Xi. When the harmonic
kT + &
from all harmonics having the required expression
Sill k7).
(28)
the frequencies
J,,, = 1 1 (C,; cos k7+ AL, sin kT) hk t, 1.. I;, Thus,
all the functions
f;“? are calculated
uniformly,
of
(29)
from the same function
f,.
4. EXAMPLE Consider equation
the oscillations
in a non-linear
conservative
“~x’~‘+x+cu~x’+a,x~=o.
system described
by the differential
(30)
40
S. S. QIU
AND
I. M.
FILANOVSKY
This is the same example which forced some authors [4] to abandon the harmonic balance techniques and was used for verification of the modification of the harmonic balance method in reference [5]. The initial conditions for the system considered are taken as the following: i-1,3,....
x,(O)= A,,xi(O)+x,+,(O)=O,
(31)
For conservative systems w, = 0. As a result the system (21) is simplified and becomes o;x:“+x,,+f”,
=o,
w;x:“+xz+fim
wW+x,+“L
= -(fo-fom), (32)
= -2x)i’w,,wz - (f, --j-i.,),
. . *.
For the example considered f= (u2x2+ (Y3x3
(33)
and to stress the possibilities of the proposed method we will carry out the calculations assuming only one harmonic in the main oscillation: i.e., we set x,, = A,” cos T.
(34)
If one substitutes equation (34) into equation (33) one obtains f0 = cxZA&cos* 7 + a3A;,, cos3 T = &,A:,
+ acy3A:, cos T+ fa2Af0 cos 2r + acu3Aiocos 37,
(35) so that j& = &,A;,
cos T
(36)
and f0 -f&, = $x,A$ + fa,A:,
cos 27 + &,A;,
cos 3 7.
(37)
Substituting equation (36) in the first equation of system (32) one obtains -w;A~,,+A,,+&~A~~=O.
Following the initial conditions
(38)
(31) A,, = A, and in this case one finds that w0= 41 +&A:.
(39)
The substitution of x = x0+ FX,+ e2x7. . . into f given by equation (33) gives f, = 2azx,x, + 3a,x;x,
)
(40)
and the substitution of x, instead of x, in this result gives sf; = (2~ax,~+3cz3x;)x,.
(41)
This is a conservative system and for it one can choose x,=
1 Ak, cos
hi E,
kr.
(42)
The substitution of the obtained x0 = A, cos T and xi given by equation (42) into equation (41) results in S,f; = (2azA, cos 7+3a,Af
cos’ T)A,, cos kT
=~(YJA~&COS(~-~)T+(Y~A,A~;COS(~-~)T
+;cu,AfAk,
cos kT+ + cuIA,Ak, cos (k+ 1)T
+&AfAki
cos (k +
2)T,
(43)
AUTONOMOUS
NON-LINEAR
SYSTEMS
41
and this simple trigonometric expansion is the key formula of the proposed method. Equation (37) gives the right side for the second equation in the system (32). As a result one has to set x, =A,,
+Az,
cos 27+A3,
cos 37.
(44)
E, = {0,2,3}
and using equation (43) (i.e., one substitutes i = 1 and then, for each the terms which give the contribution into the zeroth, or into the first, or into the third harmonic and then combines these terms together) one obtains
Thus,
k = 0, 1, 3 one selects
fin, =$,A:(A,,
+:A~,)+[&A:(A~,
+(cx2A,A2,+$r3A:A3,)
+A~,)+(Y~A,A~,]
cos 27
cos 37.
(45)
Further calculations will show that AZ, << AO, and for a more simple representation final results we will use _f ,m =&A:A~,
+ (&A:A~,
Substituting equation one can find that
+ a2~,~3,)
(46) and equation
cos 2T+(a,A,~,,
(37) into the second
+&A;)
equation
cos 3T.
of the
(46)
of the system (32)
where A = 24 +,A;
- &A;.
(48)
To solve the third equation of the system (32) we, first of all, will find f, . The substitution of x, given by equation (44) into equation (40) gives f, =;a3Af(Ao,
+A2,)
+ (2azA,A,,,
+ a>A,A,,
+&r3A:A3,)
cos T
+[~~,A:(A~,,+A~,)+(Y~A,A~,]COS~T+(((Y~A,A~,+~~~A:A~,)COS~T + (ct2A,A3, +:(Y~A;A~,)
and the subtraction
of equation
f, -f,,,,
cos
4~+~cqA~x‘t3,
cos
(45) from equation
= (2qA,A,,, +((Y,A,A,,
(49)
(49) gives
+ (Y:A,A~, +:cx3A:A3,) +:qA,AZ,)
5T,
cos
cos 7
~T+&,A;A~,
cos
57.
(50)
Ideally, the calculation (49) should not be done at all because the terms of the right side of equation (50) represent the terms discarded in the calculation ofJim (just as the terms discarded in the calculation off”,,, in equations (35)-(37) werefo-f&,,); we have to retain them only. Now one can conclude that x,=A,~cosT+A~~cos~T+A~~cos~T.
(51)
We will not carry out the calculation of the fourth and fifth harmonics (even though it is easy to do here and not so easy a task when using other methods) and instead will restrict ourselves to the calculation of the second harmonic and the frequency correction. Thus, we set E2 = (2) and, using equation (43) with k = i = 2, one can find that f?,, = &A:A,Z
COS
T+ .
+ . .
(52)
S.S.
42 Substitution equation
of equations
QIU
AND
I.
M.
FILANOVSKY
(50) and (52) in the third equation
&,A:A,~= fw,,w,A, - (2a:A,A,,, From the initial
+ a?A,A,,
of the system
+&,A;A,,).
(53)
(31) for i = 1 one obtains
conditions
A ,2=-(Ao,+Ax+A3,), and substitution
of equation
(54) into equation
w2 = (1/2w,)[(2a,-~cu,~,)A,, A brief analysis
of the results
obtained
(54)
(53) gives
+~cQ-$c+~,~A~,
-&+,A~,].
can be done by assuming
allows
(55)
that
0<(~~A:
(32) gives the
(56)
one to use the approximation
&A;
oo=l+
(57)
and to compare the results with results obtained by other methods. are of small amplitude in a Case (a): ff2” (Yeand A,<<1. In this case the oscillations system with equal influence of the non-linear terms, and one is close to the left side of the inequality (56). Thus, = -;a,A:+A,
x=x,+x,+x,
COS
d&+4;
COS27f’
. .
(58)
and w=oo+oz
~(Y~A~-~~A~.
=l+
(59)
Case (b): CY~<< 1 but (Y~> 1 and a,Af= 1. Here the oscillations are of a larger amplitude in a system with asymmetric influence of the non-linearity terms, and one is close to the right side of the inequality (56). In this case X=
-&,A:+A,
cos
7+&~3~:COS37
(60)
and w-l+
&qAf-~+x~A;.
(61)
The results (58) and (59) were obtained by other methods as well [4,5] but their modification into (60) and (61) caused by the change of system parameters and the oscillation amplitude is not so easily traced in other methods. To conclude this section we would like to stress that the corrections of the harmonics which have just been found will be corrected again not in the next equation but in the equation which follows it. This means that the calculated results have good precision.
5. DISCUSSION
The proposed method should be applied, as many other methods, in conditions when it is known that a periodic solution exists. The usual approach is to solve the first equation of the system (21) and to find the approximate solution x,). If this solution exists (i.e., one is able to find real values for A,,,", l& and wO) one assumes that the exact periodic solution exists as well. Even though this assumption requires an additional mathematical proof, it is widely used in many other methods of calculations; for example, in the describing function method [ 1,2]. As in many other methods, the calculation of the solution should be followed by verification of its stability [ 111.
AUTONOMOUS
NON-LINEAR
SYSTEMS
43
The question of the spectral content of the main oscillation should not be considered as a major defect of the harmonic balance method. Of course, the modifications providing efficient calculations of harmonics are desirable. All other methods should be considered from this point of view as well, and we think that there is no universal method which can predict all possible situations in non-linear systems. The efficient calculations of harmonics can help to discover a situation that was overlooked earlier. The correct content of the main oscillations is a desirable feature improving the convergence of the calculations. The answer is sometimes determined by the problem itself (for example, in the calculation of subharmonic oscillations the main oscillation will include a subharmonic and the oscillation with the frequency of the external force), some help can be obtained from the investigation of resonance properties of the linearized second equation in the system (21). But in any case the convergence of equation (9) and (in the limit) of its initial form (8) was proved [ 121 for the mentioned limitations of xi. In calculation practice it is reasonable to check that Si + 0, where S;= C
(A:,+@,)
(62)
k
for the correction xi. Besides, it is convenient to compare S, with Si_, to reach a conclusion about the smallness of a correction. 6. CONCLUSIONS The features of the approach presented can be summarized the following way. 1. The proposed modification of the harmonic balance method allows the efficient calculation of a small amount (up to ten) harmonics in the oscillation. This is usually enough for many engineering problems. 2. The method can be extended for systems described by differential equations of order higher than the second. In such cases, one obtains only more extended systems of linear algebraic equations. 3. The method requires that the non-linearity function be an analytical function of the variable and its derivatives. Even though many systems can be reduced to this case, this is an essential limitation of the proposed approach. 4. The method is especially efficient when the main oscillation x0 includes one harmonic only and the corrections x, and x2 are sufficient to obtain the solution.
ACKNOWLEDGMENT This work was supported in part by the National Science Foundation under Grant 6872011.
of China (PRC)
REFERENCES 1. D. P. ATHERTON 1982 Nonlinear Control Engineering. New York: Van Nostrand. 2. D. P. ATHERTON and H. T. DORRAH 1980 International Journal of Control 31(6), 1041-1105. A survey on nonlinear oscillations. 3. C. HAYASHI 1964 Nonlinear Oscillations in Physical Systems. New York: McGraw-Hill. 4. A. H. NAYFEH and D. T. MOOK 1979 Nonlinear Oscillations. New York: John Wiley. Vibrations 95,525-530. An intrinsic 5. A. S. ATADAN and K. HUSEYIN 1984 JournalofSoundand method of harmonic analysis for nonlinear oscillations (a perturbation technique). 6. SHUI-SHENG QUI and I. M. FILANOVSKY 1987 IEEE Transactions on Circuits and Systems, CAS-34, 913-918. Periodic solutions of Van der Pol equation with moderate values of damping coefficient.
44
S. S. QIU AND
I.M. FILANOVSKY
7. SHUI-SHENG QIU, I.M. FILANOVSKY and K. A. STROMSMOE 1987 ControL77wory and Advanced Technology 3(3), 189-195. On one combination of the harmonic balance method and perturbation techniques. 8. N. MINORSKY 1962 Nonlinear Oscillafions. Princeton, N.J.: .Van Nostrand. 1952 Self-oscillating Systems. Moscow: Gostechizdat. (In Russian). 9. K. F. THEODORCHIK 1954 Oscillations. Moscow: Gostechizdat. (In Russian). 10. 8. V. BULGAKOV 11. D. P. ATHERTON 1981 Stability of Nonlinear Systems. New York: John Wiley. 1942 Applied Mathematics and Mechanics (Prikladnaia Matematika i 12. B. V. BULGAKOV Mechanika, in Russian) 6(4), 263-280. On the application of Poincare’s method to free pseudolinear oscillatory systems.