Modal Filtering for Control of Flexible Aircraft

NASA/TM—2013–216597

Modal Filtering for Control of

Flexible Aircraft

Peter M. Suh Dryden Flight Research Center, Edwards, California Dimitri N. Mavris Georgia Institute of Technology, Atlanta, Georgia Click here: Press F1 key (Windows) or Help key (Mac) for help

S eptember 2013

This page is required and contains approved text that cannot be changed.

NASA STI Program ... in Profile

Since its founding, NASA has been dedicated  CONFERENCE PUBLICATION.

to the advancement of aeronautics and space Collected papers from scientific and science. The NASA scientific and technical technical conferences, symposia, seminars, information (STI) program plays a key part in or other meetings sponsored or helping NASA maintain this important role. co-sponsored by NASA.

The NASA STI program operates under the  SPECIAL PUBLICATION. Scientific, auspices of the Agency Chief Information Officer. technical, or historical information from It collects, organizes, provides for archiving, and NASA programs, projects, and missions, disseminates NASA’s STI. The NASA STI often concerned with subjects having program provides access to the NASA substantial public interest.

Aeronautics and Space Database and its public interface, the NASA Technical Reports Server,  TECHNICAL TRANSLATION.

thus providing one of the largest collections of English-language translations of foreign aeronautical and space science STI in the world. scientific and technical material pertinent to Results are published in both non-NASA channels NASA’s mission.

and by NASA in the NASA STI Report Series, which includes the following report types: Specialized services also include organizing and publishing research results, distributing  TECHNICAL PUBLICATION. Reports of specialized research announcements and completed research or a major significant feeds, providing information desk and personal phase of research that present the results of search support, and enabling data exchange NASA Programs and include extensive data services.

or theoretical analysis. Includes compila- tions of significant scientific and technical For more information about the NASA STI data and information deemed to be of program, see the following: continuing reference value. NASA counter- part of peer-reviewed formal professional  Access the NASA STI program home page papers but has less stringent limitations on at http://www.sti.nasa.gov manuscript length and extent of graphic presentations.  E-mail your question to help@sti.nasa.gov  TECHNICAL MEMORANDUM.  Fax your question to the NASA STI Scientific and technical findings that are Information Desk at 443-757-5803 preliminary or of specialized interest, e.g., quick release reports, working  Phone the NASA STI Information Desk at papers, and bibliographies that contain 443-757-5802 minimal annotation. Does not contain extensive analysis.  Write to: STI Information Desk  CONTRACTOR REPORT. Scientific and NASA Center for AeroSpace Information technical findings by NASA-sponsored 7115 Standard Drive contractors and grantees. Hanover, MD 21076-1320

NASA/TM—2013–216597

Modal Filtering for Control of

Flexible Aircraft

Peter M. Suh Dryden Flight Research Center, Edwards, California Dimitri N. Mavris Georgia Institute of Technology, Atlanta, Georgia Insert conference information, if applicable; otherwise delete Click here: Press F1 key (Windows) or Help key (Mac) for help National Aeronautics and Space Administration Dryden Flight Research Center Edwards, CA 93523-0273

S eptember 2013

Enter acknowledgments here, if applicable.

Click here: Press F1 key (Windows) or Help key (Mac) for help Available from: NASA Center for AeroSpace Information 7115 Standard Drive Hanover, MD 21076-1320 443-757-5802 Click here: Press F1 key (Windows) or Help key (Mac) for help

Abstract

Modal regulators and deformation trackers are designed for an open-loop fluttering wing model. The regulators are designed with modal coordinate and accelerometer inputs respectively. The modal coordinates are estimated with simulated fiber optics. The robust stability of the closed-loop systems is compared in a structured singular-value vector analysis.

Performance is evaluated and compared in a gust alleviation and flutter suppression simulation.

For the same wing and flight condition two wing-shape-tracking control architectures are presented, which achieve deformation control at any point on the wing.

Nomenclature

