function dAdt=ode(A)

global eps sigma N m radius U

dAdt=zeros(4*N,1);
F=zeros(4*N,1);

for i=1:4:4*N-3
    Fi=[F(i+2);F(i+3)];
    for j=1:4:4*N-3
        R1=[A(i);A(i+1)];
        R2=[A(j);A(j+1)];
        if i~=j
            dF=-4*eps*(-12*sigma^12*(norm(R1-R2))^(-13)+6*sigma^6*(norm(R1-R2))^(-7))*(R1-R2)/norm(R1-R2);
            Fi=Fi+dF;
        end
    end
    F(i+2)=Fi(1);
    F(i+3)=Fi(2);
    if norm([A(i);A(i+1)]-U)<=radius
        F(i+2)=0;
        F(i+3)=0;
    end
end

dAdt=dAdt+F/m;

for i=1:4:4*N-3
    dAdt(i)=A(i+2);
    dAdt(i+1)=A(i+3);
end
    


        
    