ScalarFields
(Electric Potential)

Version 2.3
2/18/92
John B. Schneider

Introduction

This Notebook examines electric potential for static systems consisting of point charges.  Static electric potential is one of the simplest fields in electromagnetics.  It does not depend on time and, since it is a scalar field, the value at each point is merely a number (as opposed to a vector).  Electric potential is also quite easy to interpret physically.  Most people can imagine moving a small positive test charge through a potential field.  If the charge has moved to a position where the potential is greater than where it started, then work has been done by you.  You can, of course, think of this in the opposite way.  If you move a small test charge in a potential field and you have to do work to move the charge, then the potential of that charge has increased.  Likewise, if the charge drags you around – the charge is doing work on you – then the charge will lose some of its potential.

   Despite the fact that the static electric potential is so simple to understand, it can still be difficult to visualize.  It is important that you have a good feeling for static, scalar fields before we consider time-varying or vector fields.  We will examine all the common ways that are used to depict these fields.   Several functions are used in this notebook that are not predefined Mathematica functions.  These functions are defined in the subsection entitled ™Initialization Cells,∫ at the end of this NoteBook.

Surface Plots

Introduction

We will begin by considering the electric potential about a single point charge, a monopole, located at the origin.  The potential field for this charge exists everywhere in the space about the charge.  The challenge is to represent this potential in such a way as to get a clear picture of its behavior.  Probably the best way to do this is to employ surface plots.  With a surface plot, the positions where the potential is determined are in a given plane – say the plane define by x = 0.  We then plot the potential as a function of the two variables which fix our position in that particular plane.  If we were considering the x = 0 plane, then specifing the y and z values would uniquely determine a point, and we plot our potential at that point as a height above (or below)  a plane.
       
   After the monopole, we will examine the potential about a dipole and a quadrupole.  There is a short animated sequence at the end of the dipole subsection which illustrates how the potential varies as the evaluation plane (the plane over which the potential is measured) is moved.

Monopole

The following Mathematica commands generates a surface plot of the potential about a single positive charge.  In these commands, the function Monopole[x,y,z] gives the electric potential at a point (x,y,z) when the point charge is located at the origin.  We restrict ourselves to the plane that has an x value of zero, and y and z each vary over the range ±4 to 4 (x=0, ±4 ≤ y ≤ 4, ±4 ≤ z ≤ 4).  We are not concerned with units yet so we will say that y and z vary over the range ±4 anythings to +4 anythings.  If we were to pay attention to units, that would not change the form of the output – we would just have to be more careful in our labeling.  If you want more details about the following Mathematica command, you can double click on the cell directly below the plot.  To see a diagrammatical representation of the monopole and the evaluation plane, click on the cell below the explanation cell.  Note:  If you don't know how to use the following command to generate output, you should review the Introduction NoteBook!

charge = Pi Epsilon0;
location = {0.,0.,0.};
Plot3D[ Monopole[charge,location,{0,y,z}], {y,-4,4}, {z,-4,4},
   ClipFill -> None, Lighting -> True, PlotPoints -> 31,
   PlotRange -> {Automatic, Automatic, {-.2,.2}},
   BoxRatios -> {1,1,1}, AxesLabel -> {"y","z","V"},
   PlotLabel -> "Fig. 1"]

The first line sets the variable charge to a convenient value – this gives the strength of the point charge.  The second line sets the variable location to the origin. The third line contains a plotting command that plots Monopole as a funtion of y and z.  All the entries in the command with an arrow (->) specify options.

     The options only need to be specified if you want to change their default values.  For the above plot, the ClipFill option tells Mathematica to leave the truncated values blank.  Setting the Lighting option to True shades the plot as if  it were illuminated by a real light source.  If this was omitted or set to False, the shading would be determined by the height of the surface.   The PlotPoints option tells Mathematica how many points to use along the two independent variable axes.   The PlotRange option is set so that Mathematica will use what it ™thinks∫ are the best ranges for the two independent axes and the range of the dependent variable is set to ±.2 to .2 (for reasons of comparison with plots that we will generate later).  The BoxRatios option sets all sides of the cube which contains the surface plot to have equal lengths.  The AxesLabel and PlotLabels option does the obvious – they provide a labels for the axes and the plot.

     Note that Mathematica has a little bit of difficulty plotting the truncated edge of the function.  The edge of the ™volcano∫ should be a circle, but since we are using a finite number of points to plot this function, some jaggedness is present in the plot.  The jaggedness could be further reduced by increasing the number of PlotPoints, but  this would come at the cost of increased computation time.  The number of PlotPoints chosen here more than adequately illustrates the form of the function.  On some of the following plots the jaggedness where the function is truncated becomes more pronounced, and you may want to increase the number of PlotPoints – feel free to do so, but be prepared to wait.

