- Email: [email protected]
Linear and nonlinear dynamics of reciprocating engines
International Journal of Non-Linear Mechanics 38 (2003) 723 – 738
Linear and nonlinear dynamics of reciprocating engines P. Metallidis, S. Natsiavas ∗ Department of Mechanical Engineering, Aristotle University, 54006 Thessaloniki, Greece
Abstract In the &rst part of this study, dynamic models of single- and multi-cylinder reciprocating machines are presented, which may involve torsional (exibility in the crankshaft. These models take into account the dependence of the engine moment of inertia on the crankshaft rotation. In addition, both the driving and the resisting moments are expressed as functions of the crankshaft motion. This leads to dynamic models with equations of motion appearing in a strongly nonlinear form. Initially, an approximate analytical solution is presented for a linearized version of the equations of motion, by applying a suitable asymptotic method. This provides valuable insight into some aspects of the system dynamics. In the second part, numerical results are presented for linear and nonlinear engine models by applying appropriate methodologies, leading to direct determination of complete branches of steady-state response. These results illustrate the e/ect of the system parameters on its dynamics. Finally, some results are also presented in an e/ort to investigate the e/ect of engine mis&re. ? 2002 Elsevier Science Ltd. All rights reserved.
1. Introduction The rapid improvements in the level of computing power and capabilities resulted in great advances in many technological &elds, especially during the last two decades. With the help of current computing machines and codes, large scale mechanical models can be analysed in more e/ective and reliable ways. One research area which bene&ted considerably from these advances is the area of engine dynamics and vibrations, due mainly to its practical signi&cance and the large interest of the automotive industry in reducing production costs. Previous studies on this subject have focused on the dynamic behaviour of the engine (e.g. [1– 6]), some others examined thermodynamic aspects
∗ Corresponding author. Tel.: +30-31-99-6088; fax: +30-31-996029. E-mail address: [email protected] (P. Metallidis).
mostly [7], while in some more recent studies a system approach is adopted, considering appropriate sets of coupled di/erential equations governing the dynamic, the thermodynamic as well as other aspects of the engine behavior (e.g. [8,9]). The systematic study of sophisticated engine models provides useful quantitative information for a safer and more economical design of new engines. In particular, such studies help the e/orts of reducing the amount of experimental work needed for certifying new types of engines to a minimum. However, the major disadvantage of employing large scale mechanical models is that their study can most of the times be performed by numerical integration of the corresponding equations of motion. Usually, this leads to excessively high computing times and precludes a deeper and systematic investigation of the e/ects of the important system parameters. As a consequence, results from appropriate low order models are still important, because they provide valuable insight into the
0020-7462/03/$ - see front matter ? 2002 Elsevier Science Ltd. All rights reserved. PII: S 0 0 2 0 - 7 4 6 2 ( 0 1 ) 0 0 1 2 9 - 9
724
P. Metallidis, S. Natsiavas / International Journal of Non-Linear Mechanics 38 (2003) 723 – 738
e/ect of the system parameters on the overall engine function. The emphasis of the present study is placed on improving the understanding of some important features of reciprocating machine dynamics. The engine models employed are either rigid or involve torsional (exibility in the crankshaft. An essential characteristic of the modelling approach chosen is that the equivalent mass moment of inertia of the crankshaft is a function of its rotation. Moreover, both the driving moment and the resisting moment are also considered to depend on the crankshaft rotation. In particular, the former results by the pressure developed in the engine cylinders, while the form of the latter depends on the speci&c engine load. The main objective is to develop dynamic models capturing the essential behavior of non-ideal reciprocating machines, including engine mis&re. At the same time, the relatively low order of the models examined provides the (exibility necessary to get quick qualitative estimates as well as the means to check results obtained from more involved models. The organisation of the paper is as follows. The equations of motion of the engine models examined are &rst presented in the following section. Then, the appropriate form of the driving and the resisting moments is discussed. The resulting models are strongly nonlinear and as a consequence, their study cannot be performed by employing analytical techniques. Initially, the linearized version of a two degree of freedom system is investigated in Section 4. Then, typical numerical results are presented for linear and nonlinear models, illustrating the e/ect of the engine parameters on its dynamics. Response of the system in the presence of engine mis&re is also investigated. Finally, the paper concludes with a summary of the highlights and a discussion on possible extensions. 2. Engine models This section presents sets of equations of motion for simpli&ed models of reciprocating machines. For the purposes of the present study, the machines are assumed to consist of two parts. The &rst one includes the engine, while the second part includes the load. For simplicity in the presentation, a single-cylinder engine is &rst examined. Then, multi-cylinder engines are also
Fig. 1. (a) Simpli&ed model of a single-cylinder engine and (b) Members and notation of the crank-slider engine mechanism.
considered. In both cases, the resulting equations of motion include the dependence of the engine moment of inertia on the engine rotation angle. First, Fig. 1a presents a simpli&ed model of a single-cylinder engine, while Fig. 1b shows the basic elements of the slider-crank mechanism of the engine, including a rigid connecting rod, together with the corresponding notation. For simplicity, the equations of motion for the mechanism shown in Fig. 1b are derived by application of Lagrange’s equations. For this purpose, the total kinetic energy of the mechanism is &rst expressed in the form 2
TE = 12 (Ic + m1 r 2 )˙ + 12 I2 ’˙ 2 + 12 (m2 + mp )x˙2p ; where Ic is the mass moment of inertia of the crank shaft and the elements attached to it, IG is the centroidal mass moment of inertia and mr is the mass of the connecting rod, mp is the piston mass and ‘G ‘G mr ; m2 = m1 = 1 − mr ; ‘ ‘ I2 = IG − m1 ‘G2 − m2 (‘ − ‘G )2 : Next, from the kinematics of the mechanism it is easily concluded that ˙ r sin = ‘ sin ’ ⇒ ’˙ = ();
P. Metallidis, S. Natsiavas / International Journal of Non-Linear Mechanics 38 (2003) 723 – 738
725
while the piston displacement and velocity are determined from the relations ˙ xp () = r cos + ‘ cos ’ ⇒ x˙p = −r sin [1 + ()]; respectively, with () =
cos 1−
2
and
2
sin
r = : ‘
Therefore, after making the appropriate substitutions, it turns out that the kinetic energy of the engine model examined can be put in the form TE = 12 {Ic + m1 r 2 + I2 2 () + (m2 + mp )r 2 sin2 2
2
×[1 + ()]2 }˙ ≡ 12 IE ()˙ : Consequently, application of Lagrange’s equation, after selecting as the generalized coordinate and taking into account the presence of viscous friction, yields the following equation of motion 2 IE ()K + 12 IE ()˙ + c˙ = MD () − ML ;
(1)
where MD () = Fp ()xp () is the driving torque, resulting from the force Fp () exerted by the gas on the piston, while ML represents the torque delivered from the engine to the load, including the (ywheel. In addition, the equation of motion for the load is easily set up in the form IL K = ML − MR ;
(2)
where MR is the resisting moment. In general, this moment is a function of angle and its time derivatives, as explained in the following section. The simplest case occurs when the engine is rigid, which means that ≡ . In this case, the set of equations of motion (1) and (2) is combined to yield 2
[IE () + IL ]K + 12 IE ()˙ + c˙ = MD () − MR :
(3)
Apparently, the approach chosen in the present study leads to strongly nonlinear equations of motion. This is due to the fact that the equivalent moment of inertia IE , the driving torque MD and the resisting torque MR are considered to be functions of the crank rotation, in contrast to most of the previous investigations on
Fig. 2. A four-cylinder in-line engine.
the subject, where these quantities are assumed to be constants or functions of time, instead. Among other things, this puts the system examined into the class of non-ideal dynamic systems, with important consequences that will be encountered and discussed later. The modelling diMculties increase when some of the system components exhibit signi&cant (exibility. For instance, quite frequently, the torsional (exibility of the crankshaft a/ects the dynamics in a signi&cant manner. A direct consequence of this is an increase in the number of degrees of freedom required in describing the system dynamics adequately. The simplest such case is the one where the part of the shaft connecting the engine to the load (shown in Fig. 1a) is represented by a torsional spring with sti/ness k0 and a damper with damping coeMcient c0 . Application of Lagrange’s equation to this system leads to the following set of equations of motion: 2 IE ()K + 12 IE ()˙ + c0 (˙ − ˙ ) + cE ˙ + k0 ( − )
= MD ();
(4)
˙ + cL ˙ + k0 ( − ) IL K + c0 ( ˙ − ) = − MR ( ; ˙ ; K );
(5)
which is coupled and strongly nonlinear. Finally, more complex sets of equations of motion are derived for multi-cylinder engines, like the one shown in Fig. 2, by application of similar methodologies. In particular, when the engine parts are considered to be rigid, the form of the equations of motion remains the same as that of the single-cylinder engine. However, the terms representing the equivalent moment of inertia IE () and the total driving torque MD () are determined by taking into account the
726
P. Metallidis, S. Natsiavas / International Journal of Non-Linear Mechanics 38 (2003) 723 – 738
engine set-up, including the relative position of the cylinders and their &ring sequence. For instance, the equation of motion for the system shown in Fig. 2, with a four cylinder in-line engine, takes the form 4 4 2 1 K IV + IE ( − n ) + IL + IE ( − n )˙ 2 n=1
+ (c + 4cd )˙ =
where the vector x includes all the angular coordinates of the system. The dependence of the mass moment of inertia as well as of the driving and the resisting moment on the crank rotation is re(ected in the terms of the mass matrix and the right-hand side term of (7).
n=1
4
3. Driving and resisting torques MD ( − n ) − MR ;
(6)
n=1
where IV represents the moment of inertia of the engine valve-train, while n is the phase di/erence between the angular position of the &rst and the nth cylinder, so that 1 = 4 = 0 and 2 = 3 = . Likewise, the angles n indicate the phase di/erence between the &ring angle of the engine cylinders, so that 1 = 0; 2 = 3; 3 = and 4 = 2. On the other hand, when the engine crankshaft exhibits torsional (exibility, the number of degrees of freedom of the resulting dynamic model increases. For instance, the equations of motion for the system shown in Fig. 2 appear in the following form: 2
[IE (1 ) + IV ]K1 + 12 IE (1 )˙1 + c1 (˙1 − ˙2 ) + cE ˙1 + k1 (1 − 2 ) = MD (1 ); 2 IE (2 − 2 )K2 + 12 IE (2 − 2 )˙2 + c1 (2˙2 − ˙1 − ˙3 )
+ cE ˙2 + k1 (22 − 1 − 3 ) = MD (2 − 2 ) 2 IE (3 − 3 )K3 + 12 IE (3 − 3 )˙3 + c1 (2˙3 − ˙2 − ˙4 )
First, the total driving torque is caused by the force exerted by the gas mixture on each of the engine cylinders, which in turn equals the gas pressure developed at each piston times its cross-sectional area. In general, evaluation of this pressure is a complicated task, involving the solution of coupled equations resulting by simultaneous application of thermodynamic and dynamic laws [7–9]. However, for the purposes of the present study it suMces to consider a simpli&ed ideal thermodynamic cycle, corresponding to a four stroke cycle engine [7,10]. More speci&cally, during the compression and the expansion stroke, the gas mixture inside the cylinder is assumed to follow the polytropic law pV k = const:;
(8)
where p is the pressure in the cylinder and V is the volume occupied by the gas mixture. This volume is expressed in the form V () = Vc + [r + ‘ − xp ()]
B2 ; 4
Typically, such equations can be cast in the more compact matrix form
where Vc is the clearance volume and B is the bore of the piston, while the value of the exponent k is taken to be equal to 1.3 for both the compression and the expansion phases. Moreover, during the exhaust and suction phases the gas mixture follows an isobaric process. In this way, relation (8) yields the pressure developed within a cylinder as a function of the corresponding crank angle. Also, the pressure distribution is completely determined in terms of a single parameter only, that is the mean e/ective value of the pressure, pe/ . Typical pressure diagrams are shown in Fig. 3, for two di/erent values of this parameter. On the other hand, the resisting moment MR appears usually in the general form
˙ M (x)xK + C x˙ + Kx = f(x; x);
2 MR ( ˙ ) = 0 + 0 ˙ + 0 ˙
+ cE ˙3 + k1 (23 − 2 − 4 ) = MD (3 − 3 ); 2 IE (4 − 4 )K4 + 12 IE (4 − 4 )˙4 + c0 (˙4 − ˙ )
+ c1 (˙4 − ˙3 ) + cE ˙4 + k0 (4 − ) + k1 (4 − 3 ) = MD (4 − 4 ); IL K + c0 ( ˙ − ˙4 ) + cL ˙ + k0 ( − 4 ) = −MR :
(7)
(9)
P. Metallidis, S. Natsiavas / International Journal of Non-Linear Mechanics 38 (2003) 723 – 738
of motion in the longitudinal direction is expressed in the form
80
muK + 12 AcD u˙ 2 = Fw − Rx ;
peff =10 bar
60
p (bar)
727
where m is the mass, A is the frontal area and u is the longitudinal displacement of the vehicle, is the air density, cD is the aerodynamic drag coeMcient and Rx is the total rolling resistance exerted from the ground to the wheels. Typically, this force is expressed in terms of the total vehicle weight through the relation
peff =5 bar
40
(11)
Rx = !r mg;
20
(12)
where !r is the so-called rolling resistance coeMcient. Finally, if the vehicle wheels do not slip on the ground, the following kinematic conditions are true
0 0
180
360 (deg)
540
720
u(t) ˙ = rw !w (t) ⇒ uK = rw w ;
Fig. 3. Pressure diagrams for two values of the mean e/ective pressure.
where !w =
with coeMcients depending on the nature of the engine load. In particular, the constant term may represent road slope or rolling resistance of vehicles, the linear term may originate from viscous e/ects or from resistance by an electric brake or generator, while the quadratic term may represent aerodynamic or hydraulic resistance. For instance, by applying standard procedures of vehicle dynamics [11], the equation of motion for the moving part of the vehicle can be set in the form Im K + c m ˙ = M R − M w ;
(10)
where cm represents an equivalent viscous damping coeMcient, while Im = I t +
Id Iw + 2 2 2 nt nt nf
and
Mw =
rw Fw : nt nf
In particular, the symbols It ; Id and Iw represent the moment of inertia of the transmission, the driveshaft and the wheels with the axle shafts, nt and nf are the numerical ratio of the transmission and the &nal drive, rw is the wheel radius and Fw is the tractive force exerted from the ground to the driving wheels. In addition, if the vehicle follows a straight path, its equation
(13)
1 ˙ nt nf
and
w =
1 K nt nf
is the angular speed and acceleration of the driving wheels, respectively. Therefore, substitution of Eqs. (11) – (13) in (10) yields a form of the resisting moment MR , which is similar to (9), with r w !r 0 = mg; 0 = cm and n t nf 3 1 rw 0 = AcD : 2 nt nf Moreover, the remaining term, [Im + mrw2 =(nt nf )2 ] K , is added to the mass moment of inertia of the load. The form of the driving and resisting torque examined, together with the dependence of the engine moment of inertia on the crankshaft rotation result in equations of motion which are strongly nonlinear. This is profound even in their simplest version, represented by Eq. (3), while the level of diMculty increases further when some of the system components exhibit signi&cant (exibility. In general, these complications preclude the application of analytical techniques and direct their study to development and application of numerical methodologies. However, some useful information may be obtained by studying these equations of motion with approximate analytical methods, as presented in the following section.
728
P. Metallidis, S. Natsiavas / International Journal of Non-Linear Mechanics 38 (2003) 723 – 738
4. Linear system — perturbation analysis
and
Before the application of appropriate numerical methodologies in solving the nonlinear equations of motion of the systems examined, some insight into the dynamics can be gained by studying more tractable forms of them. More speci&cally, this section presents results for the linearized version of the two degree of freedom system represented by the equations of motion (4) and (5). This is done by &rst introducing the inertia quantity I0 = Ic + m1 r 2 + 12 (m2 + mp )r 2 ; together with the dimensionless parameters (m2 + mp )r 2 &= ; = &!; 2I0 JL IL JL = ; 0 = : I0 JL + 1
f21 = ,0 (’˙ 1 − ’˙ 2 ) − − (,2 + )! − !2 : In the above expressions, the selection of the damping parameters involved is based on the following formulas cn -n = √ ; &,0 = 2 0 -0 ; &2 ,1 = 2 0 -1 ; 2 k 0 I0 &2 ,2 = 2 0 -2 ; while the normalized parameters of the resisting moment appear in the form =
I2 = &'; I0
Therefore, assuming small amplitude vibration, expressed by the conditions = 0 + &’2
and introducing the dimensionless time * = !0 t, with !0 = k0 (IL + I0 )=IL I0 , the equations of motion (4) and (5) take the linearized form (14)
with 0 = +t = !*, 1 −1 1 0 2 ; K = 0 ; M= −1 1 0 JL −!2 sin 20 f0 = 0 and elements of vector f 1 with form f11 = ’K 1 cos 20 − (,0 + 2! sin 20 )’˙ 1 + ,0 ’˙ 2 − ,1 ! − 2!2 ’1 cos 20 − !2 (3 sin 30 − sin 0 )!=2 + !fp (0 ) sin 0
&2
√ 0 0 ; k 0 I0
=
1 &2 I
0
0 :
02 ‘ Fp (0 ) &k0 ∞ (fsn sin n0 + fcn cos n0 ) ; = . f0 +
fp (0 ) =
IE () = I0 [1 − & cos 2 + O(&2 )]:
M ’K + K’ = f 0 (0 ) + &f 1 (0 ; ’; ’; ˙ ’) K
=
Moreover, the de&nition of the dimensionless periodic forcing function is as follows
Then, if & accepts small values, substitution in the expression for the equivalent moment of inertia IE () implies that
= 0 + &’1 ;
02 0 ; & 2 k0
n=1
where the parameter . represents a normalized mean e/ective pressure inside the engine cylinders. By employing standard modal analysis techniques, Eq. (14) is brought into the alternative form, uK + /u = h0 + &h1 ;
(15)
through the application of the coordinate transformation ’(*) = 1u(*); where 1 is the modal matrix of the linearized system, / is a diagonal matrix with diagonal elements equal to 0 and 1, while the elements of vectors h0 and h1 appear in the form hm0 = m0 sin 20 and hm1 = m1 u˙ 1 + m2 u˙ 2 + (m1 uK 1 + m2 uK 2 ) cos 20 + (m1 u˙ 1 + m2 u˙ 2 ) sin 20 + (,m1 u1 + ,m2 u2 ) ×cos 20 +
∞ (&mn sin n0 +-mn cos n0 )+2m : n=1
P. Metallidis, S. Natsiavas / International Journal of Non-Linear Mechanics 38 (2003) 723 – 738
729
Next, application of the multiple time scales method [12] leads to approximate analytical solutions of Eq. (15), with form
ined. For instance, for a bounded response it turns out (from the &rst of them) that the following condition must be satis&ed
u(*; &) = u 0 (*0 ; *1 ) + &u 1 (*0 ; *1 ) + O(&2 );
21 = 0:
(16)
where *0 =* and *1 =&*. More speci&cally, substituting (16) in (15) and collecting terms of order & and &2 leads to the following linear sets of equations: D02 u 0 + /u 0 = h0 (*); D02 u 1 + /u 1 = −2D0 D1 u 0 + h1 (*; u 0 ; D0 u 0 ; D02 u 0 ); (17) where Dn = @=@*n . Then, the zero order solution is determined in the form 10 u10 = c2 − sin 2!* and 4!2 20 sin 2!* + a2 (*1 ) cos[* + 2 (*1 )]; u20 = 1 − 4!2 (18) where c2 is an arbitrary constant. Substitution of u10 and u20 from (18) into (17) yields two di/erential equations for u11 and u21 . These equations have the following form: D02 u11
= 1 A2 e
i(2!+1)*
+ 1 AP 2 e
i(2!−1)*
+ ! 1 A2 e
i*
1 + '1 e4i!* + .1 e2i!t + 21 2
and D02 u21 +u21 = −2iA2 ei* + 2 A2 ei(2!+1)* + 2 AP 2 ei(2!−1)* + .2 e
which is a manifestation of an energy balance. In addition, depending on the value of the frequency !, several resonances are possible to occur at &rst order. Namely, the system exhibits such resonances at values of ! near 0, 14 ; 12 , and 1. However, the &rst and third of these resonances cannot be analyzed e/ectively for the parameter scaling chosen, as it is obvious from (18). Moreover, among the remaining two resonances, the last one presents more interest. In particular, this corresponds to the case where a fundamental parametric resonance and a primary external resonance are activated simultaneously, when the frequency ! is close to the undamped linear natural frequency of the dynamical system examined, which means that ! = 1 + &6:
(20)
In this case, the resulting solvability condition resulting from the di/erential equation for u21 can eventually be expressed by the relation
Therefore, if the solution amplitude is decomposed in the polar form
n=1
+ ! 2 A2 e + ' 2 e
(19)
+ 12 (-21 − i&21 )ei6*1 = 0:
∞
4i!*
+ ( + ,1 + ,2 )! + !2 = 12 !.fs1 ;
i 22 A2 − 2iA2 − 12 (22 + 22 − ,22 )AP 2 e2i6*1
1 + (-1n − i&1n )ein!* + cc 2
i*
Expressed in terms of the original system parameters this condition reappears in the form
2i!t
∞
1 1 + 22 + (-2n − i&2n )ein!* + cc; 2 2 n=1
where the new constant quantities n ; n ; !n ; 'n ; .n (n = 1; 2) are complicated but known functions of the system parameters, while cc denotes complex conjugate of the preceding terms. The last two di/erential equations illustrate some important characteristic features of the system exam-
A2 (*1 ) = 12 a2 (*1 )ei2 (*1 ) ; separation of the real from the imaginary component in the last solvability condition leads eventually to the following slow-(ow system for the amplitude and phase of vibration: a2 = q1 a2 + q2 a2 sin 22 + q3 sin 2 − q4 cos 2 ; a2 2 = 6a2 + q2 a2 cos 22 + q4 sin 2 + q3 cos 2
(21) (22)
with 2 =6*1 −2 . Constant solutions (a20 ; 20 ) of these equations are then obtained by imposing the conditions
730
P. Metallidis, S. Natsiavas / International Journal of Non-Linear Mechanics 38 (2003) 723 – 738
a20 = 0 = 20 . By employing relations (18), (20) and the de&nition of phase 2 it can easily be shown that these solutions correspond to 10 u1 = c2 − sin 2!* + O(&) and 4!2 20 sin 2!* + a20 cos (!* − 20 ) + O(&); u20 = 1 − 4!2 which represent periodic motions of the mechanical system examined, with fundamental frequency equal to !, whose value is determined from condition (19). Determination of these motions starts by &rst solving the resulting algebraic equations for sin 20 and cos 20 . Then, employing the trigonometric identity sin2 20 + cos2 20 = 1 and performing lengthy algebraic manipulations leads &rst to a new relation with form p1 sin 220 + p2 cos 220 = p3 + 1=a220 :
(23)
Next, solving the same equations for sin 220 and cos 220 and following a similar path leads to another relation with form p4 sin 220 + p5 cos 220 = p6 + p7 a220 − 1=4a220 : (24) Eqs. (23) and (24) constitute a linear set of algebraic equations for sin 220 and cos 220 , whose solution can be expressed in the form sin 220 = p8 + p9 a220 + p0 =a220
and
cos 220 = p10 + p11 a220 + p12 =a220 ;
(25)
where all the constant coeMcients pn are known functions of the system parameters. Elimination of the trigonometric terms leads to a fourth degree polynomial for a220 , whose solution provides the amplitude a20 . Then, back substitution in (25) determines the corresponding phase 20 . Finally, the stability properties of these solutions are determined by applying the classical method of linearizarion in (21) and (22) [12,13]. 5. Numerical results In this section, numerical results are presented in an e/ort to illustrate some characteristic features
related to the dynamic behavior of reciprocating engines mounted on road vehicles. First, a two degree of freedom system is examined, with emphasis put in comparing the results of the linearized and the fully nonlinear version of the equations of motion. Then, dynamics of a four cylinder engine nonlinear model is investigated, with attention focused on analysing the e/ect of the system parameters, including engine mis&re. 5.1. Results for two degree of freedom models First, typical numerical results are presented in Fig. 4 for two degree of freedom systems, like the one examined in the previous section. Fig. 4a presents characteristic frequency-response diagrams. In particular, the results represented by the curves denoted with 1 were obtained by solving the slow-(ow Eqs. (21) and (22) for constant solutions in frequency ranges speci&ed by condition (20). The amplitudes of the relative phase Q’ = ’1 − ’2 are plotted versus the frequency detuning parameter 6, for a combination of parameters leading to & = 0:05. For comparison purposes, the curves denoted by 2 represent branches of periodic solutions obtained by applying suitable numerical methodologies [13,14] directly to the nonlinear set of equations of motion (4) and (5). On the other hand, the response diagram denoted by 3 corresponds to branches of periodic solutions obtained for the same set of parameters, by solving Eq. (14) after setting IE () = I0 , which leads to a constant coeMcient linearized system. In all these cases, the broken parts of the curves represent branches of unstable periodic motions. Inspection of Fig. 4a, reveals that the results obtained for the nonlinear set of Eqs. (4), (5) and those obtained by the multiple time scales method appear to be quite di/erent. Besides the di/erences near the 1=3 and 1=4 superharmonic resonances, which are expected to occur, substantial di/erences appear even in the vicinity of the ! ≈ 1 resonance (the predictions of the constant coeMcient model are worse). The explanation for this discrepancy is as follows. In the linear case, the average spin speed ! of the system is not known a priori. Application of the multiple time scales method determined this speed through the condition (19). However, this condition is correct up to &rst order only and is therefore incomplete. For instance, it
P. Metallidis, S. Natsiavas / International Journal of Non-Linear Mechanics 38 (2003) 723 – 738
731
uating the loss of energy due to the torsional damping yields the new expression 1 + ( + ,1 + ,2 )! + !2 + &,0 ! 2 2 4 1 4! a = !.fs1 ; + × 20 2 2 2 (4! − 1) 2 0
Fig. 4. (a) Response diagrams for a linearized and a nonlinear engine model (1: multiple scales with (19), 2: nonlinear model, 3: constant coeMcient model, 4: multiple scales with (26)) at &=0:05 and (b) Average crankshaft speed versus the mean e/ective pressure parameter .ˆ for the same models.
represents an energy balance involving the rigid body component of the motion, but it does not involve the loss of energy due to the torsional damping present in the relative rotation of the system inertia elements. By employing the approximate expressions (18) and eval-
(26)
in place of (19). Then, if the average spin speed ! of the linear system is selected from the last expression, the results obtained are found to be virtually coincident with those obtained by employing the nonlinear equations of motion in the vicinity of the ! ≈ 1 resonance. This is illustrated by the response diagram denoted by 4 in Fig. 4a, which represents branches of periodic solutions obtained for the linearized system by application of the multiple time scales method and after selecting ! from the last expression of energy balance. The above results indicate that caution should be exercised when analyzing the linearized version of the equations of motion. Another important discrepancy between the predictions of the &rst order multiple time scales method developed in the previous section and the results obtained directly from the original nonlinear equations of motion is depicted in a better way in Fig. 4b. This &gure presents the plot of the average spin speed of the same system, as a function of the normalized pressure amplitude parameter .ˆ ≡ .fs1 . Here, the results obtained by the application of the multiple time scales method in conjunction with relation (19) are quite di/erent than those resulting for the nonlinear set of equations of motion, as in the previous &gure. However, the agreement between the results improves again when relation (26) is employed in determining the response frequency !. From Fig. 4b it is apparent that near the resonance examined in the previous section, as well as the other (superharmonic) resonances, the nonlinear system exhibits the well-known Sommerfeld phenomenon, which is expected to occur in the system considered, because it is non-ideal [12,15]. As is typical, the appearance of Sommerfeld phenomena may lead to response with considerable deviations in the mean spin speed, for the same mean e/ective pressure value. This is veri&ed by the response histories shown in Fig. 5. More speci&cally, Fig. 5a shows response histories obtained by direct integration of Eqs. (4) and (5) at
732
P. Metallidis, S. Natsiavas / International Journal of Non-Linear Mechanics 38 (2003) 723 – 738
Fig. 5. Response histories at .ˆ = 0:05, illustrating appearance of Sommerfeld phenomena (a) nonlinear model (b) linearized model with relation (26).
ˆ .=0:05 and for three di/erent sets of initial conditions. Both the amplitude and the fundamental frequency of all these motion are quite di/erent. For comparison purposes, Fig. 5b shows the corresponding response histories determined by obtaining the appropriate initial conditions from the multiple time scales method in conjunction with relation (26). The agreement between the predictions of the two models is satisfactory
Fig. 6. (a) Response diagrams for a linearized and a nonlinear engine model (1: multiple scales with (19), 2: nonlinear model, 3: constant coeMcient model, 4: multiple scales with (26)) at & = 0:24. (b) Average crankshaft speed versus the mean e/ective pressure parameter .ˆ for the same models.
even for the smaller amplitude (and frequency) motion, which occurs away from the ! ≈ 1 resonance. The accuracy of the predictions of the multiple scales method depends strongly on the value of the parameter &. Next, the response diagrams of Fig. 6a illustrate the main e/ects caused by increasing the
P. Metallidis, S. Natsiavas / International Journal of Non-Linear Mechanics 38 (2003) 723 – 738
733
the system dynamics. More speci&cally, the results of Fig. 7a were obtained by applying the multiple scales method for & = 0:05 and demonstrate the e/ect of the damping ratio -0 of the shaft. An increase in the value of this papameter results in a reduction of the nonlinearity e/ects near the ! ≈ 1 range. Likewise, the results shown in Fig. 7b illustrate the e/ect of the normalized quadratic resistance term . Here, an increase in the value of this parameter accentuates the nonlinearity e/ects around ! ≈ 1. 5.2. Results for a four cylinder engine nonlinear model
Fig. 7. (a) E/ect of the cranckshaft daping ratio -0 . (b) E/ect of the aerodynamic resistance coeMcient .
value of &. All these diagrams were obtained for a combination of parameters leading to & = 0:24, keeping the same numbering as in Fig. 4. Clearly, increasing the value of the parameter & magni&es the di/erences between the predictions of the nonlinear and the linear models. Moreover, the system exhibits higher response amplitudes around ! ≈ 1 than close to ! ≈ 1=2. Finally, the diagrams presented in Fig. 7 depict the e/ect of two other important technical parameters on
In this section, numerical results are presented for a four cylinder in-line engine, like the one shown in Fig. 2, with parameter values characterizing a typical road vehicle. In all cases examined, the engine model is nonlinear and its dynamics is governed by Eq. (6) or (7), depending on whether the system exhibits torsional sti/ness. Therefore, the numerical results presented next were obtained by applying suitable numerical methodologies [13,14]. First, Fig. 8a shows results in the form of diagrams, where the independent parameter pe/ represents the mean e/ective pressure developed within the combustion chamber of the engine cylinders [7], while the absicca indicates the resulting mean component of the crankshaft spin speed. These results were obtained by assuming that the engine is rigid and the system is in a steady-state condition. The continuous line was obtained by taking into account the exact dependence of the engine mass moment of inertia IE on the crankshaft angular position, while the broken line corresponds to the case where IE is chosen as constant, that is equal to its mean value I0 . Finally, the dots were obtained by applying the energy balance condition, expressed by (19). The variations detected in Fig. 8a on the average spin speed of the crankshaft due to the dependence of the engine moment of inertia on the crankshaft rotation are negligible. This can be justi&ed by the fact that the equivalent moment of inertia of the engine load (including the vehicle and the (ywheel) is much greater than the moment of inertia of the engine. For instance, Fig. 8b shows the history of the angular acceleration obtained for a typical engine load and for no engine load at all (thicker line). This same e/ect
734
P. Metallidis, S. Natsiavas / International Journal of Non-Linear Mechanics 38 (2003) 723 – 738
Fig. 8. (a) Average engine speed versus mean e/ective pressure. (b) Engine angular acceleration during loaded (thin line) and idle (thick line) function. (c) Angular velocity of loaded engine with variable (continuous line) and constant (broken line) moment of inertia. (d) Angular velocity of idle engine with variable (continuous line) and constant (broken line) moment of inertia.
is also apparent, though less pronounced, in the angular velocity. More speci&cally, the continuous lines in Figs. 8c and 8d represent the histories of the angular velocity obtained for the same engine load and for no load, respectively. In addition, the broken lines in Figs. 8c and 8d depict the angular velocity histories when IE is assumed to remain constant. During these calculations, it was also found that these di/erences become more pronounced for larger values of pe/ . As it is apparent from Fig. 8a, there exists a minimum e/ective pressure level required in order for the engine to start rotating. Similar trends are observed next in Fig. 9, which illustrates e/ects of the vehicle resistance and the crankshaft (exibility parameters. More speci&cally, Fig. 9a shows similar diagrams for the same rigid engine and three levels of the external damping parameter c ≡ c1 = c2 . Clearly, increasing the damping level causes a reduction in the mean
spin speed for the same level of gas pressure. Likewise, Fig. 9b depicts the e/ect of the rolling resistance parameter, while Fig. 9c depicts the e/ect of the aerodynamic resistance term. Finally, Fig. 9d shows the e/ect of the crankshaft torsional sti/ness parameter. In the last case, the di/erences in the results presented on the (pe/ ; !) plane do not seem to be signi&cant. However, it was found that a decrease in the value of the torsional sti/ness parameter or in the damping level encourages the more frequent appearance of Sommerfeld phenomena, in several frequency ranges. For the case with k0 =0:3 MN m=rad, one such area is shown magni&ed in the inset of Fig. 9d. Finally, Fig. 10 illustrates the e/ect of the internal damping parameter -0 , for a system with crankshaft torsional sti/ness parameter k0 =0:6 MN m=rad. Here, the maximum absolute value of the response parameter Q = 4 − F (representing the di/erence between
P. Metallidis, S. Natsiavas / International Journal of Non-Linear Mechanics 38 (2003) 723 – 738
735
Fig. 9. Dependence of the average engine speed on: (a) the external damping parameter; (b) the rolling resistance parameter; (c) the aerodynamic resistance term and (d) the crankshaft torsional sti/ness parameter.
the rotation angle of the (ywheel and the adjacent crank) is shown as a function of the e/ective pressure parameter pe/ , for two values of -0 . Clearly, a decrease in the internal damping level results in bigger response amplitudes. As a consequence, the nonlinear e/ects become more pronounced, which is demonstrated with the appearance of unstable solution branches for the smaller damping value (Fig. 10b). 5.3. E5ects from engine mis6re In obtaining the results presented in the previous section, it was assumed that the pressure developed within each engine cylinder takes values, which are in complete agreement with the theoretical predictions. However, due to various manufacturing or operational factors, there may appear considerable deviations between the measured and the theoretical pressure pattern inside the cylinders. With this in mind, the results
of the present section shed some light in the e/ect of unavoidable engine mis&re on the engine dynamics. Initially, Fig. 11a presents the spin speed history of the rigid engine examined in Section 5.2, when there is a complete mis&re in a single cylinder. These results were obtained for pe/ = 8 bar, which in the case of the healthy engine led to the history shown in Fig. 8c. Examination of these signals reveals &rst that there is a drop in the average spin speed, in all faulty cases. However, besides a phase shift, the recorded signals appear to be identical. Likewise, Fig. 11b presents results for the same system, when there is a partial mis&re in both the &rst and the third cylinder or in all four cylinders of the engine (the cylinder numbering is indicated in Fig. 2). More speci&cally, a 50% mis&re implies that the mean e/ective pressure developed within the combustion chamber of the faulty cylinder reaches only half of the expected value. Again, there are some changes in the signal pro&les, but the average
736
P. Metallidis, S. Natsiavas / International Journal of Non-Linear Mechanics 38 (2003) 723 – 738
Fig. 10. E/ect of the torsional damping parameter on the engine response: (a) -0 = 0:01 and (b) -0 = 0:001.
spin speed remains the same as before. In fact, the average engine spin speed is found to be a function of the total mis&re but it does not seem to depend on its exact distribution in the cylinders. This is expected to cause diMculties in the e/orts to locate faulty cylinders and determine the level of their mis&re, based on information obtained from experimental measurement of the engine spin speed (e.g., [16 –18]). A more noticeable di/erence is observed in the response, when the torsional sti/ness of the crankshaft
Fig. 11. Angular velocity of a rigid engine with: (a) a complete mis&re in a single cylinder and (b) a partial mis&re in the &rst and the third or all the cylinders.
is &nite. For instance, Fig. 12a presents results for all the cylinders and the (ywheel of the same system, but with k0 = 0:1 MN m=rad; -0 = 0:005, when all cylinders are healthy. Here, the recorded spin speed history appears to be di/erent from cylinder to cylinder, due to the torsional (exibility. Moreover, these signals are also di/erent than the corresponding signal of the rigid engine, shown in Fig. 8c. Likewise, Fig. 12b presents results for the same engine, when there is a
P. Metallidis, S. Natsiavas / International Journal of Non-Linear Mechanics 38 (2003) 723 – 738
737
6. Summary and extensions
Fig. 12. Angular velocity of a four cylinder (exible engine with: (a) all cylinders healthy and (b) a 50% mis&re in the &rst and the third cylinder.
50% mis&re in both the &rst and the third cylinder. Again, a comparison with the corresponding results of the healthy system, shown in Fig. 12a, indicates a drop in the average spin speed. From the same results it was also observed that the (ywheel spin speed of the rigid engine is slightly larger than the speed of the (exible engine, where in addition, it also depends on the fault distribution.
Dynamics of simpli&ed models of reciprocating machines was investigated. First, sets of equations of motion were set up for single- and multi-cylinder engines, taking into account the dependence of the engine moment of inertia on the crankshaft rotation. The diMculty level was increased even further, by assuming that both the driving and the resisting moments depend also on the crankshaft rotation. These led to more accurate models but at the same time they introduced strong nonlinearities in the equations of motion. As a result, an approximate analysis was &rst presented for predicting the dynamic behavior of a linearized engine model, which provided some useful insight into the system dynamics. In addition, the equations of motion of the nonlinear models were solved by applying appropriate numerical methodologies. The numerical results were obtained for typical vehicle engines in steady-state conditions. Emphasis was put on investigating the e/ect of the torsional sti/ness and damping parameters of the crankshaft as well as the characteristics of the resisting moment. Among other things, the results indicated the appearance of Sommerfeld phenomena, related to the non-ideal nature of the dynamical systems examined. Also, the dependence of the engine moment of inertia on the crankshaft rotation was found to have a weak e/ect on the average spin speed but it was shown to a/ect more signi&cantly the speed (uctuations. In the &nal part, some issues related to engine mis&re were examined, which was related to a drop in the mean e/ective pressure in the engine cylinders. The models employed in the present study can be improved in several respects. For instance, the pressure developed in the combustion chamber of the engine cylinders can be predicted in a more accurate manner by considering more involved thermodynamic models of the combustion process as well as their coupling with the engine dynamics. In a similar way, the torque resistance from a cooperating electric generator can be modelled by including details about its electric circuit and its coupling with the engine. On the other hand, suitable numerical methodologies are required before the level of mis&re in each cylinder can be predicted in a systematic and computationally eMcient way, handling both measurement and model errors [19,20]. Finally, the mechanical model of the
738
P. Metallidis, S. Natsiavas / International Journal of Non-Linear Mechanics 38 (2003) 723 – 738
engine can also be improved, by taking into account e/ects arising from the (exibility in the bearings, the connecting rod and the supports of the engine frame. Acknowledgements This research was partially supported by a grant from the General Secretariat of Research and Technology, Greek Ministry of Development, through the PENED99 program (reference number: 99EQ 580). References [1] H.H. Denman, Exact solution for the rigid slider-crank mechanism with gas pressure, Mech. Mach. Theory 26 (1991) 415–420. [2] S.-R. Hsieh, S.W. Shaw, The dynamic stability and nonlinear resonance of a (exible connecting rod: continuous parameter model, Nonlinear Dyn. 4 (1993) 573–603. [3] Y. Shiao, C.-H. Pan, J.J. Moskwa, Advanced dynamic spark ignition engine modelling for diagnostics and control, Internat. J. Veh. Design 15 (1994) 578–596. [4] H. Okamura, A. Shinno, T. Yamanaka, A. Suzuki, K. Sogabe, Simple modeling and analysis for crankshaft three-dimensional vibrations, Part 1: background and application to free vibrations, J. Vib. Acoust. 117 (1995) 70–79. [5] H. Ashra&uon, A.M. Whitman, Asymptotic analysis of the torsional vibration in reciprocating machinery, J. Vib. Acoust. 118 (1996) 485–490. [6] E. Brusa, C. Delprete, G. Genta, Torsional vibration of crankshafts: e/ects of non-constant moments of inertia, J. Sound Vib. 205 (1997) 135–150. [7] J.B. Heywood, Internal Combustion Engine Fundamentals, McGraw-Hill, New York, 1988.
[8] L. Dinca, T. Aldemir, G. Rizzoni, Fault detection and identi&cation in dynamic systems with noisy data and parameter=modelling uncertainties, Reliab. Eng. System Saf. 65 (1999) 17–28. [9] S.E. Lyshevski, Energy conversion and optimal energy management in diesel-electric drivetrains of hybrid-electric vehicles, Energy Conversion Manage. 41 (2000) 13–24. [10] J.E. Shigley, Dynamic Analysis of Machines, McGraw-Hill, New York, 1961. [11] T.D. Gillespie, Fundamentals of Vehicle Dynamics, Society of Automotive Engineers, Warrendale, PA, 1992. [12] A.H. Nayfeh, D.T. Mook, Nonlinear Oscillations, Wiley-Interscience, New York, 1979. [13] A.H. Nayfeh, B. Balachandran, Applied Nonlinear Dynamics, Wiley-Interscience, New York, 1995. [14] E. Doedel, AUTO: Software for Continuation and Bifurcation Problems in Ordinary Di/erential Equations, California Institute of Technology, Pasadena, CA, 1986. [15] J. Wauer, P. Buhrle, Dynamics of a (exible slider-crank mechanism driven by a non-ideal source of energy, Nonlinear Dyn. 13 (1997) 221–242. [16] S.J. Citron, J.E. O’Higgins, L.Y. Chen, Cylinder by cylinder engine pressure and pressure torque waveform determination utilizing speed (uctuations. SAE Paper No. 890486, 1989, pp. 933–947. [17] D. Taraza, N.A. Henein, W. Bryzik, Diesel engine diagnosis based on analysis of the crankshaft’s speed variations, SAE Paper No. 982540, 1998, pp. 2183–2195. [18] F.T. Connolly, G. Rizzoni, Real time estimation of engine torque for the detection of engine mis&res, ASME J. Dyn. Systems Meas. Control 116 (1994) 675–686. [19] J.B. Roberts, P.D. Spanos, Random Vibration and Statistical Linearization, Wiley, New York, 1990. [20] L.S. Katafygiotis, C. Papadimitriou, H.F. Lam, 1998. A Probabilistic Approach to Structural Model Updating, Internat. J. Soil Dyn. Earthquake Eng. 17 (1998) 495–507.
Linear and nonlinear dynamics of reciprocating engines P. Metallidis, S. Natsiavas ∗ Department of Mechanical Engineering, Aristotle University, 54006 Thessaloniki, Greece
Abstract In the &rst part of this study, dynamic models of single- and multi-cylinder reciprocating machines are presented, which may involve torsional (exibility in the crankshaft. These models take into account the dependence of the engine moment of inertia on the crankshaft rotation. In addition, both the driving and the resisting moments are expressed as functions of the crankshaft motion. This leads to dynamic models with equations of motion appearing in a strongly nonlinear form. Initially, an approximate analytical solution is presented for a linearized version of the equations of motion, by applying a suitable asymptotic method. This provides valuable insight into some aspects of the system dynamics. In the second part, numerical results are presented for linear and nonlinear engine models by applying appropriate methodologies, leading to direct determination of complete branches of steady-state response. These results illustrate the e/ect of the system parameters on its dynamics. Finally, some results are also presented in an e/ort to investigate the e/ect of engine mis&re. ? 2002 Elsevier Science Ltd. All rights reserved.
1. Introduction The rapid improvements in the level of computing power and capabilities resulted in great advances in many technological &elds, especially during the last two decades. With the help of current computing machines and codes, large scale mechanical models can be analysed in more e/ective and reliable ways. One research area which bene&ted considerably from these advances is the area of engine dynamics and vibrations, due mainly to its practical signi&cance and the large interest of the automotive industry in reducing production costs. Previous studies on this subject have focused on the dynamic behaviour of the engine (e.g. [1– 6]), some others examined thermodynamic aspects
∗ Corresponding author. Tel.: +30-31-99-6088; fax: +30-31-996029. E-mail address: [email protected] (P. Metallidis).
mostly [7], while in some more recent studies a system approach is adopted, considering appropriate sets of coupled di/erential equations governing the dynamic, the thermodynamic as well as other aspects of the engine behavior (e.g. [8,9]). The systematic study of sophisticated engine models provides useful quantitative information for a safer and more economical design of new engines. In particular, such studies help the e/orts of reducing the amount of experimental work needed for certifying new types of engines to a minimum. However, the major disadvantage of employing large scale mechanical models is that their study can most of the times be performed by numerical integration of the corresponding equations of motion. Usually, this leads to excessively high computing times and precludes a deeper and systematic investigation of the e/ects of the important system parameters. As a consequence, results from appropriate low order models are still important, because they provide valuable insight into the
0020-7462/03/$ - see front matter ? 2002 Elsevier Science Ltd. All rights reserved. PII: S 0 0 2 0 - 7 4 6 2 ( 0 1 ) 0 0 1 2 9 - 9
724
P. Metallidis, S. Natsiavas / International Journal of Non-Linear Mechanics 38 (2003) 723 – 738
e/ect of the system parameters on the overall engine function. The emphasis of the present study is placed on improving the understanding of some important features of reciprocating machine dynamics. The engine models employed are either rigid or involve torsional (exibility in the crankshaft. An essential characteristic of the modelling approach chosen is that the equivalent mass moment of inertia of the crankshaft is a function of its rotation. Moreover, both the driving moment and the resisting moment are also considered to depend on the crankshaft rotation. In particular, the former results by the pressure developed in the engine cylinders, while the form of the latter depends on the speci&c engine load. The main objective is to develop dynamic models capturing the essential behavior of non-ideal reciprocating machines, including engine mis&re. At the same time, the relatively low order of the models examined provides the (exibility necessary to get quick qualitative estimates as well as the means to check results obtained from more involved models. The organisation of the paper is as follows. The equations of motion of the engine models examined are &rst presented in the following section. Then, the appropriate form of the driving and the resisting moments is discussed. The resulting models are strongly nonlinear and as a consequence, their study cannot be performed by employing analytical techniques. Initially, the linearized version of a two degree of freedom system is investigated in Section 4. Then, typical numerical results are presented for linear and nonlinear models, illustrating the e/ect of the engine parameters on its dynamics. Response of the system in the presence of engine mis&re is also investigated. Finally, the paper concludes with a summary of the highlights and a discussion on possible extensions. 2. Engine models This section presents sets of equations of motion for simpli&ed models of reciprocating machines. For the purposes of the present study, the machines are assumed to consist of two parts. The &rst one includes the engine, while the second part includes the load. For simplicity in the presentation, a single-cylinder engine is &rst examined. Then, multi-cylinder engines are also
Fig. 1. (a) Simpli&ed model of a single-cylinder engine and (b) Members and notation of the crank-slider engine mechanism.
considered. In both cases, the resulting equations of motion include the dependence of the engine moment of inertia on the engine rotation angle. First, Fig. 1a presents a simpli&ed model of a single-cylinder engine, while Fig. 1b shows the basic elements of the slider-crank mechanism of the engine, including a rigid connecting rod, together with the corresponding notation. For simplicity, the equations of motion for the mechanism shown in Fig. 1b are derived by application of Lagrange’s equations. For this purpose, the total kinetic energy of the mechanism is &rst expressed in the form 2
TE = 12 (Ic + m1 r 2 )˙ + 12 I2 ’˙ 2 + 12 (m2 + mp )x˙2p ; where Ic is the mass moment of inertia of the crank shaft and the elements attached to it, IG is the centroidal mass moment of inertia and mr is the mass of the connecting rod, mp is the piston mass and ‘G ‘G mr ; m2 = m1 = 1 − mr ; ‘ ‘ I2 = IG − m1 ‘G2 − m2 (‘ − ‘G )2 : Next, from the kinematics of the mechanism it is easily concluded that ˙ r sin = ‘ sin ’ ⇒ ’˙ = ();
P. Metallidis, S. Natsiavas / International Journal of Non-Linear Mechanics 38 (2003) 723 – 738
725
while the piston displacement and velocity are determined from the relations ˙ xp () = r cos + ‘ cos ’ ⇒ x˙p = −r sin [1 + ()]; respectively, with () =
cos 1−
2
and
2
sin
r = : ‘
Therefore, after making the appropriate substitutions, it turns out that the kinetic energy of the engine model examined can be put in the form TE = 12 {Ic + m1 r 2 + I2 2 () + (m2 + mp )r 2 sin2 2
2
×[1 + ()]2 }˙ ≡ 12 IE ()˙ : Consequently, application of Lagrange’s equation, after selecting as the generalized coordinate and taking into account the presence of viscous friction, yields the following equation of motion 2 IE ()K + 12 IE ()˙ + c˙ = MD () − ML ;
(1)
where MD () = Fp ()xp () is the driving torque, resulting from the force Fp () exerted by the gas on the piston, while ML represents the torque delivered from the engine to the load, including the (ywheel. In addition, the equation of motion for the load is easily set up in the form IL K = ML − MR ;
(2)
where MR is the resisting moment. In general, this moment is a function of angle and its time derivatives, as explained in the following section. The simplest case occurs when the engine is rigid, which means that ≡ . In this case, the set of equations of motion (1) and (2) is combined to yield 2
[IE () + IL ]K + 12 IE ()˙ + c˙ = MD () − MR :
(3)
Apparently, the approach chosen in the present study leads to strongly nonlinear equations of motion. This is due to the fact that the equivalent moment of inertia IE , the driving torque MD and the resisting torque MR are considered to be functions of the crank rotation, in contrast to most of the previous investigations on
Fig. 2. A four-cylinder in-line engine.
the subject, where these quantities are assumed to be constants or functions of time, instead. Among other things, this puts the system examined into the class of non-ideal dynamic systems, with important consequences that will be encountered and discussed later. The modelling diMculties increase when some of the system components exhibit signi&cant (exibility. For instance, quite frequently, the torsional (exibility of the crankshaft a/ects the dynamics in a signi&cant manner. A direct consequence of this is an increase in the number of degrees of freedom required in describing the system dynamics adequately. The simplest such case is the one where the part of the shaft connecting the engine to the load (shown in Fig. 1a) is represented by a torsional spring with sti/ness k0 and a damper with damping coeMcient c0 . Application of Lagrange’s equation to this system leads to the following set of equations of motion: 2 IE ()K + 12 IE ()˙ + c0 (˙ − ˙ ) + cE ˙ + k0 ( − )
= MD ();
(4)
˙ + cL ˙ + k0 ( − ) IL K + c0 ( ˙ − ) = − MR ( ; ˙ ; K );
(5)
which is coupled and strongly nonlinear. Finally, more complex sets of equations of motion are derived for multi-cylinder engines, like the one shown in Fig. 2, by application of similar methodologies. In particular, when the engine parts are considered to be rigid, the form of the equations of motion remains the same as that of the single-cylinder engine. However, the terms representing the equivalent moment of inertia IE () and the total driving torque MD () are determined by taking into account the
726
P. Metallidis, S. Natsiavas / International Journal of Non-Linear Mechanics 38 (2003) 723 – 738
engine set-up, including the relative position of the cylinders and their &ring sequence. For instance, the equation of motion for the system shown in Fig. 2, with a four cylinder in-line engine, takes the form 4 4 2 1 K IV + IE ( − n ) + IL + IE ( − n )˙ 2 n=1
+ (c + 4cd )˙ =
where the vector x includes all the angular coordinates of the system. The dependence of the mass moment of inertia as well as of the driving and the resisting moment on the crank rotation is re(ected in the terms of the mass matrix and the right-hand side term of (7).
n=1
4
3. Driving and resisting torques MD ( − n ) − MR ;
(6)
n=1
where IV represents the moment of inertia of the engine valve-train, while n is the phase di/erence between the angular position of the &rst and the nth cylinder, so that 1 = 4 = 0 and 2 = 3 = . Likewise, the angles n indicate the phase di/erence between the &ring angle of the engine cylinders, so that 1 = 0; 2 = 3; 3 = and 4 = 2. On the other hand, when the engine crankshaft exhibits torsional (exibility, the number of degrees of freedom of the resulting dynamic model increases. For instance, the equations of motion for the system shown in Fig. 2 appear in the following form: 2
[IE (1 ) + IV ]K1 + 12 IE (1 )˙1 + c1 (˙1 − ˙2 ) + cE ˙1 + k1 (1 − 2 ) = MD (1 ); 2 IE (2 − 2 )K2 + 12 IE (2 − 2 )˙2 + c1 (2˙2 − ˙1 − ˙3 )
+ cE ˙2 + k1 (22 − 1 − 3 ) = MD (2 − 2 ) 2 IE (3 − 3 )K3 + 12 IE (3 − 3 )˙3 + c1 (2˙3 − ˙2 − ˙4 )
First, the total driving torque is caused by the force exerted by the gas mixture on each of the engine cylinders, which in turn equals the gas pressure developed at each piston times its cross-sectional area. In general, evaluation of this pressure is a complicated task, involving the solution of coupled equations resulting by simultaneous application of thermodynamic and dynamic laws [7–9]. However, for the purposes of the present study it suMces to consider a simpli&ed ideal thermodynamic cycle, corresponding to a four stroke cycle engine [7,10]. More speci&cally, during the compression and the expansion stroke, the gas mixture inside the cylinder is assumed to follow the polytropic law pV k = const:;
(8)
where p is the pressure in the cylinder and V is the volume occupied by the gas mixture. This volume is expressed in the form V () = Vc + [r + ‘ − xp ()]
B2 ; 4
Typically, such equations can be cast in the more compact matrix form
where Vc is the clearance volume and B is the bore of the piston, while the value of the exponent k is taken to be equal to 1.3 for both the compression and the expansion phases. Moreover, during the exhaust and suction phases the gas mixture follows an isobaric process. In this way, relation (8) yields the pressure developed within a cylinder as a function of the corresponding crank angle. Also, the pressure distribution is completely determined in terms of a single parameter only, that is the mean e/ective value of the pressure, pe/ . Typical pressure diagrams are shown in Fig. 3, for two di/erent values of this parameter. On the other hand, the resisting moment MR appears usually in the general form
˙ M (x)xK + C x˙ + Kx = f(x; x);
2 MR ( ˙ ) = 0 + 0 ˙ + 0 ˙
+ cE ˙3 + k1 (23 − 2 − 4 ) = MD (3 − 3 ); 2 IE (4 − 4 )K4 + 12 IE (4 − 4 )˙4 + c0 (˙4 − ˙ )
+ c1 (˙4 − ˙3 ) + cE ˙4 + k0 (4 − ) + k1 (4 − 3 ) = MD (4 − 4 ); IL K + c0 ( ˙ − ˙4 ) + cL ˙ + k0 ( − 4 ) = −MR :
(7)
(9)
P. Metallidis, S. Natsiavas / International Journal of Non-Linear Mechanics 38 (2003) 723 – 738
of motion in the longitudinal direction is expressed in the form
80
muK + 12 AcD u˙ 2 = Fw − Rx ;
peff =10 bar
60
p (bar)
727
where m is the mass, A is the frontal area and u is the longitudinal displacement of the vehicle, is the air density, cD is the aerodynamic drag coeMcient and Rx is the total rolling resistance exerted from the ground to the wheels. Typically, this force is expressed in terms of the total vehicle weight through the relation
peff =5 bar
40
(11)
Rx = !r mg;
20
(12)
where !r is the so-called rolling resistance coeMcient. Finally, if the vehicle wheels do not slip on the ground, the following kinematic conditions are true
0 0
180
360 (deg)
540
720
u(t) ˙ = rw !w (t) ⇒ uK = rw w ;
Fig. 3. Pressure diagrams for two values of the mean e/ective pressure.
where !w =
with coeMcients depending on the nature of the engine load. In particular, the constant term may represent road slope or rolling resistance of vehicles, the linear term may originate from viscous e/ects or from resistance by an electric brake or generator, while the quadratic term may represent aerodynamic or hydraulic resistance. For instance, by applying standard procedures of vehicle dynamics [11], the equation of motion for the moving part of the vehicle can be set in the form Im K + c m ˙ = M R − M w ;
(10)
where cm represents an equivalent viscous damping coeMcient, while Im = I t +
Id Iw + 2 2 2 nt nt nf
and
Mw =
rw Fw : nt nf
In particular, the symbols It ; Id and Iw represent the moment of inertia of the transmission, the driveshaft and the wheels with the axle shafts, nt and nf are the numerical ratio of the transmission and the &nal drive, rw is the wheel radius and Fw is the tractive force exerted from the ground to the driving wheels. In addition, if the vehicle follows a straight path, its equation
(13)
1 ˙ nt nf
and
w =
1 K nt nf
is the angular speed and acceleration of the driving wheels, respectively. Therefore, substitution of Eqs. (11) – (13) in (10) yields a form of the resisting moment MR , which is similar to (9), with r w !r 0 = mg; 0 = cm and n t nf 3 1 rw 0 = AcD : 2 nt nf Moreover, the remaining term, [Im + mrw2 =(nt nf )2 ] K , is added to the mass moment of inertia of the load. The form of the driving and resisting torque examined, together with the dependence of the engine moment of inertia on the crankshaft rotation result in equations of motion which are strongly nonlinear. This is profound even in their simplest version, represented by Eq. (3), while the level of diMculty increases further when some of the system components exhibit signi&cant (exibility. In general, these complications preclude the application of analytical techniques and direct their study to development and application of numerical methodologies. However, some useful information may be obtained by studying these equations of motion with approximate analytical methods, as presented in the following section.
728
P. Metallidis, S. Natsiavas / International Journal of Non-Linear Mechanics 38 (2003) 723 – 738
4. Linear system — perturbation analysis
and
Before the application of appropriate numerical methodologies in solving the nonlinear equations of motion of the systems examined, some insight into the dynamics can be gained by studying more tractable forms of them. More speci&cally, this section presents results for the linearized version of the two degree of freedom system represented by the equations of motion (4) and (5). This is done by &rst introducing the inertia quantity I0 = Ic + m1 r 2 + 12 (m2 + mp )r 2 ; together with the dimensionless parameters (m2 + mp )r 2 &= ; = &!; 2I0 JL IL JL = ; 0 = : I0 JL + 1
f21 = ,0 (’˙ 1 − ’˙ 2 ) − − (,2 + )! − !2 : In the above expressions, the selection of the damping parameters involved is based on the following formulas cn -n = √ ; &,0 = 2 0 -0 ; &2 ,1 = 2 0 -1 ; 2 k 0 I0 &2 ,2 = 2 0 -2 ; while the normalized parameters of the resisting moment appear in the form =
I2 = &'; I0
Therefore, assuming small amplitude vibration, expressed by the conditions = 0 + &’2
and introducing the dimensionless time * = !0 t, with !0 = k0 (IL + I0 )=IL I0 , the equations of motion (4) and (5) take the linearized form (14)
with 0 = +t = !*, 1 −1 1 0 2 ; K = 0 ; M= −1 1 0 JL −!2 sin 20 f0 = 0 and elements of vector f 1 with form f11 = ’K 1 cos 20 − (,0 + 2! sin 20 )’˙ 1 + ,0 ’˙ 2 − ,1 ! − 2!2 ’1 cos 20 − !2 (3 sin 30 − sin 0 )!=2 + !fp (0 ) sin 0
&2
√ 0 0 ; k 0 I0
=
1 &2 I
0
0 :
02 ‘ Fp (0 ) &k0 ∞ (fsn sin n0 + fcn cos n0 ) ; = . f0 +
fp (0 ) =
IE () = I0 [1 − & cos 2 + O(&2 )]:
M ’K + K’ = f 0 (0 ) + &f 1 (0 ; ’; ’; ˙ ’) K
=
Moreover, the de&nition of the dimensionless periodic forcing function is as follows
Then, if & accepts small values, substitution in the expression for the equivalent moment of inertia IE () implies that
= 0 + &’1 ;
02 0 ; & 2 k0
n=1
where the parameter . represents a normalized mean e/ective pressure inside the engine cylinders. By employing standard modal analysis techniques, Eq. (14) is brought into the alternative form, uK + /u = h0 + &h1 ;
(15)
through the application of the coordinate transformation ’(*) = 1u(*); where 1 is the modal matrix of the linearized system, / is a diagonal matrix with diagonal elements equal to 0 and 1, while the elements of vectors h0 and h1 appear in the form hm0 = m0 sin 20 and hm1 = m1 u˙ 1 + m2 u˙ 2 + (m1 uK 1 + m2 uK 2 ) cos 20 + (m1 u˙ 1 + m2 u˙ 2 ) sin 20 + (,m1 u1 + ,m2 u2 ) ×cos 20 +
∞ (&mn sin n0 +-mn cos n0 )+2m : n=1
P. Metallidis, S. Natsiavas / International Journal of Non-Linear Mechanics 38 (2003) 723 – 738
729
Next, application of the multiple time scales method [12] leads to approximate analytical solutions of Eq. (15), with form
ined. For instance, for a bounded response it turns out (from the &rst of them) that the following condition must be satis&ed
u(*; &) = u 0 (*0 ; *1 ) + &u 1 (*0 ; *1 ) + O(&2 );
21 = 0:
(16)
where *0 =* and *1 =&*. More speci&cally, substituting (16) in (15) and collecting terms of order & and &2 leads to the following linear sets of equations: D02 u 0 + /u 0 = h0 (*); D02 u 1 + /u 1 = −2D0 D1 u 0 + h1 (*; u 0 ; D0 u 0 ; D02 u 0 ); (17) where Dn = @=@*n . Then, the zero order solution is determined in the form 10 u10 = c2 − sin 2!* and 4!2 20 sin 2!* + a2 (*1 ) cos[* + 2 (*1 )]; u20 = 1 − 4!2 (18) where c2 is an arbitrary constant. Substitution of u10 and u20 from (18) into (17) yields two di/erential equations for u11 and u21 . These equations have the following form: D02 u11
= 1 A2 e
i(2!+1)*
+ 1 AP 2 e
i(2!−1)*
+ ! 1 A2 e
i*
1 + '1 e4i!* + .1 e2i!t + 21 2
and D02 u21 +u21 = −2iA2 ei* + 2 A2 ei(2!+1)* + 2 AP 2 ei(2!−1)* + .2 e
which is a manifestation of an energy balance. In addition, depending on the value of the frequency !, several resonances are possible to occur at &rst order. Namely, the system exhibits such resonances at values of ! near 0, 14 ; 12 , and 1. However, the &rst and third of these resonances cannot be analyzed e/ectively for the parameter scaling chosen, as it is obvious from (18). Moreover, among the remaining two resonances, the last one presents more interest. In particular, this corresponds to the case where a fundamental parametric resonance and a primary external resonance are activated simultaneously, when the frequency ! is close to the undamped linear natural frequency of the dynamical system examined, which means that ! = 1 + &6:
(20)
In this case, the resulting solvability condition resulting from the di/erential equation for u21 can eventually be expressed by the relation
Therefore, if the solution amplitude is decomposed in the polar form
n=1
+ ! 2 A2 e + ' 2 e
(19)
+ 12 (-21 − i&21 )ei6*1 = 0:
∞
4i!*
+ ( + ,1 + ,2 )! + !2 = 12 !.fs1 ;
i 22 A2 − 2iA2 − 12 (22 + 22 − ,22 )AP 2 e2i6*1
1 + (-1n − i&1n )ein!* + cc 2
i*
Expressed in terms of the original system parameters this condition reappears in the form
2i!t
∞
1 1 + 22 + (-2n − i&2n )ein!* + cc; 2 2 n=1
where the new constant quantities n ; n ; !n ; 'n ; .n (n = 1; 2) are complicated but known functions of the system parameters, while cc denotes complex conjugate of the preceding terms. The last two di/erential equations illustrate some important characteristic features of the system exam-
A2 (*1 ) = 12 a2 (*1 )ei2 (*1 ) ; separation of the real from the imaginary component in the last solvability condition leads eventually to the following slow-(ow system for the amplitude and phase of vibration: a2 = q1 a2 + q2 a2 sin 22 + q3 sin 2 − q4 cos 2 ; a2 2 = 6a2 + q2 a2 cos 22 + q4 sin 2 + q3 cos 2
(21) (22)
with 2 =6*1 −2 . Constant solutions (a20 ; 20 ) of these equations are then obtained by imposing the conditions
730
P. Metallidis, S. Natsiavas / International Journal of Non-Linear Mechanics 38 (2003) 723 – 738
a20 = 0 = 20 . By employing relations (18), (20) and the de&nition of phase 2 it can easily be shown that these solutions correspond to 10 u1 = c2 − sin 2!* + O(&) and 4!2 20 sin 2!* + a20 cos (!* − 20 ) + O(&); u20 = 1 − 4!2 which represent periodic motions of the mechanical system examined, with fundamental frequency equal to !, whose value is determined from condition (19). Determination of these motions starts by &rst solving the resulting algebraic equations for sin 20 and cos 20 . Then, employing the trigonometric identity sin2 20 + cos2 20 = 1 and performing lengthy algebraic manipulations leads &rst to a new relation with form p1 sin 220 + p2 cos 220 = p3 + 1=a220 :
(23)
Next, solving the same equations for sin 220 and cos 220 and following a similar path leads to another relation with form p4 sin 220 + p5 cos 220 = p6 + p7 a220 − 1=4a220 : (24) Eqs. (23) and (24) constitute a linear set of algebraic equations for sin 220 and cos 220 , whose solution can be expressed in the form sin 220 = p8 + p9 a220 + p0 =a220
and
cos 220 = p10 + p11 a220 + p12 =a220 ;
(25)
where all the constant coeMcients pn are known functions of the system parameters. Elimination of the trigonometric terms leads to a fourth degree polynomial for a220 , whose solution provides the amplitude a20 . Then, back substitution in (25) determines the corresponding phase 20 . Finally, the stability properties of these solutions are determined by applying the classical method of linearizarion in (21) and (22) [12,13]. 5. Numerical results In this section, numerical results are presented in an e/ort to illustrate some characteristic features
related to the dynamic behavior of reciprocating engines mounted on road vehicles. First, a two degree of freedom system is examined, with emphasis put in comparing the results of the linearized and the fully nonlinear version of the equations of motion. Then, dynamics of a four cylinder engine nonlinear model is investigated, with attention focused on analysing the e/ect of the system parameters, including engine mis&re. 5.1. Results for two degree of freedom models First, typical numerical results are presented in Fig. 4 for two degree of freedom systems, like the one examined in the previous section. Fig. 4a presents characteristic frequency-response diagrams. In particular, the results represented by the curves denoted with 1 were obtained by solving the slow-(ow Eqs. (21) and (22) for constant solutions in frequency ranges speci&ed by condition (20). The amplitudes of the relative phase Q’ = ’1 − ’2 are plotted versus the frequency detuning parameter 6, for a combination of parameters leading to & = 0:05. For comparison purposes, the curves denoted by 2 represent branches of periodic solutions obtained by applying suitable numerical methodologies [13,14] directly to the nonlinear set of equations of motion (4) and (5). On the other hand, the response diagram denoted by 3 corresponds to branches of periodic solutions obtained for the same set of parameters, by solving Eq. (14) after setting IE () = I0 , which leads to a constant coeMcient linearized system. In all these cases, the broken parts of the curves represent branches of unstable periodic motions. Inspection of Fig. 4a, reveals that the results obtained for the nonlinear set of Eqs. (4), (5) and those obtained by the multiple time scales method appear to be quite di/erent. Besides the di/erences near the 1=3 and 1=4 superharmonic resonances, which are expected to occur, substantial di/erences appear even in the vicinity of the ! ≈ 1 resonance (the predictions of the constant coeMcient model are worse). The explanation for this discrepancy is as follows. In the linear case, the average spin speed ! of the system is not known a priori. Application of the multiple time scales method determined this speed through the condition (19). However, this condition is correct up to &rst order only and is therefore incomplete. For instance, it
P. Metallidis, S. Natsiavas / International Journal of Non-Linear Mechanics 38 (2003) 723 – 738
731
uating the loss of energy due to the torsional damping yields the new expression 1 + ( + ,1 + ,2 )! + !2 + &,0 ! 2 2 4 1 4! a = !.fs1 ; + × 20 2 2 2 (4! − 1) 2 0
Fig. 4. (a) Response diagrams for a linearized and a nonlinear engine model (1: multiple scales with (19), 2: nonlinear model, 3: constant coeMcient model, 4: multiple scales with (26)) at &=0:05 and (b) Average crankshaft speed versus the mean e/ective pressure parameter .ˆ for the same models.
represents an energy balance involving the rigid body component of the motion, but it does not involve the loss of energy due to the torsional damping present in the relative rotation of the system inertia elements. By employing the approximate expressions (18) and eval-
(26)
in place of (19). Then, if the average spin speed ! of the linear system is selected from the last expression, the results obtained are found to be virtually coincident with those obtained by employing the nonlinear equations of motion in the vicinity of the ! ≈ 1 resonance. This is illustrated by the response diagram denoted by 4 in Fig. 4a, which represents branches of periodic solutions obtained for the linearized system by application of the multiple time scales method and after selecting ! from the last expression of energy balance. The above results indicate that caution should be exercised when analyzing the linearized version of the equations of motion. Another important discrepancy between the predictions of the &rst order multiple time scales method developed in the previous section and the results obtained directly from the original nonlinear equations of motion is depicted in a better way in Fig. 4b. This &gure presents the plot of the average spin speed of the same system, as a function of the normalized pressure amplitude parameter .ˆ ≡ .fs1 . Here, the results obtained by the application of the multiple time scales method in conjunction with relation (19) are quite di/erent than those resulting for the nonlinear set of equations of motion, as in the previous &gure. However, the agreement between the results improves again when relation (26) is employed in determining the response frequency !. From Fig. 4b it is apparent that near the resonance examined in the previous section, as well as the other (superharmonic) resonances, the nonlinear system exhibits the well-known Sommerfeld phenomenon, which is expected to occur in the system considered, because it is non-ideal [12,15]. As is typical, the appearance of Sommerfeld phenomena may lead to response with considerable deviations in the mean spin speed, for the same mean e/ective pressure value. This is veri&ed by the response histories shown in Fig. 5. More speci&cally, Fig. 5a shows response histories obtained by direct integration of Eqs. (4) and (5) at
732
P. Metallidis, S. Natsiavas / International Journal of Non-Linear Mechanics 38 (2003) 723 – 738
Fig. 5. Response histories at .ˆ = 0:05, illustrating appearance of Sommerfeld phenomena (a) nonlinear model (b) linearized model with relation (26).
ˆ .=0:05 and for three di/erent sets of initial conditions. Both the amplitude and the fundamental frequency of all these motion are quite di/erent. For comparison purposes, Fig. 5b shows the corresponding response histories determined by obtaining the appropriate initial conditions from the multiple time scales method in conjunction with relation (26). The agreement between the predictions of the two models is satisfactory
Fig. 6. (a) Response diagrams for a linearized and a nonlinear engine model (1: multiple scales with (19), 2: nonlinear model, 3: constant coeMcient model, 4: multiple scales with (26)) at & = 0:24. (b) Average crankshaft speed versus the mean e/ective pressure parameter .ˆ for the same models.
even for the smaller amplitude (and frequency) motion, which occurs away from the ! ≈ 1 resonance. The accuracy of the predictions of the multiple scales method depends strongly on the value of the parameter &. Next, the response diagrams of Fig. 6a illustrate the main e/ects caused by increasing the
P. Metallidis, S. Natsiavas / International Journal of Non-Linear Mechanics 38 (2003) 723 – 738
733
the system dynamics. More speci&cally, the results of Fig. 7a were obtained by applying the multiple scales method for & = 0:05 and demonstrate the e/ect of the damping ratio -0 of the shaft. An increase in the value of this papameter results in a reduction of the nonlinearity e/ects near the ! ≈ 1 range. Likewise, the results shown in Fig. 7b illustrate the e/ect of the normalized quadratic resistance term . Here, an increase in the value of this parameter accentuates the nonlinearity e/ects around ! ≈ 1. 5.2. Results for a four cylinder engine nonlinear model
Fig. 7. (a) E/ect of the cranckshaft daping ratio -0 . (b) E/ect of the aerodynamic resistance coeMcient .
value of &. All these diagrams were obtained for a combination of parameters leading to & = 0:24, keeping the same numbering as in Fig. 4. Clearly, increasing the value of the parameter & magni&es the di/erences between the predictions of the nonlinear and the linear models. Moreover, the system exhibits higher response amplitudes around ! ≈ 1 than close to ! ≈ 1=2. Finally, the diagrams presented in Fig. 7 depict the e/ect of two other important technical parameters on
In this section, numerical results are presented for a four cylinder in-line engine, like the one shown in Fig. 2, with parameter values characterizing a typical road vehicle. In all cases examined, the engine model is nonlinear and its dynamics is governed by Eq. (6) or (7), depending on whether the system exhibits torsional sti/ness. Therefore, the numerical results presented next were obtained by applying suitable numerical methodologies [13,14]. First, Fig. 8a shows results in the form of diagrams, where the independent parameter pe/ represents the mean e/ective pressure developed within the combustion chamber of the engine cylinders [7], while the absicca indicates the resulting mean component of the crankshaft spin speed. These results were obtained by assuming that the engine is rigid and the system is in a steady-state condition. The continuous line was obtained by taking into account the exact dependence of the engine mass moment of inertia IE on the crankshaft angular position, while the broken line corresponds to the case where IE is chosen as constant, that is equal to its mean value I0 . Finally, the dots were obtained by applying the energy balance condition, expressed by (19). The variations detected in Fig. 8a on the average spin speed of the crankshaft due to the dependence of the engine moment of inertia on the crankshaft rotation are negligible. This can be justi&ed by the fact that the equivalent moment of inertia of the engine load (including the vehicle and the (ywheel) is much greater than the moment of inertia of the engine. For instance, Fig. 8b shows the history of the angular acceleration obtained for a typical engine load and for no engine load at all (thicker line). This same e/ect
734
P. Metallidis, S. Natsiavas / International Journal of Non-Linear Mechanics 38 (2003) 723 – 738
Fig. 8. (a) Average engine speed versus mean e/ective pressure. (b) Engine angular acceleration during loaded (thin line) and idle (thick line) function. (c) Angular velocity of loaded engine with variable (continuous line) and constant (broken line) moment of inertia. (d) Angular velocity of idle engine with variable (continuous line) and constant (broken line) moment of inertia.
is also apparent, though less pronounced, in the angular velocity. More speci&cally, the continuous lines in Figs. 8c and 8d represent the histories of the angular velocity obtained for the same engine load and for no load, respectively. In addition, the broken lines in Figs. 8c and 8d depict the angular velocity histories when IE is assumed to remain constant. During these calculations, it was also found that these di/erences become more pronounced for larger values of pe/ . As it is apparent from Fig. 8a, there exists a minimum e/ective pressure level required in order for the engine to start rotating. Similar trends are observed next in Fig. 9, which illustrates e/ects of the vehicle resistance and the crankshaft (exibility parameters. More speci&cally, Fig. 9a shows similar diagrams for the same rigid engine and three levels of the external damping parameter c ≡ c1 = c2 . Clearly, increasing the damping level causes a reduction in the mean
spin speed for the same level of gas pressure. Likewise, Fig. 9b depicts the e/ect of the rolling resistance parameter, while Fig. 9c depicts the e/ect of the aerodynamic resistance term. Finally, Fig. 9d shows the e/ect of the crankshaft torsional sti/ness parameter. In the last case, the di/erences in the results presented on the (pe/ ; !) plane do not seem to be signi&cant. However, it was found that a decrease in the value of the torsional sti/ness parameter or in the damping level encourages the more frequent appearance of Sommerfeld phenomena, in several frequency ranges. For the case with k0 =0:3 MN m=rad, one such area is shown magni&ed in the inset of Fig. 9d. Finally, Fig. 10 illustrates the e/ect of the internal damping parameter -0 , for a system with crankshaft torsional sti/ness parameter k0 =0:6 MN m=rad. Here, the maximum absolute value of the response parameter Q = 4 − F (representing the di/erence between
P. Metallidis, S. Natsiavas / International Journal of Non-Linear Mechanics 38 (2003) 723 – 738
735
Fig. 9. Dependence of the average engine speed on: (a) the external damping parameter; (b) the rolling resistance parameter; (c) the aerodynamic resistance term and (d) the crankshaft torsional sti/ness parameter.
the rotation angle of the (ywheel and the adjacent crank) is shown as a function of the e/ective pressure parameter pe/ , for two values of -0 . Clearly, a decrease in the internal damping level results in bigger response amplitudes. As a consequence, the nonlinear e/ects become more pronounced, which is demonstrated with the appearance of unstable solution branches for the smaller damping value (Fig. 10b). 5.3. E5ects from engine mis6re In obtaining the results presented in the previous section, it was assumed that the pressure developed within each engine cylinder takes values, which are in complete agreement with the theoretical predictions. However, due to various manufacturing or operational factors, there may appear considerable deviations between the measured and the theoretical pressure pattern inside the cylinders. With this in mind, the results
of the present section shed some light in the e/ect of unavoidable engine mis&re on the engine dynamics. Initially, Fig. 11a presents the spin speed history of the rigid engine examined in Section 5.2, when there is a complete mis&re in a single cylinder. These results were obtained for pe/ = 8 bar, which in the case of the healthy engine led to the history shown in Fig. 8c. Examination of these signals reveals &rst that there is a drop in the average spin speed, in all faulty cases. However, besides a phase shift, the recorded signals appear to be identical. Likewise, Fig. 11b presents results for the same system, when there is a partial mis&re in both the &rst and the third cylinder or in all four cylinders of the engine (the cylinder numbering is indicated in Fig. 2). More speci&cally, a 50% mis&re implies that the mean e/ective pressure developed within the combustion chamber of the faulty cylinder reaches only half of the expected value. Again, there are some changes in the signal pro&les, but the average
736
P. Metallidis, S. Natsiavas / International Journal of Non-Linear Mechanics 38 (2003) 723 – 738
Fig. 10. E/ect of the torsional damping parameter on the engine response: (a) -0 = 0:01 and (b) -0 = 0:001.
spin speed remains the same as before. In fact, the average engine spin speed is found to be a function of the total mis&re but it does not seem to depend on its exact distribution in the cylinders. This is expected to cause diMculties in the e/orts to locate faulty cylinders and determine the level of their mis&re, based on information obtained from experimental measurement of the engine spin speed (e.g., [16 –18]). A more noticeable di/erence is observed in the response, when the torsional sti/ness of the crankshaft
Fig. 11. Angular velocity of a rigid engine with: (a) a complete mis&re in a single cylinder and (b) a partial mis&re in the &rst and the third or all the cylinders.
is &nite. For instance, Fig. 12a presents results for all the cylinders and the (ywheel of the same system, but with k0 = 0:1 MN m=rad; -0 = 0:005, when all cylinders are healthy. Here, the recorded spin speed history appears to be di/erent from cylinder to cylinder, due to the torsional (exibility. Moreover, these signals are also di/erent than the corresponding signal of the rigid engine, shown in Fig. 8c. Likewise, Fig. 12b presents results for the same engine, when there is a
P. Metallidis, S. Natsiavas / International Journal of Non-Linear Mechanics 38 (2003) 723 – 738
737
6. Summary and extensions
Fig. 12. Angular velocity of a four cylinder (exible engine with: (a) all cylinders healthy and (b) a 50% mis&re in the &rst and the third cylinder.
50% mis&re in both the &rst and the third cylinder. Again, a comparison with the corresponding results of the healthy system, shown in Fig. 12a, indicates a drop in the average spin speed. From the same results it was also observed that the (ywheel spin speed of the rigid engine is slightly larger than the speed of the (exible engine, where in addition, it also depends on the fault distribution.
Dynamics of simpli&ed models of reciprocating machines was investigated. First, sets of equations of motion were set up for single- and multi-cylinder engines, taking into account the dependence of the engine moment of inertia on the crankshaft rotation. The diMculty level was increased even further, by assuming that both the driving and the resisting moments depend also on the crankshaft rotation. These led to more accurate models but at the same time they introduced strong nonlinearities in the equations of motion. As a result, an approximate analysis was &rst presented for predicting the dynamic behavior of a linearized engine model, which provided some useful insight into the system dynamics. In addition, the equations of motion of the nonlinear models were solved by applying appropriate numerical methodologies. The numerical results were obtained for typical vehicle engines in steady-state conditions. Emphasis was put on investigating the e/ect of the torsional sti/ness and damping parameters of the crankshaft as well as the characteristics of the resisting moment. Among other things, the results indicated the appearance of Sommerfeld phenomena, related to the non-ideal nature of the dynamical systems examined. Also, the dependence of the engine moment of inertia on the crankshaft rotation was found to have a weak e/ect on the average spin speed but it was shown to a/ect more signi&cantly the speed (uctuations. In the &nal part, some issues related to engine mis&re were examined, which was related to a drop in the mean e/ective pressure in the engine cylinders. The models employed in the present study can be improved in several respects. For instance, the pressure developed in the combustion chamber of the engine cylinders can be predicted in a more accurate manner by considering more involved thermodynamic models of the combustion process as well as their coupling with the engine dynamics. In a similar way, the torque resistance from a cooperating electric generator can be modelled by including details about its electric circuit and its coupling with the engine. On the other hand, suitable numerical methodologies are required before the level of mis&re in each cylinder can be predicted in a systematic and computationally eMcient way, handling both measurement and model errors [19,20]. Finally, the mechanical model of the
738
P. Metallidis, S. Natsiavas / International Journal of Non-Linear Mechanics 38 (2003) 723 – 738
engine can also be improved, by taking into account e/ects arising from the (exibility in the bearings, the connecting rod and the supports of the engine frame. Acknowledgements This research was partially supported by a grant from the General Secretariat of Research and Technology, Greek Ministry of Development, through the PENED99 program (reference number: 99EQ 580). References [1] H.H. Denman, Exact solution for the rigid slider-crank mechanism with gas pressure, Mech. Mach. Theory 26 (1991) 415–420. [2] S.-R. Hsieh, S.W. Shaw, The dynamic stability and nonlinear resonance of a (exible connecting rod: continuous parameter model, Nonlinear Dyn. 4 (1993) 573–603. [3] Y. Shiao, C.-H. Pan, J.J. Moskwa, Advanced dynamic spark ignition engine modelling for diagnostics and control, Internat. J. Veh. Design 15 (1994) 578–596. [4] H. Okamura, A. Shinno, T. Yamanaka, A. Suzuki, K. Sogabe, Simple modeling and analysis for crankshaft three-dimensional vibrations, Part 1: background and application to free vibrations, J. Vib. Acoust. 117 (1995) 70–79. [5] H. Ashra&uon, A.M. Whitman, Asymptotic analysis of the torsional vibration in reciprocating machinery, J. Vib. Acoust. 118 (1996) 485–490. [6] E. Brusa, C. Delprete, G. Genta, Torsional vibration of crankshafts: e/ects of non-constant moments of inertia, J. Sound Vib. 205 (1997) 135–150. [7] J.B. Heywood, Internal Combustion Engine Fundamentals, McGraw-Hill, New York, 1988.
[8] L. Dinca, T. Aldemir, G. Rizzoni, Fault detection and identi&cation in dynamic systems with noisy data and parameter=modelling uncertainties, Reliab. Eng. System Saf. 65 (1999) 17–28. [9] S.E. Lyshevski, Energy conversion and optimal energy management in diesel-electric drivetrains of hybrid-electric vehicles, Energy Conversion Manage. 41 (2000) 13–24. [10] J.E. Shigley, Dynamic Analysis of Machines, McGraw-Hill, New York, 1961. [11] T.D. Gillespie, Fundamentals of Vehicle Dynamics, Society of Automotive Engineers, Warrendale, PA, 1992. [12] A.H. Nayfeh, D.T. Mook, Nonlinear Oscillations, Wiley-Interscience, New York, 1979. [13] A.H. Nayfeh, B. Balachandran, Applied Nonlinear Dynamics, Wiley-Interscience, New York, 1995. [14] E. Doedel, AUTO: Software for Continuation and Bifurcation Problems in Ordinary Di/erential Equations, California Institute of Technology, Pasadena, CA, 1986. [15] J. Wauer, P. Buhrle, Dynamics of a (exible slider-crank mechanism driven by a non-ideal source of energy, Nonlinear Dyn. 13 (1997) 221–242. [16] S.J. Citron, J.E. O’Higgins, L.Y. Chen, Cylinder by cylinder engine pressure and pressure torque waveform determination utilizing speed (uctuations. SAE Paper No. 890486, 1989, pp. 933–947. [17] D. Taraza, N.A. Henein, W. Bryzik, Diesel engine diagnosis based on analysis of the crankshaft’s speed variations, SAE Paper No. 982540, 1998, pp. 2183–2195. [18] F.T. Connolly, G. Rizzoni, Real time estimation of engine torque for the detection of engine mis&res, ASME J. Dyn. Systems Meas. Control 116 (1994) 675–686. [19] J.B. Roberts, P.D. Spanos, Random Vibration and Statistical Linearization, Wiley, New York, 1990. [20] L.S. Katafygiotis, C. Papadimitriou, H.F. Lam, 1998. A Probabilistic Approach to Structural Model Updating, Internat. J. Soil Dyn. Earthquake Eng. 17 (1998) 495–507.