01_Google_OR_Tools_demo_Build_efficient_powerful_army.ipynb

{
  "cells": [
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "view-in-github",
        "colab_type": "text"
      },
      "source": [
        "<a href=\"https://colab.research.google.com/github/gulabpatel/LinearProgramming/blob/main/01_Google_OR_Tools_demo_Build_efficient_powerful_army.ipynb\" target=\"_parent\"><img src=\"https://colab.research.google.com/assets/colab-badge.svg\" alt=\"Open In Colab\"/></a>"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "lzfCiPH4xDk2"
      },
      "source": [
        "# Introduction to Linear Programming in Python \n",
        "> Mathematical optimization with Google OR-Tools\n",
        "\n",
        "- toc: false \n",
        "- badges: true\n",
        "- comments: true\n",
        "- author: Maxime Labonne\n",
        "- categories: [linear optimization, linear programming, OR-Tools, ortools, python, tutorial]\n",
        "- image: images/ortools/thumbnail.png\n",
        "- permalink: /linearoptimization/"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "2ywEZ_RkknQC"
      },
      "source": [
        "<center><div id=\"anim_container\"><img alt=\"Illustration of linear optimization\" src=\"https://mlabonne.github.io/blog/images/ortools/units.png\" loading=\"lazy\" id=\"anim_inner\"></div></center>\n",
        "\n",
        "In real warfare and strategy games, strategists are faced with a common problem: **how to allocate their resources to be the most efficient possible**? Let's take an example where we have three resources: 🌾**food**, 🪵**wood**, and 🪙**gold**. We can produce three types of units: 🗡️**swordsmen**, 🏹**bowmen**, and 🐎**horsemen**. 🐎Horsemen are **better than 🏹bowmen**, who are in turn **better than 🗡️swordsmen**. The following table provides the **cost and power of each unit** (they're freely inspired by a real strategy game, Age of Empires IV):\n",
        "<br/>\n",
        "<br/>\n",
        "\n",
        "| Unit | 🌾Food | 🪵Wood | 🪙Gold | 💪Power |\n",
        "| :--- | :---: | :---: | :---: | :---: |\n",
        "| 🗡️Swordsman | 60 | 20 | 0 | 70 |\n",
        "| 🏹Bowman | 80 | 10 | 40 | 95 |\n",
        "| 🐎Horseman | 140 | 0 |100 | 230 |\n",
        "\n",
        "<br/>\n",
        "\n",
        "Imagine that we have 1200 🌾**food**, 800 🪵**wood**, and 600 🪙**gold**. How to **maximize the power of our army** considering these resources? We could simply find the unit with the **best power/cost ratio**, **take as many of them as possible**, and repeat the process with the other two units. But this \"guess and check\" solution **might not even be optimal**. Now imagine we have **millions of units and resources**: this greedy strategy is very likely to **completely miss the optimal solution**. It is possible to **use a machine learning** (ML) algorithm (a metaheuristic like a genetic algorithm) to solve this problem, but we have **no guarantee that the solution will be optimal either**.\n",
        "\n",
        "Fortunately for us, there is a method that can **solve our problem in an optimal way**: linear programming (or linear optimization), which is part of the field of [operations research](https://en.wikipedia.org/wiki/Operations_research) (OR). In this article, we're going to **find the best solution possible to this problem**. We'll use OR-Tools, a **powerful Python library** specialized in optimization. It will help us to **model this problem** in linear programming terms and **output an optimal solution**. Finally, we will write a **model that can take on a bigger challenge** and actually **solve a whole class of optimization problems**.\n",
        "\n",
        "## 🥇 I. Finding the optimal solution\n",
        "\n",
        "Okay, let's play with linear programming in Python. There are **different libraries** such as the multi-purposed `SciPy`, the beginner-friendly `PuLP`, the exhaustive `pyomo`, and many others. Today, we are going to use [Google OR-Tools](https://developers.google.com/optimization), which is quite **user-friendly**, comes with several **prepackaged solvers**, and has by far the **most stars** on [GitHub](https://github.com/google/or-tools). If the installation doesn't work, please **restart the kernel and try again**: it can fail sometimes. ¯\\\\\\_(ツ)\\_/¯"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 1,
      "metadata": {
        "id": "ZK0c0hhQklGb"
      },
      "outputs": [],
      "source": [
        "!python -m pip install --upgrade --user -q ortools"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "eGgO_cnFk9ce"
      },
      "source": [
        "All these libraries have a **hidden benefit**: they act as **interfaces** to **use the same model with different solvers**. Solvers like [Gurobi](https://www.gurobi.com/), [Cplex](https://www.ibm.com/analytics/cplex-optimizer), or [SCIP](https://www.scipopt.org/) have their own APIs, but the models they create are **tied to a specific solver**. OR-Tools allows us to use an abstract and quite pythonic way of **modeling our problems**, and then choose **one or several solvers** to find an optimal solution.\n",
        "\n",
        "<center class=\"dark-bg\">\n",
        "    <video autoplay loop muted width=\"800\" class=\"no-margin\">\n",
        "        <source src=\"https://mlabonne.github.io/blog/images/ortools/ortools_solvers.mp4\" type=\"video/mp4\"></source>\n",
        "    </video>\n",
        "</center>\n",
        "\n",
        "OR-Tools comes with its **own linear programming solver**, called **GLOP** for Google Linear Optimization Package. It is an open-source project created by Google's Operations Research Team and written in C++. According to Google, it is used to [stabilize YouTube videos](https://ai.googleblog.com/2014/09/sudoku-linear-optimization-and-ten-cent.html) and in various projects such as a [cooperative multi-agent reinforcement learning architecture](https://proceedings.neurips.cc/paper/2019/hash/3c3c139bd8467c1587a41081ad78045e-Abstract.html). Other solvers are available, such as **SCIP**, an excellent **non-commercial solver** created in 2005 and **updated and maintained** to this day. We could also use commercial options like Gurobi and Cplex but we would need to **install them on top of OR-Tools** and get the **appropriate licenses** (which can be quite costly). For now, let's try GLOP."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 2,
      "metadata": {
        "id": "l4RYui9hxcJG"
      },
      "outputs": [],
      "source": [
        "# Import OR-Tools wrapper for linear programming\n",
        "from ortools.linear_solver import pywraplp\n",
        "\n",
        "# Create a solver using the GLOP backend\n",
        "solver = pywraplp.Solver('Maximize army power', pywraplp.Solver.GLOP_LINEAR_PROGRAMMING)"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "GGx41ZONYF17"
      },
      "source": [
        "We created an **instance** of the OR-Tools solver using GLOP. The first thing we want to define is the **variables we want to optimize**. In our example, we have three variables: the number of 🗡️**swordsmen**, 🏹**bowmen**, and 🐎**horsemen** in the army. OR-Tools accepts three types of variables:\n",
        "\n",
        "* [`NumVar`](https://google.github.io/or-tools/python/ortools/linear_solver/pywraplp.html#Solver.NumVar) for **continuous variables**;\n",
        "* [`IntVar`](https://google.github.io/or-tools/python/ortools/linear_solver/pywraplp.html#Solver.IntVar) for **integer variables**;\n",
        "* [`BoolVar`](https://google.github.io/or-tools/python/ortools/linear_solver/pywraplp.html#Solver.BoolVar) for **boolean variables**.\n",
        "\n",
        "We're looking for **round numbers of units**, so let's choose `IntVar`. We then need to specify **lower and upper bounds** for these variables. We want at least 0 unit, but we **don't really have an upper bound**. So we can say that our **upper bound is infinity** and use `solver.infinity()`. Actually $\\infty$ might be a **little too much considering the number of resources we have**. We could also take $1000$ (for example) as an upper bound, since **we can't go that high anyway**. But let's keep $∞$ for now.\n",
        "\n",
        "Let's call the number of **🗡️swordsmen** $swordsmen$, the number of **🏹bowmen** $bowmen$, and the number of **🐎horsemen** $horsemen$. We can formally write:\n",
        "\n",
        "$$0 \\leq swordsmen < \\infty \\\\\n",
        "0 \\leq bowmen < \\infty \\\\\n",
        "0 \\leq horsemen < \\infty$$\n",
        "\n",
        "Good, let's translate it into code. The syntax is **quite straightforward** with OR-Tools."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 3,
      "metadata": {
        "id": "7Oqc6nBg0ea2"
      },
      "outputs": [],
      "source": [
        "# Create the variables we want to optimize\n",
        "swordsmen = solver.IntVar(0, solver.infinity(), 'swordsmen')\n",
        "bowmen = solver.IntVar(0, solver.infinity(), 'bowmen')\n",
        "horsemen = solver.IntVar(0, solver.infinity(), 'horsemen')"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "Yy-N4FHpY3AP"
      },
      "source": [
        "We defined our variables, but another essential parameter is the **constraints**. In our case, we have a **limited number of resources** we can use to **produce units**. In other words, **we can't spend more resources than we have**: for instance, the 🌾**food** spent to recruit units **cannot be higher than 1200**. The same is true with 🪵**wood** (800) and  🪙**gold** (600).\n",
        "\n",
        "According to our table, units have the following costs:\n",
        "\n",
        "* **1 swordsman** = 🌾60 + 🪵20;\n",
        "* **1 bowman** = 🌾80 + 🪵10 + 🪙40;\n",
        "* **1 horseman** = 🌾140 + 🪙100.\n",
        "\n",
        "We can write one constraint per resource as follows:\n",
        "\n",
        "$$60\\times swordsmen + 80\\times bowmen + 140\\times horsemen \\leq 1200 \\\\\n",
        "20\\times swordsmen + 10\\times bowmen \\leq 800 \\\\\n",
        "40\\times bowmen + 100\\times horsemen \\leq 600$$\n",
        "\n",
        "In OR-Tools, we simply add the constraints to our solver instance with [`solver.Add()`](https://google.github.io/or-tools/python/ortools/linear_solver/pywraplp.html#Solver.Add)."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 4,
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/"
        },
        "id": "jtD0ZPevY5CO",
        "outputId": "8f09d173-9ef9-4877-939e-9f1a80cd670f"
      },
      "outputs": [
        {
          "output_type": "execute_result",
          "data": {
            "text/plain": [
              "<ortools.linear_solver.pywraplp.Constraint; proxy of <Swig Object of type 'operations_research::MPConstraint *' at 0x7fd94c30bfc0> >"
            ]
          },
          "metadata": {},
          "execution_count": 4
        }
      ],
      "source": [
        "# Add constraints for each resource\n",
        "solver.Add(swordsmen*60 + bowmen*80 + horsemen*140 <= 1200) # Food\n",
        "solver.Add(swordsmen*20 + bowmen*10 <= 800) # Wood\n",
        "solver.Add(bowmen*40 + horsemen*100 <= 600) # Gold"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "PXV0Gg30Y5UD"
      },
      "source": [
        "Now that we have our variables and constraints, we want to **define our goal** (or **objective function**). In linear programming, this function **has to be linear**, so of the form $f(x, y, z) = ax + by + cz + d$. In our example, the objective is quite clear: we want to **recruit the army with the highest power**. The table gives us the following power values:\n",
        "\n",
        "* **1 swordsman** = 💪70;\n",
        "* **1 bowman** = 💪95;\n",
        "* **1 horseman** = 💪230.\n",
        "\n",
        "Maximizing the power of the army amounts to **maximizing the sum of the power of each unit**. Let's call $units$ the tuple $(swordsmen, bowmen, horsemen)$, where every element is the number of the corresponding unit. Our objective function can be written as:\n",
        "\n",
        "$$max \\ f(units) = 70\\times swordsmen + 95\\times bowmen + 230\\times horsemen$$\n",
        "\n",
        "In general, there are two types of objective functions: **minimizing** and **maximizing**. In OR-Tools, we declare this goal with [`solver.Maximize()`](https://google.github.io/or-tools/python/ortools/linear_solver/pywraplp.html#Solver.Maximize) or [`solver.Minimize()`](https://google.github.io/or-tools/python/ortools/linear_solver/pywraplp.html#Solver.Minimize)."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 5,
      "metadata": {
        "id": "HLkXypIKY8Bm"
      },
      "outputs": [],
      "source": [
        "# Maximize the objective function\n",
        "solver.Maximize(swordsmen*70 + bowmen*95 + horsemen*230)"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "wCCbr-faY8d-"
      },
      "source": [
        "And we're done! There are **three steps to model an optimization problem**:\n",
        "\n",
        "1. **Declaring the variables to optimize** with lower and upper bounds;\n",
        "2. **Adding constraints** to these variables;\n",
        "3. **Defining the objective function** to maximize or to minimize.\n",
        "\n",
        "Now that is clear, we can ask the solver to find an optimal solution for us. This is done with [`solver.Solve()`](https://google.github.io/or-tools/python/ortools/linear_solver/pywraplp.html#Solver.Solve) and returns a status. This status can then be used **to check that the solution is indeed optimal**. Let's print the **highest total power** we can get with the **best army configuration**."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 6,
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/"
        },
        "id": "ddA-aNJMY-V6",
        "outputId": "92ffacd0-4a0e-4b15-c0e7-e6265c3194d9"
      },
      "outputs": [
        {
          "output_type": "stream",
          "name": "stdout",
          "text": [
            "================= Solution =================\n",
            "Solved in 66.00 milliseconds in 2 iterations\n",
            "\n",
            "Optimal power = 1800.0 💪power\n",
            "Army:\n",
            " - 🗡️Swordsmen = 6.0000000000000036\n",
            " - 🏹Bowmen = 0.0\n",
            " - 🐎Horsemen = 5.999999999999999\n"
          ]
        }
      ],
      "source": [
        "# Solve problem\n",
        "status = solver.Solve()\n",
        "\n",
        "# If an optimal solution has been found, print results\n",
        "if status == pywraplp.Solver.OPTIMAL:\n",
        "  print('================= Solution =================')\n",
        "  print(f'Solved in {solver.wall_time():.2f} milliseconds in {solver.iterations()} iterations')\n",
        "  print()\n",
        "  print(f'Optimal power = {solver.Objective().Value()} 💪power')\n",
        "  print('Army:')\n",
        "  print(f' - 🗡️Swordsmen = {swordsmen.solution_value()}')\n",
        "  print(f' - 🏹Bowmen = {bowmen.solution_value()}')\n",
        "  print(f' - 🐎Horsemen = {horsemen.solution_value()}')\n",
        "else:\n",
        "  print('The solver could not find an optimal solution.')"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "H_QohPhkEKs6"
      },
      "source": [
        "Great! The solver found an optimal solution: our army has a **total power of 💪1800 with 6 🗡️swordsmen and 6 🐎horsemen** (sorry bowmen!).\n",
        "\n",
        "Let's unpack this result: the solver decided to take the **maximum number of 🐎horsemen** (6, since we only have 🪙600 and they each cost 🪙100). The remaining resources are **spent in 🗡️swordsmen**: we have $1200 - 6 \\times 140 = 360$ 🌾food left, which is why the solver **chose 6 🗡️swordsmen**. We can deduce that the **horsemen are the best unit** and the **bowmen are the worst one** because they haven't been chosen at all.\n",
        "\n",
        "Okay, but there's something quite weird: these numbers are **not round**, even though we specified that we wanted **integers**. So what happened? This strange behavior is due to the GLOP solver: it is a **linear programming** (LP) solver, and not an **integer linear programming** (ILP) one. In summary, it can **only use real numbers and not integers** as variables. So why did we **declare our variables as integers** if it doesn't take it into account?\n",
        "\n",
        "GLOP cannot solve ILP problems, but other solvers can. Actually, a lot of them are **mixed integer linear programming** (MILP, commonly called MIP) solvers. This means that they can consider both **continuous** (real numbers) and **discrete** (integers) variables. A particular case of discrete values is **boolean variables** to represent decisions with 0-1 values.\n",
        "\n",
        "We talked about **SCIP** earlier, which is a good example since **it can solve both MILP and MINLP** (mixed integer **nonlinear** programming) problems. Another potential candidate is [CBC](https://projects.coin-or.org/Cbc), an open-source **MIP solver directly available through OR-Tools**. Thanks to this library, we can **use the same model and just change the solver** to SCIP or CBC."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 7,
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/"
        },
        "id": "Dvxq9hyCk8uW",
        "outputId": "a2cef1c4-135e-4b94-978a-9926361e6a61"
      },
      "outputs": [
        {
          "output_type": "stream",
          "name": "stdout",
          "text": [
            "================= Solution =================\n",
            "Solved in 6.00 milliseconds in 0 iterations\n",
            "\n",
            "Optimal value = 1800.0 💪power\n",
            "Army:\n",
            " - 🗡️Swordsmen = 6.0\n",
            " - 🏹Bowmen = 0.0\n",
            " - 🐎Horsemen = 6.0\n"
          ]
        }
      ],
      "source": [
        "# Create the linear solver using the CBC backend\n",
        "solver = pywraplp.Solver('Maximize army power', pywraplp.Solver.CBC_MIXED_INTEGER_PROGRAMMING)\n",
        "\n",
        "# 1. Create the variables we want to optimize\n",
        "swordsmen = solver.IntVar(0, solver.infinity(), 'swordsmen')\n",
        "bowmen = solver.IntVar(0, solver.infinity(), 'bowmen')\n",
        "horsemen = solver.IntVar(0, solver.infinity(), 'horsemen')\n",
        "\n",
        "# 2. Add constraints for each resource\n",
        "solver.Add(swordsmen*60 + bowmen*80 + horsemen*140 <= 1200)\n",
        "solver.Add(swordsmen*20 + bowmen*10 <= 800)\n",
        "solver.Add(bowmen*40 + horsemen*100 <= 600)\n",
        "\n",
        "# 3. Maximize the objective function\n",
        "solver.Maximize(swordsmen*70 + bowmen*95 + horsemen*230)\n",
        "\n",
        "# Solve problem\n",
        "status = solver.Solve()\n",
        "\n",
        "# If an optimal solution has been found, print results\n",
        "if status == pywraplp.Solver.OPTIMAL:\n",
        "  print('================= Solution =================')\n",
        "  print(f'Solved in {solver.wall_time():.2f} milliseconds in {solver.iterations()} iterations')\n",
        "  print()\n",
        "  print(f'Optimal value = {solver.Objective().Value()} 💪power')\n",
        "  print('Army:')\n",
        "  print(f' - 🗡️Swordsmen = {swordsmen.solution_value()}')\n",
        "  print(f' - 🏹Bowmen = {bowmen.solution_value()}')\n",
        "  print(f' - 🐎Horsemen = {horsemen.solution_value()}')\n",
        "else:\n",
        "  print('The solver could not find an optimal solution.')"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "vHTmOQzUK295"
      },
      "source": [
        "Strictly speaking, our variables are still floats (`type(swordsmen.solution_value()) = float`) but we can see that **they don't have weird decimals anymore**: the CBC solver really considered them as **integers**. In general, we can just **round up these values** since the error is insignificant, but it is important to remember to choose the appropriate solver according to the **studied problem**: **LP** (continuous variables) or **MIP** (combination of continuous and discrete variables). There are other types such as **quadratic** (QP) or **nonlinear** (NLP or MINLP, with an **exponential objective function or constraints** for instance) problems.\n",
        "\n",
        "<center><div id=\"anim_container\"><img alt=\"LP, MIP, etc.\" src=\"https://mlabonne.github.io/blog/images/ortools/continuousvsdiscrete.png\" loading=\"lazy\" id=\"anim_inner\"></div></center>\n",
        "\n",
        "## 🧱 III. Building a general model\n",
        "\n",
        "But what if our **resources change**? Or if the cost of a unit **evolved**? What if we upgraded horsemen and **their power increased**? One of the **best perks** of OR-Tools is that it uses a **general-purpose programming language** like Python. Instead of **static numbers**, we can store our parameters in **objects like dictionaries or lists**. The code won't be **as readable**, but it becomes **much more flexible**: actually, it can be so flexible that we can **solve an entire class of optimization problems without changing the model** (just the parameters).\n",
        "\n",
        "Let's transform our **input parameters into Python lists** and feed them to the solver **through a function**."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 8,
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/"
        },
        "id": "qR-fbnruehR4",
        "outputId": "13b48720-feb0-4f5b-8fa6-713981621491"
      },
      "outputs": [
        {
          "output_type": "stream",
          "name": "stdout",
          "text": [
            "================= Solution =================\n",
            "Solved in 3.00 milliseconds in 0 iterations\n",
            "\n",
            "Optimal value = 1800.0 💪power\n",
            "Army:\n",
            " - 🗡️Swordsmen = 6.0\n",
            " - 🏹Bowmen = 0.0\n",
            " - 🐎Horsemen = 6.0\n"
          ]
        }
      ],
      "source": [
        "# Inputs\n",
        "UNITS = ['🗡️Swordsmen', '🏹Bowmen', '🐎Horsemen']\n",
        "\n",
        "DATA = [[60, 20, 0, 70],\n",
        "        [80, 10, 40, 95],\n",
        "        [140, 0, 100, 230]]\n",
        "\n",
        "RESOURCES = [1200, 800, 600]\n",
        "\n",
        "\n",
        "def solve_army(UNITS, DATA, RESOURCES):\n",
        "  # Create the linear solver using the CBC backend\n",
        "  solver = pywraplp.Solver('Maximize army power', pywraplp.Solver.CBC_MIXED_INTEGER_PROGRAMMING)\n",
        "\n",
        "  # 1. Create the variables we want to optimize\n",
        "  units = [solver.IntVar(0, solver.infinity(), unit) for unit in UNITS]\n",
        "\n",
        "  # 2. Add constraints for each resource\n",
        "  for r, _ in enumerate(RESOURCES):\n",
        "    solver.Add(sum(DATA[u][r] * units[u] for u, _ in enumerate(units)) <= RESOURCES[r])\n",
        "\n",
        "  # 3. Maximize the objective function\n",
        "  solver.Maximize(sum(DATA[u][-1] * units[u] for u, _ in enumerate(units)))\n",
        "\n",
        "  # Solve problem\n",
        "  status = solver.Solve()\n",
        "\n",
        "  # If an optimal solution has been found, print results\n",
        "  if status == pywraplp.Solver.OPTIMAL:\n",
        "    print('================= Solution =================')\n",
        "    print(f'Solved in {solver.wall_time():.2f} milliseconds in {solver.iterations()} iterations')\n",
        "    print()\n",
        "    print(f'Optimal value = {solver.Objective().Value()} 💪power')\n",
        "    print('Army:')\n",
        "    for u, _ in enumerate(units):\n",
        "      print(f' - {units[u].name()} = {units[u].solution_value()}')\n",
        "  else:\n",
        "    print('The solver could not find an optimal solution.')\n",
        "\n",
        "solve_army(UNITS, DATA, RESOURCES)"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "qj6lzFPsT6MS"
      },
      "source": [
        "Good, the results are the same: our code seems to work. Now let's **change the parameters** to tackle a slightly more complex problem. \n",
        "\n",
        "Imagine we have a **lot more resources**: 🌾183000, 🪵90512, and 🪙80150, so we can also **produce a lot more units**! This is the new table:\n",
        "<br/>\n",
        "<br/>\n",
        "\n",
        "| Unit | 🌾Food | 🪵Wood | 🪙Gold | 💪Attack | ❤️Health\n",
        "| :--- | :---: | :---: | :---: | :---: | :---: |\n",
        "| 🗡️Swordsman      | 60  | 20  | 0   | 6   | 70 |\n",
        "| 🛡️Man-at-arms    | 100 | 0   | 20  | 12  | 155 |\n",
        "| 🏹Bowman         | 30  | 50  | 0   | 5   | 70 |\n",
        "| ❌Crossbowman    | 80  | 0   | 40  | 12  | 80 |\n",
        "| 🔫Handcannoneer  | 120 | 0   | 120 | 35  | 150 |\n",
        "| 🐎Horseman       | 100 | 20  | 0   | 9   | 125 |\n",
        "| ♞Knight          | 140 | 0   | 100 | 24  | 230 | \n",
        "| 🐏Battering ram  | 0   | 300 | 0   | 200 | 700 |\n",
        "| 🎯Springald      | 0   | 250 | 250 | 30 | 200 |\n",
        "\n",
        "<br/>\n",
        "\n",
        "Notice that we transformed the 💪**power** into two values: 💪**attack** and ❤️ **health**, which is a little more detailed. Health values **are higher than attack values**, which is why we want to **add a weighting factor to make them more comparable**. Let's take 10 as an example, so $power = 10 \\times attack + health$. We call the tuple of variables to optimize $units = (swordsmen,  \\dots, springalds)$, where every element is the number of the corresponding unit. Our objective function becomes:\n",
        "\n",
        "$$max \\ f(units) = \\sum_{u \\in units} (10\\times power + health) \\cdot u$$\n",
        "\n",
        "Adapting our code to this new problem is actually quite simple: we just have to **change the input parameters** and update the **objective function**."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 9,
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/"
        },
        "id": "Y2EUFqCETvLZ",
        "outputId": "7b00144f-f2a1-49d3-c386-463f0b8c370f"
      },
      "outputs": [
        {
          "output_type": "stream",
          "name": "stdout",
          "text": [
            "================= Solution =================\n",
            "Solved in 136.00 milliseconds in 412 iterations\n",
            "\n",
            "Optimal value = 1393145.0 💪power\n",
            "Army:\n",
            " - 🗡️Swordsmen = 2.0\n",
            " - 🛡️Men-at-arms = 1283.0\n",
            " - 🏹Bowmen = 3.0\n",
            " - ❌Crossbowmen = 0.0\n",
            " - 🔫Handcannoneers = 454.0\n",
            " - 🐎Horsemen = 0.0\n",
            " - ♞Knights = 0.0\n",
            " - 🐏Battering rams = 301.0\n",
            " - 🎯Springalds = 0.0\n",
            " - 🪨Mangonels = 0.0\n"
          ]
        }
      ],
      "source": [
        "UNITS = [\n",
        "    '🗡️Swordsmen',\n",
        "    '🛡️Men-at-arms',\n",
        "    '🏹Bowmen',\n",
        "    '❌Crossbowmen',\n",
        "    '🔫Handcannoneers',\n",
        "    '🐎Horsemen',\n",
        "    '♞Knights',\n",
        "    '🐏Battering rams',\n",
        "    '🎯Springalds',\n",
        "    '🪨Mangonels',\n",
        "]\n",
        "\n",
        "DATA = [\n",
        "    [60, 20, 0, 6, 70],\n",
        "    [100, 0, 20, 12, 155],\n",
        "    [30, 50, 0, 5, 70],\n",
        "    [80, 0, 40, 12, 80],\n",
        "    [120, 0, 120, 35, 150],\n",
        "    [100, 20, 0, 9, 125],\n",
        "    [140, 0, 100, 24, 230],\n",
        "    [0, 300, 0, 200, 700],\n",
        "    [0, 250, 250, 30, 200],\n",
        "    [0, 400, 200, 12*3, 240]\n",
        "]\n",
        "\n",
        "RESOURCES = [183000, 90512, 80150]\n",
        "\n",
        "\n",
        "def solve_army(UNITS, DATA, RESOURCES):\n",
        "  # Create the linear solver using the CBC backend\n",
        "  solver = pywraplp.Solver('Maximize army power', pywraplp.Solver.CBC_MIXED_INTEGER_PROGRAMMING)\n",
        "\n",
        "  # 1. Create the variables we want to optimize\n",
        "  units = [solver.IntVar(0, solver.infinity(), unit) for unit in UNITS]\n",
        "\n",
        "  # 2. Add constraints for each resource\n",
        "  for r, _ in enumerate(RESOURCES):\n",
        "    solver.Add(sum(DATA[u][r] * units[u] for u, _ in enumerate(units)) <= RESOURCES[r])\n",
        "\n",
        "  # 3. Maximize the new objective function\n",
        "  solver.Maximize(sum((10*DATA[u][-2] + DATA[u][-1]) * units[u] for u, _ in enumerate(units)))\n",
        "\n",
        "  # Solve problem\n",
        "  status = solver.Solve()\n",
        "\n",
        "  # If an optimal solution has been found, print results\n",
        "  if status == pywraplp.Solver.OPTIMAL:\n",
        "    print('================= Solution =================')\n",
        "    print(f'Solved in {solver.wall_time():.2f} milliseconds in {solver.iterations()} iterations')\n",
        "    print()\n",
        "    print(f'Optimal value = {solver.Objective().Value()} 💪power')\n",
        "    print('Army:')\n",
        "    for u, _ in enumerate(units):\n",
        "      print(f' - {units[u].name()} = {units[u].solution_value()}')\n",
        "  else:\n",
        "    print('The solver could not find an optimal solution.')\n",
        "\n",
        "solve_army(UNITS, DATA, RESOURCES)"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "ID8pH6oQmKow"
      },
      "source": [
        "This problem would **take a long time for humans to address**, but the ILP solver did it **in the blink of an eye**. Better than that: it also gives us the **guarantee that our solution is optimal**, which means that our enemy cannot find a better army compositon for the same cost! We could **increase the number of units**, give **billions of resources** but you get the picture: it would **just take longer to obtain a solution**, but it **wouldn't change the problem**.\n",
        "\n",
        "## ⚔️ IV. Time for war\n",
        "\n",
        "Now, let's say we scouted our enemy and know that their **army has a 💪power of 1,000,000**. We could build a **much better army**, but our resources are precious and it **wouldn't be very efficient**: all we have to do is to build an army with a 💪power **higher than 1,000,000** (even 1,000,001 would be enough). In other words, it means we're looking for the **most cost-efficient unit** and we want the **minimum number** of them to defeat the enemy army. Knowing the **number of units to produce** will also tell us **how many resources we need**.\n",
        "\n",
        "We can reuse our **input parameters** since they didn't change. But the constraint is different this time: we want a **💪power higher than 1,000,000**. It means that the **sum of the power of the selected units** must be higher than this number.\n",
        "\n",
        "$$\\sum_{u \\in units} (10\\times attack + health) \\cdot u > 1\\,000\\,000$$\n",
        "\n",
        "In code, we can **loop through our units and resources** to design this constraint.\n",
        "\n",
        "The **objective function** also has to change. Our goal is to **minimize the sum of resources spent to build the army**.\n",
        "\n",
        "$$min \\ f(units) = \\sum_{u \\in units} (food + wood + gold) \\cdot u$$\n",
        "\n",
        "Once again, we can loop **through our resources** to implement it in OR-Tools."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 10,
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/"
        },
        "id": "nz5244z9nBYw",
        "outputId": "3291938f-e9da-4d4c-aa2f-a5e4e367e3d8"
      },
      "outputs": [
        {
          "output_type": "stream",
          "name": "stdout",
          "text": [
            "================= Solution =================\n",
            "Solved in 12.00 milliseconds in 0 iterations\n",
            "\n",
            "Optimal value = 111300.0 🌾🪵🪙resources\n",
            "Power = 💪1001700.0\n",
            "Army:\n",
            " - 🗡️Swordsmen = 0.0\n",
            " - 🛡️Men-at-arms = 0.0\n",
            " - 🏹Bowmen = 0.0\n",
            " - ❌Crossbowmen = 0.0\n",
            " - 🔫Handcannoneers = 0.0\n",
            " - 🐎Horsemen = 0.0\n",
            " - ♞Knights = 0.0\n",
            " - 🐏Battering rams = 371.0\n",
            " - 🎯Springalds = 0.0\n",
            " - 🪨Mangonels = 0.0\n",
            "\n",
            "Resources:\n",
            " - 🌾Food = 0.0\n",
            " - 🪵Wood = 111300.0\n",
            " - 🪙Gold = 0.0\n"
          ]
        }
      ],
      "source": [
        "def solve_army(UNITS, DATA, RESOURCES):\n",
        "  # Create the linear solver using the CBC backend\n",
        "  solver = pywraplp.Solver('Minimize resource consumption', pywraplp.Solver.CBC_MIXED_INTEGER_PROGRAMMING)\n",
        "\n",
        "  # 1. Create the variables we want to optimize\n",
        "  units = [solver.IntVar(0, solver.infinity(), unit) for unit in UNITS]\n",
        "\n",
        "  # 2. Add constraints for each resource\n",
        "  for r, _ in enumerate(RESOURCES):\n",
        "    solver.Add(sum((10 * DATA[u][-2] + DATA[u][-1]) * units[u] for u, _ in enumerate(units)) >= 1000001)\n",
        "\n",
        "  # 3. Minimize the objective function\n",
        "  solver.Minimize(sum((DATA[u][0] + DATA[u][1] + DATA[u][2]) * units[u] for u, _ in enumerate(units)))\n",
        "\n",
        "  # Solve problem\n",
        "  status = solver.Solve()\n",
        "\n",
        "  # If an optimal solution has been found, print results\n",
        "  if status == pywraplp.Solver.OPTIMAL:\n",
        "    print('================= Solution =================')\n",
        "    print(f'Solved in {solver.wall_time():.2f} milliseconds in {solver.iterations()} iterations')\n",
        "    print()\n",
        "\n",
        "    power = sum((10 * DATA[u][-2] + DATA[u][-1]) * units[u].solution_value() for u, _ in enumerate(units))\n",
        "    print(f'Optimal value = {solver.Objective().Value()} 🌾🪵🪙resources')\n",
        "    print(f'Power = 💪{power}')\n",
        "    print('Army:')\n",
        "    for u, _ in enumerate(units):\n",
        "      print(f' - {units[u].name()} = {units[u].solution_value()}')\n",
        "    print()\n",
        "\n",
        "    food = sum((DATA[u][0]) * units[u].solution_value() for u, _ in enumerate(units))\n",
        "    wood = sum((DATA[u][1]) * units[u].solution_value() for u, _ in enumerate(units))\n",
        "    gold = sum((DATA[u][2]) * units[u].solution_value() for u, _ in enumerate(units))\n",
        "    print('Resources:')\n",
        "    print(f' - 🌾Food = {food}')\n",
        "    print(f' - 🪵Wood = {wood}')\n",
        "    print(f' - 🪙Gold = {gold}')\n",
        "  else:\n",
        "      print('The solver could not find an optimal solution.')\n",
        "\n",
        "solve_army(UNITS, DATA, RESOURCES)"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "Qxeyd7vV3YmH"
      },
      "source": [
        "The solver found an optimal solution: we need to **build 371 🐏battering rams for a total cost of 111,300 🪵wood**. Wait, what if **we don't have that much wood**? In the previous section, we only had 🪵90512: we cannot produce 371 🐏battering rams. 😱\n",
        "\n",
        "So is it possible to **take these limited resources into account** and still try to build the **best army**? Actually, it's super easy: we just have to **copy/paste the constraints from the previous section**."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 11,
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/"
        },
        "id": "lyi0eVvu1I8Q",
        "outputId": "2569a958-13da-40d5-e0a2-c9e623bd1f7c"
      },
      "outputs": [
        {
          "output_type": "stream",
          "name": "stdout",
          "text": [
            "================= Solution =================\n",
            "Solved in 47.00 milliseconds in 1 iterations\n",
            "\n",
            "Optimal value = 172100.0 🌾🪵🪙resources\n",
            "Power = 💪1000105.0\n",
            "Army:\n",
            " - 🗡️Swordsmen = 1.0\n",
            " - 🛡️Men-at-arms = 681.0\n",
            " - 🏹Bowmen = 0.0\n",
            " - ❌Crossbowmen = 0.0\n",
            " - 🔫Handcannoneers = 0.0\n",
            " - 🐎Horsemen = 0.0\n",
            " - ♞Knights = 0.0\n",
            " - 🐏Battering rams = 301.0\n",
            " - 🎯Springalds = 0.0\n",
            " - 🪨Mangonels = 0.0\n",
            "\n",
            "Resources:\n",
            " - 🌾Food = 68160.0\n",
            " - 🪵Wood = 90320.0\n",
            " - 🪙Gold = 13620.0\n"
          ]
        }
      ],
      "source": [
        "def solve_army(UNITS, DATA, RESOURCES):\n",
        "  # Create the linear solver using the CBC backend\n",
        "  solver = pywraplp.Solver('Minimize resource consumption', pywraplp.Solver.CBC_MIXED_INTEGER_PROGRAMMING)\n",
        "\n",
        "  # 1. Create the variables we want to optimize\n",
        "  units = [solver.IntVar(0, solver.infinity(), unit) for unit in UNITS]\n",
        "\n",
        "  # 2. Add constraints for each resource\n",
        "  for r, _ in enumerate(RESOURCES):\n",
        "    solver.Add(sum((10 * DATA[u][-2] + DATA[u][-1]) * units[u] for u, _ in enumerate(units)) >= 1000001)\n",
        "\n",
        "  # Old constraints for limited resources\n",
        "  for r, _ in enumerate(RESOURCES):\n",
        "    solver.Add(sum(DATA[u][r] * units[u] for u, _ in enumerate(units)) <= RESOURCES[r])\n",
        "\n",
        "  # 3. Minimize the objective function\n",
        "  solver.Minimize(sum((DATA[u][0] + DATA[u][1] + DATA[u][2]) * units[u] for u, _ in enumerate(units)))\n",
        "\n",
        "  # Solve problem\n",
        "  status = solver.Solve()\n",
        "\n",
        "  # If an optimal solution has been found, print results\n",
        "  if status == pywraplp.Solver.OPTIMAL:\n",
        "    print('================= Solution =================')\n",
        "    print(f'Solved in {solver.wall_time():.2f} milliseconds in {solver.iterations()} iterations')\n",
        "    print()\n",
        "\n",
        "    power = sum((10 * DATA[u][-2] + DATA[u][-1]) * units[u].solution_value() for u, _ in enumerate(units))\n",
        "    print(f'Optimal value = {solver.Objective().Value()} 🌾🪵🪙resources')\n",
        "    print(f'Power = 💪{power}')\n",
        "    print('Army:')\n",
        "    for u, _ in enumerate(units):\n",
        "      print(f' - {units[u].name()} = {units[u].solution_value()}')\n",
        "    print()\n",
        "    \n",
        "    food = sum((DATA[u][0]) * units[u].solution_value() for u, _ in enumerate(units))\n",
        "    wood = sum((DATA[u][1]) * units[u].solution_value() for u, _ in enumerate(units))\n",
        "    gold = sum((DATA[u][2]) * units[u].solution_value() for u, _ in enumerate(units))\n",
        "    print('Resources:')\n",
        "    print(f' - 🌾Food = {food}')\n",
        "    print(f' - 🪵Wood = {wood}')\n",
        "    print(f' - 🪙Gold = {gold}')\n",
        "  else:\n",
        "      print('The solver could not find an optimal solution.')\n",
        "\n",
        "solve_army(UNITS, DATA, RESOURCES)"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "W2oJLwbInCIf"
      },
      "source": [
        "Since we now have a limited resource of 🪵wood, the number of 🐏battering rams **sadly dropped from 371 to 301**. In exchange, we got **681 🛡️men-at-arms and 1 lost 🗡️swordsman** (welcome to them). The total cost of the army is **172,100**, which is **much higher** than the 111,300 we previously found (+65% increase) but it **truly is the optimal solution under these constraints**. It shows that **we should produce more wood** because these 🐏battering rams are extremely cost-efficient!\n",
        "\n",
        "## 🔚 V. Conclusion\n",
        "\n",
        "Optimization is often **neglected in favor of machine learning techniques**, but both have their merits: linear programming can produce an **optimal solution in an undetermined amount of time**, while machine learning can **approximate complex functions in no time**. There is no training in LP, but an **expert is required to build a mathematical model**. Machine learning needs data, but the models can be **used as black boxes** to solve a problem. As a rule of thumb, problems that **do not have a particular time constraint and/or are not extremely complex** can be advantageously solved with **linear programming**.\n",
        "\n",
        "In this article,\n",
        "\n",
        "* We learned about **interfaces for optimization**, and especially about Google **OR-Tools**;\n",
        "* We talked about **solvers and types of optimization problems**: LP, MIP, NLP;\n",
        "* We modeled and solved an **extremely common optimization problem** in an **optimal way** and **generalized our model** through a function;\n",
        "* We **reframed this problem** and merged two sets of constraints to obtain **the best army composition for the lowest price**.\n",
        "\n",
        "There are a lot more problems where optimization can be applied: for instance, how to **create school timetables that satisfy everybody's requirements**? How to **deliver 1,000 different orders in a minimum amount of time**? Where to **create a new metro line to maximize its usefulness**? In future articles, we'll talk about **new types of applications** for these techniques, including satisfiability and nonlinear problems.\n",
        "\n",
        "I hope you enjoyed this article! Feel free to **share it and spread the knowledge about linear optimization**. Don't forget to **[follow me on Twitter](https://twitter.com/maximelabonne)** where I post summaries of these articles. Cheers!"
      ]
    }
  ],
  "metadata": {
    "colab": {
      "collapsed_sections": [],
      "name": "01: Google_OR-Tools_demo_Build_efficient_powerful_army.ipynb",
      "provenance": [],
      "include_colab_link": true
    },
    "kernelspec": {
      "display_name": "Python 3",
      "name": "python3"
    },
    "language_info": {
      "name": "python",
      "version": "3.8.12"
    }
  },
  "nbformat": 4,
  "nbformat_minor": 0
}