Integrated Approach to the Dynamics and Control of Maneuvering Flexible Aircraft

Contents Title Page Section 1 Introduction and Literature Review 1 2 Hybrid Equations of Motion in Terms of Quasi-Coordinates 9 3 Spatial Discretization of the Distributed Variables 4 The Aerodynamic and Gravity Forces 5 A Unified Theory for Flight Dynamics and Aeroelasticity 19 Control Design 25 7 Optimal State Observer 28 8 Numerical Example 33 i. Steady level flight ii. Level steady turn maneuver 46 9 Conclusions 58 10 References 60 Appendix Numerical Values of the Aircraft Parameters 63 iv 1. Introduction and Literature Review This investigation is concerned with a dynamic formulation capable of simulating the flight of flexible aircraft. It integrates seamlessly in a single and consistent mathematical formulation all the necessary material from the pertinent disciplines, namely, analytical dynamics, structural dy- namics, aerodynamics and controls. The unified formulation is based on fundamental principles and incorporates in a natural manner both the rigid body motions and the elastic deformations, and the couplings thereof, as well as the aerodynamic, propulsion, control and gravity forces. In essence, the formulation not only unifies flight dynamics and aeroelasticity, traditionally treated as separate disciplines, but, going beyond that, it permits a computer simulation of the response of flying flexible aircraft to external stimuli.

In describing the motion of rigid bodies in space, it is convenient to attach a set of axes to the body, where the axes are commonly known as body axes, and express the motions in terms of components along these body axes. It is quite common to describe the motion of rigid bodies in terms of the translation of the origin of the body axes and the rotation of the body axes; the corresponding variables, particularly the rotations, are referred to as quasi-coordinates. If the origin of this reference frame coincides with the mass center of the body, the translations and rotations are independent of each other. Moreover, if the body axes themselves coincide with the principal axes of the body, then the products of inertia are zero, so that the inertia matrix is diagonal.

The situation is more complicated forflexible bodies, in which case there are basically two types of reference frames: i. Fixed in the undeformed body - In this case, it is convenient to define the translation of the origin of the reference frame and the rotation of the reference frame as the rigid body translation and rotation of the body, and regard any displacement relative to the reference frame as elastic deformation.

ii. Moving relative to the undeformed body - In this case, it is common to choose the reference axes so that the linear momentum and angular momentum vectors due to elastic deformations vanish; axes satisfying these conditions are called mean axes. Because the elastic defor- mations depend in general on time, the mean axes are continuously moving relative to the undeformed body. Of particular interest are mean axes with the origin at the system mass center, because in this case the three types of motion, namely, the translations of the refer- ence frame, the rotations of the reference frame and the deformations measured relative to this frame, are all inertially decoupled. The price to be paid for the use of mean axes is very steep, however, as the constraints defining these axes are not easy to enforce. For this reason, it is common to invoke the use of mean axes to justify inertial decoupling, but without enforcing the constraints.

Whereas proper use of mean axes for flexible bodies in vacuum can produce inertial decoupling, in the case of aircraft any such benefits are questionable, because the equations of motion remain coupled through the aerodynamic forces. Moreover, if one insists on the use of mean axes, then the aerodynamic forces must also be expressed in terms of components along the same mean axes, which is a very tedious task at best.

The motion of force-fYee rigid bodies in space is unstable. Under certain circumstances, the motion of rigid bodies can be stabilized by imparting them some spin about the axis of either the minimum or the maximum moment of inertia. On the other hand, flexible bodies cannot be stabilized if the spin is about the axis of minimum moment of inertia. The preceding statements imply that the motion takes place in vacuum, such as in the case of a spacecraft. Matters are entirely different for flexible aircraft, which are neither force-free, nor do they operate in vacuum.

In fact, aircraft are subjected to aerodynamic and wind forces, and any stabilization is done by active means, namely, by the engine throttle and control surfaces, where the latter consist of the aileron, elevator and rudder.

The flight of aircraft tends to be quite diverse, ranging from steady level cruise to complex maneuvers, and stabilization requires the use of controls permitting the aircraft to carry out the intended maneuver while suppressing any deviations from it, whether in the form of rigid body displacements or elastic deformations. Both flight dynamics and aeroelasticity are concerned with aircraft stability. But, whereas flight dynamics is concerned mainly with rigid body motions, aeroe- lasticity is concerned with vibration and flutter.

A stability analysis tends to be limited in scope, in the sense that it can only yield a qualitative statement about the nature of motion in the neighborhood of an equilibrium state of a system.

More specifically, it can state whether the equilibrium is merely stable, asymptotically stable, or unstable. The time plays no role in a stability analysis. In fact, stability analyses tend to be limited to cases in which the time variable can be eliminated, such as when the equations of motion can be reduced to an eigenvalue problem. The stability of time-dependent maneuvers can only be evaluated numerically.

To obtain information going beyond stability statements, such as the time response of aircraft to external stimuli, it is necessary to undertake a simulation of the equations of motion, which amounts to the integration of the equations of motion. Such a dynamic simulation can be used for parametric studies in preliminary design. Moreover, it can be used to determine aircraft perfor- mance, thus reducing the time required for flight testing by “flying the aircraft on a computer.” The choice in this paper is to work with a reference frame attached to the undeformed aircraft, which has many advantages over mean axes. But, because the elastic deformations prevent the origin of a frame attached to the undeformed body from coinciding with the mass center and the axes themselves from coinciding with the principal axes for all times, there is no preferred choice of a reference frame; we base the choice on geometric considerations. In particular, we attach a set of body axes to the undeformed fuselage, where one of the axes is along the symmetry axis.

For convenience, sets of body axes are also attached to the other flexible components, such as the wing and the empennage. Ultimately, however, all motions are referred to the fuselage body axes, which act as a reference frame for the whole aircraft.

The mathematical formulation is based on equations of motion in terms of quasi-coordinates derived earlier by the first author for flexible spacecraft and later adapted by him to flexible aircraft.

The formulation is hybrid in nature, in the sense that it consists of ordinary differential equations for the rigid body translations and rotations of the aircraft as a whole and boundary-value problems for the elastic deformations of the flexible components of the aircraft, namely, the fuselage, wing and empennage. For practical reasons, the distributed variables describing the boundary-value problems for the individual components are discretized in space, obtaining a relatively large set of second-order ordinary differential equations for the whole aircraft. The discretization process amounts to representing the distributed variables as finite series of known space-dependent shape functions multiplied by time-dependent generalized coordinates. The derivation of the equations of motion in conjunction with the extended Hamilton principle requires expressions for the ki- netic energy, potential energy and virtual work, all scalar quantities. In turn, the kinetic energy requires the velocity of every point of the aircraft, which can be obtained by means of an orderly kinematical synthesis, going from the fuselage to the wing and to the empennage. The potential energy is merely equal to the strain energy. Moreover, the aerodynamic, propulsion, control and gravity forces are accounted for through the virtual work. Rather than deriving first hybrid equa- tions of motion and then discretizing them in space, it is perhaps more expeditious to carry out the discretization directly in the kinetic energy, potential energy and virtual work, thus obtaining the desired set of ordinary differential equations for the whole flexible aircraft without the need to de- rive boundary-value problems. For integration of the differential equations and for control design, it is necessary to transform the set of second-order differential equations into a set of first-order differential equations, namely, into a set of state equations. It turns out that, for the problem at hand, it is more convenient to work with momenta rather than with velocities. Although the result- ing first-order equations actually represent phase equations, we shall continue to refer to them as state equations.

The simulation of the flight of an aircraft amounts essentially to integration of the state equa- tions. Because of various nonlinearities involved, such as those due to rigid body motions and aerodynamic forces, the integration must be carried out numerically on a computer. In one way or another, computer integration can only be done in discrete time, which raises the question of the size of the sampling period, or time step. Of course, the answer depends on the desired accuracy of the simulation, and it is intimately related to the dynamic characteristics of the system. If the aircraft is to be controlled by an autopilot, then the simulation must be carried out in real time. If the dynamic characteristics are such that the time step must be relatively short, perhaps of the order of O.Ols, most aerodynamic theories in current use must be ruled out, as the computation of the dynamic pressure over the entire aircraft is sure to take considerably longer time than that. Hence, a new method for computing the dynamic pressure must be developed, one characterized by high computational speed, even if some accuracy must be sacrificed. Moreover, the method for com- puting the dynamic pressure must be compatible with the method for modeling the airframe. If the formulation is to be used for aircraft design, then real-time simulation may not really be necessary, although on-line simulation may. But, the size of the sampling period, which is determined by the system dynamic characteristics, remains the same regardless of whether the simulation is in real time or only on-line. The implication is that the computation of the dynamic pressure must still be relatively fast. Indeed, a mere 10 s simulation requires 1,000 time steps. Hence if the computation of the dynamic pressure takes one hour, the simulation requires over 40 days. This demonstrates the need for a method for computing the dynamic pressure in a very short time period. In this regard, a reasonably accurate approximate method may be acceptable.

As indicated above, the equations of motion for a flying flexible aircraft are nonlinear, where the nonlinearity is due to the rigid body motions and the aerodynamic forces. Moreover, the equations tend to be of high order, the order depending on the discretization procedure employed. Hence, one can expect difficulties both with a stability analysis and with control design. In addition, difficulties can be experienced in the integration process, because some of the variables describing the aircraft rigid body motions tend to be large and the variables describing the elastic displacements tend to be small. A perturbation approach to the solution can obviate many of these difficulties. More specifically, the solution can be represented as the sum of a zero-order part for the large rigid body variables and a first-order part for the small elastic variables and perturbations in the rigid body variables, where the zero-order quantities are larger than the first-order quantities by at least one order of magnitude. Then, the equations of motion can be separated into a zero-order problem for the rigid body motions alone and a $first-order problem for the elastic displacements and the perturbations in the rigid body motions. The state equations for the zero-order problem are of order 12 at most; they can be identified as the equations of flight dynamics and can be used to describe aircraft maneuvers. On the other hand, the state equations for the first-order problem are of order 12 + 2ne, where ne is the number of elastic degrees of freedom; they represent the extended aeroelasticity equations, where “extended” is in the sense that they include not only the elastic displacements but also perturbations in the rigid body variables, where the latter are sometimes referred to as “body freedoms.” The flight dynamics equations are in general nonlinear and describe the translations and rota- tions of the aircraft as if it were rigid. They can be used to design given maneuvers of an aircraft, which amounts to solving an “inverse” problem. In the commonly encountered direct problems in dynamics of rigid bodies, the forces are given and the equations of motion are solved for the state, i.e., for the positions and velocities. In the context of the present formulation unifying flight dynamics and aeroelasticity, however, the state representing a desired maneuver is given and the problem amounts to determining the engine thrust and the control surface forces permitting the realization of the maneuver; this represents an inverse problem. On the other hand, the extended aeroelasticity equations are linear, but they contain the state and forces from the flight dynamics problem as coefficients and as an input. Hence, there is a set o f extended aeroelasticity equations for every conceivable aircraft maneuver. If the flight dynamics problem represents steady level cruise or a steady level turn maneuver, then the zero-order state and forces are constant and the system of extended aeroelasticity equations is linear time-invariant. In this case, the state equa- tions lend themselves to a standard stability check, such as one based on the roots of the eigenvalue problem, to control design by commonly used techniques, such as the LQG method, and to ready integration for simulation of the aircraft response to external stimuli. If the flight dynamics prob- lem represents a time-dependent maneuver, such as the transition from one steady state to another, then the zero-order state and forces depend on time and the extended aeroelasticity state equations are linear time-varying, which precludes a standard stability analysis. However, the state equations still permit control design and response simulation.

The following literature review should help relate the present paper to previous investigations: Although flight dynamics and aeroelasticity have been developed traditionally as separate disci- plines, the need for considering interacting efforts was recognized quite e a r l ~ . l - ~ Still, relatively few attempts have been made to link the two disciplines, and when such attempts were made almost invariably the scope was quite limited. This lack of interest in linking flight dynamics and aeroe- lasticity can be attributed to a reluctance to increase the complexity of the problem to a significant extent at a time when powerful computers capable of solving such problems were not available.

As a result, problems combining flight dynamics and aeroelasticity effects have tended to be sub- jected to many simplifying assumptions designed to permit largely analytical solutions. In one of the first references on the subject, Bisplinghoff and Ashley derived scalar equations of motion for an unrestrained flexible vehicle. The equations consisted of three inertially decoupled sets, one for the rigid body translations, one for the rigid body rotations and one for the elastic deformations, the latter expressed in terms of aircraft structural natural modes. Although not stated explicitly, this implies the use of principal mean axes with the origin at the vehicle mass center. Moreover, the inertia matrix was assumed to be constant, which implies that the contributions from the elastic deformations to the inertia matrix were ignored. Aerodynamic forces for the case of small dis- turbed motions from steady rectilinear flight were given in terms of an influence function for an unrestrained aircraft. An integrated analytical treatment of the equilibrium and stability of flexi- ble aircraft was presented by Milne.5 In P a r t I, he derived linearized equations of motion about a steady state, assuming not only that the elastic deformations but also the rigid body translations and rotations were small. Although the constraint equations defining the mean axes were given, the formulation seems to have used body axes attached to the undeformed aircraft. The equations are expressed in a vector-dyadic form that does not permit a ready check for missing terms and, more importantly, does not lend itself to ready computer solutions. In Part 11, the general analysis was applied to the study of equilibrium and longitudinal stability about equilibrium of an aircraft having longitudinal flexibility only. A monograph by Taylor and Woodcock6 consists of two parts representing different approaches to the same problem. In Part I, Taylor presents a very lucid sum- mary of the equations of motion for deformable aircraft derived by Bisplinghoff and Ashley4 and by Milne.5 Following a reduction to scalar form, the equations are simplified to permit the study of some special cases. In Part 11, Woodcock uses an unorthodox form of Lagrange’s equations to derive scalar perturbation equations of motion about a given “datum motion,” not necessarily corresponding to steady level rectilinear flight; the equations are in terms of body-fixed axes. The question of aerodynamics receives scant attention in both parts. An extensive report by Dusto et al.,7 resulting in a computer program known as FLEXSTAB, integrates flexible body mechanics with a low frequency aerodynamics employing linear influence coefficients. The flexible aircraft mechanics uses free vibration modes superimposed on rigid body dynamics. Aerodynamic influ- ence coefficients are derived using a paneling scheme lending itself to empirical corrections. The equations are expressed in terms of steady perturbations about a reference motion to determine dy- namic stability by characteristic roots or by time histories following an initial perturbation or some gust disturbance. There are two major concerns. The first consists of the fact that the structural dynamics formulation is in terms of mean axes and the aerodynamics is in terms of a different set of axes, namely, “fluid axes,” where the latter move with a steady velocity relative to the un- deformed aircraft body axes; using two different sets of axes in the same formulation, without making the necessary transformation from one set to another, is a very questionable proposition.

The second source of concern is the time required to run FLEXSTAB (see Ref. 24). Several ana- lytical methods for mathematical modeling of aircraft active control system design are described by Roger,’ placing the emphasis on aerodynamics. Inconsistencies in control configured vehicles are highlighted by Schwanz? He suggests that familiarity of flight control specialists with a broad spectrum of pertinent disciplines, including aerodynamics, structures, modern dynamics and con- trol, can minimize and perhaps avoid altogether these inconsistencies. Free-free dynamic analyses of forward swept wing aircraft by Miller, Wykes and Brosnanl’ have shown that the static aeroe- lastic divergence exhibited by a cantilevered forward swept wing is replaced by a low-frequency flutter mode due to coupling between the wing divergent mode and the aircraft short-period mode.

This coupling is shown to have detrimental effects on flying qualities, ride qualities and gust loads, but these effects can be minimized by an active flutter control system. In the same spirit, Weisshaar and Zeiler” discuss the importance of including aircraft rigid-body modes in the aeroelastic anal- ysis of forward swept wing aircraft. They show that body-freedom flutter and aircraft aeroelastic divergence, not wing divergence, are the primary aeroelastic instabilities. Rodden and Love12point out that equations of motion derived using mean axes for the inertial terms and axes attached to the undeformed structure for the flexibility terms are likely to be incorrect; such flexibility terms are obtained when using structural influence coefficients. Cerra, Calico and No11 l3 developed a linear model of an elastic aircraft providing the capability of analyzing the coupling between the rigid body motions and the elastic motions. The model can be used for stability and control analyses.

As in Ref. 4, the rigid body translations, rigid body rotations and elastic deformations are assumed to be inertially decoupled. A framework for constructing simulation models for flexible aircraft is described by Arbuckle, Buttrill and Zeiler.14 The objectives are to increase simulation model fi- delity and to reduce the time required for developing and modifying these models. The framework has been applied to the development of an open-loop F/A- 18 simulation model. Buttrill, Zeiler and Arbuckle15 considered a mathematical model integrating nonlinear rigid body mechanics and linear aeroelasticity in conjunction with Lagrangian mechanics to derive the equations of motion for flexible aircraft. Undamped vibration modes satisfying first-order mean axes constraints were used as generalized coordinates. Considering a model of an F/A-18 aircraft, elastic modes sig- nificantly affected by inertial coupling were found to be aerodynamically decoupled from the rest of the model. Zeiler and Buttrill" used the extended Hamilton principle to derive equations of motion for a flexible body. The equations include inertial terms due to flexibility, as well as terms coupling rigid body and flexible momenta. In addition, a nonlinear strain-displacement relation permits centrifugal stiffening to be taken into account. The equations are used to simulate the motion of a structure spinning initially about an unstable principal axis in gravity-free vacuum.

Using Lagrange's equations, Waszak and Schmidt17derived the equations of motion for a flexible aircraft. The strip theory was used to obtain closed-form integral expressions for the generalized aerodynamic forces. Moreover, the use of mean axes permitted inertial decoupling of the rigid in terms body translations, rigid body rotations and elastic deformations, the latter being expressed of aircraft vibration modes. The modeling method was applied to a generic elastic aircraft, and the model was used for a parametric study of the flexibility effects. Nonlinear equations of motion for elastic panels in an aircraft executing a pull-up maneuver of given velocity and angular velocity were derived by Sipcic and Morino." The effect of the maneuver on the flutter speed and on the limit cycle amplitude was discussed for various load conditions. Accurate modeling of aeroelas- tic vehicles, with emphasis on the rigid body and elastic degrees of freedom, was discussed by Waszak, Buttrill and Schmidt.lg A comparison of the approach of Ref. 17 on the one hand and that of Refs. 15 and 16 on the other hand was presented and various model reduction techniques were reviewed. An integrated analytical framework for modeling elastic hypersonic flight vehicles was developed by Bilimoria and Schmidt.20 A Lagrangian approach was used to derive equations of motion including rigid body motions and elastic deformations, as well as effects due to fluid flow, Earth. The elastic deformations are represented rotating machinery, wind and a spherical rotating in terms of modal coordinates relative to mean axes. A paper by Olsen21 reveals new insights in the aeroelasticity and flight mechanics of flexible aircraft by obtaining and solving the equations of motion for an accelerating, rotating aircraft. General equations in terms of quasi-coordinates are first obtained and then reduced to the case of a "flat" airplane. The influence of gusts on the dy- namics of large flexible aircraft is analyzed by Teufel, Hanel and who present an integrated flight and aeroelastic control law reducing gust sensitivity. Moreover, the control laws, designed by p-synthesis, are such that flight maneuvers do not excite elastic motions. Konig and S ~ h u l e ? ~ de- scribe how an integral model for large flexible aircraft can be derived and how an integral control, covering flight control, load control and structural mode control, can be designed by multiobjec- tive parameter optimization. Samareh and Bhatia24presented a unified three-dimensional approach that reduces the number of interactions among various disciplines by using a computer-aided de- sign model. Results were presented for a blended wing body and a high-speed civil transport.

Schmidt and R a n e ~ ~ ~ considered the effects of flexibility on the flight dynamics of large flexible aircraft. In particular, when the frequencies of the lower elastic modes approach those of the rigid body modes the handling characteristics can suffer and the flight control system design tends to become significantly more complex. Expressing the motion in terms of components along mean axes, they add the flexible degrees of freedom to an existing simulation model of the vehicle's rigid body dynamics. The NASA Langley Research Center simulation facility was used to obtain the dynamic response of two different aircraft.

With some exceptions, the equations of motion in Refs. 3-23,25 were derived either by means of Newtonian mechanics or by means of standard Lagrange’s equations. These approaches are more suitable when the motions are expressed in terms of inertial axes and/or when the rotations are in terms of Euler’s angles. Yet, in the case of aircraft it is more convenient to express the motion in terms of components along body axes. In this regard, we should point out that this is common practice in flight dynamics, in which case the angular veolcities in terms of body axes are the well- known roll, pitch and yaw. Of course, equations in terms of inertial axes and/or Eulerian angles can always be transformed into equations in terms of body axes through coordinate transformations. It is appreciably simpler, however, to derive the equations of motion directly in terms of body axes, which can be done through the use of Lagrange’s equations in terms of quasi-coordinates?6 Motivated by problems in dynamics of spacecraft with flexible appendages, Meirovitch and Nelson27derived for the first time hybrid (ordinary and partial) differential equations of motion coupling rigid body rotations and elastic deformations. The elastic deformations were measured relative to a set of body axes attached to the undeformed spacecraft and the rotational motions were in terms of quasi-coordinates. The explicit formulation of Ref. 27 was extended by Meirovitch28 to a generic whole flexible body by deriving a set of hybrid equations of motion in terms of quasi- coordinates, treating for the first time translational velocities as quasi-velocities; the equations were then cast in state form. The equations of motion in terms of quasi-coordinates of Ref. 27 were used by Platus2’ to derive coupled equations of motion governing the aeroelastic stability of spinning flexible missiles. The coupling between the elastic deflections and rigid-body motions was nonlin- ear, but the equations were linearized so as to permit a stability analysis. The developments of Ref.

28 were extended by Meirovitch30 and Meirovitch and Stemple3I to flexible multibody systems.

Then, the approach of Refs. 28,30 and 3 1 was used by M e i r ~ v i t c h ~ ~ to produce a definitive unified theory for the flight dynamics and aeroelasticity of whole aircraft. Generic state equations describ- ing the flight of flexible aircraft were first derived in hybrid form and subsequently discretized in space. Then, using a perturbation approach, the discrete state equations were separated into a set of nonlinear flight dynamics equations for the rigid body variables and a set of linear extended aeroe- lasticity equations for the elastic variables and perturbations in the rigid body variables. Nydick and F ~ i e d m a n n ~ ~ applied the equations of motion in terms of quasi-coordinates derived in Ref. 28 to the analysis of a hypersonic vehicle in free flight. To this end, they simplified the equations by considering only the pitch and plunge rigid body degrees of freedom and small elastic displace- ments. The nonlinear equations were linearized about a trim state obtained by using a rigid body trim model and steady hypersonic aerodynamics. Flutter derivatives were calculated by means of piston theory. The generic formulation of Ref. 32 was used by Meirovitch and T U Z C U ~ ~ to carry out the derivation of explicit equations of motion in terms of quasi-coordinates for a flexible aircraft model and to cast the equations in a special state form suitable for simulation on a computer. Due to relative ease of integration into the unified formulation and computational speed advantages, the aerodynamic forces were derived by means of strip theory. Then, equations for flight dynamics and extended aeroelasticity were derived. An approach entirely different from that in Ref. 34 is proposed by Fornasier et al?5 Indeed, Ref. 35 is concerned essentially with the fluid-structure in- teraction in a flexible aircraft. To this end, it uses “temporal and spatial algorithms” to make two independently developed computer codes, one for the aerodynamics (CFD) and one for structural mechanics (CSM), work together. The scope of Ref. 35 is relatively limited, as the aircraft is as- sumed to follow a known preset trajectory, so that there are no rigid body degrees of freedom, and there are no controls. Worthy of notice is the fact that several 5 s simulations, including some of the wing tip displacement, took approximately 35 hrs on a 32-processor Sgi computer.

The present paper represents an extension of the developments in Ref. 34. In addition to some modeling refinements, the paper contains a numerical example for a model of a flexible aircraft containing 76 states, 12 rigid body states and 64 elastic states. Two flight dynamics problems are considered, the steady level cruise and a steady level turn maneuver. The corresponding extended aeroelasticity problems are derived and used to design feedback controls guaranteeing the vanish- ing of the rigid body perturbations and the elastic vibration, and hence the stability of the maneuver and the comfort of the flight. The control design consists of a linear quadratic regulator in conjunc- tion with a stochastic observer. The integrated process is demonstrated by means of a numerical example including a variety of rigid body and elastic displacements time simulations together with the corresponding controls time histories, all carried out on a 1 GHz PC using MATHEMATICA.

2. Hybrid Equations of Motion in Terms of Quasi-Coordinates

We regard the aircraft model shown in Fig. 1 as a flexible multibody system subjected to gravity, aerodynamic, propulsion and control forces, where the bodies can be broadly identified as the fuse- lage, wing and empennage. The motion of the aircraft can be conveniently described by attaching a reference frame x f g f z f to the undeformed fuselage (Fig. I), as well as corresponding reference frames x,yu,z, and Z e Y e Z , to the wing and empennage, where the various reference frames repre- sent respective body axes. Then, the motion can be described by six rigid body degrees of freedom of the fuselage body axes, three translations and three rotations, and by the elastic deformation of every point of each flexible component relative to the respective body axes.

From Ref. 32, and making provisions for members in torsion, as well as for structural damping, the hybrid equations of motion for the whole flexible aircraft in terms of quasi-coordinates have Figure 1. Flexible Aircraft Model the generic form d L d L - F

+Wr--cfw- dVf

i = f (fuselage), w (wing), e (empennage) where L = Lagrangian for the whole aircraft Vf, wf = vector of translational, angular quasi-velocities of xfgfzf V', W f = skew symmetric matrices derived from Vf, w Cf = matrix of direction cosines between Lcfgfzf and X Y Z (inertial axes) Rf = R f ( X f , Y f , Zf) = position vector of origin Of of xfgfzf relative to X Y Z E, = matrix relating Eulerian velocities to angular quasi-velocities 8, = symbolic vector of Eulerian angles between xfyfzf and X Y Z ui, v i = elastic displacement, velocity vectors for body i +i, ai = elastic angular displacement, velocity vectors for body i Li = Lagrangian density for body i exclusive of strain energy

Yui, FQi = Rayleigh's dissipation function37densities for body i

Cui,CQi = matrices of stiffness differential operators for body i F, M = resultant of gravity, aerodynamic, propulsion and control force, moment vectors acting on the whole aircraft in terms of fuselage body axes components Ui, !Pi = resultant of gravity, aerodynamic, propulsion and control force, moment density i vectors for body Assuming that axes x~fyfzf are obtained from axes X Y Z through the sequence of rotations $ about 2 to axes zlylzl, 0 about y1 to z2y2z2 and 4 about x2 to xfyfzf, the matrices Cf and E f have the form36 in which s = sin, c = cos. The elastic displacement vectors u , and +i are subject to given boundary conditions at the interface between bodies. Equations ( 1 ) involve the Lagrangian L = T - V, in which T is the kinetic energy and V the potential energy, the Rayleigh dissipation function densities 9 . i and .k&, which contain the information about the structural damping, and the stiffness operators Cui and CQi, which are related to the strain energy. Before more explicit equations of motion can be derived, it is necessary to produce these quantities. The kinetic energy for the whole aircraft can be written as where - T -

T , = i SV, Vidmi, i = f, w, e

(4) are kinetic energies of the individual components, in which Vi are velocity vectors of typical points in the components and dmi are corresponding mass differential elements. The velocity of a point in the fuselage can be written as

