Ecuación de ondas (Raúl, Sofía, Jaime)

Trabajo realizado por estudiantes
Título	Ecuación de ondas
Asignatura	EDP
Curso	2023-24
Autores	Raúl Ortega Sofía Gómez Jaime Sáenz de Miera
Este artículo ha sido escrito por estudiantes como parte de su evaluación en la asignatura

En este artículo vamos a estudiar la ecuación de ondas mediante el estudio de del movimiento de una cuerda. Para ello, vamos a ver cómo varía la solución en función de las condiciones iniciales...

1 Ecuación de ondas

Vamos a considerar una cuerda vibrante en el intervalo [0,1], con densidad d y tensión [math] \tau _0[/math] constante, tal que la velocidad de propagación es: [math] c=\frac{\tau_0}{d}=1[/math]. Además, vamos a suponer que la cuerda está fija en los extremos con [math] u_0(x)[/math] y [math] u_1(x)[/math] su posición e impulso iniciales respectivamente. Teniendo todo esto en cuenta y la ecuación de ondas, obtenemos el siguiente sistema:

[math] \begin{cases} u_{tt}-u_{xx}=0,\\ u(0,t)=u(1,t)=0,\\ u(x,0)=u_0(x),\\ u_t(x,0)=u_1(x). \end{cases} [/math]

Donde [math]u(x)[/math] representa la posición de cada punto de la cuerda en cada momento.

1.1 Resolución

Para resolver este sistema, vamos a utilizar el método de separación de variables, el cual consiste en suponer que [math] u(x,t)=T(t)X(x)[/math], lo cual nos da lugar a resolver dos problemas por separado, uno que depende de t y otro que depende de x. Al hacer esto, obtenemos la siguiente solución:

[math] u(x,t)=\sum_{k=1}^\infty a_k sin(k\pi x) cos(k \pi t) + b_k sin(k \pi x) sin(k \pi t), [/math]

siendo [math] a_k=\frac{\int_0^1 sin(k\pi x)u_0(x)dx}{\int_0^1 sin^2(k \pi x)dx}[/math] y [math] b_k=\frac{1}{k\pi}\frac{\int_0^1 sin(k\pi x)u_1(x)dx}{\int_0^1 sin^2(k \pi x)dx}[/math]

Ahora, para entender mejor esta solución, vamos a tomar como datos iniciales [math]u_0(x)=e^{-100(x-\frac{1}{2})^2}[/math] y [math] u_1(x)=0[/math], lo cual nos da que la solución se expresa de la siguiente manera:

[math] u(x,t)=\sum_{k=1}^\infty a_k sin(k\pi x) cos(k \pi t) .[/math]

Para representarla vamos a calcular las integrales necesarias para obtener los [math] a_k [/math] mediante el método del trapecio y a dibujar los primeros 50 términos de la serie.

Video de la evolución de la solución a lo largo del tiempo

%Definimos las funciones
u0=@(x) exp(-100*(x-0.5).^2); %Condición inicial
uk=@(x,t,k) sin(k.*pi.*x).*cos(k*pi.*t); %Solución del problema

%Aplicamos la regla del trapecio
a=0; b=1; %Extremos del intervalo                 
x=linspace(a,b,100); %Vector de espacio
t=linspace(a,2,100); %Vector de tiempo
N=length(x)-1; %Número de puntos
h=(b-a)/N;
w=ones(N+1,1); %Vector de pesos             
w(1)=1/2; w(N+1)=1/2;
f_n=zeros(length(x),length(t));

%Suma de los 10 primeros términos en todo el vector de tiempo
for i=1:length(t)
    for k=1:50
        g=u0(x)*sin(k*pi*x)'; 
        a_k=2*h*w'*g; %Resultado de la regla del trapecio (c_k)
        f_n(i,:)= f_n(i,:)+a_k.*uk(x,t(i),k);
    end
end

%Hacemos un video representando la solución para distintos tiempos
pelicula=VideoWriter("cuerda",'MPEG-4');
pelicula.FrameRate=10;
open(pelicula)
for i=1:length(t)
    f=f_n(i,:);
    figura=figure(1);
    plot(x,f,color='b',LineWidth=1.5)
    title(['t = ' num2str(t(i))]);
    axis([0 1 -1 1])
    imagen=getframe(figura);
    writeVideo(pelicula,imagen);
end
close(pelicula)

1.2 Onda que viaja en un solo sentido

En este caso tomamos como condiciones iniciales: [math]u_0(x)=e^{-100(x-\frac{1}{2})^2}[/math] y [math] u_1(x)=200(x-\frac{1}{2})e^{-100(x-\frac{1}{2})^2}[/math] que es la derivada de [math]u_0(x)[/math] cambiada de signo. Al igual que en el caso anterior para representarlo vamos a realizar las integrales correspondientes para calcular los coeficientes de la solución mediante el método del trapecio y a representar los primeros 50 términos.

Video de la evolución de la solución a lo largo del tiempo

%Definimos las funciones
u0=@(x) exp(-100*(x-0.5).^2); %Condición inicial 
u1=@(x) -(100-200.*x).*exp(-100*(x-0.5).^2); %Impulso inicial
uk_cos=@(x,t,k) sin(k.*pi.*x).*cos(k*pi.*t); %Primera parte de la solución del problema
uk_sen=@(x,t,k) sin(k.*pi.*x).*sin(k*pi.*t); %Segunda parte de la solución del problema

%Aplicamos la regla del trapecio
a=0; b=1; %Extremos del intervalo                 
x=linspace(a,b,100); %Vector de espacio
t=linspace(a,2,100); %Vector de tiempo
N=length(x)-1; %Número de puntos
h=(b-a)/N;
w=ones(N+1,1); %Vector de pesos             
w(1)=1/2; w(N+1)=1/2;
f_n=zeros(length(x),length(t));

%Suma de los 10 primeros términos en todo el vector de tiempo
for i=1:length(t)
    for k=1:50
        g1=u0(x)*sin(k*pi*x)'; 
        g2=-u1(x)*sin(k*pi*x)'/(k*pi);
        a_k=2*h*w'*g1; %Resultado de la regla del trapecio (c_k)
        b_k=2*h*w'*g2;
        f_n(i,:)= f_n(i,:)+a_k.*uk_cos(x,t(i),k)+b_k.*uk_sen(x,t(i),k);
    end
end

%Hacemos un video representando la solución para distintos tiempos
pelicula=VideoWriter("cuerda",'MPEG-4');
pelicula.FrameRate=10;
open(pelicula)
for i=1:length(t)
    f=f_n(i,:);
    figura=figure(1);
    plot(x,f)
    axis([0 1 -2 2])
    imagen=getframe(figura);
    writeVideo(pelicula,imagen);
end
close(pelicula)

En este caso vemos que efectivamente la cuerda viaja en un solo sentido y sigue manteniendo los extremos fijos en cero. Además, como el intervalo en el que viaja la cuerda tiene longitud uno y lo recorre dos veces en el intervalo de tiempo [math][0,2][/math], podemos deducir que la velocidad a la que viaja la curva es 1. Esto coincide con la constante math>c</math> que hemos elegido la cual representa esta velocidad