Ecuación de ondas en el cable de una estructura civil (Grupo 9-B)

De MateWiki
Saltar a: navegación, buscar
Trabajo realizado por estudiantes
Título Ecuación de ondas en el cable de una estructura civil (Grupo 9-B)
Asignatura Ecuaciones Diferenciales
Curso Curso 2013-14
Autores

María García Fernández

José Manuel Corbí Garrido

Estefanía Jiménez Ocampo

Diego García Vaquero

José Francisco Aguilera Corral

Este artículo ha sido escrito por estudiantes como parte de su evaluación en la asignatura


1 Introducción y modelización del problema

En este artículo analizaremos el comportamiento de un cable tenso de 10 metros de longitud en una estructura civil. El cable se encontrará sujeto por los extremos y será sometido a pequeñas vibraciones. Para estudiar el problema haremos ciertas hipótesis que simplificarán el problema, y que en muchos casos dan una buena aproximación del comportamiento real que tendría un cable como el que aquí estudiamos.

En cuanto a las características mecánicas del cable haremos las siguientes hipótesis:

  • El cable tendrá una sección despreciable frente a su longitud.
  • Será un cable homogéneo, es decir, la masa por unidad de volumen será constante.
  • Dicho cable también será perfectamente flexible, es decir, solo ofrecerá resistencia a esfuerzos tangenciales (en su dirección longitudinal), siendo por ello nula su resistencia frente a esfuerzos de flexión y corte.
  • El cable al someterlo a tracción trata de ser inextensible, es decir, que intenta ser lo suficientemente rígido (en dirección longitudinal) como para poder despreciar su extensibilidad, por el contrario al someterlo a compresión no ofrece resistencia y se arruga.

Por otra parte el cable será sometido a pequeñas vibraciones, las cuales podemos modelizar mediante la ecuación de ondas [math]u_{tt}-u_{xx}=0[/math]. Estas vibraciones, al ser pequeñas, son despreciables frente a su longitud total; preocupándonos únicamente por los desplazamientos verticales (fruto de la tracción) que sufre el cable. Estos desplazamientos verticales serán modelizados mediante la función [math]u(x,t)[/math], dependiente de dos variables; la [math]x[/math] que indica la distancia al extremo de la cuerda tomado como referencia, y la [math]t[/math] que define el tiempo.

Las condiciones de frontera en este caso nos indican que el cable no sufrirá desplazamiento vertical en sus extremos [math]u(0,t)=0; u(10,t)=0[/math]. Por otra parte las condiciones iniciales [math]u(x,0)=0; u_{t}(x,0)=0[/math] se traducen en que inicialmente el desplazamiento vertical y la velocidad del cable en todos sus puntos es cero.

El problema que se nos plantea es el siguiente: [math] \begin{cases} u_{tt}-u_{xx}=0; x∈[0,10]; t∈[0,40]\\ u(0,t)=0; u(10,t)=0\\ u(x,0)=0; u_{t}(x,0)=0\\ \end{cases} [/math]

En los sucesivos apartados estudiaremos el comportamiento del cable en diferentes situaciones, que se traducirán en diferentes condiciones para el problema (cambios en las condiciones y en la ecuación de ondas con la que modelizamos las pequeñas vibraciones).

2 Desplazamiento vertical del cable

Sujetando el cable desde el centro y desplazándolo dos metros en la dirección perpendicular, estudiaremos el comportamiento del cable una vez lo soltemos. Estará sometido a unas pequeñas vibraciones, que originarán desplazamientos verticales a lo largo del cable que aproximaremos mediante la función [math]u(x,t)[/math]. Teniendo en cuenta el desplazamiento vertical de dos metros y que el cable es soltado con una velocidad nula, el problema que vamos a estudiar es el siguiente:

[math] \begin{cases} u_{tt}-u_{xx}=0; x∈[0,10]; t∈[0,40]\\ u(0,t)=0; u(10,t)=0\\ u(x,0)=2-|2-\frac{2}{5}x|; u(5,0)=2\\ u_{t}(x,0)=0\\ \end{cases} [/math]

A continuación estudiaremos el desplazamiento vertical al que se ve cada punto del cable expuesto y para cada instante de tiempo. Para ello aproximaremos la ecuación diferencial que modeliza el comportamiento de dicho cable mediante los métodos de diferencias finitas y mediante series de Fourier, obteniendo gráficas que reflejan la posición de cada punto del cable en cada instante.

2.1 Aproximación mediante diferencias finitas

Dentro del método de diferencias finitas se nos plantea un sistema de ecuaciones diferenciales, que resolveremos aplicando diferentes métodos de los cuales luego compararemos su precisión a la hora de aproximar el comportamiento del cable y su estabilidad.

El primer método que utilizaremos para resolver el nuestro problema de ondas será el método del trapecio, el cual es incondicionalmente estable. Para implementar este método y observar el comportamiento del cable en una gráfica usaremos el siguiente código en matlab: 

%% Aproximación del desplazamiento vertical del cable (u) mediante diferencias finitas con el método del trapecio

%% Enunciado del problema

% u_{tt}-u_{xx}=0; x?[0,10]; t>0 
% u(0,t)=0; u(10,t)=0 
% u(x,0)=2-abs(2-\frac{2}{5}x); u(5,0)=2\\  u_{t}(x,0)=0\\

%% Planteamiento y discretización

clear all
L=10;
T=40;
h=0.1;
N=L/h;
x=0:h:L;
xint=h:h:L-h;