Vf(rf, t> = V,(t> + Ff + W f , t>I'[wr(t> + ' Y f h t>l + vf(rf, t )

2 v, + (Tf + iif)'Wf + + f + V f

( 5 ) where rf is the nominal position of the mass element dmf, ?f and iif are skew symmetric matrices26 corresponding to rf and uf and vf and af are elastic velocities associated with bending and tor- sion, respectively. Then, denoting by C, the matrix of direction cosines between x,yu,z, and in which rf, is the radius vector from the origin of xfpfzf to the origin of xU,pwzw, is the angular velocity of the fuselage at rfw due to bending and afw = [afw 0 0 1 Itfw is the elastic velocity of the fuselage at rfw due to torsion. The velocity Ve(re, t ) of a point on the empennage can be obtained from Eq. (6) by simply replacing w by e. The total kinetic energy can be expressed in the general form The potential energy can be expressed in terms of the operators Li, LQi (i = f, w, e), but is more conveniently expressed as the strain energy. Moreover, the Rayleigh dissipation function densities can be expressed in the form where G ~ , cQi are damping functions and Eli, GJi are flexural and torsional rigidities ( i = f, w, e ) .

and the operators .Cui and C,,,, Equation (S), the potential energy V, the functions .& and when inserted into Eqs. (l), permit the derivation of explicit hybrid equations of motion. Because for all practical purposes it is not feasible to work with hybrid equations, we do not pursue this subject any farther, and approximate instead the partial differential equations by sets of ordinary differential equations.

3. Spatial Discretization of the Distributed Variables For the most part, aircraft are modeled as discrete systems, either by regarding the inertia and stiffness properties as lumped from the onset or by describing their elastic motions in terms of aircraft structural modes. For undamped structures in vacuum, the use of structural modes yields In the case complete decoupling, Le., it yields a set of independent modal equations of motion.

of aircraft, however, complete decoupling is not possible, even when aircraft structural modes are used, as the aerodynamic forces guarantee that the equations of motion remain coupled. Spatial discretization using aircraft structural modes not only does not offer any particular advantage but also has the disadvantage that some geometric details of the aircraft are lost in the process. Hence, a discretization procedure that does not require the structural modes, which may not be readily available and/or may not be compatible with the rest of the formulation, and retains some sense of the aircraft geometry seems desirable. Consistent with this, we consider spatial discretization of the individual aircraft components separately. To this end, we use either the Galerkin method or the finite element method37and introduce the expansions where and are matrices of shape functions and qui and q@i are corersponding vectors of generalized coordinates. Note that the choice of shape functions is very important. Indeed, for accurate modeling, the shape functions must reflect the mass and stiffness characteristics of the components. Some guidelines concerning the choice of shape functions can be found in Ref. 37.

Moreover, we denote the associated generalized velocities by sui(t) = qui(t), s@i(t) = q*i(t), i = f , w , e (1 1 ) In anticipation of later needs, we write the velocity vectors for points on the individual components in the two discrete forms in which Expression analogous to QUf can be written for and @ , , , . Moreover, V, ( re,t) can be obtained from Vw(rw,t) by simply replacing w by e. Inserting the first forms of Eqs. (12) into Eq. (8) and where are damping matrices.

p = [sf pzf p z f pzu, pZe p;f

Next, we denote the momentum vector for the whole aircraft by p&, p,,] T T = [p: p r . . . p:jT, so that we can write

p = m/av = MV (20)

where the individual momenta are given by

P V f = m/avf = p1= C M l j V j

j=l j=1

Puf = m/asuf = p3 =

MSjVj j = 1 Finally, adding some obvious kinematical identities to the discretized version of Eqs. (I), the state equations can be written in the special form

Rf = Cf’Vf, bf = Ey’wf; qui = sui, q + i = s+i, i = f , w, e

Moreover, Kui = J@~iLui@uidDi, KQi = J’@;iL+i@$idDi, i = f, w, e (24) are component stiffness matrices. In practice, Kui and KQi ( i = f , w, e) can be obtained with greater ease from the strain energy directly, as shown in the Numerical Example.

The quantities F, M, Qui and QILi (i = f , w, e) appearing in the state equations, Eqs. (22), represent generalized forces. They are related to the actual forces, which consist of the distributed fi(ri, t ) over component i due to gravity, aerodynamics and controls and the engine thrust force FES(r - rE), where S(r - rE) is a spatial Dirac delta in which r E denotes the location of the engines. If some control forces are concentrated, they can also be treated as distributed, as in the case of the engine thrust. The relation between the generalized forces and the actual forces can be obtained by means of the virtual work, which can be expressed as (25) ST@ = E, J ' , [ f : + FgS(r - r,)]SRrdD, where is the virtual quasi-displacement vector of a typical point on component i (i = f, w, e).

The vector S R is related to the virtual quasi-displacement vectors corresponding to the quasi- velocities used to describe the motion of the aircraft components. Indeed, using Eqs. (12), we can write

-

SR; = SR; + (Tf + @,fq,,)'SO; + @,fSqUf + T T @ @ f S q + f

- -

SR, = C,SR; + [Cw(Tfw + @,fwq,f)T + ( T , + @uwq,wt)TCu~]60~

(26) -T

+ (~,C,a@',f, + C w ~ @ u f , ) ~ q t L f + @ u w S q u w

+ (TuTc,@@f, + cuf;u!@~ful)Sqwf + C @ @ w S q @ m

and we note that s R , * can be obtained from mu, by replacing w by e. Inserting Eqs. (26) into Eq.

(25), and collecting terms, we can write the virtual work in terms of virtual generalized displace- ments as follows:

S W = FT6R; + MT68; + C,<Q~,bqzlz + Q$,dqQz) (27)

from which, assuming that the engines are mounted on the fuselage (Fig. l), we obtain

F = S , , [ff + F E S ( ~ - rE)]dDf + C:JDuf,dD, + C,TJDefedDe

-

= JD,(Tf+',f9,,f)[ff + F E S ( r - r E ) ] d D f + S ~ _ [ ( T f w + @UfWquf)':

- - -

+~,T(r,+@,s~lw)]fiL~d~zL' + S,, [(.re + @ u j e q u f ) C , T + C,T(Fe + @ u e q z l e ) ] f e d ~ e

Q u j = JD,@:f[ff + F d ( r - rE)]dDf ,,fD,(f$CulA@~f~- C w @ ~ f w ) ~ f w d D w

(28) + JD,(pTCeA@ufe + Ce@ufe)TfedDe, &@f = JD, @:fTj[fj + F&(r - r ~ ) ] d D f + sou, ( ~ ~ C W @ Q ~ W 4- C W ~ ~ , @ @ ~ W ) ~ ~ W ~ D U

+ JD, ( T T C e Q Q j e + CeTfr,@+je)TfedDe

& , , = JD,@:,fzdD2, Q@t = ~ D , @ ~ z T z f z d D z , i w, e To complete the state equations, Eqs. (22), it is necessary to derive the stiffness matrices K,, and KQz ( i = f, w, e), the aerodynamic forces and the control forces.

4. The Aerodynamic and Gravity Forces The forces acting on the aircraft can be identified as the aerodynamic, propulsion, control and gravity forces. Equations (28) give the generalized forces in terms of actual distributed and con- centrated forces, where the first imply prescribing the forces at every point of the aircraft. Of the forces acting on an aircraft, the aerodynamic and control forces require special attention. We discuss the aerodynamic forces in this section and the control forces in a later section.

There are a number of aerodynamic theories available, some of them capable of prescribing the pressure distribution at every point of the aircraft. However, as pointed out in the Introduction and Literature Review, any aerodynamic theory to be used in a dynamic simulation such as that described in this paper must satisfy certain requirements. Indeed, one of the requirements is that with the present general dynamic for- the aerodynamic forces be expressed in a form compatible mulation, which implies that they be in terms of the same variables and be referred to the same body axes attached to the undeformed aircraft as here. Another requirement is that the aerody- namic forces lend themselves to sufficiently fast computation as to permit time simulation of the aircraft behavior. Such computational efficiency does not appear to be within the state of the art.

Hence, an aerodynamic theory capable of fitting within the framework of a computer simulation of the type envisioned here must yet be developed. Such a theory need not be unduly accurate, be- cause a feedback control design tends to be sufficiently forgiving to tolerate small deviations from the actual aerodynamic forces. Until such a theory becomes reality, it is still possible to demon- strate how aerodynamics fits in the integration process by using an existing theory satisfying the requirement described above, namely, strip even though the theory may not be entirely suitable otherwise.

(i = f, w, e) acting on the To derive the aerodynamic forces included in the distributed forces aircraft by means of strip theory, we regard the fuselage, wing and empennage as two-dimensional aerodynamic surfaces. The lift force per unit span of fuselage can be written a$8 (29) Cf = Q f C f C L a f a f = 4 f C f C L f where cf is the chord, CLaf the slope of the lift curve, CLf the lift coefficient and - ( 3 0 )

qf = ;p(v;L, + V&), af = tan-’(Vfz/Vfz>

are the dynamic pressure and angle of attack, respectively, in which p is the air density and Vfz and Vfz are components of the velocity vector Vf, Eq. (5). Similarly, the drag force per unit span of fuselage is given by where C O ~ O is the drag coefficient corresponding to af = 0 and k, is a constant. The fuselage has also vertical surfaces subjected to aerodynamic forces. The lateral force per unit span can be expressed as Sf = Q s f C s f C s P f P f = 4 s f C s f C s f (32) csf is the lateral chord, Cspf the slope of the lateral force curve, Csf the lateral force coeffi- where cient and are the dynamic pressure and the angle of attack of the lateral force, respectively, in which V j , and Vj, are components of Vj. Hence, the aerodynamic forces on the fuselage can be written in the vector form

ef sin aj - d cos af

s j sin ,f?j

0 ] 7 f s f = [ - s f r P f ] (34) f a j

[ -e, cos aj - d j sin aj

In a similar fashion, the lift and drag per unit span of wing are given by e , = qwc, (CLau&t + ~ L b a & ) , d, = quvc, ( CDWO + F;,c;, ) = QUI c, ( C o w 0 + ku!c;au>a: ) (35 where c, is the chord, Sa an aileron rotation, CLba a control effectiveness coefficient17 and 1 -

= tan-'(V,,/V,y) + +wz

(36) ( 2 , = &$/ + v : , ) , in which the velocity components V,, and VWz of V, can be obtained from Eq. (6). Moreover, ~ , z is the angular displacement of the wing about axis 2 , due to torsion. There is no meaningful so that lateral force on the wing, tu, sin a, - d, cos a, ( 3 7 ) 4, cos a , - d, sin a , The empennage has both lift and lateral surfaces. The lift, drag and lateral forces per unit span of empennage are l e = q e C e ( C L a e a e + C L b e S e ) , d e = q e C e ( C D e 0 + keCiaeaf), s e = q s e C s e ( C s p e P e + CsSrSr) (38) where c, and c,, are the chords, Se and S , are rotations of the elevator and ruder, C L ~ , and Csb, are respective control effectiveness coefficients and Hence, the aerodynamic force vectors per unit span of empennage can be written as 0 0

e, sina, - de cos a, ] , f s e = [ s, sinp, ]

-le cos a, - de sin a, -s, cos pe For a typical component, the lift, lateral force and drag per unit span are applied at the line Eqs. (28), the domain of integration for the terms involving of aerodynamic centers. Hence, in the aerodynamic forces is the line of aerodynamic centers. The gravity forces per unit volume of components are simply Note that the aerodynamic and gravity forces are in terms of respective component body axes.

5. A Unified Theory for Flight Dynamics and Aeroelasticity To complete the discussion of the forces acting on an aircraft, we turn our attention to the control forces. Aircraft control is carried out by means of control surfaces, as well as by the engine thrust.

Before we consider the problem of control design, it will prove beneficial to examine the nature of the controls. Controls are of two general types, one designed to steer the aircraft as a whole and the other to suppress the effects of any undesirable disturbances. The first type involves rigid body motions of the aircraft, which are in general large, and traditionally lies in the domain of flight dynamics. On the other hand, the second type involves elastic deformations of the aircraft, which tend to be small compared to the rigid body motions, and traditionally lies in the domain of aeroelasticity. Hence, the formulation given by Eqs. (22) can be regarded as spanning the fields of flight dynamics and aeroelasticity.

From the above discussion, it appears that control of the aircraft as a whole is likely to be different in nature from control of disturbances. In this regard, we observe that the state equations, Eqs. (22), are in general nonlinear and of high dimension, where the nonlinearity can be traced to the large rigid body variables. On the other hand, the high dimensionality can be traced to the small elastic variables. In view of this, a solution by a perturbation approach seems indicated, which amounts to a separation of the problem into a zero-orderproblem for the large variables and ajrst- order problem for the small variables, where the difference between the large and small variables is one order of magnitude, or more. Physically, in the zero-order problem the aircraft executes a given maneuver as if it were rigid, in which case the mathematical formulation consists of a maximum of six coupled, generally nonlinear ordinary differential equations, three for rigid body translations and three for rigid body rotations. They correspond to the equations commonly used in JEight dynamics. On the other hand, the first-order problem involves the elastic deformations, as well as small perturbations in the rigid body variables. In view of the inclusion of the latter, the first-order problem defined here represents an extended aeroelasticity theory, in which the rigid body degrees of freedom are included in a natural manner. This is in contrast with some in which “body freedoms” are included in an ad hoc manner. We note that occasional p r a ~ t i c e , ~ ~ ~ ~ ~ the solution of the zero-order problem enters into the first-order problem, so that this new theory provides one set of extended aeroelasticity equations for every conceivable rigid body maneuver o f the aircraft, rather than the single set of equations commonly associated with steady cruise. We express the perturbation solution as follows: where the superscripts (0) and (1) denote orders of magnitude. All the quantities related to the elastic deformations are regarded as being of first order. Then, inserting Eqs. (42) into the state equations, Eqs. (22), and separating different orders of magnitude, we obtain the zero-order prob- JEightdynamics problem lem, or the in which CF’ and E?’ can be obtained from Cf and E f , Eqs. ( 2 ) , by replacing $, 8 and 4 by $(‘I, 8(’) and 4(’), respectively. Moreover, from Eqs. (28), the zero-order generalized force and moment are given by where in which Then, from Eqs. (28), the first-order generalized forces are

’

where where hfj!) and Mi;) are the zero-order part and first-order part of M j , Eqs. (15). Equations (59) are to be solved for Vy), w (1) , s,,, sUw , . . . , Sqve in terms of pi,,, ( 1 ) p,,, ( 1 ) puf,puu,, . . . , PQe V, ( 0 ) and do) and the result inserted in Eqs. (50).

The zero-order problem, or flight dynamics problem, represents an inverse problem, which amounts to determining the controls permitting realization of a given rigid body maneuver. The first-order equations representing the extended aeroelasticity problem are linear and tend to be of high order. Moreover, they contain the zero-order variables Vy’ and w y ’ , representing a given maneuver, as coefficients and as an input. If Vy’ and wy’ are constant, then the system is time- invariant, and if Vy’ and wy’ depend on time, then the system is time-varying. In either case, controls can be designed by various methods. In the time-invariant case, a stability analysis for the closed-loop system can be carried out by solving an eigenvalue problem. Such a stability analysis is precluded in the time-varying case. Simulation of the response of the closed-loop system to external excitations, such as gusts, can be obtained in both the time-invariant case and time-varying case.

6. Control Design Flying aircraft are subjected to various disturbances tending to drive them from the intended ma- neuver and to cause vibration. If the system is controllable, these effects can be suppressed through controls, which are carried out by means of actuators; in the case of aircraft they consist of the en- gine thrust and the control surfaces. A system is said to be controllable if there exists a piecewise continuous input that will drive the initial state to any final state within a finite time interval. Sys- tem controllability can be determined by checking the rank of a so-called controllability matrix?'

which may not be feasible if the system order is large. In practice, controllability can be ascer- tained on physical grounds by making sure that the input forces, namely, the forces due to the engine thrust and control surfaces, affect all the state variables.

As indicated in earlier sections, there are two types of controls, one type designed to permit the aircraft to execute a desired maneuver as if it were rigid and the other type to reduce any deviations from the rigid body maneuver to zero, which amounts to suppressing vibration and perturbations in the rigid body motions of the aircraft. The first is associated with the flight dynamics problem and the second with the extended aeroelasticity problem. In general, the engine thrust and control surfaces are designed so as to ensure that the aircraft is able to carry out the required maneuvers, as well as to suppress any undesirable disturbances, thus addressing the needs of both the flight dynamics problem and extended aeroelasticity problem.

Using Eqs. (43), we write the flight dynamics problem, or the zero-order problem, in the com- pact state form

X ( O ) ( t ) = f [x@) ( t ) , vg (t)] + do) [V$' (t)] U(O) ( t )

(60) which, from Eqs. (49), must be considered in conjunction with (0) (0) (61) P$)(t) = Mrb vrb (t> ~ ( o ) ~ (o)T (o)= T is

where the zero-order quantities are identified as follows: x ( O ) = [RY)' f P V f Puf 1

the state vector, f is a nonlinear function of the state vector and the zero-order rigid body velocity vector ~ $ 1 = [vYlT u Y ) T ~ ~ , ~ ( 0 ) is a coefficient matrix, u(0) = [ F ~ (0) S , (0) S , ( O ) 6, (0) 1 T is the control vector, in which FE', SF', @' and 6 " ' are the engine thrust and control surfaces angles,

p$ = [ p t r p;jTlT is the rigid body momentum vector and

is the rigid-body mass matrix. In the context of the present integrated approach, Eqs. (60) and (61) represent an inverse problem, in the sense that a state vector x(O) describing a desired maneuver is postulated and a force vector do) permitting a realization of the given maneuver is determined.

Next, we assume that x(O)(t) and V $ ) ( t ) are known and use Eqs. (50) and (59) to express the extended aeroelasticity problem, or first-order problem, in the form

X ( ' ) ( t ) = A(t)xc1)(t) + B(t)u(')(t) + [0 IITFext(t)

(63) and

(1)T f p T T T (1)T p t j T *Tr T T T

respectively, in which x(l) = [Rr q : r q,, . . q,, pvr

Pu, - * P*el

is the first-order state vector, A(t) = A[x(')(t), V z ' ( t ) ] and B(t) = B [ ~ ( ~ ) ( t ) , V $ ) ( t ) ] are coefficient matrices, u(l)(t) = [ F g ) C 5 2 ) C5L1) C5:1)]T is a first-order control vector, Fext is an T i s the external disturbing force vector, such as due to gusts, p(') = [ p v f pur p u f T puw T . . . pqe] first-order momentum vector and are zero-order and first-order extended mass matrices. Note that Fext is regarded as a disturbing force of a transient nature whose effects will be eventually suppressed by the controls.

Equation (63) represents a set of linear equations, and the objective is to find a control vector u(l)(t) that drives the state vector x(') to zero. To this end, we consider a linear regulator whereby the control vector is a linear function of the state vector. In particular, we consider a linear quadratic regulator (LQR) in which the objective is to determine an optimal control vector minimizing the quadratic performance measure4' J = $ ~ ( ' ) ~ ( t f ) H x ( ' ) ( t f ) + s l t , " ~ " ' ~ ( t ) Q ( t ) ~ ( ' ) ( t ) + ~(')~(t)R(t)u(')(t)]dt (66) where H and Q are real symmetric positive semidefinite matrices R is a real symmetric positive definite matrix, t o is the initial time (commonly assumed to be zero) and t f is the final time. It is shown in Ref. 40 that the optimal feedback control vector is given by u(l)(t) = -R-l(t)BT(t)K(t)xcl)(t) = -G(t)x(')(t) (67 1 where G ( t ) = R-'(t)BT(t)K(t) (68) is a control gain matrix, in which K ( t ) is a real symmetric matrix satisfying the transient matrix Riccati equation

(69) I ? = -Q - ATK - K A + KBR-'BTK, K ( t j ) = H ( t f ) = H

a nonlinear equation which must be integrated backward in time from t f to to. Rather than inte- grating a nonlinear matrix equation, it is advisable to transform the problem into a linear one. To this end, we consider the transformation K ( t ) = E ( t ) F - ' ( t ) (70) where E ( t ) and F ( t ) can be obtained by solving the linear equation which again must be integrated backward in time. Inserting Eq. (67) into Eq. (63), we obtain the closed-loop equation

X ( ' ) ( t ) = [A(t) - B ( t ) G ( t ) ] ~ ( ' ) ( t ) + [O IITFext(t)

(72) which can be integrated to simulate the system response.

The problem is considerably simpler when the zero-order solution is constant, such as in steady level flight, as in this case the coefficient matrix A is constant. Then, if the system is controllable, H = 0 and Q and R are constant, the Riccati matrix K approaches a constant value as tf increases without bounds. In this case, Eq. (69) reduces to the steady-state matrix Riccati equation

-Q - ATK - K A + KBR-'BTK = 0

(73) a nonlinear algebraic matrix equation, which can be solved by means of Potter's algorithm?' and the gain matrix G becomes constant. The closed-loop equation reduces to one with constant coef- ficients, or

X(')(t) = ( A - BG)x(')(t) + [0 1ITFex+,(t)

(74) which can be used for response simulation. For a stability analysis, we solve the associated eigen- value problem ( A - BG - M)x = 0 (75) The closed-loop system is stable if all the eigenvalues are pure imaginary and/or complex with negative real part.

is optimal in the sense that it minimizes the performance index J , but the The control vector physical merit of this optimality is debatable. In fact, it is often necessary to adjust the otherwise arbitrary weighting matrices Q and R to achieve a desirable system performance. The real value of the LQR algorithm is that it guarantees a stable closed-loop system.

7. Optimal State Observer In Sec. 6, we carried out the control design in two stages. In the first stage, we postulated a desired aircraft maneuver x(O) and used an inverse approach to determine the control vector u(O) permitting realization of the maneuver. In the second stage, we designed a feedback control vector u(') ensuring stability of the maneuver, which amounts to driving the perturbation vector x(') to zero. The control vector is given by Eq. (67), in which G is the control gain matrix.

Implementation of the control law, Eq. (67), requires knowledge of the state vector ~ ( ' 1 , which can be obtained through measurement. This creates somewhat of a problem, as measurements represent real quantities and our state vector consists of abstract generalized coordinates rather than real coordinates. Moreover, we feed back only perturbations from the maneuver variables. In the absence of external forces, the state equations describing the extended aeroelasticity problem have the vector form

X("(t) = Ax(')(t) + Bu(')(t)

(76) where the coefficient matrices A and B were defined in Sec. 6 . Moreover, denoting the measure- ment vector by y ( t ) , we can write

y ( t ) = y(O)(t) + Y ( ' ) ( t )

(77) where y(')(t) is the contribution from the maneuver variables and y(')(t) is the contribution from the perturbations in the maneuver variables. We express the latter in the form y y t ) = Cx(l)(t) (78) and refer to y ( ' ) ( t ) as the output vector. The assumption is made here that the system is observ- able:' which implies that the initial state x(')(O) can be deduced from the outputs within a finite time period. Observability can be established by checking the rank of a so-called observability matrix,40 a matrix involving the matrices A and C. For large-order systems, working with the ob- servability matrix may not be feasible, and in practice the choice of sensors ensuring observability must be made on physical grounds. In particular, the choice must be such that the sensors signals permit reconstruction of the state at all times. The task of determining the matrix C relating the output vector to the state vector is discussed later in this section.

In reality, the state vector cannot be determined exactly from the output vector and must be estimated. A device permitting an estimate of the state vector is known as an observe#' and can be expressed in the form where KO is an observer gain matrix. Subtracting Eq. (79) from Eq. (76), we obtain e ( t ) = [A - KoC]e(t) in which KO such that the vector e ( t ) represents the observer error vector. The objective is to find a matrix approaches zero as t increases. In this case, x(t) + x ( t ) with time. For a time-invariant system, this amounts to determining KO so that all the eigenvalues of the matrix A - KoC, known as the observer poles, lie in the left half of the complex plane. In implementing feedback controls, we must use the estimated state X(')(t), because the actual state x(l)(t) is not available. Hence, Eq.

(67) must be replaced by u(')(t) = -GX(')(t) (82) One question of interest is how the choice of observer poles affect the choice of controller poles.

To answer this question, we use Eqs. (81) and (82) and rewrite Eq. (76) in the form

X ( ' ) ( t ) = [ A - BG]x(')(t) + BGe(t)

(83) Equations (80) and (83) can be combined into Because the coefficient matrix in Eq. (84) is block-triangular, the poles of the combined system consist of the sum of the poles of A - BG and the poles of A - KoC, so that the observer poles can be chosen independently of the controller poles.

An optimal observer gain matrix can be obtained by adopting a stochastic approach, leading to the so-called Kalman-Bucy filter, in contrast to a deterministic observer known as a Luenberger observer. To this end, we rewrite Eqs. (76) and (78) in the form

X ( l ) ( t ) = Ax(')(t) + Bu(l)(t) + ~ ( t )

( 8 5 ) and

y ( l ) ( t ) = Cx(')(t) + w(t)

(86) where v(t) is known as the state excitation noise and w(t) as the observation noise, or sensor noise.

It is customary to assume that v(t) and w(t) are white noise processes, with intensities V(t) and W ( t ) , respectively, so that the correlation matrices have the form (87)

E{V(tl)VT(t2)} = V(tl)W2 - h ) , E{W(t1)WT(t2)) = w(tl)W, - tl)

and that they are uncorrelated, so that E{v(tl)WT(t2)) = E(W(tl)VT(t2)} = 0 (88) in which E{ e} denotes the expected value.

The stochastic observer has the same form as that given by Eq. (79) in which the optimal gain matrix is determined by minimizing the performance measure Jo = E{eT(t)U(t)e(t)} (89) where U ( t ) is a symmetric positive definite weighting matrix. From Ref. 40, the optimal observer gain matrix is given by K,*(t) = Q(t)CTW-'(t) (90) where Q(t) is the variance matrix of e(t) satisfying the transient matrix Riccati equation

Q(t) AQ(t) + Q(t)AT + v(t) - Q(t)cTw-l(t)cQ(t>, Q(0) = Qo (91)

In the time-invariant case, Eq. (91) reduces to an algebraic matrix Riccati equation yielding a steady-state optimal observer gain matrix. The main problem in implementing a Kalman-Bucy filter lies in the selection of the noise intensities V ( t ) and W ( t ) .

At this point, we turn our attention to the determination of the matrix C defined by Eq. (78).

To this end, we must first specify the measurement vector y(’)(t), which in turn depends on the sensors used. In inertial navigation41, which represents the process of determining the position and attitude of a moving vehicle from self-contained inertial measurements made on board of the vehicle, the system consists of a platform containing accelerometers sensing translational motions, gyroscopes sensing angular motions and a computer capable of integrating the sensors signals to generate the state. Inertial navigation is widely used for aircraft, in which case there are two accelerometers aligned with the North and East directions and three rate gyroscopes with the spin axes aligned with the North, East and zenith directions. To ensure that gravity does not contaminate the accelerometers signals, the platform is made to rotate continuously so as to remain normal to the local vertical. The vertical position of the aircraft is measured by means of an altimeter. In view of this, we can assume that the system measures the vectors Rf and 8, giving the position and attitude relative to axes XYZ. In this regard, it should be mentioned that we referred to axes XYZ as inertial, but in reality they represent earth-fixed axes. If the current formulation is used to describe relatively long flights, then proper allowance must be made for the rotation of the earth.

Because the same process can be used to determine the zero-order position vectors R(O) and do), we can assume that the measurement system is capable of yielding RY) and 0:).

The above process can be used to measure the perturbations in the rigid-body motions of the aircraft and the question remains as to how to measure the balance of the variables, namely, the elastic variables. To this end, we assume that there are Ni (i = f, w, e) sensors measuring veloci- ties at the points P k ( k = 1 , 2 , . . . , Ni) of the aircraft components and express the output vector in the form where

(llT = [ V j l ) T ( ~ l , t ) Vj1’T(~2, t ) . . . vi - (1)T ( P N ~ , t)lT, i = f, w, e

(93) Yi in which Vl”(Pk, t) are vectors of velocity measurements. From Eqs. (12) and (53), we can write - - (94) where Also from Eqs. (12) and (53), we can write s ; ~ In a similar fashion, we can write where

-

block-diag Cek = [0 0 Cei;tF)@,fe 0 C e ~ ; ) @ u e ( P k ) 0 0 0 Equations (94), (96) and (98) are in terms of d(') and V(l). They can be transformed into expressions in terms of the state x(l) = [d(l)Tp(l)T]T by considering Eq. (64). To this end, we observe that = M(')(q,f, q,,, que), so that we can write M(l)V(O) = Mvd(l) (loo) where MV is a matrix depending on V(O). Hence, using Eq. (64), we obtain V(l) = (M(0))-l(p(l) - Mvd(l)) (101) so that I Finally, using Eqs. (92), (93), (94), (96), (98) and (102), we conclude that I 0 0 ... 0 0 I 0 ... 0 block-diag Cf1 block-diag Cf2 . . . . . . . . . . . . . .

block-diag Cf Nf block-diag CU,l C = block-diag Cu?2 . . . . . . . . . . . . . .

block-diag C w ~ , block-diag C,l block-diag Ce2 . . . . . . . . . . . . . .

block-diag C e ~ , 8. Numerical Example The flight of a flexible aircraft is fully described by Eqs. (43), (49), (50) and (59). A solution of these equations requires the aircraft geometry, the mass and stiffness distributions and the aerody- namic coefficients. Information pertaining to an actual aircraft was made available by an aircraft manufacturer and is listed in the Appendix.

The solution of the flight dynamics equations, Eqs. (43) and (49), requires the matrices Cy’, E?’, the total aircraft mass m, the matrix S ( O ) of the first moments of inertia of the undeformed aircraft and the inertia matrix J(O) of the undeformed aircraft. The matrices Cy’ and E?’ can be obtained from Cj and E f , Eqs. (2), by simply replacing $, 8 and 4 by $ ( O ) , t 9 ( ’ ) and @(’), respectively. The aerodynamic forces are given by Eqs. (45) and the required coefficients are given in the Appendix. Other data required is as follows: engines locations, rEl = [-108.62 37.0 -13.96IT in, r~~ = [-108.62 - 37.0 - 13.96IT in; total aircraft mass, m = 33.5896 lb - s2/in; matrix of first moments of inertia and inertia matrix -134.6827 0 0 lb * s 2 0 0 0 183183.4257 4.7745 -37624.9453

lb - in - s2

81.1583 566328.8970

704218.5794 1

Some of the above quantities involve the matrices of direction cosines between the various com- ponents body axes and the fuselage body axes. These matrices are listed in the Appendix. The flight dynamics problem essentially consists of setting the control surfaces and the engine thrust for a given aircraft maneuver, of which the steady cruise is a special case. It essentially amounts to the problem of “trimming” the aircraft. As far as this paper is concerned, it provides the input Vy’, uf (0) , pvf (0) and pUf (0) to the extended aeroelasticity problem.

The solution of the extended aeroelasticity equations, Eqs. (50) and (59), requires an explicit choice of the structural model for the aircraft of F Fig. 1. As a first approximation, the fuselage, wing, and both the horizontal and vertical sta- bilizers in the empennage are modeled as beams I C clamped at the origin of the respective body axes and undergoing bending and torsion. The fuse- lage undergoes the bending displacements ufv and u f z and the torsional displacement $fz, as shown in Fig. 2, so that uf = [0 ufv ufzIT and +f = [ $ f z 0 O]*. On the other hand, the wing and the stabilizers undergo only one bending and one torsional displacement each. Note that, as Figure 2. Aircraft Components Undergoing Bending and Torsion customary, displacements are measured relative to the elastic axis (Fig. 2). Each clamped beam is assumed to be discretized by the Galerkin method in conjunction with two shape functions per displacement component. For bending, the shape functions are chosen as the eigenfunctions of a I ~ uniform cantilever beam and for torsion, the eigenfunctions of a uniform clampe-free shaft

4Qir = sin(2r - 1)7rxi/'2Li, r = 1,2; i = f, w, e

( 106) where Li is the length of the cantilever beam. In this regard, we note from Fig. 2 that the fuselage is modeled as two cantilever beams clamped at O f , one pointing to the aircraft nose and the other to the tail. New quantities entering into the first-order equation are Cy) and E y ) , which can

be obtained from Cf and E, by using Of = O y ) + 6:) and letting the components of 0:) be

small, F(l), M(l), Qui and Q Q ~ (i = f , w, e), which are given by Eqs. (54), and Qui and which contain +uir and 4Qir as given above. Moreover, the stiffness matrices are obtained from the potential energy as follows:

(d2ufy/ax;)2 + E'fy ( d 2 u f Z / a x 2 f ) 2 + G J f Z (a$fZ/axf)2] d D f

= J D f + J D 7 L ' [ E I w ( ~ ~ w z / a x ~ ) ~ + G J w ( a @ w z / a ~ , ) ~ ] d D w + S , , [EIe ( a ' . e z / a x ; ) '

+ G J e ( ~ e z / a x e ) ~ ] dDe = ixi(qTiKutqui + qgiK+iq@i) (107)

where are the desired stiffness matrices, in which primes denote differentiations with respect to xi.

Equations (63) for the model in question are of order 76, but the equations are linear. However, in the case of certain aircraft maneuvers, the systems is time-varying.

To demonstrate the ideas, we consider two cases, steady level flight and steady level turn ma- neuver.

i. Steady level flight For steady level flight, the zero-order velocities are defined by

v:' = C?)[V(O) o 0 1 ~ = constant, wy' = o

(109) where V(O) is the aircraft forward velocity. From Eqs. (49), we conclude that the zero-order mo- l menta are pvf (0) = r n ~ y ) = constant, pWf

