Propagación de las ondas en un medio dispersivo

Movimiento ondulatorio armónico

En la página web titulada "Descripción de la propagación" vemos en una animación la propagación de un pulso, cuya forma inicial viene descrito por la función f(x), a lo largo del eje X hacia la derecha con velocidad v y sin distorsionarse (sin cambiar de forma).

Ψ =f(x-vt) describe la propagación de una perturbación representada por la función f(x), sin distorsión, a la largo del eje X, hacia la derecha, con velocidad v.

2 Ψ t 2 = v 2 2 Ψ x 2

Estudiamos un caso particular importante, aquél en el que la función f(x) es una función armónica (seno o coseno).

Ψ(x,t)=Ψ0·sin k(x-vt)

Las características de esta función de dos variables, son las siguientes:

  1. La función seno es periódica y se repite cuando el argumento se incrementa en 2π . La función Ψ(x, t) se repite cuando x se incrementa en 2π/k.

  2. k( x+ 2π k vt )=kxkvt+2π

    Se trata de una función periódica, de periodo espacial o longitud de onda λ =/k. La magnitud k se denomina número de onda.

  3. Cuando se propaga un movimiento ondulatorio armónico, un punto x del medio describe un Movimiento Armónico Simple de amplitud Ψ0 y frecuencia angular ω =kv.

  4. Ψ(x,t)=Ψ0·sin (kx-ω t)

    El periodo de la oscilación es P= y la frecuencia  f =1/P.

  5. La igualdad ω =kv, nos permite relacionar el periodo espacial o longitud de onda λ y el periodo de la oscilación P de un punto del medio.

  6. ω=kv 2π P = 2π λ vλ=vP 

    La longitud de onda λ está relacionada con la frecuencia f de la forma λ=v/f. Para una velocidad de propagación v, cuanto mayor es la longitud de onda menor es la frecuencia y viceversa.

x=-6:0.05:6;
k=2*pi/2;
w=2*pi/1;

for i=1:5
    subplot(5,1,i)
    t=(i-1)/4;
    y=sin(k*x-w*t);
    plot(x,y,'r');
    grid on
    ylim([-1.1,1.1])
    str=sprintf('t=%1.2f',t);
    title(str)
end

En la figura vemos que la longitud de onda es λ=2. Un punto del medio, por ejemplo x=0, describe un MAS cuyo periodo es P=1. Un pico (señalado con una flecha) se desplaza una longitud de onda durante un periodo de oscilación.

Superposición de ondas armónicas

Superponemos dos ondas armónicas de la misma amplitud A, la primera onda se propaga con velocidad v1=ω1/k1 y la segunda con velocidad v2=ω2/k2

Ψ 1 (x,t)=Asin( k 1 x ω 1 t) Ψ 2 (x,t)=Asin( k 2 x ω 2 t) Ψ 1 (x,t)+ Ψ 2 (x,t)=2Acos( k 1 k 2 2 x ω 1 ω 2 2 t )sin( k 1 + k 2 2 x ω 1 + ω 2 2 t )

Supongamos que k1 y k2 son números de onda cercanos, por ejemplo k1=12 y k2=10. Una fotografía de la onda resultante en el instante t=0, sería la siguiente

k1=12;
k2=10;
x=0:0.01:6;
y=2*cos((k1-k2)*x/2).*sin((k1+k2)*x/2); %superposición de dos ondas armónicas
y1=2*cos((k1-k2)*x/2); %envolvente
plot(x,y,'b',x,y1,'r',x,-y1,'r');
ylim([-2.1,2.1])
xlabel('x');
ylabel('y')
title('Superposición de ondas')
grid on

La curva en color rojo es la envolvente que varía como 2·cos(x) y la curva en color azul la superposición de las dos ondas en el instante t=0, que varía como 2·cos(x)·sin(11x). Dado que k1 k2 y ω1 ω2 la velocidad de la superposición (en color azul) de las dos ondas es v=ω1/k1 ω2/k2 y se denomina velocidad de fase. Sin embargo, la velocidad de la envolvente es Δωk=(ω1- ω2)/(k1-k2 ), se denomina velocidad de grupo.

Las definiciones de velocidad de fase vp y de grupo vg son, respectivamente