dt=h;
t=0:dt:T;

k=2*diag(ones(1,N-1))-diag(ones(1,N-2),1)-diag(ones(1,N-2),-1);
k=(1/h^2)*k;;

F=zeros(N-1,1);

u0=(2-abs(2-xint*2/5))';
v0=0*xint';

M=[zeros(N-1) eye(N-1); -k zeros(N-1) ];
G=[zeros(N-1,1); F];

%% Aproximación de u mediante el método del trapecio

w0=[u0 ; v0];
ww=w0;
sol(1,:)=[ 0 , w0(1:N-1)' , 0];

for j=1:length(t)-1
    ww=(eye(2*N-2)-dt/2*M)\((eye(2*N-2)+dt/2*M)*ww+dt*G);
    
    sol(j+1,:)=[0 , ww(1:N-1)' , 0];
end

%% Representación gráfica

figure
[xx,tt]=meshgrid(x,t);
surf(xx,tt,sol)
xlabel('Cable')
ylabel('Tiempo')
zlabel('Desplazamiento vertical')
title('u aproximado mediante diferencias finitas con el método del trapecio')


Obteniendo la gráfica:

Ecuación de ondas para el cable de una estructura civil. Aproximación mediante diferencias finitas con el método del trapecio

Usaremos también el método de Euler, el cual es un método explícito y por tanto dependemos de la ecuación para que pueda ser o no factible de ejecutar. El siguiente código matlab nos permite ejecutarlo: 

%% Aproximación del desplazamiento vertical del cable (u) mediante diferencias finitas con el método de Euler explícito

%% Enunciado del problema

% u_{tt}-u_{xx}=0; x?[0,10]; t>0 
% u(0,t)=0; u(10,t)=0 
% u(x,0)=2-abs(2-\frac{2}{5}x); u(5,0)=2\\  u_{t}(x,0)=0\\

%% Planteamiento y discretización

clear all
L=10;
T=40;
h=0.1;
N=L/h;
x=0:h:L;
xint=h:h:L-h;

dt=h;
t=0:dt:T;

k=2*diag(ones(1,N-1))-diag(ones(1,N-2),1)-diag(ones(1,N-2),-1);
k=(1/h^2)*k;

F=zeros(N-1,1);

u0=(2-abs(2-xint*2/5))';
v0=0*xint';

%% Aproximación de u mediante el método de Euler explícito

sol(1,:)=[ 0 , u0' , 0];
u=u0;
v=v0;
for j=1:length(t)-1
    u=u+dt*v+F;
    v=v-dt*k*u+F;
    sol(j+1,:)=[0 , u' , 0];
end

%% Representación gráfica

[xx,tt]=meshgrid(x,t);
surf(xx,tt,sol)
xlabel('Cable')
ylabel('Tiempo')
zlabel('Desplazamiento vertical')
title('u aproximado mediante diferencias finitas con el método de Euler explícito')


Obteniendo el siguiente gráfico del comportamiento del cable respecto al tiempo:

Ecuación de ondas para el cable de una estructura civil. Aproximación mediante diferencias finitas con el método de Euler explícito

El último método que empleamos para resolver el problema planteado mediante diferencias finitas es el método de Euler modificado, el cual tiene más precisión que el método de Euler explícito. El código que se muestra a continuación ejecuta dicho método y nos proporciona una gráfica del cable: 

%% Aproximación del desplazamiento vertical del cable (u) mediante diferencias finitas con el método de Euler modificado

%% Enunciado del problema

% u_{tt}-u_{xx}=0; x?[0,10]; t>0 
% u(0,t)=0; u(10,t)=0 
% u(x,0)=2-abs(2-\frac{2}{5}x); u(5,0)=2\\  u_{t}(x,0)=0\\

%% Planteamiento y discretización

clear all
L=10;
T=40;
h=0.1;
N=L/h;
x=0:h:L;
xint=h:h:L-h;

dt=h;
t=0:dt:T;

k=2*diag(ones(1,N-1))-diag(ones(1,N-2),1)-diag(ones(1,N-2),-1);
k=(1/h^2)*k;

F=zeros(N-1,1);

u0=(2-abs(2-xint*2/5))';
v0=0*xint';

%% Aproximación de u mediante el método de Euler modificado

sol(1,:)=[ 0 , u0' , 0];
u=u0;
v=v0;
for n=1:length(t)-1
  l=k*u;
    % calculamos k1=f(t_n,y_n)
    k1=[ v l]';
    % Calculamos k2=f(t_n+h/2,y_n+1/2k1*h)
    %xp=x(n)+dt; %argumento de tiempo
    %solp=sol+dt*k1; %argumento y 
    k2=[v l]';
    %Calculamos k3=f(t_n+h/2,y_n+1/2k2*h)
    %mismo argumento tiempo que antes
    
    
    u=u+dt/2*(k1(1,:)'+k2(1,:)')+F;
    l=+k*u;
    k1=[v l]';
    k2=[v l]';
    v=v-dt/2*(k1(2,:)'+k2(2,:)')+F;
   % v=v-dt/2*(k*u+k*u)+F;
   sol(n+1,:)=[0 , u' , 0];
     
end  

%% Representación gráfica

[xx,tt]=meshgrid(x,t);
surf(xx,tt,sol)
xlabel('Cable')
ylabel('Tiempo')
zlabel('Desplazamiento vertical')
title('u aproximado mediante diferencias finitas con el método de Euler modificado')

La gráfica que se obtiene es:

Ecuación de ondas para el cable de una estructura civil. Aproximación mediante diferencias finitas con el método de Euler modificado

2.2 Aproximación mediante el método de Fourier

2.2.1 Con un término de la serie

%% Aproximación del desplazamiento vertical del cable (u) mediante Fourier con un término

%% Enunciado del problema

% u_{tt}-u_{xx}=0; x?[0,10]; t>0 
% u(0,t)=0; u(10,t)=0 
% u(x,0)=2-abs(2-\frac{2}{5}x); u(5,0)=2\\  u_{t}(x,0)=0\\

%% Planteamiento y discretización

clear all
L=10;
T=40;
h=0.1;
N=L/h;
dt=0.1; 
M=T/dt;
x=0:h:10;
posicion=(2-abs(2-x*2/5));
velocidad=0*x;
t=0:dt:T;

%% Aproximación de u mediante Fourier con un término (Q=1)

Q=1;
for k=1:Q
p=sin(k*pi/10*x) ;
fourier1(k)=trapz(x,posicion.*p)/trapz(x, p.*p);

end
[XX,TT]=meshgrid(x,t);
U=0;

for k=1:Q
U=U+(fourier1(k)*cos(k*pi/10*TT)).*sin(k*pi/10*XX);
end

%% Representación gráfica

surf(XX,TT,U)
xlabel('Cable')
ylabel('Tiempo')
zlabel('Desplazamiento vertical')
title('u aproximado mediante Fourier con un término')

Ecuación de ondas para el cable de una estructura civil. Aproximación mediante Fourier con un término

2.2.2 Con tres términos de la serie

%% Aproximación del desplazamiento vertical del cable (u) mediante Fourier con tres términos

%% Enunciado del problema

% u_{tt}-u_{xx}=0; x?[0,10]; t>0 
% u(0,t)=0; u(10,t)=0 
% u(x,0)=2-abs(2-\frac{2}{5}x); u(5,0)=2\\  u_{t}(x,0)=0\\

%% Planteamiento y discretización

clear all
L=10;
T=40;
h=0.1;
N=L/h;
dt=0.1; 
M=T/dt;
x=0:h:10;
posicion=(2-abs(2-x*2/5));
velocidad=0*x;
t=0:dt:T;

%% Aproximación de u mediante Fourier con tres términos (Q=3)

Q=3;

for k=1:Q
    p=sin(k*pi/10*x) ;
    fourier1(k)=trapz(x,posicion.*p)/trapz(x,p.*p);
 end

[XX,TT]=meshgrid(x,t);
U=0;

for k=1:Q
    U=U+(fourier1(k)*cos(k*pi/10*TT)).*sin(k*pi/10*XX);
end

%% Representación gráfica

surf(XX,TT,U)
xlabel('Cable')
ylabel('Tiempo')
zlabel('Desplazamiento vertical')
title('u aproximado mediante Fourier con tres términos')

Ecuación de ondas para el cable de una estructura civil. Aproximación mediante Fourier con tres términos

2.2.3 Con cinco términos de la serie

%% Aproximación del desplazamiento vertical del cable (u) mediante Fourier con cinco términos

%% Enunciado del problema

% u_{tt}-u_{xx}=0; x?[0,10]; t>0 
% u(0,t)=0; u(10,t)=0 
% u(x,0)=2-abs(2-\frac{2}{5}x); u(5,0)=2\\  u_{t}(x,0)=0\\

%% Planteamiento y discretización

clear all
L=10;
T=40;
h=0.1;
N=L/h;
dt=0.1; 
M=T/dt;
x=0:h:10;
posicion=(2-abs(2-x*2/5));
velocidad=0*x;
t=0:dt:T;

%% Aproximación de u mediante Fourier con cinco términos (Q=5)

Q=5;

for k=1:Q
    p=sin(k*pi/10*x) ;
    fourier1(k)=trapz(x,posicion.*p)/trapz(x, p.*p);
end

[XX,TT]=meshgrid(x,t);
U=0;

for k=1:Q
    U=U+(fourier1(k)*cos(k*pi/10*TT)).*sin(k*pi/10*XX);
end

%% Representación gráfica

surf(XX,TT,U)
xlabel('Cable')
ylabel('Tiempo')
zlabel('Desplazamiento vertical')
title('u aproximado mediante Fourier con cinco términos')

Ecuación de ondas para el cable de una estructura civil. Aproximación mediante Fourier con cinco términos

2.2.4 Con diez términos de la serie

%% Aproximación del desplazamiento vertical del cable (u) mediante Fourier con diez términos

%% Enunciado del problema

% u_{tt}-u_{xx}=0; x?[0,10]; t>0 
% u(0,t)=0; u(10,t)=0 
% u(x,0)=2-abs(2-\frac{2}{5}x); u(5,0)=2\\  u_{t}(x,0)=0\\

%% Planteamiento y discretización

clear all
L=10;
T=40;
h=0.1;
N=L/h;
dt=0.1; 
M=T/dt;
x=0:h:10;
posicion=(2-abs(2-x*2/5));
velocidad=0*x;
t=0:dt:T;

%% Aproximación de u mediante Fourier con diez términos (Q=10)

Q=10;

for k=1:Q
    p=sin(k*pi/10*x) ;
    fourier1(k)=trapz(x,posicion.*p)/trapz(x, p.*p);
end

[XX,TT]=meshgrid(x,t);
U=0;

for k=1:Q
    U=U+(fourier1(k)*cos(k*pi/10*TT)).*sin(k*pi/10*XX);
end

%% Representación gráfica

surf(XX,TT,U)
xlabel('Cable')
ylabel('Tiempo')
zlabel('Desplazamiento vertical')
title('u aproximado mediante Fourier con diez términos')

Ecuación de ondas para el cable de una estructura civil. Aproximación mediante Fourier con diez términos

2.2.5 Con veinte términos de la serie

%% Aproximación del desplazamiento vertical del cable (u) mediante Fourier con veinte términos

%% Enunciado del problema

% u_{tt}-u_{xx}=0; x?[0,10]; t>0 
% u(0,t)=0; u(10,t)=0 
% u(x,0)=2-abs(2-\frac{2}{5}x); u(5,0)=2\\  u_{t}(x,0)=0\\

%% Planteamiento y discretización

clear all
L=10;
T=40;
h=0.1;
N=L/h;
dt=0.1; 
M=T/dt;
x=0:h:10;
posicion=(2-abs(2-x*2/5));
velocidad=0*x;
t=0:dt:T;

%% Aproximación de u mediante Fourier con diez términos (Q=20)

Q=20;

for k=1:Q
    p=sin(k*pi/10*x) ;
    fourier1(k)=trapz(x,posicion.*p)/trapz(x, p.*p);
    
end

[XX,TT]=meshgrid(x,t);
U=0;

for k=1:Q
    U=U+(fourier1(k)*cos(k*pi/10*TT)).*sin(k*pi/10*XX);
end

%% Representación gráfica

surf(XX,TT,U)
xlabel('Cable')
ylabel('Tiempo')
zlabel('Desplazamiento vertical')
title('u aproximado mediante Fourier con veinte términos')

Ecuación de ondas para el cable de una estructura civil. Aproximación mediante Fourier con veinte términos

2.3 Comparativa de diferencias finitas y Fourier. Interpretación de las gráficas anteriores

Todas las gráficas obtenidas mediante diferencias finitas (con sus distintos métodos) y Fourier representan la aproximación numérica de la solución del problema planteado, es decir, el desplazamiento vertical de todos los puntos de nuestro cable a lo largo del tiempo.

Considerando los órdenes de precisión de cada uno de los métodos empleados en los programas de diferencias finitas, para aproximar la solución de nuestro problema (la función [math]u(x,t)[/math] que representa el desplazamiento vertical): Podemos concluir que el método del trapecio (orden 2) es el más preciso, seguido del método de Euler modificado y por último del método de Euler explícito (orden 1).

En cuanto a la resolución del problema mediante Fourier, podemos observar que a medida que aumentamos el número de términos de la serie de Fourier, la forma de la solución pasa de ser más redondeada a más recta, pareciéndose cada vez más a la solución obtenida por el método de las diferencias finitas. El aumento del número de términos considerados en Fourier hace el método más estable.

3 Energía del cable

La energía mecánica del cable viene definida por la siguiente función:: [math] E(t) = \int_{0}^{L} (u_{t}^{2}(x,t)+ u_{x}^{2}(x,t)) \cdot dx [/math] Mediante el programa que se muestra a continuación procederemos a dibujar la gráfica de la función [math]E(t)[/math], que cuantifica la energía del cable. Previamente, y recurriendo al concepto de derivada, aproximaremos los valores que toman [math]u_{t}[/math] y [math]u_{x}[/math], para así después poder calcular la integral que define [math]E(t)[/math]. 

%% Aproximación del desplazamiento vertical del cable (u) mediante diferencias finitas con el método del trapecio

%% Enunciado del problema

% u_{tt}-u_{xx}=0; x?[0,10]; t>0 
% u(0,t)=0; u(10,t)=0 
% u(x,0)=2-abs(2-\frac{2}{5}x); u(5,0)=2\\  u_{t}(x,0)=0\\

%% Planteamiento y discretización

clear all
L=10;
T=40;
h=0.05;
N=L/h;
x=0:h:L;
xint=h:h:L-h;

dt=h;
t=0:dt:T;

k=2*diag(ones(1,N-1))-diag(ones(1,N-2),1)-diag(ones(1,N-2),-1);
k=(1/h^2)*k;

F=zeros(N-1,1);

u0=(2-abs(2-xint*2/5))';
v0=0*xint';

M=[zeros(N-1) eye(N-1); -k zeros(N-1) ];
G=[zeros(N-1,1); F];

%% Aproximación de u mediante el método del trapecio

w0=[u0 ; v0];
ww=w0;
sol(1,:)=[ 0 , w0(1:N-1)' , 0];

for j=1:length(t)-1
    ww=(eye(2*N-2)-dt/2*M)\((eye(2*N-2)+dt/2*M)*ww+dt*G);
    sol(j+1,:)=[0 , ww(1:N-1)' , 0];
end

[xx,tt]=meshgrid(x,t);
surf(xx,tt,sol)

%% Aproximación de u_t y u_x basándonos en el concepto de derivada

v=size(sol);
for i=2:v(1)-1
  for j=2:v(2)-1
    Ux(i,j)=(sol(i+1,j)-sol(i-1,j))/(2*h);
    Ut(i,j)=(sol(i,j+1)-sol(i,j-1))/(2*h);
  end
end

%% Cálculo de la energía

e=Ux.^2+Ut.^2;
xaux=[0,xint];

for k=1:v(1)-2
    E(k)=trapz(xaux,e(k,:));
end
t(401)=[ ];

%% Representación gráfica de la energía a lo largo del tiempo

plot(t(2:end-1),E(2:end))
xlabel('Tiempo')
ylabel('Energía')
title('Energía a lo largo del tiempo')


Energía del cable a lo largo del tiempo)


Se observa que el valor de la energía total del cable, para un paso del tiempo de 0.05 segundos, oscila aproximadamente entre 1.560 y 1.590. Este rango se reduce si elegimos un paso más pequeño para el tiempo y el espacio. Por tanto, aunque la aproximación no es perfecta, sí podemos confirmar que la energía mecánica del sistema se conserva.

4 Ejemplos de aplicación de cables en estructuras civiles

4.1 Cable sumergido en el mar (medio viscoso)

En este caso el cable está sumergido en un medio viscoso que produce amortiguamiento. El problema que se nos plantea es el siguiente: [math] \begin{cases} u_{tt}-u_{xx}+au_{t}=0; x∈[0,10]; t∈[0,40] \\ u(0,t)=0; u(10,t)=0\\ u(x,0)=0; u_{t}(x,0)=0\\ \end{cases} [/math]

Vamos a estudiar lo que ocurre con la energía en estas condiciones, dependiendo del parámetro [math]a[/math]. Damos a [math]a[/math] los valores de 0, 1 ,4, 10 y 100. Para ello utilizamos el método numérico de las diferencias finitas empleando Euler con un paso muy pequeño. En el caso de que a valga cero estamos ante el mismo caso que si no estuviese sumergido. 

%% Cable sumergido en un medio viscoso (p.ej el mar)con a=0
% Cambiando el valor de a definido más adelante obtenemos todas las gráficas
%% Enunciado del problema

% u_{tt}-u_{xx}+a*u_{t}=0; x[0,10]; t[0,40]; 
% u(0,t)=0; u(10,t)=0 
% u(x,0)=0; u(5,0)=2\\  u_{t}(x,0)=0\\

%% Planteamiento y discretización

clear all
L=10;
T=40;
h=0.1;
N=L/h;
x=0:h:L;
xint=h:h:L-h;

dt=h/2;
t=0:dt:T;

k=2*diag(ones(1,N-1))-diag(ones(1,N-2),1)-diag(ones(1,N-2),-1);

k=(1/h^2)*k;

F=zeros(N-1,1);

u0=(2-abs(2-xint*2/5))';
v0=0*xint';


%% Aproximación de u mediante el método de Euler y obtención de las derivadas

sol(1,:)=[ 0 , u0' , 0];
u=u0;
v=v0;
a=0;
for j=1:length(t)-1
    u=u+dt*v;
    v=v-dt*k*u-dt*a*v;
    sol(j+1,:)=[0 , u' , 0];
end

%% Gráficas

[xx,tt]=meshgrid(x,t);

figure

surf(xx,tt,sol);
xlabel('Cable')
ylabel('Tiempo')
zlabel('Desplazamiento vertical')
title('u aproximado mediante diferencias finitas con el método de Euler explícito')

figure

ffff=sol(1,:);
plot(x,ffff);
xlabel('Cable')
ylabel('Desplazamiento inicial del cable')
title('Desplazamiento inicial de los puntos del cable')

v=size(sol);
for i=1:v(1)-1
  for j=1:v(2)-1
    Ux(i,j)=(sol(i+1,j)-sol(i,j))/h;
    Ut(i,j)=(sol(i,j+1)-sol(i,j))/h;
  end
end

e=Ux.^2+Ut.^2;
xaux=[0,xint];
laux=length(xaux);
otro=length(e(2,:));
for k=1:v(1)-1
  E(k)=trapz(xaux,e(k,:));
 end
t(401)=[ ];
lE=length(E);
lt=length(t);

figure

plot(t,E)
xlabel('Tiempo')
ylabel('Energía')
title ('Energía con a=0')


Cable sumergido en un medio viscoso (p.ej el mar)con a=0 Cable sumergido en un medio viscoso (p.ej el mar)con a=1 Cable sumergido en un medio viscoso (p.ej el mar)con a=4 Cable sumergido en un medio viscoso (p.ej el mar)con a=10 Cable sumergido en un medio viscoso (p.ej el mar)con a=100

Observamos que, al sumergir el cable en un medio viscoso, el desplazamiento vertical de las ondas es cada vez más pequeño. Esto se debe a que la energía disminuye con el tiempo. Además, al aumentar el valor de [math]a[/math], la energía se disipa más rápidamente. Como el medio es más viscoso, la oposición al desplazamiento vertical del cable es mayor. Para el caso de a=100, el método no es estable y da problemas. Sería necesario reducir aún más el paso temporal y espacial o utilizar otro método de mayor precisión.

4.2 Cable sujeto a una estructura sometida a vibraciones periódicas

Ahora nos preguntamos qué ocurriría si el extremo derecho del cable está sujeto a una estructura que sufre vibraciones periódicas. La frecuencia de estas vibraciones es de [math]F_{0}[/math] Herzios. Tomamos por ejemplo [math]f(t) = sin(2 \pi F_{0})[/math]. Vamos a calcular el comportamiento de la energía cuando el cable parte inicialmente del reposo y tomamos un tiempo t. Repetiremos el proceso con distintos valores de [math]F_{0}[/math] para observar si se produce alguna diferencia.

Tomamos [math]F_{0}[/math] igual a 1/L+0,01 Hz, 1/L-0,01 Hz y 1/L Hz.

El problema nos quedaría planteado de la siguiente manera:

[math] \begin{cases} u_{tt}-u_{xx}=0; x∈[0,10]; t∈[0,60]\\ u(0,t)=0; u_{x}(10,t)=sin(2 \pi F_{0})\\ u(x,0)=0; u_{t}(x,0)=0\\ \end{cases} [/math]

Lo resolvemos numéricamente por diferencias finitas utilizando Euler. 

%% Cable sujeto a una estructura sometida a vibraciones periódicas

%% Enunciado del problema

% u_{tt}-u_{xx}=0; x[0,10]; t[0,60] 
% u(0,t)=0; u(10,t)=sin(2*pi*f0)
% u(x,0)=0; u(5,0)=2\\  u_{t}(x,0)=0

%% Planteamiento y discretización (f0=0.11)
 
clear all
L=10;
T=60;
h=0.1;
N=L/h;
x=0:h:L;
xint=h:h:L;
 
dt=h;
t=0:dt:T;
 
k=2*diag(ones(1,N))-diag(ones(1,N-1),1)-diag(ones(1,N-1),-1);
k(N,N-1)=-2;
 
k=(1/h^2)*k;
 
F=zeros(N-1,1);
 
u0=0*xint';
v0=0*xint';
 
%% Aproximación de u mediante el método de Euler explícito y obtención de las derivadas (f0=0.11)
 
f0=0.11;
sol(1,:)=[ 0 , u0'];
u=u0;
v=v0;
for j=1:length(t)-1
    u=u+dt*v+[F;sin(2*pi*f0*j)];
    v=v-dt*k*u+[F;sin(2*pi*f0*j)];
    sol(j+1,:)=[0 , u' ];
end
 
[xx,tt]=meshgrid(x,t);
figure(1)
surf(xx,tt,sol)
 
v=size(sol);
for i=1:v(1)-1
  for j=1:v(2)-1
    Ux(i,j)=(sol(i+1,j)-sol(i,j))/h;
    Ut(i,j)=(sol(i,j+1)-sol(i,j))/h;
  end
end

%% Energía (f0=0.11)

e=Ux.^2+Ut.^2;
xaux=[xint];
laux=length(xaux);
otro=length(e(2,:));

for k=1:v(1)-1
    E(k)=trapz(xaux,e(k,:));
end

t(401)=[ ];

%% Representación Energía-Tiempo (f0=0.11)

figure(2)
plot(t,E)
xlabel('Tiempo')
ylabel('Energía')
title('Gráfica energía-tiempo f0=0.11')

%% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 
%% Planteamiento y discretización (f0=0.09)
 
clear all
L=10;
T=60;
h=0.1;
N=L/h;
x=0:h:L;
xint=h:h:L;
 
dt=h;
t=0:dt:T;
 
k=2*diag(ones(1,N))-diag(ones(1,N-1),1)-diag(ones(1,N-1),-1);
k(N,N-1)=-2;
 
k=(1/h^2)*k;
 
F=zeros(N-1,1);
 
u0=0*xint';
v0=0*xint';
 
%% Aproximación de u mediante el método de Euler explícito y obtención de las derivadas (f0=0.09)
 
f0=0.09;
sol(1,:)=[ 0 , u0'];
u=u0;
v=v0;

for j=1:length(t)-1
    u=u+dt*v+[F;sin(2*pi*f0*j)];
    v=v-dt*k*u+[F;sin(2*pi*f0*j)];
    sol(j+1,:)=[0 , u' ];
end
 
[xx,tt]=meshgrid(x,t);
figure(3)
surf(xx,tt,sol)


v=size(sol);
for i=1:v(1)-1
  for j=1:v(2)-1
    Ux(i,j)=(sol(i+1,j)-sol(i,j))/h;
    Ut(i,j)=(sol(i,j+1)-sol(i,j))/h;
  end
end

%% Energía (f0=0.09)

e=Ux.^2+Ut.^2;
xaux=[xint];
laux=length(xaux);
otro=length(e(2,:));
for k=1:v(1)-1
    E(k)=trapz(xaux,e(k,:));
end

t(401)=[ ];

%% Representación Energía-Tiempo (f0=0.09)

figure(4)
plot(t,E)
xlabel('Tiempo')
ylabel('Energía')
title('Gráfica energía-tiempo f0=0.09')
%% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

%% Planteamiento y discretización (f0=0.1)

clear all
L=10;
T=60;
h=0.1;
N=L/h;
x=0:h:L;
xint=h:h:L;
 
dt=h;
t=0:dt:T;
 
k=2*diag(ones(1,N))-diag(ones(1,N-1),1)-diag(ones(1,N-1),-1);
k(N,N-1)=-2;
 
k=(1/h^2)*k;
 
F=zeros(N-1,1);

u0=0*xint';
v0=0*xint';
 
%% Aproximación de u mediante el método de Euler explícito y obtención de las derivadas (f0=0.1)
 
f0=0.1;
sol(1,:)=[ 0 , u0'];
u=u0;
v=v0;
for j=1:length(t)-1
    u=u+dt*v+[F;sin(2*pi*f0*j)];
    v=v-dt*k*u+[F;sin(2*pi*f0*j)];
    sol(j+1,:)=[0 , u' ];
end
 
[xx,tt]=meshgrid(x,t);
figure(5)
surf(xx,tt,sol)
v=size(sol);

for i=1:v(1)-1
  for j=1:v(2)-1
    Ux(i,j)=(sol(i+1,j)-sol(i,j))/h;
    Ut(i,j)=(sol(i,j+1)-sol(i,j))/h;
  end
end

%% Energía (f0=0.1)

e=Ux.^2+Ut.^2;
xaux=[xint];
laux=length(xaux);
otro=length(e(2,:));

for k=1:v(1)-1
    E(k)=trapz(xaux,e(k,:));
end

t(401)=[ ];

%% Representación Energía-Tiempo (f0=0.1)

figure(6)
plot(t,E)
xlabel('Tiempo')
ylabel('Energía')
title('Gráfica energía-tiempo f0=0.1')

Cable sujeto a una estructura sometida a vibraciones periódicas f0=0.11 Cable sujeto a una estructura sometida a vibraciones periódicas f0=0.09 Cable sujeto a una estructura sometida a vibraciones periódicas f0=0.1

Se observa que cuanto mayor es el valor de [math]F_{0}[/math], la energía alcanza un máximo más elevado. Esto es debido a que las vibraciones en el extremo derecho de la estructura dependen de una función senoidal que está en función de [math]F_{0}[/math]. A medida que aumenta [math]F_{0}[/math] crece el seno y, con ello, la derivada parcial en x, esto se traduce en un incremento de la energía. El valor de [math]F_{0}[/math] únicamente afecta a los máximos de la energía pero no a su carácter temporal. Se trata de una función periódica que depende del tiempo. Para cualquier [math]F_{0}[/math], el periodo es de 40 segundos.

4.3 Cable sujeto a un aparato que envía una respuesta a las vibraciones

En este apartado vamos a calcular el comportamiento de la energía para el valor b suponiendo ahora que el extremo derecho del cable está sujeto a un aparato que envía una respuesta a la vibración que recibe de manera que la condición de contorno ahora es [math]u_{x}(L,t)=bu(L,t)[/math]. Tomamos como condición inicial [math]u(x,0)=2-|2-\frac{2}{5}x|[/math] y observando en qué cambia la energía en función de distintos valores de b. estudiamos dos casos concretos: b=-2 y b=2.

[math] \begin{cases} u_{tt}-u_{xx}=0; x∈[0,10]; t∈[0,40]\\ u(0,t)=0; u_{x}(10,t)=b u(10,0) \\ u(x,0)=2-|2-\frac{2}{5}x|; u(5,0)=2\\ u_{t}(x,0)=0\\ \end{cases} [/math]

Resolvemos este problema por diferencias finitas utilizando el método de Euler, con un paso pequeño para que la aproximación sea más exacta. 


%% Cable sujeto a un aparato que envía una respuesta a las vibraciones (b=2)

%% Enunciado del problema

% u_{tt}-u_{xx}=0; x[0,10]; t[0,40]; 
% u(0,t)=0; u(10,t)=bu(10,0) 
% u(x,0)=2-abs(2-\frac{2}{5}x); u(5,0)=2\\  u_{t}(x,0)=0\\

%% Planteamiento y discretización

clear all
b=2;
L=10;
T=40;
h=0.01;
N=L/h;
x=0:h:L;
xint=h:h:L;

dt=h/2;
t=0:dt:T;

k=2*diag(ones(1,N))-diag(ones(1,N-1),1)-diag(ones(1,N-1),-1);
k(N,N-1)=-2;
k(N,N)=2+b*(2*h);

k=(1/h^2)*k;

F=zeros(N,1);

u0=(2-abs(2-xint*2/5))';
v0=0*xint';

%% Aproximación de u mediante el método de Euler y obtención de las derivadas

f0=0.1;
sol(1,:)=[ 0 , u0'];
u=u0;
v=v0;
for j=1:length(t)-1
    u=u+dt*v;
    v=v-dt*k*u;
    sol(j+1,:)=[0 , u' ];
end

[xx,tt]=meshgrid(x,t);
figure(1)
surf(xx,tt,sol)
  
v=size(sol);
for i=2:v(1)-1
  for j=2:v(2)-1
    Ux(i,j)=(sol(i+1,j)-sol(i-1,j))/(2*h);
    Ut(i,j)=(sol(i,j+1)-sol(i,j-1))/(2*h);
  end
end
e=Ux.^2+Ut.^2;
xaux=[0,xint];

%% Energía

for k=1:v(1)-2
    E(k)=trapz(xaux(2:end),e(k,:));
end
t(401)=[ ];
plot(t(2:end-1),E(2:end))
xlabel('Tiempo')
ylabel('Energía')
title('Gráfica energía-tiempo para b=2')

Cable sujeto a un aparato que envía una respuesta a las vibraciones (b=2)

No podemos apreciar para valores de b positivos (como el estudiado, b=2) cambios importantes en la energía. Esta no aumenta ni disminuye de manera significativa. 

%% Cable sujeto a un aparato que envía una respuesta a las vibraciones (b=-2)

%% Enunciado del problema

% u_{tt}-u_{xx}=0; x[0,10]; t[0,40]; 
% u(0,t)=0; u(10,t)=bu(10,0) 
% u(x,0)=2-abs(2-\frac{2}{5}x); u(5,0)=2\\  u_{t}(x,0)=0\\

%% Planteamiento y discretización

clear all
b=-2;
L=10;
T=40;
h=0.01;
N=L/h;
x=0:h:L;
xint=h:h:L;

dt=h/2;
t=0:dt:T;

k=2*diag(ones(1,N))-diag(ones(1,N-1),1)-diag(ones(1,N-1),-1);
k(N,N-1)=-2;
k(N,N)=2+b*(2*h);

k=(1/h^2)*k;

F=zeros(N,1);

u0=(2-abs(2-xint*2/5))';
v0=0*xint';

%% Aproximación de u mediante el método de Euler y obtención de las derivadas

f0=0.1;
sol(1,:)=[ 0 , u0'];
u=u0;
v=v0;
for j=1:length(t)-1
    u=u+dt*v;
    v=v-dt*k*u;
    sol(j+1,:)=[0 , u' ];
end

[xx,tt]=meshgrid(x,t);
figure(1)
surf(xx,tt,sol)
  
v=size(sol);
for i=2:v(1)-1
  for j=2:v(2)-1
    Ux(i,j)=(sol(i+1,j)-sol(i-1,j))/(2*h);
    Ut(i,j)=(sol(i,j+1)-sol(i,j-1))/(2*h);
  end
end
e=Ux.^2+Ut.^2;
xaux=[0,xint];

%% Energía

for k=1:v(1)-2
    E(k)=trapz(xaux(2:end),e(k,:));
end
t(401)=[ ];
plot(t(2:end-1),E(2:end))
xlabel('Tiempo')
ylabel('Energía')
title('Gráfica energía-tiempo para b=-2')

Cable sujeto a un aparato que envía una respuesta a las vibraciones (b=-2)

Se observa que cuanto mayor es la tensión del cable, mayor es el desplazamiento vertical. Se está otorgando energía al sistema. De está forma, para valores de b negativos (por ejemplo b=-2), la energía crece de manera exponencial y tiende a infinito. Se produce un efecto de resonancia, incrementándose la energía.

5 Visualización del movimiento del cable

Si el siguiente código basado en el método de diferencias finitas (mediante el método del trapecio) analizado anteriormente (en el apartado 2) se ejecuta en MATLAB (Octave no admite animaciones), se puede apreciar en la primera imagen que proporciona el programa la solución al problema enunciado en el apartado 2. Lo siguiente en salir por pantalla es una animación de la solución obtenida a lo largo del tiempo, de forma que cuando se agota el intervalo de tiempo en el que estudiamos el problema, queda la misma imagen que la primera que nos proporciona este programa; con este programa vemos como la onda va "viajando". La tercera imagen en salir por pantalla es la representación del desplazamiento vertical que sufre el cable a lo largo del tiempo, es decir, la vibración de todos los puntos de nuestro cable. 

%% Aproximación del desplazamiento vertical del cable (u) mediante diferencias finitas con el método del trapecio

%% Enunciado del problema

% u_{tt}-u_{xx}=0; x?[0,10]; t>0 
% u(0,t)=0; u(10,t)=0 
% u(x,0)=2-abs(2-\frac{2}{5}x); u(5,0)=2\\  u_{t}(x,0)=0\\

%% Planteamiento y discretización

clear all
L=10;
T=40;
h=0.1;
N=L/h;
x=0:h:L;
xint=h:h:L-h;

dt=h;
t=0:dt:T;

k=2*diag(ones(1,N-1))-diag(ones(1,N-2),1)-diag(ones(1,N-2),-1);
k=(1/h^2)*k;;

F=zeros(N-1,1);

u0=(2-abs(2-xint*2/5))';
v0=0*xint';

M=[zeros(N-1) eye(N-1); -k zeros(N-1) ];
G=[zeros(N-1,1); F];

%% Aproximación de u mediante el método del trapecio

w0=[u0 ; v0];
ww=w0;
sol(1,:)=[ 0 , w0(1:N-1)' , 0];

for j=1:length(t)-1
    ww=(eye(2*N-2)-dt/2*M)\((eye(2*N-2)+dt/2*M)*ww+dt*G);
    
    sol(j+1,:)=[0 , ww(1:N-1)' , 0];
end

%% Representación gráfica

figure
[xx,tt]=meshgrid(x,t);
surf(xx,tt,sol)
xlabel('Cable')
ylabel('Tiempo')
zlabel('Desplazamiento vertical')
title('u aproximado mediante diferencias finitas con el método del trapecio')

%% Animación.

% Muestra el desplazamiento vertical de cada punto del cable en cada instante

figure
[n,m]=size(sol);
aux=zeros(size(sol));
for j=1:n
    % Se va analizando el desplazamiento vertical en cada instante (j) de
    % todos los puntos (i) del cable
    for i=1:m
        aux(j,i)=sol(j,i);
    end
    surf(xx,tt,aux); % Es la representación de todos los puntos (i) del cable en un instante dado
    xlabel('Cable')
    ylabel('Tiempo')
    zlabel('Desplazamiento vertical')
    zlim([-2 2])
    MFRAMES=getframe;   % Guarda una foto del cable en cada instante
end
movie=MFRAMES;

figure
for j=1:n
    % Se va analizando el desplazamiento vertical en cada instante (j) de
    % todos los puntos (i) del cable
    for i=1:m
        aux(j,i)=sol(j,i);
    end
    plot(x,aux(j,:)); % Es la representación de todos los puntos (i) del cable en un instante dado
    xlabel('Cable')
    ylabel('Desplazamiento vertical')
    ylim([-2 2])
    MFRAMES=getframe;   % Guarda una foto del cable en cada instante
end
movie=MFRAMES;
Obtenido de «https://mat.caminos.upm.es/w/index.php?title=Ecuación_de_ondas_en_el_cable_de_una_estructura_civil_(Grupo_9-B)&oldid=13486»