线性规划问题解决开源工具(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命令界面

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

线性规划方程:

GLPK-线性规划公式

本案列就用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-计算结果

可以看到结果为


这种方法在解决简单少量的线性规划的时候很简单清晰,但是在解决大量线性规划的时候是不具备可操作性的,因此介绍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秒多:

最终计算结果1
最终计算结果2

shikanon wechat
欢迎您扫一扫,订阅我滴↑↑↑的微信公众号!