v p = ω k v g = dω dk

Ambas velocidades son funciones de k y son distintas, en general.

fichero = 'onda_1.gif';
hg=figure;
set(hg,'Position',[0,0,568,180]) %posición y tamaño de la ventana gráfica
k1=12; w1=21;
k2=10; w2=20;
x=0:0.02:20;
for t=0:0.05:2*pi
    %superposición de dos ondas armónicas
    y=2*cos((k1-k2)*x/2-(w1-w2)*t/2).*sin((k1+k2)*x/2-(w1+w2)*t/2);
    y1=2*cos((k1-k2)*x/2-(w1-w2)*t/2); %envolvente
    plot(x,y,'b',x,y1,'r',x,-y1,'r');
    ylim([-2.1,2.1])
    xlabel('x');
    ylabel('y')
    title('Superposición de dos ondas')
   
    %GIF animado
    frame=getframe;
    im = frame2im(frame);
    [imind,cm] = rgb2ind(im,256);
    if t==0
        imwrite(imind,cm,fichero,'gif','DelayTime',0,'loopcount',inf);
    else
        imwrite(imind,cm,fichero,'gif','DelayTime',0,'writemode','append');
    end
end

En esta animación, comparamos la velocidad de la envolvente (ω1- ω2)/(k1-k2 )=1/2 (en color rojo) con la velocidad de fase de las ondas (ω1+ω2)/(k1+k2 )=(21+20)/(12+10)≈2 (en color azul)

Ejemplos de medios dispersivos y no dispersivos

Velocidad de grupo

La ecuación de una onda armónica es

Ψ(x,t)= Ψ 0 sin( kxωt) )

donde k es el número de onda, ω es la frecuencia angular, y v=ω/k es la velocidad de propagación

La ecuación de una onda, es la superposición

Ψ(x,t)= a(k)exp( i(kxωt) ) dk

En un medio dispersivo la frecuencia angular ω depende del número de onda k, las ondas de diferentes longitudes de onda λ=2π/k se propagan con distinta velocidad.

Al final de la página titulada Transformada de Fourier, hemos calculado la transformada de Fourier de la función.

Ψ ( x ,0 ) = exp ( i k 0 x ) · exp ( ( x x 0 ) 2 2 σ x 2 ) a ( k ) = A exp ( i k x 0 ( k k 0 ) 2 2 σ k 2 ) σ x = 1 σ k

Como hemos visto a(k) es una función que tiene un pico agudo centrado en k0, por tanto, a(k) es distinto de cero en las proximidades de k0 y es nula en el resto.

La fase φ(k)=kx-ωt es una función de k. Aproximamos φ(k) alrededor de k0 (donde a(k) es distinto de cero) tomando los primeros términos del desarollo en serie de esta función.

φ(k)= φ 0 + dφ dk | 0 ( k k 0 )+ 1 2 d 2 φ d k 2 | 0 ( k k 0 ) 2 +...

donde

φ 0 = k 0 x ω 0 t dφ dk | 0 =x dω dk | 0 t=x v g t d 2 φ d k 2 | 0 = d 2 ω d k 2 | 0 t=αt φ(k) k 0 x ω 0 t+(x v g t)( k k 0 ) 1 2 αt ( k k 0 ) 2

vg de denomina velocidad de grupo de ondas, mientas que el cociente ω/k se denomina velocidad de fase.

Calculamos la ecuación de la onda en cualquier instante t

Ψ(x,t)=A exp( ik x 0 ( k k 0 ) 2 2 σ k 2 +i( k 0 x ω 0 t+(x v g t)(k k 0 ) 1 2 αt ( k k 0 ) 2 ) ) dk Ψ(x,t)=Aexp( i( k 0 x ω 0 t k 0 x 0 ) ) exp( i(k k 0 ) x 0 ( k k 0 ) 2 2 σ k 2 +i( (x v g t)(k k 0 ) 1 2 αt ( k k 0 ) 2 ) ) dk Ψ(x,t)=Aexp( i( k 0 x ω 0 t k 0 x 0 ) ) exp( i(x x 0 v g t)(k k 0 )(iαt σ k 2 +1) ( k k 0 ) 2 2 σ k 2 ) dk

