A Disturbance in the Force

Quick Reference

NOTE

The new sketch, 6310_PROP_22_23_SWEEP, sets offsets for the values of K_d and K_p.

The offsets are specified on lines 40

#define offsetKd 0.2

and 42

#define offsetKp 0.5

When offsets are non-zero, they are ADDED to values from the plotter send window. For example, by typing "3 1.5" in the plotter send window to set a new value of K_p, 1.5 will be added to its offset of 0.5, and therefore

By typing "1 -0.1" in the plotter send window to set a new value of K_d, -0.1 will be added to its offset of 0.2, and

Send window values DO NOT ACCUMULATE! So if you type "3 1.5" in the send window, and then type "3 2.0", the resulting K_p will just be the offset added to the latest sent value, K_p = 0.5 + 2.0.

Serial Plotter Send Window Indices and offsets

Links

• Build instructions for adding second motor are at the end of section three in here.
• Sketch for this lab is here.
• Matlab script for week B is here.

Task Summary

Week A

1. Add second motor to your arm, download the new sketch, recalibrate angle sensor offset, and level the arm. (Checkoff 1).
2. Set K_i = 4.0 and find the smallest value of K_d for which the system is stable (step response settles to its steady-state in less than 5 seconds).
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.
4. 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 (Checkoff 2).
6. Repeat 2->6 with K_i=0, K_p = 2, and for this case find a new smallest K_d for which the step response settles in less than five seconds(Checkoff 3).

Week B

Adjust \gamma_a and \beta in the provided matlab model so that 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.
Derive the disturbance to output transfer function, and update the matlab script (at the top of this page) 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.
Use your model (and matlab) 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 (clip it into the counter-weight is on the right side of the arm, see the end of section three in the assembly guide). Slide the position of the second motor on the arm so that the arm is balanced (stays level when the power is off). Connect both power supplies to the controller board and connect your laptop to the Teensy. Download the 6310_PROP_22_23_SWEEP sketch from the class web site (see link above), upload it to the Teensy, and start the serial plotter. Note that the sketch uses non-zero offset values for K_p and K_d, enough feedack that the arm should be close to level.

Transfer your angle sensor calibration (AngleOffset) BUT NOT THE NOMINAL COMMAND to the new sketch, and verify that the measured \theta_a is zero when you hold the arm level. Once the sensor is calibrated, adjust the black potentiometer (in the upper right hand corner of the board) to balance the right and left propeller forces so that the arm holds itself level.

The figure below shows the serial plotter output from running the new sketch 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 unpowered arm is balanced, and 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)?

Ask for Help Ask for Checkoff Loading...

Steps, Sinusoidals and Interference.

Step response.

Set K_i = 4.0 by typing "2 4.0" in the plotter send window, and then type "6 0.1" and "5 0.1" to generate an input that toggles between -0.1 and 0.1 every 5 seconds (two transitions in a ten second period). Reduce K_d by typing, for example, "1 -0.1" in the plotter send window. Since a nonzero offset was specificed in the Teensy sketch, offsetKp = 0.2, adding -0.1 yields K_d=0.1. You should be able to find a K_d between 0.1 and 0.3, that results in step response similar to the one in the following video.

Notice that the step response plotted in the video above is barely stable, but it is stable! Try using your hand to perturb your arm by about 45 degrees, and then let it go. If your arm persists in oscillating wildly, and does NOT eventually level out, your arm is not stable enough, please increase K_d a little (ten percent) and try again.

Measure Frequency Response

In theory, measuring the frequency response of 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 maximums. 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 microcontroller won't be able to represent the signal accurately. Reasonable sinusoidal amplitudes are values are between -0.005 to -0.02, and can be set using index 6 in the plotter send window. NOTE: USE NEGATIVE AMPITUDES, THE SKETCH WILL GENERATE STEPS INSTEAD SINES IF YOU USE POSITIVE AMPLITUDES.

Wait 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

The sketch default is to send an update to the serial plotter every 20 milliseconds, but you can decrease the time between updates and measure time more accurately. The fastest serial plotter update is once a millisecond (zero skips), to set that update rate, type "4 -19" into the plotter send window (as was done to generate the plot below).

Checkoff 2

For the second checkoff,

1. Measure the frequency of the oscillation in the step response
2. Measure the amplitude ratio (measured amplitude / input amplitude) and phase shift of sinusoidal steady-state responses for several input frequencies (at least 10 of them). Make plots of amplitude ratio vs. input frequency and phase shift vs. input frequency.
3. Find the the input frequency that maximizes the amplitude of the output, and note the associated output's phase.

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 is the relationship between the frequency of the oscillations in the step response and the frequency that maximizes the sinusoidal steady state response?

Ask for Help Ask for Checkoff Loading...

Once again, without the integral.

Set K_i=0 by typing "2 0" in the send window, set K_p = 2 by typing "3 1.5" in the send window (see note about offsets at the top of this page), and adjust K_d so that you get a "barely stable" step response like in the video below.

Try fanning your updated system, you should be able to excite a resonance, like in the video below. As you can see in the video, at the faster resonance frequency, it can be hard to "match up" to the resonance, please watch until (or skip to) the end to see what that looks like!!

Checkoff 3

For this checkoff,

1. Measure the frequency of the oscillation in the step response, and verify the relation between that frequency and the frequency that maximizes the magnitude of the sinusoidal steady-state response.
2. Try disturbing the arm by fanning it, find the "fanning resonance frequency" that maximizes the response.
3. Note the connection 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.

Checkoff 3:

• Demonstrate you can fan your arm into resonance, and show its connection to the step and frequency response of your system.

Ask for Help Ask for Checkoff Loading...

STOP HERE FOR WEEK A!!!!

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 best depicted using a more abstract picture than we used for the singly-driven case, as in the figure below.

In the above figure, we repeat the definitions of measured and desired arm angles, \theta_a and \theta_d, respectively, show the motor gravity and thrust forces, F_{nom}, f_{thrust}, and F_{gravity}^M, as well as the motor commands, C_{nom} and c. We also explicitly introduce 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 focussed on modeling disturbances, and since our example disturbances (fanning and dropping weights) are 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,

an approximation to the CT differential equation,

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

an approximation to the CT differential equation,

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,

and in continuous time where K_p is the gain of the proportional feedback.

Note that the DT and CT equations are almost identical, the only difference it the replacement of [n] by (t)? Is that generally true?

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.

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 blocks can be used directly, as shown below.

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

where \beta is assumed negative. Reorganizing the above difference equation leads to an analogous differential equation, 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 state-space, to CT, differential equations and system functions. It is those system functions, ratios of polynomials in s, that help us create blocks we can compose into block diagrams.

key identity: derivatives to s, integrals to \frac{1}{s},

for the actuator: H(s) is a ratio of polynomials in s.

for the controller: K(s) is a ratio of polynomials in s.

and for the feedback system (the well-known Black's formula): G(s) = \frac{K(s)H(s)}{1 + K(s)H(s)}.

CT Advantage

Our design approach in discrete-time has been to minimize steady-state errors, by maximizing K_p, while maintaining stability,

We could follow the same path in continuous-time, minimizing steady-state errors by maximizing K_p while maintaining stabilility, 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 z with s, we are making more than a ...symbolic... change. We will also change 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),

is given by 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,

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 (theta with a "t") system function H_{a2t}(s)? 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.

Matlab Tools for CT and Frequency Response.

You can use matlab to manipulate continuous-time transfer functions, but there are a few differences between s- and z-domain manipulation. First, you have to tell matlab to use s as the CT transform variable, by typing

>> s = tf('s')

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

then where H(s) is the CT transfer function given by

In matlab, if you type

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

Matlab will response with

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

Then if you type

>> pole(H)

matlab will respond with

ans = [0, 0, -0.1000]

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

Many of the matlab commands you've used before translate in the obvious way (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,

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 matlab, in the command window type (be sure you have defined s to matlab using s = tf('s')).

H = 1000/(s+1)

to define the system function to matlab, 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)

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

and by using K(s) = 1 in Black's formula, the closed-loop system function, G(s), is given by which matlab can determine for us using the feedback commandt

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,

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.

Checkoff 4 - Calibrate the model

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.

To get you started, we wrote a matlab script with the CT closed-loop system transfer function, nominal values for the coefficients and a weak controller, 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.

The on-line matlab for some versions of Linux produce disconnected plots for the rescaled plots, like the example below.

If you see this on your laptop 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.

Start by trying to match your measured data to the frequency plots for the K_i = 4.0, K_p = 0.5, and the value of K_d you found last week that made your arm just barely stable (you probably typed something like 1 -0.03 in the send window, which makes K_d=0.17). Use our script to plot frequency response, and compare the matlab-predicted and measured magnitude and angle of the response of your arm at a variety of frequencies. Adjust the parameters, \beta and \gamma_a of the matlab model to get a good match over a range of frequencies. Finally, check that your matlab model, WITHOUT CHANGING \gamma_a OR \beta, correctly predicts the resonance frequency (though probably won't match the amplitude at resonance) for the second case from last week, in which K_i = 0, K_p = 2.0, and K_d was set so the arm was just barely stable.

Checkoff 4:

• Calculate H_{c2a}(s) and H_{a2t}(s) as described above.
• Calibrate your matlab model as described above.
• For the two sets of controller gains you used last week, show, by measuring the response to sinusoids at a few different frequencies, that you picked good values for \gamma_a and \beta in the matlab models.
• Look particularly carefully at the frequencies corresponding to peaks in the system's responses for each gain set. How well did the model and measured magnitudes and phases match?

Ask for Help Ask for Checkoff Loading...

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 \Alpha_{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

write an additional equation for \Alpha_{thrust} in terms of \Theta_d and \Theta_a, and then solve the two equations to determine a relation between \Alpha_{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 to plot the disturbance frequency response. Do you need to recalibrate your model if you are unconcerned about the absolute scaling of the disturbance?
• Show that your matlab model can predict the frequency of maximum disturbance using your model by picking a new set of gains (for example, try K_i=2.0, K_p = 1, and K_d just large enough to make the system just barely stable). Then measure the fanning frequency that maximizes response, and see if it is the same frequency that matlab predicts.

Ask for Help Ask for Checkoff Loading...

Checkoff 6 - Design a less disturbable controller

Now that you have a calibrated matlab 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 "disturbable".

Checkoff 6:

• Demonstrate your less disturbable controller by comparing its response to fanning to one of the examples from last week.
• If K_i > 0, what can you say about the disturbance frequency response?
• 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 your arm AND SAVE IT by following the instructions for running a sweep in the resource guide. Use the default values for the controller, K_p = 0.5, K_d = 0.2 and K_i = 0.

Ask for Help Ask for Checkoff Loading...