{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "# 容量约束弧路径问题(CARP)\n", "\n", "**类别：** 路径优化\n", "\n", "使用 OptAgent 的 Python 接口描述变量、约束与目标。\n", "\n", "问题与原始示例来源：[Hexaly Code Templates](https://www.hexaly.com/templates/capacitated-arc-routing-problem-carp)。\n" ] }, { "cell_type": "markdown", "id": "2a9d7baa", "metadata": {}, "source": [ "## 问题描述\n", "\n", "**在容量约束弧路径问题(Capacitated Arc Routing Problem, CARP)**中,一组具有相同容量的配送车辆必须为具有已知需求的边提供服务。车辆从一个共同的配送中心出发并返回。每条边恰好由一辆车服务。此外,每辆卡车服务的总需求量不得超过其容量。目标是最小化所行驶的总距离。\n", "\n", "### 学习要点\n", "\n", "- 添加 列表决策变量 以建模每辆卡车的边序列\n", "- 在所有列表变量上添加 `disjoint` 约束\n", "- 使用 `contains` 保证每条必服务边的两个方向中恰好有一个被选择\n", "- 定义 lambda 函数 来计算行驶距离\n" ] }, { "cell_type": "markdown", "id": "b170f180", "metadata": {}, "source": [ "## 数据\n", "\n", "所提供的容量约束弧路径问题(CARP)算例来自 [DIMACS 网站](http://dimacs.rutgers.edu/programs/challenge/vrp/carp/)。数据文件的格式如下:\n", "\n", "- 节点数量\n", "- 必服务边的数量(具有正需求的边)\n", "- 非必服务边的数量(具有零需求的边)\n", "- 可用车辆数量\n", "- 卡车容量\n", "- 每条必服务边的需求与费用\n", "- 每条非必服务边的费用\n", "- 配送中心节点的索引\n" ] }, { "cell_type": "markdown", "id": "3656eb4f", "metadata": {}, "source": [ "## 建模思路\n", "\n", "容量约束弧路径问题(CARP)的 OptAgent 模型使用列表变量来表示分配给每辆卡车的边序列。这些边是卡车访问并提供服务的边。卡车在其路径上也可能仅经过其他边而不提供服务。\n", "\n", "边可以沿两个方向被访问,但其需求只能被满足一次。OptAgent 的 `disjoint` 约束确保每条边(沿任一方向)至多出现在一个列表变量中,而 `contains(edges_sequences, 2i) + contains(edges_sequences, 2i+1) == 1` 保证每对相反方向中恰好有一条被服务。\n", "\n", "`contains(edges_sequences, edge)` 判断某个方向的弧是否出现在任意车辆路线中。\n", "\n", "我们可以使用需求数组上的 **at** 算子来访问序列中每条边的需求。我们使用一个 lambda 函数，通过 **sum** 算子对所有被访问边的需求求和，从而计算每辆卡车的总配送量。该总配送量必须不超过卡车的容量。\n", "\n", "类似地，我们使用二维距离矩阵上的 **at** 算子来访问从一条边到下一条边所行驶的距离。\n" ] }, { "cell_type": "markdown", "id": "936716cc", "metadata": {}, "source": [ "## Python 实现\n" ] }, { "cell_type": "code", "execution_count": null, "id": "carp-model", "metadata": {}, "outputs": [], "source": [ "from pathlib import Path\n", "import re\n", "\n", "from optagent import OptModel, solve\n", "\n", "\n", "\n", "\n", "class CarpInstance:\n", " \"\"\"CARP instance parsed from the DIMACS challenge format.\"\"\"\n", "\n", " EDGE_PATTERN = re.compile(\n", " r\"\\(\\s*(\\d+)\\s*,\\s*(\\d+)\\s*\\)\\s+coste\\s+(\\d+)\"\n", " r\"(?:\\s+demanda\\s+(\\d+))?\"\n", " )\n", "\n", " def __init__(self, filename):\n", " lines = Path(filename).read_text(encoding=\"utf-8\").splitlines()\n", " nb_nodes = self._header_int(lines, \"VERTICES\")\n", " self.nb_required_edges = self._header_int(lines, \"ARISTAS_REQ\")\n", " nb_not_required_edges = self._header_int(lines, \"ARISTAS_NOREQ\")\n", " self.nb_trucks = self._header_int(lines, \"VEHICULOS\")\n", " self.truck_capacity = self._header_int(lines, \"CAPACIDAD\")\n", " depot = self._header_int(lines, \"DEPOSITO\")\n", "\n", " required_edges = self._read_edges(\n", " lines, \"LISTA_ARISTAS_REQ\", self.nb_required_edges, required=True\n", " )\n", " other_edges = self._read_edges(\n", " lines, \"LISTA_ARISTAS_NOREQ\", nb_not_required_edges, required=False\n", " )\n", "\n", " # Match the Hexaly model: an undirected adjacency matrix with 0 for\n", " # absent edges, followed by the same O(V^2) Dijkstra preprocessing.\n", " graph = [[0] * (nb_nodes + 1) for _ in range(nb_nodes + 1)]\n", " for origin, destination, cost, _ in required_edges + other_edges:\n", " graph[origin][destination] = cost\n", " graph[destination][origin] = cost\n", "\n", " self.costs_data = []\n", " self.demands_data = []\n", " self.origins_data = []\n", " self.destinations_data = []\n", " for origin, destination, cost, demand in required_edges:\n", " self.costs_data.extend((cost, cost))\n", " self.demands_data.extend((demand, demand))\n", " self.origins_data.extend((origin, destination))\n", " self.destinations_data.extend((destination, origin))\n", "\n", " relevant_nodes = set(self.origins_data + self.destinations_data + [depot])\n", " shortest_paths = {\n", " node: self._shortest_paths(node, graph) for node in relevant_nodes\n", " }\n", " self.edges_dist_data = [\n", " [\n", " shortest_paths[self.destinations_data[i]][self.origins_data[j]]\n", " for j in range(2 * self.nb_required_edges)\n", " ]\n", " for i in range(2 * self.nb_required_edges)\n", " ]\n", " self.dist_from_depot_data = [\n", " shortest_paths[depot][origin] for origin in self.origins_data\n", " ]\n", " self.dist_to_depot_data = [\n", " shortest_paths[destination][depot] for destination in self.destinations_data\n", " ]\n", "\n", " @staticmethod\n", " def _header_int(lines, key):\n", " return int(\n", " next(\n", " line.split(\":\", 1)[1] for line in lines if line.strip().startswith(key)\n", " )\n", " )\n", "\n", " @classmethod\n", " def _read_edges(cls, lines, section, count, required):\n", " if count == 0:\n", " return []\n", " start = (\n", " next(\n", " index\n", " for index, line in enumerate(lines)\n", " if line.strip().startswith(section)\n", " )\n", " + 1\n", " )\n", " edges = []\n", " for line in lines[start : start + count]:\n", " match = cls.EDGE_PATTERN.search(line)\n", " if match is None:\n", " raise ValueError(f\"Invalid edge line: {line}\")\n", " origin, destination, cost = map(int, match.group(1, 2, 3))\n", " demand = int(match.group(4)) if required else 0\n", " edges.append((origin, destination, cost, demand))\n", " return edges\n", "\n", " @staticmethod\n", " def _shortest_paths(start, graph):\n", " distances = [float(\"inf\")] * len(graph)\n", " visited = [False] * len(graph)\n", " distances[start] = 0\n", " for _ in range(len(graph)):\n", " unvisited = (\n", " index for index, distance in enumerate(distances)\n", " if not visited[index]\n", " )\n", " node = min(unvisited, key=distances.__getitem__, default=None)\n", " if node is None or distances[node] == float(\"inf\"):\n", " break\n", " visited[node] = True\n", " for neighbor, edge_cost in enumerate(graph[node]):\n", " if edge_cost and not visited[neighbor]:\n", " candidate = distances[node] + edge_cost\n", " if candidate < distances[neighbor]:\n", " distances[neighbor] = candidate\n", " return distances\n", "\n", "\n", "def build_carp_model(\n", " nb_required_edges,\n", " nb_trucks,\n", " truck_capacity,\n", " costs_data,\n", " demands_data,\n", " edges_dist_data,\n", " dist_from_depot_data,\n", " dist_to_depot_data,\n", "):\n", " \"\"\"Construct the CARP model with directed-edge list variables.\"\"\"\n", " model = OptModel()\n", "\n", " # Seed the search with a capacity-feasible assignment, as in the\n", " # original Hexaly model's initial route construction.\n", " default_sequences = [[] for _ in range(nb_trucks)]\n", " default_loads = [0] * nb_trucks\n", " edge_order = sorted(\n", " range(nb_required_edges),\n", " key=lambda edge: demands_data[2 * edge],\n", " reverse=True,\n", " )\n", " for edge in edge_order:\n", " demand = demands_data[2 * edge]\n", " candidates = [\n", " truck\n", " for truck in range(nb_trucks)\n", " if default_loads[truck] + demand <= truck_capacity\n", " ]\n", " if not candidates:\n", " raise ValueError(\"Unable to build a capacity-feasible initial assignment\")\n", " # Best-fit keeps the remaining capacity usable for tight instances.\n", " truck = min(\n", " candidates,\n", " key=lambda candidate: truck_capacity\n", " - default_loads[candidate]\n", " - demand,\n", " )\n", " default_sequences[truck].append(2 * edge)\n", " default_loads[truck] += demand\n", "\n", " # Each item is one orientation of a required edge: 2*i is direct, 2*i+1 reverse.\n", " edges_sequences_vars = [\n", " model.list(\n", " 2 * nb_required_edges,\n", " default=tuple(default_sequences[truck]),\n", " )\n", " for truck in range(nb_trucks)\n", " ]\n", " edges_sequences = model.array(edges_sequences_vars)\n", " model.constraint(model.disjoint(edges_sequences))\n", "\n", " for edge in range(nb_required_edges):\n", " model.constraint(\n", " model.contains(edges_sequences, 2 * edge)\n", " + model.contains(edges_sequences, 2 * edge + 1)\n", " == 1,\n", " )\n", "\n", " costs_array = model.array(costs_data)\n", " demands_array = model.array(demands_data)\n", " dist_from_depot_array = model.array(dist_from_depot_data)\n", " dist_to_depot_array = model.array(dist_to_depot_data)\n", " edges_dist_array = model.array(edges_dist_data)\n", "\n", " route_distances = []\n", " for k in range(nb_trucks):\n", " sequence = edges_sequences_vars[k]\n", " c = model.count(sequence)\n", "\n", " # Truck capacity.\n", " demand_lambda = model.lambda_function(lambda edge: demands_array[edge])\n", " route_quantity = model.sum(sequence, demand_lambda)\n", " model.constraint(route_quantity <= truck_capacity)\n", "\n", " # Service cost plus travel from the preceding serviced directed edge.\n", " dist_lambda = model.lambda_function(\n", " lambda i: costs_array[sequence[i]]\n", " + edges_dist_array[sequence[i - 1], sequence[i]]\n", " )\n", " route_dist = model.sum(model.range(1, c), dist_lambda) + model.iif(\n", " c > 0,\n", " costs_array[sequence[0]]\n", " + dist_from_depot_array[sequence[0]]\n", " + dist_to_depot_array[sequence[c - 1]],\n", " 0,\n", " )\n", " route_distances.append(route_dist)\n", "\n", " total_distance = model.sum(*route_distances)\n", " model.minimize(total_distance)\n", " return model, edges_sequences_vars, total_distance\n", "\n", "\n", "def main(input_file, output_file=None, time_limit=10):\n", " instance = CarpInstance(input_file)\n", " model, edges_sequences_vars, total_distance = build_carp_model(\n", " nb_required_edges=instance.nb_required_edges,\n", " nb_trucks=instance.nb_trucks,\n", " truck_capacity=instance.truck_capacity,\n", " costs_data=instance.costs_data,\n", " demands_data=instance.demands_data,\n", " edges_dist_data=instance.edges_dist_data,\n", " dist_from_depot_data=instance.dist_from_depot_data,\n", " dist_to_depot_data=instance.dist_to_depot_data,\n", " )\n", " solution = solve(model, time_limit_s=float(time_limit))\n", " if not solution.feasible:\n", " print(\n", " f\"Required edges = {instance.nb_required_edges}; Trucks = {instance.nb_trucks}; \"\n", " f\"Status = {solution.feasible}; no feasible route was found\"\n", " )\n", " return solution\n", " result_values = {'total_distance': total_distance.value, **{f'truck_{truck}': sequence.value for truck, sequence in enumerate(edges_sequences_vars)}}\n", "\n", " lines = [\n", " f\"Required edges = {instance.nb_required_edges}; Trucks = {instance.nb_trucks}; \"\n", " f\"Total distance = {result_values['total_distance']}; \"\n", " f\"Status = {solution.feasible}\"\n", " ]\n", " for truck in range(instance.nb_trucks):\n", " sequence = result_values[f\"truck_{truck}\"]\n", " if sequence:\n", " serviced_edges = [\n", " (instance.origins_data[edge], instance.destinations_data[edge])\n", " for edge in sequence\n", " ]\n", " lines.append(f\"Truck {truck + 1}: {serviced_edges}\")\n", "\n", " result_text = \"\\n\".join(lines)\n", " print(result_text)\n", " if output_file is not None:\n", " Path(output_file).write_text(result_text + \"\\n\", encoding=\"utf-8\")\n", " return solution\n" ] }, { "cell_type": "markdown", "id": "carp-runs", "metadata": {}, "source": [ "## 本地运行\n", "\n", "Notebook 直接调用 `main` 并显式传入实例路径。以下代码格相互独立，可以按需要单独运行；调整 `time_limit` 可以控制每个实例的求解时间。\n" ] }, { "cell_type": "code", "execution_count": null, "id": "carp-instance-dir", "metadata": {}, "outputs": [], "source": [ "from pathlib import Path\n", "\n", "INSTANCE_DIR = Path.cwd() / \"instances\"\n", "print(\"Instances:\", INSTANCE_DIR)\n" ] }, { "cell_type": "code", "execution_count": null, "id": "carp-run-e1-a", "metadata": {}, "outputs": [], "source": [ "solution_egl_e1_a = main(INSTANCE_DIR / \"egl-e1-A.dat\", time_limit=1)\n" ] }, { "cell_type": "code", "execution_count": null, "id": "carp-run-e1-b", "metadata": {}, "outputs": [], "source": [ "solution_egl_e1_b = main(INSTANCE_DIR / \"egl-e1-B.dat\", time_limit=1)\n" ] }, { "cell_type": "code", "execution_count": null, "id": "carp-run-e1-c", "metadata": {}, "outputs": [], "source": [ "solution_egl_e1_c = main(INSTANCE_DIR / \"egl-e1-C.dat\", time_limit=1)\n" ] } ], "metadata": { "kernelspec": { "display_name": ".venv", "language": "python", "name": "python3" }, "language_info": { "codemirror_mode": { "name": "ipython", "version": 3 }, "file_extension": ".py", "mimetype": "text/x-python", "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", "version": "3.14.6" } }, "nbformat": 4, "nbformat_minor": 5 }