Lab 3:A Disturbance in the Force

The questions below are due on Friday March 27, 2026; 10:00:00 AM.
 
You are not logged in.

Please Log In for full access to the web site.
Note that this link will take you to an external site (https://shimmer.mit.edu) to authenticate, and then you will be redirected back to this page.
You are viewing the PRELIMINARY 6.3100 Spring 2026 web site. Please note that information on this page may change throughout the semester, please check back frequently.

Quick Reference

NOTE

When the sketch is uploaded to the Teensy, it sets K_p, K_d, and K_i to 0.
To get non-zero values, you must set them via the serial plotter send window.

Serial Plotter Send Window Indices and offsets

Links

  • Sketch for this lab is here.
  • Ipython test file is here here.
  • Ipython Arm Model File here.

Task Summary

Week A

  1. Add second motor to your arm, download the new sketch, re-calibrate angle sensor offset, and level the arm. (Checkoff 1).
  2. Set K_i = 4.0, K_p = 0.5, and find the smallest value of K_d for which the system is stable (step response settles to its steady-state after at least seven visible oscillations).
  3. Measure the frequency of the oscillation in the step response, and the magnitude and phase of sinusoidal steady-state responses for several input frequencies(Checkoff 2).
  4. Set the input to zero (6 0 in send window), then disturb the arm by fanning it, and find the "fanning resonance frequency" that maximizes the disturbance response.
  5. Consider how the following are related:
    • the frequency of the oscillations in the step response,
    • the frequency that maximizes the sinusoidal steady state response
    • and the frequency that maximizes the disturbance response.
  6. Set K_i=0, K_p = 2, and for this case find a new smallest K_d for which the step response settles after at least seven visible oscillations. Turn off the input (6 0), and find the fanning resonance frequency.(Checkoff 3).
  7. Download the Ipython test file, init6310ipy.py, to a directory, start ipython in that directory, and at the python prompt type %load init6310ipy.py

Week B

  1. Adjust \gamma_a and \beta so the input-output frequency response of the model matches the steady-state magnitude and phase of your arm's input-output response at a few points. NOTE, these are not the same values from lab2, as you now have two props!
  2. Derive the disturbance to output transfer function, and so that it also computes the disturbance frequency response. Fan your arm to demonstrate that your model's disturbance frequency response accurately predicts the most disturbing frequency.
  3. Use your model (and ipython) to design a controller with better disturbance rejection by examining the disturbance frequency response. Fan the arm to demonstrate the improved disturbance rejection of your better controller.

Double the Fun

Add a second motor to the right side of the arm as shown in the picture below.

Clamp your arm to the side of the lab table, but make sure the plate sticks out beyond the end of the lab bench (to minimize air resistance). See the images below.

For this lab, please position both motors further away from the angle sensor than for last lab, a 14-hole spacing that we can measure using the "disturbing U" from last lab (see picture below).

Rotate the arm from -180 degrees to +180 degrees (-\pi to \pi radians) and check wiring, as shown below.

Check that the arm is nearly balanced. First test that the arm is nearly balanced with the power off off (the arm might or MIGHT NOT stay level when the power is off, but it should "feel" close to balanced). Then connect both power supplies (the 5V USB power supply AND the 5v 10 amp power supply) to the controller board and connect your laptop to the Teensy. Download the sketch linked at the top of this page, upload it to the Teensy, and start the serial plotter. Note that the offsets for K_p, K_d, and K_i are all zero. Below you will change these gain values to get enough feedback for the arm to level itself.

Hold the arm level and ensure the plotter shows a zero angle. If not, adjust the AngleOffset in the Teensy sketch (you may be able to use the same AngleOffset from the week 2 sketch). Once the sensor is calibrated, stop holding the arm, and try to level the arm by adjusting the black potentiometer (in the upper right hand corner of the board). The potentiometer can be used to adjust the relative thrust of right and left propellers, so that the arm holds itself level (set K_p = 0.1 first, otherwise the arm will just drift).

Tip: If your motors are not turning, here are some troubleshooting steps:

  • Is the high-current power supply (round connector) plugged in?
  • Did you download the sketch for Lab 3 and upload to your teensy?
  • Check the connection between each motor and its extension wires. It is very easy to bend the pins in the extension wire connector, THIS IS THE MOST COMMON ISSUE!!! Get new extension wires!
  • Check the connection of the extension leads to the screw terminals on the circuit board. There should be about 3mm of bare silver wire for the connector to grab, and should invisible outside the screw terminal. If uninsulated silver wire is visible outside the screw terminal, that wire could "short" and damage the board.
  • Use a multi-meter (there are about a dozen of them distributed around the lab, ask for help) to check the potential difference across the two terminals at the controller board. Touch the multi-meter voltage probes to the screws and check if there is a non-zero voltage, if not feel free to ask staff for help.

Tip: Make sure that both propellers are blowing air downwards. If they are not, there are two ways to fix this, swapping wires or editing the code. EDIT THE CODE, it is much easier. Find lines 278 and 280 and replace the order of B1L and B2L or A2R and A1R depending which motor needs to be changed, note that these ports are labeled on the PCB.

The figure below shows serial plotter output from running the this lab's sketch on a propeller arm (from prior term of 6.310!!!) with a sinusoidal input (in red) and the resulting measured angle (in blue). Notice the two intervals measured on the plot, one denoted with a black double arrow and the other with a green double arrow. How does the ratio of those two intervals relate to the phase difference between the two plots?

Checkoff 1:
  • Show that the powered arm is level and the plotted angle is zero.
  • How does the ratio of black to green intervals in figure above relate to the phase difference between the input sinusoid (red) and measured angle sinusoid (blue)?

Steps, Sinusoidal and Interference.

Step response.

Set K_p=0.5; K_d=0.3; K_i=4.0 and then type "6 0.1" and "5 0.05" to generate an input that is high for 10 seconds, then low for 10 seconds, and then repeats (Note that the plot window is 10 seconds wide, and the plot points are 20 milliseconds apart). Then try a few different values for K_d. You should be able to find a K_d between 0.1 and 0.3 that produces a step response similar to the one in the following video. Notice that in the video, the step response visibly oscillates about seven times before settling, the step response of your arm with your selected K_d should oscillate similarly (seven visible oscillations).

Again, notice the number of visible oscillations in the above figure, and ensure your arm step response oscillates similarly.

Enter the frequency of oscillation (in Hz) seen in your barely-stable step response:  

Measure Frequency Response

In theory, measuring the frequency response of a STABLE linear system is easy:

  1. Pick a frequency and apply a sinusoid of that frequency to your system's input.
  2. Wait for the system's natural responses to die out (they must die out eventually, since the system is stable). What remains will be sinusoidal at the same frequency as the input.
  3. Compare the input sinusoid to the output sinusoid and compute the amplitude ratio and the phase difference.

We can measure the frequency response for our two-propeller arm system using the above approach as long as we are careful about three issues described below.

Monitor Motor Commands

While the two-propeller arm is intended to cancel the gravity-induced nonlinearity of the one-propeller arm, the cancellation is only perfect if the arm is perfectly symmetric (so balance the arm as best you can!). In addition, there is still a "clipping" nonlinearity that arises when the controller-generated motor commands exceed power supply limitations. This nonlinearity can usually be avoided by reducing the magnitude of the sinusoidal input, provided the clipping has been detected.

Monitor the motor command plot whenever measuring frequency response, to ensure the commands are not being "clipped" by power supply limits. The clipping effect is easy to spot, as shown in the figure below. The black arrows point to cases where the motor command (the green curve) flattens, an indication that the controller's commands are being clipped because of power supply limitations.

The clipping problem can usually be eliminated by reducing the input amplitude, as shown in the picture below. Avoid reducing the input amplitude too much though, because the micro-controller won't be able to represent the signal accurately. Use index 6 in the plotter window to set a reasonable sinusoidal amplitude (do not make the sinusoidal amplitude less negative than -0.03, smaller amplitudes are not converted accurately by the controller board).

NOTE: USE NEGATIVE AMPLITUDES, THE SKETCH WILL GENERATE STEPS INSTEAD SINES IF YOU USE POSITIVE AMPLITUDES.

Wait for Steady-State

It can take several seconds to achieve sinusoidal steady-state, a typical case is shown in the picture below. Be sure to wait long enough (ten seconds or so) before pausing the plotter and taking measurements.

Spread Samples

In last week’s lab, we sent an update to the serial plotter every 10 milliseconds. For this week, the sketch default is to send an update to the serial plotter every 20 milliseconds. At higher sinusoidal frequencies, it may be helpful to decrease the time between plotter updates to measure time more accurately. The control algorithm runs every millisecond, so the fastest serial plotter update is once per millisecond (zero skips). To set the fastest update rate, type "4 -19" into the plotter send window (as was done to generate the plot below). Note that the default offset for number of skips is set to 19 in the code. As another example, to update the serial plotter every 10 milliseconds, we want to skip 9 samples so type "4 -10" (to reduce 19 skips to 9) in the plotter send window.

Checkoff 2

For the second checkoff, please measure the frequency response of your propeller arm, with K_p=0.5, K_i=4.0 and the value for K_d you found above. You should notice that for a very narrow range of input sinusoid frequencies, the arm response increases by five fold. We refer to this phenomenon as a resonance, and the center of the associated frequency range as the resonance frequency.

  1. Remember that the default time between plot points is twenty milliseconds.
  2. Measure the frequency of the oscillation in the step response.
  3. Measure the amplitude ratio (measured amplitude / input amplitude) and phase shift of sinusoidal steady-state responses for several input frequencies. Try to find the a lowest frequency for which the phase difference between input and output is at least 30 degrees (though no more than 60 degrees), and the lowest frequency for which there is a $-135$-degree phase shift. We suggest you use index 6 in the plotter window to set a very small sinusoidal amplitude (-0.03), because, at just the right frequency, the arm response can become quite large.

Checkoff 2:
  • How does phase change with frequency?
  • How much does the amplitude and phase change at frequencies that are twenty percent above and below the frequency of the peak response?
  • What frequency gives you a phase shift by -45 degrees? -135 degrees?
  • What is the relationship between the frequency of the oscillations in the step response and the frequency that maximizes the sinusoidal steady state response?

Once again, with feeling.

Set the arm input to zero (6 0 in the send window), but do not changes the feedback gains (K_i=4, K_p=0.5, and K_d is the just-barely-stable value you found in the last checkoff), and try to "fan" your system into resonance. You should be able to excite a resonance, like in the video below. Note that the resonance frequency for the arm in the video below is probably higher than yours, and you will try the higher resonance frequency case next. At higher resonance frequencies, it can be hard to "match" your fanning to the resonance, please watch until (or skip to) the end to see what that looks like!!

Now try the higher resonance frequency case, set K_i=0, K_p = 2, adjust K_d so that you get a "barely stable" step response like in the video below.

Checkoff 3

For this checkoff,

  1. Set K_i=4, K_p=0.5, and K_d the just-barely-stable value you found in the last checkoff.
  2. TURN OFF THE INPUT, and try disturbing the arm by fanning it. Find the "fanning resonance frequency" that maximizes the response.
  3. Set K_i=0, K_p=2.0, and K_d to the just-barely-stable value for this no-integral-feedback case.
  4. Measure the frequency of the oscillation in the step response, and determine the input frequency that maximizes the sinusoidal steady-state response.
  5. TURN OFF THE INPUT, and try disturbing the arm by fanning it. Find the "fanning resonance frequency" that maximizes the response.
  6. For these two sets of controller gains, what did you observe about the relationship between the frequency of oscillations in the step response, the frequency that maximizes the sinusoidal steady state response, and the fanning frequency that maximizes the disturbance response?

Finally, verify that you can run ipython with the control library by downloading our ipython test file, init6310ipy.py, to a directory. Open a terminal window in that directory (powershell on Windows, terrminal on Macs) and type ipython and return. NOT python but ipython for interactive python.

If ipython is not recognized, use pip to install ipython. Type pip install ipython followed by a return, then pip install control followed by a return. You may also need to install numpy and matplotlib using pip.

If python and pip are not installed on your laptop, but pip3 and python3 are, they may not work (or at least we don't know how to make them work). We suggest installing the anaconda package here, but you can also install python from python.org.

Once you have installed ipython and started it, you should get a prompt. At the prompt type %load init6310ipy.py. Then type return. You should see plots and a command window like in the screenshot below.

Checkoff 3:
  • Demonstrate you can fan your arm into resonance for the two sets of controller gains, and explain what you observed about the relation between the frequencies in the step response, the frequency response, and the disturbance response of your system.
  • Show that Ipython works on your laptop.

STOP HERE FOR WEEK A!!!!

PART B OF THE LAB IS UNDER CONSTRUCTION!!!! UNTIL 3/12/2026!

Continuous Time (CT) Modeling

Last lab we derived an approximate discrete-time (DT) model for the innately continuous-time (CT) propeller arm. We chose discrete-time partly to match our innately discrete-time Teensy-based controller, and partly because discrete-time models, at least the ones formulated as systems of linear difference equations, have explicit solutions in the form of easily computed sums. For this and the next few labs, we will do the reverse. That is, we will derive CT models for innately-CT physical artifacts, and create approximate CT models for our innately-DT controllers.

We are switching to CT because the concepts we are exploring next, transfer functions and frequency response, have a simpler and more intuitive description in CT. And that simplicity will allow us to gain new insights into control system design.

Doubly-Driven Arm

The larger number of forces and torques for the doubly-driven arm are depicted using a more abstract picture (below) than we used for the singly-driven case. In this picture we show the motor commands (C_{nom} and c), the measured and desired arm angles (\theta_a and \theta_d), and the motor thrust (F_{nom} and f_{thrust}) and gravity (F_{gravity}) forces.

This picture also explicitly illustrates the arm torque (\Gamma_{thrust}) due to the difference in right and left motor thrusts. In the last lab, we did not model disturbances, and implicitly assumed the net torque on the arm was due to motor thrust alone. Under that assumption, \Gamma_{thrust} was an unnecessary intermediate in the relation between motor command and arm acceleration. This lab is focused on modeling disturbances, and since our example disturbances (fanning and dropping weights) generate torques, their effect on arm acceleration will be easier to model if we can add them to an explicit \Gamma_{thrust}.

Discrete-Time (DT) to Continuous Time (CT)

In the first propeller-arm lab, we related the change in arm angle, \theta_a, to the arm angular velocity, \omega_a,

\theta_a[n] \approx \theta_a[n-1] + \Delta T \cdot \omega_a[n-1], \;\;or\;\;equivalently,\;\; \frac{\theta_a[n] - \theta_a[n-1]}{\Delta T} \approx \omega_a[n-1],
an approximation to the CT differential equation,
\frac{d\theta_a(t)}{dt} = \omega_a(t),\;\;\;or\;\;\; \theta_a(t) = \int \omega_a(\tau)) d\tau.

Similarly, we related the change in arm angular velocity to the arm angular acceleration, \alpha_a,

\omega_a[n] = \omega_a[n-1] + \Delta T \cdot\alpha_a[n-1], \;\;or\;\;equivalently,\;\; \frac{\omega_a[n] - \omega_a[n-1]}{\Delta T} \approx \alpha_a[n-1],
an approximation to the CT differential equation,
\frac{d\omega_a(t)}{dt} = \alpha_a(t),\;\;\;or\;\;\; \omega_a(t) = \int \alpha_a(\tau)) d\tau.

