VRPSolverEasy:可求解多种VRP变体(rich vehicle routing)问题的精确算法python包

当前运筹优化领域的研究者在使用 branch-cut-and-price (BCP) 算法求解车辆路径问题 (VRP)的最优解方面取得了重大进展。BCP算法是VRPSolver中常用的算法,在许多 VRP变体的求解中表现优异。然而,其复杂的底层数学模型使其对于路径规划从业者(入门者)而言难以理解。为了解决这个问题,Najib Errami等人开发了VRPSolverEasy这个python工具包(Najib Errami, Eduardo Queiroga, Ruslan Sadykov, Eduardo Uchoa (2023) VRPSolverEasy: A Python Library for the Exact Solution of a Rich Vehicle Routing Problem. INFORMS Journal on Computing),VRPSolverEasy 提供了一个 Python 接口,无需任何混合整数规划建模知识即可使用的VRPSolver。在使用时,可以直接以VRP问题常用元素定义,例如仓库、客户、链接和车辆类型。VRPSolverEasy 可以处理多种流行的 VRP 变体及其任意组合。

1. VRPSolverEasy安装


VRPSolverEasy包提供了两种使用模式,一种是免费版本,该版本中使用了开源线性规划求解器COIN-OR CLP(果然开源才是王道,COIN-OR还有其它的开源求解器,具体可参考文章开源建模框架+开源求解器 | 使用pyomo建模框架实现交通物流优化问题的求解 (shortest path p),直接使用pip安装就可获取免费版本。
python3 -m pip install VRPSolverEasy

另一个版本是学术版本(其实是使用商业版本的CPLEX求解器),该版本提供了一些基于MIP的启发式方法,在求解性能上更强,更容易获取可行解,但是安装起来稍微有点麻烦,需要以下四个步骤:
  1. 将VRPSolverEasy python包克隆到本地
  2. 安装CPLEX
  3. 下载并安装Bapcod,Bapcod是基于C++写的branch-and-price包,也是VRPSolverEasy python包的求解内核,一直在保持持续更新,可能也是感觉使用C++怕别人觉得困难,也推出了python接口

在这里插入图片描述
4. 安装完后使用cmake进行编译

cmake --build build --config Release --target bapcod-shared

编译后会生成对应不同系统的三个库文件:

"Bapcod-path/build/Bapcod/libbapcod-shared.so" // Linux
"Bapcod-path/build/Bapcod/libbapcod-shared.dylib" // MacOs
"Bapcod-path/build/Bapcod/Release/bapcod-shared.dll" // Windows

用新生成的库文件替换VRPSolverEasy中旧的库文件(lib/),然后再重新安装VRPSolverEasy包:

python -m pip install

在使用VRPSolverEasy时,需进行参数设置solver_name=‘CPLEX’,如果使用的是Windows系统,参数设置时还需要明确CPLEX的路径,cplex_path=‘ ’。如果想要用下内置的启发式算法,添加heuristic_used=True即可。

2. VRPSolverEasy模块介绍


在安装完成后就可以使用VRPSolverEasy模块进行建模了
import VRPSolverEasy
model=VRPSolverRasy.Model()

VRPSolverEasy模型主要通过车场节点(depot points)、客户节点(customer points)、弧(link)和车辆类型(vehicle types)来进行定义。

2.1 添加车场节点(depot points)


model.add_depot(id=<id>,name='',service_time=0.0,
tw_begin=0.0,tw_end=0.0))

depot点的属性及解释如下
在这里插入图片描述

2.2 添加客户节点(customer points)


每个客户与一个或多个节点关联,每个节点(车场或客户)都应该有唯一的标识

model.add_customer(id=<id>,id_customer=<id>,name='',
demand=0,penalty=0.0,service_time=0.0,tw_begin=0.0,
tw_end=0.0,incompatible_vehicles=[])

此外,VRPSolverEasy还设置了添加额外节点的功能,用于表示可替代的时间窗和/或客户位置(实现一个客户对应多个节点,以便处理不同VRP问题变体的复杂情形),以及可能具有不同兼容车辆集的情况。

model.add_point(id=<id>,id_customer=<id>,name='',
demand=0,penalty=0.0,service_time=0.0,tw_begin=0.0,
tw_end=0.0,incompatible_vehicles=[])

customer点的属性及解释如下
在这里插入图片描述

