Abstract
Energy efficiency plays important role in aeroelastic design of flying wing aircraft and may be attained by use of lightweight structures as well as solar energy. NATASHA (Nonlinear Aeroelastic Trim And Stability of HALE Aircraft) is a newly developed computer program which uses a nonlinear composite beam theory that eliminates the difficulties in aeroelastic simulations of flexible high-aspect-ratio wings which undergoes large deformation, as well as the singularities due to finite rotations. NATASHA has shown that proper engine placement could significantly increase the aeroelastic flight envelope which typically leads to more flexible and lighter aircraft. The areas of minimum kinetic energy for the lower frequency modes are in accordance with the zones with maximum flutter speed and have the potential to save computational effort. Another aspect of energy efficiency for High Altitude, Long Endurance (HALE) drones stems from needing to minimize energy consumption because of limitations on the source of energy, that is, solar power. NATASHA is capable of simulating the aeroelastic passive morphing maneuver (i.e., morphing without relying on actuators) and at as near zero energy cost as possible of the aircraft so as the solar panels installed on the wing are in maximum exposure to sun during different time of the day.
1. Introduction
Over the last decade, Hodges and coworkers [1–3] at Georgia Tech have been extensively involved in aeroelastic simulation of very light and thus highly flexible aircraft for development of the next generation of unmanned aerial vehicles (UAVs) and/or High-Altitude, Long Endurance (HALE) aircraft, including flying wings. Such aircraft typically have high-aspect-ratio wings with high flexibility, which leads to large deformation. Consequently, linear aeroelastic analyses are incapable of predicting the stability characteristics of such aircraft. They successfully proved that only nonlinear aeroelastic analysis provides correct information regarding the aeroelastic flight envelope of this class of aircraft [1, 4, 5].
Nonlinear aeroelastic trim and stability of HALE aircraft, NATASHA, is a computer program developed by the authors of references [1, 4, 5] that accommodates modeling of large deformation of high-aspect-ratio flying wings. The theory behind NATASHA is based on the geometrically exact, nonlinear, composite beam theory of Hodges [6], along with the finite-state induced flow model of Peters et al. [7].
Previous comparisons by [2] showed that results from NATASHA are in excellent agreement with well-known beam stability solutions [8, 9], the flutter problem of [10], experimental data presented by [11], and results from well-established computer codes such as DYMORE [12, 13] and RCAS [14]. Mardanpour et al. [3] considered the classical cantilever wing model of Goland and Luke [15] and verified NATASHA for the behavior of the eigenvalues as well as the effect of sweep on divergence and flutter characteristics and it was shown that gravity and load factor play an important role in aeroelasticity of high-aspect-ratio wings [16].
A principal determinant of energy consumption in aircraft is drag, which must be opposed by engine thrust for the aircraft to fly. Flying wings may achieve significant drag reduction due to a smooth outer surface and the lack of a vertical tail [3]. Consequently, the performance of such aircraft may increase significantly, relative to conventional configurations of the same size. The potential increase of performance for this class of aircraft has inspired aeroelasticians to design new generation of aircraft based on a flying wing configuration [3]. Typical aeroelastic instability of these aircraft is body-freedom flutter when the short-period mode of the aircraft couples with the elastic bending-torsion modes [3, 17–24].
In another context, a morphing solar-powered flying wing can maximize the energy absorption of solar panels on the wing surfaces by changing its configuration such that the panels have highest exposure to the sun. This change in the geometry of the flying wing is highly effective in energy absorption during times just before sunset and just after sunrise, and consequently the aircraft can sustain longer flight [25]. Use of solar energy is a novel method that eliminates one of the design constraints to a considerable extent by removing the limitation on the source of energy. The morphing concept could be based on either wing morphing systems or airfoil morphing systems, or a combination of both [26]. So far in the literature, several morphing concepts and systems have been developed based on altering various geometric parameters of the wing (such as span, chord, camber, sweep, twist, and even airfoil thickness distribution) to make the aircraft suitable for different missions and flight conditions [19, 26, 26, 26–37]. The folding wing configuration has been analyzed using linear aeroelastic models [38–40] and nonlinear aeroelastic models [41, 42]. In all the mentioned works the weight of actuators and the actuation power that the morphing mechanisms require performing their task are the problematic parts of the design [26], in particular when it comes to morphing of flying wing and/or HALE aircraft.
In this paper, after a brief outline of the theory behind NATASHA, energy efficiency in aeroelastic design and simulation for flying wing configuration will be assessed as (a) a feature of the design (i.e., engine placement), which attains instability at higher speed with lighter aircraft structure, (b) a methodology that helps to decrease computational effort required for determining favorable locations for engine placement with the potential of higher flutter speed, and (c) a scheme to passively morph a solar-powered flying wing, so that exposure to the sun of solar panels distributed on the wings is maximized for higher absorption of solar energy at specific times of the day.
2. Theory
2.1. Nonlinear Composite Beam Theory
The fully intrinsic nonlinear composite beam theory [6] is based on geometrically exact first-order partial differential equations of motion for the beam that are independent of displacement and rotation variables. They contain variables that are expressed in the bases of the reference frames of the undeformed and deformed beams, and , respectively; see Figure 1. These geometrically exact equations are written in terms of force, moment, velocity, and angular velocity, and they contain no nonlinearities higher than second degree in the unknowns. The equations of motion are where the generalized strains and velocities are related to stress resultants and moments by the structural constitutive equations and the inertial constitutive equations Finally, strain- and velocity-displacement equations are used to derive the intrinsic kinematical partial differential equations [6], which are given as
In this set of equations, and are column matrices of cross-sectional stress and moment resultant measures in the frame, respectively; and are column matrices of cross-sectional frame velocity and angular velocity measures in the frame, respectively; and are column matrices of cross-sectional linear and angular momentum measures in the frame, respectively; , , and are partitions of the cross-sectional flexibility matrix; is the identity matrix; is the cross-sectional inertia matrix; is with and being position coordinates of the cross-sectional mass center with respect to the reference line; and is the mass per unit length. Note that the tilde notation denotes the antisymmetric matrix associated with the column matrix over which the tilde is placed; denotes the partial derivative with respect to time; and denotes the partial derivative with respect to the axial coordinate, . More details about these equations can be found in [43].
This is a complete set of first-order, partial differential equations. To solve this complete set of equations, one may eliminate and using (2) and and using (3). Then, 12 boundary conditions are needed, in terms of force (), moment (), velocity (), and angular velocity (). The maximum degree of nonlinearities is only two, and because displacement and rotation variables do not appear, singularities and other concerns related to finite rotation are avoided.
If needed, the position and the orientation can be calculated as postprocessing operations by integrating where is the transformation matrix from undeformed frame of reference to the inertial frame of reference and is the transformation matrix from deformed frame of reference to the inertial frame of reference .
2.2. Finite-State Induced Flow Model of Peters et al.
The two-dimensional finite-state aerodynamic model of Peters et al. [7] is a state-space, thin-airfoil, inviscid, incompressible approximation of an infinite-state representation of the aerodynamic loads, which accounts for induced flow in the wake and apparent mass effects, using known airfoil parameters. It accommodates large motion of the airfoil as well as deflection of a small trailing edge flap (). Although the two-dimensional version of this theory does not account for three-dimensional effects associated with the wing tip, published data [2, 7] show that this theory is a suitable choice for representing the aerodynamic loads acting on high-aspect ratio wings. Expressions for lift, drag, and pitching moment at the quarter-chord are given by where and , are the measure numbers of and is the flap deflection. The effects of unsteady wake (induced flow) and apparent mass appear as and acceleration terms in the force and moment equations given asThe induced flow model of Peters et al. [7] is included to calculate as where is the column matrix of induced flow states and , , and are constant matrices derived in [7].
3. Nonlinear Aeroelastic Trim and Stability of HALE Aircraft (NATASHA)
NATASHA is a computer program that is based on a geometrically exact, composite beam formulation [43] and a finite-state induced flow aerodynamic model [7]. The governing equations for the structural model are geometrically exact, fully intrinsic, and capable of analyzing the dynamical behavior of a general, nonuniform, twisted, curved, anisotropic beam undergoing large deformation. The dependence of the partial differential equations on is approximated by spatial central differencing [5]. A discretized version of the resulting nonlinear ordinary differential equations is linearized about a state of steady motion governed by nonlinear algebraic equations which NATASHA solves to obtain the so-called trim solution, using Newton-Raphson or similar procedures [5]. This system of nonlinear aeroelastic equations, when linearized about the resulting trim state, leads to a standard eigenvalue problem, which NATASHA uses to assess stability.
4. Effect of Engine Placement on Aeroelastic Flight Envelop
4.1. High-Aspect-Ratio Wing
A lightweight, small-class thrust turboshaft engine (JetCat SP5) [44] with known thrust, mass, moments of inertia, and angular momentum, which operates at the cruise condition, is placed along the span of a highly flexible, high-aspect-ratio wing [16]. The engine mounts are displaced with an offset from the elastic axis parallel to the plane of the wing cross section, that is, the local plane, while orientations of the engines are maintained; see Figure 2.
In stability analysis with gravity ( m/s2) the wing flutters at 22.4 m/s, with a hump first bending-torsion mode and with a frequency of 12.199 rad/sec which returns to the stability region at higher speeds. At 45 m/s the second bending-torsion with a frequency of 27.406 becomes unstable [16]. NATASHA’s result for flutter characteristics captures significant increase in flutter speed, more than three times the baseline model (i.e., with engines at the root of the wing), for engines placed at 60% to 80% span, and forward of the elastic axis. Also, a significant decrease in flutter speed is observed for engine placement at 30% to 50% span and aft of the elastic axis. Here forward means in the direction, and above the elastic axis means in the direction; see Figures 3 and 4.
4.2. Two-Engine Flying Wing
The geometry of the flying wing resembles the HORTEN IV [18, 20] and the properties were tuned [3] such that the aircraft experiences body-freedom flutter with the frequencies close to those of the body-freedom flutter frequency obtained from interviews with HORTEN IV pilots [20]. Two identical engines were symmetrically moved along the span, that is, in the direction, and displaced as described in Section 4.1; see Figure 5.
The effect of engine placement showed that the maximum flutter speed occurs when the engines are just outboard of 60% span; see Figure 6. The minimum flutter speed occurs for engine placement at the wing tips. Both minima and maxima occurred when the c.g. of the engine was forward of the wing elastic axis [3]; see Figures 7 and 8.
4.3. Four-Engine Flying Wing
The four-engine configuration (see Figure 9) is an alternative because of safety and reliability considerations associated with two-engine aircraft. The engine offsets from the elastic axis are presented in polar coordinates with where is the engine number. Figures 10 and 11 show the variation in flutter speed for different engine placements along the span with different offsets from the elastic axis while one of the engines was kept fixed in a particular location and the other one moves along the span.
Former study showed that engine placement does not have a significant effect on the lift to drag ratio [24]. However, a noticeable increase in flutter speed is observed when engines are placed forward of the elastic axis. For these cases, as one of the engines is placed at the outboard portion of the span, flutter speed increases. A detailed study on this effect is presented in the work of Mardanpour et al. [24].
5. Area of Minimum Kinetic Energy Density of the Mode
Engine placement at certain location has the potential to increase the flutter speed in two ways. One is the location where lower frequency flutter mode could be relegated to a higher frequency mode and the other is the location where the fluid structure interaction is decreased. Both criteria could be met at the area of minimum kinetic energy of the mode. In other words, engine placement at the area of minimum kinetic energy of the modes has the potential to decrease fluid-structure interaction and enforce the structure to flutter at a higher mode [16]. In case of the high-aspect-ratio wing presented in Section 4.1, the area of minimum kinetic energy density of the first two bending and torsion modes of the wing, in the absence of engines, gravitational, and aerodynamic forces, is presented in Figures 12–15. For the first bending and torsion modes this area is a minimum at the root of the wing; see Figures 12–14. For the second bending mode the area of minimum kinetic energy density of the mode has a local minimum outboard of the 85% span, and for second torsion mode this minimum moves to the region between 70% and 90% span; see Figures 13–15. For engine placement forward of the elastic axis, the unstable mode contains a combination of first and second bending along with second torsion; and when the engines are placed around 50% to 70% span, there is a noticeable increase in flutter speed. This is close to the area of minimum kinetic energy of the second bending and torsion modes; see Figures 13–15.
To further explore this possibility, the same analysis is performed for the two-engine and four-engine configurations. Figure 16 depicts the kinetic energy distribution of the symmetric free-free mode of the bare wing for the case study presented in Section 4.2. One can see that the point of the minimum kinetic energy is just outboard of 60% span [3].
In case of the four-engine configuration, for engine placement forward of the elastic axis, the unstable mode contains a combination of first, second, and third bending modes [24]; and when the engines are placed around 60% to 80% span, there is a noticeable increase in flutter speed. This area is close to the area of minimum kinetic energy of the first three bending modes; see Figures 17, 18, and 19.
6. Passive Morphing of Solar Powered Flying Wing Aircraft
HALE aircraft can achieve sustained, uninterrupted flight time if they use solar power. Wing morphing of solar-powered HALE aircraft can significantly increase solar energy absorbency [25]. An example of the kind of morphing considered in this chapter requires the wings to fold so as to orient a solar panel to be hit more directly by the sun’s rays at specific times of the day [25]; see Figure 20.
In this study solar-powered HALE flying wing aircraft are modeled with three beams with lockable hinge connections. There were four identical engines with prescribed mass, inertia, and angular momentum placed at outer portion of the beams with zero offset from the elastic axis; see Figure 21. Three sets of flaps are distributed on the wings: one set (set 1) is on the side wings; the middle beam was divided into equal portions for the two other sets of flaps (sets 2 and 3); see Figures 21 and 22.
Such aircraft are shown to be capable of morphing passively, so that the solar panels are oriented toward the sun using only input from flight control surfaces and engine thrust. The analysis underlying NATASHA was extended to include the ability to simulate morphing of the aircraft into a “Z” configuration; see Figure 22. Because of the “long endurance” feature of HALE aircraft, such morphing needs to be done without relying on actuators and at as near zero energy cost as possible.
For simulation of this complete maneuver, a schedule of “trim conditions” is solved using the trim analysis in NATASHA. The resulting control history (see Figure 23) is then interpolated from these solutions and input into NATASHA simulation. The stability analysis of the trim conditions showed that this aircraft is aeroelastically stable and the instabilities are flight dynamics modes with very small nonoscillatory eigenvalues [25].
Detailed systematic processes for trim and time-marching for a passively morphing flying wing are presented in [25].
7. Concluding Remarks
Energy efficiency in aeroelastic design and simulation of flying wings may be enhanced by appropriate choice of design features, for example, engine placement, which attains instability at higher speed with lighter aircraft structure. For this purpose, a very flexible high-aspect-ratio wing, which displays geometrically nonlinear behavior in subsonic flow, is analyzed using NATASHA. The effect of engine placement along the span with various offsets from the elastic axis is studied using a small-class, lightweight thruster (e.g., JetCat SP5) operating at the cruise condition. A significant decrease in flutter speed is observed for engine placement at 30% to 50% span, aft of the elastic axis, that is, in the negative direction, and above the elastic axis, that is, in the positive direction. The same effect for a two-engine flying wing configuration, with geometry similar to the Horten IV, was studied and because of safety and reliability considerations associated with two-engine aircraft the analogy was extended to a four-engine configuration. For the case of a two-engine configuration, the maximum flutter speed occurs when the engines are just outboard of 60% span. The minimum flutter speed occurs for engine placement at the wing tips. Both minima and maxima occurred when the c.g. of the engine was forward of the wing elastic axis. In case of a four-engine configuration, a noticeable increase in flutter speed is observed when engines are placed forward of the elastic axis. For these cases, as one of the engines is placed at the outboard portion of the span, flutter speed increases. For engine placement aft of the elastic axis, flutter speed increases when both engines are close to the root.
In the absence of engines and aerodynamic and gravitational forces, NATASHA found that (a) for a high-aspect-ratio wing, the minimum kinetic energy region of the second bending mode has a local maximum around 45% to 55% span and a local minimum at 80% to 90% span and the second torsion mode has these minima at 30% to 40% and 65% to 75% span, (b) for a two-engine flying wing, it is very near the 60% span for the first symmetric elastic free-free bending mode, and (c) for a four-engine configuration, the minimum kinetic energy regions for the first and third bending modes are located at approximately 60% span and for the second mode, the kinetic energy has local minima around 20% and 80% span. The areas of minimum kinetic energy for these modes are in agreement with flutter calculations, which show a noticeable increase in flutter speed if engines are placed forward of the elastic axis at these regions. Knowledge of the minimum kinetic energy region thus has some potential to help the designer and could possibly lead to significant savings in computational effort.
The analysis underlying NATASHA was extended to include the ability to simulate morphing of the aircraft using a new set of trim and kinematical differential equations. An example of the kind of morphing considered in this study requires the wings to fold so as to orient solar panels to be hit more directly by the sun’s rays at specific times of the day. Because of the “long endurance” feature of HALE aircraft, such morphing needs to be done with as near zero energy cost as possible, that is, without relying on actuators at the hinges but instead using aerodynamic forces, controlled by flap deployments, and engine thrusts. The three-wing, solar-powered HALE aircraft morphs passively into a Z-shaped configuration while local bending moments are zeroed out at the hinge connection points, but with the hinges locked with the fold angle held at zero degrees while the aircraft morphs. The morphing motion is brought to a stop before the hinges are again locked.
8. Future Work
The present aeroelastic simulation of passive wing morphing neglects the effects of hinge stiffness and damping, but their addition is planned for a later work to enhance the realism of the model. The addition of hinge damping will allow energy dissipation to be taken into account, and it is expected that incorporation of a spring-restrained hinge may facilitate passive morphing and gust response of the aircraft during the morphing phases. Other recommendations for future studies include addition of the 3D Peters aerodynamic model to account for moderate aspect-ratio wings and axial flow. Such an aerodynamic model will extend the validity of analysis done by NATASHA to rotary wings.
In the present effort the effect of engine aerodynamics has not been considered. Explicit formulation of the same would require the use of 3D computational fluid dynamic (CFD) model. At this present moment, to the best of the author’s knowledge accurate simulation of the engine aerodynamics will retain the qualitative aspect of the current analysis but might lead to some quantitative improvements. Another aspect worthy of future investigation is to account for the flexibility of the engine nacelle and mount. It will be interesting to outline whether consideration of their flexibility and the resulting dynamics affects the aeroelastic behavior of the aircraft. Also, the inertial properties and flexibilities of control surfaces are not studied, including those in NATASHA that facilitates the simulations.
Nomenclature
| : | Deformed beam aerodynamic frame of reference |
| : | Undeformed beam cross-sectional frame of reference |
| : | Deformed beam cross-sectional frame of reference |
| : | Unit vectors in undeformed beam cross-sectional frame of reference () |
| : | Unit vectors of deformed beam cross-sectional frame of reference () |
| : | Chord |
| : | Derivative of pitch moment coefficient with respect to flap deflection |
| : | Derivative of lift coefficient with respect to angle of attack () |
| : | Derivative of lift coefficient with respect to flap deflection () |
| : | Column matrix |
| : | Offset of aerodynamic center from the origin of frame of reference along |
| : | Column matrix of distributed applied force measures in basis |
| : | Column matrix of internal force measures in basis |
| : | Gravitational vector in basis |
| : | Column matrix of cross-sectional angular momentum measures in basis |
| : | Inertial frame of reference |
| : | Unit vectors for inertial frame of reference () |
| : | Cross-sectional inertia matrix |
| : | Column matrix of undeformed beam initial curvature and twist measures in basis |
| : | Column matrix of deformed beam curvature and twist measures in basis |
| : | Column matrix of distributed applied moment measures in basis |
| : | Column matrix of internal moment measures in basis |
| : | Column matrix of cross-sectional linear momentum measures in basis |
| : | Column matrix of position vector measures in basis |
| : | partition of the cross-sectional flexibility matrix |
| : | partition of the cross-sectional flexibility matrix |
| : | partition of the cross-sectional flexibility matrix |
| : | Column matrix of displacement vector measures in basis |
| : | Column matrix of velocity measures in basis |
| : | Axial coordinate of beam |
| : | Trailing edge flap angle |
| : | Identity matrix |
| : | Column matrix of 1D generalized force strain measures |
| : | Column matrix of elastic twist and curvature measures (moment strain) |
| : | Induced flow velocity |
| : | Mass per unit length |
| : | Column matrix of center of mass offset from reference line in basis |
| : | Column matrix of cross-sectional angular velocity measures in basis |
| : | Partial derivative of with respect to |
| : | Partial derivative of with respect to time |
| : | Nodal variable |
| : | At flutter point. |
Conflict of Interests
The authors declare that there is no conflict of interests regarding the publication of this paper.
Acknowledgments
This work was supported in part by the Air Force Research Laboratory Contract FA8650-06-D-3624/000708 through Bihrle Applied Research, Inc., and Award ECCS-1101431 from the National Science Foundation.