

small (250x250 max)
medium (500x500 max)
Large
Extra Large
large ( > 500x500)
Full Resolution


DYNAMIC SIMULATION OF LOW ALTITUDE AERIAL TOW SYSTEMS By JOHNNY EARL QUISENBERRY JR. Bachelor of Science Oklahoma State University Stillwater, Oklahoma 2002 Submitted to the Faculty of the Graduate College of Oklahoma State University in partial fulfillment of the requirements for the Degree of MASTER OF SCIENCE July, 2005 ii DYNAMIC SIMULATION OF LOW ALTITUDE AERIAL TOW SYSTEMS Thesis Approved: ________________________________________________ Thesis Advisor ________________________________________________ ________________________________________________ ________________________________________________ Dean of the Graduate College iii ACKNOWLEDGEMENTS This research was conducted in cooperation with the engineers at Zivko Inc. Specifically; I would like to thank Todd C. Morse for his assistance. In addition I would like to thank the NASA Oklahoma Space Grant Consortium for their financial support for the duration of the project. I would like to express my deepest appreciation to my major advisor, Dr. Andrew S. Arena, for his guidance and support in this research. I would also like to thank the other members of my committee, Dr. R. D. Delahoussaye and Dr. P. R. Pagilla for their time and efforts. I could not have made the accomplishments I have here at OSU without the love and support of my family; mother, Susan Quisenberry, father, John Quisenberry, sister, Kiesha Quisenberry, and brother Clint Quisenberry. iv TABLE OF CONTENTS CHAPTER PAGE 1 INTRODUCTION ........................................................................................................... 1 1.1 Background ......................................................................................................... 1 1.2 Design Goals ....................................................................................................... 4 1.3 Literature Review................................................................................................ 5 1.3.1 Cable Dynamics .............................................................................................. 6 1.3.2 Lagrange’s Equations ...................................................................................... 9 1.3.3 Aircraft Dynamics......................................................................................... 13 1.3.4 Cable Towed Systems ................................................................................... 14 2 TWODIMENSIONAL DYNAMICS........................................................................... 18 2.1 Introduction....................................................................................................... 18 2.2 Cable Dynamics ................................................................................................ 19 2.2.1 Coordinate Systems ...................................................................................... 22 2.2.2 Lagrange Equations ....................................................................................... 25 2.2.3 Aerodynamic Forces ..................................................................................... 35 2.3 ATV Dynamics ................................................................................................. 37 2.3.1 Coordinate system......................................................................................... 38 2.3.2 Lagrange equations ....................................................................................... 40 2.3.3 Generalized Forces........................................................................................ 43 2.4 Host Vehicle...................................................................................................... 45 2.5 System Equations of Motion of an Alternate Segment Model ......................... 50 3 THREEDIMENSIONAL DYNAMICS ....................................................................... 54 3.1 Introduction....................................................................................................... 54 3.2 Cable Dynamics ................................................................................................ 54 3.2.1 Coordinate Systems ...................................................................................... 55 3.2.2 Lagrange Equations ....................................................................................... 60 3.2.3 Aerodynamic Modeling ................................................................................ 79 3.3 Arial Tow Vehicle............................................................................................. 82 3.3.1 Coordinate System........................................................................................ 83 3.3.2 Lagrange Equations ....................................................................................... 85 3.3.3 Generalized Forces........................................................................................ 97 3.4 Host Vehicle...................................................................................................... 99 3.4.1 Coordinate System...................................................................................... 100 3.4.2 Lagrange Equations ..................................................................................... 101 v 3.4.3 Generalized Forces...................................................................................... 106 3.5 System Equations of Motion of an Alternate Segment Model ....................... 108 4 FLIGHT PARAMETERS ............................................................................................ 111 4.1 Introduction..................................................................................................... 111 4.2 Ocean Surface Waveform............................................................................... 112 4.3 ATV Autopilot ................................................................................................ 114 5 SYSTEM SIMULATIONS.......................................................................................... 116 5.1 The Ideal Cable ............................................................................................... 116 5.2 Lumped Mass Versus Thin Rod Cable Model................................................ 119 5.3 Simulation of the CableATV System............................................................ 124 5.3.1 Dimensional Comparison............................................................................ 125 5.3.2 Altitude Effects ........................................................................................... 126 5.4 ATV Wave Tracking....................................................................................... 131 5.4.1 Surface Waveform Variation...................................................................... 132 5.4.2 Host Aircraft Cruise Altitude ...................................................................... 133 5.4.3 Lateral Motion............................................................................................. 138 5.4.4 Host Vehicle Oscillation............................................................................. 141 6 CONCLUSIONS.......................................................................................................... 144 7 REFERENCES ............................................................................................................ 147 vi LIST OF TABLES Table 5.1 Cable parameters and reference flight conditions................................... 120 Table 5.2 Aerial towed aircraft parameters............................................................. 125 vii LIST OF FIGURES Figure 1.1 The aerial tow system ................................................................................. 1 Figure 1.2 Simple pendulum shown in twodimensional a) Cartesian coordinates, and b) cylindrical coordinates.......................................................................... 11 Figure 2.1 Finite element cable segment models ....................................................... 20 Figure 2.2 Lumped mass finite element cable model................................................. 22 Figure 2.3 Twodimensional inertial and cable segment coordinate systems ............ 24 Figure 2.4 Twodimensional ATV coordinate system ............................................... 38 Figure 3.1 Threedimensional inertial and cable segment coordinate systems .......... 56 Figure 3.2 Euler angle rotations ................................................................................. 58 Figure 3.3 ATV coordinate systems ........................................................................... 83 Figure 4.1 Threedimensional ocean surface waveform .......................................... 114 Figure 5.1 Two and Threedimensional ideal cable simulation results. The figure shows, a) the height of cable base, and b) the system energy versus time. ................................................................................................................. 118 Figure 5.2 A cable in fluid flow simulation results. The cable is broken into various numbers of segments, the segment is modeled as, a) lumped masses, and b) thin rods. ............................................................................................. 121 Figure 5.3 Cable base static positions versus number of cable segments ................ 122 Figure 5.4 Potential energy versus number of cable segments. ............................... 123 Figure 5.5 Comparison between two and threedimensional cable and ATV system ................................................................................................................. 126 Figure 5.6 Host vehicle cruise altitude affect on static cable curvature and position ................................................................................................................. 127 Figure 5.7 Dynamic response to host vehicle oscillation as viewed from the host vehicle coordinate system. ...................................................................... 128 Figure 5.8 Elliptic oscillations in the x and zdirection, fo r a host vehicle altitude of a) 250 ft., b) 500 ft., c) 1000 ft., and d) 2000 ft. ..................................... 129 Figure 5.9 Tracking at with various ocean wavelengths .......................................... 133 Figure 5.10 Host vehicle cruise altitude affect on wave tracing a) 250 ft., b) 500 ft., c) 1000 ft., and d) 2000 ft........................................................................ 134 Figure 5.11 Elliptic oscillations in the x and zdirection, for a host vehicle altitude of a) 250 ft., b) 500 ft., c) 1000 ft., and d) 2000 ft. ..................................... 135 Figure 5.12 ATV motion in the x and ydirection, at a host vehicle cruise altitude of a) 250 ft., b) 500 ft, c) 1000 ft., and d) 2000 ft. ..................... 139 Figure 5.13 Threedimensional waveform tracking with a 10 knot side gust, host vehicle altitude of a) 250 ft., b) 500 ft., c) 1000 ft., and d) 2000 ft. ....... 140 Figure 5.14 Wave tracking during host vehicle oscillation, host vehicle altitude of a) 250 ft., b) 500 ft., c) 1000 ft., and d) 2000 ft. ..................................... 143 viii NOMENCLATURE A0 Þ general amplitude of oscillation A,B,C Þ augmented equation of motion matrix C Þ coordinate transformation from segment to host frame Cß,Sß, Þ shorthand notation for cosine and sine of generic angle ß cf Þ cable friction coefficient cp Þ cable pressure coefficient d Þ cable diameter Fx,Fy,Fz Þ X,Y,Z component of the aerodynamic force, in the cable segment frame f Þ generic function h Þ wave height in host vehicle frame hn Þ sea level in host vehicle frame I Þ ATV inertia matrix, in ATV frame I’ Þ ATV inertia matrix, in ATV segment frame IX,IY,IZ Þ ATV X,Y,Z inertia components, in ATV frame kd Þ derivative gain kp Þ proportional gain l Þ segment length L,M,N Þ aerodynamic moments about X,Y,Z, in ATV frame MX,MY,MZ Þ moment about X,Y,Z, in segment frame ix Mf& Þ moment matrix about the Euler rates f& Mq& Þ moment matrix about the Euler rates q& My& Þ moment matrix about the Euler rates y& m Þ mass n Þ number of cable segments p,q,r Þ angular rates about X,Y,Z, in ATV aircraft frame Qi Þ generalized force for the ith generalized coordinate q Þ generalized coordinate t Þ time T Þ kinetic energy Trot Þ rotational kinetic energy U0,V0,W0 Þ freestream flow, in inertial frame u,v,w Þ perturbation velocities, in ATV aircraft frame V Þ velocity Vx, Vy, Vz Þ velocity in the X,Y,Z direction, in segment frame Xu,…Np Þ ATV stability derivatives X,Y,Z Þ Cartesian coordinate directions, in host frame x,y,z Þ Cartesian coordinate directions, in segment frame z Þ vertical host aircraft displacement, inertial frame l Þ wavelength p Þ circular constant r Þ air density x f,q,y Þ Euler angles ? Þ Angle between the ATV segment frame and ATV aircraft frame w Þ vector of Cartesian angular rates w0 Þ general frequency Dh Þ difference in ATV and wave altitude, in host frame Diw Þ wing angle correction 1 CHAPTER 1 1INTRODUCTION 1.1 Background The aerial tow system is a relatively common flight configuration. It is a necessity for gliding aircraft to achieve altitude; it is common to have targets towed by aircraft for weapons testing and target practice; it is even being proposed by NASA as a platform for reusable launch spacevehicles. The configuration has been proven to be well suited for these uses, though they are primarily a means of transportation of the towed vehicle. Figure 1.1 The aerial tow system 2 The aerial tow system is composed of three major components (1) host vehicle, (2) aerial tow vehicle (ATV), and (3) connecting cable, figure 1.1. In general, the host vehicle is larger than that of the tow vehicle and provides the thrust for both aircrafts. The cable is connected to both the host vehicle and the ATV in such a way that both aircraft can safely maintain flight. It is not uncommon in this configuration for the connecting cable to be hundreds or thousands of feet in length. As a result, the cable orientation is subject to change due to variations in the freestream conditions, cable dynamics, and the dynamics of the two aircraft s. The thrust produced by the host aircraft’s engines is directly applied to the host aircraft and indirectly applied to the ATV thru the tension in the connecting cable. The tension in the cable produces a force on the ATV that may differ from the host vehicle thrust in magnitude and direction depending on the cable orientation. The variable cable force on the ATV causes a dynamic response resulting in further deformation of the connecting cable. Clearly, the conditions effecting these variations are coupled. Historical research suggests that, in general, the system is stable, and given sufficient time the disturbances in the system will damp out. However the degree to which the atmospheric and dynamic alterations affect the system depends principally on the system components and is, of yet, not well known. Recently the aerial tow system has been chosen as a means to carry out the measurement of atmospheric conditions in close proximity to the ocean surface. The proposed system calls for a general aviation aircraft to tow the measurement instrumentation equipped ATV near (or approximately 30 feet) above the ocean surface. The emphasis of this system is the controlled altitude of the ATV. This is a deviation of 3 the aerial tow system from its roots as a means of transportation. Using the aerial tow system simply for transportation places no strict requirement on the motion of the ATV. The aerial tow system is proposed as a safer, more robust alternative for recording oceanic atmospheric conditions. The traditional method necessitates an aircraft fitted with the appropriate sensory equipment fly immediately above the ocean surface. The aircraft’s autopilot cannot be relied upon to safely maintain such altitude in this environment. Variable terrain and unknown atmospheric conditions could, indeed, result in disaster for the aircraft. As a result, the mission puts forth a tremendous workload upon the pilot. Due to the intense flight regime and pilot stress, the testing is usually limited to short duration and is confined to calm conditions. The implementation of the aerial tow system as a method to take oceanic atmospheric measurements would allow the piloted host aircraft to fly at a safe height above sea level (such as 250 to 2000 feet) under the control of the aircraft autopilot while towing an unmanned sensory aircraft a matter of feet above the ocean. This would allow for longer tests runs, eliminate pilot stress, expand the mission envelope, as well as reduce the worstcase scenario to the mere loss of the unmanned aircraft, assuming this flight configuration is suitable for the task at hand. Given the possibility of difficult ocean and atmospheric conditions, such as swales and wind gusts, the controllability of the tow system is in question. The system component interactions are complex, and neither the system configurations for the most advantageous performance nor the circumstances at which the system fails to adequately 4 perform are obvious. As a result, a system simulation preceding construction and testing is desirable. Each component of the aerial tow system must be mathematically modeled in order to build a numerical system simulation. Simultaneous derivation of the component equations of motion will produce the governing equations for the aerial towed system from which a response to initial and applied conditions can be simulated. The utilization of a system simulation would allow for prediction of the system’s ability to perform under various conditions such as, wind gust, wave shape, and host aircraft oscillation. Previous aerial tow systems do not place such heavy emphasis on the maintenance of the altitude of the ATV. The unique mission of near ocean atmospheric measurement has driven the design of a mathematical simulationof the low altitude aerial tow system. 1.2 Design Goals The objectives of this study is to create a dynamic simulation of the aerial tow system with the main focus surrounding the derivation of the system and component equations of motion in both two and threedimensions. The study includes mathematical models for each of the system components as well as the methodology for the derivation of the equations of motion for each of the component models mentioned. The system will be adequately designed such that the effect of the following conditions on the tracking ability of the ATV may be simulated: · Variable surface oscillations (ocean swales) 5 · Variable host vehicle altitudes · Lateral gusts · Host vehicle oscillation The decisions upon the design of the system will be made with the versatility of the simulation in mind, in the hopes that this research will assist with a wide range of future studies of aerially towed systems, underwater towed systems, or cable dynamics. 1.3 Literature Review The study of the tow system and its components has been well documented. The interest of cable dynamics, alone has produced analysis on the topics of cable oscillation and damping due to uniform flow, threedimensional cable shaping, cable drag model validity, and numerical simulations of cable dynamics. Refs 14. However, the integrated effect of the ATV is necessary to simulate the aerial tow system. As various tow systems come into existence study of the systems have flourished. Papers on tow vehicle controllability, stability augmentation, and release and retrieval dynamics build upon the dynamics of the cable in a flow, Refs 5 and 6. The investigation of tethered aerostat dynamics and underwater towing illustrates the range of the application and research of the tow system, Refs 79. Each of these references provides greater understanding of the tow system and its components. 6 1.3.1 Cable Dynamics The cable has found its way into many applications involving fluid flow through out the years. As a result, countless papers have been written for discussio n on various cable related topics. In general, the research was initiated for reasons supplementary to the cable’s effect on aerial tow systems; however, each gives insight into the overall behavior of the cable in fluid flow. The past studies give a foothold for future system analysis containing cable components. Of particular interest are the studies on the physical and aerodynamic models of the cable. 1.3.1.1 Cable Models A prediction of the dynamic response of a cable to given conditions requires a mathematical model of the cable from which the equations of motion can be derived. Choo and Casarella10 produced a survey of various methods for modeling cables. The motivation behind this article was to acknowledge various methods that could be used for a basis of cable system simulations and analyzing the merits and demerits of each technique. Among the methods discussed is the finite element method. This method breaks the continuous cable into an arbitrary number of segments which themselves can be modeled in a number of ways, such as the ideal pendulum, a mass on a spring, a cylinder, or curved segment model. In any case, each segment is free to rotate about its neighboring segments. Discretizing the cable allows for the derivation of the equations of motion of the interacting cable segments. The shortcoming to this technique is that the quality of the simulation is directly dependant upon the number of segments used to model the cable; a 7 finer discretization produces a superior simulation. However, as Walton and Polachek18 discuss, when the segment length decreases a smaller step size is necessary for numerical convergence during simulation, which translates to more computation time. The benefit of using the finite element technique is that it can be used to study several types of unsteady motion in a cable or cablebody system. It is for this primary reason the Choo and Casarella10 labeled this method as the most versatile. Also discussed in the survey are the method of characteristics, the linearization method, and the equivalent lumped mass method. The method of characteristics employs constitutive laws such as Hooke’s law, to convert the partial differential equations of motion into ordinary differential equations of motion. This technique allows the study of any sort of unsteady cable motion as long as the constitutive laws are satisfied, as an example; this method cannot be used on inextensible and viscoelastic cables as there are fewer ordinary differential equations than are original partial differential equations. In addition, the method of characteristics requires a large amount of computation. The linearization method, as it suggests, linearizes the cable equations of motion. This linearization limits the technique to perturbation, small deviations from equilibrium, studies, such as stability and frequency response of the system. This technique is not applicable in cases where a limit cycle persists and additionally is not useful in the analysis of unsteady cable motion. The equivalent lumped mass technique is useful when the most important aspect of the cablebody system is the body. In this method, the dynamics of the cable are completely ignored and the equations of motion of the system are largely due to the rigidbody 8 equations of motion of the body with added terms representing the cable dynamic effects. Since the equations of motion of the cable are omitted the system is greatly simplified. The simplicity, unfortunately, comes at the cost of the robustness of the model. Each of these methods has a configuration or condition in which it gains advantage over the rest. However, the only one of the techniques with no systemtype limitations is the finite element technique. An additional detailed description of the finite element method will be given in the following chapters. 1.3.1.2 Cable Drag Model The work of Hoerner13 is directed at analyzing and modeling fluid drag about various bodies. The cable is among the bodies researched. Hoerner13 produced the crossflow method for calculating the drag forces on a cable within a fluid. This technique analyzes the friction and pressure forces on the cable. The method is limited to straight cables of circular crosssection and fluid conditions that produce subcritical Reynolds numbers with respect to the cable diameter. Despite the straight cable limitation of the cross flow principle it can predict the drag force on a flexible cable (even though the cable is rarely completely straight). A curved cable can be approximated by a series of connected straight segment. Given the orientation of the straight segment relative to the fluid flow, the drag on each segment can be predicted using the crossflow principle. Thus, the resultant fluid force on a curved cable is the accumulation of the individual segment’s drag forces. As the number of cable segments increase so to does the accuracy of the overall drag prediction on the curved cable (to an extent). 9 The crossflow principle works well with the finite element model of a cable (as it necessitates breaking the cable into straight cable segments). The applied drag forces on each of the cable segments can easily be applied to the dynamic equations of motion for each segment. The exact equation for the crossflow principle in its original twodimensional form and modified threedimensional form are presented in chapters two and three. 1.3.2 Lagrange’s Equations Newtonian mechanics states that the motion of a particle can be determined provided that there is a complete knowledge of the external forces acting on the particle and that the initial conditions of the particle are known; the same is true for groups of interacting particles. This is a fundamental step in producing the equations of motion necessary to create a simulation of a system in response to inertial and external forces. The development of the equations of motion of a system can be trivial or quite complex depending on the composition of the system and its environment. The direct application of Newton’s laws becomes difficult as the number of groups of particles needed to represent the system increase. The motion of the system components and internal and external forces acting on the system exist as vectors, with magnitude and direction. Thus, keeping track of the external forces, internal forces, and motion of the system components in vector form is often difficult. A less complex method of deriving system equations of motion is attained by implementing Lagrange’s equations, Greenwood11. 10 The Lagrange method has two main advantages over the direct application of Newton’s Laws in attaining the system equations of motion. The first is that the Lagrange approach is energy based rather than spatially based, thus it deals with scalar values which circumvents the vector bookkeeping involved with Newton’s method. The other benefit of using the Lagrange equations is that, unlike the application of Newton’s laws, it specifies an exact process to develop the equations of motion for any well defined system. After defining the system and its constraints, all one needs to do to develop the system equations of motion is to carryout the manipulation of Lagrange’s equations. An important concept in the application of Lagrange’s equations is degree of freedom of the system. The number of degrees of freedom of a system is equal to the number of coordinates necessary to define the configuration of the system minus the number equations that constrain the position of the system. For example, assume the system is that of a simple twodimensional ideal pendulum. 11 a) b) Figure 1.2 Simple pendulum shown in twodimensional a) Cartesian coordinates, and b) cylindrical coordinates In twodimensional Cartesian coordinates the location of the pendulum mass is determined by its x and z positions. However, the x and z positions are constrained to a circular path around the pendulum hinge with a radius equal to the length of the pendulum (we are assuming the length of the pendulum is rigid). The constraint equation would be the sum of the square of the x and z position is equal to the square of the length of the pendulum. Thus, if the x position is known then the application of the constraint equation will produce the z position of the pendulum mass. With two coordinates utilized to specify the orientation of the pendulum and the existence of one constraint equation the degree of freedom of the system is one. Imagine the same example in cylindrical coordinates. The position of mass of the pendulum is determined simply by the angle q of the pendulum relative to the coordinate 12 frame. There is then one coordinate necessary to define the pendulum orientation and no constraint equation resulting in one degree of freedom for the pendulum. As the pendulum example shows, the degree of freedom of the system is independent of the coordinates chosen to represent the system. The coordinates needed to represent the orientation of the system are deemed generalized coordinates. The generalized coordinates of the simple pendulum in the twodimensional Cartesian coordinate system are the x and z position of the pendulum mass. The generalized coordinate of the same pendulum in cylindrical coordinates is simply the angle q. Therefore, since there are infinite possibilities of coordinate systems that could be used to analyze the pendulum there are infinite possibilities for generalized coordinates. However, in many cases there is a set of generalized coordinates in which the number of coordinates is equal to the degree of freedom of the system. This set is referred to as the independent generalized coordinates. The independent generalized coordinates require no constraint equations; as a result it is often simpler, mathematically, to analyze the system using the set of independent generalized coordinates. In the pendulum example the angle q is the independent generalized coordinate. In a dynamic system the generalized coordinates, independent or not, vary with time and are algebraic variables in Lagrange’s equations. As stated previously, Lagrange’s equations are energy based. Prior to manipulation of Lagrange’s equations the energy of the system must be represented in terms of the generalized coordinates. The manipulation of Lagrange’s equations requires derivation of the system energies with respect to time and the generalized coordinates. The exact form of Lagrange’s equations 13 as well as the execution of Lagrange’s equations will be undertaken in detail in the next two chapters. 1.3.3 Aircraft Dynamics There has been much effort spent by researchers and scientists in the study of aerodynamics as it relates to aircraft. This work ranges from simple analysis of fluid flow over a twodimensional airfoilshape to the development of nonlinear equations of motion for high speed maneuverable aircraft. The extensive work accomplished for a wide variety of aircraft and flight conditions has produced reliable models for the prediction of aircraft response to fluid forces. The method used to represent the dynamics of the aircraft in this system is smalldisturbance theory, Nelson14. Smalldisturbance theory assumes that the motion of an aircraft consists of steadystate motion and small perturbations about its steadystate. Smalldisturbance theory linearizes the rigid body equations of motion. This is accomplished by substituting reference values plus a perturbation for the all variables in the rigid body equations of motion. Many reference values (e.g. pitch, roll, etc.,) can be set to zero by assuming a symmetric reference flight condition. Simplification of the equations of motion after substitution produces products of perturbation. These nonlinearities are ignored in the assumption that higher order small disturbances will produce negligible overall effects. Further simplifications of the equations of motion produce aerodynamic stability derivatives. The stability derivatives represent the change in forces or moments of a body due to some change in the flight conditions, such as forward velocity or angle of attack. 14 The stability derivatives are approximated based upon the size, shape, and mass distribution of an aircraft and its components. Methods for approximation of aircraft stability derivatives are presented in Nelson14 and Raymer16. The equations of motion for aircraft are written in terms of normalized stability derivatives. As a result, these general equations apply to a wide variety of general aviation aircraft over most flight conditions. However, since the equations were linearized using small disturbances, the theory breaks down and thus does not produce accurate predictions when the motion of the system involves large amplitude motions such as aircraft stall or spinning. Given the simple shape of the ATV and its proposed flight conditions the smalldisturbance theory will be adequate in modeling dynamics of the aircraft. 1.3.4 Cable Towed Systems A cable tow system in its most simple form consists of a body constrained by cable in a fluid. There are a wide range of system configurations that meet the requirements of a cable tow system, such as kites, tethered aerostats, mooring lines, towed buoys (above and below the ocean surface), gliders, and towed aircraft targets. Research on many of these systems has been accomplished in the past. Each system has its own nuances that separate one author’s work from another, and although not all are directly applicable to the aerial towed system they help to illustrate the diversity of the system. A few of the papers that do apply to the aerial tow system are discussed below. Huffman and Genin1 research the stability of a simple aerial towed system. Their research suggests that there is no direct instability in the dynamics of the towed cable 15 system. In each case tested the disturbances of the system were damped out in time. They concluded that a likelihood of instability is present, however, when the dynamic frequency of the towed body matches the natural frequency of the cable. The result would be a limit cycle oscillation. Huffman and Genin1 also show that the natural frequency of cable can vary based upon parameters such as cable length and flight speed. As a result, the limit cycle behavior could be present in tow system at particular flight conditions. However, this depends largely on many variables in the system. As a result, there was no generalized parametric study. Cochran et al.6 investigate the effect of controlling an aerial towed vehicle. The authors designed a control system that would augment the stability of the towed vehicle. Simulations were run in conditions that produced unstable dynamic motion in the uncontrolled towed vehicle (perturbed retrieval of the towed aircraft). When the control system was activated the authors showed that the system damped out quickly. This paper addresses the difficulties of designing a controlled aerially towed system and the results suggest that improved stability can be achieved during maneuvers that could otherwise result in instability. Nakagawa and Obata15 studied the longitudinal dynamic modes of the aerial tow system. They analyzed the modes of three flight configurations, a straight cable connected to a sphere, a straight cable connected to an aircraft, and a curved cable connected to an aircraft. The results of the study helped to classify the types of motion of the aerial tow system. In all of the configurations the main dynamic modes were the pitching mode, pendulum mode, and first vibration mode. The pitching mode is simply 16 the pitching of the towed body; this has little effect on the cable dynamics or the translation of the towed body. The pendulum mode oscillates the towed body about its static configuration. The motion causes slight bending in the cable and little if any towed body pitching. In both of the straight cable configurations the first vibration mode has little effect on the position or the orientation of the towed body, and simply causes oscillations in the cable. However, in the curved cable configuration the first vibration mode causes the towed aircraft to surge fore and aft in an elliptical pattern. The authors labeled this as the bowing mode. The authors complete a parametric study of the aircraft longitudinal stability derivatives and found that two of them, Zw and Mw as well as the cabletowed body hitch point played an important roll in the stability of the system. The results sho wed a strong relationship between the stability derivatives and the bowing mode. As the absolute values of Zw and Mw increase and decrease, respectively, the bowing mode becomes less stable, and eventually an unstable phenomenon known as pitching or bowing flutter occurs. Nakagawa and Obata15 noticed a similar behavior as the cable hitch point on the towed aircraft was moved more and more fore of its center of gravity. The system dynamic oscillations are of the utmost importance when attempting to fly the towed body very near the ocean surface. The research completed by Nakagawa and Obata15 suggest care should be taken when designing the ATV as its stability derivatives have a direct effect on unwanted system oscillations. Lawhon and Arena19 investigated the effects of aircraft size and shape as well as cable parameters to develop the optimum 17 system configuration for aerially towed systems tracking very near the ocean surface. The ATV used in this paper was a result of the work done by Lawhon and Arena19. 18 CHAPTER 2 2TWODIMENSIONAL DYNAMICS 2.1 Introduction It is the focus of this paper to derive and produce the mathematical equations that accurately describe the threedimensional motion of the aerially towed system. Followed by, implementing the equations of motion to produce a dynamic simulation of the aerially towed system. The simulation will then be used to draw general conclusions about the system and its ability to perform the maneuvers desired to obtain atmospheric measurements near the oceans surface. The individual dynamic behavior of the system components and the interactions between them will result in quite complex set of equations of motion. As a result of this complexity, it is easy to get bogged down by the shear volume of the equations and lose sight of the steps taken to derive the final set of equations. This is especially true in reference to the threedimensional system. Therefore, the twodimensional equations of motion will be laid out in this chapter. The twodimensional equations of motion will also act to verify that the more complex threedimensional equations of motion are derived without error. It is important to note that the same principles will be held for the derivation of both the two and threedimensional systems although they may differ in clarity. 19 The basis of any simulation is the modeling of the system components. The quality of the models directly affects the value of the simulation. The derivation of the equations that dictate the dynamics of a component is possible upon the completion of its mathematical model. To accurately model the system, the individual component dynamics are important however, the modeling of the system necessitates the equations of motion represent the dynamic interaction between the components as well. Thus, in general, the derivation of the equations of motion of a system component must be accomplished simultaneously with the dynamic equations of the remaining objects. The result will be a set of equations that represent the dynamics of the entire system, i.e. if the host vehicle oscillates the system equations of motion will reflect the oscillations effect on the dynamics of the cable and the ATV. The following sections discuss the modeling of the system components and the development of the individual dynamic equations of each component. As stated previously, the individual component dynamics are not enough to model the system; however, they may be of academic interest, and will lead to a clear understanding of the derivation of the system equations of motion. 2.2 Cable Dynamics As discussed in the literature review there has been extensive study on cable dynamics and cable modeling for a wide range of applications. Due its versatility the finite element technique will be used to model the connecting cable. The main drawback of the finite element technique is its computational expense. The great strides in computer speed over the years helps to soften this blow. There may also be situations when the operator is not 20 interested in a completely rigorous simulation but rather just a feel of the system’s behavior, this can be accomplish quickly by using fewer segments. As mentioned previously, the finite element technique requires that the cable be broken into an arbitrary number of segments which are free to rotate about their neighboring segments. The types of segments used to model the cable can vary. Notable models include the simple pendulum, spring mass, thin rod, and curved structure model. Figure 2.1 Finite element cable segment models In the simple pendulum model, the segment is straight and inextensible. The mass of the segment is lumped to a single point at the tip of the segment, resulting in a model that resembles that of a simple pendulum. Being that the segment is straight the crossflow method can be used to predict the drag force on the segment. The applied force on the segment is exerted at the point mass or node. The inextensible constraint lends this model to a system in which the strain on the cable is minimal. The benefit of this technique is its simplicity. The draw back is that placing the inertial and applied forces at the tip of the segment only becomes a good approximation when the number of segments used is large. 21 The spring mass segment also lumps the mass at tip of the segment, however the segment length acts as a spring. This allows for strain effects to be incorporated into the equations of motion of the cable. The change in length of each segment adds an additional degree of freedom to the system, which increases the complexity of the equations of motion and increase computational intensity. If the strain on the cable can be neglected the spring mass segment model expends much unneeded effort and time. The thin rod model is similar to the lumped mass model. It too is straight and inextensible. However, the model places the inertial forces and applied forces at the center of mass of the segment. The additional complexity of this model is that the mass of the rod is not lumped to a point; therefore, it has inertia that impedes segment rotation and must be taken into account the derivation of the equations of motion. The benefit of the better physical model is it will take fewer segments to accurately represent the motion of the continuous cable than the lumped mass model. The last model represents the cable segment as a curved structure. This technique models the segments as polynomials or splines which require slope continuity at the segment connections. Due to the fact that most cable configurations are curved this method reduces the number of segment needed to attain an accurate simulation. However, the evaluation of the drag on a curved cable is not well known. The more complex drag model and polynomial or spline calculations make derivation of the cable equations of motion difficult. Neglecting the strain in the cable is acceptable as we are dealing with a fairly light ATV and a high tensile steel connecting cable, thus the spring mass model unnecessary; 22 and given that there is a historically accepted method to model the drag on a straight cable the complexities of the curved structure will be avoided. The equations of motion of the cable will be derived using the simplest segment model, the pendulum model, shown in figure 2.1. However, a discussion about the changes to the cable equations of motion as a result of modeling the cable segments as thin rods will take place at the end of this chapter. Figure 2.2 Lumped mass finite element cable model 2.2.1 Coordinate Systems In this section the primary concern is the development of the twodimensional equations of motion for a flexible cable. As a result, the host aircraft and the ATV will be neglected for now. The cable is broken into n arbitrary segments in accordance with the finite element cable model. Therefore, one end of the cable is fix to a point in space and is free to rotate while the other is end is constrained only by its preceding segments. The coordinate system at the cable attach point is the inertial coordinate system. This coordinate system is orthogonal. The X axis is locally parallel to the ‘ground’ and the Z 23 axis is pointing toward and is normal to the ‘ground’. The derivation of the equations of motion of the flexible cable will be carried out in the inertial frame. As a result, it is necessary to be able to describe the orientation of the cable in this coordinate system. The location of the ith cable node in the inertial frame can be resolved using the length, li, and angle of orientation, qi, of the segment and each preceding segment, as shown below. å å = = = = i j i j i j i j j j Z l C X l S 1 1 q q (i = 1,2Kn) (2.1) where, lj is the length of the jth segment and Xi and Zi are the position of the ith segment in the inertial frame. Also note C and S are short notation for cosine and sine, this notation will be used throughout the remainder of this paper. Aside from the inertial frame, each cable segment has its own coordinate system. The cable segment coordinate system is orthogonal and is configured such that the node of that segment is located in the positive zdirection a distance equal to that of the cable segment. For example, the origin of the coordinate axis for the first cable segment is located at the exact same point as that of the origin of the inertial frame. The two frames differ only by the angle rotation of the segment, q1. When the angle is equal to zero the two frames are identical. In a similar manner, the origin of the second cable segment is located at the node of the first cable segment, and the first and second frame differ rotationally by the difference of q2 and q1. This is evident from Figure 2.2. 24 Figure 2.3 Twodimensional inertial and cable segment coordinate systems Throughout the analysis of the system it will be necessary to convert from the inertial frame to the segment frame. This is accomplished by utilizing a rotational coordinate transformation. The rotational coordina te transformation from the ith cable segment frame to the inertial frame is represented by the directional cosine matrix as follows, úû ù êë é  = i i i i i S C C S C q q q q q (2.2) The inverse of the directional cosine matrix produces a rotational coordinate transformation from the inertial frame to the ith cable segment frame. Calculating the inertial position of the node of the ith cable segment is a common implementation of the rotational transformation. Recall that in the cable segment frame the node is located at a distance equal to the length of the segment in the zdirection. Using the directional cosine matrix the position of the ith cable segment node in the inertial coordinate system is as follows, 25 å= úû ù êë é úû ù êë = é úû ù êë é i i j j i l C Z X j 1 0 q (i = 1,2,K, n) (2.3) Note that equation 2.3 and 2.1 are equivalent. 2.2.2 Lagrange Equations As discussed earlier, the connecting cable is broken into n arbitrary number of segments, and each segment has its own coordinate system. Deriving the equations of motion for the cable using the direct application of Newton’s laws of motion becomes increasingly difficult as the number of segments used to model the cable increases. This is primarily due to the fact the Newton’s laws of motion are vector based and as the system increases the number of coordinate systems increase resulting in labor intensive derivation. Another difficulty when using the Newtonian approach is that there is no general method to derive the equations of motion. As discussed in the previous chapter, an alternative to using Newton’s laws of motion is the utilization of Lagrange’s equations, which is energy based. Consequently, it circumvents many of the complexities of vectorbased derivation, as well as, gives a standard process for developing the system equations of motion. To recap, the derivation of Lagrange’s equations depends on the degree of freedom of the system, and the degree of freedom of the system is the number of coordinates needed to define the system minus the number of constraint equations of the system. In the Cartesian coordinates of the inertial frame for a cable broken into n segments there are 2n coordinates needed to define the location of the node of each segment. However, each segment has a constraint equation of the form, 26 2 1 1 2 1 1 2 ÷ ÷ø ö ç çè æ  + ÷ ÷ø ö ç çè æ = å å =  = i j i j i j i i j l X X Z Z (i = 1,2,K, n) (2.4) As a result, 2n coordinates and n constraint equations produce n degrees of freedom. In other words, the twodimensional finite element cable has a degree of freedom equal to the number of segments used to model the cable. The coordinates used to define the system are the generalized coordinates. The generalized coordinates are general due to the fact that there are an infinite number of coordinate systems that could be used to define the system all with varying number of constraint equations (although every system has n degrees of freedom). Consequently, there are no specific coordinates that need to be used to define the system. The implementation of Lagrange’s equations will be simplified if the independent generalized coordinates are used (recall, the independent generalized coordinates are the set of coordinates that define the system yet have no constraint equations). For the flexible cable system the independent generalized coordinates are the orientation or attitude angles of each of the segments, qi (i = 1, 2… n). For the duration of the chapter the attitude angles of the cable segments will be used and referred to as the generalized coordinates of the system. The generalized coordinates will be generically referred to by the variable q. As stated previously, Lagrange’s equations are energy based. The Lagrangian function, L, represents the energies of the system and is equivalent to the kinetic energy minus the potential energy of the system, 27 L = T U (2.5) where, T and U are kinetic and potential energy functions, respectively. Lagrange’s equations are as follows, i i i Q q L q L dt d = ¢ ¶ ¶  ÷ ÷ø ö ç çè æ ¶ ¶ & (i = 1,2,K, n) (2.6) where, i Q¢ represents the applied, or generalized, forces that are not derivable from a potential function, such forces include friction forces and forcing functions. However, since the potential function for this system is not velocity dependant Lagrange’s equations can be reduced to the following, i i i i i Q q U Q q T q T dt d = ¶ ¶ = ¢  ¶ ¶  ÷ ÷ø ö ç çè æ ¶ ¶ & (i = 1,2,K, n) (2.7) where, Qi are the generalized forces that include applied forces and the inertial forces derived from the potential function. Note that Equations 2.6 and 2.7 is actually a set of n equations, one equation for each generalized coordinate of the system. Thus, in the case of the twodimensional finite element modeled cable there is one equation for each attitude angle of the n cable segments. Therefore, Lagrange ’s equations can be rewritten as, i i i Q T T dt d = ¶ ¶  ÷ ÷ø ö ç çè æ ¶ ¶ q& q (i = 1,2,K, n) (2.8) 28 In order to implement Lagrange’s equations it is necessary to get the energies of the systems in terms of the generalized coordinate. Only then can the proper derivation of the equations of motion take place. The kinetic energy of the system at any given time is due to the summation of the rotational and translational motion of the inertial components of the system. For the finite element cable the inertial components are simply the segment nodes. Thus, the kinetic energy of the cable is the summation of the nodal rotational and translational kinetic energy. However, since the mass of a segment is localized to a point it has no rotational inertial and consequently has no rotational kinetic energy. As a result, the kinetic energy of the system is solely a function of the translational velocity of each node. å= = n j j j T m V 1 2 2 1 (2.9) The velocity components of a node can be found by simply taking the time derivative of its position. å å = = =  = i j i j j i j i j j j j Z l S X l C 1 1 q q q q & & & & (i = 1,2,K, n) (2.10) Equation 2.10 shows that the translational velocity of a node results from the rate of rotation of its segme nt and each preceding segment. At this moment, the velocity of each node is in vector form. However, squaring the velocity for the kinetic energy releases 29 any further vector bookkeeping. The square of velocity, V, for the kth segment is can be written in series form as follows, ( ) ( i j i j ) ( i j ) V X Z l l C C S S l l C j i j k j k i j i j i k j k i k k k i q q q q q q q q q q  = = = = = & + & =åå & & + =åå & & 1 1 1 1 2 2 2 (2.11) Similarly, the total kinetic energy for n cable segments is written in summation form as, å ååå ( ) =  = = = = = n k j i j k j k i k i n k k k i j T m V m l l C 1 1 1 1 2 2 1 2 1 q q q&q& (2.12) Being that the kinetic energy has been derived in terms of the generalized coordinates the implementation of Lagrange’s equations can commence. The derivative of kinetic energy with respect to the ith generalized coordinate and the time rate of change of the ith generalized coordinate yields the following, åå ( ) = =  = ¶ ¶ n k k j k i j i j i j i m l l S T 1 1 q q q q q & & (2.13) åå ( ) = =  = ¶ ¶ n k k j k i j j i j i m l l C T 1 1 q q q q & & (2.14) Finally evaluating the time rate of change of equation 2.14 produces, åå [ ( ) ( ) ( )] = =     = ÷ ÷ø ö ç çè æ ¶ ¶ n k k j k i j j j j i i j i j i m l l C S T dt d 1 1 q q q q q q q q q & && & & & (2.15) Combining equations 2.13 and 2.15 results in the series representation of the left hand side of the Lagrange’s equations for the twodimensional finite element modeled cable, 30 åå ( ( ) ( ) ) = =   =  ¶ ¶  ÷ ÷ø ö ç çè æ ¶ ¶ n k k j k i j j j i i j i j i m l l C S T T dt d 1 1 2 q q q q q q q q & & & (i = 1,2,K, n) (2.16) The right hand side of equation 2.8 represents the generalized forces of the system. Note that due to the fact that the generalized coordinates of the system are angles the general forces are actually moments. These moments are about the pivot point of each segment and result from the inertial and applied forces on the node of the segment as well as any applied moments. The generalized force equation is derived from the concept of virtual work. Suppose forces F1, F2, …, Fk are applied at position x1, x2, …, xk and moments M1, M2, …, Mk are applied at angles f 1, f 2, …, f k. Now let us assume that the system undergoes arbitrary small displacements dx1, dx2, …, dxk and df 1, df 2, …, df k. The work done by these applied forces and moments is known as virtual work and the small displacements are known as virtual displacements, å== + k j j j j j W F x M 1 d d dj (2.17) The virtual displacements can be rewritten in terms of the virtual displacement of the generalized coordinates, å= ¶ ¶ = n i i i j j q q x x 1 d d (2.18) å= ¶ ¶ = n i i i j j q 1 q d j dj (2.19) 31 The relationship between the generalized forces and the virtual work of the system due to the virtual displacement of the generalized coordinates is as follows, åå =å ¶ ¶ + ¶ ¶ = = = n i i i k j n i i i j i j i j j q Q q q q M q x W F d d j d d 1 1 (2.20) Simplification of equation 2.20 produces the generalized force equation for the applied forces and moments in terms of the generic variables x and f . å= ¶ ¶ + ¶ ¶ ¢ = n i i i j i j i j i j q q q M q x Q F 1 d j d (2.21) In the case of the twodimensional cable segment system the position variables are X and Z while the angles are q. Incidentally, the generalized coordinate for this system is also q. The equation for the generalized forces not produced by the potential function is as follows, å= ¶ ¶ + ¶ ¶ + ¶ ¶ ¢ = n j i j i j Z j i j i X j i M Z F X Q F 1 q q q q q (i =1,2,K,n) (2.22) The total generalized force equation is made up of the external forces and moments as well as the inertial forces, å= ¶ ¶ + ¶ ¶ + ¶ ¶ + ¶ ¶ =  n j i j j i j Z j i j X j i i M Z F X F U Q 1 q q q q q q (i =1,2,K,n) (2.23) where, X j F and Z j F are applied forces and j Mq is the applied moment of the jth segment in the inertial coordinate frame. 32 In the expansion of equations 2.23 let us look first at the term that is the derivative of potential energy with respect to the ith generalized coordinate. The potential energy of the cable is the summation of the potential energy of each segment node. åå ( ) = = =  n k k i U mk g Zi D 1 1 (2.24) where, D is the datum line from which potential energy is measured. Note that since potential energy is being differentiated the location of the datum line is inconsequential. Taking the derivative of equation 2.24 with respect to the ith generalized coordinate produces, å= ¶ ¶ = ¶ ¶ n j i j j i Z m g U q 1 q (2.25) The resulting equation allows the inertial force to be included with the zdirection applied force. Simplify the generalized force equation to that shown below, ( ) å= ¶ ¶ + ¶ ¶ + + ¶ ¶ = n j i j j i j Z j j i j i X j M Z F m g X Q F 1 q q q q q (i =1,2,K,n) (2.26) The derivative of position and angle with respect to the ith generalized coordinate is as follows, 33 ( ) ( ) ( j i) l S j i Z l C j i X i j i i j i i j i i = = ¶ ¶ =  ³ ¶ ¶ = ³ ¶ ¶ 1 q q q q q q (2.27) Note that the position of the node of a segment is dependant upon that segments attitude angle and the attitude angle of each preceding segment. Due to the nature of the dependence of position on the generalized coordinate, the evaluation of the derivative of the position of the jth segment with respect to the ith generalized coordinate results in zero when i is greater than j. The derivative of angle with respect to the ith generalized coordinate is equal to one when i is equal to j and zero otherwise. Thus, ( ) å= = +  + n j i i i X j i i Z j j i i Q M F l C F m g l Sq q q (i =1,2,K,n) (2.28) The complete derivation of Lagrange’s equations for the twodimensional finite element modeled cable is attained by combining equations 2.16 and 2.28. The resulting set of n equations of motion is presented as follows in summation form. ( ) ( ) ( ) ( ) (i n) M F l C F m g l S m l l C S n j i i X j i Z j j i n k k j k i j j j i i j i j i 1, 2, , 1 1 2 K & & = +  +  = å åå = = =   q q q q q q q q q (2.29) The set of equations of motion can be more clearly presented in matrix form, yielding, 34 [A][q&]+ [B][q&2 ]= [Q] (2.30) where, A and B are coefficient matrices and are size n by n. The remaining three matrices are column matrices of size n by 1. The values of the coefficient matrices are as follows, i j i j i j A l l C m , (q j qi ) , = (2.31) i j i j i j B l l S m , (q j qi ) , =  (2.32) where, the subscripts i and j represent the row and column of the appropriate matrix, and i j m , represents the mass distribution for row i and column j. For the lumped mass model of the n segments the mass distribution is simply a sum of segment masses. å = = n k i j i j k m m max( , ) , (2.33) For example, the equations of motion for a cable modeled using two cable segments without fluid drag would be as follows, ( ) ( ) ( ) ( ) ( ) ( ) úû ù êë é   + = úû ù êë é úû ù êë é   + úû ù êë é ú úû ù ê êë é +     2 1 1 2 2 1 1 2 2 1 2 2 1 2 1 2 2 2 1 2 1 2 2 1 2 2 1 2 2 1 2 2 2 2 1 2 2 1 2 1 0 0 q q q q q q q q q q q q q q m gl S m m gl S m l l S m l l S m l l C m l m m l m l l C & & && & (2.34) These two equations are analogous to the equations of motion of a twodimensional ideal double pendulum. 35 The resulting fully nonlinear twodimensional equations of motion describe the dynamics of the system due to inertial and applied loads, as well as, the dynamic interaction between all segments; and the set of equations place no limitations on the number of segments of the system (while the equations hold true independent of the number of segments the effects on simulation computation varies, this will be discussed in chapter 5). The angular acceleration of each segment can be calculated given appropriate initial conditions, such as segment attitude angle and angular velocity of the system. The implementation of a numerical integration routine on the set of ordinary differential equations allows the prediction of cable motion for various initial conditions, and with the addition of a drag model to predict the applied friction forces on the segments the simulation can be extended to that of a cable in fluid flow. 2.2.3 Aerodynamic Forces The generalized force matrix in the equation 2.30 is composed of the inertial forces and applied forces of each of the cable segment. These applied forces are the friction forces produced by the fluid interaction with the cable. Thus to fully derive the equations of motion for a cable in a fluid a drag model is necessary. As stated previously, the study of cable drag has a history that dates back to the early aircraft and thus has been well documented. The cable crossflowprinciple, Hoerner13, is a drag model that was presented in the sixties and has since been used as a basis for numerous research papers involving fluid immersed cables. Its limitations require the cable crosssection be circular, the cable be straight, and fluid be at a subcritical Reynolds number. In general, aircraft tow cables are 36 circular in crosssection; due to our relatively low airspeeds the Reynolds number is subcritical; and our finite element model for the cable requires that the cable be broken into straight segments. Therefore, the Hoerner crossflowprinciple applies very nicely to our model of the cable system. The crossflowprinciple calculates the tangential and normal forces on a straight cable segment produced by tangential and normal velocities about the cable, as shown below, ( ) ( 2 2 ) 2 2 2 1 2 1 j j j j j j j j z j z f x z x j x f x z p x F dl V c V V F dl V c V V c V =  + =  + + r p r p (2.35) where, r, d, cf, and cp are air density, cable diameter, frictional coefficient, and pressure coefficient, respectively, and Fx, Fz, Vx, and Vz are the aerodynamic forces and the flow field velocity in the cable segment coordinate system, respectively. The cable velocity in the cable segment frame can be attained simply by using the rotational coordinate transformation of the relative velocity of the segment node in the inertial frame, as shown below, ú úû ù ê êë é + + úû ù êë é = úû ù êë é j j T z j x j W Z U X C V V j & & 0 0 q (2.36) Presumably the freestream velocities, U0 and W0, are known and the equations for calculating the velocity of a segment node in the inertial frame were derived in equation 2.10. 37 Note that the drag forces produced by the drag model are in the cable frame yet the forces in the generalized coordinate equation are in the inertial frame. As a result, the calculated normal and tangential drag forces must go through a rotational coordinate transformation before application can be made. úû ù êë é úû ù êë = é úû ù êë é z j x j Z j X j F F C F F j q (2.37) The implementation of the Hoerner drag model produces the mathematical model of the applied loads on a cable segment. The generalized forces for each segment can be calculated based on cable motion and numerous freestream conditions and thus allows for the dynamic simulation of a cable in a fluid. 2.3 ATV Dynamics The ATV alone is an aircraft, and a very simple one at that. As a result, its rigid body dynamics could be easily modeled due to the extensive research and documentation done in the area aircraft flight dynamics. However, in the aerial towed system the dynamics of the cable and the dynamics of the ATV are intertwined, consisting of individual dynamics and the interaction of the two components. As a result, a system of equations containing both the equations of motion of the cable and that of the ATV must be derived and solved simultaneously. This is most easily accomplished by reproducing the ATV as an additional segment, one which is attached to nth cable segment. Thus the equations of motion for both the cable and ATV can be derived using Lagrange ’s equations. Since the equations of motion of the cable were derived to be independent of the number of segments this does not change the form of the equations previously derived. The 38 variation in the set of equations arises because the ATV segment model differs from the cable segment model. The set of equations for the cable and ATV and there variations from the previous cable dynamic set will be laid out in the subsequent sections. 2.3.1 Coordinate system Modeling the ATV as an additional segment adds an extra coordinate system just as adding an extra cable segment to the cable model. This system is analogous to the cable coordinate system described previously, the origin is the cableATV hitch point and its zaxis passes through the center of mass of the ATV. However, there is another coordinate system necessary to fully describe the orientation of the ATV. This coordinate system corresponds to the standard aircraft body fixed coordinate system. The aircraft body fixed coordinate system has it origin at the center of gravity of the aircraft, and for the twodimensional case the aircraft body fixed xdirection is pointing out of the nose of the aircraft and its zaxis is pointing out of the bottom of the aircraft (presumably towards the ground). The two additional axes are present in figure 2.3. Figure 2.4 Twodimensional ATV coordinate system 39 Thus, the ATV segment coordinate system and the ATV aircraft coordinate system differ directionally by a fixed rotation angle, g. This angle is geometrically defined by the shape of the ATV as well as the location of the cableATV hitch point, figure 2.3. The ATV aircraft coordinate system is necessary to calculate the aerodynamic forces and moments of the ATV using small disturbance theory. The ATV orie ntation, as seen by the inertial frame, is determined by the ATV segment angle and ATV geometric attitude angle. Thus two rotational coordinate transformations are needed to describe the ATV attitude in the inertial frame. The rotational transformation from the ATV aircraft coordinate system to the ATV segment coordinate system is as follows, úû ù êë é ¢ ¢ úû ù êë = é úû ù êë é z x C z x g (2.38) where, x’ and z’ are the x and zdirections in the aircraft fixed frame. The rotational transformation from the ATV segment coordinate system to the inertial coordinate system is as follows, úû ù êë é úû ù êë = é úû ù êë é + z x C Z X qn 1 (2.39) As a result, the equation for rotational transformation from aircraft fixed coordinate system to the inertial coordinate system is found by combining equations 2.38 and 2.39. úû ù êë é ¢ ¢ úû ù êë é úû ù êë = é úû ù êë é + z x C C Z X qn 1 g (2.40) 40 This transformation is necessary in converting the aerodynamic loads to the generalized forces of the ATV segment (due to the fact that our generalized force equations require the applied loads be in the inertial system). 2.3.2 Lagrange equations The parameters of the additional segment, which will be henceforth deemed the ATV segment, is dictated by the size and shape of the ATV. The length of the segment is actually the distance between the cableATV hitch point and the center of gravity of the ATV. The mass of the ATV segment is simply the mass of the aircraft itself. Therefore, the equations of motion for the translation of the ATV segment require no additional derivation from that accomplished in the cable dynamics section. As a result, the complex dynamic interaction between each cable segment and the ATV are inherently represented by the equations of motion. However, the equations of motion for the ATV segment are not complete. Recall that the nodes used to model the cable are point masses, thus they have no inertia to resist rotation. The same cannot be said for the ATV. To complete the equations of motion for the cableATV system the Lagrange equation must be derived for the rotational kinetic energy of the ATV segment The complete representation of Lagrange’s equations for the cable and ATV system is as follows, ( ) ( ) i i rot i rot i i i i rot i rot Q T T dt T T d dt d Q T T T T dt d = ¶ ¶  ÷ ÷ø ö ç çè æ ¶ ¶ + ¶ ¶  ÷ ÷ø ö ç çè æ ¶ ¶ = ¶ ¶ +  úû ù êë é ¶ ¶ + q q q q q q & & & (i = 1,2,K, n + 1) (2.41) 41 Due to the fact that the total kinetic energy of the system is composed of the summation of translational and rotational kinetic energy the resulting equations of motion due to the rotational kinetic energy can simply be superimposed upon the previously derived equations of motion (i.e. the derivation of the components of the equation 2.41 that represent the translational kinetic energy need not be carried out as they will yield the same results). The analysis of the rotation kinetic energy will focus on the following portion of equation 2.41, ( ) i rot i Trot T dt d q ¶q ¶  ÷ ÷ø ö ç çè æ ¶ ¶ & (i = 1,2,K, n + 1) (2.42) The equation for rotational kinetic energy of a rigid body is as follows, T {w}T [I ]{w} rot 2 1 = (2.43) where, I is the inertial matrix and w is the angular velocity vector of the rigid body. For the twodimensional case the segment rotational equation of motion can be simplified since the only allowable rotation is about the yaxis. The rotational kinetic energy of the system is the sum of the rotational kinetic energy of each segment. The rotational kinetic energy of the twodimensional system is as follows, å+ = = 1 1 2 2 n 1 i rot i yi T q& I (i = 1,2,K, n + 1) (2.44) where, Iy is the rotational inertial of the node about its yaxis. 42 The derivation of Lagrange’s equations for rotational kinetic energy term is quite straightforward, producing, = 0 ¶ ¶ i rot T q (2.45) i yi i I T q q & & = ¶ ¶ (2.46) Finally evaluating the time rate of change of equation 2.46 produces, i yi i I T dt d q q & & & = ÷ ÷ø ö ç çè æ ¶ ¶ (2.47) Combining equations 2.45 and 2.47 produces the supplementary left hand side of Lagrange’s equations due to the rotation of the each segment for the twodimensional finite element modeled cable, i yi i i I T T dt d q q q & = && ¶ ¶  ÷ ÷ø ö ç çè æ ¶ ¶ (i = 1,2,K, n + 1) (2.48) The combination of the rotational and translational components of the equations of motion and the addition of the ATV segment result in a set of equations of the form shown in equation 2.30. The inclusion of the ATV segment into the system adds an extra degree of freedom and thus an extra equation to the set of the equations of motion of the system, resulting in the coefficient matrices increasing to size n+1 by n+1 and the column matrices to size n+1 by 1. The additional column in the coefficient matrices represents the effect of the motion of the ATV segment on the dynamics of the cable segments. The 43 additional row of the coefficient matrices represents the equation of motion of the ATV segment, which contains the effects of the dynamics of each cable segment. The changes to the set of equations derived for the translation of the nodes of the segments are in the terms of the A coefficient matrix. The corrected terms of the A matrix are as follows, i j y i j i j i j A I l l C m , i , (q j qi ) , d  = + (2.49) where, the kronecker delta, di,j, is zero at all time except when i equals j at which time its value is unity. Thus the inertia term only appears on the diagonal of the A matrix, and since the segment nodes, in accordance with the chosen segment model, have a zero value for inertia the only change in A is at n+1, n+1. 2.3.3 Generalized Forces The generalized force equation derived in equation 2.28 applies to all segments including the newly added ATV segment. However, the crossflow method does not apply to aircraft. As a result, a method for modeling the aerodynamic forces and moments on the aircraft are necessary. As stated in chapter 1, there is no shortage of research on the aerodynamics of general aircraft. The technique to be used in this paper is the linearized smalldisturbance theory, Nelson14. In short, linearized smalldisturbance theory simplifies the motion of the aircraft to perturbations of the aircraft about some fixed reference flight condition, such as typical cruising conditions. Dealing with perturbations allows for the elimination of higher order terms in the aircraft equations of motion, thus linearizing them. The aerodynamics of the 44 aircraft in the equations of motion are represented by stability derivatives. Stability derivatives approximate the change in aerodynamic forces or moments due to the deviations from the reference flight condition. These stability derivatives are based upon the size, shape, and mass distribution of the aircraft, and apply to a wide range of aircraft. The method does breakdown when the aircraft enters a flight regime in which the nonlinearities are important, such as stalling and spinning. The twodimensional equations for the aerodynamic forces, X F¢ and Z F¢ , and moment, M, in the aircraft fixed frame are as follows, ( ) ( ) M M I (M u M w M w M q) F Z m Z u Z w Z w Z q F X m X u X w Y u w w q Z u w w q X u w = + + + + ¢ = + + + + ¢ = + + & & & & 0 0 0 (2.50) where, m and IY are aircraft mass and longitudinal rotational inertia. Also, u and w are the xdirection and zdirection perturbation velocities and q is the pitch rate in the ATV fixed frame. Xu, Xw,…, Mq are the stability derivatives. (Note that notation in equation 2.50 is commonly used when discussing aircraft dynamics. The pitch rate q is not to be confused with the generalized coordinate q which is itself commonly used in the notation of Lagrange’s equations.) In the aerial tow system the reference flight condition for the ATV is that of the host aircraft, and the host aircraft reference flight conditions can simply be viewed as the freestream conditions, Uo and Wo. As a result, the perturbation of the ATV about the reference conditions is simply a result of the attitude of the ATV and the motion due to segment rotation. As we have already derived the velocity of a node due to rotation in 45 the inertial frame, a coordinate transformation from the inertial frame to the ATV aircraft frame can be utilized to develop the equation for perturbation velocity, úû ù êë é  úû ù êë é + + úû ù êë é úû ù êë é = úû ù êë é + + + 1 0 1 1 W U W Z U X C C w u o o n o n T T n & & g q (2.51) With the perturbation parameters in the ATV aircraft frame the aerodynamic forces and moments on the aircraft can easily be attained. To calculate the generalized forces produced by the aerodynamics of the ATV the forces must be in the inertial frame. This requires the following rotational coordinate transformation, úû ù êë é ¢ ¢ úû ù êë é úû ù êë é = úû ù êë é + Z X Z X F F C C F F qn 1 g (2.52) No rotational transfo rmation is necessary for the twodimensional applied moment. Thus, the aerodynamic moment, M, is directly translated to the ATV segment moment Mq. The forces and moments are applied to the generalized force equation, equation 2.28, following the translation. The dynamics as well as the applied aerodynamic forces and moments of the cable ATV system have been modeled. As a result, it is possible to simulate the effect of the cable and ATV system due to various initial conditions. However, to fully analyze the aerial towed system one last component is necessary, the host aircraft. 2.4 Host Vehicle Up to this point, the model of the system has the cable, immersed in a flowing fluid, connected to a fixed point in space. In reality, the cable will be connected to the host 46 aircraft, and while the ideal case has the host vehicle flying at a constant speed and altitude it is likely that, for one reason or another, the position of the host vehicle will be disturbed. The host vehicle is an airplane, and as a result its dynamics are dictated by the aerodynamic equations presented in equation 2.50, and as it is the towing vehicle it is affected by the motion of the cable and ATV segments. On the other hand, the motion of the cable and ATV will also be affected by the host vehicle motion. In general, the host vehicle is much larger than the ATV and as a result its aerodynamic loads will overpower the loads produced by the motion of the cable and ATV. Thus, a simplifying assumption can be made that the host vehicle is unaffected by the dynamics of the cable and ATV, yet the cable and ATV are affected by the motion of the host vehicle. This assumption allows the host vehicle to be modeled as a simple point. The motion of the host aircraft is independent of the cable and the ATV dynamics and consequently the movement of this point can be chosen arbitrarily. For instance, the equations of motion for longitudinal oscillation in the zdirection of the host aircraft frame can be written simply as sin( ) cos( ) sin( ) z A 2 t z A t z A t o o o o o o o o w w w w w =  = = & & (2.53) where, z is the host vehicle vertical displacement, t is time, Ao is amplitude, and wo is frequency. The amplitude and frequency can be chosen to represent the longitudinal dynamics of a wide range of host aircraft. The effect of host aircraft oscillation on the vertical position of the ATV is of critical importance when the ATV is very near the 47 ground; however the aircraft is not limited to simply vertical motion. The equations of motion can be derived for horizontal motion, vertical motion, or any combination of the two. The technique for modeling the cable and ATV dynamics with the introduction of the varying host vehicle position remains more or less the same. In the previous models the cable connecting point always coincided with the origin of the inertial coordinate system. However, with the introduction of the movable host vehicle this is not so. Any movement of the host aircraft results in the translation of the cable connection point. The derivation of Lagrange’s equations for the cable and ATV hinges on the location and velocity of the segment nodes. Thus, to account for host vehicle movement the inertial position equation must be altered to represent the translation of the cable connection point. The resulting equation is as follows, å= úû ù êë é úû ù êë é + úû ù êë é = úû ù êë é i i j j i l C z x Z X j 1 0 q (2.54) where, x is the difference in the xdirection of the host aircraft in the inertial frame. Due to the fact that the position of the host aircraft is a function of time the derivation of Lagrange’s equations will produce a different outcome from that derived previously, the difference being the effect of the host vehicle oscillation. The derivation of Lagrange’s equations begins by taking the time rate of change of position, 48 å å = = = +  = + j i j i i j i j i i i i Z z l S X x lC 1 1 q q q q & & & & & & (2.55) The square of velocity with the incorporation of a moving host vehicle is as follows, ( ) ( ) ( ) i i i j z x x lC z l S l l C V X Z j i j k j k i i k i i i i i k k k q q q q q q q q  = = = + +å  +åå = + = & & & & & & & & & & 1 1 1 2 2 2 2 2 (2.56) Thus, the total kinetic energy for the aerial tow system is as follows, å å( ) åå ( ) å + =  = = = + = ÷ ÷ ø ö ç ç è æ + +  + = = 1 1 1 1 1 2 2 1 1 2 2 1 2 1 n k j i j k j k i i k i k i i i i n k k k i i i j m z x x l C z l S l l C T m V q q q q & & &q& &q& q&q& (2.57) Note that the last term in equation 2.57 is the same as that derived in the cable dynamics section, equation 2.12. Thus, the derivation of Lagrange equations for the last term will yield the same results. As a result, we can focus on the host vehicle oscillation terms when deriving the Lagrange equations and simply superimpose the results onto the equations of motion for the cable and ATV. The kinetic energy due to the host vehicle oscillation, Thost, is as follows, å å( ) + = = ÷ø ö çè æ = + +  1 1 1 2 2 2 2 2 1 n k k i osc k i i i i i i T m z x x l C z l Sq q & & &q& &q& (2.58) 49 Note that the first two terms in the equation 2.58 are only a function of time. Therefore, they will disappear during the differentiation of the kinetic energy with respect to the generalized coordinates. The derivative of kinetic energy with respect to the ith generalized coordinate and the time rate of change of the ith generalized coordinate yields the following, ( ) å+ = =  + ¶ ¶ n 1 k i k i i i i i osc i i m x l S z l C T q q q q q & & & & (2.59) ( ) å+ = =  ¶ ¶ n 1 k i k i i i osc i i m xl C zl S T q q q & & & (2.60) Finally, the time rate of change of equation 2.60 produces, ( ) å+ =    = ÷ ÷ø ö ç çè æ ¶ ¶ n 1 k i k i i i i i i i osc i i i i m xl C zl S x l S z l C T dt d q q q q q q q & && & & & & & (2.61) Combining equations 2.59 and 2.61, in accordance with Lagrange’s equations, results in the following, ( ) å+ = =  ¶ ¶  ÷ ÷ø ö ç çè æ ¶ ¶ n 1 k i k i i i osc i osc i i m xl C zl S T T dt d q q q q & && & (i = 1,2,K, n + 1) (2.62) Note that equation 2.62 has units of force times distance just as the generalized forces. As a result, the equation 2.62 can be thought of as generalized forces due to host vehicle acceleration, as shown below, 50 ( ) å+ = =  n 1 k i osci k i i i i Q m &x&l Cq &zl Sq (i = 1,2,K, n + 1) (2.63) The superposition of equation 2.63 into the generalized force equation 2.28 yields, ( ) ( ( )) ( 1,2, , 1) 1 = +  +   + = å+ = i n Q M F m x l C F m g z l S n j i i i X j j i i Z j j i i K && & q q q (2.64) Thus in the matrix form of the set of equations of motion, the host vehicle oscillation has no affect on the A or B coefficient matrix. The results of the derivation show that the moving host frame produces applied loads on each segment node proportional to the host vehicle acceleration. Note that the derivation due to the host vehicle motion is independent of the functions used to model it, allowing for any harmonic model of the host aircraft dynamics. 2.5 System Equations of Motion of an Alternate Segment Model At this point, the equations of motion for the twodimensional aerial tow system are complete, and given a proper numerical integration routines and a set of valid initial conditions a simulation of the twodimensional aerial tow system can be accomplished. However, the lumped mass model used for the finite element segments is quite simplistic. As discussed earlier, this model places the mass of a segment at the node at the base of the segment and all applied forces on the segment are located at the node. This simplification is particularly poor when there are few segment s that make up the model of the cable. 51 Imagine the cable is modeled as one lumped mass segment. The segment is in the presence of moving fluid, as a result there is a distributed friction load on the segment. According to the lumped mass assumption the inertial force of the segment and the accumulated drag force are applied at the base of the segment. The lumped mass segment equation of motion is as follows, ( ) ml F C F mg S q X q Z q  + && = (2.65) It would seem a much better model would have the mass and applied force applied at the midpoint of the segment. The resulting angular acceleration is twice that of the lumped mass model. ( ) ml F C F mg S q X q Z q  + && = 2 (2.66) An even better model would have the segment modeled as a thin rod or cylinder. This would have the inertial and applied loads at the center of mass of the rod as in the previous model. However, it also has its rotational inertia that resists rotation (the rotational inertia about the center of gravity of a thin rod is Iy = (ml2)/12). The resulting angular acceleration is an average of the two previous equations. ( ) ml F C F mg S q X q Z q  + = 2 && 3 (2.67) The one segment example illustrates the weakness of the lumped mass model. Presumably as the number of segment increase the difference in the dynamics of the 52 lumped mass and the thin rod models will decrease (this will be shown in chapter five). However, the as the number of segments increase so too does the computational workload. The question remains how do the equations of motion derived throughout this chapter for the lumped mass model change for a system modeled by thin rod segments. The answer can be found by rederiving Lagrange’s equations with a new position equation that represents the center of mass of each segment. Upon completion one will realized that equation 2.30 remains true, recall, [A][q&]+ [B][q&2 ]= [Q] The only changes in the equations of motion for the thin rod modeled segments are in the mass distribution in the A and B coefficient matrices and the lever arm of the generalized force equations. The mass distribution of the coefficient matrices in the lumped mass model is given in equation 2.33. The mass distribution of the coefficient matrices in the thin rod modeled system for n+1 segments is as follows, ( ) å+ = + ÷ ÷ø ö ç çè æ =  + 1 max , , , max( , ) 2 4 1 n k i j k i j i j i j m m m d (2.68) The generalized force equation remains much the same however since the forces on the ith segment are applied at the center of mass and not at the base of the segment the lever arm is halved when i is equal to j, 53 [( ) ( ( )) ] ( 1,2, , 1) 2 1 1 , = + ÷ ÷ø ö ç çè æ   +   + = å+ = i n Q M F m x l C F m g z l S n j i i j i X j j i i Z j j i i K && & d q q (2.69) Note the for equation 2.69 to be accurate for the ATV segment the ATV center of mass must coincide with the center of mass of the segment model. In other words, the length of the ATV segment will be twice the distance between the cableATV hitch point and the ATV center of mass. The set of equations of motion for the aerial towed system with segments modeled as thin rods is attained with the implementation of equations 2.32, 2.49, 2.68 and 2.69 on equation 2.30. Given a proper numerical integration routines and a set of valid initial conditions a simulation of the twodimensional aerial tow system can be accomplished. The difference in the results of the lumped mass and thin rod segment modeled system simulations with the same initial conditions will be analyzed in chapter five. It is only then a conclusion can be drawn about the pros and cons of the two methods. 54 CHAPTER 3 3THREEDIMENSIONAL DYNAMICS 3.1 Introduction Modeling the aerial tow system for the twodimensional case is moreorless straight forward. The equations of motion have been thoroughly derived, include nonlinearities, and make few assumptions. As a result, we can expect that the simulation of the system will be thorough as well; the system could be used to complete a wide variety of case studies, such as, ATV position due to wind gust or host vehicle motion, as well as ATV tracking ability under various conditions. However, it is still limited the phenomena that occur only in the twodimensional plane. For a more robust simulation it is necessary to derive the equations of motion in three dimensions. Being that the actual system exists in three dimensions, this will more accurately represent the real world dynamic occurrences. The addition of the third dimension will drastically increase the complexity of the mathematical models. The methodology for the derivation of the threedimensional equations of motion is the same as that used for the twodimensional case (chapter 2) although a few new concepts are necessary to carryout that methodology. 3.2 Cable Dynamics The threedimensional model of the cable is attained via the application of the finite element technique. As in the twodimensional case, the continuous cable is broken into an arbitrary number of segments. The lumped mass model of the cable segments will be 55 used in the initial analysis of the cable. Each segment is free to rotate about its neighbors, giving the series of rigid connections flexibility as a whole. The added dimension allows the segments two additional rotations, the rotation about the xaxis and the ability for each segment to spin (rotation about the zdirection). The equations of motion for the segment motion will include the new dynamic freedom, and as a result, the equations will increase in complexity. To focus on the dynamics of the cable we will assume that the cable is attached to the origin of the inertial frame. In other words, assume that the host vehicle is steady. Thus, the focus of the dynamics will be on the threedimensional motion of the cable. The effects of the motion of the cable connection point will be addressed in a later section. 3.2.1 Coordinate Systems The addition of the third dimension does not have an effect on the number of coordinate systems needed, although it does affect the type and complexity of the orientation of the coordinate systems. The inertial frame is orthogonal and fixed. Each cable segment has an orthogonal coordinate system which has its origin such that the location of the node of that segment is the length of the segment in the positive zdirection. For instance, imagine there is one cable segment that rotates about the inertial frame. The origin of the cable segment frame is at the point of rotation of the segment. Thus, the location of the origin of the first cable segment frame and the origin of the inertial frame are the same; the difference is in the orientation of the frames, figure 3.1. 56 Figure 3.1 Threedimensional inertial and cable segment coordinate systems The orientation of the cable segment relative to the inertial frame is due to the angle of rotation of the cable frame. The cable frame can rotate about three separate axes (x, y, and zaxes) unlike the twodimensional case in which it can rotate about only one (the yaxis). The common method for describing the threedimensional rotation requires the use of Euler angles. These Euler angles are the heading angle y, the attitude angle q, and the bank angle f. Any orientation of the cable frame relative to the inertial frame can be described by using these three angles. However, it is important to note that the order of rotation is important. Imagine that a cable segment frame is parallel to the inertial frame and their origins are coincident, figure 3.2. Allow the cable frame to rotate about its zaxis by an angle y. The resulting coordinate system is the double primed system. Thus, the coordinate transformation from the doubleprimed system to the fixed inertial frame is as follows, 57 ú ú ú û ù ê ê ê ë é ¢ ¢ ¢ úû ù êë é = ú ú ú û ù ê ê ê ë é ¢ ¢ ¢ ú ú ú û ù ê ê ê ë é =  ú ú ú û ù ê ê ê ë é ¢ ¢ ¢ = ú ú ú û ù ê ê ê ë é z y x C z y x S C C S z y x Z Y X y y y y y 0 0 1 0 0 (3.1) Next, allow the doubleprimed system to rotate about the y ¢ axis by an angle q. The resulting coordinate system is the tripleprimed system. The coordinate transformation from the tripleprimed system to the double primed system is as follows, ú ú ú û ù ê ê ê ë é ¢¢ ¢¢ ¢¢ úû ù êë é = ú ú ú û ù ê ê ê ë é ¢¢ ¢¢ ¢¢ ú ú ú û ù ê ê ê ë é  = ú ú ú û ù ê ê ê ë é ¢ ¢ ¢ z y x C z y x S C C S z y x q q q q q 0 0 1 0 0 (3.2) Finally, allow the tripleprimed system to rotate about the x¢¢ axis by an angle f. The resulting coordinate system is the unprimed cable segment frame. The coordinate transformation from the cable segment frame to the tripleprimed system is as follows, ú ú ú û ù ê ê ê ë é úû ù êë é = ú ú ú û ù ê ê ê ë é ú ú ú û ù ê ê ê ë é  = ú ú ú û ù ê ê ê ë é ¢ ¢ ¢ z y x C z y x S C C S z y x f f f f f 0 0 1 0 0 (3.3) The rotational coordinate transformation from the cable segment frame to the inertial frame is attained by combining equations 3.1, 3.2, and 3.3, ú ú ú û ù ê ê ê ë é ú ú ú û ù ê ê ê ë é  +  +  + + = ú ú ú û ù ê ê ê ë é úû ù êë é úû ù êë é úû ù êë = é ú ú ú û ù ê ê ê ë é z y x S C S C C S C C C S S S C S S S C C C S C C S S S S C S C z y x C C C Z Y X q q f q f y q y f y q f y f y q f y q y f y q f y f y q f y q f (3.4) 58 Figure 3.2 Euler angle rotations The series of three rotations of the cable frame are sufficient to attain any orientation. The same is true independent of the order of rotation. However, the order of rotation does have an effect on the final position of the cable frame. The order used throughout this paper is y, q, and f for every threedimensional segment rotation. Let us assume that the Euler angles are limited to the following rotations, f p p q p y p , 0 2 2 2 0 £ < 2 ,  £ £ £ < (3.5) When q = ± p 2 there is no distinct set of solutions for y and f , this is known as gimbal lock. To avoid complications associated with gimbal lock it is common to switch to an alternate system. In this particular aerial tow system the issue of gimbal lock will not arise as the system necessitates the ATV fly at a lower altitude than the host aircraft resulting in cable segment attitude angles between± p 2. 59 The threedimensional rotational coordinate transformation from the cable to the inertial coordinate system is represented by the following directional cosine matrix, ú ú ú û ù ê ê ê ë é  +  +  + + = úû ù êë é q q f q f y q y f y q f y f y q f y q y f y q f y f y q f f q y S C S C C S C C C S S S C S S S C C C S C C S S S S C S C C , , (3.6) The position of the ith segment’s node in its coordinate system is by definition the length of the segment in the zdirection. Since the segments are a series of connected bodies the directional cosine matrix can be used to formulate the node’s position in the inertial frame, å= ú ú ú û ù ê ê ê ë é úû ù êë = é ú ú ú û ù ê ê ê ë é i j i j i i l C Z Y X j j j 1 , , 0 0 f q y (3.7) Expansion of equation 3.9 yields a more explicit form of the ith segment position in the inertial frame, ( ) ( ) å ( ) å å = = = = =  + = + i j i j i j i j i j i j j j j j j j j j j j j j Z l C C Y l C S S S C X l S S C S C 1 1 1 q f y f y q f y f y q f (i = 1,2,K, n) (3.8) As in the towdimensional case, the location of the node of a segment is dependant upon the angles of rotation of that segment and each preceding segment. 60 3.2.2 Lagrange Equations The modeling of the cable is accomplished using the finite element method, just as with the twodimensional case. As discussed previously, the method breaks the cable into an arbitrary number of segments, and each cable segment has its own coordinate system. The accuracy of the cable model improves with an increase in the number of segments for a given cable, so too does the difficulty of the derivation of the equations of motion for the cable. Using Newton’s method to derive the equations of motion for such a complex configuration would be fraught with danger, as it would require intense vector bookkeeping and multiple coupled rotation equations for each segment. The advantages of the Lagrange technique over Newton’s method are even stronger for threedimensional analysis than in the twodimensional case. The added dimension greatly increases the system complexity and since the Lagrange technique has a standard process for the derivation of the equations of motion it is less likely that mistakes will be made in the derivation of the equations. The bounds of the Lagrange equation are set by the degree of freedom of the system (refer to section 2.2.2). In three dimensions, the complete description of a cable segment is composed of its nodal position and its orientation (i.e. angle of twist). In the Cartesian coordinate system the position of the node is defined by its x, y, and zcoordinates. Due to the fact that the segment length is fixed the coordinates are not independent. The Cartesian constraint equations are shown below, 61 2 1 1 2 1 1 2 1 1 2 ÷ ÷ø ö ç çè æ  + ÷ ÷ø ö ç çè æ  + ÷ ÷ø ö ç çè æ = å å å =  =  = i j i j i j i j i j i i j l X X Y Y Z Z (i = 1,2,K, n) (3.9) The Cartesian coordinates alone cannot specify the angle of twist of the segment. Therefore, there are in actuality four coordinates necessary to define the system, and being that there is one constraint equation each cable segment has three degrees of freedom. It is possible to completely define the position and orientation of a segment using the Euler angles. In this case, there are three angles which completely define the segment and there are no constraint equations. Therefore, the Euler angles of each segment are the independent generalized coordinates of the system. For the twodimensional case the independent generalized coordinates were the attitude angles of rotation and since there was only one angle per segment the generalized coordinate consisted of just one variable. However, for the threedimensional model each segment requires three angles for complete description of orientation. As a result, for a cable modeled by n segments there are 3n generalized coordinates made up of y, q, and f for each segment. (This assumes that each cable segment is free to fully rotate about the preceding segment. The validity of this assumption will be discussed in a later section). For the threedimensional system, Lagrange’s equations represent an equation for each of the 3n degrees of freedom, i i i Q q T q T dt d = ¶ ¶  ÷ ÷ø ö ç çè æ ¶ ¶ & (i = 1, 2,K,3n) (3.9) 62 ( ) úû ù ¶ ¶ + ¶ ¶ + ¶ ¶ êë é + ¶ ¶ + + ¶ ¶ + ¶ ¶ =å= i j j i j j i j j n j i j Z j j i j Y j i j i X j q M q M q M q Z F m g q Y F q X Q F f q y f q y 1 (i = 1, 2,K,3n) (3.10) Note that the generalized force, Q, contains applied forces, applied moments, and inertial forces as per the discussion in section 2.2.2. Equation 3.9 represents a set of 3n equations. In which the generalized coordinate is made up of three independent variables. It may be presented more clearly rewritten as three sets of n equations, in which the generalized coordinate for each is an independent variable, i i i Q T T dt d Q T T dt d Q T T dt d i i i i i i y q f y y q q f f = ¶ ¶  ÷ ÷ø ö ç çè æ ¶ ¶ = ¶ ¶  ÷ ÷ø ö ç çè æ ¶ ¶ = ¶ ¶  ÷ ÷ø ö ç çè æ ¶ ¶ & & & (i = 1, 2,K, n) (3.11) ( ) ( ) å ( ) å å = = = ¶ ¶ + ¶ ¶ + + ¶ ¶ + ¶ ¶ = ¶ ¶ + ¶ ¶ + + ¶ ¶ + ¶ ¶ = ¶ ¶ + ¶ ¶ + + ¶ ¶ + ¶ ¶ = n j i j j i j Z j j i j Y j i j i X j n j i j j i j Z j j i j Y j i j i X j n j i j j i j Z j j i j Y j i j i X j M Z F m g Y F X Q F M Z F m g Y F X Q F M Z F m g Y F X Q F 1 1 1 y y y y y q q q q q f f f f f y y q q f f (i = 1, 2,K, n) (3.12) 63 Thus, each segment will contain three equations of motion. Note that the applied moments and the moments produced by the applied forces are about the axis of rotation of their respective generalized coordinate, i f& , i q& , and i y& . The method for deriving the equations of motion is the same as that carried out in the twodimensional cable dynamic section. The only difference is the more complex kinetic energy equation and the increase in the number of generalized coordinates of the system. The kinetic energy is made up of rotational and translational kinetic energy. The nodes, being point masses, have no inertia and thus no rotational kinetic energy. As a result, the derivation of the equations of motion relies solely on the translation kinetic energy of the system as a function of the generalized coordinates. The kinetic energy can be found by calculating the square of velocity for the node of each segment as a function of the generalized coordinates, equation 2.9. The velocity components are simply the time rate of change of position [ ( ) ( ) ( )] [ ( ) ( ) ( )] [ ( ) ( )] for (i n) Z l C S S C Y l C C S S S S C C S S C S C X l S C C S S C C C C S S S C i j i j j j i j i j j j j i j i j j j j j j j j j j j j j j j j j j j j j j j j j j j j j j j j j j 1, 2 , 1 1 1 K & & & & & & & & & & & = =  +  =   + + + =  + +  å å å = = = q f q f y f y q f y q f y f y q f y f y q f y q f y f y q f f q f q y f q y (3.13) Thus, the square of velocity for the kth segment is as follows, 64 ( ) { [ ( ( ) ( ) ) ( ) ( ( ) ( )] [ ( ) ( ) ( )] [ ( ( ) ( ) ( ) ( ( ) ( )]} j i i j i j i j j i i j i i j j j i i j j j i i j i j j i i j i i j j j i i j i j i j j i i j i j i i j i j j i i j i j i j i j j i i j i j j i i j i i j i j j i i j i i j j j i i j i j i j j i i j i j j i i j j i j i j j i i j i j i j j i i j j i j i i j i j C C C S S S S S C S S S C S C S C C S C C S C C C C S S S S S C S S S S C C C S C S C C C S C C C C C C C S S C C S C S S C C C C S S C C S S S C C S S S C S S C S C S C S C S C C C S C C S C S S S S C C S C S S S C S S S C C V X Y Z j j i j j j i j j j i j k i k j k k k k y y f f q q f f y y f f q f f q y y f f q y y f f q q y y f f q f f q y y f f q q f f y y f f q y y f f q q y y f f q q f f q q y y f f q q y y f f q f f q q y y f f q f f q y y f f q q f f y y f f q q y y f f q f f q q y y f f q q f f y y f f q f f q f f q q y q y f y q q f y q f f + +  + +   +  + +  + + +  + + +  + +  +   + + + +  + + = + + =                  = = åå & & & & & & & & & & & & & & & 1 1 2 2 2 2 (3.14) Furthermore, the threedimensional kinetic energy equation is represented as a function of the independent generalized coordinates below, { [ ( ( ) ( ) ) ( ) ( ( ) ( )] [ ( ) ( ) ( )] [ ( ( ) ( ) ( ) ( ( ) ( )]} j i i j i j i j j i i j i i j j j i i j j j i i j i j j i i j i i j j j i i j i j i j j i i j i j i i j i j j i i j i j i j i j j i i j i j j i i j i i j i j j i i j i i j j j i i j i j i j j i i j i j j i i j j i j i j j i i j i j i j j i i j j i j i i j i j C C C S S S S S C S S S C S C S C C S C C S C C C C S S S S S C S S S S C C C S C S C C C S C C C C C C C S S C C S C S S C C C C S S C C S S S C C S S S C S S C S C S C S C S C C C S C C S C S S S S C C S C S S S C S S S C C T m V m j j i j j j i j j j i j n k k i k j k n k k k y y f f q q f f y y f f q f f q y y f f q y y f f q q y y f f q f f q y y f f q q f f y y f f q y y f f q q y y f f q q f f q q y y f f q q y y f f q f f q q y y f f q f f q y y f f q q f f y y f f q q y y f f q f f q q y y f f q q f f y y f f q f f q f f q q y q y f y q q f y q f f + +  + +   +  + +  + + +  + + +  + +  +   + + + +  + + = =                  = = = = å ååå & & & & & & & & & & & & 1 1 1 1 2 2 1 2 1 (3.15) Equation 3.15 is the basis from which Lagrange’s equations will be derived. All one needs to do is compare equation 3.15 and 2.12 to grasp the increase in complexity from a 65 system existing in two to threedimensions. With the cable kinetic energy equation derived the derivation of Lagrange’s equations for each set of generalized coordinates can commence. 3.2.2.1 Bank angle The implementation of the Lagrange equation in the threedimensional cable model must be accomplished for each angle necessary for cable orientation definition. For the bank angle generalized coordinate the derivation is as follows. The derivative of kinetic energy with respect to the bank angle and the time rate of change of the bank angle of the ith cable segment yield the following, { [ ( ( ) ( ) ) ( ) ( ( ) ( )] [ ( ) ( ) ( )] [ ( ( ) ( ) ( ) ( ( ) ( )]} j i i j i j i j j i i j i i j j j i i j j j i i j i j j i i j i i j j j i i j i j i j j i i j i j i i j i j j i i j i j i j i j j i i j i j j i i j i i j i j j i i j i i j j j i i j i j i j j i i j i j j i i j j i j i j j i i j i j i j j i i j j i j i i j i j C S C S S C S S S S S C C S C C C C S S C S C C S C S C S S S S S S S C C C S S C S S C C S C S C C C S C S S C S S C S S S C C S S S C C C S S S C S S C C S S S S C C C S C S S C C C C C S C C S S S S C S S S S C C S C S C C m T j j i j j j i j j j i j n k k j k i y y f f q q f f y y f f q f f q y y f f q y y f f q q y y f f q f f q y y f f q q f f y y f f q y y f f q q y y f f q q f f q q y y f f q q y y f f q f f q q y y f f q f f q y y f f q q f f y y f f q q y y f f q f f q q y y f f q q f f y y f f q f f q f f q q y q y f y q q f y q f f f  + +    +  + + +  + +   +   +   + + +  + + +  +   + + = ¶ ¶                  = = åå & & & & & & & & & & & & 1 1 (3.16) 66 { ( ( ) ( ) ) ( ) ( ( ) ( )} j i i j i i j j j i i j i j i j j i i j i j j i i j j i j i j j i i j i j i j j i i j j i j i i j i j C S S S C C S S S C S S C S C S C S C S C C C S C C S C S S S S C C S C S S S C S S S C C m T j j j n k k j k i y y f f q f f q y y f f q q f f y y f f q q y y f f q f f q q y y f f q q f f y y f f q f f q f f q q y q f f  + +    + + + +  + + = ¶ ¶       = = åå & & & & 1 1 (3.17) Finally evaluating the time rate of change of equation 3.17 produces, { ( ( ) ( ) ) ( ) ( ( ) ( ) ( ) ( ) ( ( ) ( ) ( ( ) ( ) ) ( ) ( ( ) ( ) [ ( ( ) ( ) ) ( ) ( ( ) ( )] [ ( ) ( ) ( )] [ ( ( ) ( ) ( ) ( ( ) ( )]} j i i j i j i j j i i j i i j j j i i j j j i i j i j j i i j i i j j j i i j i j i j j i i j i j i i j i j j i i j i j i j i j j i i j i j j i i j i i j i j j i i j i i j j j i i j i j i j j i i j i j j i i j j i j i j j i i j i j i j j i i j j i j i i j i j j i i j i i j j j i i j i j i j j j i i j i j j i i j j i j i j j j i i j i j i j j i i j j i j i i j i j j i i j i i j j j i i j i j i j j i i j i j j i i j j j i i j i j j i i j j i j i j j i i j i i j j j i i j i j i j j i i j i j j i i j j i j i j j i i j i j i j j i i j j i j i i j i j C S C S S C S S S S S C C S C C C C S S C S C C S C S C S S S S S S S C C C S S C S S C C S C S C C C S C S S C S S C S S S C C S S S C C C S S S C S S C C S S S S C C C S C S S C C C C C S C C S S S S C S S S S C C S C S C C S S S S C C S C S C S S C S C S C S S S C C S S C C C C S C S S C S S C C S S S S S C C C C S C S C S S S S S S S C C S S C S C C C C C C S S S C S C S C S S C S C S S S C C S S S C S S C S C S C S C S C C C S C C S C S S S S C C S C S S S C S S S C C m T dt d j j i j j j i j j j i j j j j j j j j j j j n k k j k i y y f f q q f f y y f f q f f q y y f f q y y f f q q y y f f q f f q y y f f q q f f y y f f q y y f f q q y y f f q q f f q q y y f f q q y y f f q f f q q y y f f q f f q y y f f q q f f y y f f q q y y f f q f f q q y y f f q q f f y y f f q f f q f f q q y y f f q f f q y y f f q q f f y y f f q q y y f f q f f q q y y f f q q f f y y f f q f f q f f q q y y f f q f f q y y f f q q f f y y f f q q y y f f q y y f f q q y y f f q f f q q y y f f q f f q y y f f q q f f y y f f q q y y f f q f f q q y y f f q q f f y y f f q f f q f f q q y q y f y q q f y q f f y q f y f q y f q y q f f  + +    +  + + +  + +   +   +   + + +  + + +  +   + + + +  + + + +  + + + +  + +   +  + +  +  + +  +   + + + +  + + = ÷ ÷ø ö ç çè æ ¶ ¶                                    = = åå & & & & & & & & & & & & & & & & & & & & & & && & & 2 2 1 1 2 2 2 (3.18) 67 Note that, just as in the twodimensional case, the term representing the derivative of kinetic energy with respect to generalized coordinate, equation 3.16, is completely canceled by the some of the terms in the time rate of change of the derivative of kinetic energy with respect to the time rate of change of the generalized coordinate, equation 3.18. Combining equations 3.16 and 3.18 results in the series representation of the left hand side of the Lagrange’s equations for the bank angle of the threedimensional finite element modeled cable, { ( ( ) ( ) ) ( ) ( ( ) ( ) ( ) ( ) ( ( ) ( ) ( ( ) ( ) ) ( ) ( ( ) ( )} for (i n) S S S S C C S C S C S S C S Q C S C S S S C C S S C C C C S C S S C S S C C S S S S S C C C C S C S C S S S S S S S C C S S C S C C C C C C S S S C S C S C S S C S C S S S C C S S S C S S C S C S C S C S C C C S C C S C S S S S C C S C S S S C S S S C C m T T dt d j i i j i i j j j i i j i j i j i j j i i j i j j i i j j i j i j j j i i j i j i j j i i j j i j i i j i j j i i j i i j j j i i j i j i j j i i j i j j i i j j j i i j i j j i i j j i j i j j i i j i i j j j i i j i j i j j i i j i j j i i j j i j i j j i i j i j i j j i i j j i j i i j i j j j j j j j j j j j n k k j k i i 1,2, , 2 2 2 2 2 1 1 K & & & & & & & & & & && & & = + +  = + + +  + + + +  + +   +  + +  +  + +  +   + + + +  + + = ¶ ¶  ÷ ÷ø ö ç çè æ ¶ ¶                   = = åå y y f f q f f q y y f f q q f f f y y f f q q y y f f q f f q q y y f f q q f f y y f f q f f q f f q q y y f f q f f q y y f f q q f f y y f f q q y y f f q y y f f q q y y f f q f f q q y y f f q f f q y y f f q q f f y y f f q q y y f f q f f q q y y f f q q f f y y f f q f f q f f q q y q f y f q y f q y q f f f (3.19) To fully develop the equations of motion the generalized force equation as a function of bank angle is evaluated in accordance with equation 3.12. 68 [ ( ) ( ) (F m g)l ( S C )] for (i n) Q M F l C S S S C F l C C S S S i i i i i i i i i i i i Z j j i n j i i i X j i Y j i +  =1,2,K, +   +  + = å= f q f f f y f q y f y f q y (3.20) The implementation of the Lagrange equations for the bank angle generalized coordinate is complete, however to attain the full system equations of motion Lagrange’s equations must be derived for the two remaining generalized coordinates. 3.2.2.2 Attitude Angle This section investigates the application of the Lagrange equations for the attitude angle generalized coordinate. The Lagrange equation term representing the derivative of kinetic energy with respect to the attitude angle and the rate of change of the attitude angle of the ith cable segment is as follows, { [ ( ) ( ) ( )] [ ( ) ( ) ( )] [ ( ) ( ) ( )]} j i i j i j j i i j i j i i j i j j i i j i j i i j i j j i i j i j i i j i j j i i j i j i j i j j i i j i j j i i j i i j i j j i i j i j i i j i j j i i j i j i j i j j i i j i j j i i j i i j i j C C C C S S C S C S C C C C C C C C S C S C S C C S S S C C S S C C C S C C C C S C C S S S S C C S C S C C C S S C S S C C S C S C C C S C S S C S S C S S S C C S S S C m T j j i j j j i j j j i j n k k j k i y y f f q q y y f f q y y f f q q y y f f q y y f f q q y y f f q y y f f q q y y f f q q f f q q y y f f q q y y f f q f f q q y y f f q y y f f q q y y f f q q f f q q y y f f q q y y f f q f f q q y q y f y q q f y q f f q                = = + +   +  + +  + +  + +  + +   +   + = ¶ ¶ åå & & & & & & & & & & & & 1 1 (3.21) 69 { ( ) ( ) ( j i i j i j i i j i j )} j i i j i j i j i j j i i j i j j i i j i i j i j C C S C S C C C S C C C C C C C S S C C S C S S C C C C S S C m T j j j n k k j k i y y f f q y y f f q q y y f f q q f f q q y y f f q q y y f f q f f q q y q f q      = =  + +  + + + = ¶ ¶ åå & & & & 1 1 (3.22) Taking the time rate of change of the equation 3.22 produces the first term in Lagrange’s equations for the attitude angle generalized coordinate, 70 { ( ) ( ) ( ) ( ) ( ) ( ) ( ) ( ) ( ) [ ( ) ( ) ( )] [ ( ) ( ) ( )] [ ( ) ( ) ( j i i j i j j i i j i )]} j i i j i j j i i j i j i i j i j j i i j i j i i j i j j i i j i j i j i j j i i j i j j i i j i i j i j j i i j i j i i j i j j i i j i j i j i j j i i j i j j i i j i i j i j j j i i j i j i i j i j j j i i j i j i j i j j j i i j i j j i i j i i j i j j i i j i j i i j i j j i i j i j j i i j i j i j i j j i i j i j i i j i j j i i j i j i j i j j i i j i j j i i j i i j i j C C C C S S C S C S C C C C C C C C S C S C S C C S S S C C S S C C C S C C C C S C C S S S S C C S C S C C C S S C S S C C S C S C C C S C S S C S S C S S S C C S S S C S C S C C C C C S C C C C S C C S C C C C C S S C S C C C S C C C C C S C S C S S C C C C C C S C C C S S S C C S C S C C C S C C C C C C C S S C C S C S S C C C C S S C m T dt d j j i j j j i j j j i j j j j j j j j j j n k k j k i y y f f q q y y f f q y y f f q q y y f f q y y f f q q y y f f q y y f f q q y y f f q q f f q q y y f f q q y y f f q f f q q y y f f q y y f f q q y y f f q q f f q q y y f f q q y y f f q f f q q y y f f q y y f f q q y y f f q q f f q q y y f f q q y y f f q f f q q y y f f q y y f f q q y y f f q q y y f f q q f f q q y y f f q y y f f q q y y f f q q f f q q y y f f q q y y f f q f f q q y q y f y q q f y q f f y q f y f q y f q y q f q                              = = + +   +  + +  + +  + +  + +   +   +   +  + +   + + + +  +   +  + + +  + + + = ÷ ÷ø ö ç çè æ ¶ ¶ åå & & & & & & & & & & & & & & & & & & & & & & & & & 2 2 2 1 1 2 2 2 (3.23) Again, the time rate of change of the derivative of kinetic energy with respect to the time rate of change of the generalized coordinate, equation 3.23, contains terms that are equal to the derivative of kinetic energy with respect to the generalized coordinate, equation 3.21. Combining equation 3.21 and 3.23 appropriately produces the left hand side of Lagrange’s equations for the attitude angle generalized coordinate, 71 { ( ) ( ) ( ) ( ) ( ) ( ) ( ) ( ) j ( j i i j i j i i j i j )} i j j i i j i j i j i j j j i i j i j j i i j i i j i j j i i j i j i i j i j j i i j i j j i i j i j i j i j j i i j i j i i j i j j i i j i j i j i j j i i j i j j i i j i i j i j S C S C C C C C S Q C C C C S C C S C C C C C S S C S C C C S C C C C C S C S C S S C C C C C C S C C C S S S C C S C S C C C S C C C C C C C S S C C S C S S C C C C S S C m T T dt d j j j j j j j j j n k k j k i i y y f f q y y f f q q q y y f f q q f f q q y y f f q q y y f f q f f q q y y f f q y y f f q q y y f f q q y y f f q q f f q q y y f f q y y f f q q y y f f q q f f q q y y f f q q y y f f q f f q q y q f y f q y f q y q f q q   =  + +   + + + +  +   +  + + +  + + + = ¶ ¶  ÷ ÷ø ö ç çè æ ¶ ¶               = = åå 2 2 2 1 1 2 2 2 & & & & & & & & & & && & & (3.24) The right hand side of Lagrange’s equations is the derivation of the generalized force for the ith segment. The evaluation of equation 3.12 for the attitude angle generalized coordinate results in, [ ( ) ( ) ( ) ( )] for (i n) Q M F l C C C F l C C S F m g l C S n j i i i X j i i i i Y j i i i i Z j j i i i = 1, 2,K,  + + + + = å= q q f q y f q y f q (3.25) Thus, concludes the derivation of the equations of motion for the attitude angle. However, there still remains one generalized coordinate not yet applied to the Lagrange equations. 3.2.2.3 Heading Angle The final generalized coordinate is the heading angle. The equations of motion for the heading angle are derived in the same fashion as that of the other Euler angles. The 72 derivation of the Lagrange equations for the heading angle of the ith segment as shown in equation 3.11 is presented in detail in the remainder of this section. The derivative of kinetic energy with respect to the heading angle and time rate of change of heading angle of the cable segment is evaluated below, { [ ( ( ) ( ) ( ) ( ( ) ( )] [ ( ) ( ) ( )] [ ( ( ) ( ) ( ) ( ( ) ( )]} j i i j i j i j j i i j i i j j j i i j j j i i j i j j i i j i i j j j i i j i j i j j i i j i j i i j i j j i i j i j j i i j i j j i i j i j i i j i i j j j i i j i j i j j i i j i j j i i j j j i i j i j i j j i i j j i j i S C C S S S S C C S S S C S S S C C C C C S C S C C S S S S C C S S S S C S C S C C C C C S S C C C C S C S C S C C C C S S S S C C S C S C S S C S S S C S C C C C C S S S S S C C C C S S S C S m T j j i j j j i j j j i j n k k j k i y y f f q q f f y y f f q f f q y y f f q y y f f q q y y f f q f f q y y f f q q f f y y f f q y y f f q q y y f f q q y y f f q q y y f f q y y f f q f f q y y f f q q f f y y f f q q y y f f q y y f f q q f f y y f f q f f q y q y f y q q f y q f f y +    +     + + + + +   +  +   +  + + +   + = ¶ ¶                  = = åå & & & & & & & & & & & & 1 1 (3.26) { ( ( ) ( ) ( ) ( ( ) ( )} j i i j i j i j j i i j i i j j j i i j j j i i j i j j i i j i i j j j i i j i j i j C C C S S S S S C S S S C S C S C C S C C S C C C C S S S S S C S S S S C m T j j j n k k j k i y y f f q q f f y y f f q f f q y y f f q y y f f q q y y f f q f f q y y f f q q f f y q f y + +  + +   +  + + = ¶ ¶       = = åå & & & & 1 1 (3.27) The first term of Lagrange’s equations is attained by taking the time rate of change of equation 3.27. The results are as follows, 73 { ( ( ) ( ) ( ) ( ( ) ( ) ( ) ( ) ( ( ) ( ) ( ( ) ( ) ( ) ( ( ) ( ) [ ( ( ) ( ) ( ) ( ( ) ( )] [ ( ) ( ) ( )] [ ( ( ) ( ) ( ) ( ( ) ( )]} j i i j i j i j j i i j i i j j j i i j j j i i j i j j i i j i i j j j i i j i j i j j i i j i j i i j i j j i i j i j j i i j i j j i i j i j i i j i i j j j i i j i j i j j i i j i j j i i j j j i i j i j i j j i i j j i j i j j i i j i j i j j i i j i i j j j j i i j j j i i j i j j j i i j i i j j j i i j i j i j j i i j i j i j j i i j i i j j j i i j j j i i j i j j i i j j j i i j i j j i i j i j i j j i i j i i j j j i i j j j i i j i j j i i j i i j j j i i j i j i j S C C S S S S C C S S S C S S S C C C C C S C S C C S S S S C C S S S S C S C S C C C C C S S C C C C S C S C S C C C C S S S S C C S C S C S S C S S S C S C C C C C S S S S S C C C C S S S C S S C C S S S S C C S S S C S C S C S S C C S S C C S S S C S S C C S S S S C C S S S S C S C C S S S S S S C C C C C S C C S S C S C S S C C C C S S S S S C S S S C S C S C C S C C S C C C C S S S S S C S S S S C m T dt d j j i j j j i j j j i j j j j j j j j j j n k k j k i y y f f q q f f y y f f q f f q y y f f q y y f f q q y y f f q f f q y y f f q q f f y y f f q y y f f q q y y f f q q y y f f q q y y f f q y y f f q f f q y y f f q q f f y y f f q q y y f f q y y f f q q f f y y f f q f f q y y f f q q f f y y f f q f f q y y f f q y y f f q q y y f f q f f q y y f f q q f f y y f f q q f f y y f f q f f q y y f f q y y f f q q y y f f q y y f f q q y y f f q q f f y y f f q f f q y y f f q y y f f q q y y f f q f f q y y f f q q f f y q y f y q q f y q f f y q f y f q y f q y q f y +    +     + + + + +   +  +   +  + + +   +  + +  +   +  +   +  + + + +  + +   + + +  + + +   +  + + = ÷ ÷ø ö ç çè æ ¶ ¶                                    = = åå & & & & & & & & & & & & & & & & & & & & & & && & & 2 2 2 1 1 2 2 2 (3.28) As with the previous derivations of Lagrange’s equations for the cable system the time rate of change of the derivative of kinetic energy with respect to the time rate of change of the generalized coordinate contains terms that are equal to the derivative of kinetic energy with respect to the generalized coordinate, somewhat simplifying the equations of motion. 74 Combining equation 3.26 and 3.28 appropriately produce the left hand side of the Lagrange’s equations of the heading angle generalized coordinate. { ( ( ) ( ) ( ) ( ( ) ( ) ( ) ( ) ( ( ) ( ) ( ( ) ( ) ( ) ( ( ) ( ) } j j i i j i j i j j i i j i i j j i j j i i j j j i i j i j j j i i j i i j j j i i j i j i j j i i j i j i j j i i j i i j j j i i j j j i i j i j j i i j j j i i j i j j i i j i j i j j i i j i i j j j i i j j j i i j i j j i i j i i j j j i i j i j i j S C C S S S S C C S S S C S Q C S C S S C C S S C C S S S C S S C C S S S S C C S S S S C S C C S S S S S S C C C C C S C C S S C S C S S C C C C S S S S S C S S S C S C S C C S C C S C C C C S S S S S C S S S S C m T T dt d j j j j j j j j j n k k j k i i y y f f q q f f y y f f q f f q y y y f f q y y f f q q y y f f q f f q y y f f q q f f y y f f q q f f y y f f q f f q y y f f q y y f f q q y y f f q y y f f q q y y f f q q f f y y f f q f f q y y f f q y y f f q q y y f f q f f q y y f f q q f f y q f y f q y f q y q f y y  + +  =   +  +   +  + + + +  + +   + + +  + + +   +  + + = ¶ ¶  ÷ ÷ø ö ç çè æ ¶ ¶                   = = åå 2 2 2 1 1 2 2 2 & & & & & & & & & & && & & (3.29) The final addition to the dynamics equations are the effect of the applied loads on the heading angle. The generalized force for the heading angle is derived in accordance with equation 3.12 and result in the following, [ ( ) ( )] for (i n) Q M F l S C C S S F l S S C S C n j i i i X j i i i i i i Y j i i i i i i = 1, 2,K, + +  + = å= y y f y f q y f y f q y (3.30) The equations of motion of the cable segments are dependant upon the equations of motion of the Euler angles of each segment. Simulation of the cable motion requires simultaneous solution of the angular acceleration of the generalized coordinates for all cable segments. 75 3.2.2.4 Cable Equations of Motions We have yet to discuss the threedimensional aerodynamic applied loads, but neglecting that, the equations of motion for the threedimensional cable have been derived, although they are not in a convenient form. We can rewrite the equations of motion in matrix form for an arbitrary number of cable segments as was done in the twodimensional case. The fully threedimensional nonlinear differential equations of motion for a cable of n segments can be written in the following form, ú ú ú û ù ê ê ê ë é = ú ú ú û ù ê ê ê ë é ú ú ú û ù ê ê ê ë é + ú ú ú û ù ê ê ê ë é ú ú ú û ù ê ê ê ë é + ú ú ú û ù ê ê ê ë é ú ú ú û ù ê ê ê ë é y q f yf qy fq y q f y q f Q Q Q C C C C C C C C C B B B B B B B B B A A A A A A A A A & & & & & & & & & & && && 7 8 9 4 5 6 1 2 3 2 2 2 7 8 9 4 5 6 1 2 3 7 8 9 4 5 6 1 2 3 (3.31) where A, B, C, and Q are augmented matrices. In other words, matrix A is composed of nine submatrices of size n x n. Thus A is of size 3n x 3n. The same is true for matrix B and C. A, B, and C are coefficient matrices. The generalized force column matrix, Q, is an augmented matrix of size 3n x 1. The values of the matrices components are taken directly from the results of the derivation of the Lagrange equation carried out in the previous sections, and combine to form a set of 3n equations of motion. The A, B, and C matrices in equation 3.31 are represented by the equations below, where i and j are the row and column of the matrices. [ ( ) ( ) ] j i j j i j i i j i j i j j i j i j i j S C S S S C S S S C C A m l l C S S S S C C i j i j y y f f q f f q f f q q y y f f q q f f  + = + +   i 1 , , i (3.32) ( ) i j j i j i j j i j j i j i j A m l l C S C S C S C C C S C C S i j i j y y f f q q y y f f q f f q q =   + 2 , ,  i  i (3.33) 76 [ ( ) j ( i j i j i j )] i j j i j i i j j S S C S S C S A m l l C S S S C C S i j i j y y f f q q f f y y f f q f f q  =  + +   i 3 , , i (3.34) i j j i A A 4 , 2 , = (3.35) ( ) i j j i j i j i j i j A m l l C C C C C C C S S i j i j y y f f q q f f q q = + 5 , ,  i (3.36) ( ) i j j i j i j i j i j A m l l C C S C S C C C S 6 , i, j i j y yi f f q y yi f f q q =  (3.37) i j j i A A 7 , 3 , = (3.38) i j j i A A 8 , 6 , = (3.39) [ ( ) ( )] j i j i i j j i j j i j i j i j S C S S S C S A m l l C C C S S S S i j i j y y f f q f f q y y f f q q f f  = + +   i 9 , , i (3.40) [ ( ) ( ) ] j i j j i j i i j i j i j j i j i j i j S C C S S S S S C C C B m l l C S C S S C S i j i j y y f f q f f q f f q q y y f f q q f f + + =  +   i 1 , , i (3.41) ( ) i j j i j i j j i j j i j i j B m l l C S C S S S C C S S C C C i j i j y y f f q q y y f f q f f q q = + + 2 , ,  i  i
Click tabs to swap between content that is broken into logical sections.
Rating  
Title  Dynamic Simulation of Low Altitude Aerial Tow Systems 
Date  20050701 
Author  Quisenberry, J. J. 
Department  Mechanical Engineering 
Document Type  
Full Text Type  Open Access 
Abstract  The goal of this paper is to produce a mathematical simulation of the aerial tow system with emphasis on low altitude tracking. The aerial tow system consists of a host aircraft which tows a smaller aircraft via cable. The system components are modeled and a method for deriving the system equations of motion in two and threedimensional space is presented. The implementation of the system equations of motion in conjunction with a numerical solver results in the development of system simulation. The simulation is used as a tool to inspect various flight cases to better understand the nuances of the complex system. 
Note  Thesis 
Rights  © Oklahoma Agricultural and Mechanical Board of Regents 
Transcript  DYNAMIC SIMULATION OF LOW ALTITUDE AERIAL TOW SYSTEMS By JOHNNY EARL QUISENBERRY JR. Bachelor of Science Oklahoma State University Stillwater, Oklahoma 2002 Submitted to the Faculty of the Graduate College of Oklahoma State University in partial fulfillment of the requirements for the Degree of MASTER OF SCIENCE July, 2005 ii DYNAMIC SIMULATION OF LOW ALTITUDE AERIAL TOW SYSTEMS Thesis Approved: ________________________________________________ Thesis Advisor ________________________________________________ ________________________________________________ ________________________________________________ Dean of the Graduate College iii ACKNOWLEDGEMENTS This research was conducted in cooperation with the engineers at Zivko Inc. Specifically; I would like to thank Todd C. Morse for his assistance. In addition I would like to thank the NASA Oklahoma Space Grant Consortium for their financial support for the duration of the project. I would like to express my deepest appreciation to my major advisor, Dr. Andrew S. Arena, for his guidance and support in this research. I would also like to thank the other members of my committee, Dr. R. D. Delahoussaye and Dr. P. R. Pagilla for their time and efforts. I could not have made the accomplishments I have here at OSU without the love and support of my family; mother, Susan Quisenberry, father, John Quisenberry, sister, Kiesha Quisenberry, and brother Clint Quisenberry. iv TABLE OF CONTENTS CHAPTER PAGE 1 INTRODUCTION ........................................................................................................... 1 1.1 Background ......................................................................................................... 1 1.2 Design Goals ....................................................................................................... 4 1.3 Literature Review................................................................................................ 5 1.3.1 Cable Dynamics .............................................................................................. 6 1.3.2 Lagrange’s Equations ...................................................................................... 9 1.3.3 Aircraft Dynamics......................................................................................... 13 1.3.4 Cable Towed Systems ................................................................................... 14 2 TWODIMENSIONAL DYNAMICS........................................................................... 18 2.1 Introduction....................................................................................................... 18 2.2 Cable Dynamics ................................................................................................ 19 2.2.1 Coordinate Systems ...................................................................................... 22 2.2.2 Lagrange Equations ....................................................................................... 25 2.2.3 Aerodynamic Forces ..................................................................................... 35 2.3 ATV Dynamics ................................................................................................. 37 2.3.1 Coordinate system......................................................................................... 38 2.3.2 Lagrange equations ....................................................................................... 40 2.3.3 Generalized Forces........................................................................................ 43 2.4 Host Vehicle...................................................................................................... 45 2.5 System Equations of Motion of an Alternate Segment Model ......................... 50 3 THREEDIMENSIONAL DYNAMICS ....................................................................... 54 3.1 Introduction....................................................................................................... 54 3.2 Cable Dynamics ................................................................................................ 54 3.2.1 Coordinate Systems ...................................................................................... 55 3.2.2 Lagrange Equations ....................................................................................... 60 3.2.3 Aerodynamic Modeling ................................................................................ 79 3.3 Arial Tow Vehicle............................................................................................. 82 3.3.1 Coordinate System........................................................................................ 83 3.3.2 Lagrange Equations ....................................................................................... 85 3.3.3 Generalized Forces........................................................................................ 97 3.4 Host Vehicle...................................................................................................... 99 3.4.1 Coordinate System...................................................................................... 100 3.4.2 Lagrange Equations ..................................................................................... 101 v 3.4.3 Generalized Forces...................................................................................... 106 3.5 System Equations of Motion of an Alternate Segment Model ....................... 108 4 FLIGHT PARAMETERS ............................................................................................ 111 4.1 Introduction..................................................................................................... 111 4.2 Ocean Surface Waveform............................................................................... 112 4.3 ATV Autopilot ................................................................................................ 114 5 SYSTEM SIMULATIONS.......................................................................................... 116 5.1 The Ideal Cable ............................................................................................... 116 5.2 Lumped Mass Versus Thin Rod Cable Model................................................ 119 5.3 Simulation of the CableATV System............................................................ 124 5.3.1 Dimensional Comparison............................................................................ 125 5.3.2 Altitude Effects ........................................................................................... 126 5.4 ATV Wave Tracking....................................................................................... 131 5.4.1 Surface Waveform Variation...................................................................... 132 5.4.2 Host Aircraft Cruise Altitude ...................................................................... 133 5.4.3 Lateral Motion............................................................................................. 138 5.4.4 Host Vehicle Oscillation............................................................................. 141 6 CONCLUSIONS.......................................................................................................... 144 7 REFERENCES ............................................................................................................ 147 vi LIST OF TABLES Table 5.1 Cable parameters and reference flight conditions................................... 120 Table 5.2 Aerial towed aircraft parameters............................................................. 125 vii LIST OF FIGURES Figure 1.1 The aerial tow system ................................................................................. 1 Figure 1.2 Simple pendulum shown in twodimensional a) Cartesian coordinates, and b) cylindrical coordinates.......................................................................... 11 Figure 2.1 Finite element cable segment models ....................................................... 20 Figure 2.2 Lumped mass finite element cable model................................................. 22 Figure 2.3 Twodimensional inertial and cable segment coordinate systems ............ 24 Figure 2.4 Twodimensional ATV coordinate system ............................................... 38 Figure 3.1 Threedimensional inertial and cable segment coordinate systems .......... 56 Figure 3.2 Euler angle rotations ................................................................................. 58 Figure 3.3 ATV coordinate systems ........................................................................... 83 Figure 4.1 Threedimensional ocean surface waveform .......................................... 114 Figure 5.1 Two and Threedimensional ideal cable simulation results. The figure shows, a) the height of cable base, and b) the system energy versus time. ................................................................................................................. 118 Figure 5.2 A cable in fluid flow simulation results. The cable is broken into various numbers of segments, the segment is modeled as, a) lumped masses, and b) thin rods. ............................................................................................. 121 Figure 5.3 Cable base static positions versus number of cable segments ................ 122 Figure 5.4 Potential energy versus number of cable segments. ............................... 123 Figure 5.5 Comparison between two and threedimensional cable and ATV system ................................................................................................................. 126 Figure 5.6 Host vehicle cruise altitude affect on static cable curvature and position ................................................................................................................. 127 Figure 5.7 Dynamic response to host vehicle oscillation as viewed from the host vehicle coordinate system. ...................................................................... 128 Figure 5.8 Elliptic oscillations in the x and zdirection, fo r a host vehicle altitude of a) 250 ft., b) 500 ft., c) 1000 ft., and d) 2000 ft. ..................................... 129 Figure 5.9 Tracking at with various ocean wavelengths .......................................... 133 Figure 5.10 Host vehicle cruise altitude affect on wave tracing a) 250 ft., b) 500 ft., c) 1000 ft., and d) 2000 ft........................................................................ 134 Figure 5.11 Elliptic oscillations in the x and zdirection, for a host vehicle altitude of a) 250 ft., b) 500 ft., c) 1000 ft., and d) 2000 ft. ..................................... 135 Figure 5.12 ATV motion in the x and ydirection, at a host vehicle cruise altitude of a) 250 ft., b) 500 ft, c) 1000 ft., and d) 2000 ft. ..................... 139 Figure 5.13 Threedimensional waveform tracking with a 10 knot side gust, host vehicle altitude of a) 250 ft., b) 500 ft., c) 1000 ft., and d) 2000 ft. ....... 140 Figure 5.14 Wave tracking during host vehicle oscillation, host vehicle altitude of a) 250 ft., b) 500 ft., c) 1000 ft., and d) 2000 ft. ..................................... 143 viii NOMENCLATURE A0 Þ general amplitude of oscillation A,B,C Þ augmented equation of motion matrix C Þ coordinate transformation from segment to host frame Cß,Sß, Þ shorthand notation for cosine and sine of generic angle ß cf Þ cable friction coefficient cp Þ cable pressure coefficient d Þ cable diameter Fx,Fy,Fz Þ X,Y,Z component of the aerodynamic force, in the cable segment frame f Þ generic function h Þ wave height in host vehicle frame hn Þ sea level in host vehicle frame I Þ ATV inertia matrix, in ATV frame I’ Þ ATV inertia matrix, in ATV segment frame IX,IY,IZ Þ ATV X,Y,Z inertia components, in ATV frame kd Þ derivative gain kp Þ proportional gain l Þ segment length L,M,N Þ aerodynamic moments about X,Y,Z, in ATV frame MX,MY,MZ Þ moment about X,Y,Z, in segment frame ix Mf& Þ moment matrix about the Euler rates f& Mq& Þ moment matrix about the Euler rates q& My& Þ moment matrix about the Euler rates y& m Þ mass n Þ number of cable segments p,q,r Þ angular rates about X,Y,Z, in ATV aircraft frame Qi Þ generalized force for the ith generalized coordinate q Þ generalized coordinate t Þ time T Þ kinetic energy Trot Þ rotational kinetic energy U0,V0,W0 Þ freestream flow, in inertial frame u,v,w Þ perturbation velocities, in ATV aircraft frame V Þ velocity Vx, Vy, Vz Þ velocity in the X,Y,Z direction, in segment frame Xu,…Np Þ ATV stability derivatives X,Y,Z Þ Cartesian coordinate directions, in host frame x,y,z Þ Cartesian coordinate directions, in segment frame z Þ vertical host aircraft displacement, inertial frame l Þ wavelength p Þ circular constant r Þ air density x f,q,y Þ Euler angles ? Þ Angle between the ATV segment frame and ATV aircraft frame w Þ vector of Cartesian angular rates w0 Þ general frequency Dh Þ difference in ATV and wave altitude, in host frame Diw Þ wing angle correction 1 CHAPTER 1 1INTRODUCTION 1.1 Background The aerial tow system is a relatively common flight configuration. It is a necessity for gliding aircraft to achieve altitude; it is common to have targets towed by aircraft for weapons testing and target practice; it is even being proposed by NASA as a platform for reusable launch spacevehicles. The configuration has been proven to be well suited for these uses, though they are primarily a means of transportation of the towed vehicle. Figure 1.1 The aerial tow system 2 The aerial tow system is composed of three major components (1) host vehicle, (2) aerial tow vehicle (ATV), and (3) connecting cable, figure 1.1. In general, the host vehicle is larger than that of the tow vehicle and provides the thrust for both aircrafts. The cable is connected to both the host vehicle and the ATV in such a way that both aircraft can safely maintain flight. It is not uncommon in this configuration for the connecting cable to be hundreds or thousands of feet in length. As a result, the cable orientation is subject to change due to variations in the freestream conditions, cable dynamics, and the dynamics of the two aircraft s. The thrust produced by the host aircraft’s engines is directly applied to the host aircraft and indirectly applied to the ATV thru the tension in the connecting cable. The tension in the cable produces a force on the ATV that may differ from the host vehicle thrust in magnitude and direction depending on the cable orientation. The variable cable force on the ATV causes a dynamic response resulting in further deformation of the connecting cable. Clearly, the conditions effecting these variations are coupled. Historical research suggests that, in general, the system is stable, and given sufficient time the disturbances in the system will damp out. However the degree to which the atmospheric and dynamic alterations affect the system depends principally on the system components and is, of yet, not well known. Recently the aerial tow system has been chosen as a means to carry out the measurement of atmospheric conditions in close proximity to the ocean surface. The proposed system calls for a general aviation aircraft to tow the measurement instrumentation equipped ATV near (or approximately 30 feet) above the ocean surface. The emphasis of this system is the controlled altitude of the ATV. This is a deviation of 3 the aerial tow system from its roots as a means of transportation. Using the aerial tow system simply for transportation places no strict requirement on the motion of the ATV. The aerial tow system is proposed as a safer, more robust alternative for recording oceanic atmospheric conditions. The traditional method necessitates an aircraft fitted with the appropriate sensory equipment fly immediately above the ocean surface. The aircraft’s autopilot cannot be relied upon to safely maintain such altitude in this environment. Variable terrain and unknown atmospheric conditions could, indeed, result in disaster for the aircraft. As a result, the mission puts forth a tremendous workload upon the pilot. Due to the intense flight regime and pilot stress, the testing is usually limited to short duration and is confined to calm conditions. The implementation of the aerial tow system as a method to take oceanic atmospheric measurements would allow the piloted host aircraft to fly at a safe height above sea level (such as 250 to 2000 feet) under the control of the aircraft autopilot while towing an unmanned sensory aircraft a matter of feet above the ocean. This would allow for longer tests runs, eliminate pilot stress, expand the mission envelope, as well as reduce the worstcase scenario to the mere loss of the unmanned aircraft, assuming this flight configuration is suitable for the task at hand. Given the possibility of difficult ocean and atmospheric conditions, such as swales and wind gusts, the controllability of the tow system is in question. The system component interactions are complex, and neither the system configurations for the most advantageous performance nor the circumstances at which the system fails to adequately 4 perform are obvious. As a result, a system simulation preceding construction and testing is desirable. Each component of the aerial tow system must be mathematically modeled in order to build a numerical system simulation. Simultaneous derivation of the component equations of motion will produce the governing equations for the aerial towed system from which a response to initial and applied conditions can be simulated. The utilization of a system simulation would allow for prediction of the system’s ability to perform under various conditions such as, wind gust, wave shape, and host aircraft oscillation. Previous aerial tow systems do not place such heavy emphasis on the maintenance of the altitude of the ATV. The unique mission of near ocean atmospheric measurement has driven the design of a mathematical simulationof the low altitude aerial tow system. 1.2 Design Goals The objectives of this study is to create a dynamic simulation of the aerial tow system with the main focus surrounding the derivation of the system and component equations of motion in both two and threedimensions. The study includes mathematical models for each of the system components as well as the methodology for the derivation of the equations of motion for each of the component models mentioned. The system will be adequately designed such that the effect of the following conditions on the tracking ability of the ATV may be simulated: · Variable surface oscillations (ocean swales) 5 · Variable host vehicle altitudes · Lateral gusts · Host vehicle oscillation The decisions upon the design of the system will be made with the versatility of the simulation in mind, in the hopes that this research will assist with a wide range of future studies of aerially towed systems, underwater towed systems, or cable dynamics. 1.3 Literature Review The study of the tow system and its components has been well documented. The interest of cable dynamics, alone has produced analysis on the topics of cable oscillation and damping due to uniform flow, threedimensional cable shaping, cable drag model validity, and numerical simulations of cable dynamics. Refs 14. However, the integrated effect of the ATV is necessary to simulate the aerial tow system. As various tow systems come into existence study of the systems have flourished. Papers on tow vehicle controllability, stability augmentation, and release and retrieval dynamics build upon the dynamics of the cable in a flow, Refs 5 and 6. The investigation of tethered aerostat dynamics and underwater towing illustrates the range of the application and research of the tow system, Refs 79. Each of these references provides greater understanding of the tow system and its components. 6 1.3.1 Cable Dynamics The cable has found its way into many applications involving fluid flow through out the years. As a result, countless papers have been written for discussio n on various cable related topics. In general, the research was initiated for reasons supplementary to the cable’s effect on aerial tow systems; however, each gives insight into the overall behavior of the cable in fluid flow. The past studies give a foothold for future system analysis containing cable components. Of particular interest are the studies on the physical and aerodynamic models of the cable. 1.3.1.1 Cable Models A prediction of the dynamic response of a cable to given conditions requires a mathematical model of the cable from which the equations of motion can be derived. Choo and Casarella10 produced a survey of various methods for modeling cables. The motivation behind this article was to acknowledge various methods that could be used for a basis of cable system simulations and analyzing the merits and demerits of each technique. Among the methods discussed is the finite element method. This method breaks the continuous cable into an arbitrary number of segments which themselves can be modeled in a number of ways, such as the ideal pendulum, a mass on a spring, a cylinder, or curved segment model. In any case, each segment is free to rotate about its neighboring segments. Discretizing the cable allows for the derivation of the equations of motion of the interacting cable segments. The shortcoming to this technique is that the quality of the simulation is directly dependant upon the number of segments used to model the cable; a 7 finer discretization produces a superior simulation. However, as Walton and Polachek18 discuss, when the segment length decreases a smaller step size is necessary for numerical convergence during simulation, which translates to more computation time. The benefit of using the finite element technique is that it can be used to study several types of unsteady motion in a cable or cablebody system. It is for this primary reason the Choo and Casarella10 labeled this method as the most versatile. Also discussed in the survey are the method of characteristics, the linearization method, and the equivalent lumped mass method. The method of characteristics employs constitutive laws such as Hooke’s law, to convert the partial differential equations of motion into ordinary differential equations of motion. This technique allows the study of any sort of unsteady cable motion as long as the constitutive laws are satisfied, as an example; this method cannot be used on inextensible and viscoelastic cables as there are fewer ordinary differential equations than are original partial differential equations. In addition, the method of characteristics requires a large amount of computation. The linearization method, as it suggests, linearizes the cable equations of motion. This linearization limits the technique to perturbation, small deviations from equilibrium, studies, such as stability and frequency response of the system. This technique is not applicable in cases where a limit cycle persists and additionally is not useful in the analysis of unsteady cable motion. The equivalent lumped mass technique is useful when the most important aspect of the cablebody system is the body. In this method, the dynamics of the cable are completely ignored and the equations of motion of the system are largely due to the rigidbody 8 equations of motion of the body with added terms representing the cable dynamic effects. Since the equations of motion of the cable are omitted the system is greatly simplified. The simplicity, unfortunately, comes at the cost of the robustness of the model. Each of these methods has a configuration or condition in which it gains advantage over the rest. However, the only one of the techniques with no systemtype limitations is the finite element technique. An additional detailed description of the finite element method will be given in the following chapters. 1.3.1.2 Cable Drag Model The work of Hoerner13 is directed at analyzing and modeling fluid drag about various bodies. The cable is among the bodies researched. Hoerner13 produced the crossflow method for calculating the drag forces on a cable within a fluid. This technique analyzes the friction and pressure forces on the cable. The method is limited to straight cables of circular crosssection and fluid conditions that produce subcritical Reynolds numbers with respect to the cable diameter. Despite the straight cable limitation of the cross flow principle it can predict the drag force on a flexible cable (even though the cable is rarely completely straight). A curved cable can be approximated by a series of connected straight segment. Given the orientation of the straight segment relative to the fluid flow, the drag on each segment can be predicted using the crossflow principle. Thus, the resultant fluid force on a curved cable is the accumulation of the individual segment’s drag forces. As the number of cable segments increase so to does the accuracy of the overall drag prediction on the curved cable (to an extent). 9 The crossflow principle works well with the finite element model of a cable (as it necessitates breaking the cable into straight cable segments). The applied drag forces on each of the cable segments can easily be applied to the dynamic equations of motion for each segment. The exact equation for the crossflow principle in its original twodimensional form and modified threedimensional form are presented in chapters two and three. 1.3.2 Lagrange’s Equations Newtonian mechanics states that the motion of a particle can be determined provided that there is a complete knowledge of the external forces acting on the particle and that the initial conditions of the particle are known; the same is true for groups of interacting particles. This is a fundamental step in producing the equations of motion necessary to create a simulation of a system in response to inertial and external forces. The development of the equations of motion of a system can be trivial or quite complex depending on the composition of the system and its environment. The direct application of Newton’s laws becomes difficult as the number of groups of particles needed to represent the system increase. The motion of the system components and internal and external forces acting on the system exist as vectors, with magnitude and direction. Thus, keeping track of the external forces, internal forces, and motion of the system components in vector form is often difficult. A less complex method of deriving system equations of motion is attained by implementing Lagrange’s equations, Greenwood11. 10 The Lagrange method has two main advantages over the direct application of Newton’s Laws in attaining the system equations of motion. The first is that the Lagrange approach is energy based rather than spatially based, thus it deals with scalar values which circumvents the vector bookkeeping involved with Newton’s method. The other benefit of using the Lagrange equations is that, unlike the application of Newton’s laws, it specifies an exact process to develop the equations of motion for any well defined system. After defining the system and its constraints, all one needs to do to develop the system equations of motion is to carryout the manipulation of Lagrange’s equations. An important concept in the application of Lagrange’s equations is degree of freedom of the system. The number of degrees of freedom of a system is equal to the number of coordinates necessary to define the configuration of the system minus the number equations that constrain the position of the system. For example, assume the system is that of a simple twodimensional ideal pendulum. 11 a) b) Figure 1.2 Simple pendulum shown in twodimensional a) Cartesian coordinates, and b) cylindrical coordinates In twodimensional Cartesian coordinates the location of the pendulum mass is determined by its x and z positions. However, the x and z positions are constrained to a circular path around the pendulum hinge with a radius equal to the length of the pendulum (we are assuming the length of the pendulum is rigid). The constraint equation would be the sum of the square of the x and z position is equal to the square of the length of the pendulum. Thus, if the x position is known then the application of the constraint equation will produce the z position of the pendulum mass. With two coordinates utilized to specify the orientation of the pendulum and the existence of one constraint equation the degree of freedom of the system is one. Imagine the same example in cylindrical coordinates. The position of mass of the pendulum is determined simply by the angle q of the pendulum relative to the coordinate 12 frame. There is then one coordinate necessary to define the pendulum orientation and no constraint equation resulting in one degree of freedom for the pendulum. As the pendulum example shows, the degree of freedom of the system is independent of the coordinates chosen to represent the system. The coordinates needed to represent the orientation of the system are deemed generalized coordinates. The generalized coordinates of the simple pendulum in the twodimensional Cartesian coordinate system are the x and z position of the pendulum mass. The generalized coordinate of the same pendulum in cylindrical coordinates is simply the angle q. Therefore, since there are infinite possibilities of coordinate systems that could be used to analyze the pendulum there are infinite possibilities for generalized coordinates. However, in many cases there is a set of generalized coordinates in which the number of coordinates is equal to the degree of freedom of the system. This set is referred to as the independent generalized coordinates. The independent generalized coordinates require no constraint equations; as a result it is often simpler, mathematically, to analyze the system using the set of independent generalized coordinates. In the pendulum example the angle q is the independent generalized coordinate. In a dynamic system the generalized coordinates, independent or not, vary with time and are algebraic variables in Lagrange’s equations. As stated previously, Lagrange’s equations are energy based. Prior to manipulation of Lagrange’s equations the energy of the system must be represented in terms of the generalized coordinates. The manipulation of Lagrange’s equations requires derivation of the system energies with respect to time and the generalized coordinates. The exact form of Lagrange’s equations 13 as well as the execution of Lagrange’s equations will be undertaken in detail in the next two chapters. 1.3.3 Aircraft Dynamics There has been much effort spent by researchers and scientists in the study of aerodynamics as it relates to aircraft. This work ranges from simple analysis of fluid flow over a twodimensional airfoilshape to the development of nonlinear equations of motion for high speed maneuverable aircraft. The extensive work accomplished for a wide variety of aircraft and flight conditions has produced reliable models for the prediction of aircraft response to fluid forces. The method used to represent the dynamics of the aircraft in this system is smalldisturbance theory, Nelson14. Smalldisturbance theory assumes that the motion of an aircraft consists of steadystate motion and small perturbations about its steadystate. Smalldisturbance theory linearizes the rigid body equations of motion. This is accomplished by substituting reference values plus a perturbation for the all variables in the rigid body equations of motion. Many reference values (e.g. pitch, roll, etc.,) can be set to zero by assuming a symmetric reference flight condition. Simplification of the equations of motion after substitution produces products of perturbation. These nonlinearities are ignored in the assumption that higher order small disturbances will produce negligible overall effects. Further simplifications of the equations of motion produce aerodynamic stability derivatives. The stability derivatives represent the change in forces or moments of a body due to some change in the flight conditions, such as forward velocity or angle of attack. 14 The stability derivatives are approximated based upon the size, shape, and mass distribution of an aircraft and its components. Methods for approximation of aircraft stability derivatives are presented in Nelson14 and Raymer16. The equations of motion for aircraft are written in terms of normalized stability derivatives. As a result, these general equations apply to a wide variety of general aviation aircraft over most flight conditions. However, since the equations were linearized using small disturbances, the theory breaks down and thus does not produce accurate predictions when the motion of the system involves large amplitude motions such as aircraft stall or spinning. Given the simple shape of the ATV and its proposed flight conditions the smalldisturbance theory will be adequate in modeling dynamics of the aircraft. 1.3.4 Cable Towed Systems A cable tow system in its most simple form consists of a body constrained by cable in a fluid. There are a wide range of system configurations that meet the requirements of a cable tow system, such as kites, tethered aerostats, mooring lines, towed buoys (above and below the ocean surface), gliders, and towed aircraft targets. Research on many of these systems has been accomplished in the past. Each system has its own nuances that separate one author’s work from another, and although not all are directly applicable to the aerial towed system they help to illustrate the diversity of the system. A few of the papers that do apply to the aerial tow system are discussed below. Huffman and Genin1 research the stability of a simple aerial towed system. Their research suggests that there is no direct instability in the dynamics of the towed cable 15 system. In each case tested the disturbances of the system were damped out in time. They concluded that a likelihood of instability is present, however, when the dynamic frequency of the towed body matches the natural frequency of the cable. The result would be a limit cycle oscillation. Huffman and Genin1 also show that the natural frequency of cable can vary based upon parameters such as cable length and flight speed. As a result, the limit cycle behavior could be present in tow system at particular flight conditions. However, this depends largely on many variables in the system. As a result, there was no generalized parametric study. Cochran et al.6 investigate the effect of controlling an aerial towed vehicle. The authors designed a control system that would augment the stability of the towed vehicle. Simulations were run in conditions that produced unstable dynamic motion in the uncontrolled towed vehicle (perturbed retrieval of the towed aircraft). When the control system was activated the authors showed that the system damped out quickly. This paper addresses the difficulties of designing a controlled aerially towed system and the results suggest that improved stability can be achieved during maneuvers that could otherwise result in instability. Nakagawa and Obata15 studied the longitudinal dynamic modes of the aerial tow system. They analyzed the modes of three flight configurations, a straight cable connected to a sphere, a straight cable connected to an aircraft, and a curved cable connected to an aircraft. The results of the study helped to classify the types of motion of the aerial tow system. In all of the configurations the main dynamic modes were the pitching mode, pendulum mode, and first vibration mode. The pitching mode is simply 16 the pitching of the towed body; this has little effect on the cable dynamics or the translation of the towed body. The pendulum mode oscillates the towed body about its static configuration. The motion causes slight bending in the cable and little if any towed body pitching. In both of the straight cable configurations the first vibration mode has little effect on the position or the orientation of the towed body, and simply causes oscillations in the cable. However, in the curved cable configuration the first vibration mode causes the towed aircraft to surge fore and aft in an elliptical pattern. The authors labeled this as the bowing mode. The authors complete a parametric study of the aircraft longitudinal stability derivatives and found that two of them, Zw and Mw as well as the cabletowed body hitch point played an important roll in the stability of the system. The results sho wed a strong relationship between the stability derivatives and the bowing mode. As the absolute values of Zw and Mw increase and decrease, respectively, the bowing mode becomes less stable, and eventually an unstable phenomenon known as pitching or bowing flutter occurs. Nakagawa and Obata15 noticed a similar behavior as the cable hitch point on the towed aircraft was moved more and more fore of its center of gravity. The system dynamic oscillations are of the utmost importance when attempting to fly the towed body very near the ocean surface. The research completed by Nakagawa and Obata15 suggest care should be taken when designing the ATV as its stability derivatives have a direct effect on unwanted system oscillations. Lawhon and Arena19 investigated the effects of aircraft size and shape as well as cable parameters to develop the optimum 17 system configuration for aerially towed systems tracking very near the ocean surface. The ATV used in this paper was a result of the work done by Lawhon and Arena19. 18 CHAPTER 2 2TWODIMENSIONAL DYNAMICS 2.1 Introduction It is the focus of this paper to derive and produce the mathematical equations that accurately describe the threedimensional motion of the aerially towed system. Followed by, implementing the equations of motion to produce a dynamic simulation of the aerially towed system. The simulation will then be used to draw general conclusions about the system and its ability to perform the maneuvers desired to obtain atmospheric measurements near the oceans surface. The individual dynamic behavior of the system components and the interactions between them will result in quite complex set of equations of motion. As a result of this complexity, it is easy to get bogged down by the shear volume of the equations and lose sight of the steps taken to derive the final set of equations. This is especially true in reference to the threedimensional system. Therefore, the twodimensional equations of motion will be laid out in this chapter. The twodimensional equations of motion will also act to verify that the more complex threedimensional equations of motion are derived without error. It is important to note that the same principles will be held for the derivation of both the two and threedimensional systems although they may differ in clarity. 19 The basis of any simulation is the modeling of the system components. The quality of the models directly affects the value of the simulation. The derivation of the equations that dictate the dynamics of a component is possible upon the completion of its mathematical model. To accurately model the system, the individual component dynamics are important however, the modeling of the system necessitates the equations of motion represent the dynamic interaction between the components as well. Thus, in general, the derivation of the equations of motion of a system component must be accomplished simultaneously with the dynamic equations of the remaining objects. The result will be a set of equations that represent the dynamics of the entire system, i.e. if the host vehicle oscillates the system equations of motion will reflect the oscillations effect on the dynamics of the cable and the ATV. The following sections discuss the modeling of the system components and the development of the individual dynamic equations of each component. As stated previously, the individual component dynamics are not enough to model the system; however, they may be of academic interest, and will lead to a clear understanding of the derivation of the system equations of motion. 2.2 Cable Dynamics As discussed in the literature review there has been extensive study on cable dynamics and cable modeling for a wide range of applications. Due its versatility the finite element technique will be used to model the connecting cable. The main drawback of the finite element technique is its computational expense. The great strides in computer speed over the years helps to soften this blow. There may also be situations when the operator is not 20 interested in a completely rigorous simulation but rather just a feel of the system’s behavior, this can be accomplish quickly by using fewer segments. As mentioned previously, the finite element technique requires that the cable be broken into an arbitrary number of segments which are free to rotate about their neighboring segments. The types of segments used to model the cable can vary. Notable models include the simple pendulum, spring mass, thin rod, and curved structure model. Figure 2.1 Finite element cable segment models In the simple pendulum model, the segment is straight and inextensible. The mass of the segment is lumped to a single point at the tip of the segment, resulting in a model that resembles that of a simple pendulum. Being that the segment is straight the crossflow method can be used to predict the drag force on the segment. The applied force on the segment is exerted at the point mass or node. The inextensible constraint lends this model to a system in which the strain on the cable is minimal. The benefit of this technique is its simplicity. The draw back is that placing the inertial and applied forces at the tip of the segment only becomes a good approximation when the number of segments used is large. 21 The spring mass segment also lumps the mass at tip of the segment, however the segment length acts as a spring. This allows for strain effects to be incorporated into the equations of motion of the cable. The change in length of each segment adds an additional degree of freedom to the system, which increases the complexity of the equations of motion and increase computational intensity. If the strain on the cable can be neglected the spring mass segment model expends much unneeded effort and time. The thin rod model is similar to the lumped mass model. It too is straight and inextensible. However, the model places the inertial forces and applied forces at the center of mass of the segment. The additional complexity of this model is that the mass of the rod is not lumped to a point; therefore, it has inertia that impedes segment rotation and must be taken into account the derivation of the equations of motion. The benefit of the better physical model is it will take fewer segments to accurately represent the motion of the continuous cable than the lumped mass model. The last model represents the cable segment as a curved structure. This technique models the segments as polynomials or splines which require slope continuity at the segment connections. Due to the fact that most cable configurations are curved this method reduces the number of segment needed to attain an accurate simulation. However, the evaluation of the drag on a curved cable is not well known. The more complex drag model and polynomial or spline calculations make derivation of the cable equations of motion difficult. Neglecting the strain in the cable is acceptable as we are dealing with a fairly light ATV and a high tensile steel connecting cable, thus the spring mass model unnecessary; 22 and given that there is a historically accepted method to model the drag on a straight cable the complexities of the curved structure will be avoided. The equations of motion of the cable will be derived using the simplest segment model, the pendulum model, shown in figure 2.1. However, a discussion about the changes to the cable equations of motion as a result of modeling the cable segments as thin rods will take place at the end of this chapter. Figure 2.2 Lumped mass finite element cable model 2.2.1 Coordinate Systems In this section the primary concern is the development of the twodimensional equations of motion for a flexible cable. As a result, the host aircraft and the ATV will be neglected for now. The cable is broken into n arbitrary segments in accordance with the finite element cable model. Therefore, one end of the cable is fix to a point in space and is free to rotate while the other is end is constrained only by its preceding segments. The coordinate system at the cable attach point is the inertial coordinate system. This coordinate system is orthogonal. The X axis is locally parallel to the ‘ground’ and the Z 23 axis is pointing toward and is normal to the ‘ground’. The derivation of the equations of motion of the flexible cable will be carried out in the inertial frame. As a result, it is necessary to be able to describe the orientation of the cable in this coordinate system. The location of the ith cable node in the inertial frame can be resolved using the length, li, and angle of orientation, qi, of the segment and each preceding segment, as shown below. å å = = = = i j i j i j i j j j Z l C X l S 1 1 q q (i = 1,2Kn) (2.1) where, lj is the length of the jth segment and Xi and Zi are the position of the ith segment in the inertial frame. Also note C and S are short notation for cosine and sine, this notation will be used throughout the remainder of this paper. Aside from the inertial frame, each cable segment has its own coordinate system. The cable segment coordinate system is orthogonal and is configured such that the node of that segment is located in the positive zdirection a distance equal to that of the cable segment. For example, the origin of the coordinate axis for the first cable segment is located at the exact same point as that of the origin of the inertial frame. The two frames differ only by the angle rotation of the segment, q1. When the angle is equal to zero the two frames are identical. In a similar manner, the origin of the second cable segment is located at the node of the first cable segment, and the first and second frame differ rotationally by the difference of q2 and q1. This is evident from Figure 2.2. 24 Figure 2.3 Twodimensional inertial and cable segment coordinate systems Throughout the analysis of the system it will be necessary to convert from the inertial frame to the segment frame. This is accomplished by utilizing a rotational coordinate transformation. The rotational coordina te transformation from the ith cable segment frame to the inertial frame is represented by the directional cosine matrix as follows, úû ù êë é  = i i i i i S C C S C q q q q q (2.2) The inverse of the directional cosine matrix produces a rotational coordinate transformation from the inertial frame to the ith cable segment frame. Calculating the inertial position of the node of the ith cable segment is a common implementation of the rotational transformation. Recall that in the cable segment frame the node is located at a distance equal to the length of the segment in the zdirection. Using the directional cosine matrix the position of the ith cable segment node in the inertial coordinate system is as follows, 25 å= úû ù êë é úû ù êë = é úû ù êë é i i j j i l C Z X j 1 0 q (i = 1,2,K, n) (2.3) Note that equation 2.3 and 2.1 are equivalent. 2.2.2 Lagrange Equations As discussed earlier, the connecting cable is broken into n arbitrary number of segments, and each segment has its own coordinate system. Deriving the equations of motion for the cable using the direct application of Newton’s laws of motion becomes increasingly difficult as the number of segments used to model the cable increases. This is primarily due to the fact the Newton’s laws of motion are vector based and as the system increases the number of coordinate systems increase resulting in labor intensive derivation. Another difficulty when using the Newtonian approach is that there is no general method to derive the equations of motion. As discussed in the previous chapter, an alternative to using Newton’s laws of motion is the utilization of Lagrange’s equations, which is energy based. Consequently, it circumvents many of the complexities of vectorbased derivation, as well as, gives a standard process for developing the system equations of motion. To recap, the derivation of Lagrange’s equations depends on the degree of freedom of the system, and the degree of freedom of the system is the number of coordinates needed to define the system minus the number of constraint equations of the system. In the Cartesian coordinates of the inertial frame for a cable broken into n segments there are 2n coordinates needed to define the location of the node of each segment. However, each segment has a constraint equation of the form, 26 2 1 1 2 1 1 2 ÷ ÷ø ö ç çè æ  + ÷ ÷ø ö ç çè æ = å å =  = i j i j i j i i j l X X Z Z (i = 1,2,K, n) (2.4) As a result, 2n coordinates and n constraint equations produce n degrees of freedom. In other words, the twodimensional finite element cable has a degree of freedom equal to the number of segments used to model the cable. The coordinates used to define the system are the generalized coordinates. The generalized coordinates are general due to the fact that there are an infinite number of coordinate systems that could be used to define the system all with varying number of constraint equations (although every system has n degrees of freedom). Consequently, there are no specific coordinates that need to be used to define the system. The implementation of Lagrange’s equations will be simplified if the independent generalized coordinates are used (recall, the independent generalized coordinates are the set of coordinates that define the system yet have no constraint equations). For the flexible cable system the independent generalized coordinates are the orientation or attitude angles of each of the segments, qi (i = 1, 2… n). For the duration of the chapter the attitude angles of the cable segments will be used and referred to as the generalized coordinates of the system. The generalized coordinates will be generically referred to by the variable q. As stated previously, Lagrange’s equations are energy based. The Lagrangian function, L, represents the energies of the system and is equivalent to the kinetic energy minus the potential energy of the system, 27 L = T U (2.5) where, T and U are kinetic and potential energy functions, respectively. Lagrange’s equations are as follows, i i i Q q L q L dt d = ¢ ¶ ¶  ÷ ÷ø ö ç çè æ ¶ ¶ & (i = 1,2,K, n) (2.6) where, i Q¢ represents the applied, or generalized, forces that are not derivable from a potential function, such forces include friction forces and forcing functions. However, since the potential function for this system is not velocity dependant Lagrange’s equations can be reduced to the following, i i i i i Q q U Q q T q T dt d = ¶ ¶ = ¢  ¶ ¶  ÷ ÷ø ö ç çè æ ¶ ¶ & (i = 1,2,K, n) (2.7) where, Qi are the generalized forces that include applied forces and the inertial forces derived from the potential function. Note that Equations 2.6 and 2.7 is actually a set of n equations, one equation for each generalized coordinate of the system. Thus, in the case of the twodimensional finite element modeled cable there is one equation for each attitude angle of the n cable segments. Therefore, Lagrange ’s equations can be rewritten as, i i i Q T T dt d = ¶ ¶  ÷ ÷ø ö ç çè æ ¶ ¶ q& q (i = 1,2,K, n) (2.8) 28 In order to implement Lagrange’s equations it is necessary to get the energies of the systems in terms of the generalized coordinate. Only then can the proper derivation of the equations of motion take place. The kinetic energy of the system at any given time is due to the summation of the rotational and translational motion of the inertial components of the system. For the finite element cable the inertial components are simply the segment nodes. Thus, the kinetic energy of the cable is the summation of the nodal rotational and translational kinetic energy. However, since the mass of a segment is localized to a point it has no rotational inertial and consequently has no rotational kinetic energy. As a result, the kinetic energy of the system is solely a function of the translational velocity of each node. å= = n j j j T m V 1 2 2 1 (2.9) The velocity components of a node can be found by simply taking the time derivative of its position. å å = = =  = i j i j j i j i j j j j Z l S X l C 1 1 q q q q & & & & (i = 1,2,K, n) (2.10) Equation 2.10 shows that the translational velocity of a node results from the rate of rotation of its segme nt and each preceding segment. At this moment, the velocity of each node is in vector form. However, squaring the velocity for the kinetic energy releases 29 any further vector bookkeeping. The square of velocity, V, for the kth segment is can be written in series form as follows, ( ) ( i j i j ) ( i j ) V X Z l l C C S S l l C j i j k j k i j i j i k j k i k k k i q q q q q q q q q q  = = = = = & + & =åå & & + =åå & & 1 1 1 1 2 2 2 (2.11) Similarly, the total kinetic energy for n cable segments is written in summation form as, å ååå ( ) =  = = = = = n k j i j k j k i k i n k k k i j T m V m l l C 1 1 1 1 2 2 1 2 1 q q q&q& (2.12) Being that the kinetic energy has been derived in terms of the generalized coordinates the implementation of Lagrange’s equations can commence. The derivative of kinetic energy with respect to the ith generalized coordinate and the time rate of change of the ith generalized coordinate yields the following, åå ( ) = =  = ¶ ¶ n k k j k i j i j i j i m l l S T 1 1 q q q q q & & (2.13) åå ( ) = =  = ¶ ¶ n k k j k i j j i j i m l l C T 1 1 q q q q & & (2.14) Finally evaluating the time rate of change of equation 2.14 produces, åå [ ( ) ( ) ( )] = =     = ÷ ÷ø ö ç çè æ ¶ ¶ n k k j k i j j j j i i j i j i m l l C S T dt d 1 1 q q q q q q q q q & && & & & (2.15) Combining equations 2.13 and 2.15 results in the series representation of the left hand side of the Lagrange’s equations for the twodimensional finite element modeled cable, 30 åå ( ( ) ( ) ) = =   =  ¶ ¶  ÷ ÷ø ö ç çè æ ¶ ¶ n k k j k i j j j i i j i j i m l l C S T T dt d 1 1 2 q q q q q q q q & & & (i = 1,2,K, n) (2.16) The right hand side of equation 2.8 represents the generalized forces of the system. Note that due to the fact that the generalized coordinates of the system are angles the general forces are actually moments. These moments are about the pivot point of each segment and result from the inertial and applied forces on the node of the segment as well as any applied moments. The generalized force equation is derived from the concept of virtual work. Suppose forces F1, F2, …, Fk are applied at position x1, x2, …, xk and moments M1, M2, …, Mk are applied at angles f 1, f 2, …, f k. Now let us assume that the system undergoes arbitrary small displacements dx1, dx2, …, dxk and df 1, df 2, …, df k. The work done by these applied forces and moments is known as virtual work and the small displacements are known as virtual displacements, å== + k j j j j j W F x M 1 d d dj (2.17) The virtual displacements can be rewritten in terms of the virtual displacement of the generalized coordinates, å= ¶ ¶ = n i i i j j q q x x 1 d d (2.18) å= ¶ ¶ = n i i i j j q 1 q d j dj (2.19) 31 The relationship between the generalized forces and the virtual work of the system due to the virtual displacement of the generalized coordinates is as follows, åå =å ¶ ¶ + ¶ ¶ = = = n i i i k j n i i i j i j i j j q Q q q q M q x W F d d j d d 1 1 (2.20) Simplification of equation 2.20 produces the generalized force equation for the applied forces and moments in terms of the generic variables x and f . å= ¶ ¶ + ¶ ¶ ¢ = n i i i j i j i j i j q q q M q x Q F 1 d j d (2.21) In the case of the twodimensional cable segment system the position variables are X and Z while the angles are q. Incidentally, the generalized coordinate for this system is also q. The equation for the generalized forces not produced by the potential function is as follows, å= ¶ ¶ + ¶ ¶ + ¶ ¶ ¢ = n j i j i j Z j i j i X j i M Z F X Q F 1 q q q q q (i =1,2,K,n) (2.22) The total generalized force equation is made up of the external forces and moments as well as the inertial forces, å= ¶ ¶ + ¶ ¶ + ¶ ¶ + ¶ ¶ =  n j i j j i j Z j i j X j i i M Z F X F U Q 1 q q q q q q (i =1,2,K,n) (2.23) where, X j F and Z j F are applied forces and j Mq is the applied moment of the jth segment in the inertial coordinate frame. 32 In the expansion of equations 2.23 let us look first at the term that is the derivative of potential energy with respect to the ith generalized coordinate. The potential energy of the cable is the summation of the potential energy of each segment node. åå ( ) = = =  n k k i U mk g Zi D 1 1 (2.24) where, D is the datum line from which potential energy is measured. Note that since potential energy is being differentiated the location of the datum line is inconsequential. Taking the derivative of equation 2.24 with respect to the ith generalized coordinate produces, å= ¶ ¶ = ¶ ¶ n j i j j i Z m g U q 1 q (2.25) The resulting equation allows the inertial force to be included with the zdirection applied force. Simplify the generalized force equation to that shown below, ( ) å= ¶ ¶ + ¶ ¶ + + ¶ ¶ = n j i j j i j Z j j i j i X j M Z F m g X Q F 1 q q q q q (i =1,2,K,n) (2.26) The derivative of position and angle with respect to the ith generalized coordinate is as follows, 33 ( ) ( ) ( j i) l S j i Z l C j i X i j i i j i i j i i = = ¶ ¶ =  ³ ¶ ¶ = ³ ¶ ¶ 1 q q q q q q (2.27) Note that the position of the node of a segment is dependant upon that segments attitude angle and the attitude angle of each preceding segment. Due to the nature of the dependence of position on the generalized coordinate, the evaluation of the derivative of the position of the jth segment with respect to the ith generalized coordinate results in zero when i is greater than j. The derivative of angle with respect to the ith generalized coordinate is equal to one when i is equal to j and zero otherwise. Thus, ( ) å= = +  + n j i i i X j i i Z j j i i Q M F l C F m g l Sq q q (i =1,2,K,n) (2.28) The complete derivation of Lagrange’s equations for the twodimensional finite element modeled cable is attained by combining equations 2.16 and 2.28. The resulting set of n equations of motion is presented as follows in summation form. ( ) ( ) ( ) ( ) (i n) M F l C F m g l S m l l C S n j i i X j i Z j j i n k k j k i j j j i i j i j i 1, 2, , 1 1 2 K & & = +  +  = å åå = = =   q q q q q q q q q (2.29) The set of equations of motion can be more clearly presented in matrix form, yielding, 34 [A][q&]+ [B][q&2 ]= [Q] (2.30) where, A and B are coefficient matrices and are size n by n. The remaining three matrices are column matrices of size n by 1. The values of the coefficient matrices are as follows, i j i j i j A l l C m , (q j qi ) , = (2.31) i j i j i j B l l S m , (q j qi ) , =  (2.32) where, the subscripts i and j represent the row and column of the appropriate matrix, and i j m , represents the mass distribution for row i and column j. For the lumped mass model of the n segments the mass distribution is simply a sum of segment masses. å = = n k i j i j k m m max( , ) , (2.33) For example, the equations of motion for a cable modeled using two cable segments without fluid drag would be as follows, ( ) ( ) ( ) ( ) ( ) ( ) úû ù êë é   + = úû ù êë é úû ù êë é   + úû ù êë é ú úû ù ê êë é +     2 1 1 2 2 1 1 2 2 1 2 2 1 2 1 2 2 2 1 2 1 2 2 1 2 2 1 2 2 1 2 2 2 2 1 2 2 1 2 1 0 0 q q q q q q q q q q q q q q m gl S m m gl S m l l S m l l S m l l C m l m m l m l l C & & && & (2.34) These two equations are analogous to the equations of motion of a twodimensional ideal double pendulum. 35 The resulting fully nonlinear twodimensional equations of motion describe the dynamics of the system due to inertial and applied loads, as well as, the dynamic interaction between all segments; and the set of equations place no limitations on the number of segments of the system (while the equations hold true independent of the number of segments the effects on simulation computation varies, this will be discussed in chapter 5). The angular acceleration of each segment can be calculated given appropriate initial conditions, such as segment attitude angle and angular velocity of the system. The implementation of a numerical integration routine on the set of ordinary differential equations allows the prediction of cable motion for various initial conditions, and with the addition of a drag model to predict the applied friction forces on the segments the simulation can be extended to that of a cable in fluid flow. 2.2.3 Aerodynamic Forces The generalized force matrix in the equation 2.30 is composed of the inertial forces and applied forces of each of the cable segment. These applied forces are the friction forces produced by the fluid interaction with the cable. Thus to fully derive the equations of motion for a cable in a fluid a drag model is necessary. As stated previously, the study of cable drag has a history that dates back to the early aircraft and thus has been well documented. The cable crossflowprinciple, Hoerner13, is a drag model that was presented in the sixties and has since been used as a basis for numerous research papers involving fluid immersed cables. Its limitations require the cable crosssection be circular, the cable be straight, and fluid be at a subcritical Reynolds number. In general, aircraft tow cables are 36 circular in crosssection; due to our relatively low airspeeds the Reynolds number is subcritical; and our finite element model for the cable requires that the cable be broken into straight segments. Therefore, the Hoerner crossflowprinciple applies very nicely to our model of the cable system. The crossflowprinciple calculates the tangential and normal forces on a straight cable segment produced by tangential and normal velocities about the cable, as shown below, ( ) ( 2 2 ) 2 2 2 1 2 1 j j j j j j j j z j z f x z x j x f x z p x F dl V c V V F dl V c V V c V =  + =  + + r p r p (2.35) where, r, d, cf, and cp are air density, cable diameter, frictional coefficient, and pressure coefficient, respectively, and Fx, Fz, Vx, and Vz are the aerodynamic forces and the flow field velocity in the cable segment coordinate system, respectively. The cable velocity in the cable segment frame can be attained simply by using the rotational coordinate transformation of the relative velocity of the segment node in the inertial frame, as shown below, ú úû ù ê êë é + + úû ù êë é = úû ù êë é j j T z j x j W Z U X C V V j & & 0 0 q (2.36) Presumably the freestream velocities, U0 and W0, are known and the equations for calculating the velocity of a segment node in the inertial frame were derived in equation 2.10. 37 Note that the drag forces produced by the drag model are in the cable frame yet the forces in the generalized coordinate equation are in the inertial frame. As a result, the calculated normal and tangential drag forces must go through a rotational coordinate transformation before application can be made. úû ù êë é úû ù êë = é úû ù êë é z j x j Z j X j F F C F F j q (2.37) The implementation of the Hoerner drag model produces the mathematical model of the applied loads on a cable segment. The generalized forces for each segment can be calculated based on cable motion and numerous freestream conditions and thus allows for the dynamic simulation of a cable in a fluid. 2.3 ATV Dynamics The ATV alone is an aircraft, and a very simple one at that. As a result, its rigid body dynamics could be easily modeled due to the extensive research and documentation done in the area aircraft flight dynamics. However, in the aerial towed system the dynamics of the cable and the dynamics of the ATV are intertwined, consisting of individual dynamics and the interaction of the two components. As a result, a system of equations containing both the equations of motion of the cable and that of the ATV must be derived and solved simultaneously. This is most easily accomplished by reproducing the ATV as an additional segment, one which is attached to nth cable segment. Thus the equations of motion for both the cable and ATV can be derived using Lagrange ’s equations. Since the equations of motion of the cable were derived to be independent of the number of segments this does not change the form of the equations previously derived. The 38 variation in the set of equations arises because the ATV segment model differs from the cable segment model. The set of equations for the cable and ATV and there variations from the previous cable dynamic set will be laid out in the subsequent sections. 2.3.1 Coordinate system Modeling the ATV as an additional segment adds an extra coordinate system just as adding an extra cable segment to the cable model. This system is analogous to the cable coordinate system described previously, the origin is the cableATV hitch point and its zaxis passes through the center of mass of the ATV. However, there is another coordinate system necessary to fully describe the orientation of the ATV. This coordinate system corresponds to the standard aircraft body fixed coordinate system. The aircraft body fixed coordinate system has it origin at the center of gravity of the aircraft, and for the twodimensional case the aircraft body fixed xdirection is pointing out of the nose of the aircraft and its zaxis is pointing out of the bottom of the aircraft (presumably towards the ground). The two additional axes are present in figure 2.3. Figure 2.4 Twodimensional ATV coordinate system 39 Thus, the ATV segment coordinate system and the ATV aircraft coordinate system differ directionally by a fixed rotation angle, g. This angle is geometrically defined by the shape of the ATV as well as the location of the cableATV hitch point, figure 2.3. The ATV aircraft coordinate system is necessary to calculate the aerodynamic forces and moments of the ATV using small disturbance theory. The ATV orie ntation, as seen by the inertial frame, is determined by the ATV segment angle and ATV geometric attitude angle. Thus two rotational coordinate transformations are needed to describe the ATV attitude in the inertial frame. The rotational transformation from the ATV aircraft coordinate system to the ATV segment coordinate system is as follows, úû ù êë é ¢ ¢ úû ù êë = é úû ù êë é z x C z x g (2.38) where, x’ and z’ are the x and zdirections in the aircraft fixed frame. The rotational transformation from the ATV segment coordinate system to the inertial coordinate system is as follows, úû ù êë é úû ù êë = é úû ù êë é + z x C Z X qn 1 (2.39) As a result, the equation for rotational transformation from aircraft fixed coordinate system to the inertial coordinate system is found by combining equations 2.38 and 2.39. úû ù êë é ¢ ¢ úû ù êë é úû ù êë = é úû ù êë é + z x C C Z X qn 1 g (2.40) 40 This transformation is necessary in converting the aerodynamic loads to the generalized forces of the ATV segment (due to the fact that our generalized force equations require the applied loads be in the inertial system). 2.3.2 Lagrange equations The parameters of the additional segment, which will be henceforth deemed the ATV segment, is dictated by the size and shape of the ATV. The length of the segment is actually the distance between the cableATV hitch point and the center of gravity of the ATV. The mass of the ATV segment is simply the mass of the aircraft itself. Therefore, the equations of motion for the translation of the ATV segment require no additional derivation from that accomplished in the cable dynamics section. As a result, the complex dynamic interaction between each cable segment and the ATV are inherently represented by the equations of motion. However, the equations of motion for the ATV segment are not complete. Recall that the nodes used to model the cable are point masses, thus they have no inertia to resist rotation. The same cannot be said for the ATV. To complete the equations of motion for the cableATV system the Lagrange equation must be derived for the rotational kinetic energy of the ATV segment The complete representation of Lagrange’s equations for the cable and ATV system is as follows, ( ) ( ) i i rot i rot i i i i rot i rot Q T T dt T T d dt d Q T T T T dt d = ¶ ¶  ÷ ÷ø ö ç çè æ ¶ ¶ + ¶ ¶  ÷ ÷ø ö ç çè æ ¶ ¶ = ¶ ¶ +  úû ù êë é ¶ ¶ + q q q q q q & & & (i = 1,2,K, n + 1) (2.41) 41 Due to the fact that the total kinetic energy of the system is composed of the summation of translational and rotational kinetic energy the resulting equations of motion due to the rotational kinetic energy can simply be superimposed upon the previously derived equations of motion (i.e. the derivation of the components of the equation 2.41 that represent the translational kinetic energy need not be carried out as they will yield the same results). The analysis of the rotation kinetic energy will focus on the following portion of equation 2.41, ( ) i rot i Trot T dt d q ¶q ¶  ÷ ÷ø ö ç çè æ ¶ ¶ & (i = 1,2,K, n + 1) (2.42) The equation for rotational kinetic energy of a rigid body is as follows, T {w}T [I ]{w} rot 2 1 = (2.43) where, I is the inertial matrix and w is the angular velocity vector of the rigid body. For the twodimensional case the segment rotational equation of motion can be simplified since the only allowable rotation is about the yaxis. The rotational kinetic energy of the system is the sum of the rotational kinetic energy of each segment. The rotational kinetic energy of the twodimensional system is as follows, å+ = = 1 1 2 2 n 1 i rot i yi T q& I (i = 1,2,K, n + 1) (2.44) where, Iy is the rotational inertial of the node about its yaxis. 42 The derivation of Lagrange’s equations for rotational kinetic energy term is quite straightforward, producing, = 0 ¶ ¶ i rot T q (2.45) i yi i I T q q & & = ¶ ¶ (2.46) Finally evaluating the time rate of change of equation 2.46 produces, i yi i I T dt d q q & & & = ÷ ÷ø ö ç çè æ ¶ ¶ (2.47) Combining equations 2.45 and 2.47 produces the supplementary left hand side of Lagrange’s equations due to the rotation of the each segment for the twodimensional finite element modeled cable, i yi i i I T T dt d q q q & = && ¶ ¶  ÷ ÷ø ö ç çè æ ¶ ¶ (i = 1,2,K, n + 1) (2.48) The combination of the rotational and translational components of the equations of motion and the addition of the ATV segment result in a set of equations of the form shown in equation 2.30. The inclusion of the ATV segment into the system adds an extra degree of freedom and thus an extra equation to the set of the equations of motion of the system, resulting in the coefficient matrices increasing to size n+1 by n+1 and the column matrices to size n+1 by 1. The additional column in the coefficient matrices represents the effect of the motion of the ATV segment on the dynamics of the cable segments. The 43 additional row of the coefficient matrices represents the equation of motion of the ATV segment, which contains the effects of the dynamics of each cable segment. The changes to the set of equations derived for the translation of the nodes of the segments are in the terms of the A coefficient matrix. The corrected terms of the A matrix are as follows, i j y i j i j i j A I l l C m , i , (q j qi ) , d  = + (2.49) where, the kronecker delta, di,j, is zero at all time except when i equals j at which time its value is unity. Thus the inertia term only appears on the diagonal of the A matrix, and since the segment nodes, in accordance with the chosen segment model, have a zero value for inertia the only change in A is at n+1, n+1. 2.3.3 Generalized Forces The generalized force equation derived in equation 2.28 applies to all segments including the newly added ATV segment. However, the crossflow method does not apply to aircraft. As a result, a method for modeling the aerodynamic forces and moments on the aircraft are necessary. As stated in chapter 1, there is no shortage of research on the aerodynamics of general aircraft. The technique to be used in this paper is the linearized smalldisturbance theory, Nelson14. In short, linearized smalldisturbance theory simplifies the motion of the aircraft to perturbations of the aircraft about some fixed reference flight condition, such as typical cruising conditions. Dealing with perturbations allows for the elimination of higher order terms in the aircraft equations of motion, thus linearizing them. The aerodynamics of the 44 aircraft in the equations of motion are represented by stability derivatives. Stability derivatives approximate the change in aerodynamic forces or moments due to the deviations from the reference flight condition. These stability derivatives are based upon the size, shape, and mass distribution of the aircraft, and apply to a wide range of aircraft. The method does breakdown when the aircraft enters a flight regime in which the nonlinearities are important, such as stalling and spinning. The twodimensional equations for the aerodynamic forces, X F¢ and Z F¢ , and moment, M, in the aircraft fixed frame are as follows, ( ) ( ) M M I (M u M w M w M q) F Z m Z u Z w Z w Z q F X m X u X w Y u w w q Z u w w q X u w = + + + + ¢ = + + + + ¢ = + + & & & & 0 0 0 (2.50) where, m and IY are aircraft mass and longitudinal rotational inertia. Also, u and w are the xdirection and zdirection perturbation velocities and q is the pitch rate in the ATV fixed frame. Xu, Xw,…, Mq are the stability derivatives. (Note that notation in equation 2.50 is commonly used when discussing aircraft dynamics. The pitch rate q is not to be confused with the generalized coordinate q which is itself commonly used in the notation of Lagrange’s equations.) In the aerial tow system the reference flight condition for the ATV is that of the host aircraft, and the host aircraft reference flight conditions can simply be viewed as the freestream conditions, Uo and Wo. As a result, the perturbation of the ATV about the reference conditions is simply a result of the attitude of the ATV and the motion due to segment rotation. As we have already derived the velocity of a node due to rotation in 45 the inertial frame, a coordinate transformation from the inertial frame to the ATV aircraft frame can be utilized to develop the equation for perturbation velocity, úû ù êë é  úû ù êë é + + úû ù êë é úû ù êë é = úû ù êë é + + + 1 0 1 1 W U W Z U X C C w u o o n o n T T n & & g q (2.51) With the perturbation parameters in the ATV aircraft frame the aerodynamic forces and moments on the aircraft can easily be attained. To calculate the generalized forces produced by the aerodynamics of the ATV the forces must be in the inertial frame. This requires the following rotational coordinate transformation, úû ù êë é ¢ ¢ úû ù êë é úû ù êë é = úû ù êë é + Z X Z X F F C C F F qn 1 g (2.52) No rotational transfo rmation is necessary for the twodimensional applied moment. Thus, the aerodynamic moment, M, is directly translated to the ATV segment moment Mq. The forces and moments are applied to the generalized force equation, equation 2.28, following the translation. The dynamics as well as the applied aerodynamic forces and moments of the cable ATV system have been modeled. As a result, it is possible to simulate the effect of the cable and ATV system due to various initial conditions. However, to fully analyze the aerial towed system one last component is necessary, the host aircraft. 2.4 Host Vehicle Up to this point, the model of the system has the cable, immersed in a flowing fluid, connected to a fixed point in space. In reality, the cable will be connected to the host 46 aircraft, and while the ideal case has the host vehicle flying at a constant speed and altitude it is likely that, for one reason or another, the position of the host vehicle will be disturbed. The host vehicle is an airplane, and as a result its dynamics are dictated by the aerodynamic equations presented in equation 2.50, and as it is the towing vehicle it is affected by the motion of the cable and ATV segments. On the other hand, the motion of the cable and ATV will also be affected by the host vehicle motion. In general, the host vehicle is much larger than the ATV and as a result its aerodynamic loads will overpower the loads produced by the motion of the cable and ATV. Thus, a simplifying assumption can be made that the host vehicle is unaffected by the dynamics of the cable and ATV, yet the cable and ATV are affected by the motion of the host vehicle. This assumption allows the host vehicle to be modeled as a simple point. The motion of the host aircraft is independent of the cable and the ATV dynamics and consequently the movement of this point can be chosen arbitrarily. For instance, the equations of motion for longitudinal oscillation in the zdirection of the host aircraft frame can be written simply as sin( ) cos( ) sin( ) z A 2 t z A t z A t o o o o o o o o w w w w w =  = = & & (2.53) where, z is the host vehicle vertical displacement, t is time, Ao is amplitude, and wo is frequency. The amplitude and frequency can be chosen to represent the longitudinal dynamics of a wide range of host aircraft. The effect of host aircraft oscillation on the vertical position of the ATV is of critical importance when the ATV is very near the 47 ground; however the aircraft is not limited to simply vertical motion. The equations of motion can be derived for horizontal motion, vertical motion, or any combination of the two. The technique for modeling the cable and ATV dynamics with the introduction of the varying host vehicle position remains more or less the same. In the previous models the cable connecting point always coincided with the origin of the inertial coordinate system. However, with the introduction of the movable host vehicle this is not so. Any movement of the host aircraft results in the translation of the cable connection point. The derivation of Lagrange’s equations for the cable and ATV hinges on the location and velocity of the segment nodes. Thus, to account for host vehicle movement the inertial position equation must be altered to represent the translation of the cable connection point. The resulting equation is as follows, å= úû ù êë é úû ù êë é + úû ù êë é = úû ù êë é i i j j i l C z x Z X j 1 0 q (2.54) where, x is the difference in the xdirection of the host aircraft in the inertial frame. Due to the fact that the position of the host aircraft is a function of time the derivation of Lagrange’s equations will produce a different outcome from that derived previously, the difference being the effect of the host vehicle oscillation. The derivation of Lagrange’s equations begins by taking the time rate of change of position, 48 å å = = = +  = + j i j i i j i j i i i i Z z l S X x lC 1 1 q q q q & & & & & & (2.55) The square of velocity with the incorporation of a moving host vehicle is as follows, ( ) ( ) ( ) i i i j z x x lC z l S l l C V X Z j i j k j k i i k i i i i i k k k q q q q q q q q  = = = + +å  +åå = + = & & & & & & & & & & 1 1 1 2 2 2 2 2 (2.56) Thus, the total kinetic energy for the aerial tow system is as follows, å å( ) åå ( ) å + =  = = = + = ÷ ÷ ø ö ç ç è æ + +  + = = 1 1 1 1 1 2 2 1 1 2 2 1 2 1 n k j i j k j k i i k i k i i i i n k k k i i i j m z x x l C z l S l l C T m V q q q q & & &q& &q& q&q& (2.57) Note that the last term in equation 2.57 is the same as that derived in the cable dynamics section, equation 2.12. Thus, the derivation of Lagrange equations for the last term will yield the same results. As a result, we can focus on the host vehicle oscillation terms when deriving the Lagrange equations and simply superimpose the results onto the equations of motion for the cable and ATV. The kinetic energy due to the host vehicle oscillation, Thost, is as follows, å å( ) + = = ÷ø ö çè æ = + +  1 1 1 2 2 2 2 2 1 n k k i osc k i i i i i i T m z x x l C z l Sq q & & &q& &q& (2.58) 49 Note that the first two terms in the equation 2.58 are only a function of time. Therefore, they will disappear during the differentiation of the kinetic energy with respect to the generalized coordinates. The derivative of kinetic energy with respect to the ith generalized coordinate and the time rate of change of the ith generalized coordinate yields the following, ( ) å+ = =  + ¶ ¶ n 1 k i k i i i i i osc i i m x l S z l C T q q q q q & & & & (2.59) ( ) å+ = =  ¶ ¶ n 1 k i k i i i osc i i m xl C zl S T q q q & & & (2.60) Finally, the time rate of change of equation 2.60 produces, ( ) å+ =    = ÷ ÷ø ö ç çè æ ¶ ¶ n 1 k i k i i i i i i i osc i i i i m xl C zl S x l S z l C T dt d q q q q q q q & && & & & & & (2.61) Combining equations 2.59 and 2.61, in accordance with Lagrange’s equations, results in the following, ( ) å+ = =  ¶ ¶  ÷ ÷ø ö ç çè æ ¶ ¶ n 1 k i k i i i osc i osc i i m xl C zl S T T dt d q q q q & && & (i = 1,2,K, n + 1) (2.62) Note that equation 2.62 has units of force times distance just as the generalized forces. As a result, the equation 2.62 can be thought of as generalized forces due to host vehicle acceleration, as shown below, 50 ( ) å+ = =  n 1 k i osci k i i i i Q m &x&l Cq &zl Sq (i = 1,2,K, n + 1) (2.63) The superposition of equation 2.63 into the generalized force equation 2.28 yields, ( ) ( ( )) ( 1,2, , 1) 1 = +  +   + = å+ = i n Q M F m x l C F m g z l S n j i i i X j j i i Z j j i i K && & q q q (2.64) Thus in the matrix form of the set of equations of motion, the host vehicle oscillation has no affect on the A or B coefficient matrix. The results of the derivation show that the moving host frame produces applied loads on each segment node proportional to the host vehicle acceleration. Note that the derivation due to the host vehicle motion is independent of the functions used to model it, allowing for any harmonic model of the host aircraft dynamics. 2.5 System Equations of Motion of an Alternate Segment Model At this point, the equations of motion for the twodimensional aerial tow system are complete, and given a proper numerical integration routines and a set of valid initial conditions a simulation of the twodimensional aerial tow system can be accomplished. However, the lumped mass model used for the finite element segments is quite simplistic. As discussed earlier, this model places the mass of a segment at the node at the base of the segment and all applied forces on the segment are located at the node. This simplification is particularly poor when there are few segment s that make up the model of the cable. 51 Imagine the cable is modeled as one lumped mass segment. The segment is in the presence of moving fluid, as a result there is a distributed friction load on the segment. According to the lumped mass assumption the inertial force of the segment and the accumulated drag force are applied at the base of the segment. The lumped mass segment equation of motion is as follows, ( ) ml F C F mg S q X q Z q  + && = (2.65) It would seem a much better model would have the mass and applied force applied at the midpoint of the segment. The resulting angular acceleration is twice that of the lumped mass model. ( ) ml F C F mg S q X q Z q  + && = 2 (2.66) An even better model would have the segment modeled as a thin rod or cylinder. This would have the inertial and applied loads at the center of mass of the rod as in the previous model. However, it also has its rotational inertia that resists rotation (the rotational inertia about the center of gravity of a thin rod is Iy = (ml2)/12). The resulting angular acceleration is an average of the two previous equations. ( ) ml F C F mg S q X q Z q  + = 2 && 3 (2.67) The one segment example illustrates the weakness of the lumped mass model. Presumably as the number of segment increase the difference in the dynamics of the 52 lumped mass and the thin rod models will decrease (this will be shown in chapter five). However, the as the number of segments increase so too does the computational workload. The question remains how do the equations of motion derived throughout this chapter for the lumped mass model change for a system modeled by thin rod segments. The answer can be found by rederiving Lagrange’s equations with a new position equation that represents the center of mass of each segment. Upon completion one will realized that equation 2.30 remains true, recall, [A][q&]+ [B][q&2 ]= [Q] The only changes in the equations of motion for the thin rod modeled segments are in the mass distribution in the A and B coefficient matrices and the lever arm of the generalized force equations. The mass distribution of the coefficient matrices in the lumped mass model is given in equation 2.33. The mass distribution of the coefficient matrices in the thin rod modeled system for n+1 segments is as follows, ( ) å+ = + ÷ ÷ø ö ç çè æ =  + 1 max , , , max( , ) 2 4 1 n k i j k i j i j i j m m m d (2.68) The generalized force equation remains much the same however since the forces on the ith segment are applied at the center of mass and not at the base of the segment the lever arm is halved when i is equal to j, 53 [( ) ( ( )) ] ( 1,2, , 1) 2 1 1 , = + ÷ ÷ø ö ç çè æ   +   + = å+ = i n Q M F m x l C F m g z l S n j i i j i X j j i i Z j j i i K && & d q q (2.69) Note the for equation 2.69 to be accurate for the ATV segment the ATV center of mass must coincide with the center of mass of the segment model. In other words, the length of the ATV segment will be twice the distance between the cableATV hitch point and the ATV center of mass. The set of equations of motion for the aerial towed system with segments modeled as thin rods is attained with the implementation of equations 2.32, 2.49, 2.68 and 2.69 on equation 2.30. Given a proper numerical integration routines and a set of valid initial conditions a simulation of the twodimensional aerial tow system can be accomplished. The difference in the results of the lumped mass and thin rod segment modeled system simulations with the same initial conditions will be analyzed in chapter five. It is only then a conclusion can be drawn about the pros and cons of the two methods. 54 CHAPTER 3 3THREEDIMENSIONAL DYNAMICS 3.1 Introduction Modeling the aerial tow system for the twodimensional case is moreorless straight forward. The equations of motion have been thoroughly derived, include nonlinearities, and make few assumptions. As a result, we can expect that the simulation of the system will be thorough as well; the system could be used to complete a wide variety of case studies, such as, ATV position due to wind gust or host vehicle motion, as well as ATV tracking ability under various conditions. However, it is still limited the phenomena that occur only in the twodimensional plane. For a more robust simulation it is necessary to derive the equations of motion in three dimensions. Being that the actual system exists in three dimensions, this will more accurately represent the real world dynamic occurrences. The addition of the third dimension will drastically increase the complexity of the mathematical models. The methodology for the derivation of the threedimensional equations of motion is the same as that used for the twodimensional case (chapter 2) although a few new concepts are necessary to carryout that methodology. 3.2 Cable Dynamics The threedimensional model of the cable is attained via the application of the finite element technique. As in the twodimensional case, the continuous cable is broken into an arbitrary number of segments. The lumped mass model of the cable segments will be 55 used in the initial analysis of the cable. Each segment is free to rotate about its neighbors, giving the series of rigid connections flexibility as a whole. The added dimension allows the segments two additional rotations, the rotation about the xaxis and the ability for each segment to spin (rotation about the zdirection). The equations of motion for the segment motion will include the new dynamic freedom, and as a result, the equations will increase in complexity. To focus on the dynamics of the cable we will assume that the cable is attached to the origin of the inertial frame. In other words, assume that the host vehicle is steady. Thus, the focus of the dynamics will be on the threedimensional motion of the cable. The effects of the motion of the cable connection point will be addressed in a later section. 3.2.1 Coordinate Systems The addition of the third dimension does not have an effect on the number of coordinate systems needed, although it does affect the type and complexity of the orientation of the coordinate systems. The inertial frame is orthogonal and fixed. Each cable segment has an orthogonal coordinate system which has its origin such that the location of the node of that segment is the length of the segment in the positive zdirection. For instance, imagine there is one cable segment that rotates about the inertial frame. The origin of the cable segment frame is at the point of rotation of the segment. Thus, the location of the origin of the first cable segment frame and the origin of the inertial frame are the same; the difference is in the orientation of the frames, figure 3.1. 56 Figure 3.1 Threedimensional inertial and cable segment coordinate systems The orientation of the cable segment relative to the inertial frame is due to the angle of rotation of the cable frame. The cable frame can rotate about three separate axes (x, y, and zaxes) unlike the twodimensional case in which it can rotate about only one (the yaxis). The common method for describing the threedimensional rotation requires the use of Euler angles. These Euler angles are the heading angle y, the attitude angle q, and the bank angle f. Any orientation of the cable frame relative to the inertial frame can be described by using these three angles. However, it is important to note that the order of rotation is important. Imagine that a cable segment frame is parallel to the inertial frame and their origins are coincident, figure 3.2. Allow the cable frame to rotate about its zaxis by an angle y. The resulting coordinate system is the double primed system. Thus, the coordinate transformation from the doubleprimed system to the fixed inertial frame is as follows, 57 ú ú ú û ù ê ê ê ë é ¢ ¢ ¢ úû ù êë é = ú ú ú û ù ê ê ê ë é ¢ ¢ ¢ ú ú ú û ù ê ê ê ë é =  ú ú ú û ù ê ê ê ë é ¢ ¢ ¢ = ú ú ú û ù ê ê ê ë é z y x C z y x S C C S z y x Z Y X y y y y y 0 0 1 0 0 (3.1) Next, allow the doubleprimed system to rotate about the y ¢ axis by an angle q. The resulting coordinate system is the tripleprimed system. The coordinate transformation from the tripleprimed system to the double primed system is as follows, ú ú ú û ù ê ê ê ë é ¢¢ ¢¢ ¢¢ úû ù êë é = ú ú ú û ù ê ê ê ë é ¢¢ ¢¢ ¢¢ ú ú ú û ù ê ê ê ë é  = ú ú ú û ù ê ê ê ë é ¢ ¢ ¢ z y x C z y x S C C S z y x q q q q q 0 0 1 0 0 (3.2) Finally, allow the tripleprimed system to rotate about the x¢¢ axis by an angle f. The resulting coordinate system is the unprimed cable segment frame. The coordinate transformation from the cable segment frame to the tripleprimed system is as follows, ú ú ú û ù ê ê ê ë é úû ù êë é = ú ú ú û ù ê ê ê ë é ú ú ú û ù ê ê ê ë é  = ú ú ú û ù ê ê ê ë é ¢ ¢ ¢ z y x C z y x S C C S z y x f f f f f 0 0 1 0 0 (3.3) The rotational coordinate transformation from the cable segment frame to the inertial frame is attained by combining equations 3.1, 3.2, and 3.3, ú ú ú û ù ê ê ê ë é ú ú ú û ù ê ê ê ë é  +  +  + + = ú ú ú û ù ê ê ê ë é úû ù êë é úû ù êë é úû ù êë = é ú ú ú û ù ê ê ê ë é z y x S C S C C S C C C S S S C S S S C C C S C C S S S S C S C z y x C C C Z Y X q q f q f y q y f y q f y f y q f y q y f y q f y f y q f y q f (3.4) 58 Figure 3.2 Euler angle rotations The series of three rotations of the cable frame are sufficient to attain any orientation. The same is true independent of the order of rotation. However, the order of rotation does have an effect on the final position of the cable frame. The order used throughout this paper is y, q, and f for every threedimensional segment rotation. Let us assume that the Euler angles are limited to the following rotations, f p p q p y p , 0 2 2 2 0 £ < 2 ,  £ £ £ < (3.5) When q = ± p 2 there is no distinct set of solutions for y and f , this is known as gimbal lock. To avoid complications associated with gimbal lock it is common to switch to an alternate system. In this particular aerial tow system the issue of gimbal lock will not arise as the system necessitates the ATV fly at a lower altitude than the host aircraft resulting in cable segment attitude angles between± p 2. 59 The threedimensional rotational coordinate transformation from the cable to the inertial coordinate system is represented by the following directional cosine matrix, ú ú ú û ù ê ê ê ë é  +  +  + + = úû ù êë é q q f q f y q y f y q f y f y q f y q y f y q f y f y q f f q y S C S C C S C C C S S S C S S S C C C S C C S S S S C S C C , , (3.6) The position of the ith segment’s node in its coordinate system is by definition the length of the segment in the zdirection. Since the segments are a series of connected bodies the directional cosine matrix can be used to formulate the node’s position in the inertial frame, å= ú ú ú û ù ê ê ê ë é úû ù êë = é ú ú ú û ù ê ê ê ë é i j i j i i l C Z Y X j j j 1 , , 0 0 f q y (3.7) Expansion of equation 3.9 yields a more explicit form of the ith segment position in the inertial frame, ( ) ( ) å ( ) å å = = = = =  + = + i j i j i j i j i j i j j j j j j j j j j j j j Z l C C Y l C S S S C X l S S C S C 1 1 1 q f y f y q f y f y q f (i = 1,2,K, n) (3.8) As in the towdimensional case, the location of the node of a segment is dependant upon the angles of rotation of that segment and each preceding segment. 60 3.2.2 Lagrange Equations The modeling of the cable is accomplished using the finite element method, just as with the twodimensional case. As discussed previously, the method breaks the cable into an arbitrary number of segments, and each cable segment has its own coordinate system. The accuracy of the cable model improves with an increase in the number of segments for a given cable, so too does the difficulty of the derivation of the equations of motion for the cable. Using Newton’s method to derive the equations of motion for such a complex configuration would be fraught with danger, as it would require intense vector bookkeeping and multiple coupled rotation equations for each segment. The advantages of the Lagrange technique over Newton’s method are even stronger for threedimensional analysis than in the twodimensional case. The added dimension greatly increases the system complexity and since the Lagrange technique has a standard process for the derivation of the equations of motion it is less likely that mistakes will be made in the derivation of the equations. The bounds of the Lagrange equation are set by the degree of freedom of the system (refer to section 2.2.2). In three dimensions, the complete description of a cable segment is composed of its nodal position and its orientation (i.e. angle of twist). In the Cartesian coordinate system the position of the node is defined by its x, y, and zcoordinates. Due to the fact that the segment length is fixed the coordinates are not independent. The Cartesian constraint equations are shown below, 61 2 1 1 2 1 1 2 1 1 2 ÷ ÷ø ö ç çè æ  + ÷ ÷ø ö ç çè æ  + ÷ ÷ø ö ç çè æ = å å å =  =  = i j i j i j i j i j i i j l X X Y Y Z Z (i = 1,2,K, n) (3.9) The Cartesian coordinates alone cannot specify the angle of twist of the segment. Therefore, there are in actuality four coordinates necessary to define the system, and being that there is one constraint equation each cable segment has three degrees of freedom. It is possible to completely define the position and orientation of a segment using the Euler angles. In this case, there are three angles which completely define the segment and there are no constraint equations. Therefore, the Euler angles of each segment are the independent generalized coordinates of the system. For the twodimensional case the independent generalized coordinates were the attitude angles of rotation and since there was only one angle per segment the generalized coordinate consisted of just one variable. However, for the threedimensional model each segment requires three angles for complete description of orientation. As a result, for a cable modeled by n segments there are 3n generalized coordinates made up of y, q, and f for each segment. (This assumes that each cable segment is free to fully rotate about the preceding segment. The validity of this assumption will be discussed in a later section). For the threedimensional system, Lagrange’s equations represent an equation for each of the 3n degrees of freedom, i i i Q q T q T dt d = ¶ ¶  ÷ ÷ø ö ç çè æ ¶ ¶ & (i = 1, 2,K,3n) (3.9) 62 ( ) úû ù ¶ ¶ + ¶ ¶ + ¶ ¶ êë é + ¶ ¶ + + ¶ ¶ + ¶ ¶ =å= i j j i j j i j j n j i j Z j j i j Y j i j i X j q M q M q M q Z F m g q Y F q X Q F f q y f q y 1 (i = 1, 2,K,3n) (3.10) Note that the generalized force, Q, contains applied forces, applied moments, and inertial forces as per the discussion in section 2.2.2. Equation 3.9 represents a set of 3n equations. In which the generalized coordinate is made up of three independent variables. It may be presented more clearly rewritten as three sets of n equations, in which the generalized coordinate for each is an independent variable, i i i Q T T dt d Q T T dt d Q T T dt d i i i i i i y q f y y q q f f = ¶ ¶  ÷ ÷ø ö ç çè æ ¶ ¶ = ¶ ¶  ÷ ÷ø ö ç çè æ ¶ ¶ = ¶ ¶  ÷ ÷ø ö ç çè æ ¶ ¶ & & & (i = 1, 2,K, n) (3.11) ( ) ( ) å ( ) å å = = = ¶ ¶ + ¶ ¶ + + ¶ ¶ + ¶ ¶ = ¶ ¶ + ¶ ¶ + + ¶ ¶ + ¶ ¶ = ¶ ¶ + ¶ ¶ + + ¶ ¶ + ¶ ¶ = n j i j j i j Z j j i j Y j i j i X j n j i j j i j Z j j i j Y j i j i X j n j i j j i j Z j j i j Y j i j i X j M Z F m g Y F X Q F M Z F m g Y F X Q F M Z F m g Y F X Q F 1 1 1 y y y y y q q q q q f f f f f y y q q f f (i = 1, 2,K, n) (3.12) 63 Thus, each segment will contain three equations of motion. Note that the applied moments and the moments produced by the applied forces are about the axis of rotation of their respective generalized coordinate, i f& , i q& , and i y& . The method for deriving the equations of motion is the same as that carried out in the twodimensional cable dynamic section. The only difference is the more complex kinetic energy equation and the increase in the number of generalized coordinates of the system. The kinetic energy is made up of rotational and translational kinetic energy. The nodes, being point masses, have no inertia and thus no rotational kinetic energy. As a result, the derivation of the equations of motion relies solely on the translation kinetic energy of the system as a function of the generalized coordinates. The kinetic energy can be found by calculating the square of velocity for the node of each segment as a function of the generalized coordinates, equation 2.9. The velocity components are simply the time rate of change of position [ ( ) ( ) ( )] [ ( ) ( ) ( )] [ ( ) ( )] for (i n) Z l C S S C Y l C C S S S S C C S S C S C X l S C C S S C C C C S S S C i j i j j j i j i j j j j i j i j j j j j j j j j j j j j j j j j j j j j j j j j j j j j j j j j j 1, 2 , 1 1 1 K & & & & & & & & & & & = =  +  =   + + + =  + +  å å å = = = q f q f y f y q f y q f y f y q f y f y q f y q f y f y q f f q f q y f q y (3.13) Thus, the square of velocity for the kth segment is as follows, 64 ( ) { [ ( ( ) ( ) ) ( ) ( ( ) ( )] [ ( ) ( ) ( )] [ ( ( ) ( ) ( ) ( ( ) ( )]} j i i j i j i j j i i j i i j j j i i j j j i i j i j j i i j i i j j j i i j i j i j j i i j i j i i j i j j i i j i j i j i j j i i j i j j i i j i i j i j j i i j i i j j j i i j i j i j j i i j i j j i i j j i j i j j i i j i j i j j i i j j i j i i j i j C C C S S S S S C S S S C S C S C C S C C S C C C C S S S S S C S S S S C C C S C S C C C S C C C C C C C S S C C S C S S C C C C S S C C S S S C C S S S C S S C S C S C S C S C C C S C C S C S S S S C C S C S S S C S S S C C V X Y Z j j i j j j i j j j i j k i k j k k k k y y f f q q f f y y f f q f f q y y f f q y y f f q q y y f f q f f q y y f f q q f f y y f f q y y f f q q y y f f q q f f q q y y f f q q y y f f q f f q q y y f f q f f q y y f f q q f f y y f f q q y y f f q f f q q y y f f q q f f y y f f q f f q f f q q y q y f y q q f y q f f + +  + +   +  + +  + + +  + + +  + +  +   + + + +  + + = + + =                  = = åå & & & & & & & & & & & & & & & 1 1 2 2 2 2 (3.14) Furthermore, the threedimensional kinetic energy equation is represented as a function of the independent generalized coordinates below, { [ ( ( ) ( ) ) ( ) ( ( ) ( )] [ ( ) ( ) ( )] [ ( ( ) ( ) ( ) ( ( ) ( )]} j i i j i j i j j i i j i i j j j i i j j j i i j i j j i i j i i j j j i i j i j i j j i i j i j i i j i j j i i j i j i j i j j i i j i j j i i j i i j i j j i i j i i j j j i i j i j i j j i i j i j j i i j j i j i j j i i j i j i j j i i j j i j i i j i j C C C S S S S S C S S S C S C S C C S C C S C C C C S S S S S C S S S S C C C S C S C C C S C C C C C C C S S C C S C S S C C C C S S C C S S S C C S S S C S S C S C S C S C S C C C S C C S C S S S S C C S C S S S C S S S C C T m V m j j i j j j i j j j i j n k k i k j k n k k k y y f f q q f f y y f f q f f q y y f f q y y f f q q y y f f q f f q y y f f q q f f y y f f q y y f f q q y y f f q q f f q q y y f f q q y y f f q f f q q y y f f q f f q y y f f q q f f y y f f q q y y f f q f f q q y y f f q q f f y y f f q f f q f f q q y q y f y q q f y q f f + +  + +   +  + +  + + +  + + +  + +  +   + + + +  + + = =                  = = = = å ååå & & & & & & & & & & & & 1 1 1 1 2 2 1 2 1 (3.15) Equation 3.15 is the basis from which Lagrange’s equations will be derived. All one needs to do is compare equation 3.15 and 2.12 to grasp the increase in complexity from a 65 system existing in two to threedimensions. With the cable kinetic energy equation derived the derivation of Lagrange’s equations for each set of generalized coordinates can commence. 3.2.2.1 Bank angle The implementation of the Lagrange equation in the threedimensional cable model must be accomplished for each angle necessary for cable orientation definition. For the bank angle generalized coordinate the derivation is as follows. The derivative of kinetic energy with respect to the bank angle and the time rate of change of the bank angle of the ith cable segment yield the following, { [ ( ( ) ( ) ) ( ) ( ( ) ( )] [ ( ) ( ) ( )] [ ( ( ) ( ) ( ) ( ( ) ( )]} j i i j i j i j j i i j i i j j j i i j j j i i j i j j i i j i i j j j i i j i j i j j i i j i j i i j i j j i i j i j i j i j j i i j i j j i i j i i j i j j i i j i i j j j i i j i j i j j i i j i j j i i j j i j i j j i i j i j i j j i i j j i j i i j i j C S C S S C S S S S S C C S C C C C S S C S C C S C S C S S S S S S S C C C S S C S S C C S C S C C C S C S S C S S C S S S C C S S S C C C S S S C S S C C S S S S C C C S C S S C C C C C S C C S S S S C S S S S C C S C S C C m T j j i j j j i j j j i j n k k j k i y y f f q q f f y y f f q f f q y y f f q y y f f q q y y f f q f f q y y f f q q f f y y f f q y y f f q q y y f f q q f f q q y y f f q q y y f f q f f q q y y f f q f f q y y f f q q f f y y f f q q y y f f q f f q q y y f f q q f f y y f f q f f q f f q q y q y f y q q f y q f f f  + +    +  + + +  + +   +   +   + + +  + + +  +   + + = ¶ ¶                  = = åå & & & & & & & & & & & & 1 1 (3.16) 66 { ( ( ) ( ) ) ( ) ( ( ) ( )} j i i j i i j j j i i j i j i j j i i j i j j i i j j i j i j j i i j i j i j j i i j j i j i i j i j C S S S C C S S S C S S C S C S C S C S C C C S C C S C S S S S C C S C S S S C S S S C C m T j j j n k k j k i y y f f q f f q y y f f q q f f y y f f q q y y f f q f f q q y y f f q q f f y y f f q f f q f f q q y q f f  + +    + + + +  + + = ¶ ¶       = = åå & & & & 1 1 (3.17) Finally evaluating the time rate of change of equation 3.17 produces, { ( ( ) ( ) ) ( ) ( ( ) ( ) ( ) ( ) ( ( ) ( ) ( ( ) ( ) ) ( ) ( ( ) ( ) [ ( ( ) ( ) ) ( ) ( ( ) ( )] [ ( ) ( ) ( )] [ ( ( ) ( ) ( ) ( ( ) ( )]} j i i j i j i j j i i j i i j j j i i j j j i i j i j j i i j i i j j j i i j i j i j j i i j i j i i j i j j i i j i j i j i j j i i j i j j i i j i i j i j j i i j i i j j j i i j i j i j j i i j i j j i i j j i j i j j i i j i j i j j i i j j i j i i j i j j i i j i i j j j i i j i j i j j j i i j i j j i i j j i j i j j j i i j i j i j j i i j j i j i i j i j j i i j i i j j j i i j i j i j j i i j i j j i i j j j i i j i j j i i j j i j i j j i i j i i j j j i i j i j i j j i i j i j j i i j j i j i j j i i j i j i j j i i j j i j i i j i j C S C S S C S S S S S C C S C C C C S S C S C C S C S C S S S S S S S C C C S S C S S C C S C S C C C S C S S C S S C S S S C C S S S C C C S S S C S S C C S S S S C C C S C S S C C C C C S C C S S S S C S S S S C C S C S C C S S S S C C S C S C S S C S C S C S S S C C S S C C C C S C S S C S S C C S S S S S C C C C S C S C S S S S S S S C C S S C S C C C C C C S S S C S C S C S S C S C S S S C C S S S C S S C S C S C S C S C C C S C C S C S S S S C C S C S S S C S S S C C m T dt d j j i j j j i j j j i j j j j j j j j j j j n k k j k i y y f f q q f f y y f f q f f q y y f f q y y f f q q y y f f q f f q y y f f q q f f y y f f q y y f f q q y y f f q q f f q q y y f f q q y y f f q f f q q y y f f q f f q y y f f q q f f y y f f q q y y f f q f f q q y y f f q q f f y y f f q f f q f f q q y y f f q f f q y y f f q q f f y y f f q q y y f f q f f q q y y f f q q f f y y f f q f f q f f q q y y f f q f f q y y f f q q f f y y f f q q y y f f q y y f f q q y y f f q f f q q y y f f q f f q y y f f q q f f y y f f q q y y f f q f f q q y y f f q q f f y y f f q f f q f f q q y q y f y q q f y q f f y q f y f q y f q y q f f  + +    +  + + +  + +   +   +   + + +  + + +  +   + + + +  + + + +  + + + +  + +   +  + +  +  + +  +   + + + +  + + = ÷ ÷ø ö ç çè æ ¶ ¶                                    = = åå & & & & & & & & & & & & & & & & & & & & & & && & & 2 2 1 1 2 2 2 (3.18) 67 Note that, just as in the twodimensional case, the term representing the derivative of kinetic energy with respect to generalized coordinate, equation 3.16, is completely canceled by the some of the terms in the time rate of change of the derivative of kinetic energy with respect to the time rate of change of the generalized coordinate, equation 3.18. Combining equations 3.16 and 3.18 results in the series representation of the left hand side of the Lagrange’s equations for the bank angle of the threedimensional finite element modeled cable, { ( ( ) ( ) ) ( ) ( ( ) ( ) ( ) ( ) ( ( ) ( ) ( ( ) ( ) ) ( ) ( ( ) ( )} for (i n) S S S S C C S C S C S S C S Q C S C S S S C C S S C C C C S C S S C S S C C S S S S S C C C C S C S C S S S S S S S C C S S C S C C C C C C S S S C S C S C S S C S C S S S C C S S S C S S C S C S C S C S C C C S C C S C S S S S C C S C S S S C S S S C C m T T dt d j i i j i i j j j i i j i j i j i j j i i j i j j i i j j i j i j j j i i j i j i j j i i j j i j i i j i j j i i j i i j j j i i j i j i j j i i j i j j i i j j j i i j i j j i i j j i j i j j i i j i i j j j i i j i j i j j i i j i j j i i j j i j i j j i i j i j i j j i i j j i j i i j i j j j j j j j j j j j n k k j k i i 1,2, , 2 2 2 2 2 1 1 K & & & & & & & & & & && & & = + +  = + + +  + + + +  + +   +  + +  +  + +  +   + + + +  + + = ¶ ¶  ÷ ÷ø ö ç çè æ ¶ ¶                   = = åå y y f f q f f q y y f f q q f f f y y f f q q y y f f q f f q q y y f f q q f f y y f f q f f q f f q q y y f f q f f q y y f f q q f f y y f f q q y y f f q y y f f q q y y f f q f f q q y y f f q f f q y y f f q q f f y y f f q q y y f f q f f q q y y f f q q f f y y f f q f f q f f q q y q f y f q y f q y q f f f (3.19) To fully develop the equations of motion the generalized force equation as a function of bank angle is evaluated in accordance with equation 3.12. 68 [ ( ) ( ) (F m g)l ( S C )] for (i n) Q M F l C S S S C F l C C S S S i i i i i i i i i i i i Z j j i n j i i i X j i Y j i +  =1,2,K, +   +  + = å= f q f f f y f q y f y f q y (3.20) The implementation of the Lagrange equations for the bank angle generalized coordinate is complete, however to attain the full system equations of motion Lagrange’s equations must be derived for the two remaining generalized coordinates. 3.2.2.2 Attitude Angle This section investigates the application of the Lagrange equations for the attitude angle generalized coordinate. The Lagrange equation term representing the derivative of kinetic energy with respect to the attitude angle and the rate of change of the attitude angle of the ith cable segment is as follows, { [ ( ) ( ) ( )] [ ( ) ( ) ( )] [ ( ) ( ) ( )]} j i i j i j j i i j i j i i j i j j i i j i j i i j i j j i i j i j i i j i j j i i j i j i j i j j i i j i j j i i j i i j i j j i i j i j i i j i j j i i j i j i j i j j i i j i j j i i j i i j i j C C C C S S C S C S C C C C C C C C S C S C S C C S S S C C S S C C C S C C C C S C C S S S S C C S C S C C C S S C S S C C S C S C C C S C S S C S S C S S S C C S S S C m T j j i j j j i j j j i j n k k j k i y y f f q q y y f f q y y f f q q y y f f q y y f f q q y y f f q y y f f q q y y f f q q f f q q y y f f q q y y f f q f f q q y y f f q y y f f q q y y f f q q f f q q y y f f q q y y f f q f f q q y q y f y q q f y q f f q                = = + +   +  + +  + +  + +  + +   +   + = ¶ ¶ åå & & & & & & & & & & & & 1 1 (3.21) 69 { ( ) ( ) ( j i i j i j i i j i j )} j i i j i j i j i j j i i j i j j i i j i i j i j C C S C S C C C S C C C C C C C S S C C S C S S C C C C S S C m T j j j n k k j k i y y f f q y y f f q q y y f f q q f f q q y y f f q q y y f f q f f q q y q f q      = =  + +  + + + = ¶ ¶ åå & & & & 1 1 (3.22) Taking the time rate of change of the equation 3.22 produces the first term in Lagrange’s equations for the attitude angle generalized coordinate, 70 { ( ) ( ) ( ) ( ) ( ) ( ) ( ) ( ) ( ) [ ( ) ( ) ( )] [ ( ) ( ) ( )] [ ( ) ( ) ( j i i j i j j i i j i )]} j i i j i j j i i j i j i i j i j j i i j i j i i j i j j i i j i j i j i j j i i j i j j i i j i i j i j j i i j i j i i j i j j i i j i j i j i j j i i j i j j i i j i i j i j j j i i j i j i i j i j j j i i j i j i j i j j j i i j i j j i i j i i j i j j i i j i j i i j i j j i i j i j j i i j i j i j i j j i i j i j i i j i j j i i j i j i j i j j i i j i j j i i j i i j i j C C C C S S C S C S C C C C C C C C S C S C S C C S S S C C S S C C C S C C C C S C C S S S S C C S C S C C C S S C S S C C S C S C C C S C S S C S S C S S S C C S S S C S C S C C C C C S C C C C S C C S C C C C C S S C S C C C S C C C C C S C S C S S C C C C C C S C C C S S S C C S C S C C C S C C C C C C C S S C C S C S S C C C C S S C m T dt d j j i j j j i j j j i j j j j j j j j j j n k k j k i y y f f q q y y f f q y y f f q q y y f f q y y f f q q y y f f q y y f f q q y y f f q q f f q q y y f f q q y y f f q f f q q y y f f q y y f f q q y y f f q q f f q q y y f f q q y y f f q f f q q y y f f q y y f f q q y y f f q q f f q q y y f f q q y y f f q f f q q y y f f q y y f f q q y y f f q q y y f f q q f f q q y y f f q y y f f q q y y f f q q f f q q y y f f q q y y f f q f f q q y q y f y q q f y q f f y q f y f q y f q y q f q                              = = + +   +  + +  + +  + +  + +   +   +   +  + +   + + + +  +   +  + + +  + + + = ÷ ÷ø ö ç çè æ ¶ ¶ åå & & & & & & & & & & & & & & & & & & & & & & & & & 2 2 2 1 1 2 2 2 (3.23) Again, the time rate of change of the derivative of kinetic energy with respect to the time rate of change of the generalized coordinate, equation 3.23, contains terms that are equal to the derivative of kinetic energy with respect to the generalized coordinate, equation 3.21. Combining equation 3.21 and 3.23 appropriately produces the left hand side of Lagrange’s equations for the attitude angle generalized coordinate, 71 { ( ) ( ) ( ) ( ) ( ) ( ) ( ) ( ) j ( j i i j i j i i j i j )} i j j i i j i j i j i j j j i i j i j j i i j i i j i j j i i j i j i i j i j j i i j i j j i i j i j i j i j j i i j i j i i j i j j i i j i j i j i j j i i j i j j i i j i i j i j S C S C C C C C S Q C C C C S C C S C C C C C S S C S C C C S C C C C C S C S C S S C C C C C C S C C C S S S C C S C S C C C S C C C C C C C S S C C S C S S C C C C S S C m T T dt d j j j j j j j j j n k k j k i i y y f f q y y f f q q q y y f f q q f f q q y y f f q q y y f f q f f q q y y f f q y y f f q q y y f f q q y y f f q q f f q q y y f f q y y f f q q y y f f q q f f q q y y f f q q y y f f q f f q q y q f y f q y f q y q f q q   =  + +   + + + +  +   +  + + +  + + + = ¶ ¶  ÷ ÷ø ö ç çè æ ¶ ¶               = = åå 2 2 2 1 1 2 2 2 & & & & & & & & & & && & & (3.24) The right hand side of Lagrange’s equations is the derivation of the generalized force for the ith segment. The evaluation of equation 3.12 for the attitude angle generalized coordinate results in, [ ( ) ( ) ( ) ( )] for (i n) Q M F l C C C F l C C S F m g l C S n j i i i X j i i i i Y j i i i i Z j j i i i = 1, 2,K,  + + + + = å= q q f q y f q y f q (3.25) Thus, concludes the derivation of the equations of motion for the attitude angle. However, there still remains one generalized coordinate not yet applied to the Lagrange equations. 3.2.2.3 Heading Angle The final generalized coordinate is the heading angle. The equations of motion for the heading angle are derived in the same fashion as that of the other Euler angles. The 72 derivation of the Lagrange equations for the heading angle of the ith segment as shown in equation 3.11 is presented in detail in the remainder of this section. The derivative of kinetic energy with respect to the heading angle and time rate of change of heading angle of the cable segment is evaluated below, { [ ( ( ) ( ) ( ) ( ( ) ( )] [ ( ) ( ) ( )] [ ( ( ) ( ) ( ) ( ( ) ( )]} j i i j i j i j j i i j i i j j j i i j j j i i j i j j i i j i i j j j i i j i j i j j i i j i j i i j i j j i i j i j j i i j i j j i i j i j i i j i i j j j i i j i j i j j i i j i j j i i j j j i i j i j i j j i i j j i j i S C C S S S S C C S S S C S S S C C C C C S C S C C S S S S C C S S S S C S C S C C C C C S S C C C C S C S C S C C C C S S S S C C S C S C S S C S S S C S C C C C C S S S S S C C C C S S S C S m T j j i j j j i j j j i j n k k j k i y y f f q q f f y y f f q f f q y y f f q y y f f q q y y f f q f f q y y f f q q f f y y f f q y y f f q q y y f f q q y y f f q q y y f f q y y f f q f f q y y f f q q f f y y f f q q y y f f q y y f f q q f f y y f f q f f q y q y f y q q f y q f f y +    +     + + + + +   +  +   +  + + +   + = ¶ ¶                  = = åå & & & & & & & & & & & & 1 1 (3.26) { ( ( ) ( ) ( ) ( ( ) ( )} j i i j i j i j j i i j i i j j j i i j j j i i j i j j i i j i i j j j i i j i j i j C C C S S S S S C S S S C S C S C C S C C S C C C C S S S S S C S S S S C m T j j j n k k j k i y y f f q q f f y y f f q f f q y y f f q y y f f q q y y f f q f f q y y f f q q f f y q f y + +  + +   +  + + = ¶ ¶       = = åå & & & & 1 1 (3.27) The first term of Lagrange’s equations is attained by taking the time rate of change of equation 3.27. The results are as follows, 73 { ( ( ) ( ) ( ) ( ( ) ( ) ( ) ( ) ( ( ) ( ) ( ( ) ( ) ( ) ( ( ) ( ) [ ( ( ) ( ) ( ) ( ( ) ( )] [ ( ) ( ) ( )] [ ( ( ) ( ) ( ) ( ( ) ( )]} j i i j i j i j j i i j i i j j j i i j j j i i j i j j i i j i i j j j i i j i j i j j i i j i j i i j i j j i i j i j j i i j i j j i i j i j i i j i i j j j i i j i j i j j i i j i j j i i j j j i i j i j i j j i i j j i j i j j i i j i j i j j i i j i i j j j j i i j j j i i j i j j j i i j i i j j j i i j i j i j j i i j i j i j j i i j i i j j j i i j j j i i j i j j i i j j j i i j i j j i i j i j i j j i i j i i j j j i i j j j i i j i j j i i j i i j j j i i j i j i j S C C S S S S C C S S S C S S S C C C C C S C S C C S S S S C C S S S S C S C S C C C C C S S C C C C S C S C S C C C C S S S S C C S C S C S S C S S S C S C C C C C S S S S S C C C C S S S C S S C C S S S S C C S S S C S C S C S S C C S S C C S S S C S S C C S S S S C C S S S S C S C C S S S S S S C C C C C S C C S S C S C S S C C C C S S S S S C S S S C S C S C C S C C S C C C C S S S S S C S S S S C m T dt d j j i j j j i j j j i j j j j j j j j j j n k k j k i y y f f q q f f y y f f q f f q y y f f q y y f f q q y y f f q f f q y y f f q q f f y y f f q y y f f q q y y f f q q y y f f q q y y f f q y y f f q f f q y y f f q q f f y y f f q q y y f f q y y f f q q f f y y f f q f f q y y f f q q f f y y f f q f f q y y f f q y y f f q q y y f f q f f q y y f f q q f f y y f f q q f f y y f f q f f q y y f f q y y f f q q y y f f q y y f f q q y y f f q q f f y y f f q f f q y y f f q y y f f q q y y f f q f f q y y f f q q f f y q y f y q q f y q f f y q f y f q y f q y q f y +    +     + + + + +   +  +   +  + + +   +  + +  +   +  +   +  + + + +  + +   + + +  + + +   +  + + = ÷ ÷ø ö ç çè æ ¶ ¶                                    = = åå & & & & & & & & & & & & & & & & & & & & & & && & & 2 2 2 1 1 2 2 2 (3.28) As with the previous derivations of Lagrange’s equations for the cable system the time rate of change of the derivative of kinetic energy with respect to the time rate of change of the generalized coordinate contains terms that are equal to the derivative of kinetic energy with respect to the generalized coordinate, somewhat simplifying the equations of motion. 74 Combining equation 3.26 and 3.28 appropriately produce the left hand side of the Lagrange’s equations of the heading angle generalized coordinate. { ( ( ) ( ) ( ) ( ( ) ( ) ( ) ( ) ( ( ) ( ) ( ( ) ( ) ( ) ( ( ) ( ) } j j i i j i j i j j i i j i i j j i j j i i j j j i i j i j j j i i j i i j j j i i j i j i j j i i j i j i j j i i j i i j j j i i j j j i i j i j j i i j j j i i j i j j i i j i j i j j i i j i i j j j i i j j j i i j i j j i i j i i j j j i i j i j i j S C C S S S S C C S S S C S Q C S C S S C C S S C C S S S C S S C C S S S S C C S S S S C S C C S S S S S S C C C C C S C C S S C S C S S C C C C S S S S S C S S S C S C S C C S C C S C C C C S S S S S C S S S S C m T T dt d j j j j j j j j j n k k j k i i y y f f q q f f y y f f q f f q y y y f f q y y f f q q y y f f q f f q y y f f q q f f y y f f q q f f y y f f q f f q y y f f q y y f f q q y y f f q y y f f q q y y f f q q f f y y f f q f f q y y f f q y y f f q q y y f f q f f q y y f f q q f f y q f y f q y f q y q f y y  + +  =   +  +   +  + + + +  + +   + + +  + + +   +  + + = ¶ ¶  ÷ ÷ø ö ç çè æ ¶ ¶                   = = åå 2 2 2 1 1 2 2 2 & & & & & & & & & & && & & (3.29) The final addition to the dynamics equations are the effect of the applied loads on the heading angle. The generalized force for the heading angle is derived in accordance with equation 3.12 and result in the following, [ ( ) ( )] for (i n) Q M F l S C C S S F l S S C S C n j i i i X j i i i i i i Y j i i i i i i = 1, 2,K, + +  + = å= y y f y f q y f y f q y (3.30) The equations of motion of the cable segments are dependant upon the equations of motion of the Euler angles of each segment. Simulation of the cable motion requires simultaneous solution of the angular acceleration of the generalized coordinates for all cable segments. 75 3.2.2.4 Cable Equations of Motions We have yet to discuss the threedimensional aerodynamic applied loads, but neglecting that, the equations of motion for the threedimensional cable have been derived, although they are not in a convenient form. We can rewrite the equations of motion in matrix form for an arbitrary number of cable segments as was done in the twodimensional case. The fully threedimensional nonlinear differential equations of motion for a cable of n segments can be written in the following form, ú ú ú û ù ê ê ê ë é = ú ú ú û ù ê ê ê ë é ú ú ú û ù ê ê ê ë é + ú ú ú û ù ê ê ê ë é ú ú ú û ù ê ê ê ë é + ú ú ú û ù ê ê ê ë é ú ú ú û ù ê ê ê ë é y q f yf qy fq y q f y q f Q Q Q C C C C C C C C C B B B B B B B B B A A A A A A A A A & & & & & & & & & & && && 7 8 9 4 5 6 1 2 3 2 2 2 7 8 9 4 5 6 1 2 3 7 8 9 4 5 6 1 2 3 (3.31) where A, B, C, and Q are augmented matrices. In other words, matrix A is composed of nine submatrices of size n x n. Thus A is of size 3n x 3n. The same is true for matrix B and C. A, B, and C are coefficient matrices. The generalized force column matrix, Q, is an augmented matrix of size 3n x 1. The values of the matrices components are taken directly from the results of the derivation of the Lagrange equation carried out in the previous sections, and combine to form a set of 3n equations of motion. The A, B, and C matrices in equation 3.31 are represented by the equations below, where i and j are the row and column of the matrices. [ ( ) ( ) ] j i j j i j i i j i j i j j i j i j i j S C S S S C S S S C C A m l l C S S S S C C i j i j y y f f q f f q f f q q y y f f q q f f  + = + +   i 1 , , i (3.32) ( ) i j j i j i j j i j j i j i j A m l l C S C S C S C C C S C C S i j i j y y f f q q y y f f q f f q q =   + 2 , ,  i  i (3.33) 76 [ ( ) j ( i j i j i j )] i j j i j i i j j S S C S S C S A m l l C S S S C C S i j i j y y f f q q f f y y f f q f f q  =  + +   i 3 , , i (3.34) i j j i A A 4 , 2 , = (3.35) ( ) i j j i j i j i j i j A m l l C C C C C C C S S i j i j y y f f q q f f q q = + 5 , ,  i (3.36) ( ) i j j i j i j i j i j A m l l C C S C S C C C S 6 , i, j i j y yi f f q y yi f f q q =  (3.37) i j j i A A 7 , 3 , = (3.38) i j j i A A 8 , 6 , = (3.39) [ ( ) ( )] j i j i i j j i j j i j i j i j S C S S S C S A m l l C C C S S S S i j i j y y f f q f f q y y f f q q f f  = + +   i 9 , , i (3.40) [ ( ) ( ) ] j i j j i j i i j i j i j j i j i j i j S C C S S S S S C C C B m l l C S C S S C S i j i j y y f f q f f q f f q q y y f f q q f f + + =  +   i 1 , , i (3.41) ( ) i j j i j i j j i j j i j i j B m l l C S C S S S C C S S C C C i j i j y y f f q q y y f f q f f q q = + + 2 , ,  i  i 



A 

B 

C 

D 

E 

F 

I 

J 

K 

L 

O 

P 

R 

S 

T 

U 

V 

W 