(O) = S ( O ) V ~ ' = constant

(1 10) I Hence, from the second line of Eqs. (43), we have I F(O) = 0, M(O) = 0 (1 11) The implication of Eqs. (1 1 1) is that, for steady level flight, the forces and moments due to the engine thrust, aerodynamic forces, gravitational forces and control forces balance out to zero. The angle of attack can be expressed as cy$ = tan-l(vf, (0) = cy0 = constant (1 12) For level flight, we have + ( O ) = q 5 ( O ) = 0, so that the pitch angle is equal to the angle of attack, or 0) - (0) 0( - a f (1 13) Moreover, because V j : ) = 0, the sideslip angle is zero,

pf (0) - - tan-'(v/;'/v/:') = o

(1 14) In view of this, and due to the symmetry of the gravitational and aerodynamic forces, the side force FJ0) and the roll and yaw moments, Mio) and A d : ' ) , are automatically zero. We assume that V(O) = 5000 in/s and consider a flight at a 25000 ft altitude, so that the speed of sound is 1016.1 ft/s and, hence, the Mach number is 5000/(1016.1 x 12) = 0.41. From Eq. (109), we have Vy' = 5000[cos O(O) 0 sin @ ' ) I T . Then, using Eqs. (44) in conjunction with Eqs. (45)-(48), Eqs.

( 1 11) yield Fie) = 2FE) - 695.3237 cos2 t 9 ( ' ) - 12973 sin O ( O ) + 2 1 4 9 . 5 3 5 8 ~ 5 : ' ) cos do) sin

+ 164919.9829 tan2 e(')) = 0

+ sin2 @ O ) ( 164224.6592 + 2149.53586p) tan

F,(') = - 2149.53586p) cos2 O ( ' ) + cos do)( 12973 - 202450.4921 sin e('))

+ sin2 0(0)(-2149.5358C5L0) - 239285.6777 tan@') - 36835.1856 tan3 O(O)) = 0

(1 15)

M r ) = - 27.92Ff) + (-16069.8929 - 553904.5798dp)) cos2 do)

-52019.1354 sin@') + (-2.0058 x lo6- 116986.55036p)) cos sin O ( ' )

+ sin2 19(~)(3.0203 x lo6 - 553904.57986p) + (-1.6634 x lo6

- 116986.55036p)) tan O ( ' ) + 3.0364 x lo6 tan2 + 342353.8625 tan3 @ O ) ) = 0

which can be solved for the pitch angle O ( ' ) , the engine thrust F') and the elevator angle S i ' ) . Solv- ing the nonlinear equations (115), we obtain = 0.0667 rad, FE) = 431.6465 lb and @) =

-0.2703 rad, so that the zero-order control vector is given by do) = [ F f ) 0 6;') OlT =

[431.6465 0 - 0.2703 OIT. Hence, the control force vector can be written in the matrix form 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 431.6465 0 0 0 0 2 0 143.4894 0 -0.2703 0 0 0 - 1905.7803 0 0 0 -2149.5358 0 0 803506.9168 0 - 153855.5426 -27.92 0 -561713.8624 0 0 53637.0429 0 529190.6764 = [0 0 0 0 0 0 824.5038 0 581.0801 0 139795.5057 0IT In the case of steady level flight, the aircraft experiences static deformations due to zero-order forces. Because static deformations are constant and deformations are first-order quantities, we denote all quantities involved by a subscript c and a superscript (1). Consistent with the zero-order results and using the first of Eqs. (50), we write where V(l) is the first-order forward velocity. Then, from Eqs. (59), the corresponding momenta are and, from Eqs. (50), we conclude that FP) = 0, MP) = 0 ( 1 19) e - K u i q u i c + Qui= = 0, -K$iqdric + Q ~ i c = 0, i = f, which represent algebraic equations to be solved for the first-order pitch angle, the first-order elevator incident angle, the first-order engine thrust and the static generalized displacements s i c and qli,ic. We note that FL') = [Fi;) 0 F,, (l) ] T , Mi1) = [0 Mi;) O I T , Q u i c and Q Q i c are all functions of quic and qqic ( i = f, w, e ) .

The stiffness matrices are as follows: 8687.0713 -3181.9271 0 0 -3181.9271 293633.2239 0 0 0 0 12883.4627 -4643.2903 0 0 -4643.2903 419984.1291 8891.2109 -28248.0864 0 0 -28248.0864 231517.6212 0 0 0 0 13097.4380 -36822.7531 0 0 -36822.7531 363078.6024

[ 1.0781 1.9339 ]

[ 1.5462 1.0134 ]

los 'Os' K$f = 1.9339 6.7838 K$f = 1.0134 9.0414 ( 120)

581.3939 -1488.0598 2.2974 3.6438 ]

R - K L =

] ' K ~ w = K & = [ 3.6438 13.5871 K u w - uw [ -1488.0598 13050.3375

1.5455 1.4924 406.2228 -696.8287 R - K L = K,", = K,", =

