University of Southern California Dissertations and Theses

Synthesis and analysis of high-performance controllers for high-speed autonomous aerobatic flight

(USC Thesis Other)

Synthesis and analysis of high-performance controllers for high-speed autonomous aerobatic flight

PDF

Download a page range

Download transcript

Copy asset link

Request this asset

Transcript (if available)

Content Synthesis and Analysis of High-Performance Controllers for High-Speed Autonomous Aerobatic Flight by Ying Chen A Dissertation Presented to the FACULTY OF THE USC GRADUATE SCHOOL UNIVERSITY OF SOUTHERN CALIFORNIA In Partial Fulllment of the Requirements for the Degree DOCTOR OF PHILOSOPHY in Mechanical Engineering December 2019 Synthesis and Analysis of High-Performance Controllers for High-Speed Autonomous Aerobatic Flight Copyright © 2019 by Ying Chen ii To my parents Hongde Chen & Chaonan Lyu iii Acknowledgements First and foremost, I would like to thank my doctoral advisor, Professor N estor O. P erez-Arancibia, for being an excellent mentor to me. I am grateful that N estor invited me to join the Autonomous Microrobotic Systems Laboratory (ASML) back in 2014 and provided me with all the resources necessary to complete my doctoral research work. As an advisor, N estor has taught me how to think and conduct research in a scientic way and given me enough freedom to pursue what I feel passionate about, for which I am very grateful. Finally, I would like to thank him for his tremendous patience and much advice in the past ve years to encourage me to achieve my objective. I would like to thank Professor Henryk Flashner and Professor Mihailo R. Jovanovi c for serving on my dissertation committee and providing their comments, advice and con- structive questions on my dissertation. Furthermore, I thank Professor Eva Kanso and Professor Petros Ioannou for serving on my qualifying exam committee and providing their comments on my research. I want to thank all my labmates for accompanying me in this long journey and being such wonderful friends, and I will miss the time with you guys. Specically, I thank Longlong Chang for his help with his excellent programming debugging skills; I thank Joey Ge for the joy he brings to everyone; I thank Ariel Calder on for always being willing to help everyone anytime; I thank Xiufeng Yang for being such a great research partner to work with; I thank Coco Xu for the interesting conversations about anime series; I also thank Xuan-Truc Nguyen for helping me to improve my English writing iv skills, and you are one of the nicest persons I have ever met. I would like to acknowledge my other friends: Xiaowei Ding, Hao Jiang, Dingyuan Chen, Min Zheng, Sichen Yuan, Hao Gao, Qian Wan, Wolfgang H onig, and Zaki Has- nain. Your friendship and encouragements made this journey much more enjoyable and unforgettable. I also would like to thank all the people who helped me in the past ve years. Finally, my warmest thanks belong to my parents, Mr. Hongde Chen and Mrs. Chaonan Lyu. I would not have been able to accomplish my dissertation without your consistent support and endless love. This dissertation is dedicated to my parents. v Contents Acknowledgements iv List of Figures ix List of Tables xix Abstract xx Chapter 1 Introduction 1 1.1 Motivation and Related Work . . . . . . . . . . . . . . . . . . . . . . . 1 1.2 Contributions of the Dissertation . . . . . . . . . . . . . . . . . . . . . 6 1.3 Organization of the Dissertation . . . . . . . . . . . . . . . . . . . . . . 7 1.4 Notation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 7 Chapter 2 Linear Controller Design and Analysis for Aerobatic Flight 9 2.1 System Models and Multi-Flips . . . . . . . . . . . . . . . . . . . . . . 9 2.2 Controller Synthesis and Stability Analysis . . . . . . . . . . . . . . . . 13 2.2.1 Multi-Flip Speed Planning . . . . . . . . . . . . . . . . . . . . . 13 2.2.2 Control Strategy and Closed-Loop Scheme . . . . . . . . . . . . 16 2.2.3 Stability of the Closed-Loop Angular Velocity Dynamics . . . . 19 2.2.4 Closed-Loop Attitude Analysis . . . . . . . . . . . . . . . . . . 31 2.2.5 The Eects of Controller Switching and Stability . . . . . . . . 34 2.3 Simulation and Experimental Results . . . . . . . . . . . . . . . . . . . 35 2.3.1 The Flying Vehicle and Experimental Architecture . . . . . . . 35 2.3.2 Generalized Multi-Flip and Generalized Flipping Angle . . . . . 36 vi 2.3.3 Simulation Results . . . . . . . . . . . . . . . . . . . . . . . . . 38 2.3.4 Experimental Results . . . . . . . . . . . . . . . . . . . . . . . . 38 2.3.5 Comparison with a State-of-the-Art Controller . . . . . . . . . . 42 Chapter 3 Performance Optimization of Aerobatic Flight 45 3.1 Formulation of the Optimization Problem . . . . . . . . . . . . . . . . 45 3.2 Simulation and Experimental Validations . . . . . . . . . . . . . . . . . 51 3.2.1 Optimization and Experimental Setups . . . . . . . . . . . . . . 51 3.2.2 Optimization, Simulation and Experimental Results . . . . . . . 52 Chapter 4 Nonlinear Adaptive Control for Aerobatic Flight without Torque Measurements 60 4.1 Preliminaries . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 61 4.1.1 System Modeling . . . . . . . . . . . . . . . . . . . . . . . . . . 61 4.1.2 The Process of the Aerobatic Maneuver . . . . . . . . . . . . . . 63 4.1.3 Momentum Theory Analysis for Aerobatic Maneuver . . . . . . 64 4.2 Nonlinear Adaptive Controller Design without Torque Measurements . 69 4.2.1 Backstepping Control . . . . . . . . . . . . . . . . . . . . . . . . 69 4.2.2 Nonlinear Adaptive Control . . . . . . . . . . . . . . . . . . . . 73 4.3 Experimental Results . . . . . . . . . . . . . . . . . . . . . . . . . . . . 78 4.3.1 Experiment Setup and the UAV vehicle . . . . . . . . . . . . . . 78 4.3.2 Reference Angular Speed Plannings for Aerobatic Flight . . . . 79 4.3.3 Experimental Results with Dierent Controllers . . . . . . . . . 81 4.3.4 Experimental Results for Dierent Types of Maneuvers . . . . . 86 4.3.5 Data-Based Persistent Excitation Analysis . . . . . . . . . . . . 89 Chapter 5 Nonlinear Adaptive Control for Aerobatic Flight with Torque Feedback 92 5.1 Preliminaries . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 93 5.1.1 System Dynamics with LTV Actuator Dynamics . . . . . . . . . 93 vii 5.1.2 Attitude Error Kinematics . . . . . . . . . . . . . . . . . . . . . 95 5.1.3 A Class of Aerobatic Maneuvers . . . . . . . . . . . . . . . . . . 96 5.1.4 Control Objectives . . . . . . . . . . . . . . . . . . . . . . . . . 97 5.2 Adaptive Control Design and Stability Analysis . . . . . . . . . . . . . 98 5.2.1 Lyapunov-Based Adaptive Control . . . . . . . . . . . . . . . . 98 5.2.2 Modular Adaptive Control . . . . . . . . . . . . . . . . . . . . . 103 5.2.3 Attitude Analysis . . . . . . . . . . . . . . . . . . . . . . . . . . 110 5.3 Experimental Results . . . . . . . . . . . . . . . . . . . . . . . . . . . . 112 5.3.1 Experimental Setup and Example Aerobatic Maneuvers . . . . . 112 5.3.2 Torque Feedback and Control Implementation . . . . . . . . . . 112 5.3.3 Results . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 113 Chapter 6 Flight Control of Insect-Scale Flapping-Wing Robots 118 6.1 Insect-Scale Flapping-Wing Robotic Systems . . . . . . . . . . . . . . . 119 6.1.1 Two Systems . . . . . . . . . . . . . . . . . . . . . . . . . . . . 119 6.1.2 System Dynamics . . . . . . . . . . . . . . . . . . . . . . . . . . 121 6.1.3 Actuation Command Generation . . . . . . . . . . . . . . . . . 123 6.2 Attitude and Position Control Design . . . . . . . . . . . . . . . . . . . 125 6.2.1 Attitude Control . . . . . . . . . . . . . . . . . . . . . . . . . . 125 6.2.2 Position Control . . . . . . . . . . . . . . . . . . . . . . . . . . . 126 6.3 Controlled Flight Experiments . . . . . . . . . . . . . . . . . . . . . . . 127 6.3.1 Experimental Setup . . . . . . . . . . . . . . . . . . . . . . . . . 127 6.3.2 Position Control of RoboBee . . . . . . . . . . . . . . . . . . . . 128 6.3.3 Altitude and Attitude Control of Bee + . . . . . . . . . . . . . . 130 6.3.4 Position Control of the Bee + . . . . . . . . . . . . . . . . . . . . 132 Chapter 7 Conclusions and Future Work 134 Bibliography 137 viii List of Figures 2.1 The yer and frames of reference. N =fO 0 ;n 1 ;n 2 ;n 3 g is the iner- tial frame,B =fO B ;b 1 ;b 2 ;b 3 g is the body frame, and r indicates the position of the yer's center of mass relative to the inertial originO 0 . 10 2.2 The complete process of a multi- ip maneuver is composed of three phases: the climb, the multi- ip, and descent-and-reorientation. At the beginning, the quadrotor hovers at a certain altitude and then is triggered to begin the multi- ip maneuver. The three phases are autonomously controlled employing only on-board sensors and controllers. . . . . . . 12 2.3 Generic shape of the desired ipping speed,! d (t), employed in the imple- mentation of the proposed control method and stability analysis. This signal is composed of three sections: acceleration, constant speed and deceleration. According to this denition, g;1 , g;2 and g;3 denote the amounts of rotation accumulated angle during each of the three sections. The total ipping angle is equal to P 3 i=1 g;i = 2n, where n is the number of ips. The maximum allowed ipping speed is denoted by ! max . This reference ipping speed corresponds to the multi- ip phase in Figure. 2.2. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 15 ix 2.4 Control scheme employed during the multi- ip phase of the aerobatic maneuver depicted in Fig. 2.2. In order to generate the control signal, the controller uses the angular-velocity and angular-acceleration control errors, and the reference torque, d , that is associated with the reference signals! d and _ ! d . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 16 2.5 Architecture of the experimental setup. The ground computer and joy- stick are mainly used to send references of attitude and altitude to the quadrotor and compensate drift during regular ight. The ground com- puter is also receives the real-time quadrotor data for o-line analysis. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 36 2.6 Simulation results. Plots (a), (b) & (c) show simulation results of single- , double- and triple- ip high-speed maneuvers about the body-xed axis b 1 employing the linear controller design. Plots (d), (e) & (f) show the simulations of single-, double- and triple- ip high-speed maneuvers about the body-xed axisb 0 1 employing the linear controller design. . . 39 2.7 Experimental results. Plots (a), (b) & (c) show experimental results of single-, double- and triple- ip high-speed maneuvers about the body- xed axisb 1 employing the linear controller design. Plots (d), (e) & (f) show real-time experiments of single-, double- and triple- ip high-speed maneuvers about the body-xed axisb 0 1 employing the linear controller design. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 40 2.8 Real-time ight experiments. Composite images of six high-speed ma- neuvers. Trajectories (a), (b) & (c) correspond to single-, double- and triple- ip maneuvers about the body-xed axisb 1 . Trajectories (d), (e) & (f) correspond to single-, double- and triple- ip maneuvers about the 45 -oblique axisb 0 1 = (b 1 +b 2 ) = p 2. . . . . . . . . . . . . . . . . . . . . . . 41 x 2.9 Comparison of the proposed control method with the state-of-the-art geometric controller in [1]. Both controllers are experimentally imple- mented to perform consecutive triple- ips about theb 1 axis. (a) Desired and measured angular rates of the quadrotor while own by the geometric controller. In this case, the quadrotor can nish only three consecutive triple- ips before crashing to the ground; in fact, the rst impact occurs during the DR phase of the second triple- ip. (b) Desired and mea- sured angular rates of the quadrotor while own by the control scheme in Fig. 2.4. In this case, the quadrotor can indenitely complete con- secutive triple- ips; this plot shows 13 successfully-completed aerobatic maneuvers which demonstrate both the high performance and stabil- ity robustness (against noise and disturbances) of the proposed control method when compared to that presented in [1]. . . . . . . . . . . . . 43 3.1 Optimization results for four dierent aerobatic maneuvers: (a) double ip about the b 1 axis; (b) triple ip about the b 1 axis; (c) double ip about a 45 axis as dened in Chapter 2.3.4; and (d) triple ip about a 45 axis as dened in Chapter 2.3.4. In all the cases, each data point rep- resents the solution to a quasiconvex optimization problem as specied by (3.8). The horizontal axes of the plot represent the resulting optimal values of 3 1 2 and 1 1 2 , and surface represents the values of bound as dened in (3.2). . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 53 xi 3.2 Exhaustive simulation results obtained through the performance of four types of aerobatic maneuvers. The plots labeled with (a), (c), (e) and (g) show the sim surfaces, i.e., values of obtained through simulations, as functions of 3 1 2 and 1 1 2 . The plots labeled with (b), (d), (f) and (h) show the top views of the 3D plots on the right. In all the plots, each data point represents the solution to a quasiconvex optimization problem as specied by (3.8). . . . . . . . . . . . . . . . . . . . . . . . 54 3.3 Comparison of selected experimental and simulation results. In these plots, each data point represents the relationship between a simulated sim value and an experimental exp value corresponding to the same set of controller and analysis parameters. It can be seen that, in most cases, exp decreases as sim decreases, which validates the experimental eectiveness of the proposed method for performance optimization. . . 57 3.4 Comparison of the ight performances achieved during the execution of triple- ips about the b 1 axis with two dierent controller parameter sets. (a) Experimental result obtained with the computationally-found globally-optimal parameter set. (b) Experimental result obtained with the heuristically-tuned non-optimal parameter set. In these experiments, after the triple- ip maneuver is completed,je ! (t)j 2 coincides withj!(t)j 2 becausej! d (t)j 2 becomes zero. . . . . . . . . . . . . . . . . . . . . . . 58 4.1 The 28-gram micro quadrotor and frames of reference.N =fO 0 ;n 1 ;n 2 ;n 3 g is the inertial frame,B =fO B ;b 1 ;b 2 ;b 3 g is the body frame, and r in- dicates the position of the yer's center of mass from the inertial origin O 0 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 61 xii 4.2 The complete process of an Aerobatic Maneuver. The whole process is composed of three phases, climb phase, aerobatic maneuver phase, and descent and restabilization (DR) phase. In climb and DR phases, the UAV stays in the regular ight mode with small angle variations and low angular velocity near hovering state. In the aerobatic maneuver phase, the UAV stays in aerobatic ight mode with large angle variations and high angular velocity. . . . . . . . . . . . . . . . . . . . . . . . . . . . 63 4.3 Flow model for the momentum analysis of an aerobatic quadrotor maneu- ver. The quadrotor rotates about theb 1 axis with angular speed! 1 , and the relative free-stream ow has speed v 1 resulted from the quadrotor's translational motion. . . . . . . . . . . . . . . . . . . . . . . . . . . . 65 4.4 Reference angular speed plannings of three aerobatic maneuvers. (a)-(c) show the reference angular speeds, ! d1 and ! d2 , of the three aerobatic maneuvers. ! d3 = 0 for all the maneuvers. . . . . . . . . . . . . . . . . 80 4.5 Performance comparison of the three control methods. This plot sum- marizes three sets of ten triple- ip experiments run using the LTI bench- mark controller, the backstepping controller and the nonlinear adaptive scheme. The average value of the PFM for each set is indicated with a square; each bar indicates the experimental range of for each set; the corresponding ESD for each set is written between square brackets. . . 81 4.6 Three triple- ip experiments with three dierent controllers. In (a), the linear controller in (2.12) achieves a PFM = 0:9391 rad. In (b), the backstepping controller in (4.28) compensates the time-invariant torque delay model, achieving a PFM = 0:7348 rad. In (c), the nonlinear adaptive controller (4.31) compensates the time-varying actuator model, achieving a PFM = 0:5340 rad. . . . . . . . . . . . . . . . . . . . . . 84 xiii 4.7 Filtered true and estimated torques during the triple- ip maneuver. f1 is the ltered true torque dened in (4.34). In (a), ^ 1 is estimated with the time-invariant torque model in (4.15). In (b), ^ 1 is estimated with the time-varying torque model in (4.46) and a parameter-updating law. 85 4.8 Performance comparisons of three maneuvers with linear and adaptive controllers. In (a), the average with the linear controller (2.12) is 1:2240, and the average with the adaptive controller (4.31) is 0:7233. In (b), the average with the linear controller is 1:1304, and the average with the adaptive controller is 0:7385. In (c), the average with the linear controller is 1:9715, and the average with the adaptive controller is 1:4160. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 87 4.9 Three dierent maneuvers under two dierent controllers. The reference and true angular velocities of six example maneuvers are demonstrated here. For triple ip and Pugachev's Cobra maneuvers, the reference signals onb 2 andb 3 axes are equal to zero, so are not plotted here. For mixed ipping maneuver, the reference signal onb 3 is zero. . . . . . . 88 4.10 Filtered estimated and true torques in three maneuvers. In (a)(b), the ltered estimated torque s+ ^ 1 and the ltered true torque f1 on b 1 axis are given. In (c), the ltered estimated and true torques of the ipping axes are plotted respectively. Note that the ltered estimate torques remain zero before maneuvers. . . . . . . . . . . . . . . . . . . 89 4.11 The estimated ^ a i (t) and ^ b i (t) in three maneuvers. The estimated ^ a 1 (t) and ^ b 1 (t) are plotted for the triple ipping and Pugachev's Cobra ma- neuvers. For the mixed ipping, the ipping axis is alternatively changed betweenb 1 andb 2 axes, so the estimated parameters of the two axes are plotted respectively. Note that the estimated parameters ^ a i (t) and ^ b i (t) remain constant before and after each maneuver. . . . . . . . . . . . . 90 xiv 5.1 Triple ip using three dierent controllers. (A) Triple ip using the linear controller given by (2.12) that ignores LTV actuator dynamics; the corresponding PFM = 0:9391. As! d2 =! d3 = 0 for the triple ip, they are not plotted here. (B) Triple ip using the Lyapunov-based adaptive controller (5.18); the corresponding PFM = 0:6397. (C) Triple ip using the modular adaptive controller (5.36); the corresponding PFM = 0:6286. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 114 5.2 Pugachev's Cobra using three dierent controllers. (A) Pugachev's Co- bra using the linear controller that ignores LTV actuator dynamics; the corresponding PFM = 0:7413. As! d2 =! d3 = 0 for Pugachev's Cobra, they are not plotted here. (B) Pugachev's Cobra using the Lyapunov- based adaptive controller; the corresponding PFM = 0:4625. (C) Pu- gachev's Cobra using the modular adaptive controller; the corresponding PFM = 0:4888. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 115 5.3 Mixed ip using three dierent controllers. (A) Mixed ip using the linear controller that ignores LTV actuator dynamics; the corresponding PFM = 1:1497. As ! d3 = 0 for mixed ip, it is not plotted here. (B) Mixed ip using the Lyapunov-based adaptive controller; the corre- sponding PFM = 0:6826. (C) Mixed ip using the modular adaptive controller; the corresponding PFM = 0:6920. . . . . . . . . . . . . . 116 5.4 Consecutive aerobatic maneuvers. (A) 11 consecutive Pugachev's Cobra maneuvers employing the Lyapunov-based adaptive control strategy. (B) 12 consecutive mixed ip maneuvers employing the modular adaptive control strategy. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 116 xv 6.1 Two dierent insect-scale apping-wing robotic systems. RoboBee (left) is a two-winged robot, and it is driven by two piezoelectric bimorph actuators. Bee + (right) is a four-winged robot that has a mass of 95 mg and measures 33 mm in wingspan, and it is driven by four piezoelectric unimorph actuators. A U.S. one-cent coin indicates the scale. . . . . . 119 6.2 The principle of torque generation of RoboBee (top view of the robot). In (a), the roll torque is generated by diering the magnitudes of two wings. In (b), the pitch torque is generated by shifting the apping center position. In (c), the yaw torque is generated by using an asymmetric apping pattern. All gures are top view of the robot. . . . . . . . . . 121 6.3 The principle of torque generation of Bee + . In (a), the roll torque is generated by diering the magnitudes of two opposed pairs of apping wings. In (b), the pitch torque is generated by diering the magnitudes of two pairs of apping wings that are on the same side. In (c), the stroke plane of the wing apping is pre-set or adjusted to have an angle with respect to the steering plan, and the yaw torque is generated by diering the projected forces of two diagonally opposed pairs of lift forces of two apping wings onto the steering plane. . . . . . . . . . . . . . . . . . . 122 6.4 Frames of references of two insect-scale apping-wing robots. In (a), the frames of references of the RoboBee is indicated. In (b), the frames of references of the Bee + is indicated. For both robots,fn 1 ;n 2 ;n 3 g is the inertial frame, andfb 1 ;b 2 ;b 3 g represents the body-xed frame with the origin located at the center of mass of the robot. . . . . . . . . . . . . 123 xvi 6.5 The position control of the RoboBee insect-scale apping-wing robot. In (a), then 1 position of the robot, the green line is the reference signal. In (b), then 2 position of the robot and green line is the reference. In (c), n 3 position (altitude) of the robot, and the green line is the reference signal. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 129 6.6 The roll and pitch angle of the RoboBee during the hovering. The blue line indicates the Euler Roll angle of the apping-wing robot, and the red line lines indicates the Euler Pitch angle of the apping-wing robot. 130 6.7 The altitude and attitude control experimental of the Bee + . A. It demon- strates the altitude in the experiment, in which the dash line represents the reference altitude. B. It shows the Euler roll and pitch angles during the experiment. Note that the angular oscillation remain approximately between10° and 10°, which is partially caused by the oscillation of the robot body induced by the wing apping and stays in an acceptable range. The experiment lasts for approximately 5 s, and after that the robot ies out of the specied safety space and the power is then turned o automatically. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 131 6.8 A. The time-lapse plot of the altitude and attitude control experiment. The corresponding altitude and Euler angles are demonstrated in Fig. 6.7. During the experiment, the direction of the thrust force is controlled to remain approximately perpendicular to then 1 -n 2 plane. The cable teth- ered to the robot provides the power and control signals. The robot drifts on then 1 -n 2 plane due to the lack of control of the position. After 2 s, the robots ies out of focus area of the camera. B. Image sequences of the position control experiment. The corresponding position and Euler angles of this experiment are plotted in Fig. 6.9. The images are captured with the identical view position. . . . . . . . . . . . . . . . . . . . . . 132 xvii 6.9 The position control experiment result. A. The dash lines represent the desired position, and the solid lines represent the measured posi- tion. B. The dash lines represent the desired Euler angles derived from (6.11)(6.12)(6.13). The solid lines represent the measured Euler angles. 133 xviii List of Tables 3.1 Pre-Specied Discretization of the Spacefk I ;;g. . . . . . . . . . . . 52 4.1 PE check of the regressor vector i ms i in experiments . . . . . . . . . . . 91 6.1 Comparison of parameters for Bee + and RoboBee . . . . . . . . . . . . 119 xix Abstract The capability to perform aerobatic ights is essential for aerial animals to survive in wild environments and accomplish many dierent types of tasks, such as capturing prey, escaping predators, and avoiding obstacles. Consequently, the development of fully autonomous ying robots will not become feasible until the set of basic aerobatic ying capabilities observed in nature can be automatically replicated by articial ying robots, just because that aerial robots need to survive in the wild or in unstructured environments that are full of unanticipated challenges. In the work presented in this dissertation, we aim to develop high-performance controllers and methods of analysis to enable aerial robots to achieve and even exceed the maneuverability exhibited by natural ying animals. The study of the ight control of unmanned aerial vehicles (UAVs) and micro aerial vehicles (MAVs) under dierent ight scenarios has attracted a great deal of attention in past decades. A signicant number of analytical studies have proposed numerous dierent controller designs with rigorous stability analyses of the system; these studies, however, do not present compelling experimental evidence to demonstrate the feasibil- ity of the resulting controllers during extreme aerobatic ights. On the other hand, some experimental research has demonstrated the impressive aerobatic ight capabili- ties of UAVs employing learning-based control methods; however, these works do not present rigorous analyses to explain why, how and when these solutions work. In this dissertation, we combine the analytical approach with the experimental approach to xx propose controller solutions which are examined by both rigorous theoretical analysis and challenging aerobatic ight experiments. First, we present a linear controller design that is based on the nominal dynamics of the system in order to enable a quadrotor UAV to perform aerobatic maneuvers that are characterized by rapid rotational motion in an extremely short period of time. We present a rigorous stability analysis that validates the control method and provides the stability conditions of the closed-loop system; compelling experimental results are achieved with the linear controller. Then, an optimization problem is formulated to syn- thesize the linear controller with the objective of improving the aerobatic performance. Here, an algorithm is proposed to solve the problem by nding the approximately globally optimal parameter set, employing the convex optimization. Then, both the simulation and experimental results are employed to validate the eectiveness of the optimization algorithm. From the experimental results obtained using the linear controller, the aerodynamic in uences of the varying ow elds surrounding the quadrotor UAV during rapid aero- batic ight are observed to have a noticeable negative eect on the ight performance of the system. A linear time-varying actuator dynamics model is added to the nominal plant to describe these aerodynamic eects. An adaptive control method that does not require torque measurements was developed to compensate for the negative eects of the actuator dynamics on the system. The experimental results obtained with the adap- tive control method demonstrate signicant performance improvements. Furthermore, two dierent adaptive control methods that employ torque feedback were developed to compensate for the aerodynamic eects, and rigorous stability analyses for these two approaches are provided to ensure functionality, stability and performance robustness and also explain the mechanisms that make these methods work. Consistent improve- ments in aerobatic performance are observed in the experiments executed using these two adaptive control strategies. xxi Finally, to mimic the incredible maneuverabilities and ying capabilities exhibited by aerial insects, we developed nonlinear attitude and position controllers for two dier- ent types of insect-scale apping-wing robots that are much smaller than the quadrotor UAV. This is the initially necessary step toward enabling aerobatic ights in insect-scale ying robots. The experimental results demonstrate that insect-scale ying robots are able to perform regular ights and indicate the potential capability for aerobatic ights. xxii Chapter 1 Introduction 1.1 Motivation and Related Work A great number of animal species have evolved to have superb ying capabilities that are based on apping-wing mechanisms. For some animals, the ability to execute aerobatic ight is essential for survival, as specialized high-speed maneuvers are employed to capture prey, escape predators, avoid obstacles, and perform mating rituals [2]. For instance, small animals such as fruit ies of the species Drosophila hydei, which weigh approximately 1 mg, have been reported to perform rapid banked turns at angular speeds as high as 5:3 10 3 s 1 to evade looming images on a screen [3]; similarly, fruit ies of the species Drosophila melanogaster can perform 90 banked turns in approximately 50 ms with yaw accelerations as high as 1 10 5 s 2 while freely ying within an experimental arena [4]. Larger animals such as Magnicent hummingbirds, Eugenes fulgens, have been recorded to perform an escape maneuver with a maximum roll rate of greater than 4 10 3 s 1 in 0:25 s [5]. Similarly, barn swallows (species Hirundo rustica), which typically weigh approximately 18 g, also exhibit impressive aerobatic skills; for example, they have been lmed diving into water and subsequently ying at high speeds [6] and also recorded to roll at rates as high as 5 10 3 s 1 [7], 1 thereby displaying unparalleled agility. In general, eld observations, analyses and experiments suggest that apping-wing ight is superior to other forms of aerial locomotion in terms of maneuverability and, possibly, aerodynamic eciency at low Reynolds numbers (from 10 2 to 10 5 ) and small scales (<15 cm in wingspan) [8{11]. In particular, despite the impressive progress that has been achieved in the past few years, the aerobatic capabilities exhibited by human- designed small 1 rotor-based or apping-wing-based yers remain far behind those of ying animals. We hypothesize, however, that this seemingly insurmountable gap in aerobatic ight capabilities can be reduced by employing high-performance control methods. This dissertation is motivated by the notion that the development of fully au- tonomous ying robots will not become feasible until a set of basic aerobatic ying capabilities observed in nature are reproduced automatically by articial yers, as they will require skills similar to those of birds and insects to operate and survive in the wild or other unstructured environments. We focus on the generation and control of ma- neuvers that are characterized by rapid rotational motion in an extremely short period of time (less than one second), such as multi- ip maneuver with a maximum angular speed of up to 2:0 10 3 s 1 . Moreover, ying robots should be able to recover from the high-speed aerobatic ight back to hovering ight without stalling or crashing to the ground. In this dissertation, we perform experiments employing micro quadrotors with a frame diagonal of 14 cm and an associated Reynolds number on the order of 10 4 during regular operation (four 4.5-cm blades rotate at frequencies as high as 335 Hz). This system has characteristic values that are loosely comparable to those of small y- ing birds such as hummingbirds (Lampornis clemenciae), which have typical weights of 8:4 g, average wing lengths on the order of 8.5 cm and Reynolds numbers on the order of 1:1 10 4 during regular ight (the wings ap at frequencies as high as 23 Hz) [7,12]. 1 Arbitrarily, we employ the word small to denote systems that weigh less than 1 kg and operate at low Reynolds numbers (from 10 2 to 10 5 ). 2 In the past two decades, a substantial amount of research regarding the position and attitude control of unmanned aerial vehicles (UAVs) has been published. For example, [13] presents the design and implementation of a quaternion-based controller that is analytically shown to globally stabilize the hovering attitude of a quadrotor and experimentally demonstrated to regulate the system while positionally locked; a backstepping controller that is designed to locally track reference positions and yaw- angles for a quadrotor during ight is described in [14]; in [1], a geometric tracking control method for a quadrotor is shown to approximately globally track position and attitude references; a method for achieving global position tracking control of VTOL 2 - UAVs without linear velocity measurements is introduced in [15]; and the position control of UAVs in the presence of perturbations and disturbances is studied in [16], [17] and [18]; the adaptive position control for the loss-of-thrust caused by component failure is presented in [19]. Unfortunately, the eectiveness of the methods in [1], [15], [16] and [17] is veried only through simulation and experimental evidence is not provided, and the works in [13], [14], [18] and [19] present only simple experimental verications in which either the translational motion is locked or the yer stays in a near-hovering state that is characterized by small attitude angles and slow angular velocities. It should be noted that none of these methods have been experimentally demonstrated to be eective at controlling extreme aerobatic maneuvers, and performance optimization is not explicitly considered. Furthermore, in the presence of modeling errors, the highly nonlinear terms of the control schemes in [1,13{18] might excite the high-frequency unmodeled dynamics and amplify the sensor noise that inexorably appears during extreme aerobatic ight; to substantiate this statement, we present an example in Chapter 2.3.5. In addition to analytically based approaches, such as those mentioned above, re- cent experimentally driven research has demonstrated some of the outstanding ying 2 Vertical take-o and landing. 3 capabilities of quadrotors, including their ability to perform controlled aggressive ma- neuvers. For example, the methods in [20] and [21] enable quadrotors to respectively execute multi- ips and to y through tilted windows by using iterative learning control (ILC) strategies. In [20], the maximum angular rate during the execution of triple ips is approximately 1:8 10 3 s 1 . In both cases, [20] and [21], the learning strategies provide a procedure for optimizing the control parameters; stability analyses, however, are not presented and the mathematical conditions for parameter convergence are not discussed. Thus, global stability and optimality are not guaranteed. Additionally, the ILC strategy is impractical in unknown and unstructured real scenarios in which the maneuver must be executed successfully upon the rst attempt. In general, rigorous stability analyses of the closed-loop system dened by the ight dynamics of a yer and a feedback controller are necessary to fully understand the generation of high-speed aer- obatic maneuvers and improve ight performance. These considerations are essential for creating robots that have ying capabilities similar to those exhibited by natural yers. In this dissertation, we combine the analytical approach with the experimental approach to design high-performance controllers for the high-speed aerobatic maneu- vers and validate the control solutions with both analytical proofs and challenging experimental designs. Moreover, both the analytically driven and experimentally driven research men- tioned above ignore the aerodynamic eects on the system dynamics during the extreme aerobatic ights and assume that the near-hovering ight assumption on aerodynamic eects is still valid. The rapid translational and rotational velocities in the aerobatic ight dramatically change the magnitude and direction of the local incoming ow en- countered by propellers. And the natural animals have demonstrated remarkable ma- neuverability in the presence of the surrounding time-varying ow elds. Therefore, it is important to take into account the eects of the changing surrounding ow elds when designing controllers to achieve the comparable maneuverability as natural yers. 4 Many studies concerning the dynamic modeling and control of VTOL UAVs simplify aerodynamic eects by assuming a time-invariant ow eld (constant aerodynamic co- ecients) [22{24], and this is reasonable and approximately accurate for near-hovering ights that have slow translational and angular velocities with small attitude angles. However, this assumption of the invariant aerodynamic coecients becomes question- able for aerobatic ights. In addition, the uctuation of the loading moments on the actuators (normally propellers) caused by varying aerodynamic coecients can change the torque latency of the actuators and introduce another varying component for the actuator dynamics. On the other hand, the complicated ow eld variances are dicult to model and be integrated into the controller synthesis. Therefore, it is necessary to model these aerodynamic coecients and torque latency variations in a manner that captures the essence of the complicated dynamics and is simple enough for controller design and stability analysis, which, along with the corresponding control strategies to compensate for the negative eects of actuator dynamics, is one of the main focuses of this dissertation. In addition, the design, fabrication and control of insect-scale apping-wing robots driven by piezoelectric actuators have experienced a great deal of progress in recent years [25{29]. However, articial insect-scale ying robots do not demonstrate maneu- verabilities comparable to those exhibited by natural insects, and the capability of per- forming controlled ights is still constrained by the limited thrust force and torques due to the mechanical design and hardware robustness. To perform aerobatic maneuvers, insect-scale apping-wing robots should be equipped with sucient control authority to perform rapid translational and rotational motions, which require improvements in both mechanical design and controller synthesis. 5 1.2 Contributions of the Dissertation The primary contributions of this dissertation are the following: • A linear controller synthesis method is proposed to enable a VTOL UAV to per- form the aerobatic maneuvers dened by rapid rotational motion in an extremely short period of time. This procedure includes the trajectory planning, the linear controller design and a rigorous stability analysis of the resulting closed-loop sys- tem. Compelling experimental results validate the eectiveness and performance robustness of the control strategy. Furthermore, a performance comparison with a state-of-the-art controller during an aerobatic maneuver indicates the advan- tages of the proposed controller design in terms of robustness and computational eciency. • An optimization algorithm based on the linear controller synthesis method is proposed to nd an approximation of the globally optimal parameter set that optimizes the performance of an aerobatic maneuver. Both the simulation and experimental results validate the eectiveness of the optimization approach for improving the aerobatic performance. • A linear time-varying (LTV) actuator dynamic model with time-varying unknown parameters is proposed to describe the aerodynamic-coecient and torque-latency variations during rapid aerobatic maneuvers. The model is derived from momen- tum theory analysis and experimental observations. An adaptive control strategy without torque measurements is presented to compensate for the time-varying torque dynamics. The experimental results demonstrate that the aerobatic ma- neuver performance signicantly improves with the proposed adaptive control strategy, and this also indicates that the proposed LTV actuator model captures the essence of the complicated aerodynamic eects. 6 • Two dierent types of nonlinear adaptive control strategies with torque feed- backs, a Lyapunov-based adaptive control scheme and a modular adaptive con- trol scheme, are proposed to compensate for the undesired eects produced by the LTV actuator dynamics, and rigorous stability analyses for both methods are performed to prove the stability of the closed-loop systems. Compelling exper- imental results validate the eectiveness of both adaptive control designs, and consistent signicant performance improvements are shown. • Nonlinear attitude and position controllers are developed to enable two dierent types of insect-scale apping-wing robots, with two-wing and four-wing designs, to execute attitude- and position-controlled ights in experiments. 1.3 Organization of the Dissertation The dissertation is organized as follows. The proposed linear controller synthesis method for the execution of aerobatic maneuvers is presented in Chapter 2. An opti- mization algorithm based on the linear controller synthesis method for improving the aerobatic performance is discussed in Chapter 3. An LTV actuator dynamics and the adaptive control strategy without torque measurements are presented in Chapter 4. Then, two dierent types of adaptive control approaches with torque feedbacks and the stability analyses are provided in Chapter 5. In Chapter 6, the attitude and position controller designs are provided for the two dierent kinds of insect-scale apping-wing robots. Conclusions and future work are discussed in Chapter 7. 1.4 Notation R and R + denote the sets of reals and positive reals. Scalars are denoted by regular lower-case symbols, e.g., a; vectors are denoted by bold lower-case symbols, e.g., a; 7 matrices are denoted by bold upper-case symbols, e.g.,A; and quaternions are denoted by crossed bold lower-case symbols, e.g., q. The symbol denotes the quaternion prod- uct, e.g., q p, as dened in [30]. conv(X ) denotes the convex hull of setX , as dened in [31]. The symbol (or the strict form) is used to denote an inequality. Between vectors, it represents a componentwise inequality and between Hermitian denite ma- trices, it represents a matrix inequality. S k denotes the set of symmetrickk matrices and S k + denotes the set of symmetric positive denite kk matrices. The notation jj is employed to signify the absolute value of a scalar, and the notationjj 2 denotes the Euclidean norm of a vector. For a setS, the interior of the set is denoted by S and the boundary of the set is denoted by @S. The set dierence between two set is given byAnB =fx:x2A andx = 2Bg. Symbol min () represents the extraction of the minimum eigenvalue of a symmetric matrix. kk F denotes the Frobenius norm of a matrix. The symbolr represents the del operator. `tr(A)' represents the trace of a square matrixA. 8 Chapter 2 Linear Controller Design and Analysis for Aerobatic Flight 2.1 System Models and Multi-Flips Most of Vertical Take-O and Landing (VTOL) Unmanned Aerial Vehicles (UAVs) have the property that the attitude motion is fully actuated and the translational motion is under-actuated [16,32], such as the quadrotor vehicle and the ducted-fan aircraft. In this dissertation, we utilize the quadrotor UAV as an example, and the controller designs and analyses can be seamlessly applied on other VTOL UAVs. The example quadrotor UAV is shown in Figure. 2.1. In agreement with this gure, N =fO 0 ;n 1 ;n 2 ;n 3 g denotes the inertial frame andB =fO B ;b 1 ;b 2 ;b 3 g denotes the body-xed frame with its origin and axes coinciding with the center of mass and principal axes of the quadrotor, respectively. In this conguration, the system is actuated by the four thrust forces of magnitudesff 1 ;f 2 ;f 3 ;f 4 g that are generated by rotations of four propellers at angular speedsf 1 ; 2 ; 3 ; 4 g. As explained in [33], the total eect of the propellers can be described by a single force that acts on the robot's center of mass, f , and a torque . Here, we employ two dierent dynamic models of the yer to describe 9 Figure 2.1: The yer and frames of reference. N = fO 0 ;n 1 ;n 2 ;n 3 g is the inertial frame, B = fO B ;b 1 ;b 2 ;b 3 g is the body frame, andr indicates the position of the yer's center of mass relative to the inertial originO 0 . and control the entire aerobatic maneuver process: one for aerobatic ight, in which the quadrotor generates large attitude angle variations at high angular speeds; and one for regular ight, in which the quadrotor remains in a near-hovering state that is characterized by small attitude angle variations and slow angular velocities. Assuming a perfectly rigid body for the quadrotor, the rst model (for aerobatic ight) is dened by the complete nonlinear dynamics m r =mgn 3 +f b 3 ; (2.1) _ q = 1 2 q p; (2.2) J _ ! =!J! +; (2.3) where m is the total mass of the vehicle; r denotes the position of the yer's center of mass with respect to inertial frame N ; g is the acceleration of gravity constant (9:8 m s 2 ); f is the magnitude of the total thrust force that is generated by the four propellers, f ; J is a diagonal matrix employed to represent the yer's moment of inertia, expressed with respect to the body-xed frame B; ! is the quadrotor's angular velocity with respect toN with its components expressed inB; is the total torque acting on the vehicle that is generated by the rotors; q = [q 0 q 1 q 2 q 3 ] T is the 10 quaternion, dened as in [30], used to represent the attitude ofB with respect toN ; and denotes the quaternion multiplication of q and quaternion p = 0 ! T T . Note that the model described by (2.1), (2.2) and (2.3) is not an approximation as it represents the nominal nonlinear dynamics of the system completely and globally, which is essential for analyzing and controlling high-speed aerobatic maneuvers. And, this description does not consider gyroscopic eects, motor saturation or latency, or the aerodynamic disturbances that are generated by the vehicle's translational and rota- tional motions. Additionally, the direction of the total thrust force,f , is assumed to be perpendicular to theb 1 -b 2 plane by ignoring the blade- apping eects. Experimen- tal evidence discussed later will compellingly demonstrate that this model is sucient for synthesizing controllers for aerobatic ight. During regular ight (non-aerobatic ight), the attitude of the quadrotor remains in a near-hovering state with small attitude variations and slow angular velocities; therefore, it is reasonable to linearize the model that is specied by (2.1), (2.2) and (2.3) about the hovering xed point ( = = = _ = _ = _ = 0, _ r = r = 0) and neglect the higher order terms to describe the system in the regular ight mode, which yields, m r = 2 6 6 6 6 4 f f f mg 3 7 7 7 7 5 ; (2.4) J 2 6 6 6 6 4 3 7 7 7 7 5 =; (2.5) in whichf;; g are the Euler roll, pitch and yaw angles corresponding to the Z-Y-X convention. 11 Figure 2.2: The complete process of a multi- ip maneuver is composed of three phases: the climb, the multi- ip, and descent-and-reorientation. At the beginning, the quadrotor hovers at a certain altitude and then is triggered to begin the multi- ip maneuver. The three phases are autonomously controlled employing only on-board sensors and controllers. In this chapter, multi- ip maneuver is utilized as an example of the aerobatic ight, and the denition of the generalized class of aerobatic maneuvers will be discussed in Chapter 5. To perform a multi- ip maneuver, the vehicle follows a ight process consist- ing of three phases: climb, multi- ip, and descent-and-reorientation (DR), as illustrated in Figure. 2.2. The DR phase consists of two subphases: DR-1, which lasts for 0.15 s after the multi- ip phase has been completed, and DR-2, during which the quadrotor resumes regular hovering ight. At the beginning of the process, the quadrotor hovers at an arbitrary position close to the ground; after receiving an external trigger com- mand, the yer rises up for a pre-specied amount of time to gain an initial upward speed and begins the multi- ip phase; during the multi- ip phase, the quadrotor ips about a desired axis and remains in free fall; the multi- ip phase ends when the ipping angle reaches the predened target of 2n, wheren is the number of ips; nally, in the DR phase, the vehicle stops rotating, reorients its attitude such that the b 1 -b 2 plane becomes perpendicular to the gravitational pull and ceases its free fall. The exact def- inition of ipping angle is provided in [33]. During the multi- ip phase, the direction of the total thrust forcef , which remains perpendicular to theb 1 -b 2 plane, changes 12 rapidly and periodically; therefore,f has a negligible eect on the translational mo- tion of the vehicle while freely falling under the eect of gravity. Consistently, the climb phase is required to obtain an initial upward speed in order to partially compensate the dropping speed at the end of the multi- ip phase. Since (2.3) does not depend on (2.1) or (2.2), the control strategy introduced in this work maintains the translational dynamics of the robot, specied by (2.1), in open loop during the multi- ip phase. Note that the globally dened nonlinear dynamic models (2.1)(2.2)(2.3) are employed to de- scribe the multi- ip phase due to the large attitude angle and high angular velocity, and the linearized dynamic models (2.4)(2.4) are employed to describe the climb and DR phases since the attitude of the UAV stays in near-hovering state with small angle variation and low angular velocity. In order for a quadrotor to successfully complete a high-speed multi- ip, it must not only perform a series of precisely controlled rotations about a desired axis but also recover from its rapid loss of altitude and avoid crashing to the ground. Therefore, the execution of an aerobatic maneuver of this type imposes strict performance re- quirements on the controllers that y the robot, as the entire process lasts less than one second and the maximum-thrust-to-weight ratio of the yer is only 1:7. In all the experiments presented in this dissertation, during the entire aerobatic maneuver process, the quadrotor operates completely autonomously by using onboard sensors, power, and the required computational capabilities to run the signal processing and control algorithms. 2.2 Controller Synthesis and Stability Analysis 2.2.1 Multi-Flip Speed Planning The speed planning method for the multi- ip maneuver is introduced here, and the same method can be applied for the speed planning of other types of maneuvers which 13 will be discussed in Chapter 4 and 5. We introduce a method to control multi- ip quadrotor maneuvers which are completed in less than one second. Consistently, a typical multi- ip requires the precise tracking control of angular speeds that increase from 0 s 1 to 210 3 s 1 and subsequently decrease back to 0 s 1 in fractions of a second, which resembles, to some extent, the remarkable maneuverability of some ying animals [4,7]. By design, after considering the theoretical and practical implications, a generic angular speed reference employed to generate a multi- ip,! d (t), is composed of three sections: acceleration, constant speed and deceleration, as illustrated in Figure. 2.3. The acceleration and deceleration sections are described by cubic functions such that the angular acceleration reference remains continuous during the entire maneuver and the reference signal,! d (t), is compatible with the characteristic of the robot's actuators. Similarly, the constant speed section is limited by the maximum angular speed value that the onboard sensor can measure. For a given multi- ip maneuver, the body-xed desired ipping axis is denoted by the unit vector a b f = [a 1 a 2 a 3 ] T and, consistently, the angular velocity reference, expressed inB, is given by! d (t) =! d (t)a b f . Thus, consistent with the illustration in Figure. 2.3, the rst section of the ipping speed reference is constructed as ! d;1 (t) = 1 3 (t 1 ) 3 1 2 1 t + 1 3 1 3 ; (2.6) in which 1 = ! 1 max g;1 , 1 = ( 3 =4) 3 1 ! max , and _ ! max;1 = ( 3 =4) 1 g;1 ! 2 max . In this denition, g;1 is the desired total rotation angle accumulated during the acceleration section;! max is the maximum allowed ipping speed; and _ ! max;1 is the maximum allowed ipping angular acceleration, which must remain smaller than the maximum angular acceleration that is physically attainable by the quadrotor about the ipping axis. The duration of this rst section is 1 = 2 1 ; in (2.6), therefore, 06t6 1 . The constant- 14 Figure 2.3: Generic shape of the desired ipping speed, ! d (t), employed in the implementation of the proposed control method and stability analysis. This signal is composed of three sections: acceleration, constant speed and deceleration. According to this denition, g;1 , g;2 and g;3 denote the amounts of rotation accumulated angle during each of the three sections. The total ipping angle is equal to P 3 i=1 g;i = 2n, where n is the number of ips. The maximum allowed ipping speed is denoted by ! max . This reference ipping speed corresponds to the multi- ip phase in Figure. 2.2. speed section of ! d is simply ! d;2 (t) =! max ; (2.7) with a total duration 2 =! 1 max g;2 , in which g;2 is the desired rotation accumulated during the section. Lastly, the deceleration section is constructed as ! d;3 (t) = 3 3 (t tot t 3 ) 3 3 2 3 (t tot t) + 3 3 3 3 ; (2.8) in which 3 = ! 1 max g;3 , 3 = ( 3 =4) 3 3 ! max , t tot = 1 + 2 + 3 , and _ ! max;3 = ( 3 =4) 1 g;3 ! 2 max . In this denition, g;3 is the desired total rotation angle accumulated during the deceleration section;! max is the maximum allowed ipping speed; and _ ! max;3 is the maximum allowed ipping angular acceleration, which must remain smaller than the maximum angular acceleration that is physically attainable by the quadrotor about the ipping axis. The duration of this last section is 3 = 2 3 ; in (2.8), therefore, 15 Figure 2.4: Control scheme employed during the multi- ip phase of the aerobatic maneuver depicted in Fig. 2.2. In order to generate the control signal, the controller uses the angular-velocity and angular- acceleration control errors, and the reference torque, d , that is associated with the reference signals ! d and _ ! d . 1 + 2 6t6 1 + 2 + 3 . In summary, ! d (t) = 8 > > > > < > > > > : ! d;1 (t) if 0t 1 ! d;2 (t) if 1 <t P 2 i=1 i ! d;3 (t) if P 2 i=1 i <t P 3 i=1 i : (2.9) Note that only four parameters, g;1 , g;2 , g;3 and ! max , need to be specied for constructing ! d (t); there are two main restrictions on the specied values: the angular accelerations during the rst and third sections must not exceed the maximum angular acceleration that the quadrotor can achieve and the total nominal amount of time that is needed to complete the multi- ip phase, during which the vehicle is in free fall, must be minimized. 2.2.2 Control Strategy and Closed-Loop Scheme During regular ight, the quadrotor remains in a near-hovering state that is charac- terized by small attitude angle variations and slow angular velocities. In addition, as 16 explained in Chapter 2.1, the hovering dynamics of the system are describable by the linearized model that is specied by (2.4) and (2.5), in which the third row of (2.4) and all the rows of (2.5) can be viewed as independent Single-Input{Single-Output (SISO) systems, and the rst two rows of (2.4) can be viewed as SISO systems with inputs that are generated by the dynamics described by (2.5). Consistently, the quadrotor in regular ight can be own employing six SISO proportional-integral-derivative (PID) controllers. For example, the roll control torque signal is generated as 1 =k p ( des )k d _ +k i Z ( des )dt; (2.10) in which 1 is the rst component of the torque input , des denotes the roll angle reference and k p ;k d ;k i are PID gains determined through classical control tech- niques. The control signals for the pitch and yaw degrees of freedom are computed in an identical manner. Similarly, the thrust force control signal is computed as f =k pz (z des r z )k dz _ r z +k iz Z (z des r z )dt +mg; (2.11) in which r z is the third component of the position vector r, z des denotes the altitude reference and k p z ;k dz ;k iz are PID gains determined through classical control tech- niques. The translational motion on then 1 -n 2 plane is controlled by varying the desired pitch and roll angles. In addition, many globally stabilizing controller from the literature [15,16,18] can be employed for the regular ight mode. However, they normally require more computa- tional resources than the simple PID controllers based on the linearized model proposed above. And the highly nonlinear terms in the globally stabilizing controller may excite the unmodeled dynamics and cause potential instability during the DR phase due to the complicated aerodynamic disturbances encountered by the propellers in the descending. The compelling experiment results verify the eectiveness and robustness of the simple 17 PID controllers for the climb and DR phases in the aerobatic maneuvers. During aerobatic ight, the control signal is computed as (t) = K P + s s + K I e ! (t) + d (t); (2.12) where the bracket [ ] represents a lter that operates on the signal e ! (t); d is the desired torque;e ! =!! d is the angular velocity error; s is the dierential operator; and is a real coecient. The entire control scheme that is employed by the robot during the execution of a multi- ip aerobatic maneuver is shown in Figure. 2.4. Here, K P and K I are constant diagonal matrices that are dened in terms of the vectors k P 0 and k I 0 as K P = diag k P and K I = diag k I . As illustrated in Figure. 2.4, the control law in (2.12) has a simple linear structure that employs the angular velocity and acceleration control errors as feedback to compute the control signal(t). Feedback of the angular acceleration is included in the control scheme in order to compensate for angular velocity oscillations caused by rotor response latency, which increases the experimental robustness of the closed-loop system. In addition, the simple linear time-invariant (LTI) form of the lter in (2.12) reduces the probability of exciting the high-order unmodeled dynamics of the true system; this phenomenon is more likely to occur during aerobatic ight than in regular ight. Additionally, during a multi- ip, the magnitude of the total thrust force, f , is set to half its maximum nominal physical value to avoid motor saturation. The closed-loop system switches from the regular- ight controller, composed of PID SISO lters as exemplied by (2.10) and (2.11), to the aerobatic- ight controller, speci- ed by (2.12) and illustrated by the scheme in Figure. 2.4, at the end of the climb phase; then, the closed-loop system switches back from the aerobatic- ight controller to the regular- ight controller while in the DR phase. The switching times and duration of the multi- ip phase are programmed a priori and directly depend on the parameters 18 that dene! d (t); hence, the entire aerobatic maneuver is performed automatically and autonomously by the vehicle. In general, the stability of a switching system is dicult to establish and must not be ignored in experimental implementations. In this case, straightforward extensions of, and deductions from, the stability analysis presented in Chapter 2.2.3 show that, provided that certain conditions are satised, the switches do not destabilize the closed-loop system. 2.2.3 Stability of the Closed-Loop Angular Velocity Dynamics In regular ight, each degree of freedom of the quadrotor is associated with a SISO LTI system according to (2.4) and (2.5); therefore, the stability of the closed-loop can be straightforwardly guaranteed by placing the poles of each SISO transfer function in the left-hand side of the complex plane. Since this is a well-known procedure, the details are omitted. We mainly focus on the stability analysis of the closed-loop an- gular velocity and the attitude dynamics that result when the control law specied by (2.12) is employed to control aerobatic maneuvers according to the scheme in Fig. 2.4 and the globally dened nonlinear dynamics (2.2)(2.3). Based on this analysis, we propose a constructive method for controller synthesis and performance optimization; furthermore, we demonstrate in Chapter 2.2.5 that, if several conditions are satised, the transition between the regular- ight controller and the aerobatic- ight controller, and vice versa, does not destabilize the closed-loop system. Since multi- ip maneuvers require large attitude angle variations and high angular speeds, we consider the globally-dened nonlinear dynamic model given by (2.1), (2.2) and (2.3). To begin the analysis, note that the open-loop angular velocity dynamics in (2.3) can be rewritten as J _ ! +f(!) =; (2.13) 19 where f(!) =!J! = 2 6 6 6 6 4 (j 33 j 22 )! 2 ! 3 (j 11 j 33 )! 1 ! 3 (j 22 j 11 )! 1 ! 2 3 7 7 7 7 5 ; (2.14) in whichfj 11 ;j 22 ;j 33 g are the diagonal entries ofJ andf! 1 ;! 2 ;! 3 g are the components of!. Consistently, the reference model in Figure. 2.4 is expressed as J _ ! d +f(! d ) = d : (2.15) Additionally, the nonlinear termf(!) can be exactly expressed as f(!) =f(! d ) +A f (! d )(!! d ) +g(t;!); (2.16) whereA f (! d ) is the 3 3 Jacobian matrix evaluated along! d , i.e., A f (! d ) = @f @! ! d = 2 6 6 6 6 4 @f 1 @! 1 @f 1 @! 2 @f 1 @! 3 @f 2 @! 1 @f 2 @! 2 @f 2 @! 3 @f 3 @! 1 @f 3 @! 2 @f 3 @! 3 3 7 7 7 7 5 ! d = 2 6 6 6 6 4 0 (j 33 j 22 )! 3d (j 33 j 22 )! 2d (j 11 j 33 )! 3d 0 (j 11 j 33 )! 1d (j 22 j 11 )! 2d (j 22 j 11 )! 1d 0 3 7 7 7 7 5 = 2 6 6 6 6 4 0 (j 33 j 22 )a 3 (j 33 j 22 )a 2 (j 11 j 33 )a 3 0 (j 11 j 33 )a 1 (j 22 j 11 )a 2 (j 22 j 11 )a 1 0 3 7 7 7 7 5 ! d (t) =D f ! d (t); (2.17) and g(t;!) is the high-order nonlinear residual. Consistently,f! 1d ;! 2d ;! 3d g are the 20 components of! d . By subtracting (2.15) from (2.13), we obtain J _ e ! +A f (! d )e ! +g(t;!) = d : (2.18) Then, deningv(t) = nh 2 s+ K I i e ! o (t), rewriting the control law in (2.12) as d = (K P +K I )e ! +v; (2.19) _ v =v + 2 K I e ! ; (2.20) and substituting this new form of the control law into (2.18) yields _ e ! =J 1 (A f (! d ) +K P +K I )e ! +J 1 vJ 1 g(t;!); (2.21) _ v = 2 K I e ! v: (2.22) To continue, notice thatg(t;!) can be expressed in the compact form g(t;!) =f(!)f(! d )A f (! d )(!! d ) =e ! Je ! =h(e ! ); (2.23) which allows us to rewrite the system given by (2.21) and (2.22) as _ e ! =J 1 (A f (! d ) +K P +K I )e ! +J 1 K 1 zJ 1 h(e ! ); (2.24) _ z = 2 KK I e ! z: (2.25) In this representation, the state v has been rescaled to obtain the new state z = Kv, whereK is a positive denite diagonal matrix. By doing this, the entries of the new state z take values with similar orders of magnitude as those of e ! , which is a necessary feature for the implementation of the optimization algorithm to be discussed 21 in Chapter 3. Thus, by dening x(t) = e T ! (t) z T (t) T , we obtain the state-space representation _ x =(t;x); (2.26) where (t;x) = 2 6 4 J 1 (A f (! d ) +K P +K I ) J 1 K 1 2 KK I I 33 3 7 5 x + 2 6 4 J 1 h(e ! ) 0 31 3 7 5 : (2.27) Note that this is an exact representation of the closed-loop dynamics of the angular velocity control error, and clearly (t;x) is composed of a linear part, which comes from the linear terms of f(!) in (2.16), and a nonlinear part, which comes from the nonlinear residualg(t;!) in (2.16). In order to determine the conditions in terms ofK P ,K I and under which a xed point of (2.27) is stable, we dene the proto-Lyapunov function L(x) = 1 2 x T 2 6 4 J 0 33 0 33 P z 3 7 5 x = 1 2 e T ! Je ! + 1 2 z T P z z; (2.28) whereP z is a positive denite matrix. From (2.28), it immediately follows that 1 2 minf min (J); min (P z )gjxj 2 2 L(x) 1 2 maxf max (J); max (P z )gjxj 2 2 ; (2.29) where the operators min () and max () extract the smallest and largest eigenvalues, respectively, of a symmetric matrix. Also from (2.28), after taking the time derivative of L(x) along the trajectory dened by (2.26), we obtain _ L(x) =e T ! J _ e ! +z T P z _ z 22 =e T ! (A f (! d ) +K P +K I )e ! +e T ! K 1 z e T ! h(e ! ) +z T 2 P z KK I e ! z T P z z: (2.30) It can be shown that the angular velocity error,e ! , is orthogonal to the nonlinear residual, h(e ! ), i.e., e T ! h(e ! ) = 0. Thus, an exact simplied expression for _ L(x) is given by _ L(x) =e T ! A f (! d ) +K P +K I e ! z T P z z +z T 2 P z KK I +K T e ! : (2.31) Then, by dening A e! (t) =K P +K I + 1 =2 A f (! d ) +A T f (! d ) ; (2.32) A z =P z ; (2.33) A ze! = 2 P z KK I +K T ; (2.34) we obtain the compact expression _ L(x) =e T ! A e! (t)e ! z T A z z +z T A ze! e ! : (2.35) Now, to determine the stability conditions for the closed-loop system, we establish some linear matrix inequalities (LMIs) to be used in the rest of the analysis. First, we introduce the relationship z T A ze! e ! 6z T B z z +e T ! B e! e ! ; (2.36) 23 which is equivalent to 2 6 4 B z 1 2 A ze! 1 2 A T ze! B e! 3 7 5 0; (2.37) provided that the matricesB z andB e! are real symmetric. Then, assuming that (2.36) is satised, it follows that _ L(x)6e T ! (A e! (t)B e! )e ! z T (A z B z )z; (2.38) which implies that _ L(x) is negative denite if A e! (t)B e! 0;8t>t 0 ; (2.39) A z B z 0: (2.40) In (2.39),A e! (t) is time-varying because it depends on the reference angular velocity ! d (t). This dependency on time is removed from the inequality in (2.39) and from the remainder of the analysis by employing the following propositions. Proposition 2.1. The matrixA e! (t)2 convfA e!;1 ;A e!;2 g, whereA e!;1 =K P +K I andA e!;2 =K P +K I + 1 =2 D f +D T f ! max . Additionally, (2.39) is equivalent to the condition that bothA e!;1 B e! andA e!;2 B e! are positive denite. Proof. Note that the signal ! d (t)2 convf0;! max g and, using (2.17), the relationship in (2.32) can be expressed as A e! (t) =K P +K I + 1 =2 D f +D T f ! d (t). Therefore, A e! (t)2 convfA e!;1 ;A e!;2 g. Now, we prove the second statement of the proposi- tion. Suciency (`)' direction) follows immediately from the denition of ! d (t) in Chapter 2.2.1, according to which there are times, namely, t 1 and t 2 , during the exe- cution of a multi- ip at whichA e!;1 =A e! (t 1 ) andA e!;2 =A e! (t 2 ). Necessity (`(' direction) follows from the existence of a real that satises 0 6 6 1 such that 24 A e! (t)B e! =(A e!;1 B e! ) + (1)(A e!;2 B e! ) for all t>t 0 . At this point, it is important to remark that, directly from Proposition 2.1, the sucient conditions in (2.39) and (2.40) for _ L(x) to be negative denite are equivalent to A e!;1 B e! 0; (2.41) A e!;2 B e! 0; (2.42) A z B z 0: (2.43) In addition, Proposition 2.1 enables us to derive a time-independent negative upper bound for _ L(x), which clearly enforces the negative deniteness of _ L(x), if the inequal- ities in (2.41), (2.42) and (2.43) are satised. To construct this bound, we also employ Proposition 2.2, which is discussed next. Proposition 2.2. For any ! d (t) of the form that is dened in Chapter 2.2.1, the in- equality min t min A e! (t)B e! > min min (A e!;1 B e! ); min (A e!;2 B e! ) holds. Proof. SinceA e! (t)2 convfA e!;1 ;A e!;2 g (Proposition 2.1), it follows that min t min A e! (t) B e! = min min (A e!;1 B e! ) + (1)(A e!;2 B e! ) , where 0 6 6 1. Then, by applying Weyl's Inequality Theorem [34], it immediately follows that min min A e!;1 B e! + 1 A e!;2 B e! > min min A e!;1 B e! + 1 min A e!;2 B e! > min min A e!;1 B e! ; min A e!;2 B e! ; (2.44) which completes the proof. 25 From (2.38), if the inequalities in (2.41), (2.42) and (2.43) are satised, it follows that _ L(x)6min min t min A e! (t)B e! ; min A z B z jxj 2 2 : (2.45) Then, by employing Proposition 2.2, we obtain the time-independent negative upper bound _ L(x)6min min A e!;1 B e! ; min A e!;2 B e! ; min A z B z jxj 2 2 ; (2.46) which holds true provided that the inequalities in (2.41), (2.42) and (2.43) are enforced. We have, thus far, devised a method for systematically constructing a Lyapunov function to determine the stability of an equilibrium point of the closed-loop dynamics described by (2.26). In the remainder of this chapter, taking advantage of the param- eterized and constructive nature of the followed approach, we propose a method for controller synthesis, i.e., we establish a set of criteria for selecting! d (t),K P ,K I and for the scheme in Figure. 2.4 such that the closed-loop angular velocity error dynamics remain stabilized. The next step is to employ Lyapunov's direct method [35, 36] to show that the unique equilibrium point of the system in (2.26) is exponentially stable ifP z 0 and the conditions that are specied by (2.37), (2.41), (2.42) and (2.43) are enforced simultaneously. The proof presented here is a general version of that in [37]. Formally, we state the main results relating to the stability analysis and controller synthesis using an additional set of propositions. Proposition 2.3. For a reference! d with the generic form dened in Chapter 2.2.1, the nonautonomous system _ x = (t;x), dened by (2.26) and (2.27), has a unique equilibrium point, namely,x ? = 0. 26 Proof. To begin the argument, notice that for a xed pointx ? it follows that(t;x ? ) = 0 for all t, which implies that 0 = (K P +K I +A f (! d ))e ? ! +K 1 z ? h(e ? ! ); (2.47) 0 =z ? + 2 KK I e ? ! : (2.48) Then, combining both equations, we can write (K P +A f (! d ))e ? ! +h(e ? ! ) = 0: (2.49) Because e ? ! is constant by the denition of xed point, it follows that h(e ! ) is also constant. Consequently, (2.49) is satised if and only if e ? ! = 0. The proof of the `(' direction immediately follows from noticing that ife ? ! = 0, thenh(e ? ! ) = 0. The proof of the `)' direction follows from the use of the contradiction technique by noticing that (2.49) must be satised for any feasible time-varying! d . However,e ? ! cannot be both time-varying and constant at the same time (a contradiction), which completes the argument. Note that by design, if we enforce the global exponential stability of the unique equilibrium point x ? = 0, then the angular velocity control error approaches 0 from any conceivable initial state. Proposition 2.4. For a reference ! d with the generic form that is dened in Chap- ter 2.2.1 and allr e! 2R + , the function(t;x) : [0;1)R 6 !R 6 is continuous int and uniformly Lipschitz continuous inx on [0;1)D x , whereD x =fx :je ! j 2 r e! ;z2R 3 g. Proof. Here, we prove that(t;x) is locally Lipschitz. To begin, we rewrite the non- autonomous system as _ x =(t;x) =F (t)x +G(x); (2.50) 27 where F (t) = 2 6 4 J 1 (A f (! d ) +K P +K I ) J 1 K 1 2 KK I I 33 3 7 5 ; (2.51) G(x) = 2 6 4 J 1 h(e ! ) 0 31 3 7 5 : (2.52) InF (t),A f (! d ) is a function of! d which is a polynomial oft without singularities. Therefore,(t;x) is continuous. Then, assuming two feasible statesx andy, we get j(t;x)(t;y)j 2 kF (t)k 2 jxyj 2 +jG(x)G(y)j 2 ; (2.53) where the induced 2-normkF (t)k 2 can remain bounded by proper selecting the reference ! d (t). The remaining step is to nd a Lipschitz constant ofjG(x)G(y)j 2 . Noticing that the last three entries of G(x) are zero and J 1 is constant, equivalently we can nd a Lipschitz constant forh(e ! ). To do this, we rst nd the partial derivative of J 1 h(e ! ) and get J 1 @h @e ! = 2 6 6 6 6 4 0 32 j 11 e ! 3 32 j 11 e ! 2 13 j 22 e ! 3 0 13 j 22 e ! 1 21 j 33 e ! 2 21 j 33 e ! 1 0 3 7 7 7 7 5 ; (2.54) where 32 = j 33 j 22 , 13 = j 11 j 33 and 21 = j 22 j 11 . Then, in the domain D x =fx : je ! j 2 r e! ;z2R 3 g, J 1 @h =@e! can be chosen to be continuous and the induced1-norm J 1 @h =@e! 1 to be bounded as J 1 @h @e ! 1 p 2 max 32 j 11 ; 13 j 22 ; 21 j 33 r e! ; (2.55) 28 and J 1 @h @e ! 2 p 3 J 1 @h @e ! 1 (2.56) p 6 max 32 j 11 ; 13 j 22 ; 21 j 33 r e! =L g : (2.57) This allows us to write j(t;x)(t;y)j 2 max tt 0 (kF (t)k 2 ) +L g jxyj 2 : (2.58) Therefore,(t;x) is locally Lipschitz. Remark 2.1: We employ the word uniformly to indicate that (t;x) satises the Lipschitz condition on D x with the same constant ` 0 , i.e.,j(t;x)(t;y)j 2 6 ` 0 jxyj 2 . From Theorem 2.25 in [38], we can deduce that for any initial condition x(t 0 )2 D x , a unique trajectory exists and is guaranteed to continue to exist and re- main unique as long as it stays inside D x . In this case, the bound r e! can be made arbitrarily large in order to generate a larger associated ` 0 , which implies that(t;x) is locally Lipschitz continuous on [0;1)R 6 according to the denition in [36]. Thus, even though(t;x) is not globally Lipschitz continuous, from Theorem 2.25 in [38], it can be shown that a unique trajectory starts from each nite x(t 0 )2 R 6 . Therefore, for all engineering purposes, the results that are presented here are global. Proposition 2.5. The unique equilibrium point,x ? = 0, of the nonautonomous system _ x = (t;x), that is dened by (2.26) and (2.27), is exponentially stable on R 6 if the parameter setfK P ;K I ;;! max ;P z ;B z ;B e! ;Kg satises the conditions that are specied by (2.37), (2.41), (2.42) and (2.43) andP z 0. Proof. Acknowledging Proposition 2.3, Proposition 2.4 (and the associated Remark 1), (2.29) and (2.46), the proof follows directly from the application of Lyapunov's Direct Theorem (Theorem 4.10 in [36]). 29 According to Proposition 2.5, after dening 1 = 1 2 min min J ; min P z ; (2.59) 2 = 1 2 max max J ; max P z ; (2.60) 3 = min min A e!;1 B e! ; min A e!;2 B e! ; min A z B z ; (2.61) it is straightforward to conclude that jx(t)j 2 6 r 2 1 jx(t 0 )j 2 e 3 2 2 (tt 0 ) : (2.62) Then, settingv(t 0 ) =Kz(t 0 ) = 0 yieldsjx(t 0 )j 2 =je ! (t 0 )j 2 and (2.62) becomes je ! (t)j 2 6 r 2 1 je ! (t 0 )j 2 e 3 2 2 (tt 0 ) : (2.63) Remark 2.2: For a given parameter set K P , K I , , ! max , K , the problem of nding matrix parameters P z ,B z ,B e! satisfyingP z 0 and the stability condi- tions given by (2.37), (2.41), (2.42) and (2.43) is essentially a semidenite programming problem that is solvable via interior-point methods [31]. In this case, both feasible con- trol gains and the corresponding parameters used in the stability analysis can also be selected heuristically, and empirical evidence suggests that the stabilizing parameter space is large because feasible (suboptimal) controllers can be easily identied through tuning experiments, as demonstrated in [33,37]. Remark 2.3: Because (t;x) is only locally Lipschitz continuous inx on [0;1) R 6 , and not globally Lipschitz continuous, we can not directly use a global stability argument (such as Theorem 2.33 in [38]). However, Remark 1 allows us to use Theo- rem 4.10 in [36] to extract the global conclusion in Proposition 5 (i.e.,x ? = 0 is globally exponentially stable). This conclusion is consistent with the observation that the re- gion of attraction associated with x ? = 0 can be made arbitrarily large by choosing 30 a bound r e! p 2 = 1 jx(t 0 )j 2 for any initial condition x(t 0 ). This conclusion follows from noticing that for any nite initial conditionx(t 0 ), there exists a constant such that (t;x) is uniformly Lipschitz continuous inside the corresponding D x domain. Since this constant becomes larger for larger values ofjx(t 0 )j 2 , global Lipschitz continuity is not satised. However, for all engineering purposes, the closed-loop stability of the system is global because from any nite initial statex(t 0 ), the trajectory of the system will approach x ? = 0 at an exponential rate and remain inside of D x provided that r e! p 2 = 1 jx(t 0 )j 2 . 2.2.4 Closed-Loop Attitude Analysis Thus far, we have proved that the angular-velocity control error during a multi- ip aerobatic maneuver, e ! (t), converges to zero at an exponential rate. The next step in our analysis is to determine how the attitude kinematics in (2.2) evolves according to the scheme in Fig. 2.4, in which the closed-loop !(t) is the only signal that drives the attitude dynamics of the system. In this chapter, we show that the attitude error between the true body frameB and the ideal body frameI 1 during the execution of a multi- ip maneuver remains bounded. The kinematics of the true body frameB are governed by _ q = 1 2 q p; with p = 0 ! T T ; (2.64) and the kinematics of the ideal body frameI are determined by _ q d = 1 2 q d ^ p d ; with ^ p d = 0 ^ ! T d T ; (2.65) where the quaternion q d denotes the attitude of the ideal body frameI. Because the 1 Corresponding to the desired trajectory of the system in the absence of angular-velocity control error and initial attitude error. 31 components of! d are expressed inB, and not inI, the kinematics of q d can not be directly expressed in terms of ! d in the same way as q is expressed in terms of ! in (2.64). Therefore, ^ ! d is used in (2.65), and ^ ! d has exactly the same components as! d but is expressed inI rather thanB. It follows that the attitude error quaternion q e = m e n T e T between q d and q is given by q e = q 1 d q; (2.66) as can be shown using the methods in [39]. The attitude error quaternion q e contains the information about the rotation axis and the rotation angle fromB toI. In this paper, the attitude error is described by the rotation angle fromB toI, which is computable from the scalar part of q e , namely, m e . Thus, by dierentiating the attitude error quaternion q e , specied by (2.66) along the trajectories of (2.64) and (2.65), we obtain that _ q e = 2 6 4 _ m e _ n e 3 7 5 = 1 2 2 6 4 n T e (! ^ ! d ) m e (! ^ ! d ) +n e (! + ^ ! d ) 3 7 5 : (2.67) Now, by denoting the rotation angle fromB toI as , and noticing that the corre- sponding rotation axis is specied by the unit vectora e with its components expressed in bothB andI, we obtain that m e = cos 1 2 ; n e = sin 1 2 a e ; (2.68) which directly follows from the denition of the attitude error quaternion. Thus, from (2.67) and (2.68) it follows that _ m e = sin 1 2 1 2 _ = 1 2 sin 1 2 a T e (! ^ ! d ); (2.69) 32 which, with the assumption that 6= 0, further yields that _ =a T e (! ^ ! d ) =a T e (!! d ) =a T e e ! : (2.70) Note that to obtain this expression, we use the property that ^ ! d has the same com- ponents as! d , even though they are expressed in a dierent frame of coordinates. To continue with the analysis, using (2.70) and recalling thatja e j 2 = 1, we establish the inequality: _ 6 _ = a T e e ! 6je ! j 2 : (2.71) Then, recalling from the previous chapter thate ! converges to zero at an exponential rate during a multi- ip, integrating both sides of (2.71) and substituting (2.63) into the integration result, it yields that (t)6 (t 0 ) + 2 2 p 2 3 p 1 1e 3 2 2 (tt 0 ) je ! (t 0 )j 2 ; (2.72) where (t 0 ) is the attitude error angle at the beginning of the maneuver's multi- ip phase. Note that according to the result specied by (2.72), the attitude error angle, (t), remains bounded during the entire aerobatic maneuver, provided that the xed point of the system specied by (2.26) and (2.27) is enforced to be exponentially stable by employing the control method in Chapter 2.2.3; the upper bound of the attitude error angle is also specied in (2.72). In fact, it is straightforward to see that the attitude error angle (t) is bounded uniformly by (t 0 ) + p 4 3 2= 2 3 1 je ! (t 0 )j 2 , where the initial attitude error angle and angular velocity error, (t 0 ) ande ! (t 0 ), can be regulated to zero by the regular- ight controller at the end of the climb phase. This fact allows us to conclude that the attitude error angle is bounded by a small value. Hence, theoretically, the quadrotor can perform an unlimited number of high-speed ips while the attitude 33 error remains bounded by a small value. Note, furthermore, that (t) is always positive (by denition) and that it converges to a nite values as t!1. This last statement immediately follows from applying Theorem 12.33 in [40] to the relationships specied by (2.63) and (2.71). 2.2.5 The Eects of Controller Switching and Stability The control scheme switches only once from the regular- ight controller to the aerobatic- ight controller and only once from the aerobatic- ight controller to the regular- ight controller during the entire multi- ip maneuver process. The rst switching occurs at the end of the climb phase (see Chapter 2.1 and Fig. 2.2) and the second switching occurs during the DR phase. In each case, the nal state of the pre-transition ight mode can be viewed as the initial condition of the post-transition ight mode. Thus, the theoretical stability of the rst transition is immediately guaranteed by the analyses in Chapter 2.2.3 and Chapter 2.2.4 as for any nite initial angular velocity!(t 0 ) at the start of the multi- ip phase, the angular velocity errore ! (t) will approach zero at an exponential rate and for any nonzero initial attitude and angular velocity, the attitude error will remain bounded according to (2.72). Furthermore, (2.72) establishes that (t) can be made arbitrarily small by making the initial variables,e ! (t 0 ) and (t 0 ), small. In addition, sincee ! (t 0 ) and (t 0 ) cor- respond to the nal state of the climb phase, during which the regular- ight controller ies the robot, both the theoretical stability and the experimental high performance of the closed-loop system are enforced because the experimental values of e ! (t 0 ) and (t 0 ) are always approximately zero by design. Note that a practical consequence of guaranteeing the stability and enforcing the high performance of the rst transition is that the values of the pre-transition and post-transition control signals are very close to zero and very similar to each other. 34 The stability and performance of the second switching can be enforced through theoretical and experimental considerations. According to the analysis in Chapter 2.2.3 and Chapter 2.2.4, at the end of the DR-1 phase, the angular velocity error e ! (t 1 ) is approximately zero, becausee ! (t) is attracted to zero at an exponential rate; and the attitude error remains bounded according to (2.72) and is enforced to be approximately zero, i.e., (t 1 ) 0, because bothe ! (t 0 ) and (t 0 ) are stabilized to approximately zero at the end of the climb phase by the regular- ight controller. Therefore, at the time of the second switching, the initial state of the closed-loop system that is associated with the regular- ight controller is approximately zero, and consequently, inside the region of attraction. Thus, the attitude stability of the closed-loop system after the second controller switching is guaranteed. Similar to the rst switching case, the control signal variation induced by the second transition is approximately zero by design. Hundreds of experiments performed at the USC Autonomous Microrobotic Systems Laboratory (AMSL) support the validity of the presented analysis (see [33,37] and Chapter 2.3). 2.3 Simulation and Experimental Results 2.3.1 The Flying Vehicle and Experimental Architecture The aerial robot used for the implementation of the linear controller design for aerobatic ights is Crazy ie 1.0 2 as shown in Fig. 2.1, which has a propeller-tip-to-propeller-tip distance of 13 cm and weighs 19 g including the battery. The maximum thrust-to-weight ratio is about 1.7. The micro size of the robot allows it to be more agile than a larger quadrotor [41]. This robot carries several on-board sensors used in the implementation of controllers, including a 3-axis accelerometer, a 3-axis gyroscope and a barometer. The barometer measures the altitude of the robot while the accelerometer and gyroscope are used to estimate the robot's attitude and generate the angular velocity feedback. 2 https://wiki.bitcraze.io/projects:crazy ie:index 35 The measurement range of the gyroscope is2 10 3 deg s 1 , which constrains the maximum reference angular speed in the speed planning. Figure 2.5: Architecture of the experimental setup. The ground computer and joystick are mainly used to send references of attitude and altitude to the quadrotor and compensate drift during regular ight. The ground computer is also receives the real-time quadrotor data for o-line analysis. The whole experiment architecture consists of a ground computer and the aerial robot as shown in Fig. 2.5. The ground computer is used to send the reference signal in regular ight mode and the multi- ip trigger command to the quadrotor. The algo- rithms that control the vehicle's attitude and altitude run onboard using information from the onboard sensors. During regular ight, motion on the n 1 -n 2 plane is con- trolled by manually varying the pitch and roll references. The multi- ip maneuver is triggered from the ground computer, and the quadrotor executes the multi- ip maneu- vers completely autonomously. It relies entirely on onboard sensing, signal processing and control, and no references or commands are sent from the ground computer during the maneuvers. Simultaneously, the ight information can be communicated from the yer back to the ground computer in real time at a frequency up to 100 Hz for mon- itoring and recording. More details about the robot and the experiment architecture can be found in [33] [37]. 2.3.2 Generalized Multi-Flip and Generalized Flipping Angle The multi- ip maneuver is conventionally dened as a 2n angle rotation about the b 1 (roll) or b 2 (pitch) principal axis of the yer [21] [42] [43], where n indicates the 36 number of the consecutive ips. In this paper, we generalize the concept of the multi- ip by choosing the ipping axis to be any body-xed axis passing through the yer's center of mass. The generalized ipping angle is dened as the angle with its vertex in the ipping axis, measured in the plane perpendicular to the ipping axis, generated between the horizontal plane of the inertial frame (n 1 -n 2 -plane) and the line created by the intersection of theb 1 -b 2 -plane with the plane perpendicular to the ipping axis (intersection line for short). The direction of the ipping axis with respect to the inertial frame determines the sign of the generalized ipping angle following the right- hand principle. This denition is similar to that of the roll angle in the Euler Z-Y-X convention. As a consequence of this similarity, we can rotate the axis b 1 to match the ipping axis and the axis b 2 to match the line created by the intersection of the b 1 -b 2 -plane with the plane perpendicular to the ipping axis. And the new Euler roll angle after this rotation is equal to the generalized ipping angle. The calculation of the generalized ipping angle follows the same idea. The desired multi- ip axis is denoted by a unit vector a b f = [a 1 a 2 a 3 ] T that is xed in body frame with its components expressed inB, and a 2 1 +a 2 2 +a 2 3 = 1. We use the quaternion q to represent the attitude of the vehicle, and the quaternion q u to represent the rotation of b 1 and b 2 to coincide with the generalized ipping axis and the intersection line, respectively. The components of q u depend on the components of a b f . The new attitude after the rotation q u is given by q 0 = q q u = q 0 0 q 0 1 q 0 2 q 0 3 ; (2.73) where denotes the quaternion product. Using the conversion between Euler angles and quaternion, the generalized ipping angle is g = arctan 2(q 0 0 q 0 1 +q 0 2 q 0 3 ) q 02 0 q 02 1 q 02 2 +q 02 3 : (2.74) 37 Note that the generalized ipping angle depends only on the current attitude q and the multi- ip axisa b f . 2.3.3 Simulation Results In this chapter, the simulation results that validate the linear controller design (2.12) with the system dynamics (2.2)(2.3) are presented. Six high-speed maneuvers are demonstrated including the single-, double- and triple- ip maneuvers about the body axisb 1 and a 45 -oblique axis,b 0 1 = (b 1 +b 2 ) = p 2. Fig. 2.6 shows the simulation results of six dierent multi- ip cases, and the results show that the angular velocities converge to the reference signals quickly, which is consistent with the analytical conclusion of the exponential stability. 2.3.4 Experimental Results In this chapter, we present experimental results that validate the proposed linear con- troller design (2.12). We consider six high-speed reference trajectories! d , correspond- ing to single-, double- and triple- ip maneuvers about the body axisb 1 and a 45 -oblique axis,b 0 1 = (b 1 +b 2 ) = p 2. It can be veried that there exist parameter setfP z ;B z ;B e! ;Kg satisfying the conditions stated in Proposition 2.5 for the corresponding control gains and reference parametersfK P ;K I ;;! max g used in the experiments such that the ex- ponential stability is achieved. Note that those parameters used for the experiments shown in this subsection are tuned manually through experiences. Experimental results of the six maneuver ights employing the linear controller design (2.12) are shown in Figs. 2.7. The reference ipping speed,j! d (t)j 2 , is shown in green, the measured controlled ipping speed, ! 1 (t), is shown in blue, the measured controlled ipping speed,! 2 (t), is shown in cyan, the measured ipping angle aboutb 1 and generalized ipping angle are shown in red, and the desired ipping speed about b 1 andb 2 axes is shown in magenta. The proposed linear controller design is veried in 38 (a) Simulation Single Flip (b 1 ). (d) Simulation Single Flip (b 0 1 ). (b) Simulation Double Flip (b 1 ). (e) Simulation Double Flip (b 0 1 ). (c) Simulation Triple Flip (b 1 ). (f) Simulation Triple Flip (b 0 1 ). Figure 2.6: Simulation results. Plots (a), (b) & (c) show simulation results of single-, double- and triple- ip high-speed maneuvers about the body-xed axis b 1 employing the linear controller design. Plots (d), (e) & (f) show the simulations of single-, double- and triple- ip high-speed maneuvers about the body-xed axisb 0 1 employing the linear controller design. 39 (a) Experimental Single Flip (b 1 ). (d) Experimental Single Flip (b 0 1 ). (b) Experimental Double Flip (b 1 ). (e) Experimental Double Flip (b 0 1 ). (c) Experimental Triple Flip (b 1 ). (f) Experimental Triple Flip (b 0 1 ). Figure 2.7: Experimental results. Plots (a), (b) & (c) show experimental results of single-, double- and triple- ip high-speed maneuvers about the body-xed axis b 1 employing the linear controller design. Plots (d), (e) & (f) show real-time experiments of single-, double- and triple- ip high-speed maneuvers about the body-xed axisb 0 1 employing the linear controller design. 40 (a) Single Flip (b 1 ). (b) Double Flip (b 1 ). (c) Triple Flip (b 1 ). (d) Single Flip (45°). (e) Double Flip (45°). (f) Triple Flip (45°). Figure 2.8: Real-time ight experiments. Composite images of six high-speed maneuvers. Trajectories (a), (b) & (c) correspond to single-, double- and triple- ip maneuvers about the body-xed axis b 1 . Trajectories (d), (e) & (f) correspond to single-, double- and triple- ip maneuvers about the 45 -oblique axisb 0 1 = (b1+b2) = p 2. 41 experiments to enable the quadrotor to track the high-speed angular velocity references in six dierent maneuver cases. Additionally, the vehicle is able to restabilize its attitude quickly without crashing or stalling after the maneuver is nished. Here, Figs. 2.8-(a), 2.8-(b) and 2.8-(c) show time lapse images of single-, double- and triple- ip maneuvers about theb 1 axis, respectively. Similarly, Figs. 2.8-(d), 2.8-(e) and 2.8-(f) show time lapse images of single-, double- and triple- ip maneuvers about a 45 - oblique axis,b 0 1 = (b 1 +b 2 ) = p 2, respectively. The complete video of the six experimental multi- ip cases is available at http://www.uscamsl.com/resources/MultiFlipACC.mp4 Compared with the simulation results in Fig. 2.6, the actual angular velocities in the experiments of six multi- ip maneuvers shown in Fig. 2.7 are distorted by some undesired oscillations and overshoots. This undesired phenomenon can not be explained by the nominal dynamics (2.2)(2.3), otherwise similar oscillations could be replicated in the simulations. This performance dierence between the simulation and experimental results indicate some unmodeled dynamics that result into these undesired phenomenon in the aerobatic maneuvers. In the later chapters of nonlinear adaptive control, an actuator dynamics model will be introduced to describe this unmodeled dynamics, and a nonlinear adaptive control strategy is proposed to compensate for its eects. 2.3.5 Comparison with a State-of-the-Art Controller We compare the proposed control method with a state-of-the-art geometric controller [1], in terms of their capability of performing high-speed multi- ip aerobatic maneuvers. Through ight experiments, both sets of control gains are tuned with the objective of achieving robust stability and high performance during the execution of consecutive triple- ips about the b 1 axis. As examples of typical results, Fig. 2.9-(a) shows the signalsj!(t)j 2 andj! d (t)j 2 obtained while the quadrotor is own by the geometric controller and Fig. 2.9-(b) shows the signalsj!(t)j 2 andj! d (t)j 2 obtained while the quadrotor is own by the proposed scheme in Fig. 2.4. 42 (a) Consecutive Triple Flips (b 1 ) { Ge- ometric Control Method [1]. (b) Consecutive Triple Flips (b 1 ) { Pro- posed Control Method. Figure 2.9: Comparison of the proposed control method with the state-of-the-art geometric controller in [1]. Both controllers are experimentally implemented to perform consecutive triple- ips about theb 1 axis. (a) Desired and measured angular rates of the quadrotor while own by the geometric controller. In this case, the quadrotor can nish only three consecutive triple- ips before crashing to the ground; in fact, the rst impact occurs during the DR phase of the second triple- ip. (b) Desired and measured angular rates of the quadrotor while own by the control scheme in Fig. 2.4. In this case, the quadrotor can indenitely complete consecutive triple- ips; this plot shows 13 successfully-completed aerobatic maneuvers which demonstrate both the high performance and stability robustness (against noise and disturbances) of the proposed control method when compared to that presented in [1]. In the rst case (Fig. 2.9-(a)), the yer can complete only three consecutive triple- ips; then, it crashes to the ground, which in the plot corresponds to the severe oscil- lations of the angular velocity. In contrast, in the second case (Fig. 2.9-(b)), the yer is able to perform consecutive triple- ips indenitely (thirteen are shown in the plot) until the power depletes; hence, the proposed control method exhibits both better per- formance and stability robustness than the geometric control method in [1] during the execution of triple- ip aerobatic maneuvers. To evaluate the multi- ip performance, we dene the Performance Figure of Merit (PFM) to be exp = Z t f t 0 je ! ()j 2 d; (2.75) wheret 0 andt f are the start and end time of the maneuver respectively, and this PFM denition plays an crucial role in the performance optimization that will be discussed later. 43 Consistently, the average value of exp corresponding to the experiment in Fig. 2.9- (a) is 4:4162 while the average value of exp in Fig. 2.9-(b) is 1:4372. To explain this observation, we hypothesize that the highly nonlinear terms in the geometric control method, such as those employed for coordinate conversion, may excite the unmod- eled dynamics of the system, thereby rendering it more sensitive to noise and distur- bances during extremely fast aerobatic ight. In addition, the proposed control scheme (Fig. 2.4) uses angular acceleration for feedback to inhibit oscillations caused by la- tency, which is not properly accounted for in the geometric control case. Fig. 2.9-(b) demonstrates the signicant performance robustness of the proposed linear controller design for the aerobatic maneuver even in the presence of the plant uncertainties and disturbances. 44 Chapter 3 Performance Optimization of Aerobatic Flight In Chapter 2, we proposed the linear controller synthesis for the high-speed aerobatic ights, in which the linear controller design is proved to stabilize the closed-loop system provided that several stability conditions are satised. Among all the feasible control parameters that guarantee the stability of the system, is there a way to locate the parameter set that can optimize the aerobatic ight performance? In this chapter, we will answer this question by employing the convex optimization approach. Note that the optimization strategy proposed in this chapter is based on the linear controller design and stability analysis result in Chapter 2. 3.1 Formulation of the Optimization Problem In this chapter, we introduce an optimization-based approach for controller synthesis with the objective of improving the ight performance of the quadrotor during the exe- cution of a multi- ip aerobatic maneuver. In particular, we propose a method for nding a parameter set, namely, K P ,K I , , ! max ,P z ,B z ,B e! ,K , that approximates 45 a solution that globally minimizes the performance gure of merit (PFM) = Z t f t 0 je ! ()j 2 d: (3.1) This criterion is a metric of the accumulated angular-velocity control errore ! (t), during a multi- ip that starts att 0 and ends att f , with respect to zero. To begin the formulation process, we combine (2.63) and (3.1) to obtain 6 2 2 3 r 2 1 1e 3 2 2 (t f t 0 ) je ! (t 0 )j 2 = bound ; (3.2) in which bound is a monotonically decreasing function with respect to 3 1 2 and 1 1 2 ; therefore, maximizing 3 1 2 and 1 1 2 (or a positive linear combination of the two terms) leads to the minimization (or reduction) of bound . Accordingly, if bound is a tight upper bound of , is also reduced. In this approach, however, since the minimization algorithm is not directly applied on but on bound , the optimization process may not lead to a signicant reduction of when the value of bound is excessively conservative. To prevent this from happening, the trajectory of the exponential bound p 2 = 1 je ! (t 0 )j 2 e 3 2 2 (tt 0 ) must stay suciently close to the trajectory ofje ! (t)j 2 for allt2 [t 0 :t f ]. Here, we enforce this requirement by rescalingv(t) in (2.25) to denez(t), which has a similar order of magnitude to that of e ! (t). Heuristic observations indicate that, for the system considered in this work and selected approach, state rescaling is a prerequisite for the optimization to be eective. In addition, to reduce the complexity of the optimization problem formulation, we deneK =J 1 1 ; 2R + . To continue, note that according to (2.59), (2.60) and (2.61), 1 , 2 and 3 are functions ofx op = k P ,k I ,,! max ,P z ,B z ,B e! , . Consistently, to avoid confusion with parametersf 1 ; 2 ; 3 g obtained heuristically or dened in terms of other variables, 46 we explicitly set the auxiliary identities 1 := 1 (x op ); 2 := 2 (x op ); 3 := 3 (x op ); (3.3) and dene the objective function to be minimized as f 0 (x op ) = 2 (x op ) 1 (x op ) + g 0 3 (x op ) ; (3.4) where g 0 2R + is a constant weighting factor. Unfortunately, the inequalities specied by (2.37), (2.41), (2.42), (2.43) andP z 0 do not dene a convex set in the x op space because matrices A ze! , A e!;1 and A e!;2 contain nonconvex terms. To deal with this complication, we solve the problem in two steps. First, we reduce the dimension of the optimization variable by replacingx op with z op = k P ;! max ;P z ;B z ;B e! and setting k I ;; to constant values. Then, the reformulated problem, which happens to be quasiconvex, is repeatedly solved for a large set of dierent systematically-chosen combinations of parameters k I ;; . In order to nd a solution that is guaranteed to be globally optimal for the selected parameter space, we choose the values of k I ;; to be evenly spaced within their ranges; then, we traverse through all the dened combinations of k I ;; and determine for each of them the solution that minimizes the expression specied by (3.4). Consistently with this proposed approach, for the purposes of analysis, we explicitly set the auxiliary identities 1 :=% 1 (z op ); 2 :=% 2 (z op ); 3 :=% 3 (z op ); (3.5) and dene a new objective function to be minimized, which is associated with a par- 47 ticular combination of constant parameters k I ;; , as g 0 (z op ) = % 2 (z op ) % 1 (z op ) + g 0 % 3 (z op ) : (3.6) Note that it is not dicult to show that g 0 (z op ), with the optimization variablez op , is a quasiconvex function. This fact is formally stated and proved below. Proposition 3.1. IfP z 0, and the inequalities specied by (2.37), (2.41), (2.42) and (2.43) are satised, the function% 2 (z op ) :R 3 + R + S 3 + S 3 S 3 !R + is convex, the functions% 1 (z op ) :R 3 + R + S 3 + S 3 S 3 !R + and% 3 (z op ) :R 3 + R + S 3 + S 3 S 3 ! R + are concave, and the functiong 0 (z op ): R 3 + R + S 3 + S 3 S 3 !R + is quasiconvex. Proof. First, we show that the operator that extracts the minimum eigenvalue of a symmetric matrix, min (), is concave and the operator that extracts the maxi- mum eigenvalue of a symmetric matrix, max (), is convex. The proof follows directly from Weyl's Inequality Theorem [34]. Namely, min (A + (1)B) > min (A) + (1) min (B) and max (A + (1)B)6 max (A) + (1) max (B), where 06 6 1,A2S k andB2S k . It is straightforward to see that diag J;P z ,A e!;1 B e! ,A e!;2 B e! andA z B z are ane functions of the optimization variable z op . Also, it has been shown that the pointwise minimum of a set of concave functions is concave and the pointwise maximum of a set of convex functions is convex [31]. Additionally, in this particular case, the domains of the functions % 1 ;% 2 ;% 3 and g 0 are convex. Therefore, % 1 (z op ) and % 3 (z op ) are concave functions, and % 2 (z op ) is a convex function. Now, we prove that g 0 (z op ) is a quasiconvex function. It is straightforward to see that the domain of g 0 is convex. Then, we only need to show that all the sublevel sets of g 0 , namely, S s = z op 2 domg 0 jg 0 (z op )6 s ; s 2R; (3.7) 48 are convex [31]. If s is smaller than the minimum value of g 0 , S s is empty and, therefore, convex. Next, we consider the case in which S s is non-empty. Given that P z 0, and that the inequalities specied by (2.37), (2.41), (2.42) and (2.43) are satised, it follows that% 1 ;% 2 ;% 3 andg 0 are positive. Then, for two arbitrary points x 1 ;x 2 2S s , it follows that % 2 (x 1 ) % 1 (x 1 ) + g 0 % 3 (x 1 ) 6 s ; % 2 (x 2 ) % 1 (x 2 ) + g 0 % 3 (x 2 ) 6 s : Since % 1 and % 3 are concave, and % 2 is convex, we have that % 1 (x 1 + (1)x 2 ) + g 0 % 3 (x 1 + (1)x 2 ) >% 1 (x 1 ) + (1)% 1 (x 2 ) + g 0 (% 3 (x 1 ) + (1)% 3 (x 2 )) = (% 1 (x 1 ) + g 0 % 3 (x 1 )) + (1) (% 1 (x 2 ) + g 0 % 3 (x 2 )) > 1 s (% 2 (x 1 ) + (1)% 2 (x 2 ))> 1 s % 2 (x 1 + (1)x 2 ); where 066 1. Then, it follows that % 2 (x 1 + (1)x 2 ) % 1 (x 1 + (1)x 2 ) + g 0 % 3 (x 1 + (1)x 2 ) 6 s : Thus, all the sublevel sets of g 0 are convex; hence, g 0 (z op ) is a quasiconvex function. Since the objective function g 0 (z op ) is quasiconvex, P z 0, and the constraints specied by (2.37), (2.41), (2.42) and (2.43) are ane, the optimization problem is quasiconvex with respect to the variablez op = k P ;! max ;P z ;B z ;B e! and, therefore, solvable. Formally, we formulate the problem as Minimize zop % 2 (z op ) % 1 (z op ) + g 0 % 3 (z op ) 49 subject to P z 0; (3.8) 2 6 4 B z 1 2 A ze! 1 2 A T ze! B e! 3 7 5 0; A e!;1 B e! 0; A e!;2 B e! 0; A z B z 0; % 2 b s % 1 6 0; u k P k P l k P ; u !max >! max >l !max ; e !max k P +e max k I 6 max ; in which b s is an upper bound of % 2 % 1 1 , introduced to reduce the amount of com- putation without excluding the optimal solution from the feasible set; u k P ,l k P are the selected upper and lower bounds for k P ; u !max , l !max are the selected upper and lower bounds for ! max ; e !max is an a priori empirical estimate of the maximum angular-velocity control error during the multi- ip phase; e max is an a priori empirical estimate of the maximum angular acceleration control error during the multi- ip phase; and max is the maximum torque that is allowed to be applied, introduced to ensure that the signal [K C (s)e ! ] (t) does not exceed its saturation limit. A general method for solving quasiconvex optimization problems is thoroughly dis- cussed in [31]. In that approach, a primal quasiconvex formulation is converted into a series of convex feasibility problems; then, a selection algorithm identies a solution that is proven to be globally optimal. In the case discussed here, the complete proce- dure for nding a solution numerically close to the globally-optimal parameter set that minimizes the functionf 0 dened in (3.4),x op = k P ,k I , ,! max ,P z ,B z ,B e! , , based on the sub-formulation specied in (3.8), consists of the following three steps: 50 i. Select evenly spaced points along the ranges of k I ; and . Then, generate all the combinations of the selected values fork I ; and . ii. For each selected combination of values forfk I ;;g, solve the auxiliary quasi- convex optimization problem specied by (3.8). The method in [31] guarantees that for each selected setfk I ;;g, the resulting solution to the minimization problem specied by (3.8) is globally optimal. iii. Among all the parameter sets that are identied according to (ii), select the set that yields the smallest value of f 0 dened by (3.8), f 0 . This set, x op = k P ;k I ; ;! max ;P z ;B z ;B e! ; , is globally optimal over the discretized ver- sion of the spacefk I ;;g. For a highly dense discretized space,x op is numerically close to the parameter set that globally minimizes f 0 in (3.4). Note that the numerical implementation of the proposed algorithm in cases with high-density discretizations of the spacefk I ;;g requires a large number of parameter sweepings (of order n k I n n , where n k I , n and n are the numbers of discretized points for parametersk I , and, respectively). Therefore, in the execution of the algo- rithms discussed in the next chapter, we employ the USC high performance computing cluster (HPCC). However, the resulting synthesized controller is extremely simple and computationally ecient, as demonstrated by numerous ight experimental tests. 3.2 Simulation and Experimental Validations 3.2.1 Optimization and Experimental Setups To solve the problem dened by (3.8), we used CVX with the SDPT3 solver, which is a computer package for specifying and solving convex programs [44{46]. To traverse all the selected combinations of the parameter setfk I ;;g and implement the corre- sponding quasiconvex optimization algorithms, we employ a computer cluster that has 51 Table 3.1: Pre-Specied Discretization of the Spacefk I ;;g. Parameter Range Spacing Number of Points k I i 8 10 6 to 8 10 5 4 10 6 19 30 to 330 30 11 110 to 230 10 13 parallel computing capabilities. Four multi- ip maneuver cases are selected for opti- mization: double and triple ips about theb 1 axis, and double and triple ips about the 45 oblique axis dened as p 1 =2(b 1 +b 2 ). The pre-specied discretized parameter space offk I ;;g is dened in Table 3.1, in whichk I i denotes one of the three components of k I . From Table 3.1, it follows that there are a total of 980; 837 = 19 3 11 13 combi- nations to sweep. Additionally, the quadrotor utilized for the experimental validation of the optimization strategy is also Crazy ie 1.0, same as the vehicle used in Chap- ter 2.3.4. Note that, during the aerobatic maneuver, the quadrotor operate completely autonomously. 3.2.2 Optimization, Simulation and Experimental Results First, we present the optimization results that were employed in the process of con- troller synthesis according to the method that was introduced in Chapter 3.1. For each multi- ip case, from the discretization in Table 3.1, there are 980; 837 combinations of parametersfk I ; :g to traverse. To alleviate the massive computational burden of in- vestigating this large number of possibilities, we introduce an additional constraint that bounds the maximum value of % 2 (% 1 + g 0 % 3 ) 1 for each quasiconvex sub-problem such that only a small proportion of all the formulated problems are solved, thereby re- ducing the computational cost without discarding the optimal solution. Fig. 3.1 shows the resulting bound surfaces as functions of 3 1 2 and 1 1 2 for four aerobatic maneu- vers, in which each data point represents the solution of one quasiconvex optimization 52 problem. Here, it can be seen that the value of bound decreases as 3 1 2 and 1 1 2 increase, and that the observed numerical global optima are attained when 3 1 2 and 1 1 2 take their maximum feasible values. (a) Double Flip (b 1 ). (b) Triple Flip (b 1 ). (c) Double Flip (45°). (d) Triple Flip (45°). Figure 3.1: Optimization results for four dierent aerobatic maneuvers: (a) double ip about the b 1 axis; (b) triple ip about theb 1 axis; (c) double ip about a 45 axis as dened in Chapter 2.3.4; and (d) triple ip about a 45 axis as dened in Chapter 2.3.4. In all the cases, each data point represents the solution to a quasiconvex optimization problem as specied by (3.8). The horizontal axes of the plot represent the resulting optimal values of 3 1 2 and 1 1 2 , and surface represents the values of bound as dened in (3.2). 53 (a) Double Flip (b 1 ). (b) Double Flip (b 1 ) { Top View. (c) Triple Flip (b 1 ). (d) Triple Flip (b 1 ) { Top View. (e) Double Flip (45°). (f) Double Flip (45°) { Top View. (g) Triple Flip (45°). (h) Triple Flip (45°) { Top View. Figure 3.2: Exhaustive simulation results obtained through the performance of four types of aerobatic maneuvers. The plots labeled with (a), (c), (e) and (g) show the sim surfaces, i.e., values of obtained through simulations, as functions of 3 1 2 and 1 1 2 . The plots labeled with (b), (d), (f) and (h) show the top views of the 3D plots on the right. In all the plots, each data point represents the solution to a quasiconvex optimization problem as specied by (3.8). 54 The next step is to test through simulations whether the obtained numerical global optimal parameter sets eectively improve the ight performance of the quadrotor dur- ing the execution of multi- ip maneuvers. In order to accomplish this objective, we simulate the same four aerobatic maneuvers in Fig. 3.1 by employing all the parameter sets that were obtained by solving the quasiconvex optimization problem specied by (3.8) over the discretization in Table 3.1 and calculate the corresponding according to (3.1). Fig. 3.2 shows the resulting simulated surfaces, which are labeled as sim , as functions of 3 1 2 and 1 1 2 for all the maneuvers considered here. Fig. 3.2 also shows the corresponding top views, i.e., the projection of sim on the 3 1 2 ; 1 1 2 plane. Here, it can be seen that sim decreases as both 3 1 2 and 1 1 2 increase; this char- acteristic is signicantly more pronounced along the 3 1 2 axis than along the 1 1 2 axis. Also, somehow unexpectedly, sim drastically decreases as 3 1 2 approaches the right boundary and 1 1 2 approaches the bottom boundary of the feasible region on the 3 1 2 ; 1 1 2 plane, respectively. This phenomenon indicates that along these borders, bound is not suciently close to sim for the optimization process to have an eect on sim due to the conservativeness of the exponential bound. To correct this issue, a less conservative bound can be selected to further improve the process of controller synthesis along the boundaries of the feasible region. Overall, however, the simulation results clearly demonstrate that sim decreases to its minimum value as 3 1 2 and 1 1 2 increase to their maximum values, which coincides with the shape trend of bound as a function of 3 1 2 and 1 1 2 . Thus, in conjunction, Fig. 3.1 and Fig. 3.2 provide compelling simulation-based evidence that the proposed approach for performance optimization is eective. To further validate the eectiveness of the method for improving the ight perfor- mance of the yer during the execution of multi- ips, we compare several dozens of the controllers that employ the computationally identied global optimal parameter sets with controllers that employ feasible non-optimal stabilizing parameter sets for 55 the same four aerobatic maneuvers that were used to obtain the plots in Fig. 3.1 and Fig. 3.2. For the execution of the experimental tests, we selected a relatively small number of the thousands of parameter sets that were employed in the simulations that are shown in Fig. 3.2 because it is not realistic to physically implement and run all of them. The selected parameter sets correspond to dierent magnitudes of bound , rang- ing from a feasible non-optimal value to the computationally determined global optimal value. Consistent with the notation employed in the simulation cases, the experimen- tally obtained values of are labeled as exp . In this study, the eects of sensor noise are not considered in the modeling of the system or the stability analysis. However, higher values of are expected to amplify the sensor noise that aects the feedback sig- nals. Therefore, to exclude the distortion of noise in the experimental evaluation of the optimization results, only parameter sets that have the same value of are compared. To run the experiments and calculate the values of exp , we selected three parameter sets for each value of . Then, we compared the values of exp with those of sim to determine if the experimental performances exhibit the same trends as those observed in the simulation results. The comparisons are summarized in Fig. 3.3, where it can be seen that the values of exp for most parameter sets exhibit the same trends as those of sim . Note that in these plots, the system's performance improves from right to left as sim decreases. Thus, in general, the computationally identied optimal parameter sets consistently improve the experimental ight performance of the yer, during the exe- cutions of multi- ips, compared to the performances that were obtained with feasible non-optimal parameter sets. However, it is seen in Fig. 3.3 that a few exceptions occur, as the slopes of some of the lines that connect the data points for a constant value of are negative. This phenomenon only appears at small values of , which leads us to speculate that the nominal dynamic model of the system discussed in Chapter 2.2 might deviate from the true dynamics as becomes smaller. Generally, the results in Fig. 3.3 provide concrete experimental evidence that demonstrates the eectiveness 56 (a) Double Flip (b 1 ). (b) Triple Flip (b 1 ). (c) Double Flip (45°). (d) Triple Flip (45°). Figure 3.3: Comparison of selected experimental and simulation results. In these plots, each data point represents the relationship between a simulated sim value and an experimental exp value corresponding to the same set of controller and analysis parameters. It can be seen that, in most cases, exp decreases as sim decreases, which validates the experimental eectiveness of the proposed method for performance optimization. of the proposed optimization method for improving the performance of aerobatic ma- neuvers and synthesizing high-performance stable controllers. Additionally, indirectly, these results prove that the modeling approach introduced in Chapter 2.2 captures the main characteristics of the true dynamics of the quadrotor. To visualize how a controller with an optimal parameter set improves the ight per- formance of the yer with respect to those obtained with heuristically tuned feasible controllers, we compare the angular rates that were obtained with both type of meth- ods during the execution of triple- ips about the b 1 axis, as exemplied by Fig. 3.4. Here, Fig. 3.4-(a) shows the signalsj! d (t)j 2 ,j!(t)j 2 andje ! (t)j 2 corresponding to the ight experiment that was implemented and run with the computationally identied 57 global optimal parameter set, for which the resulting performance is exp = 1:5669 rad. Similarly, Fig. 3.4-(b) shows the signalsj! d (t)j 2 ,j!(t)j 2 andje ! (t)j 2 corresponding to the ight experiment that was implemented and run with a heuristically tuned non- optimal controller, for which the resulting performance is exp = 3:2179 rad. According to hundreds of experimental observations, this dierence in performance is typical, which is consistent with the analysis summarized in Fig. 3.3. In addition, it can be clearly observed in Fig. 3.4 that the angular-velocity control error after the comple- tion of the aerobatic maneuver is signicantly smaller in the optimal case than in the non-optimal case; hence, the robot can reorient its attitude signicantly faster during the DR phase if the optimal controller, rather than a heuristic controller, is employed during the multi- ip phase. This last observation is consistent with the analysis that was presented in Chapter 2.2.5. An experimental video with performance comparisons between non-optimal and optimal parameter sets for the four multi- ip cases discussed in this paper is available at http://www.uscamsl.com/resources/TCST2019.mp4 (a) Triple Flip (b 1 ) { Optimal Controller Parameters. (b) Triple Flip (b 1 ) { Non-Optimal Con- troller Parameters. Figure 3.4: Comparison of the ight performances achieved during the execution of triple- ips about the b 1 axis with two dierent controller parameter sets. (a) Experimental result obtained with the computationally-found globally-optimal parameter set. (b) Experimental result obtained with the heuristically-tuned non-optimal parameter set. In these experiments, after the triple- ip maneuver is completed,je ! (t)j 2 coincides withj!(t)j 2 becausej! d (t)j 2 becomes zero. 58 As already mentioned, the optimization-based process for controller synthesis pre- sented in this chapter is computational intensive; and for this reason, we employ the USC HPCC to nd the optimal parameter sets that minimize the space-discretized ver- sions of the problem specied by (3.4). However, this optimization process is required and employed for performance optimization only, as it is not a necessary step for nding feasible controllers that stabilize the system. This fact is almost obvious and can be directly deduced from the stability analysis in Chapter 2.2 as the manual selection, or experimentally based identication, of stabilizing control parameters is easily accom- plished by satisfying the conditions in Proposition 2.5. Additional experimental results and analyses that explain and support this statement can be found in [33] and [37]. 59 Chapter 4 Nonlinear Adaptive Control for Aerobatic Flight without Torque Measurements The linear controller design was analytically and experimentally proved to enable the micro UAV to track the high-speed angular velocity reference signal in an extremely short period of time and accomplish the aerobatic maneuvers with a signicant perfor- mance robustness. Additionally, the proposed optimization algorithm accompanying the linear controller design can locate the optimal parameter set to further improve the performance. However, there still exists some noticeable performance dierences between the simulation and experimental results as shown in Chapter 2.3.4, such as an- gular velocity oscillations and overshoots in experimental results, which indicates some high-order dynamics play an important role here for these dierences. In this chap- ter, a momentum theory analysis for the aerobatic maneuver is conducted to identify the eects of the varying ow eld on the aerodynamic coecients and torque latency, which results into the actuator dynamics model. Then, two dierent types of nonlinear adaptive control approaches are proposed to compensate for this actuator dynamics 60 Figure 4.1: The 28-gram micro quadrotor and frames of reference.N =fO 0 ;n 1 ;n 2 ;n 3 g is the inertial frame,B =fO B ;b 1 ;b 2 ;b 3 g is the body frame, andr indicates the position of the yer's center of mass from the inertial originO 0 . in this and following chapters. In addition, the aerodynamic coecient and torque latency variations are shown to be a common phenomenon for other types of aerobatic maneuvers besides the multi- ip maneuver, and the proposed torque actuation dynam- ics and the corresponding adaptive control strategies are veried to be eective for a wide spectrum of aerobatic maneuvers in experiments. 4.1 Preliminaries 4.1.1 System Modeling The quadrotor and frames of references are shown in Fig. 4.1. The inertial frame of reference isN =fO 0 ; n 1 ; n 2 ; n 3 g, and the body frame isB =fO B ; b 1 ; b 2 ; b 3 g. The globally dened nonlinear dynamics including the actuator dynamics is given by m r =mgn 3 +fb 3 ; (4.1) _ q = 1 2 q p; (4.2) J _ ! =!J! +; (4.3) 61 _ =A(t) +B(t)u; (4.4) wherem is the mass of the vehicle,J = diagfj 11 j 22 j 33 g is the moment of inertia matrix with respect toB and assumed to be diagonal,r = [r 1 r 2 r 3 ] T denotes the displacement of the center of mass with respect toN , f is the magnitude of the total thrust force generated by the four propellers,! denotes the angular velocity with respect toN with its components expressed inB, = [ 1 2 3 ] T is the torque generated by the propellers on the vehicle,u = [u 1 u 2 u 3 ] T is the control input, symbol represents the quaternion multiplication, q is the quaternion indicating the attitude ofB with respect toN , and quaternion p = 0 ! T T . Note that (4.1)-(4.3) are exactly same as the nominal dynamics model described in (2.1)-(2.3). The added high-order model is the actuator dynamics (4.4), and the ideas and justications behind this model will be explained later. Time-varying matricesA(t) areB(t) are dened asA(t) = diagfa 1 (t) a 2 (t) a 3 (t)g and B(t) = diagfb 1 (t) b 2 (t) b 3 (t)g. Both A(t) and B(t) are positive denite, and A(t) min I; B(t) min I, whereI is identity matrix and min ; min > 0. Equation (4.4) is the LTV control torque model taking into account the varying aerodynamic co- ecients and torque latency, which will be further explained in the following subsection. Dene the reference angular velocity for the maneuver to be! d = [! d1 ! d2 ! d3 ] T , and A(t) andB(t) are modeled as follows: a i (t) = i1 + i2 ! di + i3 _ ! di ; (4.5) b i (t) = i4 + i5 ! di + i6 _ ! di ; (4.6) where ij is an unknown constant parameter, i2f1; 2; 3g; j2f1; 2; 3; 4; 5; 6g. 62 Initial Hovering Maneuver Starts Climb Phase Aerobatic Maneuver Phase Maneuver Ends Descent and Restabilization Phase Back to Hovering Figure 4.2: The complete process of an Aerobatic Maneuver. The whole process is composed of three phases, climb phase, aerobatic maneuver phase, and descent and restabilization (DR) phase. In climb and DR phases, the UAV stays in the regular ight mode with small angle variations and low angular velocity near hovering state. In the aerobatic maneuver phase, the UAV stays in aerobatic ight mode with large angle variations and high angular velocity. 4.1.2 The Process of the Aerobatic Maneuver In this and following chapters, we will extend the types of maneuvers performed by the micro UAV from the multi- ip to other high-speed aerobatic maneuvers, which will validate the feasibility of the proposed models and controllers on a large spectrum of maneuvers. Due to this, we give the process of the generalized aerobatic maneuver in this subsection, which is similar to the multi- ip process discussed in Chapter 2.1. To perform an aerobatic maneuver, the vehicle follows a ight trajectory consisting of three phases: climb, aerobatic maneuver, and descent-and-reorientation (DR), as shown in Fig.4.2. The DR phase consists of two subphases, DR-1, which lasts for 0:15 s after the aerobatic maneuver phase is completed, and DR-2, during which the UAV resumes regular hovering ight. In the climb phase and DR-2 subphase, the robot is own by simple PID controllers (2.10)(2.11) based on the linearized dynamics 63 (2.4)(2.5) for regular ight as described in Chapter 2.2.2; and, these two segments of the trajectory are not a focus of this and following chapters. The main purpose of this and following chapters is to develop new methods for synthesizing controllers for the aerobatic maneuver that can compensate for the aerodynamic eects. The robot starts the climb phase from a hovering position; then, triggered by an external command, it rises to reach an altitude high enough to perform the aerobatic maneuver. To start the maneuver, the regular- ight controller is automatically turned o and the aerobatic controller is turned on, causing the UAV to perform the maneuver while remaining in free fall. In agreement with experimental observations, due to the fast and periodical changes of direction, the thrust force is assumed to have no net eect on the translational motion of the robot, and gravity is considered to be the only direct external force acting on the vehicle's body. Once the high-speed ight maneuver is completed, the quadrotor descends, reorients and stops dropping by aligning the thrust force vertically to counteract gravity. During this last phase, the UAV is own by the aerobatic controller for 0:15 s (DR-1) and then by the regular- ight PID controller (DR-2). Note that the translational dynamics in (4.1) remains in open loop during the execution of the aerobatic maneuver; therefore, some lateral drift is expected. However, given the short duration of the complete process, possible horizontal translational drifts are presumed to be negligible. 4.1.3 Momentum Theory Analysis for Aerobatic Maneuver In the aggressive aerobatic ight with large attitude angle change and high angular and translational velocities, the rotating propellers' aerodynamic coecients will not remain constant, but change signicantly due to the variation of the direction and magnitude of the local relative incoming ow with respect to propellers. The assumption of the constant aerodynamic coecient is not longer valid in the aerobatic maneuvers. In addition, the noticeable uctuation of the moment load imposed on the propeller due 64 1 2 Figure 4.3: Flow model for the momentum analysis of an aerobatic quadrotor maneuver. The quadrotor rotates about theb 1 axis with angular speed! 1 , and the relative free-stream ow has speedv 1 resulted from the quadrotor's translational motion. to the varying aerodynamic coecients can change the dynamic response of the electric direct current (DC) motors and introduce the time-varying response latency. Momentum theory [47] is applied to model the aerodynamic thrust force and drag moment of the rotor. Even though momentum theory only provides a rst-order ap- proximation of the true dynamics, it establishes a basic principle for the complicated aerodynamics and helps reveal the time-varying essence of the aerodynamic coecients in the aerobatic ight. The ow model for the momentum analysis is plotted in Fig. 4.3, in which without loss of generality the quadrotor has an angular speed ! 1 aboutb 1 axis and encounters a relative free-stream speed, v 1 , that results from the translational ight with inclined angle relative to the rotor-disk plane. Due to the angular speed! 1 , rotor 1 encounters a downward ow with speed ! 1 l, where l is the distance between the rotation axis and the rotor center. For rotor 1, the mass ow rate is _ m f 1 =A r p (v 1 sin +! 1 l +v i1 ) 2 + (v 1 cos) 2 ; (4.7) 65 where v i1 is the induced velocity in the disk plane, is air density, and A r is the rotor-disk area. From the conservation of momentum for rotor 1, we obtain T 1 = _ m f 1 (v 1 sin +! 1 l +w 1 ) _ m f 1 (v 1 sin +! 1 l) = _ m f 1 w 1 ; (4.8) whereT 1 is thrust force,w 1 is the slipstream velocity caused by the rotor rotation, and v 1 sin+! 1 l+w 1 is the slipstream velocity in the vena contracta. Applying the energy conservation, we have T 1 (v 1 sin +! 1 l +v i1 ) = 1 2 _ m f 1 (v 1 sin +! 1 l +w 1 ) 2 1 2 _ m f 1 (v 1 sin +! 1 l) 2 = 1 2 _ m f 1 (2w 1 (v 1 sin +! 1 l) +w 2 1 ): (4.9) It follows from (4.8)(4.9) thatw 1 = 2v i1 . Therefore, from (4.7)(4.8) the thrust coecient of rotor 1 is given by C T 1 = T 1 A r 2 1 R 2 = 2v i 2 1 2 1 R 2 s v 2 1 + 2v 1 sin(! 1 l +v i1 ) v i 2 1 + ( ! 1 l v i1 + 1) 2 ; (4.10) where 1 is the rotor rotation speed and R is the rotor radius. For aerobatic maneuvers,v 1 and! 1 are time-varying and become pronounced com- pared tov i1 . Therefore, the thrust coecient, C T 1 , changes rapidly with respect tov 1 and ! 1 . The same approach can show that the drag moment coecient, C Q 1 , changes with respect to v 1 and ! 1 as well. 66 For hovering ight, v 1 =! 1 = 0, Eq. (4.10) reduces to C T 1 = 2v i 2 1 2 1 R 2 = 2v h 2 2 h R 2 =C Th ; (4.11) wherev h and h are the induced velocity and rotor speed in the hovering ight respec- tively. C Th is normally assumed to be constant and used for the quadrotor dynamic modelings and controller designs [13,14,23,41,48{50]. The ow model of rotor 2, under the conditionv 1 sin>! 1 l, is plotted in Fig. 4.10, and the same momentum analysis can show that the aerodynamic coecients of rotor 2 are time-varying relative to v 1 and ! 1 . When v 1 sin6! 1 l6 2v h +v 1 sin, rotor 2 stays in either vortex ring or turbulent wake states [47], in which the aerodynamic forces are aected by the tip vortices or turbulent wake and depend on v 1 and ! 1 as well. When ! 1 l > 2v h +v 1 sin, rotor 2 stays in windmill brake state, in which the aerodynamic coecients can be analyzed through momentum theory and depend on v 1 and! 1 l as well. In sum, the aerodynamic coecients of both rotors in Fig. 4.10 do not remain constant due to the existence of the noticeable translational speed v 1 and angular speed ! 1 . The rotor speed dynamics is given by [51] J r _ +B m + L = M ; (4.12) whereJ r is the moment of inertia of the rotor system,B m is the damping coecient, L is the load torque on the rotor system, and M is the electromagnetic torque. For rotor 1, the load torque is equal to the aerodynamic drag moment, L =Q 1 , which leads to J r _ 1 + (B m +C Q 1 A r R 3 1 ) 1 = M : (4.13) Therefore, the uctuation of the dragging moment coecient,C Q 1 , aects the dynamic 67 response of the rotor speed, which indicates a time-varying rotor speed latency. The varying aerodynamic coecients and rotor speed latency would aect the torque dy- namics. We use the LTV torque model (4.4) to take these eects into account, in which matrix B(t) re ects the uctuation of the aerodynamic coecients, and matrix A(t) re ects the varying rotor speed latency. For the three maneuvers studied in this work, the quadrotor's translational velocity has a smaller impact on the change of the ow eld than the angular velocity. Therefore, we assume v 1 0 in (4.10) which yields C T 1 = 2v i 2 1 2 1 R 2 (1 + ! 1 l v i1 ): (4.14) This equation indicates that the variation of the thrust coecient depends linearly on the angular velocity, which is why the torque dynamic parameters a i and b i depend linearly on the reference angular velocity in (4.5)(4.6). In addition, the momentum theory analysis is not able to take into account the angular acceleration's eects on the air ow. From the observations of the experiments and engineering intuition, we nd that adding the linear dependency of a i and b i on reference angular acceleration improves the modeling accuracy and ight performance, because the ow acceleration and deceleration can have dierent impacts. Note that the reference signals are imple- mented in the models (4.5)(4.6) instead of the true values which are due to the following reasons. First, the true angular acceleration is not measurable and the low-pass ltered one will inevitably introduce delay which may be problematic for such rapid maneuvers that occur in less than one second. Second, the reference signals are close enough to the true values due to the eect of feedback control, and the reference signals exhibit a prole and trend consistent with the true values, which are more important for the torque dynamics modeling than the minor discrepancies. 68 We aim to nd a torque dynamics model which is able to capture the essence of the aerodynamic eects on the aerobatic maneuvers and is also relatively simple for controller design. We believe the torque dynamic models (4.4)(4.5)(4.6) accomplish the goal in a satisfying way as shown in [52,53]. 4.2 Nonlinear Adaptive Controller Design without Torque Measurements 4.2.1 Backstepping Control At beginning, we start by assuming thatA andB matrices in (4.4) are time-invariant and known, which simplies the control problem and provides a intermediate step towards the nonlinear adaptive control design. From the assumption, we can obtain the following actuator dynamics _ =A +Bu; (4.15) whereu is the algorithmically computed control signal; is the applied control torque; andA = diagfa 1 ;a 2 ;a 3 g andB = diagfb 1 ;b 2 ;b 3 g with a i ; b i > 0, for i2f1; 2; 3g. In this work, (4.15) is mainly determined by the responses of the electric motors driving the quadrotor's propellers, whose parameters were experimentally identied using the LS-based method in [54] through ground experiments. The proposed control scheme employs a desired angular velocity, ! d , with the generic form in Chapter 2.2.1 and reference dynamics J _ ! d +! d J! d = d ; (4.16) in which d is the torque required to obtain the reference ! d according to (4.16). 69 Following the same method of analysis as in Chapter 2.2.3 and [37, 55], we dene the functionalf(!) =!J! and take its Taylor expansion about! d to obtain f (!) =f (! d ) +A f (! d ) (!! d ) +g (e ! ); (4.17) whereA f (! d ) is the 3 3 Jacobian matrix evaluated at! d , i.e., A f (! d ) = @f @! ! d = 2 6 6 6 6 6 6 6 6 4 @f 1 @! 1 @f 1 @! 2 @f 1 @! 3 @f 2 @! 1 @f 2 @! 2 @f 2 @! 3 @f 3 @! 1 @f 3 @! 2 @f 3 @! 3 3 7 7 7 7 7 7 7 7 5 ! d = 2 6 6 6 6 4 0 (j 33 j 22 )! d3 (j 33 j 22 )! d2 (j 11 j 33 )! d3 0 (j 11 j 33 )! d1 (j 22 j 11 )! d2 (j 22 j 11 )! d1 0 3 7 7 7 7 5 ; (4.18) e ! = !! d ; [! d1 ! d2 ! d3 ] T is the coordinate representation of ! d inB; and g(e ! ) is the nonlinear residual associated with the linearized expression in the right side of (4.17), which has the compact formg(e ! ) =e ! Je ! . Note that (4.17) is exact. Subtracting (4.16) from (4.3) yields J _ e ! +A f (! d )e ! +g(e ! ) = d ; (4.19) which, when combined with (4.15), denes the dynamics of the angular velocity error under the in uence of control torque latency during the execution of a high-speed multi- ip maneuver. In this approach, the model in (4.19) is the system's nominal plant studied in Chapter 2.2.3, which was mathematically proven to be stabilizable by using the LTI controller (2.12). Hence, for the synthesis and implementation of the backstepping-based controller according to the method in [56], we use the LTI control 70 scheme (2.12) as the stabilizing function a (t) for the backstepping structure, i.e., a =K 1 e ! +v + d ; _ v =v +K 2 e ! ; (4.20) in whichv is an auxiliary state variable; > 0; andK 1 andK 2 are positive denite diagonal matrices satisfying the relationship K 1 K 2 . For the purposes of stability analysis, controller synthesis, and generation of the control signalu, we dene the proto-Lyapunov function L = 1 2 e T ! Je ! + 1 2 v T Q 1 v + 1 2 ( a ) T Q 2 ( a ); (4.21) in which Q 1 and Q 2 are positive denite diagonal matrices. Next, by dierentiating both sides of (4.21) with respect to time and then substituting the dynamic models given by (4.15) and (4.19) into the resulting equation, we obtain _ L =e T ! ( a d A f (! d )e ! g (e ! )) +v T Q 1 (v +K 2 e ! ) +e T ! ( a ) + ( a ) T Q 2 (A +Bu _ a ): (4.22) By noticing thate T ! g(e ! ) = 0 and using (4.20) to replace a in the rst expression of right side of (4.22), we get the relationship _ L =e T ! (K 1 +A f (! d ))e ! +v T (I +Q 1 K 2 )e ! v T Q 1 v +e T ! ( a ) + ( a ) T Q 2 (A +Bu _ a ); (4.23) 71 which can be rewritten as _ L = e T ! v T 2 6 4 K 1 +A f (! d ) 1 2 I +K T 2 Q T 1 1 2 (I +Q 1 K 2 ) Q 1 3 7 5 2 6 4 e ! v 3 7 5 +e T ! ( a ) + ( a ) T Q 2 (A +Bu _ a ); (4.24) whereI is the identity matrix with compatible dimensions. Following the standard method of Lyapunov-based controller synthesis [35], we here impose a series of conditions in order to make L a Lyapunov function. First,K 1 ,K 2 and are chosen such that 2 6 6 4 K 1 + 1 2 A f (! d ) +A T f (! d ) 1 2 (I +K T 2 Q T 1 ) 1 2 (I +Q 1 K 2 ) Q 1 3 7 7 5 0: (4.25) Note thatA f (! d ) is time-varying; however, it can be shown that 1 2 A f (! d ) +A T f (! d ) remains in a convex hull dened by two constant symmetric matrices as explained in Proposition 2.1, which we denote here byA f 1 andA f 2 . Note that the stability condition (4.25) is equivalent to the stability conditions stated in Proposition 2.5. Furthermore, it can be shown (see Chapter 2.2.3) that (4.25) is equivalent to the conditions 2 6 6 4 K 1 +A f 1 1 2 I +K T 2 Q T 1 1 2 (I +Q 1 K 2 ) Q 1 3 7 7 5 0; (4.26) 2 6 6 4 K 1 +A f 2 1 2 I +K T 2 Q T 1 1 2 (I +Q 1 K 2 ) Q 1 3 7 7 5 0: (4.27) 72 Now, after rewriting the condition in (4.25) to this new form, the task of nding K 1 , K 2 and in order to satisfy (4.26) and (4.27) is reformulated to a semidenite programming problem as explained in Chapter 2.2.3, which we solve using the interior- point method. Then, once the matrix inequalities in (4.26) and (4.27) are satised, if we choose a control input of the form u =B 1 (A + _ a )B 1 Q 1 2 e ! B 1 K 3 ( a ) (4.28) withK 3 diagonal and positive denite, it follows that _ L = e T ! v T 2 6 4 K 1 +A f (! d ) 1 2 I +K T 2 Q T 1 1 2 (I +Q 1 K 2 ) Q 1 3 7 5 2 6 4 e ! v 3 7 5 ( a ) T Q 2 K 3 ( a ); (4.29) and consequently, L is a Lyapunov stability function for the closed-loop dynamics. Thus, by enforcing the conditions in (4.26), (4.27) and (4.28), we make the only xed point of (4.19), 0, exponentially stable. In this scheme, (4.19) and (4.20) are directly employed to obtain _ a =K 1 J 1 d A f (! d )e ! g (e ! ) v +K 2 e ! + _ d ; (4.30) and the control torque in (4.28) and (4.30) is computed employing the model in (4.15). 4.2.2 Nonlinear Adaptive Control In this chapter, we restore the assumptions on the matricesA(t) andB(t) that both matrices are time-varying and have the diagonal elements dened in (4:5) and (4:6) with unknown parameters. Through the application of the certainty equivalence (CE) prop- 73 erty [56,57], we can have the nonlinear adaptive control design from the backstepping control design to be u = ^ B(t) 1 ( ^ A(t)^ + _ ^ a Q 1 2 e ! K 3 (^ a )); (4.31) where _ ^ a =K 1 J 1 (^ d A f (! d )e ! h(e ! ))v +K 2 e ! + _ d ; (4.32) and ^ A(t); ^ B(t) and ^ are the estimations of A(t);B(t) and , respectively. Then, the next question to be answered is how to design an adaptation law to estimate the unknown parameters. A(t) andB(t) are diagonal, and the states of the torque dynamics are decoupled. To simplify the analysis, we chooseb 1 axis as an example, and the other two axes follow the same approach. The angular velocity dynamics on theb 1 axis is given by j 11 _ ! 1 + (j 33 j 22 )! 2 ! 3 = 1 ; (4.33) and the ltered true torque on theb 1 axis is given by f1 = s + 1 = s + j 11 _ ! 1 + (j 33 j 22 )! 2 ! 3 ; (4.34) where is a positive constant. As _ ! is not measurable, 1 can not be computed from (4.33). We apply low pass lter on 1 so that the ltered true torque f1 is realizable. 74 From (4.4)(4.5)(4.6), we have _ 1 =a 1 (t) 1 +b 1 (t)u 1 = 11 12 13 14 15 16 2 6 6 6 6 6 6 6 6 6 6 6 6 6 6 4 1 ! d1 1 _ ! d1 1 u 1 ! d1 u 1 _ ! d1 u 1 3 7 7 7 7 7 7 7 7 7 7 7 7 7 7 5 = T 1 1 : (4.35) Multiply both sides of the above equation by operator (s+) 2 and assume initial condi- tions to be zero and get s s + f1 = 11 12 13 14 15 16 2 6 6 6 6 6 6 6 6 6 6 6 6 6 6 4 (s+) 2 [ 1 ] (s+) 2 [! d1 1 ] (s+) 2 [ _ ! d1 1 ] (s+) 2 [u 1 ] (s+) 2 [! d1 u 1 ] (s+) 2 [ _ ! d1 u 1 ] 3 7 7 7 7 7 7 7 7 7 7 7 7 7 7 5 = T 1 (s +) 2 [ 1 ] = T 1 1 : (4.36) In this structure, the rightmost vector in the right side of (4.36) is simply composed of the signals that result from ltering 1 , ! d1 1 , _ ! d1 1 , u 1 , ! d1 u 1 and _ ! d1 u 1 through realizable transfer functions. However, since the true control torque, 1 , can not be directly measured and the model in (4.33) is not realizable, we need to rewrite the rst three entries of the column vector in (4.36) in order to obtain completely real-time computable structures. For the rst entry, this objective is attained by computing f1 using (4.34); then, employing the relationship (s +) 2 [ 1 (t)] = 1 s + [ f1 (t)]: (4.37) 75 To deal with the second and third entries of the column vector in (4.36), we employ (4.37) and invoke the swapping lemma [57,58]. Namely, to rewrite the second entry in a realizable form, we rst nd a minimal state-space realization of (s +) 2 with matrices of the form A m = 2 6 4 2 2 1 0 3 7 5 ; B m = 2 6 4 1 0 3 7 5 ; C m = 0 ; D m = 0: (4.38) Then, we dene the transfer matrices W c (s) =C m (sIA m ) 1 = (s+) 2 (s+2) (s+) 2 ; W b (s) = (sIA m ) 1 B m = 2 6 4 s (s+) 2 1 (s+) 2 3 7 5 ; (4.39) and apply the swapping lemma to obtain (s +) 2 [! d1 (t) 1 (t)] =! d1 (t) (s +) 2 [ 1 (t)] +W c (s)[(W b (s) [ 1 (t)]) _ ! d1 (t)] =! d1 (t) 1 s + [ f1 (t)] +W c (s) 2 6 4 0 B @ 2 6 4 s (s+) 1 (s+) 3 7 5 [ f1 (t)] 1 C A _ ! d1 (t) 3 7 5 : (4.40) The same method can be applied to rewrite the third entry of the column vector in (4.36) into a real-time computable form. Since the structure in (4.36) is in the standard static parametric form [57], there exists a great number of adaptive ltering algorithms suitable for the real-time es- timation of the vector 1 . Dene z 1 = s s+ [ f1 ], and let ^ 1 denote the estimation 76 of 1 . Then, ^ z 1 is the estimation of z 1 and has the form ^ z 1 = ^ T 1 1 . And dene " 1 = (z 1 ^ z 1 ) =m 2 s 1 = (z 1 ^ T 1 1 ) =m 2 s 1 , where m 2 s 1 = 1 + 1 T 1 1 , 1 > 0. We apply the modied Least-Squares algorithm with forgetting factor and normalization [57, 58] to estimate 1 and have _ ^ 1 =P 1 " 1 1 ; (4.41) _ P 1 = 8 > < > : 1 P 1 P 1 1 T 1 P 1 m 2 s 1 ; if kP 1 k F r 1 0; otherwise (4.42) where 1 is the positive forgetting factor and r 1 is the positive upper bound for the Frobenius norm ofP 1 . The parameter projection is implemented as follows ^ 1 = 8 > > > > < > > > > : ^ 1 ; if ^ 1min ^ 1 ^ 1max ; ^ 1min ; if ^ 1min ^ 1 ; ^ 1max ; if ^ 1max ^ 1 ; (4.43) where ^ 1min and ^ 1max are lower and upper bounds of ^ 1 . The estimated time-varying parameters ^ a 1 (t) and ^ b 1 (t) are ^ a 1 (t) = ^ 11 (t) + ^ 12 (t)! d1 (t) + ^ 13 (t) _ ! d1 (t); (4.44) ^ b 1 (t) = ^ 14 (t) + ^ 15 (t)! d1 (t) + ^ 16 (t) _ ! d1 (t): (4.45) Similar, the parameter projection is applied to make sure ^ a 1 and ^ b 1 remain positive as follows ^ a 1 = 8 > < > : min ; if ^ a 1 < min ; ^ a 1 ; if ^ a 1 > min ; ^ b 1 = 8 > < > : min ; if ^ b 1 < min ; ^ b 1 ; if ^ b 1 > min ; where min and min are lower bounds ofa 1 andb 1 respectively. The same approach can 77 be applied onb 2 andb 3 axes to acquire the estimated n ^ a 2 (t); ^ a 3 (t); ^ b 2 (t); ^ b 3 (t) o and construct the estimated ^ A(t) and ^ B(t), which dene the dynamic model to estimate the control torque as follows: _ ^ = ^ A(t)^ + ^ B(t)u: (4.46) The estimated ^ A(t), ^ B(t) and ^ then can be utilized to generate the nonlinear adaptive control input (4.31) and the auxiliary state (4.32). Notice that the formation of the nonlinear adaptive control input (4.31) only requires the angular velocity measurement, and the torque measurement is not needed. However, this inevitably brings the chal- lenge and diculty for the closed-loop stability analysis. Due to the application of the low pass lter to address the lack of the torque measurements, the properties of the estimation error " i can not be conveyed to the system states. 4.3 Experimental Results 4.3.1 Experiment Setup and the UAV vehicle The quadrotor UAV employed in the experiments of this chapter, shown in Fig. 4.1, is the Crazy ie 2.0, manufactured by Bitcraze AB 1 . This yer has a weight of 28 g (in- cluding the battery), a propeller-tip-to-propeller-tip distance of 14 cm and a maximum thrust-to-weight ratio of approximately 1.9. Due to its small size, this robot is more agile than larger quadrotors and capable of performing aerobatic maneuvers at much higher angular speeds [41]. Because the nonlinear adaptive control strategy requires more computation resources than the linear control strategy (2.12), the Crazy ie 2.0 UAV is employed for the nonlinear adaptive control method, and it has more memory space and more computation capacity than the Crazy ie 1.0. The control algorithms 1 https://wiki.bitcraze.io/projects:crazy ie2:index 78 run at 500 Hz, employing a 3-axis gyroscope, a 3-axis accelerometer and a barometer. In addition to the yer, the experimental setup includes a ground computer that sends the control reference signals during regular ight and the trigger command that initiates the multi- ip maneuver. However, during the execution of an aerobatic maneuver, the quadrotor relies entirely on onboard sensing, signal processing and control. Therefore, the proposed control scheme enables the robot to perform high-speed aerobatic maneu- vers completely autonomously. For monitoring, recording and analysis purposes, the measured ight information is communicated from the quadrotor to the ground com- puter in real time at a frequency of 100 Hz. More details about the experimental setup and real-time implementation of the ight controllers are in Chapter 2.3.1. 4.3.2 Reference Angular Speed Plannings for Aerobatic Flight In this chapter, three dierent aerobatic maneuvers are specied and executed in the experiments to validate the actuator dynamics and the adaptive control strategy, which are: a triple ipping maneuver with a 1080° rotation about theb 1 axis; the Pugachev's Cobra maneuver during which the quadrotor rotates 180° about theb 1 axis and then rotates 180° in the opposite direction; and a mixed ipping maneuver that involves alternating the ipping axis, during which the quadrotor rotates 140° about theb 1 axis rst, then rotates 180° about theb 2 axis, and nally rotates 40° about theb 1 axis to return to hovering state. The reference angular speeds about the b 1 and b 2 axes for each aerobatic maneuver are shown in Fig. 4.4. ! d3 = 0 for all the maneuvers, so not plotted here. A cubic function is used to describe the reference angular speed for each maneuver. For example,! d1 can be divided into three segments in Fig. 4.4a, acceleration, constant speed and deceleration, which was also explained in Chapter 2.2.1. The acceleration 79 (a) Triple Flip aboutb 1 axis. (b) Pugachev's Cobra aboutb 1 axis. (c) Mixed Flipping. Figure 4.4: Reference angular speed plannings of three aerobatic maneuvers. (a)-(c) show the reference angular speeds, ! d1 and ! d2 , of the three aerobatic maneuvers. ! d3 = 0 for all the maneuvers. segment in Fig. 4.4a is described as follows: ! d1 (t) = 1 3 (t 1 ) 3 1 2 1 t + 1 3 1 3 ; (4.47) where 1 = ! 1 max g;1 and 1 =( 3 =4) 3 1 ! max . g;1 is the total rotation angle during the acceleration segment, and ! max is the maximum angular speed. The duration of the rst segment is 1 = 2 1 , so 0 6 t 6 1 in (4.47). Therefore, the cubic function of each segment is parametrized by two parameters, the desired rotation angle and the maximum angular speed. The same planning method is applied to the segments of other two maneuvers. Through this method, the reference angular acceleration remains continuous during the whole maneuver. The more details about the speed planning can be found in Chapter 2.2.1, and the reference angular speed plannings of three dierent maneuvers are implemented in the aerobatic maneuver phase of the whole maneuver 80 Figure 4.5: Performance comparison of the three control methods. This plot summarizes three sets of ten triple- ip experiments run using the LTI benchmark controller, the backstepping controller and the nonlinear adaptive scheme. The average value of the PFM for each set is indicated with a square; each bar indicates the experimental range of for each set; the corresponding ESD for each set is written between square brackets. process discussed in Chapter 4.1.2. 4.3.3 Experimental Results with Dierent Controllers To demonstrate the capabilities and performance of both the backstepping controller and adaptive control scheme, we command the quadrotor UAV to autonomously exe- cute a triple- ip about theb 1 axis. During the performance of this aerobatic maneuver only the angular speed about the b 1 axis is required to track a non-zero time-varying ipping-speed reference while the angular speeds about the other two principal axes, b 2 and b 3 , must remain equal to zero. Therefore, in these specic test experiments, the backstepping and adaptive algorithms are only run to control the rotation about b 1 while the rotations about b 2 and b 3 are regulated using the LTI controller (2.12). This approach signicantly reduces the computational burden on the onboard micro- controller. 81 The performances of the backstepping and adaptive controllers during the execution of a multi- ip maneuver are assessed here through the PFM = Z t end t 0 je ! ()j 2 d; (4.48) in which t 0 and t end denote the start and end times of the multi- ip maneuver, re- spectively. Notice that this PFM is same as the one used in Chapter 3.1. Consistent with the denition of (4.48), it follows that lower values indicate better performance. In practice, since the algorithms are implemented digitally, we estimate using the discrete-time version of (4.48); i.e., P k je ! (k)j 2 , wherek is the discrete-time index and e ! (k) is the sampled version ofe ! (t). The eectiveness of the proposed approach for controller synthesis is compellingly demonstrated in Fig. 4.5. This plot summarizes the control results for the LTI controller (2.12) (used as a benchmark), backstepping controller (4.28) and nonlinear adaptive controller (4.31), which were each used to run ten ight experiments with the Crazy ie 2.0 micro quadrotor. For each control method, the average of the ten experimental s is indicated with a square (yellow, green and red, respectively), the corresponding measured range (MR) is depicted with a vertical bar and the value of the associated empirical standard deviation (ESD) is written between square brackets. On average, the backstepping controller decreases the PFM in (4.48) by 17.28 % and the RMS value by 17.82 % when compared to the benchmark. Furthermore, the nonlinear adaptive controller decreases the PFM in (4.48) by 38.66 % and the RMS value by 38.42 % when compared to the benchmark. As can be seen in Fig. 4.5, the ESDs associated with the three control methods are signicant; and typically, the best performances obtained with the LTI benchmark controller are better than the worst ones obtained with the backstepping controller, and the best performances obtained with the backstepping controller are very similar to the worst obtained with the nonlinear adaptive scheme. 82 However, from Fig. 4.5 it is clear that, typically, the worst performances obtained with the nonlinear adaptive scheme are signicantly better than the best ones obtained with the LTI benchmark controller. To visualize the characteristics of the control error signal that explain the dier- ences in performance obtained with the three controllers, we show the time-series of the desired and measured angular velocities in Fig. 4.6. Here, Fig. 4.6-(a) shows that in the benchmark experimental case, the angular velocity is corrupted by oscillations that are presumably caused by time-varying actuator dynamics resulted from the aero- dynamic coecients and torque latency variations. This disturbance is signicantly attenuated by the backstepping controller and almost completely eliminated by the proposed nonlinear adaptive control scheme, as shown in Fig. 4.6-(b) and Fig. 4.6-(c). The fact that the experimental results achieved with both the backstepping controller and the nonlinear adaptive scheme are signicantly better than those obtained with the LTI benchmark control method provides evidence supporting the hypothesis that the oscillatory disturbance observed in Fig. 4.6-(a) is caused by actuator dynamics. Fur- thermore, the result in Fig. 4.6-(c) indicates that the actuator dynamics aecting the system is indeed time-varying, as modeled in (4.4)(4.5)(4.6). Fig. 4.7 shows the time-series of f1 and (s +) 1 [^ 1 ] measured during the exe- cution of the triple- ip maneuver, associated with the backstepping and adaptive con- trollers. For the backstepping case in Fig. 4.7-(a), the signal ^ 1 is computed using the LTI torque model in (4.15); and for the adaptive case in Fig. 4.7-(b), ^ 1 is obtained using the linear time-varying model in (4.46) whose parameters, ^ a 1 (t); ^ b 1 (t) , are esti- mated adaptively in real time according to the algorithm underlying the computational structure in (4.44)(4.45). In Fig. 4.7-(a), the ltered estimated torque (s +) 1 ^ 1 , generated with constant identied parametersa 1 andb 1 , deviates signicantly from the ltered true torque f1 during the execution of the multi- ip maneuver. In contrast, the signal (s +) 1 ^ 1 in Fig. 4.7-(b), computed according to the adaptive control 83 (a) Linear Controller. (b) Backstepping Controller. (c) Adaptive Controller. Figure 4.6: Three triple- ip experiments with three dierent controllers. In (a), the linear controller in (2.12) achieves a PFM = 0:9391 rad. In (b), the backstepping controller in (4.28) compensates the time-invariant torque delay model, achieving a PFM = 0:7348 rad. In (c), the nonlinear adaptive controller (4.31) compensates the time-varying actuator model, achieving a PFM = 0:5340 rad. strategy employing the time-varying torque model in (4.46), matches the ltered true torque f1 considerably more accurately. This comparison explains why the adaptive scheme produces better experimental results than the backstepping controller. Addi- tionally, this result indirectly shows that the LTV actuator dynamics (4.4)(4.5)(4.6) captures the essence of the high-order torque dynamics and veries that the aerody- namic coecients and torque latency are varying in the high-speed maneuver. In all of the tests presented in this chapter, once the multi- ip maneuver is com- pleted, the robot continues to be own by the aerobatic controller for an additional 150 ms; then, it is turned o and a regular- ight LTI controller (identical in all ight tests) is turned on in order to nalize the descent and restabilization of the yer's atti- tude. We therefore expect the quadrotor's performance after the angular speed reference 84 (a) Backstepping Controller. (b) Adaptive Controller. Figure 4.7: Filtered true and estimated torques during the triple- ip maneuver. f1 is the ltered true torque dened in (4.34). In (a), ^ 1 is estimated with the time-invariant torque model in (4.15). In (b), ^ 1 is estimated with the time-varying torque model in (4.46) and a parameter-updating law. (in green in Fig. 4.6) reaches zero to be similar across all three experiments, which is indeed the case. During rapid vertical descent, the local velocities of the air ows in- teracting with the propellers can be signicantly dierent than those observed during regular hovering ight; therefore, the aerodynamic force coecients of the propellers can be drastically dierent in these two ight conditions. This phenomenon makes the restabilization of the robot nontrivial. From a computational perspective, implementing the proposed adaptive control scheme is challenging because a triple- ip maneuver lasts less than one second, thus providing an extremely short amount of time for parameter updating. The experimen- tal results compellingly demonstrate the computational eciency, stability robustness and high performance of the proposed adaptive control method. According to repeated empirical observations, employing the proposed nonlinear adaptive controller, the mi- cro quadrotor can perform 17 consecutive triple- ips until the battery is completely depleted. The video of the complete set of experiments is available at http://www.uscamsl.com/resources/IROS2018/S1.mp4 85 4.3.4 Experimental Results for Dierent Types of Maneuvers In this chapter, we test the nonlinear adaptive control strategy on the three dierent types of aerobatic maneuvers introduced in Chapter 4.3.2, and use the LTI controller (2.12) as a benchmark to compare. The aerobatic ight performance is assessed by a dierent PFM criterion as follow- ing: = Z t end t 0 je ! ()j 2 2 d 1 2 ; (4.49) wherejj 2 represents the Euclidean norm, andt 0 andt end denote the start and end time of the maneuver. It follows that lower indicates better performance. Both the linear controller (2.12) and adaptive controller (4.31) are implemented in experiments for three dierent aerobatic maneuvers: triple ipping, Pugachev's Cobra, and mixed ipping. The performance comparisons of the two controllers are demon- strated in Fig. 4.8. Each bar in Fig. 4.8 indicates the range of of at least 10 experiments, and the point in the middle denotes the average value. It follows that the adaptive controller improves the performances of the three maneuvers by 40.91%, 34.67% and 28.18% respectively. This result shows the aerobatic ight performance signicantly improves through considering the eects of time-varying aerodynamic co- ecients and torque latency for the controller design. In Fig. 4.9, six example maneuver experiments are presented to further compare the performance dierences between the linear and adaptive controllers. In Figs. (4.9a)(4.9c)(4.9e), the true angular velocities under the linear controller are distorted by undesired oscillations and overshoots in three dierent maneuvers. In Fig. (4.9b)(4.9d) (4.9f), the undesired oscillations and overshoots are substantially inhibited by the adap- tive controller. This indirectly conrms the hypothesis that the undesired angular veloc- ity oscillation results from the time-varying aerodynamic coecients and torque latency, 86 (a) Triple Flipping. (b) Pugachev's Cobra. (c) Mixed Flipping. Figure 4.8: Performance comparisons of three maneuvers with linear and adaptive controllers. In (a), the average with the linear controller (2.12) is 1:2240, and the average with the adaptive controller (4.31) is 0:7233. In (b), the average with the linear controller is 1:1304, and the average with the adaptive controller is 0:7385. In (c), the average with the linear controller is 1:9715, and the average with the adaptive controller is 1:4160. and this interpretation is general for dierent types of aerobatic maneuvers. The con- sistent performance improvements in three dierent maneuvers verify the applicability of the adaptive control strategy for a large spectrum of aerobatic maneuvers. The ltered estimated and true torques of the three maneuvers are plotted in Fig. 4.10. They are taken from experiments in Figs. (4.9b)(4.9d)(4.9f). Only the ltered estimated and true torques of the axis that utilizes the adaptive controller are plot- ted. The ltered estimated torque generally match the ltered true torque in the three maneuvers as shown in Fig. 4.10a-4.10c, which demonstrate the eectiveness of the adaptation law (4.41)(4.42) and the torque estimation model (4.4)(4.5)(4.6). The es- timated time-varying parameters ^ a i (t) and ^ b i (t) of the axis that utilizes the adaptive controller are plotted in Fig.4.11. Consistent with the hypothesis thatA(t) areB(t) of 87 (a) Triple ip with linear controller. ( = 1:0392) (b) Triple ip with adaptive con- troller. ( = 0:6886) (c) Pugachev's Cobra with linear controller. ( = 0:8984) (d) Pugachev's Cobra with adaptive controller. ( = 0:6121) (e) Mixed Flipping with linear con- troller. ( = 1:7465) (f) Mixed Flipping with adaptive controller. ( = 1:1711) Figure 4.9: Three dierent maneuvers under two dierent controllers. The reference and true angular velocities of six example maneuvers are demonstrated here. For triple ip and Pugachev's Cobra maneuvers, the reference signals on b 2 and b 3 axes are equal to zero, so are not plotted here. For mixed ipping maneuver, the reference signal onb 3 is zero. the torque dynamics are time-varying during the maneuvers, the estimated parameters ^ a 1 , ^ b 1 , ^ a 2 and ^ b 2 change dramatically during the maneuvers. 88 (a) Triple Flipping. (b) Pugachev's Cobra. (c) Mixed Flipping. Figure 4.10: Filtered estimated and true torques in three maneuvers. In (a)(b), the ltered estimated torque s+ ^ 1 and the ltered true torque f1 on b 1 axis are given. In (c), the ltered estimated and true torques of the ipping axes are plotted respectively. Note that the ltered estimate torques remain zero before maneuvers. Each maneuver takes less than one second and the transient performance of the parameter adaptations do not corrupt the aerobatic performances in experiments, which shows the robustness of the adaptive control strategy. An experimental video is available at http://www.uscamsl.com/resources/ICRA19s1.mp4 4.3.5 Data-Based Persistent Excitation Analysis In this chapter, the persistent excitation property [57{59] of the regressor vector in the adaptation law (4.41)(4.42) is analyzed via the experimental data, which provide insights for the stability of the closed-loop adaptive system. 89 (a) Triple Flipping. (b) Pugachev's Cobra. (c) Mixed Flipping. Figure 4.11: The estimated ^ a i (t) and ^ b i (t) in three maneuvers. The estimated ^ a 1 (t) and ^ b 1 (t) are plotted for the triple ipping and Pugachev's Cobra maneuvers. For the mixed ipping, the ipping axis is alternatively changed betweenb 1 andb 2 axes, so the estimated parameters of the two axes are plotted respectively. Note that the estimated parameters ^ a i (t) and ^ b i (t) remain constant before and after each maneuver. Denition 4.1. The vector2R n is persistently exciting (PE) if it satises % 2 I 1 T 0 Z t+T 0 t () T ()d% 1 I; (4.50) for some % 2 % 1 > 0, T 0 > 0 and8tt 0 . The denition (4.50) is used to justify the regressor vector i ms i in experiments is PE. The values of T 0 , % 1 and % 2 for the three dierent aerobatic maneuvers are listed in Table 4.1. % 1 and % 2 are calculated as follows % 1 = min k ( min 1 T 0 k+ T 0 t X k i (k) T i (k) m 2 s i ); (4.51) 90 % 2 = max k ( max 1 T 0 k+ T 0 t X k i (k) T i (k) m 2 s i ); (4.52) wherek 0 kk 1 T 0 t ,k 0 andk 1 denote the discrete start and end time indices of the maneuver, t is the time interval of the sampled data, and subscript i indicates the axis that utilizes the adaptive controller. It follows from Table 4.1 that the regressor vector i ms i satises the PE condition during the three aerobatic maneuvers. 2,3,4 From the property of the modied least-squares algorithm [58] [57], if i ms i is PE, the estimation error e i = i ^ i goes to zero exponentially fast, i2f1; 2; 3g, which indicates that e A =A ^ A and e B =B ^ B decay to zero. Therefore, if the state errors are assumed to grow slower than the convergence of the estimation errors, the adaptation law that has regressor vector satisfying the PE property is able to generate the stabilizing parameters before the nite time escape of the states occurs, in which the closed-loop stability is guaranteed. Table 4.1: PE check of the regressor vector i ms i in experiments T 0 (s) % 1 % 2 Triple Flipping 0.35 5:6471 10 17 9:0139 10 6 Pugachev's Cobra 0.10 2:9829 10 15 4:1598 10 5 Mixed Flipping 0.35 1:8003 10 15 2:1799 10 5 In the next chapter, we will integrate the torque measurements into the controller design that compensates for the negative aerodynamic eects during the aerobatic ma- neuver, in which the rigorous stability analysis is able to be obtained. 2 % 1 and % 2 are extremely small due to the miniature size of the regressor vectors. The matrices with minimal eigenvalue equal to % 1 are veried to be nonsingular using Matlab function `rank' with default tolerance. 3 The regressor vector i =ms i is computed from the axis that uses the adaptive controller. For the mixed ipping, 1 and 2 are computed based on the regressor vectors of bothb 1 andb 2 axes. 4 The data in Table 4.1 correspond to the experiments in Fig. (4.9b)(4.9d)(4.9f). 91 Chapter 5 Nonlinear Adaptive Control for Aerobatic Flight with Torque Feedback In the previous chapter, we presented a nonlinear adaptive control strategy to compen- sate for the linear time-varying actuator dynamics with unknown parameters during the aerobatic ight. In that approach, the adaptive control method does not require measurements or observations of the actuation torques. This brings, however, extra challenges for the closed-loop stability analysis of the proposed control method. To address this problem and provide a rigorous analysis, we include the torque feedback into the adaptive control design, and this leads to two dierent types of adaptive con- trol methods. The rst one is a Lyapunov-based adaptive control, and the second one is a modular adaptive control. In the proposed Lyapunov-based adaptive control method, the controller and the adaptation law are derived from the same Lyapunov candidate function simultaneously. In the modular adaptive control method, the adap- tive controller and the parameter adaptation law are derived separately, which relaxes the constraint on the form of the adaptation law. Rigorous stability analyses are given 92 for both adaptive control methods, and we employ experimental results to validate the eectiveness of both adaptive control schemes. 5.1 Preliminaries 5.1.1 System Dynamics with LTV Actuator Dynamics The LTV actuator dynamic model utilized in this chapter is slightly dierent from that specied by (4.4)(4.5)(4.6). To avoid any possible potential confusion and make the analysis consistent, we restate the system dynamics with the LTV actuator dynamics here. The UAV employed in the analyses and experiments, a 28-gram quadrotor aircraft, is shown in Fig. 4.1. The inertial frame of reference isN =fO 0 ; n 1 ; n 2 ; n 3 g, and the body-xed frame isB =fO B ; b 1 ; b 2 ; b 3 g. The globally dened nonlinear dynamics with the LTV actuator dynamics is m r =mgn 3 +f t b 3 ; (5.1) _ q = 1 2 q p; (5.2) J _ ! =!J! +; (5.3) _ =a(t) +b(t)u; (5.4) where m is the mass of the vehicle, J = diagfj 11 j 22 j 33 g is the moment of inertia matrix with respect toB,r denotes the displacement of the center of mass with respect toN , f t is the magnitude of the total thrust force, q is the quaternion that describes the attitude ofB with respect toN , the quaternion p = 0 ! T T , the symbol represents the quaternion multiplication,! = ! 1 ! 2 ! 3 T denotes the angular velocity with respect toN with its components expressed inB, = [ 1 2 3 ] T is the torque generated by the propellers of the vehicle, andu = [u 1 u 2 u 3 ] T is the control input. 93 Equation (5.4) describes the eects of the aerodynamic coecients and torque la- tency variations on the actuator dynamics, and it is derived employing the momentum theory analysis in Chapter 4.1.3 and engineering intuition with experimental observa- tions. In this model, the unknown time-varying parameters a(t) and b(t) represent the eects of the time-varying torque latency and time-varying aerodynamic coecients respectively, which are modeled as a(t) = 1 + 2 1 (! a ) + 3 2 ( _ ! a ) = T 1 ; (5.5) b(t) = 4 + 5 1 (! a ) + 6 2 ( _ ! a ) = T 2 : (5.6) In above equations, 1 = [ 1 2 3 ] T and 2 = [ 4 5 6 ] T are unknown constant vectors; = [1 1 (! a ) 2 (! a )] T ;! a =! T a b f is the angular speed about the rotation axis of the aerobatic maneuver; _ ! a = _ ! T a b f is the angular acceleration about the rotation axis of the aerobatic maneuver, where a b f is the body-xed rotation axis; and the functions 1 (! a ) and 2 ( _ ! a ) are given by: 1 (! a ) = ! a p 1 + 2 1 ! 2 a ; 2 ( _ ! a ) = _ ! a p 1 + 2 2 _ ! 2 a ; (5.7) where 1 ; 2 > 0. Note that it is straightforward to show that 1 < 1 (! a ) < 1 ; 2 < 2 ( _ ! a )< 2 ;8! a ; _ ! a . In addition, from the physical property of the system, we can assume that a(t)>a min > 0; b(t)>b min > 0. Note that the modeling of the time-varying parameters in (5.5)(5.6) is dierent from that dened by (4.5)(4.6). The functions 1 () and 2 () are used to bound the time-varying parameters, and the actual angular speed and acceleration,! a and _ ! a , are employed instead of the reference signals in (4.5)(4.6). The signal _ ! a is obtained by using a nonlinear derivative observer that will be discussed in Chapter 5.3.2. The model specied by (5.4) simplies the complex interaction of the propellers with the ambient air ow during the aerobatic ight; thus, simplifying the process of 94 controller synthesis while the essence of the aerodynamic eects is still captured as shown in Chapter 4. 5.1.2 Attitude Error Kinematics To make the analysis consistent with the adaptive controller synthesis to be discussed, we restate the attitude kinematics that was explained in Chapter 2.2.4. The ideal body frameI for the aerobatic maneuver is dened as the desired attitude kinematics that is driven by the reference angular velocity with a zero initial attitude error with respect to the inertial frame. According to this denition, the evolution ofI is governed by _ q d = 1 2 q d p; with p = 0 ! T d T ; (5.8) where the quaternion q d describes the attitude ofI. Since the desired angular velocity ! d is expressed inB, notI, we use the ideal angular velocity ! d in (5.8). Note that ! d has exactly the same components as! d , but is expressed inI. The attitude error quaternion q e is used to denote the error between the body frame, B, and the ideal body frame ,I, which is given by q e = q 1 d q: (5.9) Note that the attitude error quaternion q e contains the information regarding both the rotation axis and the rotation angle fromB toI; Namely, q e = e 0 e T 1 T . Thus, by dierentiating (5.9) along the trajectories of (5.2) and (5.8), we obtain _ q e = 2 6 4 _ e 0 _ e 1 3 7 5 = 1 2 2 6 4 e T 1 (! ! d ) e 0 (! ! d ) +e 1 (! + ! d ) 3 7 5 : (5.10) Hence, we dene the rotation angle fromB toI as e , which in this chapter is used to 95 measure the attitude error. And e is nonnegative from its geometric denition. From the denition of the attitude error quaternion, we have e 0 = cos 1 2 e . 5.1.3 A Class of Aerobatic Maneuvers In Chapter 4.3.2, we introduced two additional types of aerobatic maneuvers besides the multi- ip maneuver. In this chapter, we will generalize the concept of the aerobatic maneuver that can be performed by the controllers proposed in this dissertation. The generic class of aerobatic maneuvers considered is essentially a fast attitude rotational motion about unitary or multiple axes sequentially, during which the translational motion of the vehicle stays in free-fall and open-loop. The rotation axes of the sequential sub-rotations in a aerobatic maneuver are dened by a sequence of body-xed unit vectors in a order of occurrence asfa 1 ;a 2 ;:::;a m g, and the sequence of the rotation angles of the corresponding sub-rotations are given byf 1 ; 2 ;:::; m g. In addition, the sequence of the maximum angular speeds of the corresponding sub-rotations is given by n _ max1 ; _ max2 ;:::; _ maxm o . The reference angular velocity! d is composed of a sequence of reference signalsf! d1 ;! d2 ;:::;! dm g of the sub-rotations in a order of the occurrence. Then we have the following denition for the class of aerobatic maneuvers. Denition 5.1. The class of generic aerobatic maneuvers for VTOL UAVs studied in this dissertation starts from an initial attitude, in which the direction of the thrust force points vertically upward, to an ending attitude, in which the direction of the thrust force returns back to be vertically upward, and is composed of a sequence of sub-rotations. During each sub-rotation, the vehicle rotates about the body-xed axis a i accumulating the rotation angle i without surpassing the maximum rotation speed _ maxi , where i = 1; 2;:::;m. The reference signal for each sub-rotation,! di (t), starts and ends with zero, where i = 1; 2;:::;m. The reference angular velocity expressed in the body frameB,! d (t), is continuously dierentiable. In addition, the translational 96 motion stays in free-fall and follows a parabolic trajectory governed by gravity during the maneuver. In Denition 5.1, the reference angular velocity is required to be continuously dier- entiable so that the reference angular acceleration is compatible with the characteristic of yer's actuators. For example, a feasible approach to specify ! di (t) is employing cubic functions, as explained in Chapters 2.2.1 and 4.3.2. From the Denition 5.1, it follows that R ! di (t)dt = i and maxf _ ! di (t)g = _ maxi . The translational motion is assumed to stay in free-fall state, since the thrust-force direction changes extremely fast 1 between the upward and downward orientations and gravity is the most dominant external force aecting the translational motion. The translation drift on the horizon- tal plane may exist, but it remains bounded due to the extremely short duration of maneuver. The entire maneuver process can be divided into three phases: climb, aerobatic maneuver and descent-and-restabilization as discussed in Chapter 4.1.2. 5.1.4 Control Objectives The control objective is to enable the angular attitude and velocity dynamics (5.2)(5.3) to track the ideal attitude kinematics described by (5.8) and the reference angular velocity! d , respectively, in the presence of the LTV actuator dynamics (5.4), where the unknown time-varying parameters are dened according to (5.5) and (5.6). Due to the extremely short period time of the aerobatic maneuver (normally less than one second) and the extremely fast angular speeds (the maximum is about 2000 deg=s), it is more important to control the \faster" dynamics given by (5.3) than the \slower" dynamics given by (5.2). In specic, the control objectives during the aerobatic maneuver are specied as 1 In this work, the maximum angular speed of the studied maneuvers is larger than 1000 deg=s, and the maneuver takes less than one second. 97 1. !;2L 1 2. lim t!1 e ! = lim t!1 (!! d ) = 0 3. e remains bounded during the aerobatic ight Note that the translational dynamics (5.1) is under-actuated and coupled with the attitude dynamics, therefore during the maneuver the translational motion can not be controlled by the thrust force due to its extremely rapid and periodic rotation and, thus, stays in open-loop. 5.2 Adaptive Control Design and Stability Analysis In this chapter, we present two dierent adaptive control methods to achieve the control objectives along with the corresponding stability analyses. The rst is a Lyapunov- based adaptive control method, where the controller and adaptation law are generated through a single Lyapunov candidate function; the second is a modular adaptive control method, where the controller and adaptation law are designed separately. 5.2.1 Lyapunov-Based Adaptive Control The desired angular velocity dynamics is given asJ _ ! d +! d J! d = d , where d is the desired applied torque and! d = ! d1 ! d2 ! d3 T is the reference angular velocity. By design we enforce that ! d ; d 2L 1 . The reference angular velocity satises the condition in Denition 5:1. As already shown in Chapter 2.2.3, f(!) =!J! can be exactly expressed as f(!) =f(! d ) +D f (! d )e ! +g(t;!); (5.11) where D f (! d ) = @f(!) @! !=! d ; e ! = !! d ; and g(t;!) is the high-order nonlinear residue of the expansion, which has the compact form g(t;!) =h(e ! ) =e ! Je ! . 98 Then, by subtracting the desired angular velocity dynamics from (5.3), we obtain J _ e ! +D f (! d )e ! +h(e ! ) = d : (5.12) The LTI controller dened by (2.12) is used as the virtual torque input for (5:3) with the form a =K 1 e ! +v + d ; (5.13) _ v =v +K 2 e ! ; (5.14) where K 1 and K 2 are positive denite diagonal matrices such that K 1 K 2 and > 0. The estimated parameter vectors are dened to be ^ 1 and ^ 2 , and the estimation errors are ~ 1 = 1 ^ 1 and ~ 2 = 2 ^ 2 . The control Lyapunov function is dened as V 1 = 1 2 e T ! Je ! + 1 2 v T Q 1 v + 1 2 ( a ) T Q 2 ( a ) + 1 2 ~ T 1 K 1 1 ~ 1 + 1 2 ~ T 2 K 1 2 ~ 2 ; (5.15) where matricesQ 1 ;Q 2 ;K 1 ;K 2 are diagonal positive denite matrices. Thus, taking its derivative along the trajectories of (5.12)(5.14) and substituting the virtual control input specied by (5.13) into the result, we obtain _ V 1 = e T ! v T G(t) 2 6 4 e ! v 3 7 5 +e T ! ( a ) + ( a ) T Q 2 (a(t) +b(t)u _ a ) ~ T 1 K 1 1 _ ^ 1 ~ T 2 K 1 2 _ ^ 2 ; (5.16) 99 where G(t) = 2 6 4 K 1 + 1 2 (D f (! d ) +D T f (! d )) 1 2 (I +K T 2 Q T 1 ) 1 2 (I +Q 1 K 2 ) Q 1 3 7 5 : (5.17) One of the stability conditions is thatG(t) 0;8t> t 0 . In Chapters 2.2.3 and 4.2.1, G(t) was proven to stay in a convex hull dened by two constant symmetric matricesG 1 andG 2 . Then, the conditionG(t) 0 is equivalent to the conditions thatG 1 0 and G 2 0. The proof of this equivalence and more details have been given in Chapter 2.2.3. Finding K 1 , K 2 and that satisfy the conditions G 1 0 and G 2 0 is essentially semidenite programming that can be solved using the well-known methods [31]. The control inputu is designed to be u = ^ b 1 (^ a + _ a Q 1 2 e ! K 3 ( a )); (5.18) whereK 3 is a diagonal positive denite gain matrix, ^ b(t) is an estimation of b(t) with the form ^ b = ^ 4 + ^ 5 1 (! a ) + ^ 6 2 ( _ ! a ), and ^ a(t) is an estimation of a(t) with the form ^ a = ^ 1 + ^ 2 1 (! a ) + ^ 3 2 ( _ ! a ). The estimation errors ~ a = a(t) ^ a(t) = ~ T 1 (t) and ~ b = b(t) ^ b(t) = ~ T 2 (t). Note that the torque feedback is used in (5.18), and the observer for the actual torque will be given in Chapter 5.3.2. Substituting the control input given by (5.18) into (5.16), it follows that _ V 1 = e T ! v T G(t) 2 6 4 e ! v 3 7 5 ( a ) T Q 2 K 3 ( a ) ~ T 1 (( a ) T Q 2 +K 1 1 _ ^ 1 ) + ~ T 2 (( a ) T Q 2 uK 1 2 _ ^ 2 ): (5.19) The estimated ^ a and ^ b should be positive and satisfy the conditions that ^ a a min and ^ b b min , where a min ;b min > 0. To enforce this condition, n ^ 1 ; ^ 2 ; ^ 3 ; ^ 4 ; ^ 5 ; ^ 6 o 100 are bounded 1 min ^ 1 1 max ; 2 ^ 2 2 ; 3 ^ 3 3 ; (5.20) 4 min ^ 4 4 max ; 5 ^ 5 5 ; 6 ^ 6 6 ; (5.21) where 1 min ; 1 max ; 2 ; 3 ; 4 min ; 4 max ; 5 ; 6 > 0 and 1 min 1 2 2 3 a min ; (5.22) 4 min 1 5 2 6 b min ; (5.23) which guarantee that ^ a a min and ^ b b min . The parameter projections [57, 58] are applied to make sure that ^ 1 and ^ 2 satisfy the conditions (5.20)(5.21). Two convex sets are dened asS 1 = ^ 1 2R 3 j g 1 ( ^ 1 ) 0 andS 2 = ^ 2 2R 3 j g 2 ( ^ 2 ) 0 , where g 1 ( ^ 1 ) = 2 ^ 1 1 min 1 max 1 max 1 min ! 2n + ^ 2 2 ! 2n + ^ 3 3 ! 2n 1; (5.24) g 2 ( ^ 2 ) = 2 ^ 4 4 min 4 max 4 max 4 min ! 2n + ^ 5 5 ! 2n + ^ 6 6 ! 2n 1; (5.25) and n can be any positive integer. The adaptation laws for ^ 1 and ^ 2 are given by _ ^ 1 = 8 > < > : 1 ; if ^ 1 2 S 1 or if ^ 1 2@S 1 andrg T 1 1 0 1 K 1 rg 1 rg T 1 rg T 1 K 1 rg 1 1 ; otherwise (5.26) _ ^ 2 = 8 > < > : 2 ; if ^ 2 2 S 2 or if ^ 2 2@S 2 andrg T 2 2 0 2 K 2 rg 2 rg T 2 rg T 2 K 2 rg 2 2 ; otherwise (5.27) where 1 =K 1 ( a ) T Q 2 and 2 =K 2 ( a ) T Q 2 u. 101 When _ ^ = 1 and _ ^ 2 = 2 , we have that _ V 1 is negative denite given by _ V 1 = e T ! v T G(t) 2 6 4 e ! v 3 7 5 ( a ) T Q 2 K 3 ( a ): (5.28) Thus, considering that V 1 is positive denite and _ V 1 is negative denite, we have V 1 2 L 1 , which leads to thate ! ;v; a 2L 1 . Integrating both sides of (5.28) to obtain min t f min (G(t))g Z 1 0 (je ! j 2 2 +jvj 2 2 )d + min (Q 2 K 3 ) Z 1 0 j a j 2 2 d 6V 1 (0)V 1 (1): (5.29) BecauseV 1 (0) andV 1 (1) are bounded, we havee ! ;v; a 2L 2 . Additionally, from (5.4)(5.12)(5.13)(5.14)(5.18), we have _ e ! ; _ v; _ _ a 2L 1 . According to the corollary of Barbalat's lemma [56, p. 491], it can be shown that e ! ;v; a converge to zero asymptotically, which guarantees the stability of the closed-loop adaptive system. When ^ 1 2@S 1 withrg T 1 1 > 0 and ^ 2 2@S 2 withrg T 2 2 > 0, we have _ V 1 = e T ! v T G(t) 2 6 4 e ! v 3 7 5 ( a ) T Q 2 K 3 ( a ) + ~ T 1 rg 1 rg T 1 rg T 1 K 1 rg 1 1 + ~ T 2 rg 2 rg T 2 rg T 2 K 2 rg 2 2 ; (5.30) then it follows from ~ T 1 rg 1 < 0 and ~ T 2 rg 2 < 0 that _ V 1 e T ! v T G(t) 2 6 4 e ! v 3 7 5 ( a ) T Q 2 K 3 ( a ); (5.31) which indicates that _ V 1 remains negative denite. Following the same analyses as before and employing the corollary of Barbalat's lemma, it can be proved thate ! ;v; a converge to zero asymptotically. Therefore, the controller design (5.18) and adaptation 102 laws (5.26)(5.27) stabilize the closed-loop angular velocity error dynamics, and the tracking errore ! converges to zero globally asymptotically. 5.2.2 Modular Adaptive Control The adaptation laws (5.26)(5.27) in the Lyapunov-based adaptive method is constrained by the candidate Lyapunov function. In order to remove this constraint, we propose a modular adaptive control method in which the adaptation law is designed separately from the controller design. And the rst step is to propose a controller design that guarantees the closed-loop system is input-to-state stability (ISS) from the estimation errors to the system states. As described in the previous chapter, the open-loop system is J _ e ! +D f (! d )e ! +h(e ! ) = d ; (5.32) _ =a(t) +b(t)u: (5.33) To begin with, we dene the candidate Lypunov function to be V 2 = 1 2 e T ! Je ! + 1 2 v T Q 1 v + 1 2 ( a ) T Q 3 ( a ); (5.34) whereQ 1 0 andQ 3 =q 3 I 0. Then, by taking its derivative along the trajectories of (5.32)(5.33) and substituting the virtual torque specied by (5.13) (5.14) into the result, we obtain _ V 2 = e T ! v T G(t) 2 6 4 e ! v 3 7 5 +e T ! ( a ) + ( a ) T Q 3 (a(t) +b(t)u _ a ): (5.35) 103 Dene the control input to be u = ^ b(t) 1 + sgn(b(t)) b min n T Q 3 ( a ) ^ b(t) 2 T Q 3 ( a ) o ; (5.36) where = ^ a(t) + _ a K 4 ( a ); (5.37) andK 4 is the diagonal positive denite matrix. Note that b(t) b min > 1 and denez = a . Then, substituting control input (5.36) into (5.35) and deningx = e T ! v T z T T , we get _ V 2 6x T H(t)x ( T Q 3 z + ~ a(t) 2 ) 2 ( ^ b(t) 1 T Q 3 z ~ b(t) 2 ) 2 + ~ a(t) 2 4 + ~ b(t) 2 4 ; (5.38) 6x T H(t)x + ~ a(t) 2 4 + ~ b(t) 2 4 ; (5.39) where H(t) = 2 6 6 6 6 4 G(t) 1 2 I 33 0 33 1 2 I 33 0 33 Q 3 K 4 3 7 7 7 7 5 ; (5.40) ~ a(t) =a(t) ^ a(t), and ~ b(t) =b(t) ^ b(t). One of the stability conditions is thatH(t) 0;8t>t 0 . AndH(t) stays in a convex hull dened by two constant matrices because G(t) has the same property. Then, the stability condition that H(t) 0;8t > t 0 is equivalent to that the two constant matrices that dene the convex hull where H(t) stays inside are positive denite. For certain Q 1 , Q 3 and ! d (t), the problem to nd control gain matricesK 1 ;K 2 ;K 4 , and to satisfy the stability condition thatH(t) 0 104 is a semidente programming that can be solved using the method in [31]. Proposition 5.2. Given that Q 1 ;Q 3 ;H(t) 0 for8t > t 0 , the system (5.32)(5.33) with control input (5.36) is ISS with respect to ~ a; ~ b . Proof. From (5.39), it follows that _ V 2 6(1c h ) h jxj 2 (c h h jxj 2 1 4 (~ a 2 + ~ b 2 )); (5.41) where h = min t>t 0 ( min (H(t))) > 0, 0 < c h < 1. Also, it is straightforward to see that _ V 2 is negative denite wheneverjxj 2 > 1 4c h h (~ a 2 + ~ b 2 ). Therefore, from Theorem C.2 in [56] and the Claim on Page. 441 in [60], the proposition immediately follows. Note that the Proposition 5.2 indicates that the system statesfe ! ;v;zg remain bounded provided that the estimation errors ~ a and ~ b remain bounded. In addition, it follows from (5.5) and (5.6) that ~ a and ~ b remain bounded if ~ 1 and ~ 2 are bounded. Next, we proceed to design the adaptation law by substituting the control input (5.36) into (5.4) and obtain _ z =A z (t)z +W (t) T ~ (t); (5.42) where ~ = h ~ T 1 ~ T 2 i T ; A z = K 4 + ^ bsgn(b)q 3 b min ( T + ^ b 2 T ) ; W T = T ( sgn(b)q 3 b min ( T ^ b 2 T )z + ^ b 1 ) T : 105 It is straightforward that A z (t) is exponentially stable 2 because K 4 is positive denite and ^ bsgn(b)q 3 b min ( T + ^ b 2 T ) is positive semidenite since ^ b is enforced to remain positive by the parameter projection which will be discussed later in this chapter. Then, applying the Nonlinear Swapping Lemma [56] to (5.42) and obtaining the following lters _ T =A z T +W T ; (5.43) _ =A z +W T ^ ; (5.44) where ^ = h ^ T 1 ^ T 2 i T . From (5.42)(5.43)(5.44), we acquire the following static linear parametric model z + = T + ~ "; (5.45) where ~ " is governed by _ ~ " =A z ~ ", which indicates that ~ " decays to zero exponentially fast and ~ "2L 1 \L 2 because the matrixA z is exponentially stable. Then, the prediction error is given by T + ~ " T ^ = T ~ + ~ " =z + T ^ =": (5.46) Since the estimator error " is linear in ~ , this allows us to use either gradient algorithm or least-squares algorithm to estimate . In this work, we implement the Recursive Least-Squares (RLS) algorithm [58], namely, _ ^ =Proj P " 1 +tr( T P ) ; (5.47) 2 This means that there exists a positive denite symmetric matrixR such thatRA z (t)+A z (t) T R 0;8t>t 0 . 106 _ P = 8 > < > : P P T P 1+tr( T P ) ; ifkPk F 6P max 0; otherwise ; (5.48) where ;;P max > 0 and Proj() represents the projection operator in [56] given by Projf%g = 8 > > > > < > > > > : %; when ^ 2 S m orrg T m %6 0 (Ic( ^ )K m rgmrg T m rg T m Kmrgm )%; when ^ 2S m n S m andrg T m %> 0 where K m is positive denite, c( ^ ) = min 1; gm( ^ ) = ,S m = ^ 2R 6 j g m ( ^ ) 0 ,S m = ^ 2R 6 jg m ( ^ ) , and g m ( ^ ) = 2 ^ 1 1 min 1 max 1 max 1 min ! 2n + ^ 2 2 ! 2n + ^ 3 3 ! 2n + 2 ^ 4 4 min 4 max 4 max 4 min ! 2n + ^ 5 5 ! 2n + ^ 6 6 ! 2n 1; (5.49) where n can be any positive integer. The properties of the RLS adaptation law with the parameter projection can be summarized in the following proposition. Proposition 5.3. Given that Q 1 ;Q 3 ;H(t) 0 for8t > t 0 , the system dened by (5.32) and (5.33) with the control input given by (5.36) and the adaptation laws (5.47) and (5.48) guarantee the following properties: (i) ~ 2L 1 and"; _ ^ 2L 1 \L 2 (ii) e ! ;v;z2L 1 . Proof. Dene the Lyapunov candidate function to be V 3 = ~ T P 1 ~ + 1 2c 0 ~ " T ~ ", where c 0 = min (K 3 ). First, consider the case in which _ P6= 0. Take the derivative of V 3 and 107 use the property that ~ T P 1 Proj(%)6 ~ T P 1 % [56, Lemma E.1] to obtain _ V 3 6 ~ T (2" + T ~ ) 1 +tr( T P ) ~ T P 1 ~ ~ " T ~ "; (5.50) and substitute T ~ =" ~ " into above equation to obtain _ V 3 6 " T " 1 +tr( T P ) ~ T P 1 ~ : (5.51) When _ P = 0, we have _ V 3 6 2(" ~ ") T " 1 +tr( T P ) ~ " T ~ "; 6 " T " 1 +tr( T P ) " q 1 +tr( T P ) ~ " 2 2 6 " T " 1 +tr( T P ) : (5.52) Therefore, in the view of positive deniteness of P 1 , it follows that ~ 2L 1 and (1 + tr( T P )) 1 2 " 2 L 2 . Then, the boundedness of ~ implies that ~ a; ~ b 2 L 1 . Directly from Proposition 5.2, it follows thate ! ;v;z2L 1 , which, therefore, implies thatW; 2L 1 . Then, it follows from the following inequality: Z 1 0 j"j 2 2 dt6 1 +tr( T P ) 1 Z 1 0 j"j 2 2 1 +tr( T P ) <1; (5.53) that"2L 2 . Additionally, we have _ ^ 2L 2 by considering the fact that the projection algorithm in (5.47) has the property that Proj(%) T K 1 m Proj(%)6% T K 1 m %,8 ^ 2S m , and"2L 2 and the boundedness of P 1+tr( T P ) . Last, consider the above results and the relationships of " = T ~ + ~ " and (5.47), we can obtain that"; _ ^ 2L 1 . 108 Proposition 5.4. All the signals in the closed-loop adaptive system consisting of the plant given by (5.32) and (5.33), controller (5.36), lters (5.43) and (5.44), and adapta- tion laws (5.47) and (5.48) are globally uniformly bounded if the parametersK 1 ,K 2 , K 4 , Q 1 , Q 3 , and reference signal! d (t) satisfy the conditions that Q 1 ;Q 3 ;H(t) 0;8t > t 0 . Moreover, lim t!1 z(t) = lim t!1 v(t) = 0, and the global asymptotic tracking is achieved: lim t!1 e ! (t) = lim t!1 (!(t)! d (t)) = 0: (5.54) Proof. From (5.39) and (5.41), it follows thatjxj 2 6 max n jx(t 0 )j 2 ; 1 4c h h (~ a 2 + ~ b 2 ) o . Then, consider ~ a 2 + ~ b 2 6 (1 + 2 1 + 2 2 )max ~ 1 2 2 + ~ 2 2 2 , where max ~ 1 2 2 + ~ 2 2 2 can be bounded by a uniform constant due to the boundedness of ^ . From applying the Lemma F.4 in [56], it follows thatz T ~ 2L 2 . Then, we also have T ~ 2L 2 from the facts that T ~ = " ~ " and "; ~ "2L 2 . Therefore, we can obtain thatz2L 2 . Additionally, we have _ z2L 1 from the fact thatW T ~ ;A z z2L 1 in (5.42). By utilizing Barbalat's lemma, we havez! 0 as t!1. Additionally, we can have e ! ;v2L 2 through the strictly passive property of the input-to-output system fromz toe ! ;v, which can be shown using the Lyapunov func- tion V 4 = 1 2 e T ! Je ! + 1 2 v T Q 1 v. Take the derivative of V 4 and obtain _ V 4 6 g (je ! j 2 2 +jvj 2 2 ) +e T ! z 6 g jvj 2 2 (1c 1 ) g je ! j 2 2 + 1 4c 1 g jzj 2 2 ; (5.55) where g = min t>t 0 f min (G(t))g and 0<c 1 < 1. Integrating both sides of (5.55) to obtain g Z 1 0 jvj 2 2 d + (1c 1 ) g Z 1 0 je ! j 2 2 d 6V 4 (0)V 4 (1) + 1 4c 1 g Z 1 0 jzj 2 2 d: (5.56) Since z2L 2 and V 4 (1) is bounded due to the boundedness of e ! and v, we have 109 e ! ;v2L 2 from the above equation. Additionally, we have _ e ! ; _ v2L 1 , which can be seen from (5.32) and (5.14). Then, from directly applying Barbalat's lemma, we have e ! ;v! 0 as t!1. 5.2.3 Attitude Analysis In the previous chapters, we present two dierent adaptive control methods to enable VTOL UAVs to achieve global asymptotic tracking of the aerobatic reference angu- lar velocity in the presence of the LTV actuator dynamics with unknown parameters. Next, we show the attitude error remains bounded during the aerobatic maneuver when performed with the adaptive schemes. As explained in Chapter 5.1.2, e denotes the rotation angle fromB toI. The corresponding rotation axis is specied by the unit vectora e with its components ex- pressed in bothB andI. Then, from the denition of the attitude error quaternion, we have that e 0 = cos 1 2 e ande 1 = sin 1 2 e a e . Thus, it follows from (5.10) that _ e 0 = 1 2 _ e sin 1 2 e = 1 2 sin 1 2 e a T e (! ! d ); (5.57) which, under the assumption that e 6= 0, yields _ e =a T e (! ! d ) =a T e e ! ; (5.58) where we employ the property that! d and ! d have the same components. From (5.58), recalling thatja e j 2 = 1, we get the following inequality _ e 6 _ e 6je ! j 2 : (5.59) 110 Then, integrating both sides of the above equation, we obtain e (t)6 e (t 0 ) + Z t t 0 je ! j 2 d; (5.60) wheret 0 6t6t f ,t 0 andt f denote the start and end times of the aerobatic maneuver, and e (t 0 ) indicates the initial attitude error angle, which can be simply regulated to zero during the climb phase since the attitude stays near the hovering state. It is straightforward that the attitude error, e (t), remains bounded during the aerobatic maneuver considering thate ! is bounded and converges to zero asymptotically; there- fore, the attitude error kinematics does not have a nite time escape. Additionally, the total duration of the class of aerobatic maneuvers studied in this paper is normally less than one second. It is important to make sure that e (t f ) at the end of the maneuver is bounded such that the vehicle can restabilize its attitude and stop descending rapidly in the DR phase. In addition, from the experimental results shown in Chapter 5.3, the attitude error e (t f ) is reduced to a small value by employing the adaptive control strategies because the angular velocity oscillation and overshoot caused by aerodynamic eects during the aerobatic maneuver phase are inhibited. The following quantity = Z t f t 0 je ! j 2 dt (5.61) indicates the upper bound of the attitude error and is also an important performance gure of merit (PFM) to evaluate the aerobatic maneuver performance as explained in Chapter 3, and the smaller represents the better performance. 111 5.3 Experimental Results 5.3.1 Experimental Setup and Example Aerobatic Maneuvers The quadrotor UAV employed in the performance of the ight experiments is the same one used in the experiments described in Chapter 4, namely, the Crazy ie 2.0 as shown in Fig 4.1. The robot's weight is 28 g including the battery, has a propeller-tip-to- propeller-tip distance of 14 cm and exhibits a maximum thrust-to-weight ratio of 1.9. The quadrotor UAV is also equipped with an gyroscope and an accelerometer. During the execution of an aerobatic maneuver, the quadrotor UAV operates completely au- tonomously. More details about the experimental architecture can be found in Chapter 2.3.1. We implement both the Lyapunov-based adaptive controller and the modular adap- tive controller on the quadrotor UAV to perform three dierent sets of aerobatic ma- neuvers, namely, triple- ip, Pugachev's Cobra and mixed ip. In addition, the linear controller in (2.12) that does not consider the LTV actuator dynamics is implemented in experiments for the purpose of comparison. The details of the reference signals of these three maneuvers have been already described in Chapter 4.3.2. 5.3.2 Torque Feedback and Control Implementation The torques are acquired using =J _ ! +!J!, in which the angular accelerations, _ !, are estimated employing the nonlinear derivative estimation method [61] based on the measurement of angular velocity,!. The quadrotor UAV shown in Fig 4.1 has a miniature size and limited on-board computation resources. In order to reduce the computational cost and increase the execution frequency of the control loop, the parameter projections for the real-time 112 implementation of both control schemes are simplied to be ^ = 8 > > > > < > > > > : ^ ; if ^ min ^ ^ max ; ^ min ; if ^ ^ min ; ^ max ; if ^ max ^ ; where ^ min = 1 min - 2 - 3 4 min - 5 - 6 T and ^ max = 1 max 2 3 4 max 5 6 T . In addition, in order to avoid saturation of actuators, the adaptive control inputs given by (5.18) and (5.36) of both methods are implemented only on the axis that tracks the non-zero reference signals, and the linear controller (2.12) is implemented on the other two axes that track the zero reference signals. This implementation is reasonable because the axes that track zero reference signals encounter the miniature eects of the aerodynamic coecients and torque latency variations due to the approximately zero angular velocity. 5.3.3 Results In experiments, we implement both the Lyapunov-based and modular adaptive control methods on the quadrotor UAV to perform three dierent aerobatic maneuvers, and the linear controller (2.12) is implemented as comparison. The experimental results obtained with the three dierent controllers for three dierent maneuvers are demon- strated in Figs 5.1, 5.2 and 5.3, respectively. It is clear that the aerobatic performance with the linear controller displays undesired angular velocity oscillations and signi- cant overshoots, which we think are caused by the aerodynamic coecients and torque latency variations. The linear controller (2.12) is designed based on the nominal dy- namics without considering the LTV actuator dynamics as described in Chapter 2. In contrast, the aerobatic performances obtained with both the Lyapunov-based and mod- ular adaptive control schemes in all the three dierent maneuvers improve signicantly with the undesired angular velocity oscillations and overshoots greatly reduced. The 113 A: Linear control. B: Lyapunov-based adaptive control. C: Modular adaptive control. Figure 5.1: Triple ip using three dierent controllers. (A) Triple ip using the linear controller given by (2.12) that ignores LTV actuator dynamics; the corresponding PFM = 0:9391. As! d2 =! d3 = 0 for the triple ip, they are not plotted here. (B) Triple ip using the Lyapunov-based adaptive controller (5.18); the corresponding PFM = 0:6397. (C) Triple ip using the modular adaptive controller (5.36); the corresponding PFM = 0:6286. corresponding PFMs, , of three dierent maneuvers decrease signicantly compared to the PFMs obtained with the linear controller, in specic, the Lyapunov-based adaptive control reduces the PFMs of the three maneuvers by 31:88%, 37:61% and 40:63%, and the modular adaptive control by 33:06%, 34:06% and 39:81%, respectively. Because the LTV actuator dynamics is taken into account when designing the adaptive controllers, both the Lyapunov-based and modular adaptive control schemes greatly improve the aerobatic performances in the three dierent maneuvers compared to the linear con- troller's performance, which also indirectly veries the accuracy and correctness of the LTV actuator model. The aerobatic performances of the Lyapunov-based and modular adaptive control methods are similar to each other. However, there are many dierences from the im- 114 A: Linear control. B: Lyapunov-based adaptive control. C: Modular adaptive control. Figure 5.2: Pugachev's Cobra using three dierent controllers. (A) Pugachev's Cobra using the linear controller that ignores LTV actuator dynamics; the corresponding PFM = 0:7413. As! d2 =! d3 = 0 for Pugachev's Cobra, they are not plotted here. (B) Pugachev's Cobra using the Lyapunov-based adaptive controller; the corresponding PFM = 0:4625. (C) Pugachev's Cobra using the modular adaptive controller; the corresponding PFM = 0:4888. plementation perspective. In general, the Lyapunov-based adaptive control strategy is computationally much lighter than the modular adaptive control method, and this can be seen by less total number of variables and less computational complexity of the Lyapunov-based method. However, the modular adaptive control strategy has the exibility of selecting the adaptation laws; in contrast, the structure of the adaptation law of the Lyapunov-based method is xed. To demonstrate the performance robustness of the proposed adaptive strategies, we select two examples of consecutive aerobatic maneuvers. Fig 5.4a demonstrates that the quadrotor UAV is able to perform 11 consecutive Pugachev's Cobra maneuvers continuously employing the Lyapunov-based adaptive control method without crashing or stalling until the battery power is completely depleted. In Fig 5.4b, the quadrotor 115 A: Linear control. B:Lyapunov-based adaptive control. C:Modular adaptive control. Figure 5.3: Mixed ip using three dierent controllers. (A) Mixed ip using the linear controller that ignores LTV actuator dynamics; the corresponding PFM = 1:1497. As ! d3 = 0 for mixed ip, it is not plotted here. (B) Mixed ip using the Lyapunov-based adaptive controller; the corresponding PFM = 0:6826. (C) Mixed ip using the modular adaptive controller; the corresponding PFM = 0:6920. A:Lyapunov-based adaptive control. B: Modular adaptive control. Figure 5.4: Consecutive aerobatic maneuvers. (A) 11 consecutive Pugachev's Cobra maneuvers em- ploying the Lyapunov-based adaptive control strategy. (B) 12 consecutive mixed ip maneuvers em- ploying the modular adaptive control strategy. UAV executes 12 consecutive mixed ip maneuvers continuously employing the mod- ular adaptive control method without crashing or stalling until the battery power is 116 completely depleted. A video of the experimental results with both adaptive control schemes employed for aerobatic maneuvers is available at http://www.uscamsl.com/resources/XWSM3.mp4 117 Chapter 6 Flight Control of Insect-Scale Flapping-Wing Robots In this chapter, the ight control vehicle is changed from the micro quadrotor UAV to two dierent types of insect-scale apping-wing robots. The insect-scale size of the aerial robots make them capable of achieving the agility comparable to the insect aerial animals in nature, such as the fruit y (Drosophila melanogaster) [4]. Also, it is known that the maximum angular acceleration the UAV can attain is inversely proportional to the size of the vehicle [41]. Therefore, the insect-scale aerial robot is potentially more agile and maneuverable than the quadrotor UAV. However, due to the extremely small size, the manufacturing and assembling of the robots are dicult, and the mechanical robustness is also not as high as that of the quadrotor UAV. Our ultimate goal is to enable the insect-scale apping-wing robots to perform similar aerobatic maneuvers to those shown in the previous chapters. But before that, the insect-scale aerial robots should rst be able to perform regular ight in a robust and stable way. We propose a set of position and attitude controllers for two dierent types of insect-scale apping-wing robots; one is equipped with two wings, and the other one generates the lift force with four wings. We will show that the insect-scale 118 Figure 6.1: Two dierent insect-scale apping-wing robotic systems. RoboBee (left) is a two-winged robot, and it is driven by two piezoelectric bimorph actuators. Bee + (right) is a four-winged robot that has a mass of 95 mg and measures 33 mm in wingspan, and it is driven by four piezoelectric unimorph actuators. A U.S. one-cent coin indicates the scale. Table 6.1: Comparison of parameters for Bee + and RoboBee Robot Total mass (mg) Actuators mass (mg) Wingspan (mm) Flapping frequency (Hz) Lift force (mN) Wing area (mm 2 ) Bee + 95 56 33 100 1.4 200 RoboBee 75 50 35 120 1.3 104 apping-wing robots can be thought of as a single rigid body with the direction of the thrust force xed in the body frame. A similar control strategy for the quadrotor UAV can be utilized to control the insect-scale aerial robots. 6.1 Insect-Scale Flapping-Wing Robotic Systems 6.1.1 Two Systems Two insect-scale apping-wing robots are utilized in this chapter, one is called Robobee, and the other one is called Bee + as shown in Fig 6.1. Robobee is fabricated and assembled based on the design in [25] that is driven by two bimorph piezoelectric actuators. Bee + is developed at our lab, Autonomous Microrobotic Systems Laboratory (AMSL) at USC; it is driven by four unimorph piezoelectric actuators that can ap four wings independently. The parameters of the two robots are listed in Table 6.1. 119 For both robots, the direction of the thrust force is xed upwards in the body frame. However, the dierent mechanism and number of wings decide the dierence of the torque generation principles. The torque generation of the RoboBee aerial robot is demonstrated in Fig 6.2. The roll torque is generated by diering the magnitudes of two wings, the pitch torque is generated by shifting the apping center position upwards or downwards, and the yaw torque is generated by using an asymmetric apping pattern, which means that the wing aps faster in one direction than in the opposite direction. In experiments, the RoboBee robot is able to generate anticipated amounts of roll and pitch torques using this method. The yaw torque generation method requires the wing to ap much faster in one direction to generate unbalanced horizontal forces in one cycle, which inevitably requires a much higher bandwidth of the actuation system than the apping frequency. Since the apping frequency is already tunned at the resonant frequency of the system, a much higher apping frequency will decrease the apping magnitude dramatically. Therefore, no signicant yaw torque is observed in the experiments, which is consistent with the observations in [29,62,63]. The principle of the torque generation of the Bee + robot is demonstrated in Fig 6.3. The roll torque is created by diering the magnitudes of two opposed pairs of wings, and the pitch torque is created by diering the magnitudes of two pairs of wings on the same side. To generate the yaw torque, the stroke plane of the robot is pre-set and manufactured to have an inclined angle with respect to the steering plane, in which the lift force of each apping wing has a component projected on the steering plane. Then, the yaw torque can be generated by diering the projected forces of two diagonally opposed pairs of lift forces as shown in Fig 6.3c. In this way, the bandwidth constraint of the two-wing design on the yaw torque generation disappears, and the yaw torque is proportional to the lift force. Note that the four-wing design of Bee + makes its torque generation principle very similar to that of the quadrotor UAV [33]. 120 (a) Roll Torque. (b) Pitch Torque. (c) Yaw Torque. Figure 6.2: The principle of torque generation of RoboBee (top view of the robot). In (a), the roll torque is generated by diering the magnitudes of two wings. In (b), the pitch torque is generated by shifting the apping center position. In (c), the yaw torque is generated by using an asymmetric apping pattern. All gures are top view of the robot. 6.1.2 System Dynamics Both the robots can be modeled as single rigid body as the quadrotor UAV in the previous chapters. The frames of references of both insect-scale apping-wing robots are demonstrated in Fig 6.4, andfn 1 ;n 2 ;n 3 g represent the inertial frame andfb 1 ;b 2 ;b 3 g represents the body-xed frame with the origin located at the center of mass of the robot. For both robots, the direction of the total thrust force is assumed to be aligned with the b 3 axis, and the number of actuators is less than the total number of the degree of freedoms of the system. Therefore, both robots are essentially thrust-propelled underactuated systems. The globally dened nonlinear dynamics for both robots are given by m r =mgn 3 +fb 3 ; (6.1) _ q = 1 2 q p; (6.2) J _ ! =!J! +; (6.3) 121 Roll (a) Roll Torque. Pitch (b) Pitch Torque. Yaw b 3 (c) Yaw Torque. Figure 6.3: The principle of torque generation of Bee + . In (a), the roll torque is generated by diering the magnitudes of two opposed pairs of apping wings. In (b), the pitch torque is generated by diering the magnitudes of two pairs of apping wings that are on the same side. In (c), the stroke plane of the wing apping is pre-set or adjusted to have an angle with respect to the steering plan, and the yaw torque is generated by diering the projected forces of two diagonally opposed pairs of lift forces of two apping wings onto the steering plane. where m is the total mass of the robot; r = [r 1 r 2 r 3 ] T indicates the location of the robot's center of mass measured from the origin ofN ; f is the magnitude of the total thrust force generated by the apping wings; J denotes the robot's moment of the inertia, written with respect to B; ! is the yer's angular velocity with respect to N , expressed inB; = [ 1 2 3 ] T is the torque generated by the apping wings; the quaternion q describes the attitude of the robot relative toN ; p = 0 ! T T ; and the symbol denotes the standard quaternion multiplication. Note that the model specied by (6.1){(6.3) assumes that the direction of the thrust force is aligned withb 3 ; that the projection of the total aerodynamic force generated by the apping wings onto the steering plane is zero during one apping cycle, which implies that fb 3 is the only external actuation force acting on the system; that the aerodynamic disturbances aecting the yer are negligible in (6.1)(6.2); and that the gyroscopic eect resulting from the interaction of the apping wings with the rotating 122 (a) RoboBee and Frames of Reference. Roll axis Yaw axis Pitch axis 1 2 4 Center of mass 3 n3 n1 n2 b 3 b 1 b 2 r (b) Bee + and Frames of Reference. Figure 6.4: Frames of references of two insect-scale apping-wing robots. In (a), the frames of references of the RoboBee is indicated. In (b), the frames of references of the Bee + is indicated. For both robots, fn 1 ;n 2 ;n 3 g is the inertial frame, andfb 1 ;b 2 ;b 3 g represents the body-xed frame with the origin located at the center of mass of the robot. body is also negligible. Note that the model (6.1) for the Bee + robot ignores the eects of the yaw torque generation mechanism, because the experiments are conducted using the initial version of the Bee + that does not have the inclined angle between the stroking plane and steering plane. 6.1.3 Actuation Command Generation Each piezoelectric actuator of both robots is driven by a high-frequency sinusoid signal, and the mapping from the desired thrust force and torque to the magnitude of the sinusoid signal is important for the robots to have desired performance. The mapping for the RoboBee robot follows the same strategy as described in [26]. In this chapter, we focus on the mapping for the Bee + robot. By simplifying the model described in [64], we estimate the magnitude of the thrust force produced by each apping wing i =f1; 2; 3; 4g, according to f i = k f v i , where f i is the thrust force generated by the apping of the i-th wing, v i is the amplitude of the sinusoidal command signal generated by the i-th unimorph actuator, and k f is a lumped thrust force coecient. 123 Yaw torques in the steering plane can be generated by employing the strategy men- tioned in Chapter 6.1.1. Consistently, we estimate the component of the ith aerody- namic force projected on the steering plane as f s i =k s v i , with i = 1; 2; 3; 4, where k s is also a lumped force coecient. Thereby, the mapping that relates the amplitudes of the actuators' outputs, as inputs, with the total thrust force and control torques, as outputs, is given by u = v; (6.4) with u = f 1 2 3 T = 2 6 6 6 6 6 6 6 4 k f k f k f k f k f d 1 k f d 1 k f d 1 k f d 1 k f d 2 k f d 2 k f d 2 k f d 2 k s d 3 k s d 3 k s d 3 k s d 3 3 7 7 7 7 7 7 7 5 and v = v 1 v 2 v 3 v 4 T whered j , forj =f1; 2; 3g, is the equivalent lever-arm associated with the corresponding torque component j . It is straightforward that the mapping in (6.4) is similar to the mapping between the rotors' speeds and the control force and torques of the quadrotor [33]. Additionally, the four-winged design brings more control capability than the two- winged design [26] since the thrust force and the control torques are shared by four wings rather than two. Taking the inverse of (6.4) yields the mapping from the thrust force and torques to the actuator commands as v = 1 u; (6.5) 124 with 1 = 2 6 6 6 6 6 6 6 4 1 4k f 1 4d 1 k f 1 4d 2 k f 1 4d 3 ks 1 4k f 1 4d 1 k f 1 4d 2 k f 1 4d 3 ks 1 4k f 1 4d 1 k f 1 4d 2 k f 1 4d 3 ks 1 4k f 1 4d 1 k f 1 4d 2 k f 1 4d 3 ks 3 7 7 7 7 7 7 7 5 : 6.2 Attitude and Position Control Design 6.2.1 Attitude Control Here, we describe the desired attitude dynamics of the robot with the quaternion equa- tion _ q d = 1 2 q d p d ; (6.6) where q d is the quaternion employed to represent the desired attitude of the yer during ight; and p d = 0 ^ ! T d T , in which ^ ! d denotes the desired angular velocity expressed in the desired frame of reference. Consistently, it follows that the attitude error between q and q d , described by the quaternion q e = m e n T e T , is given by q e = q 1 d q: (6.7) In this case, we regulate the attitude of the yer employing the control law =K 1 sgn(m e )n e K 2 (!! d ) + d ; (6.8) whereK 1 andK 2 are positive denite diagonal gain matrices, sgn() represents the sign function,! d denotes the desired angular velocity that has the exactly same components as ^ ! d and is expressed in the body frame instead of the desired frame q d , and d 125 represents the desired torque of the reference angular velocity dynamics. Denote the axis of the rotation from q to q d by a unit vectora e , and the associated rotation angle is dened to be e , with 0 6 e < . Then, the termsgn(m e )n e is geometrically equal to sin ( 1 2 e )a e , and the multiplication of sgn(m e ) is employed to remove the ambiguity of the quaternion representation since q e and q e indicate the same rotation result. 6.2.2 Position Control As the dynamics of Bee + and RoboBee are essentially an underactuated system with the direction of thrust force aligned with b 3 axis, the position control of this system involves specifying the magnitude and the orientation of the thrust force, similar to other apping-wing robots and the quadrotors [26] [17]. Therefore, the position control design is comprised of the generations of magnitude of the thrust force and the desired attitude. In this case, the desired instantaneous total thrust force required to track a desired position of the yer's center of mass, r d , is generated according to the proportional- integral-derivative (PID) structure f d =K p (rr d )K d ( _ r _ r d )K i Z (rr d )dt +mgn 3 +m r d ; (6.9) whereK p ;K d ;K i are positive denite diagonal gain matrices. Note that the magnitude of the thrust force control input is simply given by f =f T d b 3 ; (6.10) 126 The desired direction of the b 3 axis compatible with the direction of f d can be readily computed as i 3 = f d jf d j 2 (6.11) which is chosen to be the desired yaw axis of the robot during ight. We compute the other two axes dening the desired attitude of the robot in terms of the desired instantaneous yaw rotation angle, d , andi 3 , according to i 1 = i d i 3 i d i 3 2 ; (6.12) i 2 =i 3 i 1 ; (6.13) where i d = [ sin d cos d 0] T . In the implementation of the algorithms for signal processing and control, i 1 , i 2 , i 3 and i d are expressed in the inertial frame. To im- plement the controller specied by (6.8), we compute the desired attitude quaternion q d from the desired rotation matrix S d = [i 1 i 2 i 3 ], employing standard quaternion algebra. 6.3 Controlled Flight Experiments 6.3.1 Experimental Setup The experiment architecture includes the insect-scale apping-wing robot, piezo-actuator drivers (PiezoMaster VP7206), a Vicon motion capture (VMC) system and a ground computer for processing data and generating real-time control signals. Matlab Simulink Real-Time is used for processing the sensor measurements and generating the control signal with the control algorithm embedded in the ground computer. The control al- gorithm runs at a frequency of 2 kHz and the VMC system measures the position and 127 attitude states at a frequency of 500 Hz. The angular velocity cannot be directly mea- sured with the motion capture system, thereby is estimated as 2 6 4 0 ! 3 7 5 = 2 q 1 s s + q; (6.14) where is the lter parameter, s is the complex variable of the Laplace transform, the bracket [] represents the lter that operates on the signal q. A similar low pass deriva- tive lter operates on the position states to estimate the translational velocities. Note that the low pass lter is necessary for this aerial robotic system with high apping fre- quency since the inertial force induced by the wing apping introduces the unavoidable high magnitude oscillation of the robot body. The open-loop trimming ight tests required in [26] [28] are not needed in the controlled ights of the four-winged ying robots proposed in this work, implying that the control algorithm introduced in this work does not need the tuning of the command signals for the zero oset torques. This advantage signicantly improves the eciency of experiments. 6.3.2 Position Control of RoboBee We implement the attitude and position controllers (6.8)(6.10)(6.9) in experiments to make the RoboBee insect-scale apping-wing robot hover at a desired position. In the experiments, we control the roll and pitch rotations to change the direction of the lift force so that the position control is achieved. The yaw rotation is in open-loop due to the lack of the eective yaw torque generation. The position of the RoboBee insect-scale apping-wing robot during the hovering control experiment is plotted in Fig. 6.5. In the experiment, the apping-wing robot successfully takes o and hovers about a desired position. And it is obvious that there are about 10 cm position errors in the n 1 and n 2 directions, and the altitude error 128 (a)n 1 position. (b)n 2 position. (c)n 3 position. Figure 6.5: The position control of the RoboBee insect-scale apping-wing robot. In (a), the n 1 position of the robot, the green line is the reference signal. In (b), the n 2 position of the robot and green line is the reference. In (c),n 3 position (altitude) of the robot, and the green line is the reference signal. is about 5 cm. Despite the position errors, the hovering experiment of the RoboBee demonstrates that the proposed attitude and position control strategy can achieve the control objective in a satisfying way, and the insect-scale apping-wing robot can per- form the controlled maneuver from take-o, hovering and landing without crashing or stalling. The Euler roll and pitch angles during the hovering experiment are illustrate in Fig. 6.6. A probable reason for the position errors is due to the essence of the extreme under- actuation of the system. There are only two wings in this apping robot instead of four actuators in the Bee + and the quadrotor. Each wing takes both the position and attitude control responsibilities, which makes the control margin very small. For example, the apping magnitude of one side wing has to decrease in order to generate an attitude control torque, which will decrease the total thrust force. This is another 129 Figure 6.6: The roll and pitch angle of the RoboBee during the hovering. The blue line indicates the Euler Roll angle of the apping-wing robot, and the red line lines indicates the Euler Pitch angle of the apping-wing robot. advantage of the four-wing design of Bee + , which can generate enough lift forces for both position and attitude controls with relatively larger control authority. The video of the hovering experiment of the insect-scale apping-wing robot, RoboBee, is available at http://www.uscamsl.com/resources/robotbeevideo.mp4 6.3.3 Altitude and Attitude Control of Bee + The objective of the altitude and attitude control experiment is to control the Bee + robot to y at a desired altitude with the direction of the thrust force remaining per- pendicular to the n 1 -n 2 plane. In this preliminary control experiment, yaw feedback control is not applied to reduce the control burden of the four apping wings, however the increased aerodynamic damping forces improve the open-loop stability of the yaw motion as discussed in [64]. Additionally, the direction of the thrust force (i.e. b 3 ) is theoretically irrelevant with the Euler yaw angle. Therefore, the control objective in the altitude and attitude control experiment is to regulate the altitude to a desired value, and to regulate the Euler roll and pitch angles to zero. It follows that the desired attitude quaternion in the altitude and attitude control experiment is given by q d = cos 2 0 0 sin 2 T ; (6.15) 130 The reference and measured altitudes. The Euler roll and pitch angles. Figure 6.7: The altitude and attitude control experimental of the Bee + . A. It demonstrates the altitude in the experiment, in which the dash line represents the reference altitude. B. It shows the Euler roll and pitch angles during the experiment. Note that the angular oscillation remain approximately between10° and 10°, which is partially caused by the oscillation of the robot body induced by the wing apping and stays in an acceptable range. The experiment lasts for approximately 5 s, and after that the robot ies out of the specied safety space and the power is then turned o automatically. where is the actual Euler yaw angle of Z-Y-X convention. The altitude controller can be simply derived from (6.10)(6.9) with the assumption thatb 3 b 3d as follows: f =k p (r 3 r d3 )k d _ r 3 k i Z (r 3 r d3 )dt +mg; (6.16) where r 3 is the third component ofr and r d3 is the desired altitude. The results of the altitude and attitude control experiment are illustrated in Fig. 6.7 including the desired and measured altitudes and the Euler angles. It is clear that the attitude control algorithm (6.8) enables the direction of the thrust force approxi- mately perpendicular to then 1 -n 2 plane. Also angular oscillation occurs in the range of [10 ; 10 ], which is acceptable in the sense that the magnitude of the lift force is not greatly reduced by the attitude oscillation as demonstrated in the experiment. Addi- tionally, the robot reaches the desired altitude rapidly as shown in Fig. 6.7-A despite the 131 Figure 6.8: A. The time-lapse plot of the altitude and attitude control experiment. The corresponding altitude and Euler angles are demonstrated in Fig. 6.7. During the experiment, the direction of the thrust force is controlled to remain approximately perpendicular to the n 1 -n 2 plane. The cable tethered to the robot provides the power and control signals. The robot drifts on then 1 -n 2 plane due to the lack of control of the position. After 2 s, the robots ies out of focus area of the camera. B. Image sequences of the position control experiment. The corresponding position and Euler angles of this experiment are plotted in Fig. 6.9. The images are captured with the identical view position. static altitude error. The time lapse of the corresponding experiment is demonstrated in Fig. 6.8-A, which indicates that the apping of the four wings generates enough lift force for the robot to take o and maintain its attitude upright for a long period of time. The drifting phenomenon is expected due to the lack of control on the n 1 -n 2 plane. Overall, the altitude and attitude control experiment veries the eectiveness of the design and fabrication of this four-winged insect-scale ying robot, and demonstrates that the robot is able to perform controlled ight for a long period of time. 6.3.4 Position Control of the Bee + In this experiment, the insect-scale ying robot, Bee + , is commanded to hover at a desired position in the space, in which the attitude controller (6.8) and the position controller (6.10)(6.9) are implemented. Same as the altitude and attitude control ex- periment, the desired yaw angle is set to be the true yaw angle, i.e. d = , which does not aect the position control since the direction of the thrust force is irrelevant to the yaw angle. The result of the position control experiment is demonstrated in Fig. 6.9, and the corresponding time lapse of the experiment is illustrated in Fig. 6.8-B. 132 The reference and measured positions. The Euler angles and the references. Figure 6.9: The position control experiment result. A. The dash lines represent the desired position, and the solid lines represent the measured position. B. The dash lines represent the desired Euler angles derived from (6.11)(6.12)(6.13). The solid lines represent the measured Euler angles. Fig. 6.9-A shows the controlled position along with the reference signals, from which the robot is observed to approximately track the reference signals in the rst second, then the position error on n 1 axis gradually increases. In Fig. 6.9-B, the roll and pitch angle tracking in the rst second are, to some extent, accurate, then the pitch tracking error gradually becomes large, which is consistent with the increase of the position error on n 1 axis in Fig. 6.9-A. It is hypothesized that the oscillation about pitch axis is caused by the saturation of the actuators. This problem can be solved by improving the robotic design using the method in [65] to generate more available thrust force for position regulating or trajectory following ight controls. The complete set of experimental video is available at https://www.uscamsl.com/resources/IROS2019/S1.mp4. 133 Chapter 7 Conclusions and Future Work Conclusions The capability of aerial animals to execute aerobatic ights is crucial for their survival in the wild environment. Similarly, articial aerial robots must have comparable aerobatic ying capabilities and accomplish similar tasks in order to operate in unstructured environments. The research presented in this dissertation contributes to reducing the gap in aerobatic ight capabilities between aerial animals and robots by introducing new methods for the synthesis of high-performance controllers and analytical tools to understand and explore the aerobatic potentials of ying robots. First, we presented a linear controller design to enable a quadrotor UAV to perform high-speed aerobatic maneuvers in less than one second without stalling or crashing to the ground. The proposed approach was validated through compelling experimental results. Then, rigorous stability analyses of the closed-loop systems were provided to establish the theoretical conditions for the proposed approach to work. This analytical result enabled the formulation of an optimization problem to improve the ight perfor- mance by nding the corresponding optimal parameter set. To solve the optimization problem, an algorithm based on quasiconvex programming was proposed; the result 134 is guaranteed to be a close approximation of the globally optimal. The eectiveness of the optimization approach was validated through both simulation and experimental results. From the experimental results of aerobatic ights obtained using the linear con- troller, we observed that the aerodynamic eects of the varying ow elds surrounding the yer have a noticeable in uence on the ight performance during the execution of an aerobatic maneuver. A LTV actuator dynamic model with unknown time-varying parameters is proposed to describe the eects of aerodynamic coecient and torque la- tency variations, and we developed a nonlinear adaptive control solution to compensate for the negative eects of these variations on the system. Experimental results show that the ight performance during an aerobatic maneuver signicantly improved with the use of the proposed adaptive scheme. In order to provide the necessary theoretical stability analyses, torque feedbacks are utilized in another two dierent types of adap- tive control solutions that also address the LTV actuator dynamics. The corresponding stability analyses show that the closed-loop systems are stabilized by both adaptive schemes. In addition, consistent experimental performance improvements are observed by implementing the proposed adaptive control solutions compared to that obtained with the linear controller solution. Finally, nonlinear attitude and position controllers were developed to enable two dierent types of insect-scale apping-wing robots to perform altitude and attitude and position controlled ights in experiments. Furthermore, the results obtained through controlled ight experiments serve as guidance for robotic design iterations. Overall, this research shows that the studied articial insect-scale ying robots have the potential to perform controlled aerobatic maneuvers provided that some hardware robustness and thrust force insuciency problems can be resolved. 135 Future Work The types of aerobatic maneuvers studied in this dissertation are characterized by high- speed rotational motion in an extremely short period of time. One line of future work is the development of new types of aerobatic maneuvers that combine rotational motions with high-speed translational motions, in which the translational motion is not left in open-loop. Another direction of the future research is the direct application of the aerobatic maneuvers of aerial robots in the extreme situations, such as falling recovery and obstacle avoidance. Additionally, during the development of the proposed adaptive control scheme, we realize that the aerodynamic eects associated with time-varying ow elds can sig- nicantly in uence the aerobatic ight performances, and, therefore, an appropriate high-order dynamic model is indispensable for control augmentation, stability analysis, and performance improvement. The deep neural networks with multiple layers have the potential to approximate extremely complex dynamics, and, therefore, is a suitable choice for modeling the aerodynamic eects under dierent aerobatic scenarios. Fur- thermore, the control solutions and rigorous stability analysis can be generated using the deep neural network modeling with the guaranteed modeling accuracy. Insect-scale apping-wing robots driven by piezoelectric actuators have been demon- strated to perform regular ights with small attitude angle variations and low transla- tional and rotational speeds. The next step is to enable these insect-scale ying robots to execute high-speed aerobatic maneuvers and mimic the incredible maneuverabilities exhibited by the ying insects. Theoretically speaking, the insect-scale ying robots can achieve much higher angular velocities and accelerations than most of UAVs and MAVs, which suggests a much greater maneuverability to be achieved. 136 Bibliography [1] T. Lee, M. Leoky, and N. H. McClamroch, \Geometric Tracking Control of A Quadrotor UAV on SE (3)," in Proc. 2010 IEEE Conf. Decis. Control, Atlanta, USA, Dec. 2010, pp. 5420{5425. [2] R. Dudley, The Biomechanics of Insect Flight { Form, Function, Evolution. Princeton, NJ: Princeton University Press, 2000. [3] F. T. Muijres, M. J. Elzinga, J. M. Melis, and M. H. Dickinson, \Flies Evade Looming Targets by Executing Rapid Visually Directed Banked Turns," Science, vol. 344, no. 6180, pp. 172 { 177, Apr. 2014. [4] N. F. Steven, R. Sayaman, and M. H. Dickinson, \The Aerodynamics of Free-Flight Maneuvers in Drosophila," Science, vol. 300, no. 5618, pp. 495{498, Apr. 2003. [5] B. Cheng, B. W. Tobalske, D. R. Powers, T. L. Hedrick, S. M. Wethington, G. T. Chiu, and X. Deng, \Flight mechanics and control of escape manoeuvres in hum- mingbirds. i. ight kinematics," Journal of Experimental Biology, vol. 219, no. 22, pp. 3518{3531, 2016. [6] J. Downer, P. Dalton, R. Pilley, and T. Parker (Producers), Earth ight { Africa (Episode 2) [Documentary]. BBC America, 2012. [7] W. Shyy, Y. Lian, J. Tang, D. Viieru, and H. Liu, Aerodynamics of Low Reynolds Number yers. New York, NY: Cambridge University Press, 2008. 137 [8] R. Dudley, \Mechanisms and Implications of Animal Flight Maneuverability," Integ. Comp. Biol., vol. 42, no. 1, pp. 135{140), Feb. 2002. [9] S. Ho, H. Nassef, N. Pornsinsirirak, Y.-C. Tai, and C.-M. Ho, \Unsteady Aerody- namics and Flow Control for Flapping Wing Flyers," Prog. Aerosp. Sci., vol. 39, no. 8, pp. 635{681), Nov. 2003. [10] U. Pesavento and Z. J. Wang, \Flapping Wing Flight Can Save Aerodynamic Power Compared to Steady Flight," Phys. Rev. Lett., vol. 103, no. 11, pp. 118 102 (1{4), Sep. 2009. [11] J. W. Kruyt, E. M. Quicaz an-Rubio, GJ. F. van Heijst, D. L. Altshuler, and D. Lentink, \Hummingbird Wing Ecacy Depends on Aspect Ratio and Compares with Helicopter Rotors," J. R. Soc. Interface, vol. 11, no. 99, pp. 20 140 585 (1{12), Oct. 2014. [12] W. Shyy, H. Aono, C.-K. Kang, and H. Liu, An Introduction to Flapping Wing Aerodynamics. New York, NY: Cambridge University Press, 2013. [13] A. Tayebi and S. McGilvray, \Attitude Stabilization of A VTOL Quadrotor Air- craft," IEEE Trans. Control Syst. Technol., 2006. [14] T. Madani and A. Benallegue, \Control of a Quadrotor Mini-Helicopter via Full State Backstepping Technique," in Proc. Conf. Decis. Control. IEEE, 2006. [15] A. Abdessameud and A. Tayebi, \Global trajectory tracking control of VTOL- UAVs without linear velocity measurements," Automatica, vol. 46, no. 6, pp. 1053{ 1059, 2010. [16] M. D. Hua, T. Hamel, P. Morin, and C. Samson, \A control approach for thrust- propelled underactuated vehicles and its application to VTOL drones," IEEE Transactions on Automatic Control, vol. 54, no. 8, pp. 1837{1853, 2009. 138 [17] A. Roberts and A. Tayebi, \Adaptive position tracking of VTOL UAVs," IEEE Trans. Robot., vol. 27, no. 1, pp. 129{142, 2011. [18] D. Cabecinhas, R. Cunha, and C. Silvestre, \A globally stabilizing path following controller for rotorcraft with wind disturbance rejection," IEEE Trans. Control Syst. Technol., vol. 23, no. 2, pp. 708{714, 2015. [19] Z. T. Dydek, A. M. Annaswamy, and E. Lavretsky, \Adaptive Control of Quadro- tor UAVs: A Design Trade Study with Flight Evaluations," IEEE Transactions on Control Systems Technology, vol. 21, no. 4, pp. 1400{1406, 2013. [20] S. Lupashin and R. D'Andrea, \Adaptive Fast Open-loop Maneuvers for Quadro- copters," Autonomous Robots, vol. 33, no. 1-2, pp. 89{102, 2012. [21] D. Mellinger, N. Michael, and V. Kumar, \Trajectory Generation and Control for Precise Aggressive Maneuvers with Quadrotors," Int. J. Robot. Res., vol. 31, no. 5, pp. 664{674, Apr. 2012. [22] P.-J. Bristeau, P. Martin, E. Sala un, and N. Petit, \The role of propeller aerody- namics in the model of a quadrotor uav," in 2009 European Control Conference. IEEE, 2009, pp. 683{688. [23] M. Bangura and R. Mahony, \Nonlinear Dynamic Modeling for High Performance Control of a Quadrotor," in Proc. 2012 Australasian Conf. Robot. Autom., 2012, pp. 1{10. [24] R. Mahony, V. Kumar, and P. Corke, \Multirotor aerial vehicles: Modeling, esti- mation, and control of quadrotor," IEEE robotics & automation magazine, vol. 19, no. 3, pp. 20{32, 2012. 139 [25] K. Y. Ma, S. M. Felton, and R. J. Wood, \Design, fabrication, and modeling of the split actuator microrobotic bee," in Proc. 2012 IEEE Int. Conf. Intell. Robot. Syst., 2012, pp. 1133{1140. [26] K. Y. Ma, P. Chirarattananon, S. B. Fuller, and R. J. Wood, \Controlled ight of a biologically inspired, insect-scale robot," Science, vol. 340, no. 6132, pp. 603{607, 2013. [27] N. O. P erez-Arancibia, P. J. Duhamel, K. J. Ma, and R. J. Wood, \Model-free control of a hovering apping-wing microrobot," Journal of Intelligent & Robotic Systems, vol. 77, no. 1, pp. 95{111, 2015. [28] P. Chirarattananon, K. Y. Ma, and R. J. Wood, \Adaptive control of a millimeter- scale apping-wing robot," Bioinspiration & biomimetics, vol. 9, no. 2, 2014. [29] S. B. Fuller, \Four wings: An insect-sized aerial robot with steering ability and payload capacity for autonomy," IEEE Robotics and Automation Letters, vol. 4, no. 2, pp. 570{577, 2019. [30] J. B. Kuipers, Quaternions and Rotation Sequences. Princeton, NJ: Princeton University Press, 1999, vol. 66. [31] S. Boyd and L. Vandenberghe, Convex optimization. Cambridge University Press, 2004. [32] D. Lee, \Distributed backstepping control of multiple thrust-propelled vehicles on a balanced graph," Automatica, vol. 48, no. 11, pp. 2971{2977, 2012. [33] Y. Chen and N. O. P erez-Arancibia, \Generation and Real-Time Implementation of High-Speed Controlled Maneuvers Using an Autonomous 19-Gram Quadrotor," in Proc. Int. Conf. Robot. Autom. Stockholm, Sweden: IEEE, May 2016, pp. 3204{3211. 140 [34] J. Franklin, Matrix theory. Englewood Clis, N.J.: Prentice-Hall, 1968. [35] J-J. E. Slotine and W. Li, Applied Nonlinear Control. Upper Saddle River, NJ: Prentice-Hall, Inc., 1991. [36] H. K. Khalil, Nonlinear Systems { Third Edition. Upper Saddle River, NJ: Prentice-Hall, Inc., 2002. [37] Y. Chen and N. O. P erez-Arancibia, \Lyapunov-based Controller Synthesis and Stability Analysis for the Execution of High-speed Multi- ip Quadrotor Maneu- vers," in Proc. Amer. Control Conf, Seattle, USA, May 2017, pp. 3599{3606. [38] W. M. Haddad and V. Chellaboina, Nonlinear Dynamical Systems and Control { A Lyapunov-Based Approach. Princeton, NJ: Princeton University Press, 2008. [39] R. Mahony, T. Hamel, and J.-M. P imlin, \Nonlinear Complementary Filters on the Special Orthogonal Group," IEEE Trans. Autom. Control, vol. 53, no. 5, pp. 1203{1218, Jun. 2008. [40] J. Hunter, An Introduction to Real Analysis. University of California, Davis, 2014. [41] R. Mahony, V. Kumar, and P. Corke, \Multirotor Aerial Vehicles: Modeling, Estimation, and Control of Quadrotor," IEEE Robot. Autom. Mag., vol. 19, no. 3, pp. 20{32, 2012. [42] S. Lupashin, A. Sch ollig, M. Sherback, and R. D'Andrea, \A Simple Learning Strategy for High-Speed Quadrocopter Multi-Flips," in Proc. 2010 IEEE Int. Conf. Robot. Autom. (ICRA 2010), Anchorage, USA, May. 2010, pp. 1642{1648. [43] A. Kushleyev, D. Mellinger, C. Powers, and V. Kumar, \Towards a Swarm of Agile Micro Quadrotors," Autonomous Robots, vol. 35, no. 4, pp. 287{300, 2013. 141 [44] K. Toh, M. Todd, and R. T ut unc u, \Sdpt3: A matlab software package for semidef- inite programming, version 1.3," Optimization methods and software, vol. 11, no. 1-4, pp. 545{581, 1999. [45] M. Grant and S. Boyd, \CVX: Matlab software for disciplined convex program- ming, version 2.1," http://cvxr.com/cvx, Mar. 2014. [46] ||, \Graph implementations for nonsmooth convex programs," in Recent Ad- vances in Learning and Control, V. Blondel, S. Boyd, and H. Kimura, Eds. Springer-Verlag Limited, 2008, pp. 95{110, https://web.stanford.edu/ boyd/ papers/pdf/graph dcp.pdf. [47] G. J. Leishman, Principles of Helicopter Aerodynamics, 2nd ed. Cambridge uni- versity press, 2006. [48] S. Bouabdallah, P. Murrieri, and R. Siegwart, \Design and Control of An Indoor Micro Quadrotor," in Proc. Int. Conf. Robot. Autom., 2004. [49] P. J. Bristeau, P. Martin, E. Sala un, and N. Petit, \The Role of Propeller Aero- dynamics in the Model of a Quadrotor UAV," in Proc. 2009 Eur. Control Conf. IEEE, 2009, pp. 683{688. [50] N. Guenard, T. Hamel, and V. Moreau, \Dynamic Modeling and Intuitive Control Strategy for an "X4-Flyer"," in Proc. 2005 Int. Conf. Control Autom., vol. 1. IEEE, 2005, pp. 141{146. [51] P. Krause and O. Wasynczuk and S. D. Sudho and S. Pekarek, Analysis of electric machinery and drive systems. John Wiley & Sons, 2013, vol. 75. [52] Y. Chen and N. O. P erez-Arancibia, \Nonlinear Adaptive Control of Quadrotor Multi-Flipping Maneuvers in the Presence of Time-Varying Torque Latency," in Proc. Int. Conf. Intell. Robot. Syst. IEEE, 2018, pp. 1{9. 142 [53] ||, \Adaptive control of aerobatic quadrotor maneuvers in the presence of aero- dynamic coecient and torque latency variations," in Proc. Int. Conf. Robot. Au- tom. IEEE, 2019. [54] N. O. P erez Arancibia, S. Gibson, and T.-C. Tsao, \Adaptive Control of MEMS Mirrors for Beam Steering," in Proc. ASME Int. Mech. Eng. Congr. Expo., Ana- heim, CA, Nov. 2004, pp. 1{10 (IMECE2004{60 256). [55] Y. Chen and N. O. P erez-Arancibia, \Controller Synthesis and Performance Opti- mization for Aerobatic Quadrotor Flight," IEEE Transactions on Control Systems Technology, 2019. [56] M. Krsti c, I. Kanellakopoulos, and P. V. Kokotovi c, Nonlinear and Adaptive Con- trol Design. New York, NY: John Wiley & Sons, 1995. [57] P. A. Ioannou and B. Fidan, Adaptive Control Tutorial. Philadelphia, PA: SIAM, 2006. [58] P. A. Ioannou and J. Sun, Robust Adaptive Control. Upper Saddle River, NJ: Prentice-Hall, 1996. [59] K. S. Narendra and A. M. Annaswamy, \Persistent Excitation in Adaptive Sys- tems," Int. J. Control, vol. 45, no. 1, pp. 127{160, 1987. [60] E. D. Sontag, \Smooth stabilization implies coprime factorization," IEEE trans- actions on automatic control, vol. 34, no. 4, pp. 435{443, 1989. [61] A. Levant, \Higher-order sliding modes, dierentiation and output-feedback con- trol," International journal of Control, vol. 76, no. 9-10, pp. 924{941, 2003. [62] S. B. Fuller, J. P. Whitney, and R. J. Wood, \Rotating the heading angle of underactuated apping-wing yers by wriggle-steering," in Proc. 2015 IEEE Int. Conf. Intell. Robot. Syst., 2015, pp. 1292{1299. 143 [63] N. Gravish and R. J. Wood, \Anomalous yaw torque generation from passively pitching wings," in Proc. 2016 IEEE Int. Conf. Robot. Autom., 2016, pp. 3282{ 3287. [64] X. Yang, Y. Chen, L. Chang, A. A. Calder on, and N. O. P erez-Arancibia, \Bee + : A 95-mg four-winged insect-scale ying robot driven by twinned unimorph actua- tors," arXiv preprint arXiv:1905.02253, 2019. [65] N. T. Jaeris, M. A. Graule, and R. J. Wood, \Non-linear resonance modeling and system design improvements for underactuated apping-wing vehicles," in Proc. 2016 IEEE Int. Conf. Robot. Autom., 2016, pp. 3234{3241. 144

Abstract (if available)

Abstract The capability to perform aerobatic flights is essential for aerial animals to survive in wild environments and accomplish many different types of tasks, such as capturing prey, escaping predators, and avoiding obstacles. Consequently, the development of fully autonomous flying robots will not become feasible until the set of basic aerobatic flying capabilities observed in nature can be automatically replicated by artificial flying robots, just because that aerial robots need to survive in the wild or in unstructured environments that are full of unanticipated challenges. In the work presented in this dissertation, we aim to develop high-performance controllers and methods of analysis to enable aerial robots to achieve and even exceed the maneuverability exhibited by natural flying animals. ❧ The study of the flight control of unmanned aerial vehicles (UAVs) and micro aerial vehicles (MAVs) under different flight scenarios has attracted a great deal of attention in past decades. A significant number of analytical studies have proposed numerous different controller designs with rigorous stability analyses of the system

Linked assets

University of Southern California Dissertations and Theses

Conceptually similar

PDF

Development of biologically-inspired sub-gram insect-scale autonomous robots

PDF

Novel soft and micro transducers for biologically-inspired robots

PDF

Dynamic modeling and simulation of flapping-wing micro air vehicles

PDF

Wing selection and design for flapping flight at small scale

PDF

Development and control of biologically-inspired robots driven by artificial muscles

PDF

Analysis, design, and optimization of large-scale networks of dynamical systems

PDF

Identification, control and visually-guided behavior for a model helicopter

PDF

Discrete geometric motion control of autonomous vehicles

PDF

Speeding up trajectory planning for autonomous robots operating in complex environments

PDF

BeeCAS: a collision avoidance system for micro-aerial vehicles

PDF

Optimization-based whole-body control and reactive planning for a torque controlled humanoid robot

PDF

Mathematical characterizations of microbial communities: analysis and implications

PDF

Transient modeling, dynamic analysis, and feedback control of the Inductrack Maglev system

PDF

Trajectory planning for manipulators performing complex tasks

PDF

Model based design of porous and patterned surfaces for passive turbulence control

PDF

Distributed adaptive control with application to heating, ventilation and air-conditioning systems

PDF

Optimal design, nonlinear analysis and shape control of deployable mesh reflectors

PDF

High-accuracy adaptive vibrational control for uncertain systems

PDF

Modeling and analysis of propulsion systems and components for electrified commercial aircraft

PDF

Provable reinforcement learning for constrained and multi-agent control systems

Asset Metadata

Creator Chen, Ying (author)

Core Title Synthesis and analysis of high-performance controllers for high-speed autonomous aerobatic flight

School Viterbi School of Engineering

Degree Doctor of Philosophy

Degree Program Mechanical Engineering

Publication Date 11/18/2019

Defense Date 07/03/2019

Publisher University of Southern California (original), University of Southern California. Libraries (digital)

Tag adaptive control,aerobatic maneuvers,flapping-wing MAV,Lyapunov theory,OAI-PMH Harvest,optimization,quadrotor,VTOL UAV

Language English

Contributor Electronically uploaded by the author (provenance)

Advisor Perez-Arancibia, Nestor (committee chair), Flashner, Henryk (committee member), Jovanovic, Mihailo (committee member)

Creator Email chen061@usc.edu,legencychen@gmail.com

Permanent Link (DOI) https://doi.org/10.25549/usctheses-c89-237971

Unique identifier UC11673403

Identifier etd-ChenYing-7935.pdf (filename),usctheses-c89-237971 (legacy record id)

Legacy Identifier etd-ChenYing-7935.pdf

Dmrecord 237971

Document Type Dissertation

Rights Chen, Ying

Type texts

Source University of Southern California (contributing entity), University of Southern California Dissertations and Theses (collection)

Access Conditions The author retains rights to his/her dissertation, thesis or other graduate work according to U.S. copyright law. Electronic access is being provided by the USC Libraries in agreement with the a...

Repository Name University of Southern California Digital Library

Repository Location USC Digital Library, University of Southern California, University Park Campus MC 2810, 3434 South Grand Avenue, 2nd Floor, Los Angeles, California 90089-2810, USA