Hacemos el cambio de variable y tenemos una integral del tipo descrito al final de la página titulada Transformada de Fourier o al principio del ejemplo 1 en esta página con

a= 1 2 σ ,b=i yσ 2 , c 2 = σ 2 y 2 2 Ψ(x,t)=Aexp( i( k 0 x ω 0 t k 0 x 0 ) ) π exp( c 2 ) a

El resultado es

u=k k 0 ,du=dk y=x x 0 v g t, 1 σ 2 = iαt σ k 2 +1 σ k 2 =iαt+ 1 σ k 2 Ψ(x,t)=Aexp( i( k 0 x ω 0 t k 0 x 0 ) ) exp( iyu u 2 2 σ 2 ) du Ψ(x,t)= 2π σ k 2 iαt σ k 2 +1 Aexp( i( k 0 x ω 0 t k 0 x 0 ) )exp( ( x x 0 v g t ) 2 σ k 2 2(iαt σ k 2 +1) ) Ψ(x,t)= 2π σ k 2 iαt σ k 2 +1 Aexp( i( k 0 x ω 0 t k 0 x 0 ) )exp( ( x x 0 v g t ) 2 (1iα σ k 2 t) 2( α 2 σ k 2 t 2 + 1 σ k 2 ) )

Propagación de un pulso sin distorsión, α=0, velocidad de grupo vg=2. El punto de color azul marca el centro del grupo que se mueve con velocidad constante vg.

x0=0; %posición inicial de la función de Gauss
k0=2; %número de onda
sx2=5; %extensión de la función de Gauss (cuadrado de sigma)
w0=4*k0; %velocidad de fase w0/k0
vg=2; %velocidad de grupo
alfa=0; %dispersión
sk2=1/sx2;
x=-8:0.02:30;

for i=1:5
    subplot(5,1,i)
    t=(i-1)*3;
    y=sqrt(2*pi*sk2/(1i*alfa*t*sk2+1))*exp(1i*(k0*x-w0*t-k0*x0))...
.*exp(-(x-x0-vg*t).^2*(1-1i*alfa*sk2*t)/(2*(alfa^2*sk2*t^2+1/sk2)));
    hold on
    plot(x,real(y),'r'); %movimiento ondulatorio
    plot(vg*t,0,'bo','markersize',2,'markerfacecolor','b') %centro del grupo
    hold off
    grid on
    ylim([-1.1,1.1])
    str=sprintf('t=%1.2f',t);
    title(str)
end

Propagación del pulso con distorsión α=1, velocidad de grupo vg=2. El punto de color azul marca el centro del grupo.

Animamos la segunda figura

fichero = 'onda_2.gif';
hg=figure;
%posición y tamaño de la ventana gráfica
set(hg,'Position',[0,200,690,180]) 

x0=0; %posición inicial de la función de Gauss
k0=2; %número de onda
sx2=5; %extensión de la función de Gauss (cuadrado de sigma)
w0=4*k0; %velocidad de fase w0/k0
vg=2; %velocidad de grupo
alfa=1; %dispersión
sk2=1/sx2;
x=-8:0.02:30;
for t=0:0.1:10;
    y=sqrt(2*pi*sk2/(1i*alfa*t*sk2+1))*exp(1i*(k0*x-w0*t-k0*x0))...
.*exp(-(x-x0-vg*t).^2*(1-1i*alfa*sk2*t)/(2*(alfa^2*sk2*t^2+1/sk2)));
    plot(x,real(y),'r')
    ylim([-1.1,1.1])
    grid on
    
    %GIF animado
    frame=getframe;
    im = frame2im(frame);
    [imind,cm] = rgb2ind(im,256);
    if t==0
        imwrite(imind,cm,fichero,'gif','DelayTime',0,'loopcount',inf);
    else
        imwrite(imind,cm,fichero,'gif','DelayTime',0,'writemode','append');
    end
end

Ejemplo 1

En primer lugar, obtenemos el resultado de la integral

I= exp( ( axb ) 2 c 2 )dx =exp( c 2 ) exp( ( axb ) 2 )dx u=axb,du=a·dx I= exp( c 2 ) a exp( u 2 )du = π exp( c 2 ) a

