for n=1:N-1
L(k(n,0,N),k(n,0,N))=1;
B(k(n,0,N))=sin(pi*X(n+1));
L(k(n,N,N),k(n,N,N))=1;
B(k(n,N,N))=0;
end;
Finally we define the parts of the matrix
L
and vector
B
that correspond to the equations
for the interior points.
for i=1:N-1
for j=1:N-1
L(k(i,j,N),k(i,j,N))=-4;
L(k(i,j,N),k(i+1,j,N))=1;
L(k(i,j,N),k(i-1,j,N))=1;
L(k(i,j,N),k(i,j+1,N))=1;
L(k(i,j,N),k(i,j-1,N))=1;
end
end
Now we can solve the equation for
F
and plot the result. To do this we have to put the
F
values in a two dimensional grid
FF
and use the
mesh
command to do the 3-d plot.
If
X
and
Y
are vectors of length n and
Z
is an nxn matrix then
mesh(X,Y,Z)
plots the points
[X(j),Y(i),Z(i,j)]
.
F=L\B;
FF=zeros(N+1,N+1);
for i=0:N
for j=0:N
FF(j+1,i+1)=F(k(i,j,N));
end
end
mesh(X,X,FF);
9

We can print out the resulting graph using
print laplace1.pdf
(or
print laplace1.jpg
or
print laplace1.eps
). This will produce a pdf file (or jpg or eps file) containing the graph.
Here is the result:
0
0.2
0.4
0.6
0.8
1
0
0.2
0.4
0.6
0.8
1
0
0.2
0.4
0.6
0.8
1
Run the file
laplaceeqn.m
to produce this picture. Then modify the code to solve Laplace’s
equation with boundary conditions:
f
(
x,
0) = sin(
π
x
)
0
≤
x
≤
1
f
(0
, y
) = 0
0
≤
y
≤
1
f
(
x,
1) = sin(
π
x
)
0
≤
x
≤
1
f
(1
, y
) = 0
0
≤
y
≤
1
Say what code you modified, and hand in the resulting picture. Finally modify the code
to solve Laplace’s equation with boundary conditions:
f
(
x,
0) = sin(
π
x
)
0
≤
x
≤
1
f
(0
, y
) = 0
0
≤
y
≤
1
f
y
(
x,
1) = 0
0
≤
x
≤
1
f
(1
, y
) = 0
0
≤
y
≤
1
The third boundary condition is called a Neumann boundary condition, and corresponds
to detaching the rubber membrane from the wire along the top boundary.
Again, say
what code you modified, and hand in the resulting picture.
10

The first modification is to change the code defining the top and bottom boundary conditions to
for n=1:N-1
L(k(n,0,N),k(n,0,N))=1;
B(k(n,0,N))=sin(pi*X(n+1));
L(k(n,N,N),k(n,N,N))=1;
B(k(n,N,N))=sin(pi*X(n+1));
end;
The resulting plot looks like
0
0.2
0.4
0.6
0.8
1
0
0.2
0.4
0.6
0.8
1
0
0.2
0.4
0.6
0.8
1
The second modification is to change the code defining the top and bottom boundary conditions to
for n=1:N-1
L(k(n,0,N),k(n,0,N))=1;
B(k(n,0,N))=sin(pi*X(n+1));
L(k(n,N,N),k(n,N,N))=-1;
L(k(n,N,N),k(n,N-1,N))=1;
B(k(n,N,N))=0;
end;
The resulting plot looks like
11

0
0.2
0.4
0.6
0.8
1
0
0.2
0.4
0.6
0.8
1
0
0.2
0.4
0.6
0.8
1
12