Show[ EvaluationPlane[0.], axes, mpoint, 
   Lighting -> False,
   AxesLabel -> {"x","y","z"}]

In this sketch, the evaluation plane is shown as the shaded plane which runs through the y and z axes.  The charge is represented as a sphere at the origin.  For the purpose of illustration, this sphere is given a radius which is clearly non-zero.  For a true point charge, the charge is completely contained in a sphere of vanishingly small radius.

Notice that Fig. 1 looks something like a volcano.  In fact, the surface plot should go to infinity at the origin.   Since a function obviously can't be plotted all the way to infinity, the plot is truncated at a ™reasonable∫ value.
         
   Given this suface plot, you can get a better feel for what potential represents.   Think of yourself as a positive test charge.  Your potential is given by  your height above the ™ground∫ plane.  This means that if you want to increase your potential, no matter where you are in this plane, you would have to move toward the center of the plane (meaning the origin in this case).  You would have to climb up the hill, doing work.  If you want to decrease your potential, you can let the field do the work by moving away from the center and just slide down the hill.  Sliding down the hill means that your potential decreases as you move out from the center, while in coordinate space you would be moving away from the origin and staying in the x = 0 plane.  You must keep in mind that the diagram above is just for points in the x = 0 plane.  If you move out of this plane, a new surface plot may pertain to the description of your potential.

   To illustrate what would happen if you were in the plane x = 2, go to the Mathematica command which generated Fig. 1, and change the value of x in the variable location to 2 (so that the second line reads location = {2.,0.,0.};).  After you have done this, hit the Enter key (not the Return key); wait for the output; after giving the output the once over, return here and read on.

       Notice that the volcano has been replaced by a smooth bump.  The potential no longer goes to infinity at the origin because the charge is not located at the center of this plane.  When you are at the ™center∫ of the x = 2 plane, you are at the point (2,0,0) but the charge is at the origin (0,0,0) – you are still a distance of 2 away from the charge.  You should also notice that, once again, this plot tells us that no matter where we are in this plane, if we want to increase our potential, then we must move towards the center.  But, keep in mind, that this only tells us how things behave in the plane defined by x = 2.  It doesn't necessarily tell us the best direction in general in which to move to increase our potential.  It only tells us the best direction in this particular plane.  In fact, for the simple case of a single positive charge at the origin, it should be clear that moving directly toward the charge at the origin increases potential most rapidly.  We will consider this more when we look at electric fields.

Dipole

A static electric dipole consists of two charges that have the same magnitude and opposite sign that are a fixed distance apart.  We will consider the case where there is a positive charge +Q at the point (0,0,1) and a negative charge ±Q at the point (0,0,±1).  
       
   In the following Mathematica commands, the function Dipole[x,y,z] serves the same role as the function Monopole[x,y,z] did in Fig. 1.  The dipole function gives the potential at a point (x,y,z) about a dipole.  We will start by generating a surface plot of the potential in the x = 0 plane.  If you want more information on the Mathematica commands that generate this plot, go back to the cell following Fig. 1 (in ßSurface_Plots.Monopole).  If you want the diagrammatical representation of the physical distribution of charge and the evaluation plane, double click on the cell following the Fig. 2.

charge = Pi Epsilon0;
location1 = {0,0,1};
location2 = {0,0,-1};
Plot3D[Dipole[ charge, location1, location2, {0,y,z}],
   {y,-4,4}, {z,-4,4},
   ClipFill -> None, Lighting -> True, PlotPoints -> 31,
   PlotRange -> {Automatic, Automatic, {-.2,.2}},
   BoxRatios -> {1,1,1}, AxesLabel -> {"y","z","V"},
   PlotLabel -> "Fig. 2"]

Show[ EvaluationPlane[0.], axes, dpoints, 
   Lighting -> False, AxesLabel -> {"x","y","z"}]

The evaluation plane is the shaded portion of the x = 0 plane.  The point charges are represented by spheres along the z axis.  The positive charge is black and is above the negative charge (which is ™dark∫ gray) as seen from this angle.

