## Looking for bug in a graphics program for plotting dipole fields

Greetings,

First, I apologize for postiing in the "General" area. This is the fisrt time I am asking for help in the MATLAB community and navagting this tyupe of forum. In looking of a MATAB program to plot both electric and potential fields of a dipole I came across a Book Chapter under the Academia profile of Darvin Messi on Numerical Methods. See the following link:

https://www.academia.edu/7995677/NUMERICAL_METHODS?email_work_card=view-paper&li=0

I typed in the code and worked through a majority of bugs, a couple due to typos in the text. One had to do with adding a symbol to the plot function. I got the electric fields portion to work perfectly. However, the electric potential plots still do not work. I have tried to the best of my ability to error trap. I was able to get a couple of points to plot but nothing more. I really like this approach which does not make use of MATLAB’s mesh or gradient functions because of the application I have in working with students. On the other hand, I do not know why this portion is not working. Any help would be greatly appreciated as I would not trouble the MATLAB community without exhausting the combination/permutations of what could be wrong.

I would be great to then keep the corrected version on this community as through my searches, a program like this has been requested by students very frequently.

Best wishes,

David.

%. Program below

plotit ( [-1 1], [-1.5 0; 1.5 0], 1, 1, 0.01, 0.01, 20, 20, 5)

function plotit(charges, location, ckEField, ckEq, DLE, DLV, NLE, NLV, PTS)

figure;

hold on

% Program for plotting the electric field lines

% and equipotential lines due to coplanar point charges

% the plot is to be within the range -5<x,y<5

%

% This is the correct usage:

% function plotit(charges, location,ckEField,ckEq,DLE,DLV,NLE,NLV,PTS)

%

% where,

% charges = a vector containing the charges

% location = a matrix where each row is a charge location

% ckEField = Flag set to 1 plots the Efield lines

% ckEq = Flag set to 1 plots the Equipotential lines

% DLE or DLV = the increment along E & V lines

% NLE = No. of E-Field lines per charge

% NLV = No. of Equipotential lines per charge

% PTS => Plots every PTS point (i.e. if PTS = 5 then plot every 5th point)

% note that constant Q/4*Pie*ErR is set equal to 1.0

% Determine the E-Field Lines

% For convenience, the starting points( XS,YS) are radially distributed about charge locations

Q=charges;

XQ = location(:,1);

YQ = location(:,2);

JJ=1;

NQ = length(charges);

if (ckEField)

for K=1:NQ

for I =1:NLE

THETA = 2*pi*(I-1)/(NLE);

XS=XQ(K)+0.1*cos(THETA);

YS=YQ(K)+0.1*sin(THETA);

XE=XS;

YE=YS;

JJ=JJ+1;

if (~mod(JJ,PTS))

plot(XE, YE, ‘k.’)

end

while (1)

% Find increment and new point (X,Y)

EX=0;

EY=0;

for J=1:NQ;

R =sqrt((XE-XQ(J))^2 + (YE – YQ(J))^2 );

EX = EX +Q(J)*(XE-XQ(J))/(R^3);

EY = EY +Q(J)*(YE-YQ(J))/(R^3);

end

E = sqrt(EX^2 + EY^2);

% CHECK FOR A SINGULAR POINT

if (E <=0.00005)

break;

end

DX = DLE*EX/E;

DY = DLE*EY/E;

% FOR NEGATIVE CHARGE, NEGATE DX & DY SO THAT INCREMENT IS AWAY FROM THE CHARGE

if (Q(K) < 0)

DX = -DX;

DY = -DY;

end

XE = XE + DX;

YE = YE + DY;

% CHECK WHETHER NEW POINT IS WITHIN THE GIVEN RANGE OR TOO

% CLOSE TO ANY OF THE POINT CHARGES – TO AVOID SINGULAR POINT

if ((abs(XE) >= 5) | (abs(YE) >= 5))

break;

end

if (sum(abs(XE-XQ) < 0.05 & abs(YE-YQ) < 0.05) > 0)

break;

end

