% Matlab code makemgdemo.m
% For "Applied Numerical Linear Algebra",  Question 6.16
% Written by Ivar Gustafsson, Jan 2009
%
% Initialize b, x, jac1, jac2, k, itmax and tolres to demonstrate testfmg
%
hold off, clf
% Size of problem is n=k*k+1
k=6;
n=2^k+1;
[x,y]=meshgrid(1:n,1:n);
% Continuous true solution
xtrue=sin((x-1)*3*pi/n).*sin((y-1)*pi/n);
% Right hand side computed from true solution
b=zeros(n,n);
b(2:n-1,2:n-1) = 4*xtrue(2:n-1,2:n-1) ...
                    -xtrue(1:n-2,2:n-1) ...
                    -xtrue(3:n  ,2:n-1) ...
                    -xtrue(2:n-1,1:n-2) ...
                    -xtrue(2:n-1,3:n  );
% Right hand side like a tent with two poles
%b=zeros(n,n);
%b(2^(k-1)+1,2^(k-2))=1;
%b(2^(k-1)+1,2^k+1-2^(k-2))=1;
% very smooth right hand side
%b=sin((x-1)*pi/n).*sin((y-1)*pi/n);
subplot(1,2,1); mesh(xtrue), title('true solution')
subplot(1,2,2); mesh(b), title('right hand side')
x=zeros(n,n);
jac1=10;
jac2=10;
delta=1;
itmax=80;
tolres=1e-6;