ᡓ number of accelerometers ᠧ state matrix Accel accelerometer AFS active flutter suppression ᠧ balanced and reduced-order plant state matrix 〩 ᡔ balanced ᠨ control input matrix BWB blended wing body ᠨ balanced and reduced-order input matrix 〩 ᠨ disturbance input weighting matrix for balanced and reduced-order plant ⡩ ᠨ control input matrix for balanced and reduced-order plant ⡰ ᡕ number of control surfaces ᡕᆑ mean aerodynamic chord ᠩ state space output matrix ᠩ balanced and reduced-order plant output matrix 〩 ᠩ output matrix transformation for accelerometers ぀ ᠩ state-regulated goal matrix for balanced and reduced-order plant ⡩ ᠩ output matrix for balanced and reduced-order plant ⡰ ᡖ vector of deformations at fiber optic measurement stations ᡖ vector of virtual measurement deformations ぀ ᡖ vector of reference deformations at virtual measurement locations ぅ〲〳 ᡖ(ᡶ , ᡷ , ᡸ , ᡲ) vector of deflections and rotations defined over ᡶ , ᡷ , ᡸ at time ᡲ 〄 〄 〄 〰 〰 〰 ᡖ䙣(ᡶ , ᡷ , ᡸ , ᡲ) vector of accelerations and angular accelerations defined over ᡶ , ᡷ , ᡸ at time ᡲ 〄 〄 〄 〰 〰 〰 ᠰ direct feedthrough matrix DFRC Dryden Flight Research Center ᠰ control-regulated goal matrix for balanced and reduced-order plant ⡩⡰ ᠰ output noise weighting matrix for balanced and reduced-order plant ⡰⡩ ᡗ number of goal-regulated states and control surface movements EI Effective Independence EOM equations of motion ᠲ lower linear fractional transformation 〹 ᡙ structural damping coefficient ᠳ plant GLA gust load alleviation ᠳ (ᡱ) class of all generalized plants with multiplicative input uncertainty 【 ᠳ (ᡱ) class of all generalized plants with multiplicative output uncertainty 〖 ᠳ gust model 〴えうぇ ᠳ (ᠷ) matrix transfer function from disturbance to regulated output over all controllers ᠷ ぐこ ᠳ (ᡱ) matrix transfer function from disturbance to regulated output ぐこ ぇ〵 ᡠ scalar elements of the ᡡ column and row of ᠴ 〶〶 ᠴ sensor projection matrix (also referred to as hat matrix) ᠴ h infinity control type ⦘ ⠵2 imaginary number ᠵ identity matrix ᠵᡥ(5) imaginary part ᠶ matrix defined as ᡒᡃ ᡣ reduced frequency or non-dimensional frequency ᡣ system stability margin ぀ ᠷ controller ᠷ control gain matrix 〰 ᠷ filter gain matrix 〲 ᡤ row vector index ᡤ row vector index for virtual deformation references ぅ ᠸ global stiffness matrix ᠸᠱ leading edge LFT linear fractional transformation LQG linear quadratic Gaussian LSE least-squares estimation LTI linear time-invariant ᠸᠱ reference leading edge deflection at wing tip ぅ〲〳 ᠸ elemental stiffness matrix 〲 㕠 ᠸ modal stiffness matrix ᡥ column vector index ᡥ total mass of the wing 〡 ᡥ index of reference modes ぅ ᠹ global mass matrix ᠹ elemental mass matrix 〲 ᠹ indexing of modal coordinates 〶ぁ〱〲け ᠹ multiplicative noise input matrix え 㕡 ᠹ modal mass matrix ᡦ number of states ᡦ noise vector added to plant outputs 〱 ᡦ number of states in reduced order model ぅ ᡦ noise vector added to plant inputs え N number of degrees of freedom in finite element model NASA National Aeronautics and Space Administration ᡧ number of measurements (or outputs) ᡁ lower linear fractional transformation of generalized plant and ᠷ ᡂ generalized plant ぇ〵 ᡂᡥᡥ percent modal mass of the total mass ᡥ for the ᡡ mode shape 〶 〡 ᡂ generalized plant with multiplicative input uncertainty 【 ᡂ generalized plant with multiplicative output uncertainty 〖 ᡩ vector of modal coordinates ᡩ(ᡲ) vector of modal coordinates at time ᡲ ᡩ estimated modal coordinates from modal filter in simulation 〲うぇ th (ᡲ) ᡩ ᡡ modal coordinate scalar at time ᡲ 〶 ᡩ indexed modal coordinates for input to controller ぀ ᡩ vector of reference modal coordinates ぅ〲〳 ᡩ free-stream dynamic pressure ⦘ ᡃ left singular matrix ᡉ あ ᡃ(ᡣ) generalized aerodynamic force matrix ᡰ number of reference deformations ᡰ vector of reference inputs to controller ぅ〲〳 ᡰ relative uncertainty magnitude at steady state あ ᡰ relative uncertainty magnitude at high frequency ⦘ ᡄ similarity transformation matrix ᡄᡗ(5) real part RS robust stability ᡱ(ᡶ, ᡷ, ᡸ, ᡲ) strain measurement vector at time ᡲ ᡅ strain-displacement matrix SMI structural modal interaction SVD singular-value decomposition ᡅ solution of control algebraic Ricatti equation ⦘〰 ᡲ simulation time ᡆ output complementary sensitivity function ᡆ input complementary function 【 TE trailing edge TE reference trailing edge deflection at wing tip ref ᡳ vector of control signals ᡳ(ᡲ) vector of control signals at time ᡲ ᡳ reference unit vector of deformations and rotations 〱 ᡳ control signal noise ぁ ᡳ output to uncertainty block ⥀ ᡇ right singular matrix of ᠶ ᡴ input signal to controller ᡈ left singular matrix of ᠶ Vel velocity V-g velocity versus structural damping coefficient V- velocity versus frequency ᡵ vector of disturbances ᡵ(ᡲ) vector of disturbances at time ᡲ ᡵ (ᡲ) gust velocity at time ᡲ 〴 ᡵ 䙢 (ᡲ) gust acceleration at time ᡲ 〴 ᡵ (s) weight functions on diagonals of ᡉ and ᡉ 〶 【 〖 ᡉ diagonal input frequency weighted matrix 【 ᡉ diagonal output frequency weighted matrix 〖 ᡉ controllability Gramian matrix 〰 ᡉ observability Gramian matrix あ ᡶ state vector of plant model ᡶ Cartesian coordinate 〄 ᡶ(ᡲ) state vector of a plant model at time ᡲ ᡶ (ᡲ) state vector of actuator states belonging to ᡶ(ᡲ) 〨〰ぇ ᡶ (ᡲ) state vector of modal aero lags belonging to ᡶ(ᡲ) 〨う ᡶ (ᡲ) balanced and reduced-state vector of plant model at time ᡲ 〩 ᡶ (ᡲ) state vector of modal coordinate displacements belonging to ᡶ(ᡲ) い ᡶ (ᡲ) state vector of modal coordinate velocities belonging to ᡶ(ᡲ) い䙢 ᡶ initial state at time 0 あ ᡶ䙢(ᡲ) time derivative of ᡶ(ᡲ) ᡶ䙢 (ᡲ) time derivative of ᡶ (ᡲ) 〩 〩 ᡶ 㕈(ᡲ) estimated states of balanced and reduced-order plant ᡷ output vector of plant model ᡷ Cartesian coordinate 〄 ᡷ(ᡲ) output vector at time ᡲ ᡷ input to uncertainty block ⥀ ᡸ goal-regulated outputs ᡸ Cartesian coordinate 〄 ᡸ(ᡲ) goal regulated outputs at time ᡲ ᡒ left singular matrix of ᡉ 〰 ‑ air density Z diagonal matrix of singular values of ᠶ [ structured uncertainty block [ structured input uncertainty block 【 [ structured output uncertainty block 〖 ぇ〵 ‒ ᡡ uncertainty 〶 – normal error  vector of eigenvalues ぇ〵 ’ ᡡ system eigenvalue 〶 _ diagonal matrix of eigenvalues † structured singular-value vector ‥ spectral radius ぇ〵 … ᡡ Hankel singular value of the system triple (ᠧ, ᠨ, ᠩ) 〶 … ㍤ ( 5) maximum singular values time constant th ‸ (ᡶ , ᡷ , ᡸ ) ᡡ natural mode shape vector defined over ᡶ , ᡷ , ᡸ 〶 〄 〄 〄 〄 〄 〄 f matrix of natural mode shapes, ‸ (ᡶ , ᡷ , ᡸ ) 〶 〄 〄 〄 th ′ (ᡶ , ᡷ , ᡸ ) ᡡ natural strain mode shape vector defined over ᡶ , ᡷ , ᡸ 〶 〄 〄 〄 〄 〄 〄 i matrix of strain mode shapes, ′ (ᡶ , ᡷ , ᡸ ) 〶 〄 〄 〄 ″ natural frequency — viscous damping ratio : total column vector index (that is, all columns included)

Introduction

Active flutter suppression (AFS) and gust load alleviation (GLA) may be necessary additions to the control systems of future flexible aircraft. One of the ways to reduce the weight and improve the fuel efficiency of a vehicle is to remove structure, reducing the structural stiffness of the wing. Reducing wing stiffness may lower the natural frequencies of the aircraft. As elastic frequencies shift lower, they may coalesce with rigid-body and gust frequencies. Aircraft with wing stiffness reduction modifications may require active structural control to maintain stability, improve ride quality, and alleviate gust and maneuver loads (ref. 1). These aircraft may also be more susceptible to dangerous phenomena such as buffeting or flutter.

Flutter is an aeroelastic phenomenon which is manifested when energy from the air excites the structure, and the structure moves, adding energy back to the air. This transfer of energy back and forth tends to grow, leading to an increase in structural amplitudes until structural failure occurs (ref. 2). All aircraft are subject to flutter at some particular flight condition; aircraft do not generally operate within this flutter envelope due to its destructive nature.

For this reason, AFS is generally not required in aircraft which operate today. Nevertheless, some aircraft still require suppression of marginally-stable modes and suppression of control system interactions with modes. Modal suppression is often achieved through adding mass ballast or putting notch filters in the controller (refs. 3 and 4). These methods, however, tend to add weight or reduce controller bandwidth, respectively. Some researchers are seeking AFS solutions as an alternative to passive methods (ref. 5). A recent study (ref. 6) performed by Lockheed Martin (Bethesda, Maryland) in 2010 on the SensorCraft vehicle showed that a 25-percent total weight savings could be achieved by incorporating AFS into the controller (ref. 7).

Increasing flexibility in the structure of the wing may come with other consequences. In simulation studies, drag has been shown to increase for a flexible aircraft, relative to a stiffer one, especially during maneuvers (ref. 8). As a potential solution, some researchers have suggested that the flexibility of an aircraft could be shaped to reduce this effect (ref. 9). This will naturally require adequate structural control effectors. Research has been conducted with structural-shaping technologies such as adaptive control surfaces (ref. 10), adaptive ribs (ref. 11), and piezoelectrics (ref. 12).

Research into of flexible aircraft has recently become of interest to the National Aeronautics and Space Administration (NASA). The blended wing body (BWB) aircraft X-56A (Lockheed Martin) was selected as a platform in November of 2011 (ref. 13) for experimental applications and innovation. As a first application, the X-56A aircraft control engineers at the NASA Dryden Flight Research Center (DFRC) (Edwards, California) must demonstrate AFS and GLA on the aircraft at subsonic speeds. These requirements must be met because the X-56A aircraft has several body freedom flutter modes within its flight envelope and its wings are highly flexible.

Body freedom flutter is characterized by an interaction of rigid-body and low-frequency elastic modes. Body freedom flutter occurs more often in aircraft with low tail volume or high fineness ratios (ref. 14). By learning to actively control these unstable structural modes, similar methods may be applied to full-scale aircraft in the future to realize weight savings.

The X-56A aircraft is also a distributed control and sensing platform, having discrete arrays of accelerometers placed along the wing. There are also plans for the layout of distributed strain sensors along the wings, known as fiber optic sensors with fiber Bragg gratings (ref. 15). With its many control effectors and distributed sensing the X-56A aircraft will likely make an excellent platform for wing-shape optimization research. Research has already been conducted at NASA DFRC with fiber optic sensors and fiber Bragg gratings. These sensors were implemented on the wing of the NASA Ikhana Predator B unmanned aircraft (General Atomics Aeronautical Systems, San Diego, California); the sensors feed back strain and shape data at high sample rates (ref. 15).

In this paper, it is demonstrated that distributed sensors such as the fiber optic sensors with fiber Bragg gratings can be used for the control of flexible aircraft. The strain or shape data may be passed to a modal filter which can estimate key states of the structure (ref. 16). These states may then be passed as a signal to an active controller with modal coordinate inputs (refs. 17 and 18). A controller with accelerometer inputs is also developed. The two controllers are juxtaposed and qualitatively rated on the basis of order reduction, robustness, flutter suppression, and disturbance rejection. The purpose of the comparison is to discover whether a controller using a modal filter has at least comparable performance to a controller using point-sensor feedback. Point sensors, such as accelerometers, could be considered to be the state-of-the-art sensors for experimental AFS studies (refs. 19 and 20).

After comparison of the controllers, two virtual deformation-tracking design case studies are demonstrated on the open-loop fluttering wing model. A subsonic flight condition was chosen in which the model is open-loop unstable in order to simulate what it may be like to perform the same type of deformation control on the aeroelastically-unstable X-56A aircraft. The modal filter is shown to either predict deformations at unmeasured points (virtual deformation estimation) or to estimate modes directly related to deformations, as in prior research (refs. 21 and 22).

In the first virtual deformation-tracking design case study, a reference virtual deformation signal is tracked indirectly from modal coordinate references. In the second study, a large number of modes are estimated by the modal filter, and the deflections are predicted and tracked directly. The work described in this paper is potentially of interest because to our knowledge, few if any researchers have attempted to track deformations (especially virtual deformations) on a fluttering three-dimensional wing.

Background

Since the modal filter has not yet been implemented on a flexible aircraft for control system feedback, solving the problem of integrating of the modal filter into the aeroservoelastic design may be the primary contribution of this paper. Various flutter control architectures have been developed for AFS, such as the aerodynamic energy concept (ref. 23), classical (ref. 24), linear quadratic Gaussian (LQG) (refs. 7 and 25), ᠴ (ref. 26), eigenspace (ref. 27), minimax (ref. 28), ⦘ predictive (ref. 29), neuro-adaptive (refs. 2 and 30), and many others. Most of these designs are supported by a few point sensors for feedback. Control designs (refs. 31-34) for AFS and GLA have traditionally relied on discrete arrays of single-point sensors, such as accelerometers, placed at selected locations on the aircraft wing or body. To our knowledge, there is currently no ongoing flutter research that utilizes a modal filter in the controller. Numerous control-related studies have been performed using the modal filter as an estimator in air-off structural systems (refs. 16-18, 35, and 36).

The modal filter has quite an extensive history in the structures discipline; a brief review is offered below to illustrate the significance of this filter in this field. The modal filter has also seen recent advances with regard to distributed sensing; this topic is also discussed below.

Introduction to the Modal Filter The modal filter was initially developed to solve a problem in active structural control related to the structural modal interaction (SMI) problem in aircraft (ref. 4). Aerospace engineers encounter SMI in aircraft when flexible modes are excited by control surface feedback that is meant to excite only rigid-body motion. The parallel problem in structures is known as observation spillover, which exhibits many of the same mechanisms (ref. 16). Observation spillover was observed by structural control experts utilizing observers to estimate the states for model-based controllers. They found that Luenberger (ref. 37) and Kalman (ref. 38) -type observers had the potential to feed residual modes back to the control system when reliance was placed on arrays of single-point sensors. Because single-point sensors pick up residual mode information, the controller tended to generate a control signal in response to residual modes, increasing the amplitude of the residual modes. As the amplitudes of residual modes increased, sensors picked up more residual information, and a positive feedback system resembling a flutter interaction developed, leading to instability.

Research into spillover reduction moved in several directions, including incorporating more constraints into the estimator design. One such constraint was to make the control effectors move orthogonally to the residual modes, thus, not exciting the residual modes (ref. 39).

Another direction taken to reduce spillover was to develop estimators for control which were naturally resistant to observation spillover using the concept of mutual modal orthogonality.

These estimators are known as modal filters. The first relied on the modal expansion theorem (ref. 40), which states that for structures with distinct natural frequencies, the deformation of all nodes of the structure ᡖ(ᡶ , ᡷ , ᡸ , ᡲ) can be set equal to a linear combination of the structure’s 〄 〄 〄 natural mode shapes {‸ (ᡶ , ᡷ , ᡸ ), ᡡ = 1 … ᡀ} [see eq. (1)]: 〶 〄 〄 〄 〕 (ᡲ)‸ (1) ᡖ(ᡶ , ᡷ , ᡸ , ᡲ) = 㔳 ᡩ (ᡶ , ᡷ , ᡸ ) 〄 〄 〄 〶 〶 〄 〄 〄 〶 (t) where the time-varying weights q of each mode shape may be found as shown in ⤙ equation (2) 〡 ᡩ (ᡲ) = ‸ (ᡶ , ᡷ , ᡸ ) ᠹᡖ(ᡶ , ᡷ , ᡸ , ᡲ) (2) 〶 〶 〄 〄 〄 〄 〄 〄 〕×〕 〕×⡩ ぇ〵 where ᠹ ᒈ ጷ is the global mass matrix, ‸ (ᡶ , ᡷ , ᡸ ) ᒈ ጷ is the ᡡ mode shape, and 〶 〄 〄 〄 〕×⡩ ᡖ(ᡶ , ᡷ , ᡸ , ᡲ) ᒈ ጷ is a vector of all deflections and rotations from every node in the structure.

〄 〄 〄 In 1985, Meirovitch and Baruh (ref. 17) presented the first modal filter based on the modal expansion theorem presented above. The finite-element method was used and interpolation was assumed to measure all degrees of freedom with a small set of sensors. Meirovitch also showed analytically that the filter operated well for input to a controller in a noisy environment.

Contemporary researchers thought the modal filter was impractical due to the requirement of a highly accurate analytical mass matrix (ref. 16).

In 1991 two types of modal filters were published by Shelley et al. (ref. 36): the pseudo-inverse modal filter and the modal filter with reciprocal vectors. The pseudo-inverse modal filter relies on a least-squares estimation (LSE) formulation. By forming the error equation and solving for ᡩ(ᡲ), equation (1) can be transformed into the LSE problem shown in equation (3) (ref. 36), 〡 ⡹⡩ 〡 ᡩ(ᡲ) = (f f) f ᡖ(ᡶ , ᡷ , ᡸ , ᡲ) + – (3) 〄 〄 〄 where modal matrix f = {‸ (ᡶ , ᡷ , ᡸ ), ᡡ = 1 … ᡀ} is calculated from the air-off eigenvalue 〶 〄 〄 〄 problem. The modal matrix is typically column-truncated to lower-frequency mode shapes which span the space of the deformation. The Shelley et al. reciprocal modal filter approach utilized experimentally-derived vectors. Given the measured actuator forces over a time history, the reciprocal vectors could be experimentally derived through an optimization procedure. The reciprocal vectors can be shown to be orthogonal to all modes except the desired mode. By forming a dot product of the measurements with any reciprocal vector, the corresponding modal coordinate could be calculated.

Both the pseudo-inverse and reciprocal modal filters were experimentally tested on a five-meter truss, and were shown to estimate modal states well, although the pseudo-inverse approach underperformed compared to the experiment-based reciprocal vector approach.

Rigid-body modes and very-high-frequency modes were also estimated poorly by the pseudo-inverse modal filter, possibly due to the use of a small number of accelerometers. The reciprocal approach was also utilized by Shelley et al. for structural health monitoring (ref. 41) in 1993. The approach was also shown to work in an adaptive control framework for control of the Big Darby Creek bridge structure (ref. 35) in Ohio in 1995.

Modal Filtering and Smart Sensors For good modal estimation, the structure should be fitted with many sensors (ref. 16). Since most modal filters tend to rely on an LSE calculation, and every sensor produces data containing a certain amount of error, more sensors tend to give a better estimate of the state.

Recently, researchers have favored smart sensors because they are typically lightweight and non-intrusive, facilitating distribution on the structure. Modal sensors with shaped polyvinylidene fluoride film were implemented by Lee and Moon (ref. 42) in 1990. The linear combiner was developed by Sumali et al. (ref. 43) in 1997 which relied on piezoelectric patches to estimate modal coordinates and to adaptively predict mode shapes.

Distributed strain-based modal filters with fiber optics were implemented by Kang et al.

(ref. 22) in 2007. The modal coordinates were used to estimate the deflections at any point of a cantilever beam. The strain mode filter looks very much like the pseudo-inverse modal filter given in equation (3) except the regressor matrix is formed by first multiplying the mode shapes by the strain-displacement matrix ᡅ (ref. 44). Since the strain-to-displacement relationship is linear, and equation (1) holds generally for structures with distinct natural frequencies (ref. 39), the least-squares solution of the strain mode matrix can be used to estimate the same modal coordinates calculated in equation (3), as shown in equation (4) (ref. 22), 〡 ⡹⡩ 〡 ᡩ(ᡲ) = (i i) i ᡱ(ᡶ , ᡷ , ᡸ , ᡲ) + – (4) 〄 〄 〄 where i = {ᡅ‸ (ᡶ , ᡷ , ᡸ ), ᡡ = 1 … ᡀ} is a ᡀ × ᡀ matrix with ᡀ strain modes, and 〶 〄 〄 〄 〕×⡩ ᡱ(ᡶ , ᡷ , ᡸ , ᡲ) ᒈ ጷ is the measured strain.

〄 〄 〄 The strain-mode filter clearly allows so-called virtual deflections (deflections at points which are not measured directly) to be measured. The equation can be represented in equation (5) (ref. 22).

〡 ⡹⡩ 〡 ᡖ(ᡶ , ᡷ , ᡸ , ᡲ) = f(i i) i ᡱ(ᡶ , ᡷ , ᡸ , ᡲ) (5) 〄 〄 〄 〄 〄 〄 Equation (5) demonstrates that once modal coordinates are estimated, they may be used to estimate deflections, known as virtual deformations, at any point on the structure. Others have used this concept to estimate deflections with just a few strain gauges (ref. 21). Kang et al.

found that the modal coordinates were not precisely estimated using the technique given in equation (5) due to small model errors. Treiber et al. (ref. 45) in 2008 attempted to improve upon the Kang et al. modal filtering approach with a Kalman filter. Their work showed that the filtering approaches are sensitive to structural uncertainty. Related approaches to deformation estimation such as Tessler and Spangler’s (ref. 46) strain-based least-squares variational method and Ko and Richards’s (ref. 47) Bernoulli-Euler beam bending method have also been developed. Photogrammetry (ref. 48) on the F/A-18 airplane (McDonnell Douglas, now the Boeing Company, Chicago, Illinois) Active Aeroelastic Wing has been shown to adequately capture deflections.

Wing Model

The development of the modal filter and deformation estimation was reviewed to illustrate the rich background of the modal filter and its relationship to distributed sensing. The modal filter may provide the most benefit when it is applied to a finite-element-type structure. To demonstrate the modal filtering technique, a three-dimensional wing model with 199 nodes was developed. The finite-element wing model structural characteristics were tailored to have characteristics similar to a high-aspect-ratio wing, and the model was fitted with several control surfaces and realistic actuators which could be used for AFS and GLA. A modal analysis was completed for the wing model to show some of the relevant modes which were estimated by the modal filter. The aeroservoelastic analysis is also presented to familiarize the reader with flutter and how the flight condition of the model is selected.

Geometry The wing model half-span was modeled as a 3.354- by 0.838-m aluminum rectangular plate clamped at the wing root, and made from the 6061-T6 aluminum metal alloy traditionally utilized in aircraft flying today. The dimensioned model is presented in figure 1.

Figure 1. Specifications and layout of wing test model.

Aluminum ribs and spars were added to reinforce the structure, giving it structural properties similar to a realistic aircraft wing. The leading-edge spar was made thicker than the trailing-edge spar, as it would likely carry more of the lift load. To simulate more load-bearing potential at the wing root, the spars and rib dimensions were linearly tapered toward the wing tips. The control surface panels shown outlined in red dots in figure 1 were also stiffened with leading-edge and trailing-edge spars, as well as cross-wise ribs along the dotted lines. This stiffening was done to remove the possibility of a control surface flutter interaction. Each control surface was connected to the wing structure by two six-degree-of-freedom springs. The connection joined the wing and control surface at 1/3 and 2/3 of the control surface span from the leading edge corner of a control surface. A Cartesian coordinate system was defined for the model, with origin at the lower left corner of the wing. The coordinate ᡶ is defined as increasing from the 〄 trailing edge to the leading edge of the wing. The coordinate ᡷ increases from the wing root to 〄 the wing tip. Coordinate ᡸ increases out of the page. The natural modes of the structure, 〄 presented below, were all assumed to move only in the direction of ᡸ .

〄 Modal Analysis The natural modes of a structure are shapes, which the structure will tend to take when excited at the structure’s natural frequencies (ref. 39). They are found through the eigensolution of the matrix structural equations of motion (EOM). The formulation of the EOM requires knowledge of the global mass matrix ᠹ and stiffness matrix ᠸ . The global mass and stiffness matrices may be built up from elemental mass ᠹ and stiffness ᠸ matrices using the 〲 〲 finite-element method (ref. 49). The matrix equation is presented as equation (6).

ᠹᡖ䙣(ᡶ , ᡷ , ᡸ , ᡲ) + ᠸᡖ(ᡶ , ᡷ , ᡸ , ᡲ) = 0 (6) 〄 〄 〄 〄 〄 〄 ㄔ2㄄ぇ By assuming ᡖ = ᡖ ᡗ the eigensolution of the EOM may be calculated, resulting in ⡨ eigenvalues _ = ᡖᡡᡓᡙ(’ , ᡡ = 1 … ᡀ) and eigenvectors f = {‸ (ᡶ , ᡷ , ᡸ ), ᡡ = 1 … ᡀ} . The 〶 〶 〄 〄 〄 eigenvectors can be shown to be linearly independent since they are mutually orthogonal to each other for structures with distinct natural frequencies (ref. 39). The eigenvectors are also referred to as mode shapes.

Since the ᡀ × 1 mode shapes are of same dimension as the deformation and are orthogonal to each other for structures with distinct natural frequencies, the eigenvectors also span the space of all possible deformations. Thus, they form a basis for the deformation, and the structural deformation may be expressed as a linear combination of mode shapes, as discussed above. Aircraft may have higher-order modes with frequencies near to each other due to the complexity of the structure, but orthogonality is generally still a good assumption even for structures with close frequencies or repeated frequencies (see ref. 39).

The wing model presented in figure 1 was modeled with the finite-element method. The model was designed to be a structure with distinct natural frequencies. An eigenvalue solution of equation (6) was performed to give the mode shapes and natural frequencies. The modes were then normalized with respect to the global mass matrix. When normalizing in this way, the mode shapes can be said to be orthogonal to the mass matrix. Ten modes are used in the simulations; the mode shapes and natural frequencies for the first four modes are presented in figures 2(a) through 2(d).

Figure 2. Modal representation of wing model: a) first wing bending; b) first wing torsion; c) second wing bending; and d) second wing torsion.

