## La transformada de Laplace.

### Transformada directa de Laplace.

Maxima tiene implementada una función para obtener la transformada de Laplace de la 
mayoría de las funciones que conocemos. 

Esta transformada se calcula utilizando el siguiente comando:

        laplace(f(t),t,z)

donde t es la variable de la que depende la función y z es la variable independiente de la 
que dependerá la transformada de f(t).

EJEMPLO 1: Vamos a calcular la transformada de Laplace de las funciones f(t) = t^2 y  g(t) = sen(t):

f(t):=t^2;
laplace(f(t),t,z);

g(t):=sin(t);
laplace(g(t),t,z);

Como siempre no es necesario definir las funciones f(t) y g(t) previamente para poder utilizar la función
laplace, podemos incluir la definición de la función directamente:

laplace(sin(t),t,z);

NOTA: Este mismo resultado lo encontraríamos si en el menú superior de Maxima pinchamos en 

"Análisis->Transformada de Laplace" 

y rellenamos los campos que nos piden.

Normalmente representaremos por F(z) a la transformada de Laplace de f(t).

EJEMPLO 2: Definimos F(z) como la transformada de la función

    exp(-a*t)*sin(b*t):

F(z):=laplace(exp(-a*t)*sin(b*t),t,z);

y podemos ver  el resultado escribiendo el nombre de la función.

F(z);

FUNCIÓN DE HEAVISIDE:

Podemos calcular la transformada de Laplace de funciones definidas a trozos con ayuda de la 
función escalón unitario o función de Heaviside h_0(t), que en Maxima se representa por

    unit_step(t). 

Vemos su representación gráfica:

wxplot2d(unit_step(t),[t,-2,2],[y,-2,2]);

Y su transformada de Laplace sería

F(z):=laplace(unit_step(t),t,z);
F(z);

como era de esperar.

NOTA: Realmente la función definida en MAXIMA es continua a la izquierda y no a la derecha 
como la hemos definido en clase, es decir, la función unit_step de MAXIMA está definida como:

    unit_step(t) = 0,  si 0<=t
    unit_step(t) = 1,  si t>0

aunque al estar definida mediante una integral su transformada de Laplace es la misma que para la función 
definida como :

    h_0(t) = 0 si 0<t
    h_0(t) = 1 si t>=0

Para el caso de la transformada de Laplace para 
una función de Heaviside de parámetro a>0: 

    h_a(t)

se utiliza el hecho de que h_a(t) es la versión trasladada de h_0(t) al instante "a", es decir:

    h_a(t) = h_0(t-a)

que vuelve a ser continua a la izquierda y no a la derecha como usualmente.

EJEMPLO 3: Calcula la transformada de Laplace de h_2(t).

Recordemos que
 
h_2(t) = 0 si t<2
h_2(t) = 1 si t>= 2

y gráficamente

wxplot2d([unit_step(t),unit_step(t-2)],[t,-2,4],[y,-2,2]);

Para calcular su transformada de Laplace empleamos la instrucción de MAXIMA laplace correspondiente directamente

F(z):=laplace(unit_step(t-2),t,z);
F(z);

Con esta función, unit_step(t), podemos calcular la Transformada de Laplace 
de funciones definidas a trozos,como vemos en el ejemplo siguiente:

EJEMPLO 4: Representar gráficamente y calcular la transformada de
Laplace de la función definida a  trozos

f(t) = 0,             si    t<= 0
f(t) = 3,             si  0<t<= 1; 
f(t) = 3*t^2,      si  1<t<= 3, 
f(t) = cos(2*t), si  3<t:

Utilizando unit_step, definimos la función f(t)

f(t):=3*(unit_step(t)-unit_step(t-1)) +(3*t^2)*(unit_step(t-1)-unit_step(t-3)) +
cos(2*t)*unit_step(t-3);

Representamos la función gráficamente en [0,5]

wxplot2d(f(t),[t,0,5]);

Finalmente calculamos su transformada de Laplace

F(z):=laplace(f(t),t,z);
F(z);

FUNCIÓN DELTA DE DIRAC

Recordemos que la función delta de Dirac estaba definida como

d(t) = 0 para  t ≠ 0

/ inf
| 
| delta(t) dt = 1
|  
/ -inf

y recordemos que su transformada de Laplace era la unidad.

