[Jsprit]Jsprit学习笔记-vrp问题的求解

一、整体的求解逻辑

下面是Jsprit实现的代码部分

public Collection<VehicleRoutingProblemSolution> searchSolutions() {
        logger.info("algorithm starts: [maxIterations={}]", maxIterations);
        double now = System.currentTimeMillis();
        int noIterationsThisAlgoIsRunning = maxIterations;
        counter.reset();
        //初始化解决方案集合，使用初始解决方案
        Collection<VehicleRoutingProblemSolution> solutions = new ArrayList<>(initialSolutions);
        //调用算法开始的钩子方法
        algorithmStarts(problem, solutions);
        bestEver = Solutions.bestOf(solutions);
        if (logger.isTraceEnabled()) {
            log(solutions);
        }
        logger.info("iterations start");
        for (int i = 0; i < maxIterations; i++) {
            iterationStarts(i + 1, problem, solutions);
            logger.debug("start iteration: {}", i);
            counter.incCounter();
            //从搜索策略管理器中随机选择一个搜索策略
            SearchStrategy strategy = searchStrategyManager.getRandomStrategy();
            // 运行搜索策略并获取发现的解决方案
            DiscoveredSolution discoveredSolution = strategy.run(problem, solutions);
            if (logger.isTraceEnabled()) {
                log(discoveredSolution);
            }
            // 如果发现的解决方案是迄今为止最好的，将其记住
            memorizeIfBestEver(discoveredSolution);
            // 根据发现的解决方案选择策略
            selectedStrategy(discoveredSolution, problem, solutions);
            // 如果提前终止管理器认为应该提前终止算法，则记录日志并跳出循环
            if (terminationManager.isPrematureBreak(discoveredSolution)) {
                logger.info("premature algorithm termination at iteration {}", (i + 1));
                noIterationsThisAlgoIsRunning = (i + 1);
                break;
            }
            // 调用每次迭代结束的钩子方法
            iterationEnds(i + 1, problem, solutions);
        }
        logger.info("iterations end at {} iterations", noIterationsThisAlgoIsRunning);
        // 将迄今为止最好的解决方案添加到解决方案集合中
        addBestEver(solutions);
        // 调用算法结束的钩子方法
        algorithmEnds(problem, solutions);
        // 记录算法运行所花费的时间
        logger.info("took {} seconds", ((System.currentTimeMillis() - now) / 1000.0));
        return solutions;
    }

主要步骤

日志记录和初始化：
- 记录算法开始的日志，包括最大迭代次数。
- 获取当前时间戳，用于之后计算算法运行时间。
初始化变量：
- noIterationsThisAlgoIsRunning 用于记录算法实际运行的迭代次数。
- counter 重置，可能用于跟踪某些事件的次数。
初始化解决方案集合：
- 创建一个解决方案集合 solutions，如果存在初始解决方案，则使用它们初始化。
算法开始钩子：
- 调用 algorithmStarts 方法，通知监听器算法开始。
记录最佳解决方案：
- 使用 Solutions.bestOf 方法找到当前最佳解决方案，并记录。
迭代开始：
- 进入循环，开始迭代过程。
迭代过程：
- 对于每个迭代，执行以下步骤：
  - 记录迭代开始的日志。
  - 递增计数器。
  - 从搜索策略管理器中随机选择一个搜索策略。
  - 运行搜索策略，获取可能的解决方案。
  - 如果启用了日志跟踪，记录发现的解决方案。
  - 如果发现的解决方案是迄今为止最好的，更新最佳解决方案记录。
  - 根据发现的解决方案选择策略。
  - 检查是否满足提前终止条件，如果是，则记录日志并退出循环。
迭代结束：
- 调用 iterationEnds 方法，通知监听器每次迭代结束。
记录实际迭代次数：
- 记录算法实际运行的迭代次数。
添加最佳解决方案：
- 将最佳解决方案添加到解决方案集合中。
算法结束钩子：
- 调用 algorithmEnds 方法，通知监听器算法结束。
记录运行时间：
- 计算并记录算法运行时间。
返回解决方案集合：
- 返回包含所有解决方案的集合。

这个方法体现了一个典型的优化算法结构，包括初始化、迭代、记录日志、选择策略、检查终止条件、记录最佳解决方案和最终返回结果。

二、搜索策略的选择

搜索策略的选择类似轮盘赌的思想，预先对每个搜索策略设置权重，通过遍历所有的搜索策略，计算累计的概率值与随机数进行比较，当累计概率值大于随机数时，返回此时的搜索策略，下面是实现的代码

public SearchStrategy getRandomStrategy() {
        if (random == null)
            throw new IllegalStateException("randomizer is null. make sure you set random object correctly");
        double randomFig = random.nextDouble();
        //用于累计每个策略的概率
        double sumProbabilities = 0.0;
        //如果随机生成的数 randomFig 小于 sumProbabilities，则选择当前遍历到的策略并返回
        for (int i = 0; i < weights.size(); i++) {
            sumProbabilities += weights.get(i) / sumWeights;
            if (randomFig < sumProbabilities) {
                return strategies.get(i);
            }
        }
        throw new IllegalStateException("no search-strategy found");
    }

三、搜索策略执行解

核心代码

// 运行搜索策略并获取发现的解决方案
 DiscoveredSolution discoveredSolution = strategy.run(problem, solutions);

1、解的选择

选择一个当前成本最低的解

//选择成本最低的解决方案
 VehicleRoutingProblemSolution solution = solutionSelector.selectSolution(solutions);
 public VehicleRoutingProblemSolution selectSolution(Collection<VehicleRoutingProblemSolution> solutions) {
        double minCost = Double.MAX_VALUE;
        VehicleRoutingProblemSolution bestSolution = null;
        for (VehicleRoutingProblemSolution sol : solutions) {
            if (bestSolution == null) {
                bestSolution = sol;
                minCost = sol.getCost();
            } else if (sol.getCost() < minCost) {
                bestSolution = sol;
                minCost = sol.getCost();
            }
        }
        return bestSolution;
    }

2、解的破坏

首先对选出来的解进行深度copy

//对选出来的方案进行深copy
VehicleRoutingProblemSolution lastSolution = VehicleRoutingProblemSolution.copyOf(solution);
//对copy后的解决运行搜索策略选择模块
 for (SearchStrategyModule module : searchStrategyModules) {
     lastSolution = module.runAndGetSolution(lastSolution);
 }

然后根据搜索策略模块对解进行破坏和接受

Collection<Job> ruinedJobs = ruin.ruin(previousVrpSolution.getRoutes());
public Collection<Job> ruin(Collection<VehicleRoute> vehicleRoutes) {
        ruinListeners.ruinStarts(vehicleRoutes);
        Collection<Job> unassigned = ruinRoutes(vehicleRoutes);
        logger.trace("ruin: [ruined={}]", unassigned.size());
        ruinListeners.ruinEnds(vehicleRoutes, unassigned);
        return unassigned;
 }

核心代码：

Collection<Job> unassigned = ruinRoutes(vehicleRoutes);
public Collection<Job> ruinRoutes(Collection<VehicleRoute> vehicleRoutes) {
        if (vehicleRoutes.isEmpty()) {
            return Collections.emptyList();
        }
        int nOfJobs2BeRemoved = Math.min(ruinShareFactory.createNumberToBeRemoved(), noJobsToMemorize);
        Collection<Job> jobs = vrp.getJobsWithLocation();
        if (nOfJobs2BeRemoved == 0 || jobs.isEmpty()) {
            return Collections.emptyList();
        }
        Job randomJob = RandomUtils.nextJob(jobs, random);
        return ruinRoutes(vehicleRoutes, randomJob, nOfJobs2BeRemoved);
    }

随机选择一个Job，进行ruin

ruinRoutes(vehicleRoutes, randomJob, nOfJobs2BeRemoved);
private Collection<Job> ruinRoutes(Collection<VehicleRoute> vehicleRoutes, Job targetJob, int nOfJobs2BeRemoved) {
        List<Job> unassignedJobs = new ArrayList<>();
        int nNeighbors = nOfJobs2BeRemoved - 1;
        removeJob(targetJob, vehicleRoutes);
        unassignedJobs.add(targetJob);
        Iterator<Job> neighborhoodIterator = jobNeighborhoods.getNearestNeighborsIterator(nNeighbors, targetJob);
        while (neighborhoodIterator.hasNext()) {
            Job job = neighborhoodIterator.next();
            if (removeJob(job, vehicleRoutes)) {
                unassignedJobs.add(job);
            }
        }
        return unassignedJobs;
    }

