線性規劃之單純形法【超詳解+圖解】

時間 2019-11-19

標籤線性規劃單純超詳解+圖解简体版

原文原文鏈接

1.做用

單純形法是解決線性規劃問題的一個有效的算法。線性規劃就是在一組線性約束條件下，求解目標函數最優解的問題。php

2.線性規劃的通常形式

在約束條件下，尋找目標函數z的最大值。ios

$\begin{array}{l} \max z = {x_1} + {x_2}\\ s.t.\left\{ {\begin{array}{*{20}{c}} {2{x_1} + {x_2} \le 12}\\ {{x_1} + 2{x_2} \le 9}\\ {{x_1},{x_2} \ge 0} \end{array}} \right. \end{array}$

3.線性規劃的可行域

知足線性規劃問題約束條件的全部點組成的集合就是線性規劃的可行域。若可行域有界（如下主要考慮有界可行域），線性規劃問題的目標函數最優解必然在可行域的頂點上達到最優。c++

單純形法就是經過設置不一樣的基向量，通過矩陣的線性變換，求得基可行解（可行域頂點），並判斷該解是否最優，不然繼續設置另外一組基向量，重複執行以上步驟，直到找到最優解。因此，單純形法的求解過程是一個循環迭代的過程。算法

圖1 可行域數據結構

4.線性規劃的標準形式

在說明單純形法的原理以前，須要明白線性規劃的標準形式。由於單純形算法是經過線性規劃的標準形來求解的。通常，規定線性規劃的標準形式爲：函數

$\begin{array}{l} \max z = \sum\limits_{j = 1}^n {{c_j}{x_j}} \\ s.t.\left\{ {\begin{array}{*{20}{c}} {\sum\limits_{j = 1}^n {{a_{ij}}{x_j} = {b_j}, i = 1,2,...,m} }\\ {xj \ge 0, j = 1,2,...,n} \end{array}} \right. \end{array}$

寫成矩陣形式：學習

$\begin{array}{l} \max z = CX\\ AX = b\\ X \ge 0\\ A = \left[ {\begin{array}{*{20}{c}} {{a_{11}}}&{{a_{12}}}& \cdots &{{a_{1n}}}\\ \vdots & \vdots & \ddots & \vdots \\ {{a_{m1}}}&{{a_{m2}}}& \cdots &{{a_{mn}}} \end{array}} \right] \end{array}$

標準形的形式爲：優化

1）目標函數要求maxspa

2）約束條件均爲等式.net

3）決策變量爲非負約束

普通線性規劃化爲標準形：

1）若目標函數爲最小化，能夠經過取負，求最大化

2）約束不等式爲小於等於不等式，能夠在左端加入非負鬆弛變量，轉變爲等式，好比：

${x_1} + 2{x_2} \le 9 \Rightarrow \left\{ {\begin{array}{*{20}{c}} {{x_1} + 2{x_2} + {x_3} = 9}\\ {{x_3} \ge 0} \end{array}} \right.$

同理，約束不等式爲大於等於不等式時，能夠在左端減去一個非負鬆弛變量，變爲等式。

3）若存在取值無約束的變量，可轉變爲兩個非負變量的差，好比：

$- \infty \le {x_k} \le + \infty \Rightarrow \left\{ {\begin{array}{*{20}{c}} {{x_k} = {x_m} - {x_n}}\\ {{x_m},{x_n} \ge 0} \end{array}} \right.$

本文最開始的線性規劃問題轉化爲標準形爲：

$\begin{array}{l} \max z = {x_1} + {x_2}\\ s.t.\left\{ {\begin{array}{*{20}{c}} {2{x_1} + {x_2} + {x_3} = 12}\\ {{x_1} + 2{x_2} + {x_4} = 9}\\ {{x_1},{x_2},{x_3},{x_4} \ge 0} \end{array}} \right. \end{array}$

5.單純形法

5.1幾何意義

在標準形中，有m個約束條件（不包括非負約束），n個決策變量，且（n>=m）。首先，選取m個基變量，基變量對應約束係數矩陣的列向量線性無關。經過矩陣的線性變換，基變量可由非基變量表示：

$x_i' = {C_i} + \sum\limits_{j = m + 1}^n {{m_{ij}}x_j'} (i = 1,2,...,m)$

若是令非基變量等於0，可求得基變量的值：

$x_i' = {C_i}$

若是爲可行解的話，Ci大於0。那麼它的幾何意義是什麼呢？

仍是經過上述具體的線性規劃問題來講明。

若是選擇x二、x3爲基變量，那麼令x一、x4等於0，能夠去求解基變量x二、x3的值。對係數矩陣作行變換，以下所示，x2=9/2，x3=15/2

$\left[ {\begin{array}{*{20}{c}} {\rm{X}}&{{x_1}}&{{x_2}}&{{x_3}}&{{x_4}}&b\\ {}&2&1&1&0&{12}\\ {}&1&2&0&1&9\\ {\rm{C}}&1&1&0&0&z \end{array}} \right] \to \left[ {\begin{array}{*{20}{c}} {\rm{X}}&{{x_1}}&{{x_2}}&{{x_3}}&{{x_4}}&b\\ {}&{{\textstyle{3 \over 2}}}&0&1&{ - {\textstyle{1 \over 2}}}&{{\textstyle{{15} \over 2}}}\\ {}&{{\textstyle{1 \over 2}}}&1&0&{{\textstyle{1 \over 2}}}&{{\textstyle{9 \over 2}}}\\ {\rm{C}}&{{\textstyle{1 \over 2}}}&0&0&{ - {\textstyle{1 \over 2}}}&{z - {\textstyle{9 \over 2}}} \end{array}} \right]$

X1=0表示可行解在x軸上；X4=0表示可行解在x1+2x2=9的直線上。那麼，求得的可行解即表示這兩條直線的交點，也是可行域的頂點，如圖所示：

圖2

因此，經過選擇不一樣的基變量，能夠得到不一樣的可行域的頂點。

5.2如何判斷最優

如前所述，基變量可由非基變量表示：

$x_i' = {C_i} + \sum\limits_{j = m + 1}^n {{m_{ij}}x_j'} (i = 1,2,...,m)$

目標函數z也能夠徹底由非基變量表示：

$z = {z_0} + \sum\limits_{j = m + 1}^n {{\sigma _j}x_j'}$

當達到最優解時，全部的 ${{\sigma _j}}$ 應小於等於0。當存在j， ${{\sigma _j}}$ >0時，當前解不是最優解，爲何？

當前的目標函數值爲z0，其中全部的非基變量值均取0。由以前分析可知， ${x_j'}$ =0表明可行域的某個邊界，是 ${x_j'}$ 的最小值。若是可行解逐步離開這個邊界， ${x_j'}$ 會變大，由於 ${{\sigma _j}}$ >0，顯然目標函數的取值也會變大，因此當前解不是最優解。咱們須要尋找新的基變量。

5.3如何選擇新的基變量

若是存在多個 ${{\sigma _j}}$ >0，選擇最大的 ${{\sigma _j}}$ >0對應的變量做爲基變量，這表示目標函數隨着 ${x_j'}$ 的增長，增加的最快。

5.4如何選擇被替換的基變量

假如咱們選擇非基變量做爲下一輪的基變量，那麼被替換基變量 ${x_j'}$ 在下一輪中做爲非基變量，等於0。選擇 ${x_j'}$ 的原則：替換後應該儘可能使值最大（由於上面已分析過，目標函數會隨着的增大而增大）。

繼續經過上面的例子來講明：

從最後一行能夠看到，x1的係數爲1/2>0，因此選x二、x3爲基變量並無是目標函數達到最優。下一輪選取x1做爲基變量，替換x二、x3中的某個變量。

第一行是符號

第二行：若x1替換x3做爲基變量，x3=0時，x1=(15/2)/(3/2)=5

第三行：若x1替換x2做爲基變量，x2=0時，x1=(9/2)/(1/2)=9

顯然，應該把x2做爲非基變量。

5.5終止條件

當目標函數用非基變量的線性組合表示時，全部的係數均不大於0，則表示目標函數達到最優。

若是，有一個非基變量的係數爲0，其餘的均小於0，表示目標函數的最優解有無窮多個。這是由於目標函數的梯度與某一邊界正交，在這個邊界上，目標函數的取值均相等，且爲最優。

使用單純型法來求解線性規劃，輸入單純型法的鬆弛形式，是一個大矩陣，第一行爲目標函數的係數，且最後一個數字爲當前軸值下的 z 值。下面每一行表明一個約束，數字表明係數每行最後一個數字表明 b 值。

算法和使用單純性表求解線性規劃相同。

對於線性規劃問題：

Max x1 + 14* x2 + 6*x3

s . t .　 x1 + x2 + x3 <= 4

　　　　x1<= 2

　　　　x3 <= 3

　　　　3*x2 + x3 <= 6

　　　　x1,x2,x3 >= 0

咱們能夠獲得其鬆弛形式：

Max　　x1 + 　14*x2 + 6*x3
s.t.　　 x1 + 　x2 　　+ x3 　　+ x4 = 4
　　　　x1 　　　　　　　　　　　　　　+ x5 = 2
　　　　　　　　　　　　x3 　　　　　　　　　　+ x6 = 3
　　　　　　　 3*x2 + x3 　　　　　　　　　　　　　　+ x7 = 6
　　　　x1 , x2 , x3 , x4 , x5 , x6 , x7 ≥ 0

咱們能夠構造單純性表，其中最後一行打星的列爲軸值。

單純性表
x1	x2	x3	x4	x5	x6	x7	b
c1=1	c2=14	c3=6	c4=0	c5=0	c6=0	c7=0	-z=0
1	1	1	1	0	0	0	4
1	0	0	0	1	0	0	2
0	0	1	0	0	1	0	3
0	3	1	0	0	0	1	6
			*	*	*	*

在單純性表中，咱們發現非軸值的x上的係數大於零，所以能夠經過增長這些個x的值，來使目標函數增長。咱們能夠貪心的選擇最大的c，再上面的例子中咱們選擇c2做爲新的軸，加入軸集合中，那麼誰該出軸呢？

其實咱們因爲每一個x都大於零，對於x2它的增長是有所限制的，若是x2過大，因爲其餘的限制條件，就會使得其餘的x小於零，因而咱們應該讓x2一直增大，直到有一個其餘的x恰好等於0爲止，那麼這個x就被換出軸。

咱們能夠發現，對於約束方程1，即第一行約束，x2最大能夠爲4（4/1），對於約束方程4，x2最大能夠爲3(6/3），所以x2最大隻能爲他們之間最小的那個，這樣才能保證每一個x都大於零。所以使用第4行，來對各行進行高斯行變換，使得二列第四行中的每一個x都變成零，也包括c2。這樣咱們就完成了把x2入軸，x7出軸的過程。變換後的單純性表爲：

單純性表
x1	x2	x3	x4	x5	x6	x7	b
c1=1	c2=0	c3=1.33	c4=0	c5=0	c6=0	c7=-4.67	-z=-28
1	0	0.67	1	0	0	-0.33	2
1	0	0	0	1	0	0	2
0	0	1	0	0	1	0	3
0	1	0.33	0	0	0	0.33	2
	*		*	*	*

繼續計算，咱們獲得：

單純性表
x1	x2	x3	x4	x5	x6	x7	b
c1=-1	c2=0	c3=0	c4=0	c5=-2	c6=0	c7=0	-z=-32
1.5	0	1	1.5	0	0	-0.5	3
1	0	0	0	1	0	0	2
0	0	1	0	0	1	0	3
0	1	0.33	0	0	0	0.33	2
	*		*	*	*

此時咱們發現，全部非軸的x的係數所有小於零，即增大任何非軸的x值並不能使得目標函數最大，從而獲得最優解32.

整個過程代碼以下所示：

  1 #include <bits/stdc++.h>
  2 using namespace std;
  3 vector<vector<double> > Matrix;
  4 double Z;
  5 set<int> P;
  6 size_t cn, bn;
  7 
  8 bool Pivot(pair<size_t, size_t> &p)//返回0表示全部的非軸元素都小於0
  9 {
 10     int x = 0, y = 0;
 11     double cmax = -INT_MAX;
 12     vector<double> C = Matrix[0];
 13     vector<double> B;
 14 
 15     for( size_t i = 0 ; i < bn ; i++ )
 16     {
 17         B.push_back(Matrix[i][cn-1]);
 18     }
 19 
 20     for( size_t i = 0 ; i < C.size(); i++ )//在非軸元素中找最大的c
 21     {
 22         if( cmax < C[i] && P.find(i) == P.end())
 23         {
 24             cmax = C[i];
 25             y = i;
 26         }
 27     }
 28     if( cmax < 0 )
 29     {
 30         return 0;
 31     }
 32 
 33     double bmin = INT_MAX;
 34     for( size_t i = 1 ; i < bn ; i++ )
 35     {
 36         double tmp = B[i]/Matrix[i][y];
 37        if( Matrix[i][y] != 0 && bmin > tmp )
 38        {
 39            bmin = tmp;
 40            x = i;
 41        }
 42     }
 43 
 44     p = make_pair(x, y);
 45 
 46     for( set<int>::iterator it = P.begin() ; it != P.end() ; it++)
 47     {
 48         if( Matrix[x][*it] != 0 )
 49         {
 50             //cout<<"erase "<<*it<<endl;
 51             P.erase(*it);
 52             break;
 53         }
 54     }
 55     P.insert(y);
 56     //cout<<"add "<<y<<endl;
 57     return true;
 58 }
 59 
 60 void pnt()
 61 {
 62     for( size_t i = 0 ; i < Matrix.size() ; i++ )
 63     {
 64         for( size_t j = 0 ; j < Matrix[0].size() ; j++ )
 65         {
 66             cout<<Matrix[i][j]<<"\t";
 67         }
 68     cout<<endl;
 69     }
 70     cout<<"result z:"<<-Matrix[0][cn-1]<<endl;
 71 }
 72 
 73 void Gaussian(pair<size_t, size_t> p)//行變換
 74 {
 75     size_t  x = p.first;
 76     size_t y = p.second;
 77     double norm = Matrix[x][y];
 78     for( size_t i = 0 ; i < cn ; i++ )//主行歸一化
 79     {
 80         Matrix[x][i] /= norm;
 81     }
 82     for( size_t i = 0 ; i < bn && i != x; i++ )
 83     {
 84         if( Matrix[i][y] != 0)
 85         {
 86             double tmpnorm = Matrix[i][y];
 87             for( size_t j = 0 ; j < cn ; j++ )
 88             {
 89                 Matrix[i][j] = Matrix[i][j] - tmpnorm * Matrix[x][j];
 90             }
 91         }
 92     }
 93 }
 94 
 95 void solve()
 96 {
 97     pair<size_t, size_t> t;
 98     while(1)
 99     {
100 
101         pnt();
102         if( Pivot(t) == 0 )
103         {
104             return;
105         }
106         cout<<t.first<<" "<<t.second<<endl;
107         for( set<int>::iterator it = P.begin(); it != P.end()  ; it++ )
108         {
109             cout<<*it<<" ";
110         }
111         cout<<endl;
112         Gaussian(t);
113     }
114 }
115 
116 int main(int argc, char *argv[])
117 {
118     //ifstream fin;
119     //fin.open("./test");
120     cin>>cn>>bn;
121     for( size_t i = 0 ; i < bn ; i++ )
122     {
123         vector<double> vectmp;
124         for( size_t j = 0 ; j < cn ; j++)
125         {
126             double tmp = 0;
127             cin>>tmp;
128             vectmp.push_back(tmp);
129         }
130         Matrix.push_back(vectmp);
131     }
132 
133     for( size_t i = 0 ; i < bn-1 ; i++ )
134     {
135         P.insert(cn-i-2);
136     }
137     solve();
138 }
139 /////////////////////////////////////
140 //glpk input:
141 ///* Variables */
142 //var x1 >= 0;
143 //var x2 >= 0;
144 //var x3 >= 0;
145 ///* Object function */
146 //maximize z: x1 + 14*x2 + 6*x3;
147 ///* Constrains */
148 //s.t. con1: x1 + x2 + x3 <= 4;
149 //s.t. con2: x1  <= 2;
150 //s.t. con3: x3  <= 3;
151 //s.t. con4: 3*x2 + x3  <= 6;
152 //end;
153 /////////////////////////////////////
154 //myinput:
155 /*
156 8 5
157 1 14 6 0 0 0 0 0
158 1 1 1 1 0 0 0 4
159 1 0 0 0 1 0 0 2
160 0 0 1 0 0 1 0 3
161 0 3 1 0 0 0 1 6
162 */
163 /////////////////////////////////////

【理論羅列】：

1.標準型

m個約束 n個變量用x向量表示 A是一個m*n的矩陣 c是一個n的向量 b是一個m的向量

最大化 cx

知足約束 Ax<=b x>0

2.鬆弛型

基本變量 B |B|=m 一個約束對應一個表示鬆弛量叫作鬆弛變量(基本變量)

非基變量 N |N|=n

xn+i=bi-sigma{aijxj}>=0

3.替入變量 xe(非基變量)

替出變量 xl(基本變量)

4.可行解

　基本解：全部非基變量設爲0

　基本可行解

5.單純形法的過程當中B和N不斷交換，在n維空間中不斷走

「至關於不等式上的高斯消元」

【代碼實現】：

pivot是轉動操做

基本思想就是改寫l這個約束爲xe做爲基本變量，而後把這個新xe的值帶到其餘約束和目標函數中，就消去xe了

改寫和帶入時要修改b和a 目標函數則是 c和v

轉動時l和e並無像算法導論上同樣a矩陣用了兩行分別是a[l][]和a[e][]（這樣佔用內存大），而是用了同一行，這樣a矩陣的行數=|B|，列數=|N|

也就是說，約束條件只用m個，儘管B和N不斷交換，但同一時間仍是隻有m個約束(基本變量)n個非基變量

注意改寫成鬆弛型後a矩陣實際係數爲負

（一個優化 a[i][e]爲0的約束不必帶入了

simplex是主過程

基本思想是找到一個c[e]>0的，而後找對這個e限制最緊的l，轉動這組l e

注意精度控制eps

c[e]>eps

還有找l的時候a[i][e]>eps才行

【對偶原理】：

1.原始線性規劃對偶線性規劃

2.對於

最大化 cx

知足約束 Ax<=b x>0

對偶問題爲

最小化 bx

知足約束 ATx>=c x>0 （AT爲A的轉置）

能夠轉化不少問題來避免初始解不可行

我來秀智商了……

說從前有個線性規劃
min c x^T
Ax = b
x >= 0
這裏面A是矩陣，x、b、c都是向量
x^T表示轉置

啊……咱們假設以上出現的全部元素都是整數……

那麼Ax = b要是剛好方程數等於未知數數，且解出來剛好爲正數是否是就超開心？（假設是線性無關的）
根據那個啥法則，x_i = det(A_i) / det(A)
det(A)表示A的行列式
A_i表示把A的第i列換爲b以後的矩陣
若是det(A_i)剛好是det(A)的倍數那不就超開心？這樣
可是現實是殘酷的，每每這傢伙會除出來小數，而後整數規劃就坑爹了。

可是人類的智慧是無窮的！
咱們如今仍是假定「剛好方程數等於未知數數，且解出來剛好爲正數」
咱們再加個限制：det(A) = 1或-1
就233了吧。
解必定是整數了。
因而能夠順利變爲整數規劃。咱們把det(A) = 1, -1的矩陣稱爲幺模矩陣。

可是現實是殘酷的，「剛好方程數等於未知數數，且解出來剛好爲正數」每每不知足。
可是其實不要緊。因爲每一個方程都對應着一個平面，因此解的空間是單純形，最優解必定會出如今頂點上。
何爲頂點？就是平面的交點。
何爲平面？一共m + n個：Ax = b是m個方程，x = 0是n個方程。（原本是x >= 0，咱們只靠慮切割空間的平面……）
要是頂點都是整點不是超開心？等價於從這m + n個方程中任取n個方程把它解掉獲得的解是整數解。
經過前面的分析咱們知道，若是剛好選取的這n個方程的係數矩陣的行列式值爲1，-1就超開心了。固然要是行列式值爲0對應着無解或無窮多解的狀況，它又不是頂點管它作甚……
考察係數矩陣
一個是A，好大好大
另外一個是x = 0的係數，易知就是單位矩陣I
你從I中選哪幾行……因爲行列式的性質……一行*k加到另外一行上去行列式的值不變，那麼對應的未知數就會被幹掉……
因此等價於……
從A中任意選取是一個子方陣，它的行列式的值只可能爲1, -1, 0。
這樣的矩陣咱們稱爲全幺模矩陣。

番外篇：

1. 必要不充分：只含1，-1，0。由於單個元素能夠看做行列式……
2. 充分必要：對它進行高斯消元的主元操做……（好像叫轉軸？啊反正就是消別人的未知數……），得來的仍是全幺模矩陣……這個是由於那個啥啥幺模矩陣組成了一個乘法羣？用這個性質咱們能夠不用double。
3. 您能夠手工屠矩陣用定義證它是全幺模！
4. 若是數學太差，您也能夠寫一個O(4^n * n^3)的強程序證它是全幺模！

orzorzorz

附上一道題BZOJ 1061

 1 #include<iostream>
 2 #include<cstdio>
 3 #include<cstring>
 4 #include<algorithm>
 5 #include<cmath>
 6 using namespace std;
 7 typedef long long ll;
 8 const int M=10005,N=1005,INF=1e9;
 9 const double eps=1e-6;
10 inline int read(){
11     char c=getchar();int x=0,f=1;
12     while(c<'0'||c>'9'){if(c=='-')f=-1; c=getchar();}
13     while(c>='0'&&c<='9'){x=x*10+c-'0'; c=getchar();}
14     return x*f;
15 }
16 
17 int n,m;
18 double a[M][N],b[M],c[N],v;
19 void pivot(int l,int e){
20     b[l]/=a[l][e];
21     for(int j=1;j<=n;j++) if(j!=e) a[l][j]/=a[l][e];
22     a[l][e]=1/a[l][e];
23     
24     for(int i=1;i<=m;i++) if(i!=l&&fabs(a[i][e])>0){
25         b[i]-=a[i][e]*b[l];
26         for(int j=1;j<=n;j++) if(j!=e) a[i][j]-=a[i][e]*a[l][j];
27         a[i][e]=-a[i][e]*a[l][e];
28     }
29     
30     v+=c[e]*b[l];
31     for(int j=1;j<=n;j++) if(j!=e) c[j]-=c[e]*a[l][j];
32     c[e]=-c[e]*a[l][e];
33     
34     //swap(B[l],N[e])
35 }
36 
37 double simplex(){
38     while(true){
39         int e=0,l=0;
40         for(e=1;e<=n;e++) if(c[e]>eps) break;
41         if(e==n+1) return v;
42         double mn=INF;
43         for(int i=1;i<=m;i++)
44             if(a[i][e]>eps&&mn>b[i]/a[i][e]) mn=b[i]/a[i][e],l=i;
45         if(mn==INF) return INF;//unbounded
46         pivot(l,e);
47     }
48 }
49 
50 int main(){
51     n=read();m=read();
52     for(int i=1;i<=n;i++) c[i]=read();
53     for(int i=1;i<=m;i++){
54         int s=read(),t=read();
55         for(int j=s;j<=t;j++) a[i][j]=1;
56         b[i]=read();
57     }
58     printf("%d",(int)(simplex()+0.5));
59 }

附上幾道題的題號，練習學習一下：

BZOJ 3550

BZOJ 3112

BZOJ 3265

BZOJ 1061

POJ 1275

標準型:

轉化爲標準型:

鬆弛型:

單純型算法

轉軸：

初始化：

代碼實現：

UOJ#179. 線性規劃

 1 #include<iostream>
 2 #include<cstdio>
 3 #include<cstring>
 4 #include<algorithm>
 5 #include<cmath>
 6 using namespace std;
 7 typedef long long ll;
 8 const int N=25;
 9 const double eps=1e-8,INF=1e15;
10 inline int read(){
11     char c=getchar();int x=0,f=1;
12     while(c<'0'||c>'9'){if(c=='-')f=-1; c=getchar();}
13     while(c>='0'&&c<='9'){x=x*10+c-'0'; c=getchar();}
14     return x*f;
15 }
16 int n,m,type;
17 double a[N][N],ans[N];
18 int id[N<<1];
19 int q[N];
20 void Pivot(int l,int e){
21     swap(id[n+l],id[e]);
22     double t=a[l][e]; a[l][e]=1;
23     for(int j=0;j<=n;j++) a[l][j]/=t;
24     for(int i=0;i<=m;i++) if(i!=l && abs(a[i][e])>eps){
25         t=a[i][e]; a[i][e]=0;
26         for(int j=0;j<=n;j++) a[i][j]-=a[l][j]*t;
27     }
28 }
29 bool init(){
30     while(true){
31         int e=0,l=0;
32         for(int i=1;i<=m;i++) if(a[i][0]<-eps && (!l||(rand()&1))) l=i;
33         if(!l) break;
34         for(int j=1;j<=n;j++) if(a[l][j]<-eps && (!e||(rand()&1))) e=j;
35         if(!e) {puts("Infeasible");return false;}
36         Pivot(l,e);
37     }
38     return true;
39 }
40 bool simplex(){
41     while(true){
42         int l=0,e=0; double mn=INF;
43         for(int j=1;j<=n;j++) 
44             if(a[0][j]>eps) {e=j;break;}
45         if(!e) break;
46         for(int i=1;i<=m;i++) if(a[i][e]>eps && a[i][0]/a[i][e]<mn)
47             mn=a[i][0]/a[i][e],l=i;
48         if(!l) {puts("Unbounded");return false;}
49         Pivot(l,e);
50     }
51     return true;
52 }
53 int main(){
54     freopen("in","r",stdin);
55     srand(317);
56     n=read();m=read();type=read();
57     for(int i=1;i<=n;i++) a[0][i]=read();
58     for(int i=1;i<=m;i++){
59         for(int j=1;j<=n;j++) a[i][j]=read();
60         a[i][0]=read();
61     }
62     for(int i=1;i<=n;i++) id[i]=i;
63     if(init() && simplex()){
64         printf("%.8lf\n",-a[0][0]);
65         if(type){
66             for(int i=1;i<=m;i++) ans[id[n+i]]=a[i][0];
67             for(int i=1;i<=n;i++) printf("%.8lf ",ans[i]);
68         }
69     }
70 }