Lab 13: Solving Poisson`s equation using the Relaxation Method

Here a is the spacing between the regions on the
grid. You will solve this equation for a two
dimensional box with sides kept at the constant
potential 0V. We will assume a box 1.0 m on
each side, divided into N equally spaced
intervals.
Lab 13: Solving Poisson’s equation
using the Relaxation Method
Objectives
By the end of this lab you should:
• Solve Poisson’s equation using the
relaxation method.
The charge density ρ (i, j ) with consist of two
square regions located b=0.2 m from the sides,
and c=0.2 m on each side, as shown. One region
will contain a positive charge density of 1 C/m2
and one region -1 C/m2.
Prelab question: Before lab, create a flowchart
for the program you will write for lab Activity1.
Your flowchart will be checked by the instructor
at the beginning of lab.
Programming environment: IDLE. To
implement this, do the following:
• Go to the start menu icon and click it.
• Select All Programs
• Select Anaconda 64-bit
• Select Anaconda Command Prompt
• At the prompt, type idle and then hit
return.
• This will open an IDLE window, and
give you access to a different version of
python, which has access to pylab
plotting routines.
• We will begin learning how to use these.
Activity 1: Obtaining the Potential φ
• Here are some commands you will need:
Background
The Poisson equation (PE) describes how the
electric potential φ varies in a region of space
containing a volume charge density ρ:
∇ 2φ = −
ρ
ε0
We have seen how the relaxation method may be
used to solve for φ iteratively until the
difference between successive approximations is
less than the tolerance. If φ ' is the new array
and φ is the old array, the recursion relation for
the arrays is:
φ ' [i, j] = 14 * (φ[i, j - 1] + φi[i, j + 1]
+ φ [i - 1, j] + φi[i + 1, j]) +
a2
ε0
ρ (i , j )
1
import numpy as np
import matplotlib.pyplot as
plt
Here we use the
import as command.
def rho(i,j):
your code here
We choose to define
ρ (i, j ) as a function
that checks if i and
j are within ranges
that place you in the
charge distribution.
If so, it should return
+1 or -1 depending
where you are, and 0
elsewhere. Here the
integer division // is
useful.
NECESSARY!
This reveals the
shape
of
equipotential lines.
N is the number of
intervals along the
side of the box. L
is the length of
box’s side. You
define the arrays.
Careful:
the
number of points is
different from the
number
of
intervals.
dxy
is
the
increment in x and
y (this will be used
later for plotting,
not
for
your
computation).
N=20
L=?
phi=?
phiPrime=?
tol=1e-4
delta=1.
dxy=?
• Run the program and generate a plot of
the potential. For these two regions of
opposite charge, comment on whether the
program gives you what you expect.
• You may wish to visit the PhET
simulation “Charges and Fields” and play
with it.
• Now make the regions have the same
charge density. Again, do you get what
you expect?
• Explore what happens as you increase N.
Activity 2: Circular Charge Density Region
• It is easy to create circular or elliptical
shaped region of charge density. Recall
that the formula for an circular region is:
(i − i0 )2 + ( j − j0 )2 ≤ n 2
Here i and j are integers representing
the row and column, and i0 and j0 and
constants indicating the center of the
charged region. The integer n is the
radius measured in units of grid size.
Elliptical regions are similar.
• Modify your code to create a circular
region centered L/2 from each side.
Make the radius L/8 and the charge
density +1 C/m2.
• Generate a potential map and comment
on the results.
• Create a loop that uses the relaxation
method to find the array φ . Define
δ = max φ − φ '
which is the largest element of the
absolute value of the difference array.
Use the np.max() function from the
numpy package, not the standard python
max(). The relaxation method should
run as long as
δ > tol
• The boundaries stay at the same
potential. Since the boundary conditions
apply to the rows i=0, i=N, and
columns j=0, j=N, make a conditional
that sets φ ' = φ under these conditions.
At the end of the
loop swap the
arrays before the
next iteration.
plt.imshow(phi)
These commands
plt.set_cmap(‘hsv’) create a contour
plt.colorbar()
plot: φ within a
plt.show()
certain range are
plotted using the
colormap
hsv.
CHOOSE
A
BETTER
COLORMAP IF
phi,phiPrime=phiPrime,phi
2
Activity 3: Optional: A better contour plot
method
You can try to implement the same
countourf() routine as last lab for producing
a better contour map than imshow() used in
Activities 1 and 2. Pay attention to choosing a
good colormap that most clearly illustrates your
points.
You may otherwise try a different boundary
conditions and observe how the the potential
maps change for the new conditions.
What to turn in:
• Rubric cover sheet
• Prelab flowchart
• Hardcopy of program
• Submit the python programs in dropbox
on D2L. Make sure:
o The filename is
Lab13YourName.py
o You save a copy on a flash drive
or email the program to yourself.
o Inside your program, the top line
is: #Lab13 Your name
• Your work will be evaluated and returned
by Monday. Revisions will be due in
dropbox by Thursday noon.
3