{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "# 石油の配合最適化"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##  対象とする問題\n",
    "「**原油1**」、「**原油2**」、「**原油3**」から「**スーパー**」「**レギュラー**」「**ディーゼル**」の3種類のガソリンを作る問題"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "\n",
    "gas_names = [\"super\", \"regular\", \"diesel\"]\n",
    "gas_data = np.array([[3000, 70, 10, 1], [2000, 60, 8, 2], [1000, 50, 6, 1]])\n",
    "\n",
    "oil_names = [\"crude1\", \"crude2\", \"crude3\"]\n",
    "oil_data = np.array([[5000, 45, 12, 0.5], [5000, 35, 6, 2], [5000, 25, 8, 3]])\n",
    "\n",
    "nb_gas  = len(gas_names)\n",
    "nb_oils = len(oil_names)\n",
    "range_gas = range(nb_gas)\n",
    "range_oil = range(nb_oils)\n",
    "print(\"Number of gasoline types = {0}\".format(nb_gas))\n",
    "print(\"Number of crude types = {0}\".format(nb_oils))\n",
    "\n",
    "# global data\n",
    "production_cost = 4\n",
    "production_max = 14000\n",
    "# each $1 spent on advertising increases demand by 10.\n",
    "advert_return = 10"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## データの準備"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from IPython.display import display\n",
    "import pandas as pd\n",
    "gaspd = pd.DataFrame([(gas_names[i],int(gas_data[i][0]),int(gas_data[i][1]),int(gas_data[i][2]),int(gas_data[i][3])) \n",
    "                      for i in range_gas])\n",
    "oilpd = pd.DataFrame([(oil_names[i],int(oil_data[i][0]),int(oil_data[i][1]),int(oil_data[i][2]),oil_data[i][3]) \n",
    "                      for i in range_oil])\n",
    "gaspd.columns = ['name','demand','price','octane','lead']\n",
    "oilpd.columns= ['name','capacity','price','octane','lead']"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print(\"ガソリン情報:\")\n",
    "display(gaspd)\n",
    "\n",
    "print(\"原油情報:\")\n",
    "display(oilpd)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## モデルの宣言"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from docplex.mp.model import Model\n",
    "mdl = Model(name=\"oil_blending\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 決定変数の定義"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# どの原油からどのガソリンを作るかの定義\n",
    "# 3 x 3 の行列になる\n",
    "blends = mdl.continuous_var_matrix(keys1=nb_oils, keys2=nb_gas, lb=0)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# どのガソリンにどれだけ広告を打つかの定義\n",
    "adverts = mdl.continuous_var_list(nb_gas, lb=0)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 制約条件の定義"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### 需要"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# gas_data[g][0] にガソリンごとの需要が入っている\n",
    "# 広告を打つことにより追加需要が発生することを加味する\n",
    "\n",
    "mdl.add_constraints(mdl.sum(blends[o, g] for o in range(nb_oils)) == gas_data[g][0] + advert_return * adverts[g]\n",
    "                   for g in range(nb_gas))\n",
    "mdl.print_information()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### 最大供給量"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# oil_data[o][0]には、原油ごとの供給量が入っている\n",
    "# blends[o, g]を原油単位に足した結果は、供給量以下である必要がある\n",
    "\n",
    "mdl.add_constraints(mdl.sum(blends[o,g] for g in range_gas) <= oil_data[o][0]\n",
    "                    for o in range_oil)\n",
    "mdl.print_information()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### オクタン価と鉛含有率に関する制約"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# ガソリンごとにオクタン価は最小値以上の必要がある\n",
    "# オクタン価は oil_data[o][2]、gas_data[g][2]に入っている\n",
    "\n",
    "mdl.add_constraints(mdl.sum(blends[o,g]*(oil_data[o][2] - gas_data[g][2]) for o in range_oil) >= 0\n",
    "                   for g in range_gas)\n",
    "    \n",
    "# ガソリンごとに鉛含有率は最大値以下の必要がある\n",
    "# 鉛含有率はoil_data[o][3]、gas_data[g][3]に入っている\n",
    "\n",
    "mdl.add_constraints(mdl.sum(blends[o,g]*(oil_data[o][3] - gas_data[g][3]) for o in range_oil) <= 0\n",
    "                   for g in range_gas)\n",
    "\n",
    "mdl.print_information()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### 全体生産量に関する制約"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# 全体生産量は所定の最大値以下である必要がある\n",
    "\n",
    "mdl.add_constraint(mdl.sum(blends) <= production_max)\n",
    "\n",
    "mdl.print_information()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### KPIの定義"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# KPIs\n",
    "\n",
    "# 総広告費用\n",
    "total_advert_cost = mdl.sum(adverts)\n",
    "mdl.add_kpi(total_advert_cost, \"Total advertising cost\")\n",
    "\n",
    "# 総原油コスト\n",
    "total_oil_cost = mdl.sum(blends[o,g] * oil_data[o][1] for o in range_oil for g in range_gas)\n",
    "mdl.add_kpi(total_oil_cost, \"Total Oil cost\")\n",
    "\n",
    "# 総生産コスト\n",
    "total_production_cost = production_cost * mdl.sum(blends)\n",
    "mdl.add_kpi(total_production_cost, \"Total production cost\")\n",
    "\n",
    "# 総収入\n",
    "total_revenue = mdl.sum(blends[o,g] * gas_data[g][1] for g in range(nb_gas) for o in range(nb_oils))\n",
    "mdl.add_kpi(total_revenue, \"Total revenue\")\n",
    "\n",
    "# 最終KPI = 利益の最大化\n",
    "# 利益 = 収入 - 原油コスト - 生産コスト - 広告コスト\n",
    "mdl.maximize(total_revenue - total_oil_cost - total_production_cost - total_advert_cost)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Decision Optimizerで問題を解く"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Decision Optimaizerの呼出し\n",
    "\n",
    "assert mdl.solve(), \"Solve failed\"\n",
    "mdl.report()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##  ソリューション内容の精査"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#  個別KPIの取得\n",
    "\n",
    "all_kpis = [(kp.name, kp.compute()) for kp in mdl.iter_kpis()]\n",
    "kpis_bd = pd.DataFrame(all_kpis, columns=['kpi', 'value'])\n",
    "display(kpis_bd)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "blend_values = [ [ blends[o,g].solution_value for g in range_gas] for o in range_oil]\n",
    "total_gas_prods = [sum(blend_values[o][g] for o in range_oil) for g in range_gas]\n",
    "\n",
    "prods = list(zip(gas_names, total_gas_prods))\n",
    "prods_bd = pd.DataFrame(prods)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# KPIのグラフ化\n",
    "# 総収入、個別コスト、利益\n",
    "\n",
    "%matplotlib inline\n",
    "import matplotlib.pyplot as plt\n",
    "def display_pie(pie_values, pie_labels, colors=None,title=''):\n",
    "    plt.axis(\"equal\")\n",
    "    plt.pie(pie_values, labels=pie_labels, colors=colors, autopct=\"%1.1f%%\")\n",
    "    plt.title(title)\n",
    "    plt.show()\n",
    "                                                           \n",
    "display_pie( [kpnv[1] for kpnv in all_kpis], [kpnv[0] for kpnv in all_kpis],title='KPIs: Revenue - Oil Cost - Production Cost')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### 生産量"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# ガソリン別生産量\n",
    "\n",
    "display_pie(total_gas_prods, gas_names, colors=[\"green\", \"goldenrod\", \"lightGreen\"],title='Gasoline Total Production')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# 原油・製品別生産量内訳\n",
    "\n",
    "sblends = [(gas_names[n], oil_names[o], round(blends[o,n].solution_value)) for n in range_gas for o in range_oil]\n",
    "blends_bd = pd.DataFrame(sblends)\n",
    "display(blends_bd)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# 製品別にグループ化してグラフ表示\n",
    "\n",
    "f, barplot = plt.subplots(1, figsize=(16,5))\n",
    "\n",
    "bar_width = 0.1\n",
    "offset = 0.12\n",
    "rho = 0.7\n",
    "\n",
    "# position of left-bar boundaries\n",
    "bar_l = [o for o in range_oil]\n",
    "\n",
    "mbar_w = 3*bar_width+2*max(0, offset-bar_width)\n",
    "\n",
    "tick_pos = [b*rho + mbar_w/2.0 for b in bar_l]\n",
    "\n",
    "colors = ['olive', 'lightgreen', 'cadetblue']\n",
    "\n",
    "for i in range_oil:\n",
    "    barplot.bar([b*rho + (i*offset) for b in bar_l], \n",
    "                blend_values[i], width=bar_width, color=colors[i], label=oil_names[i])\n",
    "\n",
    "plt.xticks(tick_pos, gas_names)\n",
    "barplot.set_xlabel(\"gasolines\")\n",
    "barplot.set_ylabel(\"blend\")\n",
    "plt.legend(loc=\"upper right\")\n",
    "plt.title('Blend Repartition\\n')\n",
    "                                                                          \n",
    "\n",
    "# Set a buffer around the edge\n",
    "plt.xlim([0, max(tick_pos)+mbar_w +0.5])\n",
    "\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# 原油2からディーゼルを作っていないことの確認\n",
    "\n",
    "print(\"* value of blend[crude2, diesel] is %g\" % blends[1,2])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3.6",
   "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.6.8"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 1
}