4.
Hacer el análisis de estabilidad del sistema Masa-Resorte-Amortiguador
(MRA) (Fig. 1), con los siguientes valores de los parámetros: M= 7, fv= 4, K=
10. Hacer el análisis de la respuesta transitoria del sistema MRA para la
entrada escalón unitario. Encuentre los errores en estado estacionario para
las entradas f(t)=1, f(t)=4t y f(t)=5t2 para el sistema MRA. Favor de usar
Matlab/Simulink para probar todos sus resultados. (50 puntos)
Fig. 1 Sistema Masa-Resorte-Amortiguador (MRA)
Para comenzar sacamos el modelado matemático del sistema mostrado en la
Fig.1 con ecuaciones diferenciales utilizando el método de Newton, la sumatoria
de fuerzas de este sistema es el mostrado a continuación:
M ẍ =f ( t )−fv ẋ −Kx
M ẍ + fv ẋ+ Kx=f (t )
Una vez obtenido el modelado del sistema con ecuaciones diferenciales
procederemos a obtener su función de transferencia, que se describe como la
transformada de Laplace de la salida sobre la transformada de Laplace de la
entrada cuando las condiciones iniciales son iguales a 0, ya que con la función de
transferencia podemos encontrar los polos de nuestro sistema que nos sirven para
realizar nuestro análisis de estabilidad.
L { M ẍ +fv ẋ + Kx }=L{f ( t ) }
Podemos distribuir la transformada de Laplace a cada termino para facilitar su
realización y seguidamente sacar las constantes multiplicándolas por la
transformada de la place de nuestras funciones, haciendo uso de la propiedad de
linealidad que nos dice que la transformada de Laplace es lineal, lo que significa
que la transformada de una combinación lineal de funciones es igual a la misma
combinación lineal de sus transformadas.
M L { ẍ }+ fv L{ ẋ }+ K L{x }=L{f ( t ) }
Haciendo uso de las fórmulas de la transformada de Laplace, especialmente de la
derivada de una función, nos queda:
M [ S 2 X ( s ) −SX ( 0 )−X ' ( 0 ) ] +fv [ SX ( s )− X ( 0 ) ] + KX(s)=f (s)
Ahora evaluando cuando: X ( 0 )=0 , X ' ( 0 )=0
Obtenemos:
M [ S 2 X ( s ) ]+ fv [ SX ( s ) ] + KX (s)=f (s)
Lo siguiente seria factorizar X (S)
2
X ( s ) [ M S + fvS+ K ]=f (s)
Nuestra función de transferencia es:
X (s) 1
=
f (s ) MS + fvS+ K
2
Una vez obtenida nuestra función de transferencia podemos empezar con nuestro
análisis de estabilidad. Para saber si nuestro sistema es estable o no buscaremos
los polos que tiene la función de transferencia ya que un sistema es estable si
todos sus polos se encuentran en el semiplano izquierdo del plano complejo. Para
sacar los polos de una función de transferencia se debe igualar el polinomio del
denominador a cero para encontrar los polos, después de igualarlo, se resuelven
las ecuaciones resultantes para encontrar los valores de 's' que corresponden a
los polos.
Antes de comenzar sustituiremos los valores de M=7, fv=4 y K=10 dados en el
ejercicio.
2
7 S +4 S +10=0
Haremos uso de la formula general que nos sirve para encontrar las soluciones o
raíces de cualquier ecuación cuadrática (de segundo grado).
−4 ± √ 42 −4 (7)(10)
S1 ,2 =
2 (7)
−4 ± √16−280
S1 ,2 =
14
−4 ± √−264
S1 ,2 =
14
−4 ± 2 √−66
S1 ,2 =
14
−4 ± 2 √−66
S1 ,2 =
14
Como podemos aprecias tenemos un numero imaginario ( √ −66 ¿ el cual vamos a
simplificar junto con toda la expresión y como resultado nos queda:
−2 i √ 66
S1 ,2 = ±
7 7
S1=−0.2857+1.1605 i
S2=−0.2857−1.1605 i
Podemos notar que ambos polos están en el semiplano izquierdo del plano
complejo, eso significa que el sistema es estable.
Son complejos conjugados porque tienen la misma parte real y partes imaginarias
iguales en magnitud, pero opuestas en signo.
A continuación, se muestra la implementación en Matlab de un código hace un
análisis del criterio de estabilidad (imagen1.1) y nos muestra una gráfica de los
polos en el plano complejo(imagen1.3), así como también los polos del sistema en
el command Window(imagen1.2).
Imagen1.1 Imagen1.2
Imagen1.3
El siguiente punto seria hacer el análisis de la respuesta transitoria del sistema
MRA para la entrada escalón unitario.
La forma estándar de un sistema de segundo orden, especialmente en el contexto
de la teoría de control, es una función de transferencia de la forma:
2
ωn
G(s)=K s 2 2
S +2 ζ ωn S +ω n
Donde:
𝜔𝑛 es la frecuencia natural del sistema
𝜁 es el coeficiente de amortiguación
K s es la ganancia estática
Conociendo la forma estándar de un sistema de segundo orden, tenemos que
buscar el equivalente en nuestra función de transferencia.
X (s ) 1
=
f ( s ) MS + fvS+ K
2
Sustituyendo M, fv y K tenemos:
X (s ) 1
= 2
f ( s ) 7 S + 4 S+10
Multiplicamos 1/7 en el numerador y en el denominador para dejar S2:
1
X (s ) 7
=
f (s) 4 10
S 2+ S+
7 7
2
Ahora queremos que ω n coincida tanto en el numerador como en el denominador
y para eso multiplicamos 10/10 en el denominador para no alterar nuestra
ecuación:
1
∗10
7
X (s ) 10
=
f (s) 2 4 10
S + S+
7 7
Entonces el equivalente de la forma estándar de la función de transferencia de un
sistema de segundo orden a nuestra función de transferencia es:
10
X (s ) 1 7
=
f ( s ) 10 2 4 10
S + S+
7 7
Con esta nueva expresión podemos encontrar el valor de ζ , ω n2 , ω n
Que haciendo unos despejes nos quedarían como:
2 10
ωn = =1.4285rad / s 2 2
7
ω n=
√ 10
7
=1.1952 rad/ s
4 4
7 7
ζ= = =0.2390
√
2 ωn 10
2
7
Algunas igualdades útiles son:
4
σ =ζ ω n=
7
∗
√ 10
=0.2856
√ 10 7
2
7
ω d=ωn √ 1−ζ 2=1.1605 rad/ s
−1 ωd
β=tan =1.3294
σ
Donde:
σ = es la atenuación del sistema
β = es el amortiguamiento real
Según Ogata (2010), “[Al especificar las características de la respuesta transitoria
de un sistema de control para una entrada escalón unitario, es común especificar
lo siguiente:
1. Tiempo de retardo, td
2. Tiempo de subida, tr
3. Tiempo pico, tp
4. Sobreelongación, Mp
5. Tiempo de asentamiento, ts]” (p. 170).
Las fórmulas para calcular estas especificaciones son:
1+ 0.7 ζ
Tiempo de retardo, td ≈
ωn
π −β
Tiempo de subida, tr ¿
ωd
π
Tiempo pico, tp¿
ωd
Sobreelongación, Mp ¿ e
−
( ωσ ) π
d
4
Tiempo de asentamiento (para el criterio de 2%), ts¿
σ
3
Tiempo de asentamiento (para el criterio de 5%), ts¿
σ
“Las fórmulas empleadas en este análisis se derivan siguiendo el procedimiento
mostrado en el Ejemplo 5-1 (p. 175) de Ogata (2002, 5.ª ed.).”
Ahora sustituimos las fórmulas con nuestros valores:
1+ 0.7(0.2390)
Tiempo de retardo, td ≈ ≈ 0.9766 s
1.1952
π −(1.3294 )
Tiempo de subida, tr ¿ =1.5615 s
1.1605
π
Tiempo pico, tp¿ =2.7071 s
1.1605
Sobreelongación, Mp ¿ e−( 1.1605 ) π =0.4615 , 46.15 %
0.2856
4
Tiempo de asentamiento (para el criterio de 2%), ts¿ =14.0056 s
0.2856
3
Tiempo de asentamiento (para el criterio de 5%), ts¿ =10.5042 s .
0.2856
Ahora comprobaremos nuestros resultados en Matlab y los graficaremos para su
mejor visualización.
El código de Matlab (imagen1.4):
Imagen1.4
La grafica resultante del código es la siguiente (imagen1.5):
Imagen1.5
Los resultados impresos en el command Window se muestran en la (imagen1.6)
Imagen1.6
Como puede observarse, el análisis de la respuesta transitoria obtenido mediante
MATLAB coincide con el realizado manualmente por nuestro equipo. Las ligeras
diferencias detectadas se deben al redondeo de los valores a cuatro cifras
decimales.