From visual inspection the first four modes are: first wing bending, figure 2(a); first wing torsion, figure 2(b); second wing bending, figure 2(c); and second wing torsion, figure 2(d). The maximum deformation for each mode tends to stay below 0.2 m. Each mode shape reveals that the control surfaces moved together with the rest of the structure in the static mode shape; this action is made possible by setting the stiffness in the actuators to extremely large values, as displayed in figure 1, effectively clamping the control surface movement.

The control surfaces also have mode shapes that are referred to as control modes (not shown), which are formed by enforcing a static 1-deg rotation on the actuators (see ref. 50) with actuator stiffness set to very low values. After the mode shapes (including the natural and control modes) are calculated, it is important to predict the behavior of the modes as the air interacts with them. This interaction will characterize the nature of the model during operations.

Flight Condition The aeroelastic behavior of the wing during flight is important for many reasons, including determining the flutter speed boundary due to the structural interaction with airflow. To identify the flutter boundary (the speed at which the aircraft begins to flutter) a velocity (V) versus damping (g) (V-g) and   V-) analysis is computed by solving the generalized eigenvalue problem shown in equation (7) and assuming harmonic motion, (that is, ㄔ2㄄ぇ ᡩ = ᡩ ᡗ ) of the wing modes (ref. 51): ⡨ ⡰ ‑ᡕᆑ 1 + ⠵2ᡙ 㕡 + 㕠䚃 ᡩ = 0 (7) 䚂 ᠹ ᡃ(ᡣ) ㎘ ᠸ ⡰ ⡰ 2ᡣ ″ 㕡 㕠 where ᠹ is the modal mass matrix, ᠸ is the modal stiffness matrix, ᡙ is the structural damping coefficient, ″ is the natural frequency, ‑ is the air density, ᡕᆑ is the mean aerodynamic chord, ᡣ is ㄄〰ᆑ the reduced frequency or non-dimensional frequency defined by ᡣ = , and ᡃ(ᡣ) is the ⡰〣 generalized aerodynamic force matrix. Equation (7) is a generalized eigenvalue problem solely ⡩⡸ㄔ2〴 in terms of reduced frequency, ᡣ . The eigenvalues are defined as ’ = and equation (8) ㄘ ㄄ shows the relationship of the variables in equation (7) to the eigenvalues (ref. 52).

1 ᠵᡥ(’) ″ᡕᆑ ″ = , ᡙ = 2— = , ᡈ = (8) ᡄᡗ(’) 2ᡣ 㒓ᡄᡗ(’) where — is the viscous damping ratio. By specifying a range of reduced frequencies, ᡣ , the eigenvalues of equation (7) may be calculated repeatedly and variations of frequency and structural damping coefficients with velocity may be found using equation (8). These values are then plotted against each other.

Correct interpretation of the plots is important. In the structures discipline, negative damping is defined to be good; whereas within the controls discipline, positive damping is defined to be good. This difference can cause confusion when interpreting the structural damping coefficient over the speed regime. In the case at hand, positive structural damping ( — > 0 ) corresponds to velocity of an unstable flight condition. The V-g plot is presented in figure 3.

Figure 3. Velocity versus damping plot of wing model.

One of the purposes of the V-g plot is to identify the nature of flutter in the aeroelastic model.

The characteristic of a typical two-mode flutter interaction on a V-g plot is that one mode will become more and more unstable and the mode it interacts with will become more and more stable (ref. 52). Physically, this represents a situation in which energy from the air is transferred to the structure moving it, and the movement of the structure adds energy back to the air, forming a positive feedback system as described above. As the airspeed increases past the flutter speed, the strength of the flutter tends to become stronger and more destructive, and more system damping is usually required (notice that ᡙ grows past the flutter boundary in figure 3).

