clear all;

v0 = [0,0];
e0 = 8.854e-12;
k = 1/(4*pi*e0);
q1 = 1e-9;

%axis([-10 10 -10 10]);
figure(1);
hold on

for x= -1:1
    for y=-1:1
        pos = [x,y];
        old_pos = [x,y];
        v = v0;
        plot(pos(1),pos(2),'r.');
        for t=0:10
            r = sqrt(pos(1)^2+pos(2)^2);
            Er = k*q1/r^2;
            rx = pos(1)/r;
            ry = pos(2)/r;
            Evec(1) = Er*rx;
            Evec(2) = Er*ry;
            
            plot([pos(1),pos(1)+Evec(1)],[pos(2),pos(2)+Evec(2)],'r-');
            pause(0.2);
            
            old_pos(1) = pos(1);
            old_pos(2) = pos(2);
            
            pos(1) = pos(1)+Evec(1);
            pos(2) = pos(2)+Evec(2);
            
        end
    end
end
