clear all;

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

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

plot(q1_pos(1),q1_pos(2),'bo');
plot(q3_pos(1),q3_pos(2),'ro');

[x,y] = meshgrid([-2:0.5:2]);


    r = sqrt((x-q1_pos(1)).^2+(y-q1_pos(2)).^2);
    Er = k*q1./r^2;
    rx = (x-q1_pos(1))./r;
    ry = (y-q1_pos(2))./r;
    Ex = Er.*rx;
    Ey = Er.*ry;
    r = sqrt((x-q3_pos(1)).^2+(y-q3_pos(2)).^2);
    Er = k*q3./r^2;
    rx = (x-q3_pos(1))./r;
    ry = (y-q3_pos(2))./r;
    Ex = Ex + Er.*rx;
    Ey = Ey + Er.*ry;

    %“dŠEƒxƒNƒgƒ‹‚Ì•\Ž¦
figure(1);
quiver(x,y,Ex,Ey);
