{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 1,
   "id": "7812df96",
   "metadata": {},
   "outputs": [],
   "source": [
    "from enum import Enum, auto\n",
    "from itertools import combinations\n",
    "\n",
    "import numpy as np\n",
    "from scipy.optimize import minimize"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "1d65b350",
   "metadata": {},
   "source": [
    "# Instructions\n",
    "- Ensure that you are using the `Optimization` kernel in Noto.\n",
    "- Complete all tasks and verify that your notebook is functioning correctly. To do this, restart the kernel and run all cells by clicking the double-arrow button in Jupyter notebook. Confirm that all cells execute without errors.\n",
    "- The parts of the code that you are supposed to fill are replaced by so-called \"ellipsis\", that is, three consecutive dots: `...`. If you forget to replace some of them and try to run the notebook, you will encounter errors similar to `TypeError: float() argument must be a string or a real number, not 'ellipsis'`. Make sure that all ellipses are replaced by some valid code, so that the notebook runs correctly from begin to end. \n",
    "- Do not change the specifications of the Python functions. The grading is performed by calling the functions. If the specification does not match, the result is considered wrong.\n",
    "- In principle, you should be able to complete all questions without importing extra modules. However, if you decide to import additional modules, double check that they are available, and that the notebook continues to run properly.    \n",
    "- You have access to any material you find helpful. Specifically, make sure to refer to the notebooks distributed weekly. You are also permitted to use online AI tools such as ChatGPT.\n",
    "- Communication with fellow students is not allowed."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "487dfd74",
   "metadata": {},
   "source": [
    "You will have to use some linear algebra tools, available in `numpy.linalg`. For instance, the function that\n",
    "calculates the rank of a matrix is [`matrix_rank`](\n",
    "https://numpy.org/doc/stable/reference/generated/numpy.linalg.matrix_rank.html)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "5d2ce23c",
   "metadata": {},
   "source": [
    "# Question 1\n",
    "At a position $x_0$, a projectile is launched vertically at a rate of $v_0$ meters per second in the absence of wind.\n",
    "We need to calculate the time when it will reach again the ground, that is, its **lowest** altitude.\n",
    "Your task is to write a function that returns that quantity, using the optimization routine provided by `scipy`."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "36660d4a",
   "metadata": {},
   "source": [
    "We first provide you with some useful functions that you do not need to implement yourself."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "42e3a821",
   "metadata": {},
   "source": [
    "The altitude of the projectile is given by the formula for uniformly accelerated\n",
    "movement:  $$f(t) = x_0 + v_0 t -\\frac{g}{2} t^2,$$\n",
    "where $x_0$ is the initial altitude, $v_0$ the initial speed, and $g$ the acceleration due\n",
    "to gravity."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "id": "5ff18f76",
   "metadata": {},
   "outputs": [],
   "source": [
    "GRAVITY = 9.81\n",
    "\n",
    "\n",
    "def calculate_height(\n",
    "    time: float, initial_altitude: float, initial_speed: float\n",
    ") -> float:\n",
    "    \"\"\"\n",
    "    Calculate the height of the projectile, using the formula.\n",
    "\n",
    "    :param time: time at which we need the height.\n",
    "    :param initial_altitude: initial altitude x_0.\n",
    "    :param initial_speed: initial speed v_0.\n",
    "    :return: height.\n",
    "    \"\"\"\n",
    "    return initial_altitude + initial_speed * time - GRAVITY * time * time / 2\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "b69d3f1d",
   "metadata": {},
   "source": [
    "Now, it is your task to implement a function that calculates the time it takes to hit the ground."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "0bc93d4a-eeba-4805-96b4-1e2a80f5eebd",
   "metadata": {},
   "source": [
    "<div style=\"border-left: 5px solid #f44336; background-color: #fdecea; padding: 10px;\">\n",
    "<strong>Warning:</strong> Do not include any print statement in the body of the function below. It perturbs the grading algorithm. Either add new cells, or make sure to remove any print statement before submitting your work.\n",
    "</div>"
   ]
  },
  {
   "cell_type": "raw",
   "id": "96530c30-62e4-460b-8e97-f6b5e4fccd3f",
   "metadata": {},
   "source": [
    "# BEGIN SOLUTION"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "id": "c79aa527",
   "metadata": {},
   "outputs": [],
   "source": [
    "def time_to_hit_the_ground(initial_altitude: float, initial_speed: float) -> float:\n",
    "    \"\"\"Calculates the time it takes to hit the ground\n",
    "\n",
    "    :param initial_altitude: initial altitude x_0.\n",
    "    :param initial_speed: initial speed v_0.\n",
    "    :return: time after which the projectile hit the ground.\n",
    "    \"\"\"\n",
    "\n",
    "    # Define the objective function\n",
    "    def objective_function(x: float) -> float:\n",
    "        \"\"\"Objective function of the optimization problem.\"\"\"\n",
    "        return calculate_height(\n",
    "            time=x, initial_altitude=initial_altitude, initial_speed=initial_speed\n",
    "        )\n",
    "\n",
    "    # Define the starting point\n",
    "    x0 = 200  # SOLUTION\n",
    "\n",
    "    # Define the bound constraints.\n",
    "    bounds = [(0, None)]\n",
    "\n",
    "    # Define an inequality constraint\n",
    "    def inequality_constraint(x: float) -> float:\n",
    "        \"\"\"Altitude of the projectile\"\"\"\n",
    "        return calculate_height(\n",
    "            time=x, initial_altitude=initial_altitude, initial_speed=initial_speed\n",
    "        )\n",
    "        \n",
    "\n",
    "    # Run the algorithm.\n",
    "    optimization_result = minimize(\n",
    "        fun=objective_function,\n",
    "        bounds=bounds,\n",
    "        constraints=[{'type': 'ineq', 'fun': inequality_constraint}],\n",
    "        x0=x0,\n",
    "    )\n",
    "\n",
    "    elapsed_time = optimization_result.x[0]\n",
    "    return elapsed_time\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "id": "89b23ac3",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Elapsed time:          20.4 sec. Height: 0.0\n"
     ]
    }
   ],
   "source": [
    "# Test the function with different values for $x_0$ and $v_0$. Note that $v_0$ can be negative.\n",
    "the_initial_altitude = 0\n",
    "the_initial_speed = 100\n",
    "the_elapsed_time = time_to_hit_the_ground(\n",
    "    initial_altitude=the_initial_altitude, initial_speed=the_initial_speed\n",
    ")\n",
    "optimal_height = calculate_height(\n",
    "    time=the_elapsed_time,\n",
    "    initial_altitude=the_initial_altitude,\n",
    "    initial_speed=the_initial_speed,\n",
    ")\n",
    "print(\n",
    "    f'Elapsed time:          {the_elapsed_time:.1f} sec. Height: {optimal_height:.1f}'\n",
    ")\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "id": "bd524edd",
   "metadata": {},
   "outputs": [],
   "source": [
    "\"\"\" \n",
    "success_message: Tests with the initial altitude is zero. OK.\n",
    "failure_message: Tests with the initial altitude is zero. failed. \n",
    "\"\"\"\n",
    "def _analytical_solution(initial_altitude: float, initial_speed: float) -> float:\n",
    "    delta = initial_speed * initial_speed + 2 * GRAVITY * initial_altitude\n",
    "    solution = (initial_speed + np.sqrt(delta)) / GRAVITY\n",
    "    return solution\n",
    "\n",
    "\n",
    "def _verify_solution(initial_altitude: float, initial_speed: float) -> bool:\n",
    "    the_analytical_solution = _analytical_solution(\n",
    "        initial_altitude=initial_altitude, initial_speed=initial_speed\n",
    "    )\n",
    "    the_elapsed_time = time_to_hit_the_ground(\n",
    "        initial_altitude=initial_altitude, initial_speed=initial_speed\n",
    "    )\n",
    "    return np.isclose(the_analytical_solution, the_elapsed_time)\n",
    "\n",
    "\n",
    "_tests = [\n",
    "    (0, 10),\n",
    "    (0, 20),\n",
    "    (0, 30),\n",
    "]\n",
    "\n",
    "for the_initial_altitude, the_initial_speed in _tests:\n",
    "    correct = _verify_solution(\n",
    "        initial_altitude=the_initial_altitude, initial_speed=the_initial_speed\n",
    "    )\n",
    "    assert correct == True\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "id": "eae6927e-5177-4dd3-ad9e-832a59fa7baa",
   "metadata": {},
   "outputs": [],
   "source": [
    "\"\"\" \n",
    "success_message: Tests with the initial altitude is 10. OK.\n",
    "failure_message: Tests with the initial altitude is 10. failed.\n",
    "\"\"\" \n",
    "def _analytical_solution(initial_altitude: float, initial_speed: float) -> float:\n",
    "    delta = initial_speed * initial_speed + 2 * GRAVITY * initial_altitude\n",
    "    solution = (initial_speed + np.sqrt(delta)) / GRAVITY\n",
    "    return solution\n",
    "\n",
    "\n",
    "def _verify_solution(initial_altitude: float, initial_speed: float) -> bool:\n",
    "    the_analytical_solution = _analytical_solution(\n",
    "        initial_altitude=initial_altitude, initial_speed=initial_speed\n",
    "    )\n",
    "    the_elapsed_time = time_to_hit_the_ground(\n",
    "        initial_altitude=initial_altitude, initial_speed=initial_speed\n",
    "    )\n",
    "    return np.isclose(the_analytical_solution, the_elapsed_time)\n",
    "\n",
    "\n",
    "_tests = [\n",
    "    (10, 10),\n",
    "    (10, 20),\n",
    "    (10, 30),\n",
    "]\n",
    "\n",
    "for the_initial_altitude, the_initial_speed in _tests:\n",
    "    correct = _verify_solution(\n",
    "        initial_altitude=the_initial_altitude, initial_speed=the_initial_speed\n",
    "    )\n",
    "    assert correct == True\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "id": "46587018-7923-4546-89db-75ffa35104a9",
   "metadata": {},
   "outputs": [],
   "source": [
    "\"\"\"\n",
    "success_message: Tests with the initial speed is negative. OK.\n",
    "failure_message: Tests with the initial speed is negative. failed.\n",
    "\"\"\" \n",
    "def _analytical_solution(initial_altitude: float, initial_speed: float) -> float:\n",
    "    delta = initial_speed * initial_speed + 2 * GRAVITY * initial_altitude\n",
    "    solution = (initial_speed + np.sqrt(delta)) / GRAVITY\n",
    "    return solution\n",
    "\n",
    "\n",
    "def _verify_solution(initial_altitude: float, initial_speed: float) -> bool:\n",
    "    the_analytical_solution = _analytical_solution(\n",
    "        initial_altitude=initial_altitude, initial_speed=initial_speed\n",
    "    )\n",
    "    the_elapsed_time = time_to_hit_the_ground(\n",
    "        initial_altitude=initial_altitude, initial_speed=initial_speed\n",
    "    )\n",
    "    return np.isclose(the_analytical_solution, the_elapsed_time)\n",
    "\n",
    "\n",
    "_tests = [\n",
    "    (20, -1),\n",
    "    (20, -2),\n",
    "    (20, -3),\n",
    "]\n",
    "\n",
    "for the_initial_altitude, the_initial_speed in _tests:\n",
    "    correct = _verify_solution(\n",
    "        initial_altitude=the_initial_altitude, initial_speed=the_initial_speed\n",
    "    )\n",
    "    assert correct == True\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "71ed4763-1671-4d45-a963-38c447e18adf",
   "metadata": {},
   "source": [
    "# Question 2"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "2265ecf7-0b32-444e-8650-940c74e9eee5",
   "metadata": {},
   "source": [
    "A robot is moving along a straight path on a one-dimensional line (the $x$-axis) from an initial position $x_0 = 0$ to a target position at 100 meters away: $x_T = 100$. The robot's journey is divided into three stages, and the robot moves with a constant velocity during each stage, but not faster then 4 meter/seconds. You need to program the velocity of the robot during each stage in order to minimize the total energy consumption over the journey. \n",
    "\n",
    "The time (in seconds) spent in each stage is known and given by:\n",
    "$$\n",
    "t_1 = 10, \\quad t_2 = 15, \\quad t_3 = 5.\n",
    "$$\n",
    "\n",
    "The energy consumption during each stage is proportional to the velocity, and the cost of energy per unit velocity is given by:\n",
    "$$\n",
    "c_1 = 2, \\quad c_2 = 3, \\quad c_3 = 1.\n",
    "$$\n",
    "\n",
    "Formulate this as a linear optimization problem in **standard form**. Provide the matrix $A$, the vector $b$ and the vector $c$. Do not modify the name of the variables `matrix_A`, `vector_b` and `vector_c`.\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "71f462d7-1e0e-4634-a017-7da381cb433e",
   "metadata": {},
   "source": [
    "<div style=\"border-left: 5px solid #f44336; background-color: #fdecea; padding: 10px;\">\n",
    "<strong>Warning:</strong> Do not include any print statement in the body of the function below. It perturbs the grading algorithm. Either add new cells, or make sure to remove any print statement before submitting your work.\n",
    "</div>"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "id": "4a186b99-ca90-4db4-a702-bb1e0ebf9955",
   "metadata": {},
   "outputs": [],
   "source": [
    "matrix_A = np.array([[10, 15, 5, 0, 0, 0], [1, 0, 0, 1, 0, 0], [0, 1, 0, 0, 1, 0], [0, 0, 1, 0, 0, 1]])\n",
    "vector_b = np.array([100, 4, 4, 4])\n",
    "vector_c = np.array([2, 3, 1, 0, 0, 0])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "id": "3a1e6e21-6b66-47db-8cc1-c98673be27f0",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[[10 15  5  0  0  0]\n",
      " [ 1  0  0  1  0  0]\n",
      " [ 0  1  0  0  1  0]\n",
      " [ 0  0  1  0  0  1]]\n"
     ]
    }
   ],
   "source": [
    "print(matrix_A)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "id": "606d4dc7-52ce-4607-b8c0-87fbc86f4aa6",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[100   4   4   4]\n"
     ]
    }
   ],
   "source": [
    "print(vector_b)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "id": "2259fc43-11a0-4646-b920-7a0a10480fdb",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[2 3 1 0 0 0]\n"
     ]
    }
   ],
   "source": [
    "print(vector_c)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "id": "0ecc7716-d058-40b8-b524-dcbcafc8bd9f",
   "metadata": {},
   "outputs": [],
   "source": [
    "\"\"\" \n",
    "success_message: The matrix matrix_A has the correct shape. \n",
    "failure_message: The matrix matrix_A should be 4 x 6\n",
    "\"\"\" \n",
    "_expected_number_of_variables = 6\n",
    "_expected_number_of_constraints = 4\n",
    "_number_of_constraints, _number_of_variables = matrix_A.shape\n",
    "assert((_expected_number_of_variables == _number_of_variables) and (_expected_number_of_constraints == _number_of_constraints))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "id": "b8a512b7-c198-44b0-bd39-d20e790f2422",
   "metadata": {},
   "outputs": [],
   "source": [
    "\"\"\" \n",
    "success_message: The vector vector_b has the correct length.\n",
    "failure_message: The vector vector_b should contain 4 entries. \n",
    "\"\"\" \n",
    "_expected_number_of_constraints = 4\n",
    "assert len(vector_b) == _expected_number_of_constraints"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "id": "ea13d6ff-72ed-4219-9b1d-a139bd681467",
   "metadata": {},
   "outputs": [],
   "source": [
    "\"\"\" \n",
    "success_message: The vector vector_c has the correct length. \n",
    "failure_message: The vector vector_c should contain 6 entries \n",
    "\"\"\" \n",
    "_expected_number_of_variables = 6\n",
    "assert len(vector_c) == _expected_number_of_variables"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "id": "99383fc0-4d7a-4d15-a8f0-1ec7c5fac2fc",
   "metadata": {},
   "outputs": [],
   "source": [
    "\"\"\" \n",
    "success_message: The vector vector_b is correct. \n",
    "failure_message: The vector vector_b should be [100, 4, 4, 4]\n",
    "\"\"\" \n",
    "_expected_b = np.array([100, 4, 4, 4])     \n",
    "assert np.array_equal(np.sort(_expected_b), np.sort(vector_b)) == True"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "id": "fbf627be-7484-4dd0-a6dc-2c638108a7c6",
   "metadata": {},
   "outputs": [],
   "source": [
    "\"\"\" \n",
    "success_message: The vector vector_c is correct. \n",
    "failure_message: The vector vector_c should be [2, 3, 1, 0, 0, 0]\n",
    "\"\"\" \n",
    "_expected_c = np.array([2, 3, 1, 0, 0, 0]) \n",
    "assert np.array_equal(np.sort(_expected_c), np.sort(vector_c)) == True"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "d09e669f",
   "metadata": {},
   "source": [
    "# Question 3\n",
    "Consider the linear $m$ constraints $Ax = b$, where $x\\in\\mathbb{R^n}$, $b\\in\\mathbb{R}^m$ and\n",
    "$A \\in \\mathbb{R}^{m\\times n}$.\n",
    "We want to know how many constraints are redundant and can be removed without changing the feasible set."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "748285f0-f791-4512-9977-5c59c14d5fba",
   "metadata": {},
   "source": [
    "<div style=\"border-left: 5px solid #f44336; background-color: #fdecea; padding: 10px;\">\n",
    "<strong>Warning:</strong> Do not include any print statement in the body of the function below. It perturbs the grading algorithm. Either add new cells, or make sure to remove any print statement before submitting your work.\n",
    "</div>"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "id": "2579535e",
   "metadata": {},
   "outputs": [],
   "source": [
    "def calculate_number_of_redundant_constraints(\n",
    "    matrix_a: np.array, vector_b: np.array\n",
    ") -> int:\n",
    "    number_of_rows, number_of_columns = matrix_a.shape\n",
    "    if number_of_rows >= number_of_columns:\n",
    "        raise ValueError(\n",
    "            f'The matrix must have less rows than columns. Size: ({number_of_rows}, {number_of_columns})'\n",
    "        )\n",
    "    rank = np.linalg.matrix_rank(matrix_a)\n",
    "    return number_of_rows - rank\n",
    "   \n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "id": "063132a6",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Number of redundant constraints: 1\n"
     ]
    }
   ],
   "source": [
    "# Test your function. For example, on the system defined by the following data:\n",
    "A = np.array([[1, 2, 3], [1, 2, 3]])\n",
    "b = np.array([6, 6])\n",
    "number_of_redundant_constraints = calculate_number_of_redundant_constraints(\n",
    "    matrix_a=A, vector_b=b\n",
    ")\n",
    "print(f'Number of redundant constraints: {number_of_redundant_constraints}')\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "id": "0ef6b532",
   "metadata": {},
   "outputs": [],
   "source": [
    "\"\"\" \n",
    "success_message: Testing A = np.array([[1, 2, 3], [4, 5, 6]]). No redundant constraint. OK.\n",
    "failure_message: Testing A = np.array([[1, 2, 3], [4, 5, 6]]). No redundant constraint. Failed. \n",
    "\"\"\" \n",
    "_A = np.array([[1, 2, 3], [4, 5, 6]])\n",
    "_b = np.array([0, 0])\n",
    "number_of_redundant_constraints = calculate_number_of_redundant_constraints(\n",
    "    matrix_a=_A, vector_b=_b\n",
    ")\n",
    "assert number_of_redundant_constraints == 0\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 20,
   "id": "a3356203-9b5c-4f7d-b1d9-585bc3bad4cf",
   "metadata": {},
   "outputs": [],
   "source": [
    "\"\"\" \n",
    "success_message: Testing A = np.array([[1, 2, 3], [1, 2, 3]]). 1 redundant constraint. OK. \n",
    "failure_message: Testing A = np.array([[1, 2, 3], [1, 2, 3]]). 1 redundant constraint. Failed. \n",
    "\"\"\" \n",
    "_A = np.array([[1, 2, 3], [1, 2, 3]])\n",
    "_b = np.array([0, 0])\n",
    "number_of_redundant_constraints = calculate_number_of_redundant_constraints(\n",
    "    matrix_a=_A, vector_b=_b\n",
    ")\n",
    "assert number_of_redundant_constraints == 1\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "id": "1dae303c-3266-4b78-90eb-5b4134d8e3ed",
   "metadata": {},
   "outputs": [],
   "source": [
    "\"\"\" \n",
    "success_message: Testing [[1, 2, 3, 4, 5, 6], [2, 3, 4, 5, 6, 7], [1, 2, 3, 4, 5, 6], [2, 3, 4, 5, 6, 7]]. 2 redundant constraints. OK. \n",
    "failure_message: Testing [[1, 2, 3, 4, 5, 6], [2, 3, 4, 5, 6, 7], [1, 2, 3, 4, 5, 6], [2, 3, 4, 5, 6, 7]]. 2 redundant constraints. Failed.  \n",
    "\"\"\" \n",
    "_A = np.array(\n",
    "    [[1, 2, 3, 4, 5, 6], [2, 3, 4, 5, 6, 7], [1, 2, 3, 4, 5, 6], [2, 3, 4, 5, 6, 7]]\n",
    ")\n",
    "_b = np.array([0, 0])\n",
    "number_of_redundant_constraints = calculate_number_of_redundant_constraints(\n",
    "    matrix_a=_A, vector_b=_b\n",
    ")\n",
    "assert number_of_redundant_constraints == 2\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "2ccf576a",
   "metadata": {},
   "source": [
    "# Question 4\n",
    "Consider a polyhedron written in standard form :\n",
    "$$\\mathcal{P} =\\{ x \\in \\mathbb{R}^n | Ax = b, x \\geq 0 \\}.$$\n",
    "\n",
    "Write a function that considers a list of indices and generates one of the following diagnostics:\n",
    "- 'NO_BASIS': the list of indices does not characterize a basic solution,\n",
    "- 'INFEASIBLE': the list of indices corresponds to a basic solution that is infeasible,\n",
    "- 'DEGENERATE': the list of indices corresponds to a basic solution that is feasible and degenerate,\n",
    "- 'NON_DEGENERATE': the list of indices corresponds to a basic solution that is feasible and non degenerate."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "f643feeb",
   "metadata": {},
   "source": [
    "The output of the function will be of the following type:"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "8ae87aaa",
   "metadata": {},
   "source": [
    "Now write the correct function"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "bbc5659a-e73e-4d5c-b620-77a162c58675",
   "metadata": {},
   "source": [
    "<div style=\"border-left: 5px solid #f44336; background-color: #fdecea; padding: 10px;\">\n",
    "<strong>Warning:</strong> Do not include any print statement in the body of the function below. It perturbs the grading algorithm. Either add new cells, or make sure to remove any print statement before submitting your work.\n",
    "</div>"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 22,
   "id": "239c6907",
   "metadata": {},
   "outputs": [],
   "source": [
    "def diagnostic_for_basic_solution(\n",
    "    matrix_a: np.array, vector_b: np.array, basic_indices: tuple[int, ...]\n",
    ") -> str:\n",
    "    \"\"\"\n",
    "    Identifies the type of basic solution that is defined by the list of indices.\n",
    "\n",
    "    :param matrix_a: A numpy ndarray representing the constraint matrix $A$.\n",
    "    :param vector_b: A numpy ndarray representing the right-hand side vector $b$.\n",
    "    :param basic_indices: A list of integers representing the column indices to form the basis.\n",
    "    :return: the diagnostic.\n",
    "    \"\"\"\n",
    "    n_rows, n_columns = matrix_a.shape\n",
    "    # Verify that each index is a valid column index of the matrix.\n",
    "    max_index = n_columns - 1\n",
    "    if not all(0 <= index <= max_index for index in basic_indices):\n",
    "        raise ValueError(\n",
    "            'One or more indices are out of the valid column index range of the matrix.'\n",
    "        )\n",
    "\n",
    "    if len(basic_indices) != n_rows:\n",
    "        return 'NO_BASIS'\n",
    "\n",
    "    basis_matrix = matrix_a[:, basic_indices]\n",
    "    try:\n",
    "        basic_variables = np.linalg.solve(basis_matrix, vector_b)\n",
    "    except np.linalg.LinAlgError:\n",
    "        # The system does not have a solution, as the matrix is not invertible.\n",
    "        return 'NO_BASIS'\n",
    "\n",
    "    is_feasible: bool = np.all(basic_variables >= 0)\n",
    "    is_degenerate: bool = np.any(basic_variables == 0)\n",
    "    if not is_feasible:\n",
    "        return 'INFEASIBLE'\n",
    "    if is_degenerate:\n",
    "        return 'DEGENERATE'\n",
    "    return 'NON_DEGENERATE'\n",
    "    \n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "e2ddcd72-17ca-4853-895d-2d95924b85c3",
   "metadata": {},
   "source": [
    "Test the function on the following example. Each diagnostic should be returned at least once."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 23,
   "id": "4949e1f1-1c7c-4397-b22f-04e806f94d2e",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "(0, 1, 2): NO_BASIS\n",
      "(0, 1, 3): DEGENERATE\n",
      "(0, 1, 4): INFEASIBLE\n",
      "(0, 2, 3): NON_DEGENERATE\n",
      "(0, 2, 4): NON_DEGENERATE\n",
      "(0, 3, 4): NO_BASIS\n",
      "(1, 2, 3): DEGENERATE\n",
      "(1, 2, 4): NON_DEGENERATE\n",
      "(1, 3, 4): DEGENERATE\n",
      "(2, 3, 4): INFEASIBLE\n"
     ]
    }
   ],
   "source": [
    "A = np.array([[0, 1, 1, 0, 0], [-1, -1, 0, 1, 0], [1, 1, 0, 0, 1]])\n",
    "b = np.array([2, -1, 2])\n",
    "\n",
    "valid_outputs = ['NO_BASIS', 'INFEASIBLE','DEGENERATE', 'NON_DEGENERATE']\n",
    "\n",
    "n_rows, n_columns = A.shape\n",
    "all_potential_bases: list[tuple[int, ...]] = list(\n",
    "    combinations(range(n_columns), n_rows)\n",
    ")\n",
    "for potential_basis in all_potential_bases:\n",
    "    diagnostic = diagnostic_for_basic_solution(\n",
    "        matrix_a=A, vector_b=b, basic_indices=potential_basis\n",
    "    )\n",
    "    print(f'{potential_basis}: {diagnostic}')\n",
    "    if diagnostic not in valid_outputs:\n",
    "        error_msg = f'Invalid output: {diagnostic}. It must be exactly one of the following (all upper case): {valid_outputs}'\n",
    "        raise ValueError(error_msg)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 24,
   "id": "9686b8c8-2c30-4725-a8ac-f4d02a084c7c",
   "metadata": {},
   "outputs": [],
   "source": [
    "\"\"\" \n",
    "success_message: The number of indices was insufficient. It has been correctly detected. \n",
    "failure_message: The number of indices was insufficient. It has not been correctly detected.  \n",
    "\"\"\" \n",
    "\n",
    "_A = np.array([[0, 1, 1, 0, 0], [-1, -1, 0, 1, 0], [1, 1, 0, 0, 1]])\n",
    "_b = np.array([2, -1, 2])\n",
    "\n",
    "_potential_basis = (0, 1)\n",
    "diagnostic = diagnostic_for_basic_solution(\n",
    "    matrix_a=_A, vector_b=_b, basic_indices=_potential_basis\n",
    ")\n",
    "assert diagnostic == 'NO_BASIS'"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 25,
   "id": "f9f077c2-6ddd-4023-8a18-ce7686570a68",
   "metadata": {},
   "outputs": [],
   "source": [
    "\"\"\"\n",
    "success_message: The basic matrix is singular. It has been correctly detected. \n",
    "failure_message: The basic matrix is singular. It has not been correctly detected. \n",
    "\"\"\" \n",
    "\n",
    "_A = np.array([[0, 1, 1, 0, 0], [-1, -1, 0, 1, 0], [1, 1, 0, 0, 1]])\n",
    "_b = np.array([2, -1, 2])\n",
    "\n",
    "_potential_basis = (0, 1, 2)\n",
    "diagnostic = diagnostic_for_basic_solution(\n",
    "    matrix_a=_A, vector_b=_b, basic_indices=_potential_basis\n",
    ")\n",
    "assert diagnostic == 'NO_BASIS'"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 26,
   "id": "0749968d-dbf2-41f9-a01c-14d02b5e7e69",
   "metadata": {},
   "outputs": [],
   "source": [
    "\"\"\" \n",
    "success_message: The basic solution is feasible and degenerate. It has been correctly detected. \n",
    "failure_message: The basic solution is feasible and degenerate. It has not been correctly detected. \n",
    "\"\"\" \n",
    "_A = np.array([[0, 1, 1, 0, 0], [-1, -1, 0, 1, 0], [1, 1, 0, 0, 1]])\n",
    "_b = np.array([2, -1, 2])\n",
    "\n",
    "_potential_basis = (0, 1, 3)\n",
    "diagnostic = diagnostic_for_basic_solution(\n",
    "    matrix_a=_A, vector_b=_b, basic_indices=_potential_basis\n",
    ")\n",
    "assert diagnostic == 'DEGENERATE'\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 27,
   "id": "36a20b64-9aa6-4a6b-b909-171ce590d23f",
   "metadata": {},
   "outputs": [],
   "source": [
    "\"\"\" \n",
    "success_message: The basic solution is infeasible. It has been correctly detected. \n",
    "failure_message: The basic solution is infeasible. It has not been correctly detected.\n",
    "\"\"\" \n",
    "\n",
    "_A = np.array([[0, 1, 1, 0, 0], [-1, -1, 0, 1, 0], [1, 1, 0, 0, 1]])\n",
    "_b = np.array([2, -1, 2])\n",
    "\n",
    "_potential_basis = (0, 1, 4)\n",
    "diagnostic = diagnostic_for_basic_solution(\n",
    "    matrix_a=_A, vector_b=_b, basic_indices=_potential_basis\n",
    ")\n",
    "assert diagnostic == 'INFEASIBLE'\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 28,
   "id": "b0a87f63-f078-4556-b964-241f2b766664",
   "metadata": {},
   "outputs": [],
   "source": [
    "\"\"\" \n",
    "success_message: The basic solution is feasible and non degenerate. It has been correctly detected. \n",
    "failure_message: The basic solution is feasible and non degenerate. It has not been correctly detected. \n",
    "\"\"\" \n",
    "\n",
    "_A = np.array([[0, 1, 1, 0, 0], [-1, -1, 0, 1, 0], [1, 1, 0, 0, 1]])\n",
    "_b = np.array([2, -1, 2])\n",
    "\n",
    "_potential_basis = (0, 2, 3)\n",
    "diagnostic = diagnostic_for_basic_solution(\n",
    "    matrix_a=_A, vector_b=_b, basic_indices=_potential_basis\n",
    ")\n",
    "assert diagnostic == 'NON_DEGENERATE'"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "fa09ff6c",
   "metadata": {},
   "source": [
    "# Question 5\n",
    "Consider a polyhedron written in standard form:\n",
    "$$\\mathcal{P} =\\{ x \\in \\mathbb{R}^n | Ax = b, x \\geq 0 \\},$$\n",
    "a list of indices corresponding to a feasible basic solution, and the index of a non basic variable.\n",
    "Write a function that returns one basic direction,that is, an array of length $n$."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "cae9c2d8-cc8f-4c4f-b481-418be028c322",
   "metadata": {},
   "source": [
    "<div style=\"border-left: 5px solid #f44336; background-color: #fdecea; padding: 10px;\">\n",
    "<strong>Warning:</strong> Do not include any print statement in the body of the function below. It perturbs the grading algorithm. Either add new cells, or make sure to remove any print statement before submitting your work.\n",
    "</div>"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 29,
   "id": "268db47b",
   "metadata": {},
   "outputs": [],
   "source": [
    "def identify_one_basic_direction(\n",
    "    matrix_a: np.array, basic_indices: tuple[int, ...], non_basic_index: int\n",
    ") -> np.array:\n",
    "    n_rows, n_columns = matrix_a.shape\n",
    "    basis_matrix = matrix_a[:, basic_indices]\n",
    "    non_basic_column = matrix_a[:, non_basic_index]\n",
    "    the_basic_part_of_the_direction = np.linalg.solve(basis_matrix, -non_basic_column)\n",
    "    the_basic_direction = np.zeros(n_columns)\n",
    "    the_basic_direction[non_basic_index] = 1.0\n",
    "    for variable, index in zip(the_basic_part_of_the_direction, basic_indices):\n",
    "        the_basic_direction[index] = variable\n",
    "    return the_basic_direction\n",
    "    \n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "78e13f5b-7b36-4c76-b7aa-ae5ba34b45d7",
   "metadata": {},
   "source": [
    "Test the function on the following example"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 30,
   "id": "06ddaad7-9f11-4067-87bc-b8c7b824bd2c",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Basic direction 0: [ 1. -1.  0.  0.]\n",
      "Basic direction 3: [ 0. -1. -1.  1.]\n"
     ]
    }
   ],
   "source": [
    "A = np.array([[-1, -1, 1, 0], [1, 1, 0, 1]])\n",
    "b = np.array([-1, 2])\n",
    "the_basic_indices = (1, 2)\n",
    "the_non_basic_indices = (0, 3)\n",
    "for non_basic_index in the_non_basic_indices:\n",
    "    basic_direction = identify_one_basic_direction(\n",
    "        matrix_a=A, basic_indices=the_basic_indices, non_basic_index=non_basic_index\n",
    "    )\n",
    "    print(f'Basic direction {non_basic_index}: {basic_direction}')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 31,
   "id": "a503a231",
   "metadata": {},
   "outputs": [],
   "source": [
    "\"\"\" \n",
    "success_message: The non basic part of the direction is correct. \n",
    "failure_message: The non basic part of the direction is incorrect. It should contain only zero's, except for the entry corresponding to the selected non basic variable, that should contain a 1.  \n",
    "\"\"\" \n",
    "_grading_A = np.array(\n",
    "    [\n",
    "        [1, 2, 0, 0, 0, 1, 0, 0, 0, 0],\n",
    "        [0, 1, 3, 0, 0, 0, 1, 0, 0, 0],\n",
    "        [0, 0, 0, 2, 4, 0, 0, 1, 0, 0],\n",
    "        [0, 0, 0, 0, 1, 0, 0, 0, 2, 1],\n",
    "        [1, 0, 0, 0, 0, 0, 0, 0, 1, 0],\n",
    "    ]\n",
    ")\n",
    "_grading_b = np.array([4, 6, 10, 8, 5])\n",
    "_grading_c = np.array([-1, 2, -3, 4, -5, 6, -7, 8, -9, 10])\n",
    "_n_rows, _n_columns = _grading_A.shape\n",
    "_grading_basic_indices = (0, 1, 2, 4, 8)\n",
    "_grading_non_basic_index = 3\n",
    "the_basic_direction = identify_one_basic_direction(\n",
    "    matrix_a=_grading_A,\n",
    "    basic_indices=_grading_basic_indices,\n",
    "    non_basic_index=_grading_non_basic_index,\n",
    ")\n",
    "assert the_basic_direction[_grading_non_basic_index] == 1.0\n",
    "\n",
    "_non_zero_values_indices = set(_grading_basic_indices) | {_grading_non_basic_index}\n",
    "_zero_values_indices = set(range(n_columns)) - _non_zero_values_indices\n",
    "_zero_values = the_basic_direction[list(_zero_values_indices)]\n",
    "assert np.all(np.isclose(_zero_values, 0))\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 32,
   "id": "71fce670-3c93-4e8a-bee7-b76e8126588b",
   "metadata": {},
   "outputs": [],
   "source": [
    "\"\"\" \n",
    "success_message: The basic part of the direction is correct.\n",
    "failure_message: The non bsasic part of the direction is incorrect.\n",
    "\"\"\" \n",
    "_grading_A = np.array(\n",
    "    [\n",
    "        [1, 2, 0, 0, 0, 1, 0, 0, 0, 0],\n",
    "        [0, 1, 3, 0, 0, 0, 1, 0, 0, 0],\n",
    "        [0, 0, 0, 2, 4, 0, 0, 1, 0, 0],\n",
    "        [0, 0, 0, 0, 1, 0, 0, 0, 2, 1],\n",
    "        [1, 0, 0, 0, 0, 0, 0, 0, 1, 0],\n",
    "    ]\n",
    ")\n",
    "_grading_b = np.array([4, 6, 10, 8, 5])\n",
    "_grading_c = np.array([-1, 2, -3, 4, -5, 6, -7, 8, -9, 10])\n",
    "_n_rows, _n_columns = _grading_A.shape\n",
    "_grading_basic_indices = (0, 1, 2, 4, 8)\n",
    "_grading_non_basic_index = 3\n",
    "the_basic_direction = identify_one_basic_direction(\n",
    "    matrix_a=_grading_A,\n",
    "    basic_indices=_grading_basic_indices,\n",
    "    non_basic_index=_grading_non_basic_index,\n",
    ")\n",
    "\n",
    "_basic_part = the_basic_direction[list(_grading_basic_indices)]\n",
    "_expected_result = np.array([-1 / 4, 1 / 8, -1 / 24, -1 / 2, 1 / 4])\n",
    "assert np.all(np.isclose(_basic_part, _expected_result))"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "6c3598df",
   "metadata": {},
   "source": [
    "# Question 6\n",
    "Consider a linear optimization problem:\n",
    "$$\\min_{x\\in\\mathbb{R}^n } c^T x$$\n",
    "subject to\n",
    "$$Ax = b, x \\geq 0,$$\n",
    "a list of indices corresponding to a feasible basic solution, and the index of a non basix variable.\n",
    "Write a function that calculates the corresponding reduced cost."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "e724725d-6c47-44b3-b78b-c4c941f81d03",
   "metadata": {},
   "source": [
    "<div style=\"border-left: 5px solid #f44336; background-color: #fdecea; padding: 10px;\">\n",
    "<strong>Warning:</strong> Do not include any print statement in the body of the function below. It perturbs the grading algorithm. Either add new cells, or make sure to remove any print statement before submitting your work.\n",
    "</div>"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 33,
   "id": "a7fac858",
   "metadata": {},
   "outputs": [],
   "source": [
    "def calculate_reduced_cost(\n",
    "    matrix_a: np.array,\n",
    "    vector_c: np.array,\n",
    "    basic_indices: tuple[int, ...],\n",
    "    non_basic_index: int,\n",
    ") -> float:\n",
    "    basic_direction = identify_one_basic_direction(\n",
    "        matrix_a=matrix_a, basic_indices=basic_indices, non_basic_index=non_basic_index\n",
    "    )\n",
    "    return np.inner(basic_direction, vector_c)\n",
    "    \n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "fbe6d216-c177-4b9a-892d-c06aa1f2e66b",
   "metadata": {},
   "source": [
    "Test the function on the following example:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 34,
   "id": "84151e02-2eff-4d90-8f97-2f8d197b4571",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Reduced cost for 0: -3.0\n",
      "Reduced cost for 3: -1.0\n"
     ]
    }
   ],
   "source": [
    "A = np.array([[-1, -1, 1, 0], [1, 1, 0, 1]])\n",
    "c = np.array([-1, 2, 3, 4])\n",
    "the_basic_indices = (1, 2)\n",
    "the_non_basic_indices = (0, 3)\n",
    "for non_basic_index in the_non_basic_indices:\n",
    "    reduced_cost = calculate_reduced_cost(\n",
    "        matrix_a=A,\n",
    "        vector_c=c,\n",
    "        basic_indices=the_basic_indices,\n",
    "        non_basic_index=non_basic_index,\n",
    "    )\n",
    "    print(f'Reduced cost for {non_basic_index}: {reduced_cost}')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 35,
   "id": "3ade5425",
   "metadata": {},
   "outputs": [],
   "source": [
    "\"\"\" \n",
    "success_message: Test on a m=5, n=10 problem. OK.  \n",
    "failure_message: Test on a m=5, n=10 problem. Failed.  \n",
    "\"\"\" \n",
    "_grading_A = np.array(\n",
    "    [\n",
    "        [1, 2, 0, 0, 0, 1, 0, 0, 0, 0],\n",
    "        [0, 1, 3, 0, 0, 0, 1, 0, 0, 0],\n",
    "        [0, 0, 0, 2, 4, 0, 0, 1, 0, 0],\n",
    "        [0, 0, 0, 0, 1, 0, 0, 0, 2, 1],\n",
    "        [1, 0, 0, 0, 0, 0, 0, 0, 1, 0],\n",
    "    ]\n",
    ")\n",
    "_grading_b = np.array([4, 6, 10, 8, 5])\n",
    "_grading_c = np.array([-1, 2, -3, 4, -5, 6, -7, 8, -9, 10])\n",
    "_n_rows, _n_columns = _grading_A.shape\n",
    "_grading_basic_indices = (0, 1, 2, 4, 8)\n",
    "_grading_non_basic_index = 3\n",
    "_the_reduced_cost = calculate_reduced_cost(\n",
    "    matrix_a=_grading_A,\n",
    "    vector_c=_grading_c,\n",
    "    basic_indices=_grading_basic_indices,\n",
    "    non_basic_index=_grading_non_basic_index,\n",
    ")\n",
    "assert np.isclose(_the_reduced_cost, 4.875)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "f59b7a20",
   "metadata": {},
   "source": [
    "# Question 7\n",
    "Consider a simplex tableau. During an iteration of the simplex algorithm, it is necessary to identify the row and\n",
    "the column of the pivot. Your task is to write a function that does that.\n",
    "If the row or the column cannot be identified, return None. Whenever there are several potential candidates for the row or the column, you must use the rule seen in class to identify only one."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "9bf732d7-ff12-465c-b3b5-1b36917761c8",
   "metadata": {},
   "source": [
    "**Remember that, in Python, the numbering starts at 0.** It means that, if you have $m$ rows and $n$ columns, they are numbered from $0$ to $m-1$ and from $0$ to $n-1$."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "6e1c30f5-c6c2-42c0-8086-2a494de2ca53",
   "metadata": {},
   "source": [
    "<div style=\"border-left: 5px solid #f44336; background-color: #fdecea; padding: 10px;\">\n",
    "<strong>Warning:</strong> Do not include any print statement in the body of the function below. It perturbs the grading algorithm. Either add new cells, or make sure to remove any print statement before submitting your work.\n",
    "</div>"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 36,
   "id": "fc6b46d1",
   "metadata": {},
   "outputs": [],
   "source": [
    "def identify_pivot(tableau: np.array) -> tuple[int | None, int | None]:\n",
    "    n_rows, n_columns = tableau.shape\n",
    "    n_variables = n_columns - 1\n",
    "    n_constraints = n_rows - 1\n",
    "    last_row = tableau[-1, :-1]\n",
    "    negative_indices = np.where(last_row < 0)[0]\n",
    "    if negative_indices.size == 0:\n",
    "        return None, None\n",
    "    pivot_column = int(negative_indices[0])\n",
    "    last_column = tableau[:-1, -1]\n",
    "    pivot_col_values = tableau[\n",
    "        :-1, pivot_column\n",
    "    ]  # Values in the pivot column (excluding the last row)\n",
    "\n",
    "    # Perform minimum ratio test, considering only positive pivot column values\n",
    "    count_positive_entries = np.sum(pivot_col_values > 0)\n",
    "    if count_positive_entries == 0:\n",
    "        return None, pivot_column\n",
    "    ratios = np.divide(last_column, pivot_col_values, where=pivot_col_values > 0)\n",
    "    ratios[pivot_col_values <= 0] = (\n",
    "        np.inf\n",
    "    )  # Ignore non-positive values in the pivot column for the ratio test\n",
    "\n",
    "    pivot_row = np.argmin(ratios)\n",
    "    return pivot_row, pivot_column\n",
    "    \n",
    "\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "e718d06b-aefc-4e9c-99e6-e55a6f40f0dc",
   "metadata": {},
   "source": [
    "Test the function on this example:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 37,
   "id": "c73aad4f-6faf-49d0-937a-1a80c7a119f6",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Pivot at row 2 and column 0\n"
     ]
    }
   ],
   "source": [
    "tableau = np.array(\n",
    "    [\n",
    "        [2, 3, 1, 0, 0, 30],\n",
    "        [-1, 2, 0, 1, 0, 20],\n",
    "        [3, 1, 0, 0, 1, 40],\n",
    "        [-3, -5, 0, 0, 0, 0],\n",
    "    ]\n",
    ")\n",
    "row_pivot, column_pivot = identify_pivot(tableau)\n",
    "print(f'Pivot at row {row_pivot} and column {column_pivot}')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 38,
   "id": "e28bc9a5-5c85-4dea-8695-d944167cefcb",
   "metadata": {},
   "outputs": [],
   "source": [
    "\"\"\" \n",
    "success_message: The column of the pivot corresponds to a negative reduced cost\n",
    "failure_message: The column of the pivot does not correspond to a negative reduced cost\n",
    "\"\"\" \n",
    "_grading_A = np.array(\n",
    "    [\n",
    "        [1, 2, 0, 0, 0, 1, 0, 0, 0, 0],\n",
    "        [0, 1, 3, 0, 0, 0, 1, 0, 0, 0],\n",
    "        [0, 0, 0, 2, 4, 0, 0, 1, 0, 0],\n",
    "        [0, 0, 0, 0, 1, 0, 0, 0, 2, 1],\n",
    "        [1, 0, 0, 0, 0, 0, 0, 0, 1, 0],\n",
    "    ]\n",
    ")\n",
    "_grading_b = np.array([4, 6, 10, 8, 5])\n",
    "_grading_c = np.array([-1, 2, -3, 4, -5, 6, -7, 8, -9, 10])\n",
    "_grading_basic_indices = (0, 1, 2, 4, 8)\n",
    "_grading_non_basic_index = 3\n",
    "_n_rows, _n_columns = _grading_A.shape\n",
    "_basic_matrix = _grading_A[:, _grading_basic_indices]\n",
    "_top_left = np.linalg.solve(_basic_matrix, _grading_A)\n",
    "_top_right = np.linalg.solve(_basic_matrix, _grading_b)\n",
    "_c_b = _grading_c[list(_grading_basic_indices)]\n",
    "_bottom_left = _grading_c.T - _c_b.T @ _top_left\n",
    "_bottom_right = -np.inner(_c_b, _top_right)\n",
    "_Ab = np.column_stack((_top_left, _top_right))\n",
    "_tableau = np.vstack((_Ab, np.append(_bottom_left, _bottom_right)))\n",
    "_tableau[-1][7] = -0.001\n",
    "_row_pivot, _column_pivot = identify_pivot(_tableau)\n",
    "_reduced_cost = _tableau[-1][_column_pivot]\n",
    "assert _reduced_cost < 0\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 39,
   "id": "d976342d-86ce-44b4-b357-635ccdf1d01a",
   "metadata": {},
   "outputs": [],
   "source": [
    "\"\"\" \n",
    "success_message: The column of the pivot has been selected using Bland's rule\n",
    "failure_message: The column of the pivot has not been selected using Bland's rule\n",
    "\"\"\" \n",
    "_grading_A = np.array(\n",
    "    [\n",
    "        [1, 2, 0, 0, 0, 1, 0, 0, 0, 0],\n",
    "        [0, 1, 3, 0, 0, 0, 1, 0, 0, 0],\n",
    "        [0, 0, 0, 2, 4, 0, 0, 1, 0, 0],\n",
    "        [0, 0, 0, 0, 1, 0, 0, 0, 2, 1],\n",
    "        [1, 0, 0, 0, 0, 0, 0, 0, 1, 0],\n",
    "    ]\n",
    ")\n",
    "_grading_b = np.array([4, 6, 10, 8, 5])\n",
    "_grading_c = np.array([-1, 2, -3, 4, -5, 6, -7, 8, -9, 10])\n",
    "_grading_basic_indices = (0, 1, 2, 4, 8)\n",
    "_grading_non_basic_index = 3\n",
    "_n_rows, _n_columns = _grading_A.shape\n",
    "_basic_matrix = _grading_A[:, _grading_basic_indices]\n",
    "_top_left = np.linalg.solve(_basic_matrix, _grading_A)\n",
    "_top_right = np.linalg.solve(_basic_matrix, _grading_b)\n",
    "_c_b = _grading_c[list(_grading_basic_indices)]\n",
    "_bottom_left = _grading_c.T - _c_b.T @ _top_left\n",
    "_bottom_right = -np.inner(_c_b, _top_right)\n",
    "_Ab = np.column_stack((_top_left, _top_right))\n",
    "_tableau = np.vstack((_Ab, np.append(_bottom_left, _bottom_right)))\n",
    "_tableau[-1][7] = -0.001\n",
    "_row_pivot, _column_pivot = identify_pivot(_tableau)\n",
    "assert _column_pivot == 6\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 40,
   "id": "b3676253-6ebe-4b35-b369-a143a1fc8b1a",
   "metadata": {},
   "outputs": [],
   "source": [
    "\"\"\" \n",
    "success_message: The row of the pivot corresponds to the minimum step along the direction.\n",
    "failure_message: The row of the pivot does not correspond to the minimum step along the direction.\n",
    "\"\"\" \n",
    "_grading_A = np.array(\n",
    "    [\n",
    "        [1, 2, 0, 0, 0, 1, 0, 0, 0, 0],\n",
    "        [0, 1, 3, 0, 0, 0, 1, 0, 0, 0],\n",
    "        [0, 0, 0, 2, 4, 0, 0, 1, 0, 0],\n",
    "        [0, 0, 0, 0, 1, 0, 0, 0, 2, 1],\n",
    "        [1, 0, 0, 0, 0, 0, 0, 0, 1, 0],\n",
    "    ]\n",
    ")\n",
    "_grading_b = np.array([4, 6, 10, 8, 5])\n",
    "_grading_c = np.array([-1, 2, -3, 4, -5, 6, -7, 8, -9, 10])\n",
    "_grading_basic_indices = (0, 1, 2, 4, 8)\n",
    "_grading_non_basic_index = 3\n",
    "_n_rows, _n_columns = _grading_A.shape\n",
    "_basic_matrix = _grading_A[:, _grading_basic_indices]\n",
    "_top_left = np.linalg.solve(_basic_matrix, _grading_A)\n",
    "_top_right = np.linalg.solve(_basic_matrix, _grading_b)\n",
    "_c_b = _grading_c[list(_grading_basic_indices)]\n",
    "_bottom_left = _grading_c.T - _c_b.T @ _top_left\n",
    "_bottom_right = -np.inner(_c_b, _top_right)\n",
    "_Ab = np.column_stack((_top_left, _top_right))\n",
    "_tableau = np.vstack((_Ab, np.append(_bottom_left, _bottom_right)))\n",
    "_row_pivot, column_pivot = identify_pivot(_tableau)\n",
    "assert _row_pivot == 2"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 41,
   "id": "91e56212-b135-4aed-90bf-e1813c8fe7eb",
   "metadata": {},
   "outputs": [],
   "source": [
    "\"\"\" \n",
    "success_message: All the elements in the column of the pivot are non positive. The problem is unbounded. It has correctly been detected.  \n",
    "failure_message: All the elements in the column of the pivot are non positive. The problem is unbounded. It has not correctly been detected.\n",
    "\"\"\" \n",
    "_tableau = np.array(\n",
    "    [\n",
    "        [1, 2, 0, 0, 0, 1, 0, 0, 0, 0],\n",
    "        [0, 1, 3, 0, 0, 0, -1, 0, 0, 0],\n",
    "        [0, 0, 0, 2, 4, 0, -1, 1, 0, 0],\n",
    "        [0, 0, 0, 0, 1, 0, 0, 0, 2, 1],\n",
    "        [1, 0, 0, 0, 0, 0, -1, 0, 1, 0],\n",
    "    ]\n",
    ")\n",
    "_row_pivot, _column_pivot = identify_pivot(_tableau)\n",
    "assert _row_pivot is None\n",
    "assert _column_pivot == 6\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 42,
   "id": "1ab23be4-9b69-41e1-be41-50528605a34c",
   "metadata": {},
   "outputs": [],
   "source": [
    "\"\"\" \n",
    "success_message: All the reduded costs are positive. The solution is optimal. It has correctly been detected. \n",
    "failure_message: All the reduded costs are positive. The solution is optimal. It has not correctly been detected. \n",
    "\"\"\" \n",
    "_tableau = np.array(\n",
    "    [\n",
    "        [1, 2, 0, 0, 0, 1, 0, 0, 0, 0],\n",
    "        [0, 1, 3, 0, 0, 0, 1, 0, 0, 0],\n",
    "        [0, 0, 0, 2, 4, 0, 0, 1, 0, 0],\n",
    "        [0, 0, 0, 0, 1, 0, 0, 0, 2, 1],\n",
    "        [1, 0, 0, 0, 0, 0, 0, 0, 1, 0],\n",
    "    ]\n",
    ")\n",
    "_tableau[-1, :] = 1\n",
    "_row_pivot, _column_pivot = identify_pivot(_tableau)\n",
    "assert _row_pivot is None\n",
    "assert _column_pivot is None"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "4ebfaadd-107c-4de3-b3c4-64b528f28097",
   "metadata": {},
   "source": [
    "## Before submitting"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "4f4d8e2f-3222-45fe-bdaa-55dd91e91f7a",
   "metadata": {},
   "source": [
    "Make sure you have run all cells in your notebook. To do so, click on the double arrow `Restart the kernel and run all cells`. If any exception is raised during the execution, the notebook will not be corrected. You can always insert an ellipsis (that is, three consecutive dots `...`) to represent an empty code, if you don't know what do answer."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "2430c79a-19f2-4643-8449-89bbe7a45369",
   "metadata": {},
   "source": [
    "When you are done, upload your notebook on moodle. "
   ]
  }
 ],
 "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.13.1"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}
