计算机仿真技术--第五部分.ppt

第4章 控制系统的分析与设计 4.1 控制系统的相关性能 4.1.1 稳定性分析 4.1.2 系统的相似变换 4.1.3 能控性与能观性分析 4.2 控制系统的频域分析 4.3 控制系统的时域分析 4.4 控制系统的根轨迹分析 Bode图 Nyquist 曲线 Nichols图 DateDate 1 1 4.1.1 系统的稳定性 4.1.1.1 系统的稳定性描述 从极点或特征根分布角度 设线性控制系统的传递函数为 DateDate 2 2 q对于连续时间系统，如果闭环极点全部在S 平面左半平面，则系统是稳定的。 q对于离散时间系统，如果系统全部极点都位 于Z平面的单位圆内，则系统是稳定的。 特征方程式的根就是系统传递函数的极点。 系统的特征方程为 DateDate 3 3 系统稳定：闭环系统的 极点全部分布在s平面的 左半平面； 系统不稳定：至少有一 个极点分布在s平面的右 半平面； 系统临界稳定：在s平 面上的右半平面无极点 ，至少有一个极点在虚 轴上。 DateDate 4 4 系统状态方程的稳定性 u连续线性状态方程 u解析解形式 u稳定性：A矩阵的特征根均有负实部 DateDate 5 5 离散系统的稳定性 u离散系统状态方程 u离散系统时域响应解析解 u稳定性判定：所有特征根均在单位圆内 DateDate 6 6 从能否自动恢复平衡角度 q稳定指控制系统在外作用消失后自动恢复原 有平衡状态或自动地趋向于一个新的稳定平衡 状态的能力。 q如果系统不能恢复稳定状态，则认为系统不 稳定。 DateDate 7 7 单摆系统稳定倒摆系统不稳定 各点是否稳定？ DateDate 8 8 基于特征多项式系数的解析方法 4.1.1.2 系统稳定性判别方法 间接判别法 （代数稳定判据） 霍尔维茨判据劳斯判据 DateDate 9 9 劳斯判据 判断特征方程的所有系数 是否均大于0或均小于0 劳斯(Routh)判据：系统稳定的充分必要条件是 ，劳斯表中第一列各值严格为正；如果劳斯表第 一列中出现小于零的数值，系统不稳定。 首先判断！ 然后考虑！ DateDate1010 构造劳斯阵列 由系统矩阵的特征方程 sn sn-1 sn-2 sn-3 anan-2an-4an-6 an-1an-3an-5an-7 b1b2b3b4 c1c2c3c4 其中 DateDate1111 DateDate1212 Edward John Routh (1831-1907) DateDate1313 Routh 判据的历史局限性 uRouth判据提出时，没有求多项式根的方法 u现在求解矩阵特征根、求解多项式方程的 根轻而易举，无需间接方法 uRouth判据只能得出是否稳定，进一步信息 得不出来，如系统是否振荡 u离散系统无法由Routh方法直接判定，得借 助于Jury判据，更复杂 u稳定性分析方法不统一 DateDate1414 考虑特征方程 构造霍尔维茨行列式： 霍尔维茨判据 行列式中，对角线上各元 为特征方程中自第二项开 始的各项系数。每行以对 角线上各元为准，写对角 线左方各元时，系数ai的 脚标递增；写对角线右方 各元时，系数ai的脚标递 减。当写到在特征方程中 不存在系数时，则以零来 代替。 DateDate1515 霍尔维茨(Hurwitz)判据：当且仅当由特征方程系 数构成的霍尔维茨矩阵为正定矩阵时，系统稳定。 即行列式D的各阶主子式( Di )均大于零 DateDate1616 由Lyapunov稳定性定义可知，如果能够考察系 统的全部零极点的分布情况，就可以判断系统 的稳定性。 MATLAB提供了直接求取系统所有零极点 的函数，因此可以直接根据零极点的分布情况 对系统的稳定性进行判断。 直接判别法 （特征方程求根） DateDate1717 基于 MATLAB 的稳定性判定方法 u直接判定 状态方程模型 u由 可以求出所有特征根 u离散系统： u图解判定法 连续系统： 绘制零极点图 离散系统： ，同时画出单位圆 DateDate1818 例 已知系统模型， 要求判断系统的稳 定性 程序：clear clc a=1 2 -1 2;2 6 3 0;4 7 -8 -5;7 2 1 6; b=-1 0 0 1; c=-2 5 6 1;d=7; z,p,k=ss2zp(a,b,c,d) DateDate1919 %检验极点的实部；求取极点实部大于零的个数 ii=find(real(p)0) n2=length(ii); %判断系统是否稳定 if(n20) disp(the system is unstable) disp(the unstable pole are:) disp(p(ii) else disp(the system is stable) end %绘制零极点图 pzmap(p,z) DateDate2020 ii=find(条件式) 用来求取满足条件的向量的下标向量，以列 向量表示。 条件式为real(p0)，其含义就是找出极点向量p 中满足实部的值大于0的所有元素下标，并将结 果返回到ii向量中： 1. ii的元素个数 0，则认为找到了不稳定极点 ，因而给出系统不稳定的提示 2. ii向量的元素个数 = 0，则认为没有找到不稳 定的极点，因而得出系统稳定的结论。 DateDate2121 4.1.2 线性系统的线性相似变换 u系统的状态方程表示称为系统实现 u不同状态选择下，状态方程不唯一 u相似变换 非奇异矩阵 状态变换 新状态方程模型 DateDate2222 u状态变换公式 uMATLAB 求解方法 DateDate2323 例 已知系统和转换矩阵 uMATLAB 求解 DateDate2424 u变换结果 u可见，相似变换能改变系统的结构 u引入相似变换矩阵，可以将已知系统转换 成其他的形式 DateDate2525 4.1.3 系统的可控性与可观性 u可控性定义 系统的可控性就是指系统内部的状态是否 可以由外部输入信号控制的性质！ DateDate2626 u可控性判定 1.可控性判定矩阵 u 判定依据：若矩阵Tc满秩，则系统完全可控。 u构造可控性判定矩阵 u基于 MATLAB 的判定方法 DateDate2727 例 离散状态方程的可控性 uMATLAB 求解 DateDate2828 u判定矩阵 u判定矩阵构造方法 u这样的判定方法同样适合于连续系统和离 散系统，也适用于多变量模型。 DateDate2929 2.由 Gram 矩阵判定可控性 u引入可控 Gram 矩阵 u该矩阵满足 Lyapunov 方程 uMATLAB 求解 u构造Gram矩阵 ，Gc非奇 异，则系统可控！ DateDate3030 u可观测性定义 系统的可观测性就是指系统内部的状态 是不是可以由系统输出信号重建起来的 性质！ DateDate3131 可观性判定 u判定矩阵 是否满秩？ To=obsv(A,C) rank(To) uGram 矩阵是否非奇异 Go=gram(G,o) rank(Go) DateDate3232 注意： 无论用哪种方法判别可控 性与可观性，必须判断相应矩 阵的秩！ DateDate3333 4.2 控制系统的频域分析 一、频域分析的一般方法 q频域分析法是应用频率特性研究控制系统的 一种典型方法。 1. 可直观地表达出系统的频率特性 2. 分析方法比较简单 3. 物理概念比较明确 对于诸如防止结构谐振、抑制噪声、改 善系统稳定性和暂态性能等问题，都可以从 系统的频率特性上明确地看出其物理实质和 解决途径。 DateDate3434 频率响应或频率特性: 将正弦信号加在一个线 性环节(系统)时，在稳态下环节的输出量与输 入量之比。此比值给出了不同频率下，系统传 输正弦信号的性能。这个比值称为系统的频率 响应，或频率特性。通常记为： DateDate3535 分析频率特性方法 幅相频率特性曲线 (幅相曲线,极坐标图) 对数频率特性曲线 Bode图：bode()nyquist图：nyquist() DateDate3636 q对数频率特性图包括了对数幅频特性图和 对数相频特性图。 q采用对数频率特性表示法，使得研究频率 范围很宽的频率特性时，缩小了比例尺；大 大简化绘制系统频率特性的工作，如若干环 节串联时，只需绘出各环节的对数幅频特性 ，进行加减后即可得到串联后系统的对数频 率特性。 对数频率特性曲线（Bode图） DateDate3737 qbode(num,den)：可绘制出以连续时间传 递函数表示的系统的Bode图。 qbode(a,b,c,d): 针对连续状态空间系统 a,b,c,d的每个输入的Bode图。其中频率范 围由函数自动选取，而且在响应快速变化的 位置会自动采用更多取样点。 qbode(a,b,c,d,iu)：可得到从系统第iu个输 入到所有输出的Bode图。 DateDate3838 qbode(a,b,c,d,iu,w)或bode(num,den,w)：可利 用指定的角频率矢量绘制出系统的波德图。 q当带输出变量mag,pha,w或mag,pha引 用函数时，可得到系统波德图相应的幅值 mag、相角pha及角频率w矢量或只是返回幅 值与相角。相角以度为单位，幅值可转换为 分贝单位：magdb=20log10(mag) DateDate3939 wn=6; w=0:0.1:10; num=wn2; for kosi= 0.1:0.1:10 den=1 2*kosi*wn wn2; mag,pha,w1=bode(num,den,w); % 注意mag的单位不是分贝，若需要分贝表示， 可以通过20*log10(mag)进行转换 例：求典型二阶系统自然振荡频率固定，阻尼 比变化时的Bode图 DateDate4040 subplot(121); hold on; plot(w,mag) grid on % 绘制幅值的变化规律 subplot(122); hold on; plot(w,pha) % 绘制相角的变化规律 end DateDate4141 DateDate4242 例： close all a=-2.5 -1.22 0 0;1.22 0 0 0;1 -1.14 -3.2 -2.56;. 0 0 2.56 0; b=4 1;2 0;2 0;0 0; c=0 1 0 3;0 0 0 1; d=0 -2;-2 0; figure(1) bode(a,b,c,d) figure(2) bode(a,b,c,d,1) % 绘制第一个输入到所有输出的Bode图 DateDate4343 DateDate4444 DateDate4545 将频率特性函数G(jw)表示为代数形式 幅相频率特性曲线(幅相曲线,极坐标图) 其中P(w)为频率特性的实部，称为实频特性， Q(w)为频率特性的虚部，称为虚频特性。当频 率w由0变到正无穷时，G(jw)的矢量终端将绘 出一条曲线，此曲线称为系统的幅相频率特性 ，或称为奈氏图。 DateDate4646 qnyquist(a,b,c,d)：绘制出系统的一组nyquist 曲线，每条曲线相应于连续状态空间系统 a,b,c,d的输入/输出组合对。其中频率范围由 函数自动选取，而且在响应快速变化的位置会 自动采用更多取样点。 qnyquist(a,b,c,d,iu)：可得到从系统第iu个输 入到所有输出的极坐标图。 qnyquist(num,den)：可绘制出以连续时间多 项式传递函数表示的系统的极坐标图。 DateDate4747 qnyquist(a,b,c,d,iu,w)或nyquist(num,den,w)： 可利用指定的角频率矢量绘制出系统的极坐标 图。 q当不带返回参数时，直接在屏幕上绘制出 系统的极坐标图（图上用箭头表示w的变化 方向，负无穷到正无穷）。当带输出变量 re,im,w引用函数时，可得到系统频率特性 函数的实部re和虚部im及角频率w矢量（为 正的部分）。可以用plot(re,im)绘制出对应w 从负无穷到零变化的部分。 DateDate4848 例： 已知系统的传递函数为： G(s)=K/(s3+52s2+100s)，求当K分别取 1300和5200时，系统的极坐标频率特性图 k1=1300; k2=5200; w=8:1:80; num1=k1; num2=k2; den=1 52 100 0; figure(1) subplot(211) nyquist(num1,den,w); DateDate4949 subplot(212) pzmap(num1,den); %绘制零极点图 figure(2) subplot(211) nyquist(num2,den,w); subplot(212) re,im=nyquist(num2,den); plot(re,im) xlabel(real) ylabel(image) title(w from 负无穷 to 零) DateDate5050 DateDate5151 DateDate5252 例：close all a=-2.5 -1.22 0 0;1.22 0 0 0;1 -1.14 -3.2 -2.56;. 0 0 2.56 0; b=4 1;2 0;2 0;0 0; c=0 1 0 3;0 0 0 1; d=0 -2;-2 0; figure(1) nyquist(a,b,c,d) figure(2) nyquist(a,b,c,d,2) % 绘制第二个输入到所有输出的波特图 % 注意在图中标示依然使用u(1)表示第二个输入 DateDate5353 DateDate5454 DateDate5555 二、常用频域分析函数 margin：求幅值裕度和相角裕度及对应的 转折频率 freqs： 模拟滤波器特性 nichols：求连续系统尼科尔斯频率响应曲线 （即对数幅相曲线） ngrid： 尼科尔斯方格图 DateDate5656 margin函数通常用在bode函数之后，先由 bode函数得到幅值、相角和频率矢量，然后由 margin绘制出幅值裕度和相角裕度的Bode图。 幅值裕度和相角裕度是针对开环SISO系统 而言，它指示出系统闭环时的相对稳定性。当 不带输出变量引用时，margin可在当前图形窗 口中绘制出带有裕量及相应频率显示的Bode图 ，其中幅值裕度以分贝为单位。 margin：求幅值裕度和相角裕度及对应的 转折频率 DateDate5757 q幅值裕度是在相角为-180度处使开环增 益为1的增益量，如在-180度相频处的开环 增益为g，则幅值裕度为1/g；若用分贝值 表示幅值裕度，则等于：-20*log10(g)。类 似地，相角裕度是当开环增益为1.0时，相 应的相角与180度角的和。 qmargin(num,den) ：可计算出连续系统 传递函数表示的幅值裕度和相角裕度并绘 制相应Bode图。类似，margin(a,b,c,d)可 以计算出连续状态空间系统表示的幅值裕 度和相角裕度并绘制相应的Bode图。 DateDate5858 qGm,Pm,Wcg,Wcp=margin(mag,phase,w )：由幅值mag（不是以dB为单位） 、相角 phase及角频率w矢量计算出系统幅值裕度 和相角裕度及相应的相角交界频率wcg、截 止频率wcp，而不直接绘出Bode图曲线。 DateDate5959 例： w=logspace(-1,2,39); mag,pha=bode(1,1 2 3); Gm,Pm,Wcg,Wcp=margin(mag,pha,w) % Gm幅值裕度 % Pm相角裕度 % Wcg-幅值裕度对应的角频率 % Wcp-相角裕度对应的角频率 DateDate6060 qfreqs用于计算由矢量a和b构成的模拟滤 波器H(s)=B(s)/A(s)的幅频响应。 图形上和Bode 图基本一致！ freqs：模拟滤波器特性 qh=freqs(b,a,w)用于计算模拟滤波器的幅频 响应，其中实矢量w用于指定频率值，返回 值h为一个复数行向量，要得到幅值必须对 它取绝对值，即求模。 DateDate6161 qh,w=freqs(b,a)自动设定200个频率点来 计算频率响应，这200个频率值记录在w中 。 qh,w=freqs(b,a,n)设定n个频率点计算频率 响应。 不带输出变量的freqs函数，将在当前图 形窗口中绘制出幅频和相频曲线，其中幅相 曲线对纵坐标与横坐标均为对数分度。 DateDate6262 例：有一模拟滤波器，其传递函数为： H(s)=(0.2s2+0.3s+1)/(s2+0.4s+1) 求它的幅频特性和相频特性 b=0.2 0.3 1; a=1 0.4 1; w=logspace(-1,1); h=freqs(b,a,w); %计算幅频响应 mag=abs(h); %对h求模值 semilogx(w,mag) grid on figure(2) freqs(b,a,w) %直接绘制幅频与相频曲线DateDate6363 DateDate6464 幅相曲线对纵坐标与横坐 标均为对数分度! DateDate6565 figure(3) bode(b,a,w) g