Constructor Create(_C:TExtArray; MaximiCe:boolean=false);
Constructor CreateBasis(const Simplex:TSimplex);
Constructor Copy(const Simplex:TSimplex);
Procedure AddCons(_B:extended; _A:TExtArray; Sign:TOperation);
Procedure SetAllLengths(Len:integer);
Function SimplexStep:integer;
Function CheckBasis:boolean;
Function FoundInBasis(num:integer): integer;
Function DoPrec(num:extended): extended;
Procedure NormaliCe;
Procedure MulString(Number:integer; Value:extended);
Procedure AddString(Num1,Num2:integer; Value:extended); {суммирование строки 1 со строкой 2, домноженной на коэффициент Value }
Function Solve:integer;
Function GetMin:extended;
Function GetSolution:TExtArray;
Destructor Free;
end;
TIntSimplex = class(TSimplex)
// CurX : TExtArray;
//CurL : extended;
// CurFound : boolean;
Constructor Create(_C:TExtArray; MaximiCe:boolean=false);
// Procedure DelLastCons;
Function IntSolve:integer;
Function GetIntMin:extended;
Function IsInteger(value:extended):boolean;
Function GetIntSolution:TExtArray;
// Function SearchCons(_B:extended;_A:TExtArray):integer;
end;
implementation
uses Math;
{ TSimplex }
Function TSimplex.DoPrec(num:extended): extended;
begin
if ((num < MAX_VAL) and (num > -MAX_VAL)) then
num := 0;
Result := num;
end;
procedure TSimplex.AddCons(_B: extended; _A: TExtArray; Sign: TOperation);
var
j : integer;
begin
if (Length(_A)>N) then SetAllLengths(Length(_A));
inc(M);
SetLength(Cons,M);
//if ((_B=0) and (Sign=Less)) then Sign:=Equal; //???
Cons[M-1].B:=_B;
Cons[M-1].Sign:=Sign;
SetLength(Cons[M-1].A,N);
for j:=0 to Length(_A)-1 do Cons[M-1].A[j]:=_A[j];
if Length(_A)<N then for j:=Length(_A) to N-1 do Cons[M-1].A[j]:=0;
end;
{суммирование строки 1 со строкой 2, домноженной на коэффициент Value }
procedure TSimplex.AddString(Num1, Num2: integer; Value: extended);
var
j : integer;
begin
for j:=0 to N-1 do Cons[Num1].A[j]:=Cons[Num1].A[j]+Cons[Num2].A[j]*Value;
Cons[Num1].B:=Cons[Num1].B+Cons[Num2].B*Value;
end;
function TSimplex.CheckBasis: boolean;
var
i,j,k : integer;
f : boolean;
begin
SetLength(Basis,M);
for i:=0 to M-1 do Basis[i]:=-1;
for j:=0 to N-1 do begin
f:=true;
k:=-1;
i:=0;
while (f and (i<M)) do begin
if ((Cons[i].A[j]<>0) and (Cons[i].A[j]<>1)) then f:=false;
if (Cons[i].A[j]=1) then begin
if (k=-1) then k:=i
else f:=false;
end;
inc(i);
end;
if (f and (k<>-1)) then Basis[k]:=j;
end;
f:=true;
for i:=0 to M-1 do f:=f and (Basis[i]<>-1);
Result:=f;
end;
constructor TSimplex.Create(_C: TExtArray; MaximiCe:boolean);
var
j : integer;
begin
N:=Length(_C);
RealN := N;
M:=0;
SetLength(C,N);
Max:=MaximiCe;
if (not MaximiCe) then for j:=0 to N-1 do C[j]:=-_C[j]
else for j:=0 to N-1 do C[j]:=_C[j];
Max:=MaximiCe;
L := 0;
end;
constructor TSimplex.Copy(const Simplex: TSimplex);
var
i,j : integer;
begin
M:=Simplex.M;
N:=Simplex.N;
RealN := Simplex.RealN;
SetLength(Cons,M);
SetLength(Basis,M);
SetLength(C,N);
Max:=Simplex.Max;
for i:=0 to M-1 do begin
SetLength(Cons[i].A,N);
Basis[i]:=-1;
for j:=0 to N-1 do Cons[i].A[j]:=Simplex.Cons[i].A[j];
Cons[i].B:=Simplex.Cons[i].B;
Cons[i].Sign:=Simplex.Cons[i].Sign;
end;
for i:=0 to Simplex.N-1 do C[i]:=Simplex.C[i];
L := Simplex.L;
end;
constructor TSimplex.CreateBasis(const Simplex: TSimplex);
var
i,j : integer;
begin
M:=Simplex.M;
N:=Simplex.N;
RealN := Simplex.RealN;
L := 0;
SetLength(Cons,M);
SetLength(Basis,M);
SetLength(C,N);
for i:=0 to N-1 do C[i]:=0;
for i:=0 to M-1 do begin
SetLength(Cons[i].A,N);
for j:=0 to N-1 do Cons[i].A[j]:=Simplex.Cons[i].A[j];
Cons[i].B:=Simplex.Cons[i].B;
Cons[i].Sign:=equal;
Cons[i].isT := false;
end;
for i:=0 to M-1 do begin
if (Simplex.Basis[i]<>-1) then Basis[i]:=Simplex.Basis[i]
else begin
SetAllLengths(N+1);
for j:=0 to M-1 do Cons[j].A[N-1]:=0;
Cons[i].A[N-1]:=1;
Cons[i].isT := true;
C[N-1] := 0;
for j:=0 to Simplex.N-1 do C[j] := C[j] + Simplex.Cons[i].A[j];
L := L + Cons[i].B;
end;
end;
end;
destructor TSimplex.Free;
begin
SetLength(C,0);
SetLength(Basis,0);
SetLength(Cons,0);
M:=0;
N:=0;
RealN := 0;
end;
function TSimplex.GetMin: extended;
var
i : integer;
begin
if (Max) then
Result := -L
else
Result := L;
end;
function TSimplex.GetSolution: TExtArray;
var
Solution : TExtArray;
i,j : integer;
begin
SetLength(Solution,RealN);
for j:=0 to RealN-1 do begin
Solution[j]:=0;
i:=0;
while ((i<M) and (Basis[i]<>j)) do inc(i);
if ((Basis[i]=j) and (i<M)) then Solution[j]:=Cons[i].B;
end;
Result:=Solution;
end;
procedure TSimplex.MulString(Number: integer; Value: extended);
var
j : integer;
begin
for j:=0 to N-1 do Cons[Number].A[j]:=Cons[Number].A[j]*Value;
Cons[Number].B:=Cons[Number].B*Value;
end;
procedure TSimplex.NormaliCe;
var
i : integer;
begin
for i:=0 to M-1 do if (Cons[i].Sign<>Equal) then begin
SetAllLengths(N+1);
if (Cons[i].Sign=Greater) then Cons[i].A[N-1]:=-1
else Cons[i].A[N-1]:=1;
Cons[i].Sign := Equal;
end;
end;
procedure TSimplex.SetAllLengths(Len: integer);
var
i, j : integer;
OldN : integer;
begin
OldN:=N;
N:=Len;
SetLength(C,N);
for i:=0 to M-1 do SetLength(Cons[i].A,N);
if (OldN<N) then begin
for j:=OldN to N-1 do begin
C[j]:=0;
for i:=0 to M-1 do Cons[i].A[j]:=0;
end;
end;
end;
function TSimplex.FoundInBasis(num:integer): integer;
var
i:integer;
f:boolean;
begin
f := false;
i := 0 ;
while (not f and (i<M)) do
begin
f := (Basis[i] = num);
inc(i);
end;
if (f) then
Result := i-1
else
Result := -1;
end;
function TSimplex.SimplexStep: integer;
var
i,j : integer;
f,opt : boolean;
x,y : integer; //координаты опорного элемента
CurMax : extended;
temp : array of TConstrain;
tempC : TExtArray;
begin
opt := true;
CurMax := -1;
for i := 0 to N-1 do
begin
//проверка на разрешимость
if (C[i] > 0) then
begin
opt := false; //а это попутная проверка на оптимальность
if (C[i] > CurMax) then //а это поиск ведущего столбца (максимальный элемент в C[i])
begin
CurMax := C[i];
x := i;
end;
f := true;
for j := 0 to M-1 do
f := f and (Cons[j].A[i] < 0);
if (f) then
begin
Result := SIMPLEX_NO_BOTTOM;
exit;
end;
end;
end;
if (opt) then
Result := SIMPLEX_DONE
else
begin
//зная номер ведущего столбца, ищем номер ведущей строки
CurMax := MaxExtended; //на самом деле тут будем искать минимум, а не Max
for j := 0 to M-1 do
if (Cons[j].A[x] > 0) then //идем только по положительным элементам
if (Cons[j].B/Cons[j].A[x] < CurMax) then
begin
CurMax := Cons[j].B/Cons[j].A[x];
y := j;
end
else if (DoPrec(Cons[j].B/Cons[j].A[x] - CurMax) = 0) then
if (Cons[j].isT) then
y := j;
//сохраняем текущие значения
SetLength(temp, M);
for j := 0 to M-1 do
begin
SetLength(temp[j].A, N);
for i := 0 to N-1 do
temp[j].A[i] := Cons[j].A[i];
temp[j].B := Cons[j].B;
end;
SetLength(tempC, N);
for i := 0 to N-1 do
tempC[i] := C[i];
//делаем пересчет таблицы
//строка делиться на ведущий элемент
MulString(y, 1/Cons[y].A[x]);
//преобразование остальных элементов
for j := 0 to M-1 do
begin
if (j <> y) then
begin
for i := 0 to N-1 do
begin
Cons[j].A[i] := DoPrec(temp[j].A[i] - temp[j].A[x]*temp[y].A[i]/temp[y].A[x]);
end;
Cons[j].B := DoPrec(temp[j].B - temp[j].A[x]*temp[y].B/temp[y].A[x]);
end
else
begin
for i := 0 to N-1 do
Cons[j].A[i] := DoPrec(Cons[j].A[i]);
end;
end;
//и строка с коэффициентами функции
for i := 0 to N-1 do
begin
C[i] := DoPrec(tempC[i] - tempC[x]*temp[y].A[i]/temp[y].A[x]);
end;
Basis[y] := x;
//и сама функция:
L := DoPrec(L - tempC[x]*temp[y].B/temp[y].A[x]);
for i:= 0 to M-1 do
SetLength(temp[i].A, 0);
SetLength(temp, 0);
SetLength(tempC, 0);
Result := SIMPLEX_NEXT_STEP;
end;
end;
function TSimplex.Solve: integer;
var
i,j : integer;
Simplex : TSimplex;
f : boolean;
Step : integer;
cc : extended;
begin
//oldN := N;
NormaliCe;
f:=false;
if (not CheckBasis) then begin
Simplex:=TSimplex.CreateBasis(self);
Simplex.Solve;
f:=Simplex.GetMin<>0;
if (not f) then for i:=0 to M-1 do begin
for j:=0 to N-1 do Cons[i].A[j]:=Simplex.Cons[i].A[j];
Cons[i].B:=Simplex.Cons[i].B;
Cons[i].isT := false;
Basis[i]:=Simplex.Basis[i];
cc := C[Basis[i]];
for j:=0 to N-1 do
C[j] := DoPrec(C[j] - cc*Cons[i].A[j]);
L := DoPrec(L - cc*Cons[i].B);
end;
Simplex.Free;
end;
if (f) then Step:=SIMPLEX_NO_SOLUTION
else repeat
Step:=SimplexStep;
until (Step<>SIMPLEX_NEXT_STEP);
//SetAllLengths(OldN);
Result:=Step;
end;
{ TIntSimplex }
constructor TIntSimplex.Create(_C:TExtArray; MaximiCe:boolean=false);
begin
//CurFound:=false;
inherited;
end;
function TIntSimplex.GetIntMin: extended;
begin
Result:=GetMin;
end;
function TIntSimplex.GetIntSolution: TExtArray;
begin
Result:=GetSolution;
end;
function TIntSimplex.IsInteger(Value:extended):boolean;
begin
Result:=((Value=floor(Value)) or (Value=ceil(Value)));
end;
function TIntSimplex.IntSolve: integer;
var
i : integer;
OldN : integer;
FractCol : integer;
FractRow : integer;
TmpX : TExtArray;
TmpCons : TExtArray;
NewValue : extended;
begin
if (Solve=SIMPLEX_DONE) then begin
//if (not CurFound or ((Simplex.GetMin<CurL) and not Max) or ((Simplex.GetMin>CurL) and Max)) then begin
TmpX:=GetSolution;
i:=0;
while ((i<RealN) and IsInteger(TmpX[i])) do inc(i);
FractCol:=i;
if (FractCol<>RealN) then begin // если найдена хотя бы одна нецелая переменная
OldN:=N;
SetLength(TmpCons,N);
FractRow := FoundInBasis(FractCol);
for i := 0 to N-1 do
if (FoundInBasis(i) = -1) then
TmpCons[i] := Cons[FractRow].A[i] - Floor(Cons[FractRow].A[i])
else
TmpCons[i] := 0;
NewValue := Cons[FractRow].B - Floor(Cons[FractRow].B);
//if (Max) then
AddCons(NewValue, TmpCons, Greater);
//else
// AddCons(NewValue, TmpCons, Less);
Result := IntSolve;
SetAllLengths(OldN); // удаляем пустые столбцы в конце, если они есть
end
else begin // если полученное решение - целочисленное\
Result := SIMPLEX_DONE;
end;
//end;
end
else
Result:=SIMPLEX_NO_SOLUTION;
end;
end.