There are a couple of important points about potential fields which are illustrated by this plot that could not have been easily demonstrated with the surface plots obtained from the monopole.  First, you should see that there are now two singular points (points where the function blows up).  As you approach the singularity at the positive charge, the potential goes to positive infinity; as you approach the negative charge, the potential goes to negative infinity.  This points out that ™zero potential∫ is a relative value.
              
   Again, imagine yourself as a positive test charge.  Even if you found yourself at a point of zero potential, you can still loose even more potential.  For these plots, the potential is referenced to the potential found at point infinitely far from the point charges.  If you were infinitely far away from these charges, we would say that you had zero potential.  However, there are points close to the charges that also have zero potential.  If you started infinitely far away and moved to one of these points, you would have gotten as much energy out of the system as you put into it.  Also, in some cases, there is a path along which all points are at zero potential.  This path is a line of equipotential at a potential of zero.  In Fig. 2 there is an equipotential line at zero potential that bisects (at a right angle) the line joining the two charges.  As a matter of fact, this line corresponds to the y axis – all points on the y axis are at zero potential.  This is visible in Fig. 2 as the straight line that runs from left to right, between the two singularities.
       
   The fact that the y axis is an equipotential line means that once you are on the y axis, you can move to any other point on the y axis without expending any energy.  On the y axis, the negative charge attracts you towards the origin just as much as the positive charge repels you.  Since you are being pulled in one direction just as hard as you are being pushed in the opposite direction, you can move about freely along that direction.  Be careful not to confuse the attraction and repulsion of the charges along the y axis with the forces in another direction.  If you move off the y axis, you will either gain or loose potential.
       
   If you were at the origin, you would have zero potential.  You could have gotten there in one of two ways.   You might have gone from infinitely far away, all the way to the origin, by moving just along the y axis.  In this case you never would have expended any energy.  Alternatively, you might have arrived at the origin by moving along any path which includes some segments not on along the y axis.  In this case, along some stretches you will have to work to move against the field (go up hill) and along other segments the field will do work on you (go down the hill).  But, no matter which route you take, when you arrive at the origin, you would have gotten as much energy out of the system as you put in.
       
       When you are at the origin, what will happen if you just let yourself go?  If you were a positive test charge (and you have never been anything but), you would go sailing towards the negative charge at the point (0,0,±1).  Despite the fact that you have zero potential, you still have plenty of enery to loose on your way down to a potential of negative infinity.  You would go down the potential ™hole∫ that is quite visible in Fig. 2.
       
       As in the case of the monopole, the plot in Fig. 2 is only for one particular plane, the x = 0 plane.  Feel free to change the values in the equation for Fig. 2 to answer whatever other questions you might have about the potential field for a dipole.

Moving the evaluation plane

If you felt no need to go back and change the values in the equation that generated Fig. 2, hopefully that was because you felt that you already understood everything, and not because you were being lazy.  Well, whether you are on top of things or you are lazy, we will now piece together some of concepts involved here and generate a simple animated sequence to test your understanding of these static fields.
       
   We will consider the same dipole as used in Fig. 2 except that now the evaluation plane will be moved.  The evaluation plane will be kept parallel to the x = 0 plane but the location where it intersects the x axis will vary.  Specifically, we will move the evaluation plane between the planes x = 0 and x = ±2  in small enough steps that we can observe gradual changes in electric potential with this movement.
       
   Before you jump ahead to the animation, you should try to predict how the field varies with the motion of the evaluation plane.  Also, to make sure you understand the physical layout which underlies this animation, be sure you first look at the animation in the cell for the diagrammatical representation.  Note:  If you don't know how to control an animation, review the material in the Introduction NoteBook.
          
   An explanation of the Mathematica command used to generate the animations is available below.

Double click on this cell for a diagrammatical representation of the moving of the evaluation plane

The following command generates the frames for the animation that shows the movement of the evaluation plane.  The positive charge is black and the negative charge is gray.  To animate this, select the frames and hit Command y.

For[x = 0, x >= -2., x -= 0.25,
   Show[ EvaluationPlane[x], axes, dpoints, 
   Lighting -> False,
   PlotRange -> {{-4,4}, {-4,4}, {-4,4}},
   BoxRatios -> {1,1,1}, AxesLabel -> {"x","y","z"}] ]

Double click on this cell for the electric potential as the evaluation plane is moved

Since the frames in this animation take so long to generate, they have already been generated for you.  You can view the animation by selecting the frames and then typing Command y.  Feel free to vary the command that generated the frames and generate a new set.  However, if you do this, you should be prepared to wait.

charge = Pi Epsilon0;
location1 = {0,0,1};
location2 = {0,0,-1};
For[x = 0, x >= -2., x -= 0.25,
   Plot3D[Dipole[charge, location1, location2, {x,y,z}],
   {y,-4,4}, {z,-4,4},
   ClipFill -> None, Lighting -> True, PlotPoints -> 31,
   PlotRange -> {Automatic, Automatic, {-.2,.2}},
   BoxRatios -> {1,1,1}, AxesLabel -> {"y","z","V"}]]

