↑↑↑↑↑点击上方蓝色字关注我们！

『运筹OR帷幄』原创

作者：路军

编者按

混合整数规划(MIP)，是指同时包含整数变量和连续变量的数学规划。针对二维剪板机下料问题等大规模组合问题，混合整数规划是一种速度快，精度高的求解途径，上较其他方法具有明显优势。本文介绍了二维剪板机下料问题（2D-GCSP）的混合整数规划模型并基于BWMV bin packing数据库进行模型的实现，生成排样图，为数据库中的案例提出解决方案。

二维剪板机下料问题（2D-GCSP）的混合整数规划是最优美的整数规划模型之一。以往很多人认为像2D-GCSP这样的问题由于本质上的递归性，很难建立成混合整数规划模型，但是持这种观点的人忽略了数学规划的智能特征。所谓计算智能是如果将问题描述给机器，则机器可以自行得到问题的结果。例如很多元启发算法都被称为计算智能算法，就是因为他们好像能够通用性地解决一些问题。这也给我们一个启示，就是不怕问题刁钻，只要我们把问题的特征用建模逻辑描述出来，那么它就是合理的数学规划模型，就能为我们工作。

2D-GCSP二维剪板机下料问题

2D-GCSP中的G来自法文Guillotine。约瑟夫·吉约坦（ Joseph Guillotin ）是法国一名信仰人文主义的医生，他提出认为死刑应以快速及人道的方式执行以降低死刑的痛苦。1792年，吉约坦医生与一位法国刽子手合作，动手制造了一台断头台的原型机，从此新型断头台便正式面世，并以他的姓氏Guillotine命名。为了纪念这位医生，以后凡是瞬间切断东西的方式都被命名为Guillotine。在航母、飞机、建筑等制造工业上经常会使用一种剪板机对金属板材进行切割下料，每次必须横向或纵向切断整个板材，下料图案十分独特，无论多少次切割，就是你总能找到一刀切位置。下图是剪板机。

问题数据

附件一是BWMV bin packing数据库中第40个例子。第一行表示例子编号，第二行30，30表示原始板材大小，80表示总计下料需求板材的种类数。第三行开始是需求定义，例如第三行表示需要5X8的小板材1件，第四行表示需要8X8的小板材1件。所有板材必须从WXH=30X30的大板材上G方式剪裁下来，目标是规划下料方式以极小化所需要的大板材数量。对原始数据，用以下符号表示：W-原始板材宽度，H-原始板材高度，n-需求小板材数量，ai, bi-第i类小板材宽度和高度，mi-第i类小板材的需求数量。

整数规划模型

Silva等[1]首次给出G-下料的整数规划模型，此模型描述真实下料时的逻辑。如下图所示，工人在下料时候要依照给出的需求下料，那么下料后的余料被放起来准备之后的小板材下料。于是除了原始板材外，随着下料的进展，待切割板材会堆积成板材库。那么逻辑是要保证某类待切割的板材数量大于从该类下料小板材的数量。

用xij表示第i种需求小板材从第j种余料（包括原始板材）中下料的数量。用0-1数值矩阵aijk表示从j类余料中下料i类小板擦会产生一件余料k，于是根据逻辑（即工人始终能够找到需要的余料，或者说余料的数量大于将要在该类余料进行下料的数量）：

这里的m是所有的余料种类，j=0表示原始板材（即W×H板材）。

另外，下料需求需要被满足：

目标当然是原始板材使用数量最少：

下图是根本逻辑：从余料j（j=0表示原始板材）上下料i会产生两块新生的余料备用。生成何种余料用三维0-1常数矩阵aijk。约束(1)表示一定要保证有足够的余料可用于需求小板材的下料。

至于如何生成aijk, Silva等[1]建议首先做枚举。本文以下直接将这个枚举写进模型里！

模型实现

上面的模型(1)-(3)比较抽象，现在把它具体化，并把对aijk的枚举写进模型。

首先，对余料从0到M的编码抽象性太强，这里把他们具体化为板材的长度和宽度，即使用x[w][h][i]表示从大小为wXh余料上切割第i类小板材的数量。于是目标是极小化从WxH板材（原始板材）上切割出所有小板材的数量，即：

也就是使用(4)代替上面的(3)。

其次必须满足所有需求：

也就是用(5)取代上面的(2)。

最后用下面的(6)取代上面的(1)：

这个看起来复杂，仔细分析会发现跟上面的(1)是一丘之貉，不过好处是扔掉了aijk。

整体模型看起来是这个样子的（数据区采用小案例，节省空间）：

