fopdt_tclab_widget.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import pandas as pd\n",
    "import matplotlib.pyplot as plt\n",
    "from scipy.integrate import odeint\n",
    "import ipywidgets as wg\n",
    "from IPython.display import display\n",
    "from IPython.display import clear_output "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# try to read local data file first\n",
    "try:\n",
    "    filename = 'data.csv'\n",
    "    data = pd.read_csv(filename)\n",
    "except:\n",
    "    filename = 'http://apmonitor.com/pdc/uploads/Main/tclab_data2.txt'\n",
    "    data = pd.read_csv(filename)\n",
    "\n",
    "n = 601 # time points to plot\n",
    "tf = 600.0 # final time\n",
    "\n",
    "# Use expected room temperature for initial condition\n",
    "#y0 = [23.0,23.0]\n",
    "# Use initial condition\n",
    "y0d = [data['T1'].values[0],data['T2'].values[0]]\n",
    "\n",
    "# load data\n",
    "Q1 = data['Q1'].values\n",
    "Q2 = data['Q2'].values\n",
    "T1 = data['T1'].values\n",
    "T2 = data['T2'].values\n",
    "T1p = np.ones(n)*y0d[0]\n",
    "T2p = np.ones(n)*y0d[1]\n",
    "\n",
    "def process(y,t,u1,u2,Kp,Kd,taup):\n",
    "    y1,y2 = y\n",
    "    dy1dt = (1.0/taup) * (-(y1-y0d[0]) + Kp * u1 + Kd * (y2-y1))\n",
    "    dy2dt = (1.0/taup) * (-(y2-y0d[1]) + (Kp/2.0) * u2 + Kd * (y1-y2))\n",
    "    return [dy1dt,dy2dt]\n",
    "\n",
    "def fopdtPlot(Kp,Kd,taup,thetap):\n",
    "    y0 = y0d\n",
    "    t = np.linspace(0,tf,n) # create time vector\n",
    "    iae = 0.0\n",
    "    # loop through all time steps\n",
    "    for i in range(1,n):\n",
    "        # simulate process for one time step\n",
    "        ts = [t[i-1],t[i]]         # time interval\n",
    "        inputs = (Q1[max(0,i-int(thetap))],Q2[max(0,i-int(thetap))],Kp,Kd,taup)\n",
    "        y = odeint(process,y0,ts,args=inputs)\n",
    "        y0 = y[1]                  # record new initial condition\n",
    "        T1p[i] = y0[0]\n",
    "        T2p[i] = y0[1]        \n",
    "        iae += np.abs(T1[i]-T1p[i]) + np.abs(T2[i]-T2p[i])\n",
    "\n",
    "    # plot FOPDT response\n",
    "    plt.figure(1,figsize=(15,7))\n",
    "    clear_output()\n",
    "    plt.subplot(2,1,1)\n",
    "    plt.plot(t,T1,'r.',linewidth=2,label='Temperature 1 (meas)')\n",
    "    plt.plot(t,T2,'b.',linewidth=2,label='Temperature 2 (meas)')\n",
    "    plt.plot(t,T1p,'r--',linewidth=2,label='Temperature 1 (pred)')\n",
    "    plt.plot(t,T2p,'b--',linewidth=2,label='Temperature 2 (pred)')\n",
    "    plt.ylabel(r'T $(^oC)$')\n",
    "    plt.text(200,20,'Integral Abs Error: ' + str(np.round(iae,2)))\n",
    "    plt.text(400,35,r'$K_p$: ' + str(np.round(Kp,2)))  \n",
    "    plt.text(400,30,r'$K_d$: ' + str(np.round(Kd,2)))  \n",
    "    plt.text(400,25,r'$\\tau_p$: ' + str(np.round(taup,0)) + ' sec')  \n",
    "    plt.text(400,20,r'$\\theta_p$: ' + str(np.round(thetap,0)) + ' sec')  \n",
    "    plt.legend(loc=2)\n",
    "    plt.subplot(2,1,2)\n",
    "    plt.plot(t,Q1,'b--',linewidth=2,label=r'Heater 1 ($Q_1$)')\n",
    "    plt.plot(t,Q2,'r:',linewidth=2,label=r'Heater 2 ($Q_2$)')\n",
    "    plt.legend(loc='best')\n",
    "    plt.xlabel('time (sec)')\n",
    "    plt.show()\n",
    "\n",
    "Kp_slide = wg.FloatSlider(value=0.5,min=0.1,max=1.5,step=0.05)\n",
    "Kd_slide = wg.FloatSlider(value=0.0,min=0.0,max=1.0,step=0.05)\n",
    "taup_slide = wg.FloatSlider(value=50.0,min=50.0,max=250.0,step=5.0)\n",
    "thetap_slide = wg.FloatSlider(value=0,min=0,max=30,step=1)\n",
    "wg.interact(fopdtPlot, Kp=Kp_slide, Kd=Kd_slide, taup=taup_slide,thetap=thetap_slide)\n",
    "print('FOPDT Simulator: Adjust Kp, Kd, taup, and thetap ' + \\\n",
    "      'to achieve lowest Integral Abs Error')"
   ]
  }
 ],
 "metadata": {
  "anaconda-cloud": {},
  "kernelspec": {
   "display_name": "Python 3 (ipykernel)",
   "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.12.4"
  },
  "widgets": {
   "state": {
    "b2c6e0c4193d450c8498954bcdc50530": {
     "views": [
      {
       "cell_index": 1
      }
     ]
    }
   },
   "version": "1.2.0"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}