Both the animation of the movement of the evaluation plane and the changing potential use the Mathematica command For.  The syntax of the For command is For[start_condition,test,increment,expression].  In these animations, the starting condition is x = 0 and the test is ±2.0 ≤ x– the For command will stop executing when x is no longer greater than or equal to ±2.0.  The increment command is x±=.25.  This is a shorthand notation for x = x ± 0.25 (this is the same notation used in the programming language C).  The expression which the For command executes is Plot3D.   In this case, the For command starts with x = 0;  it executes the plotting command;  it then does the increment command, setting x = ±0.25;  finally, the test is done to see if execuation of the command should stop.  Since x = ±0.25  is greater than ±2.0 the test is satisfied and the For command cycles through the instructions again.  It now increments x so that it equals ±0.5.  The plot is then redrawn with this x value and so on.

Quadrupole

Hopefully, you now have a crystal clear understanding of surface plots of electric potential about a system of point charges.  However, at the risk of overkill, we will go just one step further and discuss the quadrupole.  A quadrupole consists of two positive and two negative charges, all with the same magnitude.  These charges are placed at the corners of a square, with the positive charges on the corners linked by one diagonal, and the negative charges on the other. 

   The command NPoles is used to plot the quadrupole potential.  This command allows evaluation of the elctric potential at any point in space given any distribution of point charges.  The previous commands that were used to find potential, Monopole and Dipole, were actually special cases of the command NPoles.  NPoles is explained in greater detail below.
             
   The following is a plot of the electric potential about a quadrupole when the evaluation plane is taken to be a segment of the x = 0 plane.  Double click on the cell following the plot of the electric potential to see a diagrammatical representation of the quadrupole configuration.

charge = Pi Epsilon0;
location1 = {0,1,1};
location2 = {0,1,-1};
Plot3D[NPoles[{charge,-charge,charge,-charge},
   {location1,location2,-location1,-location2},{0,y,z}],
   {y,-4,4}, {z,-4,4},
   ClipFill -> None, Lighting -> True, PlotPoints -> 31,
   PlotRange -> {Automatic, Automatic, {-.2,.2}},
   BoxRatios -> {1,1,1}, AxesLabel -> {"y","z","V"},
   PlotLabel -> "Fig. 3"]

Diagrammatical representation of the quadrupole and the evaluation plane.

Show[ EvaluationPlane[0.], axes, qpoints, 
   Lighting -> False, AxesLabel -> {"x","y","z"}]

The positive charges are located at the upper-right and lower-left corners of the square formed by the four charges and are black.  The negative charges are gray.  The evaluation plane is the shaded portion of the x = 0 plane.  As noted earlier, the charges are represented by spheres with non-zero radius, even though we are considering point charges.

More information about the command NPoles.

The function NPoles calculates the potential at an observation point given a distribution of point charges.  The syntax is:

 NPoles[{list_of_charges},{list_of_locations},{x0,y0,z0}]
      
where list_of_charges is a comma separated list of the strength of the charges in the system, list_of_locations gives the location of the charges as a comma separated list of ordered triplets (each triplet is itself a list of comma separated numbers corresponding to the location of the charge in cartesian coordinates), and {x0,y0,z0} is the observation point.

The potential about a quadrupole has one additional feature that was not present in the potential about the monopole or the dipole.  If you look at the center of Fig. 3 (which corresponds to the origin), you see that the potential there is zero.  Unlike the line of zero potential we saw along the y axis for the dipole, if you were located exactly at the origin, you would stay there.  There is no net force trying to move you.  However, you would be in a highly unstable position – a slight movement in any direction other than along the line joining the two positive charges, and you would tend to fall down one of the holes.  The origin is known as a saddle point.  It is neither a minimun nor a maximum, but the directional derivatives there are zero.  It is as if the origin was an infinitesimally small flat spot between two mountains and two canyons.

Contour and Density Plots

Introduction

Contour and and density plots are another way in which potential fields are often displayed.  These plots are used primarily to give you a clearer picture of the equipotential lines.  This information is contained within the surface plots – they are merely the line at a constant height.  However, since the three-dimensional surface plot must be displayed in two dimensions, it is generally quite difficult to see just where the equipotential lines are.  Contour and density plots are generally evaluated over the same regions for which you would use a surface plot – that is to say a plane.  Contour plots trace several lines of constant potential (or constant height from the surface plots) in the plane being considered.
       
   Density plots assign a color or shade to each point in the evaluation plane based upon the potential at that point.  Density plots are good for giving a overall picture of the behavior of the field.  However,  there are some difficulties with density plots that can often limit their usefulness.   The problems arise when you want to consider some of the details of the behavior of the field.  It is difficult for the eye to judge differences in values by looking at change in shade or color.  Also, rates of change are hard to ascertain from changes in colors.  To get around this, density plots are usually accompanied by a key which tells you which shades or colors correspond to which values.  Nevertheless, density plots are rarely used to examine any of the fine details of a field.
              
   We will start with the dipole since the monople is just too uninteresting for these plots.  We will also examine the quadrupole using these plots.

Dipole

Contour Plot