2.3 添加弧(link)


这里link既可以是有向的arc(is_directed=True),也可以是无向的edge(is_directed=False),link的属性及解释如下
在这里插入图片描述

2.4 添加车辆类型(vehicle type)


model.add_vehicle_type(id=<id>,start_point_id=-1,end_point_id=-1,
name='',capacity=0,fixed_cost=0.0,var_cost_dist=0.0,
var_cost_time=0.0,max_number=1,tw_begin=0.0,tw_end=0.0)

vehicle type的属性及解释如下
在这里插入图片描述
此外,模型设置了总的最大车辆数(所以类型),默认取值为10000,可通过如下方式进行修改


model.set_max_total_vehicles_number(<value>)

路径成本的计算通过link成本的加和得到,link成本包括三部分:(1)link上的固定成本;(2)根据link的旅行时间和不同车型换算得到的成本;(3)根据link的距离和不同车型换算得到的成本。路径包含的属性及解释如下
在这里插入图片描述
在求解机制上,VRPSolverEasy首先将现有模型转化为VRPSolver可求解的模型(Pessoa et al.2020),然后再调用Bapcod包(Sadykov et al.2021)进行求解,求解时可设置的参数如下:

model.set_parameters(time_limit=300,upper_bound=1000000,heuristic_used=False,
time_limit_heuristic=20,solver_name='CLP',print_level=-1)

VRPSolverEasy求解的速度取决于upper_bound的取值,upper_bound越接近最优解,求解就越快,所以VRPSolverEasy实际上是先通过CLP或者CPLEX得到一个启发式的可行解作为上界,然后再进行求解,VRPSolverEasy求解参数设置的一些属性及解释如下:
在这里插入图片描述

在模型求解结束后,如果求到了解(model.solution.is_defined()==True),那么目标值可通过model.solution.value得到,路径信息可通过model.solution.routes得到。

3. VRPTW算例


设置若干客户点和车场节点构造vrptw的算例如下
在这里插入图片描述

# data
cost_per_distance = 10
begin_time = 0
end_time = 5000
nb_point = 7

# map with names and coordinates
coordinates = {"Wisconsin, USA": (44.50, -89.50),  # depot
               "West Virginia, USA": (39.000000, -80.500000),
               "Vermont, USA": (44.000000, -72.699997),
               "Texas, the USA": (31.000000, -100.000000),
               "South Dakota, the US": (44.500000, -100.000000),
               "Rhode Island, the US": (41.742325, -71.742332),
               "Oregon, the US": (44.000000, -120.500000)
               }

# demands of points
demands = [0, 500, 300, 600, 658, 741, 436]

模型构建


# Initialisation
model = solver.Model()

# Add vehicle type
model.add_vehicle_type(
    id=1,
    start_point_id=0,
    end_point_id=0,
    name="VEH1",
    capacity=1100,
    max_number=6,
    var_cost_dist=cost_per_distance,
    tw_end=5000)

# Add depot
model.add_depot(id=0, name="D1", tw_begin=0, tw_end=5000)

coordinates_keys = list(coordinates.keys())
# Add Customers
for i in range(1, nb_point):
    model.add_customer(
        id=i,
        name=coordinates_keys[i],
        demand=demands[i],
        tw_begin=begin_time,
        tw_end=end_time)

# Add links
coordinates_values = list(coordinates.values())
for i in range(0, 7):
    for j in range(i + 1, 7):
        dist = compute_euclidean_distance(coordinates_values[i][0],
                                          coordinates_values[j][0],
                                          coordinates_values[i][1],
                                          coordinates_values[j][1])
        model.add_link(
            start_point_id=i,
            end_point_id=j,
            distance=dist,
            time=dist)

求解

model.solve()

结果输出

if model.solution.is_defined():
    print(model.solution)
Route for vehicle 1:
    ID : 0 --> 2 --> 5 --> 0
    Name : DEPOT --> Vermont, USA --> Rhode Island, the US --> DEPOT
    End time : 0.0 --> 16.807 --> 19.259 --> 37.230000000000004
    Load : 0.0 --> 300.0 --> 1041.0 --> 1041.0
    Total cost : 372.29999999999995

