数学建模社区-数学中国

标题: 单纯形算法程序 [打印本页]

作者: zhyi    时间: 2005-4-18 21:45
标题: 单纯形算法程序
单纯形法程序,在VC++6.0 下测试通过! ; a9 O7 a7 N' p/ b2 b8 `9 P % B& I: S( N# G1 q1 {, ?/ ] C; v2 i* {" `! S/ I/ C

#include<iostream.h>" T. [( K4 Q8 Z' B #include<math.h>6 \" ^3 N1 H+ Z) R( J float matrix[100][100],x[100]; ( D9 V, \# L" _, @% q& {int a[100];9 g" b+ V$ b; N2 g int m,n,s,type; 6 S( |1 d& E T! j* Bint indexe,indexl,indexg; 3 p/ q: N0 `: A% @4 q///////////////////////////////// 9 r+ p3 D0 c! O$ tvoid jckxj()//基础可行解# J% i' ~& |3 }* y4 T {; _4 G$ K0 E" ^' s int i,j;0 @# D4 i3 w, j' ~& x% L, B; D6 f# L for(i=0;i<n;i++) 4 `1 X( Y/ F% d$ K% \ @$ m3 Z& S for(j=0;j<s;j++)9 Q( P' Y7 R0 \: ]! Z% r4 a if(matrix[j]==1&&a[j]==1) ( B+ [7 s: J9 M3 p# @' A { . w: ^6 L$ V3 x/ R x[j]=matrix;2 X5 U. l1 Z' ?4 G1 ~ j=s;$ O4 e. G8 y3 a$ D! U5 Y } 0 W l! a5 c) M" @- V% X/ q1 F ^ for(i=0;i<s;i++)/ V, Z, }. f. r6 H" C4 k5 v- a' q$ f if(a==0)x=0;8 ]) U7 P. V D4 ]5 S7 s/ h }

1 u9 ]+ U7 E$ s5 p( h" j0 X; o

int rj()//基解矩阵 q. ]6 j& c5 L9 S, x4 V% z$ f1 Y3 d {! x7 S% ^. R- X: ` int i;) F/ ~& P4 c4 P& A for(i=0;i<s;i++) w+ M$ q+ B m$ f& w9 G if(fabs(matrix[n])>=0.000001) 0 |4 N* ?/ X' b if(matrix[n]<0)return 0;3 f$ o( V1 r, l9 ` m return 1; % a' @1 ^2 q, Z. _) N- K. E}' Q# e- M2 i: f2 N! B int Min()//求最小的 9 G5 H+ n$ O/ X{ ' O0 G, O7 X' ^+ b int i,temp=0; ) y3 _. N4 y' J5 ~7 Y# y float min=matrix[n][0]; : q5 |; {0 U# x( z9 j; @5 l3 A4 ~7 @1 ` for(i=1;i<s;i++)! X9 p& F9 C) b if(min>matrix[n]) 7 ~% I4 ]: r* c g. U+ y l { - o6 x9 F A7 J7 T4 J2 v min=matrix[n];2 {# M& Q' h) {# w" G temp=i; # B: x2 R) L Z; Y( Q } % b1 D l! O9 G4 b" f return temp;$ N* p! R" l$ q. n# y6 r/ E/ I O8 h8 R } 9 X* _8 |* j: U' U7 Y/////////////////////////////////0 O1 S9 Q4 G! R2 Z* L/ K4 \ void JustArtificial()//人工变量$ I6 V. M7 [6 k { ' W8 W* t( |, V" V' Y9 q int i;% @- j" p, { w. K7 `* @( G for(i=m+indexe+indexl;i<s;i++)4 W5 D+ U: W P( x% } if(fabs(x)>=0.000001)$ |( B8 ` k' g" G { . d% n) o3 u! I+ A o* g cout<<"NO Answer\n";0 c7 K2 Q+ V8 b1 r9 y _ return;. o% i5 i/ a6 j. ] } " x4 J6 }, Q$ L7 S}, [2 I, b) ^( Y3 q: X# P4 Z /////////////////////7 l7 d- j! ^% k! m3 N) v% F int Check(int in)//检验5 }! M; [9 G9 ]6 Y; \( P. r, { {+ e, R. W8 t R+ T' j int i; 1 H5 `. M4 {1 Z5 { float maxl=-1; & ^/ E9 u8 n5 Y$ T3 o# { for(i=0;i<n;i++) 7 B h+ z) M+ m5 {* i( D if(fabs(matrix[in])>=0.000001&&maxl<matrix/matrix[in]). y9 A5 ?0 i' e: U. z4 m- v: E maxl=matrix/matrix[in]; - }7 T4 b$ G9 ?* Q5 ? if(maxl<0) " p- ?" I' s' |$ A% Z return 1;) h7 Y& n% ]% j3 F) p7 c return 0;: L6 Z% R( x0 z0 l- d. ? }( ]9 f9 p/ R8 p- m: W int SearchOut(int *temp,int in)//出基变量 3 e/ K" Z: o/ H0 d( R- Y{2 b8 c6 w' S8 W- o% m+ F2 V int i; % X& N, [2 o [0 F float min=10000; : z F8 ] p7 }9 o: M7 ? for(i=0;i<n;i++) 4 [* H& n- \4 T6 F; P" I' U if(fabs(matrix[in])>=0.000001&&(matrix/matrix[in]>=0) 4 t8 Z# B* z0 \% C' [9 v+ @$ S &&min>matrix/matrix[in]), V* n% |: n, O1 t) c { " s' _/ j1 F: Q, ]1 r min=matrix/matrix[in]; % h# e6 e! R0 H* w4 i. o- w; F *temp=i;, R1 F5 q% Z2 Q, b. E( h } 6 L z) j3 i1 W7 ~ for(i=0;i<s;i++) - ]9 i* E1 r' a" q, W0 Q) ]* f if(a=1&&matrix[*temp]==1) 7 f0 C1 r9 R/ M1 Q1 B0 @ return i; 9 ]/ l2 y: N. d4 y4 n}4 L1 K4 w7 [: }/ \( _ ///////////////////////////////// * L+ @: @0 Q5 m2 C9 O8 |! Z. Rvoid Mto(int in,int temp)$ _9 f6 h0 }+ T6 Q& g+ K3 b, E# n {3 C+ J( @1 o6 s+ V. ~/ y int i;% i8 F8 n' i/ M% } for(i=0;i<=s;i++)1 d( I# H# x6 f% u% a) H if(i!=in) T0 [$ P% p! J9 Y, { matrix[temp]=matrix[temp]/matrix[temp][in]; 0 E$ T4 {) ?2 K& g/ J* K6 x matrix[temp][in]=1; 3 G" E# g3 n$ w2 R+ i. t}( a5 o" c: a4 H /////////////////////////////# z4 o" K# U- I- @/ t void Be(int temp,int in)//初等变换 ' Z1 ^" ^7 J/ Y1 C{ * H2 C; a; o; k. Y5 p1 G int i,j; 5 `; O' J# v/ g float c; 8 A6 U6 f8 z. ]7 ~ for(i=0;i<=n;i++)5 o7 M7 \1 ~% g1 Y; O {3 h( { I$ O2 @+ ? c=matrix[in]/matrix[temp][in]; 3 ]. h' y9 P" m9 w& _5 f! K if(i!=temp) & O: K y: N+ k4 J( N. T for(j=0;j<=s;j++) 1 a% c+ f- K& a# ^1 y9 f* A matrix[j]=matrix[j]-matrix[temp][j]*c; / j" s% d# j6 C } " N5 i [2 f2 ^3 w# _) T: a}- ?) S( Z, p8 v- w4 { ////////////////////////// 3 e' h4 |" ?5 {0 r$ d2 h4 N2 Wvoid Achange(int in,int out)//出基入基转换 ) G: D: @5 `4 d, a{0 w( v% D( S6 H' I int temp=a[in];: O3 D' A: P. n* k a[in]=a[out]; 0 w$ o- v# F8 T7 r, T( B/ K a[out]=temp;. C* w8 |3 ^7 Q3 _, l5 w+ \ } ' E9 i! ~9 k9 `* I& x//////////////////////// * t$ K! v: C% c$ ?5 lvoid Print() , D% o" H6 D+ u{ ) ^( Y7 ~4 [- f4 e( s int i,j,k,temp=0; 5 y _0 b1 ]* D4 S/ y" P2 |1 W; k for(i=0;i<n;i++)$ |; N. ]7 d/ M. e {$ R6 V( m1 g) e for(k=temp;k<s;k++) 7 O) y/ ~) F* u2 H if(a[k]==1) . ?3 Y2 t6 `' n: o8 Q { % r7 ~4 Y3 m, ^# E" X cout<<k; 4 W7 t! O( `1 `# [1 k temp=k+1;6 ?' t+ ]6 w5 u/ l- v' [& M k=s; 9 y6 V5 b- ~/ R }8 N# [% H" L5 t9 H for(j=0;j<=s;j++) 9 b( c) d g6 M cout<<matrix[j]; 3 Z4 [ O9 @" |2 h. T0 |% U cout<<"\n"; N3 B% \& J8 D } 1 [( k; L0 R" x: P7 R cout<<"Rj"; & e4 R9 P, ^% Z. o for(j=0;j<=s;j++) 2 W4 F7 v& z0 Z$ W cout<<matrix[n][j];1 C( ^" H, x( f( w9 o# e0 @* C cout<<"\n";" l% e2 E. h9 m }% k8 m* h! T" T8 V, Z3 L @! m5 a3 q2 d. _ //////////////////////// ) G+ p! ~/ J% h) x0 B: ]4 rvoid InitPrint() 9 M. p9 m( ]) Z{ * j ^: u& B# f: f- x, o int i; 0 a" t& k! l* T' F% u5 x cout<<"X"; ' H$ ~9 `; f3 }7 S, r- X' Y for(i=0;i<s;i++)1 c h: `( {9 G+ G0 ]" f cout<<i;3 c" d. p% o. N! F+ \' X" ] cout<<"b\n"; ' z) I4 A* M! y1 o, \# g, F+ Y cout<<" "; f) S4 l% U; E* W. @" Z, N' N cout<<"\n"; ) a' c" K. \- n& s}- I q' |5 G* T- J" u% ~) R //////////////////; h- b W/ A8 {% Q void Result() 8 S3 }( ], z/ u+ a8 G5 d; @, u! K{ / x% r! S% U; Z3 s int i; 8 s( ]. a+ H: f/ h cout<<"("; * @6 ]9 Y" q: k+ I% h4 F5 N" X for(i=0;i<s;i++) , ], s n2 Z( E: ] cout<<x;' j: G) R8 |# l& S9 {) a% D cout<<")";; h4 k4 j+ P0 e$ P1 P if(type==1)8 J( Q: ?# {" ]3 h) e3 F! z( C cout<<"Zmax="<<matrix[n]; 2 S! w" C! i K$ H- a+ P: {. ]/ I8 k else cout<<"Zmin="<<matrix[n];) M# S# f' w% t9 C& m2 B$ x }8 C+ V2 S0 j4 M; S" U- E //////////////////////5 [/ c: l$ `$ k2 Z: ? void PrintResult() 6 l) O- o; I5 B+ ]4 D{ & X9 f5 ^2 o: @1 Z: ? if(type==0) 4 Z/ R& G% a* D- g8 l8 }' A cout<<"the Minimal:"<<matrix[n];" F1 `' q$ I O8 E x else cout<<"theMaximum:"<<matrix[n];7 B9 O W+ u% i0 E0 p7 G9 Z } $ P5 E- P; ^9 W/ S( S G% r////////////////////////////////7 A3 V @, m3 Y/ b% V( y3 i% d: w void Merge(float nget[][100],float nlet[][100],float net[][100],float b[])//合并 & O( B' g1 P9 s; m$ ]' X- j. e& u{ 9 J7 p' p; Y' N; Y2 b- A* z int i,j;3 i9 _6 F7 V# s | for(i=0;i<n;i++)& i4 L5 |4 [( L1 j: D \ {; u1 S: P' y0 u: }/ H for(j=m;j<m+indexe;j++)3 F% g% X5 T& n" {/ s1 I( L8 |5 @) F if(nget[j-m]!=-1)matrix[j]=0;+ _0 w0 S. B9 x5 x* v else matrix[j]=-1;& l- G; K4 c& M' n, N for(j=m+indexe;j<m+indexe+indexl;j++)8 K. a; Q) i. Y/ T M d if(nlet[j-m-indexe]!=1)matrix[j]=0;4 _, e& Y& V# }. z& P else matrix[j]=1; 4 k) D; w* t( A. ? for(j=m+indexe+indexl;j<s;j++) 0 n6 w* O# \! Y6 t if(net[j-m-indexe-indexl]!=1)matrix[j]=0;% O9 D m9 G9 p6 ^4 L* C else matrix[j]=1;0 @) s* C% r, [5 |# A } ! i$ O' }6 O c3 j" P for(i=m;i<m+indexe+indexl;i++)2 T. Q0 a+ l' w; `/ E. I matrix[n]=0; , o5 k( C1 @- s/ `, ] for(i=m+indexe+indexl;i<s;i++)1 E. X4 \ b1 _# C% k matrix[n]=100; 9 _. ~4 b: M6 ]5 n( T* j matrix[n]=0; 9 y0 Z6 O) [2 T+ K. y}

' s/ a1 O# N9 _4 d9 n5 ]

///////////////////////////# X( s+ T: J1 C; I! Y void ProcessA()//初始a[] 2 C. e# q: J. s+ a/ Y& u# C{* s# v1 p- r! r9 k4 l int i;, ?( ]+ V$ G; Y* r& e8 Y. L# o for(i=0;i<m+indexe;i++) / O! g* d1 ]- q* K' E a=0;3 U% A1 n0 `* x( _% `/ f for(i=m+indexe;i<s;i++)3 x. Q- o+ R# A6 I: h/ B a=1;' b/ Y' ]: j, z: z }

! C4 L8 _1 w5 b1 u2 ^2 x' z' u

//////////////////////////////// . {$ H5 U5 l. {% V/ bvoid Input(float b[],int code[]) * q8 ]0 V0 e( V b3 V) ]9 }{ 8 ]" b" f4 R0 ~+ c8 a+ Y* H. i int i=0;int j=0;4 D5 y" e. M' O) J cout<<"The equator variable and Restrictor\n";9 C+ E9 W2 \8 E6 u cin>>m>>n;& g/ e( w8 |! y& w6 ` for(i=0;i<n;i++)6 K9 {7 s/ j. d3 w5 B2 h' a { ; L9 q* m4 R( X3 ~0 S) y cout<<"Inputb[] and Restrictor code0:<=1 1:= 2:>=\n";6 N6 k& |9 M! Z: z cin>>b>>code;+ K1 f* y, U! V) U( ^ cout<<"The 系数 \n";' w. K* w% E A6 _. f( D4 g for(i=0;j<m;j++)& F2 g* N$ Q8 u2 L) p cin>>matrix[j];. [; b0 C3 R D& j+ ~ }3 r1 ]+ T- ^: q5 e; N( m6 A cout<<"the type 0:Min 1:max\n"; 5 c$ c+ W; Y* u E2 U0 o" M" c do{! I6 a$ B; w8 G, B6 R/ k9 v' G5 @ cin>>type;, j6 z v7 m7 l8 | if(type!=0&&type!=1): a' Y R* {4 B4 | cout<<"error,ReInput!\n"; ; W+ f: u V8 A" Z! O+ y" g ]3 Z; L }while(type!=0&&type!=1); ! R0 y4 M3 U: g- k; ^ cout<<"the Z\n"; $ \2 h4 d$ o& l7 O for(i=0;i<m;i++) 8 i! i# J) q5 }3 q9 c/ b8 R cin>>matrix[n]; ) Z% O( r+ e& b- d( J; U if(type==1)+ @9 R2 O% W5 O" V0 ? for(i=0;i<m;i++) ' K* ]( [* ?+ Z9 [/ r- H matrix[n]=-matrix[n]; 1 q2 e$ Z* g2 [! \. K( f. P* w' g- y}

$ Y# l6 h6 L8 b/ S: t" p0 ?9 l q

0 C( g0 s3 w; l1 w9 r7 ?////////////////// ! x8 r( l! {- t7 t# ]5 C; T* H void Xartificial()//消去人工变量8 T1 R( p6 t0 b+ w3 b {: i( Y8 p: O' O int i,j,k; * L& k& ?9 v" ?: W2 ?! \ if(indexg!=0)6 B8 j" S0 m, M" ^% ?5 g" n {$ a. Z- F0 y3 u5 L1 H/ V for(i=m+indexe+indexl;i<s;i++)& j& t% Y% N5 a2 @* a2 F { - V/ Y& v; ~3 F" b7 J% ^ for(j=0;j<n;j++) ) {8 l$ X- w2 H3 j if(matrix[j]==1) 3 L; h- E! ?+ E2 x5 T3 y {3 u( K/ F8 Q* Z& c6 K for(k=0;k<=s;k++)" {' S) A9 P8 d: M8 T6 q$ { matrix[n][k]=matrix[n][k]-matrix[j][k]*100;, b& p4 F2 f7 u5 r" s% ^ j=n; / _7 F; l1 Y. V. V& N } . ~( d0 Q3 p- z- q }4 e6 S ]: G5 j8 w( A S }3 }, k, B) R8 b, r9 n: ^ }

4 s; I# ]3 V v7 O1 X: O

////////////////////////////////////////////////% W1 ~, a! ^/ p3 n void Process(float c[][100],int row,int vol) & [5 Q) T# R5 P b4 u{- K1 A! f0 ]. y. C int i; 2 G8 E, Z+ N8 s8 `7 t7 @ for(i=0;i<n;i++). J( e' t7 p& x) E; { if(i!=row)c[vol]=0;) H3 a- b K* k) M9 d/ ^ }2 n# t' A+ Z: V, l4 }1 k+ Z //////////////////////5 i Z- z6 D9 b! A void Start(float b[],int code[])) N( b* L5 s% N3 E {3 a9 {4 z" U0 m5 J& T int i; & @% W1 ]4 H4 j" n float nget[100][100],nlet[100][100],net[100][100]; ' c. l9 g ~ U/ m3 m* a1 J2 H indexe=indexl=indexg=0; 6 H' a2 h; _$ N( v4 c" a0 S4 K6 X for(i=0;i<n;i++) 0 K% `) q" e" x3 j! B { ! w3 ~; n- P3 G$ v+ b, I if(code==0){nlet[indexl++]=1rocess(nlet,i,indexl-1);}! s* m4 K8 s9 i if(code==1){net[indexl++]=1rocess(net,i,indexg-1);} 5 W+ e, ^! f1 m: I$ i9 A! T if(code==2){1 L9 m3 f3 K% Q) f net[indexg++]=1;- o7 c$ K6 g E) Q$ d nget[indexe++]=-1;( G# D7 t5 \& C# {" t Process(net,i,indexg-1)rocess(nlet,i,indexe-1);/ M7 x1 s7 U7 T! e* G( \, Q& ]7 E; L }1 j) H+ [3 z. L4 L/ v } $ o( {* y3 ^2 s8 v9 ?0 I s=indexe+indexl+indexg+m; ( x. }7 S/ o3 X/ F/ ]/ A5 s Merge(nget,nlet,net,b); * X! ^* G% C0 y ProcessA(); ) R& U. H! n: j( Z9 J5 I" \, y# C InitPrint(); z! _5 H1 J- X: G Xartificial(); ; d4 |: b, ?* x7 \! Q}

! ~! f* L# t1 X( R

void Simplix()//单纯形法 ) ?( E1 n4 D$ x! f1 G0 c9 _{ " c" x! R! F3 c: d; p3 y int in,out,temp=0; f. a! H9 X, h# w! c. `' @1 W; W while(1) ; i3 @8 C: Z S! d( ?' _ {4 N% ]3 s4 A8 f) ]9 ? jckxj();& ^0 X" |6 {7 ]2 t* l9 w7 @8 D Print();. L" z" y6 {: p7 h Result(); " i6 I' I5 o! a% m- U' V, A if(!rj()) in=Min();, S/ L% t* l* _9 @) m else{ * J. U, L) |- r- [ if(indexg!=0) 9 i% E, l2 Y5 @7 t JustArtificial(); 2 G% Z4 a) S2 _6 a$ L/ d$ F" i PrintResult(); % e7 G3 z0 a7 y6 ^) U" A return;8 J J" I8 e' e6 \+ F }: G" T+ A, Q; j0 [9 l if(Check(in))* Q7 \1 {3 W5 [ {) d% K- u5 y1 N5 O F3 m cout<<"No Delimition\n"; g, U2 w0 A- _ return; 4 X4 x8 R9 r0 u( X* k, C } * z( j. r% d! x* G' ^0 Q( L; _ out=SearchOut(&temp,in);1 p* y2 v4 c6 u7 o2 u' l/ U* M+ u. d Mto(in,temp);7 K% J$ w- r1 i% s! D Be(temp,in);3 Q' p. c2 s. ^+ K Achange(in,out);$ z( u4 n+ H5 b4 _ } + ~# c/ H. s3 I3 S}

" ?9 `2 ]1 Y. q% R$ S% g

void main()) D7 X1 O2 M5 {& H8 x- o {6 L$ ^7 D" m. d( n int code[100];//输入符号标记' X$ a6 b9 E0 S5 H4 `5 Q! s float b[100];+ p6 D. ]) r5 ?, q1 C; M Input(b,code);//初始化5 j' d6 P& N- Q1 \% A Start(b,code);//标准化行: x: L# @& U0 b Simplix();4 d# f9 V# F. s& b3 [# G% ^ } # a3 A `( g9 J$ ]/ \


作者: lvming    时间: 2005-4-25 22:44
谢了
作者: M_Tramp    时间: 2005-4-29 18:10
不错啊!!
作者: pangaogao    时间: 2005-5-4 00:15
Good!
作者: 那时花开    时间: 2005-6-3 22:57
好用不?
作者: 那时花开    时间: 2005-6-3 22:58
资源分配问题可以用这个程序吗?
作者: monkeytail    时间: 2005-6-10 01:42
哇,楼主强人!
作者: monkeytail    时间: 2005-6-10 01:43
单纯行算法是指数算法啊!
作者: bnulj    时间: 2005-6-11 13:14
这么长!看不懂啊。
作者: mengfanqi    时间: 2005-9-6 19:41

so well!


作者: 天亮    时间: 2005-9-17 11:45
zhyi,你好.我有一问题请你指教阿.有一个10维的有约束的函数,现在我已经有了一些可行点,其数目可能大于或者小于10,可行点之间的关系不知道,请问怎样运用单纯行法求函数的最小值.
作者: chenbilian158    时间: 2006-3-7 19:14
试试
作者: luom    时间: 2006-3-15 13:43

还有其他算法吗关于仿生算法的谢谢


作者: jixian    时间: 2006-4-16 10:42

作者: 海云边缘    时间: 2006-5-5 21:44

有C語言的嗎?我要C語言的


作者: suolunga    时间: 2006-5-8 20:07
[em01][em01][em01]
作者: sfhnlgdx    时间: 2006-5-13 09:15

谢谢,我正在找这个程序,终于找到了!


作者: chz0829    时间: 2006-6-3 00:41
看不懂
作者: hustyangyang    时间: 2006-7-19 14:59

楼上的的确是位高手,这两天由于写论文要用到单纯形法的程序,就拿来用了,但是发现有些错误,没有结果或者得出的答案有问题,于是我仔细看了一下这个程序,将楼上的笔误,和一些其他的小错误改正如下,
注意: 我不是用纯 VC++  我的输出函数用的是 printf (),这样可以动态观察结果输出,以便于修改,这个程序是我花了一天的时间改的,太长也太烦琐,我看的也不是很清楚,只对 限制条件是 <= 的情况进行了详细的察看,和修改,其他的情况就没有看了,程序中相应的部分有比较详细的注解,如果有需要可以和我联系,我们一起共同探讨。谢谢!QQ:116490942

 

#include<stdio.h>
#include<iostream.h>
#include<math.h>
float matrix[100][100],x[100];
int a[100];
int m,n,s,type;
int indexe,indexl,indexg;
/////////////////////////////////
void jckxj()//基础可行解
{
 int i,j;                        //基础可行解即为 非基变量对应的x=b, 基变量对应的解为0
 for(i=0;i<n;i++)
  for(j=0;j<s;j++)
   if(matrix[j]==1&&a[j]==1)
   {
    x[j]=matrix;
    j=s;
   }
   for(i=0;i<s;i++)
    if(a==0)x=0;
}

int rj()//基解矩阵
{
 int i;
 for(i=0;i<s;i++)
  if(fabs(matrix[n])>=0.000001)
   if(matrix[n]<0)return 0;
   return 1;
}

int Min()//求最小的
{
    int i,temp=0;
 float min=matrix[n][0];
 for(i=1;i<s;i++)
  if(min>matrix[n])
  {
   min=matrix[n];
   temp=i;
  }
 return temp;
}
/////////////////////////////////
void JustArtificial()//人工变量
{
 int i;
 for(i=m+indexe+indexl;i<s;i++)
  if(fabs(x)>=0.000001)
  {
   cout<<"NO Answer\n";
   return;
  }
}
/////////////////////
int Check(int in)//检验
{
 int i;
 float maxl=-1;
 for(i=0;i<n;i++)
  if(fabs(matrix[in])>=0.000001&&maxl<matrix/matrix[in])
   maxl=matrix/matrix[in];
  if(maxl<0)
   return 1;
  return 0;
}

int SearchOut(int *temp,int in)//出基变量
{
 int i;
 float min=10000;
 for(i=0;i<n;i++)
  if(fabs(matrix[in])>=0.000001&&(matrix/matrix[in]>=0)
   &&min>matrix/matrix[in])
  {
   min=matrix/matrix[in];
   *temp=i;
  }  
 for(i=0;i<s;i++)
  if(a=1&&matrix[*temp]==1)
 return i;
}
/////////////////////////////////
void Mto(int in,int temp)
{
 int i;
 for(i=0;i<=s;i++)
  if(i!=in)
   matrix[temp]=matrix[temp]/matrix[temp][in];
  matrix[temp][in]=1;
}
/////////////////////////////
void Be(int temp,int in)//初等变换
{
 int i,j;
 float c;
 for(i=0;i<=n;i++)
 {
  c=matrix[in]/matrix[temp][in];
  if(i!=temp)
   for(j=0;j<=s;j++)
    matrix[j]=matrix[j]-matrix[temp][j]*c;
 }
}
//////////////////////////
void Achange(int in,int out)//出基入基转换
{
 int temp=a[in];
 a[in]=a[out];
 a[out]=temp;
}
////////////////////////
void Print()
{
 int i,j,k,temp=0;
 for(i=0;i<n;i++)
 {
  for(k=temp;k<s;k++)
   if(a[k]==1)
   {
    printf(" %d ",k);
    temp=k+1;
    k=s;
   }
  for(j=0;j<=s;j++)
   printf( " %.0f ",matrix[j]  );
  printf("\n");
 }
 printf("Rj\n");
 for(j=0;j<=s;j++)
  printf(" %.0f ",matrix[n][j] );
 printf("\n");
}
////////////////////////
void InitPrint()
{
 int i;
 printf(" X ");
 for(i=0;i<s;i++)
  printf(" %d ",i);
    printf(" b \n" );
 printf("  " );
 printf("\n" );
}
//////////////////
void Result()
{
 int i;
 printf( "(" );
 for(i=0;i<s;i++)
  printf( "%.0f",x );
 printf( ")" );
 if(type==1)
  printf("\nZmax= %.0f\n" ,matrix[n] );
 else printf("\nZmin=%.0f\n", matrix[n] );
}
//////////////////////
void PrintResult()
{
 if(type==0)
  printf("the Minimal: %.2f\n", matrix[n] );
 else
  printf("the Maximum: %.2f\n", matrix[n] );
}
////////////////////////////////
void Merge(float nget[][100],float nlet[][100],float net[][100],float b[])//合并
{                           
                                                             //置成我们需要的矩阵,最后一列置成0
                                                             // 1  2  1  0  0  0
                                                             // 4  0  0  1  0  0
                                                             // 0  4  0  0  1  0
                                                             //-2 -3  0  0  0  0
 int i,j;                                             
 for(i=0;i<n;i++)
 {
  for(j=m;j<m+indexe;j++)
   if(nget[j-m]!=-1)matrix[j]=0;
   else matrix[j]=-1;
  for(j=m+indexe;j<m+indexe+indexl;j++)
    if(nlet[j-m-indexe]!=1)matrix[j]=0;
    else matrix[j]=1;
  for(j=m+indexe+indexl;j<s;j++)
    if(net[j-m-indexe-indexl]!=1)matrix[j]=0;
    else matrix[j]=1;
 }
 for(i=m;i<m+indexe+indexl;i++)  //将目标函数中人工变量的系数补为0
  matrix[n]=0;
 for(i=m+indexe+indexl;i<s;i++)
  matrix[n]=100;              //置成M
 for(i=0;i<n;i++)                   //把b[]的值赋给matrix 
     matrix=b;
 matrix[n]=0;
}


///////////////////////////
void ProcessA()//初始a[]    a[]是标记基变量,若为基变量则为0,若为非基变量则为1
{
 int i;
 for(i=0;i<m+indexe;i++)
  a=0;
 for(i=m+indexe;i<s;i++)
  a=1;
}

////////////////////////////////
void Input(float b[],int code[])
{
 int i=0;int j=0;
 cout<<"The equator variable and Restrictor\n";
 cin>>m>>n;
 for(i=0;i<n;i++)
 {
  cout<<"Inputb[] and Restrictor code 0:<= 1:= 2:>=\n";//输入b[]和限制符号 <= ,= ,>=
  cin>>b>>code;           //分别输入到b 和  code
  cout<<"The 系数  \n";         //提示输入每个限制条件的系数
  for(j=0;j<m;j++)         
   cin>>matrix[j];
 }
 cout<<"the type 0:Min 1:max\n";//输入要求是类型极大还是极小 0:Min 1:max
 do{
  cin>>type;
  if(type!=0&&type!=1)
   cout<<"error,ReInput!\n";
 }while(type!=0&&type!=1);
 cout<<"the Z\n";
 for(i=0;i<m;i++)
  cin>>matrix[n];
 if(type==1)                    //如果是求极大,把它转化为极小来做,系数全部反号
  for(i=0;i<m;i++)
   matrix[n]=-matrix[n];
}                                           


//////////////////   
void Xartificial()//消去人工变量
{
 int i,j,k;
 
 if(indexg!=0)
 {
  for(i=m+indexe+indexl;i<s;i++)
  {
   for(j=0;j<n;j++)
    if(matrix[j]==1)
    {
     for(k=0;k<=s;k++)
      matrix[n][k]=matrix[n][k]-matrix[j][k]*100;
     j=n;
    }
  }
 }
}

////////////////////////////////////////////////
void Process(float c[][100],int row,int vol)
{
 int i;
 for(i=0;i<n;i++)     //i =0 时置第一列为 1  0  0     i=1 时置第2列为 0  1  0 
  if(i!=row)c[vol]=0;
}
//////////////////////
void Start(float b[],int code[])
{
 int i;
 float nget[100][100],nlet[100][100],net[100][100];
 indexe=indexl=indexg=0;               //indexl 表示松弛变量数   indexg 表示人工变量数, indexe表示减去的松弛变量数
 for(i=0;i<n;i++)
 {
  if(code==0)    //如果是<=
  {
   nlet[indexl++]=1;         //松弛变量数+1 且置成相应的标记
   rocess(nlet,i,indexl-1);//传 net, 行号,indexl-1 过去
  }
  if(code==1)
  {
   net[indexg++]=1;           //人工变量数+1且置成相应的标记
   rocess(net,i,indexg-1);      //将刚加入的列单位化
  }
  if(code==2)
  {
   net[indexg++]=1;         //人工变量数+1 且置成相应的标记
   nget[indexe++]=-1;       //剩余变量数+1 且置成相应的标记
            Process(net,i,indexg-1)rocess(nlet,i,indexe-1);
  }
 }
 s=indexe+indexl+indexg+m;   //变量总个数
 Merge(nget,nlet,net,b);
 rocessA();
 InitPrint();
 Xartificial();
}

void Simplix()//单纯形法
{
 int in,out,temp=0;
 while(1)
 {
  jckxj();
  rint();
        Result();
        if(!rj()) in=Min();
  else
  {
   if(indexg!=0)
    JustArtificial();
   rintResult();
   return;
  }
  if(Check(in))
  {
   cout<<"No Delimition\n";
   return;
  }
  out=SearchOut(&temp,in);
  Mto(in,temp);
  Be(temp,in);
  Achange(in,out);
 }
}

void main()
{
 int code[100];//输入符号标记
 float b[100];
 Input(b,code);//初始化
 Start(b,code);//标准化行
 Simplix();
}
/*模拟输入数据
 max z=2 x1 + 3 x2
      s.t {
        x1 + 2 x2 <=8
     4 x1      <=16
           4 x2<=12
            x1,x2>=0
   }

2 3 8 0 1 2 16 0 4 0 12 0 0  4
*/


作者: 999mmg    时间: 2006-11-3 13:55
呵呵&nbsp; 绝对的好东西&nbsp; 我喜欢
作者: hnus31    时间: 2006-11-11 09:46
不错呀
作者: pytff7    时间: 2006-11-18 15:45
hai a&nbsp;
作者: echo5183    时间: 2006-12-4 12:19


作者: jkkjk    时间: 2009-2-14 04:34
这么长!看不懂啊。
作者: aruisi    时间: 2009-4-4 16:31
我眼睛好花啊!
作者: Kadyniost    时间: 2009-8-10 00:26
。。。。。。。。。。。。。。。。




欢迎光临 数学建模社区-数学中国 (http://www.madio.net/) Powered by Discuz! X2.5