min sum{i=1,...,n}x[W][H][i] //(4)subject to sum{w=1,...,W;h=1,...,H;w>=a[i];h>=b[i]}x[w][h][i] >=m[i]|i=1,...,n //(5) sum{w=1,...,W;h=1,...,H;i=1,...,n;w>=a[i];h>=b[i];w==(a[i]+w1);b[i]==h1}x[w][h][i]+--> sum{w=1,...,W;h=1,...,H;i=1,...,n;w>=a[i];h>=b[i];h==(b[i]+h1);w1==w}x[w][h][i]--> >=sum{i=1,...,n;w1>=a[i];h1>=b[i]}x[w1][h1][i]--> | w1=1,...,W;h1=1,...H;(h1+w1)<(H+W) //(6)where n is an integer W,H are numbers dt is a set of number a[i],b[i],m[i] are numbers|i=1,...,n x[w][h][i] is a variable of nonnegative integer|w=1,...,W;h=1,...,H;i=1,...,ndata_relation n=_$(dt)/3 a[i]=dt[3(i-1)+1] | i=1,...,n b[i]=dt[3(i-1)+2] | i=1,...,n m[i]=dt[3i] | i=1,...,ndata W=10 H=10 dt={2 16 110 11 18 13 11 11 16 11 14 19 11 19 14 12 13 19 14 19 1 }

上面的模型用+Leapms解析并用cplex命令调用cplex求解即可秒解。

下料排样图的生成

上面的模型只给出x[w][h][i]的最优值，即从wxh余料中切割第i类小板材的具体数量，由此信息得到排样图是另外一个挑战。首先看结果样子（以上面80小板材数据为例）：

而后看技术。这里技术包括用+Leampms编程接口，还有细致的逻辑思考，这里不再赘述，只贴出c++代码，此源码能够一次求出文件中500个案例：

#include<iostream>#include<fstream>using namespace std;extern int W,H,n;extern int w[1000],h[1000],m[1000];extern ofstream ofcad;extern int nInstance;class c_sheet{public: int w,h; //尺寸 int n,nl; //数量，已定位的数量 int x[1000],y[1000]; //坐标 int item[1000]; //下料 bool good; c_sheet(){ w=W; h=H; n=nl=0; } c_sheet(int w,int h){ this->w=w; this->h=h; n=0; good=false; } void reset(){ w=h=n=nl=0; good=false; } void addItem(int item); void print(); void draw(ofstream &off); bool cut(); bool cut(int x, int y);};int prepare(int instance=0);//准备模型int findSheet(int w, int h);extern int nSheet;extern c_sheet sheet[1000];#include<iostream>#include<fstream>using namespace std;extern int W,H,n;extern int w[1000],h[1000],m[1000];extern ofstream ofcad;extern int nInstance;class c_sheet{public: int w,h; //尺寸 int n,nl; //数量，已定位的数量 int x[1000],y[1000]; //坐标 int item[1000]; //下料 bool good; c_sheet(){ w=W; h=H; n=nl=0; } c_sheet(int w,int h){ this->w=w; this->h=h; n=0; good=false; } void reset(){ w=h=n=nl=0; good=false; } void addItem(int item); void print(); void draw(ofstream &off); bool cut(); bool cut(int x, int y);};int prepare(int instance=0);//准备模型int findSheet(int w, int h);extern int nSheet;extern c_sheet s

#include<iostream>#include<fstream>#include "leap.h"#include "Guillotine.h"using namespace std;int W,H,n;int w[1000],h[1000],m[1000];ofstream ofcad;void addSheet(int w, int h, int item, int ns);int dealInstance();//主函数，循环计算所有500案例int main(){ while(1){ if(dealInstance()<0)break; } return 0;}int dealInstance(){ //设置时间 char time[]="600"; //准备模型 if(prepare()<0)return -1; //if(nInstance<41)return 0; //求解模型 c_leap lp; lp.loadLeap("current.leap"); lp.cplex("",time,""); //生成板材 int ns; for(int w=W;w>0;w--){ for(int h=H;h>0;h--){ for(int i=0;i<n;i++){ ns=(int)lp.getVar("x",3,w,h,i+1); if(ns>0){ addSheet(w,h,i,ns); } } } } //剪切板材 if(nSheet>1)sheet[0].cut(); //输出板材图形 char cadfname[64];

 if(lp.getState()!=IP_OPTIMAL)sprintf_s(cadfname,64,"graphic(%d)(%d_%d)!.scr",nInstance,(int)((nInstance-1)/50)+1,nInstance-(int)((nInstance-1)/50)*50);

 else sprintf_s(cadfname,64,"graphic(%d)(%d_%d).scr",nInstance,(int)((nInstance-1)/50)+1,nInstance-(int)((nInstance-1)/50)*50);

