function ShowPoisson()
% Display the approximate solution to u_{xx} + u_{yy} = -peaks(x,y)
% on the square with vertices (-3,-3), (3,-3), (3,3), (-3,3).
% u(x,y) = 0 on the boundary.

close all
a = -3; b = 3; c = -3; d = 3;
n = 500;
x = linspace(a,b,n+1);
y = linspace(d,c,n+1);

% Display the peaks function . Think of it as a heat source function.
P = zeros(n+1,n+1);
for i=1:n+1
for j=1:n+1
P(i,j) = peaks(x(j),y(i));
end
end
pcolor(x,y,P)
shading('interp')
axis equal tight
title('The Source Function peaks(x,y)', 'fontsize',14)

% Solve the Poisson equation and display...
U = FastPoisson(a,b,n,c,d,n,@peaks);
figure
pcolor(x,y,U)
shading('interp')
axis equal tight

title('The Solution u(x,y)', 'fontsize',14)
shg