El resultado de la última integral se justifica en la página titulada Integrales

Medio no dispersivo

En un medio no dispersivo, la relación entre la frecuencia angular ω y el número de onda k es ω=| k |v

Sea a(k) una función de Gauss centrada en k0. El parámetro Δk controla su anchura tal como se aprecia en la figura

a(k)=exp( ( k k 0 ) 2 2 ( Δk ) 2 )

k0=2;
hold on
for Dk=[3,5]
    g=@(k) exp(-(k-k0).^2/(2*Dk^2));
    fplot(g,[-10,14])
end
line([k0,k0],[0,1],'lineStyle','--')
hold off
grid on
xlabel('k')
legend('3','5','location','best')
ylabel('a(k)')
title('Paquete de ondas')

Tenemos que calcular la integral

Ψ(x,t)= 1 2π Δk dkexp( ( k k 0 ) 2 2 ( Δk ) 2 ) exp( i( kxωt ) ) Ψ(x,t)= 1 2π Δk dkexp( ( k k 0 ) 2 2 ( Δk ) 2 ) exp( i( kxv| k |t ) )

Esta integral es la suma de dos

Ψ(x,t)= 1 2π Δk 0 dkexp( ( k k 0 ) 2 2 ( Δk ) 2 ) exp( i( kx+vkt ) )+ 1 π Δk 0 dkexp( ( k k 0 ) 2 2 ( Δk ) 2 ) exp( i( kxvkt ) )

Resolvemos la primera integral

