- Email: [email protected]
LETTER TO THE EDITOR: ON JUMP FREQUENCIES IN THE RESPONSE OF THE DUFFING OSCILLATOR
Journal of Sound and Vibration (1996) 198(4), 522–525
ON JUMP FREQUENCIES IN THE RESPONSE OF THE DUFFING OSCILLATOR K. W Department of Mechanical and Process Engineering, University of Sheffield, Mappin Street, Sheffield S1 3JD, England (Received 8 March 1996, and in final form 18 April 1996)
1. The method of harmonic balance is a textbook approach to approximating the Frequency Response Functions (FRFs) of non-linear systems [1]. Although it is extremely simple in principle, it reproduces well the qualitative features of system FRFs. The aim of this letter is to demonstrate that the agreement can also be quantitative, at least for one widely studied system. Harmonic balance assumes a harmonic system response Y cos (vt) to a harmonic excitation X cos (vt) (Y is taken to be complex in order to encode phase information). These expressions are substituted into the equation of motion and the coefficient of the fundamental component is projected out: this leads to an equation for the FRF. In the case of the Duffing oscillator, my¨ + cy˙ + ky + k3 y 3 = X cos (vt),
(1)
harmonic balance yields the relation, X 2 = Y 2{[−mv 2 + k + 34 k3 Y 2 ]2 + c 2v 2}.
(2)
For a given amplitude of excitation X, this is a cubic in the response Y 2 (effectively Y, as one can disregard negative roots). Depending on the value of v, this equation can have one or three real roots. This leads to ‘‘jump’’ phenomena in the FRF as v passes through the bifurcation points and the amplitude switches between stable solutions. In a recent study, Friswell and Penny [2] computed the bifurcation points of the FRF, not in the first harmonic balance approximation but for a multi-harmonic series solution. They obtained excellent results using expansions up to the third and fifth harmonics for the response. Newton’s method was used to solve the equations obtained. Even up to the third harmonic, the expressions are exceedingly complex. If a trial solution of the form. y(t) = Y1 sin (vt + f1 )+Y3 sin (3vt + f3 ) (with the phases explicitly represented) is substituted in equation (1), projecting out the coefficients of sin (vt), cos (vt), sin (3vt) and cos (3vt) leads to the system of equations −mv 2Y1 cos f1 − cvY1 sin f1 + kY1cos f1 +34 k3 Y13 cos f1 + 32 k3 Y1 Y32 cos f1 − 34 k3 Y12 y3 cos f3 cos 2f1 = X, −mv 2Y1 sin f1 − cvY1 cos f1 + kY1 sin f1 +34 k3 Y13 sin f1 + 32 k3Y1 Y32 sin f1 − 34 k3 Y12 y3 sin f3 cos 2f1 = 0, −9mv2Y3 cos f3 − 3cvY3 sin f3 + kY3 cos f3 −14 k3 Y13 cos3 f1 + 34 k3 Y33 cos f3 − 34 k3Y13 cos f1 sin2 f1 + 32 k3 Y12 Y3 cos f3 = 0, 522 0022–460X/96/490522 + 04 $25.00/0
7 1996 Academic Press Limited
523
−9mv Y3 sin f3 + 3cvY3 + kY3 sin f3 2
+14 k3 Y13 sin3 f1 + 34 k3 Y33 sin f3 − 34 k3 Y13 cos2 f1 sin f1 + 32 k3 Y12 Y3 sin f3 = 0. The equations for higher order approximations are best left to computer algebra systems and in reference [2] details of the authors’ approach are given. 2. The object of the current letter is to examine the accuracy of the first order approximation, but to eschew numerical methods as far as possible. The approach taken is to compute the discriminant of the cubic equation (2) which indicates the number of real solutions [3]. In a convenient notation, equation (2) is a3 Y 6 + a2 Y 4 + a1 Y 2 + a0 = 0,
(3)
Now, dividing by a3 and making the transformation Y 2 = z − a2 /(3a3 ) yields the normal form z 3 + pz − q = 0
(4)
and the discriminant D is then given by D = −4p 3 − 27q 2.
(5)
Now, the original cubic (3) has three real solutions if D e 0 and only one if D Q 0. The bifurcation points are therefore obtained by solving the equation D = 0. For equation (1) this is an exercise in computer algebra (in this case Mathematica [4]) and the resulting discriminant is D = (256/729k36 )(−64c 2k 4v 2 − 128c 4k 2v 4 + 256c 2k 3mv 4 − 64c 6v 6 + 256c 4kmv 6 −384c 2k 2m 2v 6 − 128c 4m 2v 8 + 256c 2km 3v 8 − 64c 2m 4v 10 − 48k 3k3 X 2 −432c 2kk3 X 2v 2 + 144k 2k3 mX 2v 2 + 432c 2k3mX 2v 4 − 144kk3 m 2X 2v 4 + 48k3m 3X 2v 6 − 243k32 X 4 ). As bad as this looks, it is just a quintic in v2 and can have at most five independent solutions for v. In fact, in all the cases examined below, it had two real roots and three complex ones. The lowest real root is the bifurcation point for a downward sweep and the highest is the bifurcation point for an upward sweep. The equation D = 0 is solved T 1 ‘‘Exact’’ jump frequencies for upward sweep [2]
Forcing 0·1 0·3 1·0 3·0 10·0 30·0 100·0
Damping coefficient ZXXXXXXXXXXXXCXXXXXXXXXXXXV 0·01 0·03 0·1 0·3 3·0388 5·1622 9·3608 16·1817 29·5233 51·1258 93·3362
1·8578 3·0388 5·4362 9·3608 17·0553 29·5233 53·8909
1·2287 1·7797 3·0905 5·1624 9·3608 16·1817 29·5233
— — 1·8630 3·0411 5·4373 9·3614 17·0557
524
T 2 Estimated jump frequencies and percentage errors (bracketed) for upward sweep
Forcing
Damping coefficient ZXXXXXXXXXXXXCXXXXXXXXXXXXV 0·01 0·03 0·1 0·3
0·1
3·0290 (−0·32)
1·8520 (−0·32)
1·2271 (−0·13)
— —
0·3
5·1464 (−0·31)
3·0290 (−0·32)
1·7742 (−0·31)
— —
1·0
9·3330 (−0·30)
5·4196 (−0·31)
3·0292 (−0·33)
1·8571 (−0·32)
3·0
16·1341 (−0·29)
9·3330 (−0·30)
5·1465 (−0·31)
3·0312 (−0·33)
10·0
29·4368 (−0·29)
17·0052 (−0·29)
9·3330 (−0·30)
5·4207 (−0·31)
30·0
50·9762 (−0·29)
29·4368 (−0·29)
16·1341 (−0·29)
9·3336 (−0·30)
100·0
93·0632 (−0·29)
53·7332 (−0·29)
29·4368 (−0·29)
17·0055 (−0·29)
effortlessly by using computer algebra. However, note that an analytical solution is possible by using elliptic and hypergeometric functions [5]. 3. Following reference [2], the values m = k = k3 = 1 were chosen. This is because reference [2] presents bifurcation points for a ninth-harmonic solution of equation (1) and these can therefore be taken as reference data. A range of c and X values were examined both for the upward sweep and downward sweep (see Tables 1 and 3 respectively). The estimated bifurcation points obtained by the current method for the upward sweep are given in Table 2. Over the examples given, the percentage errors range from −0·15 to −0·33. This is in very good agreement with reference [2]. The results for the downward
T 3 ‘‘Exact’’ jump frequencies for downward sweep [2]
Forcing 0·1 0·3 1·0 3·0 10·0 30·0 100·0
Damping coefficient ZXXXXXXXXXXXXCXXXXXXXXXXXXV 0·01 0·03 0·1 0·3 1·1705 1·3307 1·6498 2·1411 3·0003 4·2009 6·1769
1·1696 1·3302 1·6495 2·1409 3·0001 4·2009 6·1768
1·1580 1·3242 1·6462 2·1388 2·9989 4·2000 6·1763
— — 1·6164 2·1208 2·9876 4·1925 6·1713
525
T 4 Estimated jump frequencies and percentage errors (bracketed) for downward sweep
Forcing
Damping coefficient ZXXXXXXXXXXXXCXXXXXXXXXXXXV 0·01 0·03 0·1 0·3
0·1
1·1703 (−0·02)
1·1694 (−0·02)
1·1577 (−0·03)
— —
0·3
1·3302 (−0·04)
1·3296 (−0·05)
1·3236 (−0·05)
— —
1·0
1·6483 (−0·09)
1·6480 (−0·09)
1·6447 (−0·09)
1·6147 (−0·11)
3·0
2·1381 (−0·14)
2·1379 (−0·14)
2·1359 (−0·14)
2·1177 (−0·15)
10·0
2·9950 (−0·18)
2·9949 (−0·17)
2·9936 (−0·18)
2·9823 (−0·18)
30·0
4·1926 (−0·20)
4·1925 (−0·20)
4·1917 (−0·20)
4·1841 (−0·20)
100·0
6·1639 (−0·21)
6·1638 (−0·21)
6·1632 (−0·21)
6·1582 (−0·21)
sweep are given in Table 4; the errors range from −0·02 to −0·21. Again, this compares well with reference [2]. 4. This report is essentially a piece of propaganda for the harmonic balance method in its simplest guise. It is shown that the bifurcation points from the method agree very well with the ‘‘exact’’ results. There is a small shift in emphasis here from reference [2]; the bifurcation points are obtained from an exact formula based on a first order approximation. In reference [2], the points are obtained from numerical analysis of a higher order approximation. The distinction is perhaps otiose, as the final resort in this case was still a computer solution of the quintic discriminant. ACKNOWLEDGMENT
The author would like to thank Dr Mike Friswell for making his data available. 1. A. H. N and D. T. M 1979 Nonlinear Oscillations. New York: Wiley Interscience. 2. M. F and J. E. T. P 1994 Journal of Sound and Vibration 169, 261–269. The accuracy of jump frequencies in series solutions of the response of a Duffing oscillator. 3. G. B and S. M 1977 A Survey of Modern Algebra. New York: Macmillan; fourth edition. 4. S. W 1991 Mathematica: A System For Doing Mathematics by Computer. Cambridge, MA: Addison Wesley; second edition. 5. F. K 1877 Lectures on the Icosahedron and the Solution of Equations of the Fifth Degree. New York: Dover.
ON JUMP FREQUENCIES IN THE RESPONSE OF THE DUFFING OSCILLATOR K. W Department of Mechanical and Process Engineering, University of Sheffield, Mappin Street, Sheffield S1 3JD, England (Received 8 March 1996, and in final form 18 April 1996)
1. The method of harmonic balance is a textbook approach to approximating the Frequency Response Functions (FRFs) of non-linear systems [1]. Although it is extremely simple in principle, it reproduces well the qualitative features of system FRFs. The aim of this letter is to demonstrate that the agreement can also be quantitative, at least for one widely studied system. Harmonic balance assumes a harmonic system response Y cos (vt) to a harmonic excitation X cos (vt) (Y is taken to be complex in order to encode phase information). These expressions are substituted into the equation of motion and the coefficient of the fundamental component is projected out: this leads to an equation for the FRF. In the case of the Duffing oscillator, my¨ + cy˙ + ky + k3 y 3 = X cos (vt),
(1)
harmonic balance yields the relation, X 2 = Y 2{[−mv 2 + k + 34 k3 Y 2 ]2 + c 2v 2}.
(2)
For a given amplitude of excitation X, this is a cubic in the response Y 2 (effectively Y, as one can disregard negative roots). Depending on the value of v, this equation can have one or three real roots. This leads to ‘‘jump’’ phenomena in the FRF as v passes through the bifurcation points and the amplitude switches between stable solutions. In a recent study, Friswell and Penny [2] computed the bifurcation points of the FRF, not in the first harmonic balance approximation but for a multi-harmonic series solution. They obtained excellent results using expansions up to the third and fifth harmonics for the response. Newton’s method was used to solve the equations obtained. Even up to the third harmonic, the expressions are exceedingly complex. If a trial solution of the form. y(t) = Y1 sin (vt + f1 )+Y3 sin (3vt + f3 ) (with the phases explicitly represented) is substituted in equation (1), projecting out the coefficients of sin (vt), cos (vt), sin (3vt) and cos (3vt) leads to the system of equations −mv 2Y1 cos f1 − cvY1 sin f1 + kY1cos f1 +34 k3 Y13 cos f1 + 32 k3 Y1 Y32 cos f1 − 34 k3 Y12 y3 cos f3 cos 2f1 = X, −mv 2Y1 sin f1 − cvY1 cos f1 + kY1 sin f1 +34 k3 Y13 sin f1 + 32 k3Y1 Y32 sin f1 − 34 k3 Y12 y3 sin f3 cos 2f1 = 0, −9mv2Y3 cos f3 − 3cvY3 sin f3 + kY3 cos f3 −14 k3 Y13 cos3 f1 + 34 k3 Y33 cos f3 − 34 k3Y13 cos f1 sin2 f1 + 32 k3 Y12 Y3 cos f3 = 0, 522 0022–460X/96/490522 + 04 $25.00/0
7 1996 Academic Press Limited
523
−9mv Y3 sin f3 + 3cvY3 + kY3 sin f3 2
+14 k3 Y13 sin3 f1 + 34 k3 Y33 sin f3 − 34 k3 Y13 cos2 f1 sin f1 + 32 k3 Y12 Y3 sin f3 = 0. The equations for higher order approximations are best left to computer algebra systems and in reference [2] details of the authors’ approach are given. 2. The object of the current letter is to examine the accuracy of the first order approximation, but to eschew numerical methods as far as possible. The approach taken is to compute the discriminant of the cubic equation (2) which indicates the number of real solutions [3]. In a convenient notation, equation (2) is a3 Y 6 + a2 Y 4 + a1 Y 2 + a0 = 0,
(3)
Now, dividing by a3 and making the transformation Y 2 = z − a2 /(3a3 ) yields the normal form z 3 + pz − q = 0
(4)
and the discriminant D is then given by D = −4p 3 − 27q 2.
(5)
Now, the original cubic (3) has three real solutions if D e 0 and only one if D Q 0. The bifurcation points are therefore obtained by solving the equation D = 0. For equation (1) this is an exercise in computer algebra (in this case Mathematica [4]) and the resulting discriminant is D = (256/729k36 )(−64c 2k 4v 2 − 128c 4k 2v 4 + 256c 2k 3mv 4 − 64c 6v 6 + 256c 4kmv 6 −384c 2k 2m 2v 6 − 128c 4m 2v 8 + 256c 2km 3v 8 − 64c 2m 4v 10 − 48k 3k3 X 2 −432c 2kk3 X 2v 2 + 144k 2k3 mX 2v 2 + 432c 2k3mX 2v 4 − 144kk3 m 2X 2v 4 + 48k3m 3X 2v 6 − 243k32 X 4 ). As bad as this looks, it is just a quintic in v2 and can have at most five independent solutions for v. In fact, in all the cases examined below, it had two real roots and three complex ones. The lowest real root is the bifurcation point for a downward sweep and the highest is the bifurcation point for an upward sweep. The equation D = 0 is solved T 1 ‘‘Exact’’ jump frequencies for upward sweep [2]
Forcing 0·1 0·3 1·0 3·0 10·0 30·0 100·0
Damping coefficient ZXXXXXXXXXXXXCXXXXXXXXXXXXV 0·01 0·03 0·1 0·3 3·0388 5·1622 9·3608 16·1817 29·5233 51·1258 93·3362
1·8578 3·0388 5·4362 9·3608 17·0553 29·5233 53·8909
1·2287 1·7797 3·0905 5·1624 9·3608 16·1817 29·5233
— — 1·8630 3·0411 5·4373 9·3614 17·0557
524
T 2 Estimated jump frequencies and percentage errors (bracketed) for upward sweep
Forcing
Damping coefficient ZXXXXXXXXXXXXCXXXXXXXXXXXXV 0·01 0·03 0·1 0·3
0·1
3·0290 (−0·32)
1·8520 (−0·32)
1·2271 (−0·13)
— —
0·3
5·1464 (−0·31)
3·0290 (−0·32)
1·7742 (−0·31)
— —
1·0
9·3330 (−0·30)
5·4196 (−0·31)
3·0292 (−0·33)
1·8571 (−0·32)
3·0
16·1341 (−0·29)
9·3330 (−0·30)
5·1465 (−0·31)
3·0312 (−0·33)
10·0
29·4368 (−0·29)
17·0052 (−0·29)
9·3330 (−0·30)
5·4207 (−0·31)
30·0
50·9762 (−0·29)
29·4368 (−0·29)
16·1341 (−0·29)
9·3336 (−0·30)
100·0
93·0632 (−0·29)
53·7332 (−0·29)
29·4368 (−0·29)
17·0055 (−0·29)
effortlessly by using computer algebra. However, note that an analytical solution is possible by using elliptic and hypergeometric functions [5]. 3. Following reference [2], the values m = k = k3 = 1 were chosen. This is because reference [2] presents bifurcation points for a ninth-harmonic solution of equation (1) and these can therefore be taken as reference data. A range of c and X values were examined both for the upward sweep and downward sweep (see Tables 1 and 3 respectively). The estimated bifurcation points obtained by the current method for the upward sweep are given in Table 2. Over the examples given, the percentage errors range from −0·15 to −0·33. This is in very good agreement with reference [2]. The results for the downward
T 3 ‘‘Exact’’ jump frequencies for downward sweep [2]
Forcing 0·1 0·3 1·0 3·0 10·0 30·0 100·0
Damping coefficient ZXXXXXXXXXXXXCXXXXXXXXXXXXV 0·01 0·03 0·1 0·3 1·1705 1·3307 1·6498 2·1411 3·0003 4·2009 6·1769
1·1696 1·3302 1·6495 2·1409 3·0001 4·2009 6·1768
1·1580 1·3242 1·6462 2·1388 2·9989 4·2000 6·1763
— — 1·6164 2·1208 2·9876 4·1925 6·1713
525
T 4 Estimated jump frequencies and percentage errors (bracketed) for downward sweep
Forcing
Damping coefficient ZXXXXXXXXXXXXCXXXXXXXXXXXXV 0·01 0·03 0·1 0·3
0·1
1·1703 (−0·02)
1·1694 (−0·02)
1·1577 (−0·03)
— —
0·3
1·3302 (−0·04)
1·3296 (−0·05)
1·3236 (−0·05)
— —
1·0
1·6483 (−0·09)
1·6480 (−0·09)
1·6447 (−0·09)
1·6147 (−0·11)
3·0
2·1381 (−0·14)
2·1379 (−0·14)
2·1359 (−0·14)
2·1177 (−0·15)
10·0
2·9950 (−0·18)
2·9949 (−0·17)
2·9936 (−0·18)
2·9823 (−0·18)
30·0
4·1926 (−0·20)
4·1925 (−0·20)
4·1917 (−0·20)
4·1841 (−0·20)
100·0
6·1639 (−0·21)
6·1638 (−0·21)
6·1632 (−0·21)
6·1582 (−0·21)
sweep are given in Table 4; the errors range from −0·02 to −0·21. Again, this compares well with reference [2]. 4. This report is essentially a piece of propaganda for the harmonic balance method in its simplest guise. It is shown that the bifurcation points from the method agree very well with the ‘‘exact’’ results. There is a small shift in emphasis here from reference [2]; the bifurcation points are obtained from an exact formula based on a first order approximation. In reference [2], the points are obtained from numerical analysis of a higher order approximation. The distinction is perhaps otiose, as the final resort in this case was still a computer solution of the quintic discriminant. ACKNOWLEDGMENT
The author would like to thank Dr Mike Friswell for making his data available. 1. A. H. N and D. T. M 1979 Nonlinear Oscillations. New York: Wiley Interscience. 2. M. F and J. E. T. P 1994 Journal of Sound and Vibration 169, 261–269. The accuracy of jump frequencies in series solutions of the response of a Duffing oscillator. 3. G. B and S. M 1977 A Survey of Modern Algebra. New York: Macmillan; fourth edition. 4. S. W 1991 Mathematica: A System For Doing Mathematics by Computer. Cambridge, MA: Addison Wesley; second edition. 5. F. K 1877 Lectures on the Icosahedron and the Solution of Equations of the Fifth Degree. New York: Dover.