Avoid root switching using roots

Question

jf 2021-9-13

0
链接

此问题的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/1452344-avoid-root-switching-using-roots

编辑： jf 2021-9-13

在 MATLAB Online 中打开

Hello, I am having following problem.

It is a very simple.

First formulate symbolically some matrix

.

Then I wish to study

. In particular, I am interested only in zeros crossings of imaginary axis

while making some sweeps of variables in the matrix, namely

gm2
gm3

So far so simple.

So I get numerator and denumerator of the resulted determinant in each iteration and then find their roots using function roots.

But the problem is, that in the last part of the code

    for pp = 1:size(real_parts_of_poles,1)
        polecross_index = zci(real_parts_of_poles(pp,:));
        
        if ~isempty(polecross_index)
        omega = imag_parts_of_poles(pp,polecross_index);
        MM = [omega, gm2, gm1vec(polecross_index), polecross_index]
        end
    end

when I am detecting these crossing events, I cant proceed like this because obviously roots and their order is mixed during almost each call of roots function.

So instead of smooth evolution like in the subplot below

A have these trajectories from which I cant detect desired crossing events.

Q: is there an alternative how to follow concrete trajectory of pole/zero?

R = 1;
L = 2;
C = 2;
gm1 = 0;
gm2 = 0;
gm3 = 0;
syms s
gm1vec =0:0.1:2.5;
gm2vec = 0:0.5:2;
Z1 = nan(length(gm2vec),length(gm1vec))
real_parts_of_poles = [];
imag_parts_of_poles = [];
for m = 1:length(gm2vec)
    for n = 1:length(gm1vec)
        gm2 = gm2vec(m);
        gm3 = gm1vec(n);
        Y = [1/R + 2/(s*L), -2/(s*L), 0,0,0,0,0,0;
            -2/(s*L), 2/(s*L) + s*C + 1/(s*L), -1/(s*L), 0,0,0,gm3,0;
            0, -1/(s*L), s*C + 2/(s*L), -1/(s*L), 0,gm2,0,0;
            0,0, -1/(s*L), s*C + 2/(s*L), -1/(s*L) + gm1, 0,0,0;
            0,0,0, -1/(s*L), s*C + 2/(s*L), -1/(s*L), 0,0;
            0,0,0,0,-1/(s*L),s*C + 2/(s*L), -1/(s*L), 0;
            0,0,0,0,0, -1/(s*L),s*C + 2/(s*L), -1/(s*L);
            0,0,0,0,0,0, -2/(s*L), 1/R + 2/(s*L);
            ];
        
        dY = det(inv(Y)); 
        
         [N,D] = numden(dY);
        % resultant(N,D,s)
        
        sys = tf(sym2poly(N),sym2poly(D));
        P = pole(sys)
        
        %%% SURFACE PLOTS 13.9.2021
        RR = imag(roots(sym2poly(D)));
        rr = real(roots(sym2poly(D)));
        real_parts_of_poles = round([real_parts_of_poles rr],2);
        imag_parts_of_poles = round([imag_parts_of_poles RR],2);
        
        OO = round(abs(RR(rr>=0)),2);
        
        if isempty(OO)
            Z1(m,n) = nan;
        else %&& (length(OO)<2)
            Z1(m,n) = min(OO(1));
%             counter_first_crossing_first_surface = 1;           
            
            %             if ~counter_first_crossing_second_surface && (length(OO)==2)
            %             Z2(m,n) = OO(2);
            %             counter_first_crossing_second_surface = 1
            %             end
            %
            %             if ~counter_first_crossing_third_surface && (length(OO)==3)
            %             Z3(m,n) = OO(3);
            %             counter_first_crossing_third_surface = 1
            %             end
        end
%         
        
        
%         %
        figure(1)
        subplot(2,1,1)
        plot(real(roots(sym2poly(D))),imag(roots(sym2poly(D))),'x','color',ps{m})
        hold on
        plot(real(roots(sym2poly(N))),imag(roots(sym2poly(N))),'o','color',ps{m})
        
        grid on
        xlabel('$\Re (s)$','interpreter','latex')
        ylabel('$\Im (s)$','interpreter','latex')
        % pause(0.002)
        hAxes = gca;
        hAxes.TickLabelInterpreter = 'latex';
        ax = gca;                        % gets the current axes
        ax.XAxisLocation = 'origin';     % sets them to zero
        ax.YAxisLocation = 'origin';     % sets them to zero
        ax.Box = 'off';                  % switches off the surrounding box
        axis equal
%         
%         
%         
        figure(1)
        subplot(2,1,2)
        plot(gm1vec(n),real(roots(sym2poly(D))),'x','color',ps{m})
        hold on
        plot(gm1vec(n),real(roots(sym2poly(N))),'o','color',ps{m})
        ylabel('$\Re \mathrm{det}\,\mathbf{Z}(s)_{\Omega_0/\Omega_\infty}$','interpreter','latex')
        xlabel('$g_{\mathrm{m}3}$','interpreter','latex')
        grid on
        ylim([-0.2, 0.1])
        
    end
    %%
    %real_parts_of_poles;
    %imag_parts_of_poles;
    
    %figure
    %plot(gm1vec,real_parts_of_poles(6,:))
    
    for pp = 1:size(real_parts_of_poles,1)
        polecross_index = zci(real_parts_of_poles(pp,:));
        
        if ~isempty(polecross_index)
        omega = imag_parts_of_poles(pp,polecross_index);
        MM = [omega, gm2, gm1vec(polecross_index), polecross_index]
        end
    end
    
    
end

0 个评论
显示 -2更早的评论隐藏 -2更早的评论

请先登录，再进行评论。

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

Avoid root switching using roots

0 个评论
显示 -2更早的评论隐藏 -2更早的评论

回答（0 个）

另请参阅

类别

标签

产品

版本

Community Treasure Hunt

Avoid root switching using roots

0 个评论 显示 -2更早的评论隐藏 -2更早的评论

回答（0 个）

另请参阅

类别

标签

产品

版本

Community Treasure Hunt

0 个评论
显示 -2更早的评论隐藏 -2更早的评论