C++实现基因遗传驱动的列生成算法

列生成算法可适用于可用列进行线性化表示的非线性问题，其核心思想为使用对偶原理寻找最有可能改善原问题的解的列，最好的情况下，原问题为非线性的但对偶模型为线性模型，因此总能求解出一个好列，例如cutstock问题（卷纸切割问题）。

那么当对偶模型也为非线性问题无法使用线性规划的方法来求解新方案时，启发或元启发算法就可以派上用场了，可以求解出局部最优的进基列。

本文要求解的问题是指派、调度问题的一种，将m个工件分配到n个机床上进行加工，并且应当具有一定的顺序，使得最终的总用时最小的类似这种问题。

class MS {
public:int m, n;vector<vector<double>>p;//p_ij加工时间vector<double>D;//D_ivector<double>w;//w_j工作j的权重
public:void read();vector<oneG> getNewSbyGenetic(IloNumArray PI, IloNumArray RHO,int i);vector<bool> getNewSbyDP(IloNumArray PI, IloNumArray RHO);double getC(vector<bool> S,int i);void cplex();void cplexSolveMILP();
};struct oneG {
public:double fitness;struct vector<bool>detail;
};class Genetic {
public:int m, n;//机器数，作业数。int i;//子问题int number;//种群规模；vector<oneG>P;//现有种群。MS* ms;IloNumArray PI;IloNumArray RHO;
public:Genetic(MS* ms,int m, int n,int i,IloNumArray PI, IloNumArray RHO, int number=100) {this->m = m;this->n = n;this->number = number;this->ms = ms;this->i = i;this->PI = PI;this->RHO=RHO;}vector<oneG> initP();//种群初始化oneG randOne();//随机产生子代void cross();//交叉遗传；
};

cpp文件：