Route for vehicle 1:
    ID : 0 --> 1 --> 3 --> 0
    Name : DEPOT --> West Virginia, USA --> Texas, the USA --> DEPOT
    End time : 0.0 --> 10.548 --> 31.625 --> 48.728
    Load : 0.0 --> 500.0 --> 1100.0 --> 1100.0
    Total cost : 487.2800000000001

Route for vehicle 1:
    ID : 0 --> 4 --> 6 --> 0
    Name : DEPOT --> South Dakota, the US --> Oregon, the US --> DEPOT
    End time : 0.0 --> 10.5 --> 31.006 --> 62.010000000000005
    Load : 0.0 --> 658.0 --> 1094.0 --> 1094.0
    Total cost : 620.1

在这里插入图片描述

4. VRPSolverEasy可求解的其它类型VRP部分算例


4.1 CVRP


容量限制车辆路径问题(CVRP)是最基础的 VRP,是指满足不同区域里的顾客需求的有容量限制的车辆路径规划问题。该问题的目标是最大限度降低总路径成本,即车辆行驶的总距离。


from VRPSolverEasy.src import solver
import os,sys,getopt
import math
import numpy as np

def read_instance(name : str):
    """ Read an instance in the folder data from a given name """
    path_project = os.path.abspath(os.getcwd())
    with open (name,
        "r",encoding="UTF-8" )as file:
        elements = [str(element) for element in file.read().split()]
    file.close()
    return elements

def compute_euclidean_distance(x_i, y_i, x_j, y_j,number_digit=3):
    """Compute the euclidean distance between 2 points from graph"""
    return round(math.sqrt((x_i - x_j)**2 +
                           (y_i - y_j)**2), number_digit)

def solve_demo(instance_path,
               solver_name="CLP",
               solver_path="",
               time_resolution=30,
               disableBuiltInHeur=False,
               upper_bound=-1):
    """return a solution from modelisation"""

    # read instance
    data = read_cvrp_instances(instance_path)

    # get data
    vehicle_type = data["VehicleTypes"]
    depot = data["Points"][0]
    customers = data["Points"][1:]
    links = data["Links"]

    # modelisation of problem
    model = solver.Model()

    # add vehicle type
    model.add_vehicle_type(id=vehicle_type["id"],
                           start_point_id=vehicle_type["start_point_id"],
                           end_point_id=vehicle_type["end_point_id"],
                           max_number=vehicle_type["max_number"],
                           capacity=vehicle_type["capacity"],
                           var_cost_dist=vehicle_type["var_cost_dist"]
                           )

    # add depot
    model.add_depot(id=depot["id"])

    # add all customers
    for customer in customers:
        model.add_customer(id=customer["id"],
                           demand=customer["demand"]
                           )
    # add all links
    for link in links:
        model.add_link(start_point_id=link["start_point_id"],
                       end_point_id=link["end_point_id"],
                       distance=link["distance"]
                       )

    # set parameters
    if disableBuiltInHeur:
        model.set_parameters(time_limit=time_resolution, solver_name=solver_name, 
                            heuristic_used=False)
        print('Built-in heuristic is disabled')
    else:
        model.set_parameters(time_limit=time_resolution, solver_name=solver_name, 
                            heuristic_used=True)
        print('Built-in heuristic is enabled')

    if upper_bound != -1:
        model.parameters.upper_bound = upper_bound

    ''' If you have cplex 22.1 installed on your laptop windows you have to specify
        solver path'''
    if (solver_name == "CPLEX" and solver_path != ""):
        model.parameters.cplex_path = solver_path

    # solve model
    model.solve()

    if model.solution.is_defined :
        path_instance_name = instance_path.split(".")[0]
        name_instance = path_instance_name.split("\\")[
            len(path_instance_name.split("\\"))-1]
        print('{0} {1} {2} {3} {4} {5} {6} {7} {8} {9}\n'.format(
                "instance_name","solver_name","ext_heuristic",
                "solution_value",
                "solution_time",
                "best_lb",
                "root_lb",
                "root_time",
                "nb_branch_and_bound_nodes",
                "status"
                ), end='')
        print('{0} {1} {2} {3} {4} {5} {6} {7} {8} {9}\n'.format(
            instance_path,solver_name,False if upper_bound == -1 else True,
            model.solution.value,
            model.statistics.solution_time,
            model.statistics.best_lb,
            model.statistics.root_lb,
            model.statistics.root_time,
            model.statistics.nb_branch_and_bound_nodes,
            model.status
            ))

    # export the result
    # model.solution.export(instance_name.split(".")[0] + "_result")

    return model.solution

