matlab用博德稳定性判据,自编劳斯判据判断系统稳定性的Matlab函数

function [routh_list,conclusion] = Routh(chara_equ)

% =======================================================

% 自编劳斯判据求解系统稳定性函数

% 输入：

% chara_equ = 特征方程向量

% 输出：

% routh_list = 劳斯表

% conclusion = 给出系统是否稳定或存在多少个不稳定的根的结论

% example:

% [routh_list,con] = Routh([1 2 3 4 5]);

% return:

% routh_list =

%

%      1     3     5

%      2     4     0

%      1     5     0

%     -6     0     0

%      5     0     0

% con =

%

% There is 2 unstable roots!

% =========================================================

n=length(chara_equ);

chara_equ=reshape(chara_equ,1,n);

if mod(n,2)==0

n1=n/2;

else

n1=(n+1)/2;

chara_equ=[chara_equ,0];

end

routh=reshape(chara_equ,2,n1);

routh_list=zeros(n,n1);

routh_list(1:2,:)=routh;

i=3;

while 1;

%  =========特殊情况1(第一列为0，其余列不为0)=====================

if routh_list(i-1,1)==0 & sum(routh_list(i-1,2:n1))~=0

chara_equ = conv(chara_equ,[1 3]);

n=length(chara_equ);

if mod(n,2)==0

n1=n/2;

else

n1=(n+1)/2;

chara_equ=[chara_equ,0];

end

routh=reshape(chara_equ,2,n1);

routh_list=zeros(n,n1);

routh_list(1:2,:)=routh;

i=3;

end

% ==========计算劳斯表===========================================

ai=routh_list(i-2,1)/routh_list(i-1,1);

for j=1:n1-1

routh_list(i,j)=routh_list(i-2,j+1)-ai*routh_list(i-1,j+1);

end

% ==========特殊情况2(全0行)======================================

if sum(routh_list(i,:))==0

k=0;

l=1;

F=zeros(1,n1);

while n-i-k>=0

F(l)=n-i+1-k;

k=k+2;

l=l+1;

end

routh_list(i,:)=routh_list(i-1,:).*F(1,:);

end

% =========更新==================================================

i=i+1;

if i>n

break;

end

end

% =============outhput===========

r=find(routh_list(:,1)<0);

if isempty(r)==1

conclusion='The system is stable!';

else

n2=length(r);

m=n2;

for i=1:n2-1

if r(i+1)-r(i)==1

m=m-1;

end

end

str1='There is ';

if r(n2)==n

str2=num2str(m*2-1);

else

str2=num2str(m*2);

end

str3=' unstable roots!';

conclusion = [str1,str2,str3];

end

有点错误     很急  高手帮个忙吧     谢谢   真心感谢