Ecuación de ondas (grupo 2B)
De MateWiki
Revisión del 13:54 17 may 2014 de Adela González (Discusión | contribuciones) (→Aproximación por el método del trapecio.)
| Trabajo realizado por estudiantes | |
|---|---|
| Título | Ecuación de ondas. Grupo 2-B |
| Asignatura | Ecuaciones Diferenciales |
| Curso | Curso 2013-14 |
| Autores |
Ignacio Díaz-Caneja Camblor Alberto Fernández Pérez Adela González Barbado Lucia López Sánchez Araceli Martín Candilejo Diego Solano López |
| Este artículo ha sido escrito por estudiantes como parte de su evaluación en la asignatura | |
1 Deslizamiento Vertical del Cable.
Consideramos un cable de 10 metros de longitud sujeto por sus extremos. Suponiendo que éste tiene una sección pequeña respecto a su longitud someteremos al cable a pequeñas vibraciones que estudiaremos con una modelización a la ecuación de ondas. Se caracteriza al cable de una masa constante por unidad de volumen, es decir, será homogéneo. Éste, flexible e inextensible a la tracción, únicamente ofrece resistencias en su dirección longitudinal, tangenciales, pero no a esfuerzos de flexión o cortes. Para iniciar el movimiento del cable lo sujetaremos por el centro subiéndolo dos metros y soltándolo.
[math] \begin{cases} u_tt- u_tx=0\\ u(0,t)=0\\ u(10,t)=0\\ u(x,0)=g(x)= 2- |\frac{2}{5}x -2|\\ u_t(x,0)=0 \end{cases} [/math]
1.1 Aproximación por el método del trapecio.
%utt-uxx=0
%u(0,t)=u(l,t)=0
%u(x,0)=2-abs(0.4x-2)
%ut(x,0)=0
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=K/h^2;
F=zeros(N-1,1);
u0=2-abs((xint*2/5-2))';
v0=zeros(N-1,1);
M=[zeros(N-1),eye(N-1);-K,zeros(N-1)];
G=[zeros(N-1,1);F];
W0=[u0;v0];
WW=W0;
sol=zeros(length(t),N+1);
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;
sol(j+1,:)=[0,WW(1:N-1)',0];
end
[xx,tt]=meshgrid(x,t);
surf(xx,tt,sol);
1.2 Aproximación por el método de Euler.
1.2.1 Euler Explícito.
%utt-uxx=0
%u(0,t)=u(l,t)=0
%u(x,0)=2-abs(0.4x-2)
%ut(x,0)=0
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=K/h^2;
F=zeros(N-1,1);
u0=2-abs((xint*2/5-2))';
v0=zeros(N-1,1);
M=[zeros(N-1),eye(N-1);-K,zeros(N-1)];
G=[zeros(N-1,1);F];
W0=[u0;v0];
WW=W0;
sol=zeros(length(t),N+1);
sol(1,:)=[0,W0(1:N-1)',0];
for j=1:length(t)-1
WW=WW+dt*M*WW;
sol(j+1,:)=[0,WW(1:N-1)',0];
end
[xx,tt]=meshgrid(x,t);
surf(xx,tt,sol);1.2.2 Euler Modificado.
%utt-uxx=0
%u(0,t)=u(l,t)=0
%u(x,0)=2-abs(0.4x-2)
%ut(x,0)=0
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=K/h^2;
F=zeros(N-1,1);
u0=2-abs((xint*2/5-2))';
v0=zeros(N-1,1);
M=[zeros(N-1),eye(N-1);-K,zeros(N-1)];
G=[zeros(N-1,1);F];
W0=[u0;v0];
WW=W0;
sol=zeros(length(t),N+1);
sol(1,:)=[0,W0(1:N-1)',0];
for j=1:length(t)-1
k1=M*WW;
WW=(eye(2*N-2)+dt*M/2)*WW+(dt/2)*k1+(dt^2)*M*k1/2;
sol(j+1,:)=[0,WW(1:N-1)',0];
end
[xx,tt]=meshgrid(x,t);
surf(xx,tt,sol);
1.3 Aproximación por Fourier con diferentes términos de series.
%utt-uxx=0
%u(0,t)=u(l,t)=0
%u(x,0)=2-abs(0.4x-2)
%ut(x,0)=0
%1 iteracion
clear all
L=10;
T=40;
Q=1;
h=0.1;
N=L/h;
x=0:h:L;
dt=h;
t=0:dt:T;
sol=zeros(length(t),N+1);
for k=1:Q
phi=sin(k/2*pi*x);
a=trapz(x,(2-abs(2/5*x-2).*phi))/(trapz(x,phi.^2));
b=0;
T=a.*cos(k*pi*t/L)/(k*pi/L);
sol=sol+T'*phi;
end
[xx,tt]=meshgrid(x,t);
figure (1)
surf(xx,tt,sol)
%3 iteracion
clear all
L=10;
T=40;
Q=3;
h=0.1;
N=L/h;
x=0:h:L;
dt=h;
t=0:dt:T;
sol=zeros(length(t),N+1);
for k=1:Q
phi=sin(k/2*pi*x);
a=trapz(x,(2-abs(2/5*x-2).*phi))/(trapz(x,phi.^2));
b=0;
T=a.*cos(k*pi*t/L)/(k*pi/L);
sol=sol+T'*phi;
end
[xx,tt]=meshgrid(x,t);
figure (2)
surf(xx,tt,sol)
%5 iteracion
clear all
L=10;
T=40;
Q=5;
h=0.1;
N=L/h;
x=0:h:L;
dt=h;
t=0:dt:T;
sol=zeros(length(t),N+1);
for k=1:Q
phi=sin(k/2*pi*x);
a=trapz(x,(2-abs(2/5*x-2).*phi))/(trapz(x,phi.^2));
b=0;
T=a.*cos(k*pi*t/L)/(k*pi/L);
sol=sol+T'*phi;
end
[xx,tt]=meshgrid(x,t);
figure (3)
surf(xx,tt,sol)
%10 iteracion
clear all
L=10;
T=40;
Q=10;
h=0.1;
N=L/h;
x=0:h:L;
dt=h;
t=0:dt:T;
sol=zeros(length(t),N+1);
for k=1:Q
phi=sin(k/2*pi*x);
a=trapz(x,(2-abs(2/5*x-2).*phi))/(trapz(x,phi.^2));
b=0;
T=a.*cos(k*pi*t/L)/(k*pi/L);
sol=sol+T'*phi;
end
[xx,tt]=meshgrid(x,t);
figure (4)
surf(xx,tt,sol)
%20 iteracion
clear all
L=10;
T=40;
Q=20;
h=0.1;
N=L/h;
x=0:h:L;
dt=h;
t=0:dt:T;
sol=zeros(length(t),N+1);
for k=1:Q
phi=sin(k/2*pi*x);
a=trapz(x,(2-abs(2/5*x-2).*phi))/(trapz(x,phi.^2));
b=0;
T=a.*cos(k*pi*t/L)/(k*pi/L);
sol=sol+T'*phi;
end
[xx,tt]=meshgrid(x,t);
figure (5)
surf(xx,tt,sol)
2 Energía del Cable.
La energía del cable que viene definida por la función: \begin{equation} E(t) = \int_{0}^{L} (u_{t}^{2}(x,t)+ u_{x}^{2}(x,t)) \cdot dx \end{equation}
FUNCIÓN:
function v = deriva(x,y,t)
n = length(x);
del=zeros(1,n);
if t==1
i=1;
del(1)=(y(i)-y(i+1))/(-x(i)+x(i+1));
for i=2:n
del(i)=(y(i)-y(i-1))/(x(i)-x(i-1));
end
elseif t==-1
for i=1:n-1
del(i)=(y(i)-y(i+1))/(x(i)-x(i+1));
end
i=n;
del(i)=(y(i)-y(i-1))/(x(i)-x(i-1));
elseif t==0
i=1;
del(1)=(y(i)-y(i+1))/(+x(i)-x(i+1));
for i=2:n-1
del(i)=(y(i+1)-y(i-1))/(2*(x(i)-x(i-1)));
end
i=n;
del(n)=(y(i)-y(i-1))/(x(i)-x(i-1));
end
v=del;
end
PROGRAMA:
clc;clear all;
L=10;dx=0.1;N=L/dx;
x=dx:dx:L-dx;
xtot=0:dx:L;
T=40;dt=0.1;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/dx^2)*K;
U=(2-abs((x*2/5-2)))';
V=(0*x)';
sol(1,:)=[0 U' 0];
for i=1:length(t)-1;
U=U+dt*V;
V=V-dt*K*U;
sol(i+1,:)=[0 U' 0];
end
%DERIVAMOS LA SOLUCION EN x Y EN t:
dsol_dx=ones(size(sol));
for i=1:length(t)
dsol_dx(i,:)=deriva(sol(i,:),xtot,0);
end
dsol_dt=ones(size(sol));
for i=1:length(xtot)
dsol_dt(:,i)=deriva(sol(:,i),t,1);
end
%SUMAMOS LAS DERIVADAS AL CUADRADO E INTEGRAMOS:
dE=zeros(1,size(t,1));
for i=1:length(xtot)
dE=dE+(dsol_dt(i,:).^2+dsol_dx(i,:).^2)*dx;
end
xtot=0:dx:L;
[xx,tt]=meshgrid(xtot,t);
figure(1)
surf(xx,tt,dsol_dx);
figure(2)
surf(xx,tt,dsol_dt);
