function [x,k,err]=PSOR(B,b,g,x0,eta,Nmax)
// Resolution de min(Bx-b,x-g)=0, algo PSOR (w=1), iteratif a partir de x0.
// En sortie: k=nbre iterations effectuees, err= norme infinie de |x^k - x^{k-1}|,
//   avec arret si err<=seuil
x=x0;
k=0; 
err= 1+eta;
n=size(x,1)
while (k< Nmax & err>eta)
  k=k+1;  
  xold=x;
  for i=1:n
    // ADAPTER ICI:
    //  x(i) <- ( b(i) - somme_{j different de i}  B(i,j) * x(j) ) / B(i,i);
    //  x(i) <- max(x(i),g(i))
    j=[1:i-1,i+1:n];
    x(i)=(b(i)-B(i,j)*x(j))/B(i,i);
    x(i)=max(x(i),g(i));
  end
  err=norm(x-xold,'inf');
  //if modulo(k,10)==0; printf('(PSOR): k=%4i, erreur |x-xold|=%16.10f\n',k,err); end 
  
end 
if k==Nmax; printf('attention augmenter Nmax !\n'); end
printf('(PSOR): sortie avec k=%4i, erreur |x-xold|=%16.10f\n',k,err); 
endfunction



