%
%   Book: MATLAB-Based Electromagnetics (Pearson Prentice Hall)
%   Author: Branislav M. Notaros
%   Instructor Resources
%   (c) 2011
%
%   This MATLAB code or any part of it may be used only for 
%   educational purposes associated with the book
%
%
%


%   Numerical proof that E-field is conservative -- movie


clear all;
close all;
EPS0 = 8.8542*10^(-12); 

Q = input('\nEnter the charge in Coulombs for point at the center of the square: ');
a = input('Enter the lenght in meters of square side : ');

% Division of square sides into N segments
N = 300; % It should be 4*k, where k is integer
n = N/4; % each side
delta = a/n;

% Coordinates of patch centers
temp = -a/2+delta/2 : delta : a/2-delta/2;
x = [-a/2*ones(1,n),temp,a/2*ones(1,n),-temp];
y = [temp,a/2*ones(1,n),-temp,-a/2*ones(1,n)];

    r = sqrt(x.^2 + y.^2);  % Distance of the center of the patch from the charge point
    ux = x./r;  % X component of unit vector of radius vector 
    uy = y./r;  % Y component of unit vector of radius vector 

% dl vector for integration
dx = [zeros(1,n),delta * ones(1,n),zeros(1,n),-delta * ones(1,n)] ; 
dy = [delta * ones(1,n),zeros(1,n),-delta * ones(1,n),zeros(1,n)];
dl = [dx;dy];

% Integration 
Int = 0; 
plot(x,y,'g');
hold on;
plot(0,0,'m*')
hold on;
for i = 1:N
    E(:,i) = Q/r(i)^2/(4*pi*EPS0)*[ux(i),uy(i)];  % Electric field vector
    if mod(i,20)==0
    vecPlot2D([x(i),y(i)],[x(i)+E(1,i),y(i)+E(2,i)],abs(1/E(1,1))/10,'b',0);
    hold on;
    end;
end;
legend('Square path' ,'Charge', 'Electric field vectors');
title('Integration in progress...');
for i = 1:N
    Int = Int + dotProduct(E(:,i),dl(:,i));
    %movie
    plot(x(i),y(i),'r*');
    hold on;
    axis equal;
    axis off;
    if i==N
        title(['Result of integral of electric field'...
            'along closed contour is 0!']);
    end;
    xlim([-1.5*a,1.5*a]);
    ylim([-1.5*a,1.5*a]);
    M(i)= getframe();
end
hold off;

Vvertex = Q/(4*pi*EPS0*a/2*sqrt(2));

fprintf(['\nPotential at square vertex due to charge Q'...
    'is %e V while circulation of E is %e V'],Vvertex,Int);
fprintf(['\nBeing aware of some numerical errors, we conclude that the'...
    'obtained result of integration is negligibly small -- zero' ]);  