def read_cvrp_instances(instance_name):
    """Read literature instances from CVRPLIB by giving the name of instance
       and returns dictionary containing all elements of model"""

    instance_iter = iter(read_instance(instance_name))
    points = []
    id_point = 0
    dimension_input = -1
    capacity_input = -1
    # Instantiate the data problem.
    data = {}

    while True:
        element = next(instance_iter)
        if element == "DIMENSION":
            next(instance_iter)  # pass ":"
            dimension_input = int(next(instance_iter))
        elif element == "CAPACITY":
            next(instance_iter)  # pass ":"
            capacity_input = int(next(instance_iter))
        elif element == "EDGE_WEIGHT_TYPE":
            next(instance_iter)  # pass ":"
            element = next(instance_iter)
            if element != "EUC_2D":
                raise Exception("EDGE_WEIGHT_TYPE : " + element + """
                is not supported (only EUC_2D)""")
        elif element == "NODE_COORD_SECTION":
            break

    # Initialize vehicle type
    vehicle_type = {"id": 1,  # we cannot have an id less than 1
                    "start_point_id": id_point,
                    "end_point_id": id_point,
                    "capacity": capacity_input,
                    "max_number": dimension_input,
                    "var_cost_dist": 1
                    }

    vehicles = []
    index = 0
    for i in range(dimension_input):
        vehicles.append(capacity_input)

    data['vehicle_capacities'] = vehicles
    data['vehicle_capacity'] = capacity_input
    data['num_vehicles'] = dimension_input
    data['depot'] = 0

    # Create points
    for current_id in range(dimension_input):
        point_id = int(next(instance_iter))
        if point_id != current_id + 1:
            raise Exception("Unexpected index")
        x_coord = float(next(instance_iter))
        y_coord = float(next(instance_iter))
        points.append({"x": x_coord,
                        "y": y_coord,
                        "demand": -1,
                        "id": id_point})
        id_point += 1

    element = next(instance_iter)
    if element != "DEMAND_SECTION":
        raise Exception("Expected line DEMAND_SECTION")
    jobs = []
    # Get the demands
    for current_id in range(dimension_input):
        point_id = int(next(instance_iter))
        if point_id != current_id + 1:
            raise Exception("Unexpected index")
        demand = int(next(instance_iter))
        points[current_id]["demand"] = demand
        jobs.append(demand)
    
    data['demands'] = jobs

    element = next(instance_iter)
    if element != "DEPOT_SECTION":
        raise Exception("Expected line DEPOT_SECTION")
    next(instance_iter) # pass id depot

    end_depot_section = int(next(instance_iter))
    if end_depot_section != -1:
        raise Exception("Expected only one depot.")

    # Compute the links of graph
    links = []
    nb_link = 0
    matrix = [[0 for i in range((len(points)))] for i in range(len(points))]
    for i, point in enumerate(points):
        for j in range(i + 1, len(points)):
            dist = compute_euclidean_distance(point["x"],
                                                    point["y"],
                                                    points[j]["x"],
                                                    points[j]["y"],
                                                    0)
            links.append({"name": "L" + str(nb_link),
                          "start_point_id": point["id"],
                          "end_point_id": points[j]["id"],
                          "distance": dist
                          })

            matrix[i][j] = dist
            matrix[j][i] = dist

            nb_link += 1
    
    data['distance_matrix'] = matrix

    return {"Points": points,
            "VehicleTypes": vehicle_type,
            "Links": links}

def main(argv):
    instance_path = ''
    solver_name = ''
    solver_path = ''
    disableBuiltInHeur = False
    time_resolution = 30
    upper_bound = -1
    opts, args = getopt.getopt(argv,"i:s:u:b:e:p:")
   
    for opt, arg in opts:
        if opt == "-i":
            instance_path = arg
        elif opt == "-s":
            solver_name = arg
        elif opt == "-u":
            upper_bound = float(arg)
        elif opt == "-b":
            disableBuiltInHeur = arg == "yes"
        elif opt == "-e":
            time_resolution = float(arg)
        if opt in ["-p"]:
            solver_path = arg

    solve_demo(instance_path, solver_name, solver_path, time_resolution, disableBuiltInHeur, upper_bound)

