function z=f(x)
    n=length(x);
    z=(x(1)-1)^2;
    for i=1:n-1
    z=z+100*(x(i+1)-x(i)^2)^2;
    end
endfunction
//
n=2;N=50;
x=-1.5:(3/N):1.5;y=x;z=zeros(N+1,N+1)
for i=1:N+1
    for j=1:N+1
      z(i,j)=f([x(i),y(j)]);
    end  
end
clf()
contour2d(x,y,z,[0:0.1:1,2:10:200]);
//
//
//
function alpha=indice(n)
    p=(n+1)*(n+2)/2;
    alpha=zeros(p,n);
    for i=1:n
       alpha(i+1,i)=1; 
    end
    q=n+1;
    for i=1:n
        for j=i:n
          q=q+1; 
          if (i==j) then
           alpha(q,i)=2;   
          else
           alpha(q,i)=1;alpha(q,j)=1;   
          end
        end
    end           
endfunction
n=2;p=(n+1)*(n+2)/2;
X=rand(p,n);
x=X(2,:);
//function y=lagrange(X,f,x)
   [p,n]=size(X);  
   A=zeros(p,p);b=zeros(p,1);c=zeros(p,1);ind=indice(n);
   for i=1:p
       for j=1:p
          val=1; 
          for k=1:n
            val=val*X(j,k)^ind(i,k)  
          end
          A(j,i)=val;
       end
         b(i)=f(X(i,:))
         val=1;
         for k=1:n
            val=val*x(k)^ind(i,k)  
          end
          c(i)=val;
   end
   d=A\b;
   y=d'*c;
//endfunction
disp('conditionnement');disp(cond(A));
disp('lagrange:');disp(y);
disp('exact:');disp(f(x));
