update CDS 110, L9

murrayrm · murrayrm · commit 354f3cc8efab · 2024-08-08T19:10:15.000-07:00
diff --git a/examples/cds110-L9_servomech-pid.ipynb b/examples/cds110-L9_servomech-pid.ipynb
@@ -42,7 +42,13 @@
    "source": [
     "## System dynamics\n",
     "\n",
-    "Consider a simple mechanism for positioning a mechanical arm whose equations of motion are given by\n",
+    "Consider a simple mechanism consisting of a spring loaded arm that is driven by a  motor, as shown below:\n",
+    "\n",
+    "<center><img src=\"https://www.cds.caltech.edu/~murray/courses/cds110/sp2024/servomech-diagram.png\" width=200 alt=\"servomech-diagram\"></center>\n",
+    "\n",
+    "The motor applies a torque that twists the arm against a linear spring and moves the end of the arm across a rotating platter. The input to the system is the motor torque $\\tau_\\text{m}$. The force exerted by the spring is a nonlinear function of the head position due to the way it is attached.\n",
+    "\n",
+    "The equations of motion for the system are given by\n",
     "\n",
     "$$\n",
     "J \\ddot \\theta = -b \\dot\\theta - k r\\sin\\theta + \\tau_\\text{m},\n",
@@ -56,20 +62,14 @@
     "  + \\begin{bmatrix} 0 \\\\ 1/J \\end{bmatrix} \\tau_\\text{m}.\n",
     "$$\n",
     "\n",
-    "The system consists of a spring loaded arm that is driven by a  motor, as shown below.\n",
-    "\n",
-    "<center><img src=\"https://www.cds.caltech.edu/~murray/courses/cds110/sp2024/servomech-diagram.png\" alt=\"servomech-diagram\" height=\"240\"></center>\n",
-    "\n",
-    "The motor applies a torque that twists the arm against a linear spring and moves the end of the arm across a rotating platter. The input to the system is the motor torque $\\tau_\\text{m}$. The force exerted by the spring is a nonlinear function of the head position due to the way it is attached.\n",
-    "\n",
     "The system parameters are given by\n",
     "\n",
     "$$\n",
     "k = 1,\\quad J = 100,\\quad b = 10,\n",
     "\\quad r = 1,\\quad l = 2,\\quad \\epsilon = 0.01.\n",
     "$$\n",
     "\n",
-    "and we assume that time is measured in msec and distance in cm.  (The constants here are made up and don't necessarily reflect a real disk drive, though the units and time constants are motivated by computer disk drives.)"
+    "and we assume that time is measured in milliseconds (ms) and distance in centimeters (cm).  (The constants here are made up and don't necessarily reflect a real disk drive, though the units and time constants are motivated by computer disk drives.)"
    ]
   },
   {
@@ -170,7 +170,7 @@
    "source": [
     "### Linearization\n",
     "\n",
-    "To study the open loop dynamicsof the system, we compute the linearization of the dynamics about the equilibrium point corresponding to $\\theta_\\text{e} = 15^\\circ$."
+    "To study the open loop dynamics of the system, we compute the linearization of the dynamics about the equilibrium point corresponding to $\\theta_\\text{e} = 15^\\circ$."
    ]
   },
   {
@@ -195,19 +195,6 @@
     "ct.bode_plot(P)"
    ]
   },
-  {
-   "cell_type": "code",
-   "execution_count": null,
-   "metadata": {},
-   "outputs": [],
-   "source": [
-    "# Colab: create a simplified form to get rid of warning messages\n",
-    "if len(P.num[0][0]) > 1:\n",
-    "  P = ct.tf(P.num[0][0][2], P.den, name='P')  # Fix up conditioning\n",
-    "print(P)\n",
-    "ct.bode_plot(P)"
-   ]
-  },
   {
    "cell_type": "markdown",
    "metadata": {
@@ -216,6 +203,10 @@
    "source": [
     "## Ziegler-Nichols tuning\n",
     "\n",
+    "Ziegler-Nichols tuning provides a method for choosing the gains of a PID controller that give reasonable closed loop response.  More information can be found in [Feedback Systems](https://fbswiki.org/wiki/index.php/Feedback_Systems:_An_Introduction_for_Scientists_and_Engineers) (FBS2e), Section 11.3.\n",
+    "\n",
+    "We show here the figures and tables that we will use (from FBS2e):\n",
+    "\n",
     "<center>\n",
     "<table>\n",
     "<tr>\n",
@@ -226,7 +217,9 @@
     "<img src=\"https://www.cds.caltech.edu/~murray/courses/cds110/sp2024/zn-step-table.png\" width=300>\n",
     "</td>\n",
     "</tr>\n",
-    "</center>"
+    "</center>\n",
+    "\n",
+    "To use the Ziegler-Nichols turning rules, we plot the step response, compute the parameters (shown in the figure), and then apply the formulas in the table:"
    ]
   },
   {
@@ -257,6 +250,13 @@
     "print(f\"{a=}, {tau=}\")"
    ]
   },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "We can then construct a controller using the parameters:"
+   ]
+  },
   {
    "cell_type": "code",
    "execution_count": null,
@@ -325,7 +325,9 @@
     "id": "6iZwB2WEeg8S"
    },
    "source": [
-    "## Loop shaping"
+    "## Loop shaping\n",
+    "\n",
+    "A better design can be obtained by looking at the loop transfer function and adjusting the controller parameters to give a loop shape that will give closed loop properties.  We show the steps for such a design here:"
    ]
   },
   {
@@ -390,13 +392,24 @@
     "ct.nyquist([ltf_shape])"
    ]
   },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "We see that the loop shaping controller has better step response (faster rise and settling time, less overshoot)."