if __name__ == "__main__":
    if len(sys.argv) > 1:
        main(sys.argv[1:])

4.2 MDVRP


多车场路径规划问题(MDVRP)是 CVRP 的一种变体,其中车辆来自于不同车场。


import os
import math
import sys
import getopt
from VRPSolverEasy.src import solver


class DataMdvrp:
    """Contains all data for MDVRP problem
    """

    def __init__(
            self,
            nb_customers: int,
            nb_depots: int,
            vehicle_capacity: int,
            cust_demands=None,
            cust_coordinates=None,
            depot_coordinates=None):
        self.nb_customers = nb_customers
        self.nb_depots = nb_depots
        self.vehicle_capacity = vehicle_capacity
        self.cust_demands = cust_demands
        self.cust_coordinates = cust_coordinates
        self.depot_coordinates = depot_coordinates


def compute_euclidean_distance(x_i, y_i, x_j, y_j, number_digit=3):
    """Compute the euclidean distance between 2 points from graph"""
    return round(math.sqrt((x_i - x_j)**2 +
                           (y_i - y_j)**2), number_digit)


def read_instance(name: str):
    """ Read an instance in the folder data from a given name """

    with open(
            os.path.normpath(name),
            "r", encoding="UTF-8") as file:
        return [str(element) for element in file.read().split()]


def solve_demo(instance_name,
               time_resolution=30,
               solver_name_input="CLP",
               solver_path=""):
    """return a solution from modelisation"""

    # read parameters given in command line
    type_instance = "MDVRP/"
    if len(sys.argv) > 1:
        print(instance_name)
        opts = getopt.getopt(instance_name, "i:t:s:p:")
        for opt, arg in opts[0]:
            if opt in ["-i"]:
                instance_name = arg
                folder_data = ""
                type_instance = ""
            if opt in ["-t"]:
                time_resolution = float(arg)
            if opt in ["-s"]:
                solver_name_input = arg
            if opt in ["-p"]:
                solver_path = arg

    # read instance
    data = read_mdvrp_instances(instance_name)
    
    
    # modelisation of problem
    model = solver.Model()

    # add vehicle types
    for i in range(data.nb_depots):
        model.add_vehicle_type(id=data.nb_customers + i,
                               start_point_id=data.nb_customers + i,
                               end_point_id=data.nb_customers + i,
                               capacity=data.vehicle_capacity,
                               max_number=data.nb_customers,
                               var_cost_dist=1
                               )

    # add all customers
    for i in range(data.nb_customers):
        model.add_customer(id=i,
                           id_customer=i+1,
                           demand=data.cust_demands[i]
                           )
    # add depots
    for i in range(data.nb_customers,data.nb_customers + data.nb_depots):
        model.add_depot(id=i)



    nb_link = 0

    # Compute the links between depots and other points
    for index,coord_depot in enumerate(data.depot_coordinates):
        for i, cust_i in enumerate(data.cust_coordinates):
            dist = compute_euclidean_distance(
                cust_i[0],
                cust_i[1],
                coord_depot[0],
                coord_depot[1])
            model.add_link(start_point_id=index+data.nb_customers,
                           end_point_id=i,
                           distance=dist
                           )

    # Compute the links between points
    for i,cust_i in enumerate(data.cust_coordinates):
        for j in range(i + 1, len(data.cust_coordinates)):
            dist = compute_euclidean_distance(cust_i[0],
                                              cust_i[1],
                                              data.cust_coordinates[j][0],
                                              data.cust_coordinates[j][1])
            model.add_link(start_point_id=i,
                           end_point_id=j,
                           distance=dist
                           )

    # set parameters
    model.set_parameters(time_limit=time_resolution,
                         solver_name=solver_name_input)

    ''' If you have cplex 22.1 installed on your laptop windows you have to specify
        solver path'''
    if (solver_name_input == "CPLEX" and solver_path != ""):
        model.parameters.cplex_path = solver_path

    #model.export(instance_name,all_elements=True)

    # solve model
    model.solve()

    # export the result
    # model.solution.export(instance_name.split(".")[0] + "_result")

    return model.solution.value