输入参数：当前的车辆路线、随机选择的Job、要移除的Job数量
遍历当前的车辆路线，从当前路线的getTourActivities中删除Job，下面是删除的代码

public boolean removeJob(Job job) {
		//先从Jobs列表中删除
        boolean jobRemoved;
        if (!jobs.contains(job)) {
            return false;
        } else {
            jobRemoved = jobs.remove(job);
        }
        //从tourActivities中删除
        boolean activityRemoved = false;
        Iterator<TourActivity> iterator = tourActivities.iterator();
        while (iterator.hasNext()) {
            TourActivity c = iterator.next();
            if (c instanceof JobActivity) {
                Job underlyingJob = ((JobActivity) c).getJob();
                if (job.equals(underlyingJob)) {
                    iterator.remove();
                    activityRemoved = true;
                }
            }
        }
        assert jobRemoved == activityRemoved : "job removed, but belonging activity not.";
        return activityRemoved;
    }

删除之后，把删除的Job添加到未分配的Job列表中，接下来删除Job的邻居
Job的邻居是存储在数组neighbors中，分析neighbors的计算逻辑
初始化

public JobNeighborhoodsOptimized(VehicleRoutingProblem vrp, JobDistance jobDistance, int capacity) {
        super();
        this.vrp = vrp;
        this.jobDistance = jobDistance;
        this.capacity = capacity;
        neighbors = new int[vrp.getJobsInclusiveInitialJobsInRoutes().size()+1][capacity];
        jobs = new Job[vrp.getJobsInclusiveInitialJobsInRoutes().size()+1];
        logger.debug("initialize {}", this);
    }

邻居的计算

private void calculateDistancesFromJob2Job() {
        logger.info("pre-process distances between locations ...");
        StopWatch stopWatch = new StopWatch();
        stopWatch.start();
        for (Job job_i : vrp.getJobsInclusiveInitialJobsInRoutes().values()) {
            if (job_i.getActivities().get(0).getLocation() == null) continue;
            jobs[job_i.getIndex()] = job_i;
            List<ReferencedJob> jobList = new ArrayList<ReferencedJob>(vrp.getJobsInclusiveInitialJobsInRoutes().values().size());
            for (Job job_j : vrp.getJobsInclusiveInitialJobsInRoutes().values()) {
                if (job_j.getActivities().get(0).getLocation() == null) continue;
                if (job_i == job_j) continue;
                double distance = jobDistance.getDistance(job_i, job_j);
                if (distance > maxDistance) maxDistance = distance;
                ReferencedJob referencedJob = new ReferencedJob(job_j, distance);
                jobList.add(referencedJob);
            }
            Collections.sort(jobList,getComparator());
            int neiborhoodSize = Math.min(capacity, jobList.size());
            int[] jobIndices = new int[neiborhoodSize];
            for (int index = 0; index < neiborhoodSize; index++) {
                jobIndices[index] = jobList.get(index).getJob().getIndex();
            }
            neighbors[job_i.getIndex()-1] = jobIndices;
        }
        stopWatch.stop();
        logger.debug("pre-processing comp-time: {}", stopWatch);
    }

分析代码，不难看出，是根据每个Job的Location计算与其他Job的Location的距离(这里用的欧几里得距离)，每个Job按照距离进行升序排列，选出来前capacity（车辆容量）个作为该Job的邻居
把Job的邻居也一起删除，添加到未分配任务列表中
把未分配的任务重新插入到当前的车辆路线中

Collection<Job> badJobs = insertUnassignedJobs(vehicleRoutes, unassignedJobs);

3、解的接受

