- Email: [email protected]
The lifting rotor inflow mode shapes and blade flapping vibration system eigen-analysis
Computer methods in applied mechanics and engineering Comput. Methods Appl. Mech. Engrg. 134 (1996) 91-105
ELSEVIER
The lifting rotor inflow mode shapes and blade flapping vibration system eigen-analysis Yi-Ren Wanga.*, David A. Petersb aDepartment of Aeronautical Engineering, Tamkang University. Tamsui, 25137, Taiwan hDepartment of Mechanical Engineering, Washington University, Sf. Louis, MO 63130, USA
Received 26 January 1995
Abstract
The three-dimensional helicopter rotor inflow vibrating modes and rigid blade flapping motion coupled with unsteady wake system have been studied through eigen-analysis. The Peters’ generalized dynamic inflow theory has been chosen as the unsteady wake model. The three-dimensional inflow velocity vibrating mode shapes for node lines have been plotted in various skew angles. The continuous motions of these mode shapes provide important information and insight into the physical phenomenon of a helicopter rotor-wake aeroelastic system. The damping of axial flight flapping modes and the flapping divergent boundary in forward flight have been compared with Su’s and Lowis’ rsults, respectively, to show the accuracy of our model. The results of this work give us stability information about a rotor aeroelastic system and also reveal that the induced flow has a profound effect on blade dynamics even at high advance ratios, The effect may cause an earlier unstable condition for a helicopter rotor system.
0. Nomenclature Transformation matrices Blade section chord length Lift per unit span of the qth blade, divided by pR2R’ Matrix of induced flow influence coefficients in non-rotating Normalized Legendre functions, (- l)“P~(u) /p’: Rotor blades number Rotor radius Induced velocity at rotor Wake skew angle function, X = tanlx/21 Induced flow coefficients Rotor solidity, QclrrR Coefficients of pressure expansion Wake skew angle Rotor rotational speed, rad./s (m)(m - 2)(m - 4). . . (2 or 1).
alat.
* Corresponding
author.
0045-7825/96/$15.00 0 1996 Elsevier Science S.A. All rights reserved PII SO0457825(96)01027-4
system
92
Y-R. Wang, D.A.
Peters I Comput. Methods Appl. Mech. Engrg.
134 (1996) 91-105
1. Introduction The analysis of a helicopter rotor is more complicated than that of a fixed-wing aircraft. Each blade of the helicopter rotor is affected by the vortex wake shed from itself and from other blades. In addition, the rotor blades cause the blade dynamics to couple strongly with the unsteady aerodynamics. The analysis of helicopter rotor dynamics involves the interactions among induced flow theory (or dynamic inflow theory), lifting theory (or airfoil aerodynamics), and rotor blade dynamics. To investigate the effects of wake dynamics on rotor blade flapping system, one of the best ways is to ‘see’ how they change during different flight conditions. Due to the complex environments of the wake-blade coupled system, even the experimental methods have difficulty in visualizing the continuous interaction of the wake and blade vibrating modes. In recent years, many studies have been done simulating the complex helicopter rotor wake configurations [l, 21. However, most of them are two-dimensional theories or are only applicable to hover. As for aeroelastic analysis, some lower order inflow approximations have been assumed in the three-dimensional analysis [3,4], and of course they are not suitable for detailed investigations. Recently, the precision of Peters’ generalized dynamic inflow theory has been proven and it is simple to use in many applications [5-81. This theory will provide us with a better choice in modeling the unsteady wake system. The objective of this work is to study the interaction of wake dynamics and blade flapping dynamics through the mode shapes and eigenvalues of coupled modes’ continuous vibrations. Our unsteady wake model is based on Peters’ inflow theory. The vibrating inflow velocity mode shapes of the unsteady wake system are studied first to see the characteristics of helicopter rotor unsteady aerodynamics in hover and forward flight. The rigid blade flapping and wake coupled aeroelastic system is also studied by eigen-analysis to explore the physical mechanisms behind the interaction of blade flapping and inflow vibrations.
2. Peters’ generalized dynamic inflow theory The theoretical background of this work is based on Peters’ generalized dynamic inflow theory [5,6]. The unsteady wake theory is developed for lifting rotors based on an acceleration potential for an actuator disk, but the pressure at the rotor disk is discretized at each rotor blade to allow for the effect of a finite number of blades. Based on Prandtl’s solution of the Laplacian equation in ellipsoidal coordinates for an actuator disk flow, Kineer applied this theory to calculate the pressure distribution on a circular wing [9]. The general form of the disturbance pressure can be written as [lo]
c
P~(u)Q~(i~)[C~(t)
cos (m$) + [D:(t)
sin (mg)]
(1)
where P,“(u) and Qr(iq) are associated Legendre functions of the first and second kind in ellipsoidal coordinates. Following Prandtl’s solution form, Peters then expressed the induced flow on the disk in terms of a Fourier series azimuthally and a polynomial distribution radially as U(T, $, f) = 2 pi(l)[,i(t) r,j
cos(r+) + p;(f) sin(@)]
(2)
where (Y; and pr are the induced flow states of the model, and the polynomial ?J’i can be chosen as Legendre functions. The matrix form of the finite state inflow equation can be expressed as [6]
(3)
Y-R. Wang, D.A.
Peters I Comput. Methods Appl. Mech. Engrg.
134 (1996) 91-105
93
(4) where
(5) HY=
(n + m - l)!!(n - m - l)!!
(n + m)!!(n -m)!!
(6)
The matrices [A] and [B] provide for coupling between radial inflow distributions of any given harmonic. They are simple integrals of the Legendre functions (PI) and of the assumed functions (PL). The [i] is the matrix of induced flow influence coefficients. The closed form of these matrices can be found in [6]. Those inflow forcing functions on the right-hand side of Eqs. (3) and (4) are defined as (7)
(8) for the cosine part, and
(9) for the sine part, and L, is the qth blade lift per unit span.
3. Inflow mode shape analysis The eigenvalues of the unsteady wake model are studied first to see the characteristics of unsteady aerodynamics in forward flight. Consider the homogeneous equation for the inflow wake model in the non-rotating system (Eqs. (3) and (4)). The cosine and sine part of the homogeneous equations can be reduced to
VIAic]-l{cl}
[K:](a)*
+
[K;](p)*
+ v[Ai”]-l(p)
=0
(10)
=0
(11)
respectively. The mass matrix [Kr] is chosen as in Eqs. (5) and (6), and the closed forms of the [AL”“] matrices are shown in [8]. We assume that the eigensolution to the homogeneous equations is of the form erv’ . Then, by solving GP>
94
Y-R. Wang, D.A.
IlO
+ !ml
Peters I Cornput. Methods Appl. Mech. Engrg.
134 (1996) 91-105
(12)
= 0
where the dynamic matrix [D] is given by [D] = [K,“]-l[Aic,s]-l
(13)
the normalized eigenvalues (5) can be found from Eq. (12). In order to normalize the dynamic matrix, we use the change of variables {a} = [K;]P”2{Xa}
= [K;]P”*{4x”}
e@
(14)
{p} = [K;J-“2{X”}
= [K;]-“*{+xfi}
e
(15)
For the cosine part, by substitution
of Eq. (14) into Eq. (lo), one obtains
{C}* + V[K’:]-“Z[A~c]-‘[K~]-1’2{X~} Then we can solve the normalized
=0
(16)
eigenvalues from Eq. (14) where the [D] matrix becomes
[D] = [K::]-“2[Ai”“]-‘[lc~]-“2 It is noted that the corresponding
(17) eigenvectors
for each normalized
eigenvalue now become
{$“} = [K;]-“*{$Y}
(18)
For an S-harmonic case, the real parts of eigenvalues (damping of this unsteady aerodynamic system) versus various skew angles from axial to edgewise flow are shown in Fig. 1. The imaginary parts (frequency of this unsteady aerodynamic system) are shown in Fig. 2. From Fig. 1, we can see that those modes become strongly coupled at higher skew angles. In order to identify each mode in axial
-22 -24 -26 -28 -30 -32 -34 0.00
0.10
0.20
0.30
0.40
0.50
0.60
0.70
0.80
0.90
Skew Angle Function X (tan(ti2)) Fig. 1. Inflow eigenvalue plots vs. skew angle function (cm, real part).
1.oo
Y-R. Wang, D.A. Peters I Comput. Methods Appl. Mech. Engrg. 134 (1996) 91-105
95
13
12 11 10 9 $8 8
7
$6 z .-EL5 w 4 3 2 1 0 0.10
0.00
0.20
0.30
0.40
0.50
0.60
0.70
0.00
0.90
1 .oo
Skew Angle Function X (tan(xl2)) Fig. 2. Inflow eigenvalue plots vs. skew angle function (cos, im. part).
flow, we further study the eigenvectors at X = 0. In such a case, once we get the eigenvalues, the corresponding eigenvectors can be found by using Eq. (18); and, thus the mode shapes of the inflow velocity for the cosine part can be obtained by using
tq(r, r/J)=
cr,i ?q(,)c#q
cos(rl))
where q(r)
=
qup
Several induced flow mode shapes contour plots for node lines are shown in Figs. 3 and 4 at some typical skew angles. It should be noted that it is not necessary for the disk edge to be a node line. The view of these contour plots is from above the rotor looking down on the disk. The free stream is from left to right for forward flight plots. Fig. 3 shows the plot of the 1st and 2nd coupled mode, and Fig. 4 is for the 24th single mode. From their eigenvalue plots (Fig. l), we can see the 1st and 2nd mode coupled at ca. a skew angle function X = 0.38 (X = tan(X/2)). Through the study of their mode shapes, we observe how the 2nd mode dominates the coupled mode. A node line occurs as X > 0.40 (Figs. 3E and F) and makes this coupled mode vibrate forward and aft. For the sine part, again, we follow the same procedure as in the cosine part. The real and imaginary parts of the eigenvalues for 8 harmonics are shwon in Figs. 5 and 6, respectively. As before, we obtain the induced flow mode shapes of the sine part by using
96
Y-R. Wang, D.A.
Peters I Comput. Methods Appl. Mech.
Engrg.
134
(1996)
91-105
Fig. 3. (A) The 1st inflow mode shape at X = 0.0 (cos part); (B) the 2nd inflow mode shape at X = 0.0 (COS part); (C) the 1st inflow mode shape at X = 0.25 (cos part); (D) the 2nd inflow mode shape at X = 0.25 (cos part); (E) the 1st and 2nd coupled inflow mode shape at X = 0.40 (cos part); (F) the 1st and 2nd coupled inflow mode shape at X = 1.00 (cos part).
Figs. 7-10 show the plots of the 17th, 18th, 19th, 20th mode, respectively, at several typical skew angle functions. Through the node line contour plots, we can see the characteristics of each inflow velocity mode and how they affect each other.
4. Coupling of wake dynamics and blade flapping dynamics The first step in deriving a coupled blade flapping-inflow equation is to formulate the lift equation (L,) in the right-hand side of the inflow equations (Eqs. (3) and (4)). We assume a small angle of attack, no stall and no reversed flow on the rotor disk. Since we desire an induced flow normal to the airfoil surface, the most logical assumption is to use lift normal to airfoil (circulatory only) as the wake forcing function [ 111,
Y-R. Wang, D.A.
Peters I Comput. Methods Appl. Mech. Engrg.
134 (1996) 91-105
Fig. 4. (A) The 24th inflow mode shape at X = 0.0 (cos part); (B) X = 0.10 (cos part); (C) X = 0.50 (cos part); (D) X = 1.00 (cos part).
-8
8 -18 1 m -18 2 al -20 iE -22 -24 -28 -28 -30 -32 -34 0.00
0.10
0.20
0.80 0.70 0.60 0.50 0.40 0.30 Skew Angle Function X (tanfxl2))
0.90
Fig. 5. Inflow eigenvalue plots vs. skew angle function (sin, real part)
1.00
98
Y-R. Wang, D.A. Peters I Comput. Methods Appl. Mech. Engrg. 134 (1996) 91-105
12 11 10 9
0.00
0.10
0.20
0.30
0.40
0.50
0.60
0.70
0.60
0.90
1 .oo
Skew Angle Function X (tan(x/2)) Fig. 6. Inflow eigenvalue
plots vs. skew angle
function
(f-3
0’)
@
@
Fig. 7. The 17th inflow mode shape part).
at (A) X = 0.0 (sin part);
(B) X = 0.05 (sin part);
(sin, im. part).
(C) X = 0.15 (sin part);
(D) X = 0.25 (sin
Y-R. Wang, D.A.
Peters
I Comput. Methods Appl. Mech. Engrg. 134 (1996) 91-105
0))
]
Fig. 8. The 18th inflow mode shape at (A) X = 0.0 (sin part); (B) X = 0.05 (sin part); (C) X = 0.20 (sin part); (D) X = 0.25 (sin part); (E) X=0.35 (sin part); (F) X=0.50 (sin part); (G) X= 0.75 (sin part); (H) X= 1.00 (sin part).
100
Y-R. Wang, D.A. Peters I Comput. Methods Appl. Mech. Engrg. 134 (19%) 91-105
Fig. 9. The 19th inflow mode shape at (A) X = 0.0 (sin part); (B) X = 0.20 (sin part); (C) X = 0.25 (sin part); (D) X = 0.40 (sin part).
Fig. 10. The 20th inflow mode shape at (A) X= 0.0 (sin part); (B) X = 0.05 (sin part); (C) X = 0.75 (sin part); (D) X = 1.00 (sin part).
Y-R. Wang, D.A.
Peters I Comput. Methods Appl. Mech. Engrg.
134 (1996) 91-105
101
(22)
L,=~acU,(w,+tw,-A”-3A,)
where a is the slope of lift curve, C is the dimensionless blade chord, U, = (? + /J, sin I++~)in forward flight, w,, and w1 is the first two terms of the Glauert series expansion. They represent the relative air speed normal to the airfoil, including both free stream and airfoil motions (positive down). A, and A, are the first two terms of a Glauert series expansion for the wake-induced velocity. We then substitute Eq. (22) into our wake forcing functions (7:’ and 7p) to obtain the wake equations. Now with the rigid blade assumption, we can write down the blade flapping equation for the qth blade by taking equilibrium of inertial, centrifugal, and aerodynamic moments about blade flapping hinge as p,** +p2p4 = O’L?dJ I
(23)
where p4 is the blade flapping angle at the qth blade, and p is the blade flapping frequency. Since we need the lift that contributes to the moment, from [ll], the total lift can be formulated the following form,
as
(24) It is noted that the first half part of this lift formulation represents a circulatory lift, and the second half part represents a non-circulatory lift. The blade flapping equation can be easily obtained by substituting Eq. (24) into the right-hand side integral in Eq. (23).
5. Rotor eigen-analysis
of the wake-blade
coupled system
The governing equations are the wake equations and blade flapping equations that were derived in the last section. We can write down the coupled equations in a state variable form and the matrix formulation is shown below:
X
1
Y-R. Wang, D.A. Peters I Comput. Methods Appl. Mech. Engrg.
102
134 (1996) 91-105
where ‘X ’ represents non-zero terms. The apparent mass matrix is a diagonal matrix and the formulation of the diagonal submatrix [K,“] was shown in Eqs. (5) and (6). It should be noted that the stiffness matrix [G] involves a time-independent term (i) which implies that the stiffness matrix is a periodic coefficient matrix. To solve the eigenvalues of this linear differential equation with periodic coefficients, we set the forcing function ({H}) on the right-hand side of this matrix formulation to be zero and apply the Floquet transition matrix [12] in the solution methodology. The eigenvalues of the Floquet transition matrix are shown to be the exponential occurring in Floquet’s form of the initial value solution. In our work, the Floquet transition matrix is computed from a fourth-order Runge-Kutta method for integration. The relation between the eigenvalues of Floquet transition matrix and the system is that if we assume the eigenvalues of Floquet transition matrix as A, = a + ib, then the eigenvalues of the system are _I b 2rm a+-
(25)
T
A general program has been written incorporating the integration and eigenvalue, eigenvector routines. This program is capable of handling any order system of periodic coefficient equations. The inputs are the harmonic numbers (up to 8 harmonics) that we need in the [AL’“] matrices, and the flight condition. The outputs of this program are the eigenvalues of this system and the eigenvectors of their corresponding eigenvalues of Floquet transition matrix.
6. Results and discussions from coupled system eigenvalues In order to compare with Su’s axial flight results [7], we choose several constant parameters to be the same as in [7], Q = 4, p = 1.15, au = 0.377, y = 5.0, V= 0.05, and p = 0.0. It should be noted that Su neglects the non-circulatory term (w,*) in his blade equation and chooses four shape functions for each harmonic in the [i] matrix. However, our shape functions numbers are based on [6]. Table 1 is the comparison of the real parts of eigenvalues for all four flapping and no inflow modes. Eight harmonics are chosen in this case. In the axial flight case, the non-circulatory effect may lessen the damping slightly, but our results still agree very well with Su’s. In forward flight, there are several interdependent variables (V, x or /.L,and AT). Therefore it is useful to fix one of them and then vary the other two in order to study their effects on blade dymamics. One possibility is to choose V constant. Then, however, the limit on advance ratio will be less than or equal to V. So we cannot obtain a CLlarger than V. Alternatively, we set x constant and vary p. However, at p = 0, V will be zero. Since the term l/V appears in the wake equation, this will diverge as P approaches zero. A third choice is fo a w thrust condition and a fixed climb rate, A, z 0 and A, constant, for each p, we have V= /--loF’ + A$, and tanx =pulh,. Thus, p can vary freely from 0 to infinity without V going to zero and with x varying from 0 to 7rl2, its entire range. Since a helicopter usually flies in the range of V= 0.05-0.10 in axial flight, we choose A, = 0.05, and let /L vary from 0 to 1.5. Figs. 11 and 12 are plots of eigenvalues vs. advance ratio with p = 1.0, y = 5.0 and A, = 0.055. Fig. 11 is for the real part plot. We can see that the ‘no-inflow’-mode becomes unstable near I_L= 1.4 which Table 1 Comparison
of flapping modes, real part eigenvalues (p = 1.15, ua = 0.337)
Now; S = 4, from [7]. With wx S based on [6]
Differential mode
Collective mode
Progressing mode
Regressing mode
No inflow
-0.2930
-0.2810
-0.3000
-0.1140
-0.3125
-0.2901
-0.2757
-0.2983
-0.1104
-0.3094
Y-R. Wang, D.A.
Peters I Comput.
Methods Appl. Mech. Engrg.
134 (1996) 91-105
103
reeressine ---_
d__
advance ratio
1
Fig. 11. Blade
flapping
eigenvalues
(coupled
with wake system)
vs. advance
ratio
(real part)
:: d
L’
I ::
‘\ ‘.
-.
N
I ii
, 1
?_
1,
.
,
P
1 \
/
I ,
”
\
EEX 30’:-- __-+-;_*
\
i.
,.
,A----____----
/
:
__--
*_-,_*-
/
I/ .: 1; . :: ,‘:
i’
8’
,;
7 -2
”
%g
... -.-....
I
Fig. 12. Blade
flapping
r
,./’
.... . ..____...
‘S &,a-
eigenvalues
(coupled
with wake system)
vs. advance
ratio
(im. part).
agrees with the result in [13] and [14]. We also see that the wake will lower the flapping damping and cause earlier instability than the no inflow mode, which is also shown in Fig. 11. Another interesting thing is that, although these flapping mode branches state at /-L= 0 as regressing, progressing, differential and collective, it is not necessary for them to maintain this character as /L increases. The legends in those figures are just marked as their initial modes. One way to see which mode is dominant is to check the eigenvectors. For a specific advance ratio, we multiply the Happing angle eigenvectors of an eigenvalue with a ‘factor matrix’ as follows:
104
Y-R. Wang, D.A.
Peters I Comput. Methods Appl. Mech. Engrg.
134 (1996) 91-105
The left-hand side matrix are complex numbers, and their absolute values (amplitudes of shape modes) can give information on the identification of flapping modes. In Fig. 11, as p = 1.3, we pick up a set of eigenvalues (0.01587, 0.00497), (0.00544, O.O), (-0.00600, O.O), (-0.51459, 0.01239), and (-0.52398, 0.0). The branches begin with collective, progressing, regrssing, differential, and the lower branch of the progressing mode, respectively. The amplitudes of the vectors of the first eigenvalue are: (1.82542, 1.88935, 0.57906, 0.44757) We can see that, although this eigenvalue begins with collective mode, it is a strongly coupled with the differential mode near p = 0.22. Thus, it becomes partly a collective and partly a differential mode. For the same reason, the 4th eigenvalue begins with differential mode, but its amplitudes at p = 0.13 are: (1.29290, 1.08725, 0.13153, 0.07587). Again, this branch has become a coupled differential-collective mode. However, the 3rd eigenvalue has the following amplitudes: (2.81553, 0.58549, 4.96667, 0.39344) which shows that it begins and ends as a regressing mode. The amplitudes of upper and lower branches of progressing mode are: (0.16252, 0.97948, 1.75680, 0.09842) and (0.27714, 0.13044, 0.18527, 1.15635) We can see that they begin and end as progressing modes.
7. Conclusions The mode shapes of an inflow wake model and the eigen-system of coupled blade-wake aeroelastic model have been studied. The results of this work give us the stability information of this aeroelastic system and the understanding of the effect of inflow dynamics on rotor blade flapping dynamics in forward flight. The conclusions from this research are as follows: (1) The visualizing of mode shape vibration provides us with a better understanding of the rotor wake mechanism. (2) The eigen-analysis exposes the complicated coupling mechanism of the blade flapping and inflow modes. Due to this inflow coupling effect, the blade flapping modes may not keep their initial characteristics, as the skew angle or advance ratio changes. (3) The induced flow has a profound effect on blade dynamics even at high advance ratios. The effect may cause an earlier unstable condition for a helicopter rotor system.
Acknowledgments When the authors were in Georgia Tech, this research was supported in part by a joint NASA-Army grant, NAG-2-462, and additional funding was provided by the Georgia Tech Center of Excellence for Rotary Wing Aircraft Technology, sponsored by the Army Research Office. The continuing part of this work was supported by the National Science Council in the Republic of China under contract No. NSC 83-0424-E-032-001. Their support is gratefully appreciated by the authors.
References [l] T. Theodorsen, General theory of aerodynamic instability [2] R.G. Loewy, A two dimensional approach to the unsteady [3] S.T. Crews, K.H. Hohenemser and R.A. Ormiston, An (1973). [4] R.A. Ormiston, Application of simplified inflow models to 34-37.
and the mechanism of flutter, NACA R.496, 1935. aerodynamics of rotary wing, J. Aerospace Sci. 24 (1957) 81-98. unsteady wake model for a hingeless rotor, J. Aircraft lO(12) rotorcraft dynamic analysis, J. Am. Helicopter Sot. 21(3) (1976)
Y-R. Wang, D.A.
Peters I Comput. Methods Appl. Mech. Engrg.
134 (1996) 91-105
105
[S] D.A. Peters and C.J. He, Comparison of measured induced velocities with results from a closed form finite state wake model in forward flight, 45th Annual National Forum of the American Helicopter Society, Boston, MA, May 22-24, 1989. [63 C.J. He, Development and application of a generalized dynamic wake theory for lifting rotors, Ph.D. Thesis, School of Aerospace Engineering, Georgia Institute of Technology, July, 1989. [7] A. Su, Application of a state space wake model to elastic blade flapping in Hover, Ph.D. Thesis, School of Aerospace Engineering, Georgia Institute of Technology, August, 1989. [S] Y.-R. Wang, The effect of wake dynamics on rotor eigenvalues in forward flight, Ph.D. Thesis, School of Aerospace Engineering, Georgia Institute of Technology, June, 1992. [9] W. Kinner, Theory of circular wing, Ingenieur Archiv, 7 (1937) 47-80, Translation No. 2345, Ministry of Aircraft Production, UK. [lo] M. Joglekar and R. Loewy, An actuator disc analysis of helicopter wake geometry and the corresponding blade response, USAAVLABS TR 69-66, December, 1970. [ 1 l] W. Johnson, Helicopter Theory (Princeton University Press, Princeton, 1980) 471-480. [12] W. Kaplan, Operational Methods of Linear Systems (Addison Wesley, 1962) 472-488. [13] O.J. Lowis, The stability of rotor blade flapping motion at high tip speed ratios, Reports and Memoranda No. 3544, January, 1963. [14] W.E. Hall, Application of Floquet theory to the analysis of rotary wing VTOL stability, Ph.D. Dissertation, Stanford University, Februrary, 1970.
ELSEVIER
The lifting rotor inflow mode shapes and blade flapping vibration system eigen-analysis Yi-Ren Wanga.*, David A. Petersb aDepartment of Aeronautical Engineering, Tamkang University. Tamsui, 25137, Taiwan hDepartment of Mechanical Engineering, Washington University, Sf. Louis, MO 63130, USA
Received 26 January 1995
Abstract
The three-dimensional helicopter rotor inflow vibrating modes and rigid blade flapping motion coupled with unsteady wake system have been studied through eigen-analysis. The Peters’ generalized dynamic inflow theory has been chosen as the unsteady wake model. The three-dimensional inflow velocity vibrating mode shapes for node lines have been plotted in various skew angles. The continuous motions of these mode shapes provide important information and insight into the physical phenomenon of a helicopter rotor-wake aeroelastic system. The damping of axial flight flapping modes and the flapping divergent boundary in forward flight have been compared with Su’s and Lowis’ rsults, respectively, to show the accuracy of our model. The results of this work give us stability information about a rotor aeroelastic system and also reveal that the induced flow has a profound effect on blade dynamics even at high advance ratios, The effect may cause an earlier unstable condition for a helicopter rotor system.
0. Nomenclature Transformation matrices Blade section chord length Lift per unit span of the qth blade, divided by pR2R’ Matrix of induced flow influence coefficients in non-rotating Normalized Legendre functions, (- l)“P~(u) /p’: Rotor blades number Rotor radius Induced velocity at rotor Wake skew angle function, X = tanlx/21 Induced flow coefficients Rotor solidity, QclrrR Coefficients of pressure expansion Wake skew angle Rotor rotational speed, rad./s (m)(m - 2)(m - 4). . . (2 or 1).
alat.
* Corresponding
author.
0045-7825/96/$15.00 0 1996 Elsevier Science S.A. All rights reserved PII SO0457825(96)01027-4
system
92
Y-R. Wang, D.A.
Peters I Comput. Methods Appl. Mech. Engrg.
134 (1996) 91-105
1. Introduction The analysis of a helicopter rotor is more complicated than that of a fixed-wing aircraft. Each blade of the helicopter rotor is affected by the vortex wake shed from itself and from other blades. In addition, the rotor blades cause the blade dynamics to couple strongly with the unsteady aerodynamics. The analysis of helicopter rotor dynamics involves the interactions among induced flow theory (or dynamic inflow theory), lifting theory (or airfoil aerodynamics), and rotor blade dynamics. To investigate the effects of wake dynamics on rotor blade flapping system, one of the best ways is to ‘see’ how they change during different flight conditions. Due to the complex environments of the wake-blade coupled system, even the experimental methods have difficulty in visualizing the continuous interaction of the wake and blade vibrating modes. In recent years, many studies have been done simulating the complex helicopter rotor wake configurations [l, 21. However, most of them are two-dimensional theories or are only applicable to hover. As for aeroelastic analysis, some lower order inflow approximations have been assumed in the three-dimensional analysis [3,4], and of course they are not suitable for detailed investigations. Recently, the precision of Peters’ generalized dynamic inflow theory has been proven and it is simple to use in many applications [5-81. This theory will provide us with a better choice in modeling the unsteady wake system. The objective of this work is to study the interaction of wake dynamics and blade flapping dynamics through the mode shapes and eigenvalues of coupled modes’ continuous vibrations. Our unsteady wake model is based on Peters’ inflow theory. The vibrating inflow velocity mode shapes of the unsteady wake system are studied first to see the characteristics of helicopter rotor unsteady aerodynamics in hover and forward flight. The rigid blade flapping and wake coupled aeroelastic system is also studied by eigen-analysis to explore the physical mechanisms behind the interaction of blade flapping and inflow vibrations.
2. Peters’ generalized dynamic inflow theory The theoretical background of this work is based on Peters’ generalized dynamic inflow theory [5,6]. The unsteady wake theory is developed for lifting rotors based on an acceleration potential for an actuator disk, but the pressure at the rotor disk is discretized at each rotor blade to allow for the effect of a finite number of blades. Based on Prandtl’s solution of the Laplacian equation in ellipsoidal coordinates for an actuator disk flow, Kineer applied this theory to calculate the pressure distribution on a circular wing [9]. The general form of the disturbance pressure can be written as [lo]
c
P~(u)Q~(i~)[C~(t)
cos (m$) + [D:(t)
sin (mg)]
(1)
where P,“(u) and Qr(iq) are associated Legendre functions of the first and second kind in ellipsoidal coordinates. Following Prandtl’s solution form, Peters then expressed the induced flow on the disk in terms of a Fourier series azimuthally and a polynomial distribution radially as U(T, $, f) = 2 pi(l)[,i(t) r,j
cos(r+) + p;(f) sin(@)]
(2)
where (Y; and pr are the induced flow states of the model, and the polynomial ?J’i can be chosen as Legendre functions. The matrix form of the finite state inflow equation can be expressed as [6]
(3)
Y-R. Wang, D.A.
Peters I Comput. Methods Appl. Mech. Engrg.
134 (1996) 91-105
93
(4) where
(5) HY=
(n + m - l)!!(n - m - l)!!
(n + m)!!(n -m)!!
(6)
The matrices [A] and [B] provide for coupling between radial inflow distributions of any given harmonic. They are simple integrals of the Legendre functions (PI) and of the assumed functions (PL). The [i] is the matrix of induced flow influence coefficients. The closed form of these matrices can be found in [6]. Those inflow forcing functions on the right-hand side of Eqs. (3) and (4) are defined as (7)
(8) for the cosine part, and
(9) for the sine part, and L, is the qth blade lift per unit span.
3. Inflow mode shape analysis The eigenvalues of the unsteady wake model are studied first to see the characteristics of unsteady aerodynamics in forward flight. Consider the homogeneous equation for the inflow wake model in the non-rotating system (Eqs. (3) and (4)). The cosine and sine part of the homogeneous equations can be reduced to
VIAic]-l{cl}
[K:](a)*
+
[K;](p)*
+ v[Ai”]-l(p)
=0
(10)
=0
(11)
respectively. The mass matrix [Kr] is chosen as in Eqs. (5) and (6), and the closed forms of the [AL”“] matrices are shown in [8]. We assume that the eigensolution to the homogeneous equations is of the form erv’ . Then, by solving GP>
94
Y-R. Wang, D.A.
IlO
+ !ml
Peters I Cornput. Methods Appl. Mech. Engrg.
134 (1996) 91-105
(12)
= 0
where the dynamic matrix [D] is given by [D] = [K,“]-l[Aic,s]-l
(13)
the normalized eigenvalues (5) can be found from Eq. (12). In order to normalize the dynamic matrix, we use the change of variables {a} = [K;]P”2{Xa}
= [K;]P”*{4x”}
e@
(14)
{p} = [K;J-“2{X”}
= [K;]-“*{+xfi}
e
(15)
For the cosine part, by substitution
of Eq. (14) into Eq. (lo), one obtains
{C}* + V[K’:]-“Z[A~c]-‘[K~]-1’2{X~} Then we can solve the normalized
=0
(16)
eigenvalues from Eq. (14) where the [D] matrix becomes
[D] = [K::]-“2[Ai”“]-‘[lc~]-“2 It is noted that the corresponding
(17) eigenvectors
for each normalized
eigenvalue now become
{$“} = [K;]-“*{$Y}
(18)
For an S-harmonic case, the real parts of eigenvalues (damping of this unsteady aerodynamic system) versus various skew angles from axial to edgewise flow are shown in Fig. 1. The imaginary parts (frequency of this unsteady aerodynamic system) are shown in Fig. 2. From Fig. 1, we can see that those modes become strongly coupled at higher skew angles. In order to identify each mode in axial
-22 -24 -26 -28 -30 -32 -34 0.00
0.10
0.20
0.30
0.40
0.50
0.60
0.70
0.80
0.90
Skew Angle Function X (tan(ti2)) Fig. 1. Inflow eigenvalue plots vs. skew angle function (cm, real part).
1.oo
Y-R. Wang, D.A. Peters I Comput. Methods Appl. Mech. Engrg. 134 (1996) 91-105
95
13
12 11 10 9 $8 8
7
$6 z .-EL5 w 4 3 2 1 0 0.10
0.00
0.20
0.30
0.40
0.50
0.60
0.70
0.00
0.90
1 .oo
Skew Angle Function X (tan(xl2)) Fig. 2. Inflow eigenvalue plots vs. skew angle function (cos, im. part).
flow, we further study the eigenvectors at X = 0. In such a case, once we get the eigenvalues, the corresponding eigenvectors can be found by using Eq. (18); and, thus the mode shapes of the inflow velocity for the cosine part can be obtained by using
tq(r, r/J)=
cr,i ?q(,)c#q
cos(rl))
where q(r)
=
qup
Several induced flow mode shapes contour plots for node lines are shown in Figs. 3 and 4 at some typical skew angles. It should be noted that it is not necessary for the disk edge to be a node line. The view of these contour plots is from above the rotor looking down on the disk. The free stream is from left to right for forward flight plots. Fig. 3 shows the plot of the 1st and 2nd coupled mode, and Fig. 4 is for the 24th single mode. From their eigenvalue plots (Fig. l), we can see the 1st and 2nd mode coupled at ca. a skew angle function X = 0.38 (X = tan(X/2)). Through the study of their mode shapes, we observe how the 2nd mode dominates the coupled mode. A node line occurs as X > 0.40 (Figs. 3E and F) and makes this coupled mode vibrate forward and aft. For the sine part, again, we follow the same procedure as in the cosine part. The real and imaginary parts of the eigenvalues for 8 harmonics are shwon in Figs. 5 and 6, respectively. As before, we obtain the induced flow mode shapes of the sine part by using
96
Y-R. Wang, D.A.
Peters I Comput. Methods Appl. Mech.
Engrg.
134
(1996)
91-105
Fig. 3. (A) The 1st inflow mode shape at X = 0.0 (cos part); (B) the 2nd inflow mode shape at X = 0.0 (COS part); (C) the 1st inflow mode shape at X = 0.25 (cos part); (D) the 2nd inflow mode shape at X = 0.25 (cos part); (E) the 1st and 2nd coupled inflow mode shape at X = 0.40 (cos part); (F) the 1st and 2nd coupled inflow mode shape at X = 1.00 (cos part).
Figs. 7-10 show the plots of the 17th, 18th, 19th, 20th mode, respectively, at several typical skew angle functions. Through the node line contour plots, we can see the characteristics of each inflow velocity mode and how they affect each other.
4. Coupling of wake dynamics and blade flapping dynamics The first step in deriving a coupled blade flapping-inflow equation is to formulate the lift equation (L,) in the right-hand side of the inflow equations (Eqs. (3) and (4)). We assume a small angle of attack, no stall and no reversed flow on the rotor disk. Since we desire an induced flow normal to the airfoil surface, the most logical assumption is to use lift normal to airfoil (circulatory only) as the wake forcing function [ 111,
Y-R. Wang, D.A.
Peters I Comput. Methods Appl. Mech. Engrg.
134 (1996) 91-105
Fig. 4. (A) The 24th inflow mode shape at X = 0.0 (cos part); (B) X = 0.10 (cos part); (C) X = 0.50 (cos part); (D) X = 1.00 (cos part).
-8
8 -18 1 m -18 2 al -20 iE -22 -24 -28 -28 -30 -32 -34 0.00
0.10
0.20
0.80 0.70 0.60 0.50 0.40 0.30 Skew Angle Function X (tanfxl2))
0.90
Fig. 5. Inflow eigenvalue plots vs. skew angle function (sin, real part)
1.00
98
Y-R. Wang, D.A. Peters I Comput. Methods Appl. Mech. Engrg. 134 (1996) 91-105
12 11 10 9
0.00
0.10
0.20
0.30
0.40
0.50
0.60
0.70
0.60
0.90
1 .oo
Skew Angle Function X (tan(x/2)) Fig. 6. Inflow eigenvalue
plots vs. skew angle
function
(f-3
0’)
@
@
Fig. 7. The 17th inflow mode shape part).
at (A) X = 0.0 (sin part);
(B) X = 0.05 (sin part);
(sin, im. part).
(C) X = 0.15 (sin part);
(D) X = 0.25 (sin
Y-R. Wang, D.A.
Peters
I Comput. Methods Appl. Mech. Engrg. 134 (1996) 91-105
0))
]
Fig. 8. The 18th inflow mode shape at (A) X = 0.0 (sin part); (B) X = 0.05 (sin part); (C) X = 0.20 (sin part); (D) X = 0.25 (sin part); (E) X=0.35 (sin part); (F) X=0.50 (sin part); (G) X= 0.75 (sin part); (H) X= 1.00 (sin part).
100
Y-R. Wang, D.A. Peters I Comput. Methods Appl. Mech. Engrg. 134 (19%) 91-105
Fig. 9. The 19th inflow mode shape at (A) X = 0.0 (sin part); (B) X = 0.20 (sin part); (C) X = 0.25 (sin part); (D) X = 0.40 (sin part).
Fig. 10. The 20th inflow mode shape at (A) X= 0.0 (sin part); (B) X = 0.05 (sin part); (C) X = 0.75 (sin part); (D) X = 1.00 (sin part).
Y-R. Wang, D.A.
Peters I Comput. Methods Appl. Mech. Engrg.
134 (1996) 91-105
101
(22)
L,=~acU,(w,+tw,-A”-3A,)
where a is the slope of lift curve, C is the dimensionless blade chord, U, = (? + /J, sin I++~)in forward flight, w,, and w1 is the first two terms of the Glauert series expansion. They represent the relative air speed normal to the airfoil, including both free stream and airfoil motions (positive down). A, and A, are the first two terms of a Glauert series expansion for the wake-induced velocity. We then substitute Eq. (22) into our wake forcing functions (7:’ and 7p) to obtain the wake equations. Now with the rigid blade assumption, we can write down the blade flapping equation for the qth blade by taking equilibrium of inertial, centrifugal, and aerodynamic moments about blade flapping hinge as p,** +p2p4 = O’L?dJ I
(23)
where p4 is the blade flapping angle at the qth blade, and p is the blade flapping frequency. Since we need the lift that contributes to the moment, from [ll], the total lift can be formulated the following form,
as
(24) It is noted that the first half part of this lift formulation represents a circulatory lift, and the second half part represents a non-circulatory lift. The blade flapping equation can be easily obtained by substituting Eq. (24) into the right-hand side integral in Eq. (23).
5. Rotor eigen-analysis
of the wake-blade
coupled system
The governing equations are the wake equations and blade flapping equations that were derived in the last section. We can write down the coupled equations in a state variable form and the matrix formulation is shown below:
X
1
Y-R. Wang, D.A. Peters I Comput. Methods Appl. Mech. Engrg.
102
134 (1996) 91-105
where ‘X ’ represents non-zero terms. The apparent mass matrix is a diagonal matrix and the formulation of the diagonal submatrix [K,“] was shown in Eqs. (5) and (6). It should be noted that the stiffness matrix [G] involves a time-independent term (i) which implies that the stiffness matrix is a periodic coefficient matrix. To solve the eigenvalues of this linear differential equation with periodic coefficients, we set the forcing function ({H}) on the right-hand side of this matrix formulation to be zero and apply the Floquet transition matrix [12] in the solution methodology. The eigenvalues of the Floquet transition matrix are shown to be the exponential occurring in Floquet’s form of the initial value solution. In our work, the Floquet transition matrix is computed from a fourth-order Runge-Kutta method for integration. The relation between the eigenvalues of Floquet transition matrix and the system is that if we assume the eigenvalues of Floquet transition matrix as A, = a + ib, then the eigenvalues of the system are _I b 2rm a+-
(25)
T
A general program has been written incorporating the integration and eigenvalue, eigenvector routines. This program is capable of handling any order system of periodic coefficient equations. The inputs are the harmonic numbers (up to 8 harmonics) that we need in the [AL’“] matrices, and the flight condition. The outputs of this program are the eigenvalues of this system and the eigenvectors of their corresponding eigenvalues of Floquet transition matrix.
6. Results and discussions from coupled system eigenvalues In order to compare with Su’s axial flight results [7], we choose several constant parameters to be the same as in [7], Q = 4, p = 1.15, au = 0.377, y = 5.0, V= 0.05, and p = 0.0. It should be noted that Su neglects the non-circulatory term (w,*) in his blade equation and chooses four shape functions for each harmonic in the [i] matrix. However, our shape functions numbers are based on [6]. Table 1 is the comparison of the real parts of eigenvalues for all four flapping and no inflow modes. Eight harmonics are chosen in this case. In the axial flight case, the non-circulatory effect may lessen the damping slightly, but our results still agree very well with Su’s. In forward flight, there are several interdependent variables (V, x or /.L,and AT). Therefore it is useful to fix one of them and then vary the other two in order to study their effects on blade dymamics. One possibility is to choose V constant. Then, however, the limit on advance ratio will be less than or equal to V. So we cannot obtain a CLlarger than V. Alternatively, we set x constant and vary p. However, at p = 0, V will be zero. Since the term l/V appears in the wake equation, this will diverge as P approaches zero. A third choice is fo a w thrust condition and a fixed climb rate, A, z 0 and A, constant, for each p, we have V= /--loF’ + A$, and tanx =pulh,. Thus, p can vary freely from 0 to infinity without V going to zero and with x varying from 0 to 7rl2, its entire range. Since a helicopter usually flies in the range of V= 0.05-0.10 in axial flight, we choose A, = 0.05, and let /L vary from 0 to 1.5. Figs. 11 and 12 are plots of eigenvalues vs. advance ratio with p = 1.0, y = 5.0 and A, = 0.055. Fig. 11 is for the real part plot. We can see that the ‘no-inflow’-mode becomes unstable near I_L= 1.4 which Table 1 Comparison
of flapping modes, real part eigenvalues (p = 1.15, ua = 0.337)
Now; S = 4, from [7]. With wx S based on [6]
Differential mode
Collective mode
Progressing mode
Regressing mode
No inflow
-0.2930
-0.2810
-0.3000
-0.1140
-0.3125
-0.2901
-0.2757
-0.2983
-0.1104
-0.3094
Y-R. Wang, D.A.
Peters I Comput.
Methods Appl. Mech. Engrg.
134 (1996) 91-105
103
reeressine ---_
d__
advance ratio
1
Fig. 11. Blade
flapping
eigenvalues
(coupled
with wake system)
vs. advance
ratio
(real part)
:: d
L’
I ::
‘\ ‘.
-.
N
I ii
, 1
?_
1,
.
,
P
1 \
/
I ,
”
\
EEX 30’:-- __-+-;_*
\
i.
,.
,A----____----
/
:
__--
*_-,_*-
/
I/ .: 1; . :: ,‘:
i’
8’
,;
7 -2
”
%g
... -.-....
I
Fig. 12. Blade
flapping
r
,./’
.... . ..____...
‘S &,a-
eigenvalues
(coupled
with wake system)
vs. advance
ratio
(im. part).
agrees with the result in [13] and [14]. We also see that the wake will lower the flapping damping and cause earlier instability than the no inflow mode, which is also shown in Fig. 11. Another interesting thing is that, although these flapping mode branches state at /-L= 0 as regressing, progressing, differential and collective, it is not necessary for them to maintain this character as /L increases. The legends in those figures are just marked as their initial modes. One way to see which mode is dominant is to check the eigenvectors. For a specific advance ratio, we multiply the Happing angle eigenvectors of an eigenvalue with a ‘factor matrix’ as follows:
104
Y-R. Wang, D.A.
Peters I Comput. Methods Appl. Mech. Engrg.
134 (1996) 91-105
The left-hand side matrix are complex numbers, and their absolute values (amplitudes of shape modes) can give information on the identification of flapping modes. In Fig. 11, as p = 1.3, we pick up a set of eigenvalues (0.01587, 0.00497), (0.00544, O.O), (-0.00600, O.O), (-0.51459, 0.01239), and (-0.52398, 0.0). The branches begin with collective, progressing, regrssing, differential, and the lower branch of the progressing mode, respectively. The amplitudes of the vectors of the first eigenvalue are: (1.82542, 1.88935, 0.57906, 0.44757) We can see that, although this eigenvalue begins with collective mode, it is a strongly coupled with the differential mode near p = 0.22. Thus, it becomes partly a collective and partly a differential mode. For the same reason, the 4th eigenvalue begins with differential mode, but its amplitudes at p = 0.13 are: (1.29290, 1.08725, 0.13153, 0.07587). Again, this branch has become a coupled differential-collective mode. However, the 3rd eigenvalue has the following amplitudes: (2.81553, 0.58549, 4.96667, 0.39344) which shows that it begins and ends as a regressing mode. The amplitudes of upper and lower branches of progressing mode are: (0.16252, 0.97948, 1.75680, 0.09842) and (0.27714, 0.13044, 0.18527, 1.15635) We can see that they begin and end as progressing modes.
7. Conclusions The mode shapes of an inflow wake model and the eigen-system of coupled blade-wake aeroelastic model have been studied. The results of this work give us the stability information of this aeroelastic system and the understanding of the effect of inflow dynamics on rotor blade flapping dynamics in forward flight. The conclusions from this research are as follows: (1) The visualizing of mode shape vibration provides us with a better understanding of the rotor wake mechanism. (2) The eigen-analysis exposes the complicated coupling mechanism of the blade flapping and inflow modes. Due to this inflow coupling effect, the blade flapping modes may not keep their initial characteristics, as the skew angle or advance ratio changes. (3) The induced flow has a profound effect on blade dynamics even at high advance ratios. The effect may cause an earlier unstable condition for a helicopter rotor system.
Acknowledgments When the authors were in Georgia Tech, this research was supported in part by a joint NASA-Army grant, NAG-2-462, and additional funding was provided by the Georgia Tech Center of Excellence for Rotary Wing Aircraft Technology, sponsored by the Army Research Office. The continuing part of this work was supported by the National Science Council in the Republic of China under contract No. NSC 83-0424-E-032-001. Their support is gratefully appreciated by the authors.
References [l] T. Theodorsen, General theory of aerodynamic instability [2] R.G. Loewy, A two dimensional approach to the unsteady [3] S.T. Crews, K.H. Hohenemser and R.A. Ormiston, An (1973). [4] R.A. Ormiston, Application of simplified inflow models to 34-37.
and the mechanism of flutter, NACA R.496, 1935. aerodynamics of rotary wing, J. Aerospace Sci. 24 (1957) 81-98. unsteady wake model for a hingeless rotor, J. Aircraft lO(12) rotorcraft dynamic analysis, J. Am. Helicopter Sot. 21(3) (1976)
Y-R. Wang, D.A.
Peters I Comput. Methods Appl. Mech. Engrg.
134 (1996) 91-105
105
[S] D.A. Peters and C.J. He, Comparison of measured induced velocities with results from a closed form finite state wake model in forward flight, 45th Annual National Forum of the American Helicopter Society, Boston, MA, May 22-24, 1989. [63 C.J. He, Development and application of a generalized dynamic wake theory for lifting rotors, Ph.D. Thesis, School of Aerospace Engineering, Georgia Institute of Technology, July, 1989. [7] A. Su, Application of a state space wake model to elastic blade flapping in Hover, Ph.D. Thesis, School of Aerospace Engineering, Georgia Institute of Technology, August, 1989. [S] Y.-R. Wang, The effect of wake dynamics on rotor eigenvalues in forward flight, Ph.D. Thesis, School of Aerospace Engineering, Georgia Institute of Technology, June, 1992. [9] W. Kinner, Theory of circular wing, Ingenieur Archiv, 7 (1937) 47-80, Translation No. 2345, Ministry of Aircraft Production, UK. [lo] M. Joglekar and R. Loewy, An actuator disc analysis of helicopter wake geometry and the corresponding blade response, USAAVLABS TR 69-66, December, 1970. [ 1 l] W. Johnson, Helicopter Theory (Princeton University Press, Princeton, 1980) 471-480. [12] W. Kaplan, Operational Methods of Linear Systems (Addison Wesley, 1962) 472-488. [13] O.J. Lowis, The stability of rotor blade flapping motion at high tip speed ratios, Reports and Memoranda No. 3544, January, 1963. [14] W.E. Hall, Application of Floquet theory to the analysis of rotary wing VTOL stability, Ph.D. Dissertation, Stanford University, Februrary, 1970.