def read_mdvrp_instances(instance_full_path):
    """Read literature instances of MDVRP by giving the name of instance
        and returns dictionary containing all elements of model"""
    instance_iter = iter(
        read_instance(instance_full_path))

    # pass type instance
    next(instance_iter)

    next(instance_iter)

    nb_customers = int(next(instance_iter))
    nb_depots = int(next(instance_iter))

    vehicle_capacities = []
    for _ in range(nb_depots):
        # pass max duration
        next(instance_iter)
        vehicle_capacities.append(int(next(instance_iter)))

    cust_demands = []
    cust_coordinates = []
    for _ in range(nb_customers):
        next(instance_iter)
        cust_x = int(next(instance_iter))
        cust_y = int(next(instance_iter))
        cust_coordinates.append([cust_x, cust_y])
        # pass duration
        next(instance_iter)

        cust_demand = int(next(instance_iter))
        cust_demands.append(cust_demand)
        for _ in range(6):
            next(instance_iter)

    depot_coordinates = []
    # add depots and vehicle types
    for _ in range(nb_depots):
        next(instance_iter)
        cust_x = int(next(instance_iter))
        cust_y = int(next(instance_iter))
        depot_coordinates.append([cust_x, cust_y])

        for _ in range(4):
            next(instance_iter)

    return DataMdvrp(nb_customers,
                     nb_depots,
                     vehicle_capacities[0],
                     cust_demands,
                     cust_coordinates,
                     depot_coordinates
                     )


if __name__ == "__main__":
    if len(sys.argv) > 1:
        solve_demo(sys.argv[1:])
    else:
        print("""Please indicates the parameters of your model like this : \n
       python MDVRP.py -i INSTANCE_PATH/NAME_INSTANCE \n
       -t TIME_RESOLUTION -s SOLVER_NAME (-p PATH_SOLVER (WINDOWS only))
       """)

4.3 RICHVRP


RichVRP 是 VRPSolverEasy 解决的特色问题之一,这里给出一个示例。它包括开放路线、多个仓库、多个时间窗口、两种类型的车辆(小型和大型)、可选客户和不兼容车辆(有些客户只可由小型车辆提供服务)。

import os
import math
import sys
import getopt
from VRPSolverEasy.src import solver

class Depot:
    def __init__(self, id, x, y, tw_begin, tw_end):
        self.id = id
        self.x = x
        self.y = y
        self.tw_begin = tw_begin
        self.tw_end = tw_end
    def __str__(self):
        return (f"Depot(id: {self.id}, x: {self.x}, y: {self.y}, "
                f"tw_begin: {self.tw_begin}, tw_end: {self.tw_end})")

class Vehicle:
    def __init__(self, capacity, fixed_cost, var_cost, max_number):
        self.capacity = capacity
        self.fixed_cost = fixed_cost
        self.var_cost = var_cost
        self.max_number = max_number
    def __str__(self):
        return (f"Vehicle(capacity: {self.capacity}, "
                f"fixed_cost: {self.fixed_cost}, var_cost: {self.var_cost}, "
                f"max_number: {self.max_number})")

class Customer:
    def __init__(self, id, x, y, demand, service_time, optional, only_small_veh, time_windows):
        self.id = id
        self.x = x
        self.y = y
        self.demand = demand
        self.service_time = service_time
        self.optional = optional # Is the customer optional?
        self.only_small_veh = only_small_veh # Can the customer served only by small vehicles? 
        self.time_windows = time_windows # List of time windows (list of pairs)
    def __str__(self):
        tw_str = ', '.join([f"({start}, {end})" for start, end in self.time_windows])
        return (f"Customer(id: {self.id}, x: {self.x}, y: {self.y}, "
                f"demand: {self.demand}, service_time: {self.service_time}, "
                f"optional: {self.optional}, only_small_veh: {self.only_small_veh}, "
                f"time_windows: [{tw_str}])")

class DataRichVrp:
    """Contains all data for RichVRP
    """
    def __init__(
            self,
            depots : None,
            big_vehicle : None,
            small_vehicle : None,
            customers : None):
        self.depots = depots
        self.big_vehicle = big_vehicle
        self.small_vehicle = small_vehicle
        self.customers = customers 
    def __str__(self):
        depots_str = '\n'.join([str(depot) for depot in self.depots])
        big_vehicle_str = str(self.big_vehicle)
        small_vehicle_str = str(self.small_vehicle)
        customers_str = '\n'.join([str(customer) for customer in self.customers])
        return (f"DataRichVrp:\n\nDepots:\n{depots_str}\n\nBig Vehicle:\n{big_vehicle_str}"
                f"\n\nSmall Vehicle:\n{small_vehicle_str}\n\nCustomers:\n{customers_str}")

