Flutter Control of a Lifting Surface via ViscoHysteretic Vibration Absorbers
 Author: Lacarbonara Walter, Cetraro Marek
 Organization: Lacarbonara Walter; Cetraro Marek
 Publish: International Journal Aeronautical and Space Sciences Volume 12, Issue4, p331~345, 30 Dec 2011

ABSTRACT
In this paper, a viscohysteretic vibration absorber (VA) is proposed to increase the flutter speed of an airfoil and enhance damping in the pre and postflutter regimes. The passive system consists of a parallel arrangement of a dashpot and a rateindependent hysteretic element, represented by the BoucWen differential model. The equations of motion are obtained and various tools of linear and nonlinear dynamics are employed to study the effects of the viscohysteretic VA in the pre and postflutter ranges.

KEYWORD
Viscohysteretic vibration absorber , Flutter control

1. Introduction
Flutter refers to a serious aeroelastic instability in both civil and military aircrafts. The flutter phenomenon is characterized by limitcycle oscillations (LCOs) of relatively high amplitude, exhibiting several harmonics. Modern fighter airplanes often carry many types and combinations of external wingmounted stores to satisfy multimission requirements. Such stores can reduce the flutter speed and thereby degrade the operational and mission effectiveness of combat airplanes. Due to the importance of flutter avoidance, considerable research has been conducted in the last decades to develop and assess the capabilities of various flutter control concepts (Marzocca et al., 2002; Njuguna, 2007). In the general context of studies on semiactive/active methodologies, an important aspect is associated with the spatial optimization of actuator and sensor parameters to facilitate control of targeted modes while providing rolloff of higherorder modes without the need for phaseinducing filters. In contrast, the complexity inherent in semiactive/ active strategies for flutter control is dramatically reduced in passive systems such as VAs. Further, these systems exhibit desirable faulttolerance characteristics.
A passive flutter control approach, known as the decoupler pylon, has been presented by Reed et al. (1980). The decoupler pylon isolates dynamically the wing from the store pitch inertia effects by means of softspring/damper elements assisted by a lowfrequency feedbackcontrol system that minimizes static pitch deflections of the store due to maneuvers and changing flight conditions. Windtunnel tests and analyses show that this relatively simple pylon suspension system provides substantial increases in flutter speed and reduces the sensitivity of flutter to changes in store inertia and center of gravity. The flutter characteristics of F16 and F17 flutter models equipped with decoupler pylon mounted stores are presented by Reed et al. (1980) and compared with results obtained on the same model configurations with active flutter control systems.
A semiactive flutter control strategy for highaspectratio wings has been presented by Hu and Zhou (2007), where multiple magnetorheological dampers are proposed. In their work, the results of the semiactive approach are contrasted with those obtained by the passive flutter control methodology based on multiple viscoelastic vibration absorbers (VAs).
A design methodology for optimized flutter control of an aeroelastic delta wing is presented by Richard and Clark (2003). Th ese experiments, employing optimized piezoelectric transducers, show substantially increased flutter control authority over nonoptimized systems. They also point out the importance of this spatial coupling as well as the transducer mass and stiff ness eff ects. Works related to semiactive techniques can be found in Schweiger et al. (1999), McGowan (1999), and Rocha et al. (2007), while Njuguna (2007) presents a good review of existing techniques.
In the last century, nonlinear VAs have been proposed by purposeful introduction of nonlinearities in either stiff ness or damping. This nonlinearity can be associated with various aspects of the VA restoring force. For example, according to Lee et al. (2007a, b) the flutter and postflutter control of a twodegreeoffreedom (twodof) system is realized by a singledof nonlinear energy sink (NES), which exhibits a nonlinear restoring force without the linear stiff ness term. The benefi t of this kind of VA lies in its ability to enhance the flutter speed or stabilize the postflutter LCOs, depending on the initial conditions, flow speed, and mass ratio. Three diff erent control phenomenologies can be seen in the NES, and these are related to the diff erent interactions of the NES with the pitch and heave modes of the twodof airfoil. In the first case, intermittent action of the NES results in intermittent oscillations of the threedof system, and thus
flutter is not completely suppressed. In the second case, flutter is completely suppressed, while in the third case, stable LCOs arise in the postflutter range. The results of this study, corroborated by experiments, showed a maximum increase of the flutter speed by 26 per cent.
With regard to material nonlinearities in the VA, Carpineto et al. (2010, 2011) and Vestroni et al. (2011) exploit the hysteresis exhibited by short wire ropes under fl exure to engineer a hysteretic VA for passive control of mechanical vibrations. In order to show the eff ectiveness of the proposed hysteretic VA, Carpineto et al. (2010, 2011) carry out a semianalytical/ numerical and experimental study of a simply supported beam equipped with the hysteretic VA, subject to harmonic and white Gaussian noise base motions.
In this work, a viscohysteretic VA is investigated for flutter control of a thin airfoil with a careful unfolding of its nonlinear features in the pre and postflutter regime.
2. The Governing Equations of Motion
Studies of the flutter condition are often carried out considering a twodof model of a wing treated as a lifting surface or a thin airfoil, accounting for the plunge and pitch degrees of freedom (Fig. 1). The flutter control strategy of this model is based on incorporation of a VA, located at a distance
ηb in the chordwise direction, which can oscillate by a displacementq normal to the chordwise direction (Fig. 1). Here,b represents half the chord length. From this, it can be seen that the augmented system is in fact a threedof system.Let (
_{1},e _{2}) be a fixed basis withe _{1} collinear with the chordwise direction in the wing’s stressfree configuration, while (e _{1},b _{2}) are bodyfixed unit vectors rotated in the current configuration by an angle denoted by α(b t ) (Fig. 1) so that the (counterclockwise) angle α denotes the pitch degree of freedom. In the stressfree configuration, the origin of the fixed frame is taken coincident with the elastic center C^{E}, while the center of mass and the aerodynamic center are denoted by C and C^{A}, respectively. The position vector of a material point that occupies position in the stressfree configuration of the surface is described in the current configuration byx where
=r ^{C}(t)h(t) ？e _{2}e so thatb _{1}(t)h(t) describes the plunge degree of freedom (see Fig. 1 for the definition of e). The position vector of the point that represents the VA is given by =r _{d}(t)h(t) ？ ηe _{2}b . In contrast,b _{1}(t)+qb _{2}is the position vector of the material point in the current configuration with respect to the center of mass.
Let
=b _{k} (α)？R , where _{k}e indicates the orthogonal tensor that describes the rotation of the surface. Given the fact that in this caseR describesR _{3} rotations, we can see thate _{1} = cos αb _{1} + sin αe _{2} ande _{2} = ？sin αb _{1} + cos αe _{2}. The velocity and acceleration of the material points of the surface and those of the VA, respectively, are given by:e where
and
Therefore, the time rates of change of linear momentum and angular momentum of the lifting surface are given by:
where
and
ρ the mass polar moment of inertia of the lifting surface about the center of mass C. Conversely, the time rate of change of the VA linear momentum is given by:J ^{C}Thus, the balance laws of linear and angular momentum of the augmented system become:
where
is the position vector of the aerodynamic center in the current configuration. The components of the elastic restoring forcer ^{A}and couple
acting on the lifting surface can be expressed as
The constitutive equation for the viscohysteretic absorber,
can be written as a direct summation of an elastic part
and a dissipative memorydependent part
according to
The dissipative part is, in turn, expressed as the summation of a linear viscous term and a purely hysteretic term z,
where the evolution of z is described by the BoucWen firstorder diff erential equation (Bouc, 1967; Lacarbonara and Vestroni, 2003; Wen, 1976).
In Eqs. (4) and (5), (
) refers to the external resultant force and moment with respect to the center of mass and (f, c ) denote the aerodynamic resultant force (lift and drag resultants) and moment, with respect to the aerodynamic centerf ^{L},c ^{L}C^{A} .Equations (46) are linearized, the incorporation of linear viscosity in the airfoil constitutive laws given by
At the same time, full nonlinearity is kept in the viscohysteretic constitutive law for the VA. The ensuing equations are:
where ρ
J ^{E} : = ρJ ^{C} + ρAe ^{2} andJ^{E}_{d} : =m_{d} (ηb )^{2} denote the polar mass moment of inertia with respect to the elastic centerC^{E} of the wing and of the VA, respectively, whilec^{E} : =df^{L} +c^{L} is the aerodynamic moment reduced to the elastic center.According to Glauert’s theory of thin airfoils (Glauert, 1947), to obtain the aerodynamic forces induced by a uniform airstream of velocity
V _{1} (with zero initial angle of attack), the eff ective angle of attack is first expressed as:e The lift force and aerodynamic moment, reduced to the aerodynamic center, are given by:
where s is the wing span, ρ^{a} is the fluid density, and
C^{α}_{M} = 0 for symmetric airfoils.Therefore, the equations of motion, as a result of neglecting other kinds of external forces and considering only the linearized parts of the inertial and aerodynamic forces, become:
The dimensional constitutive parameters of the viscohysteretic VA are denoted by (
k_{d}, k_{z}, c_{d}, β, γ, n ).The next few paragraphs illustrate the nondimensional form of the equations of motion of the augmented aeroelastic system. This form can be arrived at by dividing the vertical coordinate
h byb (i.e., the nondimensional plunge degree of freedom becomesthe VA degree of freedom
q byb (i.e.,the hysteretic variable z by z_{0} : =
k_{d}b , and introducing the characteristic timedenotes the torsional frequency of the airfoil). As a result of these steps, the following nondimensional form of the equations of motion can be obtained:
where the overbar is dropped in respect to
h and the dot over h indicates diff erentiation with respect to nondimensional timeThe most important nondimensional parameter in the design of the VA is the mass ratio μ : = ρ
A /m_{d} mμ:= md / ρA로 수식 수정해 주세요 between the VA mass and the wing mass. This parameter is very important because it scales the control force exerted by the VA on the lifting surface.Given the fact that the plunge frequency is
and the frequency of the VA by itself is
the following nondimensional frequencies ratios arise from the nondimensionalization:
(nondimensional plunge frequency) and
(nondimensional VA frequency). The nondimensional pitch frequency becomes 1 in this case. The mass and elastic distribution properties of the lifting surface are summarized by
and ε : =
e/b . The nondimensional damping coefficients of the airfoil are given by:In contrast, the remaining nondimensional constitutive parameters of the viscohysteretic VA (besides
can be expressed as:
The nondimensional velocity is
and the aerodynamic constant is given by
k_{u} : = ρsb ^{2}C ^{0}_{L} / ρA .3. Design and Performance of the Aeroelastic Viscohysteretic VA
The wellknown Den Hartog viscoelastic VA (Den Hartog, 1934; Frahm, 1911) is widely regarded as effective for vibration cosntrol of a mechanical system subject to a harmonic timevarying force. This effectiveness is attributed to the antiresonance phenomenon, which completely cancels the response when the system and the VA are undamped and the frequency of the VA is exactly tuned to the frequency of
the system driven to resonance. There is no exact resonance cancellation when the structure is damped, while the antiresonance phenomenon, accompanied by small amplitudes, is exhibited provided that the VA frequency is properly tuned to the frequency of the damped system.
Eight parameters must be determined to facilitate the design of a viscohysteretic VA for flutter and postflutter control. These parameters are the mass ratio μ, the position of the VA along the wing’s chordwise direction η, the VA frequency
which arises from the tuning condition of the VA’s overall stiffness, the VA damping ratio ζ
_{d} , and the other constitutive parameters of the hysteretic part of the VA restoring force, namely,When the hysteretic part of the restoring force z is discarded, the (dimensional) frequency of the VA can be expressed as
While ω
_{d} represents the linear frequency of the VA (i.e., at infinitesimal oscillation amplitudes), the nonlinear stiffness variations of the hysteretic VA during its finite cycles are such that its (nonlinear) frequency scales with ω^{E} , whose nondimensional counterpart, ω^{E} / ω_{d} turns out to beDuring flutter, it is expected that the VA will oscillate together with the appropriately phased airfoil so as to introduce additional damping into the augmented system.
A rational design of the VA consists of requiring ω
^{E} to be equal to the flutter frequency. This tuning condition serves as a good initial guess for the optimization process, which results in the optimum tuning. Once the flutter speed of the system without VA is calculated employing the Routh Hurwitz criterion (Meirovitch, 1970), the frequency of theflutter mode can be computed accordingly.
The mass ratio μ and the VA position described by η must satisfy certain physical restrictions. The weight limitation on the wing is such that the value of μ can be as high as 1/100 (an upper bound for the VA mass), while placing the VA at threequarters of the chord toward the trailing edge leads to
These restrictions represent limitations on the magnitude of the control force. This force is greatly enhanced by higher mass ratios and by a position closer to the tip of the wing’s trailing edge, where the torsional couple exerted by the VA on the wing profile is maximized.
In relation to the VA damping ratio ζ
_{d} , an optimum damping ratio can increase the flutter speed by orders of magnitude, while the hysteretic part of the force does not influence the onset of flutter. Hysteresis is shown to have a fundamental role in increasing the decay rates of the transients. This occurs mostly in the preflutter regime, but is also seen in the postflutter regime.By adopting the parameters of the wing model studied by Behal et al. (2006), the application of the RouthHurwitz criterion leads to calculation of the nondimensional linear flutter speed of the airfoil by itself found to be
The model in question is shown in Tables 1 and 2.
The optimal parameters for the linear viscoelastic VA can be determined by a numerical search calculating the flutter speed via the RouthHurwitz criterion on a grid that discretizes the parameter plane
in a lattice. Figure 2 shows the results of these computations, which lead to the following optimal parameters:
and ζ
_{d} = 0.36. With these optimal parameters, the flutter speed goes up to 8.3; an increase of approximately 240 percent. The sensitivity of the flutter speed with respect toand ζ
_{d} can be appreciated by computing the flutter speed at meaningful values detuned from the optimal values. For example, if the VA frequency is tuned with the frequency of the flutter mode,and if the damping ratio is set to ζ
_{d} = 0.1, flutter speed is determined to be 4.5 and the increase with respect to the case without VA is only 83 per cent. In contrast, ifand ζ
_{d} = 0.1, the flutter speed goes up to 5.9, with an increase of approximately 140 per cent.Numerical investigations via integration of the equations of motion into the effects of the viscohysteretic VA (i.e., the parallel arrangement of a dashpot and a BoucWen element) have shown that the purely hysteretic part of a system is unable to alter the flutter boundary of a system endowed with a purely viscoelastic VA. However, it does increase the effective damping of the system in the preflutter regime, as
shown in Fig. 3. The decay rate of the airfoil response to initial conditions is greatly enhanced by hysteresis. Moreover, the decay rate increases with the magnitude of the initial conditions up to a certain threshold, at which the effective damping is maximized.
It is interesting to study the effects of detuning in the optimal parameters of the viscoelastic VA while observing the effects of variations of the constitutive parameters of the hysteretic part of the VA restoring force. For example, in Fig. 4 the viscoelastic parameters are set to
and ζ
_{d} = 0.1. In such a case, flutter is reached at a flow speed of 4.5, with an increase with respect to the system without VA of only 83 per cent. However, in this detuned case there is a postflutter range in which the airfoil exhibits LCOs with relatively small amplitudes. Hence, the hysteresis has a twofold effect: it enhances damping in the preflutter stage and controls the postflutter up to a threshold value of the speed beyond which the response of the airfoil diverges if all other structural /aerodynamic nonlinearities are neglected. Some examples of such postflutter responses are shown in Figs. 57 at various flow speeds whenTo assess the effects of the purely hysteretic component of the VA restoring force, the dashpot contribution is neglected and a set of postflutter responses are shown in Figs. 810 with
at various flow speeds. While in the immediate vicinity of flutter, the LCO has one dominant frequency (the nonlinear flutter frequency); multiple harmonics of the flutter frequency are manifested at speeds further away from the flutter condition. More interestingly, when the flow speed is 3.66 times the flutter speed, the LCO undergoes a Hopf bifurcation by which the postcritical LCO becomes amplitudemodulated (i.e., quasiperiodic motion), as shown in Fig. 10.
4. Discussion and Conclusions
In this paper, a viscohysteretic VA is proposed in order to enhance the aeroelastic stability of an airfoil by increasing its flutter speed. The passive system is a parallel arrangement of a dashpot and a rateindependent hysteretic element, represented by the BoucWen differential model (Bouc, 1967; Lacarbonara and Vestroni, 2003; Wen, 1976).
This paper shows that, for a set of parameters representing an experimental airfoil (Behal et al., 2006), the optimized purely viscoelastic VA (without the hysteretic element) can result in an increase of the flutter speed by up to 240 per cent when the mass ratio is 1 per cent. Numerical investigations are carried out into the effects of the purely hysteretic VA, i.e. a VA without the dashpot. These results show that a purely
hysteretic restoring force does not have the capability to change the flutter boundary of the airfoil, since at the onset of flutter the aeroelastic destabilizing forces can only be counteracted by a VA control force in which the resistive component depends on the relative instantaneous velocity. In contrast, the hysteretic restoring force can only dissipate energy over a full cycle of oscillation. However, the presence of a purely hysteretic restoring force enhances damping in the preflutter regime and controls the postflutter behavior over a significant range of flow speeds.
The optimization described in this paper can be carried out in the presence of a viscohysteretic force by applying the RouthHurwitz criterion, accounting for only the linear viscoelastic part of the restoring force. The primary effect of the hysteretic part of the restoring force is to increase the effective damping of the system in the preflutter regime, as shown in Fig. 3. In practice, the decay rate of the airfoil response to initial conditions is greatly enhanced by hysteresis. In the optimal configuration of the absorber, decay rates can only be controlled by hysteresis. This justifies its use in conjunction with the linear viscoelasticity of the spring and dashpot. The experimental validation of the control system is part of an ongoing research project.

[Fig. 1.] Lifting surface with the coordinate frame (top) and freebody diagram of the forces (bottom). The main system is augmented by the viscohysteretic vibration absorber shown in the enlarged box (left).

[Table 1.] Dimensional wing model parameters

[Table 2.] Nondimensional wing model parameters

[Fig. 2.] Variation of the flutter velocity Vcr of the augmented threedof system, comprising the linear viscoelastic vibration absorber with ωd and ζd as obtained by the RouthHurwitz criterion.

[Fig. 3.] From left to right: (first line) time history of h, (second line) time history of α, (third line) time history of q, (bottom) loops of the total vibration absorber (VA) force, given respectively by Nd = 2ζd√δωdq and Nd = ωd2[δq + (1δ)z]+2ζd√δωdq in t∈[0,5000] for the threedof system with the linear viscoelastic VA (left) and with the viscohysteretic VA (right). The flow velocity is V = 0.999 Vcr.

[Fig. 4.] Variation of the flutter speed with the hysteretic constitutive parameter β with the vibration absorber (VA, thin solid line) and without the VA (thick solid line), when ζd= 0.1 while the other parameters are fixed to δ = 0.2, n = 1, ωd = 1, and γ = β. The shaded region represents the postflutter range, where the lifting surface exhibits limitcycle oscillations due to the nonlinear hysteretic force in the VA. The boundary is defined as the speed at which the torsional angle reaches the threshold amplitude of 0.2 rad, beyond which the response diverges.

[Fig. 5.] From left to right: (top) time history of α and phase diagrams (α, α) in the transient phase (t∈[0,600] and t∈[0,200], respectively) and at steady state (t∈[4.8？104, 5？104] and t∈[3？104, 5？104], respectively) (middle) FFT of α and q, (bottom) loops of total vibration absorber force Nd = ωd2[δq + (1δ)z]+2ζd√δωdq in the transient phase (t∈[0,600]) and at steady state (t ∈[3？104, 5？104]), when ζd = 0.1 and V = 1.02 Vcr with Vcr = 4.5. The peak in the FFT is attained by f = 0.14. f denotes the nondimensional frequency.

[Fig. 6.] From left to right: (top) time history of α and phase diagrams (α, α) in the transient phase (t ∈ [0,600] and t ∈[0,200], respectively) and at steady state (t∈[4.8？104, 5？104] and t∈[3？104, 5？104], respectively), (middle) FFT of α and q, (bottom) loops of total vibration absorber force in the transient phase (t ∈[0,600]) and at steady state (t ∈[3？104, 5？104]), when ζd = 0.1 and V = 1.33 Vcr with Vcr = 4.5. The peak in the FFT is attained by f = 0.132.

[Fig. 7.] From left to right: (top) time history of α and phase diagrams (α, α) in the transient phase (t ∈ [0,600] and t ∈[0,200], respectively) and at steady state (t∈[4.8？104, 5？104 ] and t∈[3？104, 5？104], respectively), (middle) FFT of α and q, (bottom) loops of total vibration absorber force in the transient phase (t ∈[0,600]) and at steady state (t ∈[3？104, 5？104]), when ζd = 0.1 and V = 1.64 Vcr with Vcr = 4.5. The peak in the FFT is attained by f = 0.14.

[Fig. 8.] From left to right: (top) time history of α and phase diagrams (α, α) in the transient phase (t ∈ [0,600] and t ∈[0,200], respectively) and at steady state (t∈[4.8？104, 5？104] and t∈[3？104, 5？104], respectively), (middle) FFT of α and q, (bottom) loops of total vibration absorber force in the transient phase (t ∈[0,600]) and at steady state (t ∈[3？104, 5？104]), when ζd = 0 and V = 1.05 Vcr with Vcr = 2.1. The peak in the FFT is attained by f = 0.146.

[Fig. 9.] From left to right: (top) time history of α and phase diagrams (α, α) in the transient phase (t ∈ [0,600] and t ∈[0,200], respectively) and at steady state (t∈[4.8？104, 5？104] and t∈[3？104, 5？104], respectively), (middle) FFT of α and q, (bottom) loops of total vibration absorber force Nd = ωd2[δq+(1δ)z] in the transient phase (t ∈[0,600]) and at steady state (t ∈[3？104, 5？104]), when ζd = 0 and Vcr = 2.38 Vcr with Vcr = 2.1. The peak in the FFT is attained by f = 0.136.

[Fig. 10.] From left to right: (top) time history of α and phase diagrams (α ,α) in the transient phase (t ∈[0,600] and t ∈[0,200], respectively) and at steady state (t ∈[4？104, 5？104] and t ∈[3？104, 5？104], respectively), (middle) FFT of α and q, (bottom) loops of total vibration absorber force in the transient phase (t ∈[0,600]) and at steady state (t ∈[3？104, 5？104]), when ζd = 0 and V = 3.66 Vcr with Vcr = 2.1. The peak in the FFT is attained by f= 0.122