TRABAJO DE CONTROL II 2020 1

Helicopterodedosgradosdelibertad ´
Vilca Choque Marco Antonio 110015

Resumen—En este trabajo, el modelo dinámico y el prototipo.

Se presenta el control para un helicóptero de dos grados de
libertad. El modelo del sistema es no lineal y de bucle abierto
inestable, que se obtiene utilizando las ecuaciones de movimiento
de ”Euler-Lagrange”. El modelo lineal resultante es función de
un punto operativo. El diseño de controlador consta de dos
bucles de control, uno para el angulo de elevación (Pitch) y el
otro para el movimiento de rotación (Yaw). Para el sistema de
control, se utilizan técnicas de re-alimentación de estado y se
diseña un observador para estimar todos los estados. Se describen
las ecuaciones del esquema de control de re-alimentación con
observador y, a partir de ella, se diseña un controlador central.
Este es un control basado en el observador. El controlador logra
estabilidad y buen rendimiento en diferentes puntos de operación.
Figura 1. Dinámica del helicóptero 2DOF

I. I NTRODUCCI ÓN

E n las aplicaciones de control de procesos quı́micos,

robóticos, control de aeronaves, entre otros, se requie-
re controlar diferentes variables fı́sicas, como temperatura,
presión, flujo, nivel, etc., para llevar acabo un propósito
determinado. Se utiliza el concepto de re-alimentación, que
consiste en comparar el valor deseado de una variable fı́sica
frente al valor real, y en función del error, ejecutar una acción
de control determinada. En el diseño de sistemas de control
de re-alimentación, es necesario encontrar un controlador
para mantener u obtener la estabilidad para disminuir, las
perturbaciones que son señaladas que tienden a cambiar el
valor de una señal o variable de salida. Los helicópteros DOF Figura 2. Cinemática del helicóptero
constan de dos rotores en cuadratura, que giran alrededor de
los ejes Pitch y Yaw. Para encontrar el modelo dinámico
de la planta, se utilizo el método de ”Euler-Lagrange” que se define como positivo para rotaciones en el sentido de las
describen las ecuaciones y las fuerzas implicadas. El modelo agujas del reloj.
resultante es no lineal, a partir del cual un modelo lineal esta La fuerza gravitatoria 0 fg0 genera un par en el eje Pitch que
en función del angulo de inclinación. El sistema de control en hace que la punta del helicóptero caiga por que el centro de
circuito cerrado utiliza un controlador LQR; el control LQR masa no esta en el eje sino a una distancia 0 Imc
0
a lo largo del
0 0
es una técnica de re-alimentación de estado que encuentra una fuselaje y a una altura h debajo de esto. La figura 2 muestra
solución optima a un problema de minimización. la cinemática directa del helicóptero.
El documento se divide en una parte de modelado, otra Para encontrar las ecuaciones de movimiento se utilizo el
de diseño del sistema de control, y finalmente muestra los método de ’Euler-Lagrange’, que se basa en el análisis de
resultados obtenidos en la planta. la energı́a cinética y potencial del sistema.
La posición del centro de masa del sistema en relación con
el sistema de coordenadas fijo (eje) se encuentra utilizando
II. M ODELADO
matrices de transformación homogéneas (Rotzo ,ψ )(Roty1 ,θ )

E l helicóptero con dos grados de libertad gira alrededor del

eje de inclinación (Pitch) por el angulo ’θ’ y alrededor
del eje de derrape (Yaw) por el angulo ’ψ’ como se muestra en
y se obtiene la posición cartesiana del centro de masas.

la figura 1. El angulo de inclinación se define como positivo Xmc =(Imc cos θ + h sin θ) cos ψ (1)
cuando la punta del helicóptero se eleva y el angulo de derrape Ymc =(−Imc cos θ − h sin θ) sin ψ (2)
Universidad Nacional de San Antonio Abad del Cusco. Zmc =Imc − h cos θ (3)
Ingenierı́a Electrónica.
Teorı́a de Control II
2 TRABAJO DE CONTROL II 2020

Energı́a cinética potencial.- La energı́a potencial debido a las expresiones de Euler-Lagrange da como resultado de
la gravedad es: movimiento no lineales.
P = mheli .g.Zmc = mheli .g.(Imc sin θ − h cos θ) (4) 2
(Jeqp + Mheli (lmc + h2 ))θ̈+
Donde: mheli es la masa total del helicóptero La energı́a 2
− h2 )
 
