数学模型: 目标函数:Min S=C1X1+C2X2+....+CnXn
约束条件:
a11 X1 +a12 X2+...+a1n Xn≥(=,≤)b1
a21 X1 +a22 X2+...+a2n Xn≥(=,≤)b2
am1 X1 +am2 X2+...+amn Xn≥(=,≤)bm
Xn >=0
Xn为决策变量,即各原料在配方中的比例;aij为技术系数,即各种原料相应的营养成分;bi为约束值,即配方应满足的各项营养指标或重量指标;Cj为成本系数,即每种原料的价格系数;m为约束条件的个数,n为变量个数。
目的:使得在Xn有解的情况下,得到s.问题我简单的描述了一下,我想请你帮我写一个模块函数,用单纯形法即可,请指明输入输出的参数类型(要用到多维数组) ,避免使用指针,用object pascal写,我用的开发工具是delphi,希望能直接在程序中用,输出结果有两个 ,一个是总成本s,一个是数组参数Xn, 另外需要有一个判断是否有解的函数。
约束条件:
a11 X1 +a12 X2+...+a1n Xn≥(=,≤)b1
a21 X1 +a22 X2+...+a2n Xn≥(=,≤)b2
am1 X1 +am2 X2+...+amn Xn≥(=,≤)bm
Xn >=0
Xn为决策变量,即各原料在配方中的比例;aij为技术系数,即各种原料相应的营养成分;bi为约束值,即配方应满足的各项营养指标或重量指标;Cj为成本系数,即每种原料的价格系数;m为约束条件的个数,n为变量个数。
目的:使得在Xn有解的情况下,得到s.问题我简单的描述了一下,我想请你帮我写一个模块函数,用单纯形法即可,请指明输入输出的参数类型(要用到多维数组) ,避免使用指针,用object pascal写,我用的开发工具是delphi,希望能直接在程序中用,输出结果有两个 ,一个是总成本s,一个是数组参数Xn, 另外需要有一个判断是否有解的函数。
说明
//M1 存放<=约束条件的个数 @I
//M2 存放>=约束条件的个数 @I
//M3 存放= 约束条件的个数 @I
//M 约束条件的总个数 M=M1+M2+M3 @I
//MP,NP 存储数组的物理维数,MP>=M+2,NP>=N+1
//A (M+2)*(N+1)个元素的二维数组,其中A 的前M+1行是输入,输出参数,输入时 按
// [ 0 a01 a02 a0n
b1 -a11 -a12 -a1n
b2 -a21 -a22 -a2n bm -am1 -am2 -amn bl al1 al2 aln ]IPOSV M个元素的一维整型数组,@O ,说明A的输出结果 若IZROV(J)=J,
}
/////////////////////////////////////////////////////////////////////////////////////////
//////////////////////////////////////////////////////////////////////////////////////////
matrx2 = array of array of Extended;Procedure SIMPLX(var A:matrx2; M, N, MP, NP, M1, M2, M3:integer;
var ICASE:integer;var IZROV, IPOSV:array of integer);
Label 1,2,3,4;
const
EPS =0.000001;
var
L1,L2,L3:array[0..100] of integer;
K,I,NL1,NL2,IR,KP,M12,IP,IQ,KH:integer; Q1,BMAX:real;
begin
If M <> M1 + M2 + M3 Then
begin
ShowMessage(' Bad input constraint counts');
Exit;
end;
NL1:=N;
For K:=1 To N do
begin
L1[K]:=K;
IZROV[K]:=K;
end;
NL2:=M;
For I:=1 To M do
begin
If A[I + 1, 1] < 0 Then
begin
ShowMessage('Bad input tableau.');
Exit;
End;
L2[I]:=I;
IPOSV[I]:=N + I;
end;
For I:=1 To M2 do
L3[I]:=1;
IR:=0;
If M2 + M3 = 0 Then GoTo 3;
IR:=1;
For K:=1 To N + 1 do
begin
Q1:=0;
For I:=M1 + 1 To M do
Q1:=Q1 + A[I + 1, K];
A[M + 2, K]:=-Q1;
end;
repeat
SIMP1(A, MP, NP, M + 1, L1, NL1, 0, KP, BMAX);
If (BMAX <= EPS) And (A[M + 2, 1] < -EPS) Then
begin
ICASE:=-1;
Exit;
end
Else If (BMAX <= EPS) And (A[M + 2, 1] <= EPS) Then
begin
M12:=M1 + M2 + 1;
If M12 <= M Then
begin
For IP:=M12 To M do
begin
If IPOSV[IP] = IP + N Then
begin
SIMP1(A, MP, NP, IP, L1, NL1, 1, KP, BMAX);
If BMAX > 0 Then GoTo 1;
end;
end;
end;
IR:=0;
M12:=M12 - 1;
If M1 + 1 > M12 Then Exit;
For I:=M1 + 1 To M12 do
begin
If L3[I - M1] = 1 Then
begin
For K:=1 To N + 1 do
A[I + 1, K]:=-A[I + 1, K];
end;
end;
goto 3;
end;
SIMP2(A, M, N, MP, NP, L2, NL2, IP, KP, Q1);
If IP = 0 Then
begin
ICASE:=-1;
Exit;
end;
If IPOSV[IP] >= N + M1 + M2 + 1 Then
begin
For K:=1 To NL1 do
If L1[K] = KP Then goto 4;
4: NL1:=NL1 - 1;
For IQ:=K To NL1 do
L1[IQ]:=L1[IQ + 1];
end
Else
begin
If IPOSV[IP] < N + M1 + 1 Then GoTo 2;
KH:=IPOSV[IP] - M1 - N;
If L3[KH] = 0 Then GoTo 2;
L3[KH]:=0;
end;
A[M + 2, KP + 1]:=A[M + 2, KP + 1] + 1;
For I:=1 To M + 2 do
A[I, KP + 1]:=-A[I, KP + 1];
2: IQ:=IZROV[KP];
IZROV[KP]:=IPOSV[IP];
IPOSV[IP]:=IQ;
until IR <> 0;
3: SIMP1(A, MP, NP, 0, L1, NL1, 0, KP, BMAX);
If BMAX <= 0 Then
begin
ICASE:=0;
Exit;
end;
SIMP2(A, M, N, MP, NP, L2, NL2, IP, KP, Q1);
If IP = 0 Then
begin
ICASE:=1;
Exit;
end;
SIMP3(A, MP, NP, M, N, IP, KP);
GoTo 2;
end;
Procedure SIMP1(var A:matrx2; MP, NP, MM:integer; LL:array of integer;
NLL, IABF:integer;var KP:integer; var BMAX:real);
var
K:integer; TEST:real;
begin
KP:=LL[1];
BMAX:=A[MM + 1, KP + 1];
For K:=2 To NLL do
begin
If IABF = 0 Then
TEST:=A[MM + 1, LL[K] + 1] - BMAX
Else
TEST:=Abs(A[MM + 1, LL[K] + 1]) - Abs(BMAX);
If TEST > 0 Then
begin
BMAX:=A[MM + 1, LL[K] + 1];
KP:=LL[K];
end;
end;
end;
Procedure SIMP2(var A:matrx2; M, N, MP, NP:integer; L2:array of integer;
NL2:integer;var IP, KP:integer; Q1:real);
label 1,2;
const
EPS = 0.000001;
var
I,J,II,K:integer; FLAG,Q,QP,Q0:real;
begin
IP:=0;
FLAG:=0;
For I:=1 To NL2 do
begin
If A[L2[I] + 1, KP + 1] < -EPS Then FLAG:=1;
If FLAG = 1 Then goto 1;
end;
1: If FLAG = 0 Then Exit;
Q1:=-A[L2[I] + 1, 1] / A[L2[I] + 1, KP + 1];
IP:=L2[I];
For I:=I + 1 To NL2 do
begin
II:=L2[I];
If A[II + 1, KP + 1] < -EPS Then
begin
Q:=-A[II + 1, 1] / A[II + 1, KP + 1];
If Q < Q1 Then
begin
IP:=II;
Q1:=Q;
end
Else If Q = Q1 Then
begin
For K:=1 To N do
begin
QP:=-A[IP + 1, K + 1] / A[IP + 1, KP + 1];
Q0:=-A[II + 1, K + 1] / A[II + 1, KP + 1];
If Q0 <> QP Then goto 2;
end;
2: If Q0 < QP Then IP:=II;
end;
end;
end;
end;Procedure SIMP3(var A:matrx2; MP, NP, I1, K1:integer;var IP, KP:integer);
var
PIV:real; II,KK:integer;
begin
PIV:=1 / A[IP + 1, KP + 1];
For II:=1 To I1 + 1 do
begin
If II - 1 <> IP Then
begin
A[II, KP + 1]:=A[II, KP + 1] * PIV;
For KK:=1 To K1 + 1 do
begin
If KK - 1 <> KP Then
A[II, KK]:=A[II, KK] - A[IP + 1, KK] * A[II, KP + 1];
end;
end;
end;
For KK:=1 To K1 + 1 do
If KK - 1 <> KP Then A[IP + 1, KK]:=-A[IP + 1, KK] * PIV;
A[IP + 1, KP + 1]:=PIV;
end;
所有这些为验证上述子过程而编的验证过程
implementation
//PROGRAM D9R10
//Driver for routine SIMPLX
uses
unit2;
{$R *.DFM}procedure TForm1.Button1Click(Sender: TObject);
const
s1='%8.2f';
N = 4; M = 4; NP = 5; MP = 6; M1 = 2; M2 = 1; M3 = 1; NM1M2 = N+M1+M2;
var
F:TextFile;
A:matrx2;
IZROV,IPOSV:array[0..4] of integer;
ANUM:array[0..5] of real;
TXT:array[0..7] of string;
ALPHA:array[0..5] of string;
I,J,JJ,JMAX,ICASE:integer;
STR1,STR2,STR3:string;
begin
SetLength(A,7,6);
TXT[1]:='x1'; TXT[2]:='x2'; TXT[3]:='x3'; TXT[4]:='x4';
TXT[5]:='y1'; TXT[6]:='y2'; TXT[7]:='y3';
A[1, 1]:=0; A[1, 2]:=1; A[1, 3]:=1; A[1, 4]:=3; A[1, 5]:=-0.5;
A[2, 1]:=740; A[2, 2]:=-1; A[2, 3]:=0; A[2, 4]:=-2; A[2, 5]:=0;
A[3, 1]:=0; A[3, 2]:=0; A[3, 3]:=-2; A[3, 4]:=0; A[3, 5]:=7;
A[4, 1]:=0.5; A[4, 2]:=0; A[4, 3]:=-1; A[4, 4]:=1; A[4, 5]:=-2;
A[5, 1]:=9; A[5, 2]:=-1; A[5, 3]:=-1; A[5, 4]:=-1; A[5, 5]:=-1;
A[6, 1]:=0; A[6, 2]:=0; A[6, 3]:=0; A[6, 4]:=0; A[6, 5]:=0;
SIMPLX(A, M, N, MP, NP, M1, M2, M3, ICASE, IZROV, IPOSV);
//输出计算结果到文件
AssignFile(F, 'd:\delphi_shu\p9\d9r10.dat');
Rewrite(F);
Writeln(F);
If ICASE = 1 Then
Writeln(F, ' Unbounded objective function')
Else If ICASE = -1 Then
Writeln(F, ' No solutions satisfy constraints given');
JJ:=1;
For I:=1 To N do
begin
If IZROV[I] <= N + M1 + M2 Then
begin
ALPHA[JJ]:=TXT[IZROV[I]];
JJ:=JJ + 1;
end;
end;
JMAX:=JJ - 1;
Writeln(F);
STR1:=' ';
For JJ:=1 To JMAX do
STR1:=STR1 + ALPHA[JJ]+' ';
Writeln(F, ' ',STR1);
For I:=1 To M + 1 do
begin
If I > 1 Then
ALPHA[1]:=TXT[IPOSV[I - 1]]
Else
ALPHA[1]:=' ';
ANUM[1]:=A[I, 1];
JJ:=2;
For J:=2 To N + 1 do
begin
If IZROV[J - 1] <= (N + M1 + M2) Then
begin
ANUM[JJ]:=A[I, J];
JJ:=JJ + 1;
end;
end;
JMAX:=JJ - 1;
STR2:=' ';
For JJ:=1 To JMAX do
begin
STR3:=' ';
STR3:=COPY(FloatToStr(ANUM[JJ]),1,7);
STR2:= STR2 + STR3 + ' ';
end;
Writeln(F, ALPHA[1],STR2);
end;
CloseFile(F);
//屏幕显示计算结果
memo1.Lines.LoadFromFile('d:\delphi_shu\p9\d9r10.dat');
end;