example_optimize_simple2.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Simple optimization example (using BOBYQA)\n",
    "\n",
    "The waiting times and the number of operators are each assigned a cost. The search is for the number of operators for which the costs become minimal.\n",
    "\n",
    "[BOBYQA documentation](https://numericalalgorithmsgroup.github.io/pybobyqa/build/html/index.html)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Importing modules"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Processing results arrays\n",
    "import numpy as np\n",
    "\n",
    "# Optimizer\n",
    "import pybobyqa\n",
    "\n",
    "# Simulator\n",
    "from queuesim.models import mmc_model, get_simulator_from_model"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Model parameters\n",
    "\n",
    "Since BOBYQA is slow for the considered problem compared to trying all possible solutions, we reduce the number of simulated arrivals by a factor of 10."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Arrivals to be simulated\n",
    "count = 10_000\n",
    "\n",
    "# Arrival process\n",
    "mean_I = 50\n",
    "\n",
    "# Service process\n",
    "mean_S = 600\n",
    "\n",
    "# Costs\n",
    "cost_waiting = 10  # Costs per waiting second\n",
    "cost_c = 120  # Costs per operator"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Definition of the function to be minimized\n",
    "\n",
    "* The workload (=`meanS/meanI`) is 12, which means that at least 12 operators must be used for the model to reach steady state. However, since we generally only simulate significantly less than infinite arrivals, the model will not explode even with a permanent overload.\n",
    "* Since the number of operators must be an integer, the value is rounded internally, i.e. non-integer values may also be passed.\n",
    "* The costs are returned, which include the waiting times of the customers and the operators' labor costs."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "def f(x):\n",
    "    global mean_I, mean_S, count, cost_waiting, cost_c\n",
    "    c = round(x[0])  # c must be an integer\n",
    "    model = mmc_model(mean_I, mean_S, c, count)\n",
    "    get_simulator_from_model(model).run()\n",
    "    waiting_time = model['Dispose'].statistic_client_waiting.mean\n",
    "    return waiting_time * cost_waiting + c * cost_c"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Running the optimization"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Initial solution\n",
    "x0 = np.array([15])\n",
    "\n",
    "# Search range\n",
    "a = np.array([14])\n",
    "b = np.array([24])\n",
    "\n",
    "# Running BOBYQA\n",
    "soln = pybobyqa.solve(f, x0, bounds=(a, b))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Optimization successful?"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "True"
      ]
     },
     "execution_count": 5,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "soln.flag == soln.EXIT_SUCCESS"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "****** Py-BOBYQA Results ******\n",
      "Solution xmin = [16.08110846]\n",
      "Objective value f(xmin) = 2102.671931\n",
      "Needed 27 objective evaluations (at 27 points)\n",
      "Approximate gradient = [-2.29093587e+09]\n",
      "Approximate Hessian = [[2.52304434e+17]]\n",
      "Exit flag = 0\n",
      "Success: rho has reached rhoend\n",
      "******************************\n",
      "\n"
     ]
    }
   ],
   "source": [
    "print(soln)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Result (should be 17)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Cost-optimal number of operators: 16\n"
     ]
    }
   ],
   "source": [
    "print(\"Cost-optimal number of operators:\",round(soln.x[0]))"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.11.0"
  },
  "orig_nbformat": 4
 },
 "nbformat": 4,
 "nbformat_minor": 2
}