任何系统的稳定性是系统本身的一种属性，仅仅取决于系统的结构参数，与初始条件和外作用无关。你可以设想一下，K只是增益而已，只是在一个稳定的系统上乘上放大倍数。二阶系统只要不是无阻尼系统和负阻尼系统，系统最终都会趋向于稳定。

记录一下简单的二阶理论：
1.标准型：分母是s^2+2ξωs+ω^2。注意分母第一项系数是1，对于一个任意电路的传函一定要化成1再讨论。
2.ξ是阻尼系数，ω（一般都有个下标n）是自然频率。
3.根轨迹：|ξ|<1对应圆，圆心为s-plane原点，半径为ω。其中ξ=0对应虚轴上两点（正负jωn），ξ=1对应左半轴-ωn。
4.对于稳定系统，要求ξ>0，此时根轨迹全部落在左半平面。
5.在4的条件下，ξ<1对应欠阻尼，根对称落在实轴两侧，电路振荡衰减；ξ=1临界阻尼，ξ>1对应两个左半轴实极点，
一个在-ωn左，一个在右。
6.根据根轨迹，如果分子为常数（不包含零点），可以对不同的ξ画出极零图，获得频响信息：
    1>欠阻尼情况，此时频响出现尖峰（频率靠近振荡频率处）且尖峰随ξ减小而逐渐明显（越“冲”）。
    2>临界阻尼，此时频响是完美的低通曲线
    3>过阻尼，也是低通下降的曲线
7.实际的物理意义：ξ是电路中的耗能元件，一般为电阻。
实际电路中也可能和其他因素有关（譬如锁相环中ξ反比于VCO和PD增益）。

传输矩阵法是一种常见的计算方法，例如：对于透镜系统，利用传输矩阵可以很容易计算高斯光的扩束与聚焦；激光腔的稳定性条件要利用传输矩阵计算；计算一维光子晶体（光纤光栅）的透射率。
当我用MATLAB写完for循环之后，觉得用这种老老实实一个一个算实在太慢了，想到应该尽可能向量化去解决这个问题，最后我想到了一种方法来提升传输矩阵的计算速度：
思路
我们知道对于矩阵连乘：
在其中添加括号是完全没问题的，上式可以改写为：
那么可以将一个很长很长的连乘，拼成一个方阵，形如：
MATLAB用每一列的矩阵乘以下一列的矩阵，得到每一行的结果：
最后竖着连乘，即可得到
假设原本是 
 个矩阵连乘，如果老老实实算需要计算 
 次矩阵相乘。但是用这种方法，只需要计算乘法 
 次。当然也可以连续用这种方法，计算乘法的次数可能更少一些。
对于 
 不是某个数平方的时候，只需在最后一行补足单位矩阵 
 即可。
2. 实例
随便算一个10000个2阶矩阵连乘，如果用循环一个一个做：
用之前提到改进的做法：
快了10倍有木有……而且N越大效果越明显。实际上算光纤光栅，原本半小时的程序可以缩短到几十秒，还是蛮有用的。
3. 代码
%快速二阶矩阵连乘函数
%知乎 Floron
%--------------------------------------------------------------------------
function out=SecOrderMat_Multiplication(M_11,M_12,M_21,M_22)
mat_size=ceil((length(M_11))^0.5);
Length_full=mat_size^2;
M_11=[M_11 ones(1,Length_full-length(M_11))];
M_12=[M_12 zeros(1,Length_full-length(M_12))];
M_21=[M_21 zeros(1,Length_full-length(M_21))];
M_22=[M_22 ones(1,Length_full-length(M_22))];
M_11=reshape(M_11,mat_size,mat_size).';
M_12=reshape(M_12,mat_size,mat_size).';
M_21=reshape(M_21,mat_size,mat_size).';
M_22=reshape(M_22,mat_size,mat_size).';
row_M_11_0=1;
row_M_12_0=0;
row_M_21_0=0;
row_M_22_0=1;
for lab=1:mat_size
    row_M_11=row_M_11_0.*M_11(:,lab)+row_M_12_0.*M_21(:,lab);
    row_M_12=row_M_11_0.*M_12(:,lab)+row_M_12_0.*M_22(:,lab);
    row_M_21=row_M_21_0.*M_11(:,lab)+row_M_22_0.*M_21(:,lab);
    row_M_22=row_M_21_0.*M_12(:,lab)+row_M_22_0.*M_22(:,lab);
    row_M_11_0=row_M_11;
    row_M_12_0=row_M_12;
    row_M_21_0=row_M_21;
    row_M_22_0=row_M_22;
end
M_11=row_M_11';
M_12=row_M_12';
M_21=row_M_21';
M_22=row_M_22';
row_M_11_0=1;
row_M_12_0=0;
row_M_21_0=0;
row_M_22_0=1;
for lab=1:mat_size
    row_M_11=row_M_11_0.*M_11(lab)+row_M_12_0.*M_21(lab);
    row_M_12=row_M_11_0.*M_12(lab)+row_M_12_0.*M_22(lab);
    row_M_21=row_M_21_0.*M_11(lab)+row_M_22_0.*M_21(lab);
    row_M_22=row_M_21_0.*M_12(lab)+row_M_22_0.*M_22(lab);
    row_M_11_0=row_M_11;
    row_M_12_0=row_M_12;
    row_M_21_0=row_M_21;
    row_M_22_0=row_M_22;