Block Diagrams

Simple Model, Simple Control

In our simplest arm model, angular acceleration is assumed proportional to motor command, and in the simplest controller, the motor command is assumed proportional to angle error. In discrete time,

\alpha_a[n] = \gamma_a c[n] \;\;\; c[n] = K_p \left( \theta_d[n] - \theta_a[n] \right),
and in continuous time
\alpha_a(t) = \gamma_a c(t) \;\;\; c(t) = K_p \left( \theta_d(t) - \theta_a(t) \right).
where K_p is the gain of the proportional feedback.

We can combine the equations for the arm, from the previous section, and the equations for the motor and controller (above) into a CT block diagram, shown below.

Note that in the CT diagram, we combined blocks with integral operators rather than combining blocks with derivative operators, a choice that parallels our intuition. That is, the system's controller sets the acceleration, which integrates to velocity, which in turn integrates to position, which we then measure and feed back.

Simple Model, PID control

We know that proportional feedback alone can not stabilize the propeller arm, so derivative feedback is essential, and we also know that integral feedback improves steady-state tracking. Including all three (P,I,D) in a CT block diagram is easy, because the integral and derivative relations can be represented explicitly as integral and derivative blocks (see below), rather than by approximations based on differences.

Note that even though we are representing the arm control system with a CT block diagram, the derivative and integral controller terms will still be approximated on the microcontroller using scaled angle error differences and running sums. Nevertheless, we will continue to approximate the DT behavior with a CT representation, because we can model the system with simple block diagrams and develop intuition and insight.

