UFDC Home  myUFDC Home  Help 



Full Text  
PAGE 1 1 OPTIMAL PATH PLANNING AND TRAJECTORY OPTIMIZATION FOR MULTIPLE AIRCRAFT LANDING USING RRT ALGORITHM AND PSEUDOSPECTRAL METHODS By KRITHIKA MOHAN A THESIS PRESENTED TO THE GRADUATE SCHOOL OF THE UNIVERSITY OF FLORIDA IN PARTIAL FULFILLMENT OF THE REQUI REMENTS FOR THE DEGREE OF MASTER OF SCIENCE UNIVERSITY OF FLORIDA 2011 PAGE 2 2 2011 Krithika Mohan PAGE 3 3 To my parents, Mr. D. Mohan and Mrs. Kusum Mohan; my sister Jyothika; all my friends and family members for motivati ng me and believing in me PAGE 4 4 ACKNOWLEDGMENTS I express my sincere appreciation and gratitude to my supervisory committee chair and mentor, Dr. Anil V. Rao. I thank him for the con stant support, encouragement and technical guidance that he provided me throughout my research work. I especially thank him for providing me the license to his software GPOPS which helped me accomplish my research work I also thank Dr. Warren E. Dixon for lending his knowledge and support during the course of my study at the University of Florida. All that I have learned and accomplished would not have been possible without their support and belief in me. I thank all of my colleagues for helping me with my research work motivating me towards my goal, and creating a friendly work atmosphere. I also extend my appreciation to Michael Patterson and Christopher Darby for sharing their knowledge and providing me with a better understanding of my research work. M ost importantly, I would like to express my deepest appreciation to my parents, D Mohan and Kusum Mohan and my sister, Jyothika. Their love, understanding, motivation and belief in me made this thesis possible. PAGE 5 5 TABLE OF CONTENTS page ACKNOWLEDGMENTS ................................ ................................ ................................ .. 4 LIST OF TABLES ................................ ................................ ................................ ............ 8 LIST OF FIGURES ................................ ................................ ................................ .......... 9 ABSTRACT ................................ ................................ ................................ ................... 11 CHAPTER 1 INTRODUCTION ................................ ................................ ................................ .... 12 Motivation ................................ ................................ ................................ ............... 12 Literature Review ................................ ................................ ................................ .... 12 Conflict Detection and Resolution ................................ ................................ .... 12 Trajectory Path Planning Algorithms ................................ ................................ 16 Dijkstra's algorithm ................................ ................................ ..................... 18 Bellman Ford algorithm ................................ ................................ .............. 18 A* search algorithm ................................ ................................ .................... 18 RRT algorithm ................................ ................................ ............................ 19 Aircraft Trajectory Optimization ................................ ................................ ........ 20 About the Thesis ................................ ................................ ................................ ..... 23 Thesis Outline ................................ ................................ ................................ ......... 24 2 MULTIPLE AIRCRAFT LANDING CONTROL PROBLEM ................................ ...... 26 Structure of an Optimal Control Problem ................................ ................................ 26 First Order Optimality Conditions ................................ ................................ ..... 28 Transversality Conditions ................................ ................................ ................. 28 Me thods of Solving Optimal Control Problem ................................ ......................... 29 Indirect Method ................................ ................................ ................................ 29 Direct Method ................................ ................................ ................................ ... 30 Problem Formulation ................................ ................................ ............................... 30 Aircraft Dynamics ................................ ................................ ............................. 32 States ................................ ................................ ................................ ......... 32 Control ................................ ................................ ................................ ....... 32 State space dynamic equations ................................ ................................ 33 Cost Function ................................ ................................ ................................ ... 33 Obstacle Formulation ................................ ................................ .............................. 34 Two Dimensional Obstacle Modeling ................................ .............................. 34 Three Dimensional Obstacle Modeling ................................ ............................ 36 PAGE 6 6 3 RRT AND MULTIPLE RRT ALGORITHMS ................................ ............................ 39 Original RRT Algorithm ................................ ................................ ........................... 39 Bi Dir ectional RRT Algorithm ................................ ................................ .................. 42 Multiple RRT Algorithm ................................ ................................ ........................... 45 Application of RRT in Multiple Aircraft Landing Problem ................................ ......... 49 4 INITIAL TRAJECTORY GENERATION USING RRT ................................ .............. 50 Parameters for the Problem ................................ ................................ .................... 50 Algorithm for Initial Trajectory Generation ................................ ............................... 51 Initialization ................................ ................................ ................................ ....... 51 Identification of the Nearest Aircraft ................................ ................................ 51 Specification of the Constraints ................................ ................................ ........ 51 Selection of the Way Point ................................ ................................ ............... 52 Aircraft Trajectory and Time ................................ ................................ ............. 52 Advantages and Disadvantages of the Initial Trajectory Generated ....................... 53 Trajectory Plots ................................ ................................ ................................ ....... 53 Three Dimensional Trajectory ................................ ................................ .......... 53 Check for Intersecting Trajectories ................................ ................................ ... 54 5 PSEUDOSPECTRAL METHODS AND GPOPS ................................ ..................... 58 LG, LGR, and LGL Collocation Points ................................ ................................ .... 58 Formulation of Pseudospectral Method Using LGR Points ................................ ..... 59 hp Adaptive Pseudospectral Method ................................ ................................ ...... 61 Multiple Aircraft Landing Multiple Phase Problem Formulation ............................... 62 Objective and Input Parameters ................................ ................................ ....... 62 Methodology ................................ ................................ ................................ ..... 62 Cost Function Formulation ................................ ................................ ............... 64 Constraints Formulation ................................ ................................ ................... 65 Obstacle path constraints formulation ................................ ........................ 65 Intersection constraints formulation ................................ ........................... 66 Event constraints formulation ................................ ................................ ..... 67 Bank angle constraints formulation ................................ ............................ 67 6 RESULTS AND DISCUSSION ................................ ................................ ............... 69 Anal ysis of Results ................................ ................................ ................................ 69 Validation of Results ................................ ................................ ............................... 70 Cost Function Validation ................................ ................................ ................... 70 Constraints Validation ................................ ................................ ...................... 70 Obstac le path constraints validation ................................ .......................... 70 Event constraints validation ................................ ................................ ....... 71 Intersection constraints validation ................................ .............................. 71 Bank angle constraints validation ................................ ............................... 71 7 CONCLUSION ................................ ................................ ................................ ........ 80 PAGE 7 7 LIST OF REFERENCES ................................ ................................ ............................... 82 BIOGRAPHICAL SKETCH ................................ ................................ ............................ 85 PAGE 8 8 LIST OF TABLES Table page 2 1 Notation ................................ ................................ ................................ .............. 26 2 3 Three dimensional obstacle parameters ................................ ............................ 37 5 1 User specified initial and final states of the Aircraft 1, 2, and 3 .......................... 62 6 1 Mayer cost validatio n ................................ ................................ .......................... 70 PAGE 9 9 LIST OF FIGURES Figure page 1 1 Conflict detection and resolution process ................................ ........................... 14 2 1 Figurative illustration of the multiple aircraft landing problem ............................. 31 2 2 Shapes generated for different values of p, a, and b ................................ .......... 35 2 3 Obstacle wall representation with center flyable hole. ................................ ........ 36 3 1 A visualization of a RRT graph growing denser with increasing iterations. ......... 40 3 2 Bi directional RRT algorithm solution with four obstacles and two trees ............ 43 3 3 Mu ltiple RRT algorithm solution with three wall obstacle ................................ .... 46 4 1 Initial three dimensional trajectory of the aircraft using the RRT algorithm ......... 55 4 2 X axis vs. Time plot for the three aircraft using the RRT algorithm ..................... 56 4 3 Y axis vs. Time plot for the three aircraft using the RRT algorithm ..................... 56 4 4 Z axis vs. Time plot for the three aircraft using th e RRT algorithm ..................... 57 5 1 Schematic showing LGL, LGR and LG collocation points ................................ .. 59 5 2 Schematic showing the six phases in the problem ................................ ............. 63 6 1 Three dimensional trajectory with the obst acle wall obtained using GPOPS ..... 72 6 2 Y Z axis view of the three dimensional trajectory ................................ ............... 73 6 3 X axis vs. Time plot of the three aircraft in different phases ............................... 74 6 4 Y axis vs. Time plot of the three aircraft in different phases ............................... 74 6 5 Z axis vs. Time pl ot of the three aircraft in different phases ............................... 75 6 6 Flight path angle vs. Time plot of the three aircraft in different phases ............... 75 6 7 Heading angle vs. Time plot of the three aircraft in different phases .................. 76 6 8 Load factor controls, of Aircraft 1 ................................ ................... 76 6 9 Load factor controls, of Aircraft 2 ................................ ................... 77 6 10 Load factor controls, of Aircraft 3 ................................ ................... 77 PAGE 10 10 6 11 Bank angle of Aircraft 1 in different phases ................................ ........................ 78 6 12 Bank angle of Aircraft 2 in different phases ................................ ........................ 78 6 13 Bank angle of Aircraft 3 in different phases ................................ ........................ 79 6 14 Intersection constraints validation ................................ ................................ ....... 79 PAGE 11 11 Abstract of Thesis Presented to the Graduate School of the University of Florida in Partial Fulfillment of the Requirements for the Degree of Master of Science OPTIMAL PATH PLANNING AND TRAJECTORY OPTIMIZATION FOR MULTIPLE AIRCRAFT LANDING USING RRT ALGORITHM AND PSEUDOSPECTRAL METHODS By Krithika Mohan May 2011 Chair: Anil V. Rao Major: Mechanical Engineering In this work a method is presented for solving multiple aircraft l anding problem using an RRT a lgorithm a nd p seudospectral method The multiple aircraft l anding problem involves generating a conflict free traject ory for all the aircraft waiting to land o n a single runway whil e avoiding obstacles The traject ory generated for each aircraft should satisfy time distance separation constraint s between two successive aircraft landing s and also the trajectories should not intersect The multiple aircraft landing prob lem is solved as a multi ple phase optimal control problem using the open source optimal control software GPOPS, and a nonlinear programming problem solver SNOPT In order to improve computational eff iciency a rapidly exploring random tree (RRT) algorithm is used to generate initial guesses for the optimal control software. While the RRT algorithm determines the shortest pat h between the initial position and the target position a novel algorithm has been presented which also combines other constraints like maintaining a minimum safe time distance difference and avoiding int ersecting flight paths ; thus making it suitable for initializing t he multi ple aircraft landing problem PAGE 12 12 CHAPTER 1 INTRODUCTION Motivation Currently approximately 5000 aircraft land every hour at an airport in the U nited S tates. Due to this large increase in density of aircraft arrivals, it has become important to optimize the scheduling of aircraft landings in order to reduce wait time impro ve airport efficiency, and minimize fuel consumption while maintain ing safety. As the density of aircraft arrivals increases, so does the complexity of trajectory planning and generation. Presently, t he air traffic controller (ATC) uses its radar system to maintain the p osition of each of the aircraft arriving and waiting to land at the a irport The ATC provide s inform ation and support to the pilots to prevent c ollisions with other aircraft The ATC makes use of the trajectory provided by the Instrument Landing System (ILS) which uses a combination of radio signals and high intensity lig ht arrays to guide the aircraft for safe and precision landing. All these systems generate the trajectory using the present and predicted future position of the aircraft. Alt hough these modern equipments are being used by the ATC for assistance, h uman error is always a possibility In this work aircraft landing trajectories are generated by taking into consid erati on va rious constraints such as the time difference between each aircraft la nding, obstacle collision avoidance, dynamics and minimum time for approach. Literature Review Conflic t Detection and Resolution A conflict is defined as the situation of loss of minim um safe separation between two aircraft [1] The ATC is responsible for conflict detection and its resolution. It also PAGE 13 13 tries to fulfill the desired arrival and departure time of the aircraft Conflict detection is the process of evaluating the pr obability of a conflict that will occur in the near future based on the current position and flight path of an aircraft Conflict resolution is the process of find ing the maneuver required by one or multiple aircraft to avoid the predicted c onflict. A great deal of research has been done to improve the performan ce of the air traffic c ontroller. A review of conflict detection and re solution modeling techniques has been discussed by Kuchar, J. K. and Yang, L.C. [2] The c onflict detection and resolution process consists of predicting the conflict, communicating it to the pilot an d helping resolve the conflict This process is shown in Fig ure 1 1. The surveillance is done by radars GPS and other sensors, b ased on which the current state of each aircraft is estimated. The dynamic trajectory model and the flight plan are used to determine the future stat e. The current and future state is compared with the metrics for conflict detecti on and resolution. This information is given to the air traffic c ontroller (human o perator), who communicates with the p ilot and helps resolve the conflict A number of conflict resolution algorithms that use s genetic algorithm semi definite programming, mixed integer programming (MIP) and s equential quadratic programming (SQP) method s have been proposed in literature [3 6] A real time aircraft conflict reso lution approach that uses genetic algorithm is proposed by Durand, N. and others [6] The authors solve a real variable combinatorial problem which combines both local and global optimization meth ods. The genetic algor ithm is ru n initially and is followed by a local optimization hill climbing algorithm to improve the best solution obtained for each chromosome. PAGE 14 14 Figure 1 1. Conflict detection and resolution p rocess The conflict resolution problem is formulated as a non convex quadratical l y constrained quadratic problem by Frazzoli E. and others [3] which is approximated by a convex semi definite program and solved The optimal solution to this program is used to generate feasible and locally optimal maneuvers. The authors follow the resolution methodology where they try to minimize the deviation between the desired heading proposed by th e aircraft and a conflict free heading determined by the air traffic controller. A m ulti ple airc ra ft conflict avoidance problem i s also discussed by Pallottino, L. and others [4] where a mixed integer linear program ming is used to formulate the problem with only velocit y or heading changes as admissible m aneuvers The linear PAGE 15 15 formulation of the two problems i s solved using the optimization software CPLEX to find the optimal solution to the mixed integer p rogram (MIP). A cost function which fav ors the lower priority aircraft for resolving the predicted conf licts by involving changes in heading, speed and altitude was introduced b y Hu, J. and others [5] T he vertical maneuvers are penal ized for passenger comfort. The authors developed a constrained optimization problem for a two aircraft scenario, which is extended to multiple aircraft by approximating the optimization problem with a finite dimensional convex optimization problem with linear waypoint constraints. Additional constraints on waypoint position were added to avoid sharp turns near the waypoints A t hree dimensional tr ajectory optimizati on algorithm i s developed by Menon, P.K. and others [7] to obtain a conflict free flig ht path. The algorithm use s non linear point mass model and include s realistic operational con straints on individual aircraft; making the resulting problem non convex. Also the problem of optimal cooperative three dimensional conflict resolution involvi ng multiple aircraft is addressed by Ra ghunathan, A. U. and others [8] Here the author s specify the initial and final locations along with the detailed point mass aircraft dynamics. The infinite dimensional optimal control problem is converted into a finite dimensional nonlinear program (NLP) by use of collocation on finite elements. This NLP is solved in IPOPT using an integer point algorithm which employs a line search method. The key concept in this paper is that it used the dynamic optimization by using the detailed dynamic model of the aircraft instead of using t he kinematic model; thus generating feasible trajectories that are flyable. T he Hu et al procedure is used to solve the optimal control problem and obtain the trajectory of in dividual aircraft This trajecto ry is used for initialization (initial guess) of the NLP. PAGE 16 16 In t his work a detailed point mas s dynamic model of the aircraft is used for generating flyable trajectories [9] The infinite dimensional optimal control proble m here is converted into a f inite dimen sional nonlinear program (NLP) by using the open source optimal control software, GPOPS and solved using the NLP solver, SNOPT. An additional capability of the proposed method is that every aircraft ca n specify it s in dividual dynamic model. This means that the trajectory can be genera ted simultaneously for different types of aircraft ( incl uding helicopters ) that may occupy the airspace Also fo r initi alization of the NLP, a rapidly exploring random t ree (RRT) algorithm is employed which finds the shortest path between two points avoiding any obstacles in the path. The shortest possible time required by each of the aircraft to reach the goal point maintaining the time difference constraint between successive air craft landings (without considering the aircraft dynamics) is used to sp ecify the final times in the optimal control software. The RRT algorithm is explained in detail in Chapter 3. As the speed of the aircraft is maintained constant in this problem, the o nly allowable change is the bank angle, heading angle and flight path angle of the aircraft by using the vertical and horizontal load factors as the control. Trajectory Path Planning Algorithm s The main objective of a path planning algorithm is to find a path whic h satisfies certain conditions while avoiding obstacles in the path and preventing collisions with other moving objects. Such t rajectory or motion plann ing algorithms have been primarily used in robotics and dynamics and control. In dynamics and c ontrol the guidance of vehicles depend s on the dynamic behavior of the vehicle. Both applications need motion planning algorithms to find path to reach a goal in a partially known environment with obstacle avoidance and also collision avoidanc e with othe r vehicles. However, f or PAGE 17 17 aerial vehicles additional parameters need s to be considered such as the three dimensional environment, non trivial dynamics, safety due to disturbed operating conditions, and high level of uncertainty in the state knowledge due t o sensor and radar measurements An overview of all the existing moti on planning algorithms used in unmanned aerial v ehicle, (UAV) guidance is presented by Goerzen C. and others [10] A comprehensive study and comparison of various path planning methods currently in use wa s provided by Hurni, M.A. [11] The author compa red and graded algorithms like b ug algorithm potential fields, r oadmaps, cell decomposition and optim al c ontrol based on its ability for path planning and control, o ptimality obstacle avoidance, handling v ehicle constraints, global vs. local i nformat ion, computational complexity, feasibility, c ompleteness multiple vehicles, error and u ncertainties. B ase d on his inference, the optimal c ontrol method was found to be the best in al l the criteria except computational c omplexity. The research [11 13] also discussed the obstacle modeling and avoidance techniques for use in ground v ehicle maneuver Some of the basic shortest path finding algorithms are given below : 1. Dijkstra's algorithm used to solve the single source and single destination shortest path problems with non negative edge weights 2. Bellman Ford algorithm used to solve the single source problem with negative or non negative edge weights 3. A* search algorithm used to find the least cost path using heuristics to increase the search speed. 4. RRT algorithm used to find shortest path with obstacle avoidance. It biases the search to the least explored regions, for effective and faster result. 5. Perturbation theory used to find the locally shortest path in worst case PAGE 18 18 D ijkstra's a lgorithm Dijkstra's algorithm is a single source shortest path planning algorithm which has non negative weights in its edges. It is also known as a greedy algorithm. The start point is at a source vertex, S from which a tree, T grows such tha t it ultimately spans all vertices reachable from S. The vertices that are closest to the tree are added first, i.e., first the source S is added, then the vertex closest to S, then the next closest, and so on. The worst case running time for the Dijkstra algorithm on a graph with n vertices and m edges is because it allows for directed cycles. It even finds the shortest path from a source node S to all other nodes in the graph. Bellman Ford a lgorithm The Bellman Ford algorithm is also a single so urce shortest path planning algorithm which can have negative weights in its edges. It does not follow the greedy approach used in Dij source vertex to 0 and all other ver V i s the number of vertices, then it does  V  1 passes over all edges relaxing, or updating, the distance to the destination of each edge. Thus it checks all the edges again for any negative weigh cycles in which case it returns false. No shortest path ex ists during the negative cycle. The time complexity is O ( V  E ) time, where  V  and  E  are the number of vertices and edges respectively. A* s earch a lgorithm A* use s a best first search algorithm to find the least cost path from a given initial node t o a goal node. It uses a heuristic function, to determine the order of the PAGE 19 19 search. This heuristic function is a sum of two functions, one for minimizing the distance and other for minimizing the cost. It can be written as, (1 1) w h ere is the path cost function which denotes the cost to be minimized from the starting to the goal node ; and is estimate of the shortest distance from the starting to the goal node, which is usually denoted by a straight line distance to the go al. A* can also be used as a bidirectional heuristic search algorithm like the RRT algorithm, which is explained in Chapter 3. The time complexity of A* will depend on the heuristic. In the worst case, the number of nodes added to the tree is exponential in the length of th e shortest path solution. It becomes a polynomial if the search space is a tree with a single goal state, and if the h euristic function h meets the following c ondition, where h is the optimal distance cost to get from node x to the goal node That is, the error of h will not grow faster than the logarithm h that returns the true optimal distance from node x to the goal node. RRT a lgorithm A new technique of path finding called rapidly r andom exploring t rees (RRT) wa s proposed by Lavalle, S.M. [14] The main adv antage of this algorithm over other algorithm is that it can be used for solving problems which involve s obstacles as path constraints and also kinodynamic or nonholo nomic differential constraints of the vehicle. It can be used t o generate trajectories for non linear systems with state and path constraints. It also biases the search to the least explored regions, for effective and faster result. The different types of R RT algorithm is explained in detail in Chapter 3. PAGE 20 20 The RRT algorithm is being used for various applications in recent years. It is used as a path planner for real time obstacle detection and avoidance in auton omous small scale Helicopters [15] Also Caves, A.D.J [16] use d the RRT algorithm for human supervisory control of UAV missions in environments with dynamic obstacles fields of different densities. The algorithm is u sed to assist the human operator in path planning and re planning for a group of UAVs. The RRT algorithm is also used by Aoud, G.S. [17] for designing optimal reconfiguration maneuvers for multiple spacecrafts. The author proposed a two stage approach where he used RRT algorithm in the first stage of the initializing the optimal control problem without the differential constraints and get a feasible initial guess. The RRT algorithm is used in this work for initializing the multiple aircraft landing optimal control problem, due to its ability to find faster solutions especially in an environment with obstacles. Aircraft Trajectory Optimization T he optimal flight trajectory generation in the presen ce of windshear is discussed by Miele, A. and others [18, 19] The numerical solution of this optimal control problem of the Bolza type is solved using the dual sequential gradient restoration alg orithm (DSGRA). The author optimizes the take off trajectories by minimizing the peak deviation of absolute path inclination from a reference value. Also the a bort trajectories are optimized by minimizing the peak drop of altitude from a reference value. The sur vival capability of the optimal abort landing trajectory in a severe windshear was found superior to that of the other compa rison trajectories like the constant pitch guidance trajectory and the maximum angle of attack guidance trajectory. PAGE 21 21 An a utomatic air craft landing control algorithm is also developed using both fuzzy logic and neural networks models [20 22] While the authors, Nho, K. and Agarwal, R.K. [20] use d the linear and non linear aircraft model, the authors, Juang, J. and Chio, J. [21] use d a lineariz ed aircraft model for implementing fuzzy logic technique An artificial neural network is used by Chaturvedi, D.K. and others [22] to deal with the non linearities in the aircraft model. Also Ruffier, F., and Franceschini, N. [23] d eveloped a visual gu idance based autopilot which enable s the micro air vehicle to take off, cruise and land, while reacting adequately to wind disturbances. An autonomous counter hijack control of aircraft is proposed by Patel, R.B., and Goulart, P.J. [24] where the critical building s are modeled as obstacles to be avoided in the aircraft path. These building a re modeled by employing dualization of state exclusion regions to maintain continuity. The warm start method which uses prior solutions as initial guess is used for initialization of a direct multiple shooting method to solve the optimal control problem. In recent years, p seudospectral met hod s, a type of direct trajectory optimization has been used extensively in aerospace application s as an effective and faster w ay of solving optimal control p roblems. A p seudospectral method for real time motion planning and ob stacle avoidance is discusse d by Lewis, L R. and others [12, 13] for generating minimum time trajectories for unmanned ground vehicles. The gauss p seudospectral method (GPM) is a lso used by Aoud G.S. [17] to solve the optimal control problem of reconfiguration maneuvers for multiple spacecrafts. A Gauss p seudospectral metho d using orthogonal collocation which uses the LG points for collocation i s develope d by Benson, D.A. and others [25] The orthogonal PAGE 22 22 collocation of the dynam ics is performed only at these collocation p oints, and not the boundary points. It also properly approximates the costates which leads to a set of KKT conditions that is identical to the discretized form of first order optimality conditions at the collocation points. This helps in accurate costate estimation using the KKT multipliers of the NLP. A method for direct trajecto ry optimization and costate estimation of finite horizon and infinite horizon optimal control problems using global collocation at Legendre Gauss Ra dau (LGR) points is presented by Garg, D. and others [26] It is shown that the dual multipliers for the discrete scheme correspond to a pseudospectral approximation of the adjoint equation using polynomials one degree smaller than that used for the state equation. Also it is shown that the inverse of the pseudospectr al LGR differentiation matrix is precisely the matrix associated with an implicit LGR integration scheme. Hence, the method is a global implicit integration method or a pseudospectral method. Consequently, a unified fra mework was designed for numerically s olving optimal control problems using collocation at Legendre Gauss(LG), Legendre Gauss Radau (LGR), and Legendre Gauss Loba tto (LGL) points [27] It wa s shown that the LG and LGR differentiation matrices are rectangul ar and full rank whereas the LGL differentiation matrix is square and singular. Hence, the d ifferential and integral forms are equivalent in the LG and LGR schemes while they are not equivalent in the LGL scheme. The LG and LGR discrete costate systems wer e found to be full rank while the LGL dis crete costate system i s rank defi cient. Also, t he LGL costate approximation is PAGE 23 23 found to have an error that oscillates about the true solution due to the null space in the LGL discrete costate system An hp adaptive pseudospectral m ethod with collocation at Radau points is chosen in this work [28] This method iteratively and simultaneously determines the number of segment breaks, the width of each segment, and the polyn omial degree required in each segment for approximation, until the user specified accuracy is achieved. It leads to higher accuracy solutions with less computational effort and memory than is required in global pseudospectral methods. As the problem in thi s work has multiple aircraft with an obstacle wall in the path to be avoided, the problem is broken into several overlapping phases, for each aircraft to cross the obstacle at the desired location and then reach the goal. Due to the complexity of the problem, the hp adaptive method is chosen. About the Thesis In this research work an RRT algorithm is used for simultaneous pat h planning of multiple aircraft to avoid any collision with each other, as the aircraft reach es t he goal one by one. While the RRT algorithm determines the shortest pat h between the initial position and the target position a novel algorithm has been presented which also combines various other constraints like maintaining a minimum safe distance and t ime difference between each of the aircra ft landings and also avoiding any intersecting flight path s The number of aircraft and their d ifferent start positions is to be provided initially to the algorithm which then generates a trajectory for these aircr aft using the RRT algorithm to avoid any static obstacles in the path. T hese static obstacles can be buildings or other g round or a ir vehicles or a no fly zone. PAGE 24 24 The fi nal time taken by each aircraft to reach the goal and the path genearated is then used in initializing the m ulti ple phase optimal c ontrol problem that is converted into a nonlinear p rogram (NLP) using GPOPS and solved using SNOPT. An hp adaptive p seudospectral method that uses Radau collocation poi nts is implemented to simulta neously generate the trajectory of three aircraft with overlapping phases. These trajectories satisfy the dynamics of the aircraft The research can be used to not only assist the ATC in planning the airc raft trajectory and avoiding conflict s with other aircraft but also aims at avoiding obstac les (three dimensional or two dimens ional) in its path. These obstacles can be important buildings or a no fly zone Hence it also helps in securing these areas and buildings. With future work on this research rea l time implementation may also be possible. Thesis Outline This thesis is o rganized as follows. Chapter 1 is an i ntroduction on the motivation for the research on m ultiple a ircraft l anding trajectory g eneration and a review of the previous work related to the r esearch. It also describes the c onflict detection and resolution problems, along with variou s path planning algorithms and trajectory o ptimization methods. Also the choice of using the RRT algorithm for initial path planning and p seudospectral method for optimization along with a brief description of the research work has bee n presented. In Chapter 2, the multiple aircraft l anding p roblem is formulated along with a general des cription of optimal c ontrol p r oblems with its optimality and t ransversality conditions Also the method for obstacle formulation in three d imension is described. Chapter 3 discusses in detail the three different types of rapidly random exploring (RRT) a lgorithm. The advantages and differences of one over the other are discussed. The choice of using m ulti ple RRT over the others is also PAGE 25 25 presented here In Chapter 4, a novel algorithm is developed that use s m ulti ple RRT algorithm for generating trajectories for multiple aircraft landing problem. Also examples of t rajectories generate d using three aircraft are shown. Chapter 5 gives a brief description of the Radau psedospectral and hp adaptive methods initially and then describes the formulation of the op timal control problem using overlapping phases. It also give s examples for trajec tory generated simultaneously for three aircraft using six phases. PAGE 26 26 CHAPTER 2 MULTIPLE AIRCRAFT LA NDING CONTROL PROBLE M This chapter gives a brief description of an optimal control problem and methods of solving this problem. The multiple aircraft landing problem is describ ed along with the equations of motion of the aircraft a nd the formulated cost function. The three dime nsional obstacle formulation is also described in this chapter. Structure of an Optimal Control Problem Table 2 1 Notatio n Initial Time Final Time State value at the initial time, State value at the final time, Cost Function Augmented Cost Function Mayer Cost Lagrange Cost Boundary Condition Path Constraint Hamiltonian Lagrange Multiplier for Boundary Condition Costate of the Differential Equation V Fixed Speed of the A ircraft, g Acceleration due to G ravity Time taken by Aircraft to Reach the Obstacle W all, An optimal control p roblem is formulated with the dynamic equations of motion, the objective function boundary conditions and p ath constraints. The dynamic equations should be continuous and differentiable. It is converted into the state space representation to get the individual state equations and its relation PAGE 27 27 with the control input. Hence the number of equations is equal to the number of states in the problem. These dynamic equations are represented as given in E quation ( 2 1 ) (2 1) w here is the state, is the c ontrol and t is the time. The o bj ective function also called as c ost functio n, is the parameter that has to be optimized (minimized or maximized) while solving the problem. The c ost function is a function of state variables, control variables and/ or time. That is it can be the end point states, end point time or the state and/or control over the entire period of time. The end point cost is known as Ma yer c ost and the cost function consisting of state or control variables for the entire period of time is known as Lagrange c ost. (2 2) In the above c ost function, J is the Ma yer c ost and is the Lagrange c ost. Bo undary c onditions are the initial and terminal values of states and t ime known at the beginning of the problem or to be achieved at the end. They become the initial and terminal constraints of the problem. (2 3) The path c on straints are the linear or non linear constraints to be satisfied by th e system in its path They can be equality constraint or an inequality constraint. The non linear inequality constraint is represented as (2 4) The non linear equality constraint is represented as (2 5) PAGE 28 28 First Or der O ptimality C onditions The augmented cost function is formulated as, The Hamiltonian is represented as T he augmented function in terms of Hamiltonian is given by The optimality condition for c os tate dynamics is given as The optimal c ontrol is fou nd using the condition given as Transversality C onditions T he first variation of the augmented c ost function is used to find the transversality conditions. These conditions can be used to solve the differential equations which are found from the first order optimality conditions. There are different transversality equations f or different bounda ry conditions on state and ti me of the system which are given below: For no boundary c ondition s on state and time of the syst em the boundary condition is For f ixed initial s tate, the costate at initial time is given by, PAGE 29 29 For fixed final s tate, the costate at final time is given by, For fixed initial t ime, the Hamiltonian at initial time is given by, For fixed final t ime, the Hamiltonian at the final time is given by, Methods of Solving Optimal Control Problem An optimal c ontrol problem is one where differential equations need to be sol ved subject to constraints while optimizing a performance index. The two m ain classifications of methods to solve such optimal control p robl ems with path constraints and boundary conditions are indirect methods and direct methods. Indirect M ethod This method emanates from variational c alculus. It requires solving of a t wo point boundary value problem. Here the o ptim ality conditions for solving of optimal trajectory are to be derived. The first order optimality conditions are formed using variation of calculus by tak ing the first variat ion of the c ost functional. This gives the first order diff erential equations for the states and the c ostates. The transversality conditions usually give n on linear equations which cannot be solved with Ricatti e q uation. Thus they form the nonl inear constra ints in the problem. Indirect method converts the optimal c ontrol problem into a purely di ff erential algebraic system of states, costates and dynamics. Hence it becomes a root finding problem where a set of diff erential equations has to be solved and roots are found to PAGE 30 30 satisfy the constraints. The diff erent types of indirect m ethod are Indirect shooting and Indirect multiple shooting. In Indirect multiple shooting method, the time interval is divided into several smaller intervals with indirect sho oting implemented in each interval. Direct M ethod In direct method, all the functions ( s tates and c ontrol) in the optimal control p roblem are approximated and transcribe d into a nonl inear opt imization problem. By making this approximation, the continuous time infinite dimensional problem is converted in to a nonlinear p rogra mming problem which is a finite dimensional problem. T he quadrature approximation of the i ntegral is used in the continuous c ost function as given below The path constraints are evaluated only at spec ific points of the quadrature. P seudospectral method is a direct method of solving an optimal c ont rol p roblem, which is explained in detail i n Chapter 5. P roblem F ormulation In this work three aircraft are taken into consideration to generate the landing trajectories Initially the position of all the three aircraft is recorded, assuming this time to be the start time, of the problem The objective is to make all these aircr aft land (reach a common goal p oint) on a single runway such that there is a safe time difference between each aircraft landing s The aircraft t hat is nearest to the runway should land first, then the next nearest airc ra ft, and so on. In the aircraft d ynamics, t he speed of all the aircraft is considered to be constant. Hence, only t he trajectory of these PAGE 31 31 aircraft (heading angle, flight path angle, and bank angle) can be manipulated to achieve the objective. Also the pat h has obstacles, which is formulated like a wall of some thickness with a small open ing, through which the aircraft has to pass to reach the runway This enforces additional constraints, like obstacle avoidance and also avoid ing collision with each other. A figurative il lustration of the problem with three aircraft and the obstacle wall is shown in Figure 2 1 It shows that, the trajectory of the farthest aircraft has to take a circuitous path to maintain the time difference between successive aircraft la ndings, prevent collision with other aircraft and also should avoid the obstacle wall Figure 2 1 Figurative i llustration of the m ultiple a ircraft landing problem The aircraft dynamics in three d imension is considered in this work Though all the three aircraft are considere d to have the same 3 D dynamics in this problem, the way the p seudospectral implementati on is formulated, each aircraft can have its own dynamic equations, which is explained in Chapter 5 PAGE 32 32 Aircraft Dynamics States There are five s tates considered for each aircraft. They are 1. x x position 2. y y posit i on 3. z z position 4. Flight path angle 5. Heading angle Here the first three states are the x y and z position of the aircraft in three d imension. The fourth state is the f light p ath angle, which is defined as the angle between the airc orizontal reference plane vector. The fifth state is the h eading angle, which is the yaw angle (rotation about the z axis, which is the upward/downward facing axis). It decides the left or right direction of the aircraft. Control There are two control v ariables in this problem. They are 1. Horizontal loa d factor 2. Vertical load factor The total load f actor is, This total load factor n determines the ratio of lift to w eight. The control is also related to the bank a ngle with the relation (2 6) The b ank angle is the angle between the horizontal reference plane vector towards the left and the right wing of the aircraft in the lateral plane, which is positive PAGE 33 33 (clockwise) when the right wing is facing down and negative (anti clockwise) when it is facing up, with the horizontal plane vector towards the right. State space dynamic e quations The s tate space representation for e ach of the state is given below The a cceleration due to gravity, g is taken as 0.0052952 and t he fixed speed value, V is taken as 0.04 which is about 1.5 times the stall speed of a generic medium range aircraft. This speed is near landing speed of an aircraft. Cost Function The cost function is formulated and explained in Chapter 5. The Mayer cost function is the time for the aircraft to reach the obstacle wall, which is being maximized in this problem. The Lagrange cost in this problem is minimizing the change in flight path angle and heading angle of the aircraft. The formulated cost is given as, where is the time taken by the aircraft to reach the face of the obstacle, N is the total number of aircraft taken into consideration, n is the number of aircraft present in PAGE 34 34 P hase j, and 2N is the total number of phases formulated in the problem as seen in Cha pter 5. Obstacle Formulation The o bstacles form the path constraints of the proble m. All the functions in an optim al control problem should be continuous and differentiable. Hen ce the obstacles should be represented by smooth function. It can be approxima ted by geometric shapes. In literature most of the obstacles are formulated as a circle in two d imension and a s a sphere in three d imension. But Hurni, M.A. [11] represented pol ygonal obstacles in two d imension by using the p norm method. Two Dimensional O bstacle M odeling The two dimensional obstacles can be modeled using the p nor m method [11 13, 29] Using the p norm method, obstacles can be modeled as circle, ellipse, square, rectangle of different types. These shapes can be used to fit most of the obstacles. It is very complex to model the obstacle to the exact shape and si ze. The general equation used to model these shapes is given below, The variables and are the center point of the obstacle. The value of a, b, c indicate the size of the obstacle, where distance along the y axis from the center of the obstacle and c is the radius of the obstacle. The Figure 2 2 shows the various shapes generated by using different values f or p, a, and b in the general eq uation given above with c=1. For example to represent a rectangle with center at (2, 2) and length and breadth as 8 and 4 respectively, the equation can be written as, PAGE 35 35 Assume c to be 1 here. Figure 2 2 Shapes generated for different values of p, a, and b value of the power. A l s o as the power of 50 can become a large constraint to handle, the natural logarithm is taken on both sides of the example equation PAGE 36 36 Three Dimensional O bstacle M odeling In this work the problem is defined in a three d imensional environment; hence the obsta cle should also be defined in three d imension. Thus the p norm method is used to represent the obstacle in three dimensions with x, y, and z axis The general equation for 3D obstacle repres entation is given below The center point of the obstacle is now represented as The distance along the x, y, and z axis from the center is given as a, b, c. The radius of the obstacle now is given by d. Figure 2 3 Obstacle wall r epresentation with ce nter flyable hole. PAGE 37 37 This method is used to model an obstacle, which is of t he form of a wa ll with a square hole in the center The wall has a small thickness, wh ich can be modeled in three dimension The wall with hol e is constructed by taking four obstacles such that a square flyable hole is formed by the m. This is shown in the Figure 2 3 The four obstacle walls are numbered as 1, 2, 3 and 4 in the Figure 2 3 Each of the walls is taken as an individual obstacle to be modeled. They are modeled as a function using the general equation to make it smooth a nd differentiable. These functions form the n onlinear path constraints in the problem. The limits on x, y, and z axis of the environment considered are taken as, 3 to 3, 0 to 6 and 0 to 0.3 respectively. Table 2 3 Three d imensional obstacle p arameters Obstacle Center p Distance along x axis Distance along y axis Distance along z axis 1 ( 2, 3.0005, 0.15) 50 1 0.0005 0.15 2 (2, 3.0005, 0.15) 50 1 0.0005 0.15 3 (0, 3.0005, 0.25) 50 1 0.0005 0.25 4 (0, 3.0005, 0.05) 50 1 0.0005 0.05 Table 2 2 displays the parameters of the obstacles used in the problem. The center, p value and the distance along each axis (a, b, and c value) of all the obstacles are displayed. All the distance is taken in nautical miles. The user specified parameters of the obstacles are then model ed into functions using the general equation for three dimensional obstacles. The obs tacles of the example problem solved in this work is modeled as below PAGE 38 38 Natural logarithm is taken on both the sides to make the constraint with the power of 50 smaller. The set of obstacles are modeled in the problem as below, The obstacle path constraint is then taken as, T he se path constraints are given a minimum value of 0 to prevent any collis ion of the aircraft with these obstacles. Hence the path generated should satisfy these path constraints ( i.e. its value should always be greater than 0) at all the collocation point s T he multiple aircraft landing problem is formulated along with the obstacles in the path. This problem is to be solved as an optimal control problem, to get the desired traj ectory satisfying all the constraints. Chapter 3 describes in detail about the RRT algorithm which is used to generate the initial trajectory guess that avoids th e obstacl es which acts as the non l inear path constraints of the problem. PAGE 39 39 CHAPTER 3 RRT AND MULTIPLE RRT ALGORITHMS This c hapter describes the rapidly exploring random t ree (RRT) a lgorithm, which is a type of randomized shortest path fin di ng a lgorithm. RRT is especially good for solving problems which involve obstacles as path constraints and kinodynamic or nonholonomic differential constraints of the vehicle. It can be used to generate trajectories for non linear systems with state and path constraints. It biases the search to the least explored regions for effective and faster result. That is the probability that a point is selected and added to the tree is proportional to the area of its Voronoi region. The Voronoi diagram is formed by partitioni ng an area with multiple points into convex polygons such that each polygon has only one generating point within it, and all the other point is closer to this generating point than any other point. There have been many improvements made to t he RRT a lgorithm over the years as seen in literature. Th is c hapter will discuss about these a lgorithm s including the advantage s and difference s of one over the other. There are three main algori thm s described in this c hapter : (1) Original RRT algorithm, ( 2) Bi directional RRT algorithm, and (3) Multiple RRT algorithm. The RRT algorithm that is used in this work is called multiple RRT algorithm which has multiple trees simultaneously generating towards the goal point and trying to connect with each other at every iter ation, to ultimately reach the g oal point in shortest possible time. Original RRT Algorithm The main objective of RRT a lgorithm is to find a path from point A to point B with obstacles defined in the environment. They are widely used in path planning of robot PAGE 40 40 a rm to avoid an obstacle whi le reaching another point. The r obot arm has to change its orientation and position simulta neously to make the movement. The RRT algorithm is also used in path planning of v ehicles while maneuver ing through a set of obstacles. The obstacle can be defined in 2D or 3D space. The Figure 3 1 gives a visualization of th e RRT tree growing denser as the iterations increase Figure 3 1 A visualization of a RRT graph growing denser with increasing iterations. In the RRT algorithm that was developed originally, a configurati on space, C is defined, in which specifies the position and orienta tion of the object in its path. is the set of free configuration space, w here the object can move freely without colliding with the obstacle. The start point A is taken as and final goal point is taken as The obstacle is defined by specifying its vertex points. The algorithm generates intermediate po ints randomly, biasing those points towards the Voronoi regions with larger area. It then tries to connect these poi nts with the nearest neighbor vertex/node of the tree, such that no collision occurs while making the connection. The collision is checked b y a Collision_Check function. If no collision occurs, the intermediate point becomes a part of PAGE 41 41 the main tree by adding it as its vertex/node, and creating an edge connecting it to the tree. Ho wever if a collision occurs, the algorithm rejects that intermed iate random point and selects ano ther point to continue the above mentioned steps until it reaches the final goal point. The tree is represented by T. Also n is the number of iterations for which the algorithm repeats itself. It represents the maximum number of times the algorithm should find a new intermediate point, to reach the goal point. Algorithm 3 1 RRT Original ( ) CONSTRUCT _RRT( ) 1 T.start ( ) ; 2 for j =1 to n do 3 ; 4 ADD ( T ); 5 end for 6 Return T; EXPAN D ( T ) 1 ; 2 if NEW_CONFIG ( then 3 T.add_vertex ( ; 4 T.add_edge ( ; 5 if then 6 Return Reached; 7 else 8 Return Advance; 9 end if 10 end if 11 Return Trapped; The Algorithm 3 1 i s developed by Kuffner, J.J., and LaValle, S.M. [30] The RRT grows as follows. A random configuration in the configuration space with position as is chosen initially using the RANDOM_CONFIG function. Then a vertex fr om the RRT, that is nearest to this configuration is found using the PAGE 42 42 NEAREST_NEIGHBOR function. The NEW_ CONFIG function is called to exp and the tree towards the random configuration with a small epsilon value. The new vertex and the new edge connecting with are added to the tree. If this new vertex is the goal vertex, then the goal point is reached and a path between the start and the goal vertex is formed. If the new vertex is not the goal vertex the tree keeps growing, until the goal is reached. Also, if the new vertex c annot join with the nearest vertex due to collision with an obstacle in the path, then the new configuration is rejected and the algorithm returns a trapp ed Thus a new configuration is chosen at every iteration until the RRT reaches the goal vertex Bi Directional RRT Algorithm Bi d irectional RRT is a variation of the original RRT algorithm. It has two trees growing towards each other; one startin g from the start point and the other from the goal point. This makes it more efficient and fast in finding the solution. The main difference between the function RRT_BIDIREC and RRT are listed below Bi directional RRT grows two trees originating from the source and from the destination, until they meet each other. So both the trees have to be tracked simultaneously at every iteration. Both the trees grow towards each other instead of growing towards some random configuration. This increases the probabi lity of finding a solution and also makes it faster. In original RRT algorithm the tree grows towards some random configuration with a single step of epsilo n length in each iteration. In b i directional RRT, the trees grow towards each other with multiple epsilon length steps, thus making the greediness of the algorithm stronger. The A lgorithm 3 2 is also developed by Kuffner, J.J., and LaValle, S.M. [30] which is similar to the original RRT algorithm except for the difference in the number of trees PAGE 43 43 and an additional function for connecting the two trees. Figure 3 2 depicts a b i directional RRT algorithm with its two trees growing from the start and goal node. Figure 3 2. Bi directional RRT algorithm s olutio n with four obstacles and two trees The tree has the start point as its initial vertex and the tree has the goal point as its initial vertex The RRT BIDIREC function call s an additional function called CONNECT. This CONNECT function only iterates over the EXPAN D function until the connection is made and the EXPAN message First a random configuration is chosen. The EXPAN D function is called pass ing and the goal point, as its parameters. The nearest vertex to is chosen from the tree, such that it tries to grow the tree, towards the goal point, If there are no collisions with the obstacles it goes to the next step. In this step the tree, is passed to the CONNECT function which calls the EXPAND function iteratively to grow towards the new vertex, that is added to the tree If the CONNECT function fails to reach this new vertex due to an obstacle in between ( i.e. the two trees fails to connect with each other), then the roles of the two trees are exchanged by using the SWAP function, which swaps the two trees. The SWAP function which is called in RRT_BI DIREC function takes two parameters and trees PAGE 44 44 Algorithm 3 2 RRT Bi directional ( ) RRT_BI DIREC ( ) 1 .start ( ) ; .start ( ) ; 2 f or j =1 to n do 3 ; 4 if not (ADD ( )= Trapped) then 5 if (CONNECT( then 6 Return PATH ( ; 7 end if 8 end if 9 SWAP ( ; 10 end for 11 Return Failure; EXPAND ( T ) 1 ; 2 if NEW_CONFIG ( then 3 T.add_vertex ( ; 4 T.add_edge ( ; 5 if then 6 Return Reached; 7 else 8 Return Advance; 9 end if 10 end if 11 Return Trapped; CONNECT ( T ) 1 repeat; 2 EXPAND ( T ); 3 until not ( ; 4 Return Q; Now the tree with initial vertex as tries to grow towards the goal point and the tree with initial vertex as tries to connect with the new vertex of the other tree. Hence the roles of the 2 trees are swapped in each iteration, with the first tree trying to PAGE 45 45 connect with the goal point and the second trying to connect with the first tree; until the two trees connect with each other. An impo rtant advantage of this approach is that it is very flexible with its multiple epsilon steps and the role swapping. There can be many possible variations in the solution. It explores a different solution every time and so it not a single solution approach. Also the solution is obtained much faster in a cluttered environ ment and slightly faster in a more cluttered environment. It can find solution to 2D or 3D problem with obstacles in seconds. Hence Bi directional RRT is an improved version of the original RRT algorithm which is developed for faster convergence and for use in 2D and 3D problems. However the main drawback of this approach is t hat there a lot of nearest neighbor and c ollision check function calls to be made in order to grow the trees and connect them in each iteration ; thus making it computationally more expensive. Multiple RRT Algorithm The RRT and Bi directional RRT algorithm are adept at solving very difficult and high dimensional space path planning problems. But increasing the obst acles in the environment can severely affect the convergence in these algorithms. A single or two RRTs does not affect the rate of convergence to a solution, which is largely determined by shape and position of the obstacle in the environment. For exa mple, if a start and goal point is separated by a wall with a narrow hole or passageway in it, the RRT can get lost on one side of the obstacle. This can increase the number of iterations exponentially, thus increasing the calls to the nearest neighbor and coll ision check function. This can increase the computat ion effort and time significantly PAGE 46 46 Figure 3 3 Multi ple RRT algorithm s olution with three wall obstacle The m ultiple RRT algorithm i s proposed by Clifton, M. and others [31] to address these issues. The paper introduces the method of using multiple trees and thus reducing the logarithmic complexity of the search in a high obstacle density environment. Figure 3 3 shows a m ultiple RRT algorithm solution with three wa ll obstacles. The g reen colored solution is the path found after smoothening of the initial red colored path is done. The pseudo code of an RRT and m ult iple RRT algorithm is given in A lgorithm s 3 3 and 3 4. The RRT pseudo code has only a single t ree, which grows in every iteration towards the goal point. In multiple RRT algorithm, after every c ollision c heck, a new tree is created if a collision occurs These trees try to connect with only the nearest neighbor node of the other trees at each itera tion. This way it eliminates the need to connect with PAGE 47 47 all the nodes of all the trees and hence reduces the computationally expensive collision c hecks for each connection. Also it only creates a new tree if it is required and tries to merge these trees at e ach iteration, to reduce the total number of trees. Algorithm 3 3 RRT 1 do while Path not found 2 NEW_POINT 3 NEAREST_NEIGHBOR 4 COLLISION_CHECK 5 if no Collision 6 CONNECT 7 Path Found 8 end if 9 end do Algorithm 3 4 Multiple RRT 1 do while Path not found 2 NEW_POINT 3 for each tree 4 NEAREST_NEIGHBOR 5 COLLISION_CHECK 6 if no Collision 7 CONNECT 8 else 9 NEW TREE 10 end if 11 end f or 12 end do The algorithm initially places discretely positioned se eds throughout the environment, which may or may not become the roots of new trees. It can be used for both 2D and 3D search space with obstacles. It tries to automatically connect to the goal and does path smoothing to find the shortest possible path to t he goal. The user can specify the maximum number of trees to be cr eated, with a minimum of 2. The user PAGE 48 48 can also specify the number of seeds to be placed on each axis Also the user can specify the maximum number of RRT search iterations to be done to find the path. The search space is specified by giving its limits in each axis, and also the start and goal locations have to be specified within this search space. The obstacles can be specified in the form of walls, with more than one wall allowed in the environment. The vertices of the obstacles are specified in an obstacle array file, which can be updated with new obstacles at any time. New_Point function finds a point randomly from the search space. The Nearest_Neighbor function finds the vertex of the tree that is nearest to the new point. The Collision_Check functions is called to check if there is any collisions with the obstacles while connecting the new point with the nearest neighbor vertex. If no collision occurs, the new point is added to the tree as its new vertex and along with the edge joining it with the nearest neighbor. The tree also tries to connect with the nearest neighbors of the other trees in the process. If a collision occurs, it creates a new tree with the new point as its root. The process repeats u ntil a tree reaches the goal point Once the goal point is reache d the path is smoothened to give the shortest path to the goal. The rate of convergence is greatly improved by using m ultiple tree RRT a lgorithm especially in environments with high obstacle density. Also the computational efficiency is better when compa red to the b i directio nal RRT algorithm due to fewer collision c hecks and less number of iterations. Another advantage of m ultiple RRT a lgorithm is that the search is not restricted to near the start and goal points, and hence it can be useful for path pl anning with multiple goal points. PAGE 49 49 Application of RRT in Multiple Aircraft Landing Problem Multiple RRT method has the advantage of reducing the computational complexity by using multiple trees which t ry to connect with each other in every iteration. In the application proposed, this multiple RRT is used to find the trajectory to be taken by multiple aircraft while the y are about to land. A new algorithm is develope d which takes into account the path constraints of avoiding the obstacles initially and is dis cussed in Chapter 4 The aircraft are sorted in order of their distance from the goal and the nearest aircraft lands first followed by the next nearest aircraft. This algorithm not only finds the shortest distance path of all the aircraft but also satisfie s additional constraints of keeping a minimum distance and time between aircraft landings. PAGE 50 50 CHAPTER 4 INITIAL TRAJECTORY G ENERATION USING RRT Parameters for the P roblem The problem in this work is three d imensional, so m ulti ple RRT with 3 D trajectory generation is used. The obstacles are defined in the obstacle text file of the RRT algorithm as vertices of a quadrilateral. This text file can be updated anytime with a new s et of obstacles, and the proposed algorithm can be executed again to g et a new trajectory, taking into c onsideration the new obstacles. A simple RRT algorithm is used to find the shortest path from a start point to a goal point avoiding the obstacles. However in this work the RRT algorithm is used to find the trajectory of ea ch of the aircraft, avoiding obstacles and maintaining the minimum safe distance between each landing But the problem into consideration has to include many constraints in the trajectory generated by the RRT algorithm. Also the trajectories of the aircr aft are related to each other. That is the trajectory of all the aircraft has to be generated taking into consideration the trajectory of the previous aircraft to include the minimum distance and time difference constraint in its landing Also intersectin g flight paths should be avoided So another algorithm is developed in this work which takes into consideration all these constraints for any number of aircraft However, increasing the number of aircraft would increase the computation time. In this work three aircraft and eight obstacles were taken into consideration. The user has to specif y the number of aircraft the obstacles and the initial positions of the aircraft The initial position of all the aircraft can be found using the r adar positioning system. The time at which the initial position is found, is taken as the initial time, in the problem. PAGE 51 51 Algorithm for Initial Trajectory G eneration Initialization The number of aircraft approaching to land is gi ven as n. (n=3 in th is problem ) The initial location of all the aircraft is measured using radars signals. This initial location is calculated such that it gives the position in which the aircraft would be after 60 seconds of moving in its own trajectory. This may be the comp utation time taken to pre calculate trajectory for the aircraft The goal location is specified, which is the position the aircraft should reach to start its descent to the runway. This is considered as th e same for all the aircraft ( taking into considera tion only one runway). The limit, lim specifies the airspace within which all the aircraft should reach to be considered for landing. (i.e the start positions of all the aircraft should be within this limit). Identification of the N earest A ircraft Initiall y the RRT algorithm is used to find the shortest distance path from each of the start points of the aircraft to the goal point. The RRT algorithm gives only the way po ints of the path as its output, so the total distance from the start to the goal point is found by adding up the distance between two successive points. Based on the total distan ce between start and goal point the aircraft with the nearest distance from the goal point is found. These aircraft are then sorted based on the minimum distance req uired to reach the runway while avoiding obstacles. Specif ication of the C onstraints Keeping the path of the first nearest aircraft as the same, the trajectory of rest of the aircraft has to be manipulated based on the constraints specified. For this, a waypoint is chosen randomly within the limit specified by the user initially so that the rest of the aircraft (leaving the first aircraft) reach this way point first, before reaching for the goal point. This point is chosen such that the total distance f or these aircraft to reach the goal point, should satisfy the constraints given below 1. Minimum distance between two successive aircraft landing. 2. Maximum distance between two successive aircraft landing. PAGE 52 52 3. The distance from the start point to way point sho uld be greater than half the total distance from start to goal point. The first two constraints gives a margin of safety that can be specified, while the third constraints is used to see if the way point is chosen wisely, instead of just choosing it rando mly. It helps in p reventing intersecting trajectories at th e same time, thus avoiding collisions. Selection of the W ay P oint The waypoint is chosen repeatedly to satisfy these constraints. To improve this search, it is made sure that none of the waypoints selected are the same as those selected before. For this the waypoints are stored in an array, and its elements are compared to the new random waypoint, to check that they are not equal. The RRT algorithm is used to find the trajectory from the start poi nt to the waypoint. Also the distance from the start to the way point is found as mid distance. This RRT algorithm is used again to find the trajectory from the waypoint to the goal point. The distance from waypoint to the goal point is found. This distanc e is sum med with the mid distance to obtain the full distance. If all the constraints are met, it exits out of the loop; else it keeps searching for a suitable waypoint to satisfy all the constraints. The total trajectory of the aircraft is found by adding the trajectory generated by both the RRT algorithm. This is repeated for all the aircraft except for the nearest aircraft, to find the new trajectory, which satisfies all the constraints. Aircraft T rajectory and T ime The new trajectory of the aircraft except for the one nearest to goal point is found. As the speed is assumed to be constant for all the aircraft the time at the end of each point can be found based on the distance travelled from one point to an other by using the speed, distance, time relationship However, the new trajectory may have one or more repeated points, (the waypoint will be repeated at least twice due to summing the trajectories). Such repetitions of points in the trajectory are removed and the new times associated with thes e points are found. Also the time at which the waypoint is to be reached is found, and the index value of this time, from the time array is found. PAGE 53 53 Advantages and Disadvantages of the Initial Trajectory Generated The initial trajectory of the aircraft and the time associated with it is found usi ng the developed algorithm. The time taken by aircraft is the minimum possible time it takes to reach the goal avoiding the obstacle and without considering the aircraft dynamics. It can be used to find the initial trajectories of any number of aircraft which has to be specified along with the start positions of the aircraft Also any number of 2D or 3D obstacles can be added in the trajectory. In real time, the obstacles can be updated to generate a new initial tra jectory. As the way points are selected randomly, the trajectory generated every time the algorithm is executed, will not be the same. Hence the algorithm does not restrict itself, to only one solution. Instead it searches for a new trajectory each time i t is executed. However, this also means that the computation time to find the initial trajectory may vary and is not fixed. T rajectory P lots Three Dimensional Trajectory The 3 D trajectory of all the aircraft can be plotted using the plot3 option in Matlab The plot in Figure 4 1 shows the trajectory of the three aircraft in three dimension The original trajectory and the new trajec tory that is generated for the seco nd and the third aircraft using th e a lgorithm developed are shown in Figure 4 1 It is clearly seen, that by going through the way poi nt, the trajectory of the second and thi rd aircraft becomes longer than their original trajectory. The new trajectory satisfies all the minimum and maximum distance difference constraints between each a ircraft landing. PAGE 54 54 Chec k for Intersecting Trajectories The new trajectory generated should not only avo id obstacles i n the path, but also should ensure that the trajectories do not intersect at any p oint of time, thus avoiding collision s To check this, the x, y and z positions of all the aircraft with respect to time is plotted From Figure s 4 2 4 3, and 4 4 it can be seen that though there may be overlaps in the y and z positions at different point s of time, the x position does not intersect at any point of time. This shows that the trajectory do es not intersect with the same x y and z positions at any instance of time Also even if the x position intersect s at any point of time, the instant at which the intersectio n occurs should match with the in stant of intersection in the y and z posit ions too, to have intersection s in the trajectories This co ndition almost never occurs because of th e way the constraints are formulated Hence the algorithm also avoids any intersecting trajectories and thus avoids any collisions between aircraft during landing. The initial trajectory generated can be used for initializing the optimal control problem. The aircraft can be arranged in the order of its distance from the goal point, so that the time constraints can be specified accordingly. The initial trajectory takes care of the obstacle avoidance constraint, which form s the path constraint in the optimal control problem. This helps in better convergence of the problem. PAGE 55 5 5 Figure 4 1. Initial three dimensional t rajectory of the three aircraft using the RRT and proposed algorithm PAGE 56 56 Figure 4 2. X axis vs. Time plot for the three aircraft using the RRT and proposed algorithm Figure 4 3. Y axis vs. Time plot for the three aircraft using the RRT and proposed algorithm PAGE 57 57 Figure 4 4 Z axis vs. Time plot for the three aircraft using the RRT and proposed algorithm PAGE 58 58 CHAPTER 5 PSEUDOSPECTRAL METHODS AND GPOPS Pseu dospectral methods are a direct trajectory optimization method which uses direct collocation to transcribe the optimal control problem to a nonlinear programming problem (NLP). The optimal control problems can be solved using collocation at Legendre Gauss(LG), Legendre Gauss Radau (LGR), and Legendre Gauss Loba tto (LGL) points is presented by Garg, D. and others [27] Pseudospectral methods uses global polynomials to parameterize th e state and control and then uses the nodes of the Gaussian quadrature to collocate the differential algebraic equations. It provides accurate solutions for problem with smooth solutions. For problem with solutions which are not smooth the time interval f rom [ 1,1] is broken into several intervals, so that a different global polynomial approximation is used over these intervals. LG, LGR, and LGL Collocation Points There are three most common sets of collocation point which are obtained from the roots of a Legendre polynomial and/or linear combinations of a Legendre polynomial and its derivatives; they are Legendre Gauss (LG), Legendre Gauss Radau (LGR), and Legendre Gauss Labatto (LGL) points. All the three sets are defined in the domain of [ 1,1]. The m ajor difference is that the LG points do not include any endpoints, the LGR points include one of the endpoints, and the LGL points include both of the endpoints. Also the LGR points are asymmetric relative to the origin and can be defined using the endpoi nt as either the initial point or the end point. Figure 5 1, shows a schematic of the LGL, LGR and LG collocation p oints. PAGE 59 59 Figure 5 1 Sch ematic showing LGL, LGR and LG collocation p oints Let N be the number of collocation points and be the degree Legendre polynomial, then the LG, LGR and LGL collocation points are obtained from the roots of the polynomial. LG points are the r oots obtained from LGR points are the r oots obtained from and LGL points are the r oo ts obtained from together with the points 1 and 1 Formulation of Pseudospectral Method Using LGR Points A method for direct trajectory optimization and costate estimation of finite horizon and infinite horizon optimal control problems using globa l collocation at Legendre Gauss Radau (LGR) points is shown by Garg, D. and others [26] The authors have shown that the use of LGR collocation helps in determining accurate primal and dual PAGE 60 60 solutions for both fin ite and infinite horizon optimal control problems. So the Radau Pseudospectral method is used in this work To simplify the problem, an unconstrained optimal control problem on the time interval is taken The time interval can be transformed via the affine transformation The goal is to determine the state ( ) and the control u( ) which minimize the cost functional subject to the constraints w here and is the known initial condition Let there be N LGR collocation points, described in the interval [ 1,1], where and A new non collocated point which is used to approximate the state variable is introduced [26] Each component of state x is approximated by a Lagrange polynomial, as given in equation below w ith The j th component of the state is approximated as a series, PAGE 61 61 Differentiating the above series at the collocation points gives, w here The N by N+1 matrix D is called the radau differentiation matrix. C ollocating the dynamics at the N LGR collocation points gives The finite dimensional NLP is then formulated as below, Minimize s ubject to the system dynamics, This system dynamics is then rewritten as, wh ere is the N by N differentiation matrix which is invertible hp Adaptive Pseudospectral Method An hp adaptive pseudospectral m ethod with collocation at Radau points is chosen in this work This method is proposed by Darby, C.L. and others [28] which iteratively and sim ultaneously determines the number of segment breaks, the width of each segment, and the polynomial degree required in each segment for approximation, until the user specified accuracy is achieved. It leads to higher accuracy solutions with less computation al effort and memory than is required in global pseudospectral methods. Also the Radau points are used as it is more computationally efficient and well posed. As the problem in this work has multiple aircraft and also has an obstacle wall in the path to be avoided, the problem is broken into several overlapping phases, for each PAGE 62 62 aircraft to cross the obstacle at the desired location and then reach the goal. Due to the complexity of the problem, the hp adaptive method is chosen. Multiple Aircraft Landing Multi ple Phase Problem Formulation Objective and Input Parameters The objective is to land the three aircraft one after the other, while maintaining a minimum time and distance separation between successive aircraft landings Also the initial and final de sired position and orientation specified in the problem has to be achieved The final position in this problem is the same for all the aircraft as they have to land on a common single runway. This position is [0, 0, 0], and the limits of the environment within which the aircraft has to be is [ 3 3; 0 6; 0 0.3] with the dist ance considered in n autical miles. Along x axis, it ranges from 3 to 3 n autical miles (Horizontal breadth), along y axis it ranges from 0 to 6 n autical miles (Horizontal length) and along the z axis it ranges from 0 to 0.3 n au tical miles (Vertical height). The initial and final time is taken in seconds. T h e flight path and heading angle of the aircraft are taken in radians in the problem and speci fied in degree in Table 5 1 Table 5 1. User specified initial and final s tates of the Aircraft 1, 2, and 3 Aircraft 1 0 1.2 5.8 0.28 0 1 5 155 0 0 0 0 90 2 0 1.8 6 0.25 5 10 170 0 0 0 0 9 0 3 0 3 5.5 0.22 0 10 5 185 0 0 0 0 9 0 Methodology The problem of aircraft landing is solved using three aircraft The problem is divided into 6 o verlapping phases. In P hase I all the three aircraft are hovering before the obstacle wall and are assumed to start at the same time at t=0. At the end of the PAGE 63 63 Phase I the nearest aircraft reaches t he obstacle wall. In P hase II only the second and third aircraft are hovering before the obstacle wall. At the end of P has e II second aircraft reaches the obstacle wal l. Phase III has only the third aircraft still hovering before the obstacle wa ll. Phases IV, V and VI depict the three aircraft landing and reac hing the goal (runway) In P hase IV the first aircraft is landing; the P hase V has the second aircraft landing, while the final P hase VI has the third aircraft landing. Thus t he problem is divided into the six phases. Figure 5 2 Schematic showing the six phases in the problem These phases can be seen in F igure 5 2, where and are the final landing time of Aircraft 1, 2 and 3 respectively; and is the time taken by Aircraft 1, 2 and 3 to reach the obstacle wall respectively. PAGE 64 64 These phases are linked with each o ther in GPOPS using the linkage constraints in order to make the states continuous. In this problem, as seen in Figure 5 1, at time, the first five states of Aircraft 1 in Phase I is linked with Phase IV and the last ten states of A ircraft 2 and 3 in P has e I is linked with P hase II At time, t he first five states of A ircraft 2 in P has e II is linked with P hase V and the last five states of A ircraft 3 in P hase II is linked with P has e III At time, all the states of A ircraft 3 in Phase III is linked with P hase VI. The dynamics of each of the aircraft are given separately. This gives an advantage of generating trajectories for aircraft with different dynamics such as different types of aircraft, helicopter, small UAVs and other air vehicles In the example problem, Phase I has totally fifteen states, with five states for each of the three aircraft and has total of six c ontrols (two for each aircraft). Similarly Phase II has a total of ten states and four cont rols for A ircraft 2 and 3 The P hases III, IV, V, and VI have five states and two controls for a single aircraft Cost Function Formulation The cost function for the problem with number aircraft is formulated as, where is the time taken by the aircraft to reach the face of the obstacle, N is the total number of aircraft taken into consideration, n is the number of aircraft present in P hase j and 2N is the total number of pha ses formulated in the problem. The Mayer cost function is the time for the aircraft to reach the obstacle wall, which is being maximized in this problem. The Lagrange cost in this problem is PAGE 65 65 minimizing the change in flight path angle and heading angle of t he aircraft. Maximizing the time to reach the obstacle and fixing the final time of landing of the aircraft, ensures that the aircraft satisfies the time difference constraints by hovering around before it reaches the obstacle wall, and takes an almost str aight path to the goal point, after it passes the obstacle wall. Also minimizing the change in flight path and heading angle ensure s less control effort required by the aircraft and also generates smooth aircraft trajectories Constraints Formulation Ther e are different constraints specified in each of the phases depending on the number of aircraft in that phase, the position of aircraft (reaching the obstacle wall or landing at the runway) and the aircraft dynamics There are four main constraints th at are taken into consideration (1) Obstacle path constraints, (2) Intersection constraints, (3) Event Constraints, and (4) Bank angle constraint. Obstacle path constraint s formulation The path constraint that is common to all the phases is the obstacle w all which has to be avoided in the aircraft trajectory. As each aircraft has different dynamic states, the obstacle constraints have to be modeled for every aircraft that present in each phase. In the e xample discussed, Phase I has twelve obstacle path constra ints, four for each of the three air craft. Similarly, Phase II has eight obstacle path constraints, and the rest has four obstacle path constraints. As seen before, the obstacles are modeled as below, PAGE 66 66 Intersection constraints formulation The i ntersection constraints has to be taken care of in the P hases I and II as there are more than one aircraft flying in these phases The example has all the three aircraft flying in Phase I. T he intersection constraints are modeled by finding the distance between the aircraft at e very point of time whi ch has to be always greater than the safe separation distance In this problem, the saf e distance is taken to be 0.25 n autical miles. The distance between the first and second aircraft forms the first intersection constraint; similarly distance between first and third aircraft forms the second intersection constraint ; and the distance between second and third aircraft forms the third intersection constraint. In P hase II, however there is only one intersection constraint which i s the distance between the second and the third aircraft. Let A1 be Aircraft 1, A2 be Aircraft 2 and A3 be Aircraft 3. Now let, inter11 = Minimum distance between A1 and A2 in Phase I, inter12 = Minimum distance between A1 and A3 in Phase I, inter13 = Mi nimum distance between A2 and A3 in Phase I, i nter21 = Minimum distance between A2 and A3 in Phase II. t hen the intersection constraint in Phase I is specified as, [0.25; 0.25; 0.25] PAGE 67 67 and intersection constraint in Phase II is specified as, [0.25] Event constraints formulation The event constraints are enforced at the end of Phase I, II, and III on Aircraft 1, 2 and 3 respectively These constraints are formulated to ensure that the aircraft passes only through the hole in the obstacle wall and it sho uld be orient ed parallel to the ground and perpendicular to the wall, when it reaches the face of the obstacle wall To ensure the x, y, z positions are within the limit of the hole in the wall, the constraints on these positions are specified at the time the aircraft reach es the face of the obstacle wall To satisfy the orientation constraint of the aircraft, the f light path angle has to be oriented parallel to the ground and the h eading angle should be exactly perpendicular to the obstacle wall. This is given by, Bank angle c onstraint s formulat ion Th e control is constrained by constraining the bank angle of the aircraft in its trajectory. The bank angle is related to the control as given below PAGE 68 68 where is the horizontal load factor and is the vertical load factor of th e aircraft In the example, the bank angle of all the three aircraft is constrained as in degree angle. To enforce this constraint the inverse tangent of ratio of the two controls has to be in the range As Phase I has three aircrafts, there are three bank angle constraints to be satisfied, Phase II has two aircrafts with two bank angle constraints and all the other phases has a single aircraft with a single bank angle constraint. The results are p lotted a nd discussed in Chapter 6 The three dimensional trajectory is plotted, along with all the states and controls of the three aircraft These results are analyzed and validated. PAGE 69 69 CHAPTER 6 RESULTS AND DISCUSSI ON The multiple aircraft landing problem with three aircraft is solved by converting the opt imal control problem into a non linear p rogramming (NLP) problem using GPOPS. This NLP is then solved using SNOPT. The solution obtained for the optimal control problem is plotted Also the results are analyzed and validated in this chapter. Analysis of Results The t rajectory of the three a ircraft is plotted in a three d imensional plot in Figure 6 1. It portrays the trajectory of each of the three aircraft from start to the goal point. The blue windo w shows the obsta cle wall in the environment which the aircraft has to avoid in its trajectory. Aircraft 1 has two phases it its trajectory, which are shown in different colors. Simila rly, Aircraft 2 has three phases in different colors and Aircraft 3 has four phases in different colors. Also it can be seen that P hase I which has all the three aircraft is depicted in blue color; P hase II with two aircraft is depicted in magenta color; P ha se III with one aircraft is depicted in green color, similarly P has e IV is red color; P hase V is cyan color and P hase VI is black color. Aircraft 1 reaches the obstacle at (1, 3, 0.2); Aircr aft 2 reaches the obstacle at (1, 3, 0.1); and Aircraft 3 reaches the obstacle at ( 1, 3, 0.15). Hence it satisfies the path constra ints of flying through th e obstacle hole and avoiding collisions with the obstacle wall. The x, y and z axis plot of all the three aircraft is plotted in Figures 6 3 to 6 5 They reach the goal at 165, 170, 185 seconds, which is about 5 seconds added to the solutions obtained using the in i tial trajectory guess using the RRT algorithm H ence the PAGE 70 70 software finds the shortest possible time taken by the aircraft to reach the runway while avoiding the obstacles and also taking into consideration the dynamics of the aircraft The f light path angle and the heading angle of all the aircraft is plotted in Figure s 6 6 and 6 7 The smooth curve indicates the cost of minimizing the change in the flight path angle and heading angle is satisfied Also the event c onstrai nts at the end of the P hases I, II, and III for the flight path angle and the heading angle are satisfied Validation of Results Cost Function Validation The Mayer cost of maximizing the time taken by the aircraft to reach the obstacle is satisfied and ca n be seen in the Figures 6 3 to 6 13. This can also be verified from Table 6 1, which displays the time taken by aircraft to reach the obstacle, and the final fixed time to land on the runway, for all the three aircraft The smooth curve of flight path angle and heading angle as seen in Figures 6 6 and 6 7 indicates that the Lagrange cost of minimizing the change in flight path angle and heading angle is also satisfied. Table 6 1 Mayer cost v alidation Aircraft Number in seconds in seconds 1 79.77 155 2 94.9 6 170 3 109.9 6 185 Constraints Validation The constraints formulated in Chapter 5 are validated from the plotted results Obstacle path constraints v alidation The path constraint is the obstacle wall to be avoided in the trajectory. The values of each of these obstacles are found at each phase of the aircraft. These values satisfy PAGE 71 71 the obstacle path constraints specified in Chapter 5. This can also be seen in the three dimensional trajectory of the three aircrafts shown in Figures 6 1 and 6 2, where the blue window is the obstacle wall with the hole in the center. Event c onstraints v alidation The aircraft sat isfie s the event constraints of flying only through the hole in the obstacle wall and being within the limit of [ 1, 1] along x axis, [0.1, 0.2] along z axis when it reaches the face of the obstacle wall. The obstacle wall face is placed at 3 n autical miles along the y axis. Also the state v a lue of heading path angle and f light path angl e for Aircraft 1, 2 and 3 at the end of P hase s I, II, and III respectively is found to be 0 and 1.5708 as seen in Figures 5 5 to 5 8 This s atisfies the event constraint s of flight path angle and heading ang le specified in Chapter 5. Intersection c onstraints v alidation The Figure 6 14 shows the minimum distance between the pairs of aircraft in each phase. The minimum s afe separation distance value is taken as 0.25 nautical miles, but the bar graph shows that the actual minimum distance for the four intersection constraints specified in Chapter 5 is much more than the minimum safe separation distance of 0.25 nautical miles Bank angle c onstraints v alidation The bank angles of the three aircraft are plotted in Figures 6 11 to 6 13. It is seen that the maximum and minimum bank angle values of the three aircraft is and Also it can be noted that only the first aircraft reaches these maximum and minimum bounds. However the maximum and minim um bank angle value bound of the seco nd and third aircraft is much smaller than that specified in Chapter 5 PAGE 72 72 Figure 6 1. Three d imensional t rajectory with the obstacle wall obtained using GPOPS PAGE 73 73 Figure 6 2 Y Z axis view of the t hree d imensional t rajectory PAGE 74 74 Figure 6 3 X axis vs. T ime plot of the three aircraft in different phases Figure 6 4 Y axis vs. T ime plot of the three aircraft in different phases PAGE 75 75 Figure 6 5 Z axis vs. T ime plot of the three aircraft in different phases Figure 6 6 Flight p ath a ngle vs. T ime plot of the three aircraft in different phases PAGE 76 76 Figure 6 7 Heading a ngle vs. T ime plot of the three aircraft in different phases Figure 6 8. Load f actor c ontrols of Aircraft 1 PAGE 77 77 Figure 6 9 Load f actor c ontrols, of Aircraft 2 Figure 6 1 0. Load f actor c ontrols, of Aircraft 3 PAGE 78 78 Figure 6 1 1 Bank angle of Aircraft 1 in different phases Figure 6 1 2 Bank a ngle of Aircraft 2 in different phases PAGE 79 79 Figure 6 1 3 Bank a ngle of Aircraft 3 in different phases Figure 6 1 4 Intersection c onstraints v alidation Minimum d istance between a ircraft in Phase I and II with A1 Aircraft 1, A2 Aircraft 2, and A3 Aircraft 3 0.25 0.25 0.25 0.6528 0.8565 0.6335 0.6528 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 A1&A2 A1&A3 A2&A3 Minimum safe separation distance Minimum distance in Phase I Minimum distance in Phase II PAGE 80 80 CHAPTER 7 CONCLUSION The work presented an approach to solve the multiple aircraft landing problem using a combination of RRT algorithm and p seudospectral methods. The RRT algorithm has been us ed to provide initial guess to the optimal contro l software, GPOPS, which helped in better convergence of the solution. An example problem to generate the landing trajectory of three aircraft has been solved. The problem and the obstacles has been formulated in Chapter 2. The initial trajectory genera ted using the proposed algorithm satisfied the minimum distance and time difference constraint of landing as seen in Chapter 4. The multiple phase optimal control problem ha s been formulated in Chapter 5. From the results plotted in Chapter 6, it was found that the solution satisfied all the constraints in the problem and gave an optimal solution It was also found that the aircraft final landing time was only 5 seconds more than that obtained using the proposed algorithm used for initial guess. So the aircraft takes only 5 seconds more than the short est possible time required by aircraft to reach the runway while avoiding the obstacle wall. The extra 5 seconds were used constraints, event constraints, intersection constraints and bank angle constraints in the flight path Thus the initial guess for the final landing time of the aircraft was found to be very close to the actual optimal solution to the problem. This was help ful as the proble m has been formulated with a fixed land ing time for the three aircraft and a fixed time difference between successive aircraft landings. T his research work can be used to assist the air traffic controllers in generating the optimal trajectories for aircraft which are ready to land and to improve the efficiency of the airports. Also it can help the ATC to avoid any obstacles such as buildings or no PAGE 81 81 fly zone s while generating these trajectories. Future work on this research can de crease the computation time in generating these trajectories by optimizing the number of collocation points used to solve the optimal control problem Also faster NLP solvers and better CPU processors can improve t he computation time significantly This can help in real time implementation of the work PAGE 82 82 LIST OF REFERENCES [1] Lecchini, A., Glover, W., and Lygeros, J., "Air Traffic Control in Approach Sectors: Simulation Examples and Optimisation," Hybrid Systems: Computation and Control, Series Lecture note s in Computer Science, Vol. 3414, Jan. 2005, pp. 433 448. [2] Kuchar, J.K., and Yang, L.C., "A Review of Conflict Detection and Resolution Modeling Methods," IEEE Transactions on Intelligent Transportation Systems, Vol. 1, No. 4, Dec. 2000, pp. 179 189. [3] Frazzoli, E., Mao, Z.H., and Oh, J.H., "Resolution of Conflicts Involving Many Aircraft via Semidefinite Programming," Journal of Guidance, Control, and Dynamics, Vol. 24, No. 1, Jan Feb. 2001, pp. 79 89. [4] Pallottino, L., Feron, E.M., and Bicchi, A ., "Conflict Resolution Problems for Air Traffic Management Systems Solved With Mixed Integer Programming," IEEE Transactions on Intelligent Transportation Systems, Vol. 3, No. 1, Mar. 2002, pp. 3 11. [5] Hu, J., Prandini, M., and Sastry, S., "Optimal Coo rdinated Maneuvers for Three Dimensional Aircraft Conflict Resolution," Journal of Guidance, Control, and Dynamics, Vol. 25, No. 5, Jan. 2002, pp. 888 899. [6] Durand, N., Alliot, J., and Noailles, J., "Automatic Aircraft Conflict Resolution Using Genetic Algorithm," Proceedings of the 1996 ACM symposium on Applied Computing 1996, pp. 289 298. [7] Menon, P.K., Sweriduk, G.D., and Sridhar, B., "Optimal Strategies for Free Flight Air Traffic Conflict Resolution," Journal of Guidance, Control, and Dynamics, Vol. 22, No. 2, Jan. 1999, pp. 202 212. [8] Raghunathan, A.U., Gopal, V., and Subramanian, D., "Dynamic Optimization Strategies for Three Dimensional Conflict Resolution of Multiple Aircraft," Journal of Guidance, Control, and Dynamics, Vol. 27, No. 4, J ul Aug. 2004, pp. 586 594. [9] Dai, R., and Cochran, J.E., "Three dimensional trajectory optimization in constrained airspace," Dissertation, Auburn Univ., 2007, pp. 1 162. [10] Goerzen, C., Kong, Z., and Mettler, B., "A Survey of Motion Planning Algorit hms from the Perspective of Autonomous UAV Guidance," Journal of Intelligent and Robotic Systems, Vol. 57, No. 1 4, Jan. 2010, pp. 65 100. [11] Hurni, M.A., "An information centric approach to autonomous trajectory planning utilizing optimal control techn iques," Dissertation, Naval Postgraduate School, Sept. 2009, pp. 1 296. PAGE 83 83 [12] Lewis, L.R., Ross I.M., and Gong Q., "Pseudospectral motion planning techniques for autonomous obstacle avoidance," Proceedings of the 46th IEEE Conference on Decision and Contro l, 12 14 Dec. 2007, pp. 5997 6002. [13] Lewis, L.R., and Ross, I.M., "A Pseudospectral Method for Real Time Motion Planning and Obstacle Avoidance," AVT SCI Joint Symposium on Platform Innovations and System Integration for Unmanned Air, Land and Sea Vehi cles Florence, Italy, May 2007, pp. 1 23. [14] Lavalle, S.M., "Rapidly Exploring Random Trees: A New Tool for Path Planning," Computer Science Dept., Iowa State Univ., URL: http://msl.cs.u iuc.edu/~lavalle/papers/Lav98c.pdf [cited 11 Nov. 2010] [15] Redding, J., Amin, J., and Boskovic, J., "A Real Time Obstacle Detection and Reactive Path Planning System for Autonomous Small Scale Helicopters," AIAA Guidance, Navigation and Control Confere nce and Exhibit, Aug. 2007, pp. 1 22. [16] Caves, A.D.J., "Human automation collaborative RRT for UAV mission path planning," Thesis, Massachusetts Institute of Technology, May 2010, pp. 1 111. [17] Aoud, G.S., "Two stage path planning approach for desi gning multiple spacecraft reconfiguration maneuvers and application to SPHERES onboard ISS," Thesis, Massachusetts Institute of Technology, Sep. 2007, pp. 1 149. [18] Miele, A., "Optimal Trajectories of Aircraft and Spacecraft," Rice Univ., N91109686, Hou ston, TX. Sponsor: National Aeronautics and Space Administration, Washington, DC; United States, 1990. [19] Miele, A., Wang, T., and Tzeng, C.Y., "Optimal abort landing trajectories in the presence of windshear," Journal of Optimization Theory and Applica tions, Vol. 55, No. 2, 1987, pp. 165 202. [20] Nho, K., and Agarwal, R.K., "Automatic Landing System Design Using Fuzzy Logic," Journal of Guidance, Control, and Dynamics, Vol. 23, No. 2, Mar Apr. 2000, pp. 298 310. [21] Juang, J., and Chio, J., "Fuzzy modelling control for aircraft automatic landing system," International Journal of Systems Science, Vol. 36, No. 2, Jan. 2005, pp. 77 87. [22] Chaturvedi, D.K., Chauhan, R., and Kalra, P.K., "Application of generalised neural network for aircraft landing control system," Soft Computing, Vol. 6, No. 6, 2002, pp. 441 448. [23] Ruffier, F., and Franceschini, N., "Visually guided micro aerial vehicle: automatic take off, terrain following, landing and wind reaction," Proceedings. ICRA '04. IEEE International Conference on Robotics and Automation, Vol. 3, Apr. 2004, pp. 2339 2346. PAGE 84 84 [24] Patel, R.B., and Goulart, P.J., "Trajectory Generation for Aircraft Avoidance Maneuvers Using Online Optimization," Journal of Guidance, Control, and Dynamics, Vol. 34, No. 1, J an Feb. 2011, pp. 218 230. [25] Benson, D.A., Huntington, G.T., and Thorvaldsen, T.P., "Direct Trajectory Optimization and Costate Estimation via an Orthogonal Collocation Method," Journal of Guidance, Control, and Dynamics, Vol. 29, No. 6, Nov Dec. 2006, pp. 1435 1440. [26] Garg, D., Patterson, M., and Francolin, C., "Direct trajectory optimization; costate estimation of finite horizon; infinite horizon optimal control problems using a Radau pseudospectral method," Computational Optimization and Applicat ions, Oct. 2009, pp. 1 24. [27] Garg, D., Patterson, M., and Hager, W.W., "A unified framework for the numerical solution of optimal control problems using pseudospectral methods," Automatica, Vol. 46, No. 11, Nov. 2010, pp. 1843 1851. [28] Darby, C.L., Hager, W.W., and Rao, A.V., "An hp adaptive pseudospectral method for solving optimal control problems," Optimal Control Applications and Methods, Aug. 2010, pp. 1 41. [29] Lewis, L.R., and Naval Postgraduate School (U.S.), "Rapid motion planning and autonomous obstacle avoidance for unmanned vehicles," Thesis, Naval Postgraduate School, 2006, pp. 1 161. [30] Kuffner, J.J., and LaValle, S.M., "RRT Connect: An Efficient Approach to Single Query Path Planning," Proceedings IEEE International Conference on Robotics and Automation, Vol. 2, Apr. 2000, pp. 995 1001. [31] Clifton, M., Paul, G., and Kwok, N., "Evaluating Performance of Multiple RRTs," IEEE/ASME International Co nference on Mechtronic and Embedded Systems and Applications, Jun. 2008, pp. 564 569. PAGE 85 85 BIOGRAPHICAL SKETCH degree in Electro nics and Instrumentation Engineering from Ve lammal Engineering College (Anna University), Chennai, India in April 2009. She was awarded 14 th rank Instrumentation Engineering from Anna University affiliated institutions in April, 2009. She also graduated with a Mast er of Science degr e e in Mechanical Engineering from the Department of Mechanical and Aerospace Engin eering at the University of Florida in May 2011. She completed her thesis under the supervision of Dr. Anil V. Rao and was co advised by Dr. Warren E. Dixon. Her interests lie i n the field of path planning, trajectory optimizati on, guidance and control of non linear systems and nonlinear controls 