Ejemplo de solución de SLI mediante Laplace

La transformada de Laplace se constituye en un método muy importante para solucionar ecuaciones diferenciales de orden superior con coeficientes constantes, o dicho de otra manera, Sistemas Lineales Invariantes en el tiempo (SLI).

Partiendo de la ecuación de oscilaciones: (D^2+2\zeta\omega_n D+\omega_n^2)y(t)=f(t) con condiciones iniciales: y(0)=y_0, y'(0)=p_0 y aplicando la transformada de Laplace se obtiene:

s^2 Y(s)-sy(0)-y'(0)+2\zeta\omega (sY(s)-y(0))+\omega_n^2Y(s)=F(s)

Simplificando:

(s^2+2\zeta\omega s+\omega_n^2)Y(s)=F(s)+y_0s+p_0+2\zeta\omega_n y_0
Y(s)=\dfrac{y_0s+p_0+2\zeta\omega_n y_0 }{s^2+2\zeta\omega s+\omega_n^2} + \dfrac{F(s)}{s^2+2\zeta\omega s+\omega_n^2}

Como puede verse, el primer término de Y(s) depende de las condiciones iniciales y corresponde a la transformada de Laplace de la solución homogénea o solución transitoria. Esto es:

y_c(t)=L^{-1} \lbrace \dfrac{y_0s+p_0+2\zeta\omega_n y_0 }{s^2+2\zeta\omega s+\omega_n^2} \rbrace

Esta solución dependerá de la naturaleza del sistema, es decir, si es un SLI sobre-amortiguado, criticamente-amortiguado o sub-amortiguado. Lo que es lo mismo, dependerá de las raíces del polinomio característico: p(\lambda)=\lambda^2+2\zeta\omega\lambda+\omega_n^2=0. Cuyas raíces son: \lambda_{1,2}=-(\zeta \pm \sqrt{1-\zeta^2})\omega_n= -\alpha \pm j\omega_d. Entonces tenemos las opciones:

Para movimiento sobre-amortiguado: y_c(t)=\frac{2\alpha E}{\lambda_2-\lambda_1}[e^{-\lambda_1 t} -e^{-\lambda_2 t} ]u(t)

Para movimiento sub-amortiguado: y_c(t)=\frac{2\alpha E}{\omega_d}e^{-\alpha t}\sin(\omega_d t)u(t).

Para movimiento críticamente-amortiguado: y_c(t)=2\alpha E t e^{-\alpha t}u(t).

El segundo término de la transformada Y(s), depende de la transformada de la función excitación f(t), por lo que corresponde a la solución en estado estacionario o solución particular.  Esto es:

y_{ss}(t)=L^{-1} \lbrace \dfrac{F(s)}{s^2+2\zeta\omega s+\omega_n^2} \rbrace

 De esta forma la solución al PVI es la transformada inversa de la función Y(s), así:

y(t)=y_c(t)+y_{ss}(t)=L^{-1} \lbrace \dfrac{y_0s+p_0+2\zeta\omega_n y_0 }{s^2+2\zeta\omega s+\omega_n^2} \rbrace +L^{-1} \lbrace \dfrac{F(s)}{s^2+2\zeta\omega s+\omega_n^2} \rbrace

Como se ha estudiado, la solución transitoria  de este SLI (estable) se extinguirá a lo largo del tiempo, quedando únicamente la respuesta en estado permanente o estado estacionario. Esto es: \lim \limits_{t \to \infty} y_c(t)=0  y nos queda: y(t) \approx y_{ss}(t).   Cuando la salida ha alcanzado estos intervalos temporales en los que la solución se “estabiliza” se suele usar el término “régimen de estado permanente” o “estado estacionario” de la señal.

Si hablamos del régimen permanente de la solución en el dominio del tiempo, es equivalente en el dominio de s a tener en cuenta únicamente el segundo término de la transformada de Laplace.  Así las cosas, la transformada de Laplace en régimen permanente es:

Y(s)=\dfrac{F(s)}{s^2+2\zeta\omega s+\omega_n^2}

Observe que esta también se puede obtener si las condiciones iniciales son cero, esto es: y_0=0, p_0=0.

Con esto, podemos escribir el siguiente cociente:

H(s)=\dfrac{Y(s)}{F(s)}=\dfrac{1}{s^2+2\zeta\omega s+\omega_n^2}

Que recibe el nombre de “función de transferencia del sistema“.  Observe que efectivamente, el cociente obtenido es una relación entre la salida Y(s) y la entrada del sistema F(s), o lo que es lo mismo, una relación de ganancia del sistema.

Si observamos detenidamente a H(s), esta se puede ver también como la respuesta del SLI si la excitación es el impulso unitario, esto es, si F(s)=1.  Por lo que su transformada inversa corresponderá a la respuesta del SLI al impulso unitario, es decir:

Si H(s)=1 entonces H(s)=Y(s), análogamente: y(t)=L^{-1}\lbrace H(s) \rbrace=h(t)

Recordemos que en el dominio del tiempo, la respuesta del SLI se obtiene mediante la operación de convolución entre la entrada y la respuesta al impulso unitario: y(t)=h(t)\otimes f(t) y en el dominio de la frecuencia la salida es simplemente el producto de la función de transferencia y la entrada: Y(s)=H(s)F(s).

Es evidente la equivalencia: y(t)=h(t) \otimes f(t) <=> Y(s)=H(s)F(s)

Que también podemos escribir como:

Y(s)=L\lbrace h(t)\otimes f(t)\rbrace =H(s)F(s) <===> y(t)=L^{-1}\lbrace H(s)F(s) \rbrace=h(t)\otimes f(t).

Entonces podemos ver que si en dominio del tiempo un SLI está modelado en función de su respuesta al impulso h(t), es equivalente en el dominio de la frecuencia (s) a modelarlo a través de su función de transferencia H(s).

SLI en dominio de t y s

 

En general, para un SLI de orden superior (orden n), la función de transferencia es típicamente el cociente entre dos polinomios racionales enteros de grados n, m (entero), por lo que tendrá la forma general:

H(s)=\dfrac{b_ms^m+b_{m-1}s^{m-1}+...+b_1 s+ b_0}{a_ns^n+a_{n-1}s^{n-1}+...+a_1 s +a_0}

Las raíces del denominador de la función de transferencia (dominio de s) se denominan “polos” del sistema, mientras que las raíces del numerador se llaman “ceros” del sistema, y éstos(ambos) son (en general) números complejos.  Es claro que los polos del SLI son exactamente las raíces del polinomio característico de la EDO, por lo que su interpretación tiene que ver con el tipo de movimiento natural del sistema  o “respuesta natural”.

Por otro lado, los ceros del SLI, corresponderán a los valores de frecuencia (s) de entrada para los cuales la salida del sistema, en estado estacionario, se hace nula.

De esta manera, un SLI se puede modelar en el dominio de la frecuencia (s) mediante la función de transfrencia H(s) y asociado a ella, se suele dibujar en el plano complejo, los polos y ceros del sistema.  A esto se le denomina “diagrama de polos y ceros“.  En dicho diagrama, usted podrá inferir de una manera rápida el tipo de movimiento del sistema y los valores de frecuencia de entrada para los cuales la salida es nula.

Cap3Mod21_2a

 

Ejemplo 1

Sea el SLI modelado mediante la función de transferencia: H(s)=\dfrac{s+1}{s^2+2s+2}.

Determine la respuesta y(t) del sistema ante la excitación: f(t)=e^{-t}u(t).

Solución:

Método 1: mediante convolución.

En el dominio del tiempo, la respuesta del SLI se puede determinar mediante la convolución de h(t) y la entrada o excitación, esto es: y(t)=h(t)\otimes f(t).  En cuyo caso deberemos hallar la transformada inversa de la función de transferencia:

h(t)= L^{-1} \lbrace \dfrac{s+1}{s^2+2s+2} \rbrace=L^{-1} \lbrace\dfrac{(s+1)}{(s+1)^2+1} \rbrace=e^{-t}\cos(t)u(t)

Por tanto, la respuesta es: y(t)=e^{-t}\cos(t)u(t) \otimes e^{-t}u(t) =\int{_0^t}{e^{-\tau}\cos(\tau)u(\tau) e^{-(t-\tau)}u(t-\tau)}d\tau.

Para realizar la operación de convolución, recuerde que el producto u(\tau)u(t-\tau) es unitario en el intervalo0<\tau<t, por lo que operación se resume a simplemente realizar la integral:

\int{_0^t}e^{-(t-\tau)}e^{-\tau}\cos(\tau)d\tau = e^{-t}\int{_0^t}{\cos(\tau)}=e^{-t}\sin(\tau)]_0^t = e^{-t}[\sin(t)-0].

Finalmente, como la respuesta es válida únicamente para t>0, entonces, se escribe como:

y(t)=e^{-t}\sin(t) u(t)

Método 2: mediante transformada inversa

Si la excitación es f(t)=e^{-t}u(t), entonces su transformada de Laplace es: F(s)=\dfrac{1}{s+1}.

Con esto, la salida, en el dominio de la frecuencia, nos queda:

Y(s)=H(s)F(s)=\dfrac{s+1}{(s+1)(s^2+2s+2)}=\dfrac{1}{s^2+2s+2}

Y sacando la inversa, tenemos la función y(t) que nos interesa:

y(t)=L^{-1} \lbrace \dfrac{1}{s^2+2 s+1} \rbrace=L^{-1} \lbrace \dfrac{1}{(s+1)^2+1} \rbrace = e^{-t}\sin(t)u(t)

Con ayuda de Maxima:

La secuencia de comandos para encontrar la inversa de la función es:

f(t):=exp(-t)$ 
F(s):=laplace(f(t),t,s)$
H(s):=(s+1)/(s^2+2*s+2)$Maxima
ilt(H(s)*F(s),s,t);

Cuyo resultado es:

%e^-t*sin(t)Maxima

 

Ejemplo 2

Resuelva el PVI: y''(t) + 2y'(t)+y(t)=u(t)-u(t-1) con y(0)=y'(0)=0.

Solución:

Transformando la EDO en el dominio de la frecuencia, nos queda:

s^2Y(s)-sy(0)-y'(0)+2[sY(s)-y(0)]+Y(s)=\dfrac{1}{s}-\dfrac{e^{-s}}{s}

Con los condiciones iniciales y simplificando, se obtiene: Y(s)=\dfrac{1-e^{-s}}{s(s^2+2s+1)}

Aplicando la inversa: y(t)=L^{-1} \lbrace\dfrac{1}{s(s^2+2s+1)} \rbrace -L^{-1} \lbrace\dfrac{e^{-s}}{s(s^2+2s+1)} \rbrace

Si denotamos: y_1(t)=L^{-1} \lbrace\dfrac{1}{s(s^2+2s+1)} \rbrace, entonces, usando la propiedad de traslación temporal, se tiene que el segundo término de la transformada inversa es: y_2(t)= L^{-1} \lbrace\dfrac{e^{-s}}{s(s^2+2s+1)} \rbrace = y_1(t-1)u(t-1) .  Por lo que únicamente debemos realizar la inversa para determinar y_1(t).

Para determinar a y_1(t), debemos escribir la transformada en su forma canónica mediante fracciones parciales:

\dfrac{1}{s(s^2+2s+1)}=\dfrac{1}{s(s+1)^2}= \dfrac{A}{s}+\dfrac{B}{s+1}+\dfrac{C}{(s+1)^2}

De donde: A(s+1)^2+Bs(s+1)+Cs = 1 . Haciendo s=-1 se obtiene: C=-1. Haciendo s=0 => A=1, y finalmente agrupando términos en s^2 => A+B=0 de donde B=-1.  Con esto, nos queda:

\dfrac{1}{s(s^2+2s+1)}=\dfrac{1}{s}-\dfrac{1}{s+1}-\dfrac{1}{(s+1)^2}

Entonces sacando inversa: f(t)=L^{-1} \lbrace\dfrac{1}{s}-\dfrac{1}{s+1}-\dfrac{1}{(s+1)^2} \rbrace = u(t)-e^{-t}u(t)-te^{-t}u(t)

Finalmente, la solución al PVI es:

y(t)=y_1(t)-y_1(t-1)u(t-1) con y_1(t)=[1-e^{-t}-te^{-t}]u(t)

Forma alterna para las fracciones parciales:

Tenemos que realizar la fracción parcial: \dfrac{1}{s(s+1)^2}.  Sin embargo, sabemos por la propiedad de integración que: “integrar en el tiempo equivale a dividir por s en la frecuencia”.  Usando esto, entonces solo es necesario realizar la transformada inversa del término \dfrac{1}{(s+1)^2} y luego integrar este resultado. Veamos:

y_1(t)= L^{-1}\lbrace \dfrac{1}{s(s+1)^2}\rbrace=\int{_0^t} L^{-1}\lbrace \dfrac{1}{(s+1)^2}\rbrace dt=\int{_0^t}{te^{-t}u(t)}dt = -(t+1)e^{-t}u(t)]{_0 ^t}=-(t+1)e^{-t}u(t)+u(t)=[1-e^{-t}-te^{-t}]u(t)

Con ayuda de Maxima:

Para encontrar la solución al PVI dado primero debemos definir la ecuación diferencial con sus condiciones iniciales y transformarla al dominio de Laplace, mediante la siguiente secuencia de comandos:

ode:'diff(y(t),t,2)+2*'diff(y(t),t)+y(t)=unit_step(t)-unit_step(t-1);
atvalue(y(t),t=0,0)$ atvalue('diff(y(t),t),t=0,0)$
laplace(ode,t,s);

Cuyo resultado es:

(%o1) 'diff(y(t),t,2)+2*('diff(y(t),t,1))+y(t)=unit_step(t)-unit_step(t-1)
(%o2) s^2*laplace(y(t),t,s)+2*s*laplace(y(t),t,s)+laplace(y(t),t,s)=1/s-%e^-s/s

Para despejar la transformada Y(s) resolvemos un sistema algebraico mediante el comando solve():

solve(%, 'laplace(y(t), t, s));

Observe que la comilla simple evita que Maxima intente calcular la transformada. De esta manera se obtiene:

[laplace(y(t),t,s)=(%e^-s*(%e^s-1))/(s^3+2*s^2+s)]

Ahora tomamos el lado derecho de la igualdad de la expresión en corchetes y definimos la función a transformar inversamente, mediante:

Y(s):=rhs(part(%,1))

Ya que Maxima tiene la limitación de no resolver transformadas inversas que involucren funciones exponenciales en s (traslaciones temporales), entonces, en este caso, debemos definir Y1(s) para tomar únicamente el término polinómico. Sabemos que el resultado habrá que trasladarlo en la operación indicada por la expresión 1-e^{-s}. Entonces nos queda:

Y1(s):=1/(s^3+2*s^2+s);
ilt(Y1(s),s,t);

Cuyo resultado es:

 -t*%e^-t-%e^-t+1

Finalmente, de la teoría sabemos que y(t)=y_1(t)-y_1(t-1)

 

Es de notar que el PVI en cuestión corresponde a un SLI críticamente amortiguado cuya función de transferencia es H(s)=\dfrac{1}{(s+1)^2}.  Si graficamos la salida con la correspondiente entrada (pulso unitario), se obtiene la figura:

Ej2-PVIcritica

El código en Maxima para obtener la figura anterior es:

p(t):=unit_step(t)-unit_step(t-1);
y1(t):=(1-exp(-t)-t*exp(-t))*unit_step(t);Maxima
y(t):=y1(t)-y1(t-1)$
plot2d([p(t),y(t)],[t,0,5],[y,0,1.2]);