Non-instant Motor Model

As we saw in the first propeller-arm lab, a more accurate model for the relationship between motor command and arm acceleration is

\alpha_a[n] \approx \alpha_a[n-1] + \Delta T \beta \left(\gamma_{a} c[n-1] - \alpha_a[n-1] \right) ,
where \beta and \gamma are both positive. Reorganizing the above difference equation
\frac{\alpha_a[n] - \alpha_a[n-1]}{\Delta T} = \beta \left(\gamma_{a} c[n-1] - \alpha_a[n-1] \right) ,
leads to an analogous differential equation,
\frac{d \alpha_a(t)}{dt} = \beta \left(\gamma_{a} c(\tau) - \alpha_a(\tau) \right),
or equivalently, the integral equation
\alpha_a(t) = \int \beta \left(\gamma_{a} c(\tau) - \alpha_a(\tau) \right) d \tau.
Adding the motor model to the CT block diagram leads to the figure below.

Using CT Transfer functions

In this lab we are transitioning from DT difference equations and system functions to CT differential equations and system functions. The CT system functions are represented by ratios of polynomials in s.

Key identities: derivatives \rightarrow s and integrals \rightarrow\frac{1}{s}.

The transfer function for the plant (from motor command to propeller arm angle), denoted H(s), is a ratio of polynomials in s (aka a rational function of s) given by

