# Need help finding the stability of a closed loop system and storing the values.

20 views (last 30 days)
Nimit Chovatia on 24 Apr 2022
Edited: Paul on 28 Apr 2022
I need to loop through the values of Kd (0<Kd<10) and Kp (0<kp<200), for the transfer function. and then created a new closed loop transfer function defined as GCl. Now when looping through this, and testing for the stability of the system in variable B, the new vector x or y is not storing the stable values of Kd and Kp. Basically, is there a way to loop through the values of Kd and Kp for the transfer function GCl and then output the stable values of Kd and Kp? And lastly how can I plot the stable region?
clc
clear
Gred = tf([7.385,476.7],[1,19.39,491.5,3835,0])
Kp = 0:1:200
Kd = 0:1:10
for j=1:10
for i=1:200
C = tf([Kd(j),Kp(i)],[0.005,1]);
num = Gred*C;
den=1+Gred*C;
B = isstable(num,den);
%Check if closed loop map is stable with given pair
%If it is stable, then store Kp, Kd in their vectors:
if B == 0
x = Kd;
y = Kp;%Will overwrite the values of Kp and Kd to the ones we need
end
i = i+1;
end
end
% plot(x,y)
Thank you so much in advance and looking forward to any help!
Paul on 24 Apr 2022
Before going any further can you clarify if the problem is in continuous time or discrete time? I'm asking because Gred and C are formed for continuous time transfer functions, but that use isstable() looks like it is for the numerator and denominator pair of a discrete time transfer function, but then num and den aren't polynomials as would be required for that use case.
Use the feedback() command to form the closed loop system.
Also, Kp has 201 elements and Kd as 11 elements, so the j and i loops need to be adjusted accordingly.

Paul on 24 Apr 2022
Edited: Paul on 28 Apr 2022
Assuming continuous time:
Gred = tf([7.385,476.7],[1,19.39,491.5,3835,0]);
Kp = 0:1:200;
Kd = 0:1:10;
Preallocate the matrix for the result. Preallocate to NaN so it can be checked later, if desired, to make sure all elements are filled in.
Stable = nan(numel(Kp),numel(Kd));
Loop over the compensator gains
for j=1:numel(Kd)
for i=1:numel(Kp)
C = tf([Kd(j),Kp(i)],[0.005,1]);
Clsys = feedback(C*Gred,1);
Stable(i,j) = isstable(Clsys);
end
end
Edit: There are better ways to visualize this than surf().
surf(Kd,Kp,Stable,'EdgeColor','none'),view(2),colorbar
xlabel('Kd')
ylabel('Kp')
title('Yellow = Stable, Blue = Unstable') The blue sliver for Kp = 0 and Kd >= 1 comes about because in those cases Clsys has pole and a zero at s = 0. But isstable() only only looks at the poles and considers a pole at the origin to be unstable. if that's not the desired behavior, then change the code to
Clsys = minreal(feedback(C*Gred,1));
Paul on 28 Apr 2022
Hi Nimit,
Here is the first part of your code
clc
clear
Gred = tf([7.385,476.7],[1,19.39,491.5,3835,0]);
vecKd = [];
vecKp = [];
for Kp = 0:1:200
I changed the step size of Kd to meet the execution time constraint of Answers.
for Kd = 0:10
C = tf([Kd,Kp],[0.005,1]);
Gcl = feedback(Gred*C,1);
B = isstable(Gcl);
%Check if closed loop map is stable with given pair
%If it is stable, then store Kp, Kd in their vectors:
if B == 1
vecKd(end+1) = Kd;
vecKp(end+1) = Kp;
end
end
end
plot(vecKd,vecKp,'.') Here is a problem, vecKp is being used in both arguments
% X=pid(vecKp,0,vecKp);
X = pid(vecKp,0,vecKd);
But the compensator used to derive the valid values of Kp and Kd had another term in the denominator.
C = X/(0.005*tf('s') + 1);
Closed loop systems
Z = feedback(C*Gred,1);
It takes too long on Answers to run all of them, so just run the first few.
t=0:0.01:2;
step(Z(1,1,1,1:10),t); Get the step response information
for ii = 1:10
info(ii) = stepinfo(Z(1,1,1,10));
end
info(5) % for example
ans = struct with fields:
RiseTime: 25.7201 TransientTime: 51.5012 SettlingTime: 51.5012 SettlingMin: 0.9000 SettlingMax: 0.9985 Overshoot: 0 Undershoot: 0 Peak: 0.9985 PeakTime: 93.3688

R2021b

### Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!