%Das zu lösende lin. Glsys. r:=b-Ax=0 hat hier die Besonderheit, daß
%nur r selbst einfach ausgerechnet werden kann. x ist hier formal eine Matrix,
%die auch die Randwerte enthält. Das Ergebnis r hat dieselbe Größe wie x.  
%r(j,k)=u(j,k)-0.25(u(j+1,k)+u(j-1,k)+u(j,k+1)+u(j,k-1)), (j,k) innerer Pkt,
%die Randpunkte bleiben bei r immer Null. Für x=0 ist offensichtlich r=b.  
function r=residuum(x);
  [z s]=size(x);
  r=zeros(size(x));
  hhh=x(1:z-2,2:s-1)+x(3:z,2:s-1)+x(2:z-1,1:s-2)+x(2:z-1,3:s);
  r(2:z-1,2:s-1)=0.25*hhh-x(2:z-1,2:s-1);
return;