H(s) = {\sum_{\tilde{l}=0}^{\tilde{L}_H} b^H_{\tilde{l}} s^{\tilde{l}} \over \sum_{l=0}^{L_H} a^H_l s^l} \equiv \frac{n_H(s)}{d_H(s)}.
In the above formula, \tilde{L}_H and L_H are the orders of n_H and d_H, the numerator and denominator polynomials of the H(s), respectively.

Similarly, the transfer function for the controller, K(s), is the rational function

K(s) = {\sum_{\tilde{l}=0}^{\tilde{L}_K} b^K_{\tilde{l}} s^{\tilde{l}} \over \sum_{l=0}^{L_K} a^K_l s^l} = \frac{n_K(s)}{d_K(s)},
where \tilde{L}_K and L_K are the orders of n_K and d_K, the numerator and denominator polynomials of the K(s), respectively.

The transfer function that relates the measured arm angle to the desired arm angle is denoted G(s), which can be determined by applying Black's formula to the unity gain feedback system applied to controller K(s) and plant transfer function H(s) to yield the rational function describing G(s) as in

G(s) = \frac{K(s)H(s)}{1 + K(s)H(s)} = {\sum_{\tilde{l}=0}^{\tilde{L}_G} b^G_{\tilde{l}} s^{\tilde{l}} \over \sum_{l=0}^{L_G} a^G_l s^l} = \frac{n_G(s)}{d_G(s)},
where we exploit the fact that rational functions (ratios of polynomials) are closed under addition, subtraction, multiplication and division.

