In Lecture #5, we saw how Laplace's Equation gives rise to the phenomenon
of electrostatic shielding by a conducting enclosure. We also asserted that
although this result is rigorously true for * completely enclosed *
cavities, the shielding is rather effective even for partially enclosed
regions.

As evidence of this we looked at a numerical solution to Laplace's Equation for the case of a grounded metal cylinder with a variable-sized opening in its side. You can view these solutions (along with isolated charges for comparison) by following these links:

This image shows two panels, with potential plotted as a color contour map. The top panel shows a charge alone, and thus has spherical symmetry. The lower panel shows the field with the addition of a grounded cylinder with an opening in its side. In addition to noting the shielding effect, you'll see that the influence of the conductor at V=0 is to cause the equipotentials to "shrink" about the source charge.

This pair of panels is similar, except that the cylinder now has a gaping wide 180 degree opening in it. Although the effect is far from perfect, there is a surprisingly good degree of shielding within the shell concavity.

Each year, I get an inquiry or two from students who want to know about the technique used to generate these solutions. The following is a very brief explanation, describing the numerical solution method, followed by a listing of the C program that was used:

Laplace's equation is a second order partial differential equation, and in order to solve it, you must find the unique function who derivatives satisfy (del squared) V = 0, and simultaneously satisfies the required boundary conditions. While a clever enough person might find an analytic solution in a situation with a high degree of symmetry, for an irregularly shaped configuration it is probably hopeless.

In such a case, we turn to finding an approximate solution by numerical
means. In this case, we use a method called * finite differences *. The main
idea is fairly simple. Instead of trying to find V at * every point * in
our domain, we'll seek to find V at a finite number of points that make up a
lattice or grid, which covers our region of interest. Then to calculate
"derivatives" of our function V, we'll approximate them by linear slopes
between neighboring grid points. For example, if we know V at two neighboring
points, V(i) and V(i+1), then the derivative there is approximately,
( V(i+1) - V(i) ) / delta_x . Similarly, we can use other neighbors to
approximate the second and third derivatives, and so on, as needed.

Now in the case of the two-dimensional Laplace equation, if we write it in finite difference form, it comes out something like:

4 V(i,j) - V(i-1,j) - V(i+1,j) - V(i,j-1) - V(i,j+1) = 0

A brief examination of this equation shows that for any point (i,j), the value
of V there * must * be the average of the values of its four nearest
neighbors. So we can find a solution by iteratively passing through the
points of the grid and forcing them to satisfy this condition. We begin by
fixing the values of any boundary conditions, and prohibiting them from
changing. Then we sweep through the grid one point at a time, replacing V at
each point by the average of its neighbors. By repeating this process, the
solution "relaxes" toward the correct solution and eventually stops changing,
at which point we are done.

Below is a C source code listing for the program used to accomplish this in the present example. The program is very bare bones, with no provisions for error checking or checking for convergence, but it got the job done. I hope you find it helpful.

```
/* laplace2.c - solves Laplace's eqn. in a grounded box
near a cylindrical shell with a variable sized opening
*/
#include <stdio.h>
#include <math.h>
main(argc,argv)
int argc ;
char *argv[] ;
{
FILE *infile , *outfile ;
float v[400][400] ;
int active[400][400] ;
int i , j , pass ;
float theta , fi , fj ;
infile = fopen(argv[1],"r") ;
outfile = fopen(argv[2],"w") ;
fscanf(infile,"%f",&theta) ;
printf("Shell half-angle: %f\n",theta) ;
/* initialize */
for(j=0;j<400;j++)
{
for(i=0;i<400;i++)
{
fi = (float) i ;
fj = (float) j ;
if(i==0 || i==399 || j==0 || j==399)
{
active[i][j] = 0 ;
v[i][j] = 0.0 ;
}
else if( ((fi-180.)*(fi-180.)+(fj-200.)*(fj-200.) < 675.0) &&
((fi-180.)*(fi-180.)+(fj-200.)*(fj-200.) > 575.0) &&
((180.-fi)/sqrt((fi-180.)*(fi-180.)+(fj-200.)*(fj-200.))) >
cos(theta*3.14159/180.) )
{
active[i][j] = 0 ;
v[i][j] = 0.0 ;
}
else if(i==220 && j==200)
{
active[i][j] = 0 ;
v[i][j] = 100.0 ;
}
else
{
active[i][j] = 1 ;
v[i][j] = 0.0 ;
}
}
}
/* do 10000 Gauss-Seidel updates */
printf("Starting updates...\n") ;
for(pass=0;pass<10000;pass++)
{
if((pass % 1000) == 0) printf("pass = %d\n",pass) ;
for(j=0;j<400;j++)
{
for(i=0;i<400;i++)
{
if(active[i][j] == 0) continue ;
else v[i][j] = 0.25 * (v[i][j-1]+v[i][j+1]+v[i-1][j]+v[i+1][j]) ;
}
}
}
/* write the output */
printf("Printing output...\n") ;
fprintf(outfile,"100 100\n") ;
fprintf(outfile,"5.0 0.0\n") ;
for(j=150;j<249;j++)
fprintf(outfile,"%f\n",(float)j) ;
fprintf(outfile,"\n") ;
for(i=150;i<249;i++)
fprintf(outfile,"%f\n",(float)i) ;
fprintf(outfile,"\n") ;
for(j=150;j<249;j++)
{
for(i=150;i<249;i++)
{
fprintf(outfile,"%f\n",v[i][j]) ;
}
}
fclose(infile) ;
fclose(outfile) ;
}
```

Return to Physics 51 Home Page

Last updated 9/13/95.