BBS 水木清华站 -- 精华区文章阅读
--------------------------------------------------------------------------------
发信人: FangQ (F the world), 信区: NumComp
标 题: Re: 急寻线性规划的源程序（c语言版的）
发信站: BBS 水木清华站 (Sat Jun 30 15:45:48 2001)
Message 5 in thread
寄件者：Vic Smyth ([email protected])
主旨：Re: Linear Programming
新闻群组：comp.programming
日期：1999/05/25
Tim,
Don't know about C++, but Linear Programming and Simplex are covered
in chapter 10 of Numerical Recipies in C. You can review the code
on-line at http://beta.ulib.org/webRoot/Books/Numerical_Recipes/
-- Bill Clay
A program that we use at college in the linear programming (math) class is
called Lindo. But I don't know if it's available as freeware. Math software
Maple can also do Simplex. (Though I haven't taken the class yet and haven't
tried it.)
【 在 aliali (承轩) 的大作中提到: 】
: 谢谢。
--
http://fangqq.dhs.org/
※ 来源:·BBS 水木清华站 smth.org·[FROM: 129.170.67.237]
发信人: FangQ (F the world), 信区: NumComp
标 题: Re: 急寻线性规划的源程序（c语言版的）
发信站: BBS 水木清华站 (Sat Jun 30 15:48:32 2001)
Attached below is a fairly simple Simplex Solver (written for
Turbo Pascal 3.0). Please let me know of how you use the
program, as well as any problems and any successes.
Regards,
Stephen
=====================================================================
Stephen F. Gale, B.Sc., CCP, I.S.P.
[email protected]
http://www.freenet.calgary.ab.ca/~sfgale/
=====================================================================
95/96 Program Director
Association for Systems Management - Calgary Chapter
http://www.freenet.calgary.ab.ca/populati/communit/asmcal/asmcal.html
=====================================================================
************************
*** SIMPLEX SOFTWARE ***
************************
program linearoptimization;
const rowmx = 72;
colmx = 112;
mxval = 1.0E+35;
zero = 0.0;
eqzero = 0.00000001;
var matrix : array [0..rowmx, 0..colmx] of real;
basis : array [1..rowmx] of integer;
basisp : array [1..rowmx] of integer;
minmax : real;
error : integer;
name : string[70];
filename : string[14];
ncon : integer; {number of constraints}
nvar : integer; {number of variables}
nltcon : integer; {number of less than constraints}
neqcon : integer; {number of equal to constraints}
ngtcon : integer; {number of greater than contraints}
trows1 : integer;
trows2 : integer;
tcols1 : integer;
tcols2 : integer;
tcols3 : integer;
function getmatrix( row, col : integer ) : real;
begin
getmatrix := matrix[ row, col ];
end;
procedure putmatrix( row, col : integer; value : real );
begin
matrix[ row, col ] := value;
end;
procedure initmatrix;
var row, col : integer;
begin
for row := 0 to rowmx do
for col := 0 to colmx do
putmatrix( row, col, zero );
end;
procedure price( var xcol : integer; trow : integer; var error : integer );
var quant, val : real;
col : integer;
begin
quant := -eqzero;
for col := 1 to tcols3 do
begin
val := getmatrix( trow, col );
if ( val < quant ) then
begin
xcol := col;
quant := val;
end;
end;
error := 0;
if ( quant = -eqzero ) then
error := 1;
end;
procedure leave( var xrow : integer; xcol : integer; var error : integer );
var quant, val : real;
row : integer;
begin
quant := mxval;
for row := 1 to ncon do
begin
val := getmatrix( row, xcol );
if ( val > eqzero ) then
begin
val := getmatrix( row, tcols2 ) / val;
if ( val < quant ) then
begin
xrow := row;
quant := val;
end;
end;
end;
error := 0;
if ( quant = mxval ) then
error := 2;
end;
procedure pivot( xrow, xcol : integer );
var value, val, vl : real;
row, col : integer;
begin
value := getmatrix( xrow, xcol );
for row := 1 to trows2 do
if ( row <> xrow ) then
begin
vl := getmatrix( row, xcol );
for col := 1 to tcols2 do
if ( col <> xcol ) then
begin
val := getmatrix( row, col ) - vl * getmatrix( xrow, col ) / value;
if ( abs( val ) < eqzero ) then
val := zero;
putmatrix( row, col, val );
end;
end;
for col := 1 to tcols2 do
putmatrix( xrow, col, getmatrix( xrow, col ) / value );
for row := 1 to trows2 do
putmatrix( row, xcol, zero );
putmatrix( xrow, xcol, 1.0);
basis[ xrow ] := xcol;
end;
procedure optimize( trow : integer; var error : integer);
var xrow, xcol, iterate : integer;
begin
repeat
price( xcol, trow, error );
if ( error = 0 ) then
leave( xrow, xcol, error );
if ( error = 0 ) then
pivot( xrow, xcol );
until ( error <> 0 )
end;
procedure simplex( var error : integer );
var val : real;
row, col : integer;
flag : boolean;
label 1000;
begin
if ( ncon <> nltcon ) then
begin
optimize( trows1, error );
if ( error > 1 ) then exit;
error := 3;
for row := 1 to ncon do
if ( basis[ row ] > tcols3 ) then
begin
if ( getmatrix( row, tcols2 ) > eqzero ) then exit;
flag := false;
col := 1;
repeat
if ( abs( getmatrix( row, col ) ) >= eqzero ) then
begin
pivot( row, col );
flag := true;
end;
col := col + 1;
until ( (flag) or (col > tcols3) );
end;
end;
error := 0;
optimize( trows2, error );
end;
procedure reader( var error : integer );
var row, col, column : integer;
value, amt : real;
filevar : text;
begin
error := 0;
writeln(con,'Problem file should be in the following format:');
writeln(con,' Line 1 : Up to 70 character problem description');
writeln(con,' Line 2 : +1 (for max) or -1 (for min); # of constraints; # of variables');
writeln(con,' Line 3 : # of <= constraints; # of = constraints; # of >= constraints');
writeln(con,' Next : Constraints coefficients and RHS value for each constraint');
writeln(con,' Last : Objective function coefficients');
writeln(con);
write (con,'Enter the filename containing the problem: ');
readln (con,filename);
assign(filevar,filename);
reset(filevar);
{ read the problem description }
readln( filevar, name );
{ read the minmax, number of constraints, number of variables }
readln( filevar, minmax, ncon, nvar );
minmax := -minmax;
{ read the number of less than, equal to, greater than contraints}
readln( filevar, nltcon, neqcon, ngtcon );
if ( ncon <> nltcon + neqcon + ngtcon ) then error := -1;
trows1 := ncon + 1;
trows2 := ncon + 2;
tcols1 := nvar + ncon + ngtcon;
tcols2 := tcols1 + 1;
tcols3 := nvar + nltcon + ngtcon;
{prepare matrix and basis}
for row := 1 to trows2 do
for col := 1 to tcols2 do
putmatrix( row, col, zero );
for row := 1 to ncon do
basis[ row ] := 0;
{prepare artificial and surplus variables}
for row := 1 to ncon do
if ( row <= nltcon ) then
begin
column := nvar + row;
basis[ row ] := column;
putmatrix( row, column, +1.0 );
end
else
begin
column := nvar + ngtcon + row;
basis[ row ] := column;
putmatrix( row, column, +1.0 );
if ( row > nltcon + neqcon ) then
begin
column := nvar - neqcon + row;
putmatrix( row, column, -1.0 );
putmatrix( trows1, column, +1.0 );
end
end;
{read matrix and right hand side}
for row := 1 to ncon do
begin
for col := 1 to nvar do
begin
read( filevar, value );
putmatrix( row, col, value );
end;
read( filevar, value );
putmatrix( row, 0, value );
putmatrix( row, tcols2, value );
readln( filevar );
end;
{ read the coefficients of the objective function }
for col := 1 to nvar do
begin
read( filevar, value);
putmatrix( 0, col, value * minmax );
putmatrix( trows2, col, value * minmax );
end;
readln( filevar );
{ calculate artifical variables }
for col := 1 to nvar do
begin
value := zero;
for row := nltcon+1 to ncon do
value := value - getmatrix( row, col );
putmatrix( trows1, col, value );
end;
close(filevar);
end;
procedure stats;
begin
writeln;
writeln(' * Your Variables : 1 through ', nvar);
if ( nltcon > 0 ) then
writeln(' * Slack Variables : ',nvar+1,' through ',nvar+nltcon);
if ( ngtcon > 0 ) then
writeln(' * Surplus Variables : ',nvar+nltcon+1,' through ',tcols3);
if ( nltcon <> ncon ) then
writeln(' * Artificial Variables : ',tcols3+1,' through ',tcols1);
end;
procedure setbasis;
var row, col : integer;
flag : boolean;
begin
for col := 1 to nvar+ncon do
begin
flag := false;
row := 1;
repeat
if ( basis[ row ] = col ) then
flag := true
else
row := row + 1;
until ( (flag) or (row > ncon) );
if (flag) then
basisp[ col ] := row
else
basisp[ col ] := 0;
end;
end;
procedure problem;
var row, col : integer;
begin
{filename and problem description}
writeln('Filename: ',filename);
writeln(' Problem: ',name);
writeln;
{objective function}
if minmax < 0 then
writeln('Maximize:')
else
writeln('Minimize:');
for col := 1 to nvar do
begin
write(minmax*getmatrix( trows2, col ):18:8,' * Var#',col:3);
if col <> nvar then
write(' + ');
writeln;
end;
writeln;
{constraints}
writeln('Subject to:');
for row := 1 to ncon do
begin
for col := 1 to nvar do
begin
if (col = 1) then
write(' Constraint #',row:3,'...')
else
write(' ':20);
write(getmatrix( row, col):18:8,' * Var#',col:3);
if col <> nvar then
writeln(' + ');
end;
if (row <= nltcon) then
write(' <= ')
else if (row > nltcon+neqcon) then
write(' >= ')
else
write(' = ');
writeln(getmatrix( row, tcols2 ):18:8);
end;
end;
procedure answer;
var smallpos : real;
smallptr : integer;
largeneg : real;
largeptr : integer;
value : real;
row, col, row1 : integer;
begin
setbasis;
stats;
{objective value}
writeln;
writeln('* Value of Objective Function: ',
-minmax*getmatrix( trows2, tcols2 ):18:8);
writeln;
writeln('* Variable Analysis *');
writeln('Var':3,'Value':18,'Reduced Cost':18);
writeln('---':3,'-----':18,'------------':18);
for col := 1 to nvar do
begin
write(col:3);
if (basisp[ col ] <> 0) then
write(getmatrix( basisp[ col ], tcols2 ):18:8)
else
write(0.0:18:8);
write(-minmax * getmatrix( trows2, col ):18:8);
writeln;
end;
writeln;
writeln('* Constraint Analysis *');
writeln('Row':3,'RHS Value':18,'Slack/Surplus':18,'Shadow Price':18);
writeln('---':3,'---------':18,'-------------':18,'------------':18);
for row := 1 to ncon do
begin
if (row <= nltcon) then
col := nvar + row
else if (row > nltcon+neqcon) then
col := nvar + row - neqcon
else
col := nvar + ngtcon + row;
write(row:3);
write(getmatrix( row, 0 ):18:8);
if (basisp[ col ] <> 0) then
write(getmatrix( basisp[ col ], tcols2 ):18:8)
else
write(0.0:18:8);
write(-minmax * getmatrix( trows2, col ):18:8);
writeln;
end;
writeln(' ');
writeln('* Sensitivity Analysis - Right Hand Side Ranging *');
writeln('Row':3,'Lower':18,'Current':18,'Upper':18,'OutLo':6,'OutUp':6);
writeln('---':3,'-----':18,'-------':18,'-----':18,'-----':6,'-----':6);
for row1 := 1 to ncon do
begin
if (row1 <= nltcon) then
col := nvar + row1
else if (row1 > nltcon+neqcon) then
col := nvar + row1 - neqcon
else
col := nvar + ngtcon + row1;
smallpos := +mxval;
smallptr := 0;
largeneg := -mxval;
largeptr := 0;
for row := 1 to ncon do
begin
value := getmatrix( row, col );
if ( value <> zero ) then
begin
value := getmatrix( row, tcols2 ) / value;
if (value>zero) and (value<smallpos) then
begin
smallpos := value;
smallptr := basis[ row ];
end;
if (value<zero) and (value>largeneg) then
begin
largeneg := value;
largeptr := basis[ row ];
end;
end;
end;
if (row1<=nltcon+neqcon) then
begin
write(row1:3);
if (smallpos <> mxval) then
write((getmatrix( row1, 0 )-smallpos):18:8)
else
write('No Limit':18);
write((getmatrix( row1, 0 ) ):18:8);
if (largeneg <> -mxval) then
write((getmatrix( row1, 0 )-largeneg):18:8)
else
write('No Limit':18);
if (smallptr <> 0) then
write(smallptr:6)
else
write('None':6);
if (largeptr <> 0) then
write(largeptr:6)
else
write('None':6);
end
else
begin
write(row1:3);
if (largeneg <> -mxval) then
write((getmatrix( row1, 0 )+largeneg):18:8)
else
write('No Limit':18);
write((getmatrix( row1, 0 ) ):18:8);
if (smallpos <> mxval) then
write((getmatrix( row1, 0 )+smallpos):18:8)
else
write('No Limit':18);
if (largeptr <> 0) then
write(largeptr:6)
else
write('None':6);
if (smallptr <> 0) then
write(smallptr:6)
else
write('None':6);
end;
writeln;
end;
writeln(' ');
writeln('* Sensitivity Analysis - Objective Coefficient Ranging *');
writeln('Var':3,'Lower':18,'Current':18,'Upper':18,'InLo':6,'InUp':6);
writeln('---':3,'-----':18,'-------':18,'-----':18,'----':6,'----':6);
for col := 1 to nvar do
begin
smallpos := +mxval;
smallptr := 0;
largeneg := -mxval;
largeptr := 0;
if (basisp[ col ] = 0) then
if (minmax < 0) then
begin
smallpos := -minmax*getmatrix( trows2, col );
smallptr := col;
end
else
begin
largeneg := -minmax*getmatrix( trows2, col );
largeptr := col;
end
else
for row := 1 to tcols3 do
if (basisp[ row ] = 0) then
begin
value := getmatrix( basisp[ col ], row );
if ( value <> zero ) then
begin
value := minmax*getmatrix( trows2, row ) / value;
if (value>zero) and (value<smallpos) then
begin
smallpos := value;
smallptr := row;
end;
if (value<zero) and (value>largeneg) then
begin
largeneg := value;
largeptr := row;
end;
end;
end;
write(col:3);
if (largeneg <> -mxval) then
write((minmax*getmatrix( 0, col )+largeneg):18:8)
else
write('No Limit':18);
write((minmax*getmatrix( 0, col ) ):18:8);
if (smallpos <> mxval) then
write((minmax*getmatrix( 0, col )+smallpos):18:8)
else
write('No Limit':18);
if (largeptr <> 0) then
write(largeptr:6)
else
write('None':6);
if (smallptr <> 0) then
write(smallptr:6)
else
write('None':6);
writeln;
end;
end;
procedure print;
var row, col : integer;
begin
writeln;
for row := 1 to trows2 do
begin
if (row > 0) and (row <= ncon) then
write(basis[ row ]:2,' ')
else
write(' ');
for col := 1 to tcols2 do
write(getmatrix( row, col ):9:3,' ');
writeln;
end;
writeln;
end;
begin
initmatrix;
writeln;
writeln('*** Linear Programming - Simplex Algorithm ***');
writeln;
reader( error );
problem;
if ( error = 0 ) then
simplex( error );
if ( error = 0 ) or ( error = 1 ) then
answer;
if ( error < 0 ) then
writeln('-- Inconsistent Data - Not Run --');
if ( error = 2 ) then
writeln('-- The Solution is Unbounded --');
if ( error = 3 ) then
writeln('-- The Problem is Infeasible --');
end.
********************
*** SAMPLE INPUT ***
********************
Photocopy - Sensitivity Analysis - Page 329
-1 5 4
2 0 3
1.00 0.00 0.00 0.00 2.00
0.00 1.00 0.00 0.00 1.50
300.00 1000.00 10.00 1500.00 2500.00
20.00 110.00 10.00 0.00 200.00
15.00 10.00 8.00 48.00 100.00
50.00 80.00 15.00 180.00
*********************
*** SAMPLE OUTPUT ***
*********************
*** Linear Programming - Simplex Algorithm ***
Filename: pc329.lp
Problem: Photocopy - Sensitivity Analysis - Page 329
Minimize:
50.00000000 * Var# 1 +
80.00000000 * Var# 2 +
15.00000000 * Var# 3 +
180.00000000 * Var# 4
Subject to:
Constraint # 1... 1.00000000 * Var# 1 +
0.00000000 * Var# 2 +
0.00000000 * Var# 3 +
0.00000000 * Var# 4 <= 2.00000000
Constraint # 2... 0.00000000 * Var# 1 +
1.00000000 * Var# 2 +
0.00000000 * Var# 3 +
0.00000000 * Var# 4 <= 1.50000000
Constraint # 3... 300.00000000 * Var# 1 +
1000.00000000 * Var# 2 +
10.00000000 * Var# 3 +
1500.00000000 * Var# 4 >= 2500.00000000
Constraint # 4... 20.00000000 * Var# 1 +
110.00000000 * Var# 2 +
10.00000000 * Var# 3 +
0.00000000 * Var# 4 >= 200.00000000
Constraint # 5... 15.00000000 * Var# 1 +
10.00000000 * Var# 2 +
8.00000000 * Var# 3 +
48.00000000 * Var# 4 >= 100.00000000
* Your Variables : 1 through 4
* Slack Variables : 5 through 6
* Surplus Variables : 7 through 9
* Artificial Variables : 10 through 12
* Value of Objective Function: 335.23437500
* Variable Analysis *
Var Value Reduced Cost
--- ----- ------------
1 0.00000000 -4.29687500
2 1.50000000 0.00000000
3 6.90104167 0.00000000
4 0.62065972 0.00000000
* Constraint Analysis *
Row RHS Value Slack/Surplus Shadow Price
--- --------- ------------- ------------
1 2.00000000 2.00000000 0.00000000
2 1.50000000 0.00000000 -0.46875000
3 2500.00000000 0.00000000 -0.06250000
4 200.00000000 34.01041667 0.00000000
5 100.00000000 0.00000000 -1.79687500
* Sensitivity Analysis - Right Hand Side Ranging *
Row Lower Current Upper OutLo OutUp
--- ----- ------- ----- ----- -----
1 0.00000000 2.00000000 No Limit 5 None
2 1.25469572 1.50000000 2.40506329 8 4
3 1606.25000000 2500.00000000 3316.25000000 4 8
4 No Limit 200.00000000 234.01041667 None 8
5 73.88000000 100.00000000 815.00000001 8 4
* Sensitivity Analysis - Objective Coefficient Ranging *
Var Lower Current Upper InLo InUp
--- ----- ------- ----- ---- ----
1 45.70312500 50.00000000 No Limit 1 None
2 No Limit 80.00000000 80.46875000 None 6
3 1.20000000 15.00000000 15.16363636 9 6
4 179.31645570 180.00000000 202.00000000 6 1
【 在 aliali (承轩) 的大作中提到: 】
: 谢谢。
--
http://fangqq.dhs.org/
※ 来源:·BBS 水木清华站 smth.org·[FROM: 129.170.67.237]
--------------------------------------------------------------------------------
[返回上一页]
![在DOS下运行不了ipconfig命令[图]](data:image/gif;base64,R0lGODlhAQABAIAAAP///wAAACwAAAAAAQABAAACAkQBADs=)