laplace(delta(t),t,z);

NOTA: La única función que puede utilizar la función delta(t) es laplace.

PROPIEDADES DE LA TRANSFORMADA DE LAPLACE


Podemos comprobar algunas de las propiedades de la transformada de Laplace

1) PRIMER TEOREMA DE TRASLACIÓN

        L[exp(at)*f(t)](z) = L[f(t)](z-a)

EJEMPLO 5:

g1(t):=t;
lg1(z):= laplace(g1(t),t,z)$

g2(t):=exp(2*t)*g1(t);
lg2(z):=laplace(g2(t),t,z)$

lg1(z-2);
lg2(z);

2) SEGUNDO TEOREMA DE TRASLACIÓN

    L[h_a(t)*f(t-a)](z) = exp(-az)*L[f(t)](z) para a>0

EJEMPLO 6:

g1(t):=t^2;
lg1(z):= laplace(g1(t),t,z)$

g3(t):= unit_step(t-2)*g1(t-2);
lg3(z):= laplace(g3(t),t,z)$

exp(-2*z)*lg1(z);
lg3(z);

ratsimp(lg3(z));

3) CAMBIO DE ESCALA

        L[g(a*t)](z) = (1/a)*L[f(t)](z/a) para a>0

EJEMPLO 7:

g1(t):=sin(t);
lg1(z):= laplace(g1(t),t,z)$

g4(t):=g1(3*t);
lg4(z):=laplace(g4(t),t,z)$

(1/3)*lg1(z/3);
lg4(z);

ratsimp(%o44);

4) DERIVADA TEMPORAL

Podemos calcular la transformada de Laplace de la derivada:

kill(all)$
laplace(diff(h(t),t),t,z);

EJEMPLO 8: Para calcular la transformada de Laplace de la derivada de g(t)=t^5:

g(t):=t^5;
F(z):=laplace(diff(g(t),t,1),t,z);
F(z);

expresión que coincIde al calcular la transformada de mediante la fórmula correspondiente:

z*laplace(g(t),t,z)-g(0);

Se pueden obtener la fórmula general para la derivadas siguientes, por ejemplo la derivada
segunda:

kill(all)$
laplace(diff(h(t),t,2),t,z);

tercera:

laplace(diff(h(t),t,3),t,z);

o décima

laplace(diff(h(t),t,10),t,z);

5) DERIVADA EN EL PLANO Z

Recordemos que la función transformada de Laplace es holomorfa en su espacio de definición y que esta derivada
podía utilizarse para obtener transformadas de Laplace de funciones de la forma t^n f(t), mediante la expresión

d^n                             n
--------  L(f(t),t,z) = (-1)    L(t^n f(t),t,z)
dz^n

laplace(t^5*sin(t),t,z);

(-1)^5*diff(laplace(sin(t),t,z),z,5);

6) INTEGRACIÓN

Podemos calcular la transformada de Laplace de una integral. En este caso
no se utiliza la propiedad correspondiente, sino que se calcula la integral y después
la transformada. Podemos comprobarlo

kill(all)$
laplace(integrate(h(t),t,z),t,z);

donde debería poner

(1/z)*laplace(h(t),t,z);

EJEMPLO 9: Para calcular la transformada de Laplace de la integral de g(t)=t^5:

g(t):=t^5;
F(z):=laplace(integrate(g(t),t),t,z);
F(z);

Notar que sí se cumple la propiedad, ya que si calculamos la primitiva de g(t)

G(t):=integrate(g(t),t)$
G(t);

Su transformada de Laplace

laplace(G(t),t,z);

Coincide con la que obtendríamos aplicando la fórmula correspondiente

(1/z)*laplace(g(t),t,z);

### Transformada inversa de Laplace.

Para calcular la transformada inversa de Laplace de una función, basta con teclear la orden

        ilt(F(z), z, t)

Este comando calcula bien  la transformada inversa de Laplace de funciones F(z) que sean fracciones de polinomios cuyo denominador 
tenga sólo factores lineales y cuadráticos y sus potencias, como vemos en los dos siguientes ejemplos:

EJEMPLO 10: Podemos calcular la transformada inversa de Laplace
de la función F(z)=1/(z/2+1):

F1(z) := 1/(2*z+1);
f1(t) := ilt(F1(z),z,t);
f1(t);