JJ=JJ+1;

if (~mod(JJ,PTS))

plot(XE,YE,’k.’)

end

end % while loop

end % I =1:NLE

end % K = 1:NQ

end % if

% NEXT, DETERMINE THE EQUIPOTENTIAL LINES FOR CONVENIENCE, THE STARTING POINTS (XS,YS) ARE

% CHOSEN LIKE THOSE FOR THE E-FIELD LINES

if(ckEq)

JJ=1;

DELTA = 0.2;

ANGLE = 45*pi/180;

for K =1:NQ

FACTOR = 0.5

for KK = 1:NLV

XS = XQ(K) + FACTOR*cos(ANGLE);

YS = YQ(K) + FACTOR*sin(ANGLE);

if ( abs(XS) >= 5 | abs(YS) >=5 )

break;

end

DIR = 1;

XV = XS;

YV = YS;

JJ=JJ+1;

if (~mod(JJ,PTS))

plot(XV,YV, ‘rs’)

end

% FIND INCREMENT AND NEW POINT (XV,YV)

N=1;

while(1)

EX = 0;

EY = 0;

for J = 1:NQ

R = sqrt((XV-XQ(J))^2 + (YV-YQ(J))^2);

EX = EX + Q(J)*(XV-XQ(J))/(R^3);

EY = EY + Q(J)*(YV-YQ(J))/(R^3);

end

E=sqrt(EX^2 + EY^2);

if (E <= 0.00005)

FACTOR = 2*FACTOR;

break;

end;

DX = -(DLV*EX)/E;

DY = (DLV*EY)/E;

XV = XV + DIR*DX;

YV = YV + DIR*DY;

% CHECK IF THE EQUIPOTENTIAL LINE LOOPS BACK TO (X,YS)

R0 = sqrt((XV – XS)^2 + (YV – YS)^2);

if (R0 < DELTA & N < 50)

FACTOR = 2*FACTOR;

break;

end

% CHECK WHETHER NEW POINT IS WITHIN THE GIVEN RANGE IF FOUND OUT OF RANGE, GO BACK TO THE STARTING POINT

% (S,YS) BUT INCREMENT IN THE OPPOSITE DIRECTION

if (abs(XV) > 5 | abs(YV) > 5)

DIR = DIR – 2;

XV = XS;

YV = YS;

end

if (abs(DIR) > 1)

FACTOR = 2*FACTOR;

break;

end

if ( sum( abs(XV-XQ) < 0.005 & abs(YV-YQ) < 0.005) > 0 )

break;

end

end

JJ=JJ+1;

if (~mod(JJ,PTS))

N=N+1;

plot(XV,YV,’rs’)

end

end % WHILE loop

end % KK

end % K

end % ifGreetings,

First, I apologize for postiing in the "General" area. This is the fisrt time I am asking for help in the MATLAB community and navagting this tyupe of forum. In looking of a MATAB program to plot both electric and potential fields of a dipole I came across a Book Chapter under the Academia profile of Darvin Messi on Numerical Methods. See the following link:

https://www.academia.edu/7995677/NUMERICAL_METHODS?email_work_card=view-paper&li=0

I typed in the code and worked through a majority of bugs, a couple due to typos in the text. One had to do with adding a symbol to the plot function. I got the electric fields portion to work perfectly. However, the electric potential plots still do not work. I have tried to the best of my ability to error trap. I was able to get a couple of points to plot but nothing more. I really like this approach which does not make use of MATLAB’s mesh or gradient functions because of the application I have in working with students. On the other hand, I do not know why this portion is not working. Any help would be greatly appreciated as I would not trouble the MATLAB community without exhausting the combination/permutations of what could be wrong.

I would be great to then keep the corrected version on this community as through my searches, a program like this has been requested by students very frequently.

Best wishes,

David.

%. Program below

plotit ( [-1 1], [-1.5 0; 1.5 0], 1, 1, 0.01, 0.01, 20, 20, 5)

function plotit(charges, location, ckEField, ckEq, DLE, DLV, NLE, NLV, PTS)

