How to find root locus of a polynomial

2 次查看（过去 30 天）
Supreeth D K
Supreeth D K 2023-12-9
评论： Supreeth D K 2023-12-9
I need to find the root locus of a polynomial.
m_s = 1; % Mass of the spring
m_r = 68; % Mass of the rotor
k_r = 5.45e5; % Spring constant of the rotor
c_s = 663; % Damping coefficient of the spring
c_r = 900; % Damping coefficient of the rotor
k_s = 15000; % Spring constant of the spring
% Frequency range in rad/s
omega = 0:100:60000;
% Initialize arrays to store roots
roots_real = zeros(length(omega), 5);
roots_imag = zeros(length(omega), 5);
% Initialize array to store crossing points
crossing_omega = zeros(1, 5);
crossing_from_real_to_imag = zeros(1, 5);
% Calculate and store roots for each angular frequency
for i = 1:length(omega)
a4_omega = c_s./m_s + k_r./c_r - 1i.*omega(i);
a3_omega = k_r./m_r + k_r./m_s + k_s./m_s + c_s./m_s*(k_r./c_r - 1i.*omega(i));
a2_omega = (k_r.*c_s) ./ (m_r.*m_s) + (k_r.*k_s) ./ (m_s.*c_r) - (k_r./m_r + k_r./m_s + k_r./m_s).*1i.*omega(i);
a1_omega = k_r ./m_r * (k_s./m_s - 1i.*omega(i).*c_s./m_s);
a0_omega = -1i.*omega(i).*k_s.*k_r ./(m_s.*m_r);
a = [1, a4_omega, a3_omega, a2_omega, a1_omega, a0_omega];
% Calculate roots
system_roots = roots(a);
% Store real and imaginary parts of the roots
roots_real(i, :) = real(system_roots);
roots_imag(i, :) = imag(system_roots);
I would like to plot the root locus for this system as ω varies. Could you please assist me in implementing this in MATLAB? Specifically, I am looking for guidance on expressing the characteristic equation in terms of s and ω, and how to use MATLAB functions to generate the root locus plot

请先登录，再回答此问题。

采纳的回答

Sulaymon Eshkabilov
Ran in:
Here is how to get the roots plotted.
m_s = 1; % Mass of the spring
m_r = 68; % Mass of the rotor
k_r = 5.45e5; % Spring constant of the rotor
c_s = 663; % Damping coefficient of the spring
c_r = 900; % Damping coefficient of the rotor
k_s = 15000; % Spring constant of the spring
% Frequency range in rad/s
omega = 0:10000;
% Initialize arrays to store roots
roots_real = zeros(length(omega), 5);
roots_imag = zeros(length(omega), 5);
% Initialize array to store crossing points
crossing_omega = zeros(1, 5);
crossing_from_real_to_imag = zeros(1, 5);
% Calculate and store roots for each angular frequency
for i = 1:length(omega)
a4_omega = c_s./m_s + k_r./c_r - 1i.*omega(i);
a3_omega = k_r./m_r + k_r./m_s + k_s./m_s + c_s./m_s*(k_r./c_r - 1i.*omega(i));
a2_omega = (k_r.*c_s) ./ (m_r.*m_s) + (k_r.*k_s) ./ (m_s.*c_r) - (k_r./m_r + k_r./m_s + k_r./m_s).*1i.*omega(i);
a1_omega = k_r ./m_r * (k_s./m_s - 1i.*omega(i).*c_s./m_s);
a0_omega = -1i.*omega(i).*k_s.*k_r ./(m_s.*m_r);
a = [1, a4_omega, a3_omega, a2_omega, a1_omega, a0_omega];
% Calculate roots
system_roots = roots(a);
% Store real and imaginary parts of the roots
roots_real(i, :) = real(system_roots);
roots_imag(i, :) = imag(system_roots);
end
%%
figure
plot(roots_real(:,1), roots_imag(:,1), 'ro')
hold all
plot(roots_real(:,2), roots_imag(:,2), 'bo')
plot(roots_real(:,3), roots_imag(:,3), 'mx')
plot(roots_real(:,4), roots_imag(:,4), 'k*')
plot(roots_real(:,5), roots_imag(:,5), 'kd')
ylim([-20 50])
xlim([-16 6])
ylim([-1500 1e4])
xlim([-1050 200])
xlabel('Im(s)')
ylabel('Re(s)')
figure
plot(roots_real(:,1), roots_imag(:,1), 'ro')
hold all
plot(roots_real(:,2), roots_imag(:,2), 'bo')
plot(roots_real(:,3), roots_imag(:,3), 'mx')
plot(roots_real(:,4), roots_imag(:,4), 'k*')
plot(roots_real(:,5), roots_imag(:,5), 'kd')
ylim([-20 50])
xlim([-16 6])
xlabel('Im(s)')
ylabel('Re(s)')

更多回答（1 个）

Sulaymon Eshkabilov
Ran in:
Step 1. Derive the transfer function of the system. E.g.: T = tf(1, [1 2 3])
Step 2. Use rlocus() or rlocusplot(). E.g.: rlocus(T)
m_s = 1; % Mass of the spring
m_r = 68; % Mass of the rotor
k_r = 5.45e5; % Spring constant of the rotor
c_s = 663; % Damping coefficient of the spring
c_r = 900; % Damping coefficient of the rotor
k_s = 15000; % Spring constant of the spring
% Frequency range in rad/s
omega = 0:100:6000; % Smaller range
% Initialize arrays to store roots
roots_real = zeros(length(omega), 5);
roots_imag = zeros(length(omega), 5);
% Initialize array to store crossing points
crossing_omega = zeros(1, 5);
crossing_from_real_to_imag = zeros(1, 5);
% Calculate and store roots for each angular frequency
for i = 1:length(omega)
a4_omega = c_s/m_s + k_r/c_r - omega(i);
a3_omega = k_r/m_r + k_r/m_s + k_s/m_s + c_s/m_s*(k_r/c_r - omega(i));
a2_omega = (k_r*c_s) / (m_r*m_s) + (k_r*k_s)/ (m_s*c_r) - (k_r/m_r + k_r/m_s + k_r/m_s)*omega(i);
a1_omega = k_r /m_r * (k_s./m_s - omega(i)*c_s/m_s);
a0_omega = -omega(i)*k_s*k_r/(m_s*m_r);
a = [1, a4_omega, a3_omega, a2_omega, a1_omega, a0_omega];
SYS = tf(1, a);
rlocus(SYS), hold all
% % Calculate roots
% system_roots = roots(a);
% % Store real and imaginary parts of the roots
% roots_real(i, :) = real(system_roots);
% roots_imag(i, :) = imag(system_roots);
end

类别

在 Help Center 和 File Exchange 中查找有关 Guidance, Navigation, and Control (GNC) 的更多信息

标签

Community Treasure Hunt

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

Start Hunting!

Translated by