Labwork
MECH 4310: Systems and Control Case Study, Classes 26-28
Lab 4: Controller Design Case Study
You will submit a single report for your group and should be named lab4_abc123456_def654321.pdf where you should replace your group members NetIDs in the obvious spots. This file will be submitted through eLearning. Note that you do not need to submit your Matlab code for this lab, although you are expected to use Matlab to assist with calculations and making the figures.
For this lab we will consider a satellite-attitude control problem. As satellites are flying through space (a vacuum) they need something to help stabilize their orientation so that, for example, their communications dish points towards the earth or so that the Hubble Telescope points at the distant star researchers want to observe. Aircraft in our atmosphere achieve this sort of orientation stabilization simply through the aerodynamics of flying through air - there’s no risk of a commercial airline flying through the air belly first, there’s just too much drag to stay in that orientation!
Reaction wheels are one of the standard methods used for reorienting spacecraft. Effectively a flywheel sits inside the spacecraft and by either speeding up or slowing down the wheel, opposite angular momentum is transferred to the satellite1. To keep track of the attitude movements that have been executed, craft typically employ some sort of gyroscope to provide an inertial reference point for measurement2. Because the gyroscopic sensor is a dynamical system in and of itself, the measurement equation is also dynamic (an ODE - see y equation below).
In this problem we will consider a single axis of rotation, but typically there are multiple flywheels (3 are used to get arbitrary rotations in three dimensions). The equations of motion are given as:
Satellite : Iθ̈(t) = u(t) + w(t), Flywheel : Jν̇(t) = −u(t), Measurement : ẏ(t) = θ(t) − ay(t)
where
θ = satellite position angle,
I = satellite moment of inertia (1000 kg/m2),
u = control torque,
w = disturbance torque,
ν = flywheel (angular) speed,
J = flywheel moment of inertia,
y = measurement from the sensor,
a = sensor constant (1 rad/sec),
and the control torque is generated such that it makes the measurement y(t) track a reference angle r(t).
(1) Draw the block diagram of the feedback control system.
(2) Suppose proportional control is used. Create a plot of the root locus with respect to the control gain, k, and characterize the stability of the closed loop system for different values of k.
1For more information watch: https://www.youtube.com/watch?v=EXnqTtZ5pW0. 2Another handy video: https://www.youtube.com/watch?v=XB6gTJBa1-I.
1
MECH 4310: Systems and Control Case Study, Classes 26-28
(3) Suppose now we aim to use a lead compensator which has its pole at s = −1. Our goal is to satisfy the specifications of ωBW = 0.04 rad/s and overshoot of 16.3% (Mp = 0.163) by selecting the location of the lead compensator zero and the gain. We will take an analytic approach for this problem to see another way that we can achieve specifications (in contrast to moving things around in Matlab).
We will use the approximation that ωn (the undamped natural frequency of a pair of complex poles) is close to the cross-over frequency (ωc), which is our typical conservative proxy for closed loop bandwidth. Thus we will take ωn = 0.04 rad/s. Show the magnitude Bode plot for the generic second order system with a pair of complex poles,
G(s) = ω2n
s2 + 2ζωns + ω2n ,
and annotate the location of ω = ωn. This can be a sketch or can be a Matlab plot. Use this plot to justify why ωn ≈ ωc.
Using the overshoot to calculate the damping ζ, determine the desired location of the closed loop complex poles p1,2 = −σ ± iω, where ω2n = σ2 + ω2 and ζ = σ/ωn. We need to place the zero z of the lead compensator so that the root locus passes through the desired pole locations, i.e., that the points p1,2 = −σ ± iω are on the root locus. To do this, we employ the phase constraint of the root locus to determine the angle from z to p1:
m∑ i=1
ψi − n∑
i=1
ϕi = −180◦,
which, as you recall, says that sum of the angles from each zero to the point on the root locus (ψi) minus the sum of the angles from each pole to the point on the root locus (ϕi) must sum to a multiple of −180◦ for the point to be on the root locus. Having found the angle from z to p1 you can use geometry to determine the location/value of z.
By placing the lead zero at this location (using the phase constraint of root locus), we have guaranteed that the poles we want are on the root locus. Our final task in tuning the controller is to find the gain k that positions the closed loop poles at p1,2. Here we use the magnitude constraint of root locus. Recall that a point s0 is on the root locus if kL(s0) = −1, hence |kL(s0)| = 1. This means that the gain k should be selected such that |kL(p1)| = 1. Report the final lead compensator and show the steps used to find z.
(4) Plot the root locus that corresponds to the lead compensator you just created. Zoom into the complex poles you designed and add the design requirements by right clicking the root locus and adding the constraints on ωn and ζ or Mp. Submit this plot. For what values of the gain k is the closed loop system stable?
(5) With the lead compensator and specific gain value designed above, compute the phase and gain margins. Identify the crossover frequency ωc - how does it compare with the desired ωn ≈ ωc? Now add another plot to the sisotool: add the closed loop Bode plot, i.e., r2y. From this plot identify the true closed loop band- width ωBW. Is ωBW ≈ ωn? Explain why this value is reasonable. Submit both the open loop and closed loop Bode plots to support your reported values.
(6) What is the steady-state error caused by a step in the disturbance w(t)? Provide the number and a plot of the step response. This plot can be generated by adding a new plot within the sisotool in Matlab. Matlab refers to the disturbance as du, so you are interested in the step response du2y. What is the system type with respect to rejection of the disturbance w(t)?
(7) To help build an intuition about closed loop bandwidth and the frequency content of a signal, we will now evolve the closed loop system (with the lead compensator defined above) with various sinusoids (and noise)
2
MECH 4310: Systems and Control Case Study, Classes 26-28
to see the response of the system. This information is already contained within the closed loop Bode plot, but taking the following steps will help to make it more accessible. Our goal is to visualize the response of the system when the reference changes at different rates (varying the frequency of the sinusoid). This is somewhat similar to the last section of Lab 2.
Take a look at the supplied stub code. The only changes you need to make are in the first four lines to:
1. (T) Define the closed loop transfer function (using the tuned lead compensator).
2. (freqs) Define 4-6 frequency values that exhibit notable behavior on the closed loop Bode plot (e.g., the frequency at which the magnitude plot peaks, crosses 0dB, etc). Make sure to also pick (relatively) low, medium, and high frequency values.
3. (dB) Provide the corresponding magnitudes (in dB) for the frequencies chosen above, read from the closed loop Bode magnitude plot.
4. (desc) Provide a very short description of what makes the frequency interesting. Enter these in as strings (e.g., 'peak frequency'); these are automatically included in the titles of the figures that will be generated.
The code in the for loop does two things.
• For a given frequency f, it evolves the closed loop satellite system forward in time, driven by a refer- ence r that is changing in a sinusoidal pattern with frequency f to produce the output y. It plots the reference signal in black and the output signal in blue. It overlays a dashed red line that indicates the expected amplitude of the output based on the dB magnitude from the closed loop Bode plot (the dB magnitudes you entered above). Recall that a sinusoid with amplitude M corresponds to a magnitude of 20 log10(M) decibels.
• It also plots the frequency content of the reference signal using a Fast Fourier Transform (fft). The theory behind Fourier Transform and Fourier Series says that any periodic signal can be approximated arbitrarily well by a bunch of well-selected sin and cos functions at certain frequencies. The plot we show here presents the amount of each frequency that exists in the various reference signals we are using. In the cases defined so far, there is just a single frequency present, and the fft should have a single peak representing that single frequency component.
The code after the for loop makes things more interesting.
• First, we introduce a signal that is composed (literally added) of two sinusoids at different frequencies. I have preselected these frequencies, ω = 0.1 rad/s and ω = 2 rad/s. In the plot it is easy to observe that the reference is composed of two sinusoids. Explain why the output in this scenario makes sense even though it looks like the trajectory of a single sinusoid. Identify the peak(s) associated with this reference input on the plot - explain why this signal is different from the others. In particular, dig into the code, specifically the definition of r to justify why the height of the peaks is different.
• Second, we introduce noise in the reference3. Looking at the signal r(t) (the black curve), the fast movement of the signal is easy to see. The fft plot decomposes this signal into the (many) sinusoids that are inside and the result is greatly different from the other reference inputs. Use this plot to justify the shape of the response of the system.
Submit the two figures that are generated by this analysis. On your closed loop Bode plot, annotate the frequencies that you have selected to use.
3While we typically don’t see noise in the reference signal, we do see noise present in the measurement and this mainly meant to demonstrate the idea of the frequency content of a noisy signal.
3