end
out=[row_M_11,row_M_12;row_M_21,row_M_22];
end
%--------------------------------------------------------------------------
%                               介绍
%   当有许多二阶矩阵连乘：M_1*M_2*M_3*.....*M_n。（例如传输矩阵）
%   如果采用传统的for循环，计算次数为n次。
%   这里采用向量化的办法，对于连续相乘的矩阵，在不改变顺序的前提下可以随意添加
%括号，将其改写为一个sqrt(n)*sqrt(n)的方阵，沿着行增加的方向，所有列同步连乘，
%最终得到一列向量后再沿着列增加方向连乘。最终计算次数为2*sqrt(n)次。极大提升了
%计算速度，尤其是当n足够大时。
%--------------------------------------------------------------------------
%                               用法
%   将每个单独的矩阵（1,1）位置的元素作为一个行向量M_11输入。
%   将每个单独的矩阵（1,2）位置的元素作为一个行向量M_12输入。
%   将每个单独的矩阵（2,1）位置的元素作为一个行向量M_21输入。
%   将每个单独的矩阵（2,2）位置的元素作为一个行向量M_22输入。
%--------------------------------------------------------------------------
%                               示范
%   求[1,2;3,4]*[10,20;30,40]*[100,200;300,400]*[1000,2000;3000,4000]
%   输入分别为：
%   M_11=[1,10,100,1000]
%   M_12=[2,20,200,2000]
%   M_21=[3,30,300,3000]
%   M_22=[4,40,400,4000]
%--------------------------------------------------------------------------
%           代码示范（复制同一个文件夹的另一个m文件中使用）：
% % M_11=[1,2,3,4,5,6];
% % M_12=[2,3,4,5,6,7];
% % M_21=[3,4,5,6,7,8];
% % M_22=[4,5,6,7,8,9];
% % %use function
% % A_myfun=SecOrderMat_Multiplication(M_11,M_12,M_21,M_22);
% % %use for-loop
% % A=[1 0;0 1];
% % for lab=0:5
% %     A=A*[1+lab 2+lab;3+lab 4+lab];
% % end
%--------------------------------------------------------------------------
 如果有疑问可以留言，如果觉得好用请给我口米吃。

系列文章目录
自动控制理论（1）——自动控制理论概述

自动控制理论（2）——控制系统的数学模型（微分方程、传递函数）

自动控制理论（3）——控制系统的数学模型（系统框图和信号流图）

自动控制理论（4）——系统的时域性能指标和一阶系统的时域分析

自动控制理论（5）——二阶系统的时域分析

文章目录
系列文章目录
一、高阶系统的时域分析
1.典型三阶系统的单位阶跃响应
2.高阶系统的单位阶跃响应
3.闭环主导极点
二、线性系统的稳定性分析
1.稳定的概念
2.线性系统稳定的充分必要条件
3.劳斯判据
4.稳定判据的应用（检验稳定裕量）
5.总结



一、高阶系统的时域分析
1.典型三阶系统的单位阶跃响应


$\ β=\frac{-s_3}{ζω_n}$

增加一个闭环极点，将使超调量减小，上升时间和峰值时间增大

增加的极点距离虚轴越近，上述影响越显著

$\ β<1$呈现过阻尼，响应迟缓

$\ β >5$,可忽略极点的影响
2.高阶系统的单位阶跃响应
c(t)由稳态和暂态分量组成，若极点均为负实部，则系统稳定。

各暂态分量衰减的快慢,取决于各极点负实部的绝对值大小。

各暂态分量系数的大小是F(s)零、极点共同决定的，若一对零、极点几乎重合(称偶极子),则与该极点对应的系数很小,该极点对暂态响应几乎无影响。
3.闭环主导极点
对系统的暂态响应起主导作用的极点。
满足以下两个条件：

(1)距虚轴比较近,且附近没有其它的闭环零点与极点。

(2)其实部的绝对值应比其它极点的实部绝对值小五倍以上。
二、线性系统的稳定性分析
1.稳定的概念
稳定性是指扰动消失后，系统由初始偏差状态恢复到原平衡状态的性能。若系统能恢复平衡状态，则称系统是稳定的，否则不稳定。

稳定的线性系统，大范围小范围都能稳定；

非线性系统可能小范围稳定而大范围不稳定。

2.线性系统稳定的充分必要条件
闭环特征根均为负实部
3.劳斯判据
(1)设系统特征方程$\ a_0s^n + a_1s^{n-1} + … + a_{n-1}s+ a_n = 0,     a_0>0$

(2)系统稳定的必要条件

特征方程各项系数均为正

(3)列写劳斯表



(4)系统稳定的充要条件是劳斯表中第一列元素均为正。

若第一列元素有符号变化，则符号改变的次数等于正实部根的个数。

(5)劳斯判据特殊情况

1）劳斯表某行的第一列为零，其余各项不为零

用$\ ε$（$\ ε>0$且$\ ε→0$）代替0继续计算

2）劳斯表中出现全零行

用全零行上一行的系数构成辅助方程；

将辅助方程对s求导数，得一新方程；

用新方程的系数代替全零行，按新表判稳定。

结论：

不稳定;

若第一列元素均为正,没有右根，一定有纯虚根;

若第一列元素有负数，符号改变次数等于右根个数。
4.稳定判据的应用（检验稳定裕量）
稳定裕量σ的概念：虚轴向左移动σ，系统依然稳定。

令 s = z -σ(σ>0), 将其代入特征方程，

可得关于z 的多项式，对z用劳斯判据
5.总结
线性系统的稳定性只取决于系统的结构及

参数，而与初始条件、外作用大小及形式

无关。

稳定性只取决于系统闭环极点，而与系统

零点无关。