figure;

hold on

% Program for plotting the electric field lines

% and equipotential lines due to coplanar point charges

% the plot is to be within the range -5<x,y<5

%

% This is the correct usage:

% function plotit(charges, location,ckEField,ckEq,DLE,DLV,NLE,NLV,PTS)

%

% where,

% charges = a vector containing the charges

% location = a matrix where each row is a charge location

% ckEField = Flag set to 1 plots the Efield lines

% ckEq = Flag set to 1 plots the Equipotential lines

% DLE or DLV = the increment along E & V lines

% NLE = No. of E-Field lines per charge

% NLV = No. of Equipotential lines per charge

% PTS => Plots every PTS point (i.e. if PTS = 5 then plot every 5th point)

% note that constant Q/4*Pie*ErR is set equal to 1.0

% Determine the E-Field Lines

% For convenience, the starting points( XS,YS) are radially distributed about charge locations

Q=charges;

XQ = location(:,1);

YQ = location(:,2);

JJ=1;

NQ = length(charges);

if (ckEField)

for K=1:NQ

for I =1:NLE

THETA = 2*pi*(I-1)/(NLE);

XS=XQ(K)+0.1*cos(THETA);

YS=YQ(K)+0.1*sin(THETA);

XE=XS;

YE=YS;

JJ=JJ+1;

if (~mod(JJ,PTS))

plot(XE, YE, ‘k.’)

end

while (1)

% Find increment and new point (X,Y)

EX=0;

EY=0;

for J=1:NQ;

R =sqrt((XE-XQ(J))^2 + (YE – YQ(J))^2 );

EX = EX +Q(J)*(XE-XQ(J))/(R^3);

EY = EY +Q(J)*(YE-YQ(J))/(R^3);

end

E = sqrt(EX^2 + EY^2);

% CHECK FOR A SINGULAR POINT

if (E <=0.00005)

break;

end

DX = DLE*EX/E;

DY = DLE*EY/E;

% FOR NEGATIVE CHARGE, NEGATE DX & DY SO THAT INCREMENT IS AWAY FROM THE CHARGE

if (Q(K) < 0)

DX = -DX;

DY = -DY;

end

XE = XE + DX;

YE = YE + DY;

% CHECK WHETHER NEW POINT IS WITHIN THE GIVEN RANGE OR TOO

% CLOSE TO ANY OF THE POINT CHARGES – TO AVOID SINGULAR POINT

if ((abs(XE) >= 5) | (abs(YE) >= 5))

break;

end

if (sum(abs(XE-XQ) < 0.05 & abs(YE-YQ) < 0.05) > 0)

break;

end

JJ=JJ+1;

if (~mod(JJ,PTS))

plot(XE,YE,’k.’)

end

end % while loop

end % I =1:NLE

end % K = 1:NQ

end % if

% NEXT, DETERMINE THE EQUIPOTENTIAL LINES FOR CONVENIENCE, THE STARTING POINTS (XS,YS) ARE

% CHOSEN LIKE THOSE FOR THE E-FIELD LINES

if(ckEq)

JJ=1;

DELTA = 0.2;

ANGLE = 45*pi/180;

for K =1:NQ

FACTOR = 0.5

for KK = 1:NLV

XS = XQ(K) + FACTOR*cos(ANGLE);

YS = YQ(K) + FACTOR*sin(ANGLE);

if ( abs(XS) >= 5 | abs(YS) >=5 )

break;

end

DIR = 1;

XV = XS;

YV = YS;

JJ=JJ+1;

if (~mod(JJ,PTS))

plot(XV,YV, ‘rs’)

end

% FIND INCREMENT AND NEW POINT (XV,YV)

N=1;

while(1)

EX = 0;

EY = 0;

for J = 1:NQ

R = sqrt((XV-XQ(J))^2 + (YV-YQ(J))^2);

EX = EX + Q(J)*(XV-XQ(J))/(R^3);

EY = EY + Q(J)*(YV-YQ(J))/(R^3);

end