sin 2θ(lmc
cinética total. Mheli − lmc h cos 2θ ψ˙2
2
T = Tr,p + Tr,y + Tt (5)
+Mheli g(lmc cos θ + h sin θ) + Bp θ̇
Es la adición de energı́a cinética rotacional que actúa sobre
Pitch (Tr,p ) y Yaw (Tr,y ) con energı́a cinética traslacional = (Kpp Vmp + Fcpp ) + (Kpy Vmy + Fcpy )
generada por el movimiento del centro de masa Tt . La energı́a
cinética rotacional en Pitch y Yaw es:
Jeqp + Mheli cos2 θ(lmc − h2 ) + lmc h sin 2θ + h2 ψ¨2
  
1 1
Tr,p = Jeqp θ˙2 Tr,y = Jeqy ψ̇ (6)
2 2
+Mheli sin 2θ(h2 − lmc
2
 
) + 2lmc h cos 2θ θ̇ψ̇ + By ψ̇
Donde Jeqp y Jeqy son los momentos de inercia equivalentes
respectivamente a Pitch y Yaw. La energı́a cinética traslacio- = (Kyp Vmp + Fcyp ) cos θ + (Kyy Vmy + Fcyy ) cos θ
nal es:
1 Modelo lineal de espacio de estados
2 2 2
Tt = mheli (Xmc + Ymc + Zmc ) (7) T  T
2

X = θ ψ θ̇ ψ̇ = x1 x2 x3 x4
Ecuaciones de movimiento no lineales.- Las ecuaciones
de movimiento de Euler-Lagrange se define como: y como señales de control
 T  T
∂ ∂ ∂ U = Vmp Vmy = u1 u2
L− = Q1
∂t ∂ q˙1 ∂q1
∂ ∂ ∂ Linealización alrededor del equilibrio inducido por x1 = θ
L− = Q2 (8) constante, obteniendo el modelo lineal parametrizado por θ
∂t ∂ q˙2 ∂q2
L = cos2 θ(lmc2
− h2 ) + lmc h sin 2θ + h2
 
Donde ”L” es la variable de Lagrange, que corresponde a la
diferencia entre la energı́a cinética y la potencial del sistema.
W = lmc sin θ − h cos θ
L=T −P (9) 2
N = lmc + h2
Las coordenadas generalizadas son:  
0 0 1 0
q = [q1 q2 q˙1 q˙2 ] = [θ ψ θ̇ ψ̇]  0 0 0 1 
A(θ) =  Mheli g(W ) −Bp
 
y las fuerzas generalizadas son:  Jeqp +Mheli (N ) 0 Jeqp +Mheli (N ) 0 

−By
0 0 0 Jeqy +Mheli L
Q1 = τp (Vmp , Vmy ) − Bp θ̇ Q2 = τy (Vmp , Vmy ) − By ψ̇
0 0
 
Antes de la ecuación incluye la fricción viscosa rotacional
actuando sobre los ejes Pitch y Yaw; Bp y By . Los pares  0 0 
generados en los ejes Pitch y Yaw son función de los voltajes
 Kpp Kyp 
B(θ) = 
 
 Jeqp + Mheli (N ) Jeqp + Mheli (N ) 

aplicados a los motores.  Kyp cos θ Kyy cos θ 
τp (Vmp , Vmy ) = τpp + τpy Jeqy + Mheli L Jeqy + Mheli L
τp (Vmp , Vmy ) = (Kpp Vmp + Fcpp ) + (Kpy Vmy + Fcpy )
 
1 0 0 0
C=
τy (Vmp , Vmy ) = τyp + τyy 0 1 0 0
sustituyendo valores para θ = 0 se obtiene
= (Kyp Vmp + Fcyp ) cos θ + (Kyy Vmy + Fcyy ) cos θ  
x1 x2 x3 x4
Donde ”Vmp ”es el voltaje del motor en Pitch de entrada  0 0 1 0 
y Vmy es el voltaje de entrada del motor en Yaw.los pares    
A= 0 0 0 1  ∗ θ ψ θ̇ ψ̇

que actúan sobre los ejes Pitch y Yaw están acoplados. Las −0,27 0 −0,35 0 
funciones del par están representadas por τpp , τyy , τyp y τpp
0 0 0 −22,74
Kpp mientras Kyy , Kyp y Kpy son constantes de voltaje
a torque de los actuadores de motor hélice Pitch y Yaw
 
x1 x2
se encontraron experimentalmente. Fcpp , Fcpy , Fcyp y Fcyy  0 0 
son los términos constantes que aparecen en la regrecion
   
B= 0 0  ∗ Vmp Vmy

lineal debido a la fricción de coulomb, es el voltaje necesario 0,712 0,051
para que el sistema inicie movimiento. La evaluación de 0,056 0,499
CONTROL II: HELICOPTERO DOS GRADOS DE LIBERTAD 3

D = [0];
R = ss(A, B, C, D);
stepplot(R);

III. ANALISIS EN FRECUENCIA

Para empezar a hacer el analisis en frecuencia(bode,nyquist,

nycholls) de la planta sin compesar,Debemos saber cual es la
funcion de transferencia,esta ultima se halla mediante el uso
de las matrices de estado vistas con anterioridad conocida
como matrizde funcion de transferencia.
Figura 3. Valores de prueba

Figura 4. Grafica de la planta Simulink

Ahora para hallar la matriz inversa,se debe usar la

formula siguiente.

Figura 5. Grafica de la respuesta y entrada

Figura 6. Grafica de la respuesta y entrada

 
x1 x2 x3 x4  
C =1 0 0 0  ∗ VP itch VY aw
0 1 0 0

clear all;
clc;
sistemaejemploE2GL
A = [0010; 0001; −0,270 − 0,350; 000 − 22,74, ];
B = [00; 00; 0,7120,051; 0,0560,499, ];
C = [1000; 0100];
 CODIGO MATLAB:

clc
clear all:
Determinante=d nG = [0.763];
nF = [1 0.35 0.2];
G = tf(nG,nF);
bode(G) Graficas.

Ahora reemplazamos en la ecuacion de funcion

de transferencia,para hallar su valor.

ANALISIS EN NICHOLS.-Nos permite graficar la

respuesta del sistema mediante una sola curva,donde los
ejes estan dados como magnitud vs fase
ANALISIS EN BODE.-Nos permite graficar la
respuesta del sistema en magnitud y fase frente La grafica utiliza puntos en magnitud y fase para su
a una entrda senoidal. trazado,el cual nos permite analizar la estabilidad.

G1P(S): CODIGO MATLAB:

CONSTANTE: clc
clear all:
M[db]=20log(2.82)=9.02db nG = [0.763];
nF = [1 0.35 0.2];
POLO COMPLEJO CONJUGADO: G = tf(nG,nF);
nichols(nG,nF);

Grafica.

Wn=0.52rad/s
E=0.34
 G2Y(S):

CONSTANTE:

M[db]=20log(0.024)=-32.39db

POLO eje real:

20log(1/jw)

20log(1/0.044jw+1)

CODIGO MATLAB

clc
nG1 = [0.555];
nF1 = [1 22.74 0];
G2 = tf(nG1,nF1);
bode(G2)

MF=180-154=26 grados
MG=infinito Grafica.

WcF=0.97rad/s
WcG=infinito

ANALISIS EN NYQUIST.-Tambien conocido como el

trazo polar,se emplean coordenadas polares,conforme
W varia desde 0 al infinito.

La grafica muestra el lugar geometrico de los vectores

..

CODIGO MATLAB:

clc
nG=[0.763];
dG = [1 0.35 0.27];
F=tf(nG,dG);
H=nyquistplot(F);
title('Diagrama de Nyquist')
ylabel('Eje imaginario')
xlabel('Eje real')
Grafica.

G2Y(S): NICHOLS

M[db]=20log(0.024)=-32.39db

POLO eje real:

20log(1/jw)

20log(1/0.044jw+1)

NICHOLS;
 CODIGO MATLAB: Grafica.
clc
nG1 = [0.555];
nF1 = [1 22.74 0];
G2 = tf(nG1,nF1);
bode(G2)
Grafica.

MF=180-90=90 grados IV. CONCLUSION

MG=infinito
PSe puede apreciar que el analisis en frecuencia es muy
WcF=0.015rad/s importante para determinar la estabilidad del sistema,en este
WcG=infinito
caso se vio para los 3 casos el sistema es estable por el criterio de
G2Y(S): NYQUIST estabilidad absoluta ,sin embargo se hace un analisis en
estabilidad relativa ,para ver si se procede a disenar un PID,
M[db]=20log(0.024)=-32.39db compensador ,etc,esto para hacer que el sistema sea estable ,esto
implicaria una reduccion considerable en el tiempo que aprecia
POLO eje real: el transitorio(oscilaciones).

20log(1/jw)

20log(1/0.044jw+1)

NICHOLS;

CODIGO MATLAB:

clc
nG=[0.555];
dG = [1 22.74 0];
F=tf(nG,dG);
H=nyquistplot(F);
title('Diagrama de Nyquist')
ylabel('Eje imaginario')
xlabel('Eje real')