#include"MS.h"bool sortFitness(oneG a, oneG b) { return (a.fitness < b.fitness); }
void MS::read() {//文件输入string tempString;ifstream readFile("data.txt");while (readFile >> tempString) {if (tempString == "机器数目:") {readFile >> m;//m个机器}if (tempString == "工作数目:") {readFile >> n;//n个工作}if (tempString == "加工时间:") {for (int i = 0; i < m; i++) {vector<double>pi;for (int j = 0; j < n; j++) {double pij;readFile >> pij;pi.push_back(pij);}p.push_back(pi);}}if(tempString=="中断时间:"){for (int i = 0; i < m; i++) {double di;readFile >> di;D.push_back(di);}}if (tempString == "工作权重:") {for (int j = 0; j < n; j++) {double wj;readFile >> wj;w.push_back(wj);}}}//std::cout << "数据输入完成" << endl;
}
void MS::cplexSolveMILP() {IloEnv env2;IloModel MILP(env2);IloArray<IloArray<IloNumVarArray>>X(env2);IloNumVarArray C(env2, n);for (int j = 0; j < n; j++) {string name = "C_" + to_string(j + 1);C[j] = IloNumVar(env2, 0, INFINITY,name.c_str());}IloRangeArray con1(env2);//第一组约束。IloArray<IloRangeArray> con2(env2);IloArray<IloRangeArray> con3(env2);IloArray<IloArray<IloRangeArray>> con4(env2);for (int j = 0; j < n; j++) {con1.add(IloRange(env2, 1, 1));}for (int r = 0; r < n; r++) {con2.add(IloRangeArray(env2));for (int i = 0; i < m; i++) {con2[r].add(IloRange(env2, -INFINITY, 1));}}for (int r = 0; r < n-1; r++) {con3.add(IloRangeArray(env2));for (int i = 0; i < m; i++) {con3[r].add(IloRange(env2, -INFINITY, 0));}}for (int i = 0; i < m; i++) {con4.add(IloArray<IloRangeArray>(env2));for (int j = 0; j < n; j++) {con4[i].add(IloRangeArray(env2));for (int r = 0; r < n; r++) {con4[i][j].add(IloRange(env2, -INFINITY, 0));}}}for (int i = 0; i < m; i++) {IloArray<IloNumVarArray>X_i(env2);for (int j = 0; j < n; j++) {IloNumVarArray X_ij(env2);for (int r = 0; r < n; r++) {string name = "X" + to_string(i + 1) + to_string(j + 1) + to_string(r + 1);X_ij.add(IloNumVar(env2, 0, 1, ILOBOOL,name.c_str()));}X_i.add(X_ij);}X.add(X_i);}for (int j = 0; j < n; j++) {IloExpr expr(env2);for (int i = 0; i < m; i++) {for (int r = 0; r < n; r++) {expr += X[i][j][r];}}con1[j].setExpr(expr);}for (int r = 0; r < n; r++) {for (int i = 0; i < m; i++) {IloExpr expr(env2);for (int j = 1; j < n; j++) {expr += X[i][j][r];}con2[r][i].setExpr(expr);}MILP.add(con2[r]);}for (int r = 0; r < n-1; r++) {for (int i = 0; i < m; i++) {IloExpr expr(env2);for (int j = 1; j < n; j++) {expr += X[i][j][r]-X[i][j][r+1];}con3[r][i].setExpr(expr);}MILP.add(con3[r]);}for (int i = 0; i < m; i++) {for (int j = 0; j < n; j++) {for (int r = 0; r < n; r++) {IloExpr expr(env2);expr += p[i][j] * X[i][j][r];for (int l = 0; l < r - 1; l++) {for (int u = 1; u < n; u++) {expr += X[i][u][l] * p[i][u];}}for (int v = 0; v < r; v++) {for (int q = r + 1; q < n; q++) {for (int k = 0; k < n; k++) {expr += D[i] * pow(1 - D[i], v - 1) * p[i][k] * X[i][k][q];}}}expr -= C[j]+M*(1-X[i][j][r]);con4[i][j][r].setExpr(expr);}MILP.add(con4[i][j]);}}IloObjective obj;IloExpr expr2(env2);for (int j = 0; j < n; j++) {expr2 += w[j]*C[j];}obj = IloMinimize(env2, expr2);MILP.add(obj);MILP.add(con1);IloCplex MILPSolver(MILP);MILPSolver.exportModel("MILPmodel.lp");MILPSolver.solve();cout<<"求解结果:"<< MILPSolver.getStatus()<<endl;cout << "最优解:" << MILPSolver.getObjValue()<<endl;IloNumArray val(env2);MILPSolver.getValues(C,val);cout << val<<endl;for (int i = 0; i < m; i++) {cout << "机器" + to_string(i + 1) + "上的工作:";for (int j = 0; j < n; j++) {IloNumArray val(env2);MILPSolver.getValues(X[i][j], val);for (int r = 0; r < n; r++) {if (val[r] > 0.5) {cout << j << "工作在"<<r<<"位置、";}}}cout << endl;}/*for (int i = 0; i < m; i++) {cout << "机器" + to_string(i + 1) + "上的工作:";for (int j = 0; j < n; j++) {IloNumArray val(env2);MILPSolver.getValues(X[i][j], val);cout << val<<endl;}}*/
}
void MS::cplex() {IloEnv env;try {//主问题IloModel msOpt(env);IloObjective rollsUsed = IloAdd(msOpt, IloMinimize(env));IloRangeArray constraints1(env);//第一组约束。IloRangeArray constraints2(env);for (int j = 0; j < n; j++) {constraints1.add(IloRange(env, 1, 1));}for (int i = 0; i < m; i++) {constraints2.add(IloRange(env, -INFINITY, 1));}//env.out() << constraints1.getSize();IloArray<IloNumVarArray> S(env);vector<vector<bool>>initSCs{ {0,0,0,1,0,1,0,1,0,0},{1,1,1,0,0,0,0,0,1,1},{0,0,0,0,1,0,1,0,0,0} };for (int i = 0; i < m; i++) {S.add(IloNumVarArray(env));//初始化列IloNumArray initS(env, n);//vector<bool>initSC=initSCs[i];//for (int ii = 0; ii < n; ii++) {initS[ii] = initSC[ii];}/*精确模型最优解机器1上的工作:3、5、7、机器2上的工作 : 0、1、2、8、9、机器3上的工作 : 4、6、*/double sj=getC(initSC,i);IloExpr expr(env);for (int ii = 0; ii < 1; ii++) {string name = "lambda" + to_string(i+1) + to_string(ii+1);S[i].add(IloNumVar(rollsUsed(sj) + constraints1(initS)));//cout << S[i]<<endl;S[i][ii].setName(name.c_str());expr += S[i][ii];}constraints2[i].setExpr(expr);}msOpt.add(constraints1);msOpt.add(constraints2);IloCplex msSolver(msOpt);//msSolver.exportModel("model.lp");bool Generate = 1;while (Generate) {Generate = 0;msSolver.exportModel("model.lp");msSolver.solve();IloNumArray price1(env, n);IloNumArray price2(env, m);IloNumArray newSolutions(env, n);//输出对偶值。cout << msSolver.getStatus() << endl;cout << msSolver.getObjValue() << endl;;for (int j = 0; j < n; j++) {price1[j] = msSolver.getDual(constraints1[j]);cout << price1[j] <<" ";}cout << endl;for (int i = 0; i < m; i++) {price2[i] = msSolver.getDual(constraints2[i]);cout << price2[i] << " ";}cout << endl;//子问题for (int i = 0; i < m; i++) {vector<oneG> Gs = getNewSbyGenetic(price1, price2, i);for (int l = 0; l < Gs.size(); l++) {oneG G = Gs[l];vector<bool> newSC = G.detail;for (int j = 0; j < newSC.size(); j++) {newSolutions[j] = newSC[j];}cout << newSolutions;double sj = getC(newSC, i);cout << G.fitness << endl;if (G.fitness < -RC_EPS) { Generate = 1; };IloNumArray newCons2(env, m);for (int iii = 0; iii < m; iii++) {newCons2[iii] = 1;}string name = "lambda" + to_string(i + 1) + to_string(S[i].getSize() + 1);S[i].add(IloNumVar(rollsUsed(sj) + constraints1(newSolutions) + constraints2(newCons2)));S[i][S[i].getSize() - 1].setName(name.c_str());}}}for (int i = 0; i < S.getSize(); i++) {msOpt.add(IloConversion(env, S[i], ILOBOOL));}msSolver.solve();msSolver.exportModel("model.lp");cout << "Solution status: " << msSolver.getStatus() << endl;cout << "最优解:" << msSolver.getObjValue() << endl;for (int i = 0; i < S.getSize(); i++) {IloNumArray val(env);msSolver.getValues(S[i], val);cout << val<<endl;}}catch (IloException& ex) {cerr << "Error: " << ex << endl;}catch (...) {cerr << "Error" << endl;}env.end();}
vector<oneG> MS::getNewSbyGenetic(IloNumArray PI, IloNumArray RHO,int i) {vector<oneG> solutions;Genetic genetic(this,m, n,i, PI, RHO);//cout <<"第"<<i<< "问题初始化种群" << endl;genetic.initP();//cout << "开始遗传算法" << endl;for (int k = 0; k < 50; k++) {genetic.jiaocha();}for (int k = 1;k < genetic.P.size()+1; k++) {if (genetic.P[k].fitness != genetic.P[k - 1].fitness) {solutions.push_back(genetic.P[k - 1]);//添加新方案组}if (solutions.size() > 10|| genetic.P[k].fitness>0)break;}//cout << "找到可行新方案" << endl;return solutions;
}
vector<bool> MS::getNewSbyDP(IloNumArray PI, IloNumArray RHO) {vector<bool>solution;return solution;
}
double MS:: getC(vector<bool> S,int i) {//计算方案的代价double csi = 0;double t=0;vector<int>work;for (int j = 0; j < n; j++) {if (S[j] != 0) {work.push_back(j);}}for (int j = 0; j < work.size(); j++) {t += p[i][work[j]] * pow(1 - D[i], i);for (int jj = j+1; jj < work.size(); jj++) {t += p[i][work[jj]] * D[i]*pow(1-D[i],i);}csi += t * w[work[j]];}return csi;
}
vector<oneG> Genetic::initP(){//初始化种群vector<oneG>Q;for (int k = 0; k < number; k++) {oneG np = randOne();Q.push_back(np);}P = Q;return P;
}
oneG Genetic::randOne() {oneG randO;for (int k = 0; k < n; k++) {randO.detail.push_back(0);}for (int k = 0; k < n; k++) {int temp= (rand() % (n + 1));double prob = double(rand() % 10)/10;if (prob < 0.3) {randO.detail[temp]=1;}}//计算适应度randO.fitness= ms->getC(randO.detail, i);//根据对偶值计算for (int j = 0; j < n; j++) {randO.fitness -=randO.detail[j]*PI[j];}for (int j = 0; j < m; j++) {randO.fitness -= RHO[j];}return randO;
}
void Genetic::cross() {double sum = 0;vector<oneG>Q;//新种群//cout << "交叉遗传中" << endl;for (int k = 0; k < P.size(); k++) {for (int ii = k+1; ii < P.size(); ii++) {oneG ne;int point = (rand() % (n + 1));for (int j = 0; j < n; j++) {if (j < point) {ne.detail.push_back(P[k].detail[j]);}else {ne.detail.push_back(P[ii].detail[j]);}}ne.fitness = ms->getC(ne.detail, i);for (int j = 0; j < n; j++) {ne.fitness -= ne.detail[j] * PI[j];}for (int j = 0; j < m; j++) {ne.fitness -= RHO[j];}Q.push_back(ne);}}Q.insert(Q.end(), P.begin(), P.end());//加父本sort(Q.begin(), Q.end(), sortFitness);//排序P.resize(0);for (int k = 1; k < Q.size()+1; k++) {//去重复种群if (Q[k].fitness != Q[k - 1].fitness) {P.push_back(Q[k - 1]);}if (P.size() > number)break;}//cout << "交叉完成" << endl;
}

主程序：

#include"MS.h"ILOSTLBEGINusing namespace std;void main() {MS ms;ms.read();int a = 0;cout << "请选择求解方式：1.MILP精确求解；2.基因遗传列生成求解。" << endl;cin >> a;if (a == 2) ms.cplex();else ms.cplexSolveMILP();
}

C++实现基因遗传驱动的列生成算法相关推荐

从一维cutting问题看列生成算法
从一维cutting问题看列生成算法列生成算法一维cutting问题列生成原理集合覆盖模型列生成步骤说明 python 代码实现列生成算法列生成在求解大型线性规划时往往都表现良好,在很多 ...
cutting stock 问题的列生成算法
本文以 cutting stock 问题为例,介绍列生成算法的数学模型. 1. 问题描述长板长度:10m 子板需求:2m: 100; 3m: 50; 5m: 20 - 优化目标:最小化,长板切割总数 ...
列生成算法column generation求解思路总结篇2——CW分解
目录 1.CW的思路 2.CW重写模型 3.DW分解算法步骤列生成是通过重写模型为P2,然后松弛主问题(只考虑m个变量),再通过对偶理论的思路构造子问题,即寻找最大改进主问题的非基变量进入松弛主问题 ...
列生成算法求解矩形下料问题（Matlab代码）
目录深切哀悼抗击新冠肺炎疫情斗争牺牲烈士和逝世同胞 Matlab代码算例文件运行结果最优解方案展示深切哀悼抗击新冠肺炎疫情斗争牺牲烈士和逝世同胞这学期接了<运筹学>课程,在国内 ...
下料问题的列生成算法，matlab代码，简单易懂
下料问题的列生成算法,matlab代码,简单易懂 ID:18300606016240721绿兔叽
运筹优化学习21：Java调用Cplex实现求解Cuting Stock Porblem的列生成算法详解
目录 1 CSP问题与模型 1.1 问题描述 1.2 模型构建 2 列生成方法理论 2.1 引子 2.2 单纯形法到列生成 2.3 subproblem 2.3.1 对偶理论 2.3.2 影子价格 2 ...
OR Paper Weekly(一) | 用机器学习生成列生成的列，元启发式算法=动物世界？看OR68年发文数据，哪国位居榜首？
作者:王源,徐思坤,陈贤邦 OR Paper Weekly 栏目将会从运筹学顶级期刊上选择一部分有趣的文章,对这些文章的主要研究内容进行一个概述/点评.OR Paper Weekly 的特点是不做大 ...
Talk | 东京大学博士生刘海洋：多模态驱动谈话动作生成-质量与多样性
本期为TechBeat人工智能社区第471期线上Talk! 北京时间2月1日(周三)20:00,东京大学情报理工系博士生--刘海洋的Talk将准时在TechBeat人工智能社区开播! 他与大家分享的主 ...
4列的计算机代码,干货 | 10分钟带你彻底了解column generation（列生成）算法的原理附java代码...
OUTLINE 前言预备知识预警什么是column generation 相关概念科普 Cutting Stock Problem CG求解Cutting Stock Problem 列生成代码 ...

C++实现基因遗传驱动的列生成算法

C++实现基因遗传驱动的列生成算法相关推荐

最新文章

热门文章