Here is a summary of what I am trying to do:
Use MATLAB to compute the potential $V$ at any point $(x, y, z)$ in space due to a uniform ring of charge. Use a Riemann sum to compute the integral with increments, $N$, as a variable you can change. Plot your potential and field in the plane perpendicular to the area of the ring and passing through the center.
I have two versions of code that give me the same result: correct expression for $V$ ($V_\text{tot}$) and the correct vector field. My problem is getting my contour plot to fill the 2D space. The only difference between the two codes is that one uses a for loop to perform the summation and the other uses a sum command:
Version 1 (using for loop):
%% Computing a symbolic expression for V for anywhere in space
syms x y z % phiprime is angle that an elemental dq of the circular charge is located at, x,y and z are arbitrary points in space outside the charge distribution
N = 200; % number of increments to sum
R = 2; % radius of circle is 2 meters
dphi = 2*pi/N; % discretizing the circular line of charge which spans 2pi
integrand = 0;
for phiprime = 0:dphi:2*pi
% phiprime ranges from 0 to 2pi in increments of dphi
integrand = integrand + dphi./(sqrt(((x - R.*cos(phiprime) )).^2 + ((y - R.*sin(phiprime) ).^2) + z.^2));
end
intgrl = sum(integrand);
% unnecessary but harmless step that I leave to show that I am using the
sum of the above expression for each dphi
eps0 = 8.854e-12;
kC = 1/(4*pi*eps0);
rhol = 1e-9; % linear charge density
Vtot = kC*rhol*R.*intgrl; % symbolic expression for Vtot
%% Graphing V & E in plane perpedicular to ring & passing through center
[Y1, Z1] = meshgrid(-4:.5:4, -4:.5:4);
Vcont1 = subs(Vtot, [x,y,z], {0,Y1,Z1}); % Vcont1 stands for V contour; 1 is because I do the plane of the ring next
contour(Y1,Z1,Vcont1)
xlabel('y - axis [m]')
ylabel('z - axis [m]')
title('V in a plane perpendicular to a ring of charge (N = 200)')
str = {'Red line is side view', 'of ring of charge'};
text(-1,2,str)
hold on
% visually displaying line of charge on plot
circle = rectangle('Position',[-2 0 4 .1],'Curvature',[1,1]);
set(circle,'FaceColor',[1, 0, 0],'EdgeColor',[1, 0, 0]);
% taking negative gradient of V and finding symbolic equations for Ex, Ey and Ez
g = gradient(-1.*(kC*rhol*R.*intgrl),[x,y,z]);
% now substituting all the values of the 2D coordinate system for the symbolic x and y variables to get numeric values for Ex and Ey
Ey1 = subs(g(2), [x y z], {0,Y1,Z1});
Ez1 = subs(g(3), [x y z], {0,Y1,Z1});
E1 = sqrt(Ey1.^2 + Ez1.^2); % full numeric magnitude of E in y-z plane
Eynorm1 = Ey1./E1; % This normalizes the electric field lines
Eznorm1 = Ez1./E1;
quiver(Y1,Z1,Eynorm1,Eznorm1);
hold off
Version 2 (using sum):
syms x y z
R = 2; % radius of circle is 2 meters
N=100;
dphi = 2*pi/N; % discretizing the circular line of charge which spans 2pi
phiprime = 0:dphi:2*pi; %phiprime ranges from 0 to 2pi in increments of dphi
integrand = dphi./(sqrt(((x - R.*cos(phiprime) )).^2 + ((y - R.*sin(phiprime) ).^2) + z.^2));
phiprime = 0:dphi:2*pi;
intgrl = sum(integrand); % Riemann sum performed here
eps0 = 8.854e-12;
kC = 1/(4*pi*eps0);
rhol = 1e-9; % linear charge density
Vtot = kC*rhol*R.*intgrl; % symbolic expression for Vtot
%%
[Y1, Z1] = meshgrid(-4:.5:4,-4:.5:4);
Vcont1 = subs(Vtot, [x,y,z], {0,Y1,Z1});
contour(Y1,Z1,Vcont1)
xlabel('y - axis [m]')
ylabel('z - axis [m]')
title('V in a plane perpedicular to a ring of charge (N = 100)')
str = {'Red line is side view', 'of ring of charge'};
text(-1,2,str)
hold on
circle = rectangle('Position',[-2 0 4 .1],'Curvature',[1,1]); % visually displaying ring of charge on plot
set(circle,'FaceColor',[1, 0, 0],'EdgeColor',[1, 0, 0]);
g = gradient(-1.*(kC*rhol*R.*intgrl),[x,y,z]); % taking negative gradient of V and finding symbolic equations for Ex, Ey and Ez
% substituting all the values of the 2D coordinate system for the symbolic x and y variables to get numeric values for Ex and Ey because the gradient command doesn't accept symbolic arguments
Ey1 = subs(g(2), [x y z], {0,Y1,Z1});
Ez1 = subs(g(3), [x y z], {0,Y1,Z1});
E1 = sqrt(Ey1.^2 + Ez1.^2); % full numeric magnitude of E in y-z plane
Eynorm1 = Ey1./E1; % This normalizes the electric field lines
Eznorm1 = Ez1./E1;
quiver(Y1,Z1,Eynorm1,Eznorm1);
hold off
Both versions of code produce the following graphs:
Note: the picture above this text should have the axes be $y$ and $z$ like the picture below, not $x$ and $y$. Also, the title of the picture below this text should be "$E$ in the plane..." not $V$.
As you can see, the vector field is correct while the contour plot seems to use only a few points around the ends of the ring and connect them with straight lines in a strange diamond shape. I can't get it to fill space.
As for my derivation of the formula for $V$, it is here:
- $s$ is the radial distance
- $\phi$ is the azimuthal angle
- An arbitrary point $P(x,y,z)$ has no superscript while the prime notation is used to identify a point on the ring of charge ($\phi^\prime$, $x^\prime$, etc)
- cursive $r$ indicates relative distance from a point on the ring to $P(x,y,z)$
- $R$ is the radius of the ring, I chose 2 meters
- $\mathrm{d}s = R \mathrm{d}\phi$ is a small amount of arc length