作者实现了一个阈值接收器，SchrimpfAcceptance
下面是对这个接收器的解释
阈值接受函数：
这个概念可以描述如下：大多数问题不仅仅有一个唯一的最小值（或最大值），而是有多个局部最小值（或最大值）。为了避免在搜索开始时就陷入局部最小值，这种阈值接受函数在开始时也接受较差的解（与只接受更好解的贪婪方法相反），并随着时间的推移逐渐转变为贪婪方法。
难点：
定义（i）一个合适的初始阈值和（ii）一个描述阈值如何收敛到零的相应函数，即贪婪阈值。
初始阈值的确定：
通过在搜索空间中进行随机游走来确定初始阈值。随机游走使用特定的算法，并运行直到达到预热迭代次数。在第一次迭代或游走中，算法生成一个解，这个解是下一次游走的基础，以此类推。每个解的值都被记忆，因为初始阈值本质上是这些解值的标准差的函数。更精确地说：初始阈值 = 标准差(解的值) / 2。
具体实现阈值迭代的代码：

private double getThreshold(int iteration) {
    double scheduleVariable = (double) iteration / (double) maxIterations;
    return initialThreshold * Math.exp(-1. * Math.log(2) * scheduleVariable / alpha);
 }

试着通过画图来理解这个迭代
scheduleVariable 是个0,1之间的数，表示迭代的进度
alpha表示一个迭代参数
用python实现一下这个函数的图像

import numpy as np
import matplotlib.pyplot as plt

# 定义参数alpha和scheduleVariable的范围
alpha = 0.5  # 
scheduleVariable = np.linspace(0, 1, 100)  # 从0到1的100个点

# 计算函数值
y = np.exp(-np.log(2) / alpha * scheduleVariable)

# 绘制图像
plt.plot(scheduleVariable, y)
plt.xlabel('scheduleVariable')
plt.ylabel('y')
plt.title('Plot of the function y = e^(-(log(2)/alpha) * scheduleVariable)')
plt.grid(True)
plt.show()

在这里插入图片描述
随着迭代次数的增加，阈值逐渐递减

3.1 新解的接受策略

如果当前的解的个数小于阈值，则直接把新解保留
否则与(成本最高的解+阈值)比较，如果新解小于（成本最高的解+阈值），则删除成本最高的解，添加新解，更新接受标识
下面是代码的实现

public boolean acceptSolution(Collection<VehicleRoutingProblemSolution> solutions, VehicleRoutingProblemSolution newSolution) {
        boolean solutionAccepted = false;
        if (solutions.size() < solutionMemory) {
            solutions.add(newSolution);
            solutionAccepted = true;
        } else {
            VehicleRoutingProblemSolution worst = null;
            double threshold = getThreshold(currentIteration);
            for (VehicleRoutingProblemSolution solutionInMemory : solutions) {
                if (worst == null) worst = solutionInMemory;
                else if (solutionInMemory.getCost() > worst.getCost()) worst = solutionInMemory;
            }
            if (worst == null) {
                solutions.add(newSolution);
                solutionAccepted = true;
            } else if (newSolution.getCost() < worst.getCost() + threshold) {
                solutions.remove(worst);
                solutions.add(newSolution);
                solutionAccepted = true;
            }
        }
        return solutionAccepted;
    }

这段Java代码是用于决定是否接受一个新的解决方案到一个车辆路径问题（Vehicle Routing Problem, VRP）的解决方案集合中。
代码的主要逻辑如下：

定义一个布尔变量 solutionAccepted 来标识解决方案是否被接受。
检查 solutions 集合的大小是否小于 solutionMemory（一个代表解决方案记忆容量的变量）。如果是，将 newSolution 添加到集合中，并将 solutionAccepted 设置为 true。
如果 solutions 集合已满，则执行以下步骤：
- 定义一个变量 worst 来存储当前集合中最差的解决方案（成本最高）。
- 通过 getThreshold 方法获取当前迭代的阈值 threshold。
- 遍历 solutions 集合，找到成本最高的解决方案，并将其存储在 worst 变量中。
- 检查 worst 是否为 null。如果是，这意味着集合中没有解决方案，可以将 newSolution 添加到集合中，并将 solutionAccepted 设置为 true。
- 如果 newSolution 的成本小于 worst 的成本加上阈值 threshold，则从集合中移除 worst 解决方案，将 newSolution 添加到集合中，并将 solutionAccepted 设置为 true。
返回 solutionAccepted 变量，表示是否接受新的解决方案。