{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# 选址路径问题 (LRP)\n",
    "\n",
    "**类别：** 路径优化\n",
    "\n",
    "使用 OptAgent 的 Python 接口描述变量、约束与目标。\n",
    "\n",
    "问题与原始示例来源：[Hexaly Code Templates](https://www.hexaly.com/templates/location-routing-problem-lrp)。\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 问题描述\n",
    "\n",
    "**在选址路径问题 (LRP)** 中，它是带容量车辆路径问题 ([CVRP](https://optagent.pages.dev/examples/routing/capacitated_vehicle_routing_problem_cvrp/)) 的扩展，一组具有统一容量的配送车辆必须为具有单一商品已知需求的客户提供服务。与仅有一个配送中心位置的 CVRP 不同，LRP 有多个可用的配送中心，每个配送中心都有自己的开放成本和容量。车辆从同一个配送中心出发并返回。每个客户必须恰好被一辆车辆服务，并且每辆车和每个配送中心服务的总需求不能超过其容量。目标是最小化总成本，即路径长度之和加上配送中心和路径的开放成本之和。\n",
    "\n",
    "### 学习要点\n",
    "\n",
    "- 添加 [列表决策变量](https://optagent.pages.dev/guide/modeling/) 来建模卡车的客户序列\n",
    "- 添加 [集合决策变量](https://optagent.pages.dev/guide/modeling/) 来建模配送中心关联的卡车\n",
    "- 通过 '[`find`](https://optagent.pages.dev/api/model-reference/#optagent.OptModel.find)' 算子获取每辆卡车的配送中心的索引\n",
    "- 定义 [符号 Lambda 函数](https://optagent.pages.dev/guide/modeling/) 来计算行驶距离\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 数据\n",
    "\n",
    "我们提供的选址路径问题 (LRP) 实例来自 [S. Barreto instances](http://prodhonc.free.fr/Instances/instances_us.htm)。数据文件的格式如下：\n",
    "\n",
    "- 客户数\n",
    "- 配送中心数\n",
    "- 配送中心和客户的 x、y 坐标\n",
    "- 配送车辆的容量\n",
    "- 每个配送中心的容量\n",
    "- 每个客户的需求\n",
    "- 每个配送中心的开放成本\n",
    "- 路径的开放成本\n",
    "- 一个布尔值，指示是否应对成本取整。\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 建模思路\n",
    "\n",
    "选址路径问题 (LRP) 的OptAgent模型使用 列表决策变量。对于每辆卡车，我们定义一个 list variable 表示其访问的客户的序列。在所有列表上使用 **partition** 约束，我们确保每个客户恰好被一辆卡车服务。此外，我们还使用 集合决策变量。对于每个配送中心，我们定义一个 set variable 表示在该配送中心出发和结束路径的卡车。在 set 上使用另一个 **partition** 约束，我们确保每辆卡车与一个配送中心关联。\n",
    "\n",
    "每辆卡车的总配送量通过 [符号 Lambda 函数](https://optagent.pages.dev/guide/modeling/) 对所有访问的客户应用 **sum** 算子计算得到。注意，此求和中项的数量在搜索过程中会随着列表的大小变化而变化。然后我们可以将该数量约束为小于卡车的容量。\n",
    "\n",
    "使用 **find** 算子，我们可以检索与每辆卡车关联的配送中心的索引。这使我们能够计算从每辆卡车的配送中心到其第一个客户，以及从其最后一个客户返回其配送中心的距离。为了计算每条路径的总长度，我们还需要知道客户到客户之间的距离：我们使用另一个 [符号 Lambda 函数](https://optagent.pages.dev/guide/modeling/) 来对沿途的距离求和。\n",
    "\n",
    "如果一条路径至少服务一个客户，则认为该路径是开放的；如果至少有一个非空序列使用该配送中心，则认为该配送中心是开放的。这些条件通过 **count** 算子计算。\n",
    "\n",
    "最后，我们可以通过对所有路径的总长度以及每条开放路径和每个开放配送中心的开放成本求和来计算目标函数。\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Python 实现\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "f0c91f42",
   "metadata": {},
   "outputs": [],
   "source": [
    "import math\n",
    "from pathlib import Path\n",
    "\n",
    "from optagent import OptModel, solve\n",
    "\n",
    "\n",
    "def read_elem(filename):\n",
    "    with open(filename) as f:\n",
    "        return [str(elem) for elem in f.read().split()]\n",
    "\n",
    "\n",
    "def main(instance_file, output_file=None, time_limit=20):\n",
    "    #\n",
    "    # Read instance data\n",
    "    #\n",
    "    nb_customers, nb_depots, vehicle_capacity, opening_route_cost, demands_data, \\\n",
    "        capacity_depots, opening_depots_cost, dist_matrix_data, dist_depots_data = \\\n",
    "        read_input_lrp(instance_file)\n",
    "\n",
    "    min_nb_trucks = int(math.ceil(sum(demands_data) / vehicle_capacity))\n",
    "    nb_trucks = int(math.ceil(1.5 * min_nb_trucks))\n",
    "\n",
    "    model = OptModel()\n",
    "\n",
    "    # Each route is a list of customers; partition assigns every customer once.\n",
    "    customers_sequences = [\n",
    "        model.list(nb_customers)\n",
    "        for r in range(nb_trucks)\n",
    "    ]\n",
    "    model.constraint(model.partition(customers_sequences))\n",
    "\n",
    "    # Each depot is a set of routes that start and end at that depot.\n",
    "    depots = [model.set(nb_trucks) for d in range(nb_depots)]\n",
    "    depots_array = model.array(depots)\n",
    "    model.constraint(model.partition(depots))\n",
    "\n",
    "    demands = model.array(demands_data)\n",
    "    dist_matrix = model.array(dist_matrix_data)\n",
    "    dist_depots = model.array(dist_depots_data)\n",
    "\n",
    "    sequence_used = []\n",
    "    associated_depot = []\n",
    "    quantity_served = []\n",
    "    route_distances = []\n",
    "    route_costs = []\n",
    "    for r, sequence in enumerate(customers_sequences):\n",
    "        count = model.count(sequence)\n",
    "        used = count > 0\n",
    "        depot = model.find(depots_array, r)\n",
    "        demand_lambda = model.lambda_function(lambda customer: demands[customer // 1])\n",
    "        quantity = model.sum(sequence, demand_lambda)\n",
    "        model.constraint(quantity <= vehicle_capacity)\n",
    "\n",
    "        distance_lambda = model.lambda_function(\n",
    "            lambda position: dist_matrix[sequence[position // 1], sequence[position // 1 + 1]]\n",
    "        )\n",
    "        distance = model.sum(model.range(0, count - 1), distance_lambda) + model.iif(\n",
    "            used,\n",
    "            dist_depots[sequence[0], depot] + dist_depots[sequence[count - 1], depot],\n",
    "            0,\n",
    "        )\n",
    "        sequence_used.append(used)\n",
    "        associated_depot.append(depot)\n",
    "        quantity_served.append(quantity)\n",
    "        route_distances.append(distance)\n",
    "        route_costs.append(used * opening_route_cost + distance)\n",
    "\n",
    "    quantity_served_array = model.array(quantity_served)\n",
    "    depot_costs = []\n",
    "    for d, depot in enumerate(depots):\n",
    "        depot_open = model.count(depot) > 0\n",
    "        depot_costs.append(depot_open * opening_depots_cost[d])\n",
    "        depot_lambda = model.lambda_function(lambda route: quantity_served_array[route // 1])\n",
    "        depot_quantity = model.sum(depot, depot_lambda)\n",
    "        model.constraint(depot_quantity <= capacity_depots[d])\n",
    "\n",
    "    total_cost = model.sum(route_costs) + model.sum(depot_costs)\n",
    "    model.minimize(total_cost)\n",
    "\n",
    "    solution = solve(model, time_limit_s=float(time_limit))\n",
    "    values = {\n",
    "        \"total_cost\": total_cost.value,\n",
    "        **{f\"route_{r}\": sequence.value for r, sequence in enumerate(customers_sequences)},\n",
    "        **{f\"depot_{r}\": depot.value for r, depot in enumerate(associated_depot)},\n",
    "    }\n",
    "    lines = [\n",
    "        f\"Customers = {nb_customers}; Depots = {nb_depots}; Trucks = {nb_trucks}; \"\n",
    "        f\"Total cost = {values['total_cost']}; Status = {solution.feasible}\"\n",
    "    ]\n",
    "    for r in range(nb_trucks):\n",
    "        route = values[f\"route_{r}\"]\n",
    "        if route:\n",
    "            lines.append(f\"Route {r} (depot {values[f'depot_{r}']}): {' '.join(map(str, route))}\")\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",
    "\n",
    "\n",
    "def read_input_lrp_dat(filename):\n",
    "    file_it = iter(read_elem(filename))\n",
    "\n",
    "    nb_customers = int(next(file_it))\n",
    "    nb_depots = int(next(file_it))\n",
    "\n",
    "    x_depot = [None] * nb_depots\n",
    "    y_depot = [None] * nb_depots\n",
    "    for i in range(nb_depots):\n",
    "        x_depot[i] = int(next(file_it))\n",
    "        y_depot[i] = int(next(file_it))\n",
    "\n",
    "    x_customer = [None] * nb_customers\n",
    "    y_customer = [None] * nb_customers\n",
    "    for i in range(nb_customers):\n",
    "        x_customer[i] = int(next(file_it))\n",
    "        y_customer[i] = int(next(file_it))\n",
    "\n",
    "    vehicle_capacity = int(next(file_it))\n",
    "    capacity_depots = [None] * nb_depots\n",
    "    for i in range(nb_depots):\n",
    "        capacity_depots[i] = int(next(file_it))\n",
    "\n",
    "    demands = [None] * nb_customers\n",
    "    for i in range(nb_customers):\n",
    "        demands[i] = int(next(file_it))\n",
    "\n",
    "    temp_opening_cost_depot = [None] * nb_depots\n",
    "    for i in range(nb_depots):\n",
    "        temp_opening_cost_depot[i] = float(next(file_it))\n",
    "    temp_opening_route_cost = int(next(file_it))\n",
    "    are_cost_double = int(next(file_it))\n",
    "\n",
    "    opening_depots_cost = [None] * nb_depots\n",
    "    if are_cost_double == 1:\n",
    "        opening_depots_cost = temp_opening_cost_depot\n",
    "        opening_route_cost = temp_opening_route_cost\n",
    "    else:\n",
    "        opening_route_cost = round(temp_opening_route_cost)\n",
    "        for i in range(nb_depots):\n",
    "            opening_depots_cost[i] = round(temp_opening_cost_depot[i])\n",
    "\n",
    "    distance_customers = compute_distance_matrix(x_customer, y_customer, are_cost_double)\n",
    "    distance_customers_depots = compute_distance_depot(x_customer, y_customer,\n",
    "                                                       x_depot, y_depot, are_cost_double)\n",
    "\n",
    "    return nb_customers, nb_depots, vehicle_capacity, opening_route_cost, demands, \\\n",
    "        capacity_depots, opening_depots_cost, distance_customers, distance_customers_depots\n",
    "\n",
    "# Compute the distance matrix\n",
    "def compute_distance_matrix(customers_x, customers_y, are_cost_double):\n",
    "    nb_customers = len(customers_x)\n",
    "    dist_customers = [[None for _ in range(nb_customers)] for _ in range(nb_customers)]\n",
    "    for i in range(nb_customers):\n",
    "        dist_customers[i][i] = 0\n",
    "        for j in range(nb_customers):\n",
    "            dist = compute_dist(customers_x[i], customers_x[j],\n",
    "                                customers_y[i], customers_y[j], are_cost_double)\n",
    "            dist_customers[i][j] = dist\n",
    "            dist_customers[j][i] = dist\n",
    "    return dist_customers\n",
    "\n",
    "# Compute the distance depot matrix\n",
    "def compute_distance_depot(customers_x, customers_y, depot_x, depot_y, are_cost_double):\n",
    "    nb_customers = len(customers_x)\n",
    "    nb_depots = len(depot_x)\n",
    "    distance_customers_depots = [[None for _ in range(nb_depots)] for _ in range(nb_customers)]\n",
    "    for i in range(nb_customers):\n",
    "        for d in range(nb_depots):\n",
    "            dist = compute_dist(customers_x[i], depot_x[d],\n",
    "                                customers_y[i], depot_y[d], are_cost_double)\n",
    "            distance_customers_depots[i][d] = dist\n",
    "    return distance_customers_depots\n",
    "\n",
    "\n",
    "def compute_dist(xi, xj, yi, yj, are_cost_double):\n",
    "    dist = math.sqrt(math.pow(xi - xj, 2) + math.pow(yi - yj, 2))\n",
    "    if are_cost_double == 0:\n",
    "        dist = math.ceil(100 * dist)\n",
    "    return dist\n",
    "\n",
    "\n",
    "def read_input_lrp(filename):\n",
    "    if Path(filename).suffix.lower() == \".dat\":\n",
    "        return read_input_lrp_dat(filename)\n",
    "    raise ValueError(f\"Unknown file format: {Path(filename).suffix}\")\n",
    "\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "752b7b49",
   "metadata": {},
   "source": [
    "## 本地运行\n",
    "\n",
    "以下代码格演示如何调用 OptAgent 的 LRP 模型。\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "61097c75",
   "metadata": {},
   "outputs": [],
   "source": [
    "INSTANCE_DIR = Path.cwd() / \"instances\"\n",
    "print(\"Instances:\", INSTANCE_DIR)\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "bb5a8360",
   "metadata": {},
   "outputs": [],
   "source": [
    "solution_gaspelle = main(INSTANCE_DIR / \"coordGaspelle.dat\", time_limit=1)\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "solution_christ50 = main(INSTANCE_DIR / \"coordChrist50.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
}