CT Advantage

Our design approach in discrete-time has been to minimize steady-state errors, by maximizing K_p, while maintaining stability, which in discrete time means ensuring the natural frequencies are all less than one in magnitude. We could follow the same path in continuous-time, minimizing steady-state errors by maximizing K_p while maintaining stability,

\mbox{Real}\big(\mbox{poles}(G(s)\big) < 0.
But, then switching to CT would be about trading the stability region for DT poles, the interior of the unit disk, for the stability region for CT poles, the open left-half plane. And that does not seem worth the effort.

Switching to CT is not just about replacing the unit circle with the left-half plane, we are also changing our approach to designing controllers. As we will see over the next few weeks, we can design good controllers if we know the actuator's frequency response (its sinusoidal-steady-state response to a sinusoidal inputs as a function of frequency). And in CT, if the poles of H(s) are stable, then the sinusoidal steady-state response to a sinusoidal input (in the form of a complex exponential),

e^{j\omega t}
is given by
H(j \omega) e^{j\omega t}
where H(j \omega) = H(s)|_{s = j \omega} and \omega = 2 \pi f is the frequency in radians/second.

Arm Block Diagram using System Functions

Replacing the integrals in the CT block diagram with \frac{1}{s} is the first step in using system functions, as in the diagram below.

TAKE NOTE THAT WE USE CAPITAL LETTERS when we use system functions. That is, we use capital letters when we substitute, for example,

\alpha_a(t) = A_a e^{st}.

Blocks described by system functions can be composed by multiplication, and internal feedback loops can be simplified using Black's formula, leading to the following diagram below.

What is the command to acceleration system function, H_{c2a}(s), and what is the acceleration to \theta_a system function H_{a2t}(s) (where the "t" stands for theta)? See below if you need a hint.

We can simplify the diagram further by denoting H(s) = H_{c2a}(s) H_{a2t}(s), simplifying the arm block diagram into two system functions, K(s), the controller, and H(s), the actuator, where each is a ratio of polynomials in s, as shown below.

Numerical Tools for CT and Frequency Response.

We will use numerical tools in ipython to help simulate our system. After openning ipython, but sure you load our initialization file init6310ipy.py linked at the top of the page.

You can use ipython to manipulate continuous-time transfer functions. First, you have to tell python to use s as the CT transform variable, by typing (this is in our initialization file)

>> s = ct.tf('s')

After that, you can define a CT transfer function. For example, suppose your system is described by a set of differential equations

\frac{d v(t)}{dt} = x(t), \;\;\; \frac{d w(t)}{dt} = v(t), \;\;\; \frac{d y(t)}{dt} = -0.1y(t) + w(t),
then
Y(s) = \frac{1}{s+0.1}\frac{1}{s}\frac{1}{s} X(s) = H(s) X(s),
where H(s) is the CT transfer function given by
H(s) = \frac{1}{(s+0.1)s^2}.

In ipython, if you type

>> H = (1/(s+0.1))*(1/s)*(1/s)

followed by

>> print(H)

Ipython will response with

H =
 
        1
  -------------
  s^3 + 0.1 s^2

Then if you type

>> pole(H)

ipython will respond with

[0, 0, -0.1000]

as you would expect since in CT, two of your poles are at zero.

Many (feedback, step, etc), and we will also make extensive use of the bode command

bode(H)

For example, consider a system with input x(t) and output y(t) described by a first-order differential equation,

\frac{d y(t)}{dt} = -y(t) + 1000 x(t),
where y is the output and x is the input. Suppose we know that x(t) = \cos(\omega t) = \cos(2\pi f t) where \omega is the frequency in radians per second and f is the frequency in cycles per second. Approximately what is y(t) in the limit as t \rightarrow \infty ?

The system function is easy to compute using the techniques described in the prelab, and is H(s) = \frac{1000}{s+1}. If x(t) = u(t), the unit step, then x(t)=0 for t < 0, and x(t) = 1 for t \ge 0. To plot the response to the unit step using ipython, in the command window type (be sure you have loaded our init6310ipy.py file),

H = 1000/(s+1)

to define the system function to ipython, and then

step(H)

to plot y(t) given x(t) = u(t). Note that y(t) rises to 1000 with a time constant of one second (about two-thirds of the rise in one second).

The frequency response of this first-order system can be plotted using

bode(H)

The initialization file for ipython defines two commands for analyzing feedback. The first command, fdbk takes a single arguement, H, and assuming unity feedback as in

Y = H(s) (Y_d - Y),
and a more general command, feedback for applying a more general form of Black's formula than we have seen so far. The feedback command takes two transfer function arguments, H_1 and H_2 to compute an input-output transfer function. In particular, if
Y = H_1(s) (Y_d - H_2(s)Y),
then

G = feedback(H_1, H_2)

will compute the input-output transfer function

Y = G(s)Y_d,
where
G(s) = \frac{H_1(s)}{1+H_1(s)H_2(s)}.

If we use unity-gain feedback with this first order system, then

x(t) = y_{desired}(t) - y(t)
and by using H_1(s) = K(s)H(s), H_2(s) = 1 in the generalized Black's formula above, the closed-loop system function, G(s), is given by
G(s) = \frac{K(s)H(s)}{1 + K(s)H(s)} = \frac{1000}{s+1001}
which ipython can determine for us using the feedback command

G = feedback(1000/(s+1),1)

To plot the step response of the feedback system,

step(G)

The frequency response also shows the impact of including the first-order system in a closed loop, as can be seen in the plot generated by

bode(G)

What is different about the open-loop to the closed-loop step and frequency responses?

If we want to determine y(t) given x(t) = \cos{\omega t} in the limit of large t (what we call sinusoidal steady state), we can use the frequency response. That is,

y(t) = |G(j\omega)|\cos{(\omega t + \angle G(j\omega))}
The magnitude and angle of G(j\omega) can be determined from the Bode plot, as as shown in the figure below.

Determining response at a particular frequency from the Bode plot.

Implement the model

Determine an expression for Hc2a.

Enter your expression for Hc2a, as a ratio of polynomials in s, in the box below, using Python syntax (e.g., ** for exponentiation) with Greek letters written in English (e.g., alpha for \alpha) and common constants as e, j, pi, etc.

Determine an expression for Ha2t.

Enter your expression for Ha2t, as a ratio of polynomials in s, in the box below, using Python syntax (e.g., ** for exponentiation) with Greek letters written in English (e.g., alpha for \alpha) and common constants as e, j, pi, etc.

Implement your expressions for Hc2a and Ha2t in the matlabscript provided HERE.

Running the script should produce plots like the ones below, showing the standard bode plots with the x axis in radians per second, our rescaled bode plots (with bolder lines) in which the x axis is in hertz (cycles per second) instead of radians, and the system step response.

NOTE: We use log to mean a base-10 log, but by default, log in matlab uses natural logarithm (base-e). You can get base-10 log in matlab by using the log10 command.

NOTE: For some versions of Linux, the on-line matlab produces disconnected plots for the rescaled plots, like the example below. If you see this, you can just use the non-rescaled Bode plots, but remember that the x axis in those plots is in radians per second, not cycles per second, and must be scaled by 2 \pi.

Test 1: Last week, you found the value of K_d for which your system was barely stable when K_i=0.0 and K_p=2.0. Compare those experiment results with the matlab model. Adjust the parameters, \beta and \gamma_a of the model to get a good match of the step responses of the measurement and model. Also compare frequency responses at a few frequencies (pick ones near the frequency response peak). Once you have a good match, validate your model by comparing its predictions to your arm behavior when you increase K_d enough that the step response only oscillates two or three times before settling (increasing K_d by 0.05 should do it). Does your matlab model correctly predict the step response changes you observe?

A good way to compare model and measured data (and to make checkoffs go smoothly) is to re-scale the plots so that the time axes are aligned and the amplitudes are visually similar, as in the figure below.

A good model should be predictive when you change to a higher K_d, as in the comparison below.

Test 2: Test the robustness of your model by changing the K's WITHOUT changing the system parameters \gamma_a or \beta. Try the second set of K's from last week: i.e., K_d was set so that the system was barely stable when K_i=4.0 and K_p=0.5. Compare step responses as well as the resulting resonance frequencies in the measurement and in the model. You will probably notice a poorer match between model and observed step response. What about the frequency response?

Checkoff 4:
  • Show your results for Test 1. Show your plot of the frequency response of the model, and compare results at key frequencies to those in your experiment. Did you get a good match between the experiment and the matlab model?
  • Show your results for Test 2. Do the resonance frequencies in the experiment and model match?
  • In the description of the matlab feedback command, we gave the equation for a generalized feedback system, Y = H_1(s)(Y_d - H_2(s)Y). Draw a block diagram for this generalized feedback system showing the input Y_d, the output Y , a summer, and the placement of the transfer function blocks for H_1 and H_2.

Disturbing Problems

For the propeller arm, we have examined two example disturbances, dropping weight and fanning, as shown in the image below.

Adding disturbance torques (generated by dropping or fanning) to the abstracted description of arm leads to the image below.

It is simpler to think of the disturbance torques as disturbance accelerations (which is just a rescaling), so it can be added to the CT block diagram, as shown below.

Notice that we introduced \alpha_{thrust}, the acceleration due to the thrust torque generated by the motor-propeller pairs, and added a disturbance acceleration, \alpha_{dist} (as might be generated by "dropping" or "fanning"). Then the sum of two accelerations becomes the arm angular acceleration.

What is the G_{dist}(s), the transfer function from A_{dist} to \Theta_a? To determine that transfer function, we suggest you convert the above diagram using system functions (e.g. K(s), H_1(s), and H_2(s) as we did above). Then, rather than trying to use Black's formula directly, try writing two equations and solving. That is, given

\Theta_a = H_{a2t}(s) \left( A_{dist} + A_{thrust} \right),
write an additional equation for A_{thrust} in terms of \Theta_d and \Theta_a, and then solve the two equations to determine a relation between A_{dist} and \Theta_a (assuming \Theta_d = 0).

Checkoff 5:
  • What is G_{dist}(s), the transfer function from the disturbance to arm angle?
  • Modify the matlab script to plot the disturbance frequency response.
  • Show that your model can predict the frequency of maximum disturbance by measuring the fanning frequency that maximizes response for the two test cases above, and see if they are the frequencies that your matlab model predicts.

Checkoff 6 - Design a less disturb-able controller

Now that you have a calibrated numerical model for your system, and can use the model to plot G_{dist}(j\omega) (to within a scale factor), you can pick K_p, K_i and K_d to minimize the maximum disturbance response. Please do so, and show that your new controller is less "disturb-able". What can you say about the disturbance frequency response when K_i > 0?

We also ask that you collect data from running a frequency sweep, see the instructions here. Please set Ki=0.0 and Kp=2.0, and find a value for Kd so that the arm's step response oscillates TWO OR THREE TIMES before settling. You should be able to find a good Kd value between 0.15 and 0.3. We also suggest using the offsets in the Teensy sketch (which start at approximately line 40, see the figure below) so that you do not have to enter the gains using the serial plotter window.

Checkoff 6:
  • Demonstrate your less disturb-able controller by comparing its response to fanning to one of the examples from last week.
  • If you double both K_p and K_d for your controller, how does your measured and modeled disturbance frequency response change? Pay particular attention to the peak response frequency.
  • Collect sweep data for at least 5 sweeps of the arm and SAVE it by following the instructions for running a sweep in the resource guide. Use the default values for the controller, Kp=2.0 and Ki=0.0, and Kd approximately 0.2, set so that a step response oscillates a few times. Plot the data (described in sweeper guide) and show this plot to the staff during the checkoff.