% MATH-104B SUMMER 2011 
% Hw4  Problem 5
% Jacobi method
%
clear
results=[];
for k=0:5,
    
    n=10*2^k;
    h=1./(n+1);
    %
    % Define Grid
    %
    x=zeros(n,1);
    for i=1:n,
        x(i) = i*h;
    end
    %
    % Initialization
    %
    sol = sin(pi*x);
    e=ones(n,1)/h^2;
    a=spdiags([-e 2*e -e], -1:1, n, n);
    f = a *sol;

    %
    % Iteration starts here
    %
    tol = 1.e-4;
    
    u0 = zeros(n,1);
    u1 = zeros(n,1);
    error=norm(u0-sol)/norm(sol);
    
    iter = 0;
    while error > tol
        u1(1) = (h^2*f(1)+u0(2))/2.;
        for i=2:n-1
            u1(i) = (h^2*f(i)+u0(i+1)+u0(i-1))/2.;
        end
        u1(n) = (h^2*f(n)+u0(n-1))/2.;
        
        u0 = u1;
        
        error=norm(u0-sol)/norm(sol);
        iter = iter + 1;
     %   fprintf('Iteration %d, Error: %16.10g\n', iter, norm(u1-sol));
    end
    
    results=[results; n error iter];
    
    fprintf('Iterations for n=%d: %d\n', n, iter);
end

loglog(results(:,1),results(:,3),'o-')
title('Log-Log plot of the number of iterations')
xlabel('n');
ylabel('# of Iterations')
    
        
    
    