(*   INPUT	: The associated datafile for PSIMPLEX (Revised
		  Simplex Method) is called "SimplexDatafile".

		  FIRST NUMBER in "SimplexDatafile" represents # of
		  VARIABLES in Linear  Program problem(LP).

		  SECOND NUMBER represents # of CONSTRINTS in LP.

		  Third set of number is M x N constraint matrix.

		  FOURTH set of number is 1 x M Right Hand Side(RHS)
		  matrix.

		  FIFTH set of number is 1 x N Cost Matrix in the 
		  objective function.

     Algorithm	: The Psimplex algorithm  is based on the revised
		  simplex method.

		  This algorithm solves the LP problems in STANDARD FORM
		  
			 minimize	cTx
			 s.t.
				        Ax=b
			                x>= 0
 			 where vector b is non-negative.


		  Maximum # of variables, and maximum # of constraints are
		  set to maxvar =50 and maxconstraint=40 respectively.
		  One can modify them according to their needs. But
		  maxconstraint is <= maxvar.
		  
		  EPS is small ral number such that if for ant real
		  number a, |a|< EPS, then a =0.0
		  EPS is set to 10^-16 <= EPS <= 10^-4.

    OUTPUT	: Outputs are
		  1. Check if a LP is feasible and optimal solution exists.
		  2. Determine the optimal basic variables and their 
		     corresponding values.
		  3. Determine the optimal value of the objective function

    NOTE	: This algorithm tests all nonbasic varables as candidates
		  for entering the basis.  This is very time comsuming   
								        *)
 		  


program Simplex(input,output,SimplexDatafile,SimplexOutfile);

const   maxvar  = 50;
	maxconstraint = 40;
        EPS   =  0.00001;

type	CHARFILE = file of char;
	ARRM2M2  = array[1..maxconstraint+2,1..maxconstraint+2] of real;
	ARRM2N   = array[1..maxconstraint+2,1..maxvar] of real;
        ARRM2	 = array[1..maxconstraint+2] of real;
        ARRN     = array[1..maxvar] of real;
        ARRM     = array[1..maxconstraint] of integer;

var     Nextreal : real;
        N,M	 : integer;
	A	 : ARRM2N;
        B,X	 : ARRM2;
	C	 : ARRN;
	W	 : ARRM;
	F	 : real;
	NOFEAS  : boolean;
        NOSOL    : boolean;
        SimplexDatafile, SimplexOutfile : CHARFILE;
        Nextint  : integer;


procedure Infile(var	N,M	: integer;
		 var	A	: ARRM2N;
		 var	B	: ARRM2;
		 var	C	: ARRN;
		 var	Nextreal	: real;
		 var	Nextint		: integer);

var row,column : integer;

begin
  reset(SimplexDatafile);
  readln(SimplexDatafile,Nextint);
  N := Nextint;
  readln(SimplexDatafile, Nextint);
  M := Nextint;
  for row := 1 to M do
  begin
    for column := 1 to N do
    begin
      read(SimplexDatafile,Nextreal);
      A[row,column] := Nextreal;
    end;
    readln(SimplexDatafile);
  end;

  for row := 1 to M do
  begin
    read(SimplexDatafile, Nextreal);
    B[row] := Nextreal;
  end;
  readln(SimplexDatafile);

  for column :=  1 to N do
  begin
    read(SimplexDatafile,Nextreal);
    C[column] := Nextreal;
  end;
  readln(SimplexDatafile);
end;


procedure PSIMPLEX(
       M,N          :integer;
       EPS          :real;
   var A            :ARRM2N;
   var B,X          :ARRM2;
   var C            :ARRN;
   var W            :ARRM;
   var F            :real;
   var NOFEAS,NOSOL :boolean);

   var I,J,K,L,P,Q  :integer;
       D,R,S        :real;
       U            :ARRM2M2;
       Y            :ARRM2;
       EX,PHASE,STOP:boolean;
begin
   NOFEAS:=false;  NOSOL:=false;
   P:=M+2;  Q:=M+2;
   PHASE:=true;
   K:=M+1;
   for J:=1 to N do begin
      A[K,J]:=C[J];
      S:=0.0;
      for I:=1 to M do S:=S-A[I,J];
      A[P,J]:=S
   end;  (* FOR J *)
   S:=0.0;
   for I:=1 to M do begin
      W[I]:=N+I;
      R:=B[I];  X[I]:=R;
      S:=S-R
   end;  (* FOR I *)
   X[K]:=0.0;  X[P]:=S;
   for I:=1 to P do begin
      for J:=1 to P do U[I,J]:=0.0;
      U[I,I]:=1.0
   end;
   STOP:=false;
   repeat  (* UNTIL STOP *)                              (* PHASE 1 *)
      if (X[P] >= -EPS) and PHASE then begin
         PHASE:=false;  Q:=M+1
      end;
      D:=0.0;                                             (* PHASE 2 *)
      for J:=1 to N do begin
         S:=0.0;
         for I:=1 to P do S:=S+U[Q,I]*A[I,J];
         if D > S then begin D:=S;  K:=J end
      end;  (* FOR J *)
      if D > -EPS then begin
         STOP:=true;
         if PHASE then NOFEAS:=true
         else F:=-X[Q]
      end
      else begin
         for I:=1 to Q do begin
            S:=0.0;
            for J:=1 to P do S:=S+U[I,J]*A[J,K];
            Y[I]:=S
         end;  (* FOR I *)
         EX:=true;
         for I:=1 to M do
            if Y[I] >= EPS then begin
               S:=X[I]/Y[I];
               if EX or (S < D) then begin D:=S;  L:=I end;
               EX:=false
            end;  (* IF Y[I] >= EPS *)
         if EX then begin NOSOL :=true;  STOP:=true end
         else begin
            W[L]:=K;  S:=1.0/Y[L];
            for J:=1 to M do U[L,J]:=U[L,J]*S;
            if L = 1 then I:=2 else I:=1;
            repeat
               S:=Y[I];  X[I]:=X[I]-D*S;
               for J:=1 to M do U[I,J]:=U[I,J]-U[L,J]*S;
               if I = L-1 then I:=I+2 else I:=I+1
            until I > Q;
            X[L]:=D
         end  (* ELSE: NOT EX *)
      end  (* ELSE: D <= -EPS *)
   until STOP
end;  (* PSIMPLEX *)


procedure Outfile(NOFEAS, NOSOL : boolean;
                  M	: integer;
		  F	: real;
		  W     : ARRM);

var counter : integer;

begin
  rewrite (SimplexOutfile);
  if (NOFEAS = false) and (NOSOL = false) then
  begin
  writeln(SimplexOutfile,'  NOFEAS = NOSOL  ', NOFEAS);
  writeln(SimplexOutfile,'      Basic-Var   Value');
  for counter := 1 to M do
    begin
      write(SimplexOutfile,W[counter]);
      writeln(SimplexOutfile,'    ',X[counter]:9:2);
    end;
    writeln(SimplexOutfile,'  Objective value for min problem is  ',F:9:2);
  end;
end;



begin (*  main *)
  Infile(N,M,A,B,C,Nextreal,Nextint);
  PSIMPLEX(M,N,EPS,A,B,X,C,W,F,NOFEAS,NOSOL);
  Outfile(NOFEAS,NOSOL,M,F,W);
end.