# Why is my step response going down instead of up?

10 views (last 30 days)
Justin Goh on 28 Nov 2021
Commented: Justin Goh on 29 Nov 2021
I am trying to plot the step response of a 6th order ITAE system as part of a homework assignment. One of the criteria is to use a undamped natural frequency that is 2 times of one that i used in a previous part of the problem. When I plot the step response, the plot is going down, which seems off.
I am wondering if there is a problem with my code and am hoping somebody can point out my error. My code is segmented into part i), ii) and iii). Part iii) is where I deal with the 6th order ITAE.
% i) Calculate the undamped natural frequency and damping ratio required
% for a standard generic second-order system to achieve (assuming unit step input).
% WHAT are the associated sytem eigenvalues? Plot the unit step response
% for the result and demonstrate how well the design criteria are
% met (normalize your output to ensure the final value is 1.0).
% Display the resulting rise time, peak time, settling time, and percent
% overshoot on your graph. USE Percent overshoot = 3, settling time=3
% Values used from chapter 1.
% Mass values
m1 = 1;
m2 = 2;
m3 = 3;
% Spring Coefficients
k1 = 1;
k2 = 2;
k3 = 3;
k4 = 4;
% Damping Coefficients
c1 = 1;
c2 = 2;
c3 = 3;
c4 = 4;
% A, B, C and D matrix generated from chapter 1.
A = [0 0 0 1 0 0;
0 0 0 0 1 0;
0 0 0 0 0 1;
((-1)*((k1+k2)/m1)) (k2/m1) 0 ((-1)*((c1+c2)/m1)) (c2/m1) 0;
(k2/m2) ((-1)*((k2+k3)/m2)) (k3/m3) (c2/m2) ((-1)*((c2+c3)/m2)) (c3/m2);
0 (k3/m3) ((-1)*((k3+k4)/m3)) 0 (c3/m3) ((-1)*((c3+c4)/m3))];
B = [0 0 0;
0 0 0;
0 0 0;
(1/m1) 0 0;
0 (1/m2) 0;
0 0 (1/m3)];
C = [1 0 0 0 0 0;
0 1 0 0 0 0;
0 0 1 0 0 0];
D = [0 0 0;
0 0 0;
0 0 0];
sys_1 = ss(A,B,C,D);
% Specifying percent overshoot (3%) and settling time (3s):
PO = 3;
ts = 3;
% Finding damping ratio and undamped natural frequency
term = pi^2 + log(PO/100)^2;
zeta = log(PO/100)/sqrt(term) % Damping ratio
wn = 4/(zeta*ts) % Natural frequency from settling time and zeta.
% Setting up generic desired second order system
num1 = wn^2;
den1 = [1 2*zeta*wn wn^2];
% Desired control law eignevalues
DesEig1 = roots(den1)
Des1 = tf(num1,den1) % Create desired systems from num2 and den2
% Plotting
figure();
td = [0:0.01:40.0];
step(Des1,td);
% We will utilize this to display desired info. Commented out for now
stepinfo(Des1);
%% Part ii) Augment desired second order system eigenvalues
% Need 4 additional eigenvalues since this is a 6th order system
% The first additional eigenvalue was chosen to be 10 times the
% real part of the eigenvalue obtained in part i)
DesEig3 = real(DesEig1(1))*10; % We add the real() to isolate real part from complex number
DesEig4 = DesEig3-1;
DesEig5 = DesEig4-1;
DesEig6 = DesEig5-1;
%% Part iii) Computate the ITAE 6th-order coefficients & Eigenvalues.
% Undamped natural freq is 2 times that of part i).
% Plot the 6th order ITAE and the augmented 6th order step response with
% the dominant second order step response of part i).
w_itae = wn*2
num_itae = (w_itae)^2; % NEEDS TO BE 2 TIMES w used in part i). Question is which w? or all?
% 6th order ITAE
den_itae = [1 3.25*wn 6.6*wn^2 8.6*wn^3 7.45*wn^4 3.95*wn^5 wn^6];
% Finding Eigenvalues from ITAE characteristic polynomial
DesEig_Itae = roots(den_itae)
Itae_tf = tf(num_itae,den_itae);
figure();
td = [0:0.01:60.0];
step(Itae_tf,td);

Pat Gipper on 28 Nov 2021
I would focus on your equation for zeta on line 60. It is resulting in a negative value (negative damping=bad!) which gives unstable eigenvalues in the remainder of your code.
Justin Goh on 29 Nov 2021
Thank you so much. My textbook even had absolute values around the zeta.