( k k 0 ) 2 2 ( Δk ) 2 +ik( x+vt )= ( akb ) 2 c 2 1 2 ( Δk ) 2 k 2 +( k 0 ( Δk ) 2 +i( x+vt ) )k k 0 2 2 ( Δk ) 2 = a 2 k 2 +2abk b 2 c 2 , { a= 1 Δk 2 2ab=( k 0 ( Δk ) 2 +i( x+vt ) ),b= 1 2 ( k 0 Δk +i( x+vt )Δk ) b 2 c 2 = k 0 2 2 ( Δk ) 2 1 2 ( k 0 Δk +i( x+vt )Δk ) 2 + c 2 = k 0 2 2 ( Δk ) 2 c 2 = 1 2 ( x+vt ) 2 ( Δk ) 2 i k 0 ( x+vt ) 1 2π Δk Δk 2 π 2 exp( 1 2 ( x+vt ) 2 ( Δk ) 2 +i k 0 ( x+vt ) )

De modo análogo, se calcula la segunda integral. El resultado es

Ψ(x,t)= 1 2 exp( 1 2 ( x+vt ) 2 ( Δk ) 2 )exp( i k 0 ( x+vt ) )+ 1 2 exp( 1 2 ( xvt ) 2 ( Δk ) 2 )exp( i k 0 ( xvt ) )

Se ha tenido en cuenta el resultado

0 exp( u 2 )du= 0 exp( u 2 )du = π 2

Ψ(x,t) representa una onda que se propaga sin distorsión a lo largo del eje X, hacia la derecha y hacia la izquierda, centradas en x=vt y x=-vt, modulada en amplitud

v=1; %velocidad de propagación
k0=20;
Dk=3;
hold on
f=@(x,t) (exp(-Dk^2*(x+v*t).^2/2).*exp(1i*k0*(x+v*t))+
exp(-Dk^2*(x-v*t).^2/2).*exp(1i*k0*(x-v*t)))/2;
for t=[0,4,8,12]
    g=@(x) real(f(x,t));
    fplot(g,[-14,14])
end
hold off
grid on
xlabel('x')
legend('t=0','t=4','t=8','t=12','location','best')
ylabel('\Psi(x,t)')
title('Paquete de ondas')

Otra forma de ver la figura anterior

v=1; %velocidad de propagación
k0=20;
Dk=3;
f=@(x,t) (exp(-Dk^2*(x+v*t).^2/2).*exp(1i*k0*(x+v*t))+
exp(-Dk^2*(x-v*t).^2/2).*exp(1i*k0*(x-v*t)))/2;
for i=0:3
    subplot(4,1,i+1)
    t=i*4;
    g=@(x) real(f(x,t));
    fplot(g,[-14,14])
    grid on
    xlabel('x')
    ylabel('\Psi(x,t)')
    str=sprintf('t=%1.1f',t);
    title(str)
end

Medio dispersivo

En este medio dispersivo, la relación entre la frecuencia angular ω y el número de onda k es no lineal, por ejemplo, ω=γ k 2

Tenemos que calcular la integral

Ψ(x,t)= 1 2π Δk dkexp( ( k k 0 ) 2 2 ( Δk ) 2 ) exp( i( kxωt ) ) Ψ(x,t)= 1 2π Δk dkexp( ( k k 0 ) 2 2 ( Δk ) 2 ) exp( i( kxγ k 2 t ) )

Utilizando el resultado del ejemplo 1

( k k 0 ) 2 2 ( Δk ) 2 +i( kxγ k 2 t )= ( akb ) 2 c 2 ( 1 2 ( Δk ) 2 +iγt ) k 2 +( k 0 ( Δk ) 2 +ix )k k 0 2 2 ( Δk ) 2 = a 2 k 2 +2abk b 2 c 2 , { a 2 = 1 2 ( Δk ) 2 +iγt 2ab= k 0 ( Δk ) 2 +ix,b= 1 2 Δk k 0 +ix ( Δk ) 2 1+2iγ ( Δk ) 2 t b 2 c 2 = k 0 2 2 ( Δk ) 2 1 2 ( Δk ) 2 ( k 0 +ix ( Δk ) 2 ) 2 1+2iγ ( Δk ) 2 t + c 2 = k 0 2 2 ( Δk ) 2 c 2 = k 0 2 ( k 0 +ix ( Δk ) 2 ) 2 1+2iγ ( Δk ) 2 t 2 ( Δk ) 2 = 2iγ k 0 2 ( Δk ) 2 t+ x 2 ( Δk ) 4 2i k 0 x ( Δk ) 2 2 ( Δk ) 2 ( 1+2iγ ( Δk ) 2 t )

El resultado es

Ψ(x,t)= 1 2π Δk 1 1 2 ( Δk ) 2 +iγt π exp( 2iγ k 0 2 t+ x 2 ( Δk ) 2 2i k 0 x 2( 1+2iγ ( Δk ) 2 t ) )= 1 1+2iγ ( Δk ) 2 t exp( 1 2 x 2 ( Δk ) 2 i k 0 ( xγ k 0 t ) 1+2iγ ( Δk ) 2 t ) Ψ(x,t)= 1 1+2iγ ( Δk ) 2 t exp( 1 2 ( Δk ) 2 ( x2γ k 0 t ) 2 i( k 0 ( xγ k 0 t )+ x 2 ( Δk ) 4 γt ) 1+4 γ 2 ( Δk ) 4 t 2 )

La amplitud es proporcional a

exp( 1 2 ( Δk ) 2 ( x2γ k 0 t ) 2 )

El valor máximo de la amplitud se obtiene para x=2γk0t. Marcamos con un punto de color rojo en las figuras los centros xp de los paquetes de ondas, que se mueven con la velocidad de grupo

v g = dω dk =2γk x p =2γ k 0 t

k0=10;
Dk=3;
gamma=0.5;
f=@(x,t) exp(-(x.^2*Dk^2/2-1i*k0*(x-gamma*k0*t))/(1+2*1i*gamma*Dk^2*t))
/sqrt(1+2*1i*gamma*Dk^2*t);
for i=0:3
    subplot(4,1,i+1)
    t=i*0.4;
    g=@(x) real(f(x,t));
    hold on
    fplot(g,[-2,22])
     %centro del grupo
    plot(2*gamma*k0*t,0,'ro','markersize',3,'markerfacecolor','r')
    hold off
    grid on
    xlabel('x')
    ylabel('\Psi(x,t)')
    str=sprintf('t=%1.2f',t);
    title(str)
end

Ejemplo 2

La ecuación diferencial que describe la vibración de una barra elástica

λ 2 4 y x 4 + 2 y t 2 =0

Buscamos una solución de la forma

y(x,t)=Aexp( i( kxωt ) ) ω=λ k 2

Obtenemos una relación entre la frecuencia angular ω y el número de onda k que no es lineal (como en un medio no dispersivo)

La solución de la ecuación diferencial es la superposición

y(x,t)= A(k)exp( i( kxω(k)t ) )

Consideremos el caso, en el que la amplitud A(k) es no nula en un pequeño intervalo

A(k)={ A 0 , k 0 Δk<k< k 0 +Δk 0,el resto

El desplazamiento de los puntos de la barra es

y(x,t)= A 0 k 0 Δk k 0 +Δk cos( kxλ k 2 t )dk

Esta integral no es sencilla. El resultado que proporciona Math Symbolic de MATLAB es el siguiente. Se ha sustituido el parámetro λ por la variable a.

>> syms x k a t;
>> z=int(cos(k*x-a*k^2*t),k)
z =pi^(1/2)*cos(x^2/(4*a*t))*fresnelc((2^(1/2)*(x/2 - 
a*k*t)*(-1/(a*t))^(1/2))/pi^(1/2))*(-1/(2*a*t))^(1/2) - 
pi^(1/2)*sin(x^2/(4*a*t))*fresnels((2^(1/2)*(x/2 - a*k*t)*
(-1/(a*t))^(1/2))/pi^(1/2))*(-1/(2*a*t))^(1/2)

Mediante el comando

>> latex(z)

La respuesta ans se transforma en código Latex que copiamos y pegamos en el programa MathType, que a su vez, lo transforma en código MathML

π cos( x 2 4λt )C( 2 ( x 2 λkt ) 1 λt π ) 1 2λt π sin( x 2 4λt )S( 2 ( x 2 λkt ) 1 λt π ) 1 2λt

Teniendo en cuenta que

z =i z C(ix)=iC(x) S(ix)=iS(x)

donde

S(z)= 0 z sin( π t 2 2 )dt ,C(z)= 0 z cos( π t 2 2 )dt

son las integrales de Fresnel

El desplazamiento de los puntos de la barra elástica es

y(x,t)= A 0 π 2λt { cos( x 2 4λt )C( x2λkt 2πλt )+sin( x 2 4λt )S( x2λkt 2πλt ) } k 0 Δk k 0 +Δk = A 0 π 2λt { cos( x 2 4λt ){ C( x2λ( k 0 +Δk )t 2πλt )C( x2λ( k 0 Δk )t 2πλt ) }+sin( x 2 4λt ){ S( x2λ( k 0 +Δk )t 2πλt )S( x2λ( k 0 Δk )t 2πλt ) } }

Representamos esta función en los instantes t=0.02, 5.02 y 10.02, para los valores de los parámetros

syms t x lambda k0 dk;
y=-sqrt(pi/(2*lambda*t))*(cos(x.^2/(4*lambda*t)).*
(fresnelc((x-2*lambda*t*(k0+dk))/sqrt(2*pi*lambda*t))-
fresnelc((x-2*lambda*t*(k0-dk))/sqrt(2*pi*lambda*t)))+
sin(x.^2/(4*lambda*t)).*(fresnels((x-2*lambda*t*(k0+dk))/
sqrt(2*pi*lambda*t))-fresnels((x-2*lambda*t*(k0-dk))/sqrt(2*pi*lambda*t))));

subplot(3,1,1)
y1=subs(y,{t,lambda,dk, k0},{0.02,10,0.05,1});
fplot(y1,[0,400])
grid on
xlabel('x')
ylabel('y(x)')
title('t=0.02')

subplot(3,1,2)
y1=subs(y,{t,lambda,dk, k0},{5.02,10,0.05,1});
fplot(y1,[0,400])
grid on
xlabel('x')
ylabel('y(x)')
title('t=5.02')

subplot(3,1,3)
y1=subs(y,{t,lambda,dk, k0},{10.02,10,0.05,1});
fplot(y1,[0,400])
grid on
xlabel('x')
ylabel('y(x)')
title('t=10.02')

Referencias

Masatsugu Suzuki. Lecture Notes of General Physics I and II Calculus Based. Topic 5. Waves and oscillations

David Romero-Abad, Roberto Suárez-Córdova. Energy transported by mechanical waves - two interesting examples: waves on a string and waves on a bar. No publicado