Again, Dipole[x,y,z] will be used to describe the function being plotted.  However, now we are displaying the result using a different Mathematica command.  Instead of Plot3D, ContourPlot will be used.  The configuration of charge is the same as that used in Fig. 2.

charge = Pi Epsilon0;
location1 = {0,0,1};
location2 = {0,0,-1};
ContourPlot[Dipole[ charge, location1, location2, {0,y,z} ],
   {y,-4,4}, {z,-4,4},
   PlotPoints -> 45,
   ContourShading -> False,
   Contours -> 13,
   AxesLabel -> {"y","z"},
   PlotLabel -> "Fig. 4"]

For these contour lines, we have not bothered to label the potential to which each line corresponds.  Mathematica not providing any simple means of accomplishing this.  Usually a countour plot will provide information about the potential that each line represents.  In Fig. 4, the line which cuts the plot is approximately the zero potential line which is the y axis.  The contour lines above the y axis correspond to positive voltages, and those below correspond to negative voltages.  Notice that the lines are more closely spaced around the charges where the potential changes most rapidly.
       
   The information which can be read directly from the contour plot is the behavior of the equipotential lines.  If you were a test charge and you moved along one of these lines, you would neither gain nor loose energy.  The line you move on would be at a constant height (or depth) on the surface plots.

Density Plot

Density plots can be very useful tools when examing the ™big picture∫ of the behavior of a field, but these plots generally do not give information about detailed behavior.  Limitations of the device on which a density plot is displayed can also inhibit its usefulness.  
       
   Shown below is the density plot for the dipole of Figs. 2 and 4.

charge = Pi Epsilon0;
location1 = {0,0,1};
location2 = {0,0,-1};
DensityPlot[Dipole[ charge, location1, location2, {0,y,z} ], 
   {y,-4,4}, {z,-4,4},
   PlotPoints -> 45,
   ColorFunction -> Hue,
   Mesh -> False,
   AxesLabel -> {"y","z"},
   PlotLabel -> "Fig. 5"]

In Fig. 5, you can see some of the difficulties determining detail about a field from its density plot.  Even if you were provided with a key that described the relationship between shade and field strength, it might be rather difficult to say were the zero potential line was.  However,  you can easily see approximately where the charges are.

  Figure 5 was generated from the same dipole function that we have used for all the other plots of dipole potential.  The only thing that has changed is that the Mathematica function DensityPlot is now being used.  This is significant in that one underlying equation describes all the dipole plots you are seeing – despite the immense differences in the displayed form.   (Feel free to look at the last section of this NoteBook, where all the Mathematica definitions of the electromagnetic equations are stored.  The equations had to be translated into forms that Mathematica understands and they loose some of there elegence in the translation.  Therefore, you are probably better off just looking in an introductory electromagnetics text book if you need a refresher on the form of some of these equations.)

Combination Plots

Of course, it is possible to generate a combination contour-density plot as shown here.

charge = Pi Epsilon0;
location1 = {0,0,1};
location2 = {0,0,-1};
ContourPlot[Dipole[ charge, location1, location2, {0,y,z} ],
   {y,-4,4}, {z,-4,4},
   PlotPoints -> 45,
   ColorFunction -> Hue,
   Contours -> 13,
   AxesLabel -> {"y","z"}]

Quadrupole

Contour Plot

Figure 6 shows the contour plot for the quadrupole.  The system of charge is the same as used in Fig. 3.

charge = Pi Epsilon0;
location1 = {0,1,1};
location2 = {0,1,-1};
ContourPlot[NPoles[{charge,-charge,charge,-charge},
   {location1,location2,-location1,-location2},{0,y,z}],
   {y,-4,4}, {z,-4,4},
   PlotPoints -> 45,
   ContourLevels -> 11,
   AxesLabel -> {"y","z"},
   PlotLabel -> "Fig. 6"]

Density Plot

The density plot for a quadrupole is shown in Fig. 7.

charge = Pi Epsilon0;
location1 = {0,1,1};
location2 = {0,1,-1};
DensityPlot[NPoles[{charge,-charge,charge,-charge},
   {location1,location2,-location1,-location2},{0,y,z}],
   {y,-4,4}, {z,-4,4},
   PlotPoints -> 45,
   Mesh -> False,
   PlotLabel -> "Fig. 7"]

Parametric Plots (Spherical Coordinates)

Introduction

