authors/gomfy/dev/5_benders_parallel_stoch_floc.ipynb

{
  "cells": [
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "# AMPL Development Tutorial 5/6 -- Parallelizing Subproblem Solves in Benders Decomposition\n",
        "[![5_benders_parallel_stoch_floc.ipynb](https://img.shields.io/badge/github-%23121011.svg?logo=github)](https://github.com/ampl/colab.ampl.com/blob/master/authors/gomfy/dev/5_benders_parallel_stoch_floc.ipynb) [![Open In Colab](https://colab.research.google.com/assets/colab-badge.svg)](https://colab.research.google.com/github/ampl/colab.ampl.com/blob/master/authors/gomfy/dev/5_benders_parallel_stoch_floc.ipynb) [![Kaggle](https://kaggle.com/static/images/open-in-kaggle.svg)](https://kaggle.com/kernels/welcome?src=https://github.com/ampl/colab.ampl.com/blob/master/authors/gomfy/dev/5_benders_parallel_stoch_floc.ipynb) [![Gradient](https://assets.paperspace.io/img/gradient-badge.svg)](https://console.paperspace.com/github/ampl/colab.ampl.com/blob/master/authors/gomfy/dev/5_benders_parallel_stoch_floc.ipynb) [![Open In SageMaker Studio Lab](https://studiolab.sagemaker.aws/studiolab.svg)](https://studiolab.sagemaker.aws/import/github/ampl/colab.ampl.com/blob/master/authors/gomfy/dev/5_benders_parallel_stoch_floc.ipynb) [![Hits](https://h.ampl.com/https://github.com/ampl/colab.ampl.com/blob/master/authors/gomfy/dev/5_benders_parallel_stoch_floc.ipynb)](https://colab.ampl.com)\n",
        "\n",
        "Description: In the fifth installment of our six-part series, we delve deeper by showing how to evolve our Benders decomposition Python script from a serial execution to one that solves subproblems in parallel. \n",
        "This transition has the potential to significantly reduce solution times, especially when dealing with complex subproblems or a large number of scenarios that generate many subproblems.\n",
        "In our concluding notebook, we will introduce the AMPLS library, a robust resource designed to further enhance Benders decomposition's efficiency. \n",
        "It achieves this by enabling the export of AMPL model instances as persistent solver representations, thus streamlining the solution process.\n",
        "\n",
        "Tags: amplpy, ampl, mip, stochastic, facility location, benders, decomposition, parallel solves\n",
        "\n",
        "Notebook author: Gyorgy Matyasfalvi <<gyorgy@ampl.com>>\n",
        "\n",
        "References:\n",
        "\n",
        "* AMPL a Modeling Language for Mathematical Programming -- Robert Fourer et al.\n",
        "* Intro to Linear Optimization -- Dimitris Bertsimas et al.\n",
        "* SCIP Optimization Suite example projects -- Stephen J. Maher (ZIB)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 15,
      "metadata": {},
      "outputs": [],
      "source": [
        "# Install dependencies\n",
        "%pip install -q amplpy pandas"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 16,
      "metadata": {},
      "outputs": [],
      "source": [
        "# Google Colab & Kaggle integration\n",
        "from amplpy import AMPL, ampl_notebook\n",
        "\n",
        "ampl = ampl_notebook(\n",
        "    modules=[\"open\"],  # modules to install\n",
        "    license_uuid=\"default\",  # license to use\n",
        ")  # instantiate AMPL object and register magics"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 17,
      "metadata": {},
      "outputs": [],
      "source": [
        "# Import all necessary libraries\n",
        "import os\n",
        "import warnings\n",
        "import concurrent.futures\n",
        "\n",
        "# import queue\n",
        "import pandas as pd\n",
        "from typing import Tuple, Dict"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "## Introduction\n",
        "\n",
        "To ensure a cohesive understanding and for ease of reference, we will revisit the foundational concepts of stochastic programming and the theory behind Benders decomposition in this notebook. \n",
        "This repetition will help reinforce these critical concepts and provide a seamless transition into the PYTHON script implementing Benders decomposition of our stochastic facility location problem."
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "## Two-stage stochastic programming problem formulation\n",
        "\n",
        "Consider a decision-making process occurring over two sequential stages. \n",
        "Initially, a decision vector $x$ is selected. \n",
        "Following this, new information becomes available, prompting the selection of a second-stage decision vector $y$. \n",
        "We introduce $S$ possible scenarios for the newly acquired information, with the actual scenario being revealed only after the selection of $x$. \n",
        "Each scenario is indexed by $s$, and its occurrence probability is denoted by $\\alpha^s$, which we assume to be positive. Given the sequential nature of decision-making, the choice of $y$ can be scenario-dependent, denoted as $y^s$ to reflect this.\n",
        "\n",
        "The decision vectors $x$ and $y^s$ incur costs represented by vectors $f$ and $q$, respectively. \n",
        "The first-stage decisions are constrained by:\n",
        "\n",
        "\\begin{equation}\n",
        "\\begin{array}{crlllllll}\n",
        "\t& A x \t\t&  =  \t& {b} \\\\\n",
        "\t& x\t\t\t& \\geq\t& 0.\n",
        "\\end{array}\n",
        "\\tag{1}\n",
        "\\end{equation}\n",
        "\n",
        "Furthermore, the combination of first and second-stage decisions must satisfy the scenario-specific constraints:\n",
        "\n",
        "\\begin{equation}\n",
        "\\begin{array}{crcllllll}\n",
        "\t& {B}^s x + {D}^s y^s \t\t& =\t\t& {d}^s \\\\\n",
        "\t& y^s\t\t\t\t\t\t\t\t\t\t& \\geq \t& \t0,\n",
        "\\end{array}\n",
        "\\tag{2}\n",
        "\\end{equation}\n",
        "\n",
        "for every scenario $s$. \n",
        "The uncertainty in the model is encapsulated by $\\xi := (q, B, D, d)$, where a constraint is deemed part of the second stage if it involves any nonzero coefficients of variables ${y}$. \n",
        "The matrices ${B}^s$ are referred to as *technology matrices*, which determine the influence of the decision at the first stage on the second stage. \n",
        "The objective is to determine the choices of $x$ and $y^s$ for all $s=1, \\ldots, S$ that minimize the expected total cost:\n",
        "\n",
        "\\begin{equation}\n",
        "\\begin{array}{cllllllll}\n",
        "\t& f^T x & + & \\alpha_1 {(q^1)^T} y^1 & + & \\alpha^2 {(q^2)^T} y^2 & \\ldots & + & \\alpha^S {(q^S)^T} y^S \n",
        "\\end{array}\n",
        "\\tag{3}\n",
        "\\end{equation}\n",
        "\n",
        "This formulation leads to the following two-stage stochastic optimization problem:\n",
        "\n",
        "\n",
        "\\begin{equation}\n",
        "\\begin{array}{crcrcrccrcl}\n",
        "\\min\t\t\t& f^T x \t\t\t& + & \\alpha^1 {(q^1)^T} y_1 & + \t& \\alpha^2 {(q^1)^T} y^2 \t& \\ldots  \t& + & \\alpha^s {(q^s)^T} y^S\t&   & \\\\\n",
        "\\textrm{s.t.}\t& A x \t\t\t\t&   &\t\t \t\t\t   &\t& \t\t\t\t \t\t& \t \t\t& \t&\t\t\t\t\t\t&   &= {b} \\\\\n",
        "\t\t\t\t& {B^1} x \t\t\t& + & {D^1} y^1 \t \t   &    &\t\t \t\t\t \t&\t \t\t&\t&\t\t\t\t\t\t& \t&= {d^1} \\\\\n",
        "\t\t\t\t& {B^2} x \t\t\t& + &\t\t\t\t\t\t&   & {D^2} y^2 \t \t\t&\t\t\t&\t&\t\t\t\t\t\t& \t&= {d^2} \\\\\n",
        "\t\t\t\t& \\vdots \t\t\t&   &  \t\t\t  \t\t  \t&   & \t\t\t\t\t\t & \\ddots \t&\t&\t\t\t\t\t\t& \t&\\vdots \\\\\n",
        "\t\t\t\t& {B^S} x \t\t\t& + & \t  \t\t\t\t  \t&   & \t\t\t\t\t\t &\t \t  \t& \t& {D^S} y^S\t\t\t\t&\t&= {d^S} \\\\\n",
        "\t\t\t\t& x \\geq 0, \t\t& \t& y^1 \\geq 0,\t\t  & \\ldots\t&\t\t\t \t\t&\t\t\t& \t& y^S \\geq 0\t\t\t&   &\n",
        "\\end{array}\n",
        "\\tag{4}\n",
        "\\end{equation}\n",
        "\n",
        "we will refer to the above problem as the **master problem**. \n",
        "Notice that even if the number of possible scenarios $S$ is moderate, this formulation can be a large linear programming problem. \n",
        "However, a decomposition method can help."
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "## Benders decomposition for two-stage problems\n",
        "\n",
        "1. Given a vector $\\hat{x}$, using the dual simplex method solve the subproblem (for all $s \\in S$):\n",
        "\n",
        "    \\begin{equation}\n",
        "    \\begin{array}{crl}\n",
        "        \\min\t    &  (q^s)^T y^s\t& \\\\\n",
        "    \\textrm{s.t.}   &  D^s y^s \t\t& = d^s - B^s \\hat{x} \\\\\n",
        "    \t            &  y^s       \t& \\geq 0 \n",
        "    \\end{array}\n",
        "    \\tag{5}\n",
        "    \\end{equation}\n",
        "\n",
        "2. The dual of (5) takes the following form:\n",
        "\n",
        "    \\begin{equation}\n",
        "    \\begin{array}{crl}\n",
        "        \\max\t    & (p^s)^T (d^s - B^s \\hat{x}) & \\\\\n",
        "    \\textrm{s.t.}   & (p^s)^T D^s           \t\t& <=  (q^s)^T\n",
        "    \\end{array}\n",
        "    \\tag{6}\n",
        "    \\end{equation}\n",
        "\n",
        "    Let  \n",
        "    \n",
        "    \\begin{equation}\n",
        "    \\begin{array}{crrl}\n",
        "        \t    & P \\;\\; = & \\{ p & | \\;\\; p^T D \\leq q^T \\}  & \\\\\n",
        "    \\end{array}\n",
        "    \\tag{7}\n",
        "    \\end{equation}\n",
        "\n",
        "    We assume that $P$ is nonempty and has at least one extreme point (so we can apply duality theory). Let $p_i$, $i = 1, \\ldots, I$, be the extreme points, and let $w^j$, $j = 1, \\ldots, J$, be a complete set of extreme rays of $P$.\n",
        "\n",
        "3. Define $z^s(\\hat{x})^*$ as the optimal cost obtained from problem (5) described above. Returning to the optimization of $x$, we encounter the following problem:\n",
        "\n",
        "    \\begin{equation}\n",
        "    \\begin{array}{crcrcrccrcl}\n",
        "    \\min\t\t\t& f^T x \t    \t& + & \\sum_{s=1}^S \\alpha^s z^s(x) &    \\\\\n",
        "    \\textrm{s.t.}\t& A x \t\t\t\t&   &\t\t \t\t\t\t\t& = {b} \\\\\n",
        "    \t\t\t\t& x          \t    &\t\\geq 0 &              \t\t&\n",
        "    \\end{array}\n",
        "    \\tag{8}\n",
        "    \\end{equation}\n",
        "\n",
        "    Under our assumption that the set $P$ is nonempty, either the dual problem (6) has an optimal solution and $z^s(\\hat{x})$ is finite, or the optimal dual cost is infinite (indicating that the \n",
        "    primal (5) is infeasible).\n",
        "\n",
        "4. In particular, $z^s(x) \\leq \\infty \\leftrightarrow (w^j)^T (d^s - B^s x) \\leq 0, \\; \\forall j$ and whenever $z^s(x)$ is finite the optimum of (6) must be attained at an extreme point of the set $P$ meaning $z^s(x)$ is the smallest number $z^s$ such that $(p^i)^T (d^s - B^s x) \\leq z^s, \\; \\forall i$, leading us to the following master-problem formulation:\n",
        "\n",
        "\n",
        "    \\begin{equation}\n",
        "    \\begin{array}{crcrcrccrcl}\n",
        "    \\min\t\t\t& f^T x \t    \t    & + & \\sum_{s=1}^S \\alpha^s z^s &           & \\\\\n",
        "    \\textrm{s.t.}\t& A x \t\t\t\t    &   &\t\t \t\t\t\t\t& = b       & \\\\\n",
        "                \t& (p^i)^T (d^s - B^s x)\t&   &\t\t \t\t\t\t\t& \\leq z^s  & \\forall \\; i, s \\\\\n",
        "                \t& (w^j)^T (d^s - B^s x)\t&   &\t\t \t\t\t\t\t& \\leq 0    & \\forall \\; j, s \\\\\n",
        "    \t\t\t\t& x          \t        &\t\\geq 0 &              \t\t&           & \n",
        "    \\end{array}\n",
        "    \\tag{9}\n",
        "    \\end{equation}\n",
        "\n",
        "5. **In our algorithm we will only generate constraints that we find to be violated by the current solution.**\n"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "## The algorithm\n",
        "\n",
        "0. **Initialization step**: initialize $(\\hat{x}, \\hat{z}^1, \\ldots, \\hat{z}^S)$ to zero.\n",
        "\n",
        "1. **Sub-problem step**: using the dual-simplex method solve the sub-problems individually for every $s \\in S$:\n",
        "\n",
        "\\begin{equation}\n",
        "\\begin{array}{crl}\n",
        "    \\min\t    &  (q^s)^T y^s\t& \\\\\n",
        "\\textrm{s.t.}   &  D^s y^s \t\t& = d^s - B^s \\hat{x} \\\\\n",
        "\t            &  y^s       \t& \\geq 0 \n",
        "\\end{array}\n",
        "\\tag{5}\n",
        "\\end{equation}\n",
        "\n",
        "2. **Check optimality of master-problem**: If for every $s$, the corresponding sub-problem is feasible and **$(q^s)^T \\hat{y}^s \\leq \\hat{z}^s$** then all constraints are satisfied and we have an optimal solution to the master problem.\n",
        "\n",
        "3. **Add optimality cut**: If the sub-problem corresponding to some $s$ has an optimal solution such that **$(q^s)^T \\hat{y}^s > \\hat{z}^s$** the following cut is added to the master problem:  \n",
        "\\begin{equation}\n",
        "(p^s)^T (d^s - B^s x) \\leq z^s,\n",
        "\\tag{10}\n",
        "\\end{equation}\n",
        "where $p^s$ is an optimal basic feasible solution to the dual of the sub-problem. \n",
        "\n",
        "4. **Add feasibility cut**: If the sub-problem corresponding to some $s$ is infeasible, its dual has infinite cost and the following cut is added to the master problem:  \n",
        "\\begin{equation}\n",
        "(w^s)^T (d^s - B^s x) \\leq 0, \n",
        "\\tag{11}\n",
        "\\end{equation}\n",
        "where $w^s$ is a positive cost extreme ray. "
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "## Stochastic facility location a concrete example\n",
        "\n",
        "**Given:** \n",
        "* A set of facilities: $I$.\n",
        "* A set of customers: $J$.\n",
        "* Set of scenarios: $S$ (representing different customer demands).\n",
        "\n",
        "**Task:** \n",
        "* Find the minimum cost facilities to open such that the customer demand can be satisfied in all scenarios.\n",
        "\n",
        "### Variables\n",
        "* $x_i \\in \\{0, 1\\} \\quad \\forall i \\in I$\n",
        "    * $x_i = 1$ if facility $i$ is opened.\n",
        "* $y_{ij}^s \\geq 0 \\quad \\forall i \\in I, \\forall j \\in J, \\forall s \\in S$\n",
        "    * $y_{ij}^s$ is the level of demand for customer $j$ satisfied by facility $i$ in scenario $s$.\n",
        "\n",
        "### Parameters:\n",
        "* $\\alpha^s$: the probability of scenario $s$.\n",
        "* $f_i$: the fixed cost for opening facility $i$,\n",
        "* $q_{ij}$: the cost of servicing customer $j$ from facility $i$,\n",
        "* $\\lambda_j^s$: the demand of customer $j$ in scenario $s$,\n",
        "* $k_i:$ the capacity of facility $i$.\n",
        "\n",
        "\n",
        "### The extensive form\n",
        "The extensive form of our stochastic program can be formulated as follows:\n",
        "\n",
        "$\n",
        "\\begin{equation}\n",
        "\\begin{array}{rrlrrll}\n",
        " \\min \\quad & \\sum_{i \\in I} f_i x_i & + & \\sum_{s \\in S} \\sum_{i \\in I} \\sum_{j \\in J} \\alpha^s q_{ij} y_{ij}^s & & & \\\\\n",
        " & & & & & & \\\\\n",
        " \\textrm{s.t.} \\quad & & & \\sum_{i \\in I} y_{ij}^s & \\geq & \\lambda_j^s & \\forall j \\in J, \\forall s \\in S \\\\\n",
        " & & & \\sum_{j \\in J} y_{ij}^s & \\leq & k_i x_i & \\forall i \\in I, \\forall s \\in S \\\\\n",
        " & \\sum_{i \\in I} k_i x_i & & & \\geq & \\max_{s \\in S} \\sum_{j \\in J} \\lambda_j^s & \\\\\n",
        " & & & & & & \\\\\n",
        " & x_i  & \\in \\{0, 1\\} & & & & \\forall i \\in I \\\\\n",
        " & y_{ij}^s & \\geq 0 & & & & \\forall i \\in I, \\forall j \\in J, \\forall s \\in S  \n",
        "\\end{array} \\tag{12}\n",
        "\\end{equation}\n",
        "$\n",
        "\n",
        "### The sub-problem\n",
        "Given the above formulation for the extensive form (12) we can express our scenario specific sub-problems as follows\n",
        "\n",
        "$\n",
        "\\begin{equation}\n",
        "\\begin{array}{lcrllll}\n",
        " z_s(\\hat{x}) & = & \\min \\quad &  \\sum_{i \\in I} \\sum_{j \\in J} q_{ij} y_{ij}^s & & &  \\\\\n",
        " & & & & & & \\\\\n",
        " & & \\textrm{s.t.} \\quad & \\sum_{i \\in I} y_{ij}^s & \\geq  & \\lambda_j^s & \\forall j \\in J \\\\\n",
        " & & & \\sum_{j \\in J} y_{ij}^s & \\leq & k_i \\hat{x_i} & \\forall i \\in I \\\\\n",
        " & & & & & &  \\\\\n",
        " & & & y_{ij}^s \\geq 0 & & & \\forall \\; i \\in I, \\forall j \\in J\n",
        "\\end{array} \\tag{13}\n",
        "\\end{equation}\n",
        "$\n",
        "\n",
        "### The master-problem\n",
        "The master-problem takes the following form\n",
        "\n",
        "$\n",
        "\\begin{equation}\n",
        "\\begin{array}{rllll}\n",
        " \\min \\quad          & \\sum_{i \\in I} f_i x_i + \\sum_{s \\in S} \\alpha^s z^s & & & \\\\\n",
        " & & & & \\\\\n",
        " \\textrm{s.t.} \\quad & \\sum_{i \\in I} k_i x_i & \\geq & \\max_{s \\in S} \\sum_{j \\in J} \\lambda_j^s & \\\\\n",
        "                     & \\sum_{j \\in J} p_j^s \\lambda_j^s + \\sum_{i \\in I} p_i^s k_i x_i & \\leq & z^s & \\forall s \\in S, \\forall p \\in P^s \\\\\n",
        "                     & \\sum_{j \\in J} w_j^s \\lambda_j^s  + \\sum_{i \\in I} w_i^s k_i x_i &  \\leq & 0 & \\forall s \\in S, \\forall w \\in P^s \\\\\n",
        " & & & & \\\\\n",
        "                     & x_i \\in \\{0, 1\\} & & & \\forall i \\in I \\\\\n",
        "\\end{array} \\tag{12}\n",
        "\\end{equation}\n",
        "$"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "## AMPL implementation\n",
        "\n",
        "Converting the mathematical formulation of our decomposition scheme into an AMPL model is straightforward. \n",
        "Similar to the extensive form, AMPL enables the modeling of problems in a way that closely aligns with their conceptual structure. \n",
        "In AMPL, sets and parameters are globally accessible, allowing for seamless integration across different problem components.\n",
        "\n",
        "AMPL facilitates the assignment of variables, objective functions, and constraints to specific named problems using the syntax `problem <problemName>: <list of model entities>`. \n",
        "In the model file presented below, we've organized the content to clearly indicate which elements pertain to the master problem and which to the subproblems. \n",
        "This organization is primarily for clarity and ease of understanding; the AMPL interpreter treats the entire file as a single problem until it encounters explicit problem statements.\n",
        "\n",
        "Beyond the entities familiar from the extensive form, we introduce additional parameters like `sub_facility_open`, which convey the master problem's solutions to the subproblems. \n",
        "We also incorporate parameters for the dual values from the subproblems, namely `customer_price` and `facility_price`. \n",
        "New constraints, often referred to as cuts, are added to the master problem based on the subproblems' feasibility and optimality relative to the master problem's objectives. \n",
        "\n",
        "Furthermore, a new set of variables, `sub_variable_cost`, is added to the master problem. These variables account for the variable costs identified in each scenario's subproblem, contributing to the overall total cost calculation."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 18,
      "metadata": {},
      "outputs": [
        {
          "name": "stdout",
          "output_type": "stream",
          "text": [
            "Overwriting floc_bend.mod\n"
          ]
        }
      ],
      "source": [
        "%%writefile floc_bend.mod\n",
        "\n",
        "# Global sets\n",
        "set FACILITIES; # set of facilities\n",
        "set CUSTOMERS;  # set of customers\n",
        "set SCENARIOS;  # set of scenarios\n",
        "\n",
        "# Global parameters\n",
        "param customer_demand{CUSTOMERS, SCENARIOS} >= 0;   # customer_demand of customer j in scenario s\n",
        "param facility_capacity{FACILITIES} >= 0;           # facility_capacity of facility_open i\n",
        "\n",
        "\n",
        "### SUBPROBLEM ###\n",
        "\n",
        "# Subproblem parameters\n",
        "param variable_cost{FACILITIES, CUSTOMERS} >= 0;    # variable cost of satisfying customer_demand of customer j from facility_open i\n",
        "\n",
        "# Bender's parameters\n",
        "param sub_facility_open{FACILITIES} binary, default 1; # 1 if facility i is open, 0 otherwise\n",
        "param sub_scenario symbolic in SCENARIOS;\n",
        "\n",
        "# Subproblem variables\n",
        "var production{FACILITIES, CUSTOMERS, SCENARIOS} >= 0;  # production from facility i to satisfy customer demand j in scenario s\n",
        "\n",
        "# Subproblem objective\n",
        "minimize operating_cost: \n",
        "    sum{i in FACILITIES, j in CUSTOMERS} variable_cost[i,j]*production[i,j,sub_scenario]; \n",
        "\n",
        "# Subproblem constraints\n",
        "s.t. satisfying_customer_demand{j in CUSTOMERS}:\n",
        "    sum{i in FACILITIES} production[i,j,sub_scenario] >= customer_demand[j,sub_scenario];\n",
        "\n",
        "s.t. facility_capacity_limits{i in FACILITIES}:\n",
        "    sum{j in CUSTOMERS} production[i,j,sub_scenario] <= facility_capacity[i] * sub_facility_open[i];\n",
        "\n",
        "\n",
        "\n",
        "### MASTER PROBLEM ### \n",
        "\n",
        "# Master parameters\n",
        "param fixed_cost{FACILITIES} >= 0;               # fixed cost of opening facility_open i\n",
        "param prob{SCENARIOS} default 1/card(SCENARIOS); # probability of scenario s\n",
        "\n",
        "# Bender's parameters\n",
        "param nITER >= 0 integer;\n",
        "param cut_type {1..nITER, SCENARIOS} symbolic within {\"opt\", \"feas\", \"none\"};\n",
        "param customer_price{CUSTOMERS, SCENARIOS, 1..nITER} >= 0;  # dual variable for the demand constraint\n",
        "param facility_price{FACILITIES, SCENARIOS, 1..nITER} <= 0; # dual variable for the capacity constraint\n",
        "\n",
        "# Master variables\n",
        "var sub_variable_cost{SCENARIOS} >= 0; # subproblem objective bound\n",
        "var facility_open{FACILITIES} binary;  # 1 if facility i is open, 0 otherwise\n",
        "\n",
        "# Master objective\n",
        "minimize total_cost: \n",
        "    sum{i in FACILITIES} fixed_cost[i]*facility_open[i] + sum{s in SCENARIOS} prob[s]*sub_variable_cost[s];\n",
        "\n",
        "# Master constraints\n",
        "s.t. sufficient_production_capacity:\n",
        "    sum{i in FACILITIES} facility_capacity[i]*facility_open[i] >= max{s in SCENARIOS} sum{j in CUSTOMERS} customer_demand[j,s];\n",
        "\n",
        "s.t. optimality_cut {k in 1..nITER, s in SCENARIOS: cut_type[k,s] == \"opt\"}:\n",
        "        sum{j in CUSTOMERS} customer_price[j,s,k]*customer_demand[j,s] + \n",
        "        sum{i in FACILITIES} facility_price[i,s,k]*facility_capacity[i]*facility_open[i] <= sub_variable_cost[s]; \n",
        "\n",
        "s.t. feasibility_cut {k in 1..nITER, s in SCENARIOS: cut_type[k,s] == \"feas\"}: \n",
        "        sum{j in CUSTOMERS} customer_price[j,s,k]*customer_demand[j,s] + \n",
        "        sum{i in FACILITIES} facility_price[i,s,k]*facility_capacity[i]*facility_open[i] <= 0; "
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "### Data as Python objects\n",
        "\n",
        "In this notebook, we continue to work with the same dataset as in the previous installments. \n",
        "However, the format in which we store our data differs; we now utilize Python lists and Pandas data frames. Initially, we create these Python objects and subsequently employ `amplpy` methods to seamlessly transfer the data into AMPL."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 19,
      "metadata": {},
      "outputs": [],
      "source": [
        "# Define the data\n",
        "\n",
        "# Define sets\n",
        "FACILITIES = [\"Baytown_TX\", \"Beaumont_TX\", \"Baton_Rouge_LA\"]\n",
        "CUSTOMERS = [\"San_Antonio_TX\", \"Dallas_TX\", \"Jackson_MS\", \"Birmingham_AL\"]\n",
        "SCENARIOS = [\"Low\", \"Medium\", \"High\"]\n",
        "\n",
        "# Define parameters indexed over FACILITIES\n",
        "fixed_cost_df = pd.DataFrame(\n",
        "    [400000, 200000, 600000], index=FACILITIES, columns=[\"fixed_cost\"]\n",
        ")\n",
        "\n",
        "facility_capacity_df = pd.DataFrame(\n",
        "    [1550, 650, 1750], index=FACILITIES, columns=[\"facility_capacity\"]\n",
        ")\n",
        "\n",
        "# Define parameters indexed over SCENARIOS\n",
        "prob_df = pd.DataFrame([0.25, 0.5, 0.25], index=SCENARIOS, columns=[\"prob\"])\n",
        "\n",
        "# Define parameters indexed over CUSTOMERS and SCENARIOS\n",
        "customer_demand_df = pd.DataFrame(\n",
        "    [450, 650, 887, 910, 1134, 1456, 379, 416, 673, 91, 113, 207],\n",
        "    index=pd.MultiIndex.from_product(\n",
        "        [CUSTOMERS, SCENARIOS], names=[\"CUSTOMERS\", \"SCENARIOS\"]\n",
        "    ),\n",
        "    columns=[\"customer_demand\"],\n",
        ")\n",
        "\n",
        "# Define parameters indexed over FACILITIES and CUSTOMERS\n",
        "variable_cost_df = pd.DataFrame(\n",
        "    [\n",
        "        5739.725,\n",
        "        6539.725,\n",
        "        8650.40,\n",
        "        22372.1125,\n",
        "        6055.05,\n",
        "        6739.055,\n",
        "        8050.40,\n",
        "        21014.225,\n",
        "        8650.40,\n",
        "        7539.055,\n",
        "        4539.72,\n",
        "        15024.325,\n",
        "    ],\n",
        "    index=pd.MultiIndex.from_product(\n",
        "        [FACILITIES, CUSTOMERS], names=[\"FACILITIES\", \"CUSTOMERS\"]\n",
        "    ),\n",
        "    columns=[\"variable_cost\"],\n",
        ")"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "## Solving subproblems in parallel\n",
        "\n",
        "Leveraging multiple processor cores enables us to solve subproblems in parallel, enhancing efficiency by allowing each subproblem to be addressed independently. \n",
        "To facilitate this, we'll develop a dedicated worker function specifically tailored for solving these subproblems. \n",
        "Each subproblem will be processed in its own AMPL process, with the solutions integrated back into the master problem to update its parameters.\n",
        "\n",
        "### **Noteworthy Changes**\n",
        "\n",
        "In contrast to our approach in previous notebooks, where a single AMPL process managed both the master problem and all subproblems, our new setup diverges significantly. \n",
        "**We will initiate separate AMPL processes for each subproblem, in addition to one for the master problem, resulting in a total of \"number of scenarios plus one\" AMPL processes.**\n",
        "This change precludes the possibility of direct data transfers between subproblems and the master problem using AMPL's `let` commands, such as:\n",
        "\n",
        "```ampl\n",
        "let {j in CUSTOMERS} customer_price[j, sub_scenario, nITER] := satisfying_customer_demand[j].dual;\n",
        "let {i in FACILITIES} facility_price[i, sub_scenario, nITER] := facility_capacity_limits[i].dual;\n",
        "```\n",
        "\n",
        "Instead, data transfer between the subproblems and the master problem will be managed through Python. \n",
        "This requires the master problem's Python process to extract data from the subproblems and then load it into the master problem's AMPL process, and vice versa.\n",
        "\n",
        "To manage these data transformations efficiently, additional Python code is necessary. \n",
        "A practical solution is to develop a Python module comprised of utility and worker functions tailored for these tasks. \n",
        "Given our Jupyter Notebook environment, we will define these functions in a separate code cell for ease of access and organization."
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "### Utility and worker Functions\n",
        "\n",
        "To streamline our approach before tackling the core of our task, we will establish a series of utility functions designed to handle data transformations and to reduce redundancy in our code. \n",
        "Alongside these, we'll introduce a pivotal worker function, `solve_sub()`, designated to be run by individual threads in our parallel setup.\n",
        "\n",
        "Worker functions typically have access to global symbols and data. \n",
        "However, in the realm of parallel computing, it's advisable to directly supply data to the worker function. \n",
        "This practice promotes clarity and ensures that each operation is distinct and self-contained. \n",
        "Adhering to this principle, the `solve_sub()` function will use a specific AMPL process for each scenario, handling the necessary computations for each subproblem independently."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 20,
      "metadata": {},
      "outputs": [],
      "source": [
        "# Define data loading functions\n",
        "def load_set_data(list_dict: dict, ampl: AMPL) -> None:\n",
        "    for key, value in list_dict.items():\n",
        "        ampl.set[key] = value\n",
        "\n",
        "\n",
        "def load_param_data(df_dict: dict, ampl: AMPL) -> None:\n",
        "    # Since we have named the columns of our data frame with the same names as the AMPL parameters, we don't need the keys\n",
        "    for _, value in df_dict.items():\n",
        "        ampl.set_data(value)\n",
        "\n",
        "\n",
        "# Define functions for setting options\n",
        "def set_benders_options(ampl: AMPL) -> None:\n",
        "    # Set options\n",
        "    ampl.option[\"solver\"] = \"highs\"  # Select the solver\n",
        "    # Ensure that the solver uses the dual simplex method and returns dual rays that's why presolve is turned off:\n",
        "    ampl.option[\"highs_options\"] = \"pre:solve off alg:simplex 1\"\n",
        "    # Also turn off presolve in AMPL so we get info from the solver:\n",
        "    ampl.option[\"presolve\"] = 0\n",
        "    ampl.option[\"display_eps\"] = 0.000001  # Set the display precision\n",
        "    ampl.option[\"omit_zero_rows\"] = 1  # Omit zero rows in the output\n",
        "\n",
        "\n",
        "# Define the subproblem initialization function\n",
        "def init_sub(\n",
        "    model_file_path: str, list_dict: dict, df_dict: dict, s: str\n",
        ") -> Tuple[str, AMPL]:\n",
        "    ampl = AMPL()  # Instantiate an AMPL object\n",
        "    ampl.read(model_file_path)  # Read the model from file\n",
        "    # Load sets\n",
        "    load_set_data(list_dict, ampl)\n",
        "    # Load parameters\n",
        "    load_param_data(df_dict, ampl)\n",
        "    # Set options\n",
        "    set_benders_options(ampl)\n",
        "    ampl.param[\"sub_scenario\"] = s\n",
        "    # Define subproblem entities\n",
        "    ampl.eval(\n",
        "        r\"\"\"\n",
        "        problem Sub: production, operating_cost, satisfying_customer_demand, facility_capacity_limits;\n",
        "        \"\"\"\n",
        "    )\n",
        "    return s, ampl\n",
        "\n",
        "\n",
        "# Define the subproblem worker function\n",
        "# def solve_sub(s: str, ampl: AMPL) -> None:\n",
        "def solve_sub(s: str, ampl: AMPL) -> None:\n",
        "    # Solve the subproblem\n",
        "    ampl.param[\"sub_scenario\"] = s\n",
        "    ampl.solve(problem=\"Sub\")\n",
        "    # After solving put the scenario and the AMPL object into the queue and make it available to the master problem\n",
        "    # results_queue.put((s, ampl))\n",
        "\n",
        "\n",
        "# Define the master initialization function\n",
        "def init_master(model_file_path: str, list_dict: dict, df_dict: dict) -> AMPL:\n",
        "    ampl = AMPL()  # Instantiate an AMPL object\n",
        "    ampl.read(model_file_path)  # Read the model from file\n",
        "    # Define suffix for storing dual rays\n",
        "    ampl.eval(r\"suffix dunbdd;\")\n",
        "    # Load sets\n",
        "    load_set_data(list_dict, ampl)\n",
        "    # Load parameters\n",
        "    load_param_data(df_dict, ampl)\n",
        "    # Set options\n",
        "    set_benders_options(ampl)\n",
        "    # Define master problem entities\n",
        "    ampl.eval(\n",
        "        r\"\"\"\n",
        "        problem Master: facility_open, sub_variable_cost, total_cost, optimality_cut, feasibility_cut, sufficient_production_capacity; \n",
        "        problem Sub: production, operating_cost, satisfying_customer_demand, facility_capacity_limits;\n",
        "        \"\"\"\n",
        "    )\n",
        "    return ampl\n",
        "\n",
        "\n",
        "# Define add feasibility cut function for the master problem\n",
        "def add_feas_cut(master: AMPL, sub: AMPL, n_iter: int, s: str) -> None:\n",
        "    # Extract dual rays for constraints related to customers\n",
        "    # from subproblem and store them in such a way that we can readily load them into master\n",
        "    c_rays_df = (\n",
        "        sub.con[\"satisfying_customer_demand\"]\n",
        "        .get_values(suffixes=[\"dunbdd\"])\n",
        "        .to_pandas()\n",
        "    )\n",
        "    c_rays_df.rename(\n",
        "        columns={\"satisfying_customer_demand.dunbdd\": \"customer_price\"}, inplace=True\n",
        "    )\n",
        "    c_rays_df.index = pd.MultiIndex.from_tuples(\n",
        "        [(i, s, n_iter) for i in c_rays_df.index]\n",
        "    )\n",
        "    # Extract dual rays for constraints related to facilities\n",
        "    f_rays_df = (\n",
        "        sub.con[\"facility_capacity_limits\"].get_values(suffixes=[\"dunbdd\"]).to_pandas()\n",
        "    )\n",
        "    f_rays_df.rename(\n",
        "        columns={\"facility_capacity_limits.dunbdd\": \"facility_price\"}, inplace=True\n",
        "    )\n",
        "    f_rays_df.index = pd.MultiIndex.from_tuples(\n",
        "        [(i, s, n_iter) for i in f_rays_df.index]\n",
        "    )\n",
        "    # Update cut_type, customer_price, and facility_price inside the master\n",
        "    master.param[\"cut_type\"][(n_iter, s)] = \"feas\"\n",
        "    master.set_data(c_rays_df)\n",
        "    master.set_data(f_rays_df)\n",
        "\n",
        "\n",
        "# Define add optimality cut function for the master problem\n",
        "def add_opt_cut(master: AMPL, sub: AMPL, n_iter: int, s: str) -> None:\n",
        "    # Extract dual rays for constraints related to customers\n",
        "    # from subproblem and store them in such a way that we can readily load them into master\n",
        "    c_duals_df = (\n",
        "        sub.con[\"satisfying_customer_demand\"].get_values(suffixes=[\"dual\"]).to_pandas()\n",
        "    )\n",
        "    c_duals_df.rename(\n",
        "        columns={\"satisfying_customer_demand.dual\": \"customer_price\"}, inplace=True\n",
        "    )\n",
        "    c_duals_df.index = pd.MultiIndex.from_tuples(\n",
        "        [(i, s, n_iter) for i in c_duals_df.index]\n",
        "    )\n",
        "    # Extract dual rays for constraints related to facilities\n",
        "    f_duals_df = (\n",
        "        sub.con[\"facility_capacity_limits\"].get_values(suffixes=[\"dual\"]).to_pandas()\n",
        "    )\n",
        "    f_duals_df.rename(\n",
        "        columns={\"facility_capacity_limits.dual\": \"facility_price\"}, inplace=True\n",
        "    )\n",
        "    f_duals_df.index = pd.MultiIndex.from_tuples(\n",
        "        [(i, s, n_iter) for i in f_duals_df.index]\n",
        "    )\n",
        "    # Update cut_type, customer_price, and facility_price inside the master\n",
        "    master.param[\"cut_type\"][(n_iter, s)] = \"opt\"\n",
        "    master.set_data(c_duals_df)\n",
        "    master.set_data(f_duals_df)\n",
        "\n",
        "\n",
        "# Define update subproblem function\n",
        "def update_sub(master: AMPL, sub_dict: Dict[str, AMPL]) -> None:\n",
        "    # Pass master solution to subproblem as parameters\n",
        "    sub_facility_open_df = master.var[\"facility_open\"].get_values().to_pandas()\n",
        "    sub_facility_open_df.rename(\n",
        "        columns={\"facility_open.val\": \"sub_facility_open\"}, inplace=True\n",
        "    )\n",
        "    for _, ampl in sub_dict.items():\n",
        "        ampl.set_data(sub_facility_open_df)\n",
        "\n",
        "\n",
        "# Close subproblem AMPL objects\n",
        "def close_sub(sub_dict: Dict[str, AMPL]) -> None:\n",
        "    for _, ampl in sub_dict.items():\n",
        "        ampl.close()"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "### **Parallel Benders**\n",
        "Key elements of our parallel Benders implementation include:\n",
        "- **Data Management**: We employ dictionaries for each data type (Python lists, Pandas data frames) to streamline and clarify data loading operations, reducing code redundancy.\n",
        "- **Core Availability Check**: We assess the number of available processor cores. If the number of scenarios exceeds the available cores, we issue a warning to the user.\n",
        "- **AMPL Object Initialization**: For each scenario, we initialize an `ampl_dict` object, each containing an AMPL object configured to run as an independent process.\n",
        "- **Subproblem Executor**: We establish `sub_executor`, a dedicated object for managing the execution of subproblems.\n",
        "- **While Loop Execution**: The iterative process begins similarly to previous implementations, but with the solve operations are conducted outside the for loop in parallel.\n",
        "- **Violation Checks and Cuts**: Within the for loop, we simply identify violations and append the necessary cuts to the master problem.\n",
        "- **Cleanup Operations**: Upon exiting the while loop, we undertake cleanup operations to prevent the persistence of unused processes, ensuring system resources are efficiently managed.\n",
        "\n",
        "This structured approach allows us to leverage parallel computing capabilities, potentially enhancing the efficiency and speed of the Benders decomposition process."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 21,
      "metadata": {},
      "outputs": [
        {
          "name": "stdout",
          "output_type": "stream",
          "text": [
            "option solver highs;\n",
            "\n",
            "\n",
            "ITERATION 1\n",
            "\n",
            "HiGHS 1.6.0:   pre:solve = off\n",
            "  alg:simplex = 1\n",
            "HiGHS 1.6.0: HiGHS 1.6.0:  pre:solve = off\n",
            "  pre:solve = off\n",
            "  alg:simplex = 1\n",
            "HiGHS 1.6.0: optimal solution; objective 11621793.46\n",
            "4 simplex iterations\n",
            "0 barrier iterations\n",
            "  alg:simplex = 1\n",
            "HiGHS 1.6.0: optimal solution; objective 21050711.2\n",
            "6 simplex iterations\n",
            "0 barrier iterations\n",
            "HiGHS 1.6.0: optimal solution; objective 14779784.86\n",
            "5 simplex iterations\n",
            "0 barrier iterations\n",
            "1: Optimality cut added for scenario Low\n",
            "1: Optimality cut added for scenario Medium\n",
            "1: Optimality cut added for scenario High\n",
            "\n",
            "SOLVING THE MASTER PROBLEM\n",
            "\n",
            "HiGHS 1.6.0:   pre:solve = off\n",
            "  alg:simplex = 1\n",
            "\bHiGHS 1.6.0: optimal solution; objective 16688018.6\n",
            "0 simplex iterations\n",
            "1 branching nodes\n",
            "\n",
            "\n",
            "ITERATION 2\n",
            "\n",
            "HiGHS 1.6.0: HiGHS 1.6.0: HiGHS 1.6.0:   pre:solve = off\n",
            "  alg:simplex = 1\n",
            "  pre:solve = off\n",
            "  pre:solve = off\n",
            "  alg:simplex = 1\n",
            "\b  alg:simplex = 1\n",
            "HiGHS 1.6.0: optimal solution; objective 11621793.46\n",
            "0 simplex iterations\n",
            "0 barrier iterations\n",
            "\bHiGHS 1.6.0: optimal solution; objective 14966984.87\n",
            "1 simplex iterations\n",
            "0 barrier iterations\n",
            "\bHiGHS 1.6.0: optimal solution; objective 21570711.2\n",
            "0 simplex iterations\n",
            "0 barrier iterations\n",
            "2: No cut needed for scenario Low\n",
            "2: Optimality cut added for scenario Medium\n",
            "2: No cut needed for scenario High\n",
            "\n",
            "SOLVING THE MASTER PROBLEM\n",
            "\n",
            "HiGHS 1.6.0:   pre:solve = off\n",
            "  alg:simplex = 1\n",
            "\bHiGHS 1.6.0: optimal solution; objective 16758018.6\n",
            "1 simplex iterations\n",
            "1 branching nodes\n",
            "\n",
            "\n",
            "ITERATION 3\n",
            "\n",
            "HiGHS 1.6.0:   pre:solve = off\n",
            "HiGHS 1.6.0:   alg:simplex = 1\n",
            "HiGHS 1.6.0:  pre:solve = off\n",
            "  pre:solve = off\n",
            "  alg:simplex = 1\n",
            "  alg:simplex = 1\n",
            "\bHiGHS 1.6.0: optimal solution; objective 11621793.46\n",
            "0 simplex iterations\n",
            "0 barrier iterations\n",
            "HiGHS 1.6.0: optimal solution; objective 14779784.87\n",
            "1 simplex iterations\n",
            "0 barrier iterations\n",
            "\bHiGHS 1.6.0: optimal solution; objective 21050711.2\n",
            "0 simplex iterations\n",
            "0 barrier iterations\n",
            "3: No cut needed for scenario Low\n",
            "3: No cut needed for scenario Medium\n",
            "3: No cut needed for scenario High\n",
            "\n",
            "\n",
            "OPTIMAL SOLUTION FOUND\n",
            "Expected cost = 16758018.596250001\n",
            "\n",
            "facility_open [*] :=\n",
            "Baton_Rouge_LA  1\n",
            "    Baytown_TX  1\n",
            "   Beaumont_TX  1\n",
            ";\n",
            "\n"
          ]
        }
      ],
      "source": [
        "# Define data dictionaries for easier loading\n",
        "list_dict = {\"FACILITIES\": FACILITIES, \"CUSTOMERS\": CUSTOMERS, \"SCENARIOS\": SCENARIOS}\n",
        "df_dict = {\n",
        "    \"fixed_cost\": fixed_cost_df,\n",
        "    \"variable_cost\": variable_cost_df,\n",
        "    \"customer_demand\": customer_demand_df,\n",
        "    \"facility_capacity\": facility_capacity_df,\n",
        "    \"prob\": prob_df,\n",
        "}\n",
        "\n",
        "# Check available threads if there aren't enough then warn the user\n",
        "num_cores = os.cpu_count()\n",
        "num_cores_needed = len(SCENARIOS)\n",
        "if num_cores < num_cores_needed:\n",
        "    warnings.warn(\n",
        "        f\"Number of threads ({num_cores}) is less than the number of scenarios ({num_cores_needed}).\"\n",
        "    )\n",
        "\n",
        "# Instantiate an AMPL object for the master problem, load the model, and set options\n",
        "file_path = \"floc_bend.mod\"  # Define the file path\n",
        "master = init_master(file_path, list_dict, df_dict)\n",
        "master.eval(\"option solver;\")\n",
        "\n",
        "# Set up AMPL and Python for Benders decomposition\n",
        "# Loop logic variables\n",
        "epsilon = 0.00001\n",
        "n_iter = 0\n",
        "master.param[\"nITER\"] = n_iter\n",
        "violation_ctr = 0\n",
        "\n",
        "# Create a queue to hold the results of the subproblem\n",
        "# results_queue = queue.Queue()\n",
        "# Create AMPL objects for each scenario\n",
        "sub_ampl_dict = dict(init_sub(file_path, list_dict, df_dict, s) for s in SCENARIOS)\n",
        "# Create thread handlers for each scenario\n",
        "sub_executor = concurrent.futures.ThreadPoolExecutor(max_workers=len(SCENARIOS))\n",
        "\n",
        "while True:\n",
        "    violation_ctr = 0\n",
        "    n_iter += 1\n",
        "    master.param[\"nITER\"] = n_iter\n",
        "    print(f\"\\n\\nITERATION {n_iter}\\n\")\n",
        "\n",
        "    # Solve each subproblem in parallel\n",
        "    ### BEGIN PARALLEL COMPUTATION ###\n",
        "    # futures = {sub_executor.submit(solve_sub, s, sub_ampl_dict[s], results_queue): s for s in SCENARIOS}\n",
        "    futures = {\n",
        "        sub_executor.submit(solve_sub, s, sub_ampl_dict[s]): s for s in SCENARIOS\n",
        "    }\n",
        "    concurrent.futures.wait(futures)\n",
        "    ### END PARALLEL COMPUTATION ###\n",
        "\n",
        "    # Check results of each sub problem\n",
        "    for s, sub in sub_ampl_dict.items():\n",
        "        # Check for feasibility\n",
        "        if sub.solve_result == \"infeasible\":\n",
        "            print(f\"{n_iter}: Feasibility cut added for scenario {s}\")\n",
        "            add_feas_cut(master, sub, n_iter, s)\n",
        "        # Check for optimality\n",
        "        elif (\n",
        "            sub.obj[\"operating_cost\"].value()\n",
        "            > master.var[\"sub_variable_cost\"][s].value() + epsilon\n",
        "        ):\n",
        "            print(f\"{n_iter}: Optimality cut added for scenario {s}\")\n",
        "            add_opt_cut(master, sub, n_iter, s)\n",
        "        else:\n",
        "            print(f\"{n_iter}: No cut needed for scenario {s}\")\n",
        "            violation_ctr += 1\n",
        "            master.param[\"cut_type\"][(n_iter, s)] = \"none\"\n",
        "\n",
        "    # Exit loop if no cuts are needed indicating optimality\n",
        "    if violation_ctr == len(SCENARIOS):\n",
        "        break\n",
        "\n",
        "    # Solve the master problem\n",
        "    print(f\"\\nSOLVING THE MASTER PROBLEM\\n\")\n",
        "    master.solve(problem=\"Master\")\n",
        "    update_sub(master, sub_ampl_dict)\n",
        "\n",
        "# Clean up threads and AMPL objects\n",
        "sub_executor.shutdown()\n",
        "close_sub(sub_ampl_dict)\n",
        "\n",
        "print(\n",
        "    f\"\\n\\nOPTIMAL SOLUTION FOUND\\nExpected cost = {master.obj['total_cost'].value()}\\n\"\n",
        ")\n",
        "master.eval(\"display facility_open;\")"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "## Conclusion\n",
        "\n",
        "In this section, we've illustrated the process of integrating parallel computations into the Benders decomposition algorithm, with only minimal modifications required from our original serial Python implementation. \n",
        "As the complexity of our code increases, we've adopted certain practices to enhance the modularity and robustness of our application.\n",
        "\n",
        "Typically, this stage would involve the development of Python modules and the implementation of comprehensive exception handling to prepare the application for deployment. \n",
        "However, a detailed exploration of these advanced topics falls beyond the scope of this series.\n",
        "\n",
        "While our current implementation efficiently addresses subproblems, each solver instance still begins its computation from scratch. \n",
        "In the concluding notebook of this series, we will explore the use of the AMPLS library to maintain persistent solver instances, thereby further improving the efficiency of our solution process."
      ]
    }
  ],
  "metadata": {
    "kernelspec": {
      "display_name": "Python 3 (ipykernel)",
      "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.11.7"
    },
    "vscode": {
      "interpreter": {
        "hash": "ee4198aac13dae9b6fccd443041680d99c2643c4956ef7c80ff6dcc0057cb523"
      }
    }
  },
  "nbformat": 4,
  "nbformat_minor": 2
}