Tank_Level_Control.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Tank Level Control\n",
    "\n",
    "Consider a cylindrical tank with no outlet flow and an adjustable inlet flow that is controlled by a valve. The inlet flow rate is not measured but there is a level measurement that shows how much fluid has been added to the tank. The objective of this exercise is to develop a controller that can maintain a certain water level by automatically adjusting the inlet flow rate.\n",
    "\n",
    "See https://apmonitor.com/pdc/index.php/Main/TankLevel"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from IPython import display\n",
    "import numpy as np\n",
    "import matplotlib.pyplot as plt\n",
    "from scipy.integrate import odeint\n",
    "\n",
    "# animate plots?\n",
    "animate=True # True / False\n",
    "\n",
    "# define tank model\n",
    "def tank(Level,time,c,valve):\n",
    "    rho = 1000.0 # water density (kg/m^3)\n",
    "    A = 1.0      # tank area (m^2)\n",
    "    # calculate derivative of the Level\n",
    "    dLevel_dt = (c/(rho*A)) * valve\n",
    "    return dLevel_dt\n",
    "\n",
    "# time span for the simulation for 10 sec, every 0.1 sec\n",
    "ts = np.linspace(0,10,101)\n",
    "\n",
    "# valve operation\n",
    "c = 50.0          # valve coefficient (kg/s / %open)\n",
    "u = np.zeros(101) # u = valve % open\n",
    "\n",
    "# level initial condition\n",
    "Level0 = 0\n",
    "\n",
    "# initial valve position\n",
    "valve = 10\n",
    "\n",
    "# for storing the results\n",
    "z = np.zeros(101)\n",
    "es = np.zeros(101)\n",
    "sps = np.zeros(101)\n",
    "\n",
    "# TO DO: what is the value for ubias?\n",
    "# ubias = ?\n",
    "\n",
    "# TO DO: decide on a tuning value for Kc\n",
    "# Kc = ?\n",
    "\n",
    "# TO DO: record the desired level (set point)\n",
    "SP = 1\n",
    "\n",
    "plt.figure(1,figsize=(12,5))\n",
    "if animate:\n",
    "    plt.ion()\n",
    "    plt.show()\n",
    "\n",
    "# simulate with ODEINT\n",
    "for i in range(100):\n",
    "    # calculate the error\n",
    "    error = SP - Level0\n",
    "\n",
    "    # TO DO: put P-only controller here\n",
    "    # valve = 1\n",
    "\n",
    "    u[i+1] = valve   # store the valve position\n",
    "    es[i+1] = error  # store the error \n",
    "    y = odeint(tank,Level0,[0,0.1],args=(c,valve))\n",
    "    Level0 = y[-1] # take the last point\n",
    "    z[i+1] = Level0 # store the level for plotting\n",
    "    sps[i+1] = SP\n",
    "\n",
    "    if animate:\n",
    "        # update plot\n",
    "        display.clear_output(wait=True)\n",
    "        plt.figure(1,figsize=(12,5))\n",
    "        plt.subplot(3,1,1)\n",
    "        plt.plot(ts[0:i+1],z[0:i+1],'r-',linewidth=3,label='level PV')\n",
    "        plt.plot(ts[0:i+1],sps[0:i+1],'k:',linewidth=3,label='level SP')\n",
    "        plt.ylabel('Tank Level')\n",
    "        plt.legend(loc='best')\n",
    "        plt.subplot(3,1,2)\n",
    "        plt.plot(ts[0:i+1],u[0:i+1],'b--',linewidth=3,label='valve')\n",
    "        plt.ylabel('Valve')    \n",
    "        plt.legend(loc='best')\n",
    "        plt.subplot(3,1,3)\n",
    "        plt.plot(ts[0:i+1],es[0:i+1],'g-',linewidth=3,label='error')\n",
    "        plt.ylabel('Error = SP-PV')\n",
    "        plt.xlabel('Time (sec)')\n",
    "        plt.legend(loc='best')\n",
    "        plt.pause(0.1)\n",
    "\n",
    "if not animate:\n",
    "    # plot results\n",
    "    plt.subplot(3,1,1)\n",
    "    plt.plot(ts,z,'r-',linewidth=3,label='level PV')\n",
    "    plt.plot(ts,sps,'k:',linewidth=3,label='level SP')\n",
    "    plt.ylabel('Tank Level')\n",
    "    plt.legend(loc='best')\n",
    "    plt.subplot(3,1,2)\n",
    "    plt.plot(ts,u,'b--',linewidth=3,label='valve')\n",
    "    plt.ylabel('Valve')    \n",
    "    plt.legend(loc='best')\n",
    "    plt.subplot(3,1,3)\n",
    "    plt.plot(ts,es,'g-',linewidth=3,label='error')\n",
    "    plt.ylabel('Error = SP-PV')    \n",
    "    plt.xlabel('Time (sec)')\n",
    "    plt.legend(loc='best')\n",
    "    plt.show()"
   ]
  }
 ],
 "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.5"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}