Dispersión por un potencial V(r)

Dispersión por un potencial repulsivo V(r)=k/r2

La energía E y el momento angular L en coordenadas polares de una partícula de masa m y carga q son

E= 1 2 m ( dr dt ) 2 + 1 2 m r 2 ( dθ dt ) 2 +qV(r) L=m r 2 dθ dt

Estas dos magnitudes son constantes en todos los puntos de la trayectoria e iguales a su valor inicial

L=m v 0 b E= 1 2 m v 0 2

b se denomina parámetro de impacto y v0 es la velocidad de la partícula en el infinito. La ecuación diferencial de la trayectoria es

dθ dr = L m r 2 2 m ( E L 2 2m r 2 qV(r) )

Sustituimos E, L y V(r)

dθ dr = m v 0 b m r 2 2 m ( 1 2 m v 0 2 m 2 v 0 2 b 2 2m r 2 q Q 4π ε 0 r 2 ) dθ= b r 2 ( 1 b 2 r 2 qQ 2π ε 0 m v 0 2 1 r 2 ) dr

Hacemos el cambio de variable u=1/r, du=-dr/r2

dθ= bdu 1 k 2 u 2 = b k du 1 k 2 u 2 , k 2 = b 2 + qQ 2π ε 0 m v 0 2

La ecuación de la trayectoria en coordenadas polares es

θ= b k arcsin( u 1 k )+C c k b θ=arcsin( k r ) r= k sin( c k b θ )

La constante c se determina sabiendo que la partícula parte del infinito, r→∞ con un parámetro de impacto b, la posición angular θ→0 (véase la figura más abajo). Por tanto, c=π. La ecuacion de la trayectoria en coordenadas polares es

r= k sin( π k b θ )

La otra asíntota se obtiene para el ángulo θ1 tal que π-1/b=0, es decir, cuando r→∞

El ángulo de dispersión es

Φ=π θ 1 =π b k π=( 1 b k )π

Para un parámetro de impacto dado b, la distancia mínima a entre la partícula y el centro fijo de fuerzas es

E= 1 2 m v 0 2 = L 2 2m a 2 + Qq 4π ε 0 a 2 1 2 m v 0 2 = m 2 v 0 2 b 2 2m a 2 + Qq 4π ε 0 a 2 a 2 = b 2 + Qq 2π ε 0 m v 0 2 a=k

En la posición de mínima distancia se cumple que la componente radial de la velocidad es nula, dr/dt=0

Cuando el parámero de impacto b=0, la partícula se mueve hacia el centro de fuerzas hasta una distancia h y luego, retrocede

1 2 m v 0 2 = Qq 4π ε 0 h 2 h 2 = Qq 2π ε 0 m v 0 2 a 2 = b 2 + h 2

La ecuación de la trayectoria en términos de los parámetros a y h se escribe

r= a sin( π a a 2 h 2 θ ) = a sin( π θ 1 h 2 a 2 )

Para representar la trayectoria se ha empleado el código MATLAB

h2=0.2; %se fija la energía de la partícula
hold on
b=0.3; % parámetro de impacto
k2=b^2+h2;
r=@(x) sqrt(k2)./sin(pi-sqrt(k2)*x/b);
Th=pi*b/sqrt(k2); %asíntota
fplot(@(x) r(x).*cos(x),@(x) r(x).*sin(x),[0.05,Th-0.05])
plot(0,0,'ro', 'markersize',8,'markerfacecolor','y')
plot(sqrt(k2)*cos(Th/2),sqrt(k2)*sin(Th/2),'bo', 'markersize',3
,'markerfacecolor','b')
plot(sqrt(h2),0,'ko', 'markersize',3,'markerfacecolor','k')
line([0,-0.5],[0,0.5*tan(pi-Th)],'lineStyle','--') %asíntota
hold off
axis equal
xlim([-0.5,2])
ylim([0,1.5])
xlabel('x')
ylabel('y')
title('Dispersión')

Representamos las trayectorias para una energía fijada por el parámetro h2, para varios parámeros de imapacto b

h2=0.2; %se fija la energía de la partícula
hold on
for b=0.1:0.1:1 % parámetro de impacto
    k2=b^2+h2;
    r=@(x) sqrt(k2)./sin(pi-sqrt(k2)*x/b);
    fplot(@(x) r(x).*cos(x),@(x) r(x).*sin(x),[0.05, pi*b/sqrt(k2)-0.05])
    %simétrica
    fplot(@(x) r(x).*cos(x),@(x) -r(x).*sin(x),[0.05, pi*b/sqrt(k2)-0.05])
end
plot(0,0,'ro', 'markersize',8,'markerfacecolor','y') %centro de fuerzas
hold off
grid on
axis equal
xlim([-2,2])
ylim([-1.5,1.5])
xlabel('x')
ylabel('y')
title('Dispersión')

Envolvente