Podemos comprobar como la transformada de Laplace de f(t) es F(z)

laplace(f1(t),t,z);

EJEMPLO 11:  Calcular la transformada inversa de Laplace de la función 

    F(z) = z^2/(z^2+1)^2

F2(z):= z^2/(z^2+1)^2;
f2(t):= ilt(F2(z),z,t);
f2(t);

EJEMPLO 12:  Calcular la transformada inversa de Laplace de la función
 
F(z)=1/(z^4-4*z^3+5*z^2-4*z+4):

F3(z):=1/(z^4-4*z^3+5*z^2-4*z+4);
f3(t):= ilt(F3(z),z,t);
f3(t);

Si fuera necesario, podríamos realizar primero la descomposición de dicho cociente en fracciones simples y con posterioridad 
calcular la transformada inversa de Laplace de cada una de estas fracciones.

Para reducir un polinomio a fracciones simples podemos usar la sentencia

    partfrac(f(s),s);

Lo vemos en el siguiente ejemplo:

EJEMPLO 13: Descomponemos la función en fracciones simples

F(z):=1/(z^4-4*z^3+5*z^2-4*z+4);
partfrac(F(z),z);

Utilizando las propiedades de linealidad de la transformada inversa de Laplace, obtendremos la 
transformada inversa de Laplace de F(z) como suma de las transformadas inversas de cada una 
de las fracciones simples:

ilt((3+4*z)/(25*(z^2+1)),z,t)-
ilt(4/(25*(z-2)),z,t)+
ilt(1/(5*(z-2)^2),z,t);

Aunque para esta función se obtiene el mismo resultado de forma directa

ilt(F(z),z,t);

PROPIEDADES DE LA TRANSFORMADA DE LAPLACE INVERSA

1) TEOREMA DE TRASLACIÓN

Para obtener la transformada inversa de Laplace de funciones de la forma 

    exp(-a*z)*F(z), para a>0; 

donde F(z) es un cociente de polinomios, necesariamente hemos de usar el segundo teorema de traslación, 
que esquemáticamente nos dice: 

    TIL(exp(-a*z)*F(z))=h_a(t)*TIL(F(z))(t-a),

siendo h_a(t) la funcion de Heaviside de parámetro a.

Si lo hacemos directamente tendremos un error:

EJEMPLO 14: Calcula la transformada inversa de

e^(-2z)
---------
   z^2

kill(all)$;
F(z) := %e^(-2*z)/z^2;

Si ahora utilizamos ilt para calcular la inversa obtendremos:

ilt(F(z),z,t);

¡MAXIMA no sabe cómo calcular el resultado!

Para obtener la transformada inversa hay que utilizar el segundo teorema de traslación.

Tomamos la función a la que queremos calcular su transformada inversa

    F(z)

F(z);

Prescindimos de la exponencial

Ftemp(z) := 1/z^2;

Y calculamos la transformada inversa de esta función, pero en lugar de evaluarla en el punto t, lo hacemos en el punto (t-2)
ya que estábamos multiplicando por %e^(-2z) y después el resultado se multiplica por la función de heaviside correspondiente, unit_step(t-2):

f(t):= ilt(Ftemp(z),z,t-2)*unit_step(t-2)$
f(t);

Comprobamos que el resultado es correcto, calculando la transformada de Laplace de la función f(t)

F2(z):= laplace(f(t),t,z)$
F2(z);

Y simplificando

ratsimp(F2(z));

Comprobamos que coincide con la función F(z) inicial.

EJEMPLO 15:  Calcula la transformada inversa de Laplace de la función 

F(z)= exp(-%pi*z)*(1/(z^2+1))

F(z):=exp(-%pi*z)*(1/(z^2+1));

Definimos la función Ftemp(z) , prescindiendo de la función exponencial

Ftemp(z):= (1/(z^2+1));

Como antes, usamos el 2º teorema de traslación y la transformada inversa pedida vendrá dada por

TIL(G(z)) = TIL(exp(-%pi*z)*F(z)) = h_%pi(t)*TIL(F(z))(t-%pi)

f(t):=unit_step(t-%pi)*ilt(Ftemp(z),z,t-%pi);
f(t);

wxplot2d(f(t),[t,0,4*%pi]);

