optimization/01_linear_program.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Simple LP (Linear Program) Problem\n",
    "\n",
    "Reference: https://developers.google.com/optimization/lp/lp_example\n",
    "\n",
    "\\begin{align*}\n",
    "\\text{max} \\quad & 3x + 4y \\\\\n",
    "\\text{s.t} \\quad & x + 2y \\leq 14 \\\\\n",
    "                        & 3x - y \\geq 0 \\\\\n",
    "                        & x - y \\leq 2 \\\\\n",
    "                        & x \\geq 0 \\\\\n",
    "                        & y \\geq 0\n",
    "\\end{align*}\n",
    "\n",
    "<div style=\"text-align: center;\">\n",
    "  <img src=\"https://developers.google.com/static/optimization/images/mip/feasible_region_sol.png\">\n",
    "</div>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## SciPy"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Optimal value: 34.0\n",
      "Optimal x: 6.0\n",
      "Optimal y: 4.0\n"
     ]
    }
   ],
   "source": [
    "import numpy as np\n",
    "from scipy.optimize import linprog\n",
    "\n",
    "# Objective Function:\n",
    "c = [-3, -4]  # Coefficients for the objective function (negated for maximization)\n",
    "\n",
    "# Constraints:\n",
    "A = [\n",
    "    [1, 2],\n",
    "    [-3, 1],\n",
    "    [1, -1]\n",
    "]\n",
    "b = [14, 0, 2]\n",
    "\n",
    "# Bounds for variables x and y\n",
    "x_bounds = (0, None)\n",
    "y_bounds = (0, None)\n",
    "\n",
    "# Solve:\n",
    "result = linprog(c, A_ub=A, b_ub=b, bounds=[x_bounds, y_bounds], method='highs')\n",
    "\n",
    "if result.success:\n",
    "    print(f\"Optimal value: {-result.fun}\")  # Negate the result.fun to get the maximized value\n",
    "    print(f\"Optimal x: {result.x[0]}\")\n",
    "    print(f\"Optimal y: {result.x[1]}\")\n",
    "else:\n",
    "    print(\"The problem does not have an optimal solution.\")\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## CVXPY"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Optimal value: 33.999999439816996\n",
      "Optimal x: 5.999999994270371\n",
      "Optimal y: 3.9999998642514707\n"
     ]
    }
   ],
   "source": [
    "import cvxpy as cp\n",
    "\n",
    "# Variable Definition:\n",
    "x = cp.Variable()\n",
    "y = cp.Variable()\n",
    "\n",
    "# Objective Function:\n",
    "objective = cp.Maximize(3*x + 4*y)\n",
    "\n",
    "# Constraints:\n",
    "constraints = [\n",
    "    x + 2*y <= 14,\n",
    "    3*x - y >= 0,\n",
    "    x - y <= 2,\n",
    "    x >= 0,\n",
    "    y >= 0\n",
    "]\n",
    "\n",
    "# Define the problem:\n",
    "problem = cp.Problem(objective, constraints)\n",
    "\n",
    "# Solve:\n",
    "problem.solve(solver=cp.SCS)\n",
    "\n",
    "if problem.status == cp.OPTIMAL:\n",
    "    print(f\"Optimal value: {problem.value}\")\n",
    "    print(f\"Optimal x: {x.value}\")\n",
    "    print(f\"Optimal y: {y.value}\")\n",
    "else:\n",
    "    print(\"The problem does not have an optimal solution.\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## PuLP"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Optimal value: 34.0\n",
      "Optimal x: 6.0\n",
      "Optimal y: 4.0\n"
     ]
    }
   ],
   "source": [
    "import pulp\n",
    "\n",
    "# Define the problem:\n",
    "problem = pulp.LpProblem(\"MaximizeObjective\", pulp.LpMaximize)\n",
    "\n",
    "# Variable Definition:\n",
    "x = pulp.LpVariable(\"x\", lowBound=0)\n",
    "y = pulp.LpVariable(\"y\", lowBound=0)\n",
    "\n",
    "# Objective Function:\n",
    "problem += 3*x + 4*y, \"Objective\"\n",
    "\n",
    "# Constraints:\n",
    "problem += x + 2*y <= 14, \"Constraint 1\"\n",
    "problem += 3*x - y >= 0, \"Constraint 2\"\n",
    "problem += x - y <= 2, \"Constraint 3\"\n",
    "\n",
    "# Solve:\n",
    "solver = pulp.PULP_CBC_CMD()\n",
    "status = problem.solve()\n",
    "\n",
    "if pulp.LpStatus[status] == \"Optimal\":\n",
    "    print(f\"Optimal value: {pulp.value(problem.objective)}\")\n",
    "    print(f\"Optimal x: {pulp.value(x)}\")\n",
    "    print(f\"Optimal y: {pulp.value(y)}\")\n",
    "else:\n",
    "    print(\"The problem does not have an optimal solution.\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Pyomo"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Optimal value: 34.0\n",
      "Optimal x: 6.0\n",
      "Optimal y: 4.0\n"
     ]
    }
   ],
   "source": [
    "import pyomo.environ as pyo\n",
    "\n",
    "# Ceate a model:\n",
    "model = pyo.ConcreteModel()\n",
    "\n",
    "# Variable Definition:\n",
    "model.x = pyo.Var(within=pyo.NonNegativeReals)\n",
    "model.y = pyo.Var(within=pyo.NonNegativeReals)\n",
    "\n",
    "# Objective Function:\n",
    "model.obj = pyo.Objective(expr=3*model.x + 4*model.y, sense=pyo.maximize)\n",
    "\n",
    "# Constraints:\n",
    "model.constraint1 = pyo.Constraint(expr=model.x + 2*model.y <= 14)\n",
    "model.constraint2 = pyo.Constraint(expr=3*model.x - model.y >= 0)\n",
    "model.constraint3 = pyo.Constraint(expr=model.x - model.y <= 2)\n",
    "\n",
    "# Solve:\n",
    "solver = pyo.SolverFactory('glpk')\n",
    "result = solver.solve(model)\n",
    "\n",
    "if (result.solver.status == pyo.SolverStatus.ok) and (result.solver.termination_condition == pyo.TerminationCondition.optimal):\n",
    "    print(f\"Optimal value: {pyo.value(model.obj)}\")\n",
    "    print(f\"Optimal x: {pyo.value(model.x)}\")\n",
    "    print(f\"Optimal y: {pyo.value(model.y)}\")\n",
    "else:\n",
    "    print(\"The problem does not have an optimal solution.\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Google OR-Tools"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Optimal value: 33.99999999999999\n",
      "Optimal x: 5.999999999999998\n",
      "Optimal y: 3.9999999999999996\n"
     ]
    }
   ],
   "source": [
    "from ortools.linear_solver import pywraplp\n",
    "\n",
    "# Create a solver\n",
    "solver = pywraplp.Solver.CreateSolver('GLOP')\n",
    "\n",
    "# Variable Definition:\n",
    "x = solver.NumVar(0, solver.infinity(), 'x')\n",
    "y = solver.NumVar(0, solver.infinity(), 'y')\n",
    "\n",
    "# Objective Function:\n",
    "solver.Maximize(3 * x + 4 * y)\n",
    "\n",
    "# Constraints:\n",
    "solver.Add(x + 2 * y <= 14)\n",
    "solver.Add(3 * x - y >= 0)\n",
    "solver.Add(x - y <= 2)\n",
    "\n",
    "# Solve:\n",
    "status = solver.Solve()\n",
    "\n",
    "# Print results\n",
    "if status == pywraplp.Solver.OPTIMAL:\n",
    "    print(f\"Optimal value: {solver.Objective().Value()}\")\n",
    "    print(f\"Optimal x: {x.solution_value()}\")\n",
    "    print(f\"Optimal y: {y.solution_value()}\")\n",
    "else:\n",
    "    print(\"The problem does not have an optimal solution.\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Gekko"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Optimal value: 33.999999122\n",
      "Optimal x: 5.9999995199\n",
      "Optimal y: 4.0000001405\n"
     ]
    }
   ],
   "source": [
    "from gekko import GEKKO\n",
    "\n",
    "# Create a model:\n",
    "m = GEKKO(remote=False)\n",
    "\n",
    "# Variable Definition:\n",
    "x = m.Var(lb=0, value=0)\n",
    "y = m.Var(lb=0, value=0)\n",
    "\n",
    "# Objective Function:\n",
    "m.Maximize(3 * x + 4 * y)\n",
    "\n",
    "# Constraints:\n",
    "m.Equation(x + 2 * y <= 14)\n",
    "m.Equation(3 * x - y >= 0)\n",
    "m.Equation(x - y <= 2)\n",
    "\n",
    "# Solve:\n",
    "m.solve(disp=False)\n",
    "\n",
    "# Check solution:\n",
    "if m.options.APPSTATUS == 1:\n",
    "    print(f\"Optimal value: {m.options.OBJFCNVAL * -1}\") # Negate the result.fun to get the maximized value\n",
    "    print(f\"Optimal x: {x.value[0]}\")\n",
    "    print(f\"Optimal y: {y.value[0]}\")\n",
    "else:\n",
    "    print(\"The problem does not have an optimal solution.\")"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": ".venv",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.10.11"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}