Power transmission line inspection with UAVs
Abstract
El uso de vehículos aéreos no tripulados (UAVs) para realizar tareas de inspección o transporte de mercancias está en alza, y es por ello que se necesitan algoritmos que sean capaces de planificar la ruta que siguen los UAVs de la mejor manera posible. En este trabajo se presentan distintos heurísticos para el problema de inspección óptima de las aristas de una red eléctrica con un UAV de batería limitada. La red está modelada por un grafo tipo árbol y se plantean dos funciones a minimizar: el tiempo total de inspección o el número de viajes que realiza el dron.
Full text
Proyecto Fin de Carrera Ingeniería de Telecomunicación Formato de Publicación de la Escuela Técnica Superior de Ingeniería Autor: F. Javier Payán Somet Tutor: Juan José Murillo Fuentes Dep. Teoría de la Señal y Comunicaciones Escuela Técnica Superior de Ingeniería Universidad de Sevilla Sevilla, 2013 Trabajo Fin de Grado Grado en Ingeniería Aeroespacial Power transmission line inspection with UAVs Autor: Manuel Francisco Fernández Melgosa Tutor: José Miguel Díaz-Bañez Dpto. Matemática Aplicada II Escuela Técnica Superior de Ingeniería Universidad de Sevilla Sevilla, 2022
Trabajo Fin de Grado Grado en Ingeniería Aeroespacial Power transmission line inspection with UAVs Autor: Manuel Francisco Fernández Melgosa Tutor: José Miguel Díaz-Bañez Catedrático Dpto. Matemática Aplicada II Escuela Técnica Superior de Ingeniería Universidad de Sevilla Sevilla, 2022
Trabajo Fin de Grado: Power transmission line inspection with UAVs Autor: Manuel Francisco Fernández Melgosa Tutor: José Miguel Díaz-Bañez El tribunal nombrado para juzgar el trabajo arriba indicado, compuesto por los siguientes profesores: Presidente: Vocal/es: Secretario: acuerdan otorgarle la calificación de: El Secretario del Tribunal Fecha:
Agradecimientos En primer lugar, quiero agradecer a mi familia, a mis amigos y a mi pareja por haber estado siempre presente en los momentos difíciles. Estos cuatro años de carrera han sido muy complicados y sin ellos no hubiera podido seguir adelante. También quiero agradecer a aquellos profesores que con su enseñanza me han motivado a estudiar y a sacar la mejor versión de mi mismo. Por último, agradecer al profesor Diaz-Bañez por toda la dedicación que me ha brindado en la realización de este proyecto. Manuel Francisco Fernández Melgosa Estudiante de ingeniería aeroespacial Sevilla, 2022 I
Resumen El uso de vehículos aéreos no tripulados (UAVs) para realizar tareas de inspección o transporte de mercancias está en alza, y es por ello que se necesitan algoritmos que sean capaces de planificar la ruta que siguen los UAVs de la mejor manera posible. En este trabajo se presentan distintos heurísticos para el problema de inspección óptima de las aristas de una red eléctrica con un UAV de batería limitada. La red está modelada por un grafo tipo árbol y se plantean dos funciones a minimizar: el tiempo total de inspección o el número de viajes que realiza el dron. III
1 Introduction 1.1 Motivation UAV is the acronym of Unmanned Aerial Vehicles, an aircraft without on-boar crew. The definition of UAV includes both radio-controlled vehicles and autonomous drones. Many terms are used for vehicles without on-board crew. The term drone has been used from the early days of aviation, but this type of autonomous vehicle can be aerial, ground-based or aquatic, whereas UAVs use aerodynamic forces to provide lift. The term Unmanned Aircraft System (UAS) was adopted by the United States department and FAA in 2005, and involves not only the aircraft itself but also other important elements such as control stations, data links and other support equipment. Under the recent regulations in FAA and EASA, the term RPAS (Remotely Piloted Aircraft System) and UAS (Unmanned Aerial System) came up to describe a set of configurable elements consisting of remotely piloted aircraft, its control station, the command and control links and any other system elements required during flight operator [ 3 ]. UAVs can be classified based on weight, degree of autonomy, altitude, endurance and materials. These type of aircrafts were firstly developed at the end of First World War as a target training for anti-aircraft guns, and the production of those increased during Second World War. Also, in this period a different type of UAV can be found, the aerial torpedoes. This missile was the first military aircraft drone, but was not exactly an UAV, as these are intended to be recovered after a mission, and the aerial torpedo is a munition. Recognition UAVs appeared during the Cold War in order to prevent the looses of experience pilots, and at the beginning of the Gulf War the first UAVs with military capability were used. Nowadays, UAVs have many other uses than the military, as they can be implemented in recognition and inspection task, cargo transport, surveillance, filming images and domestic use. For a few years, automation has increased dramatically as well as the use of drones and UAVs. This use is favored by the following characteristics •Dangerous or uncomfortable work environments. •Repetitive work cycles. •Ability to connect with multiple computer equipment. •Capability for reprogramming. The inspection of power transmission lines is one of the perfect tasks to use UAVs. Transmission lines are essential in the distribution process of electricity. Its importance entails an exhaustive 1
2Chapter 1. Introduction inspection to be able to detect any damage. These inspections are currently inspected manly (figure 1.1), which implies some disadvantages: • Low accuracy. The distribution of power lines involves numerous towers connected with kilometers of electric fedders in complex geographical conditions. Therefore, manual inspection is very difficult and the reliability of the inspection cannot be guaranteed. • Danger. Power lines transport around 750 kV. At these levels of voltage, the effects on the human body of the radiation would be fatal. Also, the task must be performed in an aerial vehicle due to the geographical complexity. These things results in a high risk environment for the inspecting personnel. • Economic cost. The inspection is carried out by helicopters with vision monitoring systems. The operation implies that not only a pilot helicopter crew is needed, but also qualified personnel capable of manage the inspection properly. These problems make the use of UAVs and remote sensing sensors to inspect power lines a necessary topic of research and development. UAVs have light weight, small size, low cost and high accuracy, ideal characteristics to perform the inspection. But there are also inconveniences in its use. One of the biggest limiting factors is the battery capacity and the maximum path the UAV can cover before it returns to base. Also, due to the high-voltage environment, electromagnetic interference can appear and the electronic systems of the UAV could be affected and even the UAV itself may be damaged. As a result, it is important to plan and control appropriate flight paths of the UAV for the corresponding missions, maximizing the efficiency and accomplishing the task in a autonomously, safely and reliably way. Figure 1.1 Chopper doing a visual inspection of electric transmission lines in Minnesota, U.S..
1.2 UASs in power transmission line inspection 3 1.2 UASs in power transmission line inspection In power transmission line inspection are used different types of system that can be divide in the UAVs used to survey the installations and the ground stations used to receive all the information that the UAV is collecting, monitor movements and, if necessary, the possibility for a pilot to take control of the aircraft. Furthermore, there can be ground charge stations spread all over the UAV path and data-link towers. Regarding to UAVs, there can be separated in two categories whether are used to inspect tower transmissions or the power lines. The inspection of towers involves local state monitoring such as explosion of insulator string or deformations in the tower structure. This type of work requires stable flight and accurate positioning of the UAV, therefore, a multi-rotor UAV is the better choice to accomplish this task. Multirotors are UAVs with a slow inspection speed, but they can maintain a very stable hover flight, which allows them to capture images with high resolution. The low endurance is compensated by its great maneuverability and easy mounting. This inspection of transmission towers is not studied in this project. For more information about this topic, see [7]. Surveying power transmission lines is the main objective of this work. This type of mission implies the recording of images to detect defects in cables, high trees along the corridors or spot birds nests in all the fedders lines, involving paths of various kilometers. As a high-endurance mission, fixed-wing UAV is the perfect aircraft due to its large autonomous flight capacity and faster speed. The main problem with fixed wing UAVs is its aerodynamic complexity. It is necessary carry out a CFD analysis before the construction to select the correct airfoils and load position. Also, it is needed to compare aerodynamic features with flight data to get control design parameters. Regarding to ground stations, it can be divided in fixed stations and stations mounted in vehicles. The hardware in both types of stations includes typically an human-machine interface, telemetry, video capture cards and a computer to process the information obtained from the UAV. The software is used for planning and flying a mission. Figure 1.2 Hexarotor with a camera, UAV ELIMCO and a portable ground station.
4Chapter 1. Introduction 1.3 State of the art and related work As previously commented, power lines inspection with UAVs is a hot topic and numerous papers and research works have been published over the years in order to give an optimal solution to the path planning of the UAV. With the development of UAV sensor technology, multisensor cooperation can be realized, which improves the accuracy and automation of equipment defect judgment. UAVs carry on infrared and ultraviolet cameras alongside other equipment capable of repairing and carry out maintenance operations. For all different electric companies and investigation studies, the UAV power line inspections can be divided into three load categories: visible imaging load for a global view of the area inspected; infrared imaging load for accurate and efficient detection; and ultraviolet load fro abnormal corona discharge and surface pollution. To maintain a safe distance between the UAV and power lines infrastructure or other constructions, it is used LiDAR technology. It can obtain reliable and accurate three-dimensional topography and feature information of the corridor, mainly including power towers, trees and buildings [ 13 ]. Also it is worth to mention the works [ 14 ], [5] and [8]. Figure 1.3 (a) 3D test platform, (b) Diagram of UAV flight, (c) Infrared images, (d) Computergenerated red dots identifying tower before filtering and (e) after filtering. There are different research projects that explore the possibility of using the power line itself as a charging element for the UAV [ 4 ], either through a kind of charging pole hooked to the line or through different chargers distributed throughout the electric grid. The possibility of connecting the UAV to the electric wire not only has endurance advantages, but also can drive into maintenance operations with the UAV alongside the cables. Furthermore, it has been tested with small solar panels connected on the UAV and with charging UAVs, drones that can refill the UAV battery, similar to air-to-air refueling. Finally, a set of research works related to routing optimization of UAVs is presented. In [ 7 ], genetic simulated annealing algorithms (GSA) are simulated to solve the power inspection line problem optimally alongside genetic algorithms to give an optimal path for the tower inspection problem. To continue, in [ 15 ], a physics model for the energy consumption based on propulsion, anti-collision and data delivery is presented beside an algorithm to minimize the energy consumption on an industrial internet of inspecting and relay UAVs, based on discretization-time points over multi-state generations. About arc routing problems, the work presented in [ 1 ] shows a model to minimize the total cost of travel in a given graph for the single robot line coverage problem, focused on windy and asymmetric problems, both eulerian and not eulerian graphs and digraphs. Different
1.4 Project structure 5 types of genetic algorithms, its simulation and its real implementation in a drone are shown in [ 10 ]. Finally, the work [ 6 ] studies some drone arc routing problems as well as their application to infrastructure inspection. They propose an algorithm that iteratively solves the Rural Postman Problem, useful in vehicles with long endurance capacity. 1.4 Project structure The structure of this project is defined as follows: • Chapter 1. In this chapter, a definition of UAV is presented alongside its historical background. Then, the current situation of power line inspection is analyzed, not only with investigation works but also with new technology that is being used. • Chapter 2. As an introduction, the different parts involved in the inspection of electrical networks are described. In addition, the importance of obtaining an appropriate model of the problem is explained. • Chapter 3. The description of the Power Transmission Line Inspection Problem (PTLIP) is carried out with simplification hypotheses and an energy model. Also, it is presented two models of velocity that are considered in this project. Furthermore, the complexity of the problem is studied. Finally, a mathematical formulation of the problem is presented alongside the implementation with GUROBI. • Chapter 4. Three heuristics are given to solve the PTLIP: the Local Algorithm, the Semilocal Algorithm and the Global Algorithm. For the three heuristics, two types of criteria are presented to make decisions: the Down Wind criterion and the Minimum Energy criterion. To test the performance, the algorithms are used to inspect a given tree. • Chapter 5. In order to draw conclusions about the capabilities of the algorithms, a set of computational results for 1000 random trees are presented. A comparison between the heuristics as well as the optimal solution is given. • Chapter 6. The conclusions obtained in chapter 5 are presented and some ideas and improvements to continue the research started in this project are proposed.
2 Power Transmission Line Inspection Problem 2.1 First approach to the problem In this section a first approach to the Power Transmission Line Inspection Problem (PTLIP) is developed. This approach is going to be a description of all the factors and variables that come into play when carrying out inspections of electrical networks. The following description of the problem is based on the work presented in [15]. Figure 2.1 Development of the electricity transmission grid in Andalusia until 2014. In figure 2.1, is shown a map of the electricity transmission grid in Andalusia. As it can be seen, the grid has multiples power transmission lines of hundreds of kilometers in multiple directions, specifically, it contains 436 substations, 60524 km of power transmission lines and 59384 MVA of transformation power. If a review of the entire grid wants to be carried out manually, the costs of personnel, vehicles and time would be too high. Therefore, using UAVs for inspection is a better choice, but it also implies some logistical challenges. First of all, due to the grid‘s dimension, it is obviously that it would be needed more than one UAV. Next, not only multiples UAVs but also 7
8Chapter 2. Power Transmission Line Inspection Problem various recharge stations distributed throughout the grid are required, in order to refill the battery of the UAVs. To continue, it would be necessary some communication stations in order to transfer the information recollected by the UAVs or in case a pilot-on-ground wanted to take control of the vehicle if a problem occurs. Let I be the number of inspection UAVs (denoted as IUAVs) that are in charge of travel alongside the grid seeking defects on the power lines. The communication devices mounted on each IUAV have two main functions. First of all, it have to ensure a correct data transmission. Secondly, it have to prevent for two IUAV to collide. To prevent this to happen, a minimum safe distance must be adopted. Related to data transmission, when the data is firstly collected it is stored in IUAVs and later on, sent to other UAVs responsible for linking the IUAVs with the base of operations. Let R be the number of relay UAVs (denoted as RUAVs) that its mission is to serve as a data transistors for the IAUVs. It is interesting to have RUAVs instead of information transmission towers due to the displacement capacity that the RUAVs have, being able to approach the IUAVs at very short distances to ensure a better transmission of information. Also, both IUAVs and RUAVs can charge their batteries in different base stations (denoted as BSs) all over the grid, reducing the reducing the environmental impact of placing multiple transmission towers throughout the area. So, in order to inspect a power transmission line a set of IIUAVs, RRUAVs and SBSs are required. Due to the fact that in the same space of time the UAVs are traveling, collecting and sending data a energy model based on three types of consumption is required. First of all, as the UAV are moving, a energy consumption based on propulsion must be formulated. This energy depends on UAV velocity, air density, wind direction and geometrical characteristics of the UAV. Secondly, it is necessary to develop an energy consumption model based on collision avoidance. This energy model shelter three parts: distance estimation, warning collision and control actions. The consumption depends on the frequency an IUAV transmit its coordinates and velocity to the other UAVs and an error on the probabilities of the distance estimation between two IUAVs and failure in the warning system and control action system. To finish, it is necessary a consumption model for the data delivery. This model is based on the fact that not all the RUAVs are available for each IUAV at any time, so a priority list must be developed to connect an IUAV to an RUAV, in order to send the data to the central server. The connection is established based on the SNR (signal-to-noise ratio) between the IUAV and the RUAV, and also the RUAV and the central server. The IUAVs can rove the power lines with any velocity between a minimum and a maximum, but the velocity selected to inspect a transmission line is constant alongside the line. The altitude the IUAVs have to travel depends on a minimum safety altitude to prevent collisions with the towers and minimize EM interference and also to adapt to the geographical environment of the power lines that are being inspected. Also, every IUAV has to inspect zones of power lines so that not more than one IUAV roves the same power line. Every RUAV has to maintain a minimum distance with every IUAV and with RUAVs to avoid collisions and interference. Furthermore, the frequency of communications for each IAUV to a set of RUAVs must be appropriate to avoid saturation of the channel used for data delivery. Not only UAVs, but also the BSs must be positioned in a certain manner so that every UAV can select to which BS travel to recharge the battery, and furthermore, if the base stations are not well positioned, there will be points of the grid to inspect inaccessible for the IUAVs. More than one UAV can charge its battery in the same base station at the same time, but it has a maximum number of UAVs that can be recharged. In addition, the IUAV sensors can detect multiple failures and classify them by its severity. If the IUAV detects a failure with a more severity than one detected by other IUAV, the data delivery of the first one must have priority over the second one. Even the possibility for a pilot-on-ground to take control of each IUAV must exist if a problem occurs during the inspection. Finally, for each power line or any corridor that connects to power towers, there are cost, associated to energy, depending on which way the line is traveled or inspected. Indeed, a power line can be inspected but also use as a travel corridor from one point to
2.2 Problem modeling 9 another. Summarizing, a set of IUAVs have to collect data from power transmission lines and send the information into a central base through a link with RUAVs. The grid of power lines also shelter some base stations to charge the battery of the UAVs, battery that can be consumed by propulsion, collision avoidance and data delivery. With all these things presented, the idea now is to optimize the path of the IUAVs trying to minimize, for instance, the energy consumption, the traveling time or the number of tours the IUAVs have to make. 2.2 Problem modeling Each time a solution to a problem is given, what is really happening is that only a solution for one of the possibilities the problem can be modeled is found. In every scientific area, all the models are a simplification of the real world, and these simplifications are necessary in order to find a solution, otherwise the problems would be as difficult to tackle as the universe itself. The problem resolution process comprise two parts: generate a model of the problem and using that model to obtain a solution. If the model is a great approximation of reality, the solution obtained is more useful than a solution obtained for a bad-description model. Indeed, if the model has to many simplifications and no-sense assumptions, the solution to that model has no utilities at all. Between this two possibilities, it can be found two different approaches. The first one consist on having an approximation model ( Ma ) of a problem and give a precise solution to that model ( Sp ). On the other hand an approximate solution ( Sa ) can be found for a precise model of the reality ( Mp ). The first approach can be a good starting point in order to, once a solution is given, complicates the model and use the previous solution as an start point for the resolution of the next model. Nevertheless, finding a precise solution to a problem can last long time. The second approach is a good technique in case the problem is itself not complex, so a good modeling can be carried out and, for that model, iterates over the solution the solver can found. The problem with this approach is that the majority of the problems in physics, chemistry or any other scientific area are not easy to model in a precisely way, so it can take a long time to find a precise model for the concrete problem, time that could be employ into finding a solution. Problem −→ Ma−→ Sp(Ma) Problem −→ Mp−→ Sa(Mp) The description given of the PTLIP suggests that its resolution will not be easy. This problem involves multiple UAVs with different tasks, connection and data transfer as well as various energy models. That is why it is necessary to make some simplifying hypotheses in order to address the problem in a first approximation and, as solutions are obtained, increase the complexity of the models that are proposed. Dynamic scenarios Not only a good-described model must be used in order to have an useful solution of a problem, but also the model must consider time-changing conditions and how these changes affects the model. In real world, problem characteristics change before a model is given, while the solution to that model is being calculated and after the solution is obtained. In the PTLIP, one of the biggest challenges are whether conditions. It is known the whether behaves chaotically over time, which means that UAVs have to continually take parameters such as ambient pressure, wind direction or rainfall to adapt its route to environmental conditions that, in many cases, may be adverse. Furthermore, periods of sunlight must be consider in order to record images properly, but also to protect electronic
16 Chapter 3. Problem formulation model is used in which the value of ρ decreases with altitude. In this case, as the altitude is constant (100 m), the density remains constant. Table 3.2 UAV characteristics to solve the problem. W (N) 30 S (m2) 0,4332 Battery (mAh) 6000 Battery (V) 22,5 e 0,85 Λ3 CD00.018 In equation 3.1, the drag produced by the wind is one of the multiple effects that could rise the energy to travel an edge. To explain this graphically, just look at the three graphs in figure 3.2. It compares the energy of two edges, one of 50 m and the other of 52 m, for different CDw and for 2, 5 and 10 m/s. It can be clearly seen how the trend for the three speeds is linear, although as the speed increases, the slope of the line is greater. At the same time, it can be seen how as the speed increases, the 52 and 50 m edge lines approach each other, especially for very small resistance coefficients. Furthermore, the energy consumed to travel along the edge decrease with an increase on the speed, due to the weight of the c1 and c2 terms in 3.1. With this, it has wanted to make clear that, although the wind direction affects the energy that is consumed, the speed with which an edge is traveled and its length have a great importance in the SSFL equation.
3.3 Energy model and UAV characteristics 17 Figure 3.2 Energy to travel an edge depending on CDwfor different speeds.
18 Chapter 3. Problem formulation 3.4 Formulation In order to analyze the results of the heuristics, it is presented now a mathematical formulation for the PTLIP based on Capacitated Arc Routing Problems (CARPs) like the one in [ 2 ]. This formulation is made in Integer Lineal Programming (ILP) to be able of computing in the solver GUROBI [ 11 ]. Let TOW be the set of N power transmission towers, Let V the set of nodes in the grid, V={n} ∪ TOW , where n is the node of the base. Let G(V,E) be the graph with V vertices and E edges, so E={(i,j)∈V2:i=j} . The tree to inspect T(TOW,Er) is part of the graph G , with Er the set of required edges, that is, the edges the UAV has to inspect, Er⊂E and Er={(i,j)∈N:i=j} . Let ti j , esi j and edi j be the time, energy of service (i,j)∈Er and energy of traveling through (i,j)∈E . Let B be the capacity of the UAV. Let M be the number of tours. Let xi j be a binary variable that is 1 if edge i j in Er is serviced and 0 if it is not. Let yi j be a integer variable that represents the number of times the UAV travels along the edge i j in deadheading. min{ M ∑ k=1 ∑ i,j∈Er xk i,jti,j+ M ∑ k=1 ∑ i,j∈E yk i,jti,j+(M−1)tc}(3.4) Subject to ∑ i∈N yk ni =1k={1,2,...M}(3.5) ∑ j:(i,j)∈Er xk i j +∑ j:(i,j)∈E yk i j =∑ j:(i,j)∈Er xk ji +∑ j:(i,j)∈E yk ji ∀i∈V k ={1,2,...M}(3.6) M ∑ k=1 (xk i j +xk ji) = 1∀i,j∈Er(3.7) ∑ i,j∈E (esi j xk i j +edi j yk i j)<=B k ={1,2,...M}(3.8) In 1v model, the energy to travel through an edge is the same as servicing it, due to the fact that it is only one speed. Constrain 3.5 assure that in every tour there exists a link between the ground station and a node. To fix the continuity of the paths followed by the UAV it is written the constrain 3.6 and constrain 3.7 ensure that every required edge is serviced exactly one. The constrain 3.8 assure the energy requirements. In this mathematical formulation the number of tours must be known. That is why in 3.1 an iterative algorithm is presented by which the solution of the GUROBI solver can be obtained through the gurobipy python module. Through this code, it is tried to find the value of M that makes the problem solvable (model is feasible). There is the possibility that, for a value of M , gurobipy will not be able to find the solution and therefore the problem would be unsolvable (model is infeasible) for that determined number of tours. With this algorithm, it would be automatically searched the solution with M+1 tours or higher values. If for M=N−1 the problem cannot be solved, it means that the UAV cannot service the tree with the full battery available. This algorithm is going to be called OPT, as it tries to find the optimal solution.
3.4 Formulation 19 Code 3.1 OPT Algorithm. Input: mathematical formulation Output: GUROBI solution set M = 1 while M < N-1: get GUROBI solution if model == feasible: return GUROBI solution break else: M = M+1
4 Heuristics 4.1 Decision criteria Before introducing the heuristics, we establish some decision criteria (CRIT) that can be considered when the UAV arrives at a node or is at the base and there are several edges to be serviced. In this project two criteria will be followed: • Down Wind criterion. Using this criteria, if the UAV is at the base, it has to select the node that is most against the wind direction. If the UAV is on a node, it travels to the node that is most aligned with the wind direction. • Minimum Energy criterion. Using this criteria, if the UAV is at the base, it has to select the node that minimize the traveling energy. If the UAV is on a node, it travels to the node that using the edge that minimizes the traveling energy. 4.2 H1: Local Algorithm The decisions at the nodes of the local algorithm is as follows: • When the UAV is at the base, it has to travel to nodes that have required edges without servicing using CRIT. • When the UAV is at a node with required edges without servicing, it has to select the node to travel using CRIT. If there is no enough battery to serve that edge, it has to return to base and come back to that node. •When the UAV is at a node without required edges to give service, it has to return to base. The motivation behind the creation of this heuristic is to obtain quick solutions due to its local behavior. This algorithm is very useful when you want to obtain immediate solutions or as initial estimates to find the optimal solution. 21
22 Chapter 4. Heuristics Code 4.1 Local Algorithm. Input: T(V,E), UAV characteristics, p, energy model Output: Tc, M while e in E not serviced: if UAV in BASE: travel to node using CRIT update Tc, M, battery if UAV in node: if neighbours(node) not serviced: if enough battery travel to node using CRIT else: return BASE and come back to node else: return BASE 4.3 H2: Semilocal Algorithm The decisions at the nodes of the semilocal algorithm is as follows: • When the UAV is at the base, it has to travel to nodes that have required edges without servicing using CRIT. • When the UAV is at a node with required edges without servicing, it has to select the node to travel using CRIT. If there is no enough battery to serve that edge, it has to return to base and come back to that node. • When the UAV is at a node without required edges to give service, it has to look up a node to travel in the neighbors of the node in the graph using CRIT and select the first option. If the battery is enough to inspect at least another edge, it has to travel, but if it is not, the UAV has to return to base. This algorithm incorporates a small improvement over the local algorithm. Instead of returning to base when he reaches a node with no edges to visit, the opportunity to search for a node to travel to in deadheading is given. The edge is chosen using CRIT. If there is not enough battery to inspect at least another edge when that selection is made, the UAV returns to base.
4.4 H3: Global Algorithm 23 Code 4.2 Semilocal Algorithm. Input: T(V,E), UAV characteristics, p, energy model Output: Tc, M while e in E not serviced: if UAV in BASE: travel to node using CRIT if UAV in node: if neighbours(node) not serviced: if enough battery travel to node using CRIT else: return BASE and come back to node else: look up node with required edges not serviced if enough battery: travel to node using CRIT else: return BASE 4.4 H3: Global Algorithm The decisions at the nodes of the global algorithm is as follows: • When the UAV is at the base, it has to travel to nodes that have required edges without servicing using CRIT. • When the UAV is at a node with required edges without servicing, it has to select the node to travel using CRIT. If there is no enough battery to serve that edge, it has to return to base and come back to that node. • When the UAV is at a node without required edges to give service, it has to look up a node to travel in the neighbors of the node in the graph using CRIT and select the first option. If the battery is enough to inspect at least another edge, it has to travel. If it is not, it has to keep looking if there is another node to travel in deadheading and then inspect an edge with the battery available using CRIT. If there is none, the UAV has to return to base. The Global Algorithm is presented as the improvement of the two previous algorithms. When it reaches a node with no edges to visit, it has to search for edges to deadhead through using CRIT. However, this time the algorithm is allowed to explore more than one edge in case there is no battery available in the first selection. With this, the UAV can, for instance, explore edges without much energy requirement that are close to the base.
24 Chapter 4. Heuristics Code 4.3 Global Algorithm. Input: T(V,E), UAV characteristics, p, energy model Output: Tc, M while e in E not serviced: if UAV in BASE: travel to node using CRIT if UAV in node: if neighbours(node) not serviced: if enough battery travel to node using CRIT else: return BASE and come back to node else: look up all edges to travel in deadheading using CRIT if enough battery: travel to node using CRIT else: return BASE
4.5 Comparison with a particular tree 25 4.5 Comparison with a particular tree In this section we compare the performance of the heuristics using a tree generated at random. For this we used a code existing in the graph module NetworkX [12] in Python 4.4. The tree is shown in figure 4.1. This example has 15 nodes, that represent the power transmission towers, and 14 required edges, representing the power lines. The UAV base is located in point (0,0) (shown as a blue square), the wind direction is 1, the minimum velocity is 2 m/s and the maximum velocity is 10 m/s. We consider p=0and p=0.5as the probabilities of wind change. Code 4.4 Aleatory Tree Algorithm. Input: N Output: T, Lr, pos import networkx as nx T = nx.random_tree(N, seed = None, create_using = nx.MultiGraph) E = set() for i in nx.nodes(T): for j in nx.nodes(T): edge1 = (i,j) edge2 = (j,i) if edge1 in nx.edges(T) and edge2 not in nx.edges(T): E.add(edge2) T.add_edges_from(E) Lr = nx.adjacency_matrix(T) pos = nx.spring_layout(T, scale=1000) As the heuristics can be implemented using different criteria, for this example we consider the Down Wind CRIT to H2 (4.2) and the Minimum Energy CRIT to both H1 (4.1) and H3 (4.3). It is obvious that, using the same CRIT for the heuristics, we have: CT(H3)≤CT (H2)≤CT(H1)(4.1) NR(H3)≤NR(H2)≤NR(H1)(4.2) where CT stands for the coverage time and NR stands for the number of refueling.
32 Chapter 4. Heuristics 4.5.3 Model 2V with p = 0 The results for the different heuristics using a service speed of 2 m/s and a deadheading speed of 10 m/s are presented below. Coverage_time_H1 = 90.39 Coverage_time_H2 = 50.815 Coverage_time_H3 = 47.59 Path_H1: [11, 8, 10, 4, 7] [11, 2, 1, 13] [3, 5, 14, 9] [3, 2] [4, 0] [4, 12] [14, 6] Path_H2: [13, 1, 2, 3, 5, 14, 9, 6, 14, 0, 4, 7] [12, 4, 10, 8, 11, 2] Path_H3: [6, 14, 9, 14, 5, 3, 2, 1, 13, 2, 11, 8, 10, 4, 7] [4, 0, 4, 12] Computation_time_H1 = 0.06236 Computation_time_H2 = 0.07813219 Computation_time_H3 = 0.096877 The main difference with respect to the 1V model is that both the coverage time and the number of tours are reduced, or at least maintained, as in the case of the H1(ME). In the case of the number of tours, using high speeds in the SSFL model reduces the energy needed to traverse an edge until 10 m/s approximately, as can be seen in figure 4.5. This makes it possible for more unneeded edges to be used before the UAV runs out of battery. In turn, the coverage time is reduced as the speed is increased when making deadheading travels Figure 4.5 Energy to travel an edge with different CDw.
4.5 Comparison with a particular tree 33 The tours done by the UAV following the three heuristics are illustrated in figure 4.6. Figure 4.6 Path followed by the UAV with H1(ME), H2(DW) and H3(ME) for 2V and p = 0. For the 2V model, the optimal results alongside the optimal path followed by the UAV (figure 4.7) are the following: Coverage_time = 37.7121 Path [7, 4, 0, 12, 4, 10, 8, 11, 2, 1, 13, 2, 3, 5, 14, 6, 9, 14] Computation_time = 0.15625
34 Chapter 4. Heuristics Notice that the UAV can cover the tree only using one tour due to the reduction of energy in deadheading using 10 m/s. This reduction of energy is possible due to the selection of velocities, the wind direction and the weight of c1 and c2 in SSFL formula (3.1). Figure 4.7 Path followed by the UAV in optimal solution for 2V and p = 0.
4.5 Comparison with a particular tree 35 4.5.4 Model 2V with p = 0.5 To finish this section, we present the solutions obtained considering a probability of change in wind direction of 50 % in the 2V model (figure 4.8). The results shown a better coverage time for the three algorithms and for the H1(ME) a less number of tours compared to p = 0. At higher speeds, the CDW has a higher impact in the selection of nodes for the H1(ME) and H3(ME), not starting at the same points in p = 0 and p = 0.5. Coverage_time_H1 = 83.279166 Coverage_time_H2 = 48.86278 Coverage_time_H3 = 47.62 Path_H1: [3, 5, 14, 9] [3, 2, 1, 13] [11, 8, 10, 4, 7] [6, 14] [0, 4, 12] [11, 2] Path_H2: [13, 1, 2, 11, 8, 10, 4, 12] [0, 4, 7, 2, 3, 5, 14, 6, 9, 14] Path_H3: [7, 4, 0, 4, 12, 4, 10, 8, 11, 2, 1, 13, 2, 3, 5, 14, 9] [14, 6] Computation_time_H1 = 0.078056 Computation_time_H2 = 0.07844 Computation_time_H3 = 0.18753
36 Chapter 4. Heuristics Figure 4.8 Path followed by the UAV with H1(ME), H2(DW) and H3(ME) for 2V and p = 0.5.
5 Experiments 5.1 Introduction The results of the previous chapter have been obtained for a generic tree of 15 nodes in order to show the performance of the heuristics. However, if the objective is to know which heuristic offers the best solutions, they should be tested for a large number of trees. In this section, a statistical study will be carried out testing the H1, H2 and H3 as well as the optimal solution for 1000 random trees. The number of nodes of the trees varies between 10 and 100 and the probability of wind change (p) varies between 0 and 1. An UAV will be simulated with the characteristics of table 3.2 with a minimum speed of 2 m/s and a maximum of 10 m/s. It has to be said that the optimal solution is obtained only for p=0 , due to the fact that it is needed to work with deterministic models. The computer used to obtain the results is a HP LAPTOP-99RQA79B with an 11th Gen Intel Core i5-1135G7 2.4GHz and a 16GB RAM. To begin with, it is important to sum up some ideas drawn from the previous chapter: • It seems that, for minimizing coverage time, H3(ME) shows better results than H1(ME) and H2(DW). Also, for minimizing the number of tours, H3(ME) has slightly better results than H2(DW) and H1(ME). • For small speeds, wind direction is determinant only for the Down Wind criterion. At high speeds, the both ME and DW are hardly affected by wind direction and its changes. • The optimal solution for 1V and 2V model minimize both coverage time and the number of tours, so it seems that the solution that minimizes the coverage time is correlated with the solution that minimizes the number of tours. • Regarding to computational time, the H2(DW) and H3(ME) require more time to obtain the solution than the H1(ME). Obviously, the optimal solution will require more time to be calculated that the ones obtained with the heuristics. • In 1V model, the results between the heuristics are much different between them than in 2V model. In the next section, relevance information extracted from the UAV path will be exposed for the generated random trees. This information will be presented in the form of means (5.1), medians (5.2) and population deviations (5.3). ¯x= N ∑ i=1 xi N(5.1) ˜x=xN/2(5.2) 37
38 Chapter 5. Experiments sp(x) = s1 N−1 N ∑ i=1 (xi−¯x)(5.3) Also, we introduce the concept of relative deviation of x to y as it appears in 5.4. In the next section it is specified whether the numerator is in absolute value or not Ratio =x−y y(5.4) The coverage time is going to be referenced as Tc , the number of tours as M and the computational time as Te . The number of not required edges the UAV uses is going to be referenced as NRE , and the number of required edges traveled in deadheading as RDE . As done in the previous chapter, the experiments are going to be carried out with both H1 and H3 using CRIT = ME and H2 using CRIT = DW. This can be referenced as H1(ME), H2(DW) and H3(ME) or as H1-ME, H2-DW and H3-ME.
5.2 Statistics 39 5.2 Statistics 5.2.1 1V model First of all, it is going to be studied the type of distribution followed by the coverage time, the number of tours and the computation time for the heuristics and the optimal solution. In figures 5.1, 5.2 and 5.3 the coefficient of variation is presented. This coefficient, denoted as cv , is the ratio between the sample deviation and the mean of a distribution, in this case, the coverage time. The data is presented for different N and p. It can be seen how this ratio is under 25% for H1(ME), H2(DW) and H3(ME), so the data collected is relatively close to the mean coverage time. Also, the H2(DW) is the algorithm with more variations in cvas p increases. Figure 5.1 Ratio between sp(tc) and ¯ tc for 10 to 50 nodes (left) and 60 to 100 nodes (right) in H1(ME) 1V. Figure 5.2 Ratio between sp(tc) and ¯ tc for 10 to 50 nodes (left) and 60 to 100 nodes (right) in H2(DW) 1V.
40 Chapter 5. Experiments Figure 5.3 Ratio between sp(tc) and ¯ tc for 10 to 50 nodes (left) and 60 to 100 nodes (right) in H3(ME) 1V. Furthermore, in figures 5.4, 5.5 and 5.6 the relative deviation of median coverage time to mean coverage time is presented. This relative deviation is presented in absolute value. It can be seen how this deviation is greater for N = 10 (still under 8 % ) but it quickly decreases as N increases. As the number of nodes increases, the distribution of coverage time for all p trends to be symmetrical. For H1(ME), this relative deviation seems to be more stable with the increase of p than the other algorithms. All together, the mean coverage time can be used to draw conclusions about the results obtained. Another thing observed is that there is no a correlation between an increase in p and a more symmetrical distribution. Figure 5.4 Relative deviation of median coverage time to mean coverage time H1(ME) 1V. Figure 5.5 Relative deviation of median coverage time to mean coverage time H2(DW) 1V.
5.2 Statistics 41 Figure 5.6 Relative deviation of median coverage time to mean coverage time H3(ME) 1V. In figures 5.7, 5.8, 5.9 and 5.10, 5.11, 5.12 the same analysis for the distribution of number of tours is carried out, obtaining the same conclusions: the distribution is close to the mean and as N increases, the distribution is more symmetrical. Nevertheless, it is appreciated how the distributions are more asymmetrical than the ones for coverage time. Figure 5.7 Ratio between sp(M) and ¯ M for 10 to 50 nodes (left) and 60 to 100 nodes (right) in H1(ME) 1V. Figure 5.8 Ratio between sp(M) and ¯ M for 10 to 50 nodes (left) and 60 to 100 nodes (right) in H2(DW) 1V.
48 Chapter 5. Experiments Secondly, it is going to be studied how the mean coverage time, number of tours and computational time changes as p increases for the three heuristics. Furthermore, it is studied how the increase of N affects the mean. It is presented below the variations of mean coverage time when p changes for different number of nodes (figures 5.25, 5.26, 5.27). Figure 5.25 Mean coverage time H1(ME) for 10 to 50 nodes (left) and 60 to 100 nodes (right) 1V. Figure 5.26 Mean coverage time H2(DW) for 10 to 50 nodes (left) and 60 to 100 nodes (right) 1V. Figure 5.27 Mean coverage time H3(ME) for 10 to 50 nodes (left) and 60 to 100 nodes (right) 1V. As can be seen, the H2(DW) presents more variations in mean coverage time when the wind direction changes. Both H1(ME) and H2(ME) present almost no variations in coverage time when the wind direction changes. The same data for number of tours is presented in figures 5.28, 5.29 and 5.30. The same conclusions can be drawn: H2(DW) is more sensitive to changes on p than H1(ME) and H3(ME).
5.2 Statistics 49 Figure 5.28 Mean number of tours H1(ME) for 10 to 50 nodes (left) and 60 to 100 nodes (right) 1V. Figure 5.29 Mean number of tours H2(DW) for 10 to 50 nodes (left) and 60 to 100 nodes (right) 1V. Figure 5.30 Mean number of tours H3(ME) for 10 to 50 nodes (left) and 60 to 100 nodes (right) 1V. To finish this part, the variations of mean computational are presented in figures 5.31, 5.32 and 5.33. It can be seen how the grater the p, the higher the time the computer needs to find the solution. This can be explained by how the wind is changed in the algorithms. This change is made through a function that changes the wind direction to one of the neighbors wind directions when UAV reaches a node or the base. As p increases, this function will be called more times, increasing the computational time. Furthermore, the increase in Te is different depending on the heuristic. For H1(ME) and H2(DW), there is a linear trend between the increase in p and the increase in Te for the same N. On the other hand, for H3(ME), the increase of Teslows down as p increases.
50 Chapter 5. Experiments Figure 5.31 Mean computational time H1(ME) for 10 to 50 nodes (left) and 60 to 100 nodes (right) 1V. Figure 5.32 Mean computational time H2(DW) for 10 to 50 nodes (left) and 60 to 100 nodes (right) 1V. Figure 5.33 Mean computational time H3(ME) for 10 to 50 nodes (left) and 60 to 100 nodes (right) 1V.
5.2 Statistics 51 Thirdly, it is going to be compared the results of the three heuristics and the optimal solution between them. To do this, it has to be used the data obtained for p = 0. As previously commented, for the coverage time and number of tours we use the mean and for computational time we use the median. In table 5.1 and figure 5.34, the variations of mean coverage time with the increase in N are presented. It can be seen in figure 5.34 how the coverage time follows a linear trend for the the heuristics and the optimal solution. In terms of coverage time, H3(ME) is the best heuristics followed by H2(DW) and H1(ME). The trend fits very well the experimental data, so it is possible to calculate the coverage time for the heuristics in short time and compare the result with the optimal solution without having to calculate it, it can be simply looked at the value it would have following the linear trend. In figure 5.35, the mean relative deviation between the heuristic and the optimal solution for Table 5.1 Mean coverage time for the heuristics and the optimal solution 1V. Mean coverage time (minutes) N H2(DW) H3(ME) H1(ME) OPT 10 72.15 55.84 76.83 55.8 20 128.69 80.44 132.42 70.357 30 182.73 104.83 187.69 82.41 40 236.49 127.83 242.03 97.36 50 281.5 144.41 288.43 108.83 60 332.35 164.29 339.26 121.61 70 382.13 183.65 389.96 132.15 80 429.89 201.83 439.15 142.85 90 471.02 215.8 483.19 155.01 100 517.34 233.84 532.03 165.67 Figure 5.34 Mean coverage time for the heuristics and the optimal solution 1V.
52 Chapter 5. Experiments coverage time is presented. Again, it can be seen that with H3(ME) are obtained better results than with the other heuristics. Figure 5.35 Mean relative deviation of the heuristics to optimal solution for coverage time 1V.
5.2 Statistics 53 For the number of tours, the same conclusions can be extracted in figures 5.36, 5.37 and in table 5.2. H3(ME) offers better results on minimizing the number of tours and the trend is linear for the three heuristics and the optimal solution. But there is more to be talked about. Table 5.2 Mean number of tours for the heuristics and the optimal solution 1V. Mean number of tours N H2(DW) H3(ME) H1(ME) OPT 10 2.718 1.792 4.3 1.786 20 4.928 2.53 7.97 2.006 30 7.056 3.23 11.564 2.114 40 9.208 3.882 15.288 2.466 50 11.126 4.338 19.038 2.772 60 13.024 4.874 22.532 3.002 70 15.082 5.436 26.356 3.136 80 17.002 5.934 30.402 3.344 90 18.756 6.262 34.008 3.654 100 20.576 6.79 37.62 3.866 Figure 5.36 Mean number of tours for the heuristics and the optimal solution 1V.
54 Chapter 5. Experiments Figure 5.37 Mean relative deviation of the heuristics to optimal solution for number of tours 1V. Using H1(ME), the UAV has to do much more tours than in H2(DW), but the results obtained in coverage time shown that there is a small difference between H1(ME) and H2(DW). In figure 5.38, the relative deviation of CCA to MDA for ¯ Tc and ¯ M is presented. It can be seen that the deviation of coverage time is less than the deviation of the number of tours. This leads to the idea that, the reduction of one of them may not lead to the reduction of the other, although there is a connection between this two variables. For instance, with H1(ME) the UAV has to travel more to base because of its local behavior, not because a need in recharge the battery. So, when the UAV arrives to the base in H1(ME), it spends less time than in H2(DW). Indeed, looking at the data obtained from the trees, it can be notice how for some aleatory trees, the UAV has to do more tours using H1(ME) than using H2(DW) but the coverage time is minor in H1(ME). In table 5.3, the number of trees where this happens for every N for p = 0 is presented. Also, as wind changes affects in a high way the Down Wind criterion, the proportion of tress where this happens for three different probabilities of change in wind direction is presented in figure 5.38. The decrease in the proportion of the number of trees is expected because the coverage time decreases as p increases for the H2(DW). Figure 5.38 Relative deviation of mean coverage time and mean number of tours of H2(DW) to H1(ME) 1V.
5.2 Statistics 55 Table 5.3 Proportion of trees where UAV spends less time inspecting in H1(ME) than in H2(DW) with more number of tours. Mean coverage time and mean number of tours for those trees for p = 0 1V. H1(ME) H2(DW) N Ratio (%) Tc(minutes) M Tc(minutes) M 10 19.4 82.18 4.48 86.3 3.25 20 39.8 123 7.52 133.4 5.32 30 40.2 177.6 11.21 191.99 7.58 40 43.4 221.39 14.49 238.56 9.52 50 40.8 274.28 18.64 292.98 11.73 60 40.2 313.63 21.65 335.47 13.44 70 39.8 366.11 25.62 390.89 15.69 80 38.8 413.04 29.37 439.68 17.76 90 41.8 454.01 33 482.75 19.51 100 38 498.6 36.81 528.51 21.55 Figure 5.39 Proportion of trees where Tc in H1(ME) is less than Tc in H2(DW) but the number of tours is greater for p = 0, p = 0.4 and p = 0.7 1V.
56 Chapter 5. Experiments Table 5.4 and figure 5.40 show the results of the median computational time. It can be seen how the heuristics H1 and H2 has a Te of O(N) , H3 has a computational time of O(N2) and the optimal solution has a computational time of O(exp(N)) . The order of the H1 an H2 are expected because this algorithms have a more local behavior than H3. Table 5.4 Median computational time for the heuristics and the optimal solution 1V. Median computational time (s) N H3(ME) H2(DW) H1(ME) OPT 10 0.005982 0.00459 0.00448 0.07564 20 0.03125 0.05062 0.013623 0.21865 30 0.07825 0.05946 0.018638 0.48517 40 0.17298 0.09075 0.03124 1.7347 50 0.33005 0.11593 0.04686 4.56 60 0.5349 0.19622 0.06248 8.54 70 0.8498 0.27304 0.06924 16.27 80 1.2751 0.30947 0.09872 42.021 90 1.80539 0.33951 0.11362 107.615 100 2.40819 0.4381 0.1425 208.16 Figure 5.40 Median computational time for the heuristics (top) and for the optimal solution (bottom) 1V.
5.2 Statistics 57 However, at the beginning of this section was presented the variations for the mean computational time when p is modified, so the trend of the heuristics must be inspected for some values of p. As have been explained in the previous point, the relative deviation between the median and mean of computational time is not significant, so to look at the variations of the trend both mean and median can be used. Figures 5.41, 5.42 and 5.43 shown the trend of the median computational time for p = 0.2, p = 0.5 and p = 1. It can be seen how the linear trend of the H1 and H2 are lost but are still smaller than H3 in every N. H3 is still following a O(N2) trend. Looking a the value of R2 , as p increases, it tends to increase a little or maintain for the polynomial trend of H3 and decrease in the linear trend of both H2 and H1. In statistics, R-squared ( R2 ) or the coefficient of determination is the proportion of the variation in the dependent variable (in this case ˜ Te ) that is predictable from the independent variable (in this case N ). If R2 approaches to one, it can be predicted with minor errors how the model would behave for the data, whereas if R2 approaches to 0, the dependent variable has small correlation to the independent variable. Figure 5.41 Median computational time for p = 0.2 1V. Figure 5.42 Median computational time for p = 0.5 1V.
64 Chapter 5. Experiments Figure 5.56 Distribution of coverage time for N = 70 and N = 80 H3(ME) 2V. Figure 5.57 Distribution of coverage time for N = 90 and N = 100 H3(ME) 2V.
5.2 Statistics 65 Figure 5.58 Distribution of coverage time for N = 10 and N = 20 H1(ME) 2V. Figure 5.59 Distribution of coverage time for N = 30 and N = 40 H1(ME) 2V. Figure 5.60 Distribution of coverage time for N = 50 and N = 60 H1(ME) 2V.
66 Chapter 5. Experiments Figure 5.61 Distribution of coverage time for N = 70 and N = 80 H1(ME) 2V. Figure 5.62 Distribution of coverage time for N = 90 and N = 100 H1(ME) 2V.
5.2 Statistics 67 In figures 5.63, 5.64, 5.65 and 5.66 the distribution of the number of tours for different values of p and different number of nodes is shown. It is not going to be presented all N possibilities because the objective is to show how the distribution behaves as the number of nodes increase. The nodes selected to plot are N = 20, N = 40, N = 60 and N = 90. It can be clearly appreciated how H2(DW) is more sensitive to wind changes than H1(ME) and H3(ME), being this last one the most stable of the heuristics. Also, for H2(DW) and H3(ME) are several non-negligible peaks associated with the different values of number of tours. There are peaks higher than the others, denoting the most likely M it will take to inspect a random tree. Nevertheless, in H1(ME) it is only one peak for each p but much wider. That means that there are more values of M dispersed than in H2(DW) and H3(ME). Figure 5.63 Distribution of number of tours for N = 20 2V.
68 Chapter 5. Experiments Figure 5.64 Distribution of number of tours for N = 40 2V. Figure 5.65 Distribution of number of tours for N = 60 2V.
5.2 Statistics 69 Figure 5.66 Distribution of number of tours for N = 90 2V.
70 Chapter 5. Experiments To finish with the analysis of the algorithm data distribution, it is going to be looked at the computational time. The distribution of Te is presented in figures 5.67, 5.68, 5.69, 5.70 and 5.71 for N = 10, 30, 50, 70 and 100. The results shown below are quite interesting. As mentioned in 1V model, an increase on p implies an increase on computational time, due to the fact that the function that changes the wind direction is called more times. This effect is superior as N increases. For the three heuristics, as N increases the differences between the distribution of p = 0 and p = 1 are greater. It is also noticed how the distribution is more symmetrical as N becomes greater, being N = 10 the distribution more asymmetrical for the three heuristics, as in 1V model. There are also difference between the heuristic distributions. In H3(ME) the distributions from p = 0 to p = 0.9 are very close each other, but for p = 1 there is an important gap, being this effect greater as N increases. The distributions of Te in H2(DW) show how, as p increases, the peak of the distribution moves progressively to the right (mean Te becomes greater). Furthermore, in H2(DW), it can be seen how as N increases, the distributions tend to have only one peak, meaning there are less deviation of data. In H1(ME) happens the same as in H2(DW), but it can be appreciated how the peaks are smaller as p increases. Also, in H1(ME) for p = 0, there are several non-negligible peaks in the distribution until N reaches 90. In general terms, there is a significant difference in time to obtain the solution as p increases. Furthermore, the distributions are quite symmetrical and do not usually have a high deviation to the mean. Nevertheless, for small N, much more dispersion can be observed. Figure 5.67 Distribution of computational time for N = 10 2V.
5.2 Statistics 71 Figure 5.68 Distribution of computational time for N = 30 2V. Figure 5.69 Distribution of computational time for N = 50 2V.
72 Chapter 5. Experiments Figure 5.70 Distribution of computational time for N = 70 2V. Figure 5.71 Distribution of computational time for N = 100 2V.
5.2 Statistics 73 Before choosing an appropriated statistic to drawn conclusions about the algorithms performances, it has to be looked at the distributions of the optimal solution (OPT) for Tc , M and Te . In figures 5.72, 5.73, 5.74, 5.75 and 5.76 the distribution of coverage time for different number of nodes can be seen. The coverage time follows a symmetrical distribution for almost every N; the bigger difference between mean and median is in N = 60. Figure 5.72 Distribution of coverage time for N = 10 and N = 20 OPT 2V. Figure 5.73 Distribution of coverage time for N = 30 and N = 40 OPT 2V. Figure 5.74 Distribution of coverage time for N = 50 and N = 60 OPT 2V.
80 Chapter 5. Experiments Figure 5.86 Distribution of coverage time for N = 50 and N = 60 2V. Figure 5.87 Distribution of coverage time for N = 70 and N = 80 2V. Figure 5.88 Distribution of coverage time for N = 90 and N = 100 2V. In table 5.8 it is presented the relative deviation of the mean coverage time of H3(ME) to OPT. The coverage time values for the H3(ME) and the OPT that have the biggest and the smallest difference between them are in the table 5.9. In that table, it can be seen how the minimum difference of inspecting a tree using H3(ME) and OPT is less than a minute for each N. Furthermore, if the optimal coverage time is divided by the coverage time of H3(ME), the maximum value is between 1.21%.
5.2 Statistics 81 Table 5.8 Relative deviation of mean coverage time of H3(ME) to OPT 2V. Ratio = ¯ Tc(H3)−¯ Tc(OPT ) ¯ Tc(OPT ) N10 20 30 40 50 60 70 80 90 100 Ratio (%) 0,2 0,523 1,1 6 9,68 11,61 11,1 9,7 9,4 10,3 Table 5.9 Values of coverage time with minimum and maximum difference between H3(ME) and OPT 2V. Minimum difference (minutes) Maximum difference (minutes) N H3(ME) OPT H3(ME) OPT 10 52.328 52.291 40.056 37.651 20 53.265 53.258 46.501 40.559 30 51.862 51.859 82.379 69.851 40 74.949 74.945 109.72 91.976 50 67.1013 67.1011 114.82 92.95 60 76.826 76.561 115.39 95.04 70 106.922 106.832 112.829 94.02 80 107.78 107.5 115.37 93.29 90 116.304 116.14 117.55 96.315 100 112.35 111.99 161.09 133.55 Finishing with Tc , in figure 5.89 the trend of mean coverage time of the heuristics and OPT is shown in figure 5.89. It can be appreciated how, as in 1V model, the trend is linear. Also, the difference between the mean coverage time of H1(ME) and H2(DW) has increased. Figure 5.89 Mean coverage time for the heuristics and the optimal solution 2V.
82 Chapter 5. Experiments For the number of tours, in figures 5.90, 5.91, 5.92, 5.93 and 5.94 a comparison between the distributions of the three heuristics and OPT is presented. As mentioned before, H1(ME) has more spread dispersion of data than H2(DW) and H3(ME). As can be seen, both the H3(ME) and OPT have distributions with values that are very similar to each other, even for high values of N. The highest values for the number of tours are in H1(ME) followed by H2(DW). Obviously, the OPT minimize the number of tours that the UAV needs to do to inspect the tree. Figure 5.90 Distribution of number of tours of the heuristics and the optimal solution for N = 10 and N = 20 2V. Figure 5.91 Distribution of number of tours of the heuristics and the optimal solution for N = 30 and N = 40 2V.
5.2 Statistics 83 Figure 5.92 Distribution of number of tours of the heuristics and the optimal solution for N = 50 and N = 60 2V. Figure 5.93 Distribution of number of tours of the heuristics and the optimal solution for N = 70 and N = 80 2V. Figure 5.94 Distribution of number of tours of the heuristics and the optimal solution for N = 90 and N = 100 2V. Furthermore, in figure 5.95 it is shown the trend of mean number of tours. The trend is linear for the heuristics and the optimal solution
84 Chapter 5. Experiments Figure 5.95 Mean number of tours for the heuristics and the optimal solution 2V. Finally, in figures 5.96, 5.97, 5.98, 5.99, and 5.100 the comparison between the relative deviation of the heuristics to OPT of Te is presented. This is done due to the big difference between the computational time of the heuristics and the computational time of the optimal solution. It can be clearly seen how the deviation is very close to −100% for all the heuristics, but the data dispersion is bigger in H3(ME). This happens due to the fact that is the heuristic with the biggest computational time. Also, it can be seen how H1(ME) has biggest (in absolute value) deviation, meaning that is the heuristic with minor computational time, followed pretty close by the H2(DW). It should be mentioned that even for small values of N, there are no trees in which the heuristic solutions take longer time to obtain than the optimal solution. Figure 5.96 Relative deviation of heuristics to optimal solution of Tefor N = 10 and N = 20 2V.
5.2 Statistics 85 Figure 5.97 Relative deviation of heuristics to optimal solution of Tefor N = 30 and N = 40 2V. Figure 5.98 Relative deviation of heuristics to optimal solution of Tefor N = 50 and N = 60 2V. Figure 5.99 Relative deviation of heuristics to optimal solution of Tefor N = 70 and N = 80 2V.
86 Chapter 5. Experiments Figure 5.100 Relative deviation of heuristics to optimal solution of Tefor N = 90 and N = 100 2V.
5.2 Statistics 87 The trends of the median computational time of the heuristics and the optimal solution are in figure 5.101. We can extract the same conclusions than in 1V model: Both H1 and H2 follow a linear trend, H3 follows a second-order polynomial trend and OPT follows an exponential trend. Figure 5.101 Median coverage time for H1(ME), H2(DW), H3(ME) and OPT 2V.
88 Chapter 5. Experiments As done with 1V model, we must look at the trend of median computational time for different values of p. In figures 5.102, 5.103, 5.104 and 5.105 it is presented the polynomial trend for H3 and the linear trend for both H1 and H2 for p = 0.1, p = 0.4, p = 0.7 and p = 1. As mentioned in 1V model, the linear trends of the H1 and H2 are lost with the increase of p. Again, although the linear trend is lost, ˜ Tevalues of H1 are under H2 and H3. It is worth to mention that, although the linearity is lost, R2values are still above 90%. Figure 5.102 Median computational time for the heuristics p = 0.1 2V. Figure 5.103 Median computational time for the heuristics p = 0.4 2V.
5.2 Statistics 89 Figure 5.104 Median computational time for the heuristics p = 0.7 2V. Figure 5.105 Median computational time for the heuristics p = 1 2V.
96 Chapter 6. Conclusions and future work 6.2 Future work To finish this project, a guide about how to continue this work is presented below: • Population. The statistics presented in chapter 5 have been obtained with a population of 1000 random trees. In order to get better results, the same statistics can be extracted from a larger population and to test if the results are affected increasing the sample size. • Energy model. The energy model presented in 3.1 is a good model related to the propulsive and aerodynamic characteristics of an UAV. As commented in this project, the UAVs used for power line inspections may have other devices that are capable of detecting failures in the network, such as infrared cameras, LiDAR technology or anti-collision sensors. This devices also have a energy consumption that must be taken into account. Then, it would be interesting to work with an energy model that includes energy consumption due to the operation of multiple electronic devices. • Sensitivity study. In table 3.2 the different characteristics of the UAV that affects the energy consumption are presented. A repetition of the results obtained in this project should be carried out but with different UAVs to show how the solutions change when the type of UAV used for the PTLIP varies. Thus, not only it can be found the solution that optimizes the path of the UAV, but also what type of UAV should be used to carry out that path. The same applies to battery capacity and the chosen speeds to inspect the power lines. • Multiple UAVs. For large electrical networks, the battery of a single UAV is not enough to inspect all the transmission lines. A quick solution would be to combine in a single UAV the batteries of multiple UAVs, obtain the solution for any of the three heuristics and, subsequently, distribute the number of tours among the available UAVs. However, a more detailed study would need to be done on how multiple UAVs can inspect a tree efficiently so that the workload of the UAVs is balanced. A suitable partition of the tree should be applied to address this problem. • Tree structure. The power lines are modeled as required edges of a graph that forms a tree. It should be interesting using heuristics that take advantage of the topology of tree-like structures. • Uncertainty. As stated in chapter 2, the assumption of a dynamic scenario is usual in the applications. Then, it should be interesting to consider a model with uncertainty and study the inspection problem in these scenarios.
List of Figures 1.1 Chopper doing a visual inspection of electric transmission lines in Minnesota, U.S. 2 1.2 Hexarotor with a camera, UAV ELIMCO and a portable ground station 3 1.3 (a) 3D test platform, (b) Diagram of UAV flight, (c) Infrared images, (d) Computer-generated red dots identifying tower before filtering and (e) after filtering 4 2.1 Development of the electricity transmission grid in Andalusia until 2014 7 3.1 Euler diagram for P, NP, NP-complete, and NP-hard set of problems for both assumptions P=NP and P=NP 14 3.2 Energy to travel an edge depending on CDwfor different speeds 17 4.1 Example of a power line distribution 26 4.2 Path followed by the UAV with H1(ME), H2(DW) and H3(ME) for 1V and p = 0 28 4.3 Path followed by the UAV in optimal solution for 1V and p = 0 29 4.4 Path followed by the UAV with H2(DW) with 1V and p = 0.5 31 4.5 Energy to travel an edge with different CDw32 4.6 Path followed by the UAV with H1(ME), H2(DW) and H3(ME) for 2V and p = 0 33 4.7 Path followed by the UAV in optimal solution for 2V and p = 0 34 4.8 Path followed by the UAV with H1(ME), H2(DW) and H3(ME) for 2V and p = 0.5 36 5.1 Ratio between sp(tc)and ¯ tcfor 10 to 50 nodes (left) and 60 to 100 nodes (right) in H1(ME) 1V 39 5.2 Ratio between sp(tc)and ¯ tcfor 10 to 50 nodes (left) and 60 to 100 nodes (right) in H2(DW) 1V 39 5.3 Ratio between sp(tc)and ¯ tcfor 10 to 50 nodes (left) and 60 to 100 nodes (right) in H3(ME) 1V 40 5.4 Relative deviation of median coverage time to mean coverage time H1(ME) 1V 40 5.5 Relative deviation of median coverage time to mean coverage time H2(DW) 1V 40 5.6 Relative deviation of median coverage time to mean coverage time H3(ME) 1V 41 5.7 Ratio between sp(M) and ¯ M for 10 to 50 nodes (left) and 60 to 100 nodes (right) in H1(ME) 1V 41 5.8 Ratio between sp(M) and ¯ M for 10 to 50 nodes (left) and 60 to 100 nodes (right) in H2(DW) 1V 41 5.9 Ratio between sp(M) and ¯ M for 10 to 50 nodes (left) and 60 to 100 nodes (right) in H3(ME) 1V 42 5.10 Relative deviation of median number of tours to mean number of tours H1(ME) 1V 42 5.11 Relative deviation of median number of tours to mean number of tours H2(DW) 1V 42 5.12 Relative deviation of median number of tours to mean number of tours H3(ME) 1V 42 5.13 Ratio between sp(te)and ¯ tefor 10 to 50 nodes (left) and 60 to 100 nodes (right) in H1(ME) 1V 43 5.14 Ratio between sp(te)and ¯ tefor 10 to 50 nodes (left) and 60 to 100 nodes (right) in H2(DW) 1V 43 5.15 Ratio between sp(te)and ¯ tefor 10 to 50 nodes (left) and 60 to 100 nodes (right) in H3(ME) 1V 44 5.16 Relative deviation of median computational time to mean computational time H1(ME) 1V 44 5.17 Relative deviation of median computational time to mean computational time H2(DW) 1V 44 97
98 List of Figures 5.18 Relative deviation of median computational time to mean computational time H3(ME) 1V 45 5.19 Ratio between sp(Tc)and ¯ TcOPT 1V 45 5.20 Ratio between sp(M)and ¯ Min %OPT 1V 46 5.21 Ratio between sp(Te)and ¯ TeOPT 1V 46 5.22 Relative deviation of median coverage time to mean coverage time OPT 1V 47 5.23 Relative deviation of median number of tours to mean number of tours OPT 1V 47 5.24 Relative deviation of median computational time to mean computational time OPT 1V 47 5.25 Mean coverage time H1(ME) for 10 to 50 nodes (left) and 60 to 100 nodes (right) 1V 48 5.26 Mean coverage time H2(DW) for 10 to 50 nodes (left) and 60 to 100 nodes (right) 1V 48 5.27 Mean coverage time H3(ME) for 10 to 50 nodes (left) and 60 to 100 nodes (right) 1V 48 5.28 Mean number of tours H1(ME) for 10 to 50 nodes (left) and 60 to 100 nodes (right) 1V 49 5.29 Mean number of tours H2(DW) for 10 to 50 nodes (left) and 60 to 100 nodes (right) 1V 49 5.30 Mean number of tours H3(ME) for 10 to 50 nodes (left) and 60 to 100 nodes (right) 1V 49 5.31 Mean computational time H1(ME) for 10 to 50 nodes (left) and 60 to 100 nodes (right) 1V 50 5.32 Mean computational time H2(DW) for 10 to 50 nodes (left) and 60 to 100 nodes (right) 1V 50 5.33 Mean computational time H3(ME) for 10 to 50 nodes (left) and 60 to 100 nodes (right) 1V 50 5.34 Mean coverage time for the heuristics and the optimal solution 1V 51 5.35 Mean relative deviation of the heuristics to optimal solution for coverage time 1V 52 5.36 Mean number of tours for the heuristics and the optimal solution 1V 53 5.37 Mean relative deviation of the heuristics to optimal solution for number of tours 1V 54 5.38 Relative deviation of mean coverage time and mean number of tours of H2(DW) to H1(ME) 1V 54 5.39 Proportion of trees where Tc in H1(ME) is less than Tc in H2(DW) but the number of tours is greater for p = 0, p = 0.4 and p = 0.7 1V 55 5.40 Median computational time for the heuristics (top) and for the optimal solution (bottom) 1V 56 5.41 Median computational time for p = 0.2 1V 57 5.42 Median computational time for p = 0.5 1V 57 5.43 Median computational time for p = 1 1V 58 5.44 Polynomial 2 trend (P2) and linear trend (L) for H1(ME) for different values of p 1V 59 5.45 Median relative deviation of the heuristics to the optimal solution of computational time 1V 59 5.46 RDE/(RDE + NRE) for H2(DW), H3(ME) and OPT 1V 60 5.47 Not required edges for H1(ME), H2(DW), H3(ME) and OPT 1V 60 5.48 Distribution of coverage time for N = 10 and N = 20 H2(DW) 2V 61 5.49 Distribution of coverage time for N = 30 and N = 40 H2(DW) 2V 61 5.50 Distribution of coverage time for N = 50 and N = 60 H2(DW) 2V 62 5.51 Distribution of coverage time for N = 70 and N = 80 H2(DW) 2V 62 5.52 Distribution of coverage time for N = 90 and N = 100 H2(DW) 2V 62 5.53 Distribution of coverage time for N = 10 and N = 20 H3(ME) 2V 63 5.54 Distribution of coverage time for N = 30 and N = 40 H3(ME) 2V 63 5.55 Distribution of coverage time for N = 50 and N = 60 H3(ME) 2V 63 5.56 Distribution of coverage time for N = 70 and N = 80 H3(ME) 2V 64 5.57 Distribution of coverage time for N = 90 and N = 100 H3(ME) 2V 64 5.58 Distribution of coverage time for N = 10 and N = 20 H1(ME) 2V 65 5.59 Distribution of coverage time for N = 30 and N = 40 H1(ME) 2V 65 5.60 Distribution of coverage time for N = 50 and N = 60 H1(ME) 2V 65 5.61 Distribution of coverage time for N = 70 and N = 80 H1(ME) 2V 66 5.62 Distribution of coverage time for N = 90 and N = 100 H1(ME) 2V 66 5.63 Distribution of number of tours for N = 20 2V 67 5.64 Distribution of number of tours for N = 40 2V 68 5.65 Distribution of number of tours for N = 60 2V 68 5.66 Distribution of number of tours for N = 90 2V 69
List of Figures 99 5.67 Distribution of computational time for N = 10 2V 70 5.68 Distribution of computational time for N = 30 2V 71 5.69 Distribution of computational time for N = 50 2V 71 5.70 Distribution of computational time for N = 70 2V 72 5.71 Distribution of computational time for N = 100 2V 72 5.72 Distribution of coverage time for N = 10 and N = 20 OPT 2V 73 5.73 Distribution of coverage time for N = 30 and N = 40 OPT 2V 73 5.74 Distribution of coverage time for N = 50 and N = 60 OPT 2V 73 5.75 Distribution of coverage time for N = 70 and N = 80 OPT 2V 74 5.76 Distribution of coverage time for N = 90 and N = 100 OPT 2V 74 5.77 Number of tours in optimal solution from N = 10 to N = 50 2V 75 5.78 Number of tours in optimal solution from N = 60 to N = 100 2V 75 5.79 Distribution of computational time for N = 10 and N = 20 2V 76 5.80 Distribution of computational time for N = 30 and N = 40 2V 76 5.81 Distribution of computational time for N = 50 and N = 60 2V 77 5.82 Distribution of computational time for N = 70 and N = 80 2V 77 5.83 Distribution of computational time for N = 90 and N = 100 2V 77 5.84 Distribution of coverage time for N = 10 and N = 20 2V 79 5.85 Distribution of coverage time for N = 30 and N = 40 2V 79 5.86 Distribution of coverage time for N = 50 and N = 60 2V 80 5.87 Distribution of coverage time for N = 70 and N = 80 2V 80 5.88 Distribution of coverage time for N = 90 and N = 100 2V 80 5.89 Mean coverage time for the heuristics and the optimal solution 2V 81 5.90 Distribution of number of tours of the heuristics and the optimal solution for N = 10 and N = 20 2V 82 5.91 Distribution of number of tours of the heuristics and the optimal solution for N = 30 and N = 40 2V 82 5.92 Distribution of number of tours of the heuristics and the optimal solution for N = 50 and N = 60 2V 83 5.93 Distribution of number of tours of the heuristics and the optimal solution for N = 70 and N = 80 2V 83 5.94 Distribution of number of tours of the heuristics and the optimal solution for N = 90 and N = 100 2V 83 5.95 Mean number of tours for the heuristics and the optimal solution 2V 84 5.96 Relative deviation of heuristics to optimal solution of Tefor N = 10 and N = 20 2V 84 5.97 Relative deviation of heuristics to optimal solution of Tefor N = 30 and N = 40 2V 85 5.98 Relative deviation of heuristics to optimal solution of Tefor N = 50 and N = 60 2V 85 5.99 Relative deviation of heuristics to optimal solution of Tefor N = 70 and N = 80 2V 85 5.100 Relative deviation of heuristics to optimal solution of Tefor N = 90 and N = 100 2V 86 5.101 Median coverage time for H1(ME), H2(DW), H3(ME) and OPT 2V 87 5.102 Median computational time for the heuristics p = 0.1 2V 88 5.103 Median computational time for the heuristics p = 0.4 2V 88 5.104 Median computational time for the heuristics p = 0.7 2V 89 5.105 Median computational time for the heuristics p = 1 2V 89 5.106 Required edges traveled in deadheading by the UAV 2V 90 5.107 Not required edges traveled by the UAV 2V 90 5.108 Mean coverage time trend for 1V and 2V model 91 5.109 Mean number of tours trend for 1V and 2V model 92 5.110 Median computational time trend for 1V and 2V model 93
List of Tables 3.1 Complexity of the RPP. S stands for symmetric and A for asymmetric 13 3.2 UAV characteristics to solve the problem 16 5.1 Mean coverage time for the heuristics and the optimal solution 1V 51 5.2 Mean number of tours for the heuristics and the optimal solution 1V 53 5.3 Proportion of trees where UAV spends less time inspecting in H1(ME) than in H2(DW) with more number of tours. Mean coverage time and mean number of tours for those trees for p = 0 1V 55 5.4 Median computational time for the heuristics and the optimal solution 1V 56 5.5 Polynomial trend of order 2 (P2) for H3(ME) and linear trend (L) for both H1(ME) and H2(DW) for different values of p 1V 58 5.6 Mean and median computational time for the optimal solution and relative deviation of median to mean 2V 78 5.7 Maximum and minimum computational time of the optimal solution 2V 78 5.8 Relative deviation of mean coverage time of H3(ME) to OPT 2V 81 5.9 Values of coverage time with minimum and maximum difference between H3(ME) and OPT 2V 81 101
List of Codes 3.1 OPT Algorithm 19 4.1 Local Algorithm 22 4.2 Semilocal Algorithm 23 4.3 Global Algorithm 24 4.4 Aleatory Tree Algorithm 25 103
Bibliography [1] Saurav Agarwal and Srinivas Akella, Approximation algorithms for the single robot line coverage problem, International Workshop on the Algorithmic Foundations of Robotics, Springer, 2020, pp. 534–550. [2] ,Line coverage with multiple robots, 2020 IEEE International Conference on Robotics and Automation (ICRA), IEEE, 2020, pp. 3248–3254. [3] European Aviation Safety Agency, Normativa europea de UAS/drones, 2020, Last accessed 4 April 2022. [4] Boaz Ben-Moshe, Power line charging mechanism for drones, Drones 5(2021), no. 4, 108. [5] Alvaro Calvo, Giuseppe Silano, and Jesús Capitán, Mission planning and execution in heterogeneous teams of aerial robots supporting power line inspection operations, International Conference on Unmanned Aircraft Systems, 2022. [6] James F Campbell, Ángel Corberán, Isaac Plana, and José M Sanchis, Drone arc routing problems, Networks 72 (2018), no. 4, 543–559. [7] Jingkui Cui, Youmin Zhang, Sha Ma, Yingmin Yi, Jing Xin, and Ding Liu, Path planning algorithms for power transmission line inspection using unmanned aerial vehicles, 2017 29th Chinese Control And Decision Conference (CCDC), IEEE, 2017, pp. 2304–2309. [8] Lyubomyr Demkiv, Massimiliano Ruffo, Giuseppe Silano, Jan Bednar, and Martin Saska, An application of stereo thermal vision for preliminary inspection of electrical power lines by mavs, 2021 Aerial Robotic Systems Physically Interacting with the Environment (AIRPHARO), IEEE, 2021, pp. 1–8. [9] Javier Galve Francés et al., Algorítmica: diseño y análisis de algoritmos funcionales e imperativos, Paradigma; (1993). [10] Américo Gaudín, Gabriel Madruga, and Carlos Rodríguez, Planificación de vuelo utilizando algoritmos evolutivos, (2020). [11] Gurobi Optimization, Inc., Gurobi optimizer reference manual, 2012. [12] Aric Hagberg, Pieter Swart, and Daniel S Chult, Exploring network structure, dynamics, and function using networkx, Tech. report, Los Alamos National Lab.(LANL), Los Alamos, NM (United States), 2008. 105