EJEMPLO 16: Calcula la transformada inversa de la función

F(z) = exp(-2*z)*(z/(z^2+1)) + exp(-4*z)*(1/z)

Definimos la función F(z)

kill(all)$
F(z):= exp(-2*z)*(z/(z^2+1)) + exp(-4*z)*(1/z);

Notar que hay dos funciones exponenciales y que podemos usar la propiedad de linealidad

Si definimos F1(z) y F2(z) como

F1(z):= z/(z^2+1);
F2(z):= 1/z;

entonces

F(z) : = exp(-2z)*F1(z) + exp(-4z)*F2(z)

Por linealidad

TIL(F(z)) = TIL(exp(-2z)*F1(z)) + TIL(exp(-4z)*F2(z))

y empleando el segundo teorema de traslación para la transformada inversa

TILF(F(z)) = h_2(t) * TIL(F1)(t-2) + h_4(t) * TIL(F2)(t-4)

f(t):= unit_step(t-2)*ilt(F1(z),z,t-2) + unit_step(t-4)*ilt(F2(z),z,t-4)$
f(t);

Cuya representación es

wxplot2d(f(t),[t,1,10]);

Notar las discontinuidades en t=2 y en t=4.

Y cuya transformada de Laplace

laplace(f(t),t,z);

coincide con la función original F(z).

Fórmula de INVERSIÓN de BROMWICH

Cuando la función F(z) es una función racional podremos aplicarl Fórmula de Inversión mediante 
números complejos usando los residuos en las singularidades de la función exp(zt)*F(z)



EJEMPLO 17:  Calcular la transformada inversa de Laplace de la función  

    F(z) = z^2/(z^2+1)^2

mediante residuos.

F(z):=z^2/(z^2+1)^2;
polos : solve(denom(F(z))=0,z);
suma:0;
for k in polos do (ldisp(k),temp:residue(exp(z*t)*F(z),z,rhs(k)),ldisp(temp),suma:suma+temp);
suma;
ratsimp(rectform(suma));

Notar que es el resultado que obtenemos para la función ilt

ilt(F(z),z,t);

suma;

Podemos calcular los residuos mediante la definición de límite. Para la función anterior

F(z) = z^2/(z^2+1)^2

Calculamos las singularidades y sus multiplicidades

F(z):=z^2/(z^2+1)^2;
polos : solve(denom(F(z))=0,z);
multiplicities;

Luego son polos dobles, así que los residuos se calcularán segújn la fórmula general

lim          1/k-1! (d/dz)^k-1 ( (z-z0)^k G(z))
z->z_0 

En este caso G(z) = exp(zt)*F(z)

z0:%i$
m0:2$
factor_numerico: (1/(m0-1)!)$
factor_funcion(z) := (z-z0)^(m0)$
G(z) := diff(factor_funcion(z)*exp(z*t)*F(z),z,m0-1)$
residuo1: factor_numerico*limit(G(z),z,z0);

z0:-%i$
m0:2$
factor_numerico: (1/(m0-1)!)$
factor_funcion(z) := (z-z0)^(m0)$
G(z) := diff(factor_funcion(z)*exp(z*t)*F(z),z,m0-1)$
residuo2: factor_numerico*limit(G(z),z,z0);

ratsimp(rectform(residuo1+residuo2));

### Resolución de problemas de valor inicial.

Veremos como podemos resolver una EDO, concretamente los llamados problemas de valor inicial mediante el uso de la transformada
de Laplace (recordamos que usar la transformada de Laplace para resolver ecuaciones diferenciales hace que el problema a resolver se
transforme en un problema de resolver una ecuación lineal):

kill(all)$

Vamos a seguir el procedimiento visto en clase para resolver el siguiente problema de valor inicial:

    y'' + y = t

    y(0) = 1
    y'(0) = 2.

Definimos la ecuación: 

ecuacion : diff(y(t),t,2)+ y(t) = t;

Calculamos la transformada de Laplace de toda la ecuación:

lecuacion : laplace(ecuacion,t,z);

Notar que la instrucción laplace() de MAXIMA utiliza las propiedades de linealidad y de derivación.

Utilizamos las condiciones iniciales, sustituyendo sus valores en la ecuación mediante la instrucción correspondiente

