%
% This is an example to show how a vector based on a 2D mesh
% is mapped into a 1D array, sent to a solver, and mapped back.
%
% It will also show you how to use my fish2d.m Poisson solver.
%
% We'll use a 33 x 33 mesh with h=1/32 and 31^2 uknowns
%
% WATCH OUT! fish2d does not make the boundary conditions part of the
% solution or the right hand side. This is just like the direct solvers
% we used in 1D. This means that you need to figure out how to extract what
% fish2d needs from your vector.
%
m=31; 
h=1/(m+1);
h2=(m+1)*(m+1);
h1=(m+1);
tol=1/h2;
%
% set up the equation
%
x=1:m;
x=x'*h;
%
% Form the solution: sol(x,y)= 10 x y (1-x) (1-y) exp(x^4.5)
% Notice how I map the x variable to the first coordinate and the y
% variable to the second.
%
zsoly=x.*(1-x);
zsolx=zsoly.*exp((1-x).^4.5);
solt=10*zsolx*zsoly';
%
% Now it's time to map the whole works on a 1d array with a simple matlab
% command.
%
sol=solt(:);
%
% fix the right side so that the solution is known. lampmf.m is the 
% discrete laplacian
%
b = lapmf(sol); 
%
% send b to the solver and the difference from sol should be near machine
% roundoff
% 
v=fish2d(b);
norm(v-sol,inf)