毕业论文-连续系统simulink状态空间建模分析方法程序设计.doc

陕西理工学院毕业设计 毕业设计题 目 连续系统simulink状态空间建模 分析方法程序设计 学 连续系统simulink状态空间建模分析方法程序设计摘要 本课题基于对信号与信息处理课程中用matlab/simulink建模及应用分析滤波器问题的深入研究。通过自身掌握的理论知识，主要以高阶连续系统（模拟滤波器）为例，并将其离散化，转化为离散系统，从而对离散系统处理。用simulink状态空间函数模块建模，观察并分析波形。其次，用matlab中的M文件编程，求解系统，绘制波形并进行频谱分析。在本课题中，主要将连续系统转化为离散系统，再用计算机和matlab软件进行研究，用simulink对高阶离散系统建模，并设置模块参数，自定义函数为正弦波的叠加，传输函数的相关参数后运行并进行频谱分析，使信号的性态都能得到处理和研究。通过编程，求解高阶离散系统的零输入响应，零状态响应和完全响应，求解实际生活中的各种问题，改变参数并对信号进行适当的频谱分析。关键字 连续系统;离散化;simulink ;M文件 Continuous system simulink state spacemodeling analysis program design Yang lijuan(Grade10,Class2,Major Electronic Information Science and TechnologyDepartment of Physics,Shannxi University of Technology,Hanzhong,723001) Tutor: Long shumingAbstract: In-depth study of the signal and information processing program，using matlab / simulink modeling and application analysis on this topic filter problem. Through my own master theoretical knowledge, mainly in the high-end continuous system (analog filters), for example, and discrete, into a discrete system, and thus the discrete system processing. Simulink modeling with spatial function module status, to observe and analyze waveforms. Secondly, using the M-file matlab programming, solving system, drawing the waveform and spectral analysis. In this paper, the main continuous system into a discrete system, then the computer and matlab software research for high-end discrete system with simulink modeling, and set the module parameters, a custom function as the sine wave, the transfer function after running parameters and spectral analysis of the behavior of the signal can be processed and research. By programming, solving high-end discrete system response zero input, zero state response and complete response, solving real-life problems, changing the parameters and signals for proper spectral analysis. Keywords: Continuous system ;Discrete; Simulink; M file目录引言 11 状态空间分析方法的概述12 快速创建LTI连续系统状态空间模型的方法22.1 创建LTI连续系统传递函数的方法22.2 构造描述LTI连续系统的状态空间模型矩阵33 利用simulink状态空间建模求解LTI系统数值解的思路43.1 MATLAB编程设计并描述低通数字滤波器43.2 创建系统的simulink状态空间模型63.3 模块内部参数设置及数据存储74 利用simulink状态空间建模求解LTI系统的优缺点85 连续系统simulink状态空间建模分析方法程序设计的思路85.1 调用模型文件及编程求解系统响应85.2 分析系统的频谱与相位106 状态空间分析方法的应用实例116.1 分析求解低阶电路系统116.2 设计分析滤波器系统127 结束语12致谢 12参考文献13引言随着科学技术的发展，系统的组成也日益复杂。连续系统Simulink状态空间建模分析方法连续系统的各种解析解法虽然便于理论分析系统响应的变化趋势和系统特性，但实际系统总是多输入多输出的高阶系统，它们的解析形式的响应求解极为困难、即便较低阶系统的解析响应能够得到，其函数表示也比较复杂。而连续系统的区间数值解法本质上用的是迭代解法，总是能够方便、快速地的得到，之后如果企图观察其响应随时间演化的趋势，可用数值解画出其波形来观察，甚至必要时做数据拟合寻找区间解的拟合函数也是人可能的，而且数值解法还可以求一定区间上的非线性问题。在许多情况下，人们不禁关心系统输出的变化情况，而且还要研究与系统内部一些变量有关的问题，比如，系统的客观性和可控性行。系统的最优控制与设计等问题。为适应这一变化，引入了状态变量法，也称内部法。在科学研究和实际生活中，离散时间信号具有高精度、可靠性好、便于集成等优点，渐渐在很多领域得以运用,，而离散信号与系统的分析在通信与信息系统，信号处理，自动控制，检测，监控等领域也都有十分重要的作用。人员伤亡。因此，系统的性能必须在设计、制造和使用环节应用先进的可靠性工程技术来提高，在多状态、多阶段理论框架下开展对其可靠性理论与方法的研究就显得尤为重要1。所以可以将本课题的连续系统离散化，为离散系统，方便求解。本课题基于对信号与信息处理课程中用matlab/simulink建模及应用分析实际问题的深入研究。由于simulink是可以对动态系统进行建模，仿真及分析的软件包，依托于matlab的丰富的仿真资源，可以使用其快速对系统进行仿真，准确分析所描述动态系统。除此之外，还可用其中的M文件编程以改变系统中参数，更加方便快捷地研究参数对系统的影响，包括求解系统，绘制波形并进行频谱分析等等。1 状态空间分析方法的概述状态空间法是一种基于解答空间的问题表示和求解方法，它是以状态和操作符为基础的。在利用状态空间图表示时，从某个初始状态开始，每次加一个操作符，递增地建立起操作符的试验序列，直到达到目标状态为止。其中，可由微分方程求状态空间表达式。系统的实现:根据系统的外部描述构造一个内部结构，要求既保持外部描述的输入输出关系，又要将系统的内部结构确定下来。这是一个复杂的问题，但也是一个非常重要的问题。一方面，描述系统输入输出关系的微分方程或传递函数可以用实验的方法得到，我们可以从输入输出关系描述建立状态空间描述，这是建立状态空间描述的一条途径。一般描述为: 【2 （式1.1）设状态变量选为: （式1.2）则 （式1.3）由微分方程有所以：因此，系统的状态方程为：1 （式1.4）输出方程为 由微分方程表达为矩阵形式，如下： （式1.5） （式1.6）2 快速创建LTI连续系统状态空间模型的方法2.1 创建LTI连续系统传递函数的方法对于一般的低阶电路系统，例如RLC三阶电路（电路），要得到它的系统函数，可以先由把三阶电路系统映射到S域，可以计算出其等效阻抗，再在复频域进行求解，这样比较方便简单。具体做法如下：将电路映射到S域，在频域的电路如下图所示：图2.1 三阶 RLC复频域电路利用拉普拉斯变换进行分析：电阻不变，电感，电容，由复频域电路图建立复频域代数方程: （式2.1）其等效阻抗： （式2.2）输出像函数： （式2.3）代入R,L,C的值，最后的输出象函数： （式2.4）由输出象函数得系统函数： （式2.5）在本课题里，主要以高阶连续系统离散化为例，并选取合适典型离散系统（低通数字滤波器系统），在这里，由于滤波器内部的结构不明确，所以无法用上述方法求解，但可以利用MATLAB编程求解传递函数。举例如下：若给出滤波器的技术指标如下：wp=0.2rad,p=1dB,ws=0.35rad, s=25dB。具体MATLAB编程求解传递函数如下：T=1;%T=1swp=0.2*pi/T;ws=0.35*pi/T;rp=1;rs=25;%T=1s的模拟滤波器的指标N,wc=buttord(wp,ws,rp,rs,s);%计算响应的模拟滤波器阶数N和3dB截止频率wcB,A=butter(N,wc,s);%计算该响应的模拟滤波器的系统函数程序运行结果如下：B = 0 0 0 0 0 0 0 0.1094A =1.0000 3.2762 5.3668 5.6538 4.1218 2.0795 0.6747 0.1094说明传递函数的分母的各个系数为A，分子的系数为B。分子分母都是降序排列。由此便可写出系统的传递函数H（S）。2.2 构造描述LTI连续系统的状态空间模型矩阵要由传递函数求解状态空间矩阵的方法，有两种。一种是通过迭代法演算，先将传递函数写成微分（差分）方程，再根据定义或画出流程图进行计算，这种方法只适合低阶系统，高阶系统会变的复杂，繁琐，所以应用不广泛，但不抽象，有利于理解。另一种是通过MATLAB编程直接求解，它简单，快捷，但具体算法不明确。例如上述滤波器已求得H（S）的分子，分母系数，分别为B,A,在MATLAB里可以调用函数tf2ss进行转换，具体做法如下：设Az=1.0000, -3.8255, 6.7082 , -6.8507,4.3564, -1.7136,0.3843 ,-0.0378;Bz= -0.0000, 0.0001,0.0032,0.0105, 0.0066, 0.0008,0.0000, 0;A,B,C,D=tf2ss(Bz,Az)运行结果如下：A = 3.8255 -6.7082 6.8507 -4.3564 1.7136 -0.3843 0.0378 1.0000 0 0 0 0 0 0 0 1.0000 0 0 0 0 0 0 0 1.0000 0 0 0 0 0 0 0 1.0000 0 0 0 0 0 0 0 1.0000 0 0 0 0 0 0 0 1.0000 0B = 1 0 0 0 0 0 0C = 0.0001 0.0032 0.0105 0.0066 0.0008 0 0D = 0由此，便可由传递函数得到状态空间矩阵。3 利用simulink状态空间建模求解LTI系统数值解的思路3.1 MATLAB编程设计并描述低通数字滤波器在现代通信系统中，由于信号中经常混有各种复杂成分，所以很多信号分析都是基于滤波器而进行的，而数字滤波器是通过数值运算实现滤波，具有处理精度高、稳定、灵活、不存在阻抗匹配问题，可以实现模拟滤波器无法实现的特殊滤波功能。数字滤波器根据其冲激响应函数的时域特性，可分为两种，即无限长冲激响应(IIR)数字滤波器和有限长冲激响应(FIR)数字滤波器。实现IIR滤波器的阶次较低，所用的存储单元较少，效率高，精度高，而且能够保留一些模拟滤波器的优良特性，因此应用很广。Matlab软件以矩阵运算为基础，把计算、可视化及程序设计有机融合到交互式工作环境中，并且为数字滤波的研究和应用提供了一个直观、高效、便捷的利器3。在数字信号处理中起着重要的作用并已获得广泛应用的是数字滤波器（DF，Digital Filter），根据其单位冲激响应函数的时域特性可分为两类：无限冲激响应IIR（Infinite Impulse Response）滤波器和有限冲激响应FIR（Finite Impulse Response）滤波器。与FIR滤波器相比，IIR的实现采用的是递归结构，极点须在单位圆内，在相同设计指标下，实现IIR滤波器的阶次较低，即所用的存储单元少，从而经济效率高。其差分方程为： （式3.1）系统函数为： （式3.2）设计IIR滤波器的任务就是寻求一个物理上可实现的系统函数H(z)，使其频率响应H(z)满足所希望得到的频域指标，即符合给定的通带截止频率、阻带截止频率、通带衰减系数和阻带衰减系数。IIR滤波器指标参数如下图所示。图中，p和s分别为通带边界频率和阻带边界频率；1和2分别为通带波纹和阻带波纹；允许的衰减一般用dB数表示，通带内所允许的最大衰减（dB）和阻带内允许的最小衰减（dB）分别为p和s表示： （式3.3） （式3.4）图3.1 低通滤波器的技术要求 IIR数字滤波器的设计一般有两种方法：一个是借助模拟滤波器的设计方进行。其设计步骤是，先设计模拟滤波器，再按照某种方法转换成数字滤波器。另外一种直接在频率或者时域内进行，由于需要解联立方程，设计时需要计算机做辅助设计。利用模拟滤波器成熟的理论设计IIR数字滤波器的过程是：(1)确定数字低通滤波器的技术指标：通带边界频率、通带最大衰减、 阻带截止频率、阻带最小衰减。(2)将数字低通滤波器的技术指标转换成相应的模拟低通滤波器的技术指标。(3)按照模拟低通滤波器的技术指标设计过渡模拟低通滤波器。(4)用所选的转换方法，将模拟滤波器转换成数字低通滤波器系统函数4。先在模拟（连续）域根据技术指标设计模拟低通滤波器，后转化到数字域，用脉冲响应不变法设计数字滤波器，要求通带和阻带具有单调下降的特性，指标参数如下：wp=0.2rad,p=1dB,ws=0.35rad, s=25dB。在matlab里求解。设计程序如下：%用脉冲响应不变法设计数字低通滤波器T=1;%T=1swp=0.2*pi/T;ws=0.35*pi/T;rp=1;rs=25;%T=1s的模拟滤波器的指标N,wc=buttord(wp,ws,rp,rs,s);%计算响应的模拟滤波器阶数N和3dB截止频率wcB,A=butter(N,wc,s);%计算该响应的模拟滤波器的系统函数Bz,Az=impinvar(B,A)%用脉冲响应不变法将模拟转为数字滤波器Az=1.0000, -3.8255, 6.7082 , -6.8507,4.3564, -1.7136,0.3843 ,-0.0378;Bz= -0.0000, 0.0001,0.0032,0.0105, 0.0066, 0.0008,0.0000, 0;A,B,C,D=tf2ss(Bz,Az)在以上程序中，先根据模拟滤波器的指标，求出阶数和截止频率wc，然后调用语句计算此滤波器的系统函数，再用语句A,B,C,D=tf2ss(Bz,Az)将系统函数转化为状态空间，并用传递函数的系数，依据相关知识，写出此离散系统的差分方程便于理论研究，差分方程如下：y(k)+a*y(k-1)+b*y(k-2)+c*y(k-3)+d*y(k-4)+l*y(k-5)+m*y(k-6)+n*y(k-7)=a1*f(k-1)+b1*f(k-2)+c1*f(k-3)+d1*f(k-4)+l1*f(k-5)。 y(-1)=y(-2)=q,y(-3)=p1,y(-4)=q1,y(-5)=h,y(-6)=h1,y(-7)=j以上方程系数变量及初始状态完全描述此滤波器系统。3.2 创建系统的simulink状态空间模型 在仿真过程中，开始也学习认识到MATLAB/SIMULINK系统建模设计与仿真基础知识，例如MATLAB在系统建模中的应用有求取系统传递函数（矩阵），求取幅频特性和相频特性，系统输出响应，求系统状态空间模型，数值分析与函数逼近等等。数学仿真采用数学模型，用数学语言对系统的特性进行描述，其工作过程是：1.建立系统的数学模型；2. 建立系统仿真模型，即设计算法，并转化为计算机程序，使系统的数学模型能为计算机所接受并能在计算机上运行；3.运行仿真模型，进行仿真试验，再根据仿真试验的结果进一步修正系统的数学模型和仿真模型5。模型基本结构：一个典型的Simulink 模型包括如下三种类型的元素：信号源模块，被模拟的系统模块，输出显示模块。两种Simulink 运行仿真的方法使用窗口运行仿真；使用MATLAB 命令运行仿真。使用窗口运行仿真优点：人机交互性强，不必记住繁琐的命令语句即可进行操作。使用窗口运行仿真主要可以完成以下一些操作。1. 设置仿真参数 2. 应用仿真参数 3. 启动仿真4. 停止仿真 5. 中断仿真 6. 仿真诊断。状态空间模型是动态时域模型，以隐含着的时间为自变量。状态空间模型在经济时间序列分析中的应用正在迅速增加6。根据以上准备工作，先在simulink中建模如下，信源选择正弦波的叠加。如下： 图3.2 用状态空间创建仿真模块3.3 模块内部参数设置及数据存储 在模块创建中，各个模块的内部参数设置如下。对于信号源模块，是用两个不同的正弦信号源的叠加，内部设置参数如下图3.3。 其中，一个信号源的频率设置为0.2*pi，另一个设置为0.5*pi，这使得只有一个频率成分可以通过，正好有滤波效果。状态空间内部参数设置为矩阵A,B,C,D，初始条件为x0.如图3.4。 数据存储即工作空间如图3.5，其中时间为t为结算时间，x表示states，是指输入变量的存储，y为output,是输出变量的存储。图3.3 正弦信号源内部参数设置图3.4 状态空间模型内部参数设置图3.5 数据存储4 利用simulink状态空间建模求解LTI系统的优缺点用状态空间建模设计程序求解高阶离散系统系统。其优点为：可应用范围广，实用性强，可刻画系统的内部特征，用“一阶微分方程组”来描述系统的动态特性。它的主要特点是： （1）利用描述系统内部特性的状态变量替代了仅能描述系统外部特性的系统函数，能完整地揭示系统的内部特性，从而使得控制系统的分析和设计产生根本性的变革。 （2）便于处理多输入多输出系统。 （3）一阶微分（或差分）方程组便于计算机数值计算。 （4）容易推广用于时变系统和非线性系统7。 缺点为：（1）、不直观，几何、物理意义不明显：不象经典法那样，能用Bode 图及根轨迹进行直观的描述。对于简单问题，显得有点烦琐。 （2）对数学模型要求很高：而实际中往往难以获得高精度的模型，这妨碍了它推广和应用。5 连续系统simulink状态空间建模分析方法程序设计的思路5.1 调用模型文件及编程求解系统响应求解响应MATLAB程序如下：% 高阶离散系统用simulink状态空间函数模块建模，仿真求解系统的方法%y(k)+a*y(k-1)+b*y(k-2)+c*y(k-3)+d*y(k-4)+l*y(k-5)+m*y(k-6)+n*y(k-7)=a1*f(k-1)+b1*f(k-2)+c1*f(k-3)+d1*f(k-4)+l1*f(k-5)% y(-1)=p, y(-2)=q,y(-3)=p1,y(-4)=q1,y(-5)=h,y(-6)=h1,y(-7)=j% 以上系统完全描述了此高阶离散系统clc;clear all;close all;a=-3.8255; b=6.7082; c= -6.8507;d=4.3564;l=-1.7136;m=0.3843;n=-0.0378;a1=0.0001;b1=0.0032;c1=0.0105;d1=0.0066;l1=0.0008;p=0.3; q=0.4;p1=0.6;q1=0.8;h=0.25;h1=0.4;j=0.2;A=3.8255 -6.7082 6.8507 -4.3564 1.7136 -0.3843 0.0378;1 0 0 0 0 0 0;0 1 0 0 0 0 0;0 0 1 0 0 0 0;0 0 0 1 0 0 0;0 0 0 0 1 0 0;0 0 0 0 0 1 0;u=5;%正弦波的幅度为5% 定步长求解，步长为Ts, 求解时区为 Tspan=0,100T=1; Ts=0.01; %信号源采样时间，即定步长解算器的解算步长tf=20; L=round(tf/Ts); %解算时区取为0,tfTspan=0,(L-1)*Ts;% =yzi=u=0;B=0;0;0;0;0;0;0; C=0.0001 0.0032 0.0105 0.0066 0.0008 0 0;D=0;x0=j;h1;h;q1;p1;q;p;t,x,y=sim(ylj1,Tspan);%状态变量x的第一列是x(k-7),第二列是x(k-6)，以此类推yzi=-a*x(:,7)-b*x(:,6)-c*x(:,5)-d*x(:,4)-l*x(:,3)-m*x(:,2)-n*x(:,1);%=yzs=u=5;B=1;0;0;0;0;0;0;C=0.0001 0.0032 0.0105 0.0066 0.0008 0 0; D=0;x0=0;t,x,y=sim(ylj1,Tspan);%状态变量x的第一列是x(k-2),第二列是x(k-1);yzs=y(:,1);yt=yzi+yzs;% = 绘制波形 =subplot(311);plot(t,yzi,r-);xlabel(k);ylabel(yzi(k);title(高阶离散系统在正弦波激励下的零输入响应);subplot(312);plot(t,yzs,g-);xlabel(k);ylabel(yzs(k);title(高阶离散系统在正弦波激励下的零状态响应);subplot(313);plot(t,yt,k-);xlabel(k);ylabel(y(k); title(高阶离散系统在正弦波激励下的全响应);程序运行结果与分析： 图5.1高阶系统在激励下的各种响应在运行结果第一幅图中，主要绘制滤波器系统在激励下的各种响应。在零输入响应中，由于系统本身的初始状态存储一定的能量，所以输出波形在大约0-0.4s内发生谐振现象，之后能量损耗完，波形迅速衰减为零。其中，此时系统的零状态响应与其的初始状态成线性关系。在零状态响应中，由于系统的响应与输入信号呈线性关系，与初始状态无关，而系统的信号源是两个不同频率正弦波的叠加，所以输出波形类似于正弦波，它是以20s为周期的，且具有一定的滤波作用，使波形平滑，而全响应是两者的波形的叠加。5.2 分析系统的频谱与相位程序如下： fs=1/Ts;tao=Tspan(2)+Ts;df=1/tao;L=length(t);L2=round(L/2);Y1=fft(yzi)*(2/L);Y1(1)=Y1(1)/2;AY1=abs(Y1(1:L2);Am1=max(AY1);BY1=angle(Y1(1:L2);Bm1=max(BY1);Y2=fft(yzs)*(2/L);Y2(1)=Y2(1)/2;AY2=abs(Y2(1:L2);Am2=max(AY2);BY2=angle(Y2(1:L2);Bm2=max(BY2);F=(0:(L2-1)*df;figure(2);subplot(211);stem(F,AY1,r.);xlabel(f/Hz);ylabel(A1(f);title(系统零输入响应频谱);subplot(212);stem(F,AY2,b.);xlabel(f/Hz);ylabel(A2(f);title(系统零状态响应频谱);figure(3);subplot(211);stem(F,BY1,r.);xlabel(f/Hz);ylabel(B1(f);title(系统零输入响应相位);subplot(212);stem(F,BY2,g.);xlabel(f/Hz);ylabel(B1(f);title(系统零状态响应相位); 图5.2 高阶系统的频谱分析 图5.3高阶系统的相位分析在第一幅图中，是分析两种响应的频谱。其中，对于零状态响应，由于输入不为零，又因为此系统为低通滤波器系统，所以从运行结果来看，其幅频特性大约在0.2，即0.628HZ处截止，并设置了一定的过渡带，在技术指标中阻带最小衰减由于比较大，所以滤波效果比较理想。振幅最大值为5.在零输入响应中，其频谱主要是由初始状态决定，而初始状态设置较小，所以振幅较小，最大为0.037，大约在25HZ左右逐渐趋于0，反应了滤波器系统本身的属性。6 状态空间分析方法的应用实例6.1 分析求解低阶电路系统 例如分析二阶RLC电路， 图6.1 二阶RLC电路其中各元件的伏安关系为 , , （式6.1）以为响应,则得到以为响应的微分方程为： （式6.2）将电阻电路映射到S域，利用拉普拉斯变换进行分析8：电阻不变，电感，电容,则求出系统函数 （式6.3）之后根据H(S),求出状态空间矩阵，编程求解此电路系统。编程思想大致与本课题相同。6.2 设计分析滤波器系统例如合计数字低通滤波器，要求指标参数如下：wp=0.2rad, p=1dB,ws=0.35rad, s=25dB。具体程序如下：%用双线性变换法设计DFT=1;Fs=1/T;wpz=0.2;wsz=0.35;wp=2*tan(wpz*pi/2);ws=2*tan(wsz*pi/2);rp=1;rs=25；%预畸变矫正转换指标N,wc=buttord(wp,ws,rp,rs,s); %计算响应的模拟滤波器阶数N和3dB截止频率wcB,A=butter(N,wc,s); %计算该响应的模拟滤波器的系统函数Bz,Az=bilinear(B,A,Fs); %用双线性变换法转换成数字滤波器Nd,wdc=buttord(wpz,wsz,rp,rs); %调用butter和buttord直接设计数字滤波器Bdz,Adz=butter(Nd,wdc)9得出状态空间矩阵后，可以写出微分方程，后编程求解系统。7 结束语现代的社会是信息的社会，很多与人们生活紧密相关的都是通信领域，比如决定现代社会发展方向的电脑以及几乎每个同学都拥有的收音机，手机等等无不说明通信系统的重要性。学会的软件如Matlab，并对通信系统应用和开发的设计思想有了更进一步的了解和掌握。通过这次做毕业设计发现自己的看法有点太片面。毕业设计不仅是对前面所学知识的一种检验，而且也是对自己能力的一种提高。通过这次毕业设计使我明白了自己原来知识还比较欠缺。自己要学习的东西还太多，以前老是觉得自己什么东西都会，什么东西都懂，有点眼高手低。通过这次毕业设计，我才明白学习是一个长期积累的过程，无论是在以后的工作还是生活中都应该不断的学习，努力提高自己知识和综合素质。在设计的过程中我总会遇到这样那样的问题，我有失落过，烦恼过，悲伤过，但我明白这又是我人生中的一大挑战，角色的转换，这除了有较强的适应力和乐观