def compute_euclidean_distance(x_i, y_i, x_j, y_j, number_digit=3):
    """Compute the euclidean distance between 2 points from graph"""
    return round(math.sqrt((x_i - x_j)**2 +
                           (y_i - y_j)**2), number_digit)

def read_instance(name: str):
    """ Read an instance in the folder data from a given name """
    with open(
            os.path.normpath(name),
            "r", encoding="UTF-8") as file:
        return [str(element) for element in file.read().split()]

def solve_demo(instance_name,
               time_resolution=30,
               solver_name_input="CLP",
               solver_path=""):
    """return a solution from modelisation"""

    # read parameters given in command line
    if len(sys.argv) > 1:
        print(instance_name)
        opts = getopt.getopt(instance_name, "i:t:s:p:")
        for opt, arg in opts[0]:
            if opt in ["-i"]:
                instance_name = arg
            if opt in ["-t"]:
                time_resolution = float(arg)
            if opt in ["-s"]:
                solver_name_input = arg
            if opt in ["-p"]:
                solver_path = arg

    # read instance
    data = read_richvrp_instance(instance_name)
    print(data)

    # modelisation of problem
    model = solver.Model()

    big_vehicle_id = 1 # Big vehicles have odd ids
    big_vehicle_ids = [] # Store all the ids for big vehicles
    small_vehicle_id = 2 # Small vehicles have even ids
    # Add depots and vehicles
    for depot in data.depots:
        # Add the depot
        model.add_depot(id = depot.id,
                        tw_begin = depot.tw_begin,
                        tw_end = depot.tw_end)
        # Add the big vehicle to this depot
        model.add_vehicle_type(id = big_vehicle_id,
                               start_point_id = depot.id,
                               end_point_id = -1, # A vehicle may end anywhere 
                               capacity = data.big_vehicle.capacity,
                               max_number = data.big_vehicle.max_number,
                               fixed_cost = data.big_vehicle.fixed_cost,
                               var_cost_dist = data.big_vehicle.var_cost,
                               tw_begin = depot.tw_begin,
                               tw_end = depot.tw_end)
        # Add the small vehicle to this depot
        model.add_vehicle_type(id = small_vehicle_id,
                               start_point_id = depot.id,
                               end_point_id = -1, # A vehicle may end anywhere 
                               capacity = data.small_vehicle.capacity,
                               max_number = data.small_vehicle.max_number,
                               fixed_cost = data.small_vehicle.fixed_cost,
                               var_cost_dist = data.small_vehicle.var_cost,
                               tw_begin = depot.tw_begin,
                               tw_end = depot.tw_end)
        big_vehicle_ids.append(big_vehicle_id)
        big_vehicle_id += 2
        small_vehicle_id += 2

    # ID for alternative points of a customer with multiple time windows
    next_id = max([customer.id for customer in data.customers]) + 1
    # IDs of a customer, including the alternative ones
    customer_ids = {customer.id : [customer.id] for customer in data.customers}
    # Add customers
    for customer in data.customers:
        # Add a point for each time window of a customer
        for i, tw in enumerate(customer.time_windows):
            point_id = customer.id
            if i > 0: # Is it an alternative point?
                point_id = next_id
                customer_ids[customer.id].append(next_id)
                next_id += 1 # Get the next alternative ID
            model.add_customer(id = point_id,
                                id_customer = customer.id,
                                demand = customer.demand,
                                service_time = customer.service_time,
                                tw_begin = tw[0],
                                tw_end = tw[1],
                                penalty = 1.0 if customer.optional else 0.0,
                                incompatible_vehicles = big_vehicle_ids if customer.only_small_veh else [])

    # Compute the links between depots and other points
    for depot in data.depots:
        for customer in data.customers:
            dist = compute_euclidean_distance(customer.x, customer.y, depot.x, depot.y)
            for point_id in customer_ids[customer.id]:
                model.add_link(start_point_id = depot.id,
                                end_point_id = point_id,
                                distance = dist,
                                time = dist)
                
    # Compute the links between customer points
    for i, c1 in enumerate(data.customers):
        for j, c2 in enumerate(data.customers):
            if j <= i:
                continue
            dist = compute_euclidean_distance(c1.x, c1.y, c2.x, c2.y)
            # Add a link for each pair of points from c1 to c2
            for point_id_c1 in customer_ids[c1.id]:
                for point_id_c2 in customer_ids[c2.id]:
                    model.add_link(start_point_id = point_id_c1,
                                    end_point_id = point_id_c2,
                                    distance = dist,
                                    time = dist)

    # set parameters
    model.set_parameters(time_limit=time_resolution,
                         solver_name=solver_name_input)

    ''' If you have cplex 22.1 installed on your laptop windows you have to specify
        solver path'''
    if (solver_name_input == "CPLEX" and solver_path != ""):
        model.parameters.cplex_path = solver_path

    # model.export(instance_name)

    # Solve model
    model.solve()

    print("\nCustomer IDs and time windows:")
    for customer in data.customers:
        ids_and_tws = []
        for i, point_id in enumerate(customer_ids[customer.id]):
            ids_and_tws.append(f"(id: {point_id}, tw: {list(customer.time_windows[i])})")
        print(f"Customer {customer.id}: {', '.join(ids_and_tws)}")

    if model.solution.is_defined():
        print(model.solution)

    # export the result
    # model.solution.export(instance_name.split(".")[0] + "_result")

    return model.solution.value