E=sqrt(EX^2 + EY^2);

if (E <= 0.00005)

FACTOR = 2*FACTOR;

break;

end;

DX = -(DLV*EX)/E;

DY = (DLV*EY)/E;

XV = XV + DIR*DX;

YV = YV + DIR*DY;

% CHECK IF THE EQUIPOTENTIAL LINE LOOPS BACK TO (X,YS)

R0 = sqrt((XV – XS)^2 + (YV – YS)^2);

if (R0 < DELTA & N < 50)

FACTOR = 2*FACTOR;

break;

end

% CHECK WHETHER NEW POINT IS WITHIN THE GIVEN RANGE IF FOUND OUT OF RANGE, GO BACK TO THE STARTING POINT

% (S,YS) BUT INCREMENT IN THE OPPOSITE DIRECTION

if (abs(XV) > 5 | abs(YV) > 5)

DIR = DIR – 2;

XV = XS;

YV = YS;

end

if (abs(DIR) > 1)

FACTOR = 2*FACTOR;

break;

end

if ( sum( abs(XV-XQ) < 0.005 & abs(YV-YQ) < 0.005) > 0 )

break;

end

end

JJ=JJ+1;

if (~mod(JJ,PTS))

N=N+1;

plot(XV,YV,’rs’)

end

end % WHILE loop

end % KK

end % K

end % if Greetings,

First, I apologize for postiing in the "General" area. This is the fisrt time I am asking for help in the MATLAB community and navagting this tyupe of forum. In looking of a MATAB program to plot both electric and potential fields of a dipole I came across a Book Chapter under the Academia profile of Darvin Messi on Numerical Methods. See the following link:

https://www.academia.edu/7995677/NUMERICAL_METHODS?email_work_card=view-paper&li=0

I typed in the code and worked through a majority of bugs, a couple due to typos in the text. One had to do with adding a symbol to the plot function. I got the electric fields portion to work perfectly. However, the electric potential plots still do not work. I have tried to the best of my ability to error trap. I was able to get a couple of points to plot but nothing more. I really like this approach which does not make use of MATLAB’s mesh or gradient functions because of the application I have in working with students. On the other hand, I do not know why this portion is not working. Any help would be greatly appreciated as I would not trouble the MATLAB community without exhausting the combination/permutations of what could be wrong.

I would be great to then keep the corrected version on this community as through my searches, a program like this has been requested by students very frequently.

Best wishes,

David.

%. Program below

plotit ( [-1 1], [-1.5 0; 1.5 0], 1, 1, 0.01, 0.01, 20, 20, 5)

function plotit(charges, location, ckEField, ckEq, DLE, DLV, NLE, NLV, PTS)

figure;

hold on

% Program for plotting the electric field lines

% and equipotential lines due to coplanar point charges

% the plot is to be within the range -5<x,y<5

%

% This is the correct usage:

% function plotit(charges, location,ckEField,ckEq,DLE,DLV,NLE,NLV,PTS)

%

% where,

% charges = a vector containing the charges

% location = a matrix where each row is a charge location

% ckEField = Flag set to 1 plots the Efield lines

% ckEq = Flag set to 1 plots the Equipotential lines

% DLE or DLV = the increment along E & V lines

% NLE = No. of E-Field lines per charge

% NLV = No. of Equipotential lines per charge

% PTS => Plots every PTS point (i.e. if PTS = 5 then plot every 5th point)

% note that constant Q/4*Pie*ErR is set equal to 1.0

% Determine the E-Field Lines

% For convenience, the starting points( XS,YS) are radially distributed about charge locations

Q=charges;

XQ = location(:,1);

YQ = location(:,2);

JJ=1;

NQ = length(charges);

if (ckEField)

for K=1:NQ

for I =1:NLE

THETA = 2*pi*(I-1)/(NLE);

XS=XQ(K)+0.1*cos(THETA);

YS=YQ(K)+0.1*sin(THETA);

XE=XS;

YE=YS;

JJ=JJ+1;

if (~mod(JJ,PTS))

plot(XE, YE, ‘k.’)

end

while (1)

% Find increment and new point (X,Y)

EX=0;

EY=0;

for J=1:NQ;

R =sqrt((XE-XQ(J))^2 + (YE – YQ(J))^2 );

EX = EX +Q(J)*(XE-XQ(J))/(R^3);

EY = EY +Q(J)*(YE-YQ(J))/(R^3);

end

E = sqrt(EX^2 + EY^2);

% CHECK FOR A SINGULAR POINT

if (E <=0.00005)

break;

end

DX = DLE*EX/E;

DY = DLE*EY/E;

% FOR NEGATIVE CHARGE, NEGATE DX & DY SO THAT INCREMENT IS AWAY FROM THE CHARGE

if (Q(K) < 0)

DX = -DX;

DY = -DY;

end

XE = XE + DX;

YE = YE + DY;

% CHECK WHETHER NEW POINT IS WITHIN THE GIVEN RANGE OR TOO

% CLOSE TO ANY OF THE POINT CHARGES – TO AVOID SINGULAR POINT

if ((abs(XE) >= 5) | (abs(YE) >= 5))

break;

end

if (sum(abs(XE-XQ) < 0.05 & abs(YE-YQ) < 0.05) > 0)

break;

end

JJ=JJ+1;

if (~mod(JJ,PTS))

plot(XE,YE,’k.’)

end

end % while loop

end % I =1:NLE

end % K = 1:NQ

end % if

% NEXT, DETERMINE THE EQUIPOTENTIAL LINES FOR CONVENIENCE, THE STARTING POINTS (XS,YS) ARE

% CHOSEN LIKE THOSE FOR THE E-FIELD LINES

if(ckEq)

JJ=1;

DELTA = 0.2;

ANGLE = 45*pi/180;

for K =1:NQ

FACTOR = 0.5

for KK = 1:NLV

XS = XQ(K) + FACTOR*cos(ANGLE);

YS = YQ(K) + FACTOR*sin(ANGLE);

if ( abs(XS) >= 5 | abs(YS) >=5 )

break;

end

DIR = 1;

XV = XS;

YV = YS;

JJ=JJ+1;

if (~mod(JJ,PTS))

plot(XV,YV, ‘rs’)

end

% FIND INCREMENT AND NEW POINT (XV,YV)

N=1;

while(1)

EX = 0;

EY = 0;

for J = 1:NQ

R = sqrt((XV-XQ(J))^2 + (YV-YQ(J))^2);

EX = EX + Q(J)*(XV-XQ(J))/(R^3);

EY = EY + Q(J)*(YV-YQ(J))/(R^3);

end

E=sqrt(EX^2 + EY^2);

if (E <= 0.00005)

FACTOR = 2*FACTOR;

break;

end;

DX = -(DLV*EX)/E;

DY = (DLV*EY)/E;

XV = XV + DIR*DX;

YV = YV + DIR*DY;

% CHECK IF THE EQUIPOTENTIAL LINE LOOPS BACK TO (X,YS)

R0 = sqrt((XV – XS)^2 + (YV – YS)^2);

if (R0 < DELTA & N < 50)

FACTOR = 2*FACTOR;

break;

end

% CHECK WHETHER NEW POINT IS WITHIN THE GIVEN RANGE IF FOUND OUT OF RANGE, GO BACK TO THE STARTING POINT

% (S,YS) BUT INCREMENT IN THE OPPOSITE DIRECTION

if (abs(XV) > 5 | abs(YV) > 5)

DIR = DIR – 2;

XV = XS;

YV = YS;

end

if (abs(DIR) > 1)

FACTOR = 2*FACTOR;

break;

end

if ( sum( abs(XV-XQ) < 0.005 & abs(YV-YQ) < 0.005) > 0 )

break;

end

end

JJ=JJ+1;

if (~mod(JJ,PTS))

N=N+1;

plot(XV,YV,’rs’)

end

end % WHILE loop

end % KK

end % K

end % if graphics and visualization pl MATLAB Answers — New Questions