An alisis Din amico del Fallo de Rotores en un Hexac...
Transcript of An alisis Din amico del Fallo de Rotores en un Hexac...
Analisis Dinamico del Fallo de Rotores en un Hexacoptero
Daniel Felipe Torres Cardozo
Hugo Manuel Romero Pineda
Universidad Distrital Francisco Jose de Caldas
Facultad de Ingenierıa
Proyecto Curricular de Ingenierıa Electronica
Bogota D.C.
2018
Analisis Dinamico del Fallo de Rotores en un Hexacoptero
Daniel Felipe Torres Cardozo Codigo: 20082005107
Hugo Manuel Romero Pineda Codigo: 20082005113
Trabajo de grado para optar al tıtulo de:
Ingeniero Electronico en la modalidad de investigacion
Directora:
Diana Marcela Ovalle Martınez. PhD.
Lınea de Investigacion:
Senales y Control
Grupo de Investigacion:
IDEAS
Universidad Distrital Francisco Jose de Caldas
Facultad de Ingenierıa
Proyecto Curricular de Ingenierıa Electronica
Bogota D.C.
2018
Agradecimientos
Agradezco Dios y a mi familia por el apoyo durante todo es tiempo en el que deasarrolle
mi carrera, por prestarme incondicional apoyo en cada momento dificıl y hacer posible
lograr este gran objetivo como es el de ser un profesional.
Hugo Romero
Doy gracias a mi familia por el apoyo incondicional que me brindaron durante la
carrera. A todos los companeros que me acompanaron durante esta aventura, con los
cuales crecimos tanto academica como personalmente.
Daniel Torres
Especial agradecimiento a la profesora Diana Ovalle por el apoyo brindado como direc-
tora de grado de este proyecto la cual fue de vital importancia en este arduo proceso.
Resumen
Los Vehıculos Aereos no Tripulados (UAV,Unmanned Aerial Vehicles) han sido objeto
de estudio en multiples trabajos durante la ultima decada. Su popularidad ha crecido
gracias a su precio asequible, el cual es debido principalmente a los avances en minia-
turizacion electronica y a los aumentos en velocidad de procesamiento de dispositivos
digitales de bajo costo. Teniendo en cuenta su capacidad de maniobra, estos vehıcu-
los son utilizados en muchos campos como puede ser el reconocimiento en areas de
desastre, fotografıa, vigilancia y hasta entrega de mercancıa. Muchos esfuerzos se han
enfocado a modelar su dinamica y controlarla. Sin embargo, estudios sobre falla en uno
o mas rotores y analisis de su comportamiento no es comun encontrar, siendo este el
objetivo de este trabajo.
Este trabajo busca estudiar el comportamiento de un hexacoptero al momento de
presentar falla o dano en una o varias de sus helices durante el tiempo de vuelo o antes
de que inicie el vuelo y bajo un control no lineal. Para lograr el proposito, se definira
un modelo no lineal que represente de manera adecuada la dinamica del hexacoptero.
Asimismo, se planteara un problema de control optimo no lineal que permita encontrar
los controles necesarios para que el hexacoptero defina una trayectoria determinada en
un solo plano del espacio tridimensional. Finalmente, se simula la falla de uno o varios
rotores del hexacoptero, con el fin de evidenciar su implicacion en las capacidades
dinamicas del mismo.
Abstract
Unmanned Aerial Vehicles (UAVs) have been studied in multiple works during the last
decade. Its popularity has grown thanks to its affordable price, which is mainly due to
advances in electronic miniaturization and increases in the speed of processing of low-
cost digital devices. Considering their ability to maneuver, these vehicles are used in
many fields such as recognition in areas of disaster, photography, surveillance and even
delivery of goods. Many efforts have focused on modeling its dynamics and controlling
it. However, studies on failure in one or more rotors and analysis of their behavior is
not common to find, this is the objective of this work.
This work aims to study the behavior of a hexacopter at the moment of presenting
failure or damage in one or more of its propellers during the flight time or before the
flight starts and under a non-linear control. To achieve the purpose, a nonlinear model
that adequately represents the dynamics of the hexacopter will be defined. Likewise,
a problem of optimal non-linear control will be raised, which will allow to find the
necessary controls so that the hexacopter defines a certain trajectory in a single plane
of three-dimensional space. Finally, the failure of one or several rotors of the hexacopter
is simulated, in order to demonstrate its involvement in the dynamic capacities of the
same.
Contenido
Resumen VII
Abstract IX
Lista de Figuras XIII
Lista de Tablas XV
Lista de Sımbolos XV
Introduccion 1
1. Generalidades 3
1.1. Planteamiento del problema . . . . . . . . . . . . . . . . . . . . . . . . 3
1.2. Justificacion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 4
1.3. Antecedentes . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 4
1.4. Objetivos . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 9
1.4.1. Objetivo general . . . . . . . . . . . . . . . . . . . . . . . . . . 9
1.4.2. Objetivos especıficos . . . . . . . . . . . . . . . . . . . . . . . . 9
2. Analisis Dinamico del Hexacoptero 11
2.1. Aproximaciones a la modelizacion de hexacopteros . . . . . . . . . . . . 11
2.1.1. Modelamiento con aproximacion Euler-Newton . . . . . . . . . . 11
2.1.2. Modelamiento con aproximacion Euler-Lagrange . . . . . . . . . 14
2.1.3. Cuaterniones . . . . . . . . . . . . . . . . . . . . . . . . . . . . 15
2.2. Modelamiento del Hexacoptero . . . . . . . . . . . . . . . . . . . . . . 16
2.2.1. Analisis Cinematico del Hexacoptero . . . . . . . . . . . . . . . 17
2.2.2. Ecuaciones de Movimiento . . . . . . . . . . . . . . . . . . . . . 20
2.3. Verificacion del comportamiento del modelo . . . . . . . . . . . . . . . 23
xii CONTENIDO
3. Problema de Maniobrabilidad 28
3.1. Introduccion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 28
3.2. Problema de Maniobrabilidad . . . . . . . . . . . . . . . . . . . . . . . 28
3.3. Una Aproximacion de Control Optimo . . . . . . . . . . . . . . . . . . 29
3.4. Solucion Numerica del Problema (PPP t) . . . . . . . . . . . . . . . . . . . 30
3.4.1. Calculo del Gradiente de la Funcion de Costo . . . . . . . . . . 30
3.4.2. Las restricciones del problema . . . . . . . . . . . . . . . . . . . 31
3.5. Ejemplos numericos de la maniorabilidad del Hexacoptero . . . . . . . 31
4. Fallos en los rotores
(simulaciones y analisis) 39
4.1. Modelamiento de hexacopteros con falla en rotores . . . . . . . . . . . . 39
4.2. Simulacion Numerica . . . . . . . . . . . . . . . . . . . . . . . . . . . . 41
4.3. Resultados . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 42
4.3.1. Resultados con falla en el rotor 1 desde el inicio de simulacion . 42
4.3.2. Resultados con falla en 2 rotores contiguos desde el inicio de
simulacion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 48
4.3.3. Resultados con fallas en los rotores durante el vuelo . . . . . . . 54
5. Conclusiones y observaciones 65
6. Anexos 69
6.1. Codigos del Matlab . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 69
6.1.1. Codigo General . . . . . . . . . . . . . . . . . . . . . . . . . . . 69
6.1.2. Dinamica del hexacoptero . . . . . . . . . . . . . . . . . . . . . 71
6.1.3. Dinamica del coestado . . . . . . . . . . . . . . . . . . . . . . . 73
6.1.4. Calculo del Gradiente . . . . . . . . . . . . . . . . . . . . . . . . 79
Lista de Figuras
1-1. Rotacion motores. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 7
2-1. Sistema Coordinado fijo del Hexacoptero. . . . . . . . . . . . . . . . . . 18
2-2. Verificacion de la condicion de Hover. . . . . . . . . . . . . . . . . . . . 25
2-3. Desplazamiento sobre el eje z con velocidad en los rotores mayor a la
velocidad de hover. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 26
2-4. Comportamiento con τK > 0. . . . . . . . . . . . . . . . . . . . . . . . 27
3-1. Movimiento en el eje Z: (a) Avance en x, (b) Avance en y, (c)Avance
en z, (d) Angulo de balanceo (roll) φ, (e) Angulo de cabeceo (pitch) θ,
(f) Angulo de viraje (yaw) ψ. . . . . . . . . . . . . . . . . . . . . . . . 33
3-2. Velocidad de los rotores para el movimiento sobre el eje Z. . . . . . . . 34
3-3. Movimiento en el eje y: (a) Avance en x, (b) Avance en y, (c)Avance en
z, (d) Angulo de balanceo (roll) φ, (e) Angulo de cabeceo (pitch) θ, (f)
Angulo de viraje (yaw) ψ. . . . . . . . . . . . . . . . . . . . . . . . . . 35
3-4. Velocidad de los rotores para movimiento sobre el eje Y . . . . . . . . . 36
3-5. Movimiento en el eje x: (a) Avance en x, (b) Avance en y, (c)Avance en
z, (d) Angulo de balanceo (roll) φ, (e) Angulo de cabeceo (pitch) θ, (f)
Angulo de viraje (yaw) ψ. . . . . . . . . . . . . . . . . . . . . . . . . . 37
3-6. Velocidad de los rotores para movimiento sobre el eje X. . . . . . . . . 38
4-1. Movimiento en el eje x con falla en el rotor 1: (a) Avance en x, (b)
Avance en y, (c)Avance en z, (d) Angulo de balanceo (roll) φ, (e) Angulo
de cabeceo (pitch) θ, (f) Angulo de viraje (yaw) ψ. . . . . . . . . . . . 43
4-2. Velocidad de los rotores en r.p.m para movimieno en eje x con falla en
el rotor 1. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 44
4-3. Movimiento en el eje z con falla en el rotor 1: (a) Avance en x, (b)
Avance en y, (c)Avance en z, (d) Angulo de balanceo (roll) φ, (e) Angulo
de cabeceo (pitch) θ, (f) Angulo de viraje (yaw) ψ. . . . . . . . . . . . 45
xiv LISTA DE FIGURAS
4-4. Velocidad de los rotores en r.p.m para movimiento en eje z con falla en
el rotor 1. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 46
4-5. Movimiento en el eje y con falla en el rotor 1: (a) Avance en x, (b)
Avance en y, (c)Avance en z, (d) Angulo de balanceo (roll) φ, (e) Angulo
de cabeceo (pitch) θ, (f) Angulo de viraje (yaw) ψ. . . . . . . . . . . . 47
4-6. Velocidad de los rotores en r.p.m para movimiento en eje y con falla en
el rotor 1. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 48
4-7. Movimiento en el eje x con falla en el rotor 1: (a) Avance en x, (b)
Avance en y, (c)Avance en z, (d) Angulo de balanceo (roll) φ, (e) Angulo
de cabeceo (pitch) θ, (f) Angulo de viraje (yaw) ψ. . . . . . . . . . . . 49
4-8. Velocidad de los rotores en r.p.m para movimiento en eje x con falla en
el rotor 1 y 2. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 50
4-9. Movimiento en el eje y con falla en el rotor 1: (a) Avance en x, (b)
Avance en y, (c)Avance en z, (d) Angulo de balanceo (roll) φ, (e) Angulo
de cabeceo (pitch) θ, (f) Angulo de viraje (yaw) ψ. . . . . . . . . . . . 51
4-10.Velocidad de los rotores en r.p.m para movimiento en eje y con falla en
el rotor 1. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 52
4-11.Movimiento en el eje z con falla en 2 rotores consecutivos: (a) Avance
en x, (b) Avance en y, (c)Avance en z, (d) Angulo de balanceo (roll) φ,
(e) Angulo de cabeceo (pitch) θ, (f) Angulo de viraje (yaw) ψ. . . . . . 53
4-12.Velocidad de los rotores en r.p.m para movimiento en eje z con falla en
2 rotores consecutivos. . . . . . . . . . . . . . . . . . . . . . . . . . . . 54
4-13.Lazos PD para el control del hexacoptero. . . . . . . . . . . . . . . . . 56
4-14.Respuesta del control multilazos PD para el control del hexacoptero. . . 57
4-15.Velocidades del control multilazos PD para el control del hexacoptero. . 58
4-16.Esquema de realimentacion de estados para el hexacoptero. . . . . . . . 60
4-17.Esquema de realimentacion de estados para el hexacoptero. . . . . . . . 61
4-18.Respuesta del control por realimentacion de estados e integrador para
el hexacoptero. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 62
4-19.Comortamiento de las velocidades de los rotores para el control por
realimentacion de estados e integrador para el hexacoptero. . . . . . . . 63
Lista de Tablas
1-1. Clasificacion por peso. . . . . . . . . . . . . . . . . . . . . . . . . . . . 5
1-2. Clasificacion por sistema de propulsion. . . . . . . . . . . . . . . . . . . 5
1-3. Clasificacion por tipo de despegue. . . . . . . . . . . . . . . . . . . . . 6
1-4. Movimiento del hexacoptero. . . . . . . . . . . . . . . . . . . . . . . . . 7
2-1. Tabla de Componentes de Movimiento. . . . . . . . . . . . . . . . . . . 17
2-2. Tabla de Matrices de Movimiento. . . . . . . . . . . . . . . . . . . . . . 17
4-1. Definicion de constantes para los bloques PD de la Figura 4-13. . . . . 56
Lista de sımbolos
A Matrix que relaciona las velocidades angulares de las helices con las fuerzas y
los movimientos que producen sobre el hexacoptero.
k Vector de funciones no lineales asociadas a la cinematica del hexacoptero.
d Vector de funciones no lineales asociadas a la dinamica del hexacoptero.
e fZ Fuerza generada por las helices en el eje z del eje a bordo del hexacoptero.
h Funcion de coste del problema de control optimo asociado a la maniobrabilidad
del hexacoptero.
I Momento de inercia del hexacoptero
Ix Momento de inercia del hexacoptero alrededor del eje x.
Iy Momento de inercia del hexacoptero alrededor del eje y.
Iz Momento de inercia del hexacoptero alrededor del eje z.
K Conjunto acotado que se refiere a las limitantes fısicas de las variables de
control u.
l Porcion incremental del coste.
L Longitud de cada brazo del hexacoptero.
m Masa del hexacoptero.
O Origen del plano fijo en el hexacoptero.
P1 Matriz de ponderacion para el estado final deseado.
p vector de coestado del hexacoptero.
p1i i−esimo elemento de P1.
p Vector de primeras derivadas de las variables de coestado del hexacoptero.
p Velocidad angular de rotacion alrededor del eje x, medida desde el eje a bordo
del hexacopeto.
Q Matriz de ponderacion para el estado deseado.
q Velocidad angular de rotacion alrededor del eje y, medida desde el eje a bordo
del hexacopeto.
q Unidad cuaternion.
q∗ Conjugado cuaternion.
R Matriz de ponderacion para las variables de control.
Rζ Matriz de velocidades angulares
LISTA DE SIMBOLOS xvii
r Velocidad angular de rotacion alrededor del eje z, medida desde el eje a bordo
del hexacopeto.
t Vector de tiempo para la ejecuacion de maniobras del hexacoptero.
τττ Duracion maxima de ejecucion de una maniobra.
ue Vector de variables de control en un punto de equilibrio.
u Velocidad de avance longitudinal, medida desde el eje a bordo del hexacopeto.
v Velocidad de avance lateral, medida desde el eje a bordo del hexacopeto.
w Velocidad de descenso, medida desde el eje a bordo del hexacopeto.
x Vector de variables de estado del hexacoptero.
x Vector de primeras derivadas de las variables de estado del hexacoptero.
x0 Vector de condiciones iniciales de las variables de estado.
xd(t) Vector de comportamiento deseado de las variables de estado.
xe Vector de variables de estado en un punto de equilibrio.
xT Vector de estado final deseado.
x Eje longitudinal, dirigido hacia la helice #1.
y Eje transversal, dirigido al punto medio entre las helices #2 y #3.
z Eje normal, dirigido hacia abajo.
ε Tolerancia para el criterio de parada del algoritmo de gradiente.
λ Paso fijo del algoritmo de gradiente.
φ Angulo de rotacion alrededor del eje x.
θ Angulo de rotacion alrededor del eje y.
ψ Angulo de rotacion alrededor del eje z.
ηηη Vector de posiciones lineales y angulares del hexacoptero, medidas en el marco
de referencia inercial.
ηηη1 Vector de posiciones lineales del hexacoptero, medidas en el marco de referencia
inercial.
ηηη2 Vector de posiciones angulares del hexacoptero, medidas en el marco de refe-
rencia inercial.
ννν Vector de velocidades lineales y angulares del hexacoptero, medidas en el marco
de referencia abordo del mismo.
ννν1 Vector de velocidades lineales del hexacoptero, medidas en el marco de refe-
rencia abordo del mismo.
ννν2 Vector de velocidades angulares del hexacoptero, medidas en el marco de re-
ferencia abordo del mismo.
∇ Gradiente de la funcion que lo acompana.
∇x Gradiente de la funcion que lo acompana, respecto a las variables de estado.
∇u Gradiente de la funcion que lo acompana, respecto a las variables de control.
Ω Vector de las velocidades angulares de las helices del hexacoptero.
xviii LISTA DE SIMBOLOS
Ωi Velocidad angular de rotacion de la helice #i.
T Vector de fuerzas y momentos que actuan sobre el hexacoptero.
τ1 Vector de fuerzas que actuan sobre el hexacoptero.
τ2 Vector de momentos que actuan sobre el hexacoptero.
τK momento angular alrededor del eje x.
τM momento angular alrededor del eje y.
τN momento angular alrededor del eje z.
Introduccion
El presente trabajo de investigacion se refiere al estudio de las fallas de los rotores en los
hexacopteros ya sea al inicio de un vuelo o durante su trayecto, esto con el fin de poner
en estudio el comportamiento del sistema al momento de estar bajo el control no lineal.
La motivacion principal para realizar este estudio es el creciente uso de estos UAV con
fines que van desde recreativo hasta militar [1], y con la evolucion de la tecnologıa se
busca tener mejores prestaciones del vehıculo en su estabilidad, maniobrabilidad entre
otras.
El hexacoptero se encuentra propulsado por un grupo de actuadores los cuales se en-
cuentran compuestos por motor, rotor y propulsor (helice), cualquiera de estos com-
ponentes puede presentar falla ya sea por desgaste del motor, fractura del propulsor,
etc..., sin embargo para este estudio no se tiene presente que parte del actuador pre-
senta falla, sino que se considera como un sistema completo que puede funcionar o
presentar falla.
La finalidad de este trabajo de investigacion es analizar el comportamiento simulado
del hexacoptero con un control no lineal aplicado al momento de presentar falla ya sea
al comienzo de realizar un desplazamiento o durante el mismo, las fallas comprenden
el dano ya sea en uno o en varios de sus actuadores, esto con el fin de verificar la
estabilidad del sistema. Al igual se realizaron las consultas de los diferentes modela-
mientos que se aplican al sistema como lo son los mas representativos Euler-Newton,
Euler-Lagrange y cuaterniones ademas de herramientas para el diseno como lo son las
rotaciones matriciales en las cuales se abarca el modelamiento del sistema desde dife-
rentes puntos de vista.
El trabajo esta estructurado como se detalla a continuacion. En el primer capıtulo del
trabajo se muestra las generalidades del trabajo tal como lo es el planteamiento del
problema su respectiva justificacion donde se valida la motivacion de la investigacion y
antecedentes de los estudios realizados. En el segundo capıtulo se realiza el analisis sobre
diferentes trabajos realizados usando metodos de diseno para hexacoptero como lo son
2 INTRODUCCION
las aproximaciones de Euler-Newton, Euler-Lagrange y cuaterniones, donde finalmente
se tiende a un analisis newtoniano, que se detalla hasta la obtencion del modelo no
lineal del hexacoptero y su validacion. En un tercer capıtulo se plantea el problema
de maniobrabilidad del hexacoptero desde el punto de vista del control optimo, que
es resuelto para obtener las senales de control que le permiten al hexacoptero ir de
un punto a a un punto b. En el capıtulo cuarto, se presenta la simulacion de fallos
y el analisis de los resultados obtenidos. Finalmente en el capıtulo quito se presentan
algunas conclusiones y observaciones.
Capıtulo 1
Generalidades
1.1. Planteamiento del problema
Los vehıculos aereos no tripulados (UAV, Unmanned Aerial Vehicles) han sido objeto
de estudio en los ultimos anos. El auge y popularidad que han tenido estos vehıculos
se debe, principalmente, a los avances tecnologicos en materia de semiconductores de
estado solido, en particular a los avances en la minuaturizacion de estos dispositivos
y a los incrementos en la velocidad de procesamiento de los mismos; por otro lado,
aunque muy ligado a lo anterior, se encuentra la dismunicion sustancial en los precios
de los UAVs. Todo lo anterior incentiva la versatilidad de sus aplicaciones entre las
que se pueden mencionar las imagenes aereas, para uso cientıfico y de entretenimiento,
el monitoreo de cultivos, zonas de desastres, obras civiles y cualquier otro entorno de
observacion e incluso la entrega de mercancıa.
Teniendo en cuenta la popularidad de estos vehıculos y su diversidad de aplicaciones,
en la ultima decada se han dedicado muchos esfuerzos a entender la dinamica de es-
te tipo de vehıculos y las diversas posibilidades de controlar su movimiento. En los
ultimos dos o tres anos se han empezado a desarrollar esfuerzos en torno a analizar
fallas en rotores de hexacopteros, en muchos de los casos basandose en la linealizacion
del modelo dinamico del mismo y aplicando tecnicas de control adaptivo o de control
robusto [2], [3], [4].
De acuerdo a lo anterior, en este trabajo se hace un analisis dinamico de los fallos en
rotores de un hexacoptero, utilizando herramientas de simulacion, basicamente a partir
del modelo no lineal y la aplicacion de un control optimo no lineal al movimiento del
hexacoptero, con el fin de conocer sus implicaciones en el movimiento tridimensional
del mismo y, ademas, calcular los controles necesarios para mantener la trayectoria que
se tenıa antes de presentarse la falla, en caso de ser posible.
4 1 Generalidades
1.2. Justificacion
Justificacion Tecnica
Cuando se realiza el diseno de cualquier dispositivo, ademas de pensar en su rendimien-
to, que tan potente o eficiente sera, es necesario tambien pensar en como se comportara
el dispositivo en condiciones adversas, debido a que en el ambiente las condiciones no
son las ideales y estara expuesto a lugares que puedan desgastar rapidamente elementos
claves y causar estados de falla.
Para el caso del proyecto, el desarrollo del modelamiento del hexacoptero es el principio
para el diseno de los controladores y simulaciones, pensando en ambientes y condiciones
ideales, en que normalmente solo se toman en cuenta las fuerzas esenciales para acercar
el comportamiento de la simulacion a un entorno real. Se busca en el trabajo que se
realizara es analizar el caso extremo cuando el hexacoptero pierde uno o mas rotores
durante su vuelo o justo antes de iniciar el vuelo.
Justificacion Academica
Para el desarrollo de este trabajo, se utilizaran conceptos de sistemas dinamicos, de
sistemas de control y el uso de herramientas de simulacion, ampliando las competen-
cias que se adquieren luego de cursar el plan de estudios del Proyecto Curricular de
Ingenierıa Electronica hacia el control optimo no lineal, el modelamiento dinamico de
cuerpos en el espacio tridimensional y el uso mas especializado de herramientas como
Matlab®. Se espera fortalecer la formacion de los autores en el area de los sistemas
de control y la automatica.
1.3. Antecedentes
La aviacion no tripulada tuvo sus comienzos en los modelos construidos y volados por
inventores como Cayley, Stringfellow, Du Temple y otros pioneros de la aviacion, que
fueron previos a sus propios intentos de desarrollar aeronaves tripuladas a lo largo de
la primera mitad del siglo XIX. Estos modelos sirvieron como bancos de pruebas tec-
nologicos para el posterior desarrollo de modelos de mayor tamano con piloto a bordo
y, en este sentido, fueron los precursores de la aviacion tripulada [5].
El uso de los sistemas no tripulados se ha venido incrementando con el avance de la
tecnologıa en la cual se han realizado mejoras en los tamanos de los componentes, su
1.3 Antecedentes 5
consumo, rendimiento, etc. . . , esto junto con las necesidades que se presentan a dia-
rio con el fin de solventar problemas cada dıa mas complejos, los cuales han llegado
a comprometer la vida humana en dichas labores y disminuir costos[6]; tareas como
rescate, reconocimiento, entregas de paquetes, entre otras, se han mejorado para que
no sea necesario exponer a una persona con la labor directamente sino pueda estar en
una distancia prudente y realizar la accion mediante un sistema no tripulado [7].
Los UAV son vehıculos aereos no tripulados, los cuales son usados para aplicaciones
tanto militares como civiles, actualmente estos estan equipados con diferentes sensores
y dispositivos electronicos con el fin de dar un uso especializado [8]. De los sistemas
no tripulados son los que mayor cantidad abarcan ya que su uso es mayor ya sea en
labores especıficas o recreativas.
El desarrollo de estos dispositivos ya puntualmente para aeronaves no tripuladas se
presento desde sus inicios en el campo militar durante la primera y segunda guerra
mundial, con el fin de realizar labores de espionaje en terreno enemigo y evitar bajas
humanas en estos casos. Actualmente el uso en el campo civil se ha incrementado, y ha
aumentado el interes entre las universidades con el fin de realizar innovacion y optimizar
procesos de automatizacion en estos dispositivos a pequena escala y para aplicaciones
puntuales. De igual manera a traves de la historia se han dado clasificaciones a los
UAV ya sea por su peso, sistema de propulsion o tipo de despegue como se muestra en
las tablas 1-1, 1-2 y 1-3.
Tabla 1-1: Clasificacion por peso.
Peso
Micro <1Kg
Mini 1Kg - 10Kg
Pequeno 10Kg - 50Kg
Mediano 50Kg - 100Kg
Grande >100kg
Tabla 1-2: Clasificacion por sistema de propulsion.
Sisitema de prpulsion
Motores de helice Gasolina, disel, otros
Turbina de gas comprimido
Motores electricos Baterıas, energıa solar
Motores de hidrogeno
Al tratarse de una tecnologıa que podemos llamar novedosa muchas areas se encuentran
realizando desarrollos, en particular ya sea dando mayor soporte para los impactos a
gran velocidad como lo es el dron desarrollado por NCCR Robotics, o como se men-
ciono al inicio de este capıtulo, usados por AMAZON [9] para entrega de mercancıa,
6 1 Generalidades
Tabla 1-3: Clasificacion por tipo de despegue.
Tipo de despegue
Vertical
Ala RotativaHelicopteros
Quad-rotor
Auto-sustentadosDirigible
Globo aeroestatico
No verticalAla Flexible
Parapente
Ala delta
Ala fija Aereoplano
practicamente se puede observar que los productos existentes, buscan la manera de
acoplarse a estas nuevas tecnologıas con el fin de minimizar los costos y aumentar su
versatilidad. Se debe conocer que al igual que existen reglas para los vuelos tripulados
y se tiene un espacio virtual de division del espacio aereo, para los UAV se debe regir
la misma segmentacion y se debe dar paso a una integracion en las reglas que permita
una armonıa en los vuelos de cada equipo.
Visto de esa manera se han realizado muchos estudios sobre la ampliacion en el campo
de aplicacion para estos sistemas y a su vez se ha mejorado la tecnologıa de su fabri-
cacion, aunque tristemente se han realizado mejoras mas relevantes en el campo belico
ya que es en la guerra donde se realiza un mayor uso de estos equipos para minimizar
las bajas y obtener informacion de una manera mas rapida.
Teniendo presente lo anterior se debe tener una plataforma de interaccion Humano –
Dron sobre la cual se puedan realizar labores de una manera sencilla, esto con el fin
de brindar las ordenes correctas para el manejo del sistema, puede que a diferencia
con un piloto de avion las ordenes que se toman en caso de ser las erroneas estrıan
poniendo en riesgo la vida de muchas personas, en caso de los UAV se estarıa poniendo
en riesgo solo el equipo o la mercancıa que se lleve, de tal manera la interfaz debe ser lo
mas sencilla posible y a su vez lo mas completa para tener la informacion vital para el
manejo del sistema. Por tal motivo en caso de presentarse el dano sobre alguno de sus
rotores se debe informar por esta interfaz para que se realice la respectiva correccion
por parte del piloto [10]. Pero en caso de no tener una respuesta rapida o en caso de
estar en estado de desconexion con el piloto, el sistema de control implementado sobre
el UAV deberıa ser capaz de sortear la emergencia y evitar el siniestro de la aeronave.
Por tal motivo se realiza el estudio sobre los avances que se tienen en cuanto a mode-
lamiento, control y estudio de fallas en hexacopteros, ya que estos son los mas usados
debido a que tienen una mayor velocidad, maneorabilidad y capaciad de carga util.
En la figura 1-1, se puede observar como es la configuracion del exacorptero como un
1.3 Antecedentes 7
Y0
Z0
X0
Ω1Ω2
Ω3
Ω4Ω5
Ω6
Figura 1-1: Rotacion motores.
cuerpo rigıdo sobre el marco de referencia X0, Y0 y Z0, la distribucion de los rotores
Ωi y el setido de giro.
Del analisis de la figura 1-1, se obtiene la tabla 1-4 la cual indica que rotores deben
interactuar para realizar movimientos de avance, retroceso, ascenso, descenso y giros
sobre su propio eje. Sobre la tabla el signo = se toma como la velocidad de sustento
para el sistema, el signo + se tomarıa como un aumento sobre la velocidad angular ω
y el signo - se tomarıa como una disminucion sobre la velocidad obteniendo.
Tabla 1-4: Movimiento del hexacoptero.Ω1 Ω2 Ω3 Ω4 Ω5 Ω6
Sustento = = = = = =
Ascenso + + + + + +
Descenso - - - - - -
Avance = = - = = +
Retroceso = = + = = -
Giro a la derecha = + = + = +
Giro a la izquierda + = + = + =
La maniobrabilidad es uno de los temas de estudio academico que se presentan en
relacion al control de sistemas de vuelo, por tal motivo para realizar dichos estudios
es requerido tener una ecuacion caracterıstica del sistema la cual permita realizar los
estudios matematicamente, esto serıa el modelamiento para el sistema.
Tal como se ha venido tratando la velocidad de las helices de cada rotor son las que
controlaran la trayectoria del vehıculo, teniendo esto presente se toman como premisas
para el modelamiento que:
8 1 Generalidades
El hexacoptero es un cuerpo solido en tres dimensiones.
Su centro de masa esta localizado en el centro.
Los efectos giroscopicos se cancelan.
Los efectos externos por el roce con el aire son despreciables
Los efectos giroscopicos se cancelan por la disposicion de las helices, estas se ubican de
manera equidistantes con respecto al centro y tomando la ubicacion al rededor como
si fueran los vertices de un hexagono, lo cual simplifica el estudio del mismo.
Teniendo en cuenta las presentes premisas basicas sobre la dinamica del hexacoptero,
de igual manera conociendo el movimiento de los rotores para tener presente como serıa
el comportamiento del sistema, se proceden a generar las acciones de control sobre el
sistema las cuales estan definidas a mantener la estabilidad del sistema y a monitorear
la posicion del UAV en el espacio con el fin de determinar la accion a tomar para llegar
al punto deseado.
Las arquitecturas de control de UAVs integran diversos sensores de navegacion inclu-
yendo giroscopos, acelerometros, magnetometros, GPS y sensores barometricos, entre
otros. Normalmente existen diferentes modos de control para distintas condiciones de
despegue, aterrizaje y vuelo, las estrategias basicas de control incluyen aterrizaje, se-
guimiento de secuencias de puntos definidos por sus coordenadas GPS, con eventual
vuelo estacionario, y aterrizaje.
La robotica en sus diferentes estudios ha llevado al desarrollo de arquitecturas de con-
trol con percepcion del entorno, seguimiento de objetos en movimiento, programas de
planificacion de trayectorias, entre otros, un claro ejemplo es el LaTrax Alias Quad-
Rotor, disenado por la companıa LaTrax, una empresa fundada en 1974 por Jim Jenkins
la cual fue pionera e innovadora en los sistemas de radio-control durante la decada de
1970. Esta empresa se dedica a la construccion especializada de vehıculos radio-control
de altas prestaciones [11].
El modelo en mencion tiene en ejecucion rutinas que no afectan la maniobrabilidad
del UAV y la experiencia de vuelo para el usuario, pero que permiten recuperar la
estabilidad del UAV en cualquier momento o situacion en la que se encuentre.
Al igual que LaTrax [11] implementa en sus sistemas rutinas para mantener su posi-
cion existen companıas dedicadas a realizar “vuelos seguros” con estos UAV en ejemplo
1.4 Objetivos 9
existe Horizon Hobby, una empresa fundada en Octubre de 1985 por Rick Stephens es-
pecializada en el hobby de los vehıculos radio control, con sucursales alrededor del
mundo USA, UK, Francia, Alemania y China entre otros. Esta companıa ha desarro-
llado diferentes tecnologıas en las cuales se potencializa la destreza al manejar estos
vehıculos, una de ellas es llamada SAFE (Sensor Assisted Flight Envelope) [12], la cual
maneja varios niveles de “libertad” en el vuelo, llegando al punto en el que si el sistema
se encuentra fuera de control con solo presionar una palanca, retoma el curso del UAV
llevandolo a una altura segura y manteniendo un curso estable, algunos de las funciones
se mencionan a continuacion.
Control de altitud.
Limitacion de los angulos de cabeceo y balanceo.
Recuperacion de la posicion de vuelo estable.
Aterrizaje automatico.
1.4. Objetivos
1.4.1. Objetivo general
Estudiar por medio de simulacion el comportamiento de un hexacoptero cuando este
presenta averıa en uno o mas rotores, ya sea durante pleno vuelo o antes de emprender
vuelo.
1.4.2. Objetivos especıficos
1.4.2.1. Proponer un modelo matematico no lineal para un hexacoptero que se mueve de
forma tridimensional en un espacio abierto y validarlo a partir de simulacion.
1.4.2.2. Verificar la capacidad de maniobra de un hexacoptero a partir de la solucion
numerica de un problema de control optimo no lineal, que permita encontrar los
controles necesarios para llegar a un estado final a partir de un estado inicial,
ambos determinados.
1.4.2.3. Simular la falla de uno de los rotores del hexacoptero, con el fin de concluir
respecto a sus implicaciones en el movimiento tridimensional del mismo.
1.4.2.4. Simular la falla de dos rotores aledanos del hexacoptero, con el fin de concluir
respecto a sus implicaciones en el movimiento tridimensional del mismo.
10 1 Generalidades
1.4.2.5. Simular la falla de dos de los rotores del hexacoptero ubicados de forma simetri-
ca (ya sea respecto al eje x o al eje y), con el fin de concluir respecto a sus
implicaciones en el movimiento tridimensional del mismo.
Capıtulo 2
Analisis Dinamico del Hexacoptero
2.1. Aproximaciones a la modelizacion de hexacopte-
ros
Durante la busqueda bibliografica sobre el modelamiento dinamico de UAV se eviden-
cia que sobresalen tres metodos para la realizacion de este modelamiento matematico
y en donde se parte del concepto de estudiar el UAV como un cuerpo rıgido. En primer
lugar, en contramos el modelamiento basado netamente en conceptos Newtonianos,
donde se tienen en cuenta las fuerzas que actuan sobre el objeto. Una segunda aproxi-
macion, desde el punto de vista lagrangiano, donde se realiza un modelamiento a partir
de la suma de energıas. La tercera aproximacion, desde el concepto de cuaterniones,
el cual esta basado en el teorema de la rotacion de Euler, en el cual se establece que
cualquier desplazamiento de un cuerpo rıgido en un punto fijo equivale a una rotacion.
A continuacion se detallan las aproximaciones mencionadas, para finalmente decidir la
que se seguira para modelar el hexacoptero.
2.1.1. Modelamiento con aproximacion Euler-Newton
Esta aproximacion se tomarıa como las mas basica e intuitiva de realizar, debido que se
vale de las leyes de Newton, con el fin de modelar la dinamica del hexacoptero a partir
de las leyes de movimiento y a las fuerzas que interactuan. Se presentaran ejemplos de
como diferentes autores colocaron en practica esta aproximacion para el modelamiento
del hexacoptero.
En este ejemplo el autor busca ademas de realizar un modelamiento matematico del
hexacoptero si no tambien un sistema de control [13].
12 2 Analisis Dinamico del Hexacoptero
El esquema cinematico del modelo desde el que se realizara el modelamiento es el si-
guiente.
Xe
Ye
Ze
η1 =(x, y, z)
φ
ψ
θ
X0Y0
Z0
El marco de referencia se usara el basico que es el (x, y, z) donde el centro de gravedad
del hexacoptero esta sobre el eje z, ademas ωi corresponde a la direccion de rotacion
respecto al empuje Ti estos a su vez con respecto a los torques Mi de cada uno de los
rotores.
La leyes de Newton pueden ser escritas como vectores de la forma∑~F =
d
dt(m~V )
∣∣∣∣0
, (2-1)
∑~M =
d ~H
dt
∣∣∣∣0
, (2-2)
donde la parte de la izquierda representa la sumatoria de todas las fuerzas y momentos
que actuan sobre el heacoptero. La ecuacion 2-1 puede ser expresada sobre el marco de
referencia del hexacoptero de la siguiente manera∑Fxr = m(vxr + vzrωyr − vyrωzr),∑Fyr = m(vyr + vxrωzr − vzrωxr),∑Fzr = m(vzr + vyrωxr − vxrωyr).
En las ecuaciones anteriores se muestran de manera individual las fuerzas que actuan
sobre cada uno de los ejes. Adicionalmente, vir y ωir son componetes de velocidad y
velocidad angular respectivamente.
2.1 Aproximaciones a la modelizacion de hexacopteros 13
La ecuacion 2-2 puede ser expresada tambien sobre el marco de referencia del he-
xacoptero, pero se debe conocer la distribucion de masas. Se realiza el supuesto de
que la distribucion real sobre el hexacoptero esta dada de la siguiente manera: (x, y)r,
(x, z)r y (y, z)r, de manera simetrica. Entonces los productos de la inercia pueden ser
omitidos y 2-2 puede ser escrito como∑Mxr = Jxωxr + ωyrωzr(Jz − Jy),∑Myr = Jyωyr + ωzrωxr(Jx − Jz),∑Mzr = Jzωzr + ωxrωyr(Jy − Jx).
Al lado izquierdo se encuentra la suma de todos los momentos que actuan sobre cada
uno de los ejes y J representa los momentos de inercia que actuan sobre los ejes.
Se realiza ahora el modelo matematico del movimiento en el hexacoptero, teniendo en
cuenta que este se toma sobre el marco de referencia (x, y, z)0. La relacion entre la
velocidad angular de un cuerpo y la inercia en un marco de referencia fijo esta dada
por los movimientos angulares al rededor de los ejes x, y y z, llamados roll, pitch y yaw,
respectivamente. Por lo que las velocidades angulares pueden ser expresadas como
θ = ωyr cosφ− ωzr cosφ,
φ = ωzr + ωyr sinφ tan θ + ωzr cosφ tan θ,
ψ = ωyrsinφ
cos θ+ ωzr
cosφ
cos θ.
Finalmente, las fuerzas y los momentos que actuan sobre el hexacopero deben ser
especificadas. Las principales fuerzas son las de empuje de los rotores y la gravitacional,
por lo que sobre cada eje tendrıamos las sumatorias de fuerzas a continuacion:∑Fxr = −mg sin θ,∑
Fyr = mg cos θ sinφ,∑Fzr = T +mg cos θ cosφ,
donde fZ representa el empuje total de los rotores.
Lo momentos principalmente son reactivos y giroscopicos debido al empuje de los ro-
tores, por lo que los podrıamos expresar como
14 2 Analisis Dinamico del Hexacoptero
∑Mxr = Mx + Jpωyr
(−∑i=1,3,5
ωi +∑i=2,4,6
ωi
),
∑Myr = My + Jpωyr
( ∑i=1,3,5
ωi −∑i=2,4,6
ωi
),
∑Mzr = Mz,
donde Jp es el momento de inercia de un motor con su helice, ωi representa la velocidad
angular de cada uno los rotores. Mx,My y Mz son momentos causados por el empuje
de los rotores y los restantes son los momentos giroscopicos.
2.1.2. Modelamiento con aproximacion Euler-Lagrange
Esta aproximacion varios autores la usan para realizar el modelamiento dinamico del
UAV, independientemente que este sea un hexacoptero, un cuadricoptero u otros tipos
de UAV con diferentes disenos del que es objeto este trabajo. A continuacion se pre-
sentan un par de ejemplos en los que se realiza el modelamiento dinamico a partir de
la aproximacion Eluer-Lagrange (analisis a partir de energıas).
Este tipo de aproximacion se utiliza por ejemplo en [3] y [4] . Iniciamente, se definen dos
marcos de referencia, uno respecto a un punto fijo en la tierra nombrado por sus siglas
en ingles (EF-Earth Frame) y otro con respecto al cuerpo rıigido como tal (BF-Body
Frame). La fuerza de la masa actua sobre el centro de gravedad y siempre a lo largo del
eje negativo z. La energıa cinetica del hexacoptero debido al movimiento transalacional
es expresada como
Ttrans ,1
2−mξT , ξ
donde m es la masa del hexacoptero, ξ es el vector de posicion ξ = [x y z]T y η es el
vector de posicion angular η = [φ θ ψ]T .
La cinetica de hexacoptero se puede expresar debido al movimiento rotacional, el cual
puede ser expresado como
Trot ,1
2ηTJη,
donde J es la matriz de inercia expresada en terminos de la coordenadas η:
J = W Tη IW η,
2.1 Aproximaciones a la modelizacion de hexacopteros 15
con
Wη =
1 0 −Sθ0 Cφ SφCθ0 −Sφ CφCθ
,donde I es la matriz de inercia, la cual es una matriz diagonal, debido a la simetrıa del
hexacoptero.
La energıa potencial es expresada en funcion del campo gravitacional de la tierra
U = mgz.
Se definen de manera general las coordenadas como el vector q = [ξT ηT ]T = [x, y,
z,φ,θ, ψ]T ∈ R6, entonces el Lagrangiano se puede expresar como
L(q, q) = Ttrans + Trot − U =1
2mξT ξ +
1
2ηTJη −mgz (2-3)
El modelo total del hexacoptero es obtenido desde las ecuaciones de Euler-Lagrange
con la fuerza externa generalizada
d
dt
∂L∂q− ∂L∂q
= F. (2-4)
Donde F = [F Tξ τ
Tη ]T . Fxi contiene todas las fuerzas trasacionales debido a las entradas
de control y expresadas en EF x, y, z y τη en el vector que generaliza los torques.
2.1.3. Cuaterniones
Los angulos de Euler son ampliamente utilizados para describir la dinamica de un
cuadrotor. A pesar de que son muy intuitivos y faciles de interpretar y visualizar,
adolecen de singularidades. Ademas, la representacion de los angulos de Euler va de
la mano con el calculo del seno y el coseno, lo que aumenta el costo computacional.
La forma alternativa, mas eficiente y sin singularidades para describir la dinamica del
cuadrotor es usar un cuaternion. La orientacion se puede caracterizar por una sola
rotacion α alrededor de un eje [14], [16], [17]:
q[
cos α2
aT sin α2
]=[q0 q13
]=[q0 q1 q2 q3
]. (2-5)
La rotacion tridimensional de cualquier vector es descrita como una multiplicacion a la
izquierda por la unidad cuaternion q y a la derecha por su conjugado q∗, que se puede
escribir como una multiplicacion de la matriz Rq y el vector mencionado anteriormente,
donde [q13] es una matriz simetrica [14], [17]:
16 2 Analisis Dinamico del Hexacoptero
Rq = (q0I + [q13])2+q13q13T =
q20 + q2
1 − q22 − q2
3 2(q1q2 − q3q0) 2(q1q3 + q2q0)
2(q1q2 + q3q0) q20 − q2
1 + q22 − q2
3 2(q2q3 − q1q0)
2(q1q3 − q2q0) 2(q2q3 + q1q0) q20 − q2
1 − q22 + q2
3
.(2-6)
La relacion entre la velocidad angular y el cuaternion esta dada por
q =1
2q ·[
0
mη
]=
1
2Smη =
1
2
−q1 −q2 −q3
q0 −q3 q2
q3 q0 −q1
−q2 q1 q0
mη. (2-7)
Para derivar el modelo de la dinamca del cuadrotor en forma de un cuaternion, el vector
de velocidades angulares es reemplazado en el Lagraniano usando la expresion 2-3.
Posterior a esto se utiliza la ecuacion 2-4 para obtener la ecuacion de Euler-Lagrange
de la rotacion y traslacion del cuadrotor:
mτ − 2Jr
4∑i=1
(STmq
)×
0
0
(−1)iωi
− 4Kηdiag(STmq)(STmq
)=
2
(2Jmq + 2
d
dt(J)mq − ∂
∂t(mqTJmq)
)(2-8)
mmζ = Rq
0
0∑4i=1 Fi
− 0
0
mg
−Kvdiag( ˙Rqζ)( ˙Rqζ) (2-9)
2.2. Modelamiento del Hexacoptero
En la seccion anterior se mostraron de manera general las diferentes aproximaciones
para modelizar el comportamiento dinamico de un hexacoptero. Pese a que los modelos
obtenidos con las tres son equivalentes, a continuacion se presentara de forma detallada
el modelamiento desde la aproximacion newtoniana con el fin de obtener un modelo
claro que describa el comportamiento no lineal del hexacoptero.
Es facil notar que el hexacoptero puede describir libremente movimientos rotacionales
y traslacionales a lo largo de las tres dimensiones, es por esta razon que un hexacoptero
tiene seis grados de libertad. Por conveniencia, las variables que tomaremos en conside-
racion para el analisis dinamico del hexacoptero seran las posiciones lineales y angulares
del mismo y sus derivadas, como se muestra en la Tabla 2-1. Teniendo en cuenta el
2.2 Modelamiento del Hexacoptero 17
numero de variables, es conveniente utilizar una notacion vectorial como se representa
en la Tabla 2-2. Donde ηηη es el vector de posicion y orientacion del hexacoptero en el
marco de referencia fijo en la tierra, ννν es el vector de velocidades lineales y angulares en
el marco de referencia a bordo del hexacoptero y τττ es el vector de fuerzas y momentos
actuando sobre el hexacoptero, tambien referenciado respecto al hexacoptero.
Tabla 2-1: Tabla de Componentes de Movimiento.Hexacoptero
Componentes del Movimiento
Fuerzas y
Momentos
Velocidades
Lineales y
Angulares
Posiciones
Lineales y
Angulares
Movimiento en el eje x - u x
Movimiento en el eje y - v y
Movimiento en el eje z fZ w z
Rotaciones al rededor de x τK p φ
Rotaciones al rededor de y τM q θ
Rotaciones al rededor de z τN r ψ
Tabla 2-2: Tabla de Matrices de Movimiento.
ηηη1 = [x y z]T ηηη2 = [φ θ ψ]T ηηη = [ηηη1 ηηη2]T
ννν1 = [u v w]T ννν2 = [p q r]T ννν = [ννν1 ννν2]T
τττ 1 = [X Y Z]T τττ 2 = [K M N ]T τττ = [τττ 1 τττ 2]T
2.2.1. Analisis Cinematico del Hexacoptero
Para poder entender el analisis cinematico del hexacoptero, iniciamos definiendo un
sistema coordinado inercial Xe, Ye, Ze en la convencion NED (North, East, Down)
tambien denominado marco de referencia espacial. Los angulos de rotacion alrededor
de los ejes de este marco de referencia se definen como:
Roll: Rotacion sobre el eje Xe, positiva a la derecha.
Pitch: Rotacion sobre el eje Ye, positiva hacia abajo.
Yaw: Rotacion sobre el eje Ze , positiva desde el Norte al Este.
Tambien se define un sistema coordenado fijo centrado en el centro de gravedad del
hexacoptero X0 , Y0, Z0, como el mostrado en la Figura 2-1.
18 2 Analisis Dinamico del Hexacoptero
Ω5
Ω6
Ω1
Ω2Ω3
Ω4
X0
Z0
Y0
Figura 2-1: Sistema Coordinado fijo del Hexacoptero.
Con el fin de desarrollar un analisis cinematico adecuado para el hexacoptero, se hacen
las siguientes suposiciones:
La dinamica de los motores es mucho mas rapida que la del cuerpo rıgido y por
ende se desprecia.
El eje de rotacion de los motores es paralelo a Z0.
las fuerzas y momentos de propulsion dentro del plano del motor son pequenos y
son perpendiculares al plano del rotor.
A partir de las definiciones vectoriales en la tabla 2-2 y la definicion de los dos sistemas
coordenados, sabemos que los angulos ηηη2 = [φ θ ψ]T representa la orientacion del marco
de referencia fijo al hexacoptero respecto al marco de referencia inercial (espacial). Por
lo tanto, se relaciona la velocidad lineal espacial ηηη1 = [x y z]T y el vector de velocidades
angulares del hexacoptero ννν1 = [u v w]T a traves de una matriz de rotacion, ası
ηηη1 = R(ηηη2)ννν1, (2-10)
donde R(ηηη2) representa la rotacion alrededor de los tres ejes, definida como
R(ηηη2) = RTz (ψ)RT
y (θ)RTx (φ),
siendo RTx(φ) la traspuesta de la matriz de rotacion respecto al eje x en el angulo φ,
que puede ser escrita como
Rx(φ) =
1 0 0
0 cφ sφ0 −sφ cφ
,
2.2 Modelamiento del Hexacoptero 19
con cφ = cos(φ) y sφ = sin(φ). De manera similar,
Ry(θ) =
cθ 0 sθ0 1 0
−sθ 0 cθ
y Rz(ψ) =
cψ −sψ 0
sψ cψ 0
0 0 1
.Por lo que R(ηηη2) puede ser escrito como
R(ηηη2) =
cθcψ −cθsψ + sφsθcψ sφsψ + cφsθcψcθsψ cφcψ + sφsθsψ −sφcψ + cφsθsψ−sθ sφcθ cφcθ
El orden de las rotaciones se hace siempre considerando que primero ocurre la rotacion
respecto a x, luego a y y, finalmente, con respecto a z.
En este momento es importante precisar que para el modelo del hexacoptero conside-
ramos
−π2≤ φ ≤ π
2,−π2≤ θ ≤ π
2y 0 ≤ ψ ≤ 2π.
En este punto, buscamos una transformacion entre las velocidades angulares en el
marco de referencia abordo del hexacoptero y las variaciones de los angulos de Euler
definidas en el marco de referencia inercial (espacial). Para lo anterior, consideramos
ννν2 = [p, q, r]T como las velocidades angulares respecto al marco de referencia del he-
xacoptero, tal que:
R(ηηη2) = R(ηηη2)ννν2 y ννν2 = RT (ηηη2)R(ηηη2),
donde
ννν2 =
0 −r q
r 0 −p−q p 0
es una representacion matricial (anti-simetrica) del mapeo lineal
ννν1 → ννν2 × ννν1 = ννν 2ννν1.
Un calculo directo muestra que la transformacion J(.) de ηηη2 a ννν2 satisface
ννν2 = J(ηηη2)ηηη2 =
1 0 −sθ0 cφ sφcθ0 −sφ cφcθ
φθψ
.
20 2 Analisis Dinamico del Hexacoptero
J(ηηη2) es invertible siempre que cos(θ) 6= 0, entonces en la region de interes,
ηηη2 = J−1(ηηη2)ννν2 (2-11)
con
J−1(ηηη2) =
1 sφtθ cφtθ0 cφ −sφ0 −sφ/cθ cφ/cθ
.El conjunto de tranformaciones (2-10) y (2-11) constituyen las ecuaciones cinematicas
del hexacoptero, haciendo uso de estas, podemos definir el mapeo de ννν a ˙(ηηη) como:
T (ηηη) =
[R(ηηη2) 03x3
03x3 J−1(ηηη2)
],
para expresar las ecuaciones cinematicas como
ηηη2 = T (ηηη)ννν. (2-12)
2.2.2. Ecuaciones de Movimiento
El movimiento de un cuerpo rıgido puede ser descompuesto en sus componentes rota-
cionales y traslacionales. Por lo tanto, con el fin de describir la dinamica de hexacoptero
asumido como un cuerpo rıgido, utilizaremos las ecuaciones de Newton-Euler, que go-
biernan el movimiento linear y angular. Primero que todo, la fuerza que actua sobre el
hexacoptero esta dada por
τττ 1 = maaa = m(ννν1 + ννν2 × ννν1), (2-13)
donde m es la masa del hexacoptero y es asumida como constante y , aaa la aceleracion
del hexacoptero.
Para determinar las fuerzas que generan las helices sobre el hexacoptero, es importante
considerar que cada rotor i tiene una velocidad angular Ωi, y estan ubicados paralelos
al eje Z0. Por lo que la rotacion de los rotores genera una fuerza de empuje, F , paralela
al eje Z0, con componentes iguales a cero en los otros ejes, es decir
Fi =[0 0 FZ
]T,=[0 0 kΩi
2]T, (2-14)
siendo k la constante de elevacion.
El empuje total junto con la fuerza gravitacional representa la fuerza total que actua
sobre el hexacoptero. Dado que la fuerza de empuje esta expresada en el marco de
2.2 Modelamiento del Hexacoptero 21
referencia del hexacoptero, al igual que la fuerza total actuando sobre este y la gravedad
actua en el marco de referencia inercial, podemos escribir la fuerza total actuando sobre
el hexacoptero como:
τττ 1 = F −mgRT (ηηη2)~eee3,
con ~eee3 =[0 0 Fi
]T. Reemplazando la expresion de τττ 1 en (2-13) y despejando para
ννν1, tendremos
ννν1 =1
mF − g RT (ηηη2)~eee3 − ννν2 × ννν1. (2-15)
Con respecto al movimiento rotacional del hexacoptero, se inicia llamando III a la matriz
de inercia. El hexacoptero tiene una estructura simetrica respecto a todos los ejes del
marco de referencia centrado en el centro de masa del hexacoptero, por lo tanto la
matriz de inercia se puede escribir como
I =
Ix −Ixy −Ixz−Ixy Iy Iyz−Ixz −Iyz Iz
=
Ix 0 0
0 Iy 0
0 0 Iz
Lo anterior se deriva de las simetrıas con respecto a los ejes X0, Y0 y Z0, aunque la
simetrıa con respecto a este ultimo no es completa, asumiendo que no cambie mucho el
comportamiento del modelo. Podemos escribir el momento de inercia total que actua
sobre el hexacoptero como:
τττ 2 = III(ννν2 + ννν2 × ννν2) (2-16)
De la estructura geometrica del hexacoptero, es posible obtener informacion respecto
a los momentos que causan roll, pitch y yaw, ası
τττ 2 =
τKτMτN
=
3
4k l (Ω2
2 + Ω32 − Ω5
2 − Ω62)
k l (−Ω12 − Ω2
2
4+
Ω32
4+ Ω4
2 +Ω5
2
4− Ω6
2
4)
b (−Ω12 + Ω2
2 − Ω32 + Ω4
2 − Ω52 + Ω6
2)
, (2-17)
donde l es la distancia entre el eje del rotor y el centro de gravedad del hexacoptero y
b es la constante de arrastre de los motores. Por lo que para el movimiento rotacional
tendremos
I ννν2 + ννν2 × (Iννν2) = τττ 2,
entonces,
ννν2 = I−1(τττ 2 − ννν2 × (Iννν2)). (2-18)
22 2 Analisis Dinamico del Hexacoptero
De acuerdo a lo anterior, el modelo del hexacoptero se reduce a (2-12) y (2-18), de
forma explıcita se puede escribir como:
x =u cosψ cos θ − v(cosφ sinψ − sinφ cosψ sin θ) + w(sinφ sinψ + cosψ cosφ sin θ),
y =u cos θ sinψ + v(cosφ cosψ + sinφ sinψ sin θ)− w(cosψ sinφ− cosφ sinψ sin θ),
z =− u sin θ + v cos θ sinφ+ w cosφ cos θ,
φ =p+ q sinφ tan θ + r cosφ tan θ
θ =q cosφ− r sinφ,
ψ =rcosφ
cos θ+ q
sinφ
cos θ,
u =rv − qw − g sin θ,
v =pw − ru+ g cos θ sinφ,
w =qu− pv − fZm
+ g cosφ cos θ,
p =Iy − IzIx
qr +τKIx,
q =Iz − IxIy
pr − τMIy,
r =Ix − IyIz
pq − τNIz.
(2-19)
Tambien podrıamos escribir el modelo de una forma mas compacta, como:
xxx = fff(xxx,uuu), (2-20)
donde xxx =[ηηη′ ννν ′
]′sera el vector de variables de estado del hexacoptero y uuu =[
Ω1 Ω2 Ω3 Ω4 Ω5 Ω6
]sera el vector de variables de control, que se puede obtener
como funcion de τττ . Note que τττ =[0 0 fZ τK τM τN
]y que
1
kfZ
2√3 k l
τK
2
k lτM
1
bτN
︸ ︷︷ ︸
τττ
=
1 1 1 1 1 1
0 1 1 0 −1 −1
−2 −1 1 2 1 −1
−1 1 −1 1 −1 1
︸ ︷︷ ︸
AAA
Ω21
Ω22
Ω23
Ω24
Ω25
Ω26
︸ ︷︷ ︸
ΩΩΩ2
. (2-21)
Por lo que,
ΩΩΩ2 = A′ (AA′)−1A′ (AA′)−1A′ (AA′)−1 τττ ,
2.3 Verificacion del comportamiento del modelo 23
o, lo que es equivalente
Ω21
Ω22
Ω23
Ω24
Ω25
Ω26
=
1
60 −1
6−1
61
6
1
4− 1
12
1
61
6
1
4
1
12−1
61
60
1
6
1
61
6−1
4
1
12−1
61
6−1
4− 1
12
1
6
1
kfZ
2√3 k l
τK
2
k lτM
1
bτN
.
Por lo que el problema se puede ver en dos sentidos, o encontrar las velocidades de
las helices que generan las maniobras deseadas en el hexacoptero o buscar los valores
de las 4 componentes de τττ y, a partir de la transformacion anterior, encontrar las 6
velocidades ΩΩΩ de los rotores.
2.3. Verificacion del comportamiento del modelo
Con el fin de verificar que el modelo representa la dinamica de un hexacoptero, se im-
plemento un script de Matlab® (vease seccion 6.1.1), que permite observar la evolucion
de los estados del hexacoptero para valores de las velocidades de los rotores.
Es necesario definir son los parametros del hexacoptero, que incluyen sus propiedades
fısicas generales y las de los rotores, [14]:
m = 0,468 Kg,
g = 9,81 m/s2,
l = 0,225 m,
k = 0,00000298 Kg-m/rad2,
b = 0,000000114 N-m/ rad,
Ix = Iy = 0,004856 Kg-m2,
Iz = 0,008801 Kg-m2.
Para verificar el comportamiento del modelo, vamos a hallar su punto de equilibrio,
que para hexacopteros se denomina “hover”. En la posicion del equilibrio, como es
bien sabido, todas las derivadas del modelo son iguales a cero. Es decir, el punto de
24 2 Analisis Dinamico del Hexacoptero
equilibrio son las soluciones de (2-20) cuando x = 0. En el caso de un hexacoptero [15],
el equilibrio consiste en
xe = [xe ye ze 0 0 ψe 0 0 0 0 0 0],
con xe, ye, ze puntos arbitrarios del espacio, debido a que las variables x, y y z no
influyen en la dinamica del sistema y φe una orientacion arbitraria, porque esta variable
tampoco influye en el comportamiento dinamico del hexacoptero. Por su parte para los
controles, se encuentra que
ue =mg
6 k[1 1 1 1 1 1],
que genera la fuerza que contraresta el peso del cuerpo y hace que el hexacoptero “le-
vite”.
De acuerdo a lo anterior, las pruebas que haremos sobre el programa de simulacion
seran:
Verificacion de la condicion de hover, para Ωi = 506,7266 r.p.m., i = 1, · · · , 6.
El resultado se muestra en la Figura 2-2, donde se evidencia que el hexacoptero
mantiene indefinidamente la condicion inicial.
Simulacion para velocidades de rotores superiores a la velocidad de hover. Es-
pecıficamente se probo para Ωi = 520, i = 1, · · · , 6., donde lo que se espera
que suceda es que se generara mayor fZ , que generara un movimiento sobre el
eje z hacia arriba, es decir sobre el eje real negativo, ya que z se toma positivo
hacia abajo. En Figura 2-3 muestra se demuestra que efectivamente se presenta
el desplazamiento sobre el eje z.
Simulacion para velocidades tal que τK > 0. Si buscamos que el momento al
rededor del eje x sea mayor que cero, se necesita que Ω22+Ω2
3 > Ω25+Ω2
6, por lo que
dejaremos Ω2 = Ω3 = 520 y Ωi = 506,7266, i = 1, 4, 5, 6, esperando obtener un
movimiento en φ y por ende en y. La Figura 2-4 muestra los resultados esperados,
donde ademas se observa ademas un desplazamiento en el eje z, debido a que fZes mayor al valor de equilibrio.
2.3 Verificacion del comportamiento del modelo 25
0 2 4 6 8 10
t[x]
0
0.2
0.4
0.6
0.8
1
[m]
(a)
0 2 4 6 8 10
t[x]
0
0.2
0.4
0.6
0.8
1
[m]
(b)
0 2 4 6 8 10
t[x]
-10
-8
-6
-4
-2
0
[m]
(c)
0 2 4 6 8 10
t[x]
-10
-5
0
5
10
[gra
dos]
(d)
0 2 4 6 8 10
t[x]
-10
-5
0
5
10
[gra
dos]
(e)
0 2 4 6 8 10
t[x]
-10
-5
0
5
10
[gra
dos]
(f)
Figura 2-2: Verificacion de la condicion de Hover.
26 2 Analisis Dinamico del Hexacoptero
t (s)
0 2 4 6 8 10
m
-1
-0.5
0
0.5
1(a)
t (s)
0 2 4 6 8 10
m
-1
-0.5
0
0.5
1(b)
t (s)
0 2 4 6 8 10
m
-40
-30
-20
-10(c)
t (s)
0 2 4 6 8 10
Gra
do
s
-1
-0.5
0
0.5
1(d)
t (s)
0 2 4 6 8 10
Gra
do
s
-1
-0.5
0
0.5
1(e)
t (s)
0 2 4 6 8 10
Gra
do
s
-1
-0.5
0
0.5
1(f)
Figura 2-3: Desplazamiento sobre el eje z con velocidad en los rotores mayor a la
velocidad de hover.
2.3 Verificacion del comportamiento del modelo 27
0 2 4 6 8 10
t[x]
0
0.2
0.4
0.6
0.8
1
[m]
(a)
0 2 4 6 8 10
t[x]
-0.5
0
0.5
1
1.5
2
[m]
(b)
0 2 4 6 8 10
t[x]
-10
-8
-6
-4
-2
0
[m]
(c)
0 2 4 6 8 10
t[x]
0
50
100
150
[gra
dos]
(d)
0 2 4 6 8 10
t[x]
-10
-5
0
5
10
[gra
dos]
(e)
0 2 4 6 8 10
t[x]
-10
-5
0
5
10
[gra
dos]
(f)
Figura 2-4: Comportamiento con τK > 0.
Capıtulo 3
Problema de Maniobrabilidad
3.1. Introduccion
En este capıtulo, primero definimos el problema de maniobrabilidad del hexacoptero.
Luego, con el fin de aproximarnos al problema de maniobrabilidad del hexacoptero,
formularemos un problema de control optimo no lineal con restricciones. Para la so-
lucion numerica del problema de control optimo usaremos el metodo del descenso por
gradiente. Al final del capıtulo mostraremos la solucion numerica al problema con sus
restricciones.
3.2. Problema de Maniobrabilidad
Fısicamente, la trayectoria de un hexacoptero esta compuesta particularmente por
(xxx(t),uuu(t)) tal que ˙xxx(t) = fff(xxx(t),uuu(t)) y xxx(0) = xxx0, nosotros podemos llamar maniobra
a la porcion de trayectoria xxx(t). Para hexacopteros en general, las principales maniobras
de interes son:
Movimiento en lınea recta: Manteniendo el rumbo y la altura, avanzando tanto
sobre el eje x como sobre el eje y.
Movimiento de cambio de Altura: Manteniendo el rumbo y ascendiendo o des-
cendiendo a un punto determinado, o sencillamente manteniendo la posicion en
x y/o y.
El problema de maniobrabilidad para el hexacoptero se puede abordar de la siguiente
manera. Encontrar (si es posible ) un funcion de (xxx(t),uuu(t)) o trayectoria, que representa
la ley de estados del hexacoptero como ˙xxx(t) = fff(xxx(t),uuu(t)), tal que xxx(t) defina una
maniobra especifica de interes, la cual es fijada comenzando. Adicionalmente, podemos
requerir algun comportamiento optimo para alguna variable de estado o variable de
control durante la maniobra.
3.3 Una Aproximacion de Control Optimo 29
3.3. Una Aproximacion de Control Optimo
En terminos matematicos, el problema de maniobrabilidad puede ser formulado para
nuestro caso en concreto como se presenta a continuacion:
Dado un estado inicial xxx(0) = xxx0 ∈ Ω, y un estado final deseado xxxT ∈ Ω, buscamos un
control uuu(t) tal que resuelve:
(PT ) =
Minimizar uuu: h(uuu) = m(xxx(T )) +∫ T
0l(xxx(τ),uuu(τ), τ)δτ
Sujeto a:
xxx = fff(xxx(t),uuu(t));
xxx(0) = xxx0 ∈ Ω;
xxx(t) ∈ Ω y uuu(t) ∈ K, 0 ≤ t ≤ T.
(3-1)
Donde:
h(uuu) es la funcion de coste,
l(xxx(τ),uuu(τ), τ) es el coste incremental,
y
m(xxx(T )) es el coste terminal, todas estas son funciones escalares.
Tıpicamente,
m(xxx(T )) =1
2[xxx(T )− xxxT (t)]T PPP 1 [xxx(T )− xxxT (t)] (3-2)
y
l(xxx(t),uuu(t), t) =1
2[xxx(t)− xxxd(t)]T QQQ [xxx(t)− xxxd(t)] +
1
2uuu(t)T RRRuuu(t) (3-3)
con PPP 1, QQQ, RRR, matrices diagonales reales, y xxxd el estado deseado. PPP 1, QQQ son conside-
radas matrices positivas semi definidas, mientras que RRR es considerada extrictamente
positiva definida.
Ajustando los valores de PPP 1 y QQQ nosotros podemos pesar la importancia relativa de
la desviacion de cada estado de su valor deseado en el tiempo final TTT y durante todo
el tiempo de la maniobra, respectivamente. Para el momento en que se incrementan
los elementos de la diagonal de PPP 1, p1i nosotros le damos mas significancia que la
desviacion de xxxi(T ) para el valor deseado, por que haciendo P1i cero indicamos que el
valor final de xxxi no es importante para cualquier valor. De igual manera, podemos darle
un peso de relativa importancia a la desviacion de las variables de estado para algun
valor en particular a lo largo de toda la maniobra y el esfuerzo del control, ajustando
los elementos de las diagonales de QQQ y RRR, respectivamente.
30 3 Problema de Maniobrabilidad
3.4. Solucion Numerica del Problema (PPP t)
Hay varios metodos de optimizacion que pueden ser aplicados para solucionar (3-1).
Debido a la complejidad de la ley de estado y el numero de variables del problema, es
bastante razonable usar el metodo del descenso del gradiente. Brevemente, el esquema
de este metodo consiste en seguir los siguientes pasos:
1. Calcular el gradiente de la funcion de costo, 5h(uuu).
2. Escoger un tamano de paso pequeno para el paso tj > 0.
3. Actualizar para j ≥ 0 ∈ Z
uuu(j+1) = u(j) − t(j)5 h(u(j)). (3-4)
Hasta la convergencia. El criterio de parada puede ser escogido como:
|h(uuu(j+1)) − h(uuu(j))| ≤ ε|h(uuu(0))|, (3-5)
con ε > 0, una tolerancia conveniente.
3.4.1. Calculo del Gradiente de la Funcion de Costo
El paso crucial es el computo del gradiente 5h(u(k)). Esto puede ser obtenido usando
el metodo del estado adjunto que es descrito a continuacion:
1. Dado el control u(j), j ≥ 0, soluciona la ecuacion de estado:˙xxx(t) = fff(xxx(t),uuu(j)(t))
xxx(0) = xxx(j)(t)(3-6)
para obtener el estado xxx(j+1)(t).
2. Con el par (uuu(j)(t),xxx(j+1)(t)), resolvemos la ecuacion del sistema de ecuaciones
diferenciales para el estado adjunto ppp(t), hacia atras en el tiempo,ppp(t) = −5x l(xxx
(j+1)(t),uuu(j)(t))− [5xfff(xxx(j+1)(t),uuu(j)(t))]Tppp(t),
ppp(T ) = −5x m(xxx(j+1)(T )xxxT ),(3-7)
donde 5x es el gradiente con respecto a las variables de estado xxx. Ası, obtenemos
ppp(j+1)(t).
3. Finalmente, el gradiente de la funcion de coste puede ser encontrado como:
5h(uuu(j)) = 5ul(xxx(j+1)(t),uuu(j)(t)) + [5ufff(xxx(j+1)(t),uuu(j)(t))]Tppp(j+1)(t), (3-8)
donde 5u es el gradiente con respecto a las variables de control uuu.
3.5 Ejemplos numericos de la maniorabilidad del Hexacoptero 31
3.4.2. Las restricciones del problema
Como es visto en la formulacion (PT ) hay lımites tanto en el estado como en las variables
de control. En este trabajo, tratamos con estos lımites a seguir:
Si durante el proceso de optimizacion, el comportamiento de variables de esta-
do como son los angulos de roll ψ(t) o elevacion θ(t) estan fuera de las fron-
teras realısticas (−45o < φ < 45o ,−45o < θ < 45o , 0o < ψ < 360o), se
trata de mantener los lımites agregando o incrementando los parametros q4 y
q5 (de la matriz QQQ), respectivamente. Lo mismo es aplicable a las variables de
control (K = [0, 1000 r.p.s.] × [0, 1000 r.p.s.] × [0, 1000 r.p.s.] × [0, 1000 r.p.s.] ×[0, 1000 r.p.s.]× [0, 1000 r.p.s.]) en cuanto a los lımites fısicos asociados a los roto-
res. No obstante, incluso si despues del aumento de los valores de los parametros
de peso de la variable(s) de interes, estamos fuera de los lımites, se opta por
incrementar el tiempo final, ya que la razon de tal comportamiento del algoritmo
esta socidado a que el sistema no es fısicamente capaz de ejecutar la maniobra
deseada, dentro de los lımites establecidos. Obviamente, a nivel matematico, la
forma correcta de tratar con lımites tanto en el estado como en las variables de
control debe introducir funciones de barrera asociadas con las restricciones en la
funcion de costo; aunque no seguimos este acercamiento porque nos beneficiamos
del fısico del problema que hace innecesario el poner en practica un algoritmo
mas sofisticado.
Como es de esperarse, el hecho de incrementar los valores de q4 y q5, normalmente
repercutira en una necesidad en la reduccion del paso fijo del algoritmo, aumentando
el numero de iteraciones, y por ende el tiempo, para lograr la convergencia.
3.5. Ejemplos numericos de la maniorabilidad del
Hexacoptero
A continuacion se realizan pruebas de maniobrabilidad del hexacoptero con el fin de
validar el modelo matematico propuesto anteriormente y verificar que podemos encon-
trar los controles necesarios para describir una maniobra en particular y, en lo posible,
que la maniobra termine en la posicion de hover. Se realizan pruebas realizando movi-
mientos sencillos en cada uno de los 3 ejes, x, y y z.
Se busca que en las pruebas de desplazamiento sobre los ejes x y y que se realice el
movimiento sin modificar la posicion en el eje z y en la pruebas sobre el eje z que no
se tenga variacion sobre los ejes x y y, todo esto durante un tiempo definido en que se
32 3 Problema de Maniobrabilidad
evidencie estabilidad, para esta oportunidad se realizaron pruebas para un tiempo de
10 segundos.
Partiendo de la formulacion definida en (3-1), (3-2) y (3-3), se definen:
PPP 1 = diag(α),
QQQ = diag(γ),
RRR = diag(β),
y para cada una de las pruebas mostradas a continuacion, los valores para dichos valores
se detallan. Ası mismo, el paso del algoritmo es definido como λ y se decidio dejar fijo
en todas las iteraciones.
Se realizaron pruebas sobre el eje Z realizando movimientos desde z0 = −10
a z0 = −5, que teniendo en cuenta que el eje Z es negativo hacia arriba sera
un movimiento en ascenso. Los parametros del algoritmo se ajustaron como se
muestra a continuacion:
γ = [ 0 0 1 0 0 0 0 0 3 0 0 0]
α = [ 0 0 1 0 0 0 0 0 0 0 0 0]
xxxd = [ 0 0 −5 0 0 0 0 0 0 0 0 0]
xxx0 = [ 0 0 −10 0 0 0 0 0 0 0 0 0]
xxxf = [ 0 0 −5 0 0 0 0 0 0 0 0 0]
(3-9)
con λ = 1× 10−4.
Se resalta que, en este caso en particular, se decidio tanto penalizar el estado
final como el estado a lo largo de todo el tiempo de la maniobra, ya que si solo
se penaliza el estado final, el algoritmo no garantiza una velocidad en z de cero
al final de la maniobra, lo que podrıa implicar que el estado final no serıa el
estado de hover, que a su vez implicarıa perder control del vehıculo una vez
terminada la maniobra. Ası mismo, se decidio penalizar la velocidad w durante
todo el tiempo, que tiene cierta equivalencia con z, como se observa a partir de las
relaciones cinematicas, con el fin de que la velocidad sobre el eje Z se garantice
en cero lo mas pronto posible. Los resultados se muestran en las Figuras 3-1 y
3-2, donde se observa que para lograr el movimiento en el eje Z no se producen
desplazamientos ni lineales ni angulares en ninguna otra variable y que se logra
la altura deseada luego de 6 s. Respecto a la velocidad de todos los rotores se
observa una completa simetrıa, ya que cada uno de ellos aporta la sexta parte del
movimiento en z, cuando la maniobra se desarrolla exclusivamente en dicho eje.
3.5 Ejemplos numericos de la maniorabilidad del Hexacoptero 33
0 2 4 6 8 10
t (s)
-1
-0.5
0
0.5
1
m
(a)
0 2 4 6 8 10
t (s)
-1
-0.5
0
0.5
1
m
(b)
0 2 4 6 8 10
t (s)
-10
-8
-6
-4
m
(c)
0 2 4 6 8 10
t (s)
-1
-0.5
0
0.5
1
Gra
dos
(d)
0 2 4 6 8 10
t (s)
-1
-0.5
0
0.5
1
Gra
dos
(e)
0 2 4 6 8 10
t (s)
-1
-0.5
0
0.5
1
Gra
dos
(f)
Figura 3-1: Movimiento en el eje Z: (a) Avance en x, (b) Avance en y, (c)Avance
en z, (d) Angulo de balanceo (roll) φ, (e) Angulo de cabeceo (pitch) θ, (f) Angulo de
viraje (yaw) ψ.
34 3 Problema de Maniobrabilidad
0 2 4 6 8 10
t(s)
400
450
500
550
r.p.m
1
0 2 4 6 8 10
t(s)
400
450
500
550
r.p.m
2
0 2 4 6 8 10
t(s)
400
450
500
550
r.p.m
3
0 2 4 6 8 10
t(s)
400
450
500
550r.
p.m
4
0 2 4 6 8 10
t(s)
400
450
500
550
r.p.m
5
0 2 4 6 8 10
t(s)
400
450
500
550
r.p.m
6
Figura 3-2: Velocidad de los rotores para el movimiento sobre el eje Z.
Se realizaron pruebas sobre el eje Y realizando un movimiento desde y0 = 0 a yT =
1. Los parametros del algoritmo se ajustaron como se muestra a continuacion:
γ = [ 0 10 10 0 0 0 0 10 10 0 0 0]
α = [ 0 5 5 0 0 0 0 0 1 0 0 0]
xxxd = [ 0 1 −10 0 0 0 0 0 0 0 0 0]
xxx0 = [ 0 0 −10 0 0 0 0 0 0 0 0 0]
xxxf = [ 0 1 −10 0 0 0 0 0 0 0 0 0]
(3-10)
3.5 Ejemplos numericos de la maniorabilidad del Hexacoptero 35
con λ = 1× 10−7.
Se resalta que, en este caso en particular, ademas de penalizar tanto el estado
final como el estado a lo largo de todo el tiempo de la maniobra en la componente
y, tambien fue necesario penalizar en z, ya que los movimientos sobre x y y no
son posibles de realizar sin movimiento en z, pero este ultimo desea ser penali-
zado para que sea lo mınimo posible. De la misma forma, no solo se penalizo la
velocidad v, sino tambien w, de nuevo, para lograr que se alcance rapidamente el
estado final y garantizar que la maniobra es finalizada en la condicion de hover.
Los resultados se muestran en las Figuras 3-3 y 3-4, donde se observa que para
lograr el movimiento en el eje Y es necesario desplazamiento sobre el eje Z y al
rededor del eje X. Respecto a la velocidad de todos los rotores se observa una
simetrıa entre las parejas de rotores 1 y 4, 2 y 3 y, 5 y 6.
0 2 4 6 8 10
t (s)
-1
-0.5
0
0.5
1
m
(a)
0 2 4 6 8 10
t (s)
0
0.5
1
1.5m(b)
0 2 4 6 8 10
t (s)
-10.2
-10.1
-10
-9.9
-9.8
m
(c)
0 2 4 6 8 10
t (s)
-0.5
0
0.5
1
Gra
dos
(d)
0 2 4 6 8 10
t (s)
-1
-0.5
0
0.5
1
Gra
dos
(e)
0 2 4 6 8 10
t (s)
-1
-0.5
0
0.5
1
Gra
dos
(f)
Figura 3-3: Movimiento en el eje y: (a) Avance en x, (b) Avance en y, (c)Avance en z,
(d) Angulo de balanceo (roll) φ, (e) Angulo de cabeceo (pitch) θ, (f) Angulo de viraje
(yaw) ψ.
36 3 Problema de Maniobrabilidad
0 2 4 6 8 10
t(s)
506.5
507
507.5
508
r.p.m
1
0 2 4 6 8 10
t(s)
400
500
600
700
r.p.m
2
0 2 4 6 8 10
t(s)
400
500
600
700
r.p.m
3
0 2 4 6 8 10
t(s)
506.5
507
507.5
508
r.p.m
4
0 2 4 6 8 10
t(s)
200
400
600
r.p.m
5
0 2 4 6 8 10
t(s)
200
400
600
r.p.m
6
Figura 3-4: Velocidad de los rotores para movimiento sobre el eje Y .
Se realizaron pruebas sobre el eje X realizando un movimiento desde x0 = 0 a
xT = 1. Los parametros del algoritmo se ajustaron como se muestra a continua-
cion:γ = [ 10 0 10 0 0 0 10 0 10 0 0 0]
α = [ 5 0 5 0 0 0 0 0 0 0 0 0]
xxxd = [ 1 0 −10 0 0 0 0 0 0 0 0 0]
xxx0 = [ 0 0 −10 0 0 0 0 0 0 0 0 0]
xxxf = [ 1 0 −10 0 0 0 0 0 0 0 0 0]
(3-11)
3.5 Ejemplos numericos de la maniorabilidad del Hexacoptero 37
con λ = 1× 10−8.
Se resalta que, al igual que en el movimiento sobre el eje Y , ademas de penalizar
tanto el estado final como el estado a lo largo de todo el tiempo de la maniobra
en la componente x, tambien fue necesario penalizar en z. De la misma forma,
no solo se penalizo la velocidad u, sino tambien w, de nuevo, para lograr que se
alcance rapidamente el estado final y garantizar que la maniobra es finalizada en
la condicion de hover. Los resultados se muestran en las Figuras 3-5 y 3-6, donde
se observa que para lograr el movimiento en el eje X es necesario desplazamiento
sobre el eje Z y al rededor del eje Y . Respecto a la velocidad de todos los rotores
se observa una simetrıa entre los rotores 1, 2 y 6, y 2, 3, y 5.
0 2 4 6 8 10
t (s)
0
0.5
1
1.5
m
(a)
0 2 4 6 8 10
t (s)
-1
-0.5
0
0.5
1
m
(b)
0 2 4 6 8 10
t (s)
-10.2
-10.1
-10
-9.9
-9.8
m
(c)
0 2 4 6 8 10
t (s)
-1
-0.5
0
0.5
1
Gra
dos
(d)
0 2 4 6 8 10
t (s)
-1
-0.5
0
0.5
Gra
dos
(e)
0 2 4 6 8 10
t (s)
-1
-0.5
0
0.5
1
Gra
dos
(f)
Figura 3-5: Movimiento en el eje x: (a) Avance en x, (b) Avance en y, (c)Avance en z,
(d) Angulo de balanceo (roll) φ, (e) Angulo de cabeceo (pitch) θ, (f) Angulo de viraje
(yaw) ψ.
38 3 Problema de Maniobrabilidad
0 2 4 6 8 10
t(s)
300
400
500
600
r.p.m
1
0 2 4 6 8 10
t(s)
450
500
550
r.p.m
2
0 2 4 6 8 10
t(s)
450
500
550
600
r.p.m
3
0 2 4 6 8 10
t(s)
400
500
600
700r.
p.m
4
0 2 4 6 8 10
t(s)
450
500
550
600
r.p.m
5
0 2 4 6 8 10
t(s)
450
500
550
r.p.m
6
Figura 3-6: Velocidad de los rotores para movimiento sobre el eje X.
Capıtulo 4
Fallos en los rotores
(simulaciones y analisis)
4.1. Modelamiento de hexacopteros con falla en ro-
tores
En los ultimos anos la investigacion en el modelamiento de soluciones sobre fallas en
UAVs en general se ha disparado, esto es debido a que cada dıa estos dispositivos ganan
importancia en el campo militar, civil o de vigilancia y aun mas cuando los vehıculos
aereos son equipados con camaras y sensores muy costosos.
Varios estudios se han dedicado a estudiar el como salvaguardar la integridad del UAV
en caso de falla, principalmente a fallas como lo son degradaciones del sistema de pro-
pulsion, de manera mas puntual se puede hablar de la falla en sus rotores, ya sea de
solo un motor o de varios motores lo cual causa inestabilidad o perdida total del control
del equipo causando su siniestro. Se ha propuesto como solucion disenando un control
en lazo cerrado que compense la perdida o degradacion de los rotores y ası evitar la
perdida total del UAV [20], [21].
Los enfoques de control tolerante a fallas pueden clasificarse en dos enfoques, uno ac-
tivo y pasivo. La idea del enfoque pasivo es en encontrar un control que estabilice el
sistema sin la informacion sobre la falla. Este enfoque a un diseno de control robusto.
El enfoque activo por el contrario trabaja con la reconfiguracion del control basado
en la identificacion de la falla usando controladores LPV (Linear Parameter-Varaying)
basados en un modelamiento no lineal del hexacoptero. Para el caso de [20] se realiza
el diseno de un controlador tolerante a falla para compensar una degradacion de los
propulsores usando una tecnica con un vector de parametros parcialmente medidos.
Esto significa que un parametro adicional representa una la falla, esta es introducida a
404 Fallos en los rotores
(simulaciones y analisis)
al modelo y un controlador sin dependencia en este adicional parametro es disenado.
Basado en el ultimo hecho, este puede ser clasificado como un enfoque pasivo.
Otro enfoque de trabajo realizado con otra perspectiva es el presentado en [21], este
artıculo presenta un control en el que se combina un control de modo deslizante adap-
tativo con un control de asignacion el cual puede acomodar simultaneamente fallas en
los actuador en el mismo grupo de actuadores y mantener la estabilidad del sistema
en lazo cerrado, esto se hace con el fin que el control tolere la falla de un rotor como
de varios, donde no solo el control de re-asignacion necesita ser ajustado para redistri-
buir las senales de control. El esquema del control incluye dos esquemas separados de
control, un esquema de bajo nivel y otro de alto nivel, el de bajo nivel es el modulo
de asignacion o reasignacion usado para distribuir las senales de control en la cantidad
requerida a los propulsores el cual solo puede reconfigurar la distribucion de las senales
de control en presencia de falla de los propulsores, el de alto nivel, es construido por
un control de deslizamiento adaptativo el cual es empleado para mantener en general
el rendimiento del sistema de seguimiento.
Otra manera de manejar la falla en el sistema de propulsion se propone en [3] un siste-
ma de deteccion y asilamiento de la falla por sus siglas en ingles (FDI, Fault dectection
and isolation) para fallas en los actuadores. Esta solucion es basada en un observador
dedicado a diagnosticar aplicado en el modelo del hexacoptero no lineal. Para el desa-
rrollo del controlador se debe primero obtener el modelo matematico del hexacoptero,
depues de esto se disena los observadores para sistemas no lineales, par este caso se
uso un observador para sistemas no lineales de Lipschitz. El sistema de deteccion de
falla (FD, Fault Detection) es hecho de un modulo de generacion residual y un modulo
de evaluacion residual. El modulo de generacion residual es obtenido por el calculo
entre la diferencia de la salida estimada y la salida medida por el observador Thau
(non-linear observer). El resultado obtenido alimenta el modulo de evaluacion residual
el cual entrega resultados binarios 1 o 0 indicando la presencia o no de una falla. El
modulo de aislamiento de la falla (FI, Fault Isolation), es activado solo cuando la sali-
da del FD es alta, en este caso este calcula la diferencia entre las referencias de roll y
pitch establecidas en el controlador de lazo cerrado y estas medias desde el sistema, el
asilamiento es realizado por la comparacion de estas diferencias.
Por ultimo tomamos como referencia el punto de vista explicado en [4], este artıculo
explica un metodo basado en el analisis redundante o por sus siglas en ingles (RA) el
cual explota el espacio nulo del espacio estado, permitiendo a una matriz observable la
creacion de un con conjunto de pruebas residuales para aprender sobre comportamiento
4.2 Simulacion Numerica 41
del sistema. El RA lineal es extendido dentro el campo no lineal, llamado como NLAR
(Non-linear Analytical Redundancy). Tanto como SMO (Sliding Mode Observer) y
NLAR tecnicas usadas para generacion de residuales, sin embargo en esta ocasion
debido a que la NLAR muestra que la suavidad y la reconstruccion de la senal es
mejor que la SMO se uso la tecnica NLAR. Utilizando el algoritmo NLAR se deriva
los residuales de la dinamica del hexacoptero basado en el modelo de espacio estado.
4.2. Simulacion Numerica
El programa en el que se implemento la solucion de (3-1) fue hecho usando Matlab®,
para las pruebas con fallo en 1 y 2 rotores, para todo el tiempo de simulacion. Las ecua-
ciones diferenciales ordinarias de primer orden, resultantes del modelo del hexacoptero,
son resueltas usando la funcion de Matlab ode45, que es un solucionador de paso fijo
basado el metodo Runge-Kutta explıcito, de orden entre 4 y 5 (Shampine 1997). No
fue posible lograr resultados satisfactorios en las simulaciones de fallos en pleno vuelo,
debido a que la forma de simular el fallo es directamente modificar el valor de la veloci-
dad de una o dos helices de forma arbitraria y, dado que se busca resolver un problema
de minimizacion por el metodo de descenso por gradiente, el modificar directamente
las variables de control implica alterar el metodo de solucion, genera la no disminucion
del valor de la funcion de coste y, por ende, la no solucion del problema de optimizacion.
Teniendo en cuenta la situacion anterior, con el fin de observar las implicaciones de los
fallos en rotores a partir de algun instante durante el vuelo, se decidio disenar un par
de controladores multivariables con el fin de verificar si estos son capaces de soportar o
no el fallo. Especıficamente, se diseno una realimentacion de estado multivariable con
seguimiento por integrador, donde las variables de salida se seleccionaron como x, y, z
y ψ y se mantuvieron como variables de control las velocidades de los rotores; la reali-
mentacion de estados se diseno para el modelo linealizado del hexacoptero. Tambien
se disenaron lazos de control PD, para las mismas variables de la realimentacion de
estados, pero considerando como entradas la fuerza vertical fZ y los torque al rede-
dor de cada uno de los ejes coordenados τK , τM y τN , que resultan directamente de
las velocidades de los rotores. La simulacion del hexacoptero con los controladores se
hizo embebiendo el modelo en una S-function que pudiese ser ejecutada desde Simulink.
Como es de esperarse, los resultados de la solucion del problema de optimizacion en
(3-1) dependen de la inicializacion del algoritmo. Debido a lo anterior, se decidio inicia-
lizar cada una de las maniobras del hexacoptero en el punto de equilibrio o posicion de
“hover”, donde el hexacoptero logra mantenerse fijo en un punto del espacio en cual-
424 Fallos en los rotores
(simulaciones y analisis)
quier valor positivo del eje z. Para hallar el punto de equilibrio, siguiendo la definicion
clasica en [22], hallamos las soluciones de (2-20) cuando todos las derivadas son iguales
a cero, que resulta en todas las variables de estado iguales a cero y
fZe = mg,
y por ende
Ωie =
√mg
6 k= 506,7266 r.p.m.
Debido a lo anterior, la totalidad de las pruebas se realizan partiendo desde el mis-
mo punto inicial en el espacio tridimensional para los tres ejes y el movimiento se
completara en un tiempo de 10 segundos, la posicion inicial para el movimiento sera:
x0 = 0, φ0 = 0, u0 = 0, p0 = 0,
y0 = 0, θ0 = 0, v0 = 0, q0 = 0,
z0 = −10, ψ0 = 0, ω0 = 0, r0 = 0.
Ası mismo, los valores de los controles seran Ωie = 506,7266.
4.3. Resultados
A continuacion se muestran los diferentes resultados obtenidos, donde se observa la
respuesta del hexacoptero en situacion de falla sobre sus rotores. Se plantea un primer
escenario donde se tiene falla sobre el rotor 1, despues sobre los rotores contiguos 1 y
2, sobre los rotores opuestos 1 y 4 y, finalmente, pruebas de fallas durante el vuelo.
4.3.1. Resultados con falla en el rotor 1 desde el inicio de
simulacion
Se realizaron simulaciones con falla en el rotor 1, para un desplazamiento deseado sobre
el eje x, obteniendo una configuracion final para los vectores, como se muestra en la
ecuacion 3-11.
La posicion final para el eje x sera con xd = 1, yd = 0 y zd = −10. Para lograr el re-
sultado esperado se deben modificar los vectores α y γ, siendo α la ponderacion sobre
el estado final deseado, y γ la ponderacion sobre estado deseado a lo largo de todo el
tiempo de simulacion y la variable λ es el paso fijo de las iteraciones.
La figura 4-1 muestra el resultado de la simulacion cuando se desea alcanzar un des-
plazamiento de un metro sobre el eje x con falla en el rotor 1, en cuanto a los despla-
zamientos sobre los tres ejes x, y y z, donde se observa que se logra el desplazamiento
4.3 Resultados 43
en x, sin movimiento en y y con alguna pequena perturbacion en z. En la figura 4-2
se muestra la velocidad de cada uno de los rotores en r.p.m (revoluciones por minuto)
durante el tiempo de simulacion, donde se aprecia una simetrıa en el movimiento de
los rotores 2 y 3, por un lado, y 4 y 6, por el otro; ademas, se observa que el algortimo
tiende a no usar el rotor 5.
0 2 4 6 8 10
t (s)
0
0.5
1
1.5
m
(a)
0 2 4 6 8 10
t (s)
-1
-0.5
0
0.5
1
m
(b)
0 2 4 6 8 10
t (s)
-10.2
-10.1
-10
-9.9
-9.8
m
(c)
0 2 4 6 8 10
t (s)
-1
-0.5
0
0.5
1
Gra
dos
(d)
0 2 4 6 8 10
t (s)
-1
-0.5
0
0.5
Gra
dos
(e)
0 2 4 6 8 10
t (s)
-1
-0.5
0
0.5
1
Gra
dos
(f)
Figura 4-1: Movimiento en el eje x con falla en el rotor 1: (a) Avance en x, (b) Avance
en y, (c)Avance en z, (d) Angulo de balanceo (roll) φ, (e) Angulo de cabeceo (pitch)
θ, (f) Angulo de viraje (yaw) ψ.
444 Fallos en los rotores
(simulaciones y analisis)
0 2 4 6 8 10
t(s)
-1
0
1
r.p.m
1
0 2 4 6 8 10
t(s)
500
600
700
r.p.m
2
0 2 4 6 8 10
t(s)
620.5
621
621.5
622
r.p.m
3
0 2 4 6 8 10
t(s)
0
200
400
600
r.p.m
4
0 2 4 6 8 10
t(s)
620.5
621
621.5
622
r.p.m
5
0 2 4 6 8 10
t(s)
500
600
700
r.p.m
6
Figura 4-2: Velocidad de los rotores en r.p.m para movimieno en eje x con falla en el
rotor 1.
Continuando con las simulaciones de falla en el rotor 1, ahora para un desplazamiento
deseado sobre el eje z, obteniendo una configuracion final para los vectores, como se
muestra en la ecuacion 3-9.
La figura 4-3 muestra el resultado de la simulacion cuando se desea alcanzar un des-
plazamiento de 5 metros en descenso sobre el eje z con falla en el rotor 1, en cuanto
a los desplazamientos sobre los tres ejes x, y y z, donde se observa que se logra el
4.3 Resultados 45
desplazamiento en z, sin movimiento en x y y, pese al fallo, gracias a que los rotores
actuan directamente sobre este eje. En la figura 4-4 se muestra la velocidad de cada
uno de los rotores en r.p.m (revoluciones por minuto) durante el tiempo de simulacion,
donde se aprecia una simetrıa en el movimiento de los rotores 2, 3, 5 y 6; ademas, se
observa que el algortimo no usa en lo absoluto el rotor 4.
0 2 4 6 8 10
t (s)
-1
-0.5
0
0.5
1
m
(a)
0 2 4 6 8 10
t (s)
-1
-0.5
0
0.5
1
m
(b)
0 2 4 6 8 10
t (s)
-10
-8
-6
-4
m
(c)
0 2 4 6 8 10
t (s)
-1
-0.5
0
0.5
1
Gra
dos
(d)
0 2 4 6 8 10
t (s)
-1
-0.5
0
0.5
1
Gra
dos
(e)
0 2 4 6 8 10
t (s)
-1
-0.5
0
0.5
1
Gra
dos
(f)
Figura 4-3: Movimiento en el eje z con falla en el rotor 1: (a) Avance en x, (b) Avance
en y, (c)Avance en z, (d) Angulo de balanceo (roll) φ, (e) Angulo de cabeceo (pitch)
θ, (f) Angulo de viraje (yaw) ψ.
464 Fallos en los rotores
(simulaciones y analisis)
0 2 4 6 8 10
t(s)
-1
0
1
r.p.m
1
0 2 4 6 8 10
t(s)
500
550
600
650
r.p.m
2
0 2 4 6 8 10
t(s)
500
550
600
650
r.p.m
3
0 2 4 6 8 10
t(s)
-1
0
1r.
p.m
4
0 2 4 6 8 10
t(s)
500
550
600
650
r.p.m
5
0 2 4 6 8 10
t(s)
500
550
600
650
r.p.m
6
Figura 4-4: Velocidad de los rotores en r.p.m para movimiento en eje z con falla en
el rotor 1.
Continuando con las simulaciones de falla en el rotor 1, ahora para un desplazamiento
deseado sobre el eje y, obteniendo una configuracion final para los vectores, como se
muestra en la ecuacion 3-10
La figura 4-5 muestra el resultado de la simulacion cuando se desea alcanzar un des-
plazamiento de un metro sobre el eje y, de nuevo con falla en el rotor 1, en cuanto
a los desplazamientos sobre los tres ejes x, y y z, donde se observa que se logra el
desplazamiento en y, sin movimiento en x y con un ligero desplazamiento sobre el eje
z, lo que confirma la simetrıa respecto a los ejes x y y que presenta el hexacoptero. En
4.3 Resultados 47
la figura 4-6 se muestra la velocidad de cada uno de los rotores en r.p.m (revoluciones
por minuto) durante el tiempo de simulacion, donde se aprecia una simetrıa en el mo-
vimiento de los rotores 2 y 4, por un lado, y 3 y 6, por el otro; ademas, se observa que
el algortimo no usa en lo absoluto el rotor 5.
0 2 4 6 8 10
t (s)
-1
-0.5
0
0.5
1
m
(a)
0 2 4 6 8 10
t (s)
0
0.5
1
1.5
m
(b)
0 2 4 6 8 10
t (s)
-10.2
-10.1
-10
-9.9
-9.8
m
(c)
0 2 4 6 8 10
t (s)
-0.5
0
0.5
1
Gra
dos
(d)
0 2 4 6 8 10
t (s)
-1
-0.5
0
0.5
1
Gra
dos
(e)
0 2 4 6 8 10
t (s)
-1
-0.5
0
0.5
1
Gra
dos
(f)
Figura 4-5: Movimiento en el eje y con falla en el rotor 1: (a) Avance en x, (b) Avance
en y, (c)Avance en z, (d) Angulo de balanceo (roll) φ, (e) Angulo de cabeceo (pitch)
θ, (f) Angulo de viraje (yaw) ψ.
484 Fallos en los rotores
(simulaciones y analisis)
0 2 4 6 8 10
t(s)
-1
0
1
r.p.m
1
0 2 4 6 8 10
t(s)
500
600
700
800
r.p.m
2
0 2 4 6 8 10
t(s)
500
600
700
800
r.p.m
3
0 2 4 6 8 10
t(s)
-1
0
1
r.p.m
4
0 2 4 6 8 10
t(s)
400
500
600
700
r.p.m
5
0 2 4 6 8 10
t(s)
400
500
600
700
r.p.m
6
Figura 4-6: Velocidad de los rotores en r.p.m para movimiento en eje y con falla en
el rotor 1.
4.3.2. Resultados con falla en 2 rotores contiguos desde el
inicio de simulacion
Se realizan simulaciones con falla en los rotores 1 y 2, para movimientos independientes
sobre los ejes x, y y z. Inicialmente, buscamos un movimiento de un metro sobre el eje
x, para lograrlo se configuro el algoritmo con los siguientes parametros mostrados en
la ecuacion 3-11
4.3 Resultados 49
La figura 4-7 muestra el resultado de la simulacion cuando se desea alcanzar un des-
plazamiento de un metro sobre el eje y, ahora en la presencia de fallas simultaneas en
los rotores 1 y 2, en los desplazamientos sobre los tres ejes x, y y z, se observa que
se logra el desplazamiento en x, sin movimiento en y y con un ligero desplazamiento
sobre el eje z. En la figura 4-8 se muestra la velocidad de cada uno de los rotores
en r.p.m (revoluciones por minuto) durante el tiempo de simulacion, donde ya no se
aprecia simetrıa de movimiento en ninguno de los rotores funcionales respecto a otros;
lo que se mantiene, respecto a cuando se considero unicmente la falla del rotor 1, es el
uso despreciable del rotor 5.
0 2 4 6 8 10
t (s)
0
0.5
1
1.5
m
(a)
0 2 4 6 8 10
t (s)
-1
-0.5
0
0.5
1
m
(b)
0 2 4 6 8 10
t (s)
-10.2
-10.1
-10
-9.9
-9.8
m
(c)
0 2 4 6 8 10
t (s)
-1
-0.5
0
0.5
1
Gra
dos
(d)
0 2 4 6 8 10
t (s)
-1
-0.5
0
0.5
Gra
dos
(e)
0 2 4 6 8 10
t (s)
-1
-0.5
0
0.5
1
Gra
dos
(f)
Figura 4-7: Movimiento en el eje x con falla en el rotor 1: (a) Avance en x, (b) Avance
en y, (c)Avance en z, (d) Angulo de balanceo (roll) φ, (e) Angulo de cabeceo (pitch)
θ, (f) Angulo de viraje (yaw) ψ.
504 Fallos en los rotores
(simulaciones y analisis)
0 2 4 6 8 10
t(s)
-1
0
1
r.p.m
1
0 2 4 6 8 10
t(s)
-1
0
1
r.p.m
2
0 2 4 6 8 10
t(s)
800
850
900
r.p.m
3
0 2 4 6 8 10
t(s)
0
200
400
600r.
p.m
4
0 2 4 6 8 10
t(s)
0
200
400
r.p.m
5
0 2 4 6 8 10
t(s)
700
800
900
1000
r.p.m
6
Figura 4-8: Velocidad de los rotores en r.p.m para movimiento en eje x con falla en
el rotor 1 y 2.
Continuando con las fallas de los rotores 1 y 2, ahora con el fin de lograr un movimiento
de un metro sobre el eje y, para lograrlo se configuro el algoritmo con los siguientes
parametros en la ecuacion 3-10.
La figura 4-9 muestra el resultado de la simulacion cuando se desea alcanzar un des-
plazamiento de un metro sobre el eje y, de nuevo en la presencia de fallas simultaneas
en los rotores 1 y 2, en los desplazamientos sobre los tres ejes x, y y z, se observa que se
4.3 Resultados 51
logra el desplazamiento en y, sin movimiento en x y con un ligero desplazamiento sobre
el eje z. En este caso, los resultados no son tan similares al movimiento sobre el eje
x mostrado previamente porque intencionalmente se selecciono un paso del algoritmo
mas grande, con el fin de ilustrar la dependencia de los resultados con este parametro.
0 2 4 6 8 10
t (s)
-1
-0.5
0
0.5
1
m
(a)
0 2 4 6 8 10
t (s)
0
0.5
1
1.5
m
(b)
0 2 4 6 8 10
t (s)
-10.2
-10.1
-10
-9.9
-9.8
m
(c)
0 2 4 6 8 10
t (s)
-0.5
0
0.5
1
Gra
dos
(d)
0 2 4 6 8 10
t (s)
-1
-0.5
0
0.5
1
Gra
dos
(e)
0 2 4 6 8 10
t (s)
-1
-0.5
0
0.5
1
Gra
dos
(f)
Figura 4-9: Movimiento en el eje y con falla en el rotor 1: (a) Avance en x, (b) Avance
en y, (c)Avance en z, (d) Angulo de balanceo (roll) φ, (e) Angulo de cabeceo (pitch)
θ, (f) Angulo de viraje (yaw) ψ.
En la figura 4-9 se muestra la velocidad de cada uno de los rotores en r.p.m (revolucio-
524 Fallos en los rotores
(simulaciones y analisis)
nes por minuto) durante el tiempo de simulacion, donde de nuevo es imposible apreciar
simetrıa de movimiento en ninguno de los rotores funcionales respecto a otros; lo que
se mantiene, respecto a cuando se considero unicmente la falla del rotor 1, es el uso
despreciable del rotor 4 (notese la escala vertical en 10−5).
0 2 4 6 8 10
t(s)
-1
0
1
r.p.m
1
0 2 4 6 8 10
t(s)
-1
0
1
r.p.m
2
0 2 4 6 8 10
t(s)
800
900
1000
1100
r.p.m
3
0 2 4 6 8 10
t(s)
9.24
9.25
9.26
9.27
r.p.m
10-6 4
0 2 4 6 8 10
t(s)
0
200
400
r.p.m
5
0 2 4 6 8 10
t(s)
877
878
879
880
r.p.m
6
Figura 4-10: Velocidad de los rotores en r.p.m para movimiento en eje y con falla en
el rotor 1.
Finalmente, mostramos los resultados para un desplazamiendo sobre el eje z de 5
metros en descenso, de nuevo con 2 rotores fallando simulataneamente desde el inicio
de la simulacion. Para lograrlo se configuro el problema con los siguientes parametros
que se muestran en la ecuacion 3-9.
4.3 Resultados 53
La figura 4-11 muestra el resultado de la simulacion cuando se desea alcanzar un
desplazamiento de 5 metros en descenso sobre el eje z, de nuevo en la presencia de
fallas simultaneas en los rotores 1 y 2. En los desplazamientos sobre los tres ejes x,
y y z, se observa que se logra el desplazamiento en z, sin movimiento en x y y, tal y
como ocurrio en la presencia de falla de un solo rotor. En la figura 4-12 se muestra la
velocidad de cada uno de los rotores en r.p.m. durante el tiempo de simulacion, donde
se aprecia simetrıa entre los rotores 3 y 6, unicos en funcionamiento. El algoritmo ante
la falla de los rotores 1 y 2, “decide” no utilizar sus simetricos, los rotores 4 y 5, (notese
la escala vertical en 10−5 para el rotor 4).
0 2 4 6 8 10
t (s)
-1
-0.5
0
0.5
1
m
(a)
0 2 4 6 8 10
t (s)
-1
-0.5
0
0.5
1
m
(b)
0 2 4 6 8 10
t (s)
-10
-8
-6
-4
m
(c)
0 2 4 6 8 10
t (s)
-1
-0.5
0
0.5
1
Gra
dos
(d)
0 2 4 6 8 10
t (s)
-1
-0.5
0
0.5
1
Gra
dos
(e)
0 2 4 6 8 10
t (s)
-1
-0.5
0
0.5
1
Gra
dos
(f)
Figura 4-11: Movimiento en el eje z con falla en 2 rotores consecutivos: (a) Avance
en x, (b) Avance en y, (c)Avance en z, (d) Angulo de balanceo (roll) φ, (e) Angulo de
cabeceo (pitch) θ, (f) Angulo de viraje (yaw) ψ.
544 Fallos en los rotores
(simulaciones y analisis)
0 2 4 6 8 10
t(s)
-1
0
1
r.p.m
1
0 2 4 6 8 10
t(s)
-1
0
1
r.p.m
2
0 2 4 6 8 10
t(s)
750
800
850
900
r.p.m
3
0 2 4 6 8 10
t(s)
8
9
10
r.p.m
10-6 4
0 2 4 6 8 10
t(s)
-1
0
1
r.p.m
5
0 2 4 6 8 10
t(s)
750
800
850
900
r.p.m
6
Figura 4-12: Velocidad de los rotores en r.p.m para movimiento en eje z con falla en
2 rotores consecutivos.
4.3.3. Resultados con fallas en los rotores durante el vuelo
A partir de la solucion del problema de maniobrabilidad definido en (3-3) no fue posi-
ble lograr resultados aceptables para la simulacion de fallas a partir de algun instante
durante el vuelo. Lo anterior es debido a que la simulacion de las fallas durante todo
el vuelo se obtuvieron al omitir de la optimizacion la velocidad del rotor o rotores
4.3 Resultados 55
particulares y mandandola(s) a cero, por lo que de alguna forma es como si la falla se
tuviese en cuenta directamente en la dinamica del sistema y el metodo de descenso por
gradiente es capaz de encontrar la solucion del problema de optimizacion. En el caso
de las fallas que ocurren, no durante todo el vuelo, sino a partir de cierto instante de
tiempo, lo que se intento hacer fue modificar el valor de la velocidad del rotor a partir
de cierto punto, pero como esto es imposible que los “vea” el algoritmo, el metodo de
descenso por gradiente no puede encontrar la solucion, debido a que en cada iteracion
nosotros modificamos la solucion que se encuentra.
Con el fin de estudiar las implicaciones de las fallas en los rotores a partir de cierto
instante de tiempo, se decidio disenar un controlador multivariable que permitieran
el control de las variables x, y, z y ψ y ver si el controlador en cuestion es capaz de
recuperar el sistema, o definitivamente sera mas viable implementar tecnicas de control
tolerante a fallos.
Lazos de control PD
Para el diseno de lazos de control PD que logren estabilizar el sistema no lineal, lo
primero que se decidio hacer fue linealizar el sistema y tratar de encontrar funciones de
transferencia entre la fuerza y los torques producidos por las velocidades de las helices
y alguna variable de estado. Es ası como se llega a que
Z(s)
fZ(s)= − 1
m
1
s2,
Φ(s)
τK(s)=
1
Ix
1
s2,
Θ(s)
τK(s)= − 1
Iy
1
s2,
Ψ(s)
τK(s)= − 1
Iz
1
s2.
De donde es evidente el doble integrador para cada lazo de control, lo que nos sugiere la
implementacion de controles PD para cada uno de los lazos. Sin embargo, para proposi-
tos de comparacion de los resultados obtenidos en esta seccion, resulta mas interesante
controlar x y y en lugar de φ y θ. Con el fin de lograr este ultimo control, se genera una
senal de error entre la posicion deseada en x y su valor real (al igual que para y), ese
error pasa por un lazo PD que se vuelve una referencia variable para el lazo del control
de θ (φ en el caso de x), lo que permite el control indirecto de las variables de interes.
El esquema de simulacion de este tipo de control multilazo con referencia variable se
muestra en la Figura 4-13.
564 Fallos en los rotores
(simulaciones y analisis)
Los controladores PD en cada lazo, definidos por el nombre debajo del bloque, fueron
disenados para lograr alcanzar el valor deseado en un tiempo de establecimiento cercano
a los 10 segundos. Las definiciones de las constantes de cada controlador se presentan
en la Tabla 4-1.
Figura 4-13: Lazos PD para el control del hexacoptero.
Tabla 4-1: Definicion de constantes para los bloques PD de la Figura 4-13.Bloque Kp Kd
PD x -0.0524 -0.1047
PD y 0.0524 0.1047
PD z -1.25 -1.25
PD phi 0.014 0.014
PD theta -0.014 -0.014
PD psi -0.02 -0.02
El resultado del controlador es bastante eficiente, logra en 10 segundos desplazar el
hexacoptero en los tres ejes y con una orientacion particular, como se muestra en la
Figura 4-14. Allı se evidencia que en las variables z y ψ se llega al valor de estado
estacionario a los 5 segundos, mientras que en x y y se oscila al rededor de el valor
deseado. En la Figura 4-18 se muestra el comportamiento de la velocidad de los roto-
res, obtenidas con ayuda de la matriz A definida en (2-21).
4.3 Resultados 57
Con el fin de verificar las implicaciones del fallo de un rotor durante el vuelo y la
capacidad del control para superar el fallo, se procedio a hacer modificaciones en el rotor
1. Sin embargo, con este esquema, para ningun fallo el control es capaz de recuperar el
sistema.
0 2 4 6 8 10 12
t[x]
0
0.2
0.4
0.6
0.8
1
1.2
[m]
(a)
0 2 4 6 8 10 12
t[x]
0
0.2
0.4
0.6
0.8
1
1.2
[m]
(b)
0 2 4 6 8 10 12
t[x]
-1.2
-1
-0.8
-0.6
-0.4
-0.2
0
[m]
(c)
0 2 4 6 8 10 12
t[x]
-1
-0.5
0
0.5
1
[gra
do
s]
(d)
0 2 4 6 8 10 12
t[x]
-1
-0.5
0
0.5
1
[gra
do
s]
(e)
0 2 4 6 8 10 12
t[x]
0
0.2
0.4
0.6
0.8
1
[gra
do
s]
(f)
Figura 4-14: Respuesta del control multilazos PD para el control del hexacoptero.
584 Fallos en los rotores
(simulaciones y analisis)
0 2 4 6 8 10 12
t[x]
450
500
550
600
[r.p
.m.]
1
0 2 4 6 8 10 12
t[x]
480
500
520
540
560
[r.p
.m.]
4
0 2 4 6 8 10 12
t[x]
480
500
520
540
560
580
600
[r.p
.m.]
2
0 2 4 6 8 10 12
t[x]
480
500
520
540
560
[r.p
.m.]
5
0 2 4 6 8 10 12
t[x]
450
500
550
600
[r.p
.m.]
3
0 2 4 6 8 10 12
t[x]
480
500
520
540
560
[r.p
.m.]
6
Figura 4-15: Velocidades del control multilazos PD para el control del hexacoptero.
Control por realimentacion de estados, con seguimiento por integrador
Para el diseno de este control se utilizo un modelo linealizado del hexacoptero definido
por (2-20) en el punto de hover x = 0 y u = 506,7266. El modelo lineal obtenido es
una representacion en variables de estado contınua de la forma:
4.3 Resultados 59
x = A x + B u,
y = C x,
donde para el hexacoptero en el punto de hover especıficamente tenemos
A =
0 0 0 0 0 0 1 0 0 0 0 0
0 0 0 0 0 0 0 1 0 0 0 0
0 0 0 0 0 0 0 0 1 0 0 0
0 0 0 0 0 0 0 0 0 1 0 0
0 0 0 0 0 0 0 0 0 0 1 0
0 0 0 0 0 0 0 0 0 0 0 1
0 0 0 0 −0,1712 0 0 0 0 0 0 0
0 0 0 0,1712 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0 0 0
,
B =
0 0 0 0 0 0
0 0 0 0 0 0
0 0 0 0 0 0
0 0 0 0 0 0
0 0 0 0 0 0
0 0 0 0 0 0
0 0 0 0 0 0
0 0 0 0 0 0
−0,0065 −0,0065 −0,0065 −0,0065 −0,0065 −0,0065
0 0,1050 0,1050 0 −0,1050 −0,1050
0,1399 0,0350 −0,0350 −0,1399 −0,0350 0,0350
0,0131 −0,0131 0,0131 −0,0131 0,0131 −0,0131
,
y
C =
1 0 0 0 0 0 0 0 0 0 0 0
0 1 0 0 0 0 0 0 0 0 0 0
0 0 1 0 0 0 0 0 0 0 0 0
0 0 0 0 0 1 0 0 0 0 0 0
,Para este sistema se diseno una realimentacion de estados con un esquema de segui-
miento por integrador, como se muestra en la Figura 4-16. Para el diseno de esta
realimentacion de estados, debe utilizarse un sistema ampliado, considerando como va-
riables de estado adicionales las variables de error [29]. Ası las cosas, la ecuacion de
604 Fallos en los rotores
(simulaciones y analisis)
estado del sistema ampliado tendra la forma[x
0
]=
[A 0
−C 0
] [x
0
]+
[B
0
]u,
para el que finalmente se disena la realimentacion de estados, ubicando todos los polosde lazo cerrado en el semiplano derecho, de forma tal que el sistema logre cada maniobraal rededor de los 10 segundos, para que los resultados sean comparables con las seccionesanteriores. Ası las cosas, la matriz de realimentacion K de estados es obtenida conayuda de Matlab y queda de la forma
K =
0,2652 0,0835 −0,0252 0,0433 −0,1345 0,0074 0,6897 0,2183 −0,0443 0,0188 −0,0546 0,0130
−0,0937 −0,2492 −0,0255 −0,1244 0,0498 −0,0157 −0,2460 −0,6425 −0,0445 −0,0504 0,0234 −0,0268
0,3999 0,2045 −0,0251 0,1098 −0,2117 0,0151 1,0540 0,5404 −0,0442 0,0510 −0,0943 0,0263
0,0408 −0,1432 −0,0254 −0,0737 −0,0273 −0,0080 0,1181 −0,3728 −0,0444 −0,0318 −0,0162 −0,0135
0,3998 0,1895 −0,0251 0,0941 −0,2116 0,0151 1,0538 0,4880 −0,0442 0,0373 −0,0942 0,0263
−0,0938 −0,2642 −0,0255 −0,1401 0,0499 −0,0157 −0,2463 −0,6949 −0,0445 −0,0640 0,0235 −0,0268
,
por 1000. Mientras que
Ki =
−37,9851 −11,9640 4,7804 −1,3955
13,4491 36,0134 4,8649 3,0584
−56,8409 −29,1139 4,7524 −2,8807
−5,3907 20,5741 4,8369 1,5734
−56,8249 −27,4034 4,7524 −2,8805
13,4651 37,7240 4,8649 3,0586
.
Figura 4-16: Esquema de realimentacion de estados para el hexacoptero.
Teniendo en cuenta que el modelo lineal del hexacoptero es un modelo incremental,
para que la realimentacion de estado previamente disenada funcione con el modelo
no lineal del hexacoptero, es necesario adicionarles a las variables de control el valor
constante de la velocidad de hover. Dado que el punto de equilibrio de las variables de
estado es cero, no hace falta modificar las variables de estado. El esquema de simulacion
con el modelo no lineal del hexacoptero se muestra en la figura 4-17.
4.3 Resultados 61
Figura 4-17: Esquema de realimentacion de estados para el hexacoptero.
Luego de probar el control no lineal se evidencia que funciona perfecto para z y ψ,
pero no funciona bien para x y y. Lo anterior muy seguramente debe tener que ver
con que pese a que x y y no intervienen en la dinamica del sistema directamente, al
igual que z y ψ, dependen del comportamiento de φ y θ, dado que la unica forma en la
que el hexacoptero puede lograr desplazarse en x y y es rotando al rededor de dichos
ejes y como la linealizacion se debe hacer para valores fijos de φ y θ, para ninguna
linealizacion se logra que la realimentacion de estados controle x y y del modelo no
lineal, que es el que vale la pena estudiar.
En vista de lo anterior probaremos lo que sucede con fallos para movimientos sobre el
eje z y alrededor del mismo. La Figura 4-18 muestra el hexacoptero describiendo un
ascenso con este esquema de control y la Figura 4-19 muestra el comportamiento de las
velocidades de los rotores. En ambos casos parece pasar lo mismo, dado que el control
se disena para cierto numero de entradas, si una de ellas falla en cualquier momento,
no es posible que el control retome su tarea, dado que el diseno no esta pensado para eso.
Situacion distinta cuando se disena el control para una o dos entradas menos en cual-
quier configuracion. En este ultimo caso, el control es capaz de cumplir su objetivo
con uno o dos rotores detenidos, porque se diseno para esas condiciones. Curiosamente
cuando se hacen estudios de controlabilidad sobre el modelo linealizado para el punto
de hover, se encuentra que tal y como se encontro a traves de la solucion del problema
de optimizacion, el sistema es completamente controlable a falta de uno o incluso 2
rotores, solo que los rotores no pueden tener ubicacion simetrica, pero incluso pueden
estar uno al lado del otro.
624 Fallos en los rotores
(simulaciones y analisis)
0 2 4 6 8 10 12
t[x]
0
0.2
0.4
0.6
0.8
1
[m]
(a)
0 2 4 6 8 10 12
t[x]
0
0.2
0.4
0.6
0.8
1
[m]
(b)
0 2 4 6 8 10 12
t[x]
-5
-4
-3
-2
-1
0
[m]
(c)
0 2 4 6 8 10 12
t[x]
-10
-5
0
5
10
[gra
do
s]
(d)
0 2 4 6 8 10 12
t[x]
-10
-5
0
5
10
[gra
do
s]
(e)
0 2 4 6 8 10 12
t[x]
-10
-5
0
5
10
[gra
do
s]
(f)
Figura 4-18: Respuesta del control por realimentacion de estados e integrador para el
hexacoptero.
4.3 Resultados 63
0 5 10 15 20
t[x]
500
505
510
515
520
525[r
.p.m
.]1
0 5 10 15 20
t[x]
500
505
510
515
520
525
[r.p
.m.]
4
0 5 10 15 20
t[x]
500
505
510
515
520
525
[r.p
.m.]
2
0 5 10 15 20
t[x]
500
505
510
515
520
525
[r.p
.m.]
5
0 5 10 15 20
t[x]
500
505
510
515
520
525
[r.p
.m.]
3
0 5 10 15 20
t[x]
500
505
510
515
520
525
[r.p
.m.]
6
Figura 4-19: Comortamiento de las velocidades de los rotores para el control por
realimentacion de estados e integrador para el hexacoptero.
Capıtulo 5
Conclusiones y observaciones
Los resultados muestran lo facil que es que el hexacoptero se mueva verticalmente,
porque los rotores estan actuando directamente sobre el eje z, sin embargo, los
movimientos longitudinales y laterales son mas difıciles de conseguir debido a la
falta de actuadores en los ejes x y y.
Se observa que para cada uno de los puntos finales deseados en las simulaciones,
los demas pesos varıan de una manera no proporcional, es decir en caso de requerir
un movimiento al doble de la distancia obtenida no basta con duplicar los valores
sobre los vectores γ y α debido que el sistema no se comporta de manera lineal.
Se observa en las simulaciones que al realizar los desplazamientos sobre los ejes
x y y se presenta una leve variacion en la posicion de z esto indica que hay
una necesidad implıcita de realizar un movimiento sobre el eje z para realizar
movimientos en los ejes x y y.
La falla puede en un rotor puede ser simulada poniendo en cero la columna del
rotor perteneciente sobre la matriz A lo que significarıa fısicamente que se realizan
todos los calculos como si existieran 6 rotores pero en su implementacion fueran
solo 5 o la cantidad de rotores activos.
Al realizar la prueba apagando un motor el resultado en posicion no varia, sin
embargo las velocidades calculadas varıan y como primer punto se observa un
apagado automatico al rotor opuesto en posicion y variacion de las r.p.m. para
mantener la posicion de hover.
Se observa que cuando se realizan pruebas con falla en rotores 1 y 2 el control
apagaba el 4 y 5, esto para guardar la simetrıa en la propulsion y mantener
la altura y estabilidad, por tanto aumentaba la velocidad de las helices para
compensar fZ al igual que en el caso de apagado del rotor 1.
66 5 Conclusiones y observaciones
Los resultados arrojados en las pruebas con falla en rotor 1 son iguales a las
que se realizaron apagando los rotores 1 y 4, debido a lo explicado anteriormente
donde el control apaga motores simetricos a la falla para guardar simetrıa de
propulsion.
Debido que con el control de optimizacion disenado a partir problema de ma-
niorabilidad no fue posible realizar pruebas con fallas en los rotores durante el
vuelo, por tanto se disenaron otras dos estrategias de control en lazo cerrado (PD
y seguimiento por integradores) las cuales no fueron exitosas en x y y por su
dependencia cruzada con los angulos φ y θ respectivamente.
Capıtulo 6
Anexos
6.1. Codigos del Matlab
6.1.1. Codigo General
1 %function costegeneralf
2
3 clear all
4 close all
5 clc
6
7 t=cputime;
8
9 global x t1 alpha xf fz TK TM TN p ci h
10
11 load hexacopter params
12
13 Tf=10; % tiempo final
14 h=Tf/100; % paso de tiempo
15
16 % El coste esta definido como
17 % J = int (F(x,u,t)) dt + Phi(x(f))
18 % donde
19 % F = gamma*(x(t)−xd)ˆ2 + beta*u(t)ˆ2
20 % y
21 % Phi = alpha*(x(Tf)−xf)ˆ222
23 alpha = [10 10 10 0 0 0 100 100 100 0 0 0]; % estado ...
final
24 beta = 0*[0 1 1 0]; % control
25 lambda = 1e−9; % paso de las ...
iteraciones
70 6 Anexos
26 gamma = [1 1 1 0 0 0 0 0 0 0 0 0]; % estado deseado
27 xd = [1 1 1 0 0 0 0 0 0 0 0 0]; % valor ...
estacionario deseado de x(t)
28 x0 = [0 0 0 0 0 0 0 0 0 0 0 0]; % estado inicial
29 xf = [1 1 1 0 0 0 0 0 0 0 0 0]; % estado final
30
31 t1=0:h:Tf; % vector de tiempo
32
33 fz=mass*gg+0*t1';
34 TK=0+0*t1';
35 TM=0+0*t1';
36 TN=0+0*t1';
37
38
39 % ...
−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−40 % NONLINEAR SYSTEM resolucion con ODE
41
42 tspan = [0,Tf];
43
44 [T,x] = ode45(@(t,x) hex state(t,x,t1,fz,TK,TM,TN),tspan,x0');
45
46 x = interp1(T,x,t1);
47
48
49 % ...
−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−50
51 % c?lculo del primer coste
52
53 a = length(x);
54 F = beta*[(fz−gg*mass).ˆ2 (TK).ˆ2 (TM).ˆ2 (TN).ˆ2]' + ...
gamma*(((x−ones(a,1)*xd)').ˆ2);55 q = simps v(F,h);
56 ci(2) = q + alpha*(((x(a,:)−xf)').ˆ2)57 ci(1) = 1.1*ci(2)
58
59 j=2;
60 %% iteraciones
61 while (ci(j−1)−ci(j))>(1e−06*ci(2))62 %
63 j=j+1;
64
65 % solucion ecuacion coestado a partir de x y u
66
67 pt2 = 2*alpha.*(x(length(x),:)−xf);
6.1 Codigos del Matlab 71
68 tspan = [Tf,0];
69
70 [t2,p] = ode45(@(t2,p) hex costate(t2,p,t1,x,gamma,xd),tspan,pt2);
71
72 p = interp1(t2,p,t1);
73
74 % calculo del gradiente a partir de x u y p
75
76 gradJ = gradientefinalp2(fz,TK,TM,TN,p,t1,beta);
77
78 % calculo del nuevo control
79
80 uk = [fz TK TM TN];
81
82 ukp1 = uk − lambda*gradJ;
83
84 fz = ukp1(:,1);
85 TK = ukp1(:,2);
86 TM = ukp1(:,3);
87 TN = ukp1(:,4);
88
89 % volvemos a resolver la dinamica para conocer el nuevo coste y ...
decidir si seguimos iterando o encontramos la solucion
90
91 tspan = [0,Tf];
92
93 [T,x] = ode45(@(t,x) hex state(t,x,t1,fz,TK,TM,TN),tspan,x0');
94
95 x = interp1(T,x,t1);
96
97 % calculo del coste con el nuevo control debido a lambda
98
99 a = length(x);
100 F = beta*[(fz−mass*gg).ˆ2 (TK).ˆ2 (TM).ˆ2 (TN).ˆ2]' + ...
gamma*(((x−ones(a,1)*xd)').ˆ2);101 q1 = simps v(F,h);
102 ci(j) = q1 + alpha*(((x(a,:)−xf)').ˆ2)103
104 end
105 cputime−t106 %% −−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−
6.1.2. Dinamica del hexacoptero
72 6 Anexos
1 function out = hex state(t,x,t1,fz,TK,TM,TN)
2
3 % aqu? se define la funcion f(x(t),u(t)), que es el segundo ...
miembro de
4 % la EDO dx(t)=f(x(t),u(t))
5
6 fz = interp1(t1,fz,t);
7 TK = interp1(t1,TK,t);
8 TM = interp1(t1,TM,t);
9 TN = interp1(t1,TN,t);
10
11 phi = pi/180*x(4);
12 th = pi/180*x(5);
13 psi = pi/180*x(6);
14 uu = x(7);
15 vv = x(8);
16 ww = x(9);
17 pp = pi/180*x(10);
18 qq = pi/180*x(11);
19 rr = pi/180*x(12);
20
21 load hexacopter params
22
23 out = [ %f1
24 uu*cos(th)*cos(psi) + vv*sin(phi)*sin(th)*cos(psi) + ...
25 − vv*cos(phi)*sin(psi) + ww*sin(phi)*sin(psi) + ...
26 ww*cos(phi)*sin(th)*cos(psi)
27 %f2
28 uu*cos(th)*sin(psi) + vv*cos(phi)*cos(psi) + ...
29 vv*sin(phi)*sin(th)*sin(psi) − ww*sin(phi)*cos(psi) + ...
30 ww*cos(phi)*sin(th)*sin(psi)
31 %f3
32 − uu*sin(th) + vv*sin(phi)*cos(th) + ww*cos(phi)*cos(th)
33 %f4
34 pp + qq*sin(phi)*tan(th) + rr*cos(phi)*tan(th)
35 %f5
36 qq*cos(phi) − rr*sin(phi)
37 %f6
38 qq*sin(phi)/cos(th) + rr*cos(phi)/cos(th)
39 %f7
40 vv*rr − ww*qq − gg*sin(th)
41 %f8
42 ww*pp − uu*rr + gg*sin(phi)*cos(th)
43 %f9
44 uu*qq − vv*pp + gg*cos(phi)*cos(th) + ...
6.1 Codigos del Matlab 73
45 − fz/mass
46 %f10
47 ((Iy−Iz)*qq*rr + TK)/Ix;
48 %f11
49 ((Iz−Ix)*pp*rr − TM)/Iy;
50 %f12
51 (Ix−Iy)/Iz*pp*qq − TN/Iz];
52
53 % xxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxx
6.1.3. Dinamica del coestado
1 function out1 = hex costate(t2,p,t1,x,gamma,xd)
2
3 % g1 define el segundo miembro de la ecuacion del coestado:
4 % dp(x)=−A'*p(x)−2*0.001*b;5 % donde A esta formada por las derivadas respecto de x(t) de cada ...
funcion del segundo miembro de la ecuacion de estado.
6
7 x = interp1(t1,x,t2);
8
9 load hexacopter params
10
11 f1x1 = 0;
12 f1x2 = 0;
13 f1x3 = 0;
14 f1x4 = ...
pi/180*(((x(:,8).*cos(1/180*pi*x(:,4))).*sin(pi/180*x(:,5)))....
15 *cos(pi/180*x(:,6)))+...
16 pi/180*((x(:,8).*sin(pi/180*x(:,4))).*sin(pi/180*x(:,6)))+...
17 pi/180*((x(:,9).*cos(pi/180*x(:,4))).*sin(pi/180*x(:,6)))+...
18 −pi/180*(((x(:,9).*sin(pi/180*x(:,4))).*sin(pi/180*x(:,5)))....19 *cos(pi/180*x(:,6)));
20 %plot(T,f1x4)
21 f1x5 = −(pi/180*x(:,7).*sin(pi/180*x(:,5))).*cos(pi/180*x(:,6))+...22 ((pi/180*x(:,8).*sin(pi/180*x(:,4))).*cos(pi/180*x(:,5)))....
23 *cos(pi/180*x(:,6))+...
24 ((pi/180*x(:,9).*cos(pi/180*x(:,4))).*cos(pi/180*x(:,5)))....
25 *cos(pi/180*x(:,6));
26 %plot(T,f1x5)
27 f1x6 = −(pi/180*x(:,7).*cos(pi/180*x(:,5))).*sin(pi/180*x(:,6))+...28 −((pi/180*x(:,8).*sin(pi/180*x(:,4))).*sin(pi/180*x(:,5)))....29 *sin(pi/180*x(:,6))+...
74 6 Anexos
30 −(pi/180*x(:,8).*cos(pi/180*x(:,4))).*cos(pi/180*x(:,6))+...31 (pi/180*x(:,9).*sin(pi/180*x(:,4))).*cos(pi/180*x(:,6))+...
32 −((pi/180*x(:,9).*cos(pi/180*x(:,4))).*sin(pi/180*x(:,5)))....33 *sin(pi/180*x(:,6));
34 %plot(T,f1x6)
35 f1x7 = cos(pi/180*x(:,5)).*cos(pi/180*x(:,6));
36 %plot(T,f1x7)
37 f1x8 = ...
((sin(pi/180*x(:,4))).*sin(pi/180*x(:,5))).*cos(pi/180*x(:,6))−...38 (cos(pi/180*x(:,4))).*sin(pi/180*x(:,6));
39 %%plot(T,f1x8)
40 f1x9 = ...
(sin(pi/180*x(:,4))).*sin(pi/180*x(:,6))+((cos(pi/180*x(:,4)))....
41 *sin(pi/180*x(:,5))).*cos(pi/180*x(:,6));
42 %plot(T,f1x9)
43 f1x10 = 0;
44 f1x11 = 0;
45 f1x12 = 0;
46
47 f2x1 = 0;
48 f2x2 = 0;
49 f2x3 = 0;
50 f2x4 = −(pi/180*x(:,8).*sin(pi/180*x(:,4))).*cos(pi/180*x(:,6))+...51 ((pi/180*x(:,8).*cos(pi/180*x(:,4))).*sin(pi/180*x(:,5)))....
52 *sin(pi/180*x(:,6))+...
53 −((pi/180*x(:,9).*sin(pi/180*x(:,4))).*sin(pi/180*x(:,5)))....54 *sin(pi/180*x(:,6))+...
55 −(pi/180*x(:,9).*cos(pi/180*x(:,4))).*cos(pi/180*x(:,6));56 %plot(T,f2x4)
57 f2x5 = −(pi/180*x(:,7).*sin(pi/180*x(:,5))).*sin(pi/180*x(:,6))+...58 ((pi/180*x(:,8).*sin(pi/180*x(:,4))).*cos(pi/180*x(:,5)))....
59 *sin(pi/180*x(:,6))+...
60 ((pi/180*x(:,9).*cos(pi/180*x(:,4))).*cos(pi/180*x(:,5)))....
61 *sin(pi/180*x(:,6));
62 %plot(T,f2x5)
63 f2x6 = (pi/180*x(:,7).*cos(pi/180*x(:,5))).*cos(pi/180*x(:,6))+...
64 −(pi/180*x(:,8).*cos(pi/180*x(:,4))).*sin(pi/180*x(:,6))+...65 ((pi/180*x(:,8).*sin(pi/180*x(:,4))).*sin(pi/180*x(:,5)))....
66 *cos(pi/180*x(:,6))+...
67 ((pi/180*x(:,9).*cos(pi/180*x(:,4))).*sin(pi/180*x(:,5)))....
68 *cos(pi/180*x(:,6))+...
69 (pi/180*x(:,9).*sin(pi/180*x(:,4))).*sin(pi/180*x(:,6));
70 %plot(T,f2x6)
71 f2x7 = (cos(pi/180*x(:,5))).*sin(pi/180*x(:,6));
72 %plot(T,f2x7)
6.1 Codigos del Matlab 75
73 f2x8 = ...
(cos(pi/180*x(:,4))).*cos(pi/180*x(:,6))+((sin(pi/180*x(:,4)))....
74 *sin(pi/180*x(:,5))).*sin(pi/180*x(:,6));
75 %plot(T,f2x8)
76 f2x9 = ...
((cos(pi/180*x(:,4))).*sin(pi/180*x(:,5))).*sin(pi/180*x(:,6))...
77 −(sin(pi/180*x(:,4))).*cos(pi/180*x(:,6));78 %plot(T,f2x9)
79 f2x10 = 0;
80 f2x11 = 0;
81 f2x12 = 0;
82
83 f3x1 = 0;
84 f3x2 = 0;
85 f3x3 = 0;
86 f3x4 = (pi/180*x(:,8).*cos(pi/180*x(:,4))).*cos(pi/180*x(:,5))+...
87 −(pi/180*x(:,9).*sin(pi/180*x(:,4))).*cos(pi/180*x(:,5));88 %plot(T,f3x4)
89 f3x5 = −pi/180*x(:,7).*cos(pi/180*x(:,5))+...90 −(pi/180*x(:,8).*sin(pi/180*x(:,4))).*sin(pi/180*x(:,5))+...91 −(pi/180*x(:,9).*cos(pi/180*x(:,4))).*sin(pi/180*x(:,5));92 %plot(T,f3x5)
93 f3x6 = 0;
94 f3x7 = −sin(pi/180*x(:,5));95 %plot(T,f3x7)
96 f3x8 = (sin(pi/180*x(:,4))).*cos(pi/180*x(:,5));
97 %plot(T,f3x8)
98 f3x9 = (cos(pi/180*x(:,4))).*cos(pi/180*x(:,5));
99 %plot(T,f3x9)
100 f3x10 = 0;
101 f3x11 = 0;
102 f3x12 = 0;
103
104 f4x1 = 0;
105 f4x2 = 0;
106 f4x3 = 0;
107 f4x4 = (pi/180*x(:,11).*cos(pi/180*x(:,4))).*tan(pi/180*x(:,5))+...
108 −(pi/180*x(:,12).*sin(pi/180*x(:,4))).*tan(pi/180*x(:,5));109 %plot(T,f4x4)
110 f4x5 = ...
(pi/180*x(:,11).*sin(pi/180*x(:,4))).*(1+(tan(pi/180*x(:,5))).ˆ2)+...
111 (pi/180*x(:,12).*cos(pi/180*x(:,4))).*(1+(tan(pi/180*x(:,5))).ˆ2);
112 %plot(T,f4x5)
113 f4x6 = 0;
114 f4x7 = 0;
115 f4x8 = 0;
76 6 Anexos
116 f4x9 = 0;
117 f4x10 = 1;
118 f4x11 = (sin(pi/180*x(:,4))).*tan(pi/180*x(:,5));
119 %plot(T,f4x11)
120 f4x12 = (cos(pi/180*x(:,4))).*tan(pi/180*x(:,5));
121 %plot(T,f4x12)
122
123 f5x1 = 0;
124 f5x2 = 0;
125 f5x3 = 0;
126 f5x4 = −pi/180*x(:,11).*sin(pi/180*x(:,4)) − ...
pi/180*x(:,12).*cos(pi/180*x(:,4));
127 %plot(T,f5x4)
128 f5x5 = 0;
129 f5x6 = 0;
130 f5x7 = 0;
131 f5x8 = 0;
132 f5x9 = 0;
133 f5x10 = 0;
134 f5x11 = cos(pi/180*x(:,4));
135 %plot(T,f5x11)
136 f5x12 = −sin(pi/180*x(:,4));137 %plot(T,f5x12)
138
139 f6x1 = 0;
140 f6x2 = 0;
141 f6x3 = 0;
142 f6x4 = pi/180*x(:,11).*cos(pi/180*x(:,4))./cos(pi/180*x(:,5))+...
143 −pi/180*x(:,12).*sin(pi/180*x(:,4))./cos(pi/180*x(:,5));144 %plot(T,f6x4)
145 f6x5 = pi/180*x(:,11).*sin(pi/180*x(:,4))./(cos(pi/180*x(:,5))).ˆ2....
146 *sin(pi/180*x(:,5))+...
147 pi/180*x(:,12).*cos(pi/180*x(:,4))./(cos(pi/180*x(:,5))).ˆ2....
148 *sin(pi/180*x(:,5));
149 %plot(T,f6x5)
150 f6x6 = 0;
151 f6x7 = 0;
152 f6x8 = 0;
153 f6x9 = 0;
154 f6x10 = 0;
155 f6x11 = (sin(pi/180*x(:,4)))./cos(pi/180*x(:,5));
156 %plot(T,f6x11)
157 f6x12 = (cos(pi/180*x(:,4)))./cos(pi/180*x(:,5));
158 %plot(T,f6x12)
159
160 f7x1 = 0;
6.1 Codigos del Matlab 77
161 f7x2 = 0;
162 f7x3 = 0;
163 f7x4 = 0;
164 f7x5 = −pi/180*gg*cos(pi/180*x(:,5));165 f7x6 = 0;
166 f7x7 = 0;
167 f7x8 = (pi/180)*x(:,12);
168 f7x9 = −(pi/180)*x(:,11);169 f7x10 = 0;
170 f7x11 = −(pi/180)*x(:,9);171 f7x12 = (pi/180)*x(:,8);
172
173 f8x1 = 0;
174 f8x2 = 0;
175 f8x3 = 0;
176 f8x4 = (pi/180*gg*cos(pi/180*x(:,4))).*cos(pi/180*x(:,5));
177 f8x5 = −(pi/180*gg*sin(pi/180*x(:,4))).*sin(pi/180*x(:,5));178 f8x6 = 0;
179 f8x7 = −(pi/180)*x(:,12);180 f8x8 = 0;
181 f8x9 = (pi/180)*x(:,10);
182 f8x10 = (pi/180)*x(:,9);
183 f8x11 = 0;
184 f8x12 = −(pi/180)*x(:,7);185
186 f9x1 = 0;
187 f9x2 = 0;
188 f9x3 = 0;
189 f9x4 = −(pi/180*gg*sin(pi/180*x(:,4))).*cos(pi/180*x(:,5));190 f9x5 = −(pi/180*gg*cos(pi/180*x(:,4))).*sin(pi/180*x(:,5));191 f9x6 = 0;
192 f9x7 = (pi/180)*x(:,11);
193 f9x8 = −(pi/180)*x(:,10);194 f9x9 = 0;
195 f9x10 = −(pi/180)*x(:,8);196 f9x11 = (pi/180)*x(:,7);
197 f9x12 = 0;
198
199
200 f10x1 = 0;
201 f10x2 = 0;
202 f10x3 = 0;
203 f10x4 = 0;
204 f10x5 = 0;
205 f10x6 = 0;
206 f10x7 = 0;
78 6 Anexos
207 f10x8 = 0;
208 f10x9 = 0;
209 f10x10 = 0;
210 f10x11 = (pi/180)ˆ2*x(:,12)*(Iy−Iz)/Ix;211 f10x12 = (pi/180)ˆ2*x(:,11)*(Iy−Iz)/Ix;212
213 f11x1 = 0;
214 f11x2 = 0;
215 f11x3 = 0;
216 f11x4 = 0;
217 f11x5 = 0;
218 f11x6 = 0;
219 f11x7 = 0;
220 f11x8 = 0;
221 f11x9 = 0;
222 f11x10 = (pi/180)ˆ2*x(:,12)*(Iz−Ix)/Iy;223 f11x11 = 0;
224 f11x12 = (pi/180)ˆ2*x(:,10)*(Iz−Ix)/Iy;225
226 f12x1 = 0;
227 f12x2 = 0;
228 f12x3 = 0;
229 f12x4 = 0;
230 f12x5 = 0;
231 f12x6 = 0;
232 f12x7 = 0;
233 f12x8 = 0;
234 f12x9 = 0;
235 f12x10 = (pi/180)ˆ2*x(:,11)*(Ix−Iy)/Iz;236 f12x11 = (pi/180)ˆ2*x(:,10)*(Ix−Iy)/Iz;237 f12x12 = 0;
238
239
240 A = [f1x1 f1x2 f1x3 f1x4 f1x5 f1x6 f1x7 f1x8 f1x9 f1x10 ...
f1x11 f1x12
241 f2x1 f2x2 f2x3 f2x4 f2x5 f2x6 f2x7 f2x8 f2x9 f2x10 ...
f2x11 f2x12
242 f3x1 f3x2 f3x3 f3x4 f3x5 f3x6 f3x7 f3x8 f3x9 f3x10 ...
f3x11 f3x12
243 f4x1 f4x2 f4x3 f4x4 f4x5 f4x6 f4x7 f4x8 f4x9 f4x10 ...
f4x11 f4x12
244 f5x1 f5x2 f5x3 f5x4 f5x5 f5x6 f5x7 f5x8 f5x9 f5x10 ...
f5x11 f5x12
245 f6x1 f6x2 f6x3 f6x4 f6x5 f6x6 f6x7 f6x8 f6x9 f6x10 ...
f6x11 f6x12
6.1 Codigos del Matlab 79
246 f7x1 f7x2 f7x3 f7x4 f7x5 f7x6 f7x7 f7x8 f7x9 f7x10 ...
f7x11 f7x12
247 f8x1 f8x2 f8x3 f8x4 f8x5 f8x6 f8x7 f8x8 f8x9 f8x10 ...
f8x11 f8x12
248 f9x1 f9x2 f9x3 f9x4 f9x5 f9x6 f9x7 f9x8 f9x9 f9x10 ...
f9x11 f9x12
249 f10x1 f10x2 f10x3 f10x4 f10x5 f10x6 f10x7 f10x8 f10x9 f10x10 ...
f10x11 f10x12
250 f11x1 f11x2 f11x3 f11x4 f11x5 f11x6 f11x7 f11x8 f11x9 f11x10 ...
f11x11 f11x12
251 f12x1 f12x2 f12x3 f12x4 f12x5 f12x6 f12x7 f12x8 f12x9 f12x10 ...
f12x11 f12x12];
252
253 out1 = −A'*p−(2*diag(gamma)*(x−ones(length(t2),1)*xd)');%−2*0.01*b;254
255 % xxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxx
6.1.4. Calculo del Gradiente
1 function gradJ = gradientefinalp2(fz,TK,TM,TN,p,t1,beta)
2
3 load hexacopter params
4
5 gradJ = zeros(length(t1),4);
6 for i = 1:length(p);
7 A2 = [zeros(8,4)
8 −1/mass 0 0 0
9 0 1/Ix 0 0
10 0 0 −1/Iy 0
11 0 0 0 −1/Iz];12
13 %aux=inv(M)*(A2)
14 A3 = 2*beta*[fz(i)−mass*gg TK(i) TM(i) TN(i)]';
15 gradJ(i,:) = A3 + p(i,:)*(A2);
16 end
17
18 % xxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxx
Bibliografıa
[1] G. Calvo, Perfiles IDS Sistemas no tripulados, Febrero 2014, Madrid.
[2] M. F. Santos, L. M. Honorio, E. B. Costa, E. J. Oliveiray, J. P. Portella Guedes
Visconti, Active Fault-Tolerant Control Applied to a Hexacopter under Propulsion
System Failures, 2015 19th International Conference on System Theory, Control
and Computing (ICSTCC), October 14-16, Cheile Gradistei, Romania.
[3] A. Freddi, S. Longhi, A. Monteriu and M. Prist, Actuator Fault Detection and
Isolation System for an Hexacopter, Mechatronic and Embedded Systems and Ap-
plications (MESA), 2014 IEEE/ASME 10th International Conference.
[4] J. Su, J. He, P. Cheng y J. Chen, Actuator Fault Diagnosis of a Hexacopter: A
Nonlinear Analytical Redundancy Approach, 2017 25th Mediterranean Conference
on Control and Automation (MED) July 3-6, 2017. Valletta, Malta.
[5] J. Beltran, Origen y desarrollo de los drones [online], Espana, Universidad de Va-
lencia, 2015, Disponible en http://drones.uv.es/origen-y-desarrollo-de-los-drones/
[6] J. Sullivan, Revolution or Evolution? The Rise of UAVs, Sibley School of Mecha-
nical and Aerospace Engineering, Cornell University, Ithaca, NY 14853.
[7] H. Lim, J. Park, D. Lee, H. J. Kim, Build Your Own Quadrotor, Open-Source
Projects on Unmanned Aerial Vehicles, IEE Robotics & Automation Magazine,
September 2012.
[8] Fahlstrom y T. Gleason, Introduction to UAV systems, Wiley, 2012.
[9] Amazon Prime Air [online], Amazon, 2016, Disponible en
https://www.amazon.com/Amazon-Prime-Air/b?node=8037720011
[10] P. Lara, Diseno de interfaces graficas y de algoritmos de localizacion de un qua-
drotor, Dep. de Ingenierıa de Sistemas y Automatica, Escuela Tecnica Superior
de Ingenierıa Universidad de Sevilla, September 2015.
82 BIBLIOGRAFIA
[11] J. Jenkins and M. Jenkins, LaTrax Alias Quad-Rotor [online], LATRAX, Dispo-
nible en https://latrax.com/products/alias?t=overview
[12] R. Stephens, Horizon Hobby SAFETM [online], disponible en
http://www.flysaferc.com/#homeSection
[13] R. Baranek, F. Solc, Modelling and Control of a Hexa-copter, Carpathian Control
Conference (ICCC), 2012 13th International.
[14] A. Alaimo, V. Artale, C. Milazzo, A. Ricciardello, L. Trefiletti, Mathematical Mo-
deling and Control of a Hexacopter, 2013 International Conference on Unmanned
Aircraft Systems (ICUAS).
[15] Camilo Morales; Diana Ovalle; Alain Gauthier, Hexacopter maneuverability ca-
pability: An optimal control approach. 7th International Conference on Modeling,
Simulation, and Applied Optimization (ICMSAO), United Arab Emirates, 2017.
[16] E. Reyes-Valeria, R. Enriquez-Caldera, S. Camacho-Lara, J. Guichard, LQR con-
trol for a quadrotor using unit quaternions: Modeling and simulation, 23rd inter-
national conference on electronics, communications and computing, (CONIELE-
COMP 2013), holula, Puebla, Mexico (2013), pp. 172-178.
[17] E. Fresk, G. Nikolakopoulos, Full quaternion based attitude control for a quadrotor,
2013 12th European Control Conference (ECC 2013), Zurich, Switzerland (2013),
pp. 3864-3869.
[18] M. Ranjbaran, K. Khorasani, Fault Recovery of an Under-Actuated Quadrotor
Aerial Vehicle, in: 49th Proc. Conference on Decision and Control, Atlanta, 2010,
pp. 4385-4392.
[19] Z.T. Dydel, A.M. Annaswamy, E. Lavretsky, Adaptive Control of Quadrotor UAVs:
A Design Trade Study With Flight Evaluations, IEEE Transactions on Control
Systems Technology. 21(4) (2013) 1400-1406.
[20] F. Goßmann, G. P. Falconı and F. Svaricek, Fault Tolerant Control Approach for
Passive Handling of Actuator Degradation of a Hexacopter using an LPV State-
Feedback Controller with Partly-Measured Parameters, 2017 25th Mediterranean
Conference on Control and Automation (MED) July 3-6, 2017. Valletta, Malta.
[21] B. Wang, and Y. Zhang, An Adaptive Fault-Tolerant Sliding Mode Control Allo-
cation Scheme for Multirotor Helicopter Subject to Simultaneous Actuator Faults,
IEEE Transactions on industrial electronics.
BIBLIOGRAFIA 83
[22] H. K. Khalil, Nonlinear Systems 2nd Edition, Prentice Hall,Michigan State Uni-
versity, 1996.
[23] K. Shimizu, K. Otsuka y T. Yamaguchi, Direct Gradient Descent Control Nonlie-
nar Systems, Procedings of the American Control Conference, May 8-2002.
[24] K. Shimizu, K. Otsuka y J. Naiborhu, Improved Direct Gradient Descent Control
Nonlienar Systems, 1999 European Control Conference (ECC), 31 August - 3
September 1999, Karlsruhe, Germany.
[25] K. Shimizu, K. Otsuka, Performance Improvement of Direct Gradient Descent
Control Nonlienar Systems, International Conference on control applications
Kohala Coast-Island of Hawai 22-27 August 1999.
[26] D. Ballesteros, A. Gaona, L. Pedraza. Comparacion por simplicidad de metodos de
aprendizaje en estimacion de funciones. Vision Electronica. Universidad Distrital
Francisco Jose de Caldas. Octubre 2010.
[27] C. Vega, R. Castano, Control optimo inverso para sistemas no lineales en tiempo
continuo. Respuestas. 2014; 19(1):13-18.
[28] D.M. Ovalle, J. Garcıa, and F. Periago. Analysis and numerical simulation of
a nonlinear mathematical model for testing the manoeu-vrability capabilities of a
submarine. NonLinear Analysis: Real world Applications. Elsevier, 2010.
[29] K. Ogata, Ingenierıa de Control Moderna, 3ra Ed, Universidad de Minnesota,
Editorial Prentice Hall 1998.