干货|十分钟快速掌握CPLEX求解VRPTW数学模型(附JAVA代码及CPLEX安装流程)
内容提要:
*什么是VRPTW
*CPLEX求解VRPTW实例
*CPLEX操作补充说明
1.什么是VRPTW
提到 带时间窗车辆路径问题 (vehicle routing problems with time windows,VRPTW),就不得不先说说 车辆路径问题 (VRP)。
- 什么是VRP?
车辆路径问题 (VRP)最早是由 Dantzig 和 Ramser 于1959年首次提出,它是指一定数量的客户,各自有不同数量的货物需求,配送中心向客户提供货物,由一个车队负责分送货物,组织适当的行车路线, 目标 是使得客户的需求得到满足,并能在一定的约束下,达到诸如 路程最短、成本最小、耗费时间最少 等目的。
VRP图示
- 什么是VRPTW?
由于VRP问题的持续发展,考虑需求点对于 车辆到达的时间 有所要求之下,在车辆途程问题之中加入时窗的限制,便成为 带时间窗车辆路径问题 (VRP with Time Windows, VRPTW)。带时间窗车辆路径问题(VRPTW)是在VRP上加上了客户的被访问的时间窗约束。在VRPTW问题中,除了行驶成本之外, 成本函数还要包括由于早到某个客户而引起的 等待时间 和 客户需要的服务时间 。
在VRPTW中,车辆除了要满足VRP问题的限制之外,还必须要满足需求点的时窗限制,而需求点的时窗限制可以分为两种,一种是 硬时窗 (Hard Time Window),硬时窗要求车辆必须要在时窗内到达,早到必须等待,而迟到则拒收;另一种是 软时窗 (Soft Time Window),不一定要在时窗内到达,但是在时窗之外到达必须要处罚,以处罚替代等待与拒收是软时窗与硬时窗最大的不同。
2.CPLEX求解VRPTW实例
解决 带时间窗车辆路径问题 (vehicle routing problems with time windows,VRPTW)的常用 求解方法 :
- 1.精确解算法(Exact methods)
精确解算法解VRPTW问题主要有三个策略, 拉格朗日松弛、列生成 和 动态规划 ,但是可以求解的算例规模非常小。
- 2.途程构建启发式算法(Route-building heuristics)
在问题中以某节点选择原则或是路线安排原则,将需求点一一纳入途程路线的解法。
- 3.途程改善启发式算法(Route-improving heuristics)
先决定一个可行途程,也就是一个起始解,之后对这个起始解一直做改善,直到不能改善为止。
- 4.通用启发式算法(Metaheuristics)
传统区域搜寻方法的最佳解常因起始解的特性或搜寻方法的限制,而只能获得局部最佳解,为了改善此一缺点,近年来在此领域有重大发展,是新一代的 启发式解法 ,包含禁忌搜索算法(Tabu Search)、模拟退火法(Simulated Annealing)、遗传算法(Genetic Algorithm)和门坎接受法(Threshold Accepting)等,可以有效解决陷入局部最优的困扰。
VRPTW问题建模实例
接下来分享一波代码和算例 ↓ ↓ ↓
package vrptw;
import java.io.BufferedReader;
import java.io.FileReader;
import java.util.ArrayList;
import java.util.Scanner;
import ilog.concert.IloException;
import ilog.concert.IloNumExpr;
import ilog.concert.IloNumVar;
import ilog.concert.IloNumVarType;
import ilog.cplex.IloCplex;
* @author:huangnan
* @School: HuaZhong University of science and technology
* @操作说明:读入不同的文件前要手动修改vetexnum参数,参数值为所有点个数,包括配送中心0和n+1,
* 代码算例截取于Solomon测试算例
//定义参数
class Data{
int vetexnum; //所有点集合n(包括配送中心和客户点,首尾(0和n)为配送中心)
double E; //配送中心时间窗开始时间
double L; //配送中心时间窗结束时间
int vecnum; //车辆数
double cap; //车辆载荷
int[][] vertexs; //所有点的坐标x,y
int[] demands; //需求量
int[] vehicles; //车辆编号
double[] a; //时间窗开始时间【a[i],b[i]】
double[] b; //时间窗结束时间【a[i],b[i]】
double[] s; //客户点的服务时间
int[][] arcs; //arcs[i][j]表示i到j点的弧
double[][] dist; //距离矩阵,满足三角关系,暂用距离表示花费 C[i][j]=dist[i][j]
//截断小数3.26434-->3.2
public double double_truncate(double v){
int iv = (int) v;
if(iv+1 - v <= 0.000000000001)
return iv+1;
double dv = (v - iv) * 10;
int idv = (int) dv;
double rv = iv + idv / 10.0;
return rv;
//类功能:解的可行性判断
class Solution{
double epsilon = 0.0001;
Data data = new Data();
ArrayList<ArrayList<Integer>> routes = new ArrayList<>();
ArrayList<ArrayList<Double>> servetimes = new ArrayList<>();
public Solution(Data data, ArrayList<ArrayList<Integer>> routes, ArrayList<ArrayList<Double>> servetimes) {
super();
this.data = data;
this.routes = routes;
this.servetimes = servetimes;
//函数功能:比较两个数的大小
public int double_compare(double v1,double v2) {
if (v1 < v2 - epsilon) {
return -1;
if (v1 > v2 + epsilon) {
return 1;
return 0;
//函数功能:解的可行性判断
public void fesible() throws IloException {
//车辆数量可行性判断
if (routes.size() > data.vecnum) {
System.out.println("error: vecnum!!!");
System.exit(0);
//车辆载荷可行性判断
for (int k = 0; k < routes.size(); k++) {
ArrayList<Integer> route = routes.get(k);
double capasity = 0;
for (int i = 0; i < route.size(); i++) {
capasity += data.demands[route.get(i)];
if (capasity > data.cap) {
System.out.println("error: cap!!!");
System.exit(0);
//时间窗、车容量可行性判断
for (int k = 0; k < routes.size(); k++) {
ArrayList<Integer> route = routes.get(k);
ArrayList<Double> servertime = servetimes.get(k);
double capasity = 0;
for (int i = 0; i < route.size()-1; i++) {
int origin = route.get(i);
int destination = route.get(i+1);
double si = servertime.get(i);
double sj = servertime.get(i+1);
if (si < data.a[origin] && si > data.b[origin]) {
System.out.println("error: servertime!");
System.exit(0);
if (double_compare(si + data.dist[origin][destination],data.b[destination]) > 0) {
System.out.println(origin + ": [" + data.a[origin] + ","+data.b[origin]+"]"+ " "+ si);
System.out.println(destination + ": [" + data.a[destination] + ","+data.b[destination]+"]"+ " "+ sj);
System.out.println(data.dist[origin][destination]);
System.out.println(destination + ":" );
System.out.println("error: forward servertime!");
System.exit(0);
if (double_compare(sj - data.dist[origin][destination],data.a[origin]) < 0) {
System.out.println(origin + ": [" + data.a[origin] + ","+data.b[origin]+"]"+ " "+ si);
System.out.println(destination + ": [" + data.a[destination] + ","+data.b[destination]+"]"+ " "+ sj);
System.out.println(data.dist[origin][destination]);
System.out.println(destination + ":" );
System.out.println("error: backward servertime!");
System.exit(0);
if (capasity > data.cap) {
System.out.println("error: cap!!!");
System.exit(0);
//类功能:建立模型并求解
public class Vrptw {
Data data; //定义类Data的对象
IloCplex model; //定义cplex内部类的对象
public IloNumVar[][][] x; //x[i][j][k]表示弧arcs[i][j]被车辆k访问
public IloNumVar[][] w; //车辆访问所有点的时间矩阵
double cost; //目标值object
Solution solution;
public Vrptw(Data data) {
this.data = data;
//函数功能:解模型,并生成车辆路径和得到目标值
public void solve() throws IloException {
ArrayList<ArrayList<Integer>> routes = new ArrayList<>(); //定义车辆路径链表
ArrayList<ArrayList<Double>> servetimes = new ArrayList<>(); //定义花费时间链表
//初始化车辆路径和花费时间链表,链表长度为车辆数k
for (int k = 0; k < data.vecnum; k++) {
ArrayList<Integer> r = new ArrayList<>(); //定义一个对象为int型的链表
ArrayList<Double> t = new ArrayList<>(); //定义一个对象为double型的链表
routes.add(r); //将上述定义的链表加入到链表routes中
servetimes.add(t); //同上
//判断建立的模型是否可解
if(model.solve() == false){
//模型不可解
System.out.println("problem should not solve false!!!");
return; //若不可解,则直接跳出solve函数
else{
//模型可解,生成车辆路径
for(int k = 0; k < data.vecnum; k++){
boolean terminate = true;
int i = 0;
routes.get(k).add(0);
servetimes.get(k).add(0.0);
while(terminate){
for (int j = 0; j < data.vetexnum; j++) {
if (data.arcs[i][j]>=0.5 && model.getValue(x[i][j][k])>=0.5) {
routes.get(k).add(j);
servetimes.get(k).add(model.getValue(w[j][k]));
i = j;
break;
if (i == data.vetexnum-1) {
terminate = false;
solution = new Solution(data,routes,servetimes);
cost = model.getObjValue();
System.out.println("routes="+solution.routes);
//函数功能:根据VRPTW数学模型建立VRPTW的cplex模型
private void build_model() throws IloException {
//model
model = new IloCplex();
model.setOut(null);
//variables
x = new IloNumVar[data.vetexnum][data.vetexnum][data.vecnum];
w = new IloNumVar[data.vetexnum][data.vecnum]; //车辆访问点的时间
//定义cplex变量x和w的数据类型及取值范围
for (int i = 0; i < data.vetexnum; i++) {
for (int k = 0; k < data.vecnum; k++) {
w[i][k] = model.numVar(0, 1e15, IloNumVarType.Float, "w" + i + "," + k);
for (int j = 0; j < data.vetexnum; j++) {
if (data.arcs[i][j]==0) {
x[i][j] = null;
else{
//Xijk,公式(10)-(11)
for (int k = 0; k < data.vecnum; k++) {
x[i][j][k] = model.numVar(0, 1, IloNumVarType.Int, "x" + i + "," + j + "," + k);
//加入目标函数
//公式(1)
IloNumExpr obj = model.numExpr();
for(int i = 0; i < data.vetexnum; i++){
for(int j = 0; j < data.vetexnum; j++){
if (data.arcs[i][j]==0) {
continue;
for(int k = 0; k < data.vecnum; k++){
obj = model.sum(obj, model.prod(data.dist[i][j], x[i][j][k]));
model.addMinimize(obj);
//加入约束
//公式(2)
for(int i= 1; i < data.vetexnum-1;i++){
IloNumExpr expr1 = model.numExpr();
for (int k = 0; k < data.vecnum; k++) {
for (int j = 1; j < data.vetexnum; j++) {
if (data.arcs[i][j]==1) {
expr1 = model.sum(expr1, x[i][j][k]);
model.addEq(expr1, 1);
//公式(3)
for (int k = 0; k < data.vecnum; k++) {
IloNumExpr expr2 = model.numExpr();
for (int j = 1; j < data.vetexnum; j++) {
if (data.arcs[0][j]==1) {
expr2 = model.sum(expr2, x[0][j][k]);
model.addEq(expr2, 1);
//公式(4)
for (int k = 0; k < data.vecnum; k++) {
for (int j = 1; j < data.vetexnum-1; j++) {
IloNumExpr expr3 = model.numExpr();
IloNumExpr subExpr1 = model.numExpr();
IloNumExpr subExpr2 = model.numExpr();
for (int i = 0; i < data.vetexnum; i++) {
if (data.arcs[i][j]==1) {
subExpr1 = model.sum(subExpr1,x[i][j][k]);
if (data.arcs[j][i]==1) {
subExpr2 = model.sum(subExpr2,x[j][i][k]);
expr3 = model.sum(subExpr1,model.prod(-1, subExpr2));
model.addEq(expr3, 0);
//公式(5)
for (int k = 0; k < data.vecnum; k++) {
IloNumExpr expr4 = model.numExpr();
for (int i = 0; i < data.vetexnum-1; i++) {
if (data.arcs[i][data.vetexnum-1]==1) {
expr4 = model.sum(expr4,x[i][data.vetexnum-1][k]);
model.addEq(expr4, 1);
//公式(6)
double M = 1e5;
for (int k = 0; k < data.vecnum; k++) {
for (int i = 0; i < data.vetexnum; i++) {
for (int j = 0; j < data.vetexnum; j++) {
if (data.arcs[i][j] == 1) {
IloNumExpr expr5 = model.numExpr();
IloNumExpr expr6 = model.numExpr();
expr5 = model.sum(w[i][k], data.s[i]+data.dist[i][j]);
expr5 = model.sum(expr5,model.prod(-1, w[j][k]));
expr6 = model.prod(M,model.sum(1,model.prod(-1, x[i][j][k])));
model.addLe(expr5, expr6);
//公式(7)
for (int k = 0; k < data.vecnum; k++) {
for (int i = 1; i < data.vetexnum-1; i++) {
IloNumExpr expr7 = model.numExpr();
for (int j = 0; j < data.vetexnum; j++) {
if (data.arcs[i][j] == 1) {
expr7 = model.sum(expr7,x[i][j][k]);
model.addLe(model.prod(data.a[i], expr7), w[i][k]);
model.addLe(w[i][k], model.prod(data.b[i], expr7));
//公式(8)
for (int k = 0; k < data.vecnum; k++) {
model.addLe(data.E, w[0][k]);
model.addLe(data.E, w[data.vetexnum-1][k]);
model.addLe(w[0][k], data.L);
model.addLe(w[data.vetexnum-1][k], data.L);
//公式(9)
for (int k = 0; k < data.vecnum; k++) {
IloNumExpr expr8 = model.numExpr();
for (int i = 1; i < data.vetexnum-1; i++) {
IloNumExpr expr9 = model.numExpr();
for (int j = 0; j < data.vetexnum; j++) {
if (data.arcs[i][j] == 1) {
expr9=model.sum(expr9,x[i][j][k]);
expr8 = model.sum(expr8,model.prod(data.demands[i],expr9));
model.addLe(expr8, data.cap);
//函数功能:从txt文件中读取数据并初始化参数
public static void process_solomon(String path,Data data,int vetexnum) throws Exception{
String line = null;
String[] substr = null;
Scanner cin = new Scanner(new BufferedReader(new FileReader(path))); //读取文件
for(int i =0; i < 4;i++){
line = cin.nextLine(); //读取一行
line = cin.nextLine();
line.trim(); //返回调用字符串对象的一个副本,删除起始和结尾的空格
substr = line.split(("\\s+")); //以空格为标志将字符串拆分
//初始化参数
data.vetexnum = vetexnum;
data.vecnum = Integer.parseInt(substr[1]);
data.cap = Integer.parseInt(substr[2]);
data.vertexs =new int[data.vetexnum][2]; //所有点的坐标x,y
data.demands = new int[data.vetexnum]; //需求量
data.vehicles = new int[data.vecnum]; //车辆编号
data.a = new double[data.vetexnum]; //时间窗开始时间
data.b = new double[data.vetexnum]; //时间窗结束时间
data.s = new double[data.vetexnum]; //服务时间
data.arcs = new int[data.vetexnum][data.vetexnum];
data.dist = new double[data.vetexnum][data.vetexnum]; //距离矩阵,满足三角关系,用距离表示cost
for(int i =0; i < 4;i++){
line = cin.nextLine();
//读取vetexnum-1行数据
for (int i = 0; i < data.vetexnum - 1; i++) {
line = cin.nextLine();
line.trim();
substr = line.split("\\s+");
data.vertexs[i][0] = Integer.parseInt(substr[2]);
data.vertexs[i][1] = Integer.parseInt(substr[3]);
data.demands[i] = Integer.parseInt(substr[4]);
data.a[i] = Integer.parseInt(substr[5]);
data.b[i] = Integer.parseInt(substr[6]);
data.s[i] = Integer.parseInt(substr[7]);
cin.close();//关闭流
//初始化配送中心参数
data.vertexs[data.vetexnum-1] = data.vertexs[0];
data.demands[data.vetexnum-1] = 0;
data.a[data.vetexnum-1] = data.a[0];
data.b[data.vetexnum-1] = data.b[0];
data.E = data.a[0];
data.L = data.b[0];
data.s[data.vetexnum-1] = 0;
double min1 = 1e15;
double min2 = 1e15;
//距离矩阵初始化
for (int i = 0; i < data.vetexnum; i++) {
for (int j = 0; j < data.vetexnum; j++) {
if (i == j) {
data.dist[i][j] = 0;
continue;
data.dist[i][j] = Math.sqrt((data.vertexs[i][0]-data.vertexs[j][0])*(data.vertexs[i][0]-data.vertexs[j][0])+
(data.vertexs[i][1]-data.vertexs[j][1])*(data.vertexs[i][1]-data.vertexs[j][1]));
data.dist[i][j]=data.double_truncate(data.dist[i][j]);
data.dist[0][data.vetexnum-1] = 0;
data.dist[data.vetexnum-1][0] = 0;
//距离矩阵满足三角关系
for (int k = 0; k < data.vetexnum; k++) {
for (int i = 0; i < data.vetexnum; i++) {
for (int j = 0; j < data.vetexnum; j++) {
if (data.dist[i][j] > data.dist[i][k] + data.dist[k][j]) {
data.dist[i][j] = data.dist[i][k] + data.dist[k][j];
//初始化为完全图
for (int i = 0; i < data.vetexnum; i++) {
for (int j = 0; j < data.vetexnum; j++) {
if (i != j) {
data.arcs[i][j] = 1;
else {
data.arcs[i][j] = 0;
//除去不符合时间窗和容量约束的边
for (int i = 0; i < data.vetexnum; i++) {
for (int j = 0; j < data.vetexnum; j++) {
if (i == j) {
continue;
if (data.a[i]+data.s[i]+data.dist[i][j]>data.b[j] || data.demands[i]+data.demands[j]>data.cap) {
data.arcs[i][j] = 0;
if (data.a[0]+data.s[i]+data.dist[0][i]+data.dist[i][data.vetexnum-1]>data.b[data.vetexnum-1]) {
System.out.println("the calculating example is false");
for (int i = 1; i < data.vetexnum-1; i++) {
if (data.b[i] - data.dist[0][i] < min1) {
min1 = data.b[i] - data.dist[0][i];
if (data.a[i] + data.s[i] + data.dist[i][data.vetexnum-1] < min2) {
min2 = data.a[i] + data.s[i] + data.dist[i][data.vetexnum-1];
if (data.E > min1 || data.L < min2) {
System.out.println("Duration false!");
System.exit(0);//终止程序
//初始化配送中心0,n+1两点的参数
data.arcs[data.vetexnum-1][0] = 0;
data.arcs[0][data.vetexnum-1] = 1;
for (int i = 1; i < data.vetexnum-1; i++) {
data.arcs[data.vetexnum-1][i] = 0;
for (int i = 1; i < data.vetexnum-1; i++) {
data.arcs[i][0] = 0;
public static void main(String[] args) throws Exception {
Data data = new Data();
int vetexnum = 102;//所有点个数,包括0,n+1两个配送中心点
//读入不同的文件前要手动修改vetexnum参数,参数值等于所有点个数,包括配送中心
String path = "data/c101.txt";//算例地址
process_solomon(path,data,vetexnum);
System.out.println("input succesfully");
System.out.println("cplex procedure###########################");
Vrptw cplex = new Vrptw(data);
cplex.build_model();
double cplex_time1 = System.nanoTime();
cplex.solve();