For simulation studies it is desired that the flight condition be selected substantially past the flutter speed to demonstrate the use of modal filtering in flutter suppression. The V-g plot indicates that the first wing bending (mode # 1) passes the stability boundary at 72 m/s. It interacts with the first wing torsion (mode # 2) pushing it into the highly-stabilizing region, which is typical of a clamped plate bending/torsion flutter mode. Because only these two modes interact, only their frequencies are plotted in the V- in figure 4. The variation of frequency over the speed range is also typical for rectangular plate aeroelastic analysis (ref. 52). For the studies to be conducted on this wing model, the flight condition is selected to be at a speed of 80 m/s at an altitude of 305 m, where at least 2.8 percent structural damping is required and flutter frequency (see fig. 4) is approximately 4.42 Hz, which is the frequency shown for the first mode at 80 m/s. The damping ratio predicted with this frequency matching method tends to be incorrect at conditions above and below flutter because aerodynamic behavior is typically not harmonic.

Figure 4. Velocity versus frequency plot of wing model.

Methodology

The methodology for incorporating the modal filter into the control design is brief, as the primary requirement is to update the state space model so that it reflects the use of the modal filter. Fulfilling this requirement is actually as simple as choosing which modal coordinates one wishes to measure and updating the output matrix to reflect this choice. The update of the output matrix when accelerometers are utilized is given for comparison and because the accelerometers are used later for a benchmark of the controller with modal filtering. A guideline for the incorporation of the modal filter into the controller is also given with a bias toward control of flexible aircraft. This guideline also serves as a roadmap for the work described in this paper.

State Space Model – Modal Filter or Accelerometers The aeroservoelastic wing model described above may be represented by a linear time-invariant (LTI) system of finite dimensions. The LTI system may be cast into linear constant coefficient differential equations which are valid only at the flight condition selected in figure 3.

The Mach number and free-stream dynamic pressure ᡩ are required to be known to formulate ⦘ the state space equations. Assuming known quantities, the aeroservoelastic state space model can assume the matrix time derivative form for simulation and control design shown in equation (9) (ref. 50): ᡶ䙢(ᡲ) = ᠧᡶ(ᡲ) + ᠨᡳ(ᡲ) ᡶ(0) = ᡶ (9) あ ᡷ(ᡲ) = ᠩᡶ(ᡲ) + ᠰᡳ(ᡲ) with the initial state ᡶ(0) = ᡶ . The ᡦ -dimensional vector ᡶ(ᡲ) is referred to as a state vector and あ at any time during a simulation can be accessed to give the current “state” of the system. The ᡧ -dimensional vector ᡷ is the system measurements. According to equation (9), only the current state and the ᡕ -dimensional input ᡳ is required to know the state in the next time step. The ᠧ , ᠨ , ᠩ and ᠰ matrices are real constant matrices with ᡦ × ᡦ , ᡦ × ᡕ , ᡧ × ᡦ and ᡧ × ᡕ dimensions.

The exact formulation of these matrices is not the focus of this study. The matrices are derived from the finite-element nodal model, linear aero-paneling and frequency-based methods, rational function approximations, and modal transformations. An excellent development of the state space formulation is given in reference 50. The states, however, are the focus of this study. The state space vector is composed of vectors of states with different units. The state vector may be defined as shown in equation (10), 〡 ᡶ(ᡲ) ᕤ {ᡶ (ᡲ), ᡶ (ᡲ), ᡶ (ᡲ), ᡶ (ᡲ)} (10) い い䙢 〨う 〨〰ぇ ぀×⡩ ぀×⡩ where ᡶ (ᡲ) ᒈ ጷ is a vector of modal coordinate displacements, ᡶ (ᡲ) ᒈ ጷ is a vector of い い䙢 〲×⡩ modal coordinate velocities, ᡶ (ᡲ) ᒈ ጷ is a vector of aerodynamic lag states and 〨う こ×⡩ ᡶ (ᡲ) ᒈ ጷ may be formed from a vector of actuator accelerations, velocities, and 〨〰ぇ displacements. The modal filter may be used to sense some or all of the states in the ᡶ (ᡲ) い vector, as shown in equation (3). Any of these states (or combinations of these states) can be used to form the output matrix if the appropriate sensor is utilized. In turn, the output matrix is used during control design.

Modal Coordinate Sensor Output Matrix The output matrix must be modified to incorporate the modal filter in the control design. The output matrix ᠩ is a matrix of row vectors relating the output sensor to the state vector ᡶ . Since the modal filter directly measures some or all of the ᡶ state, the form of the output matrix may い be cast into the form shown in equation (11), assuming all modal coordinates are measured.

ᠩ ᕤ [ᠵ 0 0 0 ] (11) ぀×぀ ぀×぀ ぀×〲 ぀×こ Measuring all modes is not always required. Higher-order modes which do not significantly contribute to the overall modal deformation in the system may sometimes be ignored and the row dimension of the output matrix may be reduced.

Accelerometer Sensor Output Matrix The output matrix may also be modified to account for accelerometer measurements. Rather than measuring directly any part of the state vector, accelerometers measure linear combinations of modal coordinate accelerations. The relationship between the accelerometers and the state of the system is often modeled as shown in equation (12).

ᡷ(ᡲ) ᕤ ᠩ ᡶ䙢(ᡲ) = ᠩ (ᠧᡶ(ᡲ) + ᠨᡳ(ᡲ)) (12) ぀ ぀ 〨×ぁ The matrix ᠩ ᒈ ጷ has as many rows as accelerometers, ᡓ, and as many columns as the ぀ state space vector size. The matrix ᠩ is defined as shown in equation (13) (ref. 53).

぀ ᠩ ᕤ [0 f (ᡤ, : ) 0 0 ] (13) 〨×぀ 〨×぀ 〨×〲 〨×こ ぀ where f (ᡤ, : ) is the displacement modal matrix row indexed by ᡤ .

〨×぀ Modal Filtering Control Design The above section clarified the way in which the modal filter could be implemented into the state space model. This takes care of the modeling side of things. But questions linger over which modes should be selected to update the output matrix in equation (11). It is implied through the use of the modal filter that distributed sensing will be utilized. But it is not clear if some sensors should be excluded or not. There are multiple types of modal filtering techniques, some of which cater to particular types of sensors, and the determination of which sensor platform to use should be based on the application. There are many options for design of the controller with the modal filter. The following ten steps were loosely followed in this study: 1. Identify which modes are significant to the response of the model based on a percent modal mass calculation and the V-g analysis.

2. Update the output matrix with the identified modes as in equation (11).

3. Place or select sensors, either using an optimal sensor placement technique or an intuition of which modes are being measured.

4. Set up the modal matrix and the pseudo-inverse of the modal matrix corresponding to sensor locations ᡤ , in order to use equations (3)-(5) for modal coordinate or displacement estimation.

5. Reduce the order of the aeroservoelastic plant.

6. Specify the performance requirements for the control design.

7. Design a controller for the plant and verify that performance requirements are met.

If they are not, iterate steps 1-6.

Optional steps for virtual deformation tracking: 8. Upgrade the plant with an integral tracking state and design a modal tracking or deformation tracking controller.

9. Define optimal or desired reference deflections for virtual deformation tracking.

10. Verify that the selected modes for tracking produce an indirect deformation tracking controller, with minimal deformation error due to residual modes. If the error is too high, track more residual modes or redesign the controller or the actuator placement so that the residual modes are not strongly excited.

A few of these steps are discussed because they may help to illuminate some of the paths taken below and they may also provide space for complementary work. Liberties were taken during this study, and the steps were not precisely implemented as listed. Step 1 has many options. It is probable that a modal contribution percentage (see ref. 54) to flutter is the best way to define the significant modes, although this specific technique was not used directly in this study. Step 3 was followed for both accelerometer and fiber optic sensor placement; however, the fiber was placed using intuition and the accelerometers were placed using optimal sensor placement (ref. 55). Step 9 was not implemented in the present study. Finding desirable shapes of the aircraft wing/body may be an entire research topic in itself, likely requiring computational fluid dynamics studies. Step 10 was not implemented because the number of control surfaces was limited to two; thus, only two modes could be tracked at one time.

Design

This section focuses on the design of the modal-filter-type controller, its comparison to a controller with accelerometer inputs, and virtual deformation case studies.

To ensure that the instability of the model is indeed flutter, the wing model is simulated in an open-loop study when perturbed from equilibrium. The wing model is upgraded with either accelerometers or a modal filter (with fiber optic inputs), and plants [see Eqs. (9) and (11)-(13)] are designed for both systems. Both plants are different only in that one is using a modal filter and the other is using accelerometers for output. The plants are reduced in order and differences are noted for each.

An ᠴ optimal regulator is designed for each plant with the objective of being robustly stable ⦘ to input and output multiplicative uncertainty and having good disturbance rejection properties.

The simulation is presented for time-history analysis and the performance during a disturbance is demonstrated, showing that controllers with either modal coordinate input or accelerometer input can be used for GLA. The same controllers are also shown to be capable of AFS.

The modal coordinate controller is updated with integral tracking states corresponding to its measured modal coordinates. Two methods of virtual deformation tracking are presented. The first virtual deformation tracking case study uses a modal approximation of the deformation tracking signal. The second virtual deformation tracking case study uses the modal filter to predict the deformations [see eq. (1)] of the plant model, in order to achieve improved deformation tracking. The open-loop flutter of the wing model is presented first.

Open-Loop Flutter The aeroelastic wing model is in an open-loop flutter condition, based on the prediction made by the V-g analysis (see fig. 3) at the freestream velocity of 80 m/s and altitude of 305 m. To observe the characteristics of the flutter instability, the model was perturbed from equilibrium at time 0 by a unit deflection command to the control surfaces which lasted for 0.01 s. The time history of the modal coordinates is shown in figure 5.

Figure 5. Wing model in open-loop flutter.

The modal amplitudes oscillate at 4.49 Hz, which was almost predicted by the V-   (see fig. 4). The damping ratio for this mode (from a controls perspective) was -5.3 percent, making it more unstable than what was predicted in the V-g analysis (see fig. 3). The fluttering mode is predominantly a contribution from the first wing bending and first wing torsion mode.

The frequencies of the two modes have coalesced into that of a single frequency, which is typical of a flutter mode. The time to double amplitude is 0.46 s. The bending mode reaches higher amplitudes than the torsion mode, which is expected due to its lower natural frequency.

The torsion mode is slightly out of phase with the bending mode. From the above results (see figs. 3 and 5), it is apparent that the wing will experience strong flutter at the chosen flight condition.

Sensor System Design Two regulators and corresponding estimators for the fluttering wing model are now described. The regulators are given accelerometer inputs and modal coordinate inputs, respectively. Several objectives were defined to guide the controller designs so that the controllers could be qualitatively compared. The first objective was to stabilize the plant. As indicated above, the wing experiences strong open-loop flutter at the flight condition. The second objective was to reduce the controller order as far as possible. This has the effect of later improving numerical stability of the estimator and reducing the computational burden in the flight computer. A third objective was for the controller to have good disturbance rejection properties, especially from low-frequency turbulence. Passenger comfort can be improved by actively rejecting gust disturbances (ref. 56). The fourth objective of the controller design was that it be robust to modeled uncertainty. A controller designed about a linearized model rarely performs the same way in practice as it does in the laboratory (ref. 57). Below, the control design is described in detail, beginning with accelerometer placement and selection.

Accelerometer Placement The first controller was designed with two accelerometers. The accelerometers were placed using Kammer’s Effective Independence (EI) procedure (ref. 55), in which only the first wing bending and first wing torsion modes were maintained in the modal matrix. The EI procedure 〡 ⡹⡩ 〡 begins by forming the sensor projection matrix, ᠴ = f(f f) f from the modal matrices.

Large diagonals, ᡠ of the projection matrix ᠴ correspond to relatively important sensor 〶〶 locations ᡡ on the finite-element wing model. The EI sensor reduction procedure begins by removing the row of the modal matrix, f corresponding to the smallest ᡠ . The reduced 〶〶 projection matrix is then recalculated and the row corresponding to the smallest ᡠ is again 〶〶 removed. This reduction process can be iterated to the desired number of sensors, which are then used to form the index ᡤ used in equation (13). The EI procedure tends to select sensors that carry the highest amount of the desired modal information. The resulting sensor placements are shown in figure 6.

Figure 6. Accelerometer placement on wing model.

By down-selecting the rows of the modal matrix to the wing tips, the EI procedure indicates that the wing tips (both leading edge and trailing edge) carry the most modal information.

“The most modal information” means that the first wing bending and wing torsion modes had the highest deformation at these points; also see figures 2(a) and 2(b).

The EI procedure is only a first step and was used for its simplicity. The EI procedure contains flaws in that it does not attempt to optimize sensor layout for residual modes. It has been shown that residual modes tend to corrupt the sensor signal, leading to observation spillover (ref. 16). Thus, optimization procedures, such as the modified EI procedure (ref. 58), have been developed to place sensors on a BWB-type aircraft while minimizing residual mode information. The EI procedure is also sensitive to structure with many nodes, because nodes tend to be very close together and many good locations may fall in the same spot. The severity of this problem can be reduced by using correlation matrices to select sensor locations that maximize modal information without redundancy (ref. 59).

Modal Coordinate Selection The second controller is designed with modal coordinate feedback. The modes which dominate the deformation of the structure and can be easily controlled and observed should be selected for feedback. The modal mass, which may be calculated for each mode, is an indicator for whether a mode dominates the structural deformation. The percent of the total modal mass ぇ〵 for the ᡡ mode may be calculated as shown in equation (14), 〡 1 ‸ (ᡶ, ᡷ, ᡸ) ᠹᡳ 〶 〱 (14) ᡂᡥᡥ = 〶 〡 ᡥ ‸ (ᡶ, ᡷ, ᡸ) ᠹ‸ (ᡶ, ᡷ, ᡸ) 〡 〶 〶 where ᡥ is the total mass of the structure, and ᡳ is an ᡀ × 1 reference vector of unit 〡 〱 deflections and rotations. Modes with the highest percent modal mass should be selected for feedback, since they may contribute heavily to the response (ref. 39). A good rule of thumb is to include modes which sum to approximately 90 percent of the total mass of the structure. The percent of the modal mass for each mode as calculated by equation (14) is given in figure 7.

Figure 7. Percent modal mass per mode shape.

The percent modal mass of the first bending mode is clearly the highest at 85 percent, indicating that it will be a mode which both contributes significantly to the modal response and can likely be easily controlled and observed (ref. 52). It is difficult to state an exact measure of observability and controllability of modes, as the measure of observability and controllability will be determined by the placement of the control effectors and sensors (see refs. 30 and 60).

Thus, another way must be found to decide which modes should be included for measurement.

Recall that the airflow interacts with the structure (see fig. 3) near flutter so that some modes tend to dominate the structural deformation more than others. The flutter analysis at the selected flight condition indicates that the torsion mode will also be highly mobile near flutter. As such, based on the observation of the high modal mass of mode 1 and the flutter interaction due to mode 2, both first wing-bending and first wing-torsion modal coordinates were selected for feedback. For more complicated structures it might be beneficial to utilize the percent contribution of each mode to flutter as a criterion. Pak shows how to expand the flutter mode as a linear combination of the natural modes and calculate the percent contribution of each mode to flutter (ref. 54).

Fiber Optic Sensor Placement The means of estimating the modal coordinates is probably one of the most important aspects of the modal coordinate feedback controller design. The sensors could come in the form of piezoelectric materials (ref. 43), fiber optics (ref. 22), strain gages (ref. 21), or, potentially, photogrammetry (ref. 61). Each sensor type has different characteristics, which may make some sensors more appropriate for certain systems than others. Since fiber optic sensors had already been tested on aircraft (see ref. 15) this type of sensor was selected for the true sensor feedback to the modal filter. The modeled layout on the wing model is shown in highlighted green in figure 8.

Figure 8. Fiber optic sensor placement on wing.

At the locations where the fiber is placed, the deflections are assumed to be directly measured during the simulations. Several methods have been developed for the purpose of estimating deflection from strain. Equation (5) shows how strain can be directly used to achieve this objective. Intuition was used to place the fibers so that necessary modal information could be estimated. The selected modes for feedback include strong bending and torsion effects. To capture sufficient bending information, the sensors are placed span-wise along the entire wing.

To capture torsional effects, three fibers are placed chord-wise. The spacing between each location was set at one-half-inch intervals, which is the same spacing used on the NASA Ikhana Predator B unmanned aircraft (ref. 15). The use of a strain-shape algorithm need not be required to use the fiber optic sensors for modal coordinate estimation. Work by Kang et al.

(ref. 22) has shown that fiber optic strain measurements can be utilized to estimate modal coordinates. But instead of a modal matrix f , a strain mode matrix i is formed, which can then be utilized as shown in equation (4).

Controller Design After sensor locations are selected for both sensor types: accelerometers and fiber optics, the process of control design may begin. Here this is demonstrated by reducing the order of the plant, respectively for each sensor type. The controllers are then designed around each reduced order plant. The robust stability is analyzed for each closed-loop system with two structured uncertainty cases. The model reduction approach is presented first.

Model Reduction There is a multitude of literature on the topic of reducing the order of a controller; two references are given here (refs. 62 and 63). Some researchers reduce the plant and design the controller around the reduced-order plant, others reduce the controller after it has been designed around the full plant. It is typical to preserve the dominant eigenvalues in the reduced-order model. A common method of order reduction is to first balance the plant and then reduce it based on the Hankel singular values (ref. 59). This method was selected for the current research; a short description is provided here. The Hankel singular values {’ 㐙… = 㒓’ (ᡉ ᡉ ), ᡡ = 1 … ᡦ㐣 are derived from the eigenvalues , ᡡ = 1 … ᡦ} of the square root of 〶 〶 〰 あ 〶 the product of the controllability Gramian ᡉ and observability Gramian ᡉ . The Gramians are 〰 あ found from a solution of algebraic equations known as Lyapunov equations, shown in equation (15) (ref. 59): 〡 〡 ᠧᡉ + ᡉ ᠧ + ᠨᠨ = 0 〰 〰 (15) 〡 〡 ᠧ ᡉ + ᡉ ᠧ + ᠩ ᠩ = 0 あ あ The Gramians give a degree of relative observability and controllability if the plant is internally balanced. To balance the plant a transformation on the states ᡶ = ᡄᡶ may be found so that the 〩 controllability and observability Gramians are both diagonal and equal. The diagonality means that each state has its own independent measure of controllability and observability. The equality of the Gramians indicates that each state is equally controllable and observable (is excited to the same degree to which it is sensed). The transformation is found by decomposing the solutions of equation 15: ᡉ and ᡉ of the unbalanced system using a 〰 あ singular-value decomposition (SVD). The left singular matrix ᡒ of ᡉ may be multiplied with left 〰 singular matrix ᡃ of ᡉ , as in reference 59: あ ᠶ = ᡃᡒ (16) 〡 Another SVD of ᠶ may be performed to arrive at ᠶ = ᡈZᡇ . The transformation matrix may then be found as discussed in reference 59 and presented in equation (17).

⡹⡩ 㐕 ⡰ (17) ᡄ = ᡒᡇZ The transformed state space matrices are given by the method shown in equation (18).

⡹⡩ ᠧ = ᡄ ᠧᡄ 〩 ⡹⡩ (18) ᠨ = ᡄ ᠨ 〩 ᠩ = ᠩᡄ 〩 To reduce the order of the system, states (rows or columns) of the balanced system may be removed which correspond to relatively low Hankel singular values. From the reduced-order model, the reduced-order controller and estimator may be designed as it would be from the original plant matrices.

The discussion above on the topic of control order reduction illustrates the manner in which sensor selection can affect the plant order reduction and resulting controller order reduction.

The transformation matrix ᡄ is directly dependent on the SVD of ᠶ . The matrix ᠶ is, in turn, directly dependent on the SVD of ᡉ . Equation (15) shows a direct relationship between the あ output matrix ᠩ and ᡉ . The use of either modal coordinates or accelerometers affects the form あ of the output matrix and thus will affect the relative Hankel singular values through the eigensolution. To illustrate this effect more clearly, the relative unit-normalized Hankel singular values of the balanced systems with accelerometer outputs and modal coordinate outputs are given in figure 9.

Figure 9. Comparison of normalized Hankel singular values for two plants.

The first plant with accelerometer outputs has significant unit-normalized Hankel singular values out to state 14. The Hankel singular values for the second plant with modal coordinate outputs show a very steep drop-off after the fifth state. It is not precisely known why this occurs.

The steep drop-off may be accounted for by the fact that no relationship in the output matrix is given for modes past the first two modal states. Since the higher modes are not as observable to the system, their input-output contribution may be less, which seems to indicate that modal filtering may focus the observability in the system. This assumption has merit, since the modal filter was first introduced to reduce observation spillover.

The presentation of the Hankel singular values in this form may indicate that the controller order could be reduced based only on the relative magnitude of the singular values. It was found, however, that proceeding thus blindly could lead to an unstable controller. A more rigorous approach was taken, by reducing the order of the controller by removing states corresponding to the lowest Hankel singular values until the reduced-order controller performance diverged significantly from the performance of the original full controller (see ref. 59). The pole comparison of the two plants is given in figures 10(a) and 10(b).

Figure 10. Order reduction: (a) plant with accelerometer output; and (b) plant with modal coordinate output.

The plant with accelerometer outputs was not reduced as far as was the plant with modal coordinate outputs. The final order of the plant with accelerometer outputs was 27. The final order of the plant with modal coordinate outputs was 5. The plant poles shown in figure 10 qualitatively compare well to many structures, in which lightly-damped modes are very near the imaginary axis (ref. 39). Sometimes these lightly-damped modes can become unstable due to interaction with a flight control system, described previously as SMI. For the present case, aerodynamic coupling (see fig. 3), is the cause of the pole migration to the right half-plane (RHP). Simulation results in figure 5 indicated that two modes are unstable in the model. The poles in figure 10 which have moved into the RHP characterize this instability as a flutter instability. One of the primary objectives of the present work is to actively suppress this flutter mode. That process is described below.

Ⅴ ↑↖ↈ↑↖↑∂∇ Optimal Control Many choices for control design are available once a state space model has been defined as in equation (9). The ᠴ optimal controller (ref. 57) was chosen for this study because it has a ⦘ wide range of applicability. One excellent feature is that it is not assumed that the disturbances are collocated with the control inputs. The locations of controlled outputs are not necessarily collocated with the location of system performance as they are in the LQG controller design (ref. 59). The ᠴ optimal controller was designed for both reduced-order plants shown in figures ⦘ 10(a) and 10(b). The reduced-order plants were cast into the state-space realization shown in equation (19) (ref. 59): ᡶ䙢 (ᡲ) = ᠧ ᡶ (ᡲ) + ᠨ ᡵ(ᡲ) + ᠨ ᡳ(ᡲ) 〩 〩 〩 ⡩ ⡰ ᡸ(ᡲ) = ᠩ ᡶ (ᡲ) + ᠰ ᡳ(ᡲ) (19) ⡩ 〩 ⡩⡰ ᡷ(ᡲ) = ᠩ ᡶ (ᡲ) + ᠰ ᡵ(ᡲ) ⡰ 〩 ⡰⡩ ぁ ×ぁ ぁ ×ぐ ㊓ ㊓ ㊓ where ᠧ ᒈ ጷ is the reduced-order state matrix, ᠨ ᒈ ጷ is the disturbance matrix, 〩 ⡩ ぁ ×〰 〲×ぁ 〲×え ㊓ ㊓ ᠨ ᒈ ጷ is the control input matrix, ᠩ ᒈ ጷ is the state-regulated goal matrix, ᠰ ᒈ ጷ ⡰ ⡩ ⡩⡰ ぀×ぁ ぀×ぐ ㊓ is the control-regulated goal matrix, ᠩ ᒈ ጷ is the measurement matrix, and ᠰ ᒈ ጷ ⡰ ⡰⡩ ぁ ×⡩ ㊓ is the measurement noise matrix. The reduced-order states ᡶ (t) ᒈ ጷ are driven by the 〩 ぐ×⡩ 〰×⡩ disturbances ᡵ(ᡲ) ᒈ ጷ and the control inputs ᡳ(ᡲ) ᒈ ጷ .

The goal of the ᠴ optimal control methodology is to find the controller ᠷ which minimizes ⦘ the ᠴ norm of the transfer function ᠳ (ᠷ) from disturbance ᡵ to regulated output ᡸ , over all ⦘ ぐこ possible controllers (ref. 59). Recall that the ᠴ norm of a transfer function is defined as the ⦘ supremum or least upper bound of the transfer function over all frequencies. So, more simply stated, the ᠴ synthesis routine results in a controller which best suppresses the peak of ⦘ ᠳ (ᡱ) , where ᡱ is the Laplace variable.

ぐこ The suboptimal solution of the ᠴ problem requires that two algebraic Ricatti equations be ⦘ solved in which the observer and the controller matrices are coupled by an inequality constraint on the spectral radius ‥ . Another constraint is also included to assure that the Hamiltonian matrices do not have eigenvalues on the imaginary axis, which may cancel poles or zeros on the imaginary axis and lead to instability. Typically, the objective function is minimized with a local optimization technique such as the bi-section method, which generally performs well since the objective function is convex (ref. 57). The resulting ᠴ controller may then be represented ⦘ by equation (20) (ref. 59), ⡹⡰ 〡 ᡶ 㕈䙢(ᡲ) = (ᠧ + ‥ ᠨ ᠨ ᡅ ㎘ ᠨ ᠷ ㎘ ᠷ ᠩ )ᡶ 㕈(ᡲ) + ᠷ ᡷ(ᡲ) ⡩ ⡩ ⦘〰 ⡰ 〰 〲 ⡰ 〲 (20) ᡳ(ᡲ) = ㎘ᠷ ᡶ 㕈(ᡲ) 〰 ぁ×ぁ where, ᡅ ᒈ ጷ is the solution to the control algebraic Ricatti equation, ᠷ is the filter gain ⦘〰 〲 ぁ ×⡩ ㊓ matrix, and ᡶ 㕈(t) ᒈ ጷ is the estimated state vector of the reduced-order plant. For more information on the ᠴ problem formulation and solution, refer to references 57 and 59.

⦘ For control design the two reduced-order plants were subjected to the same disturbance input matrix, ᠨ . The matrix was formed from a gust model and process noise weights. The gust ⡩ model representing the interaction of a wind gust with the modal velocities was derived from a sinusoidal gust column (ref. 64). The sinusoidal gust column is initialized from 0 m/s wind speed at the leading edge of the wing model building in strength toward the trailing edge of the wing.

The derived gust basis represents the physical effect of the gust in modal space and is used in the weighting scheme as well as for simulation later on. The first column of disturbance weighting matrix ᠨ was represented as the velocity basis of the gust weighted with a sustained ⡩ gust velocity ᡵ (ᡲ) of 5 m/s. The second column of ᠨ was represented by the acceleration 〴 ⡩ basis of the gust weighted with a sustained gust acceleration ᡵ 䙢 (ᡲ) of 9.81 m/s . The values of 〴 the gust weights were selected arbitrarily, but it is speculated that the values represent what a passenger aircraft could experience when passing through rough air. The weighting matrix ᠨ ⡩ was also augmented with unit vectors characterizing the presence of process noise.

For control design on both reduced-order plants, the goal state matrix ᠩ was modified so ⡩ that the first two modes would receive highest weights. The first two modes were weighted highest since figure 3 indicates that these two modes will have a flutter interaction. By giving the first two modes higher weightings, the optimization technique emphasizes the reduction of the peaks of the first two modes due to the gust disturbance across the frequency range of ᠳ (ᡱ) ぐこ as much as possible. The control-regulated goal matrix ᠰ was given equal weightings which ⡩⡰ penalized high control surface movement.

The measurement matrix ᠩ was set for each controller respectively to either equation (11) ⡰ for the controller with modal coordinate input or equation (13) for the controller with accelerometer input. The controller’s sensitivity to measurement noise was also reduced by modifying the weighting matrix ᠰ with moderately high weights corresponding to each sensor.

⡰⡩ Robust Stability The ᠴ synthesis does not always result in a controller which meets performance ⦘ specifications; thus the controller design process is an iterative one. To bound the design process, objectives were defined for both control designs. Assuming an initially-stable control design, two objectives for the controllers were that each should be robust to structured uncertainty and also have good disturbance rejection properties. More specifically, the first goal for both controllers was that they be robust to at least 5 percent multiplicative uncertainty on the inputs or outputs of the plant at a low frequency of 1 Hz. They should also be robust to at least 25 percent multiplicative uncertainty near the higher flutter frequency of 4.49 Hz. Secondly, each controller was designed to mitigate the modeled sinusoidal gust disturbance described above.

A quantitative measure of the robust stability margin of both controllers is the structured singular values (SSV) or μ. The SSV is defined as shown in equation (21) (ref. 57), †(ᡁ) ᕤ (21) ᡥᡡᡦ{ᡣ | ᡖᡗᡲ(ᠵ ㎘ ᡣ ᡁ[) = 0 for structured [, … ㍤([) 㐉 1} ぀ ぀ where ᡣ is the stability margin defined as ᡣ = 1/†(ᡁ) , ᡁ is the lower linear fractional ぀ ぀ transformation (LFT) of the generalized plant ᡂ , [ represents a structured uncertainty block, and the maximum singular value of [ is defined by … ㍤([) . The value of † = 1 occurs when there is a perturbation with … ㍤([) 㐉 1 , which is just large enough to make ᠵ ㎘ ᡁ[ singular. A larger value of † is undesirable, as it means that a smaller perturbation makes ᠵ ㎘ ᡁ[ singular. The “generalized small gain theorem” states the robust stability (RS) condition, as shown in equation (22) (ref. 57).

ᡄᡅ 㙃 †㐵ᡁ(⠵2″)㐹 < 1 ᡓᡦᡖ … ㍤㐵[(⠵2″)㐹 < 1, ᒂ″ (22) To verify RS for both plants, both reduced-order plants were subjected to multiplicative input or output uncertainty. The class of all generalized plants for multiplicative input uncertainty is given as shown in equation (23) (ref. 57).

ᠳ (ᡱ) = ᠳ(ᠵ + ᡉ [ ) (23) 【 【 【 The class of all generalized plants for multiplicative output uncertainty is given as shown in equation (24) (ref. 57).

ᠳ (ᡱ) = (ᠵ + ᡉ [ )ᠳ (24) 〖 〖 〖 The structure of the input and output weights is defined to be diagonal; that is, ᡉ ᕤ ᡖᡡᡓᡙ{ᡵ (ᡱ), ᡡ = 1 … ᡔ} . For ᡉ , b is a scalar equal to the dimension of the inputs.

〶 【 For ᡉ , ᡔ is set equal to the scalar dimension of the outputs. The uncertainty block [ for both 〖 input and outputs is defined to be diagonal, [ ᕤ ᡖᡡᡓᡙ{‒ , ᡡ = 1 … ᡔ} , where ᡔ is set accordingly 〶 for inputs and outputs. The uncertainties on the diagonal ‒ are also defined to be < 1 , so that 〶 the RS condition shown in equation (22) always holds. Generally, uncertainty is greater at higher frequencies, so the uncertainty is made to vary with frequency by weight functions, as shown in equation (25) (ref. 57).

ᡱ + ᡰ あ ᡵ (ᡱ) = 〶 (25) 䙲 䙳 ᡱ + 1 ᡰ ⦘ The constant ᡰ is the relative uncertainty magnitude at steady state; and 1/ is approximately あ the frequency where the relative uncertainty reaches 100 percent. The constant ᡰ is the ⦘ magnitude of the weight at higher frequencies. The constants were selected so that the uncertainty would be greater than 5 percent at 1 Hz and greater than 25 percent near the flutter frequency of 4.49 Hz.

The plant with input multiplicative uncertainty is modeled as that given in figure 11.

Figure 11. Generalized plant with structured input uncertainty.

〡 〡 For the inputs [ᡷ , ᡳ] and the outputs [ᡳ , ᡴ] , the generalized plant with multiplicative input ⥀ ⥀ uncertainty defined above may be shown to be as presented in equation (26).

0 ᡉ 【 ᡂ = 䙴 䙵 (26) 【 ㎘ᠳ ㎘ᠳ The ᡁ structure may be formed from a lower LFT of ᡂ and ᠷ represented as ᠲ (ᡂ , ᠷ) . By 【 〹 【 ⡹⡩ carrying out the matrix operations it can be shown that ᡁ = ㎘ᡉ ᠷᠳ(ᠵ + ᠷᠳ) = ㎘ᡉ ᡆ , 【 【 【 where ᡆ is the input complementary function. The †(ᡁ) for the ᠴ controllers with multiplicative 【 ⦘ control input uncertainty is given in figure 12.

Figure 12. Structured singular-value vector analysis with input multiplicative uncertainty.

The stability boundary is defined to be 1 for † as is a requirement for RS from equation (22).

The maximum † for the controller with modal coordinate inputs is approximately 0.36 at 2.14 Hz. The maximum † for the controller with accelerometer inputs was approximately 0.35 at 2.23 Hz. The † for both controllers was bell-shaped across the frequency range. The ᠴ controller with accelerometer inputs resulted in a closed-loop system which was slightly ⦘ more robust to input uncertainty. The difference is not substantial, however. Both controllers seemed to experience a peak in † near 2.2 Hz. One might expect the peak to occur at the open-loop flutter frequency of 4.49 Hz; however, the frequencies correspond to the closed-loop pole locations. From this analysis, it was determined that both controllers meet and exceed expectations with respect to input multiplicative uncertainty.

Output uncertainty is also of interest, since different measurement systems are being utilized.

The plant with output multiplicative uncertainty is modeled as that given in figure 13.

Figure 13. Generalized plant with multiplicative output uncertainty.

〡 〡 For inputs, [ᡷ [ᡳ , ᡳ] and outputs , ᡴ] , the generalized plant with multiplicative output ⥀ ⥀ uncertainty defined above may be shown to be as presented in equation (27).

0 ᡉ ᠳ 〖 ᡂ = 䙴 䙵 (27) 〖 ㎘ᠵ ㎘ᠳ The ᡁ structure may be formed from a lower LFT of ᡂ and ᠷ represented by ᠲ (ᡂ , ᠷ) . By 〖 〹 〖 ⡹⡩ (ᠵ + ᠳᠷ) carrying out the matrix operations, it can be shown that ᡁ = ㎘ᡉ ᠳᠷ = ㎘ᡉ ᡆ , where 〖 〖 ᡆ is the output complementary sensitivity function. The †(ᡁ) calculated for both ᠴ controllers ⦘ with multiplicative measurement output uncertainty is given in figure 14.

Figure 14. Structured singular-value vector analysis with output multiplicative uncertainty.

The maximum † for the controller with modal coordinate inputs is approximately 0.47 at 2.10 Hz. The † for the controller with modal coordinate inputs experienced a second peak near a frequency of 3.38 Hz. The maximum † for the controller with accelerometer inputs was approximately 0.46 at 2.08 Hz. The † for the controller with accelerometer inputs was nearly bell-shaped across most of the frequency range, and descended until a frequency of 3.54 Hz, at which point it climbed for a short time. Overall, the characteristics of the † for both controllers indicate that the ᠴ controller with accelerometers resulted in a closed-loop system which was ⦘ slightly more robust to output uncertainty. The relative stability margin difference between the two controllers, however, is negligible. Both controllers meet and exceed expectations with respect to output multiplicative uncertainty.

Simulation

The time simulation of controllers is a reliable way to diagnose performance and make comparisons. Several case studies were selected to be performed with different objectives. The first two case studies pertain to GLA and AFS. The next two case studies focused on the use of modal filtering for virtual deformation estimation and tracking on the wing. The original plants were upgraded with “integral of modal position” states and reduced. New controllers were then derived using the same methodology presented for the regulators. The first controller tracks virtual deformations from a modal command; the second controller tracks the same virtual deformations with a virtual deformation modal approximation.

Gust Disturbance The regulators described above were tested in a simulation environment, with a gust disturbance input. The simulation structure that was used to model the gust disturbance for the controller with modal coordinate inputs is given in figure 15.

Figure 15. Control simulation with modal filter for gust modeling.

The simulation structure for the accelerometer inputs is the same except that the connection after the gust inputs is input directly into the summing block with the reference input. The 〡 exogenous inputs to the system are [ᡦ , ᡵ, ᡦ ] corresponding to input multiplicative noise, gust え 〱 disturbance states, and additive measurement noise. Zero mean multiplicative Gaussian noise, ᡦ , with a standard deviation of 0.1, is modeled on each control input for both control systems, え making ᠹ in figure 15 the identity matrix of size ᡕ .

え The standard 1 ㎘ cos(ᡶ) gust profile (see ref. 64) is modeled to characterize the transient 〡 shape of the gust disturbance ᡵ = 㐙ᡵ (ᡲ) ᡵ 䙢 (ᡲ)㐣 which lasts for 1.6 s. The gust velocity and 〴 〴 acceleration at time 0 are both initialized to zero. The gust is shaped to achieve a maximum velocity, ᡵ (ᡲ) of 5 m/s and a maximum acceleration ᡵ 䙢 (ᡲ) of 9.81 m/s . The signals, ᡵ (ᡲ) and 〴 〴 〴 ᡵ 䙢 (ᡲ) are input to the LTI gust system ᠳ which acts as a disturbance on the output of the 〴 〴えうぇ plant ᠳ .

For the controller with modal coordinate inputs, zero mean Gaussian noise ᡦ with a 〱 standard deviation of 1 cm was added to the measurement signals. Deformations are used instead of strain because strain was not available directly in the model. When accelerometers were used in place of the deformation measurements, it was assumed that the additive noise had a standard deviation of 1.0 m/s . For the controller with modal coordinate inputs, the noise was added to the simulated deformations, shown in figure 8. Accelerometer measurements tend to be somewhat noisy, whereas fiber optic measurement systems are expected to produce measurements with a very high signal-to-noise ratio (ref. 65).

This model being a simulation model, true displacement measurements were not available.

The controller with modal coordinate inputs makes use of a deformation simulation, by multiplying the modal matrix f indexed at measurement index stations ᡤ with the true modal coordinates ᡩ . The modeled displacement information and additive displacement noise is input into the least-squares modal filter introduced in equation (3). The estimated modal coordinates are then indexed ( ᠹ = {1,2}) to obtain ᡩ . This signal is then used to form the control 〶ぁ〱〲け ぀ signal input to the ᠴ controller ᠷ . Since gust alleviation is the objective of the controller, the ⦘ reference on each measurement is set to 0. The gust simulation results are presented in figure 16 for the controllers with acceleration inputs and modal coordinate inputs respectively.

Figure 16. Controller performance in a gust: a) controller with accelerometer inputs, modal amplitude time history; b) controller with accelerometer inputs, control surface time history; c) controller with modal coordinate inputs, modal amplitude time history; and d) controller with modal coordinate inputs, control surface time history.

The modal response, figure 16(a), to the gust for the controller with accelerometer inputs indicates that the first two modes (first wing torsion and first wing bending) responded most to the gust as expected. The peak amplitude of the bending mode was approximately -2.2, which corresponds to roughly a 44-cm-upward bending deflection at the wing tip. A negative bending modal coordinate corresponds to a positive wing tip deflection, figure 2(a). The twist mode moves to a maximum amplitude of 0.4, which corresponds to approximately 10 deg of positive wing twist, leading-edge up. Over the gust time history, the control surfaces move to counteract the effect of the gust, figure 16(b). As the wing experiences a lift increase, the control surfaces rotate upward to reduce the angle of attack of the wing and reduce lift. The rotation of the control surfaces stayed well within the bounds of reason for wing control surface rotations.

The gust disturbance rejection performance of the controller with modal coordinate inputs was comparable to that of the controller with accelerometer inputs, shown in figure 16(c); the peak amplitude of the bending mode was slightly higher at -2.4. The torsion angle was nearly the same at approximately 10 deg. Little can be said as to which controller has better disturbance rejection. The differences were negligible. The controllers rejected the specified gust disturbance adequately.

Flutter Suppression A major theme in this study has been to demonstrate that the controller with the modal coordinate inputs may be used for flutter suppression. Figure 5 shows that the model is open-loop unstable, resulting in modal motions characteristic of flutter. The same simulation structure (see fig. 15) that was used for the gust disturbance modeling was used for flutter suppression with very small changes. The gust model inputs ᡵ were set to zero. A small control input at time 0 was introduced to perturb the wing model from its trim state. The modal amplitudes of the model are allowed to increase without control input until 3.5 s. The controller was linearly phased-in from 3.5 s to 4.5 s. The controller was not turned on to full instantly at 3.5 s to avoid large oscillations due to the output magnitudes being far from the reference condition of zero. The simulation results are presented in figure 17 for the controllers with acceleration inputs and modal coordinate inputs respectively.

Figure 17. Controller performance in a flutter suppression: a) controller with accelerometer inputs, modal amplitude time history; b) controller with accelerometer inputs, control surface time history; c) controller with modal coordinate inputs, modal amplitude time history; and d) controller with modal coordinate inputs, control surface time history.

The time history for the controller with accelerometer inputs is examined first in figure 17(a).

The modal amplitudes oscillate with a frequency of 4.49 Hz increasing in amplitude until approximately 4.2 s. At this time, the controller force begins to remove a sufficient amount of energy from the flutter mode to begin to reduce the amplitudes of the modes. The oscillations die out quickly at approximately 5 s.

The time history of the controller matches what would physically be required to reduce flutter in the wing. As the wing bends upward, figure 17(a), the control surfaces rotate upward, figure 17(b), to reduce the angle of attack of the wing and reduce lift. The net aerodynamic force has the effect of moving the wing downward. As the wing moves down, the control surfaces rotate downward to increase the lift on the wing. This counterbalancing effect performs work and removes energy from the flutter mode. The overall effect asymptotically stabilizes the structure.

After the flutter mode stabilizes, the control surface movements appear to oscillate at low frequency and the modal coordinates remain near zero. The movement from equilibrium is in response to the additive noise on the accelerometers.

The controller with modal coordinate inputs performed similarly to the controller with accelerometer inputs, as before with the gust inputs. Figure 17(c) shows that the modal coordinates begin to flutter up to 3.5 s and are slowly damped out once the controller is enabled. As before, the control surface movements worked to extract energy from the flutter mode as seen in figure 17(d). The modal coordinates also stay near zero as expected. From time analyses it was clear that both controllers performed well in meeting the primary objective to suppress the flutter mode at the selected flight condition.

So far, the actual modal coordinate time histories are presented in each plot, representing the motion of the model in either a gust or fluttering condition. To satisfy curiosity, the modal coordinates which were given to the controller for feedback during the flutter suppression ᡩ are ぀ presented in figure 18. Recall that zero mean Gaussian noise with a standard deviation of 1 cm was added to all deflections that were used to estimate the modal coordinates. This means that the error can likely go up approximately 3 cm for the deformations some of the time, assuming that some data points will fall roughly 3 standard deviations away. All of the deflections with error were put through a least-squares modal filter [see fig. 15 and equation (3)]. The measurement error of the modal coordinates indicates that a typical least-squares smoothing has taken place, as the modal amplitudes measurement error tended to stay near 0.01.

Although this smoothing is not substantial enough to raise eyebrows, it does show that the errors tend to average out when many sensors are utilized in forming the least-squares estimates. It is also clear that the increasing amplitude of the modal coordinates did not affect the modal coordinate estimation error, even during open-loop flutter.

Figure 18. Modal coordinate measurement error.

Virtual Deformation Control (Case Study: Modal Tracking) Until this point, much of the focus of this study has been to compare regulators with accelerometer inputs or modal coordinate inputs. It was observed that the controllers had similar performance during GLA and AFS. In this section, shape control is demonstrated. This type of research has been pursued in other works as well; see references 11 and 66.

Shape control on the aeroelastic model is implemented here making use of the modal filter.

Further work with accelerometers is not continued. While the controller with accelerometers might be able to track deformations by double-integrating the accelerations, it is not a natural fit; the deformations predicted with accelerometers may start to drift and require deformation updates.

The modal controller may be a suitable match for shape control since modes are linearly-related to deflections [see eq. (1)]. If a set of reference deflections ᡖ are known at ぅ〲〳 specific locations ᡤ , perhaps from a flight computer, these deformations may be transformed to ぅ modal coordinate reference values, ᡩ , indicating that the reference deformations need not be ぅ〲〳 at the points of measurement on the model. These deflections can be referred to as virtual deformations, and the control of these deformations as virtual deformation tracking. The simulation structure for deformation control through modal reference tracking is presented in figure 19.

Figure 19. Control simulation for virtual deformation control: modal control.

To achieve zero steady-state tracking error, “integrator of modal position” states are added to the controller. All modes are estimated by the modal filter and indexed (that is, ᠹ = {1,2} ) 〶ぁ〱〲け at the first two locations to give, ᡩ . The modal measurements are input to a differencing ぀ junction with the modal reference and sent through a single continuous-time integrator.

〡 To simulate a small torsion angle command, ᡖ = {ᡆᠱ , ᠸᠱ } was set to {1 ᡥᡥ, ㎘1 ᡥᡥ} , ぅ〲〳 ぅ〲〳 ぅ〲〳 which represents a leading-edge-down rotation. Since small deformations are used as references, the noise was adjusted accordingly, so that the standard deviation of the multiplicative control noise was set to 0.001. The standard deviation on the deformation measurements was assumed to be 1 mm. These noise settings allows one to better see what the controller is doing in the resulting plots; this was not done to simulate actual sensor noise characteristics. The simulation results are presented in figures 20(a) through 20(c).

Figure 20. Virtual deformation control: modal commands: a) modal amplitudes; b) deformations and predicted deformations at wing tip; c) control surface rotations.

Figure 20(a) shows that the first torsion modal coordinate overlays the reference torsion modal coordinate within 50 s. The bending modal coordinate reference is near zero and is also tracked within 50 s. The other eight modal coordinate time histories are also plotted so that the effects of residual modes may be observed. The second bending mode becomes highly excited.

The deformations achieved through modal tracking are presented in figure 20(b). The deflections achieve what would be a torsional angle with the leading edge down, but the net deformations of both are up approximately 2 mm. There is significant error between the desired deformations of {1 mm, -1 mm} and what is achieved {3.6 mm, 2.2 mm}. The prominent second bending modal coordinate is clearly to blame for this error. From desktop simulations it was observed that to reduce this error, the second bending mode, figure 2(c), could be tracked if more actuators were available. Note that the deformations to be tracked were at the same locations as where the accelerometers were placed. They were not measured by the fiber optics themselves. The action of tracking a deformation which is not directly measured is observed here.

Figure 20(c) indicates that to achieve small deformations, very large control surface rotations were required, almost up to 40 deg. The large rotations are a result of either high stiffness in the wing or potentially low control surface steady-state effectiveness. For this reason, the deflection references were kept small to ensure the control surfaces rotated within reasonable limits.

These results indicate that wings with low torsional modal mass (see fig. 7) may be difficult to structurally morph using aerodynamic effectors.

Virtual Deformation Control (Case Study: Predicted Deformation Tracking) In the previous case study, the deformation command was transformed to a reduced modal command and the modes were tracked. Due to the effect of residual modes, the wing deformation reference command was tracked poorly. To reduce this effect, the modal filter can also be used to form a predicted estimate of the deformation of the structure at any point by including residual modes into the estimate. To prepare the controller, the output matrix in the state space may be defined to have the form shown in equation (28), , : ) 0 0 0 ]ᡶ ᡷ ᕤ [f(ᡤ (28) ぅ ぅ×぀ ぅ×〲 ぅ×こ where ᡰ is the number of deformations (or virtual deformations) desired to be tracked. This method of definition has the effect of making the outputs of the plant equivalent to the deformations. A similar transformation is used to model the accelerometers, where the sensors are assumed to measure linear combinations of the modal states [see eqs. (12) and (13)].

The simulation scheme used for virtual deformation tracking is presented in figure 21.

The simulation for tracking deformations estimated by modal coordinates is similar to the modal tracking simulation presented above (see fig. 19). The difference is that after the modal coordinates are estimated and indexed, the modal matrix is again multiplied by the modal matrix indexed at desired virtual measurement locations, ᡤ . Since all modes are utilized to ぅ estimate the virtual deflections, ᡩ = ᡩ . Alternatively, it might be stated that ぀ 〲うぇ ᠹ = {1,2,3,4,5,6,7,8,9,10} . The noise levels in the previous simulation are also used, as 〶ぁ〱〲け were used previously for modal tracking.

Figure 21. Control simulation for virtual deformation tracking: predicted deformation control.

The same reference values of {1mm,-1mm} were used for the deformation reference and the simulation results are presented in figures 22(a) through 22(c).

Figure 22. Virtual deformation control: predicted deformation: a) modal amplitudes; b) deformations and predicted deformations at wing tip; c) control surface rotations.

The time history of the modal coordinates in figure 22(a) shows that the torsion modal coordinate moved comparative to the way it did previously, seen in figure 20(a). The second wing bending modal coordinate also moved positively; however, this time, the first wing bending moves from zero to a large positive value, which has the effect of offsetting the second wing-bending effects. The net effect of this was that the actual virtual deformations and reference virtual deformations were overlaid, seen in figure 22(b).

Figure 22(c) indicates that the control surface movements were lower than they were previously, figure 20(c). The outboard control surface moved to approximately 30 deg and the inboard surface moved to -25 deg. The use of either the tracking strategy presented here or that presented above may depend on the application. If a multitude of points on the wing are required to be tracked or moved to a particular shape, then the strategy first presented may be more useful. The reduction of the reference signals to a few modal coordinates may alleviate the control design effort. If only a few virtual deflections at a few points on the structure (that is, 2-4) are required to be tracked, then the strategy presented here may be more applicable. The selection of the appropriate tracking strategy will be application-dependent.

Conclusions and Future Work

The findings described herein suggest that the modal filter may be suitable for integration into flexible aircraft control systems. A comparison of regulators with modal filters and regulators with accelerometers indicated that the stability and performance of the regulators during gust load alleviation and active flutter suppression changes very little with use of the modal filter.

During the control design, an unexpected finding was that when using balanced truncation, the order of the plant with modal filter outputs, could be reduced substantially. A lower order plant (and thus controller) may facilitate the development of a digital controller.

Modal filters were shown to support the capability for wing-shape control. It was demonstrated in two case studies, using two different architectures, how deformations (a few or many) at any point on the wing could be controlled. The choice of architectures, it is believed, will be application-dependent. The virtual deformation control architecture with modal control will be the most generally applicable and further research into limiting the impact of residual modes is needed. During the virtual deformation control simulation studies it was discovered that high control rotations may be required if wing torsional modal mass is very low.

The model used for the study had a few shortcomings which limited the analysis. The first shortcoming was that the structure was notional and not based on an operating or experimental aircraft. The second shortcoming was that the deformation measurements were essentially perfect before noise was added. The accuracy of the modal estimation may also be called into question. All of these modes in the state space simulation model were incorporated into the basis of the modal filter, so the bias of residual modes was not present in the modal coordinate estimates. Additionally, aircraft require guidance systems and rigid-body mode control; the interacting effects of the two control systems with the modal controller were not seen here.

Further work will likely be to revisit the implementation of the modal filter and potentially its variations. A deformation tracking scheme may also be incorporated onto a model that is more representative of an aircraft.

Acknowledgments

This research was approved by Steve Jacobson, Branch Chief of the controls branch at the Dryden Flight Research Center, Edwards, California, and funded by the National Aeronautics and Space Administration Aeronautics Research Mission Directorate Fixed Wing program.

Steve Jacobson; John Bosworth, the Dryden X-56A chief engineer; Chris Regan, the former Dryden X-56A controls lead; John Ryan, Dryden X-56A controls lead; and Dryden aeroelasticity expert Marty Brenner; have been excellent advisors in the development of the results of this research.

References

1. O'Connell, R. F., and A.F. Messina, “Development of an Active Flutter Margin Augmentation System for a Commercial Transport,” Journal of Guidance and Control , Vol. 3, No. 4, 1979, pp. 352-360.

2. Bernelli-Zazzera, Franco, Paolo Mantegazza, Giovanni Mazzoni, and Matteo Rendina, “Active Flutter Suppression Using Recurrent Neural Networks,” Journal of Guidance, Control, and Dynamics , Vol. 23, No. 6, 2000, pp. 1030-1036.

3. North Atlantic Treaty Organization Science and Technology Organization, Structural Aspects of Flexible Aircraft Control, RTO-MP-36, NATO, 2000.

4. Levin, Jason, Néstor O. Pérez-Arancibia, and Petros A. Ioannou, “Adaptive Notch Filter using Real-Time Parameter Estimation,” IEEE Transactions on Control Systems Technology , Vol. 19, No. 3, 2011, pp. 673-681.

5. Zeng, Jie, Jiang Wang, Raymond de Callafon, and Martin Brenner, “Suppression of the Aeroelastic/Aeroservoelastic Interaction Using Adaptive Feedback Control Instead of Notching Filters,” AIAA-2011-6459, 2011.

6. Nicolai, Leland, Keith Hunten, Scott Zink, and Pete Flick, “System Benefits of Active Flutter Suppression for a SensorCraft-Type Vehicle,” AIAA-2010-9349, 2010.

7. Beranek, Jeff, Lee Nicolai, Mike Buonanno, Edward Burnett, Christopher Atkinson, Brian Holm-Hansen, et al., “Conceptual Design of a Multi-utility Aeroelastic Demonstrator,” AIAA-2010-9350, 2010.

8. Kusmina, S., F. Ishmuratov, and V. Kuzmin, “Minimization of induced drag of elastic airplane,” AIAA-2004-4611, 2004.

9. Heinze, Sebastian, Aeroelastic Concepts for Flexible Aircraft Structures, PhD Thesis; Royal Institute of Technology, Stockholm, Sweden, 2007.

10. Weisshaar, Terrence A., and David K. Duke, “Induced Drag Reduction Using Aeroelastic Tailoring with Adaptive Control Surfaces,” AIAA-2000-1619, 2006.

11. Austin, Fred, Michael J. Rossi, William Van Nostrand, Gareth Knowles, and Antony Jameson, “Static Shape Control for Adaptive Wings,” AIAA Journal , Vol. 32, No. 9, 1994, pp.

1895-1901.

12. Ehlers, S. M., and T. A. Weisshaar, “Static Aeroelastic Control of an Adaptive Lifting Surface,” AIAA-90-1078, 1993.

13. Norris, Guy, “Breaking the Flutter Barrier,” Aviation Week & Space Technology , February 6, 2012, pp. 28-29.

14. Su, Weihua, and Carlos E. S. Cesnik, “Nonlinear Aeroelasticity of a Very Flexible Blended-Wing-Body Aircraft,” AIAA-2009-2402, 2009.

15. Richards, Lance, Allen R. Parker, William L. Ko, and Anthony Piazza, “Real-time In-Flight Strain and Deflection Monitoring with Fiber Optic Sensors,” Space Sensors and Measurements Techniques Workshop , Nashville, Tennessee, 2008.

16. Preumont, A., A. Francois, P. De Man, and V. Piefort, “Spatial filters in structural control,” Journal of Sound and Vibration , Vol. 265, 2003, pp. 61-79.

17. Meirovitch, L., and H. Baruh, “The Implementation of Modal Filters for Control of Structures,” Journal of Guidance , Vol. 8, No. 6, 1985, pp. 707-716.

18. Hwang, Jai-Hyuk, Jung-Soo Kim, and Joon-Soo Kim, “Effect of Modal Filter Errors on Vibration Control Characteristics,” KSME International Journal , Vol. 12, No. 5, 1998, pp. 836-846.

19. Liebst, Bradley S., “Accelerometer Placement in Active Flutter Suppression Systems,” Journal of Guidance , Vol. 10, No. 5, 1987, pp. 441-446.

20. Mukhopadhyay, Vivek, “Flutter Suppression Control Law Design and Testing for the Active Flexible Wing,” AIAA-92-2095, 1995.

21. Pisoni, Attilio C., Claudio Santolini, Dagmar E. Hauf, and Steven Dubowsky, “Displacements in a Vibrating Body by Strain Gauge Measurements,” Proceedings of the 13th International Conference on Modal Analysis , Nashville, Tennessee, 1995.

22. Kang, Lae-Hyong, Dae-Kwan Kim, and Jae-Hung Han, “Estimation of dynamic structural displacements using fiber bragg grating strains sensors,” Journal of Sound and Vibration , Vol. 305, No. 3, 2007, pp. 534-542.

23. Nissim, E., Flutter Suppression Using Active Controls Based on the Concept of Aerodynamic Energy , NASA TN D-6199, 1971.

24. Abel, Irving, Jerry R. Newsom, and Henry J. Dunn, “Application of Two Synthesis Methods for Active Flutter Suppression of an Aeroelastic Wind Tunnel Model,” AIAA 79-1633, 1979.

25. Mukhopadhyay, V., J. R. Newsom, and I. Abel, “Reduced-Order Optimal Feedback Control Law Synthesis for Flutter Suppression,” AIAA-80-1613, 1982.

26. Waszak, Martin R., “Robust Multivariable Flutter Suppression for Benchmark Active Control Technology Wind-Tunnel Model,” Journal of Guidance, Control and Dynamics , Vol. 24, No.

1, 2001, pp. 147-153.

27. Garrard, William L., and Bradley S. Liebst, “Active Flutter Suppression Using Eigenspace and Linear Quadratic Design Techniques,” Journal of Guidance , Vol. 8, No. 3, 1985, pp. 304-311.

28. Mukhopadhyay, Vivek “Transonic Flutter Suppression Control Law Design Using Classical and Optimal Techniques with Wind-Tunnel Results,” AIAA-99-1396, 1999.

29. Haley, Pam, and Don Soloway, “Generalized Predictive Control for Active Flutter Suppression,” Journal of Guidance, Control, and Dynamics , Vol. 24, No. 1, 2001, pp. 154-159.

30. Kim, Do-Hyung, Jae-Hung Han, and In Lee, “Application of Fiber Optic Sensor and Piezoelectric Actuator to Flutter Suppression,” Journal of Aircraft , Vol. 41, No. 2, 2004, pp. 409-411.

31. Abel, Irving, and Maynard C. Sandford, Status of Two Studies on Active Control of Aeroelastic Response , NASA TM X-2909, 1973.

32. Roger, Kenneth L., Garold E. Hodges, and Larry Felt, “Active Flutter Suppression - A Flight Test Demonstration,” AIAA-74-402, 1974.

33. Sensburg, O., H. Honlinger, T. E. Noll, and L. J. Huttsell, “Active Flutter Suppression on an F-4F Aircraft,” AIAA-80-0770, 1982.

34. Peloubet Jr., R. P., R. L. Haller, and R. M. Bolding, “F-16 Flutter Suppression System Investigation Feasibility Study and Wind Tunnel Tests,” AIAA-80-0786, 1982.

35. Shelley, S. J., K. L. Lee, T. Aksel, and A. E. Aktan, “Active-Control and Forced-Vibration Studies on Highway Bridge,” Journal of Structural Engineering , Vol. 121, No. 9, 1995, pp. 1306-1312.

36. Shelley, S., L. Freudinger, R. J. Allemang, and Q. Zhang, “Implementation of a Modal Filter on a Five Meter Truss Structure,” Proceedings of the 9th International Modal Analysis Conference , Schenectady, New York, 1991, pp. 1036-1044.

37. Luenberger, David G., “Observing the State of a Linear System,” IEEE Transactions on Military Electronics , Vol. 8, No. 4, 1964, pp. 74-80.

38. Kalman, R. E., and R. S. Bucy, “New Results in Linear Filtering and Prediction Theory,” Journal of Basic Engineering , March, 1961, pp. 95-108.

39. Preumont, Andre, Vibration Control of Active Structures An Introduction 2nd Edition, Kluwer Academic Publishers, Dodrecht, The Netherlands, 2002.

40. Meirovitch, L. Analytical Methods in Vibrations, The Macmillan Company, New York, 1967.

41. Shelley, S. J., R. J. Allemang, G. L. Slater, and J. F. Schultze, “Active Vibration Control Utilizing an Adaptive Modal Filter Based Modal Control Method,” 11th International Modal Analysis Conference , Kissimmee, Florida, 1993.

42. Lee, C.-K., and F. C. Moon, “Modal Sensors/Actuators,” Journal of Applied Mechanics , Vol. 57, No. 2, 1990, pp. 434-441.

43. Sumali, Hartono, Karsten Meissner, and Harley H. Cudney, “A piezoelectric array for sensing vibration modal coordinates,” Sensors and Actuators A: Physical , Vol. 93, No. 2, 2001, pp. 123-131.

44. Foss, G. C., and E. D. Haugse, “Using Modal Test Results to Develop Strain to Displacement Transformations,” 13th International Modal Analysis Conference , Nashville, Tennessee, 1995, pp. 112-118.

45. Treiber, Johannes, Uwe C. Mueller, Jae-Hung Han, and Horst Baier, “Filtering techniques in the dynamic deformation estimation using multiple strains measured by FBGs,” Proceedings of SPIE, Vol. 6932 69322A-1, 2008.

46. Tessler, Alexander, and Jan L. Spangler, “A least-squares variational method for full-field reconstruction of elastic deformations in shear-deformable plates and shells,” Computer Methods in Applied Mechanics and Engineering , Vol. 194, 2005, pp. 327-339.

47. Ko, William L., and William Lance Richards, Method for real-time structure shape-sensing, Patent No. 7,520,176 B1, April 21, 2009.

48. Burner, Alpheus W., William A. Lokos, and Danny A. Barrows, “In-flight aeroelastic measurement technique development,” Proc. SPIE 5191 Optical Diagnostics for Fluids, Solids, and Combustion II, 186 , 2003.

49. Szilard, Rudolph, Theories and Applications of Plate Analysis: Classical, Numerical and Engineering Methods , John Wiley & Sons, New Jersey, 2004.

50. Pototzky, Anthony S., “Enhanced Modeling of First-Order Plant Equations of Motion for Aeroelastic and Aeroservoelastic Applications,” AIAA-2010-7801, 2010.

51. Conyers, Howard J., Earl H. Dowell, and Kenneth C. Hall, “Aeroelastic Studies of a Rectangular Wing with a Hole: Correlation of Theory and Experiment,” 2010 Aerospace Systems Conference, February 5-10, 2010, Los Angeles, California, 2010.

52. Wright, Jan R., and Jonathan E. Cooper, Introduction to Aircraft Aeroelasticity and Loads, John Wiley & Sons, Ltd, West Sussex, England, 2007.

53. Gupta, K. K., M. J. Brenner, and L. S. Voelker, Development of an Integrated Aeroservoelastic Analysis Program and Correlation with Test Data , NASA TP-3120, 1991.

54. Pak, Chan-gi, Unsteady Aerodynamic Model Tuning for Precise Flutter Prediction, NASA/TM-2011-215963, 2011.

55. Kammer, D. C., “A Hybrid Approach to Test-Analysis-Model Development for Large Space Structures,” Journal of Vibration and Acoustics , Vol. 113, 1991, pp. 325-332.

56. Zeng, Jie, Boris Moulin, Raymond de Callafon, and Martin J. Brenner, “Adaptive Feedforward Control for Gust Load Alleviation,” AIAA-2008-6373, 2010.

57. Skogestad, Sigurd, and Ian Postlethwaite, Multivariable Feedback Control: Analysis and Design, John Wiley & Sons, Inc., West Sussex, England, 2005.

58. Hanis, Tomas, and Martin Hromcik, “Optimal sensors placement and spillover suppression,” Mechanical Systems and Signal Processing 28 , 2012, pp. 367-378.

59. Gawronski, Wodek, Advanced Structural Dynamics and Active Control of Structures, Springer-Verlag, New York, 2004.

60. Sjövall, Per, Thomas Abrahamsson, and Tomas McKelvey, “Indirect Vibration Sensing and Optimal Sensor Placement,” 2008 IMAC-XXVI: Conference & Exposition on Structural Dynamics , Orlando, Florida, 2008.

61. Ryall, T. G., and C. S. Fraser, “Determination of Structural Modes of Vibration Using Digital Photogrammetry,” Journal of Aircraft , Vol. 39, No. 1, 2002, pp. 114-119.

62. Vishwakarma, C. B., and R. Prasad, “Clustering Method for Reducing Order of Linear System using Padé Approximation,” IETE Journal of Research , Vol. 54, No. 5, 2008, pp. 326-330.

63. Mukhopadhyay, Vivek, Jerry R. Newsom, and Irving Abel, A Method for Obtaining Reduced- Order Control Laws for High-Order Systems Using Optimization Techniques , NASA Technical Paper 1876, 1981.

64. Kroll, R. I., and R. E. Clemmons, “ A Computer Program to Generate Equations of Motion Matrices, L217 (EOM) Volume I: Engineering and Usage,” NASA Contractor Report 2851, 1979.

65. Chuang, Kuo-Chih, Heng-Tseng Liao, and Chien-Ching Ma, “Dynamic Sensing Performance of a Point-Wise Fiber Bragg Grating Displacement Measurement System Integrated in an Active Structural Control System,” Sensors , Vol. 11, 2011, pp. 11605-11628.

66. Hu, Yan-Ru, and George Vukovich, “Active robust shape control of flexible structures,” Mechatronics , Vol. 15, 2005, pp. 807-820.