%Exercise 2: Find the solution of a 2-D Laplace Equation, on unit square, with
%boundary conditions equal to 1 on left and right boundaries and 0 on 
%upper and lower boundaries. Use h=0.2 in both dimensions.
%Use any system solver available. Plot the solution.

%"exact" solution (taken from Verena Horak's paper; she had this one assigned)
U =[1.0000,         0,         0,         0,         0,    1.0000;
    1.0000,    0.5000,    0.3333,    0.3333,    0.5000,    1.0000;
    1.0000,    0.6667,    0.5000,    0.5000,    0.6667,    1.0000;
    1.0000,    0.6667,    0.5000,    0.5000,    0.6667,    1.0000;
    1.0000,    0.5000,    0.3333,    0.3333,    0.5000,    1.0000;
    1.0000,         0,         0,         0,         0,    1.0000];

%Exercise 4: Implement 20 steps of Jacobi iteration for the system matrix and boundary 
%conditions from exercise 2. Estimate the maximal error in the iterative solution.
%Plot the solution.

D=[1 0 0 0 0 1;1 0 0 0 0 1;1 0 0 0 0 1;1 0 0 0 0 1;1 0 0 0 0 1;1 0 0 0 0 1];

Dold = D;
for k = 1:20 %we do 20 iteration steps
    for i=2:5 %we also skip the values for the boundary conditions at i=1,6
        for j = 2:5 %skip the boundary values again
            %now follows the literal implementation of system of equations
            D(i,j) = (Dold(i-1,j)+Dold(i,j-1)+Dold(i+1,j)+Dold(i,j+1))/4;
        end
    end
    Dold=D;
end
D
%note that the axes labes of the plane of the definition domain are the indices of the
%matrix entries and NOT the real coordinates (i.e. the arguments of u)
surfc(D);

%to see the error, we compare our result with the "exact" solution from exercise 2:
%the largest error is 1/100 at the innermost points of the definition domain
U-D