[ -696.8287 10685.4900 ' K G e - qe

1 1

5.0679 5.3834 2126.9677 -4303.4586

K,", = [ -4303.4586 52929.3173 ] ' K G e - 5.3834 31.9145

1 - [

where the superscripts F, A , R, L and V denote the fore part, aft part, right half, left half and verti- cal (stabilizer), respectively. Assuming that V(l) = - 5 in/s, so that Vk) = V(l) [cos O ( O ) 0 sin +V(o)[-O(l)sO(o) 0 I ~ ( ~ ) c O ( O ) ] ~ = [-4.6953 - 333.0272O(l) 0 - 4.7316 + 4988.89700(l)IT, we can solve Eqs. ( I 19) and obtain O ( l ) = - 0.00088 rad, F,!$) = -3.1169 lb, 6:') = -0.0047 rad S u f c - - [0 0 0.1206 - 0.0015ITin, q;fc = [0 OIT rad

S u f c A - - [0 0 0.0585 0.0038IT in, qif, = [0 OIT rad

q,,, R - - qtwc = [-3.4537 - 0.2918IT in, q&, = q$wc = [0.0021 - 0.0002IT rad

que, R - - qiec = [-0.0279 0.0006IT in, qGeC = q$,, - - [0.0024 - 0.0002]T rad

q : e c = [0 OlT in, q:ec = [0 OIT rad Assuming that the first-order state x(l) is measured from the static elastic displacements posi- tion, the first-order state equations for steady level flight can be written in the customary form, Eq.

(63), where A and B are constant coefficient matrices. This requires the mass matrices, defined by Eqs. (15). They can be computed as follows: Mll = mI, 131.1002 0 0 0 -4,2385 0 0 0 - 131.1002 0 M12 = S T = 0 4.2385 0 0 0 4.2385 0

+ [ 8 8 3'20026 ] qcfyz + [ -4.2385 0

4 u f z l + * - -3.2026 0 0 0 0 0

I F

0 -0.0807 -0.0128 0 0 -0.1252 0 0 0.0128 0 0 0.1252 0 0 0 0 0.0219 -0.0219 0 0 0 0 8 1 , . . . , 0 4.2385 -3.2026 0 -863.2814 0 182862.0586 4.7745 -37732.0967 F -863.2814 0 -21.0393 4.7745 566132.2704 81.1583 M22 = J = Quf Yl 0 -21.0393 0 -37732.0967 81.1583 704093.8389

I +

367.2116 0 0 -863.2814 0 7.3021 qUfy2+ 42.0786 0 qufzl+. . . (122)

7.3021 0 ] [ I C 1 4 0 0

I F

0 35.9911 0 6.9600 3.2162 -20.30651

3.2162 7.5428 -1.2425 QUez + 35.9911 0 12-8989 QUel

-20.3065 -1.2425 -0.5828 0 12.8989 0

1 ' - [

0 -4.0598 0

4.0183 qL2, + [ - 4 . y 0

4.0183 0 0 0 2.4862 -0.7090 -21.0393 7.3021 0 M& = 0 0 0 0 0 863.2814 -367.2116 0 68.4129 72.7595

+ [ 0 . 7 : -2.3655 0 0

0 0 0 ... ) M; = -344.8056 -402.1912 0 0 0 2.2363 -1.5318 0 -1.5318 8.4942 0 M A = 0 -1.5310 7.9257 17.3669 10.8317 M\J - 34.1002 32.0129

Mk3 = [ 10.8317

32.5338 ’ 88 - [ 32.0129 65.4667

1 1

Before we can compute A and B, we must still generate the damping matrices Cui and Cqi (i = f, w, e ) . To this end, we assume that the damping functions G~ and clLi (i = f , w, e) are all constant, so that the damping matrices, Eqs. (19), reduce to (123) Cui = cuiKui, Cqi = cwiK*i, i = f, w, e Equations (123) state that the damping matrices are proportional to the stiffness matrices, which permits us to write the relations37 cui = 2c/A;{’, cVi = 2</AQi 112 , i = f, W , e ( 124) where 5 is a structural damping factor and gi2 and Ai? (i = f , w , e ) are the lowest natural

frequencies of the respective components. We assume that < = 0.03. Moreover, we obtain the

component natural frequencies by solving the eigenvalue problems det[Kcf - AffM&] = 0, det[K,’lf - AtfMi] = 0, det[K,R, - AfwMz] = 0 (1 25) det[K,L, - A$~kf&] = 0, . . . , det[K:e - AteM&] = 0, det[K:e - AzeM&] = 0 with the results

& = 59.1046 rad/s, & = 54.1788 rad/s, a = a = 36.2911 rad/s

a = a = 72.1276 rad/s, a = 131.9557 rad/s,

& = 221.8479 rad/s, & = 61.9940 rad/s, @ = & = 253.9293 rad/s

,/Age = ,/ASe = 294.1938 rad/s,/A:e = 384.8014 rad/s Hence, using Eqs. (120), (123), (124) and (126), the damping matrices are 8.8187 -3.2301 0 298.0817 0

c ; f = [ -3.i301

0 13.0786 -4.7136 0 0 -4.7136 426.3467 9.8465 -31.2832 0 0 -31.2832 256.3929 0 0 0 14.5047 -40.7792 0 0 -40.7792 402.0894 ( 127) 418 18.3661 27408.5814 104342.4068 187172.4172

[ 187172.4172 656562.6583 ‘ $ f = [ 27408.5814 244530.1083 1’ ‘ $ f

5428.3762 8609.7859 0.9612 -2.4602 C~ =

C R =c= =

uw

uw [ -2.4602 21.5761 ’ iW 8609.7859 32104.4936

1 1

315.1931 304.3728 0.3379 -0.5797

C R = c= =

] CR =c;e= [

304.3728 2056.9682 ue ue [ -0.5797 8.8888 ’ *e

0.9671 -1.9568 CV = 790.2110 839.4058 ]

c v =

ue [ -1.9568 24.06691’ @e [ 839.4058 4976.2592

Then. the coefficient matrices in Eq. (63) can be shown to be 0 0 0 0 0 0 ... 0 0 0 0 0 -333.0272 0 5000 ... -0.0014 -0.0001 0 0 0 0 -5000 0 ... 0 0 0 0 0 0 0 0 ... 0 0 A = . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

0 0 0 36.3883 19.0770 0 ... -0.0212 -0.0020 0 0 0 40.8390 41.4025 0 ... 0.1414 0.0257 0 0 0 455.5429 0 0 . .. -31.2223 -9.6060 0 0 0 538.2891 0 0 ... 63.7274 -136.8998 0 0 0 0 0 0 0 0 ........................

0 0 0 0 B = 2 0 143.4894 0 0 0 0 - 1905.7803 ........................

0 0 0 -3097.7273 Moreover, the feedback control vector is in which it is assumed that the right and the left ailerons rotate by angles of the same magnitude 6?)

but of opposite sense and that the right and left elevators rotate by angles of the same magnitude and sense.

Choosing the weighting matrices Q and R in the performance index, Eq. (66), as follows 1 0 0 0 0 0 0 ...

0 1 0 0 0 0 0 ... O l 0 0 0 1 0 0 0 0 ... 0 0 0 0 1 0 0 0 ... 0 0 0 0 0 1 0 0 ... 0 R =

1 ; 5 x ; 0 8 : : (130)

Q =

108 0 0 0 0 0 0 1 0 ... 0 0 0 108 0 0 0 0 0 0 0 ... 0 .........................

1 0 0 0 0 0 0 0 I . . O J solving the steady-state Riccati equation, Eq. (73), and using Eq. (68), we obtain the gain matrix T 0.2994 0 0 0 0 0 0 0 -0.1016 0 0.0001 0 0.0003 0.5222 0 0.1024 -1112.0150 0 -1.8943 0 -0.0187 1.9350 0 0.5103 G = ....................................

-0.0006 0 0 0 0.0001 0 0 0 0 0 0 0 0 0 0 0 Then, solving the closed-loop eigenvalue problem, Eq. ( 7 3 , we obtain the closed-loop eigenvalues A1,2 = -0.1038 f 0.0868 i, A3,4 = -0.1174 f 0.2033 i, As = -0.2349, x6,; = -0.2918 f 0.3064 i , . . . , &9,70 = -40.5617 f 599.4532 i , (132) x71,72 = -40.7342 f 604.4685 i, X73,74 = -79.0030 f 791.6011i, A75.76 = -211.4842 f 1063.7436 i Clearly, all the eigenvalues are real and negative or complex with negative real part, so that the closed-loop first-order system is asymptotically stable. The implication is that any disturbances from the steady level flight are driven to zero. This is in contrast with the open-loop eigenvalues, the eigenvalues of A, the first four of which are zero and the fifth is real and positive.

Finally we consider the response of a closed-loop system to a gust acting on the wing and having the linearly distributed form R R f , (zw, t ) = [0 0 - 0.5(3 + z : / L , ) sinxih(1 - t)]’, 0 < x : < L , (133) f,(z,,t) L L = [0 0 - 0.5(3 - zf;,/L,) sinrTTt((1 - t)]’, 0 < xf;, < L , i where sin T&( 1 - t ) represents a half-sine pulse, in which M ( 1 - t ) is a rectangular function of unit amplitude and unit length. Inserting Eq. (133) into Eqs. (54), we obtain the generalized force components of the disturbance vector Fext entering into E $ (74), which can be integrated to obtain the system response. Figures 3-5 show the response for the rigid body variables and a selected number of elastic variables, and Fig. 6 shows the control inputs.

y/”’ [ i n ]

-10 -20 -30 -40 -50 -60 Figure 3. Rigid Body Displacements 2 . 5 0 - r - .

-2.5 -5 -7.5 -10 y 6 g [rad/s] 0.003' I\ 0.002 0 . 0 1 0 . 0 0 1 A / \ 0 O'J

v Yo "lYY 0 - 5

-0.001 t [SI -0.01 -0.002 -0.003 -0.02 -0.004 ' J Figure 4. Rigid Body Velocities O.OOO4 0.0oO2 -0.0002 -0.0004 94 [in1 0.00002 0.025 4 6 8 1 0 12 1 4 0.02 I -0.oooo2 t [SI 0.015 -0.00004 0 . 0 1 -0.00006 0.005 -0.oooO8 ' 0 - -0.o001 g / 4 6 8 1 0 12 1 4 t [SI Figure 5. Generalized Elastic Displacements 4 1 G!)[rad] 0.0125 0.0 1 0.0075 0.005 0.0025

0 -’

15 20 25 Right L e f t Fore Aft Right Left Horizontal Horizontal Vertical Stabilizer Fuselage Fuselage Wing Wing Stabilizer Stabilizer 111.21,36.25 111.92.35.0 120.35.58.23 120.35, -58.23 50.54.32.0 50.54. -32.0 50.0.42.0 PI : 1 1 ~ ~ 1

Velocities vj:) and vji) vjt) and vji) v g VA;’ v ; : ) vi; ) v ; : )

111.21, -36.25 111.92, -35.0 120.35, -28.0 120.35, 28.0 46.6, -20.0 46.6.20.0 41.0, -28.0 Pz 1 2 , y 2 Velocities VU) vi;) v w z - ( 1 ) vg ) vi;) v ; ; fz f z 271.21.24.0 279.79.25.0 327.66,25.07 327.66, -25.07 127.0, 19.0 127.0. -19.0 112.0. 27.0 P 3 : 1 3 . ~ 3

Velocities vj:) and vji) vjt) and vjtl v w z - ( 1 ) v w z - ( 1 ) v ; : v : ; vd,‘ )

277.21, -24.0 279.79, -25.0 327.66, -10.0 327.66. 10.0 127.0. -12.0 127.0, 12.0 113.0, -16.0 P 4 : 1 4 ~ ~ 4 - ( 1 ) - (11 VU) V ( ’ ) Velocities VWZ v w z

f z f z v ; : ) vi; ) v ; ; )

To compute the matrix C, Eiq. (103), relating the output vector to the state vector, we use Eq. (100) and write 0 ... 0 0 0 0 0 ... 0 0 0 0 0 ... 0 0 0 0 0 ... 0 1411.5402 -1066.5377 0 0 ... 0 0 0 21145.5039 M v = 0 ... 0 -21145.5039 15977.2118 0 0 ... 0 0 0 0 .............................................

1 0 ... 0 0 0 0 0 0 0 0 ... 0 0 0 0 ... 0 0 0 0 0 ... 0 ... 4.2544 41.7049 -7.3013 0 ... 0 ... -402.3892 0 0 0 ... 0 -63.7322 - -624.7574 109.3773 0 ... 0 0 0 0 0 ... 0 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . f . . . . . . . . .

... 0 0 0 0 ... 0 Then, inserting Eq. (134) into Eq. (103), we obtain 0 ...

- 1 0 0 0 0 0 0 0 ........................... . . . . . . . . . . . . . . . . . . . . . . . . . .

0 0 0 0 0 1 0 0 0 0 0 0 0 0 0 0 0.6890 -0.5206 0.0026 0.0019 0 0 0 0 0 0 -1.8548 1.4015 4.1411 -3.1290 C = 0 0 0 0 0 0 1.8590 -1.4046 4.1328 -3.1227 ........................... . . . . . . . . . . . . . . . . . . . . . . . . . .

0 0 0 0 0 0 1.5244 -1.1518 -0.0066 0.0050 0 0 0 0 0 0 -2.4135 1.8236 -0.1611 0.1217 0 0 0 0 0 0 1.8189 - 1.3970 0.1101 -0.0832 0 0 ... 0 0 . . . . . . . . . . . . . . . . . . .....................

0 0 ... 0 0 -0.00003 0.0123 ... -0.0002 0 0.0008 0.0320 ... -0.o004 0 0.0011 -0.0373 . . . 0.0004 0 . . . . . . . . . , . . . . . . . . .....................

-0.0015 0.0871 ... -0.4412 -0.6612 0.0002 -0.0260 ... 2.1761 -1.4006 0.0009 -0.0385 . . . -1.4522 0.9932 We assume that the excitation noise v(t) and observation noise w(t) represent white noise pro- cesses with intensities V and W, respectively, and choose

V = diag[100000 100000 100000 100000 100000 100000 lo-’ lo-’

1 0 - ~ 1 0 - ~ 1 0 - ~ 1 0 - ~ 10-~ 1 0 - ~ 1 0 - ~ 10-9 10-9 10-9 10-9 1 0 - ~ 1 0 - ~ 1 0 - ~ 1 0 - ~ 1 0 - ~ 1 0 - ~ 1 0 - ~ 1 0 - ~ 1 0 - ~ 1 0 - ~ o . . . 0 1 (136) and so that, from Eq. (91), the steady-state Riccati equation

AQ + QAT + V - QCTW-lCQ = 0

yields . ..

Q = diag [ 0.0003 0.0003 0.0003 0.3162 0.3162 0.3162

(1 39)

0 0 ]

Hence, from Eq. (90), the optimal observer gain matrix is - 3.1623 x lo8 0 0 0 0 0 ... 0 0 3.1623 x lo8 0 -0.0003 0 0.0050 . .. 0 0 0 3.1623 x lo8 0 -0.0050 0 ... 0 0 -332.6945 0 316227.7660 0 0 ... 0 0 0 -4995.005 1 0 316227.7661 0 . . . 0 0 4995.0050 0 0 0 316227.7661 . . . 0 0.0018 -0.6163 0.0023 -0.0471 -0.0003 0.1358 . . . 0 - - -0.0013 0.4656 -0.0018 0.0348 0.0002 -0.1029 ... 0 3.4719 0.0042 4.5563 0.0001 -0.5724 0.0003 ... 0 -2.6233 -0.0032 -3.4428 -0.o001 0.4346 -0.0002 . .. 0 0.0016 -0.5422 0.0021 -0.0414 - 0.0003 0.1193 ... 0 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

-0.0014 0 0 0 0 ... 0 -0.0001 0 0 0 0 ... 0 .

( 140) Finally, solving the eigenvalue problem for A - KoC, we obtain the observer eigenvalues A1 = -0.0103, A2,3 = -0.7807 f 1.66182, X4,5 = -0.1874 f 2.21662, = -5.7278, x7,8 = -5.7237 f 36.97982, xg,lo = -3.7696 f 45.60802, x11,12 = -4.6367 f 62.63482, x13,14 = -4.2692 f 68.79292, . . . , (141) &3,64 = -40.5621 599.44942', &5,66 = -40.7343 f 604.46872, x67,68 = -79.0066 f 791.58212, &9,70 = -211.4807 f 1063.74322, A71,72,73 = -316228, A74,75,76 = -3.1623 X 10' The performance of the observer design can be demonstrated by simulating the response of the combined system, defined by Eq. (84), to an initial state and initial observer error. To this end, we choose the values x(l)(O) = [5 5 5 0.005 0.005 0.005 0.2 0.1 0.2 0.1 0.2 0.1 -0.2 -0.1 - 0 . 2 - 0.1 - 0.2 - 0.1 - 0.2 - 0.1 -0.2 - 0.1 -0.2 - 0.1 0.05 0.01 - 0.05 ( 142) - 0.01 0.05 0.01 - 0.05 -0.01 0.05 0.01 - 0.05 - 0.01 0.05 0.01 O . . . OlT e(0) = - 0.15d1)(0) Figures 7-9 show the response for a selected number of rigid body and elastic variables and their observer estimates. Figure 10 shows the control inputs as given by Eq. (82).

4 4 L 5 10 1 5 - 0.0 1 -: 15 20 25 0.005 -0.005 -0.005 -0.01 - -0.01 -0.015 -0.015 -0.02

- Actual - Estimated

Figure 7. Rigid Body Displacements - 4 . 2 100 -4.4 -4.6 -4.8 1 0 15 2 0 2 5 -5.2 -50 -5.4 0 . 2 0.5 0.1 0 0 -0.5 -0.1 -0.2 -1 -0.3 -1.5

- Actual - Estimated

Figure 8. Rigid Body Velocities -2 -4 0.04 0.04 0.02 0.02 -0.02 -0.02 -0.04 -0.04

-0.06 [I

I '

- Actual

- Estimated

Figure 9. Generalized Elastic Displacements &rad] ,Fill [ l b ] 0 . 2 0.005 -0.2 -0.4 -0.005 - 0 . 6 -0.01 - 0 . 8 0.001 ' 0.003

-

0 - 1

20 25 0.002 -0.001 b 0.001 -0.002 p r e 10. Control Inputs ii. Level steady turn maneuver We consider the case in which in the zero-order problem the aircraft flies at a constant velocity around a circular path of radius R in the horizontal X, Y-plane. In this case, it is convenient to refer the rigid body motions to a set of axes zlylzl obtained through a rotation $ ( O ) about 2, where $ O ) = R = constant. It is not difficult to see that axes zl, y1 and z1 represent a set of cylindrical (0) axes t , n and 2, where t is tangent to the circle and n is normal to it. Denoting by R, the velocity of O f in terms of cylindrical components, the kinematical relation corresponding to the first of Eqs. (43) can be written as where = is the matrix of direction cosines between tnZ and the fuselage body axes zf yf z f , obtained from C f , the first of Eqs. (2), by letting $ = 0 and replacing O and 4 by and $('I, respectively.

Similarly, the second of Eqs. (43) can be written as

6 , ( O ) = [4 '(0) $0) p'] = [O 0 Q]T = ( E p ) ) - q )

(145) where Ep' can be obtained from E f , the second of Eqs. (2), by replacing O and 4 by O(O) and $ ( O ) , respectively. Equations (143) and (145) yield Vf (0) - - [ ~ ( 0 ) f. ~ ( 0 ) f y ~ ( 0 ) f z I T - - C(o)R(o) f f = RR[cO(') S O ( ~ ) S C $ ( ~ ) S O ( ~ ) C ~ ( ~ ) ] ~ = constant (146) and uy) = [ w : ) wfy (0) w,,] (0) T - - E f (0) 6 , * ( O ) = R[-sO(') cd0)s4(O) C O ( ~ ) S ~ ( ~ ) ] ~ = constant (147) so that, from Eqs. (49), we have P V f ( O ) - - mv?) + S ( 0 ) T W y ) = constant, p : ; = S(o)v(O) f + J ( O ) W ~ ) = constant (148) It follows that the equations of motion, the last two of Eqs. (43), reduce to -(o) (0) + F(0) = 0 -V(O) pvf ( O ) - Wf (0) pWf (0) +M(O) = 0 ? ( 1 49) -Wf P V f which are independent of time.

To determine the parameters defining the steady level turn maneuver, we choose the turn radius and control R and angular velocity R and solve Eqs. (149) for the bank angle qh(O), pitch angle (0) ~ ( 0 ) a(0) (0) T S , ] . Hence, assuming the values R = 1.5mi = 95037in and vector do) = [Fe a e R = 0.0526 rad/s, so that RR = 5000 in/s, and using Eqs. (146) and (147), we have (1 50) We consider a flight at a 25000 ft altitude, so that the speed of sound is 1016.1 ft/s and, hence, the Mach number for the flight is 5000/(1016.1 x 12) = 0.41. Inserting Eqs. (150) into Eqs. (149) in conjunction with Eqs. (148) and solving the resulting transcendental equations, we obtain O(O) = 0.0986 rad, 4 ( O ) = 0.6160 rad (151)

Ff) = 468.7429 lb, SF) = -0.0028 rad, Sp) = -0.3233 rad, a?) = -0.3445 rad

Then, the zero-order control force vector can be written as 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 178.2611 77.3023 844.5606 2 1.0370 0 -1916.1116 654.9960 0 0 -2875.1495 0 692.9931 0 -13.7394 1.4415 - 154689.6019 50641.1455 0 1.0746 x lo6 -6240.6944 171482.9784 -27.92 - 164.5690 - 750586.5476 -182061.4664 0 65027.5713 11.7592 532059.4480 As in the case of steady level flight, the aircraft experiences static deformations in the steady level turn maneuver as well. Denoting the corresponding constant quantities by the superscript (1) and subscript c, using the first of Eqs. (50) and letting V(') be the cylindrical coordinates L ( 1 ) counterpart of Rf , we can write

(1) - c(0) v(1) - pya)] =

= 0 Vfc - f [ f (153) sUic = 0, sqic = 0, i = f, w , e The momenta are all constant and can be expressed in'terms of the zero-order and first-order (59).

velocities and the static deformations using Eqs.

A constant solution of the first-order equations, Eqs. (50), can be obtained by letting the left sides be equal to zero and solving for &), e:'), u b . " = [F$; 6:) 6::' 6si)]T and the static deformations quic and qQic ( i = f, w, e). Assuming that V(') = [-5 0 OIT in/s and using the first of Eqs. (153), we obtain -4.9754 - 492.1217eP

vk) =

-0.2843 + 2874.81219L1) + 486.4652&'

(1 54)

[ -0.4017 + 4061.1909&) - 344.35624;''

so that Eqs. (50) yield q$p) = - 0.0011 rad, 0:') = -0.0013 rad

