- num and den must be of the same length, which is probably not in your case. So, den+num.*k will not work if length(den) is NOT EQUAL to length(num). You can avoid it by using newnum instead of num. See the code for the definition of newnum.
- You must use the matrix as in qs(ii,jj) where ii and jj must be integers. So, if k=0.1, qs(k,:) will result in error.
- roots command may lead to another bug. See the code below.
How to plot Root Locus without using special Matlab functions?
11 views (last 30 days)
Show older comments
Hello,
I need to plot the Root Locus with a chaging "k" of a given transfer function without using any special Matlab functions i.e., "rlocus", "tf". I'm allow to use roots. The code bellow displays an error/warning message (Subscript indices must either be real positive integers or logicals.) that I have not been able to figure out .
See my code.
num = input('Enter the coefficients of numerator of J(s): '); %In vector form
den = input('Enter the coefficients of denominator of J(s): '); %In vector form
qs = 0;
for k = 0:0.1:1000;
qs(k,:) = roots(den + num.*k);
end
plot(qs,'+'), xlabel('\sigma'), ylabel('j\omega'),
title ('Root-Locus'), grid
Thank you
1 Comment
Vivek
on 24 Oct 2019
There are two major bugs in your code.
I suggets the following code that still contains Bug 3 above. (See comments)
%% Plot ROOT LOCI without using the control Systems Tool Box command rlocus()
% my_loci(num, den, kend, usrtitle) plots the Root Loci of the open loop transfer function
% G(s)H(s)=N(s)/D(s)
% num are the coefficients of N(s)
% den are the coefficients of D(s)
%
% Written by P Vivekananda Shanmuganathan.
% Email: Vivek.Shanmuganathan@gmail.com
%
% Oct 24, 2019
% ========================
% CAUTION:
%
% There is a possible error in algorithm. The roots command may possibly
% change the order of the roots as the value of K is changed in each loop.
% This will lead to jump in the loci, i.e. roots of one locus will be
% mistaken for another.:-)
% For example, you can observe this error if you try with num = [3 10] and
% den = 1 1 1 1 ]
% ========================
function locus = my_loci(num, den, kend, usrtitle)
%% Create the matrix locus to hold root loci data.
locus = zeros(kend, length(den)-1);
%% Make dimensions of numerator same as the denominator and fill with zeros.
newnum = [zeros(1,length(den)-length(num)) num];
% If num=[3 1] and den=[1 1 1 1], newnum will be newnum=[0 0 3 1].
% This will make obtaining the coefficients of the characteristic
% polynomial as den+K*num .
kdata=0:kend; % Range of values for the gain K
for K=kdata % K is the Gain
% den+K*newnum is the characteristic polynomial.
% Roots of the chatacteristic polynimal are the poles for given K.
locus(K+1,:)=roots(den+K*newnum)';
end
%% Poles and zeros of the open loop
if length(den>1)
op = roots(den); % Open loop poles
opRe = real(op); % Real part
opIm = imag(op); % Imaginary part
end
if length(num>1)
oz = roots(num); % Open loop zeros
ozRe = real(oz); % Real part
ozIm = imag(oz); % Imaginary part
end
% Mark the open-loop poles and zeros and hold the graph.
clf; % Clear figure
plot(opRe, opIm,'x','LineWidth', 2, 'MarkerSize',15); hold;
plot(ozRe, ozIm,'o','LineWidth', 2, 'MarkerSize',8);
%% Plot the root loci.
% Each element of the locus contains a complex number.
% Plot command plots the imag parts against the real parts
plot(locus); % Plot the root loci
sgrid; % Display the special grid Radial lines indicate zeta and circles indicate w_n
title(usrtitle);
xlabel('Real Axis: Attenuation \sigma=\zeta\omega_n (seconds^-^1)')
ylabel('Imaginary Axis: Natural Frequency j\omega_n (seconds^-^1)')
axis equal % Scale same for both the axes
hold off
end
Answers (1)
manikanta dharmana
on 22 Jul 2017
How to find rootlocus of a transfer fuction wihtout using matlab inbuilt functions??
2 Comments
John D'Errico
on 22 Jul 2017
Please don't answer a question with your own question. You can add a comment, as I just did here.
See Also
Categories
Find more on Classical Control Design in Help Center and File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!