{ "cells": [ { "cell_type": "markdown", "id": "b1015ddb", "metadata": {}, "source": [ "\n", "*This notebook contains material from [CBE60499](https://ndcbe.github.io/CBE60499);\n", "content is available [on Github](git@github.com:ndcbe/CBE60499.git).*\n" ] }, { "cell_type": "markdown", "id": "36fefbea", "metadata": {}, "source": [ "\n", "< [2.10 Pyomo Homework 1](https://ndcbe.github.io/CBE60499/02.10-Pyomo1.html) | [Contents](toc.html) | [Tag Index](tag_index.html) | [2.12 Pyomo Homework 3](https://ndcbe.github.io/CBE60499/02.12-Pyomo3.html) >

\"Open

\"Download\"" ] }, { "cell_type": "markdown", "metadata": { "nbpages": { "level": 1, "link": "[2.11 Pyomo Homework 2](https://ndcbe.github.io/CBE60499/02.11-Pyomo2.html#2.11-Pyomo-Homework-2)", "section": "2.11 Pyomo Homework 2" } }, "source": [ "# 2.11 Pyomo Homework 2\n", "\n", "**Due Date:** 2/23/2021" ] }, { "cell_type": "code", "execution_count": 1, "metadata": { "nbpages": { "level": 1, "link": "[2.11 Pyomo Homework 2](https://ndcbe.github.io/CBE60499/02.11-Pyomo2.html#2.11-Pyomo-Homework-2)", "section": "2.11 Pyomo Homework 2" } }, "outputs": [], "source": [ "# This code cell installs packages on Colab\n", "\n", "import sys\n", "if \"google.colab\" in sys.modules:\n", " !wget \"https://raw.githubusercontent.com/ndcbe/CBE60499/main/notebooks/helper.py\"\n", " import helper\n", " helper.install_idaes()\n", " helper.install_ipopt()\n", " helper.install_glpk()" ] }, { "cell_type": "code", "execution_count": 2, "metadata": { "nbpages": { "level": 1, "link": "[2.11 Pyomo Homework 2](https://ndcbe.github.io/CBE60499/02.11-Pyomo2.html#2.11-Pyomo-Homework-2)", "section": "2.11 Pyomo Homework 2" } }, "outputs": [], "source": [ "## IMPORT LIBRARIES\n", "from pyomo.environ import *\n", "import pandas as pd" ] }, { "cell_type": "markdown", "metadata": { "nbpages": { "level": 1, "link": "[2.11 Pyomo Homework 2](https://ndcbe.github.io/CBE60499/02.11-Pyomo2.html#2.11-Pyomo-Homework-2)", "section": "2.11 Pyomo Homework 2" } }, "source": [ "Special thanks to the Pyomo team for create these excercises as part of their excellent PyomoFest workshop." ] }, { "cell_type": "markdown", "metadata": { "nbpages": { "level": 2, "link": "[2.11.1 Some Advanced Pyomo Tricks](https://ndcbe.github.io/CBE60499/02.11-Pyomo2.html#2.11.1-Some-Advanced-Pyomo-Tricks)", "section": "2.11.1 Some Advanced Pyomo Tricks" } }, "source": [ "## 2.11.1 Some Advanced Pyomo Tricks" ] }, { "cell_type": "markdown", "metadata": { "nbpages": { "level": 3, "link": "[2.11.1.1 Using the decorator notation for rules](https://ndcbe.github.io/CBE60499/02.11-Pyomo2.html#2.11.1.1-Using-the-decorator-notation-for-rules)", "section": "2.11.1.1 Using the decorator notation for rules" } }, "source": [ "### 2.11.1.1 Using the decorator notation for rules\n", "\n", "Alternative notation for declaring and defining Pyomo components using decorators exists. Starting with the warehouse location problem code below, change the model to use the decorator notation." ] }, { "cell_type": "code", "execution_count": null, "metadata": { "nbpages": { "level": 3, "link": "[2.11.1.1 Using the decorator notation for rules](https://ndcbe.github.io/CBE60499/02.11-Pyomo2.html#2.11.1.1-Using-the-decorator-notation-for-rules)", "section": "2.11.1.1 Using the decorator notation for rules" } }, "outputs": [], "source": [ "# warehouse_location.py: Warehouse location determination problem\n", "model = ConcreteModel(name=\"(WL)\")\n", "\n", "W = ['Harlingen', 'Memphis', 'Ashland']\n", "C = ['NYC', 'LA', 'Chicago', 'Houston']\n", "d = {('Harlingen', 'NYC'): 1956, \\\n", " ('Harlingen', 'LA'): 1606, \\\n", " ('Harlingen', 'Chicago'): 1410, \\\n", " ('Harlingen', 'Houston'): 330, \\\n", " ('Memphis', 'NYC'): 1096, \\\n", " ('Memphis', 'LA'): 1792, \\\n", " ('Memphis', 'Chicago'): 531, \\\n", " ('Memphis', 'Houston'): 567, \\\n", " ('Ashland', 'NYC'): 485, \\\n", " ('Ashland', 'LA'): 2322, \\\n", " ('Ashland', 'Chicago'): 324, \\\n", " ('Ashland', 'Houston'): 1236 }\n", "P = 2\n", "\n", "model.x = Var(W, C, bounds=(0,1))\n", "model.y = Var(W, within=Binary)\n", "\n", "@model.Objective()\n", "def obj(m):\n", " return sum(d[w,c]*m.x[w,c] for w in W for c in C)\n", "\n", "@model.Constraint(C)\n", "def one_per_cust(m, c):\n", " return sum(m.x[w,c] for w in W) == 1\n", "\n", "# YOUR SOLUTION HERE\n", "def warehouse_active(m, w, c):\n", " return m.x[w,c] <= m.y[w]\n", "\n", "# Note: This is only split across cells because of a bug in nbpages (notebook/website software).\n", "# There is no other reason to split your code across cells." ] }, { "cell_type": "code", "execution_count": null, "metadata": { "nbpages": { "level": 3, "link": "[2.11.1.1 Using the decorator notation for rules](https://ndcbe.github.io/CBE60499/02.11-Pyomo2.html#2.11.1.1-Using-the-decorator-notation-for-rules)", "section": "2.11.1.1 Using the decorator notation for rules" } }, "outputs": [], "source": [ "# YOUR SOLUTION HERE\n", "def num_warehouses(m):\n", " return sum(m.y[w] for w in W) <= P\n", "\n", "SolverFactory('glpk').solve(model)\n", "\n", "model.y.pprint()\n", "model.x.pprint()" ] }, { "cell_type": "markdown", "metadata": { "nbpages": { "level": 3, "link": "[2.11.1.2 Changing parameter values](https://ndcbe.github.io/CBE60499/02.11-Pyomo2.html#2.11.1.2-Changing-parameter-values)", "section": "2.11.1.2 Changing parameter values" } }, "source": [ "### 2.11.1.2 Changing parameter values\n", "\n", "A parameter can be specified to be mutable. This tells Pyomo that the value of the parameter may change in the future, and allows the user to change the parameter value and resolve the problem without the need to rebuild the entire model each time. We will use this functionality to find a better solution to the knapsack problem. We would like to find when the wrench becomes valuable enough to be a part of the optimal solution. Create a Pyomo Parameter for the value of the items, make it mutable, and then write a loop that prints the solution for different wrench values." ] }, { "cell_type": "code", "execution_count": null, "metadata": { "nbpages": { "level": 3, "link": "[2.11.1.2 Changing parameter values](https://ndcbe.github.io/CBE60499/02.11-Pyomo2.html#2.11.1.2-Changing-parameter-values)", "section": "2.11.1.2 Changing parameter values" } }, "outputs": [], "source": [ "A = ['hammer', 'wrench', 'screwdriver', 'towel']\n", "b = {'hammer':8, 'wrench':3, 'screwdriver':6, 'towel':11}\n", "w = {'hammer':5, 'wrench':7, 'screwdriver':4, 'towel':3}\n", "W_max = 14\n", "\n", "model = ConcreteModel()\n", "model.x = Var( A, within=Binary )\n", "# YOUR SOLUTION HERE" ] }, { "cell_type": "markdown", "metadata": { "nbpages": { "level": 3, "link": "[2.11.1.3 Integer cuts](https://ndcbe.github.io/CBE60499/02.11-Pyomo2.html#2.11.1.3-Integer-cuts)", "section": "2.11.1.3 Integer cuts" } }, "source": [ "### 2.11.1.3 Integer cuts\n", "\n", "Often, it can be important to find not only the \"best\" solution, but a number of solutions that are equally optimal, or close to optimal. For discrete optimization problems, this can be done using something known as an integer cut. Consider again the knapsack problem where the choice of which items to select is a discrete variable $x_i \\forall i \\in A$. Let $x_i^*$ be a particular set of $x$ values we want to remove from the feasible solution space. We define an integer cut using two sets. The first set $S_0$ contains the indices for those variables whose current solution is 0, and the second set $S_1$ consists of indices for those variables whose current solution is 1. Given these two sets, an integer cut constraint that would prevent such a solution from appearing again is defined by,\n", "\n", "$\\sum_{i \\in S_0}x[i] + \\sum_{i \\in S_1}(1-x_i) \\geq 1$\n", "\n", "Write a loop that solves the problem 5 times, adding an integer cut to remove the previous solution, and printing the value of the objective function and the solution at each iteration of the loop." ] }, { "cell_type": "code", "execution_count": null, "metadata": { "nbpages": { "level": 3, "link": "[2.11.1.3 Integer cuts](https://ndcbe.github.io/CBE60499/02.11-Pyomo2.html#2.11.1.3-Integer-cuts)", "section": "2.11.1.3 Integer cuts" } }, "outputs": [], "source": [ "from pyomo.environ import *\n", "\n", "A = ['hammer', 'wrench', 'screwdriver', 'towel']\n", "b = {'hammer':8, 'wrench':3, 'screwdriver':6, 'towel':11}\n", "w = {'hammer':5, 'wrench':7, 'screwdriver':4, 'towel':3}\n", "W_max = 14\n", "\n", "model = ConcreteModel()\n", "model.x = Var( A, within=Binary )\n", "\n", "def obj_rule(m):\n", " return sum( b[i]*m.x[i] for i in A )\n", "model.obj = Objective(rule=obj_rule, sense = maximize )\n", "\n", "def weight_con_rule(m):\n", " return sum( w[i]*m.x[i] for i in A ) <= W_max\n", "model.weight_con = Constraint(rule=weight_con_rule)\n", "\n", "opt = SolverFactory('glpk')\n", "\n", "# YOUR SOLUTION HERE" ] }, { "cell_type": "markdown", "metadata": { "nbpages": { "level": 3, "link": "[2.11.1.4 Putting it all together: Lot sizing example (Hart et al., 2017)](https://ndcbe.github.io/CBE60499/02.11-Pyomo2.html#2.11.1.4-Putting-it-all-together:-Lot-sizing-example-(Hart-et-al.,-2017))", "section": "2.11.1.4 Putting it all together: Lot sizing example (Hart et al., 2017)" } }, "source": [ "### 2.11.1.4 Putting it all together: Lot sizing example (Hart et al., 2017)\n", "\n", "We will now write a complete model from scratch using a well-known multi-period optimization problem for optimal lot-sizing adapted from Hagen et al. (2001) shown below.\n", "\n", "$\\min \\sum_{t \\in T}c_ty_t+h_t^+I_t^+ +h_t^-I_t^-$\n", "\n", "s.t. $I_t=I_{t-1}+X_t-d_t, \\forall t \\in T$\n", "\n", "$I_t=I_t^+-I_t^-, \\forall t \\in T$\n", "\n", "$X_t \\leq Py_t, \\forall t \\in T$\n", "\n", "$X_t, I_t^+, I_t^- \\geq 0, \\forall t \\in T$\n", "\n", "$y_t \\in \\{0,1\\}, \\forall t \\in T$\n", "\n", "Our goal is to finnd the optimal production $X_t$ given known demands $d_t$, fixed cost $c_t$ associated with active production in a particular time period, an inventory holding cost $h_t^+$ and a shortage cost $h_t^-$ (cost of keeping a backlog) of orders. The variable $y_t$ (binary) determines if we produce in time $t$ or not, and $I_t^+$ represents inventory that we are storing across time period $t$, while $I_t^-$ represents the magnitude of the backlog. Note that $X_t \\leq Py_t$ is a constraint that only allows production in time period $t$ if the indicator variable $y_t$=1.\n", "\n", "Write a Pyomo model for this problem and solve it using glpk using the data provided below.\n", "\n", "|Parameter|Description|Value|\n", "|---|---|---|\n", "|$c$|fixed cost of production|4.6|\n", "|$I_0^+$|initial value of positive inventory|5.0|\n", "|$I_0^-$|initial value of backlogged orders|0.0|\n", "|$h^+$|cost (per unit) of holding inventory|0.7|\n", "|$h^-$|shortage cost (per unit)|1.2|\n", "|$P$|maximum production amoung (big-M value)|5|\n", "|$d$|demand|[5,7,6.2,3.1,1,7]|\n", "\n", "**Reference**: Hart, W. E., Laird, C. D., Watson, J. P., Woodruff, D. L., Hackebeil, G. A., Nicholson, B. L., and Siirola, J. D. Pyomo: Optimization Modeling in Python (Second Edition), Vol (67), Springer Verlag, 2017." ] }, { "cell_type": "code", "execution_count": null, "metadata": { "nbpages": { "level": 3, "link": "[2.11.1.4 Putting it all together: Lot sizing example (Hart et al., 2017)](https://ndcbe.github.io/CBE60499/02.11-Pyomo2.html#2.11.1.4-Putting-it-all-together:-Lot-sizing-example-(Hart-et-al.,-2017))", "section": "2.11.1.4 Putting it all together: Lot sizing example (Hart et al., 2017)" } }, "outputs": [], "source": [ "model = ConcreteModel()\n", "model.T = RangeSet(5) # time periods\n", "\n", "i0 = 5.0 # initial inventory\n", "c = 4.6 # setup cost\n", "h_pos = 0.7 # inventory holding cost\n", "h_neg = 1.2 # shortage cost\n", "P = 5.0 # maximum production amount\n", "\n", "# demand during period t\n", "d = {1: 5.0, 2:7.0, 3:6.2, 4:3.1, 5:1.7}\n", "\n", "# YOUR SOLUTION HERE\n", "\n", "# solve the problem\n", "solver = SolverFactory('glpk')\n", "solver.solve(model)\n", "\n", "# print the results\n", "for t in model.T:\n", " print('Period: {0}, Prod. Amount: {1}'.format(t, value(model.x[t]))) " ] }, { "cell_type": "markdown", "metadata": { "nbpages": { "level": 2, "link": "[2.11.2 Nonlinear programs: initial and problem formulation are very important!](https://ndcbe.github.io/CBE60499/02.11-Pyomo2.html#2.11.2-Nonlinear-programs:-initial-and-problem-formulation-are-very-important!)", "section": "2.11.2 Nonlinear programs: initial and problem formulation are very important!" } }, "source": [ "## 2.11.2 Nonlinear programs: initial and problem formulation are very important!" ] }, { "cell_type": "markdown", "metadata": { "nbpages": { "level": 3, "link": "[2.11.2.1 Alternative initialization](https://ndcbe.github.io/CBE60499/02.11-Pyomo2.html#2.11.2.1-Alternative-initialization)", "section": "2.11.2.1 Alternative initialization" } }, "source": [ "### 2.11.2.1 Alternative initialization\n", "\n", "Effective initialization can be critical for solving nonlinear problems, since they can have several local solutions and numerical diffculties. Solve the [Rosenbrock problem](https://en.wikipedia.org/wiki/Rosenbrock_function) using different initial values for the x variables. Write a loop that varies the initial value from 2.0 to 6.0, solves the problem, and prints the solution for each iteration of the loop." ] }, { "cell_type": "code", "execution_count": null, "metadata": { "nbpages": { "level": 3, "link": "[2.11.2.1 Alternative initialization](https://ndcbe.github.io/CBE60499/02.11-Pyomo2.html#2.11.2.1-Alternative-initialization)", "section": "2.11.2.1 Alternative initialization" } }, "outputs": [], "source": [ "model = ConcreteModel()\n", "model.x = Var()\n", "model.y = Var()\n", "\n", "def rosenbrock(m):\n", " return (1.0-m.x)**2 + 10000.0*(m.y - m.x**2)**2\n", "model.obj = Objective(rule=rosenbrock, sense=minimize)\n", "\n", "\n", "solver = SolverFactory('ipopt')\n", "\n", "print('x_init, y_init, x_soln, y_soln')\n", "\n", "# YOUR SOLUTION HERE \n" ] }, { "cell_type": "markdown", "metadata": { "nbpages": { "level": 3, "link": "[2.11.2.1 Alternative initialization](https://ndcbe.github.io/CBE60499/02.11-Pyomo2.html#2.11.2.1-Alternative-initialization)", "section": "2.11.2.1 Alternative initialization" } }, "source": [ "As elaborated here, the [Rosenbrock problem](https://en.wikipedia.org/wiki/Rosenbrock_function) is a classic \"hard\" test case for optimization algorithms. Your results may surprise you (and show the effectiveness of Pyomo and Ipopt!)." ] }, { "cell_type": "markdown", "metadata": { "nbpages": { "level": 3, "link": "[2.11.2.2 Evaluation errors](https://ndcbe.github.io/CBE60499/02.11-Pyomo2.html#2.11.2.2-Evaluation-errors)", "section": "2.11.2.2 Evaluation errors" } }, "source": [ "### 2.11.2.2 Evaluation errors\n", "\n", "Consider the following problem with **initial values** $x$=5, $y$=5.\n", "\n", "$\\min_{x,y} f(x,y)=(x-1.01)^2+y^2$\n", "\n", "s.t. $y=\\sqrt{x-1.0}$\n" ] }, { "cell_type": "markdown", "metadata": { "nbpages": { "level": 3, "link": "[2.11.2.2 Evaluation errors](https://ndcbe.github.io/CBE60499/02.11-Pyomo2.html#2.11.2.2-Evaluation-errors)", "section": "2.11.2.2 Evaluation errors" } }, "source": [ "1. Formulate this Pyomo model and solve using IPOPT. You should get a list of errors from the solver. Add the IPOPT solver option solver.options['halt_on_ampl_error']='yes' to find the problem. *Hint*: the error output might be ordered strangely, look this up in the console output. What did you discover? How might you fix this?" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "nbpages": { "level": 3, "link": "[2.11.2.2 Evaluation errors](https://ndcbe.github.io/CBE60499/02.11-Pyomo2.html#2.11.2.2-Evaluation-errors)", "section": "2.11.2.2 Evaluation errors" } }, "outputs": [], "source": [ "# YOUR SOLUTION HERE" ] }, { "cell_type": "markdown", "metadata": { "nbpages": { "level": 3, "link": "[2.11.2.2 Evaluation errors](https://ndcbe.github.io/CBE60499/02.11-Pyomo2.html#2.11.2.2-Evaluation-errors)", "section": "2.11.2.2 Evaluation errors" } }, "source": [ "**Question Answers**\n", "\n", "*Fill in here*" ] }, { "cell_type": "markdown", "metadata": { "nbpages": { "level": 3, "link": "[2.11.2.2 Evaluation errors](https://ndcbe.github.io/CBE60499/02.11-Pyomo2.html#2.11.2.2-Evaluation-errors)", "section": "2.11.2.2 Evaluation errors" } }, "source": [ "2. Add bounds $x \\geq 1$ to fix this problem and resolve. Comment on the number of iterations and the quality of solution. Note, the problem still occurs because $x \\geq 1$ is not enforced exactly, and small numerical values still cause the error." ] }, { "cell_type": "code", "execution_count": null, "metadata": { "nbpages": { "level": 3, "link": "[2.11.2.2 Evaluation errors](https://ndcbe.github.io/CBE60499/02.11-Pyomo2.html#2.11.2.2-Evaluation-errors)", "section": "2.11.2.2 Evaluation errors" } }, "outputs": [], "source": [ "# YOUR SOLUTION HERE" ] }, { "cell_type": "markdown", "metadata": { "nbpages": { "level": 3, "link": "[2.11.2.2 Evaluation errors](https://ndcbe.github.io/CBE60499/02.11-Pyomo2.html#2.11.2.2-Evaluation-errors)", "section": "2.11.2.2 Evaluation errors" } }, "source": [ "**Discussion**\n", "\n", "*Fill in here*" ] }, { "cell_type": "markdown", "metadata": { "nbpages": { "level": 3, "link": "[2.11.2.2 Evaluation errors](https://ndcbe.github.io/CBE60499/02.11-Pyomo2.html#2.11.2.2-Evaluation-errors)", "section": "2.11.2.2 Evaluation errors" } }, "source": [ "3. Think about other solutions for this problem and attempt to implement one of these solutions. *Hint*: $x \\geq 1.001$." ] }, { "cell_type": "code", "execution_count": null, "metadata": { "nbpages": { "level": 3, "link": "[2.11.2.2 Evaluation errors](https://ndcbe.github.io/CBE60499/02.11-Pyomo2.html#2.11.2.2-Evaluation-errors)", "section": "2.11.2.2 Evaluation errors" } }, "outputs": [], "source": [ "model = ConcreteModel()\n", "\n", "model.x = Var(initialize=5.0, bounds=(1.001,None))\n", "model.y = Var(initialize=5.0)\n", "\n", "def obj_rule(m):\n", " return (m.x-1.01)**2 + m.y**2\n", "model.obj = Objective(rule=obj_rule)\n", "\n", "def con_rule(m):\n", " return m.y == sqrt(m.x - 1.0)\n", "model.con = Constraint(rule=con_rule)\n", "\n", "solver = SolverFactory('ipopt')\n", "solver.options['halt_on_ampl_error'] = 'yes'\n", "solver.solve(model, tee=True)\n", "\n", "print( value(model.x) )\n", "print( value(model.y) )" ] }, { "cell_type": "markdown", "metadata": { "nbpages": { "level": 3, "link": "[2.11.2.3 Alternative formulations](https://ndcbe.github.io/CBE60499/02.11-Pyomo2.html#2.11.2.3-Alternative-formulations)", "section": "2.11.2.3 Alternative formulations" } }, "source": [ "### 2.11.2.3 Alternative formulations\n", "\n", "Consider the following problem with initial values $x$=5, $y$=5.\n", "\n", "$\\min_{x,y} f(x,y)=(x-1.01)^2+y^2$\n", "\n", "s.t. $\\frac{x-1}{y}=1$\n", "\n", "Note, the solution to this problem is $x$=1.005 and $y$=0.005. There are several ways that the problem above can be reformulated. Some examples are shown below. Which ones do you expect to be better? Why?\n", "\n", "1. $\\min_{x,y} f(x,y)=(x-1.01)^2+y^2$\n", "\n", "s.t. $\\frac{x-1}{y}=1$\n", "\n", "2. $\\min_{x,y} f(x,y)=(x-1.01)^2+y^2$\n", "\n", "s.t. $\\frac{x}{y+1}=1$\n", "\n", "3. $\\min_{x,y} f(x,y)=(x-1.01)^2+y^2$\n", "\n", "s.t. $y=x-1$\n", "\n", "Implement Pyomo models for each formulation and solve with IPOPT." ] }, { "cell_type": "markdown", "metadata": { "nbpages": { "level": 4, "link": "[2.11.2.3.1 Formulation 1](https://ndcbe.github.io/CBE60499/02.11-Pyomo2.html#2.11.2.3.1-Formulation-1)", "section": "2.11.2.3.1 Formulation 1" } }, "source": [ "#### 2.11.2.3.1 Formulation 1" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "nbpages": { "level": 4, "link": "[2.11.2.3.1 Formulation 1](https://ndcbe.github.io/CBE60499/02.11-Pyomo2.html#2.11.2.3.1-Formulation-1)", "section": "2.11.2.3.1 Formulation 1" } }, "outputs": [], "source": [ "# YOUR SOLUTION HERE" ] }, { "cell_type": "markdown", "metadata": { "nbpages": { "level": 4, "link": "[2.11.2.3.2 Formulation 2](https://ndcbe.github.io/CBE60499/02.11-Pyomo2.html#2.11.2.3.2-Formulation-2)", "section": "2.11.2.3.2 Formulation 2" } }, "source": [ "#### 2.11.2.3.2 Formulation 2" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "nbpages": { "level": 4, "link": "[2.11.2.3.2 Formulation 2](https://ndcbe.github.io/CBE60499/02.11-Pyomo2.html#2.11.2.3.2-Formulation-2)", "section": "2.11.2.3.2 Formulation 2" } }, "outputs": [], "source": [ "# YOUR SOLUTION HERE" ] }, { "cell_type": "markdown", "metadata": { "nbpages": { "level": 4, "link": "[2.11.2.3.3 Formulation 3](https://ndcbe.github.io/CBE60499/02.11-Pyomo2.html#2.11.2.3.3-Formulation-3)", "section": "2.11.2.3.3 Formulation 3" } }, "source": [ "#### 2.11.2.3.3 Formulation 3" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "nbpages": { "level": 4, "link": "[2.11.2.3.3 Formulation 3](https://ndcbe.github.io/CBE60499/02.11-Pyomo2.html#2.11.2.3.3-Formulation-3)", "section": "2.11.2.3.3 Formulation 3" } }, "outputs": [], "source": [ "# YOUR SOLUTION HERE" ] }, { "cell_type": "markdown", "metadata": { "nbpages": { "level": 4, "link": "[2.11.2.3.3 Formulation 3](https://ndcbe.github.io/CBE60499/02.11-Pyomo2.html#2.11.2.3.3-Formulation-3)", "section": "2.11.2.3.3 Formulation 3" } }, "source": [ "Note the nuber of iterations and quality of solutions. What can you learn about the problem formulation from these examples?\n", "\n", "**Discussion**\n", "\n", "*Fill in here*" ] }, { "cell_type": "markdown", "metadata": { "nbpages": { "level": 4, "link": "[2.11.2.3.3 Formulation 3](https://ndcbe.github.io/CBE60499/02.11-Pyomo2.html#2.11.2.3.3-Formulation-3)", "section": "2.11.2.3.3 Formulation 3" } }, "source": [ "Bounds and initialization can be very helpful when solving nonlinear optimization problems. Resolve the original problem below, but add bounds, $y \\geq 0$. Note the number of iterations and quality of solution, and compare with what you found for Formulation 1." ] }, { "cell_type": "code", "execution_count": null, "metadata": { "nbpages": { "level": 4, "link": "[2.11.2.3.3 Formulation 3](https://ndcbe.github.io/CBE60499/02.11-Pyomo2.html#2.11.2.3.3-Formulation-3)", "section": "2.11.2.3.3 Formulation 3" } }, "outputs": [], "source": [ "# YOUR SOLUTION HERE" ] }, { "cell_type": "markdown", "metadata": { "nbpages": { "level": 4, "link": "[2.11.2.3.3 Formulation 3](https://ndcbe.github.io/CBE60499/02.11-Pyomo2.html#2.11.2.3.3-Formulation-3)", "section": "2.11.2.3.3 Formulation 3" } }, "source": [ "**Discussion**\n", "\n", "*Fill in here*" ] }, { "cell_type": "markdown", "metadata": { "nbpages": { "level": 3, "link": "[2.11.2.4 Reactor design problem (Hart et al., 2017; Bequette, 2003)](https://ndcbe.github.io/CBE60499/02.11-Pyomo2.html#2.11.2.4-Reactor-design-problem-(Hart-et-al.,-2017;-Bequette,-2003))", "section": "2.11.2.4 Reactor design problem (Hart et al., 2017; Bequette, 2003)" } }, "source": [ "### 2.11.2.4 Reactor design problem (Hart et al., 2017; Bequette, 2003)\n", "\n", "Here we will consider a chemical reactor designed to produce product B from reactant A using a reaction scheme known as the Van de Vusse reaction:\n", "\n", "$A \\overset{k_1}{\\rightarrow} B \\overset{k_2}{\\rightarrow} C$\n", "\n", "$2A \\overset{k_3}{\\rightarrow} D$\n", "\n", "Under appropriate assumptions, $F$ is the volumetric flowrate through the tank. The concentation of component A in the feed is $c_{Af}$, and the concentrations in the reactor are equivalent to the concentrations of each component \n", "flowing out of the reactor, given by $c_A$, $c_B$, $c_C$, and $c_D$.\n", "\n", "If the reactor is too small, we will not produce sufficient quantity of B, and if the reactor is too large, much of B will be further reacted to form the undesired product C. Therefore, our goal is to solve for the reactor volume that maximizes the outlet concentration for product B.\n", "\n", "The steady-state mole balances for each of the four components are given by\n", "\n", "$0=\\frac{F}{V}c_{Af}-\\frac{F}{V}c_A-k_1C_A-2k_3c^2_A$\n", "\n", "$0=-\\frac{F}{V}c_{B}+k_1C_A-k_2c_B$\n", "\n", "$0=-\\frac{F}{V}c_{C}+k_2c_B$\n", "\n", "$0=-\\frac{F}{V}c_{D}+k_3c^2_A$\n", "\n", "The known parameters for the system are:\n", "\n", "$c_{Af}=10\\frac{\\mathrm{gmol}}{\\mathrm{m}^3}$\n", "\n", "$k_1=\\frac{5}{6}\\mathrm{min}^{-1}$\n", "\n", "$k_2=\\frac{5}{3}\\mathrm{min}^{-1}$\n", "\n", "$k_3=\\frac{1}{6000}\\frac{\\mathrm{m}^3}{\\mathrm{gmol}~\\mathrm{min}^{-1}}$\n", "\n", "Formulate and solve this optimization problem using Pyomo. Since the volumetric flowrate $F$ always appears as the numerator over the reactor volume $V$, it is common to consider this ratio as a single variable, called the space-velocity $SV$.\n", "\n", "**References**:\n", "Hart, W. E., Laird, C. D., Watson, J. P., Woodruff, D. L., Hackebeil, G. A., Nicholson, B. L., and Siirola, J. D. Pyomo: Optimization Modeling in Python (Second Edition), Vol (67), Springer Verlag, 2017.\n", "\n", "B.W. Bequette. Process control: modeling, design, and simulation. Prentice Hall, 2003." ] }, { "cell_type": "code", "execution_count": 3, "metadata": { "nbpages": { "level": 3, "link": "[2.11.2.4 Reactor design problem (Hart et al., 2017; Bequette, 2003)](https://ndcbe.github.io/CBE60499/02.11-Pyomo2.html#2.11.2.4-Reactor-design-problem-(Hart-et-al.,-2017;-Bequette,-2003))", "section": "2.11.2.4 Reactor design problem (Hart et al., 2017; Bequette, 2003)" } }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "5 Var Declarations\n", " ca : Size=1, Index=None\n", " Key : Lower : Value : Upper : Fixed : Stale : Domain\n", " None : 0 : 3874.2588672317133 : None : False : False : PositiveReals\n", " cb : Size=1, Index=None\n", " Key : Lower : Value : Upper : Fixed : Stale : Domain\n", " None : 0 : 1072.437200108632 : None : False : False : PositiveReals\n", " cc : Size=1, Index=None\n", " Key : Lower : Value : Upper : Fixed : Stale : Domain\n", " None : 0 : 1330.093533408881 : None : False : False : PositiveReals\n", " cd : Size=1, Index=None\n", " Key : Lower : Value : Upper : Fixed : Stale : Domain\n", " None : 0 : 1861.6051996253875 : None : False : False : PositiveReals\n", " sv : Size=1, Index=None\n", " Key : Lower : Value : Upper : Fixed : Stale : Domain\n", " None : 0 : 1.343811761067278 : None : False : False : PositiveReals\n", "\n", "1 Objective Declarations\n", " obj : Size=1, Index=None, Active=True\n", " Key : Active : Sense : Expression\n", " None : True : maximize : cb\n", "\n", "4 Constraint Declarations\n", " ca_bal : Size=1, Index=None, Active=True\n", " Key : Lower : Body : Upper : Active\n", " None : 0.0 : 10000.0*sv - sv*ca - 0.8333333333333334*ca - 0.0003333333333333333*ca**2.0 : 0.0 : True\n", " cb_bal : Size=1, Index=None, Active=True\n", " Key : Lower : Body : Upper : Active\n", " None : 0.0 : - sv*cb + 0.8333333333333334*ca - 1.6666666666666667*cb : 0.0 : True\n", " cc_bal : Size=1, Index=None, Active=True\n", " Key : Lower : Body : Upper : Active\n", " None : 0.0 : - sv*cc + 1.6666666666666667*cb : 0.0 : True\n", " cd_bal : Size=1, Index=None, Active=True\n", " Key : Lower : Body : Upper : Active\n", " None : 0.0 : - sv*cd + 0.00016666666666666666*ca**2.0 : 0.0 : True\n", "\n", "10 Declarations: sv ca cb cc cd obj ca_bal cb_bal cc_bal cd_bal\n" ] } ], "source": [ "# YOUR SOLUTION HERE" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "nbpages": { "level": 3, "link": "[2.11.2.4 Reactor design problem (Hart et al., 2017; Bequette, 2003)](https://ndcbe.github.io/CBE60499/02.11-Pyomo2.html#2.11.2.4-Reactor-design-problem-(Hart-et-al.,-2017;-Bequette,-2003))", "section": "2.11.2.4 Reactor design problem (Hart et al., 2017; Bequette, 2003)" } }, "outputs": [], "source": [] }, { "cell_type": "markdown", "id": "2802fa53", "metadata": {}, "source": [ "\n", "< [2.10 Pyomo Homework 1](https://ndcbe.github.io/CBE60499/02.10-Pyomo1.html) | [Contents](toc.html) | [Tag Index](tag_index.html) | [2.12 Pyomo Homework 3](https://ndcbe.github.io/CBE60499/02.12-Pyomo3.html) >

\"Open

\"Download\"" ] } ], "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.8.6" } }, "nbformat": 4, "nbformat_minor": 2 }