UP) = [-4.0727 - 0.0002 - 0.0063 O.OOOIIT

S u f c - - [-0.0005 0.0001 0.1458 - 0.00191T, qcfc = [-0.0001 OIT Sufc A - - [-0.0408 - 0.0008 0.0827 0.00391T, q$fc = [0.0003 - O.OOOIIT - [-4.1523 - 0.35111T, qgwc = [0.0025 - 0.00031T q u w c - (155) q u w c - - [-4.1765 - 0.35261T, qGWc = [-0.0025 0.0003]T que, = [-0.0473 O.OOO1]T, q$ec = [0.0029 - 0.0002]T que, - - [-0.0489 0.000081T, qGeC = [-0.0029 0.00021T qYec = [0.0032 0.00141T, qTec = [0.0017 - O.OOO1]T The mass matrices, Eqs. (19, for the steady turn case are as follows: 0 0 -4.2385 130.4114 0.1477 0

0 0 k f 1 2 = sT = [ -130.4114 0

4,s Yl

0 4.2385 0 0 O I F -0.1477 0

0 4.2385 0

+ [ x x QffY2 + [ 0

-3.2026 0 0 -4.2385 0 0 0.1095 0.0173 0 0 -0.1252 0 0 0 0.1252 0 0 0 0 0 0 + [ 0 8 O'Oilg] qL2, MA = [ 4.20385 -3.2026 0 x 1, ... , -0.0219 0 0 0 0 4.2385 -3.2026 182795.8456 - 18.6695 -37747.5473 -0.0025 -863.2814

-18.6695 566091.8052 74.4042 ] + [ -86%2814 -21.1022 quful 0

-37747.5473 74.4042 704068.0911 -21.1022 -0.0025

O I F

0.0010 367.2116 0 42.2044 0 -863.2814 0 42.2044 0.0012 qufzl + . . .

7.3208 0.0010 -863.2814 0.0012 0

I F ( 156)

6.9557 3.2162 -20.3065 ] -0.0075 35.9911

12.8965 quel 3.2162 7.5394 -1.2456 35.9911 0

que2 +

-20.3065 -1.2456 -0.5837 0 12.8965 -0.0075

O I v

[

0.0022 -4.0598 0 - 4 . y 0 4.0187 0.0022

+ [

0 0 2.4362 -0.7090 -21.1022 7.3208 -0.0012 0 0 863.2814 -367.2116 68.3900 72.7325 0.7090 -2.3655 0 0 0 0

x 0 0 0 -344.8056 -402.1912

+ [

2.4862 -0.7090 0 -0.7090 2.3655 0 0 -0.7090 2.3655 34.1002 32.0129 17.3669 10.8317

Mk= [ 10.8317 32.5338 1 7 M L = [ 32.0129 65.4667

Next, we consider the time-varying part of the first-order problem, Eqs. (50) with all quantities measured from the constant static solution. To this end, we use Eq. (63) in which the coefficient - 0 0 0 0 0 0 ... 0.0004 -0.0001 0 -596.0114 0 0 0 0 5000 ... -0.0008 -0.0012 - 0 .ooo 1 0 -5000 0 0 0 0 ... 0.0010 -0.0008 0 0 0.0529 0 ... 0 0 0 0 0 0 -0.0524 0 0 0 0 0 ... 0 0 0 A = , . . . . . . . . . . . . . . . . . . . . . . . . . . . ............................. . . . . . . .

174.902 1 42.3914 0 . .. -76.1004 0.1415 0.0257 0 0 0 0 0 0 370.8328 -25.9627 0 ... 0.0604 -31.2428 -9.6184 0 0 0 438.1920 -30.6787 0 ... 0.1026 63.7171 -136.9109

' (157)

0 0 0

r o

. . . . . . . .......................... . . . . . . . . . .

0 0 0 0 2 1.0326 176.1982 104.3594 0 0 0 - 1901.4388 0 -10.2357 -2143.1016 0 B = 0 800910.394 1 1.8988 -153505.0524 -27.92 - 116.9527 -561836.0209 - 8425.04 11 0 64654.3829 11.7819 527985.1550 ...........................................

0 0 0 -7461.0801 0 0 0 -3090.6705 and in which the feedback control vector has the form u(l) = [F;) 6 : ) 6;') 6?)IT. In the case in which u(') = 0 and F,, = 0, the state equations admit an exponential solution yielding an eigenvalue problem. Solving the eigenvalue problem, we conclude that the system is unstable, with four eigenvalues being equal to zero and one being real and positive. Using a linear quadratic regulator in conjunction with the weighting matrices 1 0 0 0 0 0 ... 0 0 1 0 0 0 0 ... 0 0 0 1 0 0 0 ... 0

[ 8 1

0 0 0 1 0 0 ... 0 , R = & = 0 0 0 0 1 0 ... 0 108 0 0 0 0 0 0 1 ... 0 0 0 108 .........................

0 0 0 0 0 0 ... 0 solving the corresponding steady-state Riccati equation, Eq. (73), and using Eq. (68), we obtain the gain matrix T 0.2974 0 0 0 -0.0163 0 -0.0001 0 -0.1062 0 0.0001 0 12.4530 0.5355 0.0531 0.1015 -1022.6944 -1.0086 -1.2660 -0.2651 -71.9234 1.3797 -1.4307 0.3965 G = (159) . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

-0.0006 0 0 0 0.0001 0 0 0 0 0 0 0 0 0 0 0 Then, solving the closed-loop eigenvalue problem, Eq. (75), we obtain the closed-loop eigenvalues I X1,2 = -0.1031 f 0.08702, A3 = -0.2342, X4,5 = -0.1194 f 0.21722, X6,7 = -0.2998 f 0.30892, Xs,g = -0.7842 f 1.62262, . . . , ( 160) X69,70 = -40.5623 f 599.45582, X71,72 = -40.7352 f 604.47132 ~ A73,74 = -79.0053 f 791.6059i, X75,76 = -211.4850 f 1063.74652 Clearly, all the closed-loop eigenvalues are either real and negative or complex with negative real part, so that the closed-loop first-order system is asymptotically stable. Hence, any disturbances from the steady level turn maneuver will be driven to zero.

~ Finally, we compute the response of the closed-loop system to the gust given by Eqs. (133).

Figures 11-13 show a selected number of rigid body and elastic variables, and Fig. 14 shows the control inputs.

-5 -7.5 -30 -10 -12.5 -40 -15 -50 -17.5 0.001 0.005 -0.005 -0.001 -0.01 -0.002 I -0.015 -0.003 Figure 1 1. Rigid Body Displacements 5 1 -3 -4 -5 -5 -10 -6 -15

- I -

-0.01 -0.002 -0.02 ' -0.004 9; m1 0.025 0.1 0.02 0.015 1 7 15 20 25 0.01 -0.1 t [SI 0.005 -0.2 c 0 5 - 5 20 25 -0.005 t [SI 0.00006 0.004 0.00004 0.002 0 . m 2 -o,00002 -0.002 -0.00004 -0.004 -0.00006 -0.oooO8 Figure 13. Generalized Elastic Displacements 0.0125 0.01 0.0075 0.005 0.0025 S 3 * a d ] 0.0025 0.001 0.002 0.0015 0.001 -0.001 0.0005 -0.002 -0.0005 Figure 14. Control Inputs Next, we turn our attention to the observer. To this end, we must compute the matrix C, Eq.

(103), relating the output vector to the state vector. First, we use Eq. (100) and write 0.1369 0.1282 ...

[ ::: -0.1811

0 0.0219 ...

, 0.0166 0 ...

0 ... 0 -0.0219 0 . . . 0 1676.3711 -1275.2629 -1242.2537 . . .

Mr. = 0 ... 0 3.5685 -1.5887 21090.9413 . . .

15935.2420 0 ...

0 ... 0 -21090.3030 .................................................

1 0 ... 0 0 0 0 ...

-0.0030 -0.0054 0.0009 0 ... 0 0 0 0 ... 0 - 0.0004 0 ... 0 -0.OOO1 -0.0006 0.0001 27.2581 51.3897 -8.9291 0 ... 0 -401.1684 0.3647 0.1928 0 ... 0 -63.5214 -622.7178 109.2101 0 ... 0 ........................................

0 0 0 0 ... 0 1

Then, inserting Eq. (161) into Eq. (103), we obtain 0 0 0 ... 0 ... 0 0 ‘ 1 0 0 0 0 0 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

0 0 0 0 0 0 ... 0 0 0 0 0 0 1 0.6400 -0.4820 0.2213 -0.1651 0 ... 0 0 0 0 0 0 0 0 0 0 0 0 0 0 -2.2285 1.6958 5.8787 -4.4245 0 ... 0 0.0010 C = 2.2323 - 1.6989 2.3810 -1.8161 0 ... 0 0.0011 ‘ 0 0 0 0 0 0 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

. . . . . . . . . . . . . . . . . .

0.3731 -0.2782 0 ... 0 -0.0015 0 0 0 0 0 0 1.4370 -1.0831 1.7578 -0.5259 0.3938 0 . . . 0 0.0002 0 0 0 0 0 0 -2.3298 0.2414 -0.1811 0 ... 0 0.0009 ~ 0 0 0 0 0 0 1.8150 - 1.3705 0 ... 0 0 ....................... . . . . .

0 ... 0 0 0.0122 ... -0.0002 0 0.0374 . . . -0.0004 0 -0.0372 ... 0.0004 0 ....................... I . . . . . . .

0.0871 . . . -0.4412 -0.6612 -0.0259 ... 2.1761 - 1.4006 -0.0384 ... -1.4522 0.9932 We assume that the excitation noise v(t) and observation noise w(t) represent white noise pro- cesses with intensities V and W , respectively, and choose V=diag[lOOOOO 1OOOOO 1OoooO 1OOOOO 1OOOOO 1OOOOO lo-’ lo-’ lo-’ lo-’ lo-’ lo-’ lo-’ lo-’ lo-’ lo-’ lo-’ 1 0 - ~ lo-’ lo-’ lo-’ lo-’ lo-’ lo-’ lo-’ io-’ io-’ (163) io-’ 1 0 - ~ io-’ io-’ io-’ lo-’ lo-’ lo-’ lo-’ o ... 0 1 and so that, from Eq. (91), the steady-state Riccati equation

AQ + QAT + V - QCTW-lCQ = 0

yields 0.00032 0 0 0 0 0 ... 0 0 ... 0 0 0.00032 0 0 0 0 ... 0 0 ... 0 0.00032 0 0 0 ...

0 0 0 0 ... 0 0 0 0 0.3162 0 0 ... 0 0 ... 0 0 0 0.3162 0 . ..

0 0 0 0 ... 0 0 0 0 0 0.3162 ... 0 0 ... 0 0 0 0 0 0 0 ... 0.0001 0 ... 0 0 0 0 0 0 0 ... -0.0001 -0.0001 ... 0 & = . . . . . . . . . . . . . . . . . . . . . . . . ............................. . . . . . . . . . . . . .

. . . . . . . .

0 0 0 0 0 0 ... 0 0 ... 0 0 0 0 0 0 0 ... 0.0705 0.0288 ... 0 0 0 0 0 0 0 ... 0.0288 0.0133 ... 0 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . ............................. . . . . . . . . . . . . . . .

0 0 0 0 0 0 ... 0 0 ... 0 0 0 0 0 0 0 ... 0 0 ... 0 Hence, from Eq. (90), the optimal observer gain matrix is KO = QCTW-' 3.1623 x 10' 0 0 0 0 ... 0 0 3.1623 x 10' 0 -0.0006 0 0.0050 ... 0 0 0 3.1623 x lo8 0 0 ... 0 0.0050 0 -595.4161 0 3 16227.7660 0.0003 0 1 . . 0 0 0.0001 -4995.0050 316227.7661 0.0003 0.0026 ... 0 0 4995.0050 0 0 0.0026 316227.7661 ... 0 -0.0427 -0.7057 -0.4951 -0.0554 -0.0783 0.1107 ... 0 0.0326 0.5398 0.3786 0.0414 0.0594 -0.0839 ... 0 3.5739 -1.8039 4.2355 0.0150 -0.4665 -0.3316 ... 0 -2.6997 1.3723 -3.1937 -0.0107 0.3542 0.2518 ... 0 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

0.0003 -0.0009 0.0010 0 0 0 ... 0 0.0004 -0.0008 0.0010 0 0 0 ... 0 -0 .ooo 1 -0.0012 -0.0008 0 0 0 ... 0 0 -0.0001 0 0 0 0 ... 0 - ( 1 Finally, solving the eigenvalue problem for A - K,C, we obtain the observer eigenvalues A1 = -0.0111, A2,3 = -0.7858 f 1.62912, A4,5 -0.1935 f 2.22312, A6 = -5.7234, A7,8 = -5.7294 f 36.97852, Ag,lo = -3.7702 f 45.60962, A11,12 = -4.6312 f 62.62802,. . . , &3,64 = -40.5646 f 599.45142, &5,66 = -40.7353 f 604.47192, &7,68 = -79.0086 f 791.58922, &9,70 = -211.4819 f 1063.74632, A71,72 = -316227.8058 f 0.03462, A73 = -316227.8459, A74,75,76 = -3.1623 X I O s To check the performance of the observer just designed, we simulate the response of the com- bined system, Eq. (84) to the initial conditions ~ ( ~ ' ( 0 ) = [5 5 5 0.02 0.005 0.001 0.2 0.1 0.2 0.1 0.2 0.1 -0.2 -0.1 -0.2 - 0.1 - 0.2 - 0.1 -0.2 -0.1 -0.2 - 0.1 -0.2 - 0.1 0.05 0.01 - 0.05 - 0.01 0.05 0.01 - 0.05 - 0.01 (169) 0.05 0.01 -0.05 -0.01 0.05 0.01 0 ... 0 1 e(0) = 0.15d1)(0) Figures 15-17 show a selected number of rigid body and elastic variables and Fig. 18 shows the control inputs as given by Eq. (82).

-10 -20 -30 -40 0 . 0 2 0 . 0 1 0.005 0.015 0 . 0 1 -0.005 -0.01 0.005 -0.015 n 5 Y O - 25 -0.02

- Actual - Estimated

Figure 15. Rigid Body Displacements 0 100 -2 -4 -6 -50 -8 -10 0.2 0 . 5 0 . 1 0 0 -0.5 -0.1 -0.2 -1 -0.3 -1.5

- Actual - Estimated

Figure 16. Rigid Body Velocities 1.5 0.5 -0.5 -1 -1.5 -2 -4 -6

- Actual --- Estimated

Figure 17. Generalized Elastic Displacements F i l l [lb] 0.0075 0.005 0.0025

J

-0.005 -0.0025 -0.0075 -2 0.01 0.01 0.005 0.005 -0.005 -0.0 1 -0.005 Figure 18. Control Inputs 9. Conclusions A problem of increasing interest in aeronautics is the simulation of the flight of flexible aircraft on a computer. Note that on-line computer simulations can help expedite the design process, and possibly reduce the time devoted to flight testing. Moreover, the new class of autonomous aerial ve- hicles requires autopilots, which in turn requires real-time simulations. This points to the need for a in the treatment of flying flexible aircraft, one using the system concept to produce new paradigm a rigorous formulation of the dynamics and control problem and one using potent methodology permitting an efficient solution of the problem. This work develops such a paradigm in the form of an on-line formulation integrating pertinent material from the disciplines of analytical dynamics, structural dynamics, aerodynamics and controls. Moreover, the formulation is cast in a certain matrix form ideally suited for computer processing. A perturbation approach permits division of the problem into a nonlinear flight dynamics problem for maneuvering quasi-rigid aircraft and a linear “extended aeroservoelasticity” problem for the elastic deformations and small perturbations in the rigid body translations and rotations, where the solution of the first problem enters as an input into the second problem. As a result, there is a different extended aeroservoelasticity prob- lem corresponding to each aircraft maneuver. The controls for the flight dynamics problem are obtained by an inverse process, which amounts to prescribing a given maneuver and determining the controls permitting realization of the given maneuver. On the other hand, the elastic deforma- tions and perturbations in the rigid body motions characterizing the extended aeroservoelasticity problem are driven to zero by means of feedback controls designed using linear quadratic theory in conjunction with a stochastic observer, thus ensuring the stability of the aircraft maneuver. A numerical example presents a variety of time simulations of rigid body perturbations and elastic deformations about 1) a steady level flight and 2) a level steady turn maneuver.

The integration of the aerodynamics into the unified process is particularly challenging, and requires elaboration. For seamless integration, the aerodynamic forces must be referred to the same generally noninertial reference frame as that used for all other forces and they must be expressed in terms of variables compatible with the variables used throughout the entire formulation. Moreover, because the simulation of the system response on a computer is carried out in discrete time, the size of the time step is of vital importance. Indeed, on line and real time simulations require that the time step be quite small, the time step for the latter being of the order of a minute fraction of a second. The implication is that, to be ready to compute the state at the next sampling time, it is necessary to be able to compute the aerodynamic forces within the time step; the computation of the other forces within the same time step presents no problem. These two requirements are quite difficult to satisfy as most aerodynamic theories have been developed for purposes other than time response simulations, and the computation of the aerodynamic forces are notorious for consuming a great deal of time. In view of this, the development of a new aerodynamic method seems highly desirable. Such a method need not be unduly accurate, as robust feedback controls should be able to tolerate small errors in the aerodynamic forces. An aerodynamic theory satisfying the two requirements outlined above is strip theory. Even though strip theory may not be entirely satisfactory for describing the aerodynamic forces acting on whole aircraft, the method is often used in aircraft design. In fact it is being used by the same company that provided the data for the Numerical Example. The development of an aerodynamic technique suitable for on line, or real time response simulations is likely to require a great deal of time and effort. Although the use of a more suitable aerodynamic theory is desirable, it is not really necessary at this time. Indeed, at this time it is more important to demonstrate how the unified formulation works and how an eventual 5 8 aerodynamic theory is to be integrated into the overall process. Certainly, such a demonstration can provide valuable guidance in the development of an appropriate aerodynamic method to be used in conjunction with aircraft time response simulations, thus helping reduce the time and effort required for the development in question.

Finally, it should be pointed out that, with appropriate modifications, the present formulation is eminently suited for UAVs, and in particular for autonomous UAVs. The fact that it was applied here to an executive jet is due entirely to the ready availability of data from an actual flying aircraft.

It should also be pointed out that all the time response simulations presented here were carried out on a 1 GHz PC using MATHEMATICA. This is particularly important for autonomous UAVs, as the required onboard computer is likely to be much which must be controlled by autopilots, closer to a PC than to a multiprocessor supercomputer.

10. References 1. Collar, A. R., “The Expanding Domain of Aeroelasticity,” Journal o f the Royal Aeronautical Society, Vol. 50, 1946, p. 428.

2. J. B. Rea Co., “Aeroelasticity in Stability and Control,” Wright Air Development Centel; Technical Report 55-173, 1957.

3. Taylor, A. S., “The Present Status of Aircraft Stability Problems in the Aeroelastic Domain,” Journal o f the Royal Aeronautical Society, Vol. 63, 1959, p. 580.

4. Bisphinghoff, R. L. and Ashley, H., Principles o f Aeroelasticity (Ch. 9), John Wiley & Sons, Inc., New York, 1962. (Reprinted by Dover Publications, Inc., New York, 1975).

5. Milne, R. D., “Dynamics of the Deformable Aeroplane, Parts I and I I , ” Her Majesty’s Sta- tionery OfJice, Reports and Memoranda No. 3345, London, 1962.

6. Taylor, A. S. and Woodcock, D. L., “Mathematical Approaches to the Dynamics of De-

formable Aircraft, Parts I and n, “Her Majesty’s Stationery OfJice, Reports and Memoranda

No. 3776, London, 1971.

7. Dusto, A.R. et al., “A Method for Predicting the Stability Characteristics of an Elastic Air- plane, Vol. 1-FLEXSTAB Theoretical Description,” NASA CR-114712, October 1974.

8. Roger, K. L., “Airplane Math Modeling Methods for Active Control Design,” Structural As- pects o f Active Control, AGARD-CP-228, 1977, pp. 4-1 - 4-1 1.

9. Schwanz, R. C., “Consistency in Aircraft Structural and Flight Control Analysis,” Structural Aspects o f Active Control, AGARD-CP-228, 1977, pp. 5-1 - 5- 16.

10. Miller, G . D., Wykes, J.H. and Brosnan, M. J., “Rigid-Body Structural Mode Coupling on a Forward Swept Wing Aircraft,” Journal o f Aircrafr, Vol. 20, No. 8, 1983, pp. 696-702.

11. Weisshaar, T. A. and Zeiler, T. A., “Dynamic Stability of Flexible Forward Swept Wing I Aircraft,” Journal o f Aircrafr, Vol. 20, No. 12, 1983, pp. 1014-1020.

I 12. Rodden, W. P. and Love, J. R., “Equations of Motion of a Quasi-Steady Flight Vehicle Utiliz- ing Restrained Static Aeroelastic Characteristics,” Journal o f Aircraft, Vol. 22, No. 9, 1985, pp. 802-809.

I 13. Cerra, J. J.II, Calico, R. A. and Noll, T. E., “Modelling of Rigid-Body and Elastic Aircraft Dy- Proceedings o f the Atmospheric Flight Mechanics namics for Flight Control Development,” Conference, Williamsburg, VA, Aug. 18-20, 1986, pp. 386-395.

14. Arbuckle, P. D., Buttrill, C.S. and Zeiler, T. A., “A New Simulation Model Building Process , for Use in Dynamic Systems Integration Research,” Paper 87-2498, AIAA Flight Simulation Technologies Conference, Monterey, CA, August 17-19, 1987.

15. Buttrill, C. S., Zeiler, T. A. and Arbuckle, P.D., “Nonlinear Simulation of a Flexible Aircraft in Maneuvering Flight, Paper 87-2501- C e AIAA Flight Simulation Technologies Conference, Monterey, CA, Aug. 17-19, 1987.

16. Zeiler, T. A. and Butrill, C. S., “Dynamic Analysis of an Unrestrained, Rotating Structure Through Nonlinear Simulation,” Proceedings of the AIAALASME/ASCE/AHS29th Structures, Structural Dynamics and Materials Corzference, Williamsburg, VA, April 18-20, 1988, pp.

167- 1 77.

17. Waszak, M. R. and Schmidt, D. K., “Flight Dynamics of Aeroelastic Vehicles,” Journal o f Aircraft, Vol. 25, No. 6, 1988, pp. 563-571.

18. Sipcic, S. R. and Morino, L., “Dynamic Behavior of Fluttering Two-Dimensional Panels on an Airplane in Pull-up Maneuver,”AIAA Journal, Vol. 29, No. 8, 1991, pp. 1304-1312.

19. Waszak, M. R., Buttrill, C. S. and Schmidt, D. K., “Modeling and Model Simplification of Aeroelastic Vehicles: An Overview,” NASA TM 107691, September 1992.

20. Bilimoria, K. D. and Schmidt, D. K., “Integrated Development of the Equations of Motion for Elastic Hypersonic Flight Vehicles,” Journal o f Guidance, Control, and Dynamics, Vol.

18, NO. 1, 1995, pp. 73-81.

2 1. Olsen, J. J., “Unified Flight Mechanics and Aeroelasticity for Accelerating, Maneuvering, Flexible Aircraft,” RTO AVT Specialists ’ Meeting on “Structural Aspects o f Flexible Aircraft Control ”, Ottawa, Canada, Oct. 18-20, 1999, pp. 6- 1 - 6- 12.

22. Teufel, P., Hanel, M. and Well, K. H., “Integrated Flight Mechanics and Aeroelastic Mod- elling and Control of a Flexible Aircraft Considering Multidimensional Gust Input,” RTO AVT Specialists ’ Meeting on “Structural Aspects o f Flexible Aircraft Control ”, Ottawa, Canada, Oct. 18-20, 1999,21, pp. 18-1 - 18-9.

23. Konig, K. and Schuler, J., “Integral Control of Large Flexible Aircraft,” RTO AVTSpeciulists ’ Meeting on “Structural Aspects of Flexible Aircrajl Control ”, Ottawa, Canada, Oct. 18-20, 1999,21, pp. 19-1 - 19-12.

24. Samareh, J. A. and Bhatia, K. G., “A Unified Approach to Modeling Multidisciplinary In- teraction s ,” Paper 2000-4 704, 8th AI M N A SA/USA F/ISSMO Mu1 tidisciplina ry A nalysis and Optimization Symposium, Long Beach, CA, September 6-8,2000.

25. Schmidt, D. K. and Raney, D. L., “Modeling and Simulation of Flexible Flight Vehicles,” Journal of Guidance, Control, and Dynamics, Vol. 24, No. 3,2001, pp. 539-546.

26. Meirovitch, L., Methods of Analytical Dynamics, McGraw-Hill Book Co., New York, 1970.

(A Dover Publications edition is to appear in 2003).

27. Meirovitch, L. and Nelson, H. D., “On the High-Spin Motion of a Satellite Containing Elastic Parts,” Journal of Spacecraft and Rockets, Vol. 3, No. 11, 1966, pp. 1597-1602.

28. Meirovitch, L., “Hybrid State Equations of Motion for Flexible Bodies in Terms of Quasi- Coordinates,” Journal of Guidance, Control, and Dynamics, Vol. 14, No. 5, 1991, pp. 1008- 1013.

29. Platus, D. H., “Aeorelastic Stability of Slender, Spinning Missiles,” Journal of Guidance, Control, and Dynamics, Vol. 15, No. 1, 1992, pp. 144-151.

30. Meirovitch, L., “Derivation of Equations for Flexible Multibody Systems in Terms of Quasi- Coordinates from the Extended Hamilton’s Principle,” Journal o f Shock and Vibration, Vol.

1, Issue 2, 1993/94, pp. 107-119.

3 1. Meirovitch, L. and Stemple, T. J., “Hybrid Equations of Motion for Flexible Multibody Sys- tems Using Quasi-Coordinates,” Journal o f Guidance, Control, and Dynamics, Vol. 18, No.

4, 1995, pp. 678-688.

32. Meirovitch, L., “A Unified Theory for the Flight Dynamics and Aeroelasticity of Whole Air- craft,” Proceedings o f the Eleventh Symposium on Structural Dynamics and Control, Blacks- burg, VA, May 12-14, 1997, pp. 461-468.

33. Nydick, I. and Friedmann, P. P., “Aeroelastic Analysis of a Trimmed Generic Hypersonic Ve- hicle,” Proceedings o f the CEAS/AIAA/ICASE/NASA International Forum on Aeroelasticity, and Structural Dynamics, Langley Research Center, Hampton, VA, 1999, pp. 777-809.

34. Meirovitch, L. and Tuzcu, I., “Multidisciplinary Approach to the Modeling of Flexible Air- craft,” Proceedings o f the CEAS/AIAA/AIAE International Forum on Aeroelasticity and Struc- tural Dynamics, Madrid, Spain, June 5-7,2001, pp. 435-448.

I 35. Fornasier, L., Rieger, H. Tremel, U. and van der Weide, E., “Time-Dependent Aeroelastic I Simulation of Rapid Manoeuvering Aircraft,” Paper AIAA-2002-0949,40rh Aerospace Sci- ences Meeting and Exhibit, Reno, NV, January 14-17,2002.

~ 36. Etkin, B., Dynamics o f Flight (2nd Ed.), John Wiley & Sons, New York, 1982.

37. Meirovitch, L., Principles and Techniques of Vibrations, Prentice-Hall, Inc., Upper Saddle River, NJ, 1997.

~ 38. Dowell, E. H. (Editor), A Modern Course in Aeroelasticity (3rd Ed.), Kluwer Academic Pub- lishers, DordrechtBostonLondon, 1995.

39. Meirovitch, L. and Seitz, T. J., “Structural Modeling for Optimization of Low-Aspect-Ratio Wings, ” Journal o f Aircraft, Vol. 32, No. 5, 1995, pp. 1114-1123.

40. Meirovitch, L., Dynamics and Control of Structures, John Wiley & Sons, New York, 1990.

41. Biezad, D. J., Integrated Navigation and Guidance Systems, AIAA Education Series, Reston, VA, 1999.

Appendix - Numerical Values of the Aircraft Parameters

The values of the pertinent aircraft parameters were provided by an aircraft manufacturer in lumped form, and the current formulation assumes distributed parameters. This presents no problem, how- ever, as lumped parameters can always be treated as distributed by means of spatial Dirac delta

functions.37

In the first place, we consider the aircraft geometry. To this end, we regard the fuselage as consisting of two cantilever beams, a fore part and an aft part, with the origin of both sets of body axes at point 0, and with axes x ? and x ; ' collinear (Fig. 2), where superscripts denote the fore part and aft part. The wing is also divided into two parts, the right half-wing and the left half- wing, both with the origin of the respective body axes at 0; and 0; and with the longitudinal axes coinciding with the respective elastic axes. Moreover, the empenage consists of the horizontal stabilizer, divided into a right half and a left half, both with the origin of the respective body axes at 0: and O:, and a vertical stabilizer with the origin of the body axes at 0:. The radius vectors

from 0, to the corresponding origins are

rfw R = rFw = [-5.04 0 38.33]* in

( A I )

r7e = rFe = [-244.75 0 - 43.13IT in, rye = [-238.97 0 - 24.01]* in

The formulation calls for matrices of direction cosines between the various component body axes and the fuselage body axes, and in particular the body axes of the fore part of the fuselage denoted by xfyfzf. The component body axes can be obtained from x:fyfzf through a sequence of rota- tions. For example, axes xEyEzf for the right half-wing can be obtained through a rotation 31 about xfyfzf to an intermediate set of axes xLyk,zL and a rotation 7 2 about yL to xwyu,zu,, R R R where in the case of the wing y 2 is known as the dihedral angle. Table 2 gives the rotation angles for the individual components.

Table 2 - Rotation Angles for Component Body Axes

Component Rotations Body Axes

71 I Y2 I Y3

x A A A

I f Yf Z f I I I I

Righi half-wing 900 40 0 Righi half-wing 900 40 0 . . .

L L L ~ W Y U J Z W - L e f t half-wing -goo 4 O 0 L e f t half-wing -goo 4 O 0

- - -

R R R Xe Y e ze Vertical Stabilizer 118.74 -900 Using analogies with the first of Eqs. (2), the various matrices of direction cosines are as follows: 0 0.9976 -0.0698 0 0 0 0.0698 0.9976 0 -0.9976 -0.0698 -0.0750 0.9849 -0.1560 -0.0741 0.0117 0 -0.0698 0.9976 0.1564 0.9877 -0.0750 -0.9849 -0.1560 -0.4808 0 -0.8768 0 -0.1564 0.9877 0 The inertia properties of the aircraft model are given in lumped form. To this end, the flexible components are divided into certain numbers of bays and the mass corresponding to each bay is lumped at the mass center of the bay. Moreover, the manner in which the mass is distributed over the associated cross-sectional area is represented by mass moments and mass products of inertia about axes with the origin at the mass center and parallel to the body axes of the respective com- ponent. Table 3 lists the coordinates of the mass centers, the mass values and the mass moments and mass products of inertia, in which the symmetry of the inertia matrices is implied. The masses have units lb - s2/in and the mass moments and mass products of inertia have units lb s2 - in.

The cantilever beams lengths are LF = 295.86 in, L;? = 279.79 in, LE = L i = 328.83 in, LF = L t = 127.46 in, Ly = 113.48 in.

Table 3 - Inertia Properties

Fore Part of the Fuselage -

JL f Jry JLf J : y f JA f JLf

No. z ; yfF 2 ;

A

- 5 14.6542 230.8500 313.3728 -0.1010 -7.9229 0 7.41 0 5.20 0.7237 I 1.4603 0 1.98 0.2095 206.4986 115.1881 97.2035 -0.3495 2 24.74 0 -1.9160 0 0 - 1.54 0.3014 333.7084 182.4735 160.9599 0.1424 3 38.95 -0.6913 0 0 3.25 0.7558 354.9372 129.905 1 25 1.8536 1.2350 4 48.49 3.0682 0 0 8.60 0.3555 35 1 S220 18 I .4974 177.6291 1.3775 5 55.08 2.9413 0 0 0.5562 559.965 I 274.9544 306.9516 -2.2111 6 69.88 7.88 7.32 0.3167 298.6947 182.3156 124.5506 -0.3573 -2.2086 0 7 84.47 0 0.7053 418.5278 480.2490 0.1657 5.3648 0 8 101.30 0 5.29 501.1153 0.1546 137.191 1 84.7780 54.2436 -0.1605 0.3754 0 9 11 1.10 0 5.40 1.2842 0 1.62 1.0289 934.8339 556.8193 535.0830 -10.4733 10 129.67 0 17.3295 0 6.60 1.2664 919.97 19 477.2896 595.6675 3.1044 11 160.23 0 0 -0.30 1.0170 673.8353 340.958 1 437.3046 2.3821 -20.5815 0 12 189.38 1.0789 354.2640 386.3648 -0.4686 -1.1859 0 13 204.35 0 1.98 623.35 12 1.16 1.0497 394.733 1 245.32 10 192.5843 -0.8907 9.7354 0 14 223.81 0 0.2496 86.3807 41.6213 48.6587 -0.5308 0.0829 0 15 237.7 1 0 6.55 6.57 0.7234 135.0187 96.6805 69.4034 0.9425 -9.6525 0 16 251.30 0 7.6562 0 14.02 0.4702 65.4238 49.4018 57.8918 1.1289 17 266.50 0 0 17.71 0.0979 14.8309 8.2310 10.2739 -0.0104 -0.0647 0 18 28 1.47 o nio4 0 0.0640 3.0164 2.1283 2.1801 -0.0673 19 292.08 0 17.66 - Aft Part of the Fuselage -

J,", f J,", J,", f Jlc", f J,", f J k

No. X ; ' y j zj

m;3

- 160.7476 76.2673 84.8272 -0.0052 0.1320 0 1 1.67 0 0.81 0.1434 0.7172 8.8913 0 17.72 0.6 147 2 1 14.0956 253.7152 1887.3572 0 2 10.66 1.4282 3571.1313 848.8 1 3 1 2879.2026 7.9514 -102.4285 0 3 22.18 0 19.00 0 2.00 0.1 183 133.0820 72.5725 60.7063 0.0233 0.0777 0 4 40.41 0 8.04 0.5808 354.6912 22 1.2985 163.5387 -0.1709 -4.8263 0 5 46.82 0 9.5 1 0.2506 259.3829 186.0932 84.9230 0.4246 -0.7664 0 6 66.1 1 0.7843 470.6690 478.1052 8.3087 -14.4348 0 7 83.74 0 8.39 900.7860 222.4144 138.4934 -3,6275 2.4908 0 8 98.2 1 0 3.87 0.5821 317.6683 0.5870 217.1299 530.1 143 -0.3728 I 1 S375 0 9 117.50 0 -8.33 634.89 13 0.6048 3 17.3757 146.8824 2 12.2078 0.5593 -4.3913 0 I O 142.62 0 -6.45 0.4945 168.9475 127.7664 71.1511 3.7414 -6.0380 0 11 156.28 0 -10.77 0.1 165 62.1588 48.9953 21.8166 0.4609 1.1496 0 12 173.22 0 -11.04 0.2131 49.4821 37.8 178 19.5795 0.5644 1.3852 0 189.26 0 -12.24 28.0876 15.0225 -0.1605 -3.7207 0 14 203.4 1 0 -6.82 0.1220 36.0701 37.9395 11.4727 -0.0129 3.5368 0 218.31 0 -19.66 0.0880 40.5364 5.9785 0.0078 -0.8596 0 238.30 0 -10.17 0.0585 14.7040 11 SO64 6.2555 3.7440 0.1036 -0.4220 0 253.16 0 -14.94 0.0580 6.649 I 0.8415 -0.0052 -0.1087 0 278.12 0 -17.35 0.0098 0.4143 0.9580 -

Table 3 - Continued

Wing R R R R R No.

X C YUP 2.w m W J X X W

J,"r 'u' J,",w J E W Jll", w J x z w

- 1 17.34 -5.67 -0.41 0.7932 592.7715 135.2155 675.5300 -11.4074 -20.6970 23.8870 2 41.44 -0.36 0.64 0.3380 259.1733 13.6761 253.5778 2.9790 0.1402 14.4486 3 55.80 7.21 0.72 0.548 1 203.2508 22.5351 202.0929 -2.4412 -0.1781 16.5185 4 68.45 0.92 -1.55 0.2900 145.9962 10.7568 144.9995 1.1106 0.4620 4.8288 5 80.75 4.29 -0.63 109.3290 0.3517 19.5652 117.3550 0.4698 -1.1316 4.4349 6 95.1 1 2.34 -0.89 87.6687 0.1960 7.3973 86.707 1 -3.5075 -0.1369 3.391 1 7 109.2 1 5.70 -0.49 0.1405 68.3377 5.1512 67.9640 0.3682 -0.0707 3.2896 8 123.58 10.93 -0.20 0.2014 98.035 1 8.6246 100.8629 - 5.1787 -0.2981 4.7601 9 141.58 6.17 -0.36 0.1441 62.7718 6.3306 64.6489 -0.3631 -0.0224 2.5948 10 159.69 5.62 -0.59 0.1550 58.9409 7.3054 62.3486 -2.4040 -0.0168 1.5852 1 1 180.67 2.46 -0.60 0.1481 57.806 1 6.8510 61.632 1 -1.7038 0.0205 1.76 13 12 199.99 5.62 -0.49 0.1504 44.7968 8.9780 5 1.3853 -4.0862 -0.1240 1.2184 13 223.12 3.28 -0.46 0.1 146 26.24 12 4.9868 29.6906 0.2463 0.0071 0.5415 14 242.88 1.87 -0.49 0.0825 16.9014 4.1 168 19.9998 -0.2915 -0.0464 0.3456 15 260.53 -0.02 -0.71 0.0700 1 1.9896 2.0960 13.3314 0.1592 0.068 1 0.32 13 16 278.05 1.04 -0.46 0.0628 9.2044 1.9160 10.6262 0.1316 0.0019 0.0663 295.21 -0.06 -0.38 0.047 1 6.9537 1.6338 8.2733 0.2092 -0.0066 -0.0035 18 31 1.97 -3.95 -0.38 0.0334 4.6844 0.6266 5.1908 0. I694 -0.0063 -0.0284 19 328.86 0.04 -0.40 0.0287 2.3874 0.7263 3.0305 -0.0109 -0.0045 -0.0293 - Horizontal Stabilizer R R R R R R X e Y e z e m,R J x x e JyRye J E e JEy e J x z e J y z e 3.4154 -2.0199 4.8360 0.1298 152.7019 24.5226 24.9826 9.4614 -20.3966 0.6193 14.6498 4.1648 -0.0644 0.06 18 13.3528 2.0288 14.9172 I .0676 - 0.0060 0.0244 32.5615 6.2772 -0.0954 0.0625 12.1552 1.9324 13.7307 0.2775 0.01 12 0.0409 50.7350 3.5463 -0.1694 0.0574 8.3296 2.7220 10.7676 0.5878 0.06 I2 -0.0078 66.0100 3.8814 -0.0389 0.0273 3.6287 0.3475 3.8652 0.2049 - 0.0006 0.0 143 76.5595 2.0048 0.0296 0.0173 2.3491 0.2014 2.4673 0.1894 0.0006 0.0055 88.2275 2.79 18 -0.1213 0.0348 3.6527 0.9998 4.5598 0.3189 0.0 164 0.0305 105.1312 2.4328 -0.0898 0.0259 2.271 1 0.6068 2.8 183 -0.0003 -0.0048 0.01 87

g 10

121.3090 3.3819 -0.0468 0.0199 1.6192 0.4246 2.0069 0.1574 -0.0007 0.0040 128.5213 -2.0644 -0.0276 0.01 35 0.3576 0.01 14 0.3586 0.0197 -0.0001 -0.0010 Vertical Stahilizer .~.~...~ - V V V V v

X e Y e I my I J x z e JyVye J L e Z e J Z e J,c; e J x z e

-3.0130 -4.0388 0 1 0.0867 I 26.0215 26.8758 50.3184 18.1691 0.1131 0.1417

17.3464

-5.0753 1 0.0419 8.2448 2.9690 10.2273 3.9298 0.0052 0.0149

32.7117 -3.0312 0.0347 6.6885 2.6766 8.6842 3.2180 0.0097 0.0124 48.0499 -0.6415 0.0298 7.0021 2.2154 8.7929 3.2766 0.0395 0.0935 49.3577 25.8074 0.074 1 5.0234 42.2968 46.9395 -2.2970 -0.0025 -0.0045 90.2614 -3.2365 0.0176 1.7168 1.1 158 2.6177 1.0828 -0.1246 -0.1088 10 I .949 1 19.9499 0.0062 0.21 11 0.1462 0.3470 0.0991 0 0 103.6312 -3.8959 0.0150 1.3610 0.8632 2.1464 0.8286 -0.0010 0.0035 114.7248 - 10.6405 0.0065 0.3567 0.2077 0.5619 0.2227 0 0

10 I 116.1657 2.4489 0 0.0523 4.5419 I .0689 5.5228 1.0377 0 0

2 1 108.62 -37.00 -13.96 I 0.2647 1 35.2545 137.1807 162.3658 7.4957 11.3756 5.7131

The stiffness properties consist of the flexural rigidity and torsional rigidity of the cross-sectional area at certain locations of the elastic components. Both have units lb - in2. Table 4 gives the loca- tions of the cross-sectional areas and the values of the corresponding rigidities.

Table 4 - Stiffness ProDerties

Fuselage Fore Part Aft Part

No. E Ir', EIfF, EIfAy

E l f z ," G J ~

- L 22-

1 0.00 1.10 x 10" 7.40 x 10'" 1.29 1.1 x 10" 7.44 x 10 +

2 10.21 1.10 x 10" 7.44 x 1 o ' O 5.04 1.09 x loll 7.40 x 10" 4.07 x 10"

3 25.21 1.10 x 10" 7.47 x 10'0 8.79 1.04 x 10" 7.37 x 1 o ' O 4.04 x lo1'

4 38.2 1 1.10 x 10" 7.47 x 1 o ' O 26.29 9.44 x 1 o 1 O 6.87 x 10'' 3.67 x lolo 5 47.2 1 1.10 x 10" 7.47 x 10'0 40.29 7.98 x 10" 6.17 x 10" 3.23 x 10" 6 55.2 1 1.10 x 10" 7.47 x 10'0 49.79 6.72 x 10'' 5.32 x 10" 2.87 x 10" 7 70.2 I 1.10 x 10" 7.45 x 1 o ' O 65.79 6.72 x 10'' 4.15 x 10" 2.27 x 10" 8 84.7 1 1.10 x 10" 7.36 x 10" 82.79 6.45 x 10" 3.21 x 10" 1.72 x lolo 9 99.7 1 1.10 x 10" 7.20 x 10" 01.29 5.15 x 10" 2.43 x 10" 1.37 x 10" 10 111.21 1.10 x 10" 7.03 x lolo 20.29 4.07 x 10" 1.94 x 10'' 1.12 x lo1'

1 1 129.21 1.08 x 10" 6.80 x 10" 38.79 2.30 x 10" 1.37 x 10'' 8.39 x io9

3.87 x lo1' ~ 12 144.21 1.02 x 10" 6.54 x 10" 156.29 1.09 x 10'0 9.68 x io9 5.79 x io9

3.14 x 10" 13 :: 160.21 9.06 x 10" 6.14 x 10'' 172.79 1.63 x 10'' 6.62 x io9 3.92 x io9

1.71 x 10" ~ 14 14 185.7 1 6.5 x 10" 5.25 x 10" 188.79 2.10 x 10'0 4.39 x io9 2.58 x io9

8.67 x io9 1 15 15 205.71 4.37 x 10" 4.23 x 10" 204.79 1.38 x 10'' 2.80 x 109 1.63 x io9

2.98 x lo9 1 16 16 224.7 1 2.68 x 10" 3.97 x 10" 220.79 7.48 x io9 1.70 x io9 9.67 x 10'

2.19 x lo9 17 17 238.1 1 1.76 x 10'' 2.04 x 10" 237.29 4.67 x io9 9.22 x 10' 5.46 x 10'

1.27 x lo9 18 18 251.11 1.11 x 10'' 1.28 x lolo 256.29 3.02 x io9 4.88 x 10' 3.08 x 10'

6.43 x 10' 19 19 267.06 6.01 x io9 6.42 x io9 279.79 1.30 x io9 2.40 x 10' x 10'

1.17

20 282.06 3.13 x io9 2.87 x io9

7.72 x io7 21 295.86 1.16 x lo9 6.15 x 10'

I 1 1 ; Wing Horizontal Stabilizer - R R No. EIE GJ: No. EIF GJ: 2 , 2, - - 1 0.00 1.09 x 10'" 1.07 x 10'" 1 0.00 3.92 x 10' 2.43 x 10'

2 17.10 9.70 x io9 1.07 x 10" 2 2.41 3.78 x 10'

2.36 x 10'

3 34.19 8.65 x io9 1.04 x lolo 3 13.84 3.16 x 10' x 10'

2.09

4 40.95 8.10 x io9 9.95 x io9 4 3 1.89 2.32 x 10' 1.73 x 10'

5 54.45 6.95 x io9 8.80 x io9 5 49.94 1.66 x 10'

1.35 x 10' c . L A 11 i 3'2 i n 8

6 67.85 5.70 x io9 7.60 x io9 l . & U h I"

1.03 x 10'

7 81.11 4.84 x io9 6.40 x io9 9.70 x io7 8.05 x io7

8 94.45 4.28 x io9 5.30 x io9 7.10 x io7

5.85 x io7

9 107.95 3.49 x io9 4.50 x io9 103.98 4.82 x io7 4.14 x io7

10 122.95 2.81 x io9 3.92 x io9 3.27 x io7 2.98 x io7

11 140.45 2.18 x io9 3.19 x io9 2.74 x io7

2.19 x io7

12 159.45 1.68 x io9 2.43 x io9 Vertical Stabilizer

V

13 179.45 1.35 x io9 EIY

GJ: X e

14 200.70 1.02 x io9 0.00 1.56 x 10'

7.30 x 10'

15 22 1.95 7.35 x 10' 7.40 1.33 x io9 6.70 x 10'

16 24 1.70 5.65 x 10' 21.15 9.66 x 10' 5.57 x 10' 17 3.56 x 10' 4 259.95 4.28 x 10' 34.90 6.86 x 10' 4.45 x 10' 18 2.53 x 10' 5 276.70 3.54 x 10' 49.73 4.61 x 10' 3.55 x 10' 1.71 x 10' 6 293.95 2.82 x 108 63.98 3.06 x 10' 19 2.41 x 10'

8.55 x io7 7 2.02 x 10' 77.63

20 31 1.10 2.01 x 10' 1.38 x 10'

6.50 x lo6 8 1.27 x 10' 91.31 1.26 x 10'

21 327.01 8.69 x io7

9 105.82 6.88 x lo7

3.56 x io7

113.48

10 4.25 x io7 8.51 x lo6

The aerodynamic forces involve the slope C L ~ ~ of the lift curve, the drag coefficient CDio, the slope Csai of the lateral force curve and the control effectiveness coefficients CLba,Csbe and CshT. To determine the aerodynamic forces, the aircraft components are divided into a given num- ber of sections and the coefficients are given for each of the sections. A typical section has a trapezoidal shape defined by two points Z i a Y i a Z i a and XibYibZib and by respective chords qa and G b , where the chords are parallel to the longitudinal axis xf of the fuselage (Fig. A). For all R . R R XwbYwbZwb

Wing

/

Elevator

Horizontal

Stabilizer 5o

P l

n -250 -300 I f * 150 n 0 ' I -150 -200 -250 -300 -350 50 0 -50 -100

Fore Fuselage (Horizontal)

Fore Fuselage (Vertical)

&I2 50 50 --L-....-

I j

r i 1

0 0 I I

I I

p 4 250 200 150 p2100 50 250 200 150 100 50 0 -50 I 0.2 0 . 4 0 . 6 0.8 Mach Number

Figure A. Aerodynamic Sections for the Model

sections, CDfo = C D ~ , O = 0.016 and k f = k , = 0.04. The lines of the aerodynamic centers are also shown in Fig. A; they are located at one quarter of the chord of the wing and horizontal and vertical stabilizers. The line of aerodynamic centers for the fuselage are located at one half of the chord of each aerodynamic section, as shown in Fig. A. The aerodynamic forces for the fore and aft parts of the fuselage are acting at the line of the aerodynamic centers. The control forces, on the other hand, are acting at 0.55 of the chord. The aerodynamic coefficients corresponding to the sections just mentioned are listed in Table 5. The slope of the lift coefficients and the control effectiveness for the aileron and the elevator listed in Table 5 are given for a single Mach number of 0.75. They must be corrected for different Mach numbers by the compressibility factor given in Fig. A. Note that there are only seven sections for the wing with control effectiveness coefficients; they correspond to the aileron.

Table 5

Fore Part of the Fuselage - Horizontal Lifting Surface

-36.25 0.17 40.00 151.21 36.25 0.17 40.00 0.065 151.21 -36.25 0.17 40.00 191.21 36.25 0.17 40.00 0.20 191.21 40.00 231.21 36.25 0.17 40.00 0.50 -36.25 0.17 231.21 46.00 0.50 -36.25 0.17 277.21 46.00 277.21 36.25 0.17 f X T a 13.21 0.142 69.2 1 0.142 126.23 0.142 0.242 182.21 248.21 0.305 Wing - n R C L a w

d a Z a .wb ?db zwb Cwb

0 -62.49 0 129.49 34.00 -62.49 0 129.49 4.9675 34.00 -62.49 0 129.49 58.36 -52.40 0 115.51 6.1879 58.36 -52.40 0 115.51 80.7 I -43.15 0 102.68 6.5604 80.7 1 -43.15 0 102.68 101.07 -34.72 0 90.99 6.7609 101.07 -34.72 0 90.99 120.35 -33.06 0 86.25 6.9328 -33.06 0 86.25 138.62 -31.49 0 81.75 7.1620 120.35 138.62 -31.49 0 81.75 155.90 -30.00 0 77.50 7.4485 155.90 -30.00 0 77.50 172.17 -28.60 0 73.50 7.7349 172.17 -28.60 0 73.50 187.43 -27.29 0 69.75 7.9641 187.43 -27.29 0 69.75 203.59 -25.89 0 65.77 8.1647 2.4494 -25.89 0 65.77 219.20 -24.55 0 6 1.93 8.2792 2.4838 203.59 -24.55 0 6 1.93 234.30 -23.25 0 58.21 8.3365 2.5010 2 19.20 -23.25 58.21 248.91 -21.99 0 54.62 8.3365 2.5010 234.30 0 -21.99 0 54.62 263.01 -20.78 0 51.15 8.2506 2.4752 248.91 -20.78 0 51.15 276.62 -19.61 0 47.80 8.0214 2.4064 263.01 -19.61 0 47.80 289.72 -18.48 0 44.58 7.7063 2.31 19 276.62 - 18.48 44.58 30 1.25 -17.49 0 41.74 7.2766 289.72 0 311.27 -16.62 0 39.28 6.6463 30 1.25 -17.49 0 41.74 - 15.85 0 37.06 5.7869 31 1.27 - 16.62 0 39.28 320.29 327.66 - 15.21 0 35.25 4.5837 320.29 - 15.85 0 37.06 335.02 -14.58 0 33.43 2.5783 327.66 -15.21 0 35.25 able 5 (Continued) ~~~~~~ orizontal Stabilizer ~~ R R R R R R R C L a e c16e Ceb X e a Y e a z e a C e a x e b Y e b zeb 0 -25.40 0 65.3 1 10.98 -24.06 0 61.62 1 S783 0.7891 58.67 1.6844 0.8422 10.98 -24.06 0 61.92 23.37 -22.77 0 35.24 -21.54 0 55.55 1.8104 0.9052 23.37 -22.77 0 58.67 -21.54 46.60 -20.36 0 52.57 1.9828 0.9414 35.24 0 55.55 46.60 -20.36 0 52.57 57.44 -19.23 0 49.72 2.2282 1.1 141 57.44 - 19.23 0 49.72 67.88 -18.14 0 46.98 2.5398 1.2699 67.88 -18.14 0 46.98 77.91 -17.10 0 44.34 2.89 13 1 .4457 77.9 1 -17.10 0 44.34 87.53 -16.10 0 41.82 3.1632 1.5816 87.53 -16.10 0 41.82 96.74 -15.15 0 39.40 3.3025 1.6512 96.74 -15.15 0 39.40 105.54 -14.23 0 37.09 3.2958 1.6479 105.54 - 14.23 0 37.09 113.93 -13.36 0 34.88 3.1300 1 S650 1.3727 1 13.93 -13.36 0 34.88 121.91 -12.53 0 32.79 2.7454 121.91 -12.53 0 32.79 129.08 -11.79 0 30.90 1.9828 0.9914 iertical Stabilizer v V V V \' 1 ' 1' V ceb C G e c S v 6 ? .

X e a Y e a z e a C e a x e b Y e b zeb -31.15 -39.58 0 102.00 -20.81 -37.96 0 97.85 0.9482 0.4267 -20.81 -37.96 0 97.85 -4.81 -35.44 0 91.43 1.4740 0.6633 -4.81 -35.44 0 9 I .43 3.05 -34.20 0 88.28 1.7290 0.778 1 3.05 -34.20 0 88.28 8.92 -33.28 0 85.92 2.6550 1.1948 8.92 -33.28 0 85.92 25.14 -30.73 0 79.42 2.8540 1.2843 25.14 -30.73 0 79.42 41.37 -28.18 0 72.91 3.1030 I .3964 41.37 -28.18 0 72.91 56.34 -25.82 0 66.90 3.3020 1.4859 56.34 -25.82 0 66.90 71.32 -23.47 0 60.89 3.4280 1 S426 71.32 -23.47 0 60.89 86.30 -21.11 0 54.88 3.4350 1 S458 86.30 -21.11 0 54.88 103.37 -18.43 0 48.03 3.2020 1.4409 103.37 -18.43 0 48.03 113.44 -16.84 0 43.99 2.3430 1 .os44 Engine Pylon - A

.fb Y f b zfb I CL", f

$a xC;a Y f a cC;a

50.79 0.00 -5.13 126.00 50.79 20.00 -11.83 I 126.00 1.2376

50.79 20.00 -11.83 126.00 50.79 33.00 -16.16 117.33 1.2892 50.79 33.00 -16.16 117.33 50.79 47.00 -20.83 108.00 1.2900 42.79 47.00 -20.83 1 16.00 42.79 56.00 -23.84 116.00 1.1250

42.79 56.00 -23.84 1 16.00 42.79 65.00 -26.85 116.00 I .os00

42.79 65.00 -26.85 1 16.00 42.79 74.00 -29.86 116.00 0.8500 42.79 74.00 -29.86 116.00 42.79 82.00 -32.53 116.00 0.5500

I REPORT DOCUMENTATION PAGE I Form Approved

Public reporting burden lor this collection 01 inlormation is estimated to average 1 hour per response including the time loneviewing instructions searching existing data scurces gathering and maintainin the data needed and completing and reviewing the collection of information Send comments regarding this burden estimate or any other aspect 01 this collection of in?ormation including suggestions lor reducing this burden to Washington Headquarters Services Directorate lor Information Operations and R e m s 1215 Jellerson Davis Highway Suite 1204 Arlington VA 22202-4302 and to the Ollice of Management and Budget Paperwork Reduction Project (0704 0188) Washington DC 20503 1. AGENCY USE ONLY (Leave blank) 2. REPORT DATE 3. REPORT TYPE AND DATES COVERED July 2003 Contract0 Report 4. TITLE AND SUBTITLE 5. FUNDING NUMBERS Integrated Approach to the Dynamics and Control of Maneuvering Flexible Aircraft NCC 1-346 6. AUTHOR(S) 706-3 1-4 1-02 Leonard Meirovitch and Ilhan Tuzcu 7. PERFORMINGORGANIZATION NAME(S) AND ADDRESS(ES) 8. PERFORMING ORGANlZATlON REPORT NUMBER Virginia Polytechnic Institute and State University Engineering Science and Mechanics Department, MC 02 19 Blacksburg, VA 2406 1 10. SPONSORING/MONITORING 9. SPONSORlNGlMONITORlNG AGENCY NAME(S) AND ADDRESS(ES) AGENCY REPORT NUMBER National Aeronautics and Space Administration NASA/CR-2003-2 1 1748 Langley Research Center Hampton, VA 2368 1-2 199 12b. DISTRIBUTION CODE 12a. DISTRIBUTION/AVAILABlLlTY STATEMENT Unclassified-Unlimited Subject Category 08 Distribution: Standard Availability: NASA CAS1 (301) 621-0390

I 14. SUBJECT TERMS I 15. NUMBER OF PAGES

Flexible Aircraft, Integrated Approach, Flight Dynamics,

F 16. PRICE CODE

Extended Aeroservoelasticity, Time Simulations

I

20. LlFJlTmml 17. SECURITY CLASSIFICATION 18. s t - 19. S t - o N OF THIS PAGE OF ABSTRACT OF ABSTRACT OF REPORT Unclassified Unclassified Unclassified UL I I 1 I standera Form 298 R ev. 2-B9 N - Prescribed by ANSI St .2-39-18 298-102