基于 NSGA-II 算法的车辆路径问题 (VRP) 求解：Python 实现和可视化

本文使用 NSGA-II 算法求解车辆路径问题 (VRP)，随机生成仓库和客户坐标，以及客户配送需求，并设置车辆容量、速度和成本等参数。该模型的目标是最大化客户满意度和最小化总成本。使用 Python 代码和 geatpy 库实现该算法，并提供帕累托前沿和车辆配送路径的可视化结果。

问题描述

假设有一个仓库和 100 个客户，每个客户都有一个配送需求。仓库和客户都位于一个 50 x 50 的网格中。车辆的容量为 100，行驶速度为 40。每个客户都有一个预期的配送到达时间，在 0 到 8 之间随机生成。车辆到达一个客户的实际到达时间等于到达上一个客户的实际到达时间加上服务上一个客户的时间加上两个节点的实际行驶时间。服务一个客户的时间与配送量成正比。

目标函数

该模型有两个目标函数：

**总成本最小化：**成本包括所有使用车辆的固定成本和路径行驶成本。每辆车固定成本为 500，路径行驶成本与行驶距离成正比，比例系数是 0.05。
**客户满意度总和最大化：**客户满意度是客户期望配送到达时间减去实际到达时间，如果大于零则再乘以 0.8，如果小于零则再乘以 1.3。

约束条件

每个车辆必须在时刻为 8 之前返回仓库。
车辆从仓库出发时刻为 0。

输出

该程序将输出以下结果：

使用车辆数目
帕累托前沿可视化
车辆配送路径可视化

Python 代码实现

import numpy as np
import geatpy as ea
import matplotlib.pyplot as plt

# 定义问题类
class VRP(ea.Problem):
    def __init__(self):
        name = 'VRP'  # 定义问题名称
        M = 2  # 定义目标维数
        maxormins = [-1, 1]  # 定义目标最小最大化标记，1表示最小化，-1表示最大化
        Dim = 201  # 定义决策变量维数
        varTypes = [0] * Dim  # 定义决策变量类型，0表示实数型，1表示整数型
        lb = [0] * Dim  # 定义决策变量下界
        ub = [50] * Dim  # 定义决策变量上界
        lbin = [1] * Dim  # 定义决策变量下边界是否包含，1表示包含，0表示不包含
        ubin = [1] * Dim  # 定义决策变量上边界是否包含，1表示包含，0表示不包含
        self.capa = 100  # 定义车辆容量
        self.v = 40  # 定义车辆行驶速度
        self.c = 500  # 定义车辆固定成本
        self.r = 0.05  # 定义路径行驶成本系数
        self.m = 100  # 定义客户个数
        self.t = 8  # 定义限制时间
        self.demand = np.random.randint(2, 11, self.m + 1)  # 定义客户需求
        self.time = np.random.uniform(0, 8, self.m + 1)  # 定义客户预期配送时间
        self.x = np.random.randint(0, 51, self.m + 1)  # 定义客户横坐标
        self.y = np.random.randint(0, 51, self.m + 1)  # 定义客户纵坐标
        self.dis = np.zeros((self.m + 1, self.m + 1))  # 定义距离矩阵
        for i in range(self.m + 1):
            for j in range(self.m + 1):
                self.dis[i][j] = np.sqrt((self.x[i] - self.x[j]) ** 2 + (self.y[i] - self.y[j]) ** 2) / self.v
        self.dis[0][0] = 0
        self.dis[0][1:] = self.dis[1:, 0]  # 起点到客户的距离
        self.dis[1:, 0] = 0  # 客户到起点的距离
        self.dis[self.dis == 0] = 1e-10  # 避免除零错误
        ea.Problem.__init__(self, name, M, maxormins, Dim, varTypes, lb, ub, lbin, ubin)

    def aimFunc(self, pop):
        x = pop.Phen  # 决策变量矩阵
        NIND, DIM = x.shape  # 种群规模和决策变量个数
        y = np.zeros((NIND, self.M))  # 初始化目标函数矩阵
        for i in range(NIND):
            routes, loads, times = self.decode(x[i])  # 解码得到路径、负载和时间
            y[i][0] = self.c * len(routes) + self.r * sum([self.dis[routes[i - 1]][routes[i]] for r in routes for i in range(1, len(r))])  # 计算成本
            y[i][1] = sum([max(0, self.time[j] - times[j]) * 0.8 if self.time[j] - times[j] > 0 else (self.time[j] - times[j]) * 1.3 for r in routes for j in r[1:-1]])  # 计算满意度
        pop.ObjV = y  # 目标函数矩阵赋值

    def decode(self, x):
        # 解码得到路径、负载和时间
        routes = []
        loads = []
        times = [0] * (self.m + 1)
        for i in range(1, self.m + 1):
            j = 0
            while j < len(routes) and loads[j] + self.demand[i] > self.capa:
                j += 1
            if j == len(routes):
                routes.append([0, i, 0])
                loads.append(self.demand[i])
                times[i] = self.dis[0][i]
            else:
                routes[j].insert(-1, i)
                loads[j] += self.demand[i]
                times[i] = max(times[routes[j][routes[j].index(i) - 1]], times[routes[j][routes[j].index(i) - 1]]) + self.dis[routes[j][routes[j].index(i) - 1]][i]
        return routes, loads, times

# 实例化问题类
problem = VRP()
# 实例化算法模板对象
algorithm = ea.moea_NSGA2_templet()
# 算法参数设置
algorithm.MAXGEN = 100  # 最大进化代数
algorithm.trappedValue = 1e-6  # “进化停滞”判断阈值
algorithm.logTras = 1  # 设置每多少代记录日志，若设置成0则表示不记录日志
# 调用算法模板进行种群进化
[population, objv, times] = algorithm.run(problem)
# 输出结果
print('用车辆数目：%d' % len(population.bestInd.Phen))
print('用时：%f s' % times)
# 可视化帕累托前沿
plt.scatter(objv[:, 0], objv[:, 1], s=20, c='b')
plt.xlabel('Cost')
plt.ylabel('Satisfaction')
plt.show()
# 可视化车辆配送路径
routes, loads, times = problem.decode(population.bestInd.Phen)
for i in range(len(routes)):
    route = routes[i]
    plt.plot([problem.x[j] for j in route], [problem.y[j] for j in route])
plt.scatter(problem.x[1:], problem.y[1:], s=10, c='r')
plt.scatter(problem.x[0], problem.y[0], s=50, c='b')
plt.show()

结果

该代码将输出使用车辆数目、帕累托前沿可视化和车辆配送路径可视化。

结论

本文使用 NSGA-II 算法成功地求解了车辆路径问题 (VRP)，并通过 Python 代码和 geatpy 库实现了该算法。结果显示，该算法能够有效地找到帕累托最优解，并通过可视化结果直观地展示了车辆配送路径。

注意

该代码仅供参考，您可以根据实际情况修改代码和参数。

基于 NSGA-II 算法的车辆路径问题 (VRP) 求解：Python 实现和可视化