Digital Control Systems: Difference between revisions
Line 18: | Line 18: | ||
===MATLAB/Octave=== | ===MATLAB/Octave=== | ||
*[http://laris.fesb.hr/digitalno_vodjenje/octave/doc/octave_30.html Octave Control Systems Toolbox] This is not the same thing that is on [http://octave.sourceforge.net/control/overview.html Octave Forge here]. | *[http://laris.fesb.hr/digitalno_vodjenje/octave/doc/octave_30.html Octave Control Systems Toolbox] This is not the same thing that is on [http://octave.sourceforge.net/control/overview.html Octave Forge here]. | ||
====Pendulum Specific Scripts==== | |||
=====Single Pendulum===== | |||
<nowiki> | |||
% Double Pendulum Parameters (Tentative: There are two pendulum setups, each with different parameters. I'm not sure which these go to.) | |||
% This script is for balancing only the long rod. | |||
% Run parameters | |||
%f = input('Control Frequency (Hz) = '); | |||
%Trun = input('Run Time (s) = '); | |||
%f=130; | |||
f=1000; | |||
Trun=30; | |||
kmax = round(f*Trun); | |||
T = 1/f; | |||
Maxpos = 0.25; % Max carriage travel +- 0.25 m | |||
Maxangle = 0.175; % Max rod angle -- 10 deg | |||
Maxvoltage = 20; % Max motor voltage, V | |||
pstart = 0.005; % Carriage position starting limit, m | |||
astart = 1*pi/180; % Angle starting limit, rad | |||
g = 9.81; % m/s^2 Gravitational constant | |||
% SYSTEM PARAMETERS | |||
% Measured Mechanical Parameters | |||
d1 = 0.323; % m Length of pendulum 1 (long) | |||
d2 = 0.079; % m Length of pendulum 2 (short) | |||
%mp1 = 0.0208; % kg Mass of pendulum 1 | |||
mp1 = 0.0318; | |||
%mp2 = 0.0050; % kg Mass of pendulum 2 | |||
mp2 = 0.0085; | |||
%m=.5; | |||
m = 0.3163; % kg Mass of carriage | |||
rd = 0.0254/2; % m Drive pulley radius | |||
md = 0.0375; % kg Mass of drive pulley (cylinder) | |||
%mc1 = 0.0036; % kg Mass of clamp 1* | |||
%mc2 = 0.0036; % kg Mass of clamp 2* | |||
mc1 = 0.0085; | |||
mc2 = mc1; | |||
% *Clamp Dimensions | |||
% Rectangular 0.0254 x 0.01143 m | |||
% The pivot shaft is 0.00714 m from the end | |||
% Motor Parameters (Data Sheet) | |||
Im = 43e-7; % kg m^2/rad Rotor moment of inertia | |||
R = 4.09; % ohms Resistance | |||
kt = 0.0351; % Nm/A Torque constant | |||
ke = 0.0351; % Vs/rad Back emf constant | |||
% Derived Mechanical Parameters | |||
% kg m^2/rad Moment of inertia, clamp 1 | |||
%Ic1 = mc1*(0.01143^2 + 0.0254^2)/12 + mc1*(0.0127-0.00714)^2; | |||
Ic1 = mc1*(0.0098^2 + 0.0379^2)/12; | |||
Ic2 = Ic1; % kg m^2/rad Moment of inertia, clamp 2 | |||
Id = md*(rd^2)/2; % kg m^2/rad Moment of inertia, drive pulley | |||
Imd = Im + Id; % kg m^2/rad Moment of inertia, combined | |||
J1 = Ic1 + mp1*(d1^2)/3; % Total moment of inertia, pendulum 1 (long) | |||
J2 = Ic2 + mp2*(d2^2)/3; % Total moment of inertia, pendulum 2 (short) | |||
Jd = Im + Id; % Total moment of inertia, motor drive | |||
Mc = m + mc1 + mc2; % Total carriage mass | |||
% Friction Test Data | |||
% Carriage Slope = 19 deg; Terminal Velocity xdotss = 0.312 m/s; From | |||
% twincarriage.m; formula b = m g sin(theta)/xdotss | |||
% Pendulum 1 (long) Exponent a1 = 0.0756 1/s; From longfit.m | |||
% Pendulum 2 (short) Exponent a2 = 0.2922 1/s; From shortfit.m | |||
% formula b = 2 a J | |||
%alpha = 19; | |||
alpha = 12.2; | |||
%xdotss = 0.312; | |||
xdotss = 0.4852; | |||
%a1 = 0.0756; | |||
%a2 = 0.2922; | |||
a1 = 0.0185; | |||
a2 = 0.012; | |||
% Ns/m Viscous friction of carriage system | |||
b = (Mc + mp1 + mp2)*g*sin(alpha*pi/180)/xdotss; | |||
b1 = 2*a1*J1; % Nms/rad Viscous friction of pendulum 1 (rotational) | |||
b2 = 2*a2*J2; % Nms/rad Viscous friction of pendulum 2 (rotational) | |||
scale = [rd*2*pi/4096 2*pi/4096 -0.05/250]; | |||
T = 1/f; | |||
% The data above comes from the fweb wiki. | |||
M=Mc; %mass of cart | |||
m=mp1; %mass of pendulum 1 | |||
b=b; %friction | |||
l=d1/2; %length of pendulum | |||
I=J1; %inertia of pendulum | |||
%q=(M+m)*(l+m*l^2)-(m*l)^2; | |||
%num=[m*l,0]; %numerator for transfer function | |||
%den=[q,b*(l+m*l^2),-m*g*l*(M+m),-b*m*g*l]; %denominator for transfer function | |||
%[A,B,C,D]=tf2ss(num,den) | |||
%A,B,C,D matricies for the state space model | |||
% x_vec is [x, x_dot, theta, theta_dot]' | |||
% See the web site: https://www.library.cmu.edu/ctms/ctms/examples/pend/invpen.htm | |||
A=[ 0 1 0 0; | |||
0 ((-(I+m*l^2)*b)/(I*(M+m)+M*m*l^2)) ((m^2*g*l^2)/(I*(M+m)+M*m*l^2)) 0; | |||
0 0 0 1; | |||
0 ((-m*l*b)/(I*(M+m)+M*m*l^2)) ((m*g*l*(M+m))/(I*(M+m)+M*m*l^2)) 0]; | |||
B=[ 0; | |||
((I+m*l^2)/(I*(M+m)+M*m*l^2)); | |||
0; | |||
((m*l)/(I*(M+m)+M*m*l^2))]; | |||
C=[ 1 0 0 0; | |||
0 0 1 0]; | |||
D=[ 0; | |||
0]; | |||
cont_sys = ss(A,B,C,D) | |||
rank_ctrb = rank(ctrb(A,B)) | |||
original_poles = eig(A) %poles for our system | |||
</nowiki> | |||
===Z Transforms=== | ===Z Transforms=== |
Revision as of 10:51, 19 May 2014
Links
- Textbooks
- Introduction to Applied Digital Control, Greg Starr, University of New Mexico
- Greg's ME481/ME581 web pages contain solutions to the book problems and other things.
- Optimal Sampled-Data Control Systems, Chen & Francis
- Digital Control System Analysis and Design, 3rd Ed., Philips
- Control Systems and Control Engineering
- Introduction to Signal Processing
- Zoran Gajic's Book on Control Systems (Continuous & Discrete) and some other documents. Look for the documents entitled Chapter N.
- Introduction to Applied Digital Control, Greg Starr, University of New Mexico
- Lectures
- Sampled Time Control
- Zoran Gajic's notes on Observers. This is a very short set of lecture notes that explains the why's of observers, with several interesting derivations.
- Papers
MATLAB/Octave
- Octave Control Systems Toolbox This is not the same thing that is on Octave Forge here.
Pendulum Specific Scripts
Single Pendulum
% Double Pendulum Parameters (Tentative: There are two pendulum setups, each with different parameters. I'm not sure which these go to.) % This script is for balancing only the long rod. % Run parameters %f = input('Control Frequency (Hz) = '); %Trun = input('Run Time (s) = '); %f=130; f=1000; Trun=30; kmax = round(f*Trun); T = 1/f; Maxpos = 0.25; % Max carriage travel +- 0.25 m Maxangle = 0.175; % Max rod angle -- 10 deg Maxvoltage = 20; % Max motor voltage, V pstart = 0.005; % Carriage position starting limit, m astart = 1*pi/180; % Angle starting limit, rad g = 9.81; % m/s^2 Gravitational constant % SYSTEM PARAMETERS % Measured Mechanical Parameters d1 = 0.323; % m Length of pendulum 1 (long) d2 = 0.079; % m Length of pendulum 2 (short) %mp1 = 0.0208; % kg Mass of pendulum 1 mp1 = 0.0318; %mp2 = 0.0050; % kg Mass of pendulum 2 mp2 = 0.0085; %m=.5; m = 0.3163; % kg Mass of carriage rd = 0.0254/2; % m Drive pulley radius md = 0.0375; % kg Mass of drive pulley (cylinder) %mc1 = 0.0036; % kg Mass of clamp 1* %mc2 = 0.0036; % kg Mass of clamp 2* mc1 = 0.0085; mc2 = mc1; % *Clamp Dimensions % Rectangular 0.0254 x 0.01143 m % The pivot shaft is 0.00714 m from the end % Motor Parameters (Data Sheet) Im = 43e-7; % kg m^2/rad Rotor moment of inertia R = 4.09; % ohms Resistance kt = 0.0351; % Nm/A Torque constant ke = 0.0351; % Vs/rad Back emf constant % Derived Mechanical Parameters % kg m^2/rad Moment of inertia, clamp 1 %Ic1 = mc1*(0.01143^2 + 0.0254^2)/12 + mc1*(0.0127-0.00714)^2; Ic1 = mc1*(0.0098^2 + 0.0379^2)/12; Ic2 = Ic1; % kg m^2/rad Moment of inertia, clamp 2 Id = md*(rd^2)/2; % kg m^2/rad Moment of inertia, drive pulley Imd = Im + Id; % kg m^2/rad Moment of inertia, combined J1 = Ic1 + mp1*(d1^2)/3; % Total moment of inertia, pendulum 1 (long) J2 = Ic2 + mp2*(d2^2)/3; % Total moment of inertia, pendulum 2 (short) Jd = Im + Id; % Total moment of inertia, motor drive Mc = m + mc1 + mc2; % Total carriage mass % Friction Test Data % Carriage Slope = 19 deg; Terminal Velocity xdotss = 0.312 m/s; From % twincarriage.m; formula b = m g sin(theta)/xdotss % Pendulum 1 (long) Exponent a1 = 0.0756 1/s; From longfit.m % Pendulum 2 (short) Exponent a2 = 0.2922 1/s; From shortfit.m % formula b = 2 a J %alpha = 19; alpha = 12.2; %xdotss = 0.312; xdotss = 0.4852; %a1 = 0.0756; %a2 = 0.2922; a1 = 0.0185; a2 = 0.012; % Ns/m Viscous friction of carriage system b = (Mc + mp1 + mp2)*g*sin(alpha*pi/180)/xdotss; b1 = 2*a1*J1; % Nms/rad Viscous friction of pendulum 1 (rotational) b2 = 2*a2*J2; % Nms/rad Viscous friction of pendulum 2 (rotational) scale = [rd*2*pi/4096 2*pi/4096 -0.05/250]; T = 1/f; % The data above comes from the fweb wiki. M=Mc; %mass of cart m=mp1; %mass of pendulum 1 b=b; %friction l=d1/2; %length of pendulum I=J1; %inertia of pendulum %q=(M+m)*(l+m*l^2)-(m*l)^2; %num=[m*l,0]; %numerator for transfer function %den=[q,b*(l+m*l^2),-m*g*l*(M+m),-b*m*g*l]; %denominator for transfer function %[A,B,C,D]=tf2ss(num,den) %A,B,C,D matricies for the state space model % x_vec is [x, x_dot, theta, theta_dot]' % See the web site: https://www.library.cmu.edu/ctms/ctms/examples/pend/invpen.htm A=[ 0 1 0 0; 0 ((-(I+m*l^2)*b)/(I*(M+m)+M*m*l^2)) ((m^2*g*l^2)/(I*(M+m)+M*m*l^2)) 0; 0 0 0 1; 0 ((-m*l*b)/(I*(M+m)+M*m*l^2)) ((m*g*l*(M+m))/(I*(M+m)+M*m*l^2)) 0]; B=[ 0; ((I+m*l^2)/(I*(M+m)+M*m*l^2)); 0; ((m*l)/(I*(M+m)+M*m*l^2))]; C=[ 1 0 0 0; 0 0 1 0]; D=[ 0; 0]; cont_sys = ss(A,B,C,D) rank_ctrb = rank(ctrb(A,B)) original_poles = eig(A) %poles for our system
Z Transforms
- Relationship between the Laplace and Z transforms
- Convolution and Z Transforms
- Z Transforms and Convolution
- Here is an animation of convolution with continuous signals. Look at it so you understand what is happening with the mechanics of convolution. To convolve with , you flip shift into on the axis, then you multiply it by to get , then you integrate with respect to , so that the convolution is: . The first animation shows this happening with sampled waveforms: and .
- Notice that this is really the same as Polynomial Multiplication.
Discretization
- This is how the c2d zero order hold works in MATLAB/octave. It uses the "exact" solution to the discretization problem.
- Also check Chapter 7.3 of Greg Starr's book.