ofcad.open(cadfname); ofcad<<"erase all "<<endl; for(int i=0;i<nSheet;i++){ sheet[i].draw(ofcad); } ofcad<<"zoom e"<<endl; ofcad.close(); //清空sheet for(int i=0;i<nSheet;i++){ sheet[i].reset(); } return 0;}void addSheet(int w, int h, int item, int ns){ int sht=findSheet(w,h); if(sht<0){ sht=nSheet++; sheet[sht].w=w; sheet[sht].h=h; } for(int k=0;k<ns;k++){ sheet[sht] }}

#include<iostream>#include<fstream>#include "leap.h"#include "Guillotine.h"using namespace std;ifstream iffd;//原始数据文件句柄bool iffdOpen=false;int nInstance=0;//案例编号int prepare(int instance){ //准备模型 if(nInstance>500)return -1; //如果没有打开原始数据文件，打开之====================== const char iffdname[]="BWMV bin packing instances.txt"; if(!iffdOpen){ int ntmp; iffd.open(iffdname); iffd>>ntmp; cout<<"ntmp="<<ntmp<<endl; iffdOpen=true; } if(!iffd){ cout<<"Can not open the file"<<iffdname<<endl; exit(0); } //打开模板文件并复制==================================== ifstream iffm; ofstream off; iffm.open("Guillotine.leap"); if(!iffm){ cout<<"\tCannot open file 'Guillotine.leap'"<<endl; return 0; } off.open("current.leap"); if(!off){ cout<<"\tCannot create file 'current.txt'"<<endl; return 0; } char line[4096]; while(!iffm.eof()){ iffm.getline(line,4096); off<<line<<endl; } //生成Leap模型数据区===================================== nInstance++; int dummy; iffd>>dummy; iffd>>W>>H>>n; off<<"\tW="<<W<<endl; off<<"\tH="<<H<<endl; off<<"\tdt={"<<endl; for(int i=0;i<n;i++){ iffd>>w[i]>>h[i]>>m[i]; off<<"\t\t"<<w[i]<<"\t"<<h[i]<<"\t"<<m[i]<<endl; } off<<"\t}"<<endl; iffm.close();off.close(); //===================================================== nSheet=0; return nInstance;}#include<iostream>#include<fstream>#include "Guillotine.h"int nSheet=0;c_sheet sheet[1000];int findSheet(int w,int h){ for(int i=0;i<nSheet;i++){ if(w==sheet[i].w&&h==sheet[i].h)return i; } return -1;}void c_sheet::addItem(int item){ this->item[n]=item; n++;}void c_sheet::print(){ cout<<"Sheet"<<w<<"X"<<h; cout<<"\titems:"; for(int i=0;i<n;i++){ cout<<item[i]+1<<","; } cout<<endl;}bool c_sheet::cut(){ bool rt=true; int w1,h1,w2,h2,s1,s2; if(w==W&&h==H){ nl=n; for(int i=0;i<=n;i++){ x[i]=(int)(i*1.2*W+.5); y[i]=0; w1=w-::w[item[i]]; h1=::h[item[i]]; s1=findSheet(w1,h1); if(s1>=0)sheet[s1].cut(x[i]+::w[item[i]],y[i]); w2=w; h2=h-::h[item[i]]; s2=findSheet(w2,h2); if(s2>=0)sheet[s2].cut(x[i],y[i]+::h[item[i]]); } } return rt;}bool c_sheet::cut(int xx, int yy){ if(nl==n)return false; int w1,h1,w2,h2,s1,s2; x[nl]=xx;y[nl]=yy; w1=w-::w[item[nl]]; h1=::h[item[nl]]; s1=findSheet(w1,h1); if(s1>=0)sheet[s1].cut(xx+::w[item[nl]],yy); w2=w; h2=h-::h[item[nl]]; s2=findSheet(w2,h2); if(s2>=0)sheet[s2].cut(xx,yy+::h[item[nl]]); nl++; return true;}void c_sheet::draw(ofstream &off){ bool hangle=true; for(int i=0;i<nl;i++){ if(w==W&&h==H)off<<"rectangle "<<x[i]<<","<<y[i]<<" @"<<w<<","<<h<<endl; off<<"rectangle "<<x[i]<<","<<y[i]<<" @"<<::w[item[i]]<<","<<::h[item[i]]<<endl; off<<"-hatch p ansi31 0.25 "<<hangle*90<<" "<<endl; hangle=!hangle; off<<"-hatch "<<x[i]+::w[item[i]]/2.0<<","<<y[i]+::h[item[i]]/2.0<<" "<<endl; off<<"line "<<x[i]<<","<<y[i] }}

原始数据文件见附件二

附件二.txt

参考文献：