{
    "cells": [
        {
            "cell_type": "markdown",
            "metadata": {
                "colab_type": "text",
                "id": "9O_tHTY_PL_s"
            },
            "source": [
                "```{index} disjunctive programming\n",
                "```\n",
                "```{index} single: solver; cbc\n",
                "```\n",
                "```{index} single: application; maintenance planning\n",
                "```\n",
                "\n",
                "# Extra material: Maintenance planning"
            ]
        },
        {
            "cell_type": "markdown",
            "metadata": {
                "colab_type": "text",
                "id": "Y_yzzUJQDJtr"
            },
            "source": [
                "## Problem statement\n",
                "\n",
                "\n",
                "A process unit is operating over a maintenance planning horizon from $1$ to $T$ days.  On day $t$ the unit makes a profit $c[t]$ which is known in advance. The unit needs to shut down for $P$ maintenance periods during the planning period.  Once started, a maintenance period takes $M$ days to finish.\n",
                "\n",
                "Find a maintenance schedule that allows the maximum profit to be produced."
            ]
        },
        {
            "cell_type": "markdown",
            "metadata": {},
            "source": [
                "## Imports"
            ]
        },
        {
            "cell_type": "code",
            "execution_count": null,
            "metadata": {},
            "outputs": [],
            "source": [
                "# install dependencies and select solver\n",
                "%pip install -q amplpy matplotlib\n",
                "\n",
                "SOLVER = \"cbc\"\n",
                "\n",
                "from amplpy import AMPL, ampl_notebook\n",
                "\n",
                "ampl = ampl_notebook(\n",
                "    modules=[\"cbc\"],  # modules to install\n",
                "    license_uuid=\"default\",  # license to use\n",
                ")  # instantiate AMPL object and register magics"
            ]
        },
        {
            "cell_type": "code",
            "execution_count": 2,
            "metadata": {},
            "outputs": [],
            "source": [
                "import matplotlib.pyplot as plt\n",
                "import numpy as np\n",
                "import pandas as pd"
            ]
        },
        {
            "cell_type": "markdown",
            "metadata": {
                "colab_type": "text",
                "id": "RrSgeD0ewWpE"
            },
            "source": [
                "## Modeling with disjunctive constraints\n",
                "\n",
                "The model is comprised of two sets of the binary variables indexed 1 to $T$. Binary variables $x_t$ correspond to the operating mode of the process unit, with $x_t=1$ indicating the unit is operating on day $t$ and able to earn a profit $c_t$. Binary variable $y_t=1$ indicates the first day of a maintenance period during which the unit is not operating and earning $0$ profit."
            ]
        },
        {
            "cell_type": "markdown",
            "metadata": {
                "colab_type": "text",
                "id": "RrSgeD0ewWpE"
            },
            "source": [
                "### Objective\n",
                "\n",
                "The planning objective is to maximize profit realized during the days the plant is operational. \n",
                "\n",
                "$$\n",
                "\\begin{align*}\n",
                "\\mbox{Profit} & = \\max_{x, y} \\sum_{t=1}^T c_t x_t\n",
                "\\end{align*}\n",
                "$$\n",
                "\n",
                "subject to completing $P$ maintenance periods.  "
            ]
        },
        {
            "cell_type": "markdown",
            "metadata": {
                "colab_type": "text",
                "id": "RrSgeD0ewWpE"
            },
            "source": [
                "### Constraints\n",
                "\n",
                "**Number of maintenance periods is equal to P.**\n",
                "\n",
                "Completing $P$ maintenance periods requires a total of $P$ starts.\n",
                "\n",
                "$$\n",
                "\\begin{align*}\n",
                "\\sum_{t=1}^T y_t &  = P \\\\\n",
                "\\end{align*}\n",
                "$$\n",
                "\n",
                "**Maintenance periods do not overlap.**\n",
                "\n",
                "No more than one maintenance period can start in any consecutive set of M days.\n",
                "\n",
                "$$\n",
                "\\begin{align*}\n",
                "\\sum_{s=0}^{M-1}y_{t+s} & \\leq 1 \\qquad \\forall t = 1, 2, \\ldots, T-M+1\n",
                "\\end{align*}\n",
                "$$\n",
                "\n",
                "This last requirement could be modified if some period of time should occur between maintenance periods.\n",
                "\n",
                "**The unit must shut down for M days following a maintenance start.**\n",
                "\n",
                "The final requirement is a disjunctive constraint that says either $y_t = 0$ or the sum $\\sum_{s}^{M-1}x_{t+s} = 0$, but not both. Mathematically, this forms a set of constraints reading\n",
                "\n",
                "$$\n",
                "\\begin{align*}\n",
                "\\left(y_t = 0\\right) \\veebar \\left(\\sum_{s=0}^{M-1}x_{t+s} = 0\\right)\\qquad \\forall t = 1, 2, \\ldots, T-M+1\n",
                "\\end{align*}\n",
                "$$\n",
                "\n",
                "where $\\veebar$ denotes an exclusive or condition."
            ]
        },
        {
            "cell_type": "markdown",
            "metadata": {},
            "source": [
                "## AMPL solution"
            ]
        },
        {
            "cell_type": "markdown",
            "metadata": {
                "colab_type": "text",
                "id": "Bae6-dR_lYkm"
            },
            "source": [
                "### Parameter values"
            ]
        },
        {
            "cell_type": "code",
            "execution_count": 3,
            "metadata": {
                "colab": {},
                "colab_type": "code",
                "id": "L5TaVPbkEPZ0"
            },
            "outputs": [],
            "source": [
                "# problem parameters\n",
                "T = 90  # planning period from 1..T\n",
                "M = 3  # length of maintenance period\n",
                "P = 4  # number of maintenance periods\n",
                "\n",
                "# daily profits\n",
                "np.random.seed(0)\n",
                "c = {k: np.random.uniform() for k in range(1, T + 1)}"
            ]
        },
        {
            "cell_type": "markdown",
            "metadata": {
                "colab_type": "text",
                "id": "e28FWA1GlNyQ"
            },
            "source": [
                "### AMPL model\n",
                "\n",
                "The disjunctive constraints can be represented directly in AMPL using [Logic, Nonlinear & Constraint Programming Extensions](https://ampl.com/products/ampl/logic-and-constraint-programming-extensions/). The extension transforms the disjunctive constraints to a mixed integer linear optimization (MILO) problem using convex hull and cutting plane methods."
            ]
        },
        {
            "cell_type": "code",
            "execution_count": 4,
            "metadata": {},
            "outputs": [
                {
                    "name": "stdout",
                    "output_type": "stream",
                    "text": [
                        "Writing maintenance_planning.mod\n"
                    ]
                }
            ],
            "source": [
                "%%writefile maintenance_planning.mod\n",
                "\n",
                "param T; # number of planning periods\n",
                "param P; # number of maintenance periods\n",
                "param M; # number of days for each maintenance period\n",
                "\n",
                "set PH := 1..T; # set of planning horizon time periods\n",
                "\n",
                "set Y := 1..(T - M + 1);\n",
                "set S := 0..(M - 1);\n",
                "\n",
                "param c{PH}; # profit\n",
                "\n",
                "var x{PH} binary;\n",
                "var y{PH} binary;\n",
                "\n",
                "maximize profit: sum{t in PH} c[t] * x[t];\n",
                "\n",
                "s.t. required_maintenance: P == sum{t in Y} y[t];\n",
                "s.t. no_overlap {t in Y}: sum{s in S} y[t + s] <= 1;\n",
                "s.t. required_shutdown {t in Y}:\n",
                "    (y[t] == 0 or sum{s in S} x[t + s] == 0)\n",
                "    and not\n",
                "    (y[t] == 0 and sum{s in S} x[t + s] == 0);"
            ]
        },
        {
            "cell_type": "code",
            "execution_count": 5,
            "metadata": {
                "colab": {
                    "base_uri": "https://localhost:8080/",
                    "height": 1014
                },
                "colab_type": "code",
                "executionInfo": {
                    "elapsed": 1543,
                    "status": "ok",
                    "timestamp": 1557947668999,
                    "user": {
                        "displayName": "Jeffrey Kantor",
                        "photoUrl": "https://lh5.googleusercontent.com/-8zK5aAW5RMQ/AAAAAAAAAAI/AAAAAAAAKB0/kssUQyz8DTQ/s64/photo.jpg",
                        "userId": "09038942003589296665"
                    },
                    "user_tz": 240
                },
                "id": "cUARDFzP9fla",
                "outputId": "ba0f2867-3a14-4314-c09f-bb4a76ebd1a7",
                "scrolled": false
            },
            "outputs": [
                {
                    "name": "stdout",
                    "output_type": "stream",
                    "text": [
                        "cbc 2.10.7: \b\b\b\b\b\b\b\b\b\b\b\bcbc 2.10.7: optimal solution; objective 41.92584964\n",
                        "21 simplex iterations\n",
                        "21 barrier iterations\n"
                    ]
                },
                {
                    "data": {
                        "image/png": "iVBORw0KGgoAAAANSUhEUgAAA2wAAAGGCAYAAAAHLAKcAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjcuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/bCgiHAAAACXBIWXMAAA9hAAAPYQGoP6dpAABDSklEQVR4nO3deXwV9b3/8ffJdpKQDUIWAgkEpLIpRYKIQKWFirkUb0TBKGgEqrckICEtFuoFQhEQuFoLWpZ6xVrBWmhBbYuy1KaXS8AIIlL2GiElLAEJiRAC5nx/f/jjXA5ZTwg5c5LX8/HI48F85zszn5l8Sc47s9mMMUYAAAAAAMvx8XQBAAAAAICqEdgAAAAAwKIIbAAAAABgUQQ2AAAAALAoAhsAAAAAWBSBDQAAAAAsisAGAAAAABZFYAMAAAAAiyKwAQAAAIBFEdgAAAAAwKIIbAAAAABgUQQ2AEC9ZGdny2azub3c66+/LpvNpi+++KLGNm+Sl5enu+++Wy1atJDNZtPu3bu9fp8AANZAYAMA4AZcuXJFI0eO1Jdffqlf/OIX+u1vf6v27dtX2Xfbtm3Kzs5WcXFx4xYJAPBaBDYAQKN67LHHVFZWVm2o8Tb//Oc/dfToUf3kJz/RU089pTFjxqhly5ZV7ue2bds0e/ZsAhsAoM78PF0AAKB58fX1la+vr6fLqNaFCxfUokWLOvc/ffq0JCkiIsKl3er7CQDwDpxhAwDUauvWrerTp48CAwPVqVMnLV++vFKfo0ePKj09XbfeequCgoIUGRmpkSNHVrqHq7Z7uz788EPZbDatW7eu0rzVq1fLZrMpNze3xnqv3l934MABjRo1SmFhYYqMjNTkyZN16dKlSv327dunRx99VC1bttSAAQOc8z/55BMlJycrLCxMISEhGjx4sLZv3+6c/8QTT+iee+6RJI0cOVI2m02DBg2qcj+zs7M1depUSVJiYqJsNlud7nHbuXOnHnroIUVHRyswMFBdunTRz3/+8xqXAQA0HZxhAwDU6LPPPtO9996rqKgoZWdn6+uvv9asWbMUExPj0i8vL0/btm1Tamqq2rVrpy+++EJLly7VoEGDtG/fPgUHB9dpe4MGDVJ8fLxWrVqlBx54wGXeqlWr1KlTJ/Xr169O6xo1apQ6dOig+fPna/v27Vq8eLHOnTunN954w6XfyJEj1blzZ82bN0/GGEnSP/7xDw0cOFBhYWF65pln5O/vr+XLl2vQoEHKyclR37599R//8R9q27at5s2bp6efflp9+vSpdFyuGjFihA4dOqS33npLv/jFL9S6dWtJUlRUVLX1r1u3TqmpqerYsaOmTp2qkJAQ53EGADQTBgCAGqSkpJjAwEBz9OhRZ9u+ffuMr6+vufbXyMWLFystm5ubaySZN954w9m2cuVKI8nk5+dX2zZ9+nRjt9tNcXGxs8/p06eNn5+fmTVrVq01z5o1y0gy999/v0t7enq6kWQ+/fRTl36PPPJIlfsdEBBg/vnPfzrbCgsLTWhoqPnOd77jbPvwww+NJLNmzRqX5avaz0WLFlVqq87BgwdNcHCweeCBB0xZWZnLvEuXLtW6PACgaeCSSABAtSoqKvTBBx8oJSVFCQkJzvauXbtq6NChLn2DgoKc/75y5YrOnj2rW265RREREdq1a5db23388cdVXl6utWvXOtvefvttff311xozZkyd15ORkeEyPWnSJEnSX/7yF5f2H/3oRy7TFRUV2rhxo1JSUtSxY0dne5s2bfToo49q69atKikpqXMd9TFz5kwFBgbq9ddfV2BgoMs8u91+U7cNALAOAhsAoFpFRUUqKytT586dK8279dZbXabLyso0c+ZMxcfHy263q3Xr1oqKilJxcbHOnz/v1na7dOmiPn36aNWqVc62VatW6a677tItt9xS5/VcX3enTp3k4+NT6b6xxMREl+mioiJdvHix0j5K34RVh8OhgoKCOtfhrvLycr333nsaM2aMwsLCauxrjFFISIjz4ScAgKaFe9gAAA1i0qRJWrlypTIzM9WvXz+Fh4fLZrMpNTVVDofD7fU9/vjjmjx5sv71r3+pvLxc27dv18svv3xDNVb3ou9rzw5aweeff66LFy+qd+/etfbNz89XcHCwoqOjG6EyAEBjI7ABAKoVFRWloKAgHT58uNK8gwcPukyvXbtWaWlpeuGFF5xtly5dqvc7x1JTU5WVlaW33npLZWVl8vf318MPP+zWOg4fPuxy9uzIkSNyOBzq0KFDjctFRUUpODi40j5K0oEDB+Tj46P4+Hi3apGqD4zXKysrq1P//fv3q1evXqqoqFBISIi6du2qvLw8t+sCAFgXl0QCAKrl6+uroUOHav369Tp27Jizff/+/frggw8q9TX//wmLVy1ZskQVFRX12nbr1q2VnJysN998U6tWrdJ9993nfLJiXb3yyiuV6pGk5OTkGpfz9fXVvffeq3feecfl8slTp05p9erVGjBgQK2XKlbl6vvdaguxVy/73Lx5c6V5V65ccf67a9eumjVrliZMmKCvvvqKsAYATRBn2AAANZo9e7bef/99DRw4UOnp6fr666+1ZMkSde/eXXv27HH2+8EPfqDf/va3Cg8PV7du3ZSbm6vNmzcrMjKy3tt+/PHH9dBDD0mS5syZ4/by+fn5uv/++3XfffcpNzdXb775ph599FH17Nmz1mWfe+45bdq0SQMGDFB6err8/Py0fPlylZeXa+HChW7XIsl5ieOzzz6r1NRU+fv7a/jw4ZVe1B0WFqYnnnhCr7/+usrLyzVo0CCVlpbqww8/1LBhw1weprJnzx5973vfq1c9AADrI7ABAGp0++2364MPPlBWVpZmzpypdu3aafbs2Tpx4oRLYPvlL38pX19frVq1SpcuXVL//v21efPmSk+TdMfw4cPVsmVLORwO3X///W4v//bbb2vmzJmaNm2a/Pz8NHHiRC1atKhOy3bv3l3/8z//o+nTp2v+/PlyOBzq27ev3nzzTfXt29ftWiSpT58+mjNnjpYtW6b3339fDodD+fn5lQKb9M3ZwbZt2+r3v/+91q9fr1atWunuu+/Wfffd59Jvz549yszMrFc9AADrs5nrr18BAMAivv76a8XFxWn48OH67//+7zovl52drdmzZ6uoqMjtyyi9yeXLlxUSEqJz585VGfoAAN6Pe9gAAJa1fv16FRUV6fHHH/d0KZZUWloq6ZvgBgBomrgkEgBgOTt27NCePXs0Z84c9erVS/fcc4+nS7KkyMhIPfLII0pISFD37t21fft2T5cEAGhgnGEDAFjO0qVLNWHCBEVHR+uNN97wdDmW9pvf/EalpaWENQBooriHDQAAAAAsijNsAAAAAGBRbge2v//97xo+fLji4uJks9m0fv36m1AWAAAAAMDth45cuHBBPXv21Lhx4zRixAi3N+hwOFRYWKjQ0FDZbDa3lwcAAAAAb2eMUWlpqeLi4uTjU/15NLcDW3JyspKTk+tdWGFhoeLj4+u9PAAAAAA0FQUFBWrXrl218xv9sf6hoaGSviksLCyssTcPAAAAAB5XUlKi+Ph4Zz6qzk0PbOXl5SovL3dOX33JZ1hYGIENAAAAQLNW221iNz2wzZ8/X7Nnz77ZmwEAAG7qMO3PNc7/4vlhjVQJAKA6Nz2wTZ8+XVlZWc7pq6f+gOakpg9FfCACAOD/NOQfEvj9i6bgpgc2u90uu91+szfTLPBDBwAAAGhe3A5sX331lY4cOeKczs/P1+7du9WqVSslJCQ0aHEAmgYuu7Im/ggEAID1uR3YPv74Y333u991Tl+93DEtLU2vv/56gxWGpoMP6wAAAGgObsYfQ90ObIMGDZIxpl4bAwAAAADUXfWv1AYAAAAAeFSjvzgbVeNeEgAAAADXI7ABN4iwDQAAgJuFSyIBAAAAwKI4w4ZmiSdXAgDQsLjiBLg5CGwAAADNGEELsDYCG2AR/MIEAADA9QhszRThAAAAALA+AhtQDW++z41AXjtv/v4CAIDmg6dEAgAAAIBFcYYNAAAAqAFXZcCTOMMGAAAAABblNWfY+MsGAAAAgObGawIbAKDp44E5ANC0cRLGfQQ2AEC1mvovVgIiAMDqCGwAgCanqQdNeC/+SADAXQS2G8AHAjR1fLAA0Nj43QoArnhKJAAAAABYFGfYAMAi6npmoS5nPjk7CgBA09DkAhsfUgAAAFBXXIYLq2tygQ0AAAAAGpqnTgwR2HBDOKMJ8P8AjY8xB6Cp4+fc/7FEYGvsU9EMAAAAAADewBKBDdZEsAUAAAA8i8AGwDL4IwEAAGhITeGhMryHDQAAAAAsisAGAAAAABbFJZEAAKBJ4jLrpo3vL5oLAhvgRZrCddjehg8EAADAkwhsjYAPfADAz8LGxvEGgKaBwFYNzmQAAAAA8DQCGwCvwlkDAAAaH79/PYfA1sRwZhAAAABoOghsAAAATRB/xAWaBgIbAOCG8KGwYXHZEQDgWrw4GwAAAAAsisAGAAAAABbFJZEAADQCLh0FYFVcim1tBDYAAADAy/BHoOaDwAYAAAA0Es5mwV0ENgBNDn91BACg6Wsu4ZfABgAAvEpD/lGmuXzgA+C9eEokAAAAAFgUZ9gAAKgBZ2DqhuMEADcHgQ0AAMDLEJCB5oNLIgEAAADAojjDBgAAGgVnhQDAfQQ2AACAGvCqEACeRGADcEP4IAMADYszkQCuRWADmik+EAANhz9cAGhIDfk7mt/33o/ABgAAasQHvtoR2gHcLAQ2eBU+NAAA0Pga6vcvwRZWZPVxWa/A9sorr2jRokU6efKkevbsqSVLlujOO+9s6NrQjFj9PwoAAADgCW4HtrfffltZWVlatmyZ+vbtq5deeklDhw7VwYMHFR0dfTNqBADAiT/wAIA18fP55nA7sL344ot68sknNXbsWEnSsmXL9Oc//1mvvfaapk2b1uAFAgAAAN6AWzdwM7gV2C5fvqydO3dq+vTpzjYfHx8NGTJEubm5DV4cUB/8sIQVMS4BAEB9uBXYzpw5o4qKCsXExLi0x8TE6MCBA1UuU15ervLycuf0+fPnJUklJSXONkf5xRq3W1JSUqc+ta2rofpYsSZ36rZiTVaq24o18f1t2jXx/aWmxq6JMde0a+L727Rr4vvr/TVdP22MqXH7NlNbj2sUFhaqbdu22rZtm/r16+dsf+aZZ5STk6MdO3ZUWiY7O1uzZ8+u6yYAAAAAoNkoKChQu3btqp3v1hm21q1by9fXV6dOnXJpP3XqlGJjY6tcZvr06crKynJOOxwOffnll4qMjJTNZnNn84CllZSUKD4+XgUFBQoLC/N0OUCDY4yjqWOMoyljfFuPMUalpaWKi4ursZ9bgS0gIEC9e/fWli1blJKSIumbALZlyxZNnDixymXsdrvsdrtLW0REhDubBbxKWFgYPwjRpDHG0dQxxtGUMb6tJTw8vNY+bj8lMisrS2lpaUpKStKdd96pl156SRcuXHA+NRIAAAAA0DDcDmwPP/ywioqKNHPmTJ08eVLf/va39f7771d6EAkAAAAA4Ma4HdgkaeLEidVeAgk0V3a7XbNmzap0CTDQVDDG0dQxxtGUMb69l1tPiQQAAAAANB4fTxcAAAAAAKgagQ0AAAAALIrABgAAAAAWRWAD3DR//nz16dNHoaGhio6OVkpKig4ePOjS59KlS8rIyFBkZKRCQkL04IMPVnrhPOANnn/+edlsNmVmZjrbGN/wdsePH9eYMWMUGRmpoKAg3Xbbbfr444+d840xmjlzptq0aaOgoCANGTJEhw8f9mDFQN1UVFRoxowZSkxMVFBQkDp16qQ5c+bo2kdWML69D4ENcFNOTo4yMjK0fft2bdq0SVeuXNG9996rCxcuOPtMmTJF7733ntasWaOcnBwVFhZqxIgRHqwacF9eXp6WL1+u22+/3aWd8Q1vdu7cOfXv31/+/v7asGGD9u3bpxdeeEEtW7Z09lm4cKEWL16sZcuWaceOHWrRooWGDh2qS5cuebByoHYLFizQ0qVL9fLLL2v//v1asGCBFi5cqCVLljj7ML69kAFwQ06fPm0kmZycHGOMMcXFxcbf39+sWbPG2Wf//v1GksnNzfVUmYBbSktLTefOnc2mTZvMPffcYyZPnmyMYXzD+/30pz81AwYMqHa+w+EwsbGxZtGiRc624uJiY7fbzVtvvdUYJQL1NmzYMDNu3DiXthEjRpjRo0cbYxjf3oozbMANOn/+vCSpVatWkqSdO3fqypUrGjJkiLNPly5dlJCQoNzcXI/UCLgrIyNDw4YNcxnHEuMb3u/dd99VUlKSRo4cqejoaPXq1Uu//vWvnfPz8/N18uRJlzEeHh6uvn37MsZheXfffbe2bNmiQ4cOSZI+/fRTbd26VcnJyZIY396qXi/OBvANh8OhzMxM9e/fXz169JAknTx5UgEBAYqIiHDpGxMTo5MnT3qgSsA9v/vd77Rr1y7l5eVVmsf4hrf7/PPPtXTpUmVlZelnP/uZ8vLy9PTTTysgIEBpaWnOcRwTE+OyHGMc3mDatGkqKSlRly5d5Ovrq4qKCs2dO1ejR4+WJMa3lyKwATcgIyNDe/fu1datWz1dCtAgCgoKNHnyZG3atEmBgYGeLgdocA6HQ0lJSZo3b54kqVevXtq7d6+WLVumtLQ0D1cH3Jjf//73WrVqlVavXq3u3btr9+7dyszMVFxcHOPbi3FJJFBPEydO1J/+9Cd9+OGHateunbM9NjZWly9fVnFxsUv/U6dOKTY2tpGrBNyzc+dOnT59WnfccYf8/Pzk5+ennJwcLV68WH5+foqJiWF8w6u1adNG3bp1c2nr2rWrjh07JknOcXz9k08Z4/AGU6dO1bRp05SamqrbbrtNjz32mKZMmaL58+dLYnx7KwIb4CZjjCZOnKh169bpr3/9qxITE13m9+7dW/7+/tqyZYuz7eDBgzp27Jj69evX2OUCbhk8eLA+++wz7d692/mVlJSk0aNHO//N+IY369+/f6VXsRw6dEjt27eXJCUmJio2NtZljJeUlGjHjh2McVjexYsX5ePj+vHe19dXDodDEuPbW3FJJOCmjIwMrV69Wu+8845CQ0Od13yHh4crKChI4eHhGj9+vLKystSqVSuFhYVp0qRJ6tevn+666y4PVw/ULDQ01Hk/5lUtWrRQZGSks53xDW82ZcoU3X333Zo3b55GjRqljz76SCtWrNCKFSskyfneweeee06dO3dWYmKiZsyYobi4OKWkpHi2eKAWw4cP19y5c5WQkKDu3bvrk08+0Ysvvqhx48ZJYnx7LU8/phLwNpKq/Fq5cqWzT1lZmUlPTzctW7Y0wcHB5oEHHjAnTpzwXNHADbj2sf7GML7h/d577z3To0cPY7fbTZcuXcyKFStc5jscDjNjxgwTExNj7Ha7GTx4sDl48KCHqgXqrqSkxEyePNkkJCSYwMBA07FjR/Pss8+a8vJyZx/Gt/exGXPNq88BAAAAAJbBPWwAAAAAYFEENgAAAACwKAIbAAAAAFgUgQ0AAAAALIrABgAAAAAWRWADAAAAAIsisAEAAACARRHYAAAAAMCiCGwAAAAAYFEENgAAAACwKAIbAAAAAFgUgQ0AAAAALIrABgAAAAAWRWADAAAAAIsisAEAAACARRHYAAAAAMCiCGwAAAAAYFEENgAAAACwKAIbAAAAAFgUgQ0ALOz111+XzWbTF1984elSLKUpHJfs7GzZbDadOXOm0bZVH03hWAOANyOwAYAX2bZtm7Kzs1VcXOzpUhpFc9tfAACuR2ADAAt77LHHVFZWpvbt20v6JsDMnj272QSY6vb3+uMCAEBT5efpAgAA1fP19ZWvr6+ny2gwFy5cUIsWLW54PU3tuAAAUB3OsAFAA3niiSfUoUOHSu3X3z90dfrIkSN64oknFBERofDwcI0dO1YXL150Wfba+4eys7M1depUSVJiYqJsNlud7i365JNPlJycrLCwMIWEhGjw4MHavn17lTUeOHBAo0aNUlhYmCIjIzV58mRdunSp0jqPHz+ucePGKSYmRna7Xd27d9drr71W5Tr37dunRx99VC1bttSAAQMkSUePHlV6erpuvfVWBQUFKTIyUiNHjnTZl5r2t6r7qtw5rpL0t7/9TUlJSQoMDFSnTp20fPnyOt3rVVpaqszMTHXo0EF2u13R0dH6/ve/r127dlU6RuPHj1dcXJzsdrsSExM1YcIEXb582aVfcXFxneqtyzGXpK1bt6pPnz4u+1WVuo7X6tS1nqqWCwwM1Lhx41zaN2/eLH9/f02ZMqXWdQBAc8IZNgDwkFGjRikxMVHz58/Xrl279Oqrryo6OloLFiyosv+IESN06NAhvfXWW/rFL36h1q1bS5KioqKq3cY//vEPDRw4UGFhYXrmmWfk7++v5cuXa9CgQcrJyVHfvn0r1dShQwfNnz9f27dv1+LFi3Xu3Dm98cYbzj6nTp3SXXfdJZvNpokTJyoqKkobNmzQ+PHjVVJSoszMTJd1jhw5Up07d9a8efNkjJEk5eXladu2bUpNTVW7du30xRdfaOnSpRo0aJD27dun4ODgeu1vXY/rJ598ovvuu09t2rTR7NmzVVFRoZ///Oe1rluSfvSjH2nt2rWaOHGiunXrprNnz2rr1q3av3+/7rjjDklSYWGh7rzzThUXF+upp55Sly5ddPz4ca1du1YXL15UQECAW/XW9Zh/9tlnuvfeexUVFaXs7Gx9/fXXmjVrlmJiYmrdL3e4Owau1bZtW/3whz/UihUrNGvWLLVv314HDhzQyJEjlZycrBdeeKFBawUAr2cAAA0iLS3NtG/fvlL7rFmzzLU/bq9Ojxs3zqXfAw88YCIjI13aVq5caSSZ/Px8Y4wxixYtcpmuTUpKigkICDD//Oc/nW2FhYUmNDTUfOc736lU0/333++yfHp6upFkPv30U2fb+PHjTZs2bcyZM2dc+qampprw8HBz8eJFl3U+8sgjleq62udaubm5RpJ54403nG3V7e/1x+Xa7dXluA4fPtwEBweb48ePO9sOHz5s/Pz8TG2/GsPDw01GRkaNfR5//HHj4+Nj8vLyKs1zOBxu11vXY56SkmICAwPN0aNHnX327dtnfH19K+1XXcdrVce6rvVU51//+pex2+1mwoQJ5syZM6ZTp07m29/+tvnqq69qXA4AmiMuiQQAD/nRj37kMj1w4ECdPXtWJSUlDbL+iooKbdy4USkpKerYsaOzvU2bNnr00Ue1devWStvKyMhwmZ40aZIk6S9/+YskyRijP/zhDxo+fLiMMTpz5ozza+jQoTp//nylSwOv309JCgoKcv77ypUrOnv2rG655RZFRERUWt5dtR3XiooKbd68WSkpKYqLi3P2u+WWW5ScnFzr+iMiIrRjxw4VFhZWOd/hcGj9+vUaPny4kpKSKs2//nLD2uqt6zGvqKjQBx98oJSUFCUkJDjX17VrVw0dOrTW/aqr+oyB67Vt21ZPPvmkXnvtNQ0bNkxlZWX605/+1CD3NwJAU0NgAwAPufZDtSS1bNlSknTu3LkGWX9RUZEuXryoW2+9tdK8rl27yuFwqKCgwKW9c+fOLtOdOnWSj4+P816xoqIiFRcXa8WKFYqKinL5Gjt2rCTp9OnTLutITEystP2ysjLNnDlT8fHxstvtat26taKiolRcXKzz58/fyG7XelxPnz6tsrIy3XLLLZWWrartegsXLtTevXsVHx+vO++8U9nZ2fr888+d84uKilRSUqIePXo0SL11PeZFRUUqKyur9D2UVOUYqK/6jIGq/OQnP1F5ebn27Nmjd999V23btnWZb4xRSEhIndYFAE0Z97ABQAOp7kENFRUVVbZX95RD8//v87KC6/fJ4XBIksaMGaO0tLQql7n99ttdpq89m3bVpEmTtHLlSmVmZqpfv34KDw+XzWZTamqqcxv1dbOP66hRozRw4ECtW7dOGzdu1KJFi7RgwQL98Y9/rNMZuuvVVm9dj7m7x83d8XpVfcZAVebOnStJ+vrrr9WqVatK8/Pz8xUcHKzo6Oha1wUATRmBDQAaSMuWLat8P9rRo0cbbBt1eXrfVVFRUQoODtbBgwcrzTtw4IB8fHwUHx/v0n748GGXM2JHjhyRw+FwPk0wKipKoaGhqqio0JAhQ+q3E5LWrl2rtLQ0lwdMXLp0qdLxc2d/6yo6OlqBgYE6cuRIpXlVtVWlTZs2Sk9PV3p6uk6fPq077rhDc+fOVXJysqKiohQWFqa9e/c2SL11PeYVFRUKCgrS4cOHK82ragzUd7w2xBhYtGiRXn31Vb388suaOnWq5s6dq1dffdU5f//+/erVq5cqKioUEhKirl27Ki8vr17bAgBvxyWRANBAOnXqpPPnz2vPnj3OthMnTmjdunUNto2r9/jU5cXZvr6+uvfee/XOO++4PP7+1KlTWr16tQYMGKCwsDCXZV555RWX6SVLlkiS88yRr6+vHnzwQf3hD3+oMpAUFRXVaT98fX0rnfFasmRJpbM77uxvXfn6+mrIkCFav369y31oR44c0YYNG2pctqKiotIlm9HR0YqLi1N5ebkkycfHRykpKXrvvff08ccfV1qHu2f66nrMfX19NXToUK1fv17Hjh1zzt+/f78++OCDSsvVd7ze6BhYv369pk2bpjlz5igjI0NPPfWU3njjDeXn5zv7dO3aVbNmzdKECRP01VdfEdYANGucYQOABpKamqqf/vSneuCBB/T000/r4sWLWrp0qb71rW/d8IM0rurdu7ck6dlnn1Vqaqr8/f01fPjwah/W8Nxzz2nTpk0aMGCA0tPT5efnp+XLl6u8vFwLFy6s1D8/P1/333+/7rvvPuXm5urNN9/Uo48+qp49ezr7PP/88/rwww/Vt29fPfnkk+rWrZu+/PJL7dq1S5s3b9aXX35Z63784Ac/0G9/+1uFh4erW7duys3N1ebNmxUZGVmn/b1R2dnZ2rhxo/r3768JEyaooqJCL7/8snr06KHdu3dXu1xpaanatWunhx56SD179lRISIg2b96svLw8l7OF8+bN08aNG3XPPffoqaeeUteuXXXixAmtWbNGW7duVUREhFv11vWYz549W++//74GDhyo9PR0ff3111qyZIm6d+/uEsykGxuv9R0DO3fu1OjRozV69Gg9++yzkqRnnnlGy5Ytq3SWbc+ePfre977n1nECgCbJMw+nBICmaePGjaZHjx4mICDA3HrrrebNN9+s9rH+RUVFLstW9fj0qtrmzJlj2rZta3x8fOr0iP9du3aZoUOHmpCQEBMcHGy++93vmm3btrn0uVrTvn37zEMPPWRCQ0NNy5YtzcSJE01ZWVmldZ46dcpkZGSY+Ph44+/vb2JjY83gwYPNihUrat1PY4w5d+6cGTt2rGndurUJCQkxQ4cONQcOHDDt27c3aWlpLn2r2t+aHutfl+NqjDFbtmwxvXr1MgEBAaZTp07m1VdfNT/+8Y9NYGBgtceyvLzcTJ061fTs2dOEhoaaFi1amJ49e5pf/epXlfoePXrUPP744yYqKsrY7XbTsWNHk5GRYcrLy+tVb12OuTHG5OTkmN69e5uAgADTsWNHs2zZskpj8Kq6jNcbreeqgoIC06ZNG9O/f39z6dIll3kTJkww/v7+5vPPP3e2devWzWzfvr3KdQFAc2IzxkJ3twMAPCI7O1uzZ89WUVGR8wXVzVFKSor+8Y9/VHkfGBrP5cuXFRISonPnzvGofwDNHvewAQCapbKyMpfpw4cP6y9/+YsGDRrkmYLgVFpaKumb4AYAzR33sAEAmqWOHTvqiSeeUMeOHXX06FEtXbpUAQEBeuaZZzxdWrMXGRmpRx55RAkJCerevbu2b9/u6ZIAwGMIbACAZum+++7TW2+9pZMnT8put6tfv36aN29elS+eRuP7zW9+o9/85jeeLgMAPI572AAAAADAoriHDQAAAAAsisAGAAAAABbldmD7+9//ruHDhysuLk42m03r16+/CWUBAAAAANx+6MiFCxfUs2dPjRs3TiNGjHB7gw6HQ4WFhQoNDZXNZnN7eQAAAADwdsYYlZaWKi4uTj4+1Z9HczuwJScnKzk5ud6FFRYWKj4+vt7LAwAAAEBTUVBQoHbt2lU7v9Ef6x8aGirpm8LCwsIae/MAAAAA4HElJSWKj4935qPq3PTAVl5ervLycud0aWmpJCksLIzABgAAAKBZq+02sZse2ObPn6/Zs2fX2KfDtD/XOP+L54fVqU9t62qoPlasyZ26rVpTY+L723z/r1ixJr6/1NTYNfHz2TP7Jlnz+2u1mqxUt7s1NSYrHafm8v1tqJrcddMf6z99+nSdP3/e+VVQUHCzNwkAAAAATcJNP8Nmt9tlt9tv9mYAAAAAoMlxO7B99dVXOnLkiHM6Pz9fu3fvVqtWrZSQkNCgxQEAAABAc+Z2YPv444/13e9+1zmdlZUlSUpLS9Prr7/eYIUBAAAAQHPndmAbNGiQjDE3oxYAAAAAwDVu+kNHAAAAAAD1Q2ADAAAAAIsisAEAAACARRHYAAAAAMCiCGwAAAAAYFEENgAAAACwKAIbAAAAAFgUgQ0AAAAALIrABgAAAAAWRWADAAAAAIsisAEAAACARRHYAAAAAMCiCGwAAAAAYFEENgAAAACwKAIbAAAAAFgUgQ0AAAAALIrABgAAAAAWRWADAAAAAIsisAEAAACARRHYAAAAAMCiCGwAAAAAYFEENgAAAACwKAIbAAAAAFgUgQ0AAAAALIrABgAAAAAWRWADAAAAAIsisAEAAACARRHYAAAAAMCiCGwAAAAAYFEENgAAAACwKAIbAAAAAFgUgQ0AAAAALIrABgAAAAAWRWADAAAAAIsisAEAAACARRHYAAAAAMCiCGwAAAAAYFEENgAAAACwKAIbAAAAAFgUgQ0AAAAALIrABgAAAAAWRWADAAAAAIsisAEAAACARRHYAAAAAMCiCGwAAAAAYFEENgAAAACwKAIbAAAAAFgUgQ0AAAAALIrABgAAAAAWRWADAAAAAIsisAEAAACARdUrsL3yyivq0KGDAgMD1bdvX3300UcNXRcAAAAANHtuB7a3335bWVlZmjVrlnbt2qWePXtq6NChOn369M2oDwAAAACaLbcD24svvqgnn3xSY8eOVbdu3bRs2TIFBwfrtddeuxn1AQAAAECz5VZgu3z5snbu3KkhQ4b83wp8fDRkyBDl5uY2eHEAAAAA0Jz5udP5zJkzqqioUExMjEt7TEyMDhw4UOUy5eXlKi8vd06fP39eklRSUuJsc5RfrHG7JSUldepT27oaqo8Va3KnbqvW1Jj4/jbf/ytWrInvLzU1dk38fPbMvknW/P5arSYr1e1uTY3JSsepuXx/G6qm66eNMTVu32Zq63GNwsJCtW3bVtu2bVO/fv2c7c8884xycnK0Y8eOSstkZ2dr9uzZdd0EAAAAADQbBQUFateuXbXz3TrD1rp1a/n6+urUqVMu7adOnVJsbGyVy0yfPl1ZWVnOaYfDoS+//FKRkZGy2WzubB6wtJKSEsXHx6ugoEBhYWGeLgdocIxxNHWMcTRljG/rMcaotLRUcXFxNfZzK7AFBASod+/e2rJli1JSUiR9E8C2bNmiiRMnVrmM3W6X3W53aYuIiHBns4BXCQsL4wchmjTGOJo6xjiaMsa3tYSHh9fax63AJklZWVlKS0tTUlKS7rzzTr300ku6cOGCxo4dW68iAQAAAABVczuwPfzwwyoqKtLMmTN18uRJffvb39b7779f6UEkAAAAAIAb43Zgk6SJEydWewkk0FzZ7XbNmjWr0iXAQFPBGEdTxxhHU8b49l5uPSUSAAAAANB43HpxNgAAAACg8RDYAAAAAMCiCGwAAAAAYFEENgAAAACwKAIb4Kb58+erT58+Cg0NVXR0tFJSUnTw4EGXPpcuXVJGRoYiIyMVEhKiBx98UKdOnfJQxUD9Pf/887LZbMrMzHS2Mb7h7Y4fP64xY8YoMjJSQUFBuu222/Txxx875xtjNHPmTLVp00ZBQUEaMmSIDh8+7MGKgbqpqKjQjBkzlJiYqKCgIHXq1Elz5szRtc8YZHx7HwIb4KacnBxlZGRo+/bt2rRpk65cuaJ7771XFy5ccPaZMmWK3nvvPa1Zs0Y5OTkqLCzUiBEjPFg14L68vDwtX75ct99+u0s74xve7Ny5c+rfv7/8/f21YcMG7du3Ty+88IJatmzp7LNw4UItXrxYy5Yt044dO9SiRQsNHTpUly5d8mDlQO0WLFigpUuX6uWXX9b+/fu1YMECLVy4UEuWLHH2YXx7IQPghpw+fdpIMjk5OcYYY4qLi42/v79Zs2aNs8/+/fuNJJObm+upMgG3lJaWms6dO5tNmzaZe+65x0yePNkYw/iG9/vpT39qBgwYUO18h8NhYmNjzaJFi5xtxcXFxm63m7feeqsxSgTqbdiwYWbcuHEubSNGjDCjR482xjC+vRVn2IAbdP78eUlSq1atJEk7d+7UlStXNGTIEGefLl26KCEhQbm5uR6pEXBXRkaGhg0b5jKOJcY3vN+7776rpKQkjRw5UtHR0erVq5d+/etfO+fn5+fr5MmTLmM8PDxcffv2ZYzD8u6++25t2bJFhw4dkiR9+umn2rp1q5KTkyUxvr2Vn6cLALyZw+FQZmam+vfvrx49ekiSTp48qYCAAEVERLj0jYmJ0cmTJz1QJeCe3/3ud9q1a5fy8vIqzWN8w9t9/vnnWrp0qbKysvSzn/1MeXl5evrppxUQEKC0tDTnOI6JiXFZjjEObzBt2jSVlJSoS5cu8vX1VUVFhebOnavRo0dLEuPbSxHYgBuQkZGhvXv3auvWrZ4uBWgQBQUFmjx5sjZt2qTAwEBPlwM0OIfDoaSkJM2bN0+S1KtXL+3du1fLli1TWlqah6sDbszvf/97rVq1SqtXr1b37t21e/duZWZmKi4ujvHtxbgkEqiniRMn6k9/+pM+/PBDtWvXztkeGxury5cvq7i42KX/qVOnFBsb28hVAu7ZuXOnTp8+rTvuuEN+fn7y8/NTTk6OFi9eLD8/P8XExDC+4dXatGmjbt26ubR17dpVx44dkyTnOL7+yaeMcXiDqVOnatq0aUpNTdVtt92mxx57TFOmTNH8+fMlMb69FYENcJMxRhMnTtS6dev017/+VYmJiS7ze/fuLX9/f23ZssXZdvDgQR07dkz9+vVr7HIBtwwePFifffaZdu/e7fxKSkrS6NGjnf9mfMOb9e/fv9KrWA4dOqT27dtLkhITExUbG+syxktKSrRjxw7GOCzv4sWL8vFx/Xjv6+srh8MhifHtrbgkEnBTRkaGVq9erXfeeUehoaHOa77Dw8MVFBSk8PBwjR8/XllZWWrVqpXCwsI0adIk9evXT3fddZeHqwdqFhoa6rwf86oWLVooMjLS2c74hjebMmWK7r77bs2bN0+jRo3SRx99pBUrVmjFihWS5Hzv4HPPPafOnTsrMTFRM2bMUFxcnFJSUjxbPFCL4cOHa+7cuUpISFD37t31ySef6MUXX9S4ceMkMb69lqcfUwl4G0lVfq1cudLZp6yszKSnp5uWLVua4OBg88ADD5gTJ054rmjgBlz7WH9jGN/wfu+9957p0aOHsdvtpkuXLmbFihUu8x0Oh5kxY4aJiYkxdrvdDB482Bw8eNBD1QJ1V1JSYiZPnmwSEhJMYGCg6dixo3n22WdNeXm5sw/j2/vYjLnm1ecAAAAAAMvgHjYAAAAAsCgCGwAAAABYFIENAAAAACyKwAYAAAAAFkVgAwAAAACLIrABAAAAgEUR2AAAAADAoghsAAAAAGBRBDYAAAAAsCgCGwAAAABYFIENAAAAACyKwAYAAAAAFkVgAwAAAACLIrABAAAAgEUR2AAAAADAoghsAAAAAGBRBDYAAAAAsCgCGwAAAABYFIENAAAAACyKwAYAqNLrr78um82mL774wtOlAADQbBHYAMDD5s6dK5vNph49eni6lAazbds2ZWdnq7i42NOleK2bfQz5HgGAdyCwAYAH/etf/9K8efPUokULT5dSyWOPPaaysjK1b9/e7WW3bdum2bNnEwZuwM0+hnyPAMA7+Hm6AABozn7yk5/orrvuUkVFhc6cOePpclz4+vrK19fX02WggV24cMGSfyAAAFSNM2wA4CF///vftXbtWr300kt1XiY7O1s2m02HDh3SmDFjFB4erqioKM2YMUPGGBUUFOjf//3fFRYWptjYWL3wwguV1nH06FGlp6fr1ltvVVBQkCIjIzVy5MhK96pdfw/b1W0fOXJETzzxhCIiIhQeHq6xY8fq4sWLLjVOnTpVkpSYmCibzeaynuPHj2vcuHGKiYmR3W5X9+7d9dprr1W7r7Vtr677U9f1XXX8+HGNHz9ecXFxstvtSkxM1IQJE3T58mWXPnXZl6qUlpYqMzNTHTp0kN1uV3R0tL7//e9r165dtR5Dd/d53759evTRR9WyZUsNGDCgxvXXVFdN/u3f/k0dOnSo1G6M0R133KGBAwfW6bgAAFxxhg0APKCiokKTJk3SD3/4Q912221uL//www+ra9euev755/XnP/9Zzz33nFq1aqXly5fre9/7nhYsWKBVq1bpJz/5ifr06aPvfOc7zmXz8vK0bds2paamql27dvriiy+0dOlSDRo0SPv27VNwcHCN2x41apQSExM1f/587dq1S6+++qqio6O1YMECSdKIESN06NAhvfXWW/rFL36h1q1bS5KioqJ06tQp3XXXXbLZbJo4caKioqK0YcMGjR8/XiUlJcrMzHR7e+7uT23rk6TCwkLdeeedKi4u1lNPPaUuXbro+PHjWrt2rS5evKiAgIB67cu1fvSjH2nt2rWaOHGiunXrprNnz2rr1q3av39/jcewPvs8cuRIde7cWfPmzZMxRgMGDKh2/U899VS1dd1xxx3V7k+fPn20YcMGnTt3Ti1btnS2/+53v9Mnn3yirVu31ng8AADVMACARvfyyy+b8PBwc/r0aWOMMffcc4/p3r17rcvNmjXLSDJPPfWUs+3rr7827dq1MzabzTz//PPO9nPnzpmgoCCTlpbmso6LFy9WWm9ubq6RZN544w1n28qVK40kk5+f77LtcePGuSz7wAMPmMjISJe2RYsWuSx71fjx402bNm3MmTNnXNpTU1NNeHi4S2113V5d98ed+h9//HHj4+Nj8vLyKq3b4XC4vS9VCQ8PNxkZGdXOr+4YGuP+Pj/yyCN1Xn9tdVXn3XffNZLMli1bnG2XL182nTp1MsOHD3d7fQCAb3BJJAA0srNnz2rmzJmaMWOG84yJu374wx86/+3r66ukpCQZYzR+/Hhne0REhG699VZ9/vnnLssGBQU5/33lyhWdPXtWt9xyiyIiImq97E365szQtQYOHKizZ8+qpKSkxuWMMfrDH/6g4cOHyxijM2fOOL+GDh2q8+fPV7n92rbn7v7Utj6Hw6H169dr+PDhSkpKqrS8zWar975cKyIiQjt27FBhYWGN/apyo/t8M+rq06ePJLlsf8WKFcrPz9e8efPcWhcA4P8Q2ACgkf3nf/6nWrVqpUmTJtV7HQkJCS7T4eHhCgwMdF7adm37uXPnXNrKyso0c+ZMxcfHy263q3Xr1oqKilJxcbHOnz/v9ravXv52/XauV1RUpOLiYq1YsUJRUVEuX2PHjpUknT592u3tubs/ta2vqKhIJSUlNb5mob77cq2FCxdq7969io+P15133qns7OxK4bo67u5zYmJindZ7I3XFxsaqbdu2+uSTTyR983CTOXPmaMyYMS7H0hijkJCQWo8PAOAb3MMGAI3o8OHDWrFihV566SWXMxiXLl3SlStX9MUXXygsLEytWrWqcT1VPb2xuic6GmNcpidNmqSVK1cqMzNT/fr1U3h4uGw2m1JTU+VwOGrdh7pu53pX1z1mzBilpaVV2ef22293e3vu7k99679WffflWqNGjdLAgQO1bt06bdy4UYsWLdKCBQv0xz/+UcnJyTUu6+4+X3tGrjY3UlefPn2cge3FF1/UuXPn9POf/9ylT35+voKDgxUdHV3nmgCgOSOwAUAjOn78uBwOh55++mk9/fTTleYnJiZq8uTJbj050l1r165VWlqayxMkL1261KDv47LZbJXaoqKiFBoaqoqKCg0ZMqTBttXQ+xMVFaWwsDDt3bu3xj4NsS9t2rRRenq60tPTdfr0ad1xxx2aO3eukpOTqzyGVzXEPte0/prqqkmfPn307rvv6tixY/qv//ovTZgwweU9fvv371evXr1UUVGhkJAQde3aVXl5eXWuGQCaIwIbADSiHj16aN26dZXa//M//1OlpaX65S9/qU6dOt3UGnx9fSudTVqyZIkqKioabBtX3/N1bYDw9fXVgw8+qNWrV2vv3r2VLjksKiqq1z19Db0/Pj4+SklJ0ZtvvqmPP/640n1sxpgb3peKigp99dVXCg8Pd7ZFR0crLi5O5eXlkqo+hlc1xD5Xtf661FWTpKQkORwOPfroozLG6Nlnn3WZ37VrV82aNUsnTpzQ4sWL61wrADRnBDYAaEStW7dWSkpKpfarZ9SqmtfQfvCDH+i3v/2twsPD1a1bN+Xm5mrz5s2KjIxssG307t1bkvTss88qNTVV/v7+Gj58uJ5//nl9+OGH6tu3r5588kl169ZNX375pXbt2qXNmzfryy+/tMT+zJs3Txs3btQ999yjp556Sl27dtWJEye0Zs0abd26VRERETe0L6WlpWrXrp0eeugh9ezZUyEhIdq8ebPy8vKcZ82qO4YtWrRokH2uav0DBw7UrbfeWmNdNbkabv/3f/9X2dnZVYbWPXv26Hvf+16d6wSA5o7ABgDNzC9/+Uv5+vpq1apVunTpkvr376/Nmzdr6NChDbaNPn36aM6cOVq2bJnef/99ORwO5efnq0OHDvroo4/085//XH/84x/1q1/9SpGRkerevbvLe9A8vT9t27bVjh07NGPGDK1atUolJSVq27atkpOTne84i4mJqfe+BAcHKz09XRs3btQf//hHORwO3XLLLfrVr36lCRMmSKr+GLZo0aJB9rmq9R88eLDWumrSqlUrdejQQRcuXNCPf/zjKvvs2bOn1nfUAQD+j824c5c1AABANT7//HN961vf0osvvljlPZqXL19WSEiIzp0757wkEwBQMwIbAABoEA8//LB27typffv2KSAgoNL8s2fPqk2bNjp16pTzdQoAgJpxSSQAAKi34uJibdiwQX/729+0Zs0abdiwocqwJkmRkZF65JFHlJCQoO7du2v79u2NXC0AeB/OsAEAgHr7wx/+oIceekjt2rXTzJkz9eSTT3q6JABoUghsAAAAAGBRPp4uAAAAAABQNQIbAAAAAFiU2w8d+fvf/65FixZp586dOnHihNatW+fWi14dDocKCwsVGhoqm83m7uYBAAAAwOsZY1RaWqq4uDj5+FR/Hs3twHbhwgX17NlT48aN04gRI9wurLCwUPHx8W4vBwAAAABNTUFBgdq1a1ftfLcDW3JyspKTk+tdUGhoqLOwsLCweq8HAAAAALxVSUmJ4uPjnfmoOjf9PWzl5eUqLy93TpeWlkqSwsLCCGwAAAAAmrXabhO76YFt/vz5mj179s3eDLxch2l/rnH+F88Pa6RKAADXq+lntLf/fG7K+4bGx3jCzXDTnxI5ffp0nT9/3vlVUFBwszcJAAAAAE3CTT/DZrfbZbfbb/ZmAAAAAKDJ4T1sAAAAAGBRbp9h++qrr3TkyBHndH5+vnbv3q1WrVopISGhQYsDAAAAgObM7cD28ccf67vf/a5zOisrS5KUlpam119/vcEKAwAAAIDmzu3ANmjQIBljbkYtAAAAAIBrcA8bAAAAAFgUgQ0AAAAALIrABgAAAAAWRWADAAAAAIsisAEAAACARRHYAAAAAMCiCGwAAAAAYFEENgAAAACwKAIbAAAAAFgUgQ0AAAAALIrABgAAAAAWRWADAAAAAIsisAEAAACARRHYAAAAAMCiCGwAAAAAYFEENgAAAACwKAIbAAAAAFgUgQ0AAAAALIrABgAAAAAWRWADAAAAAIsisAEAAACARRHYAAAAAMCiCGwAAAAAYFEENgAAAACwKAIbAAAAAFgUgQ0AAAAALIrABgAAAAAWRWADAAAAAIsisAEAAACARRHYAAAAAMCiCGwAAAAAYFEENgAAAACwKAIbAAAAAFgUgQ0AAAAALIrABgAAAAAWRWADAAAAAIsisAEAAACARRHYAAAAAMCiCGwAAAAAYFEENgAAAACwKAIbAAAAAFgUgQ0AAAAALIrABgAAAAAWRWADAAAAAIsisAEAAACARRHYAAAAAMCiCGwAAAAAYFEENgAAAACwKAIbAAAAAFgUgQ0AAAAALIrABgAAAAAWVa/A9sorr6hDhw4KDAxU37599dFHHzV0XQAAAADQ7Lkd2N5++21lZWVp1qxZ2rVrl3r27KmhQ4fq9OnTN6M+AAAAAGi23A5sL774op588kmNHTtW3bp107JlyxQcHKzXXnvtZtQHAAAAAM2WnzudL1++rJ07d2r69OnONh8fHw0ZMkS5ublVLlNeXq7y8nLn9Pnz5yVJJSUl9akXTZSj/GKN8xkvAOA5Nf2M9vafz01539D4GE9wx9UxYYypsZ9bge3MmTOqqKhQTEyMS3tMTIwOHDhQ5TLz58/X7NmzK7XHx8e7s2k0c+EveboCAEBVmvLP56a8b2h8jCdUp7S0VOHh4dXOdyuw1cf06dOVlZXlnHY4HPryyy8VGRkpm812szcPNJqSkhLFx8eroKBAYWFhni4HaHCMcTR1jHE0ZYxv6zHGqLS0VHFxcTX2cyuwtW7dWr6+vjp16pRL+6lTpxQbG1vlMna7XXa73aUtIiLCnc0CXiUsLIwfhGjSGONo6hjjaMoY39ZS05m1q9x66EhAQIB69+6tLVu2ONscDoe2bNmifv36uV8hAAAAAKBabl8SmZWVpbS0NCUlJenOO+/USy+9pAsXLmjs2LE3oz4AAAAAaLbcDmwPP/ywioqKNHPmTJ08eVLf/va39f7771d6EAnQ3Njtds2aNavSJcBAU8EYR1PHGEdTxvj2XjZT23MkAQAAAAAe4faLswEAAAAAjYPABgAAAAAWRWADAAAAAIsisAEAAACARRHYADfNnz9fffr0UWhoqKKjo5WSkqKDBw+69Ll06ZIyMjIUGRmpkJAQPfjgg5VeOA94g+eff142m02ZmZnONsY3vN3x48c1ZswYRUZGKigoSLfddps+/vhj53xjjGbOnKk2bdooKChIQ4YM0eHDhz1YMVA3FRUVmjFjhhITExUUFKROnTppzpw5uvYZg4xv70NgA9yUk5OjjIwMbd++XZs2bdKVK1d077336sKFC84+U6ZM0Xvvvac1a9YoJydHhYWFGjFihAerBtyXl5en5cuX6/bbb3dpZ3zDm507d079+/eXv7+/NmzYoH379umFF15Qy5YtnX0WLlyoxYsXa9myZdqxY4datGihoUOH6tKlSx6sHKjdggULtHTpUr388svav3+/FixYoIULF2rJkiXOPoxvL2QA3JDTp08bSSYnJ8cYY0xxcbHx9/c3a9ascfbZv3+/kWRyc3M9VSbgltLSUtO5c2ezadMmc88995jJkycbYxjf8H4//elPzYABA6qd73A4TGxsrFm0aJGzrbi42NjtdvPWW281RolAvQ0bNsyMGzfOpW3EiBFm9OjRxhjGt7fiDBtwg86fPy9JatWqlSRp586dunLlioYMGeLs06VLFyUkJCg3N9cjNQLuysjI0LBhw1zGscT4hvd79913lZSUpJEjRyo6Olq9evXSr3/9a+f8/Px8nTx50mWMh4eHq2/fvoxxWN7dd9+tLVu26NChQ5KkTz/9VFu3blVycrIkxre38vN0AYA3czgcyszMVP/+/dWjRw9J0smTJxUQEKCIiAiXvjExMTp58qQHqgTc87vf/U67du1SXl5epXmMb3i7zz//XEuXLlVWVpZ+9rOfKS8vT08//bQCAgKUlpbmHMcxMTEuyzHG4Q2mTZumkpISdenSRb6+vqqoqNDcuXM1evRoSWJ8eykCG3ADMjIytHfvXm3dutXTpQANoqCgQJMnT9amTZsUGBjo6XKABudwOJSUlKR58+ZJknr16qW9e/dq2bJlSktL83B1wI35/e9/r1WrVmn16tXq3r27du/erczMTMXFxTG+vRiXRAL1NHHiRP3pT3/Shx9+qHbt2jnbY2NjdfnyZRUXF7v0P3XqlGJjYxu5SsA9O3fu1OnTp3XHHXfIz89Pfn5+ysnJ0eLFi+Xn56eYmBjGN7xamzZt1K1bN5e2rl276tixY5LkHMfXP/mUMQ5vMHXqVE2bNk2pqam67bbb9Nhjj2nKlCmaP3++JMa3tyKwAW4yxmjixIlat26d/vrXvyoxMdFlfu/eveXv768tW7Y42w4ePKhjx46pX79+jV0u4JbBgwfrs88+0+7du51fSUlJGj16tPPfjG94s/79+1d6FcuhQ4fUvn17SVJiYqJiY2NdxnhJSYl27NjBGIflXbx4UT4+rh/vfX195XA4JDG+vRWXRAJuysjI0OrVq/XOO+8oNDTUec13eHi4goKCFB4ervHjxysrK0utWrVSWFiYJk2apH79+umuu+7ycPVAzUJDQ533Y17VokULRUZGOtsZ3/BmU6ZM0d1336158+Zp1KhR+uijj7RixQqtWLFCkpzvHXzuuefUuXNnJSYmasaMGYqLi1NKSopniwdqMXz4cM2dO1cJCQnq3r27PvnkE7344osaN26cJMa31/L0YyoBbyOpyq+VK1c6+5SVlZn09HTTsmVLExwcbB544AFz4sQJzxUN3IBrH+tvDOMb3u+9994zPXr0MHa73XTp0sWsWLHCZb7D4TAzZswwMTExxm63m8GDB5uDBw96qFqg7kpKSszkyZNNQkKCCQwMNB07djTPPvusKS8vd/ZhfHsfmzHXvPocAAAAAGAZ3MMGAAAAABZFYAMAAAAAiyKwAQAAAIBFEdgAAAAAwKIIbAAAAABgUQQ2AAAAALAoAhsAAAAAWBSBDQAAAAAsisAGAAAAABZFYAMAAAAAiyKwAQAAAIBFEdgAAAAAwKL+H1jO7iHhaboHAAAAAElFTkSuQmCC",
                        "text/plain": [
                            "<Figure size 900x400 with 3 Axes>"
                        ]
                    },
                    "metadata": {},
                    "output_type": "display_data"
                }
            ],
            "source": [
                "class MP(object):\n",
                "    def __init__(self, c, M, P):\n",
                "        self.T = len(c)\n",
                "        self.M = M\n",
                "        self.P = P\n",
                "        self.c = c\n",
                "\n",
                "        self.modfile = \"maintenance_planning.mod\"\n",
                "\n",
                "        self.ampl = AMPL()\n",
                "        self.ampl.option[\"solver\"] = SOLVER\n",
                "        self.solved = False\n",
                "\n",
                "    def load_data(self):\n",
                "        self.ampl.param[\"T\"] = self.T\n",
                "        self.ampl.param[\"P\"] = self.P\n",
                "        self.ampl.param[\"M\"] = self.M\n",
                "        self.ampl.param[\"c\"] = self.c\n",
                "\n",
                "    def set_solver(self, solver):\n",
                "        self.ampl.option[\"solver\"] = solver\n",
                "\n",
                "    def set_modfile(self, modfile):\n",
                "        self.modfile = modfile\n",
                "\n",
                "    def solve(self):\n",
                "        self.ampl.read(self.modfile)\n",
                "        self.load_data()\n",
                "        self.ampl.solve()\n",
                "        self.solved = True\n",
                "\n",
                "    def plot_schedule(self):\n",
                "        if not self.solved:\n",
                "            self.solve()\n",
                "\n",
                "        T = self.T\n",
                "        M = self.M\n",
                "        P = self.P\n",
                "        PH = list(range(1, T + 1))\n",
                "        Y = list(range(1, T - M + 2))\n",
                "        S = list(range(M))\n",
                "\n",
                "        fig, ax = plt.subplots(3, 1, figsize=(9, 4))\n",
                "\n",
                "        ax[0].bar(PH, [self.c[t] for t in PH])\n",
                "        ax[0].set_title(\"daily profit $c_t$\")\n",
                "\n",
                "        x = self.ampl.var[\"x\"].to_dict()\n",
                "\n",
                "        ax[1].bar(PH, [x[t] for t in PH], label=\"normal operation\")\n",
                "        ax[1].set_title(\"unit operating schedule $x_t$\")\n",
                "\n",
                "        y = self.ampl.var[\"y\"].to_dict()\n",
                "\n",
                "        ax[2].bar(Y, [y[t] for t in Y])\n",
                "        ax[2].set_title(str(P) + \" maintenance starts $y_t$\")\n",
                "        for a in ax:\n",
                "            a.set_xlim(0.1, T + 0.9)\n",
                "\n",
                "        plt.tight_layout()\n",
                "\n",
                "\n",
                "model = MP(c, M, P)\n",
                "model.plot_schedule()"
            ]
        },
        {
            "cell_type": "markdown",
            "metadata": {
                "colab_type": "text",
                "id": "gcy3yD4J_6OC"
            },
            "source": [
                "## Ramping constraints\n",
                "\n",
                "Prior to maintenance shutdown, a large processing unit may take some time to safely ramp down from full production. And then require more time to safely ramp back up to full production following maintenance. To provide for ramp-down and ramp-up periods, we modify the problem formation in the following ways.\n",
                "\n",
                "* The variable denoting unit operation, $x_t$ is changed from a binary variable to a continuous variable $0 \\leq x_t \\leq 1$ denoting the fraction of total capacity at which the unit is operating on day $t$.\n",
                "\n",
                "* Two new variable sequences, $0 \\leq u_t^+ \\leq u_t^{+,\\max}$ and $0\\leq u_t^- \\leq u_t^{-,\\max}$, are introduced which denote the fraction increase or decrease in unit capacity to completed on day $t$.\n",
                "\n",
                "* An additional sequence of equality constraints is introduced relating $x_t$ to $u_t^+$ and $u_t^-$.\n",
                "\n",
                "$$\n",
                "\\begin{align*}\n",
                "x_{t} & = x_{t-1} + u^+_t - u^-_t\n",
                "\\end{align*}\n",
                "$$\n",
                "\n",
                "We begin the AMPL model by specifying the constraints, then modifying the Big-M formulation to add the features described above."
            ]
        },
        {
            "cell_type": "code",
            "execution_count": 6,
            "metadata": {
                "colab": {},
                "colab_type": "code",
                "id": "9fRtZc5nCInh"
            },
            "outputs": [
                {
                    "name": "stdout",
                    "output_type": "stream",
                    "text": [
                        "Writing maintenance_planning_ramp.mod\n"
                    ]
                }
            ],
            "source": [
                "%%writefile maintenance_planning_ramp.mod\n",
                "\n",
                "param T; # number of planning periods\n",
                "param P; # number of maintenance periods\n",
                "param M; # number of days for each maintenance period\n",
                "\n",
                "param upos_max;\n",
                "param uneg_max;\n",
                "\n",
                "set PH := 1..T; # set of planning horizon time periods\n",
                "\n",
                "set Y := 1..(T - M + 1);\n",
                "set S := 0..(M - 1);\n",
                "\n",
                "param c{PH}; # profit\n",
                "\n",
                "var x{PH} >= 0, <= 1;\n",
                "var y{PH} binary;\n",
                "var upos{PH} >= 0, <= upos_max;\n",
                "var uneg{PH} >= 0, <= uneg_max;\n",
                "\n",
                "maximize profit: sum{t in PH} c[t] * x[t];\n",
                "\n",
                "s.t. required_maintenance: P == sum{t in Y} y[t];\n",
                "s.t. no_overlap {t in Y}: sum{s in S} y[t + s] <= 1;\n",
                "s.t. required_shutdown {t in Y}:\n",
                "    y[t] == 1 ==> sum{s in S} x[t + s] == 0;\n",
                "s.t. ramp {t in PH: t != 1}:\n",
                "    x[t] == x[t-1] + upos[t] - uneg[t];"
            ]
        },
        {
            "cell_type": "code",
            "execution_count": 7,
            "metadata": {
                "colab": {
                    "base_uri": "https://localhost:8080/",
                    "height": 297
                },
                "colab_type": "code",
                "executionInfo": {
                    "elapsed": 3307,
                    "status": "ok",
                    "timestamp": 1557947674808,
                    "user": {
                        "displayName": "Jeffrey Kantor",
                        "photoUrl": "https://lh5.googleusercontent.com/-8zK5aAW5RMQ/AAAAAAAAAAI/AAAAAAAAKB0/kssUQyz8DTQ/s64/photo.jpg",
                        "userId": "09038942003589296665"
                    },
                    "user_tz": 240
                },
                "id": "ESHZWcEaA9gP",
                "outputId": "030e96a8-1ae1-47f5-f4c6-d83378818838"
            },
            "outputs": [
                {
                    "name": "stdout",
                    "output_type": "stream",
                    "text": [
                        "cbc 2.10.7: \b\b\b\b\b\b\b\b\b\b\b\bcbc 2.10.7: optimal solution; objective 39.53508979\n",
                        "4839 simplex iterations\n",
                        "4839 barrier iterations\n",
                        "12 branching nodes\n"
                    ]
                },
                {
                    "data": {
                        "image/png": "iVBORw0KGgoAAAANSUhEUgAAA2wAAAGGCAYAAAAHLAKcAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjcuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/bCgiHAAAACXBIWXMAAA9hAAAPYQGoP6dpAABDXklEQVR4nO3deXxV9Z3/8ffNdpOQDUIWAgkEZGRTigQRgUoLI2YoTkTBKGgEqlMSkJApFuoAocg+WgtaljpirWAttKC2RVnqpMMQMIKIlL0GSAlLQEIihIC5398f/rjDNesNIffc5PV8PPJ4cL7ne875nJMvyX3nbDZjjBEAAAAAwHJ8PF0AAAAAAKBqBDYAAAAAsCgCGwAAAABYFIENAAAAACyKwAYAAAAAFkVgAwAAAACLIrABAAAAgEUR2AAAAADAoghsAAAAAGBRBDYAAAAAsCgCGwAAAABYFIENAFAv2dnZstlsbi/3xhtvyGaz6dixYzW2eZO8vDzde++9atGihWw2m/bs2eP1+wQAsAYCGwAAN+HatWsaOXKkvvzyS/385z/Xb37zG7Vv377Kvtu3b1d2draKi4sbt0gAgNcisAEAGtUTTzyhsrKyakONt/n73/+u48eP68c//rGeeeYZjRkzRi1btqxyP7dv367Zs2cT2AAAdebn6QIAAM2Lr6+vfH19PV1GtS5duqQWLVrUuf/Zs2clSRERES7tVt9PAIB34AwbAKBW27ZtU58+fRQYGKhOnTppxYoVlfocP35c6enpuv322xUUFKTIyEiNHDmy0j1ctd3b9dFHH8lms2n9+vWV5q1Zs0Y2m025ubk11nv9/rqDBw9q1KhRCgsLU2RkpCZPnqwrV65U6rd//349/vjjatmypQYMGOCc/+mnnyo5OVlhYWEKCQnR4MGDtWPHDuf8p556Svfdd58kaeTIkbLZbBo0aFCV+5mdna2pU6dKkhITE2Wz2ep0j9uuXbv0yCOPKDo6WoGBgerSpYt+9rOf1bgMAKDp4AwbAKBGn3/+ue6//35FRUUpOztbX3/9tWbNmqWYmBiXfnl5edq+fbtSU1PVrl07HTt2TMuWLdOgQYO0f/9+BQcH12l7gwYNUnx8vFavXq2HHnrIZd7q1avVqVMn9evXr07rGjVqlDp06KD58+drx44dWrJkiS5cuKA333zTpd/IkSPVuXNnzZs3T8YYSdLf/vY3DRw4UGFhYXruuefk7++vFStWaNCgQcrJyVHfvn31b//2b2rbtq3mzZunZ599Vn369Kl0XK4bMWKEDh8+rLfffls///nP1bp1a0lSVFRUtfWvX79eqamp6tixo6ZOnaqQkBDncQYANBMGAIAapKSkmMDAQHP8+HFn2/79+42vr6+58dfI5cuXKy2bm5trJJk333zT2bZq1SojyeTn51fbNn36dGO3201xcbGzz9mzZ42fn5+ZNWtWrTXPmjXLSDIPPvigS3t6erqRZD777DOXfo899liV+x0QEGD+/ve/O9sKCwtNaGio+e53v+ts++ijj4wks3btWpflq9rPxYsXV2qrzqFDh0xwcLB56KGHTFlZmcu8K1eu1Lo8AKBp4JJIAEC1Kioq9OGHHyolJUUJCQnO9q5du2ro0KEufYOCgpz/vnbtms6fP6/bbrtNERER2r17t1vbffLJJ1VeXq5169Y529555x19/fXXGjNmTJ3Xk5GR4TI9adIkSdKf//xnl/Yf/ehHLtMVFRXatGmTUlJS1LFjR2d7mzZt9Pjjj2vbtm0qKSmpcx31MXPmTAUGBuqNN95QYGCgyzy73X5Ltw0AsA4CGwCgWkVFRSorK1Pnzp0rzbv99ttdpsvKyjRz5kzFx8fLbrerdevWioqKUnFxsS5evOjWdrt06aI+ffpo9erVzrbVq1frnnvu0W233Vbn9Xy77k6dOsnHx6fSfWOJiYku00VFRbp8+XKlfZS+CasOh0MFBQV1rsNd5eXlev/99zVmzBiFhYXV2NcYo5CQEOfDTwAATQv3sAEAGsSkSZO0atUqZWZmql+/fgoPD5fNZlNqaqocDofb63vyySc1efJk/eMf/1B5ebl27NihV1555aZqrO5F3zeeHbSCL774QpcvX1bv3r1r7Zufn6/g4GBFR0c3QmUAgMZGYAMAVCsqKkpBQUE6cuRIpXmHDh1ymV63bp3S0tL04osvOtuuXLlS73eOpaamKisrS2+//bbKysrk7++vRx991K11HDlyxOXs2dGjR+VwONShQ4cal4uKilJwcHClfZSkgwcPysfHR/Hx8W7VIlUfGL+trKysTv0PHDigXr16qaKiQiEhIeratavy8vLcrgsAYF1cEgkAqJavr6+GDh2qDRs26MSJE872AwcO6MMPP6zU1/z/Jyxet3TpUlVUVNRr261bt1ZycrLeeustrV69Wg888IDzyYp19eqrr1aqR5KSk5NrXM7X11f333+/3n33XZfLJ8+cOaM1a9ZowIABtV6qWJXr73erLcRev+xzy5YtleZdu3bN+e+uXbtq1qxZmjBhgr766ivCGgA0QZxhAwDUaPbs2frggw80cOBApaen6+uvv9bSpUvVvXt37d2719nvBz/4gX7zm98oPDxc3bp1U25urrZs2aLIyMh6b/vJJ5/UI488IkmaM2eO28vn5+frwQcf1AMPPKDc3Fy99dZbevzxx9WzZ89al33hhRe0efNmDRgwQOnp6fLz89OKFStUXl6uRYsWuV2LJOcljs8//7xSU1Pl7++v4cOHV3pRd1hYmJ566im98cYbKi8v16BBg1RaWqqPPvpIw4YNc3mYyt69e/X973+/XvUAAKyPwAYAqNGdd96pDz/8UFlZWZo5c6batWun2bNn69SpUy6B7Re/+IV8fX21evVqXblyRf3799eWLVsqPU3SHcOHD1fLli3lcDj04IMPur38O++8o5kzZ2ratGny8/PTxIkTtXjx4jot2717d/3P//yPpk+frvnz58vhcKhv375666231LdvX7drkaQ+ffpozpw5Wr58uT744AM5HA7l5+dXCmzSN2cH27Ztq9/97nfasGGDWrVqpXvvvVcPPPCAS7+9e/cqMzOzXvUAAKzPZr59/QoAABbx9ddfKy4uTsOHD9d//dd/1Xm57OxszZ49W0VFRW5fRulNrl69qpCQEF24cKHK0AcA8H7cwwYAsKwNGzaoqKhITz75pKdLsaTS0lJJ3wQ3AEDTxCWRAADL2blzp/bu3as5c+aoV69euu+++zxdkiVFRkbqscceU0JCgrp3764dO3Z4uiQAQAPjDBsAwHKWLVumCRMmKDo6Wm+++aany7G0X//61yotLSWsAUATxT1sAAAAAGBRnGEDAAAAAItyO7D99a9/1fDhwxUXFyebzaYNGzbcgrIAAAAAAG4/dOTSpUvq2bOnxo0bpxEjRri9QYfDocLCQoWGhspms7m9PAAAAAB4O2OMSktLFRcXJx+f6s+juR3YkpOTlZycXO/CCgsLFR8fX+/lAQAAAKCpKCgoULt27aqd3+iP9Q8NDZX0TWFhYWGNvXkAAAAA8LiSkhLFx8c781F1bnlgKy8vV3l5uXP6+ks+w8LCCGwAAAAAmrXabhO75YFt/vz5mj179q3eDAAAcFOHaX+qcf6xBcMaqRIAQHVueWCbPn26srKynNPXT/0BzUlNH4r4QAQAwP9pyD8k8PsXTcEtD2x2u112u/1Wb6ZZ4IcOAAAA0Ly4Hdi++uorHT161Dmdn5+vPXv2qFWrVkpISGjQ4gA0DVx2ZU38EQgAAOtzO7B98skn+t73vuecvn65Y1pamt54440GKwxNBx/WAQAA0Bzcij+Guh3YBg0aJGNMvTYGAAAAAKi76l+pDQAAAADwqEZ/cTaqxr0kAAAAAL6NwAbcJMI2AAAAbhUuiQQAAAAAi+IMG5olnlwJAEDD4ooT4NYgsAEAADRjBC3A2ghsgEXwCxMAAADfRmBrpggHAAAAgPUR2IBqePN9bgTy2nnz9xcAADQfPCUSAAAAACyKM2wAAABADbgqA57EGTYAAAAAsCivOcPGXzYAAAAANDdeE9gAAE0fD8wBgKaNkzDuI7ABAKrV1H+xEhABAFZHYAMANDlNPWjCe/FHAgDuIrDdBD4QoKnjgwWAxsbvVgBwxVMiAQAAAMCiOMMGABZR1zMLdTnzydlRAACahiYX2PiQAgAAgLriMlxYXZMLbAAAAADQ0Dx1YojAhpvCGU2A/wdofIw5AE0dP+f+jyUCW2OfimYAAAAAAPAGlghssCaCLQAAAOBZBDYAlsEfCQAAQENqCg+V4T1sAAAAAGBRBDYAAAAAsCguiQQAAE0Sl1k3bXx/0VwQ2AAv0hSuw/Y2fCAAAACeRGBrBHzgAwB+FjY2jjcANA0EtmpwJgMAAACApxHYAHgVzhoAAND4+P3rOQS2JoYzgwAAAEDTQWADAABogvgjLtA0ENgAADeFD4UNi8uOAAA34sXZAAAAAGBRBDYAAAAAsCguiQQAoBFw6SgAq+JSbGsjsAEAAABehj8CNR8ENgAAAKCRcDYL7iKwAWhy+KsjAABNX3MJvwQ2AADgVRryjzLN5QMfAO/FUyIBAAAAwKI4wwYAQA04A1M3HCcAuDUIbAAAAF6GgAw0H1wSCQAAAAAWxRk2AADQKDgrBADuI7ABAADUgFeFAPAkAhuAm8IHGQBoWJyJBHAjAhvQTPGBAGg4/OECQENqyN/R/L73fgQ2AABQIz7w1Y7QDuBWIbDBq/ChAQCAxtdQv38JtrAiq4/LegW2V199VYsXL9bp06fVs2dPLV26VHfffXdD14ZmxOr/UQAAAABPcDuwvfPOO8rKytLy5cvVt29fvfzyyxo6dKgOHTqk6OjoW1EjAABO/IEHAKyJn8+3htuB7aWXXtLTTz+tsWPHSpKWL1+uP/3pT3r99dc1bdq0Bi8QAAAA8AbcuoFbwa3AdvXqVe3atUvTp093tvn4+GjIkCHKzc1t8OKA+uCHJayIcQkAAOrDrcB27tw5VVRUKCYmxqU9JiZGBw8erHKZ8vJylZeXO6cvXrwoSSopKXG2Ocov17jdkpKSOvWpbV0N1ceKNblTtxVrslLdVqyJ72/TronvLzU1dk2MuaZdE9/fpl0T31/vr+nb08aYGrdvM7X1uEFhYaHatm2r7du3q1+/fs725557Tjk5Odq5c2elZbKzszV79uy6bgIAAAAAmo2CggK1a9eu2vlunWFr3bq1fH19debMGZf2M2fOKDY2tsplpk+frqysLOe0w+HQl19+qcjISNlsNnc2D1haSUmJ4uPjVVBQoLCwME+XAzQ4xjiaOsY4mjLGt/UYY1RaWqq4uLga+7kV2AICAtS7d29t3bpVKSkpkr4JYFu3btXEiROrXMZut8tut7u0RUREuLNZwKuEhYXxgxBNGmMcTR1jHE0Z49tawsPDa+3j9lMis7KylJaWpqSkJN199916+eWXdenSJedTIwEAAAAADcPtwPboo4+qqKhIM2fO1OnTp/Wd73xHH3zwQaUHkQAAAAAAbo7bgU2SJk6cWO0lkEBzZbfbNWvWrEqXAANNBWMcTR1jHE0Z49t7ufWUSAAAAABA4/HxdAEAAAAAgKoR2AAAAADAoghsAAAAAGBRBDbATfPnz1efPn0UGhqq6OhopaSk6NChQy59rly5ooyMDEVGRiokJEQPP/xwpRfOA95gwYIFstlsyszMdLYxvuHtTp48qTFjxigyMlJBQUG644479MknnzjnG2M0c+ZMtWnTRkFBQRoyZIiOHDniwYqBuqmoqNCMGTOUmJiooKAgderUSXPmzNGNj6xgfHsfAhvgppycHGVkZGjHjh3avHmzrl27pvvvv1+XLl1y9pkyZYref/99rV27Vjk5OSosLNSIESM8WDXgvry8PK1YsUJ33nmnSzvjG97swoUL6t+/v/z9/bVx40bt379fL774olq2bOnss2jRIi1ZskTLly/Xzp071aJFCw0dOlRXrlzxYOVA7RYuXKhly5bplVde0YEDB7Rw4UItWrRIS5cudfZhfHshA+CmnD171kgyOTk5xhhjiouLjb+/v1m7dq2zz4EDB4wkk5ub66kyAbeUlpaazp07m82bN5v77rvPTJ482RjD+Ib3+8lPfmIGDBhQ7XyHw2FiY2PN4sWLnW3FxcXGbrebt99+uzFKBOpt2LBhZty4cS5tI0aMMKNHjzbGML69FWfYgJt08eJFSVKrVq0kSbt27dK1a9c0ZMgQZ58uXbooISFBubm5HqkRcFdGRoaGDRvmMo4lxje833vvvaekpCSNHDlS0dHR6tWrl371q1855+fn5+v06dMuYzw8PFx9+/ZljMPy7r33Xm3dulWHDx+WJH322Wfatm2bkpOTJTG+vVW9XpwN4BsOh0OZmZnq37+/evToIUk6ffq0AgICFBER4dI3JiZGp0+f9kCVgHt++9vfavfu3crLy6s0j/ENb/fFF19o2bJlysrK0k9/+lPl5eXp2WefVUBAgNLS0pzjOCYmxmU5xji8wbRp01RSUqIuXbrI19dXFRUVmjt3rkaPHi1JjG8vRWADbkJGRob27dunbdu2eboUoEEUFBRo8uTJ2rx5swIDAz1dDtDgHA6HkpKSNG/ePElSr169tG/fPi1fvlxpaWkerg64Ob/73e+0evVqrVmzRt27d9eePXuUmZmpuLg4xrcX45JIoJ4mTpyoP/7xj/roo4/Url07Z3tsbKyuXr2q4uJil/5nzpxRbGxsI1cJuGfXrl06e/as7rrrLvn5+cnPz085OTlasmSJ/Pz8FBMTw/iGV2vTpo26devm0ta1a1edOHFCkpzj+NtPPmWMwxtMnTpV06ZNU2pqqu644w498cQTmjJliubPny+J8e2tCGyAm4wxmjhxotavX6+//OUvSkxMdJnfu3dv+fv7a+vWrc62Q4cO6cSJE+rXr19jlwu4ZfDgwfr888+1Z88e51dSUpJGjx7t/DfjG96sf//+lV7FcvjwYbVv316SlJiYqNjYWJcxXlJSop07dzLGYXmXL1+Wj4/rx3tfX185HA5JjG9vxSWRgJsyMjK0Zs0avfvuuwoNDXVe8x0eHq6goCCFh4dr/PjxysrKUqtWrRQWFqZJkyapX79+uueeezxcPVCz0NBQ5/2Y17Vo0UKRkZHOdsY3vNmUKVN07733at68eRo1apQ+/vhjrVy5UitXrpQk53sHX3jhBXXu3FmJiYmaMWOG4uLilJKS4tnigVoMHz5cc+fOVUJCgrp3765PP/1UL730ksaNGyeJ8e21PP2YSsDbSKrya9WqVc4+ZWVlJj093bRs2dIEBwebhx56yJw6dcpzRQM34cbH+hvD+Ib3e//9902PHj2M3W43Xbp0MStXrnSZ73A4zIwZM0xMTIyx2+1m8ODB5tChQx6qFqi7kpISM3nyZJOQkGACAwNNx44dzfPPP2/Ky8udfRjf3sdmzA2vPgcAAAAAWAb3sAEAAACARRHYAAAAAMCiCGwAAAAAYFEENgAAAACwKAIbAAAAAFgUgQ0AAAAALIrABgAAAAAWRWADAAAAAIsisAEAAACARRHYAAAAAMCiCGwAAAAAYFEENgAAAACwKAIbAAAAAFgUgQ0AAAAALIrABgAAAAAWRWADAAAAAIsisAEAAACARRHYAAAAAMCiCGwAAAAAYFEENgCwsDfeeEM2m03Hjh3zdCmW0hSOS3Z2tmw2m86dO9do26qPpnCsAcCbEdgAwIts375d2dnZKi4u9nQpjaK57S8AAN9GYAMAC3viiSdUVlam9u3bS/omwMyePbvZBJjq9vfbxwUAgKbKz9MFAACq5+vrK19fX0+X0WAuXbqkFi1a3PR6mtpxAQCgOpxhA4AG8tRTT6lDhw6V2r99/9D16aNHj+qpp55SRESEwsPDNXbsWF2+fNll2RvvH8rOztbUqVMlSYmJibLZbHW6t+jTTz9VcnKywsLCFBISosGDB2vHjh1V1njw4EGNGjVKYWFhioyM1OTJk3XlypVK6zx58qTGjRunmJgY2e12de/eXa+//nqV69y/f78ef/xxtWzZUgMGDJAkHT9+XOnp6br99tsVFBSkyMhIjRw50mVfatrfqu6rcue4StJ///d/KykpSYGBgerUqZNWrFhRp3u9SktLlZmZqQ4dOshutys6Olr//M//rN27d1c6RuPHj1dcXJzsdrsSExM1YcIEXb161aVfcXFxneqtyzGXpG3btqlPnz4u+1WVuo7X6tS1nqqWCwwM1Lhx41zat2zZIn9/f02ZMqXWdQBAc8IZNgDwkFGjRikxMVHz58/X7t279dprryk6OloLFy6ssv+IESN0+PBhvf322/r5z3+u1q1bS5KioqKq3cbf/vY3DRw4UGFhYXruuefk7++vFStWaNCgQcrJyVHfvn0r1dShQwfNnz9fO3bs0JIlS3ThwgW9+eabzj5nzpzRPffcI5vNpokTJyoqKkobN27U+PHjVVJSoszMTJd1jhw5Up07d9a8efNkjJEk5eXlafv27UpNTVW7du107NgxLVu2TIMGDdL+/fsVHBxcr/2t63H99NNP9cADD6hNmzaaPXu2Kioq9LOf/azWdUvSj370I61bt04TJ05Ut27ddP78eW3btk0HDhzQXXfdJUkqLCzU3XffreLiYj3zzDPq0qWLTp48qXXr1uny5csKCAhwq966HvPPP/9c999/v6KiopSdna2vv/5as2bNUkxMTK375Q53x8CN2rZtqx/+8IdauXKlZs2apfbt2+vgwYMaOXKkkpOT9eKLLzZorQDg9QwAoEGkpaWZ9u3bV2qfNWuWufHH7fXpcePGufR76KGHTGRkpEvbqlWrjCSTn59vjDFm8eLFLtO1SUlJMQEBAebvf/+7s62wsNCEhoaa7373u5VqevDBB12WT09PN5LMZ5995mwbP368adOmjTl37pxL39TUVBMeHm4uX77sss7HHnusUl3X+9woNzfXSDJvvvmms626/f32cblxe3U5rsOHDzfBwcHm5MmTzrYjR44YPz8/U9uvxvDwcJORkVFjnyeffNL4+PiYvLy8SvMcDofb9db1mKekpJjAwEBz/PhxZ5/9+/cbX1/fSvtV1/Fa1bGuaz3V+cc//mHsdruZMGGCOXfunOnUqZP5zne+Y7766qsalwOA5ohLIgHAQ370ox+5TA8cOFDnz59XSUlJg6y/oqJCmzZtUkpKijp27Ohsb9OmjR5//HFt27at0rYyMjJcpidNmiRJ+vOf/yxJMsbo97//vYYPHy5jjM6dO+f8Gjp0qC5evFjp0sBv76ckBQUFOf997do1nT9/XrfddpsiIiIqLe+u2o5rRUWFtmzZopSUFMXFxTn73XbbbUpOTq51/REREdq5c6cKCwurnO9wOLRhwwYNHz5cSUlJleZ/+3LD2uqt6zGvqKjQhx9+qJSUFCUkJDjX17VrVw0dOrTW/aqr+oyBb2vbtq2efvppvf766xo2bJjKysr0xz/+sUHubwSApobABgAecuOHaklq2bKlJOnChQsNsv6ioiJdvnxZt99+e6V5Xbt2lcPhUEFBgUt7586dXaY7deokHx8f571iRUVFKi4u1sqVKxUVFeXyNXbsWEnS2bNnXdaRmJhYaftlZWWaOXOm4uPjZbfb1bp1a0VFRam4uFgXL168md2u9biePXtWZWVluu222yotW1Xbty1atEj79u1TfHy87r77bmVnZ+uLL75wzi8qKlJJSYl69OjRIPXW9ZgXFRWprKys0vdQUpVjoL7qMwaq8uMf/1jl5eXau3ev3nvvPbVt29ZlvjFGISEhdVoXADRl3MMGAA2kugc1VFRUVNle3VMOzf+/z8sKvr1PDodDkjRmzBilpaVVucydd97pMn3j2bTrJk2apFWrVikzM1P9+vVTeHi4bDabUlNTnduor1t9XEeNGqWBAwdq/fr12rRpkxYvXqyFCxfqD3/4Q53O0H1bbfXW9Zi7e9zcHa/X1WcMVGXu3LmSpK+//lqtWrWqND8/P1/BwcGKjo6udV0A0JQR2ACggbRs2bLK96MdP368wbZRl6f3XRcVFaXg4GAdOnSo0ryDBw/Kx8dH8fHxLu1HjhxxOSN29OhRORwO59MEo6KiFBoaqoqKCg0ZMqR+OyFp3bp1SktLc3nAxJUrVyodP3f2t66io6MVGBioo0ePVppXVVtV2rRpo/T0dKWnp+vs2bO66667NHfuXCUnJysqKkphYWHat29fg9Rb12NeUVGhoKAgHTlypNK8qsZAfcdrQ4yBxYsX67XXXtMrr7yiqVOnau7cuXrttdec8w8cOKBevXqpoqJCISEh6tq1q/Ly8uq1LQDwdlwSCQANpFOnTrp48aL27t3rbDt16pTWr1/fYNu4fo9PXV6c7evrq/vvv1/vvvuuy+Pvz5w5ozVr1mjAgAEKCwtzWebVV191mV66dKkkOc8c+fr66uGHH9bvf//7KgNJUVFRnfbD19e30hmvpUuXVjq7487+1pWvr6+GDBmiDRs2uNyHdvToUW3cuLHGZSsqKipdshkdHa24uDiVl5dLknx8fJSSkqL3339fn3zySaV1uHumr67H3NfXV0OHDtWGDRt04sQJ5/wDBw7oww8/rLRcfcfrzY6BDRs2aNq0aZozZ44yMjL0zDPP6M0331R+fr6zT9euXTVr1ixNmDBBX331FWENQLPGGTYAaCCpqan6yU9+ooceekjPPvusLl++rGXLlumf/umfbvpBGtf17t1bkvT8888rNTVV/v7+Gj58eLUPa3jhhRe0efNmDRgwQOnp6fLz89OKFStUXl6uRYsWVeqfn5+vBx98UA888IByc3P11ltv6fHHH1fPnj2dfRYsWKCPPvpIffv21dNPP61u3brpyy+/1O7du7VlyxZ9+eWXte7HD37wA/3mN79ReHi4unXrptzcXG3ZskWRkZF12t+blZ2drU2bNql///6aMGGCKioq9Morr6hHjx7as2dPtcuVlpaqXbt2euSRR9SzZ0+FhIRoy5YtysvLczlbOG/ePG3atEn33XefnnnmGXXt2lWnTp3S2rVrtW3bNkVERLhVb12P+ezZs/XBBx9o4MCBSk9P19dff62lS5eqe/fuLsFMurnxWt8xsGvXLo0ePVqjR4/W888/L0l67rnntHz58kpn2fbu3avvf//7bh0nAGiSPPNwSgBomjZt2mR69OhhAgICzO23327eeuutah/rX1RU5LJsVY9Pr6ptzpw5pm3btsbHx6dOj/jfvXu3GTp0qAkJCTHBwcHme9/7ntm+fbtLn+s17d+/3zzyyCMmNDTUtGzZ0kycONGUlZVVWueZM2dMRkaGiY+PN/7+/iY2NtYMHjzYrFy5stb9NMaYCxcumLFjx5rWrVubkJAQM3ToUHPw4EHTvn17k5aW5tK3qv2t6bH+dTmuxhizdetW06tXLxMQEGA6depkXnvtNfPv//7vJjAwsNpjWV5ebqZOnWp69uxpQkNDTYsWLUzPnj3NL3/5y0p9jx8/bp588kkTFRVl7Ha76dixo8nIyDDl5eX1qrcux9wYY3Jyckzv3r1NQECA6dixo1m+fHmlMXhdXcbrzdZzXUFBgWnTpo3p37+/uXLlisu8CRMmGH9/f/PFF18427p162Z27NhR5boAoDmxGWOhu9sBAB6RnZ2t2bNnq6ioyPmC6uYoJSVFf/vb36q8DwyN5+rVqwoJCdGFCxd41D+AZo972AAAzVJZWZnL9JEjR/TnP/9ZgwYN8kxBcCotLZX0TXADgOaOe9gAAM1Sx44d9dRTT6ljx446fvy4li1bpoCAAD333HOeLq3Zi4yM1GOPPaaEhAR1795dO3bs8HRJAOAxBDYAQLP0wAMP6O2339bp06dlt9vVr18/zZs3r8oXT6Px/frXv9avf/1rT5cBAB7HPWwAAAAAYFHcwwYAAAAAFkVgAwAAAACLcjuw/fWvf9Xw4cMVFxcnm82mDRs23IKyAAAAAABuP3Tk0qVL6tmzp8aNG6cRI0a4vUGHw6HCwkKFhobKZrO5vTwAAAAAeDtjjEpLSxUXFycfn+rPo7kd2JKTk5WcnFzvwgoLCxUfH1/v5QEAAACgqSgoKFC7du2qnd/oj/UPDQ2V9E1hYWFhjb15AAAAAPC4kpISxcfHO/NRdW55YCsvL1d5eblzurS0VJIUFhZGYAMAAADQrNV2m9gtD2zz58/X7Nmza+zTYdqfapx/bMGwOvWpbV0N1ceKNblTtxVrslLdVqyJ72/TronvLzU1dk3ePubQsKz2/bXSuLyxHxqGFb+/jV2Tu275Y/2nT5+uixcvOr8KCgpu9SYBAAAAoEm45WfY7Ha77Hb7rd4MAAAAADQ5bge2r776SkePHnVO5+fna8+ePWrVqpUSEhIatDgAAAAAaM7cDmyffPKJvve97zmns7KyJElpaWl64403GqwwAAAAAGju3A5sgwYNkjHmVtQCAAAAALjBLX/oCAAAAACgfghsAAAAAGBRBDYAAAAAsCgCGwAAAABYFIENAAAAACyKwAYAAAAAFkVgAwAAAACLIrABAAAAgEUR2AAAAADAoghsAAAAAGBRBDYAAAAAsCgCGwAAAABYFIENAAAAACyKwAYAAAAAFuXn6QIAAADQ8DpM+1ON848tGNZIlaC5YMzdGpxhAwAAAACLIrABAAAAgEUR2AAAAADAoghsAAAAAGBRBDYAAAAAsCgCGwAAAABYFIENAAAAACyKwAYAAAAAFkVgAwAAAACLIrABAAAAgEX5eboAAAAAd3SY9qca5x9bMKyRKkFzwZiDJ3GGDQAAAAAsisAGAAAAABZFYAMAAAAAiyKwAQAAAIBFEdgAAAAAwKIIbAAAAABgUQQ2AAAAALAoAhsAAAAAWBSBDQAAAAAsys/TBQAAAMA9Hab9qdp5xxYMa8RK0Fww5jyHM2wAAAAAYFEENgAAAACwKAIbAAAAAFgUgQ0AAAAALIrABgAAAAAWRWADAAAAAIsisAEAAACARRHYAAAAAMCiCGwAAAAAYFEENgAAAACwKAIbAAAAAFgUgQ0AAAAALIrABgAAAAAWRWADAAAAAIsisAEAAACARdUrsL366qvq0KGDAgMD1bdvX3388ccNXRcAAAAANHtuB7Z33nlHWVlZmjVrlnbv3q2ePXtq6NChOnv27K2oDwAAAACaLbcD20svvaSnn35aY8eOVbdu3bR8+XIFBwfr9ddfvxX1AQAAAECz5VZgu3r1qnbt2qUhQ4b83wp8fDRkyBDl5uY2eHEAAAAA0Jz5udP53LlzqqioUExMjEt7TEyMDh48WOUy5eXlKi8vd05fvHhRklRSUuJsc5RfrnG7JSUldepT27oaqo8Va3KnbivWZKW6rVgT39+mXRPfX2pq7Jqay5hryvj+WrOmpozvb8PV9O1pY0yN27eZ2nrcoLCwUG3bttX27dvVr18/Z/tzzz2nnJwc7dy5s9Iy2dnZmj17dl03AQAAAADNRkFBgdq1a1ftfLfOsLVu3Vq+vr46c+aMS/uZM2cUGxtb5TLTp09XVlaWc9rhcOjLL79UZGSkbDabO5sHLK2kpETx8fEqKChQWFiYp8sBGhxjHE0dYxxNGePbeowxKi0tVVxcXI393ApsAQEB6t27t7Zu3aqUlBRJ3wSwrVu3auLEiVUuY7fbZbfbXdoiIiLc2SzgVcLCwvhBiCaNMY6mjjGOpozxbS3h4eG19nErsElSVlaW0tLSlJSUpLvvvlsvv/yyLl26pLFjx9arSAAAAABA1dwObI8++qiKioo0c+ZMnT59Wt/5znf0wQcfVHoQCQAAAADg5rgd2CRp4sSJ1V4CCTRXdrtds2bNqnQJMNBUMMbR1DHG0ZQxvr2XW0+JBAAAAAA0HrdenA0AAAAAaDwENgAAAACwKAIbAAAAAFgUgQ0AAAAALIrABrhp/vz56tOnj0JDQxUdHa2UlBQdOnTIpc+VK1eUkZGhyMhIhYSE6OGHH9aZM2c8VDFQfwsWLJDNZlNmZqazjfENb3fy5EmNGTNGkZGRCgoK0h133KFPPvnEOd8Yo5kzZ6pNmzYKCgrSkCFDdOTIEQ9WDNRNRUWFZsyYocTERAUFBalTp06aM2eObnzGIOPb+xDYADfl5OQoIyNDO3bs0ObNm3Xt2jXdf//9unTpkrPPlClT9P7772vt2rXKyclRYWGhRowY4cGqAffl5eVpxYoVuvPOO13aGd/wZhcuXFD//v3l7++vjRs3av/+/XrxxRfVsmVLZ59FixZpyZIlWr58uXbu3KkWLVpo6NChunLligcrB2q3cOFCLVu2TK+88ooOHDighQsXatGiRVq6dKmzD+PbCxkAN+Xs2bNGksnJyTHGGFNcXGz8/f3N2rVrnX0OHDhgJJnc3FxPlQm4pbS01HTu3Nls3rzZ3HfffWby5MnGGMY3vN9PfvITM2DAgGrnOxwOExsbaxYvXuxsKy4uNna73bz99tuNUSJQb8OGDTPjxo1zaRsxYoQZPXq0MYbx7a04wwbcpIsXL0qSWrVqJUnatWuXrl27piFDhjj7dOnSRQkJCcrNzfVIjYC7MjIyNGzYMJdxLDG+4f3ee+89JSUlaeTIkYqOjlavXr30q1/9yjk/Pz9fp0+fdhnj4eHh6tu3L2Mclnfvvfdq69atOnz4sCTps88+07Zt25ScnCyJ8e2t/DxdAODNHA6HMjMz1b9/f/Xo0UOSdPr0aQUEBCgiIsKlb0xMjE6fPu2BKgH3/Pa3v9Xu3buVl5dXaR7jG97uiy++0LJly5SVlaWf/vSnysvL07PPPquAgAClpaU5x3FMTIzLcoxxeINp06appKREXbp0ka+vryoqKjR37lyNHj1akhjfXorABtyEjIwM7du3T9u2bfN0KUCDKCgo0OTJk7V582YFBgZ6uhygwTkcDiUlJWnevHmSpF69emnfvn1avny50tLSPFwdcHN+97vfafXq1VqzZo26d++uPXv2KDMzU3FxcYxvL8YlkUA9TZw4UX/84x/10UcfqV27ds722NhYXb16VcXFxS79z5w5o9jY2EauEnDPrl27dPbsWd11113y8/OTn5+fcnJytGTJEvn5+SkmJobxDa/Wpk0bdevWzaWta9euOnHihCQ5x/G3n3zKGIc3mDp1qqZNm6bU1FTdcccdeuKJJzRlyhTNnz9fEuPbWxHYADcZYzRx4kStX79ef/nLX5SYmOgyv3fv3vL399fWrVudbYcOHdKJEyfUr1+/xi4XcMvgwYP1+eefa8+ePc6vpKQkjR492vlvxje8Wf/+/Su9iuXw4cNq3769JCkxMVGxsbEuY7ykpEQ7d+5kjMPyLl++LB8f14/3vr6+cjgckhjf3opLIgE3ZWRkaM2aNXr33XcVGhrqvOY7PDxcQUFBCg8P1/jx45WVlaVWrVopLCxMkyZNUr9+/XTPPfd4uHqgZqGhoc77Ma9r0aKFIiMjne2Mb3izKVOm6N5779W8efM0atQoffzxx1q5cqVWrlwpSc73Dr7wwgvq3LmzEhMTNWPGDMXFxSklJcWzxQO1GD58uObOnauEhAR1795dn376qV566SWNGzdOEuPba3n6MZWAt5FU5deqVaucfcrKykx6erpp2bKlCQ4ONg899JA5deqU54oGbsKNj/U3hvEN7/f++++bHj16GLvdbrp06WJWrlzpMt/hcJgZM2aYmJgYY7fbzeDBg82hQ4c8VC1QdyUlJWby5MkmISHBBAYGmo4dO5rnn3/elJeXO/swvr2PzZgbXn0OAAAAALAM7mEDAAAAAIsisAEAAACARRHYAAAAAMCiCGwAAAAAYFEENgAAAACwKAIbAAAAAFgUgQ0AAAAALIrABgAAAAAWRWADAAAAAIsisAEAAACARRHYAAAAAMCiCGwAAAAAYFEENgAAAACwKAIbAAAAAFgUgQ0AAAAALIrABgAAAAAWRWADAAAAAIsisAEAAACARRHYAAAAAMCiCGwAgCq98cYbstlsOnbsmKdLAQCg2SKwAYCHzZ07VzabTT169PB0KQ1m+/btys7OVnFxsadL8Vq3+hjyPQIA70BgAwAP+sc//qF58+apRYsWni6lkieeeEJlZWVq376928tu375ds2fPJgzchFt9DPkeAYB38PN0AQDQnP34xz/WPffco4qKCp07d87T5bjw9fWVr6+vp8tAA7t06ZIl/0AAAKgaZ9gAwEP++te/at26dXr55ZfrvEx2drZsNpsOHz6sMWPGKDw8XFFRUZoxY4aMMSooKNC//uu/KiwsTLGxsXrxxRcrreP48eNKT0/X7bffrqCgIEVGRmrkyJGV7lX79j1s17d99OhRPfXUU4qIiFB4eLjGjh2ry5cvu9Q4depUSVJiYqJsNpvLek6ePKlx48YpJiZGdrtd3bt31+uvv17tvta2vbruT13Xd93Jkyc1fvx4xcXFyW63KzExURMmTNDVq1dd+tRlX6pSWlqqzMxMdejQQXa7XdHR0frnf/5n7d69u9Zj6O4+79+/X48//rhatmypAQMG1Lj+muqqyb/8y7+oQ4cOldqNMbrrrrs0cODAOh0XAIArzrABgAdUVFRo0qRJ+uEPf6g77rjD7eUfffRRde3aVQsWLNCf/vQnvfDCC2rVqpVWrFih73//+1q4cKFWr16tH//4x+rTp4+++93vOpfNy8vT9u3blZqaqnbt2unYsWNatmyZBg0apP379ys4OLjGbY8aNUqJiYmaP3++du/erddee03R0dFauHChJGnEiBE6fPiw3n77bf385z9X69atJUlRUVE6c+aM7rnnHtlsNk2cOFFRUVHauHGjxo8fr5KSEmVmZrq9PXf3p7b1SVJhYaHuvvtuFRcX65lnnlGXLl108uRJrVu3TpcvX1ZAQEC99uVGP/rRj7Ru3TpNnDhR3bp10/nz57Vt2zYdOHCgxmNYn30eOXKkOnfurHnz5skYowEDBlS7/meeeabauu66665q96dPnz7auHGjLly4oJYtWzrbf/vb3+rTTz/Vtm3bajweAIBqGABAo3vllVdMeHi4OXv2rDHGmPvuu89079691uVmzZplJJlnnnnG2fb111+bdu3aGZvNZhYsWOBsv3DhggkKCjJpaWku67h8+XKl9ebm5hpJ5s0333S2rVq1ykgy+fn5LtseN26cy7IPPfSQiYyMdGlbvHixy7LXjR8/3rRp08acO3fOpT01NdWEh4e71FbX7dV1f9yp/8knnzQ+Pj4mLy+v0rodDofb+1KV8PBwk5GRUe386o6hMe7v82OPPVbn9ddWV3Xee+89I8ls3brV2Xb16lXTqVMnM3z4cLfXBwD4BpdEAkAjO3/+vGbOnKkZM2Y4z5i464c//KHz376+vkpKSpIxRuPHj3e2R0RE6Pbbb9cXX3zhsmxQUJDz39euXdP58+d12223KSIiotbL3qRvzgzdaODAgTp//rxKSkpqXM4Yo9///vcaPny4jDE6d+6c82vo0KG6ePFilduvbXvu7k9t63M4HNqwYYOGDx+upKSkSsvbbLZ678uNIiIitHPnThUWFtbYryo3u8+3oq4+ffpIksv2V65cqfz8fM2bN8+tdQEA/g+BDQAa2X/8x3+oVatWmjRpUr3XkZCQ4DIdHh6uwMBA56VtN7ZfuHDBpa2srEwzZ85UfHy87Ha7WrduraioKBUXF+vixYtub/v65W/f3s63FRUVqbi4WCtXrlRUVJTL19ixYyVJZ8+edXt77u5PbesrKipSSUlJja9ZqO++3GjRokXat2+f4uPjdffddys7O7tSuK6Ou/ucmJhYp/XeTF2xsbFq27atPv30U0nfPNxkzpw5GjNmjMuxNMYoJCSk1uMDAPgG97ABQCM6cuSIVq5cqZdfftnlDMaVK1d07do1HTt2TGFhYWrVqlWN66nq6Y3VPdHRGOMyPWnSJK1atUqZmZnq16+fwsPDZbPZlJqaKofDUes+1HU733Z93WPGjFFaWlqVfe688063t+fu/tS3/hvVd19uNGrUKA0cOFDr16/Xpk2btHjxYi1cuFB/+MMflJycXOOy7u7zjWfkanMzdfXp08cZ2F566SVduHBBP/vZz1z65OfnKzg4WNHR0XWuCQCaMwIbADSikydPyuFw6Nlnn9Wzzz5baX5iYqImT57s1pMj3bVu3TqlpaW5PEHyypUrDfo+LpvNVqktKipKoaGhqqio0JAhQxpsWw29P1FRUQoLC9O+fftq7NMQ+9KmTRulp6crPT1dZ8+e1V133aW5c+cqOTm5ymN4XUPsc03rr6mumvTp00fvvfeeTpw4of/8z//UhAkTXN7jd+DAAfXq1UsVFRUKCQlR165dlZeXV+eaAaA5IrABQCPq0aOH1q9fX6n9P/7jP1RaWqpf/OIX6tSp0y2twdfXt9LZpKVLl6qioqLBtnH9PV83BghfX189/PDDWrNmjfbt21fpksOioqJ63dPX0Pvj4+OjlJQUvfXWW/rkk08q3cdmjLnpfamoqNBXX32l8PBwZ1t0dLTi4uJUXl4uqepjeF1D7HNV669LXTVJSkqSw+HQ448/LmOMnn/+eZf5Xbt21axZs3Tq1CktWbKkzrUCQHNGYAOARtS6dWulpKRUar9+Rq2qeQ3tBz/4gX7zm98oPDxc3bp1U25urrZs2aLIyMgG20bv3r0lSc8//7xSU1Pl7++v4cOHa8GCBfroo4/Ut29fPf300+rWrZu+/PJL7d69W1u2bNGXX35pif2ZN2+eNm3apPvuu0/PPPOMunbtqlOnTmnt2rXatm2bIiIibmpfSktL1a5dOz3yyCPq2bOnQkJCtGXLFuXl5TnPmlV3DFu0aNEg+1zV+gcOHKjbb7+9xrpqcj3c/u///q+ys7OrDK179+7V97///TrXCQDNHYENAJqZX/ziF/L19dXq1at15coV9e/fX1u2bNHQoUMbbBt9+vTRnDlztHz5cn3wwQdyOBzKz89Xhw4d9PHHH+tnP/uZ/vCHP+iXv/ylIiMj1b17d5f3oHl6f9q2baudO3dqxowZWr16tUpKStS2bVslJyc733EWExNT730JDg5Wenq6Nm3apD/84Q9yOBy67bbb9Mtf/lITJkyQVP0xbNGiRYPsc1XrP3ToUK111aRVq1bq0KGDLl26pH//93+vss/evXtrfUcdAOD/2Iw7d1kDAABU44svvtA//dM/6aWXXqryHs2rV68qJCREFy5ccF6SCQCoGYENAAA0iEcffVS7du3S/v37FRAQUGn++fPn1aZNG505c8b5OgUAQM24JBIAANRbcXGxNm7cqP/+7//W2rVrtXHjxirDmiRFRkbqscceU0JCgrp3764dO3Y0crUA4H04wwYAAOrt97//vR555BG1a9dOM2fO1NNPP+3pkgCgSSGwAQAAAIBF+Xi6AAAAAABA1QhsAAAAAGBRbj905K9//asWL16sXbt26dSpU1q/fr1bL3p1OBwqLCxUaGiobDabu5sHAAAAAK9njFFpaani4uLk41P9eTS3A9ulS5fUs2dPjRs3TiNGjHC7sMLCQsXHx7u9HAAAAAA0NQUFBWrXrl21890ObMnJyUpOTq53QaGhoc7CwsLC6r0eAAAAAPBWJSUlio+Pd+aj6tzy97CVl5ervLzcOV1aWipJCgsLI7ABAAAAaNZqu03slge2+fPna/bs2bd6MwAAAC46TPtTtfOOLRjWaNu6FduDNTHmcCvc8qdETp8+XRcvXnR+FRQU3OpNAgAAAECTcMvPsNntdtnt9lu9GQAAAABocngPGwAAAABYlNtn2L766isdPXrUOZ2fn689e/aoVatWSkhIaNDiAAAAAKA5czuwffLJJ/re977nnM7KypIkpaWl6Y033miwwgAAAACguXM7sA0aNEjGmFtRCwAAAADgBtzDBgAAAAAWRWADAAAAAIsisAEAAACARRHYAAAAAMCiCGwAAAAAYFEENgAAAACwKAIbAAAAAFgUgQ0AAAAALIrABgAAAAAWRWADAAAAAIsisAEAAACARRHYAAAAAMCiCGwAAAAAYFEENgAAAACwKAIbAAAAAFgUgQ0AAAAALIrABgAAAAAWRWADAAAAAIsisAEAAACARRHYAAAAAMCiCGwAAAAAYFEENgAAAACwKAIbAAAAAFgUgQ0AAAAALIrABgAAAAAWRWADAAAAAIsisAEAAACARRHYAAAAAMCiCGwAAAAAYFEENgAAAACwKAIbAAAAAFgUgQ0AAAAALIrABgAAAAAWRWADAAAAAIsisAEAAACARRHYAAAAAMCiCGwAAAAAYFEENgAAAACwKAIbAAAAAFgUgQ0AAAAALIrABgAAAAAWRWADAAAAAIsisAEAAACARRHYAAAAAMCiCGwAAAAAYFEENgAAAACwKAIbAAAAAFgUgQ0AAAAALIrABgAAAAAWRWADAAAAAIsisAEAAACARdUrsL366qvq0KGDAgMD1bdvX3388ccNXRcAAAAANHtuB7Z33nlHWVlZmjVrlnbv3q2ePXtq6NChOnv27K2oDwAAAACaLbcD20svvaSnn35aY8eOVbdu3bR8+XIFBwfr9ddfvxX1AQAAAECz5edO56tXr2rXrl2aPn26s83Hx0dDhgxRbm5ulcuUl5ervLzcOX3x4kVJUklJSX3qBQAAqBNH+eVq5zX055CatnUrtgdrYszBHde/R8aYGvu5FdjOnTuniooKxcTEuLTHxMTo4MGDVS4zf/58zZ49u1J7fHy8O5sGAABoMOEvN+3twXoYc6hOaWmpwsPDq53vVmCrj+nTpysrK8s57XA49OWXXyoyMlI2m+1Wbx5oNCUlJYqPj1dBQYHCwsI8XQ7Q4BjjaOoY42jKGN/WY4xRaWmp4uLiauznVmBr3bq1fH19debMGZf2M2fOKDY2tspl7Ha77Ha7S1tERIQ7mwW8SlhYGD8I0aQxxtHUMcbRlDG+raWmM2vXufXQkYCAAPXu3Vtbt251tjkcDm3dulX9+vVzv0IAAAAAQLXcviQyKytLaWlpSkpK0t13362XX35Zly5d0tixY29FfQAAAADQbLkd2B599FEVFRVp5syZOn36tL7zne/ogw8+qPQgEqC5sdvtmjVrVqVLgIGmgjGOpo4xjqaM8e29bKa250gCAAAAADzC7RdnAwAAAAAaB4ENAAAAACyKwAYAAAAAFkVgAwAAAACLIrABbpo/f7769Omj0NBQRUdHKyUlRYcOHXLpc+XKFWVkZCgyMlIhISF6+OGHK71wHvAGCxYskM1mU2ZmprON8Q1vd/LkSY0ZM0aRkZEKCgrSHXfcoU8++cQ53xijmTNnqk2bNgoKCtKQIUN05MgRD1YM1E1FRYVmzJihxMREBQUFqVOnTpozZ45ufMYg49v7ENgAN+Xk5CgjI0M7duzQ5s2bde3aNd1///26dOmSs8+UKVP0/vvva+3atcrJyVFhYaFGjBjhwaoB9+Xl5WnFihW68847XdoZ3/BmFy5cUP/+/eXv76+NGzdq//79evHFF9WyZUtnn0WLFmnJkiVavny5du7cqRYtWmjo0KG6cuWKBysHardw4UItW7ZMr7zyig4cOKCFCxdq0aJFWrp0qbMP49sLGQA35ezZs0aSycnJMcYYU1xcbPz9/c3atWudfQ4cOGAkmdzcXE+VCbiltLTUdO7c2WzevNncd999ZvLkycYYxje8309+8hMzYMCAauc7HA4TGxtrFi9e7GwrLi42drvdvP32241RIlBvw4YNM+PGjXNpGzFihBk9erQxhvHtrTjDBtykixcvSpJatWolSdq1a5euXbumIUOGOPt06dJFCQkJys3N9UiNgLsyMjI0bNgwl3EsMb7h/d577z0lJSVp5MiRio6OVq9evfSrX/3KOT8/P1+nT592GePh4eHq27cvYxyWd++992rr1q06fPiwJOmzzz7Ttm3blJycLInx7a38PF0A4M0cDocyMzPVv39/9ejRQ5J0+vRpBQQEKCIiwqVvTEyMTp8+7YEqAff89re/1e7du5WXl1dpHuMb3u6LL77QsmXLlJWVpZ/+9KfKy8vTs88+q4CAAKWlpTnHcUxMjMtyjHF4g2nTpqmkpERdunSRr6+vKioqNHfuXI0ePVqSGN9eisAG3ISMjAzt27dP27Zt83QpQIMoKCjQ5MmTtXnzZgUGBnq6HKDBORwOJSUlad68eZKkXr16ad++fVq+fLnS0tI8XB1wc373u99p9erVWrNmjbp37649e/YoMzNTcXFxjG8vxiWRQD1NnDhRf/zjH/XRRx+pXbt2zvbY2FhdvXpVxcXFLv3PnDmj2NjYRq4ScM+uXbt09uxZ3XXXXfLz85Ofn59ycnK0ZMkS+fn5KSYmhvENr9amTRt169bNpa1r1646ceKEJDnH8beffMoYhzeYOnWqpk2bptTUVN1xxx164oknNGXKFM2fP18S49tbEdgANxljNHHiRK1fv15/+ctflJiY6DK/d+/e8vf319atW51thw4d0okTJ9SvX7/GLhdwy+DBg/X5559rz549zq+kpCSNHj3a+W/GN7xZ//79K72K5fDhw2rfvr0kKTExUbGxsS5jvKSkRDt37mSMw/IuX74sHx/Xj/e+vr5yOBySGN/eiksiATdlZGRozZo1evfddxUaGuq85js8PFxBQUEKDw/X+PHjlZWVpVatWiksLEyTJk1Sv379dM8993i4eqBmoaGhzvsxr2vRooUiIyOd7YxveLMpU6bo3nvv1bx58zRq1Ch9/PHHWrlypVauXClJzvcOvvDCC+rcubMSExM1Y8YMxcXFKSUlxbPFA7UYPny45s6dq4SEBHXv3l2ffvqpXnrpJY0bN04S49trefoxlYC3kVTl16pVq5x9ysrKTHp6umnZsqUJDg42Dz30kDl16pTnigZuwo2P9TeG8Q3v9/7775sePXoYu91uunTpYlauXOky3+FwmBkzZpiYmBhjt9vN4MGDzaFDhzxULVB3JSUlZvLkySYhIcEEBgaajh07mueff96Ul5c7+zC+vY/NmBtefQ4AAAAAsAzuYQMAAAAAiyKwAQAAAIBFEdgAAAAAwKIIbAAAAABgUQQ2AAAAALAoAhsAAAAAWBSBDQAAAAAsisAGAAAAABZFYAMAAAAAiyKwAQAAAIBFEdgAAAAAwKIIbAAAAABgUf8PO07uJyNDaFUAAAAASUVORK5CYII=",
                        "text/plain": [
                            "<Figure size 900x400 with 3 Axes>"
                        ]
                    },
                    "metadata": {},
                    "output_type": "display_data"
                }
            ],
            "source": [
                "class MPRamp(MP):\n",
                "    def __init__(self, c, M, P, upos_max, uneg_max):\n",
                "        super(MPRamp, self).__init__(c, M, P)\n",
                "        super(MPRamp, self).set_modfile(\"maintenance_planning_ramp.mod\")\n",
                "        self.upos_max = upos_max\n",
                "        self.uneg_max = uneg_max\n",
                "\n",
                "    def load_data(self):\n",
                "        super(MPRamp, self).load_data()\n",
                "        self.ampl.param[\"upos_max\"] = self.upos_max\n",
                "        self.ampl.param[\"uneg_max\"] = self.uneg_max\n",
                "\n",
                "\n",
                "upos_max = 0.3334\n",
                "uneg_max = 0.5000\n",
                "\n",
                "model = MPRamp(c, M, P, upos_max, uneg_max)\n",
                "model.plot_schedule()"
            ]
        },
        {
            "cell_type": "markdown",
            "metadata": {
                "colab_type": "text",
                "id": "HM_k7HkbThAe"
            },
            "source": [
                "## Specifying the minimum number of operational days between maintenance periods\n",
                "\n",
                "Up to this point we have imposed no constraints on the frequency of maintenance periods. Without such constraints, particularly when ramping constraints are imposed, is that maintenance periods will be scheduled back-to-back, which is clearly not a useful result for most situations.\n",
                "\n",
                "The next revision of the model is to incorporate a requirement that $N$ operational days be scheduled between any maintenance periods. This does allow for maintenance to be postponed until the very end of the planning period. The disjunctive constraints read\n",
                "\n",
                "$$\n",
                "\\begin{align*}\n",
                "\\left(y_t = 0\\right) \\veebar \\left(\\sum_{s=0}^{(M + N -1) \\wedge (t + s \\leq T)}x_{t+s} = 0\\right)\\qquad \\forall t = 1, 2, \\ldots, T-M+1\n",
                "\\end{align*}\n",
                "$$\n",
                "\n",
                "where the upper bound on the summation is needed to handle the terminal condition. \n",
                "\n",
                "Paradoxically, this is an example where the Big-M method provides a much faster solution.\n",
                "\n",
                "$$\n",
                "\\begin{align*}\n",
                "\\sum_{s=0}^{(M + N -1) \\wedge (t + s \\leq T)}x_{t+s} \\leq (M+N)(1-y_t) \\qquad \\forall t = 1, 2, \\ldots, T-M+1\n",
                "\\end{align*}\n",
                "$$\n",
                "\n",
                "The following cell implements both sets of constraints. "
            ]
        },
        {
            "cell_type": "code",
            "execution_count": 8,
            "metadata": {
                "colab": {},
                "colab_type": "code",
                "id": "fEOMGHPvUJzJ"
            },
            "outputs": [
                {
                    "name": "stdout",
                    "output_type": "stream",
                    "text": [
                        "Writing maintenance_planning_ramp_operational.mod\n"
                    ]
                }
            ],
            "source": [
                "%%writefile maintenance_planning_ramp_operational.mod\n",
                "\n",
                "param T; # number of planning periods\n",
                "param P; # number of maintenance periods\n",
                "param M; # number of days for each maintenance period\n",
                "param N;\n",
                "\n",
                "param upos_max;\n",
                "param uneg_max;\n",
                "\n",
                "set PH := 1..T; # set of planning horizon time periods\n",
                "\n",
                "set Y := 1..(T - M + 1);\n",
                "set S := 0..(M - 1);\n",
                "set W := 0..(M + N - 1);\n",
                "    \n",
                "param c{PH}; # profit\n",
                "\n",
                "var x{PH} >= 0, <= 1;\n",
                "var y{PH} binary;\n",
                "var upos{PH} >= 0, <= upos_max;\n",
                "var uneg{PH} >= 0, <= uneg_max;\n",
                "\n",
                "maximize profit: sum{t in PH} c[t] * x[t];\n",
                "\n",
                "# ramp constraint\n",
                "s.t. ramp {t in PH: t > 1}:\n",
                "    x[t] == x[t-1] + upos[t] - uneg[t];\n",
                "# required number P of maintenance starts\n",
                "s.t. sumy: P == sum{t in Y} y[t];\n",
                "# no more than one maintenance start in the period of length M\n",
                "s.t. sprd {t in Y}:\n",
                "    sum{s in W: t + s <= T} y[t+s] <= 1;\n",
                "\n",
                "# Choose one of the following constraints (disj or bigm).\n",
                "# Comment out the constraint not used.\n",
                "#s.t. disj {t in Y}:\n",
                "#    y[t] == 1 ==> sum{s in W: t + s <= T} x[t+s] == 0;\n",
                "\n",
                "s.t. bigm {t in Y}:\n",
                "    sum{s in W: t + s <= T} x[t+s] <= (M+N)*(1 - y[t]);"
            ]
        },
        {
            "cell_type": "code",
            "execution_count": 9,
            "metadata": {
                "colab": {
                    "base_uri": "https://localhost:8080/",
                    "height": 297
                },
                "colab_type": "code",
                "executionInfo": {
                    "elapsed": 3249,
                    "status": "ok",
                    "timestamp": 1557947681109,
                    "user": {
                        "displayName": "Jeffrey Kantor",
                        "photoUrl": "https://lh5.googleusercontent.com/-8zK5aAW5RMQ/AAAAAAAAAAI/AAAAAAAAKB0/kssUQyz8DTQ/s64/photo.jpg",
                        "userId": "09038942003589296665"
                    },
                    "user_tz": 240
                },
                "id": "m93Na6TkUqDb",
                "outputId": "cae87fae-7eba-46c3-ea3b-7e2427c4f529"
            },
            "outputs": [
                {
                    "name": "stdout",
                    "output_type": "stream",
                    "text": [
                        "cbc 2.10.7: \b\b\b\b\b\b\b\b\b\b\b\bcbc 2.10.7: optimal solution; objective 26.64390078\n",
                        "90873 simplex iterations\n",
                        "90873 barrier iterations\n",
                        "492 branching nodes\n"
                    ]
                },
                {
                    "data": {
                        "image/png": "iVBORw0KGgoAAAANSUhEUgAAA2wAAAGGCAYAAAAHLAKcAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjcuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/bCgiHAAAACXBIWXMAAA9hAAAPYQGoP6dpAABDUUlEQVR4nO3de3xU9Z3/8ffkNknIDUIuBBIIyMpNKRJEBCotrJiluBEFo6ARqG5JQEK2WKgLhCL31VrQcqkr1grWQgtqW5RL3XRZAkYQkXKvAVLCJSAhEULAzPf3hz9mGXOdEDJnktfz8cjjwXzP95zzOWe+TOadc7MZY4wAAAAAAJbj4+kCAAAAAABVI7ABAAAAgEUR2AAAAADAoghsAAAAAGBRBDYAAAAAsCgCGwAAAABYFIENAAAAACyKwAYAAAAAFkVgAwAAAACLIrABAAAAgEUR2AAAAADAoghsAIB6yc7Ols1mc3u+N954QzabTceOHauxzZvk5eXp3nvvVYsWLWSz2bRnzx6v3yYAgDUQ2AAAuAnXrl3TyJEj9eWXX+rnP/+5fvOb36h9+/ZV9t2+fbuys7NVXFzcuEUCALwWgQ0A0KieeOIJlZWVVRtqvM3f//53HT9+XD/+8Y/1zDPPaMyYMWrZsmWV27l9+3bNnj2bwAYAqDM/TxcAAGhefH195evr6+kyqnXp0iW1aNGizv3Pnj0rSYqIiHBpt/p2AgC8A0fYAAC12rZtm/r06aPAwEB16tRJK1asqNTn+PHjSk9P1+23366goCBFRkZq5MiRla7hqu3aro8++kg2m03r16+vNG3NmjWy2WzKzc2tsd7r19cdPHhQo0aNUlhYmCIjIzV58mRduXKlUr/9+/fr8ccfV8uWLTVgwADn9E8//VTJyckKCwtTSEiIBg8erB07djinP/XUU7rvvvskSSNHjpTNZtOgQYOq3M7s7GxNnTpVkpSYmCibzVana9x27dqlRx55RNHR0QoMDFSXLl30s5/9rMZ5AABNB0fYAAA1+vzzz3X//fcrKipK2dnZ+vrrrzVr1izFxMS49MvLy9P27duVmpqqdu3a6dixY1q2bJkGDRqk/fv3Kzg4uE7rGzRokOLj47V69Wo99NBDLtNWr16tTp06qV+/fnVa1qhRo9ShQwfNnz9fO3bs0JIlS3ThwgW9+eabLv1Gjhypzp07a968eTLGSJL+9re/aeDAgQoLC9Nzzz0nf39/rVixQoMGDVJOTo769u2rf/u3f1Pbtm01b948Pfvss+rTp0+l/XLdiBEjdPjwYb399tv6+c9/rtatW0uSoqKiqq1//fr1Sk1NVceOHTV16lSFhIQ49zMAoJkwAADUICUlxQQGBprjx4872/bv3298fX3Njb9GLl++XGne3NxcI8m8+eabzrZVq1YZSSY/P7/atunTpxu73W6Ki4udfc6ePWv8/PzMrFmzaq151qxZRpJ58MEHXdrT09ONJPPZZ5+59Hvssceq3O6AgADz97//3dlWWFhoQkNDzXe/+11n20cffWQkmbVr17rMX9V2Ll68uFJbdQ4dOmSCg4PNQw89ZMrKylymXblypdb5AQBNA6dEAgCqVVFRoQ8//FApKSlKSEhwtnft2lVDhw516RsUFOT897Vr13T+/HnddtttioiI0O7du91a75NPPqny8nKtW7fO2fbOO+/o66+/1pgxY+q8nIyMDJfXkyZNkiT9+c9/dmn/0Y9+5PK6oqJCmzZtUkpKijp27Ohsb9OmjR5//HFt27ZNJSUlda6jPmbOnKnAwEC98cYbCgwMdJlmt9tv6boBANZBYAMAVKuoqEhlZWXq3LlzpWm33367y+uysjLNnDlT8fHxstvtat26taKiolRcXKyLFy+6td4uXbqoT58+Wr16tbNt9erVuueee3TbbbfVeTnfrrtTp07y8fGpdN1YYmKiy+uioiJdvny50jZK34RVh8OhgoKCOtfhrvLycr3//vsaM2aMwsLCauxrjFFISIjz5icAgKaFa9gAAA1i0qRJWrVqlTIzM9WvXz+Fh4fLZrMpNTVVDofD7eU9+eSTmjx5sv7xj3+ovLxcO3bs0CuvvHJTNVb3oO8bjw5awRdffKHLly+rd+/etfbNz89XcHCwoqOjG6EyAEBjI7ABAKoVFRWloKAgHTlypNK0Q4cOubxet26d0tLS9OKLLzrbrly5Uu9njqWmpiorK0tvv/22ysrK5O/vr0cffdStZRw5csTl6NnRo0flcDjUoUOHGueLiopScHBwpW2UpIMHD8rHx0fx8fFu1SJVHxi/raysrE79Dxw4oF69eqmiokIhISHq2rWr8vLy3K4LAGBdnBIJAKiWr6+vhg4dqg0bNujEiRPO9gMHDujDDz+s1Nf8/zssXrd06VJVVFTUa92tW7dWcnKy3nrrLa1evVoPPPCA886KdfXqq69WqkeSkpOTa5zP19dX999/v959912X0yfPnDmjNWvWaMCAAbWeqliV6893qy3EXj/tc8uWLZWmXbt2zfnvrl27atasWZowYYK++uorwhoANEEcYQMA1Gj27Nn64IMPNHDgQKWnp+vrr7/W0qVL1b17d+3du9fZ7wc/+IF+85vfKDw8XN26dVNubq62bNmiyMjIeq/7ySef1COPPCJJmjNnjtvz5+fn68EHH9QDDzyg3NxcvfXWW3r88cfVs2fPWud94YUXtHnzZg0YMEDp6eny8/PTihUrVF5erkWLFrldiyTnKY7PP/+8UlNT5e/vr+HDh1d6UHdYWJieeuopvfHGGyovL9egQYNUWlqqjz76SMOGDXO5mcrevXv1/e9/v171AACsj8AGAKjRnXfeqQ8//FBZWVmaOXOm2rVrp9mzZ+vUqVMuge0Xv/iFfH19tXr1al25ckX9+/fXli1bKt1N0h3Dhw9Xy5Yt5XA49OCDD7o9/zvvvKOZM2dq2rRp8vPz08SJE7V48eI6zdu9e3f9z//8j6ZPn6758+fL4XCob9++euutt9S3b1+3a5GkPn36aM6cOVq+fLk++OADORwO5efnVwps0jdHB9u2bavf/e532rBhg1q1aqV7771XDzzwgEu/vXv3KjMzs171AACsz2a+ff4KAAAW8fXXXysuLk7Dhw/Xf/3Xf9V5vuzsbM2ePVtFRUVun0bpTa5evaqQkBBduHChytAHAPB+XMMGALCsDRs2qKioSE8++aSnS7Gk0tJSSd8ENwBA08QpkQAAy9m5c6f27t2rOXPmqFevXrrvvvs8XZIlRUZG6rHHHlNCQoK6d++uHTt2eLokAEAD4wgbAMByli1bpgkTJig6Olpvvvmmp8uxtF//+tcqLS0lrAFAE8U1bAAAAABgURxhAwAAAACLcjuw/fWvf9Xw4cMVFxcnm82mDRs23IKyAAAAAABu33Tk0qVL6tmzp8aNG6cRI0a4vUKHw6HCwkKFhobKZrO5PT8AAAAAeDtjjEpLSxUXFycfn+qPo7kd2JKTk5WcnFzvwgoLCxUfH1/v+QEAAACgqSgoKFC7du2qnd7ot/UPDQ2V9E1hYWFhjb16AAAAAPC4kpISxcfHO/NRdW55YCsvL1d5ebnz9fWHfIaFhRHYAAAAADRrtV0mdssD2/z58zV79uxbvRoAAOCmDtP+VOP0YwuGNVIlAIDq3PLANn36dGVlZTlfXz/0BzQnNX0p4gsRAAD/pyH/kMDvXzQFtzyw2e122e32W72aZoEPHQAAAKB5cTuwffXVVzp69KjzdX5+vvbs2aNWrVopISGhQYsD0DRw2pU18UcgAACsz+3A9sknn+h73/ue8/X10x3T0tL0xhtvNFhhaDr4sg4AAIDm4Fb8MdTtwDZo0CAZY+q1MgAAAABA3VX/SG0AAAAAgEc1+oOzUTWuJQEAAADwbQQ24CYRtgEAAHCrcEokAAAAAFgUR9jQLHHnSgAAGhZnnAC3BoENAACgGSNoAdZGYAMsgl+YAAAA+DYCWzNFOAAAAACsj8AGVMObr3MjkNfOm99fAADQfHCXSAAAAACwKI6wAQAAADXgrAx4EkfYAAAAAMCivOYIG3/ZAAAAANDceE1gAwA0fdwwBwCaNg7CuI/ABgCoVlP/xUpABABYHYENANDkNPWgCe/FHwkAuIvAdhP4QoCmji8WABobv1sBwBV3iQQAAAAAi+IIGwBYRF2PLNTlyCdHRwEAaBqaXGDjSwoAAADqitNwYXVNLrABAAAAQEPz1IEhAhtuCkc0Af4foPEx5gA0dXzO/R9LBLbGPhTNAAAAAADgDSwR2GBNBFsAAADAswhsACyDPxIAAICG1BRuKsNz2AAAAADAoghsAAAAAGBRnBIJAACaJE6zbtp4f9FcENgAL9IUzsP2NnwhAAAAnkRgawR84QMAPgsbG/sbAJoGAls1OJIBAAAAwNMIbAC8CkcNAABofPz+9RwCWxPDkUEAAACg6SCwAQAANEH8ERdoGghsAICbwpfChsVpRwCAG/HgbAAAAACwKAIbAAAAAFgUp0QCANAIOHUUgFVxKra1EdgAAAAAL8MfgZoPAhsAAADQSDiaBXcR2AA0OfzVEQCApq+5hF8CGwAA8CoN+UeZ5vKFD4D34i6RAAAAAGBRHGEDAKAGHIGpG/YTANwaBDYAAAAvQ0AGmg9OiQQAAAAAi+IIGwAAaBQcFQIA9xHYAAAAasCjQgB4EoENwE3hiwwANCyORAK4EYENaKb4QgA0HP5wAaAhNeTvaH7fez8CGwAAqBFf+GpHaAdwqxDY4FX40gAAQONrqN+/BFtYkdXHZb0C26uvvqrFixfr9OnT6tmzp5YuXaq77767oWtDM2L1/ygAAACAJ7gd2N555x1lZWVp+fLl6tu3r15++WUNHTpUhw4dUnR09K2oEQAAJ/7AAwDWxOfzreF2YHvppZf09NNPa+zYsZKk5cuX609/+pNef/11TZs2rcELBAAAALwBl27gVnArsF29elW7du3S9OnTnW0+Pj4aMmSIcnNzG7w4oD74sIQVMS4BAEB9uBXYzp07p4qKCsXExLi0x8TE6ODBg1XOU15ervLycufrixcvSpJKSkqcbY7yyzWut6SkpE59altWQ/WxYk3u1G3FmqxUtxVr4v1t2jXx/lJTY9fEmGvaNfH+Nu2aeH+9v6ZvvzbG1Lh+m6mtxw0KCwvVtm1bbd++Xf369XO2P/fcc8rJydHOnTsrzZOdna3Zs2fXdRUAAAAA0GwUFBSoXbt21U536whb69at5evrqzNnzri0nzlzRrGxsVXOM336dGVlZTlfOxwOffnll4qMjJTNZnNn9YCllZSUKD4+XgUFBQoLC/N0OUCDY4yjqWOMoyljfFuPMUalpaWKi4ursZ9bgS0gIEC9e/fW1q1blZKSIumbALZ161ZNnDixynnsdrvsdrtLW0REhDurBbxKWFgYH4Ro0hjjaOoY42jKGN/WEh4eXmsft+8SmZWVpbS0NCUlJenuu+/Wyy+/rEuXLjnvGgkAAAAAaBhuB7ZHH31URUVFmjlzpk6fPq3vfOc7+uCDDyrdiAQAAAAAcHPcDmySNHHixGpPgQSaK7vdrlmzZlU6BRhoKhjjaOoY42jKGN/ey627RAIAAAAAGo+PpwsAAAAAAFSNwAYAAAAAFkVgAwAAAACLIrABbpo/f7769Omj0NBQRUdHKyUlRYcOHXLpc+XKFWVkZCgyMlIhISF6+OGHKz1wHvAGCxYskM1mU2ZmprON8Q1vd/LkSY0ZM0aRkZEKCgrSHXfcoU8++cQ53RijmTNnqk2bNgoKCtKQIUN05MgRD1YM1E1FRYVmzJihxMREBQUFqVOnTpozZ45uvGUF49v7ENgAN+Xk5CgjI0M7duzQ5s2bde3aNd1///26dOmSs8+UKVP0/vvva+3atcrJyVFhYaFGjBjhwaoB9+Xl5WnFihW68847XdoZ3/BmFy5cUP/+/eXv76+NGzdq//79evHFF9WyZUtnn0WLFmnJkiVavny5du7cqRYtWmjo0KG6cuWKBysHardw4UItW7ZMr7zyig4cOKCFCxdq0aJFWrp0qbMP49sLGQA35ezZs0aSycnJMcYYU1xcbPz9/c3atWudfQ4cOGAkmdzcXE+VCbiltLTUdO7c2WzevNncd999ZvLkycYYxje8309+8hMzYMCAaqc7HA4TGxtrFi9e7GwrLi42drvdvP32241RIlBvw4YNM+PGjXNpGzFihBk9erQxhvHtrTjCBtykixcvSpJatWolSdq1a5euXbumIUOGOPt06dJFCQkJys3N9UiNgLsyMjI0bNgwl3EsMb7h/d577z0lJSVp5MiRio6OVq9evfSrX/3KOT0/P1+nT592GePh4eHq27cvYxyWd++992rr1q06fPiwJOmzzz7Ttm3blJycLInx7a3q9eBsAN9wOBzKzMxU//791aNHD0nS6dOnFRAQoIiICJe+MTExOn36tAeqBNzz29/+Vrt371ZeXl6laYxveLsvvvhCy5YtU1ZWln76058qLy9Pzz77rAICApSWluYcxzExMS7zMcbhDaZNm6aSkhJ16dJFvr6+qqio0Ny5czV69GhJYnx7KQIbcBMyMjK0b98+bdu2zdOlAA2ioKBAkydP1ubNmxUYGOjpcoAG53A4lJSUpHnz5kmSevXqpX379mn58uVKS0vzcHXAzfnd736n1atXa82aNerevbv27NmjzMxMxcXFMb69GKdEAvU0ceJE/fGPf9RHH32kdu3aOdtjY2N19epVFRcXu/Q/c+aMYmNjG7lKwD27du3S2bNnddddd8nPz09+fn7KycnRkiVL5Ofnp5iYGMY3vFqbNm3UrVs3l7auXbvqxIkTkuQcx9++8yljHN5g6tSpmjZtmlJTU3XHHXfoiSee0JQpUzR//nxJjG9vRWAD3GSM0cSJE7V+/Xr95S9/UWJiosv03r17y9/fX1u3bnW2HTp0SCdOnFC/fv0au1zALYMHD9bnn3+uPXv2OH+SkpI0evRo578Z3/Bm/fv3r/QolsOHD6t9+/aSpMTERMXGxrqM8ZKSEu3cuZMxDsu7fPmyfHxcv977+vrK4XBIYnx7K06JBNyUkZGhNWvW6N1331VoaKjznO/w8HAFBQUpPDxc48ePV1ZWllq1aqWwsDBNmjRJ/fr10z333OPh6oGahYaGOq/HvK5FixaKjIx0tjO+4c2mTJmie++9V/PmzdOoUaP08ccfa+XKlVq5cqUkOZ87+MILL6hz585KTEzUjBkzFBcXp5SUFM8WD9Ri+PDhmjt3rhISEtS9e3d9+umneumllzRu3DhJjG+v5enbVALeRlKVP6tWrXL2KSsrM+np6aZly5YmODjYPPTQQ+bUqVOeKxq4CTfe1t8Yxje83/vvv2969Ohh7Ha76dKli1m5cqXLdIfDYWbMmGFiYmKM3W43gwcPNocOHfJQtUDdlZSUmMmTJ5uEhAQTGBhoOnbsaJ5//nlTXl7u7MP49j42Y2549DkAAAAAwDK4hg0AAAAALIrABgAAAAAWRWADAAAAAIsisAEAAACARRHYAAAAAMCiCGwAAAAAYFEENgAAAACwKAIbAAAAAFgUgQ0AAAAALIrABgAAAAAWRWADAAAAAIsisAEAAACARRHYAAAAAMCiCGwAAAAAYFEENgAAAACwKAIbAAAAAFgUgQ0AAAAALIrABgAAAAAWRWADAAAAAIsisAGAhb3xxhuy2Ww6duyYp0uxlKawX7Kzs2Wz2XTu3LlGW1d9NIV9DQDejMAGAF5k+/btys7OVnFxsadLaRTNbXsBAPg2AhsAWNgTTzyhsrIytW/fXtI3AWb27NnNJsBUt73f3i8AADRVfp4uAABQPV9fX/n6+nq6jAZz6dIltWjR4qaX09T2CwAA1eEIGwA0kKeeekodOnSo1P7t64euvz569KieeuopRUREKDw8XGPHjtXly5dd5r3x+qHs7GxNnTpVkpSYmCibzVana4s+/fRTJScnKywsTCEhIRo8eLB27NhRZY0HDx7UqFGjFBYWpsjISE2ePFlXrlyptMyTJ09q3LhxiomJkd1uV/fu3fX6669Xucz9+/fr8ccfV8uWLTVgwABJ0vHjx5Wenq7bb79dQUFBioyM1MiRI122pabtreq6Knf2qyT993//t5KSkhQYGKhOnTppxYoVdbrWq7S0VJmZmerQoYPsdruio6P1z//8z9q9e3elfTR+/HjFxcXJbrcrMTFREyZM0NWrV136FRcX16neuuxzSdq2bZv69Onjsl1Vqet4rU5d66lqvsDAQI0bN86lfcuWLfL399eUKVNqXQYANCccYQMADxk1apQSExM1f/587d69W6+99pqio6O1cOHCKvuPGDFChw8f1ttvv62f//znat26tSQpKiqq2nX87W9/08CBAxUWFqbnnntO/v7+WrFihQYNGqScnBz17du3Uk0dOnTQ/PnztWPHDi1ZskQXLlzQm2++6exz5swZ3XPPPbLZbJo4caKioqK0ceNGjR8/XiUlJcrMzHRZ5siRI9W5c2fNmzdPxhhJUl5enrZv367U1FS1a9dOx44d07JlyzRo0CDt379fwcHB9dreuu7XTz/9VA888IDatGmj2bNnq6KiQj/72c9qXbYk/ehHP9K6des0ceJEdevWTefPn9e2bdt04MAB3XXXXZKkwsJC3X333SouLtYzzzyjLl266OTJk1q3bp0uX76sgIAAt+qt6z7//PPPdf/99ysqKkrZ2dn6+uuvNWvWLMXExNS6Xe5wdwzcqG3btvrhD3+olStXatasWWrfvr0OHjyokSNHKjk5WS+++GKD1goAXs8AABpEWlqaad++faX2WbNmmRs/bq+/HjdunEu/hx56yERGRrq0rVq1ykgy+fn5xhhjFi9e7PK6NikpKSYgIMD8/e9/d7YVFhaa0NBQ893vfrdSTQ8++KDL/Onp6UaS+eyzz5xt48ePN23atDHnzp1z6ZuammrCw8PN5cuXXZb52GOPVarrep8b5ebmGknmzTffdLZVt73f3i83rq8u+3X48OEmODjYnDx50tl25MgR4+fnZ2r71RgeHm4yMjJq7PPkk08aHx8fk5eXV2maw+Fwu9667vOUlBQTGBhojh8/7uyzf/9+4+vrW2m76jpeq9rXda2nOv/4xz+M3W43EyZMMOfOnTOdOnUy3/nOd8xXX31V43wA0BxxSiQAeMiPfvQjl9cDBw7U+fPnVVJS0iDLr6io0KZNm5SSkqKOHTs629u0aaPHH39c27Ztq7SujIwMl9eTJk2SJP35z3+WJBlj9Pvf/17Dhw+XMUbnzp1z/gwdOlQXL16sdGrgt7dTkoKCgpz/vnbtms6fP6/bbrtNERERleZ3V237taKiQlu2bFFKSori4uKc/W677TYlJyfXuvyIiAjt3LlThYWFVU53OBzasGGDhg8frqSkpErTv326YW311nWfV1RU6MMPP1RKSooSEhKcy+vatauGDh1a63bVVX3GwLe1bdtWTz/9tF5//XUNGzZMZWVl+uMf/9gg1zcCQFNDYAMAD7nxS7UktWzZUpJ04cKFBll+UVGRLl++rNtvv73StK5du8rhcKigoMClvXPnzi6vO3XqJB8fH+e1YkVFRSouLtbKlSsVFRXl8jN27FhJ0tmzZ12WkZiYWGn9ZWVlmjlzpuLj42W329W6dWtFRUWpuLhYFy9evJnNrnW/nj17VmVlZbrtttsqzVtV27ctWrRI+/btU3x8vO6++25lZ2friy++cE4vKipSSUmJevTo0SD11nWfFxUVqaysrNJ7KKnKMVBf9RkDVfnxj3+s8vJy7d27V++9957atm3rMt0Yo5CQkDotCwCaMq5hA4AGUt2NGioqKqpsr+4uh+b/X+dlBd/eJofDIUkaM2aM0tLSqpznzjvvdHl949G06yZNmqRVq1YpMzNT/fr1U3h4uGw2m1JTU53rqK9bvV9HjRqlgQMHav369dq0aZMWL16shQsX6g9/+EOdjtB9W2311nWfu7vf3B2v19VnDFRl7ty5kqSvv/5arVq1qjQ9Pz9fwcHBio6OrnVZANCUEdgAoIG0bNmyyuejHT9+vMHWUZe7910XFRWl4OBgHTp0qNK0gwcPysfHR/Hx8S7tR44ccTkidvToUTkcDufdBKOiohQaGqqKigoNGTKkfhshad26dUpLS3O5wcSVK1cq7T93treuoqOjFRgYqKNHj1aaVlVbVdq0aaP09HSlp6fr7NmzuuuuuzR37lwlJycrKipKYWFh2rdvX4PUW9d9XlFRoaCgIB05cqTStKrGQH3Ha0OMgcWLF+u1117TK6+8oqlTp2ru3Ll67bXXnNMPHDigXr16qaKiQiEhIeratavy8vLqtS4A8HacEgkADaRTp066ePGi9u7d62w7deqU1q9f32DruH6NT10enO3r66v7779f7777rsvt78+cOaM1a9ZowIABCgsLc5nn1VdfdXm9dOlSSXIeOfL19dXDDz+s3//+91UGkqKiojpth6+vb6UjXkuXLq10dMed7a0rX19fDRkyRBs2bHC5Du3o0aPauHFjjfNWVFRUOmUzOjpacXFxKi8vlyT5+PgoJSVF77//vj755JNKy3D3SF9d97mvr6+GDh2qDRs26MSJE87pBw4c0IcfflhpvvqO15sdAxs2bNC0adM0Z84cZWRk6JlnntGbb76p/Px8Z5+uXbtq1qxZmjBhgr766ivCGoBmjSNsANBAUlNT9ZOf/EQPPfSQnn32WV2+fFnLli3TP/3TP930jTSu6927tyTp+eefV2pqqvz9/TV8+PBqb9bwwgsvaPPmzRowYIDS09Pl5+enFStWqLy8XIsWLarUPz8/Xw8++KAeeOAB5ebm6q233tLjjz+unj17OvssWLBAH330kfr27aunn35a3bp105dffqndu3dry5Yt+vLLL2vdjh/84Af6zW9+o/DwcHXr1k25ubnasmWLIiMj67S9Nys7O1ubNm1S//79NWHCBFVUVOiVV15Rjx49tGfPnmrnKy0tVbt27fTII4+oZ8+eCgkJ0ZYtW5SXl+dytHDevHnatGmT7rvvPj3zzDPq2rWrTp06pbVr12rbtm2KiIhwq9667vPZs2frgw8+0MCBA5Wenq6vv/5aS5cuVffu3V2CmXRz47W+Y2DXrl0aPXq0Ro8ereeff16S9Nxzz2n58uWVjrLt3btX3//+993aTwDQJHnm5pQA0DRt2rTJ9OjRwwQEBJjbb7/dvPXWW9Xe1r+oqMhl3qpun15V25w5c0zbtm2Nj49PnW7xv3v3bjN06FATEhJigoODzfe+9z2zfft2lz7Xa9q/f7955JFHTGhoqGnZsqWZOHGiKSsrq7TMM2fOmIyMDBMfH2/8/f1NbGysGTx4sFm5cmWt22mMMRcuXDBjx441rVu3NiEhIWbo0KHm4MGDpn379iYtLc2lb1XbW9Nt/euyX40xZuvWraZXr14mICDAdOrUybz22mvm3//9301gYGC1+7K8vNxMnTrV9OzZ04SGhpoWLVqYnj17ml/+8peV+h4/ftw8+eSTJioqytjtdtOxY0eTkZFhysvL61VvXfa5Mcbk5OSY3r17m4CAANOxY0ezfPnySmPwurqM15ut57qCggLTpk0b079/f3PlyhWXaRMmTDD+/v7miy++cLZ169bN7Nixo8plAUBzYjPGQle3AwA8Ijs7W7Nnz1ZRUZHzAdXNUUpKiv72t79VeR0YGs/Vq1cVEhKiCxcucKt/AM0e17ABAJqlsrIyl9dHjhzRn//8Zw0aNMgzBcGptLRU0jfBDQCaO65hAwA0Sx07dtRTTz2ljh076vjx41q2bJkCAgL03HPPebq0Zi8yMlKPPfaYEhIS1L17d+3YscPTJQGAxxDYAADN0gMPPKC3335bp0+flt1uV79+/TRv3rwqHzyNxvfrX/9av/71rz1dBgB4HNewAQAAAIBFcQ0bAAAAAFgUgQ0AAAAALMrtwPbXv/5Vw4cPV1xcnGw2mzZs2HALygIAAAAAuH3TkUuXLqlnz54aN26cRowY4fYKHQ6HCgsLFRoaKpvN5vb8AAAAAODtjDEqLS1VXFycfHyqP47mdmBLTk5WcnJyvQsrLCxUfHx8vecHAAAAgKaioKBA7dq1q3Z6o9/WPzQ0VNI3hYWFhTX26gEAAADA40pKShQfH+/MR9W55YGtvLxc5eXlztelpaWSpLCwMAIbAAAAgGattsvEbnlgmz9/vmbPnl1jnw7T/lTj9GMLhtWpT23Laqg+VqzJnbqtWJOV6rZiTby/TbsmK7+/AFBX3vr5DFjdLb+t//Tp03Xx4kXnT0FBwa1eJQAAAAA0Cbf8CJvdbpfdbr/VqwEAAACAJsftwPbVV1/p6NGjztf5+fnas2ePWrVqpYSEhAYtDgAAAACaM7cD2yeffKLvfe97ztdZWVmSpLS0NL3xxhsNVhgAAAAANHduB7ZBgwbJGHMragEAAAAA3OCW33QEAAAAAFA/BDYAAAAAsCgCGwAAAABYFIENAAAAACyKwAYAAAAAFkVgAwAAAACLIrABAAAAgEUR2AAAAADAoghsAAAAAGBRBDYAAAAAsCgCGwAAAABYFIENAAAAACyKwAYAAAAAFkVgAwAAAACL8vN0AQAAALC2DtP+VO20YwuGNWIlQPPDETYAAAAAsCgCGwAAAABYFIENAAAAACyKwAYAAAAAFkVgAwAAAACLIrABAAAAgEUR2AAAAADAoghsAAAAAGBRBDYAAAAAsCgCGwAAAABYlJ+nCwAAWFeHaX+qcfqxBcMaqRIAuDX4nIPVcYQNAAAAACyKwAYAAAAAFkVgAwAAAACLIrABAAAAgEUR2AAAAADAoghsAAAAAGBRBDYAAAAAsCgCGwAAAABYFIENAAAAACzKz9MFAAAAwHM6TPtTtdOOLRjWiJUAqApH2AAAAADAoghsAAAAAGBRBDYAAAAAsCgCGwAAAABYFIENAAAAACyKwAYAAAAAFkVgAwAAAACLIrABAAAAgEUR2AAAAADAoghsAAAAAGBRBDYAAAAAsCgCGwAAAABYFIENAAAAACyKwAYAAAAAFkVgAwAAAACLqldge/XVV9WhQwcFBgaqb9+++vjjjxu6LgAAAABo9twObO+8846ysrI0a9Ys7d69Wz179tTQoUN19uzZW1EfAAAAADRbbge2l156SU8//bTGjh2rbt26afny5QoODtbrr79+K+oDAAAAgGbLrcB29epV7dq1S0OGDPm/Bfj4aMiQIcrNzW3w4gAAAACgOfNzp/O5c+dUUVGhmJgYl/aYmBgdPHiwynnKy8tVXl7ufH3x4kVJUklJibPNUX65xvWWlJTUqU9ty2qoPlasyZ26rViTleq2Yk28v027Jm9/fwF4N6t9pljps/DGfkBDuz62jDE19rOZ2nrcoLCwUG3bttX27dvVr18/Z/tzzz2nnJwc7dy5s9I82dnZmj17dl1XAQAAAADNRkFBgdq1a1ftdLeOsLVu3Vq+vr46c+aMS/uZM2cUGxtb5TzTp09XVlaW87XD4dCXX36pyMhI2Ww2d1YPWFpJSYni4+NVUFCgsLAwT5cDNDjGOJo6xjiaMsa39RhjVFpaqri4uBr7uRXYAgIC1Lt3b23dulUpKSmSvglgW7du1cSJE6ucx263y263u7RFRES4s1rAq4SFhfFBiCaNMY6mjjGOpozxbS3h4eG19nErsElSVlaW0tLSlJSUpLvvvlsvv/yyLl26pLFjx9arSAAAAABA1dwObI8++qiKioo0c+ZMnT59Wt/5znf0wQcfVLoRCQAAAADg5rgd2CRp4sSJ1Z4CCTRXdrtds2bNqnQKMNBUMMbR1DHG0ZQxvr2XW3eJBAAAAAA0HrcenA0AAAAAaDwENgAAAACwKAIbAAAAAFgUgQ0AAAAALIrABrhp/vz56tOnj0JDQxUdHa2UlBQdOnTIpc+VK1eUkZGhyMhIhYSE6OGHH9aZM2c8VDFQfwsWLJDNZlNmZqazjfENb3fy5EmNGTNGkZGRCgoK0h133KFPPvnEOd0Yo5kzZ6pNmzYKCgrSkCFDdOTIEQ9WDNRNRUWFZsyYocTERAUFBalTp06aM2eObrzHIOPb+xDYADfl5OQoIyNDO3bs0ObNm3Xt2jXdf//9unTpkrPPlClT9P7772vt2rXKyclRYWGhRowY4cGqAffl5eVpxYoVuvPOO13aGd/wZhcuXFD//v3l7++vjRs3av/+/XrxxRfVsmVLZ59FixZpyZIlWr58uXbu3KkWLVpo6NChunLligcrB2q3cOFCLVu2TK+88ooOHDighQsXatGiRVq6dKmzD+PbCxkAN+Xs2bNGksnJyTHGGFNcXGz8/f3N2rVrnX0OHDhgJJnc3FxPlQm4pbS01HTu3Nls3rzZ3HfffWby5MnGGMY3vN9PfvITM2DAgGqnOxwOExsbaxYvXuxsKy4uNna73bz99tuNUSJQb8OGDTPjxo1zaRsxYoQZPXq0MYbx7a04wgbcpIsXL0qSWrVqJUnatWuXrl27piFDhjj7dOnSRQkJCcrNzfVIjYC7MjIyNGzYMJdxLDG+4f3ee+89JSUlaeTIkYqOjlavXr30q1/9yjk9Pz9fp0+fdhnj4eHh6tu3L2Mclnfvvfdq69atOnz4sCTps88+07Zt25ScnCyJ8e2t/DxdAODNHA6HMjMz1b9/f/Xo0UOSdPr0aQUEBCgiIsKlb0xMjE6fPu2BKgH3/Pa3v9Xu3buVl5dXaRrjG97uiy++0LJly5SVlaWf/vSnysvL07PPPquAgAClpaU5x3FMTIzLfIxxeINp06appKREXbp0ka+vryoqKjR37lyNHj1akhjfXorABtyEjIwM7du3T9u2bfN0KUCDKCgo0OTJk7V582YFBgZ6uhygwTkcDiUlJWnevHmSpF69emnfvn1avny50tLSPFwdcHN+97vfafXq1VqzZo26d++uPXv2KDMzU3FxcYxvL8YpkUA9TZw4UX/84x/10UcfqV27ds722NhYXb16VcXFxS79z5w5o9jY2EauEnDPrl27dPbsWd11113y8/OTn5+fcnJytGTJEvn5+SkmJobxDa/Wpk0bdevWzaWta9euOnHihCQ5x/G373zKGIc3mDp1qqZNm6bU1FTdcccdeuKJJzRlyhTNnz9fEuPbWxHYADcZYzRx4kStX79ef/nLX5SYmOgyvXfv3vL399fWrVudbYcOHdKJEyfUr1+/xi4XcMvgwYP1+eefa8+ePc6fpKQkjR492vlvxje8Wf/+/Ss9iuXw4cNq3769JCkxMVGxsbEuY7ykpEQ7d+5kjMPyLl++LB8f16/3vr6+cjgckhjf3opTIgE3ZWRkaM2aNXr33XcVGhrqPOc7PDxcQUFBCg8P1/jx45WVlaVWrVopLCxMkyZNUr9+/XTPPfd4uHqgZqGhoc7rMa9r0aKFIiMjne2Mb3izKVOm6N5779W8efM0atQoffzxx1q5cqVWrlwpSc7nDr7wwgvq3LmzEhMTNWPGDMXFxSklJcWzxQO1GD58uObOnauEhAR1795dn376qV566SWNGzdOEuPba3n6NpWAt5FU5c+qVaucfcrKykx6erpp2bKlCQ4ONg899JA5deqU54oGbsKNt/U3hvEN7/f++++bHj16GLvdbrp06WJWrlzpMt3hcJgZM2aYmJgYY7fbzeDBg82hQ4c8VC1QdyUlJWby5MkmISHBBAYGmo4dO5rnn3/elJeXO/swvr2PzZgbHn0OAAAAALAMrmEDAAAAAIsisAEAAACARRHYAAAAAMCiCGwAAAAAYFEENgAAAACwKAIbAAAAAFgUgQ0AAAAALIrABgAAAAAWRWADAAAAAIsisAEAAACARRHYAAAAAMCiCGwAAAAAYFEENgAAAACwKAIbAAAAAFgUgQ0AAAAALIrABgAAAAAWRWADAAAAAIsisAEAAACARRHYAAAAAMCiCGwAgCq98cYbstlsOnbsmKdLAQCg2SKwAYCHzZ07VzabTT169PB0KQ1m+/btys7OVnFxsadL8Vq3eh/yHgGAdyCwAYAH/eMf/9C8efPUokULT5dSyRNPPKGysjK1b9/e7Xm3b9+u2bNnEwZuwq3eh7xHAOAd/DxdAAA0Zz/+8Y91zz33qKKiQufOnfN0OS58fX3l6+vr6TLQwC5dumTJPxAAAKrGETYA8JC//vWvWrdunV5++eU6z5OdnS2bzabDhw9rzJgxCg8PV1RUlGbMmCFjjAoKCvSv//qvCgsLU2xsrF588cVKyzh+/LjS09N1++23KygoSJGRkRo5cmSla9W+fQ3b9XUfPXpUTz31lCIiIhQeHq6xY8fq8uXLLjVOnTpVkpSYmCibzeaynJMnT2rcuHGKiYmR3W5X9+7d9frrr1e7rbWtr67bU9flXXfy5EmNHz9ecXFxstvtSkxM1IQJE3T16lWXPnXZlqqUlpYqMzNTHTp0kN1uV3R0tP75n/9Zu3fvrnUfurvN+/fv1+OPP66WLVtqwIABNS6/prpq8i//8i/q0KFDpXZjjO666y4NHDiwTvsFAOCKI2wA4AEVFRWaNGmSfvjDH+qOO+5we/5HH31UXbt21YIFC/SnP/1JL7zwglq1aqUVK1bo+9//vhYuXKjVq1frxz/+sfr06aPvfve7znnz8vK0fft2paamql27djp27JiWLVumQYMGaf/+/QoODq5x3aNGjVJiYqLmz5+v3bt367XXXlN0dLQWLlwoSRoxYoQOHz6st99+Wz//+c/VunVrSVJUVJTOnDmje+65RzabTRMnTlRUVJQ2btyo8ePHq6SkRJmZmW6vz93tqW15klRYWKi7775bxcXFeuaZZ9SlSxedPHlS69at0+XLlxUQEFCvbbnRj370I61bt04TJ05Ut27ddP78eW3btk0HDhyocR/WZ5tHjhypzp07a968eTLGaMCAAdUu/5lnnqm2rrvuuqva7enTp482btyoCxcuqGXLls723/72t/r000+1bdu2GvcHAKAaBgDQ6F555RUTHh5uzp49a4wx5r777jPdu3evdb5Zs2YZSeaZZ55xtn399demXbt2xmazmQULFjjbL1y4YIKCgkxaWprLMi5fvlxpubm5uUaSefPNN51tq1atMpJMfn6+y7rHjRvnMu9DDz1kIiMjXdoWL17sMu9148ePN23atDHnzp1zaU9NTTXh4eEutdV1fXXdHnfqf/LJJ42Pj4/Jy8urtGyHw+H2tlQlPDzcZGRkVDu9un1ojPvb/Nhjj9V5+bXVVZ333nvPSDJbt251tl29etV06tTJDB8+3O3lAQC+wSmRANDIzp8/r5kzZ2rGjBnOIybu+uEPf+j8t6+vr5KSkmSM0fjx453tERERuv322/XFF1+4zBsUFOT897Vr13T+/HnddtttioiIqPW0N+mbI0M3GjhwoM6fP6+SkpIa5zPG6Pe//72GDx8uY4zOnTvn/Bk6dKguXrxY5fprW5+721Pb8hwOhzZs2KDhw4crKSmp0vw2m63e23KjiIgI7dy5U4WFhTX2q8rNbvOtqKtPnz6S5LL+lStXKj8/X/PmzXNrWQCA/0NgA4BG9h//8R9q1aqVJk2aVO9lJCQkuLwODw9XYGCg89S2G9svXLjg0lZWVqaZM2cqPj5edrtdrVu3VlRUlIqLi3Xx4kW313399Ldvr+fbioqKVFxcrJUrVyoqKsrlZ+zYsZKks2fPur0+d7entuUVFRWppKSkxscs1HdbbrRo0SLt27dP8fHxuvvuu5WdnV0pXFfH3W1OTEys03Jvpq7Y2Fi1bdtWn376qaRvbm4yZ84cjRkzxmVfGmMUEhJS6/4BAHyDa9gAoBEdOXJEK1eu1Msvv+xyBOPKlSu6du2ajh07prCwMLVq1arG5VR198bq7uhojHF5PWnSJK1atUqZmZnq16+fwsPDZbPZlJqaKofDUes21HU933Z92WPGjFFaWlqVfe6880631+fu9tS3/hvVd1tuNGrUKA0cOFDr16/Xpk2btHjxYi1cuFB/+MMflJycXOO87m7zjUfkanMzdfXp08cZ2F566SVduHBBP/vZz1z65OfnKzg4WNHR0XWuCQCaMwIbADSikydPyuFw6Nlnn9Wzzz5baXpiYqImT57s1p0j3bVu3TqlpaW53EHyypUrDfo8LpvNVqktKipKoaGhqqio0JAhQxpsXQ29PVFRUQoLC9O+fftq7NMQ29KmTRulp6crPT1dZ8+e1V133aW5c+cqOTm5yn14XUNsc03Lr6mumvTp00fvvfeeTpw4of/8z//UhAkTXJ7jd+DAAfXq1UsVFRUKCQlR165dlZeXV+eaAaA5IrABQCPq0aOH1q9fX6n9P/7jP1RaWqpf/OIX6tSp0y2twdfXt9LRpKVLl6qioqLB1nH9OV83BghfX189/PDDWrNmjfbt21fplMOioqJ6XdPX0Nvj4+OjlJQUvfXWW/rkk08qXcdmjLnpbamoqNBXX32l8PBwZ1t0dLTi4uJUXl4uqep9eF1DbHNVy69LXTVJSkqSw+HQ448/LmOMnn/+eZfpXbt21axZs3Tq1CktWbKkzrUCQHNGYAOARtS6dWulpKRUar9+RK2qaQ3tBz/4gX7zm98oPDxc3bp1U25urrZs2aLIyMgGW0fv3r0lSc8//7xSU1Pl7++v4cOHa8GCBfroo4/Ut29fPf300+rWrZu+/PJL7d69W1u2bNGXX35pie2ZN2+eNm3apPvuu0/PPPOMunbtqlOnTmnt2rXatm2bIiIibmpbSktL1a5dOz3yyCPq2bOnQkJCtGXLFuXl5TmPmlW3D1u0aNEg21zV8gcOHKjbb7+9xrpqcj3c/u///q+ys7OrDK179+7V97///TrXCQDNHYENAJqZX/ziF/L19dXq1at15coV9e/fX1u2bNHQoUMbbB19+vTRnDlztHz5cn3wwQdyOBzKz89Xhw4d9PHHH+tnP/uZ/vCHP+iXv/ylIiMj1b17d5fnoHl6e9q2baudO3dqxowZWr16tUpKStS2bVslJyc7n3EWExNT720JDg5Wenq6Nm3apD/84Q9yOBy67bbb9Mtf/lITJkyQVP0+bNGiRYNsc1XLP3ToUK111aRVq1bq0KGDLl26pH//93+vss/evXtrfUYdAOD/2Iw7V1kDAABU44svvtA//dM/6aWXXqryGs2rV68qJCREFy5ccJ6SCQCoGYENAAA0iEcffVS7du3S/v37FRAQUGn6+fPn1aZNG505c8b5OAUAQM04JRIAANRbcXGxNm7cqP/+7//W2rVrtXHjxirDmiRFRkbqscceU0JCgrp3764dO3Y0crUA4H04wgYAAOrt97//vR555BG1a9dOM2fO1NNPP+3pkgCgSSGwAQAAAIBF+Xi6AAAAAABA1QhsAAAAAGBRbt905K9//asWL16sXbt26dSpU1q/fr1bD3p1OBwqLCxUaGiobDabu6sHAAAAAK9njFFpaani4uLk41P9cTS3A9ulS5fUs2dPjRs3TiNGjHC7sMLCQsXHx7s9HwAAAAA0NQUFBWrXrl21090ObMnJyUpOTq53QaGhoc7CwsLC6r0cAAAAAPBWJSUlio+Pd+aj6tzy57CVl5ervLzc+bq0tFSSFBYWRmADAAAA0KzVdpnYLQ9s8+fP1+zZs2/1agAAHtJh2p9qnH5swbBGqgS3Sk3vMe8vmgM+55o2q7+/t/wukdOnT9fFixedPwUFBbd6lQAAAADQJNzyI2x2u112u/1WrwYAAAAAmhyewwYAAAAAFuX2EbavvvpKR48edb7Oz8/Xnj171KpVKyUkJDRocQAAAADQnLkd2D755BN973vfc77OysqSJKWlpemNN95osMIAAAAAoLlzO7ANGjRIxphbUQsAAAAA4AZcwwYAAAAAFkVgAwAAAACLIrABAAAAgEUR2AAAAADAoghsAAAAAGBRBDYAAAAAsCgCGwAAAABYFIENAAAAACyKwAYAAAAAFkVgAwAAAACLIrABAAAAgEUR2AAAAADAoghsAAAAAGBRBDYAAAAAsCgCGwAAAABYFIENAAAAACyKwAYAAAAAFkVgAwAAAACLIrABAAAAgEUR2AAAAADAoghsAAAAAGBRBDYAAAAAsCgCGwAAAABYFIENAAAAACyKwAYAAAAAFkVgAwAAAACLIrABAAAAgEUR2AAAAADAoghsAAAAAGBRBDYAAAAAsCgCGwAAAABYFIENAAAAACyKwAYAAAAAFkVgAwAAAACLIrABAAAAgEUR2AAAAADAoghsAAAAAGBRBDYAAAAAsCgCGwAAAABYFIENAAAAACyKwAYAAAAAFkVgAwAAAACLIrABAAAAgEUR2AAAAADAoghsAAAAAGBRBDYAAAAAsCgCGwAAAABYFIENAAAAACyKwAYAAAAAFkVgAwAAAACLIrABAAAAgEXVK7C9+uqr6tChgwIDA9W3b199/PHHDV0XAAAAADR7bge2d955R1lZWZo1a5Z2796tnj17aujQoTp79uytqA8AAAAAmi23A9tLL72kp59+WmPHjlW3bt20fPlyBQcH6/XXX78V9QEAAABAs+XnTuerV69q165dmj59urPNx8dHQ4YMUW5ubpXzlJeXq7y83Pn64sWLkqSSkpL61AsAsBhH+eUap/N57/1qeo95f9Ec8DnXtHnq/b2+XGNMjf3cCmznzp1TRUWFYmJiXNpjYmJ08ODBKueZP3++Zs+eXak9Pj7enVUDALxU+MuergC3Eu8vwP+Dpu5Wv7+lpaUKDw+vdrpbga0+pk+frqysLOdrh8OhL7/8UpGRkbLZbLd69UCjKSkpUXx8vAoKChQWFubpcoAGxxhHU8cYR1PG+LYeY4xKS0sVFxdXYz+3Alvr1q3l6+urM2fOuLSfOXNGsbGxVc5jt9tlt9td2iIiItxZLeBVwsLC+CBEk8YYR1PHGEdTxvi2lpqOrF3n1k1HAgIC1Lt3b23dutXZ5nA4tHXrVvXr18/9CgEAAAAA1XL7lMisrCylpaUpKSlJd999t15++WVdunRJY8eOvRX1AQAAAECz5XZge/TRR1VUVKSZM2fq9OnT+s53vqMPPvig0o1IgObGbrdr1qxZlU4BBpoKxjiaOsY4mjLGt/eymdruIwkAAAAA8Ai3H5wNAAAAAGgcBDYAAAAAsCgCGwAAAABYFIENAAAAACyKwAa4af78+erTp49CQ0MVHR2tlJQUHTp0yKXPlStXlJGRocjISIWEhOjhhx+u9MB5wBssWLBANptNmZmZzjbGN7zdyZMnNWbMGEVGRiooKEh33HGHPvnkE+d0Y4xmzpypNm3aKCgoSEOGDNGRI0c8WDFQNxUVFZoxY4YSExMVFBSkTp06ac6cObrxHoOMb+9DYAPclJOTo4yMDO3YsUObN2/WtWvXdP/99+vSpUvOPlOmTNH777+vtWvXKicnR4WFhRoxYoQHqwbcl5eXpxUrVujOO+90aWd8w5tduHBB/fv3l7+/vzZu3Kj9+/frxRdfVMuWLZ19Fi1apCVLlmj58uXauXOnWrRooaFDh+rKlSserByo3cKFC7Vs2TK98sorOnDggBYuXKhFixZp6dKlzj6Mby9kANyUs2fPGkkmJyfHGGNMcXGx8ff3N2vXrnX2OXDggJFkcnNzPVUm4JbS0lLTuXNns3nzZnPfffeZyZMnG2MY3/B+P/nJT8yAAQOqne5wOExsbKxZvHixs624uNjY7Xbz9ttvN0aJQL0NGzbMjBs3zqVtxIgRZvTo0cYYxre34ggbcJMuXrwoSWrVqpUkadeuXbp27ZqGDBni7NOlSxclJCQoNzfXIzUC7srIyNCwYcNcxrHE+Ib3e++995SUlKSRI0cqOjpavXr10q9+9Svn9Pz8fJ0+fdpljIeHh6tv376McVjevffeq61bt+rw4cOSpM8++0zbtm1TcnKyJMa3t/LzdAGAN3M4HMrMzFT//v3Vo0cPSdLp06cVEBCgiIgIl74xMTE6ffq0B6oE3PPb3/5Wu3fvVl5eXqVpjG94uy+++ELLli1TVlaWfvrTnyovL0/PPvusAgIClJaW5hzHMTExLvMxxuENpk2bppKSEnXp0kW+vr6qqKjQ3LlzNXr0aElifHspAhtwEzIyMrRv3z5t27bN06UADaKgoECTJ0/W5s2bFRgY6OlygAbncDiUlJSkefPmSZJ69eqlffv2afny5UpLS/NwdcDN+d3vfqfVq1drzZo16t69u/bs2aPMzEzFxcUxvr0Yp0QC9TRx4kT98Y9/1EcffaR27do522NjY3X16lUVFxe79D9z5oxiY2MbuUrAPbt27dLZs2d11113yc/PT35+fsrJydGSJUvk5+enmJgYxje8Wps2bdStWzeXtq5du+rEiROS5BzH377zKWMc3mDq1KmaNm2aUlNTdccdd+iJJ57QlClTNH/+fEmMb29FYAPcZIzRxIkTtX79ev3lL39RYmKiy/TevXvL399fW7dudbYdOnRIJ06cUL9+/Rq7XMAtgwcP1ueff649e/Y4f5KSkjR69Gjnvxnf8Gb9+/ev9CiWw4cPq3379pKkxMRExcbGuozxkpIS7dy5kzEOy7t8+bJ8fFy/3vv6+srhcEhifHsrTokE3JSRkaE1a9bo3XffVWhoqPOc7/DwcAUFBSk8PFzjx49XVlaWWrVqpbCwME2aNEn9+vXTPffc4+HqgZqFhoY6r8e8rkWLFoqMjHS2M77hzaZMmaJ7771X8+bN06hRo/Txxx9r5cqVWrlypSQ5nzv4wgsvqHPnzkpMTNSMGTMUFxenlJQUzxYP1GL48OGaO3euEhIS1L17d3366ad66aWXNG7cOEmMb6/l6dtUAt5GUpU/q1atcvYpKysz6enppmXLliY4ONg89NBD5tSpU54rGrgJN97W3xjGN7zf+++/b3r06GHsdrvp0qWLWblypct0h8NhZsyYYWJiYozdbjeDBw82hw4d8lC1QN2VlJSYyZMnm4SEBBMYGGg6duxonn/+eVNeXu7sw/j2PjZjbnj0OQAAAADAMriGDQAAAAAsisAGAAAAABZFYAMAAAAAiyKwAQAAAIBFEdgAAAAAwKIIbAAAAABgUQQ2AAAAALAoAhsAAAAAWBSBDQAAAAAsisAGAAAAABZFYAMAAAAAiyKwAQAAAIBF/T9Ir+c+zd7+GgAAAABJRU5ErkJggg==",
                        "text/plain": [
                            "<Figure size 900x400 with 3 Axes>"
                        ]
                    },
                    "metadata": {},
                    "output_type": "display_data"
                }
            ],
            "source": [
                "class MPRampOperational(MPRamp):\n",
                "    def __init__(self, c, M, P, upos_max, uneg_max, N):\n",
                "        super(MPRampOperational, self).__init__(c, M, P, upos_max, uneg_max)\n",
                "        super(MPRampOperational, self).set_modfile(\n",
                "            \"maintenance_planning_ramp_operational.mod\"\n",
                "        )\n",
                "        self.N = N\n",
                "\n",
                "    def load_data(self):\n",
                "        super(MPRampOperational, self).load_data()\n",
                "        self.ampl.param[\"N\"] = self.N\n",
                "\n",
                "\n",
                "N = 10  # minimum number of operational days between maintenance periods\n",
                "model = MPRampOperational(c, M, P, upos_max, uneg_max, N)\n",
                "model.plot_schedule()"
            ]
        },
        {
            "cell_type": "markdown",
            "metadata": {
                "colab_type": "text",
                "id": "3j0fLxBOwSch"
            },
            "source": [
                "## Exercises\n",
                "\n",
                "1. Rather than specify how many maintenance periods must be accommodated, modify the model so that the process unit can operate no more than $N$ days without a maintenance shutdown.  (Hint. You may to introduce an additional set of binary variables, $z_t$ to denote the start of an operational period.)\n",
                "\n",
                "2. Do a systematic comparison of the Big-M, Convex Hull, and Cutting Plane techniques for implementing the disjunctive constraints. Your comparison should include a measure of complexity (such as the number of decision variables and constraints in the resulting transformed problems), computational effort, and the effect of solver (such as HiGHS vs cbc)."
            ]
        }
    ],
    "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.10.6"
        }
    },
    "nbformat": 4,
    "nbformat_minor": 4
}