def read_richvrp_instance(instance_full_path):
    """Read instance of RichVRP by giving the name of instance
        and returns dictionary containing all elements of model"""
    instance_iter = iter(
        read_instance(instance_full_path))

    nb_depots = int(next(instance_iter))

    # Read depots
    depots = []
    for _ in range(nb_depots):
        id = int(next(instance_iter))
        x = int(next(instance_iter))
        y = int(next(instance_iter))
        tw_begin = int(next(instance_iter))
        tw_end = int(next(instance_iter))
        depots.append(Depot(id, x, y, tw_begin, tw_end))

    # Read big vehicle
    capacity_big_veh = int(next(instance_iter))
    fixed_cost_big_veh = int(next(instance_iter))
    var_cost_big_veh = float(next(instance_iter))
    max_numb_big_veh = int(next(instance_iter))
    big_vehicle = Vehicle(capacity_big_veh, fixed_cost_big_veh, var_cost_big_veh, max_numb_big_veh)

    # Read small vehicle
    capacity_small_veh = int(next(instance_iter))
    fixed_cost_small_veh = int(next(instance_iter))
    var_cost_small_veh = float(next(instance_iter))
    max_numb_small_veh = int(next(instance_iter))
    small_vehicle = Vehicle(capacity_small_veh, fixed_cost_small_veh, var_cost_small_veh, max_numb_small_veh)

    # Read customers
    nb_customers = int(next(instance_iter))
    customers = []
    for _ in range(nb_customers):
        id = int(next(instance_iter))
        x = int(next(instance_iter))
        y = int(next(instance_iter))
        demand = int(next(instance_iter))
        service_time = int(next(instance_iter))
        optional = int(next(instance_iter)) == 1 if True else False
        small_vehicle_flag = int(next(instance_iter)) == 1 if True else False
        nb_time_windows = int(next(instance_iter))
        time_windows = []
        for _ in range(nb_time_windows):
            tw_begin = int(next(instance_iter))
            tw_end = int(next(instance_iter))
            time_windows.append((tw_begin, tw_end))
        customers.append(Customer(id, x, y, demand, service_time, optional, small_vehicle_flag, time_windows))

    return DataRichVrp(depots, big_vehicle, small_vehicle, customers)

if __name__ == "__main__":
    if len(sys.argv) > 1:
        solve_demo(sys.argv[1:])
    else:
        print("""Please indicates the parameters of your model like this : \n
       python RichVRP.py -i INSTANCE_PATH/NAME_INSTANCE \n
       -t TIME_RESOLUTION -s SOLVER_NAME (-p PATH_SOLVER (WINDOWS only))
       """)


5. 性能


毕竟是精确算法求解器,在大规模问题上还是拼不过启发式算法,对于200 个以上客户的算例得到最优解可能需要很长时间(数小时或数天),但少于 100 个客户的算例一般可以在几分钟内得到解决。

评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值