+   ]
+  },
   {
    "cell_type": "markdown",
    "metadata": {
     "id": "GyXQXykafzWs"
    },
    "source": [
-    "### Gang of Four"
+    "### Gang of Four\n",
+    "\n",
+    "When designing a controller, it is important to look at all of the input/output responses, not just the response from reference to output (which is what the step response above focuses on). \n",
+    "\n",
+    "In the frequency domain, the Gang of 4 plots provide useful information on all (important) input/output pairs:"
    ]
   },
   {
@@ -408,15 +421,24 @@
     "ct.gangof4(P, ctrl_shape)"
    ]
   },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "These all look pretty resonable, except that the transfer function from the reference $r$ to the system input $u$ is getting large at high frequency.  This occurs because we did not filter the derivative on the PID controller, so high frequency components of the reference signal (or the measurement noise!) get amplified.  We will fix this in the more advanced controller below."
+   ]
+  },
   {
    "cell_type": "markdown",
    "metadata": {
     "id": "uFO3wiWXhBAK"
    },
    "source": [
-    "## Anti-windup\n",
+    "## Anti-windup + derivative filtering\n",
+    "\n",
+    "In addition to the amplification of high frequency signals due to the derivative term, another practical consideration in the use of PID controllers is integrator windup.  Integrator windup occurs when there are limits on the control inputs so that the error signal may not descrease quickly.  This causes the integral term in the PID controller to see an error for a long period of time, and the resulting integration of the error must be offset by making the error have opposite sign for some period of time.  This is often undesireable.\n",
     "\n",
-    "We now implement the full PID controller with anti-windup and derivative filtering:\n",
+    "To see how to address both amplification of noise due to the derivative term and integrator windup effects in the presence of input constraints, we now implement PID controller with anti-windup and derivative filtering, as shown in the following figure (see also Figure 11.11 in [FBS2e](https://fbswiki.org/wiki/index.php/Feedback_Systems:_An_Introduction_for_Scientists_and_Engineers)):\n",
     "\n",
     "<center>\n",
     "<img src=\"https://www.cds.caltech.edu/~murray/courses/cds110/sp2024/pid-aw-diagram.png\"</img>\n",
@@ -529,17 +551,9 @@
     "    clsys, timepts, stepsize, params={'umax': 5, 'kaw': 2*ki, 'a': 100})\n",
     "\n",
     "# Plot the time responses in a single plot\n",
-    "out = ct.time_response_plot(resp_ln, color='b', plot_inputs=False)\n",
-    "ct.time_response_plot(resp_cl, color='r', plot_inputs=False)\n",
-    "ct.time_response_plot(resp_aw, color='g', plot_inputs=False)\n",
-    "\n",
-    "# Annotations\n",
-    "axs = ct.get_plot_axes(out)\n",
-    "axs[0, 0].legend(['linear', 'clipped', 'anti-windup'])\n",
-    "axs[0, 0].plot([0, timepts[-1]], [stepsize, stepsize], 'k--')\n",
-    "axs[1, 0].plot([0, timepts[-1]], [0, 0], 'k--')\n",
-    "axs[1, 0].set_ylim([-5, 15])\n",
-    "axs[1, 0].legend([]);"
+    "ct.time_response_plot(resp_ln, color='b', plot_inputs=False, label=\"linear\")\n",
+    "ct.time_response_plot(resp_cl, color='r', plot_inputs=False, label=\"clipped\")\n",
+    "ct.time_response_plot(resp_aw, color='g', plot_inputs=False, label=\"anti-windup\");"
    ]
   },
   {
@@ -558,21 +572,12 @@
    "outputs": [],
    "source": [
     "resp_aw = ct.input_output_response(\n",
-    "#    clsys, timepts, stepsize, params={'umax': 5, 'kaw': 1, 'a': 100})\n",
     "    clsys, timepts, stepsize, params={'umax': 5, 'kaw': 0.05 * ki, 'a': 100})\n",
     "\n",
     "# Plot the time responses in a single plot\n",
-    "out = ct.time_response_plot(resp_ln, color='b', plot_inputs=False)\n",
-    "ct.time_response_plot(resp_cl, color='r', plot_inputs=False)\n",
-    "ct.time_response_plot(resp_aw, color='g', plot_inputs=False)\n",
-    "\n",
-    "# Annotations\n",
-    "axs = ct.get_plot_axes(out)\n",
-    "axs[0, 0].legend(['linear', 'clipped', 'anti-windup'])\n",
-    "axs[0, 0].plot([0, timepts[-1]], [stepsize, stepsize], 'k--')\n",
-    "axs[1, 0].plot([0, timepts[-1]], [0, 0], 'k--')\n",
-    "axs[1, 0].set_ylim([-5, 15])\n",
-    "axs[1, 0].legend([]);"
+    "ct.time_response_plot(resp_ln, color='b', plot_inputs=False, label=\"linear\")\n",
+    "ct.time_response_plot(resp_cl, color='r', plot_inputs=False, label=\"clipped\")\n",
+    "ct.time_response_plot(resp_aw, color='g', plot_inputs=False, label=\"anti-windup\");"
    ]
   },
   {
