线性规划问题解决开源工具(GNU Linear Programming Kit)

GNU Linear Programming Kit (GLPK)一个开源的线性规划工具，用了一下感觉语法还挺简单了（有点像python的感觉，但没python清晰）向大家介绍一下

入门实践

最近在做一个叫交通最小通勤计算问题，需要用到线性规划来解决，因此在网上搜了一下啊线性规划工具，因为不想装MATLAB，(实在是太大了，电脑c盘剩下不到4g了)就找了一个开源的线性规划小工具，感觉还蛮实用的，(GNU Linear Programming Kit, GLPK)[http://gnu.april.org/software/glpk/] 一个开源的线性规划工具，再这里给大家介绍介绍。

GLPK文件目录

glpsol.exe就是主程序了，glpsol.exe主要是通过命令行运行，可以通过 –help 命令了解下他的主要命令：

GLPK所使用的编译语言主要是 GNU MathProg language，我主要尝试了glpsol的两个命令–math 和 –model，分别介绍下：

线性规划方程：

本案列就用Sriram在Coursera公开课的上讲的案例直接进行介绍了，math方法是最简单的方法，就是直接把线性方程写下来，找一个txt记事本：

var x1>=0;
var x2>=0;
maximize obj:x1+2*x2;
c1:-3*x1+x2<=2;
c2:x2<=11;
c3:x1-x2<=3;
c4:x1<=6;
solve;
display x1;
display x2;
end;

可以看出MathProg language很简单，定义变量范围var,定义目标maximize obj:和约束条件就可以了，最后求解solve和显示display 然后保持为first.ampl

在CMD命令行直接输入glpsol –math fitst.ampl就可以了

可以看到结果为 这种方法在解决简单少量的线性规划的时候很简单清晰，但是在解决大量线性规划的时候是不具备可操作性的，因此介绍GLPK的第二种命令--model,这种命令可以用两个文件存储一个为MODEL文件，一个为DATA文件，MODEL文件主要通过构建矩阵进行线性规划计算，同样以上面的线性规划为例，可以得出其实上面的约束方程可以看出两个矩阵相乘，分别为一个系数矩阵A和所求矩阵X相乘小于等于b矩阵(A*x<=b)：

param m; param n;

param c{i in 1..n}; param A{i in 1..m,j in 1..n}; param b{i in 1..m};

var x{i in 1..n}>=0;

maximize obj:sum{i in 1..n} c[i]x[i]; s.t. e{j in 1..m}:sum{i in 1..n} A[j,i]x[i]<=b[j];

solve; display x; end;

写完model文件还需要写一个赋值的data文件对model中的参数赋值：

param n:=2;

param m:=4;

param c:=1 2： 1 2;

param A:=1 2： 1 -3 1 2 0 1 3 1 -1 4 1 0;

param b:= 1： 1 2 2 11 3 3 4 6;

具体data文件的写作格式可以参考我的百度云盘上的gmpl文件（ http://pan.baidu.com/s/1i3CDq8t  ）

里面有详细说明，

然后只要再在cmd中用命令--model执行就可以了：

![GLPK-计算结果](./img/GLPK-计算结果.png)

结果和math命令一样，不过内存使用稍微大了点。


## 过剩通勤应用

> 本文继续对 GNU Linear Programming Kit进行解释，本次介绍引进实际应用，计算最大过剩通勤和最小过剩通勤。

---------------

通过上面对GLPK的建模计算有了大概了解，本章完成BOSS下达任务，完成一个过剩通勤计算。

首先当然是写好model文件,其中最大通勤为：

param n; param m;

param Population{i in 1..n,j in 1..m}; param Distance{i in 1..n,j in 1..m};

var x{i in 1..n,j in 1..m}>=0,integer;

maximize obj:sum{i in 1..n}sum{j in 1..m} Distance[i,j]*x[i,j]; s.t. e{i in 1..n}:sum{j in 1..m}x[i,j]=sum{j in 1..m}Population[i,j]; f{j in 1..m}:sum{i in 1..n}x[i,j]=sum{i in 1..n}Population[i,j];

solve; printf “min sum:%d”,sum{i in 1..n}sum{j in 1..m} Distance[i,j]*x[i,j]; end;

最小通勤为：

param n; param m;

param Population{i in 1..n,j in 1..m}; param Distance{i in 1..n,j in 1..m};

var x{i in 1..n,j in 1..m}>=0,integer;

minimize obj:sum{i in 1..n}sum{j in 1..m} Distance[i,j]*x[i,j]; s.t. e{i in 1..n}:sum{j in 1..m}x[i,j]=sum{j in 1..m}Population[i,j]; f{j in 1..m}:sum{i in 1..n}x[i,j]=sum{i in 1..n}Population[i,j];

solve; printf “min sum:%d”,sum{i in 1..n}sum{j in 1..m} Distance[i,j]*x[i,j]; end;

就是将maximize改为miniminze，写了个run.bat文件方便输出：

glpsol –model MaxTrafficCommuting.model –data Data.data –output max.solve glpsol –model MinTrafficCommuting.model –data Data.data –output min.solve

接下来是data文件部分了，由于原始数据是excel数据，需要先进行格式整理，主要就用pandas进行整理，操作方便，直接附上python代码：（data文件的格式可以参照上一篇GPLK解释的文章）

-- coding: cp936 --

import pandas as pd

def toformat():

#data is big table, and data1 is small table
#make the big equal to the small
data=pd.read_excel('Distance.xlsx',index_col=None)#big
data1=pd.read_excel('Population.xlsx',index_col=None)#small
#change its column
data=data[data1.columns]
if len(data1.columns)!=len(data.columns):
    print "row exist Duplicate items"
#change its row's index
data=data.loc[data1.index]
if len(data1.index)!=len(data.index):
    print "row exist Duplicate items"
data.to_excel('Distance.xlsx',sheet_name='Sheet1',engine='xlsxwriter')
#to creat txt file
data.to_csv('Distance.txt',sep='\t')
data1.to_csv('Population.txt',sep='\t')
pop,popnum=changeformat('Population.txt')
dis,disnum=changeformat('Distance.txt')
if popnum!=disnum:
    print '两个矩阵大小不相等，请检查数据'
print '实际通勤为%d'%(sum(data*data1))
return pop,dis,str(len(data.index)),str(len(data.columns))

def todealmiss(): data=pd.read_excel(‘Population.xlsx’,’Sheet1’,index_col=None,na_values=[‘0’]) data=data.fillna(0) data.to_excel(‘Population.xlsx’,sheet_name=’Sheet1’,engine=’xlsxwriter’)

def changeformat(filename=’Population.txt’): with open(filename,’r’) as datafile: data=datafile.read().split(‘\n’) n=len(data)

data[0]='param '+filename[:-4]+':'+data[0][1:]+':='
#从excel转换为txt时最后一个空格键多出了一行，
data[n-2]=data[n-2]+';'
outputname='new_'+filename
return data,n

todealmiss() pop,dis,row,col=toformat() with open(‘Data.data’,’w’) as datafile: datafile.write(‘param n:=’+row+’;’+’\n’) datafile.write(‘param m:=’+col+’;’+’\n’) for p in pop: datafile.write(p+’\n’) for d in dis: datafile.write(d+’\n’) ```

代码写的很乱，没优化了，直接运行结果，总体来说速度还是挺快的，381*381（145161）个数据大概用了5秒多：