Parametric plots display the potential as a function of one or two parameters.  These parameters indicate the evaluation point.  The surface plots that we have already seen can be thought of as parametric plots.  The parameters were the y and z values of the evaluation point – the x value was fixed.  Points in space can also be described by spherical coordinates.  The location is be determined by a distance from the origin and two angles (see the NoteBook CoordinateSystems).  If the distance from the origin is constant, the set of points described by varying the two angles is a sphere.   We might be interested in the intensity of the fields over this sphere (note intensity is a magnitude – it implies absolute value).  The intensity of the field as a function of the two angles can be displayed as a distance from the origin.  In cartesian coordinates we fixed the value of x and varied y and z.  The field at any point was shown as a height above a plane.  In spherical coordinates, the radius from the origin to the evaulation point is fixed and the two angles are varied.  We display the intensity at the evaluation point as a distance from the origin (note that this distance is unrelated to the distance from the origin to the evaluation point).
      
   In this section we will consider the spherical parametric plots of the potential field about a dipole, quadrupole and octupole.  In each of these cases the distance from the origin to the evaluation point is held fixed at 4 m.  The monopole will not be shown because it is not very interesting.  If we were to measure the potential a fixed distance away from a single charge (a monopole), what would we find?  The potential would not vary as the angles specifying the evaulation point changed.  The result would just be a sphere.
      
      Unfortunately the following plots take some time to generate.  Unless you are patient, you probably don't want to tweak a parameter and run these yourself.

Dipole

The parametric plot of the electric potential for a dipole is shown below.  Remeber that the evaluation point is a fixed distance from the origin and the evaluation point is changed by changing the angles θ and ϕ.  If a line were drawn from the origin to the evaluation point, it would form an angle θ with the positive z axis.  If this line were projected into the xy plane, it would form an angle ϕ with the positive x axis. (See the NoteBook CoordinateSystems for a further description of spherical coordinates.)
      
   In the plot below, the origin is at the center of the cube.  The labels on the side of the cube show you the orientation of the axes but they are not the axes themselves – the ™real∫ axes intersect at the center of the cube.  The distance from the origin to the surface of the object in a particular direction indicates how strong the field is in that same direction, but keep in mind that we are evaluating the strength of this field 4m from the origin.  You should be able to see that the magnitude of the potential is greatest on the z axis (the z axis contains the points on the surface most distant from the origin).  In the xy plane the potential goes to zero and hence the distance from the origin is zero in that plane.  The reason the magnitude of the potential is greatest on the z axis is because those points are as close to one charge as possible while being as distant from the opposite charge as possible.  
      
   The reason the absolute value of the potential is used is because it is confusing to plot negative distances from the origin in spherical coordinates.  If you are curious about why the absolute value is needed, feel free to delete the Abs[ from the following command (you will also have to remove the ] associate with this command so that the ™]],∫ at the end of the seventh line reads ™],∫).  In the two minutes it takes to get the output, think about what the output should look like.

charge = Pi Epsilon0;
location1 = {0,0,1};
location2 = {0,0,-1};
radius = 4.;
SphericalPlot3D[
  Abs[Dipole[ charge, location1, location2,
     {radius*xloc[t,p], radius*yloc[t,p], radius*zloc[t,p]} ]],
  {t,0,Pi,Pi/25},{p,0,2 Pi,Pi/15},
  PlotRange -> {{-.035,.035},{-.035,.035},{-.035,.035}},
  AxesLabel -> {"x","y","z"},
  PlotLabel -> "Fig. 8"]

Mathematica command used to generate the  parametric plot.

SphericalPlot3D is a Mathematica command that plots functions that are of the form r(θ,ϕ) where r is the distance from the origin, and θ and ϕ specify a direction from the origin.  The syntax of the command is 
      
SphericalPlot3D[
   r[theta,phi],{theta,start,stop,step},{phi,start,stop,step}
   ]
      
      The command Abs takes the absolute value of its argument.  This is done to insure the distance from the origin is non-negative.

Even though θ and ϕ are spherical coordinates that determine where the potential will be evaluated,  we are not restricted to plotting the output only in spherical coordinates.  We can create a surface plot where the two independent variables are θ and ϕ, and the dependent variable is the potential.  Since θ varies between 0 and π and ϕ varies between 0 and 2 π, a segment of a plane with these dimensions is all that is needed.  The plot shown below is a surface plot that corresponds to the parametric plot of the dipole potential shown in Fig. 8.  The distance from the origin is still fixed at 4 m.  Notice that it is harder to get a feel for the distribution of potential over the sphere in which we are interested.

charge = Pi Epsilon0;
location1 = {0,0,1};
location2 = {0,0,-1};
radius = 4.;
Plot3D[
  Abs[Dipole[ charge, location1, location2,
     {radius*xloc[t,p], radius*yloc[t,p], radius*zloc[t,p]} ]],
  {t,0,Pi},{p,0,2 Pi},
  AxesLabel -> {"theta","phi","Potential"},
  PlotLabel -> "Fig. 9"]

