- Email: [email protected]
An approximate analysis of non-linear, non-conservative systems using orthogonal polynomials
Journal of Sound and Vibration (1973) 29(4), 463--474
AN APPROXIMATE NON-CONSERVATIVE
ANALYSIS OF NON-LINEAR, SYSTEMS USING
ORTHOGONAL
POLYNOMIALS G. L. ANDERSON Bdndt Weapons Laboratory, Watervliet Arsenal, Watervliet, New York 12189, U.S.A.
(Received 10 February 1973, and fll revised form 10 April 1973) Approximate solutions of second-order non-linear ordinary differential equations containing a linear velocity dependent damping term and a small parameter of non-linearity are obtained with the aid of sets of orthogonal polynomials. The method of solution treats the effect of linear damping in such a way that the approximate solutions of the non-linear equations reduce to the exact solution of the corresponding linear problems as the small parameter of non-linearity tends to zero. Ultraspherical polynomials, which depend upon a parameter 2, the value of which is selected at the time of performing numerical calculations, are found to be efficacious for solving problems of engineering significance. Numerical calculations, based upon the approximate method for several values of 2 and upon the fourth-order Runge-Kutta method of numerical integration, for three non-linear equations, were performed and results were compared. It was found that, in these cases, the calculations effected with 2 = 3/2 most nearly agreed with the Runge-Kutta calculations.
1. INTRODUCTION In a recent investigation, Sinha and Srinivasan [1] developed a method for obtaining approximate solutions of the general second-order differential equation for a non-linear autonomous system 5~+ x +~f(x,2)=O, 5~=d2x/dt 2, (1)
subject to the initial conditions x(0) = Ao,
~(0) = 0,
(2)
where e is a small non-linearity factor and A-o is the initial displacement of the system. Their method, based upon the expansions of certain functions in terms of ultraspherical polynomials, represents relative to the well-known method of Krylov and Bogoliubov [2] an improvement in accuracy in determining the response of such systems. However, in the event that equation (1) can be written as 5~+ 2~.~ + x + eg (x, 2) = 0,
(3)
which include~ explicitly the effect of linear damping in the term 2~5c, it is well known that the Krylov-Bogoliubov method and, consequently, the Sinha-Srinivasan method when applied to equation (3) do not reduce to the limiting case of the solution of the corresponding linear problem, namely, + 2 ~ + :~ = 0, (4) as the non-linearity parameter e tends to zero, and hence, according to Brunelle [3], these methods do not properly account for the phase lag introduced by the linear damping term ;~. 463
464
o.L. ANDERSON
In reference [3] the author presented a modification of the Krylov-Bogoliubov method which does assure that the approximate solutions of equation (3) do reduce to the correct limiting cases as e tends to zero. Furthermore, it was shown by means of examples and graphical results that this modified procedure led to approximate solutions that were more accurate than those obtained through the Krylov-Bogoliubov scheme. Murty [4] considered nonlinear second-order differential equations including linear damping in still another way. Several authors [5-15] have used various sets of orthogonal polynomials for the purpose of obtaining approximate solutions to several non-linear oscillation problems. With the exception of reference [15], these studies have been considered only with steady-state oscillations and the direct, equivalent linearization of the non-linear differential equation of motion. In the present study, the objective consists of employing orthogonal polynomials, in the spirit of reference [1], to determine approximately the transient response of non-linear systems in which the effect of linear damping is included. The approximate solution is constructed, as in reference [3], in such a manner that it reduces to the proper solution of the corresponding linear problem as e approaches zero. Moreover, if), = 0, the method presented herein reduces to that discussed in reference [I ]. 2. METHOD OF APPROXIMATION Consider the non-linear ordinary differential equation
5?(t) + 2)'.~(t) + 12Zx(t) + ef ( t,x, YO = 0
(5)
subject to the initial conditions in equation (2), where 12 is a constant and e and 7 are "small" parameters. If e = 0, equation (5) reduces to the linear equation in equation (4), which has a solution of the form ~(t) = e "f, where m must satisfy the auxiliary equation m 2
+ 2ym + 0 2 = O.
But it follows that m = - y +/co,
o9 = (O" -- y2)t/2,
i = (--1) 1/2.
(6)
Throughout this discussion, it is assumed that y < f2. Thus, by virtue of equation (6), the solution of equation (4) may be represented as
s
= ,,l(t)cos ~,
(7)
where .4(t) = A l e - ' ,
~ = cot + ~bl,
with Aa and ~bl denoting constants of integration. Also, for this representation, .~(t) --=-t/(t)()' cos ~ + o9sin ~).
(8)
Ife # 0, the solution of equation (5) can be expected to be not too different from equation (7). Consequently, one is motivated to assume that
x(t) = A(t)cos ~p,
(9)
~9(t) = tot + ~b(t),
(10)
where and the functions A(t) and q~(t) are now to be determined. Upon differentiation of equation (9) one obtains ~(t) = / / c o s Ip - (co + q~)A sin ft.
(1 I)
NON-LINEAR, NON-CONSERVATIVE SYSTEMS
465
If next the identity ?A cos ~k - ~A cos ~k= 0 is inserted into equation (11), then the equation :~(t) ---A(Tcos ~ + tn sin ~) + [(.4 + TA)cos ~k-- ~A sinG]
(12)
is obtained. The first set of terms in equation (12) is identical in form to the right side of equation (8). Therefore, one can require that (.4 + ~A) cos ~ - q~Asin ~k= 0,
(13)
which is, in essence, merely an application of the method of variation of parameters. Consequently, equation (12) reduces to :~(t) = --d (},cos ~k + o9sin tp).
(14)
Furthermore, differentiation of equation (14) and use of equation (13) yield dr(t) = -o9(.4 sin ~ + ~A cos ~) + A [~o9sin r + (2~2 - f22) cos ~].
(15)
Finally, substituting equations (9), (14), and (15) into equation (5) gives (.4 + ~A) sin ~, + Ag~eos ~k= (e/og)f(t, A e o s ~ , - ? A c o s r
(16)
Equations (13) and (16), namely,
(.4 + yA) cos ~ - q~a sin ~ = 0, (.4 + yA) sin ff + q~A cos ~k= (e/og)f(t, A cos ~k,-yA cos (t - ogA sin if), represent a pair of coupled equations in the unknowns .4 + yA and Aq~. The solution of this system for .4 + yA and ~ is
.4 + yA = (e/o~)f(t, A cos ~,-~A cos ~k- ogA sin ~k) sin ~p,
(17)
= (e/ogA)f(t, A cos ~k,-?A cos @- ogA sin @) cos ft.
(18)
To obtain the associated initial conditions, equations (9) and (14) are inserted into equation (2), to give A(0) cos ~b(0) = Ao,
?A(0) cos ~b(0) + ogA(0) sin qS(0) = 0,
whence
A(0) = aAo/o9,
cos qb(o) = o91f2,
sin qs(o) = '~/s'2.
(19)
The right sides of equations (17) and (18) are clearly periodic functions in r with period 2re, and the quantities .4 + ?A and q~ are directly proportional to e, a small parameter, so it is evident that these quantities are slowly varying functions of time. In reference [3], equations (17) and (18) were simplified by expanding the right sides as Fourier series in ~ over the interval 0 < ~ < 2n. Because .4 + ~A and q~ are slowly varying functions oftime, the Fourier coefficients were treated as being approximately constant over the interval 0 < ~k < 2rr, and the reduced equations were found to be
d + 7A -~ (e/o9)Ko(t,A),
~ = (e/coA)Po(t,A),
(20)
where 2~
Ko(t, A) = (l[2n) f f(t, A cos ff,-),A cos ~k-- OgAsin ~9)sin ~d~,,
(21)
o 2I
Po(t, A) = (l/2rc) f f(t, A cos ~ , - ? A cos Ip - o9A sin if) cos Ipd~,. o
(22)
466
G.'L. ANDERSON
It is the objective here to proceed in a somewhat more general fashion. Instead of the right sides of equations (17) and (18) being expanded in Fourier series, they will be expanded in terms of Jacobi polynomials, P~,"P)(z). The ultraspherical (Gegenbauer), Legendre, and Chebyshev polynomials are special cases of Jacobi polynomials. In order to expedite subsequent calculations, it is convenient to introduce the change of variable z = ( r - rr)/~,
- 1 < z < 1,
and the notation
F~(t, A, ~) = f ( t , A cos ~k, - y A cos ~k- toA sin ~) sin ~, F2(t, A, tp) = f ( t , A cos ip,--yA cos ~p-- o~A sin ~) cos Ip. Next the definitions
Fk(t, A, ~) = Fk[t, A, rr(z + 1)1 -- Gk(t, A, z),
k = 1,2,
are introduced, where
Gl(t, A, z) = --f(t,--A cos nz, vA cos rrz + toA sin ~rz)sin ~rz,
(23)
G2(t, A, z) = - f ( t , - - A cos nz, yA cos nz + toA sin rcz) cos r~z.
(24)
Hence, insertion of equations (23) and (24) into equations (17) and (18) yields
A + rA = (~1o~)G,(t, A, z),
dp = (e/toA) G2(t, A, z).
(25)
The expansion of Gk(t,A, z) in terms of Jacobi polynomials shall be denoted by
Gk(t, A, z) = ~ a ~' c~" k,nB~tt \ , A~P I n aJt..~ k~"l,
(26)
where Re (co), Re (fl) > - 1 . Since the Jacobi polynomials form an orthogonal set of functions with respect to the weight function (1 - z ) ' ( 1 + z) B on the interval - 1 < z < 1, it is easily demonstrated that 1
f Gk(t,A,z)p~.t~(z)(1 - z)~(l + z)Bdz a~.'f'(t,A) = -1 1
(27)
f [P~," a'(z)lZ(1 -- z)'(1 + z) adz --1
where, according to Rainville [16, p. 260], 2x+~+~F(I + cc+ n) F(I + fl + n) -1 [P~'B~(z)]2(l -- z)~(l + z)Bdz = n!(1 + o c + f l + 2 n ) F ( 1 + ~ + f l + n ) '
(28)
where F(...) denotes the gamma function. If the expansions in equation (26) are restricted to the first term only,
Gk(t, A, z) ~_ a~,'oB'(t, A) =- gk(t, A; ~, fl)
(29)
since P~o~' B)(z) = 1. But using equations (27) and (28) with n = 0 gives 1
gk(t,
F(2 + C~ + /~) | a~(t,A, z)(1 - z)'(1 + z ) B d z . d ; cq fl) = 2z+~+B F(I + ~) F(I + ,6) -13
(30)
N O N - L I N E A R , N O N - C O N S E R V A T I V E SYSTEMS
467
Therefore, in view of equation (29), one can approximate equation (25) by
X + ?A = (e/og)g~(t,A;~z, fl),
$ = (e/oJA)g,(t,A;o~,fl),
(31)
where the functions gk(t, A ;cz,fl), k -- I, 2, are obtained from equation (30). Equations (31) are the counterparts of equations (20), and indeed it will be shown later that the latter are merely a special case of the more general form given in equations (31). The procedure outlined here implies a more general averaging technique because the free choice of the parameters 0~and fl is still available, subject only to the conditions Re (c0 > - 1 , Re (fl) > -1. The merit of the present technique is that for suitable choices of ~ and fl the solutions of equations (31) will lead to somewhat more accurate approximations than do those of equations (20). It frequently happens that the integration of the first-order differential equations in equations (31) is a laborious process, but, of course, the degree of complication is dependent upon the form off(t, x, :~) in equation (5) and upon the choices of the parameters cc and fl in equation (30). However, for a suitable choice of c~ and fl, namely, c~= f l - - - 2 - 1/2, where Re (2) > - 1 / 2 , which implies only a modest sacrifice of generality, more efficacious forms of equation (31), at least for certain practical problems, can be obtained. In this case the Jacobi polynomials P(nx-ll2" ;t-ll2)(z), which can be obtained from (2 + I/2)n C~(z) (2).), '
p~x-l/2, ~-1/2)(z ) =
where Ca(z) denotes the Gegenbauer polynomial of degree n and (...), the factorial function [16, p. 22], are known as ultraspherical polynomials. Thus, equation (30) becomes 1
F(22 + 1) f g,(t, A; 2 - 1/2, ;. - 1/2) = 22x[r(;. + l/2)]" ,J Ga(t, A, z)(l - z2) x-xn dz.
(32)
--1
However, by using Legendre's duplication formula ([16], p. 24) in the form F(2). + 1) =
22x F(2 + 1/2) F(2 + 1) ~112
equation (32) can be expressed as follows: 1
gk(t,A;2-- 1/2,2-- I/2) --=h,(t,A;2)= r d ntO. F 0 . ++ I)I/2)_f"t Gk(t,A,z)(1 -z2)X-ll2dz.
(33)
Consequently, equation (31) can be replaced by
d + TA = (e/o) ht(t, A ; 2),
$ = (e/o~A) h2(t, A; 2),
(34)
which are somewhat more general than are the corresponding expressions reported by Sinha and Srinivasan [1]. Furthermore, if). = 1/2, then equation (33) reduces to 1
hk(t,A;l/2)=~
G(t,A,z)dz --! 2n
0
Ko(t, A) Po(t, A)
if k = 1, if k = 2,
(35)
468
O.L. ANDERSON
upon comparing equations (35) and (21)-(22). Therefore, the procedure proposed by Brunelle [3] to include explicitly the effect o f linear damping in the amplitude and phase differential equations is a special case o f equation (31) with c~= ,8 = 0 and, of course, equation (34) with ). = 1/2. In applications o f equation (33) to specific problems, one often encounters the definite integrals 1
-J sin (nnz)(l -- z2)a-,/2dz = O,
(36)
--1
i
_ cos (nrcz)(l
-- z2) 1-1/2
dz =
nt/2F(2 + 1/2) F(2 + 1) r/~,
(37)
where tlo = 1,
r/~ =
F(2 + 1)j~(nn), (nx/2) x
n = 1 , 2 , 3 .....
(38)
Jx denoting the Bessel function of the first kind o f order 2. The integral in equation (37) is found in Gradshteyn and Ryzhik [17, p. 403, formula 21]. Selected numerical values oft/, are given in the Appendix. 3. AN EXAMPLE Next, as an example o f a p r o b i e m t o which the theory derived above can be applied, one may consider the non-linear differential equation 2 + (2~ + ~ l x 2) 2 + [[22 + tp(t) cos O(t)] x + e(~t x 3 + ct2x s) = O,
(39)
where/l(t) and 0(0 are functions o f t, which, for the present, shall remain unspecified, and )'1, cq, and ct2 are constants. Comparing equations (5) and (39) indicates that
f(t,X,5:)=~tx2yc + Iz(t)cosO(t)x +~tX3 + ~2xS.
(40)
Substituting equation (40) into equations (23) and (24) leads to
G~(t, A, z) = (A/2)/t cos 0 sin 27rz + (~1 -- Y~',)(A[2)3( 2 sin 27rz + sin 47rz) - co~1(A/2)3(1 - cos 4nz) + %(A/2)5(5 sin 27rz + 4 sin 47rz + sin 6nz),
(41)
G2(t, A, z) = (A/2) p cos 0(1 + cos 2nz) + (cq -- ~,)q) (A/2)3(3 + 4 cos 2rrz + cos 47rz) + + Ocz(A/2)s(10+ 15 cos 2nz + 6 cos 4zz + cos 6rrz) - mvI(A/2)3(2 sin 2nz + sin 4zrz).
(42)
Insertion o f equations (41) and (42) into equation (33) yields
hi(t, A ;2) = --~oyt(A/2)3(l -- r/,),
(43)
h2(t, A ; ).) = (A/2)/t cos 0(1 + r/2) + (r -- 7Yl)(A/2)3( 3 + 4t12 + 1/4) + + =2(Z/2)s(10 + 15r/2 + 6r/, + r/6),
(44)
where equations (36)-(38) have been used. Combining equations (43) and (44) with equation (34) yields the amplitude and phase differential equations: A + yA = --~?,1(A/2)3(1 -- n,),
(45)
= (e/2w)[(1 + r/2)/a(t)cos O(t) + (oq - ~,~h)(A/2)2(3 + 4r/2 + r/4) + + ~2(A/2)4(10 + 15r/2 + O/4 + r/6)].
(46)
NON-LINEAR, NON-CONSERVATIVE SYSTEMS
469
In view o f the initial conditions in equation (19), integration o f equation (45) leads to A(t) = f2Ao e-"[co 2 + B (12Ao)'(l -- e-2~')] -1/2,
(47)
where a = (e~q[8y)(1 - t/,). With the objective of integrating equation (46), this equation may now be written as
t~(t) --- ~(t) + r, A 2 + r2 A',
(48)
~(t) = (e]2co)(l + q2) It(t) cos O(t)
(49)
where and rl = (e/8w)(cq - ~yl)(3 + 4q2 + q4),
r2 = (e~2/32co)(10 + 151/2 + 6r/4 + r/6).
(50)
But, according to equation (45) and the chain rule of calculus, d dA d d ~ ( 4 , - 0 = d-'7d---J(4' - ~) = - r A ( 1 + BA 2) ~-j (4, - ~), and consequently equation (48) becomes d
(rt A + r2 A 3)
T~ (4, - O =
~(l + B a ~) "
(51)
The solution of equation (51) subject to the initial conditions in equation (19) is 4,(0 = cos-l(co/f2) + ~(t) - ~(0) + (r2/2~B)[(f2Ao/o)) 2 -- A 2] +
1 + 2~
[ I+BA z ] (r2 - rl B ) I n 1 +-B(---~oflO)i]"
(52)
However, if 7~ = 0, then B = 0, and equation (52) is not suitable for computer calculations. In this case equation (47) reduces to A( t ) = (~?Ao/w) e - ' ,
(53)
whereas equation (48) must be replaced by = ~ +r*a 2+r2A*,
where r* = rt [~=o- The counterpart of equation (51) is d d---A(4, - ~) = - - A ( r * + r2 A2)/~, whence 4,(0 = cos-l(co/O) + ~(t) - ~(0) + 7(f2Ao/27co)'(l - e-2")[2r] ' + rz(f2AofiO)2(l + e-Z')].
(54)
4. NUMERICAL RESULTS Linear expressions in t for the functions It(t) and O(t), namely, It(t) = Ito + It1 t,
O(t) = Oo + 01 t,
(55)
where Ito, Iti, 0o, and 01 are constants, are ofconsiderable importance in the theory of dynamic stability--see, for example, reference [18]. Inserting equation (55) into equation (49) and integrating the result gives ~(t) = (e/2coO~)(1 + r/2)~ o Ot sin (0o + 01 t) + It1 01 t sin (0o + Ot t) + Itt cos (0o + 01 t)],
470
G.L. ANDERSON
and consequently ~(t) -- ~(0) = (e/2to02)(1 + q2){Po 0~ [sin (0o + 0~ 0 - sin 0o] + tq 0~ t sin (0o + 0~ t) + + tq[cos (0o + 0~ t) - cos 0o]}. (56) Therefore, the Solution of 5/+ (2~ + ~),ax 2) ~ + [f22 + ~(Po + P~ t) cos (0o + 0~ t)] x + ~(cq x 3 + cq x 5) = 0
(57)
is, approximately,
x(t) = A(t)cos r
(58)
with, by virtue of equation (56), r2 2-~1 [1 I + B A 2 O(t) =- cot + cos-'(oJfl2) + 2-~-~[(QAo]co)2 --A 2] + (rz - rl B)ln +~w)2j
"]
+
+ (e/2co0~)(l + rh)(po 0 x [sin (0o + Ot t) -- sin 0o] + Pt 01 t sin (0o + 0~ t) + tq[cos (0o + 0~ t) - cos 0o]},
(59)
and
A(t) = I2Aoe-~t[to 2 + B(I2Ao)2(I -- e-2~')] -~/2,
(60)
provided that 71 :A 0, and
x(t) = (t2Ao/co) e-" cos ~b(t),
(61)
with ~,(t) = cot + cos-X(to/f2) + y(t2Ao/2),co)X(l -- e -2~') [2r* + r2(f2Ao/co)2(l + e-a")] + + (e/2co0~)(1 + q2){Po 01[sin (0o + 0, t) -- sin 0o] + lq 0, tsin (0o + 0x t) + + ,tq [cos (0o + 01 t) -- cos 0o]},
(62)
whenever ~,, = 0. At this point it should be remarked that the approximate solution of equation (57) as embodied in equations (59) and (60) does not reveal the possibility of the existence of unstable motions that are known to occur in certain systems whose motions are described by linear or non-linear differential equations of the Mathieu type (if F~ = 0, equation (57) falls in this category). This phenomenon, frequently called parametric resonance, has been discussed in some detail by Bolotin [19] and Bogoliubov and Mitropolsky [20]. In this study only a firstorder approximate solution was obtained; however, in reference [20] it is shown that the Krylov-Bogoliubov method predicts the possibility of resonance only in the second and higher order approximations. Consequently, the method of solution described herein can be used only in those situations in which the ratio of the frequency of the driving force and the natural frequency of the associated linear system is not a rational number, i.e., in non-resonance cases. In the event that the driving frequency is equal to or approximately equal to the natural frequency, a different method of analysis must be applied, and the interested reader is referred to pages 196-236 and 267-284 of reference [20]. Finally, it is worth mentioning that Huston and Doty [21] have applied the Krylov-Bogoliubov method in the first-order approximation to linear differential equations of the type 5/ + ~(t)5/ + ~ ( t ) x = o,
where ~(t) and ~j(t) are known functions, and they reported an error analysis which indicates the domain of validity of their method, which may be considered a special case of the technique described here.
NON-LINEAR, NON-CONSERVATIVESYSTEMS
471
The expression for the velocity :~(t) in equation (14) can be simplified somewhat. Set ~(t) = ~b(t) -- COS-l(to/12). Then, in view of equation (19), one can show that equation (14) assumes the form 2 ( 0 = - f 2 A ( t ) sin 7'(0.
Numerical calculations were performed with Ao=e=f2=cq
= 1,
ct2 = 0 ,
It(t) = 0,-0"02 (1 + t),
1,
yl =0,0"001
y = 0"05,
O(t) = 1.8t
and the values 2 = 0, I/2, 1, 3/2, 2, 5/2, 7/2, and 9/2. Specifically, the following three nonlinear differential equations, all possessing a linear damping term, were solved approximately by means of equations (58)-(62), as appropriate: 2 + O. I Ec+ x + xa + xS = O,
5( + 0"1 .~ + [I -- 0"02(1 + t ) c o s ( l . 8 t ) ] x + x 3 + x 5 = O, .~ + 0"1(1 + 0"01 x2).~ + [1 -- 0"02(1 + t)cos (l'8t)]x + x 3 = 0,
(63)
subject to the initial conditions x(O) = l,
.~(0) = O.
In order to judge the accuracy of the approximate formulas in equations (58)-(62), these same equations were also solved numerically with the fourth-order Runge-Kutta method. The variations o f x ( t ) as a function of time t for these three differential equations are plotted in Figures 1, 2, and 3 for only the first few cycles ofmotion. The value of the initial displacement Ao was selected to be a fairly large number in order to emphasize the effect of the value of the parameter 2. If 1, --- 0, the approximate expression for the amplitude is independent o f any of the coefficients ~I,, which depend upon 2, and if ~i :~ 0, the approximate expression for A ( t ) in equation (60) depends upon only 114,which, generally, is a small number relative to unity. For the three differential equations given in equation (63), the numerical computations and the graphs in Figures 1-3 reveal that the choice of 2 = 3/2 leads to somewhat more accurate values of the phase than did the choice of 2 = 1/2, which corresponds to the method of reference [3]. The approximations obtained for 2 = 1 and 2 were only slightly less accurate than those obtained for ). = 3/2, but for the other values of 2 that were tried the results were I.O
8
,
,
i
'
~
J
'
I
l
t
t
;
I
'
t
~
i
I
'
i
~
i
I
I
I
l
,
i
t
~
I
J
l-
/ E x a c t
05
\\
/;/I,/
\\
/ ./
\,\
v
',,\ ,;'/
-0.5
-I.0
i
,
i
I
~/2
,
i
t
J
I
T
i
,
T
i
I
i
i
i
3r/2
i
t
2~"
,
,
,
,
" l
5~/2
,
,
i
i
3t-
t
,
'
i
r 7r/:
t
Figure 1. Variation of x(t) with t for ~ + 0.12 + x + x 3 + x s = O, x(O)= 1, 2(0) = O.
472
G.L. ANDERSON I'O
,
~ i l i t l l l l ~ l l t l
\
,~
iILi
l l ~ l
j l l l l t
/,r ///
i L l l l
i
tl
iii8
ii
!1 'x\
//
\\
e
o
-oV V --IC
I l t l
III
I!11
~'/2
I
~"
111
I l t l
37r/2
I III
2zr
I I I I I I
5~/2
tit
t
3~"
III
?z~/2
l l l l l
4~"
II
9~'/2
t
Figure 2. Variation of x(t) with t for .~ + 0.1.~ + [1 - 0.02(I + t)cos(l.8t)]x + x a + x s = 0, x(0) = 1, ~(o)
=
o.
I'O
it
1111t
iii
l l l l l l l l l
~-
--I O
--
IIII
t t l l l l
z-/2
~r
l i l t
3~r/2
i l l l l l l l l l
i1111
i l l l l l
i i i i
~ExQct
iIII
iIi
i I i 1 1 t l
2~"
5r/2
i
37r
ii~
I , ~ , l ! ! t
71r/2
4r
pp
9rr/2
l
Figure 3. Variation of x(t) with t for $ + 0"l(l + 0'01x2)~ + [I - 0.02(1 + t)cos(l.8t)]x + x 3 = 0, x(0) = 1, ~(o) = 0.
n o t i c e a b l y less a c c u r a t e . Therefore, for an a p p o s i t e choice o f the p a r a m e t e r 2, the present m e t h o d a p p e a r s to be s o m e w h a t m o r e accurate t h a n the m e t h o d suggested in reference [3] a n d has the a d d i t i o n a l feature that as e tends to zero the s o l u t i o n o f e q u a t i o n (5) tends to the correct solution o f the c o r r e s p o n d i n g linear p r o b l e m , whereas the m e t h o d o f s o l u t i o n given in reference [1], as a p p l i e d to e q u a t i o n (5), does not.
REFERENCES 1. S. C. SINHA and P. SRINIVASAN1971 Journal of Sound and Vibration 18, 55-60. Application of ultraspherical polynomials to non-linear autonomous systems. 2. N. KRYLOV and N. BOGOLIUaOV 1947 Introduction to Non-Lhtear Mechanics. Princeton, New Jersey: Princeton University Press. See pp. 8-14. 3. E.J. BRONELLE1967 InternationalJournalofNon-Lhzear Mechanics 2, 404-415. Transient response of second-order non-linear equations. 4. I. N. S. MORAY 1971 lnternationalJournalofNon-Linear Mechanics 6, 45-53. A unified K r y l o v Bogoliubov method for solving second-order non-linear systems. 5. H . H . DENt,IAN 1959 AmericanJournalofPhysics 27, 524-525. Amplitude dependence of frequency in a linear approximation to the simple pendulum equation. 6. H. H. DENMAN 1964 Journal of the Industrial Mathematics Society 14, 9-20. Application of ultraspherieal polynomials to asymmetric non-linear oscillations.
NON-LINEAR, NON-CONSERVATIVESYSTEMS
473
7. H. H. DENMANand J. E. HOWARD 1964 Quarterly of Applied Mathematics 21,325-330. Application ofultraspherical polynomials to non-linear oscillations. I. Free oscillations of the pendulum. 8. H. H. DENMAN and Y. K. LIU 1965 Journal of the Industrial Mathematics Society 15, 19-35. Application of ultraspherical polynomials to non-linear systems with step function excitation. 9. R. M. GARDE 1965 Indian Journal of Mathematics 7, I 1 I - I 17. Application of Gegenbauer polynomials to non-linear oscillations--forced and free oscillations without damping. 10. H. H. DENMANand Y. K. LIU 1965 Quarterly of Applied Mathematics 22, 273-292. Application of ultraspherical polynomials to non-linear oscillations. II. Free oscillations. 11. R. K. SAXENAand K. S. KUSHAWA 1970 Proceedings ofthe National Academy of Sciences, India A40, 65-72. Application of Jacobi polynomials to non-linear oscillations. 12. R. M. GARDE 1967 Journal of Scientific and Engineering Research 11, 157-166. Application of Gegenbauer polynomials to non-linear damped oscillations. 13. Y. K. LIu 1967 Journal of Applied Mechanics 34, 223-226. Application of ultraspherical polynomials to nonlinear forced oscillations. 14. R. M. GARDE 1969 Proceedings of the National Academy of Sciences, India A39, 445-450. Application of Jacobi polynomials to some nonlinear oscillations. 15. S. C. SINHAand P. SRINIVASAN1972JournalofSoundand Vibration 22, 211-219. An approximate analysis of non-linear non-conservative systems subjected to step function excitation. 16. E. D. RAINVILLE1960 SpecialFunctions. New York: The Macmillan Company. 17. I. S. GRADSHTEYNand I. M. RYZHIK 1965 Tables of Integrals, Series, and Products. New York: Academic Press. 18. R. M. EVAN-IWANOWSKI,W. F. SANFORDand T. KEHAGIOGLOU1971 DecelopmentsofTheoretical and Applied Mechanics, Vol. 5, 714-743. Chapel Hill, North Carolina: University of North Carolina Press. Nonstationary parametric response of a nonlinear column. 19. V.V. BOLOrIN 1964 The Dynamic Stability of Elastic Systems. San Francisco: Holden-Day, Inc. See pp. 18-27 and 63-74. 20. N. N. BOGOLIUBOVand Y. A. MrrROPOLSKY 1961 Asymptotic Methods in the Theory of NonLhtear Oscillations. Delhi, India: Hindustan Publishing Corporation. See pp. 196-236 and 267-284. 21. R. L. HUSTON and L. T. DOTY 1966 SIAM Journal on Applied Mathematics 1,1, 424-428. Note on the Krylov-Bogoliubov method applied to linear differential equations.
APPENDIX I n e q u a t i o n (38), the coefficients q.(2) were defined as
r(;. §
(Al)
=
F o r 2 = 0, 1, 2, the first six values o f !1.(2) are given in T a b l e .I. 9 TABLE 1
Vahtes of q.(2) for some integer vahtes of 2 n 1 2 3 4 5 6
tl.(O) -0"304 0"220 -0"181 0"157 -0"141 0"129
242 18 276 91 211 46 507 39 182 05 063 52
t/.(1) 0"181 -0.676 0"375 -0"245 0"177 -0"135
191 76 034 59 022 53 943 43 012 25 185 03
t/.(2)
x x x x x
10 - l 10 -1 10 - l 10 -1 10 -x
0-393 477 -0"583 367 0"196 980 -0"922 538 0"515 143 -0"321 035
93 59 72 17 73 10
x x x • x
I0 -1 10 -1 10 -2 10 -2 10 -2
474
G . L . ANDERSON
If 2 = k + 1/2, k = 0, I, 2, ..., then the Bessel functions Jk+l/2 can be expressed in terms of elementary functions:
Jl/2(x)= \~x]
=
J3/2(x)=\~x]
sinx,
[ 2 \"2 ) [[ 3
I)sin x- 3cosx]
Jw2(x)=(~x) ',2[{156]sinx+( l-:Oc~ [2\ 1/2
45
105/sinx + (~_-- ~)
cosx].
In these cases, equation (A1) assumes rather simple forms: namely, 0.(1/2) = 0,
0.(7/2)
q.(3/2)
3 ( - 1 ) "§ (n~) 2 ,
105(-I)" (nr0-----~ [(n~) 2 - 15],
~.(5/2)
45(--1) "§ (n~) 4 ,
4725 th(9/2) = (-1)" ~ - ~ [2(n~) 2 - 21].
In Table 2, the first six values oft/.().) are given for 2 = k + 1[2, k = 0, 1, 2, 3, and 4. TABLE 2
Vahtes of lh(2)for certain vahtes of 2 n r/.(112)
I 0 2 3 4 5 6
0 0 0 0 0
r/.(3/2)
0"30396355 -0"75990888 • 10 -1 0"33773729 x 10 -I -0"18997722 x 10 -1 0"12158542x10 -1 -0"84434320x10 -2
r/.(5/2)
1/.(7/2)
t1.(912)
0"46196920
0"56032618
0"62783569
-0.28873075 x 10 -1 0"57033235 x 10 -2 -0"18045672 • 10 -2 0"73915072x10 -3 -0"35645772x10 -3
0"41772785 • 10 -t -0"11060492 x 10 -1 0"38106924 x 10 -2 -0"16198367x10 -2 0"79662120x10 -3
0"11273727 -0"118 89707 x 10 -I 0-22402199 x 10 -2 -0"60231881• -3 0.20445481x10 -3
AN APPROXIMATE NON-CONSERVATIVE
ANALYSIS OF NON-LINEAR, SYSTEMS USING
ORTHOGONAL
POLYNOMIALS G. L. ANDERSON Bdndt Weapons Laboratory, Watervliet Arsenal, Watervliet, New York 12189, U.S.A.
(Received 10 February 1973, and fll revised form 10 April 1973) Approximate solutions of second-order non-linear ordinary differential equations containing a linear velocity dependent damping term and a small parameter of non-linearity are obtained with the aid of sets of orthogonal polynomials. The method of solution treats the effect of linear damping in such a way that the approximate solutions of the non-linear equations reduce to the exact solution of the corresponding linear problems as the small parameter of non-linearity tends to zero. Ultraspherical polynomials, which depend upon a parameter 2, the value of which is selected at the time of performing numerical calculations, are found to be efficacious for solving problems of engineering significance. Numerical calculations, based upon the approximate method for several values of 2 and upon the fourth-order Runge-Kutta method of numerical integration, for three non-linear equations, were performed and results were compared. It was found that, in these cases, the calculations effected with 2 = 3/2 most nearly agreed with the Runge-Kutta calculations.
1. INTRODUCTION In a recent investigation, Sinha and Srinivasan [1] developed a method for obtaining approximate solutions of the general second-order differential equation for a non-linear autonomous system 5~+ x +~f(x,2)=O, 5~=d2x/dt 2, (1)
subject to the initial conditions x(0) = Ao,
~(0) = 0,
(2)
where e is a small non-linearity factor and A-o is the initial displacement of the system. Their method, based upon the expansions of certain functions in terms of ultraspherical polynomials, represents relative to the well-known method of Krylov and Bogoliubov [2] an improvement in accuracy in determining the response of such systems. However, in the event that equation (1) can be written as 5~+ 2~.~ + x + eg (x, 2) = 0,
(3)
which include~ explicitly the effect of linear damping in the term 2~5c, it is well known that the Krylov-Bogoliubov method and, consequently, the Sinha-Srinivasan method when applied to equation (3) do not reduce to the limiting case of the solution of the corresponding linear problem, namely, + 2 ~ + :~ = 0, (4) as the non-linearity parameter e tends to zero, and hence, according to Brunelle [3], these methods do not properly account for the phase lag introduced by the linear damping term ;~. 463
464
o.L. ANDERSON
In reference [3] the author presented a modification of the Krylov-Bogoliubov method which does assure that the approximate solutions of equation (3) do reduce to the correct limiting cases as e tends to zero. Furthermore, it was shown by means of examples and graphical results that this modified procedure led to approximate solutions that were more accurate than those obtained through the Krylov-Bogoliubov scheme. Murty [4] considered nonlinear second-order differential equations including linear damping in still another way. Several authors [5-15] have used various sets of orthogonal polynomials for the purpose of obtaining approximate solutions to several non-linear oscillation problems. With the exception of reference [15], these studies have been considered only with steady-state oscillations and the direct, equivalent linearization of the non-linear differential equation of motion. In the present study, the objective consists of employing orthogonal polynomials, in the spirit of reference [1], to determine approximately the transient response of non-linear systems in which the effect of linear damping is included. The approximate solution is constructed, as in reference [3], in such a manner that it reduces to the proper solution of the corresponding linear problem as e approaches zero. Moreover, if), = 0, the method presented herein reduces to that discussed in reference [I ]. 2. METHOD OF APPROXIMATION Consider the non-linear ordinary differential equation
5?(t) + 2)'.~(t) + 12Zx(t) + ef ( t,x, YO = 0
(5)
subject to the initial conditions in equation (2), where 12 is a constant and e and 7 are "small" parameters. If e = 0, equation (5) reduces to the linear equation in equation (4), which has a solution of the form ~(t) = e "f, where m must satisfy the auxiliary equation m 2
+ 2ym + 0 2 = O.
But it follows that m = - y +/co,
o9 = (O" -- y2)t/2,
i = (--1) 1/2.
(6)
Throughout this discussion, it is assumed that y < f2. Thus, by virtue of equation (6), the solution of equation (4) may be represented as
s
= ,,l(t)cos ~,
(7)
where .4(t) = A l e - ' ,
~ = cot + ~bl,
with Aa and ~bl denoting constants of integration. Also, for this representation, .~(t) --=-t/(t)()' cos ~ + o9sin ~).
(8)
Ife # 0, the solution of equation (5) can be expected to be not too different from equation (7). Consequently, one is motivated to assume that
x(t) = A(t)cos ~p,
(9)
~9(t) = tot + ~b(t),
(10)
where and the functions A(t) and q~(t) are now to be determined. Upon differentiation of equation (9) one obtains ~(t) = / / c o s Ip - (co + q~)A sin ft.
(1 I)
NON-LINEAR, NON-CONSERVATIVE SYSTEMS
465
If next the identity ?A cos ~k - ~A cos ~k= 0 is inserted into equation (11), then the equation :~(t) ---A(Tcos ~ + tn sin ~) + [(.4 + TA)cos ~k-- ~A sinG]
(12)
is obtained. The first set of terms in equation (12) is identical in form to the right side of equation (8). Therefore, one can require that (.4 + ~A) cos ~ - q~Asin ~k= 0,
(13)
which is, in essence, merely an application of the method of variation of parameters. Consequently, equation (12) reduces to :~(t) = --d (},cos ~k + o9sin tp).
(14)
Furthermore, differentiation of equation (14) and use of equation (13) yield dr(t) = -o9(.4 sin ~ + ~A cos ~) + A [~o9sin r + (2~2 - f22) cos ~].
(15)
Finally, substituting equations (9), (14), and (15) into equation (5) gives (.4 + ~A) sin ~, + Ag~eos ~k= (e/og)f(t, A e o s ~ , - ? A c o s r
(16)
Equations (13) and (16), namely,
(.4 + yA) cos ~ - q~a sin ~ = 0, (.4 + yA) sin ff + q~A cos ~k= (e/og)f(t, A cos ~k,-yA cos (t - ogA sin if), represent a pair of coupled equations in the unknowns .4 + yA and Aq~. The solution of this system for .4 + yA and ~ is
.4 + yA = (e/o~)f(t, A cos ~,-~A cos ~k- ogA sin ~k) sin ~p,
(17)
= (e/ogA)f(t, A cos ~k,-?A cos @- ogA sin @) cos ft.
(18)
To obtain the associated initial conditions, equations (9) and (14) are inserted into equation (2), to give A(0) cos ~b(0) = Ao,
?A(0) cos ~b(0) + ogA(0) sin qS(0) = 0,
whence
A(0) = aAo/o9,
cos qb(o) = o91f2,
sin qs(o) = '~/s'2.
(19)
The right sides of equations (17) and (18) are clearly periodic functions in r with period 2re, and the quantities .4 + ?A and q~ are directly proportional to e, a small parameter, so it is evident that these quantities are slowly varying functions of time. In reference [3], equations (17) and (18) were simplified by expanding the right sides as Fourier series in ~ over the interval 0 < ~ < 2n. Because .4 + ~A and q~ are slowly varying functions oftime, the Fourier coefficients were treated as being approximately constant over the interval 0 < ~k < 2rr, and the reduced equations were found to be
d + 7A -~ (e/o9)Ko(t,A),
~ = (e/coA)Po(t,A),
(20)
where 2~
Ko(t, A) = (l[2n) f f(t, A cos ff,-),A cos ~k-- OgAsin ~9)sin ~d~,,
(21)
o 2I
Po(t, A) = (l/2rc) f f(t, A cos ~ , - ? A cos Ip - o9A sin if) cos Ipd~,. o
(22)
466
G.'L. ANDERSON
It is the objective here to proceed in a somewhat more general fashion. Instead of the right sides of equations (17) and (18) being expanded in Fourier series, they will be expanded in terms of Jacobi polynomials, P~,"P)(z). The ultraspherical (Gegenbauer), Legendre, and Chebyshev polynomials are special cases of Jacobi polynomials. In order to expedite subsequent calculations, it is convenient to introduce the change of variable z = ( r - rr)/~,
- 1 < z < 1,
and the notation
F~(t, A, ~) = f ( t , A cos ~k, - y A cos ~k- toA sin ~) sin ~, F2(t, A, tp) = f ( t , A cos ip,--yA cos ~p-- o~A sin ~) cos Ip. Next the definitions
Fk(t, A, ~) = Fk[t, A, rr(z + 1)1 -- Gk(t, A, z),
k = 1,2,
are introduced, where
Gl(t, A, z) = --f(t,--A cos nz, vA cos rrz + toA sin ~rz)sin ~rz,
(23)
G2(t, A, z) = - f ( t , - - A cos nz, yA cos nz + toA sin rcz) cos r~z.
(24)
Hence, insertion of equations (23) and (24) into equations (17) and (18) yields
A + rA = (~1o~)G,(t, A, z),
dp = (e/toA) G2(t, A, z).
(25)
The expansion of Gk(t,A, z) in terms of Jacobi polynomials shall be denoted by
Gk(t, A, z) = ~ a ~' c~" k,nB~tt \ , A~P I n aJt..~ k~"l,
(26)
where Re (co), Re (fl) > - 1 . Since the Jacobi polynomials form an orthogonal set of functions with respect to the weight function (1 - z ) ' ( 1 + z) B on the interval - 1 < z < 1, it is easily demonstrated that 1
f Gk(t,A,z)p~.t~(z)(1 - z)~(l + z)Bdz a~.'f'(t,A) = -1 1
(27)
f [P~," a'(z)lZ(1 -- z)'(1 + z) adz --1
where, according to Rainville [16, p. 260], 2x+~+~F(I + cc+ n) F(I + fl + n) -1 [P~'B~(z)]2(l -- z)~(l + z)Bdz = n!(1 + o c + f l + 2 n ) F ( 1 + ~ + f l + n ) '
(28)
where F(...) denotes the gamma function. If the expansions in equation (26) are restricted to the first term only,
Gk(t, A, z) ~_ a~,'oB'(t, A) =- gk(t, A; ~, fl)
(29)
since P~o~' B)(z) = 1. But using equations (27) and (28) with n = 0 gives 1
gk(t,
F(2 + C~ + /~) | a~(t,A, z)(1 - z)'(1 + z ) B d z . d ; cq fl) = 2z+~+B F(I + ~) F(I + ,6) -13
(30)
N O N - L I N E A R , N O N - C O N S E R V A T I V E SYSTEMS
467
Therefore, in view of equation (29), one can approximate equation (25) by
X + ?A = (e/og)g~(t,A;~z, fl),
$ = (e/oJA)g,(t,A;o~,fl),
(31)
where the functions gk(t, A ;cz,fl), k -- I, 2, are obtained from equation (30). Equations (31) are the counterparts of equations (20), and indeed it will be shown later that the latter are merely a special case of the more general form given in equations (31). The procedure outlined here implies a more general averaging technique because the free choice of the parameters 0~and fl is still available, subject only to the conditions Re (c0 > - 1 , Re (fl) > -1. The merit of the present technique is that for suitable choices of ~ and fl the solutions of equations (31) will lead to somewhat more accurate approximations than do those of equations (20). It frequently happens that the integration of the first-order differential equations in equations (31) is a laborious process, but, of course, the degree of complication is dependent upon the form off(t, x, :~) in equation (5) and upon the choices of the parameters cc and fl in equation (30). However, for a suitable choice of c~ and fl, namely, c~= f l - - - 2 - 1/2, where Re (2) > - 1 / 2 , which implies only a modest sacrifice of generality, more efficacious forms of equation (31), at least for certain practical problems, can be obtained. In this case the Jacobi polynomials P(nx-ll2" ;t-ll2)(z), which can be obtained from (2 + I/2)n C~(z) (2).), '
p~x-l/2, ~-1/2)(z ) =
where Ca(z) denotes the Gegenbauer polynomial of degree n and (...), the factorial function [16, p. 22], are known as ultraspherical polynomials. Thus, equation (30) becomes 1
F(22 + 1) f g,(t, A; 2 - 1/2, ;. - 1/2) = 22x[r(;. + l/2)]" ,J Ga(t, A, z)(l - z2) x-xn dz.
(32)
--1
However, by using Legendre's duplication formula ([16], p. 24) in the form F(2). + 1) =
22x F(2 + 1/2) F(2 + 1) ~112
equation (32) can be expressed as follows: 1
gk(t,A;2-- 1/2,2-- I/2) --=h,(t,A;2)= r d ntO. F 0 . ++ I)I/2)_f"t Gk(t,A,z)(1 -z2)X-ll2dz.
(33)
Consequently, equation (31) can be replaced by
d + TA = (e/o) ht(t, A ; 2),
$ = (e/o~A) h2(t, A; 2),
(34)
which are somewhat more general than are the corresponding expressions reported by Sinha and Srinivasan [1]. Furthermore, if). = 1/2, then equation (33) reduces to 1
hk(t,A;l/2)=~
G(t,A,z)dz --! 2n
0
Ko(t, A) Po(t, A)
if k = 1, if k = 2,
(35)
468
O.L. ANDERSON
upon comparing equations (35) and (21)-(22). Therefore, the procedure proposed by Brunelle [3] to include explicitly the effect o f linear damping in the amplitude and phase differential equations is a special case o f equation (31) with c~= ,8 = 0 and, of course, equation (34) with ). = 1/2. In applications o f equation (33) to specific problems, one often encounters the definite integrals 1
-J sin (nnz)(l -- z2)a-,/2dz = O,
(36)
--1
i
_ cos (nrcz)(l
-- z2) 1-1/2
dz =
nt/2F(2 + 1/2) F(2 + 1) r/~,
(37)
where tlo = 1,
r/~ =
F(2 + 1)j~(nn), (nx/2) x
n = 1 , 2 , 3 .....
(38)
Jx denoting the Bessel function of the first kind o f order 2. The integral in equation (37) is found in Gradshteyn and Ryzhik [17, p. 403, formula 21]. Selected numerical values oft/, are given in the Appendix. 3. AN EXAMPLE Next, as an example o f a p r o b i e m t o which the theory derived above can be applied, one may consider the non-linear differential equation 2 + (2~ + ~ l x 2) 2 + [[22 + tp(t) cos O(t)] x + e(~t x 3 + ct2x s) = O,
(39)
where/l(t) and 0(0 are functions o f t, which, for the present, shall remain unspecified, and )'1, cq, and ct2 are constants. Comparing equations (5) and (39) indicates that
f(t,X,5:)=~tx2yc + Iz(t)cosO(t)x +~tX3 + ~2xS.
(40)
Substituting equation (40) into equations (23) and (24) leads to
G~(t, A, z) = (A/2)/t cos 0 sin 27rz + (~1 -- Y~',)(A[2)3( 2 sin 27rz + sin 47rz) - co~1(A/2)3(1 - cos 4nz) + %(A/2)5(5 sin 27rz + 4 sin 47rz + sin 6nz),
(41)
G2(t, A, z) = (A/2) p cos 0(1 + cos 2nz) + (cq -- ~,)q) (A/2)3(3 + 4 cos 2rrz + cos 47rz) + + Ocz(A/2)s(10+ 15 cos 2nz + 6 cos 4zz + cos 6rrz) - mvI(A/2)3(2 sin 2nz + sin 4zrz).
(42)
Insertion o f equations (41) and (42) into equation (33) yields
hi(t, A ;2) = --~oyt(A/2)3(l -- r/,),
(43)
h2(t, A ; ).) = (A/2)/t cos 0(1 + r/2) + (r -- 7Yl)(A/2)3( 3 + 4t12 + 1/4) + + =2(Z/2)s(10 + 15r/2 + 6r/, + r/6),
(44)
where equations (36)-(38) have been used. Combining equations (43) and (44) with equation (34) yields the amplitude and phase differential equations: A + yA = --~?,1(A/2)3(1 -- n,),
(45)
= (e/2w)[(1 + r/2)/a(t)cos O(t) + (oq - ~,~h)(A/2)2(3 + 4r/2 + r/4) + + ~2(A/2)4(10 + 15r/2 + O/4 + r/6)].
(46)
NON-LINEAR, NON-CONSERVATIVE SYSTEMS
469
In view o f the initial conditions in equation (19), integration o f equation (45) leads to A(t) = f2Ao e-"[co 2 + B (12Ao)'(l -- e-2~')] -1/2,
(47)
where a = (e~q[8y)(1 - t/,). With the objective of integrating equation (46), this equation may now be written as
t~(t) --- ~(t) + r, A 2 + r2 A',
(48)
~(t) = (e]2co)(l + q2) It(t) cos O(t)
(49)
where and rl = (e/8w)(cq - ~yl)(3 + 4q2 + q4),
r2 = (e~2/32co)(10 + 151/2 + 6r/4 + r/6).
(50)
But, according to equation (45) and the chain rule of calculus, d dA d d ~ ( 4 , - 0 = d-'7d---J(4' - ~) = - r A ( 1 + BA 2) ~-j (4, - ~), and consequently equation (48) becomes d
(rt A + r2 A 3)
T~ (4, - O =
~(l + B a ~) "
(51)
The solution of equation (51) subject to the initial conditions in equation (19) is 4,(0 = cos-l(co/f2) + ~(t) - ~(0) + (r2/2~B)[(f2Ao/o)) 2 -- A 2] +
1 + 2~
[ I+BA z ] (r2 - rl B ) I n 1 +-B(---~oflO)i]"
(52)
However, if 7~ = 0, then B = 0, and equation (52) is not suitable for computer calculations. In this case equation (47) reduces to A( t ) = (~?Ao/w) e - ' ,
(53)
whereas equation (48) must be replaced by = ~ +r*a 2+r2A*,
where r* = rt [~=o- The counterpart of equation (51) is d d---A(4, - ~) = - - A ( r * + r2 A2)/~, whence 4,(0 = cos-l(co/O) + ~(t) - ~(0) + 7(f2Ao/27co)'(l - e-2")[2r] ' + rz(f2AofiO)2(l + e-Z')].
(54)
4. NUMERICAL RESULTS Linear expressions in t for the functions It(t) and O(t), namely, It(t) = Ito + It1 t,
O(t) = Oo + 01 t,
(55)
where Ito, Iti, 0o, and 01 are constants, are ofconsiderable importance in the theory of dynamic stability--see, for example, reference [18]. Inserting equation (55) into equation (49) and integrating the result gives ~(t) = (e/2coO~)(1 + r/2)~ o Ot sin (0o + 01 t) + It1 01 t sin (0o + Ot t) + Itt cos (0o + 01 t)],
470
G.L. ANDERSON
and consequently ~(t) -- ~(0) = (e/2to02)(1 + q2){Po 0~ [sin (0o + 0~ 0 - sin 0o] + tq 0~ t sin (0o + 0~ t) + + tq[cos (0o + 0~ t) - cos 0o]}. (56) Therefore, the Solution of 5/+ (2~ + ~),ax 2) ~ + [f22 + ~(Po + P~ t) cos (0o + 0~ t)] x + ~(cq x 3 + cq x 5) = 0
(57)
is, approximately,
x(t) = A(t)cos r
(58)
with, by virtue of equation (56), r2 2-~1 [1 I + B A 2 O(t) =- cot + cos-'(oJfl2) + 2-~-~[(QAo]co)2 --A 2] + (rz - rl B)ln +~w)2j
"]
+
+ (e/2co0~)(l + rh)(po 0 x [sin (0o + Ot t) -- sin 0o] + Pt 01 t sin (0o + 0~ t) + tq[cos (0o + 0~ t) - cos 0o]},
(59)
and
A(t) = I2Aoe-~t[to 2 + B(I2Ao)2(I -- e-2~')] -~/2,
(60)
provided that 71 :A 0, and
x(t) = (t2Ao/co) e-" cos ~b(t),
(61)
with ~,(t) = cot + cos-X(to/f2) + y(t2Ao/2),co)X(l -- e -2~') [2r* + r2(f2Ao/co)2(l + e-a")] + + (e/2co0~)(1 + q2){Po 01[sin (0o + 0, t) -- sin 0o] + lq 0, tsin (0o + 0x t) + + ,tq [cos (0o + 01 t) -- cos 0o]},
(62)
whenever ~,, = 0. At this point it should be remarked that the approximate solution of equation (57) as embodied in equations (59) and (60) does not reveal the possibility of the existence of unstable motions that are known to occur in certain systems whose motions are described by linear or non-linear differential equations of the Mathieu type (if F~ = 0, equation (57) falls in this category). This phenomenon, frequently called parametric resonance, has been discussed in some detail by Bolotin [19] and Bogoliubov and Mitropolsky [20]. In this study only a firstorder approximate solution was obtained; however, in reference [20] it is shown that the Krylov-Bogoliubov method predicts the possibility of resonance only in the second and higher order approximations. Consequently, the method of solution described herein can be used only in those situations in which the ratio of the frequency of the driving force and the natural frequency of the associated linear system is not a rational number, i.e., in non-resonance cases. In the event that the driving frequency is equal to or approximately equal to the natural frequency, a different method of analysis must be applied, and the interested reader is referred to pages 196-236 and 267-284 of reference [20]. Finally, it is worth mentioning that Huston and Doty [21] have applied the Krylov-Bogoliubov method in the first-order approximation to linear differential equations of the type 5/ + ~(t)5/ + ~ ( t ) x = o,
where ~(t) and ~j(t) are known functions, and they reported an error analysis which indicates the domain of validity of their method, which may be considered a special case of the technique described here.
NON-LINEAR, NON-CONSERVATIVESYSTEMS
471
The expression for the velocity :~(t) in equation (14) can be simplified somewhat. Set ~(t) = ~b(t) -- COS-l(to/12). Then, in view of equation (19), one can show that equation (14) assumes the form 2 ( 0 = - f 2 A ( t ) sin 7'(0.
Numerical calculations were performed with Ao=e=f2=cq
= 1,
ct2 = 0 ,
It(t) = 0,-0"02 (1 + t),
1,
yl =0,0"001
y = 0"05,
O(t) = 1.8t
and the values 2 = 0, I/2, 1, 3/2, 2, 5/2, 7/2, and 9/2. Specifically, the following three nonlinear differential equations, all possessing a linear damping term, were solved approximately by means of equations (58)-(62), as appropriate: 2 + O. I Ec+ x + xa + xS = O,
5( + 0"1 .~ + [I -- 0"02(1 + t ) c o s ( l . 8 t ) ] x + x 3 + x 5 = O, .~ + 0"1(1 + 0"01 x2).~ + [1 -- 0"02(1 + t)cos (l'8t)]x + x 3 = 0,
(63)
subject to the initial conditions x(O) = l,
.~(0) = O.
In order to judge the accuracy of the approximate formulas in equations (58)-(62), these same equations were also solved numerically with the fourth-order Runge-Kutta method. The variations o f x ( t ) as a function of time t for these three differential equations are plotted in Figures 1, 2, and 3 for only the first few cycles ofmotion. The value of the initial displacement Ao was selected to be a fairly large number in order to emphasize the effect of the value of the parameter 2. If 1, --- 0, the approximate expression for the amplitude is independent o f any of the coefficients ~I,, which depend upon 2, and if ~i :~ 0, the approximate expression for A ( t ) in equation (60) depends upon only 114,which, generally, is a small number relative to unity. For the three differential equations given in equation (63), the numerical computations and the graphs in Figures 1-3 reveal that the choice of 2 = 3/2 leads to somewhat more accurate values of the phase than did the choice of 2 = 1/2, which corresponds to the method of reference [3]. The approximations obtained for 2 = 1 and 2 were only slightly less accurate than those obtained for ). = 3/2, but for the other values of 2 that were tried the results were I.O
8
,
,
i
'
~
J
'
I
l
t
t
;
I
'
t
~
i
I
'
i
~
i
I
I
I
l
,
i
t
~
I
J
l-
/ E x a c t
05
\\
/;/I,/
\\
/ ./
\,\
v
',,\ ,;'/
-0.5
-I.0
i
,
i
I
~/2
,
i
t
J
I
T
i
,
T
i
I
i
i
i
3r/2
i
t
2~"
,
,
,
,
" l
5~/2
,
,
i
i
3t-
t
,
'
i
r 7r/:
t
Figure 1. Variation of x(t) with t for ~ + 0.12 + x + x 3 + x s = O, x(O)= 1, 2(0) = O.
472
G.L. ANDERSON I'O
,
~ i l i t l l l l ~ l l t l
\
,~
iILi
l l ~ l
j l l l l t
/,r ///
i L l l l
i
tl
iii8
ii
!1 'x\
//
\\
e
o
-oV V --IC
I l t l
III
I!11
~'/2
I
~"
111
I l t l
37r/2
I III
2zr
I I I I I I
5~/2
tit
t
3~"
III
?z~/2
l l l l l
4~"
II
9~'/2
t
Figure 2. Variation of x(t) with t for .~ + 0.1.~ + [1 - 0.02(I + t)cos(l.8t)]x + x a + x s = 0, x(0) = 1, ~(o)
=
o.
I'O
it
1111t
iii
l l l l l l l l l
~-
--I O
--
IIII
t t l l l l
z-/2
~r
l i l t
3~r/2
i l l l l l l l l l
i1111
i l l l l l
i i i i
~ExQct
iIII
iIi
i I i 1 1 t l
2~"
5r/2
i
37r
ii~
I , ~ , l ! ! t
71r/2
4r
pp
9rr/2
l
Figure 3. Variation of x(t) with t for $ + 0"l(l + 0'01x2)~ + [I - 0.02(1 + t)cos(l.8t)]x + x 3 = 0, x(0) = 1, ~(o) = 0.
n o t i c e a b l y less a c c u r a t e . Therefore, for an a p p o s i t e choice o f the p a r a m e t e r 2, the present m e t h o d a p p e a r s to be s o m e w h a t m o r e accurate t h a n the m e t h o d suggested in reference [3] a n d has the a d d i t i o n a l feature that as e tends to zero the s o l u t i o n o f e q u a t i o n (5) tends to the correct solution o f the c o r r e s p o n d i n g linear p r o b l e m , whereas the m e t h o d o f s o l u t i o n given in reference [1], as a p p l i e d to e q u a t i o n (5), does not.
REFERENCES 1. S. C. SINHA and P. SRINIVASAN1971 Journal of Sound and Vibration 18, 55-60. Application of ultraspherical polynomials to non-linear autonomous systems. 2. N. KRYLOV and N. BOGOLIUaOV 1947 Introduction to Non-Lhtear Mechanics. Princeton, New Jersey: Princeton University Press. See pp. 8-14. 3. E.J. BRONELLE1967 InternationalJournalofNon-Lhzear Mechanics 2, 404-415. Transient response of second-order non-linear equations. 4. I. N. S. MORAY 1971 lnternationalJournalofNon-Linear Mechanics 6, 45-53. A unified K r y l o v Bogoliubov method for solving second-order non-linear systems. 5. H . H . DENt,IAN 1959 AmericanJournalofPhysics 27, 524-525. Amplitude dependence of frequency in a linear approximation to the simple pendulum equation. 6. H. H. DENMAN 1964 Journal of the Industrial Mathematics Society 14, 9-20. Application of ultraspherieal polynomials to asymmetric non-linear oscillations.
NON-LINEAR, NON-CONSERVATIVESYSTEMS
473
7. H. H. DENMANand J. E. HOWARD 1964 Quarterly of Applied Mathematics 21,325-330. Application ofultraspherical polynomials to non-linear oscillations. I. Free oscillations of the pendulum. 8. H. H. DENMAN and Y. K. LIU 1965 Journal of the Industrial Mathematics Society 15, 19-35. Application of ultraspherical polynomials to non-linear systems with step function excitation. 9. R. M. GARDE 1965 Indian Journal of Mathematics 7, I 1 I - I 17. Application of Gegenbauer polynomials to non-linear oscillations--forced and free oscillations without damping. 10. H. H. DENMANand Y. K. LIU 1965 Quarterly of Applied Mathematics 22, 273-292. Application of ultraspherical polynomials to non-linear oscillations. II. Free oscillations. 11. R. K. SAXENAand K. S. KUSHAWA 1970 Proceedings ofthe National Academy of Sciences, India A40, 65-72. Application of Jacobi polynomials to non-linear oscillations. 12. R. M. GARDE 1967 Journal of Scientific and Engineering Research 11, 157-166. Application of Gegenbauer polynomials to non-linear damped oscillations. 13. Y. K. LIu 1967 Journal of Applied Mechanics 34, 223-226. Application of ultraspherical polynomials to nonlinear forced oscillations. 14. R. M. GARDE 1969 Proceedings of the National Academy of Sciences, India A39, 445-450. Application of Jacobi polynomials to some nonlinear oscillations. 15. S. C. SINHAand P. SRINIVASAN1972JournalofSoundand Vibration 22, 211-219. An approximate analysis of non-linear non-conservative systems subjected to step function excitation. 16. E. D. RAINVILLE1960 SpecialFunctions. New York: The Macmillan Company. 17. I. S. GRADSHTEYNand I. M. RYZHIK 1965 Tables of Integrals, Series, and Products. New York: Academic Press. 18. R. M. EVAN-IWANOWSKI,W. F. SANFORDand T. KEHAGIOGLOU1971 DecelopmentsofTheoretical and Applied Mechanics, Vol. 5, 714-743. Chapel Hill, North Carolina: University of North Carolina Press. Nonstationary parametric response of a nonlinear column. 19. V.V. BOLOrIN 1964 The Dynamic Stability of Elastic Systems. San Francisco: Holden-Day, Inc. See pp. 18-27 and 63-74. 20. N. N. BOGOLIUBOVand Y. A. MrrROPOLSKY 1961 Asymptotic Methods in the Theory of NonLhtear Oscillations. Delhi, India: Hindustan Publishing Corporation. See pp. 196-236 and 267-284. 21. R. L. HUSTON and L. T. DOTY 1966 SIAM Journal on Applied Mathematics 1,1, 424-428. Note on the Krylov-Bogoliubov method applied to linear differential equations.
APPENDIX I n e q u a t i o n (38), the coefficients q.(2) were defined as
r(;. §
(Al)
=
F o r 2 = 0, 1, 2, the first six values o f !1.(2) are given in T a b l e .I. 9 TABLE 1
Vahtes of q.(2) for some integer vahtes of 2 n 1 2 3 4 5 6
tl.(O) -0"304 0"220 -0"181 0"157 -0"141 0"129
242 18 276 91 211 46 507 39 182 05 063 52
t/.(1) 0"181 -0.676 0"375 -0"245 0"177 -0"135
191 76 034 59 022 53 943 43 012 25 185 03
t/.(2)
x x x x x
10 - l 10 -1 10 - l 10 -1 10 -x
0-393 477 -0"583 367 0"196 980 -0"922 538 0"515 143 -0"321 035
93 59 72 17 73 10
x x x • x
I0 -1 10 -1 10 -2 10 -2 10 -2
474
G . L . ANDERSON
If 2 = k + 1/2, k = 0, I, 2, ..., then the Bessel functions Jk+l/2 can be expressed in terms of elementary functions:
Jl/2(x)= \~x]
=
J3/2(x)=\~x]
sinx,
[ 2 \"2 ) [[ 3
I)sin x- 3cosx]
Jw2(x)=(~x) ',2[{156]sinx+( l-:Oc~ [2\ 1/2
45
105/sinx + (~_-- ~)
cosx].
In these cases, equation (A1) assumes rather simple forms: namely, 0.(1/2) = 0,
0.(7/2)
q.(3/2)
3 ( - 1 ) "§ (n~) 2 ,
105(-I)" (nr0-----~ [(n~) 2 - 15],
~.(5/2)
45(--1) "§ (n~) 4 ,
4725 th(9/2) = (-1)" ~ - ~ [2(n~) 2 - 21].
In Table 2, the first six values oft/.().) are given for 2 = k + 1[2, k = 0, 1, 2, 3, and 4. TABLE 2
Vahtes of lh(2)for certain vahtes of 2 n r/.(112)
I 0 2 3 4 5 6
0 0 0 0 0
r/.(3/2)
0"30396355 -0"75990888 • 10 -1 0"33773729 x 10 -I -0"18997722 x 10 -1 0"12158542x10 -1 -0"84434320x10 -2
r/.(5/2)
1/.(7/2)
t1.(912)
0"46196920
0"56032618
0"62783569
-0.28873075 x 10 -1 0"57033235 x 10 -2 -0"18045672 • 10 -2 0"73915072x10 -3 -0"35645772x10 -3
0"41772785 • 10 -t -0"11060492 x 10 -1 0"38106924 x 10 -2 -0"16198367x10 -2 0"79662120x10 -3
0"11273727 -0"118 89707 x 10 -I 0-22402199 x 10 -2 -0"60231881• -3 0.20445481x10 -3