Respuesta de un oscilador a una fuerza impulsiva (II)

Pulso triangular

Sea la fuerza

F ( t ) = { F 0 a t   0 ≤ t ≤ a 0   t > a

Su transformada de Laplace es

f(t)= F(t) m =( u(t)−u(t−a) ) F 0 m t a = f 0 a ( t−u(t−a)(t−a)−u(t−a)a ) g(s)= f 0 a ( 1 s 2 − exp(−as) s 2 −a exp(−as) s )

Para el caso de que las condiciones iniciales sean: posición inicial x0=0 y velocidad inicial v0=0 en el instante t=0. La transformada de Laplace de la ecuación diferencial que describe el comportamiento del oscilador de frecuencia propia ω0 y coeficiente de la fuerza de rozamiento γ, es F(s)

d 2 x d t 2 +2γ dx dt + ω 0 2 x=f(t) F(s)= f 0 a ( 1 s 2 − exp(−as) s 2 −a exp(−as) s ) 1 s 2 +2γs+ ω 0 2

Descomponemos la fracción de la izquierda en suma de fracciones más simples

1 s( s 2 +2γs+ ω 0 2 ) = A 1 s + C 1 s+ D 1 s 2 +2γs+ ω 0 2 A 1 = 1 ω 0 2   C 1 =− 1 ω 0 2   D 1 =− 2γ ω 0 2 1 s 2 ( s 2 +2γs+ ω 0 2 ) = A 2 s + B 2 s 2 + C 2 s+ D 2 s 2 +2γs+ ω 0 2 A 2 =− 2γ ω 0 4   B 2 = 1 ω 0 2   C 2 = 2γ ω 0 4   D 2 = 4 γ 2 − ω 0 2 ω 0 4

Escribimos F(s) en forma apropiada para obtener la transformada inversa de Laplace, mirando a las tablas.

F(s)= f 0 a ( A 2 s + B 2 s 2 + C 2 s+ D 2 s 2 +2γs+ ω 0 2 ) − f 0 a exp(−as)( A 2 s + B 2 s 2 + C 2 s+ D 2 s 2 +2γs+ ω 0 2 ) − f 0 exp(−as)( A 1 s + C 1 s+ D 1 s 2 +2γs+ ω 0 2 )

La transformada inversa de Laplace es

x(t)= f 0 a ( g(t)−u(t−a)g(t−a)+a·u(t−a)h(t−a) ) h(t)= A 1 + C 1 exp(−γt)cos(ωt)+( D 1 − C 1 γ ω )exp(−γt)sin(ωt) g(t)= A 2 + B 2 t+ C 2 exp(−γt)cos(ωt)+( D 2 − C 2 γ ω )exp(−γt)sin(ωt) ω= ω 0 2 − γ 2

Aplicamos a un oscilador de frecuencia propia ω0=4, coeficiente de la fuerza de rozamiento proporcional a la velocidad γ=0. La fuerza dura un tiempo a=0.4 s y su valor máximo es f0=F0/m=40.

g=0; %rozamiento
w0=4; %frecuencia propia del oscilador
F0=40; %máximo valor de la fuerza
a=0.4; %tiempo de la fuerza

A1=1/w0^2;
C1=-1/w0^2;
D1=-2*g/w0^2;

A2=-2*g/w0^4;
B2=1/w0^2;
C2=2*g/w0^4;
D2=(4*g^2-w0^2)/w0^4;
w=sqrt(w0^2-g^2);

h=@(t) A1+C1*exp(-g*t).*cos(w*t)+(D1-C1*g)*exp(-g*t).*sin(w*t)/w;
g=@(t) A2+B2*t+C2*exp(-g*t).*cos(w*t)+(D2-C2*g)*exp(-g*t).*sin(w*t)/w;

t=linspace(0,2.4,200);
x=F0*(g(t)-g(t-a).*heaviside(t-a)-a*h(t-a).*heaviside(t-a))/a;
plot(t,x)
ylabel('x(t)')
xlabel('t')
title('Pulso F_0·(t/a)')
grid on

Obtenemos un resultado similar empleando Math Symbolic, salvo que aparece un error cuando se describe el oscilador sin rozamiento γ=0.

>> syms w0 g s t;
>> ft=100*(1-heaviside(t-0.4))*t;
>> gs=laplace(ft);
>> Fs=gs/(s^2+2*g*s+w0^2);
>> x=ilaplace(Fs);
>> xx=subs(x,{g,w0},{0.001,4}); %rozamiento muy pequeño
>> ezplot(xx,[0 2.4])
>> grid on
>> ylabel('x(t)')
>> xlabel('t')
>> title('Pulso F_0·(t/a)')

Función rampa

Aplicamos una fuerza de la forma F0·t/a durante un tiempo a≤2, la fuerza es F0 para t>a.

La función f(t) se define

f(t)= F(t) m =( u(t)−u(t−a) ) f 0 t a +u(t−a) f 0 = f 0 a ( t−u(t−a)( t−a ) )

Su transformada de Laplace g(s) es

g(s)= f 0 a ( 1 s 2 − exp(−as) s 2 )

Para un oscilador de frecuencia propia ω0 y coeficiente de la fuerza de rozamiento γ, y para las condiciones iniciales: posición inicial x0=0 y velocidad inicial v0=0 en el instante t=0, calculamos la transformada inversa de Laplace de la función.

F(s)= f 0 a ( 1 s 2 − exp(−as) s 2 ) 1 s 2 +2γs+ ω 0 2

Expresamos la siguiente fracción como suma de fracciones más simples

1 s 2 ( s 2 +2γs+ ω 0 2 ) = A s + B s 2 + Cs+D s 2 +2γs+ ω 0 2 A=− 2γ ω 0 4  B= 1 ω 0 2  C= 2γ ω 0 4  D= 4 γ 2 − ω 0 2 ω 0 4

Escribimos F(s) en forma apropiada para obtener la transformada inversa de Laplace, mirando a las tablas.

F(s)= f 0 ω f ( As+B s 2 + ω f 2 + Cs+D (s+γ) 2 +( ω 0 2 − γ 2 ) )+ f 0 π a exp(−as)( As+B s 2 + ω f 2 + Cs+D (s+γ) 2 +( ω 0 2 − γ 2 ) ) x(t)= f 0 ω f ( g(t)+u(t−a)g(t−a) ) g(t)=Acos( ω f t )+Bsin( ω f t )+Cexp(−γt)cos(ωt)+( D−Cγ ω )exp(−γt)sin(ωt) ω= ω 0 2 − γ 2   ω f = π a

Representamos gráficamente x(t) para f0=10, γ=0.5, ω0=3 rad/s

g=0.5; %rozamiento
w0=3; %frecuencia propia
f0=10; %valor máximo de la fuerza
a=2; %anchura del pulso

A=-2*g/w0^4;
B=1/w0^2;
C=2*g/w0^4;
D=(4*g^2-w0^2)/w0^4;
w=sqrt(w0^2-g^2);
g=@(t) A+B*t+C*exp(-g*t).*cos(w*t)+(D-C*g)*exp(-g*t).*sin(w*t)/w;

t=linspace(0,10,400);
x=f0*(g(t)-heaviside(t-a).*g(t-a))/a;
plot(t,x)
ylabel('x(t)')
xlabel('t')
title('Pulso rampa')
grid on

Obtenemos un resultado similar empleando Math Symbolic.

>> syms a t s g w0 f0;
>> ft=f0*(t-heaviside(t-a)*(t-a))/a;
>> gs=laplace(ft);
>> Fs=gs/(s^2+2*g*s+w0^2);
>> x=ilaplace(Fs);
>> xx=subs(x,{g w0 f0 a},{0.5 3 10 2});
>> ezplot(xx,[0 10])
>> grid on
>> ylabel('x(t)')
>> xlabel('t')
>> title('Pulso rampa')

Pulso media onda

Aplicamos una fuerza de la forma F(t)=F0·sin(ωft) durante un tiempo a=π/ωf. La función f(t) se define

f(t)= F(t) m ={ f 0 sin( ω f t) 0≤t≤a 0 t>a

La transformada de Laplace de esta función se calcula a partir de su definición, integrando dos veces por partes.

g(s)= ∫ 0 ∞ e −st f(t)·dt= f 0 ∫ 0 a e −st sin( π a t )·dt= f 0 ω f 1+exp(−as) ω f 2 + s 2

Para un oscilador de frecuencia propia ω0 y coeficiente de la fuerza de rozamiento γ, y para las condiciones iniciales: posición inicial x0=0 y velocidad inicial v0=0 en el instante t=0, calculamos la transformada inversa de Laplace de la función.

F(s)= f 0 ω f 1+exp(−as) ω f 2 + s 2 1 s 2 +2γs+ ω 0 2

Expresamos la siguiente fracción como suma de fracciones más simples

1 ( s 2 + ω f 2 )( s 2 +2γs+ ω 0 2 ) = As+B s 2 + ω f 2 + Cs+D s 2 +2γs+ ω 0 2 A= −2γ ( ω 0 2 − ω f 2 ) 2 +4 γ 2  B= ω 0 2 − ω f 2 ( ω 0 2 − ω f 2 ) 2 +4 γ 2   C= 2γ ( ω 0 2 − ω f 2 ) 2 +4 γ 2  D= 4 γ 2 − ω 0 2 + ω f 2 ( ω 0 2 − ω f 2 ) 2 +4 γ 2

Escribimos F(s) en forma apropiada para obtener la transformada inversa de Laplace, mirando a las tablas.

F(s)= f 0 ω f ( As+B s 2 + ω f 2 + Cs+D (s+γ) 2 +( ω 0 2 − γ 2 ) )+ f 0 π a exp(−as)( As+B s 2 + ω f 2 + Cs+D (s+γ) 2 +( ω 0 2 − γ 2 ) ) x(t)= f 0 ω f ( g(t)+u(t−a)g(t−a) ) g(t)=Acos( ω f t )+Bsin( ω f t )+Cexp(−γt)cos(ωt)+( D−Cγ ω )exp(−γt)sin(ωt) ω= ω 0 2 − γ 2   ω f = π a

Representamos gráficamente x(t) para f0=1, γ=1, a=π, ω 0 = 3 rad/s

g=1; %rozamiento
f0=1; %amplitud de la fuerza
w0=sqrt(3); %frecuencia propia del oscilador
a=pi;       %semiperiodo
wf=pi/a;    %frecuencia de la fuerza oscilante

A=-2*g/((w0^2-wf^2)^2+4*g^2);
C=-A;
B=(w0^2-wf^2)/((w0^2-wf^2)^2+4*g^2);
D=(wf^2-w0^2+4*g^2)/((w0^2-wf^2)^2+4*g^2);
w=sqrt(w0^2-g^2);
g=@(t) A*cos(wf*t)+(B/wf)*sin(wf*t)+
C*exp(-g*t).*cos(w*t)+(D-C*g)*exp(-g*t).*sin(w*t)/w;

t=linspace(0,10,100);
x=f0*wf*(g(t)+heaviside(t-a).*g(t-a));
plot(t,x)
ylabel('x(t)')
xlabel('t')
title('Pulso sin(t)')
grid on

El máximo valor del desplazamiento es 0.3315 en el instante 2.424

Obtenemos un resultado similar empleando Math Symbolic.

>> syms a t s g w0 f0;
>> ft=f0*(1-heaviside(t-a))*sin(pi*t/a);
>> gs=laplace(ft);
>> Fs=gs/(s^2+2*g*s+w0^2);
>> x=ilaplace(Fs);
>> xx=subs(x,{g w0 f0 a},{1 sqrt(3) 1 pi});
>> ezplot(xx,[0 10])
>> grid on
>> ylabel('x(t)')
>> xlabel('t')
>> title('Pulso seno')

Nota: La versión 2014 de MATLAB no ejecuta este código, da un mensaje de error 'Division by zero'. Sin embargo, es ejecutado sin problemas por la versión 2007