% % Written by: % -- % John L. Weatherwax 2009-04-21 % % email: wax@alum.mit.edu % % Please send comments and especially bug reports to the % above email address. % %----- clc; close all; clear; M_q = -0.568; M_alpha = -17.98; % negated from books number L_alpha_over_V = 1.237; M_dE = -0.175; % negated from books number L_dE_over_V = 0.001; w_n = sqrt( -M_alpha ) xi = ( 0.5/(w_n) ) * ( L_alpha_over_V - M_q ) % the poles of the original transfer function (roots when k=0) s1 = -xi*w_n - w_n * sqrt( xi^2 - 1 ) s2 = -xi*w_n + w_n * sqrt( xi^2 - 1 ) Nk = 1000; kArray = linspace(-100,+100,Nk); S1 = zeros(1,Nk); S2 = zeros(1,Nk); for ki = 1:Nk, % form the coefficients of the polynomial we want the roots of: k = kArray(ki); A = 1; B = 2 * xi * w_n + k * M_dE; C = w_n^2 + k * ( L_alpha_over_V * M_dE - L_dE_over_V * M_alpha ); coef = [ A, B, C ]; rt = flipud( roots( coef ) ); % order these roots so that they are the same as S1(ki) = rt(1); S2(ki) = rt(2); end indsN = find( kArray<0 ); indsP = find( kArray>=0 ); figure; phs1 = plot( real(S1(indsN)), imag(S1(indsN)), '.b' ); hold on; phs1 = plot( real(S1(indsP)), imag(S1(indsP)), 'xb' ); plot( real(s1), imag(s1), '.k', 'markersize', 30 ); phs2 = plot( real(S2(indsN)), imag(S2(indsN)), '.r' ); phs2 = plot( real(S2(indsP)), imag(S2(indsP)), 'xr' ); plot( real(s2), imag(s2), '.k', 'markersize', 30 ); legend( [phs1,phs2], {'root 1(k)','root 2(k)'} ); xlabel('real part'); ylabel('imag part'); grid on; %saveas( gcf, '../../WriteUp/Graphics/Chapter2/sect_6_prob_1_part_a_root_locus.eps', 'epsc' );