This surface plot differs from those shown previously in that the magnitude of the potential is plotted.  You can see that the potential is independent of phi and it obtains its maximum value at theta equal 0 or π (these points correspond to points on the postive and negative z axis.

Quadrupole

The parametric plot for a quadrupole is shown in Fig. 10.  Owing to the four charges in the quadrupole, there are four ™lobes∫ in the distribution of potential.  The maximum in these lobes occurs when the evaluation point is as close as possible to a charge while being as distant as possible from the oposite charges.

charge = Pi Epsilon0;
location1 = {0,1,1};
location2 = {0,1,-1};
radius = 4.;
SphericalPlot3D[
  Abs[NPoles[{charge,-charge,charge,-charge},
     {location1,location2,-location1,-location2},
     {radius*xloc[t,p], radius*yloc[t,p], radius*zloc[t,p]}]],
  {t,0,Pi,Pi/30},{p,0,2Pi,Pi/30},
  PlotRange -> {{-.022,.022},{-.022,.022},{-.022,.022}},
  AxesLabel -> {"x","y","z"},
  PlotLabel -> "Fig. 10"]

Figure 11 is a surface plot presentation of the information contained in Fig. 10.  Notice that it is easier to ascertain the physical distribution of potential from Fig. 10 than it is in Fig. 11.

charge = Pi Epsilon0;
location1 = {0,1,1};
location2 = {0,1,-1};
radius = 4.;
Plot3D[
  Abs[NPoles[{charge,-charge,charge,-charge},
     {location1,location2,-location1,-location2},
     {radius*xloc[t,p], radius*yloc[t,p], radius*zloc[t,p]}]],
  {t,0,Pi},{p,0,2 Pi},
  PlotPoints -> 25,
  AxesLabel -> {"theta","phi","Potential"},
  PlotLabel -> "Fig. 11"]

Octupole

An octupole consists of eight charges that are located at the corners of a cube.  Like charges appear on the diagonals of any face and the other diagonal will have charges of the opposite sign.  A diagrammatical view of the distribution of charges about the origin is available in a cell below.  The spherical parametric plot for a octupole is shown in Fig. 12.

charge = Pi Epsilon0;
radius = 4.;
SphericalPlot3D[
  Abs[
     NPoles[
        { charge, -charge,  charge, -charge,
         -charge,  charge, -charge,  charge},
        {{1,1,1},{1,-1,1},{-1,-1,1},{-1,1,1},
         {1,1,-1},{1,-1,-1},{-1,-1,-1},{-1,1,-1}},
     {radius*xloc[t,p], radius*yloc[t,p], radius*zloc[t,p]}]],
  {t,0,Pi,Pi/40},{p,0,2 Pi,Pi/40},
  PlotRange -> {{-.022,.022},{-.022,.022},{-.022,.022}},
  AxesLabel -> {"x","y","z"},
  PlotLabel -> "Fig. 12"]

Diagrammatical representation of the octupole.

Show[ opoints, axes, AxesLabel -> {"x","y","z"}]

The positive charges are black and the negative charges are gray.  There is a positive charge at (1,1,1).  The points that are one segment away from this charge are negative – the charges at (1,±1,1), (±1,1,1), and (1,1,±1) are negative.  You are viewing this diagramatical representation at the same angle used in the display of the parametric plot.

Figure 13 shows the surface plot representation of Fig. 12.  However, in this case the absolute value is not taken so you can see which angular directions are associated with which types of charge.  The positive ™bumps∫ correspond to the angular direction closest to a positive charge whereas the negative ™bumps∫ indicate the direction closest to a negative charge.

charge = Pi Epsilon0;
radius = 4.;
Plot3D[
  NPoles[
     { charge, -charge,  charge, -charge,
      -charge,  charge, -charge,  charge},
     {{1,1,1},{1,-1,1},{-1,-1,1},{-1,1,1},
      {1,1,-1},{1,-1,-1},{-1,-1,-1},{-1,1,-1}},
  {radius*xloc[t,p], radius*yloc[t,p], radius*zloc[t,p]}],
  {t,0,Pi},{p,0,2 Pi},
  PlotPoints -> 30,
  AxesLabel -> {"theta","phi","Potential"},
  PlotLabel -> "Fig. 13"]

Initialization Cells

The initialization cells do the following:
	* Define values for constants and give default values to some variable used in the
	   equations.
	* Define EM functions used in this section.
	* Define some graphics primitives used in the illustrative graphics.

The ElectricPotential package.

This is a Mathematica package.  It defines commands in a manner that allows the user to receive information about the command (via ?Command), and it keeps local variables hidden.

BeginPackage["ElectricPotential`"]

Dipole::usage =  "Dipole[q,p01,p02,p1] give the potential
at the point p1 = {x1,y1,z1} due to a charge of q at point p01 =
{x01,y01,z01} and a charge of -q at the point p02 = {x01,y02,z02}."

Distance::usage = "Distance[p0,p1] gives the distance between the
point p1 = {x0,y0,z0} and p2 = {x1,y1,z1}."

Epsilon0::usage = "Epsilon0 is the permittivity of freespace =
8.5837 10^-12."

Potential::usage = "Potential[q,r] gives the potential a distance r
away from a charge of q (in Coulombs).  The permittivity of free
space is assumed."

Monopole::usage =  "Monopole[q,p0,p1] give the potential
at the point p1 = {x1,y1,z1} due to a charge of q at point p0 =
{x0,y0,z0}."

NPoles::usage =  "NPoles[{q1,q2,...,qn},{p01,p02,...,p0n},{x1,y1,z1}]
give the potential at the point p1 = {x1,y1,z1} due to discrete charges
q1, q2, ... qn located at points p01, p02, ... p0n, respectively."


Begin["`Private`"]

Epsilon0 = 8.8537 10^-12

Potential[q_Real,r_Real] := 
   N[q/(4 Pi Epsilon0 r)]

Distance[{x_,y_,z_},{x0_,y0_,z0_}] := 
    Sqrt[(x-x0)^2 + (y-y0)^2 + (z-z0)^2 ]

(* Monopole calls NPoles with the proper syntax for a single
   charge.  Monopole provides a slightly more intuitive syntax
   than NPoles. *)
Monopole[q_Real, p0_List, p1_List] := 
    NPoles[{q},{p0},p1]

(* Dipole calls NPoles with the proper syntax for two charges of
   equal magnitude and opposite signs.  Dipole provides a slightly
   more intuitive syntax than NPoles. *)
Dipole[q_Real,p01_List, p02_List, p1_List] := 
    NPoles[{q,-q},{p01,p02},p1]
     	
NPoles[qList_List, p0List_List, p1_List] :=
   Block[{i},
      N[Apply[Plus,
         Table[Potential[qList[[i]], Distance[ p0List[[i]], p1 ]],
            {i,1,Length[qList]}]
         ]
      ]
   ]
 	
End[(* "`Private`" *)]

EndPackage[(* "ElectricPotential`" *)]

Some definitions needed for the illustrative graphics.

(* Evaluation plane used in the dipole surface plot. *)
EvaluationPlane[d_] :=
   Graphics3D[ 
      {GrayLevel[0.666667],
      Polygon[{{d,4,4},{d,4,-4},{d,-4,-4},{d,-4,4}}]} ];

(* Axes used in several of the plots. *)
axes = Graphics3D[ {Thickness[.006], GrayLevel[0],
   (*z axis*)  Line[{{0,0,-4},{0,0,4}}],
   (*y axis*)  Line[{{0,-4,0},{0,4,0}}],
   (*x axis*)  Line[{{-4,0,0},{4,0,0}}]}];

(* Large point drawn at the origin.  Used in the monopole 
   discussion. *)
mpoint  = Graphics3D[ {PointSize[0.03], Point[{0,0,0}]} ];

(* Two large points drawn on the z axis.  Used in the dipole
   discussion. *)
dpoints = Graphics3D[
   {PointSize[0.03],
   GrayLevel[0],  Point[{0,0,1}],
   GrayLevel[.3333], Point[{0,0,-1}]} ];

(* Four large points drawn in the yz plane.  Used in the 
   quadrupole discussion. *)
qpoints = Graphics3D[ 
   {PointSize[0.03],
   GrayLevel[0],  Point[{0,1,1}],
   GrayLevel[.3333], Point[{0,1,-1}],
   GrayLevel[.3333], Point[{0,-1,1}],
   GrayLevel[0],  Point[{0,-1,-1}]} ];

(* Eight large points drawn at the corner of a cube
   centered at the origin.  Used in the octopole discussion. *)
opoints = Graphics3D[ 
   {PointSize[0.03],
   GrayLevel[0],  Point[{1,1,1}],
   GrayLevel[.3333], Point[{1,-1,1}],
   GrayLevel[0],  Point[{-1,-1,1}],
   GrayLevel[.3333], Point[{-1,1,1}],
   GrayLevel[.3333], Point[{1,1,-1}],
   GrayLevel[0],  Point[{1,-1,-1}],
   GrayLevel[.3333], Point[{-1,-1,-1}],
   GrayLevel[0],  Point[{-1,1,-1}]} ];

Miscellaneous commands

Needs["Graphics`ParametricPlot3D`"]

(* Define functions to project the unit vectors in the
   spherical coordinate system to the cartesian system *)
xloc[theta_,phi_] := Sin[theta]Cos[phi];
yloc[theta_,phi_] := Sin[theta]Sin[phi];
zloc[theta_,phi_] := Cos[theta];