function [xo,Ot,nS,lambda]=karmarkar(c,x0,A,b,tol,mxit,ip)
//   Solves the linear programming problem
//      min S(x)=c'x   subject to: A x = b, x >= 0
//   using the Karmarkar's algorithm.
//
//   [xo,Ot,nS,lambda]=karmarkar(c,x0,A,b,tol,mxit)
//
//   c: coefficients of the objective function
//   x0: initial point (must be feasible)
//   A: constraint matrix
//   b: constraint vector
//   tol: tolerance (default = 1e-4)
//   mxit: maximum number of iterations (default = 100)
//   ip: (0): no plot (default), (>0) plot figure ip with pause, (<0) plot figure ip
//   xo: optimal point
//   Ot: optimal value of S
//   nS: number of objective function evaluations
//   lambda: Lagrange multipliers

//   Copyright (c) 2001 by LASIM-DEQUI-UFRGS
//   $Revision: 1.0 $  $Date: 2001/08/15 23:05:15 $
//   Argimiro R. Secchi (arge@enq.ufrgs.br)

 if nargin < 4,
   error('karmarkar requires four input arguments');
 end
 if nargin < 5 || isempty(tol),
   tol=1e-4;
 end
 if nargin < 6 || isempty(mxit),
   mxit=100;
 end
 if nargin < 7 || isempty(mxit),
   ip=0;
 end

[m,n]=size(A);
erro=norm(c);
I=eye(n);
um=ones(n,1);
//gama=0.9995;
gama=0.8;

nS=0;
Q=I;
x1=x0;
x0=um;
c0=c;

while erro > tol && nS < mxit,
 nS=nS+1;
 D=diag(x1);
 Q=Q*D;
 A=A*D;
 c=D*c;

 B=inv(A*A')*A;
 P=I-A'*B;
 dp=-P*c;
 lambda=-B*c;
 z=-x0./dp;
 j=find(z>0);
 alfa=gama*min(z(j));
 x1=x0+alfa*dp;
 xo=Q*x1;
 x01=Q*x0;
 if ip,
   plot3([x01(1);xo(1)],[x01(2);xo(2)],[x01(3);xo(3)],'-ro'); hold on
   if ip > 0,
     disp('Pause: hit any key to continue...'); pause;
   end
 end
 erro=norm(dp);
end

Ot=c0'*xo;