La ecuación de la trayectoria depende del parámetro de impacto b, f(r, θ, b)=0

rsin( π b 2 + h 2 b θ ) b 2 + h 2 =0

La ecuación de la envolvente de las trayectorias se obtiene derivando con respecto a b e igualando a cero.

f b =0 r( b 2 b 2 + h 2 b 2 + h 2 b 2 )cos( π b 2 + h 2 b θ ) b b 2 + h 2 =0 r h 2 b 2 cos( π b 2 + h 2 b θ )b=0

No podemos despejar b de la ecuación de la trayectoria y sustituirlo en esta última. La ecuación de la envolvente viene descrita por el sistema de dos ecuaciones

{ rsin( π b 2 + h 2 b θ ) b 2 + h 2 =0 r h 2 cos( π b 2 + h 2 b θ ) b 3 =0

Dado el parámetro b, se resuelve el sistema de dos ecuaciones no lineales mediante fsolve para obtener r y θ que representa un punto, en coordenadas polares, de la envolvente (curva de color rojo)

h2=0.8; %se fija la energía de la partícula
hold on
for b=0.1:0.1:1 % parámetro de impacto
    k2=b^2+h2;
    r=@(x) sqrt(k2)./sin(pi-sqrt(k2)*x/b);
    fplot(@(x) r(x).*cos(x),@(x) r(x).*sin(x),[0.05, pi*b/sqrt(k2)-0.05])
    fplot(@(x) r(x).*cos(x),@(x) -r(x).*sin(x),[0.05, pi*b/sqrt(k2)-0.05])
end
plot(0,0,'ro', 'markersize',8,'markerfacecolor','y')
%envolvente
bb=0.1:0.05:1;
z=zeros(length(bb),2);
j=1;
for b=bb
%r es x(1) y th es x(2)
    k2=b^2+h2;
    F=@(x) [x(1)*sin(pi-sqrt(k2)*x(2)./b)-sqrt(k2); x(1)*h2*
cos(pi-sqrt(k2)*x(2)./b)-b.^3];
    zz=fsolve(F,[sqrt(k2),pi*b/(2*sqrt(k2))]);
        z(j,1)=zz(1);
        z(j,2)=zz(2);
        j=j+1;
end
plot(z(:,1).*cos(z(:,2)),z(:,1).*sin(z(:,2)),'r')
plot(z(:,1).*cos(z(:,2)),-z(:,1).*sin(z(:,2)),'r')
hold off
grid on
axis equal
xlim([-2,2])
ylim([-1.5,1.5])
xlabel('x')
ylabel('y')
title('Dispersión')

Para h2=0.2, el error en la envolvente es importante para valores del parámetro de impacto b grandes

Alternativamente, utilizamos el procedimiento de Newton-Raphson, para representar la evolvente, para ello necesitamos determinar la matriz J del siguiente modo

{ f( r,θ )=rsin( π b 2 + h 2 b θ ) b 2 + h 2 =0 g( r,θ )=r h 2 cos( π b 2 + h 2 b θ ) b 3 =0 J=( f r f θ g r g θ )=( sin( π b 2 + h 2 b θ ) b 2 + h 2 b rcos( π b 2 + h 2 b θ ) h 2 cos( π b 2 + h 2 b θ ) b 2 + h 2 b r h 2 sin( π b 2 + h 2 b θ ) )

function hiperbola_8
    h2=0.8; %se fija la energía de la partícula
    hold on
    for b=0.1:0.05:1 % parámetro de impacto
        k2=b^2+h2;
        r=@(x) sqrt(k2)./sin(pi-sqrt(k2)*x/b);
        fplot(@(x) r(x).*cos(x),@(x) r(x).*sin(x),[0.05, pi*b/sqrt(k2)-0.05])
        fplot(@(x) r(x).*cos(x),@(x) -r(x).*sin(x),[0.05, pi*b/sqrt(k2)-0.05])
    end
    plot(0,0,'ro', 'markersize',8,'markerfacecolor','y')
%envolvente
    bb=0.1:0.05:1;
     z=zeros(length(bb),2);
     j=1;
    for b=bb
        zz=envolvente(b);
        z(j,1)=zz(1);
        z(j,2)=zz(2);
        %disp([zz(1),zz(2)])
        j=j+1;
     end
    plot(z(:,1).*cos(z(:,2)),z(:,1).*sin(z(:,2)),'r')
    hold off
    grid on
    axis equal
    xlim([-2,2])
    ylim([-1.5,1.5])
    xlabel('x')
    ylabel('y')
    title('Dispersión')

    function x=envolvente(b)
        k2=b^2+h2;
        F=@(x) [x(1)*sin(pi-sqrt(k2)*x(2)./b)-sqrt(k2); 
x(1)*h2*cos(pi-sqrt(k2)*x(2)./b)-b.^3];
        J=@(x) [sin(pi-sqrt(k2)*x(2)./b), -x(1)*sqrt(k2)*cos(pi-sqrt(k2)*x(2)./b)/b
; h2*cos(pi-sqrt(k2)*x(2)./b), x(1)*h2*sin(pi-sqrt(k2)*x(2)./b)*sqrt(k2)/b];
        x=[sqrt(k2),pi*b/(2*sqrt(k2))]; %valor inicial
        for i=1:100
            dx=-J(x)\F(x);
            if sqrt(norm(dx)/norm(x))<0.001
%                   disp(x+dx)
                break;
            end
            x=x+dx;
        end
        if i==100
            disp('Se ha soprepasado el número de iteracciones');
        end
    end
end

Para los primeros valores del parámetro de impacto b, el procedimiento no realiza bien los cálculos tal como se puede ver en la siguiente tabla

La primera columna b es el parámetro de impacto. La segunda y tercera columna, es la posición (r,θ) de un punto de la envolvente en coordendas polares, calculado mediante fsolve de MATLAB. La cuarta y quinta columna es la posición de dicho punto calculada medinte el método de Newton-Raphson

brθrθ
0.10.90000.17420.90000.8728
0.20.91660.3452-0.91661.0307
0.30.94400.5109-0.9440-0.4881
0.40.98310.67450.98310.6745
0.51.03650.84031.03650.8403
0.61.11041.01191.11041.0119
0.71.21401.19061.21401.1906
0.81.36001.37381.36001.3738
0.91.56221.55591.56221.5559
1.01.83371.72991.83371.7299

A partir de b=0.4 ambos procedimientos coinciden en los resultados

Dispersión de un ión por un átomo neutro

Consideremos un átomo inmóvil, situado en el origen y un ión de msa m y carga q situado a una distancia r

El campo eléctrico producido por el ión polariza el átomo que podemos considerar ahora como un dipolo de cargas +Q y -Q separadas una distancia 2a.

El campo eléctrico producido por el dipolo en la posición del ión es

E p = 1 4π ε 0 Q (r+a) 2 1 4π ε 0 Q (ra) 2 = Q 4π ε 0 r 2 ( ( 1+ a r ) 2 ( 1 a r ) 2 )

Como a<<r, aproximamos a

E p Q 4π ε 0 r 2 ( ( 1 2a r )( 1+ 2a r ) )= 4Qa 4π ε 0 r 3

El campo eléctrico Ei producido por el ión polariza el átomo, p=α·Ei, donde p es el momento dipolar p=Q(2a). Calculamos la carga Q del dipolo

p=α 1 4π ε 0 q r 2 Q= α 2a 1 4π ε 0 q r 2

El campo producido por el átomo neutro situado en el origen, en la posición del ión es

E p = 2α ( 4π ε 0 ) 2 q r 5

La fuerza sobre el ión es atractiva

f p = 2α ( 4π ε 0 ) 2 q 2 r 5

El potencial V(r) correspondiente al campo Ep es

V(r)= r 2α ( 4π ε 0 ) 2 q r 5 dr= α 2 ( 4π ε 0 ) 2 q r 4

La posición de máximo acercamiento es la raíz real positiva de la ecuación bicuadrada

r 2 b 2 + k r 2 E =0k= α q 2 2 ( 4π ε 0 ) 2

Las raíces reales se obtienen siempre que

b 4 4 k E b ( α q 2 4 ( π ε 0 ) 2 m v 0 2 ) 1/4

Trayectorias

Las ecuaciones del movimiento del ión de masa m son

m d 2 x d t 2 = 4k r 5 x r m d 2 y d t 2 = 4k r 5 y r

Representamos las trayectorias seguidas por los iones, resolviendo el sistema de ecuaciones diferenciales por procedimientos numéricos. Superponemos en líneas a trazos, la curva de energía potencial Ep(r)=-k/r4

k=1;
energia=5;
N=20; %partículas
tspan=[0,10];
hold on
potencial=@(r) -k./r.^4;
for i=1:N
    b=0.7+i/N; %parámetro de impacto
    if b>(4*k/energia)^(1/4)
        x0=[-10,sqrt(2*energia),b,0];
        fg=@(t,x)[x(2);-4*k*x(1)/(x(1)^2+x(3)^2)^3; x(4); 
-4*k*x(3)/(x(1)^2+x(3)^2)^3];
        [t,x]=ode45(fg,tspan,x0);
        plot(x(:,1),x(:,3))    
    end
end
fp=fplot(potencial,[0.05,5],'--k');
plot(-fp.XData,fp.YData,'--k')
ylim([-3,3]);
xlim([-4,4])
hold off
grid on
xlabel('x')
ylabel('y');
title('dispersión')

Referencias

Thomas Curtright, Gaurav Verma. Scattering shadows. Am. J. Phys. 92 (12), December 2024, pp. 924-930

Problema propuesto en la XLII Olimpiada Internacional de Física. Bangkok (2011)