lecuacion : lecuacion,y(0)=1,at('diff(y(t),t,1),t=0)=2;

Se despeja el valor de la transformada de Laplace

solucion:solve(lecuacion,laplace(y(t),t,z));

Y se calcula la inversa de Laplace de la función resultante.

ilt(rhs(solucion[1]),z,t);

La función rhs toma el miembro de la derecha de la ecuación descrita en solucion[1], 
que es el primer elemento de la variable solucion.

NOTA: MAXIMA tiene ya definida una función que emplea la transformadade Laplace para resolver este tipo de problemas. 

EJEMPLO 16: Utiliza la función desolve para resolver el PVI anterior.

Definimos la ecuación diferencial

ec:'diff(y(t),t,2)+ y(t) = t;

Incluimos las condiciones iniciales

atvalue(y(t),t=0,1);
atvalue('diff(y(t),t),t=0,2);

Usamos la instrucción desolve para resolver el PVI:

desolve([ec],[y(t)]);

Finalmente, este mismo problema también se puede resolver a través del menú superior de Maxima: 

Ecuaciones -> Condición inicial: Para cada condición del PVI debemos incluir la expresión correspondiente a la derivada, 
                                                              el punto donde se debe evaluar dicha derivada y el valor de la misma.

Ecuaciones -> Resolver EDO con Laplace: Se abrirá un ventana de captura de datos donde debemos incluir la EDO con la misma 
                                                                                  sintáxis que si estuviera escrita en línea de comandos y la función respecto de la cual 
                                                                                  se construye la ecuación.

kill(all)$
atvalue(y(t), t=0, 0);
atvalue('diff(y(t),t,1), t=0, 1);
desolve([diff(y(t),t,2)+diff(y(t),t,1)=0],[y(t)]);

EJEMPLO 17: Resuelve el siguiente Problema de Valor Inicial

2*y''(t)+y'(t)-1=sin(t),
y(0) = 0;
y'(0)= 0;

Definimos la ecuación:

kill(all)$

ec2:2*'diff(y(t),t,2)+'diff(y(t),t)-1=sin(t);

Aplicamos la transformada de Laplace

lec2: laplace(ec2,t,z);

Incluimos las condiciones iniciales

lec2: lec2,y(0)=0,at('diff(y(t),t,1),t=0)=0;

Despejamos laplace(y(t),t,z)

solucion:solve(lec2,laplace(y(t),t,z));

Y calculamos la transforma inversa para obtener la solución:

ilt(rhs(solucion[1]),z,t);

Utilizando la instrucción dsolve el problema se resolvería de la siguiente forma:

kill(all);
ec2:2*'diff(y(t),t,2)+'diff(y(t),t)-1=sin(t);
atvalue(y(t),t=0,0);
atvalue('diff(y(t),t),t=0,0);
desolve([ec2],[y(t)]);

SISTEMAS DE EDO:

También es posible resolver un sistema de ecuaciones diferenciales, siguiendo el mismo procedimiento

EJEMPLO 18: Vamos a resolver el sistema 

x'(t) + y(t) = 1
y'(t) - x(t) = 0

junto con las condiciones iniciales 

x(0)=1
y(0)=0

Definimos las 2 ecuaciones del sistema:

ec_1 : diff(x(t),t,1) + y(t) = 1;
ec_2 : diff(y(t),t,1) - x(t) = 0;

Aplicamos la transformada de Laplace a ambas ecuaciones:

lec_1 : laplace(ec_1,t,z);
lec_2 : laplace(ec_2,t,z);

Sustituimos las condiciones iniciales en ambas ecuaciones:

lec_1 : lec_1,x(0)=1,y(0)=0;
lec_2 : lec_2,x(0)=1,y(0)=0;

Despejamos las transformadas correspondientes. En este caso como es un sistema lineal utilizamos 
la instrucción linsolve que es más eficiente:

solucion:linsolve([lec_1,lec_2],[laplace(x(t),t,z),laplace(y(t),t,z)]);

Finalmente utilizamos la transformada de Laplace inversa con cada elemento de la solución:

x(t):= ilt(rhs(solucion[1]),z,t);
y(t):= ilt(rhs(solucion[2]),z,t);

x(t);
y(t);

También podemos resolver el sistema de forma directa con la instrucción desolve que se ha empleado para una 
sola ecuación.

Primero borramos la memoria de MAXIMA, para evitar errores

kill(all)$

Definimos las ecuaciones:

ec_1: 'diff(x(t),t)+y(t)=1;
ec_2: 'diff(y(t),t)-x(t)=0;

Las condiciones iniciales:

atvalue(x(t),t=0,1);
atvalue(y(t),t=0,0);

Utilizamos la función desolve:

desolve([ec_1,ec_2],[x(t),y(t)]);

PROBLEMAS DE VALOR INICIAL CON ENTRADAS CONTINUAS A TROZOS

Si la función de entrada es continua a trozos, tendremos que tener en cuenta esta situación a la hora de calcular la transformada inversa, ya que como
se ha comprobado en la sección anterior no es posible calcular esta inversa de forma directa.


EJEMPLO 19:

Por ejemplo, vamos a resolver el siguiente problema de valor inicial:

2*y''(t)+y'(t)-1=f(t), con las cond. iniciales y(0)=y'(0)=0:

donde la función f(t) está definida como

f(t) = t si  0<=t<1
f(t) = 1 si  1<=t<2
f(t) = 0 si  t>=2

Borramos memoria

kill(all)$

Definimos la función de entrada f(t) utilizando la función unit_step()

f(t):= t*(unit_step(t)-unit_step(t-1))+1*(unit_step(t-1)-unit_step(t-2));

Vamos a representar f(t):

wxplot2d(f(t),[t,0,4],[y,-1,4]);

Definimos la ecuación:

ecuacion:2*'diff(y(t),t,2)+'diff(y(t),t)-1=f(t);

Aplicamos la transformada de Laplace

lecuacion: laplace(ecuacion,t,z);

Incluimos las condiciones iniciales

lecuacion: lecuacion,y(0)=0,at('diff(y(t),t,1),t=0)=0;

Despejamos laplace(y(t),t,z)

solucion:solve(lecuacion,laplace(y(t),t,z));

Calcularemos la transforma inversa para obtener la solución, pero ahora hay que tener en cuenta la propiedad de traslación.

La función Y(z) se expresa como

Y1(z):=rhs(solucion[1])$
Y1(z);

O mejor en fracciones como:

Y(z) := (z+1)/(2*z^4+z^3)-%e^(-z)*(1/(2*z^4+z^3))-%e^(-2*z)*(z/(2*z^4+z^3));

Y para la inversa se utilizará la propiedad de linealidad para expresar como suma de inversas y después el segundo teorema de traslación a cada sumando 
para obtener la inversa.

sol: ilt((z+1)/(2*z^4+z^3),z,t)-ilt(z/(2*z^4+z^3),z,t-1)*unit_step(t-1)-ilt(1/(2*z^4+z^3),z,t-2)*unit_step(t-2);

wxplot2d(sol,[t,0,5]);

## Ejercicios para practicar.

PROBLEMA 1. Calcular la transformada de Laplace de las funciones

a) f(t)=-1.25+3.25*exp(-2*t)+1.25*exp(-2*t).

b) f(t)=exp(5*t)*cos(3*t).

c) f(t)=exp(t)*cos(t)*sin(2*t)


PROBLEMA 2. Representar gráficamente y calcular las transformadas
de Laplace de la funciones definidas por:

a) f(t)=1, si 0<=t<=1; f(t)=t^3, si 1<t<=4; f(t)=sin(2*t), si t>4.

b) f(t)=t^2-1, si 0<=t<=2; f(t)=2, si t>2.


PROBLEMA 3. Calcular la transformada inversa de Laplace para las funciones 

a) F(z)=(z+4)/(5*z^2+7*z+2).

b) F(z)=1/((z+1)*(z+2)*(z^2+2*z+10)).

c) F(z)=exp(-z)*((z-1)/(z^3+2)).


PROBLEMA 4. Resolver, usando la TL, los siguientes problemas de
condiciones iniciales:

a) La ecuación y'''(t)+5*y''(t)+17*y'(t)+13*y(t)=1, 
con y(0)=y'(0)=1; y''(0)=0.

b) La ecuación y''(t)+y(t)=t, con y(0)=1, y'(0)=-2.

c) El sistema de ecuaciones diferenciales
x'=x-z, y'=2*y, z'=x+z, con x(0)=-2, y(0)=2, z(0)=-1.