tic
format long

%discretization and initialization%%%%%%%%%%
N=101; %number of nodes in x- direction
       %%number of nodes in y- direction
h=1/(N-1);
ur=0.3;
T=zeros(N,N);
T_new=T;
x=0:h:(N-1)*h;y=0:h:(N-1)*h;
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
k=0;
nr_it=0;    
stop=0;
while (stop~=1)
    nr_it=nr_it+1;
    
     err_u=0;err_T=0;  
    for i=2:N-1
        for j=2:N-1
            T_new(i,j)=0.25*(T(i+1,j)+T(i-1,j)+T(i,j+1)+T(i,j-1)+h*h);
            if abs(T(i,j)-T_new(i,j))>err_T 
                err_T=abs(T(i,j)-T_new(i,j));
            end
        
            T(i,j)=T_new(i,j);
        end
    end;
     if err_T<1e-6
            stop=1;
        end
        
    T=T_new; 
    if mod(nr_it,250)==0
    fprintf('nr_it=%g err=%g\n', nr_it,err_T);
%     k=k+1;
%     [C,H]=contourf(x,y,T',20);
%     axis off 
%     %axis([0,1,0,1])
%     m(k)=getframe;
    end
end


[C,H]=contour(x,y,T',20);
axis equal
axis([0,1,0,1])
 
T_max=max(max(T))
nr_it
%movie(m)