信号与系统仿真实验报告（精选8篇）
信号与系统仿真实验报告 篇1
仿真
实 验 报 告
班级： 学号： 姓名： 学院：
实验一
一、实验者姓名：
二、实验时间：
三、实验地点：
四、实验题目：
5(s25s6)求三阶系统H(s)3的单位阶跃响应，并绘制响应波形图。
s6s210s8
五、解题分析：要知道求单位阶跃响应需知道所用函数，以及产生波形图所需要用到的函数。
六、试验程序：
num=[5 25 30];den=[1 6 10 8];step(num,den,10);title(‘Step response’)
七、实验结果：
实验所得波形图如下：
Step response4.543.53Amplitude2.521.510.50012345Time(sec)678910
八、实验心得体会：通过本次试验了解学会了一些新的函数的应用。了解到了N阶系统的单位阶跃响应的计算方法，和系统的响应波形图的函数应用和绘制方法。为后面的实验打下基础，并对信号仿真和《信号与系统》这门课程之间的联系有所增加，对《信号与系统》这门课里的问题也有了更加深入地了解。
九、实验改进想法：无。
实验二
一、实验者姓名：
二、实验时间：
三、实验地点：
四、实验题目:
一个因果线性移不变系统y(n)0.81y(n2)x(n)x(n2)，求：（1）H(z)；（2）冲激响应h(n)；（3）单位阶跃响应u(n)；（4）H(ej)，并绘出幅频和相频特性。
五、解题分析：离散卷积是数字信号处理中的一个基本运算，MTLAB提供的计算两个离散序列卷积的函数是conv，其调用方式为 y=conv(x,h)。其中调用参数x,h为卷积运算所需的两个序列，返回值y是卷积结果。
MATLAB函数conv的返回值y中只有卷积的结果，没有y的取值范围。由离散序列卷积的性质可知，当序列x和h的起始点都为k=0时，y的取值范围为k=0至length(x)+length(h)-2。
许多离散LTI都可用如下的线性常系数的差分方程描述
ay[kn]bx[kn]
nnn0n0NN其中x[k]、y[k]分别系统的输入和输出。在已知差分方程的N个初始状态y[k]，和输入x[k]，就可由下式迭代计算出系统的输出
y[k](an/a0)y[kn](bn/b0)x[kn]
n1n0NM利用MATLAB提供的filter函数，可方便地计算出上述差分方程的零状态响应。filter函数调用形式为 y=filter(b,a,x)。其中 a[a0,a1,...,aN]，b[b0,b1,...,bM]，分别表示差分方程系数。X表示输入序列，y表示输出序列。输出序列的长度和序列相同。
当序列的DTFT可写成ej的有理多项式时，可用MATLAB信号处理工具箱提供的freqz函数计算DTFT的抽样值。另外，可用MATLAB提供的abs、angle、real、imag等基本函数计算 DTFT的幅度、相位、实部、虚部。若X（ej）可表示为
b0b1ej...bMejMB(ej)X(e)jjjNA(e)a0a1e...aNe则freqz的调用形式为 X=freqz(b,a,w)，其中的b和 a分别是表示前一个
j式子中分子多项式和分母多项式系数的向量，即a[a0,a1,...,aN]，w为抽样的频率点，向量w的长度至少为2。返回值X就是DTFTb[b0,b1,...,bM]。在抽样点w上的值。注意一般情况下，函数freqz的返回值X是复数。
六、实验程序:
clc;clear;close;b=[1 0-1];a=[1 0-0.81];figure(1);subplot(2,1,1);dimpulse(b,a,20)subplot(2,1,2);dstep(b,a,50)w=[0:1:512]*pi/512;figure(2);freqz(b,a,w)
七、实验结果：
冲击响应图及阶跃响应图:
Impulse Response1Amplitude0.50-0.50246810Time(sec)Step Response12141618201Amplitude0.500510152025Time(sec)3035404550 100Magnitude(dB)0-100-200-30000.10.20.30.40.50.60.70.8Normalized Frequency( rad/sample)0.91100Phase(degrees)500-50-10000.10.20.30.40.50.60.70.8Normalized Frequency( rad/sample)0.91
八、实验心得体会：通过实验我们知道了使用Matlab来绘出出一个线性移不变系统的幅频和相频曲线。并知道了在《信号与系统》中得一些差分方程和各种响应，譬如零输入相应、零状态响应、全响应、自由响应、强迫响应、冲击响应、单位阶跃响应等等各种响应在Matlab中的函数表达方式和他们的求法，以及系统的幅频和相频曲线的绘制都有了一定深刻的认识。
九、实验改进想法：无。
实验三
一、实验者姓名：
二、实验时间：
三、实验地点：
四、实验题目：
模拟信号x(t)2sin(4t)5cos(8t)，求N64的DFT的幅值谱和相位谱。
五、解题分析：在MATLAB信号处理工具箱中，MATLAB提供了4个内部函数用于计算DFT和IDFT，它们分别是：fft(x)，fft(x,N)，ifft(X)，ifft(X,N)。
fft(x)计算M点的DFT。M是序列x的长度，即M=length(x)。
fft(x,N)计算N点的DFT。若M>N，则将原序列截短为N点序列，再计算其N点DFT；若M ifft(X)计算M点的IDFT。M是序列X的长度。 ifft(X,N)计算N点IDFT。若M>N，则将原序列截短为N点序列，再计算其N点IDFT；若M 六、实验程序： clc;clear;close;N=64;n=0:63;t=d*n;q=n*2*pi/N;x=2*sin(4*pi*t)+5*cos(8*pi*t);y=fft(x,N);subplot(3,1,1);plot(t,x);title(‘source signal’);subplot(3,1,2);plot(q,abs(y));title(‘magnitude’);subplot(3,1,3);plot(q,angle(y));title(‘phase’); 七、实验结果： ***0100806040200|F(k)|05101520Frequency253035 ***0100806040200|F(k)|05101520Frequency253035 4321|jW|0-1-2-3-405101520Frequency253035Step Response400020000-2000 Amplitude-4000-6000-8000-10000-12000-1400001234n(samples)5678 八、实验心得体会：通过本次试验我知道了求取模拟信号在N等于一定值时的的DFT的幅值谱和相位谱的求法。通过本次实验，对幅值谱和相位谱有了更深的了解，并与课程《信号与系统》里的一些相关知识连接到了一起，使得学到的只是更加深刻、有意义。 九、实验改进想法：无。 实验四 一、实验者姓名： 二、实验时间： 三、实验地点： 四、实验题目： 将信号x(t)sin(240t)做离散序列，比较原序列与经过FFT和IFFT变换后的序列，并做出说明。 五、解题分析：此题需要对信号做离散序列，还要做FFT和IFFT变换，然后得到图像进行比较。连续时间函数与离散时间函数在编程中的区别主要体现在如下两个方面：第一，自变量的取值范围不同，离散时间函数的自变量是整数，而连续时间函数的自变量为一定范围内的实数；第二，绘图所用的函数不同，连续函数图形的绘制不止一个。本实验中要求绘制离散时间信号图，可以应用MATLAB中的函数来实现。用MATLAB表示一离散序列，可用两个向量来表示。其中一个向量表示自变量的取值范围，另一个向量表示序列的值。之后画出序列波形。当序列是从0开始时，可以只用一个向量x来表示序列。由于计算机内寸的限制，MATLAB无法表示一个无穷长的序列。对于典型的离散时间信号，可用逻辑表达式来实现不同自变量时的取值。 六、实验程序： t=0:1/255:1;x=sin(2*pi*120*t);y=real(ifft(fft(x)));subplot(2,1,1);plot(t,x);title(‘原波形’);subplot(2,1,2);plot(t,y); 七、实验结果： 原波形10.50-0.5-100.10.20.30.40.50.60.70.80.91恢复的波形10.50-0.5-100.10.20.30.40.50.60.70.80.91 八、实验心得体会：通过对做信号的离散序列以及经FFT和IFFT的变换，了解了相关特性。通过计算机做出的信号波形图，我们能够很直白的看出原波形和经过变换后的波形的差别。 九、实验改进想法：无。 实验五 一、实验者姓名： 二、实验时间： 三、实验地点： 四、实验题目： 2s，激励信号22(s1)100x(t)(1cot)sco1s0(t)0，求（1）带通滤波器的频率响应；（2）输出稳态响应并绘制图形。已知带通滤波器的系统函数为H(s) 五、解题分析：需要知道求频率响应的方法，并绘制图形。 六、实验程序： clear;t=linspace(0,2*pi,1001);w=[99,100,101];U=[0.5,1,0.5];b=[2,0];a=[1,2,10001];u1=U*cos(w’*t+angle(U’)*ones(1,1001));H=polyval(b,j*w)./polyval(a,j*w);H=freqs(b,a,w);subplot(2,1,1),plot(w,abs(H)),grid;subplot(2,1,2),plot(w,angle(H)),grid;u21=abs(U(1)*H(1))*cos(99*t+angle(U(1)*H(1)));u22=abs(U(2)*H(2))*cos(100*t+angle(U(2)*H(2)));u23=abs(U(3)*H(3))*cos(101*t+angle(U(3)*H(3)));u2=u21+u22+23;figure(2);subplot(2,1,1),plot(t,u1);subplot(2,1,2),plot(t,u2); 七、实验结果： 10.90.80.79910.50-0.5-19999.299.499.699.8100100.2100.4100.6100.810199.299.499.699.8100100.2100.4100.6100.8101 210-1-***222101234567 八、实验心得体会：通过本次试验，了解了频率响应求法，加深了对输出稳态响应的印象。 1 系统设计的原理和方法 本设计描述了信号与系统仿真系统的作用、要求和构成。详细介绍了利用MATLAB软件来实现信号与系统仿真系统的基本原理及功能, 利用MATLAB软件提供的图形用户界面 (Graphical User Inte rface, GUI) 设计具有界面友好、人机交互便利等特点的用户界面。 在仿真实验中, 包括了许多信号与系统的子实验, 涉及范围广泛表现良好。本系统是一种演示型的软件, 可以利用可视化的仿真工具, 以及图形和动态仿真的方式来演示部分基本信号的波形和变换, 可以使学习人员清楚、明了的明白和掌握信号与系统的一些基础知识。 2 系统设计的总体方案和实现 通过添加GUI中的功能按钮, 可以实现各种各样功能转换, 如图1所示。 3 基础信号设计 基础信号包括:信号选择、正弦信号、指数信号、阶跃信号、脉冲信号、矩形脉冲信号、抽样信号。Axis1用于显示图形, 同时使图形可以放大缩小。将text1~6的string改为“自定义参数”, “A=”, “B=”, “C=”, “D=”, “简介”。Edit1~4用于改变信号的参数, 并将更改后的参数传递到popupmenu1中参与计算。Edit5的作用是显示各个函数的简介, 具体显示的内容由popupmenu1控制。pushbutton1转换界面。 4 语音信号设计 4.1 设计思路 本设计可以将本地语音信号读取到软件中, 同时也可以录制语音并保存。录制时用户可以根据自己的要求改变录制时长和采样频率。读取信号之后能对信号进行各种处理, 如变速播放、显示信号时域图形、快速傅里叶变换、显示信号频域图形。 pus hbutton1用于读取本地语音信号。pus hbutton2的作用是播放读取在软件中的变量。pushbutton4的作用是通过电脑自带的麦克风录制声音。pushbutton6的作用是将软件中的数据保存在本地文件中。pushbutton10的作用是将数据进行FFT变换。pushbutton3的作用为把数据以图形的形式显示在axes1中。 4.2 快速傅里叶变换程序设计 有限长序列可以通过离散傅里叶变换 (DFT) 将其频域也离散化成有限长序列。但其计算量太大, 很难实时地处理问题, 因此引出了快速傅里叶变换 (FFT) 。1965年, Cooley和Tukey提出了计算离散傅里叶变换 (DFT) 的快速算法, 将DFT的运算量减少了几个数量级。从此, 对快速傅里叶变换 (FFT) 算法的研究便不断深入, 数字信号处理这门新兴学科也随FFT的出现和发展而迅速发展。根据对序列分解与选取方法的不同而产生了FFT的多种算法, 基本算法是基2DIT和基2DIF。FFT在离散傅里叶反变换、线性卷积和线性相关等方面也有重要应用。 快速傅里叶变换 (FFT) 是计算离散傅里叶变换 (DFT) 的快速算法, 实现的流程图如图2所示。 5 设计结果与分析 整个系统设计编写已经完成并调试成功, 能够明确显示各种信号的仿真。经测试软件的全部内容与《信号与系统》一书内容符合, 完全可用于帮助该课程的图形显示。 本设计思路明确, 全部内容已用MATLABGUI来实现, 对于课程的初学者具有加深理解的现实意义。 参考文献 [1]奥本海姆AV, 威斯基AS.信号与系统[M].刘树棠, 译.西安:西安交通大学出版社, 1998. [2]宗孔德, 胡广书.数字信号处理[M].北京:清华大学出版社, 1988. [3]李强, 明艳, 陈前斌, 等.基于Matlab的数字信号处理实验仿真系统的实现[J].实验技术与管理, 2006. 【关键词】优质课程 信号与系统 实验教学 体系 【中图分类号】G642 【文献标识码】A 【文章编号】1006-9682（2012）10-0010-02 高等教育“十一五”规划纲要中明确指出，21世纪的高等教育要以质量工程建设为核心。在高等教育质量工程建设中，优质课程及精品课程的建设占有极其重要的地位。“信号与系统”是工科院校开设的一门重要专业基础课，在许多高校都列为重点建设课程。[1～2]我学院针对两个一级学科（仪器科学与技术、光学工程）开设“信号与系统”课程，不仅要求充分讲授该课程的经典内容，还必须兼顾学科及专业需求，进行有特色的理论教学及实验教学。作为校级优质课程，仪器专业“信号与系统”课程长期以来持续建设，已经形成了理论教学扎实、实验教学强化的特点。但由于目前本科培养计划中学时数的限制，无法兼顾本课程理论及实验多方面的要求，往往只能采用软件仿真作为实验教学的主要内容。考虑到我校的办学定位及仪器专业应用型人才培养的需求，我们在原有理论及实验教学的基础上，开展全方位实验教学体系建设，进一步提高本课程的教学效果。 一、“信号与系统”教学内容及实验内容的演变 1.教学内容的演变 信号是信息学科研究的基本内容，信号与系统是两个用得极为广泛且密切相关的基本概念。在许多实际应用中，尤其是在信号提取、信号恢复、信号增强、语音识别等信号处理的问题中，以及在大规模集成电路的整体设计中，信号、系统、处理往往是有机结合在一起的。因此，教学内容从单纯讲信号、系统演变为信号、系统与数字信号处理融合，以信号分析为基础，以系统分析为桥梁，以处理技术为手段，形成新的教学体系，适应新的科技条件下对专业基础课的教学要求。[3～4]因此，仪器专业“信号与系统”课程的教学内容主要包括：信号与系统的基本概念、采样定理、连续及离散线性时不变系统的时域分析、连续及离散信号的频域分析（含离散傅里叶变换、快速傅里叶变换）、连续及离散信号与系统的复频域分析、数字滤波器设计等。 2.实验内容的演变 实验教学始终是为理论教学服务的，其最终的目的是为了强化理论学习，全面提高教学质量。早期“信号与系统”的实验以硬件为主，尤其是以电路实验为主。随着信息技术的不断发展和信息技术应用领域的不断扩展，这门课程已从电子信息工程类专业的专业基础课演变为众多工科专业（如计算机技术、自动控制、测控技术及仪器、生物医学工程等）的专业基础课，其实验领域也获得拓展。特别是伴随着计算机软硬件技术的快速发展，本课程实验由早期的硬件实验演变为软件仿真实验。[5]由于本课程仍然是一门实践性很强的课程，无论技术如何发展演变，信号也脱离不了实际应用中的物理系统，因此现在又强调从软件仿真实验 演变为软硬件实验结合。[6] 二、“信号与系统”实验教学计算机软硬件应用方案 “信号与系统”课程理论性强，实践性强，实验对于理论具有巩固和强化的作用。由于学时的限制，课内实验具有较大的局限性。一方面是实验数量受到局限，只能开设最基本、最重要的实验，无法进行多个理论的验证，而本课程又是一门逻辑性很强的课程，基本概念环环紧扣，实验数量少非常妨碍学生对课程基本内容的理解和掌握，从而使学生的知识面受到局限；另一方面，在实验方法上也受到局限，不能尝试多种解决问题的方法，只能用常规方法去做，这样限制了学生的思维，不利于学生综合能力的培养。目前很多学校以软件仿真实验为主，由于学时的限制以及软件仿真实验的内容十分丰富且复杂的，所以就舍弃了硬件实验。无论软件仿真怎样逼近实际，毕竟与实际物理系统存在差异。如果没有硬件实验，学生便不知道其结果应该怎样去应用，在什么情况下采用软件计算的结果为好，在什么情况下用硬件实现较好。从课程教学质量和应用型人才培养要求来看，这是一个亟待解决的严重问题。因此，计算机软硬件实验对于“信号与系统”实验体系都是不可缺少的，要两条腿走路，其关键是要完成硬件实验和软件实验的全方位比较，才能深入理解信号与系统的理论实质。 为了更好地完成辅助理论教学的任务，我们需要对本课程的实验教学体系进行精心的设计。实验教学体系的建立要从整个课程教学体系来考虑，在课堂讲授内容、学生作业内容、教学目标等方面综合权衡。在一个典型的数字信号处理系统中，输入通道中的传感器输出信号需要进行调理，属于物理系统部分，硬件实现比较好；数字信号处理部分可以采用单片机、数字信号处理器、计算机等硬件设备，同时利用硬件平台上的软件来完成数字信号处理任务，这部分以软件仿真实验为好，而在输出通道需要进行必要的信号变换并输出模拟信号，还是硬件实验较好。本课程采用的实验教学体系，见图1。 在课内实验安排方面：实验内容强调经典、重要、基本，少而精，始终抓住信号产生→获取→处理→使用这条主线，不仅有利于教学安排，而且保证学生基础牢固，知识更系统，理解更全面。课内实验以软件仿真为核心，其编程软件采用公认的优秀软件，即MATLAB。MATLAB是优秀的科学计算和仿真软件，研究设计单位和工业部门同样公认它的重要价值。如美国NI公司的信号测量与分析软件LabVIEW、Cadence公司的信号和通信分析设计软件、TI公司的DSP等都和MATLAB具有良好的接口。现在的计算机硬件（PC机）配置很高，计算精度很高，用于进行数字信号处理和数字图像处理十分理想，结果显示也很直观。 MATLAB 综合实验 项目二 连续系统的频域阐发 目的： 周期信号输入连续系统的响应可用傅里叶级数阐发。由于盘算历程啰嗦，最适适用MATLAB 盘算。通过编程实现对输入信号、输出信号的频谱和时域响应的盘算，认识盘算机在系统阐发中的作用。 任务： 线性连续系统的系统函数为11)(jj H，输入信号为周期矩形波如图 1 所示，用MATLAB 阐发系统的输入频谱、输出频谱以及系统的时域响应。 -3-2-1 0 1 2 300.511.52Time(sec) 图 1 要领： 1、确定周期信号 f(t)的频谱nF。基波频率 Ω。 2、确定系统函数 )( jn H。 3、盘算输出信号的频谱 n nF jn H Y )( 4、系统的时域响应 nt jnn eY t y)（MATLAB 盘算为 y=Y_n*exp(j*w0*n“*t); 要求（画出 3 幅图）： 1、在一幅图中画输入信号 f(t)和输入信号幅度频谱|F(j)|。用两个子图画出。 2、画出系统函数的幅度频谱|H(j)|。 3、在一幅图中画输出信号 y(t)和输出信号幅度频谱|Y(j)|。用两个子图画出。 解:(1)阐发盘算： 输入信号的频谱为 (n)输入信号最小周期为 =2，脉冲宽度，基波频率Ω=2π/ =π，所以 (n)系统函数为 因此 输出信号的频谱为 系统响应为 (2)步伐： t=linspace(-3,3,300); tau_T=1/4; % n0=-20;n1=20; n=n0:n1; %盘算谐波次数20 F_n=tau_T*Sa(tau_T*pi*n); f=2*(rectpuls(t+1.75,0.5)+rectpuls(t-0.25,0.5)+rectpuls(t-2.25,0.5)); figure(1),subplot(2,1,1),line(t,f,”linewidth“,2); %输入信号的波形 axis([-3,3,-0.1,2.1]);grid on xlabel(”Time(sec)“,”fontsize“,8),title(”输入信号“,”fontweight“,”bold“)%设定字体巨细，文本字符的粗细 text(-0.4,0.8,”f(t)“) subplot(2,1,2),stem(n,abs(F_n),”.“); %输入信号的幅度频谱 xlabel(”n“,”fontsize“,8),title(”输入信号的幅度频谱“,”fontweight“,”bold“) text(-4.0,0.2,”|Fn|“) H_n=1./(i*n*pi+1); figure(2),stem(n,abs(H_n),”.“); %系统函数的幅度频谱 xlabel(”n“,”fontsize“,8),title(”系统函数的幅度频谱“,”fontweight“,”bold“) text(-2.5,0.5,”|Hn|“) Y_n=H_n.*F_n;y=Y_n*exp(i*pi*n”*t); figure(3),subplot(2,1,1),line(t,y,“linewidth”,2); %输出信号的波形 axis([-3,3,0,0.5]);grid on xlabel(“Time(sec)”,“fontsize”,8),title(“输出信号”,“fontweight”,“bold”) text(-0.4,0.3,“y(t)”) subplot(2,1,2),stem(n,abs(Y_n),“.”); %输出信号的幅度频谱 xlabel(“n”,“fontsize”,8),title(“输出信号的幅度频谱”,“fontweight”,“bold”) text(-4.0,0.2,“|Yn|”) (3)波形： -3-2-1 0 1 2 300.511.52Time(sec)输 入 信 号f(t)-20-15-10-5 0 5 10 15 2000.10.20.30.4n输 入 信 号 的 幅 度 频 谱|Fn|-20-15-10-5 0 5 10 15 2000.10.20.30.40.50.60.70.80.91n系 统 函 数 的 幅 度 频 谱|Hn| -3-2-1 0 1 2 300.10.20.30.4Time(sec)输 出 信 号y(t)-20-15-10-5 0 5 10 15 2000.10.20.30.4n输 出 信 号 的 幅 度 频 谱|Yn| 项目三 连续系统的复频域阐发 目的： 周期信号输入连续系统的响应也可用拉氏变更阐发。用 MATLAB 的标记盘算成果，通过编程实现对系统瞬态响应和稳态响应的阐发，加深理解拉氏变更在阐发系统中的作用。 任务： 线性连续系统的系统函数为11)(ss H，输入信号为周期矩形波如图2所示，用MATLAB阐发系统的响应和稳态响应。 0 1 2 3 4 5 6 700.511.52Time(sec) 图 2 要领： 1、确定第一个周期拉氏变更)(0s F。 2、确定前 6 个周期的拉氏变更)(s F。 3、盘算输出信号的拉氏变更)()()(s F s H s Y 4、系统的时域响应)()(s Y t y MATLAB 盘算为 y=ilaplace(Y);5、系统的稳态响应和稳态值，即经过 4 个周期后，系统响应趋于稳态，两个稳态值可取为 t=8s 和 t=8.5s 要求： 1、画出输入信号 f(t)波形。 2、画出系统输出信号 y(t)的波形。 3、画出系统稳态响应 yss(t)的波形，4 个周期后。并盘算出稳态值。 解:(1)步伐 syms s; H=1/(s+1); F0=1/s*(1-exp(-0.5*s)); %输入信号第一个周期的laplace变更 F=F0+F0*exp(-2*s)+F0*exp(-4*s)+F0*exp(-6*s); Y=H.*F; Y0=H.*F0; y=ilaplace(Y); y=simple(y); t=linspace(0,12,300); f=2*(rectpuls(t-0.25,0.5)+rectpuls(t-2.25,0.5)+rectpuls(t-4.25,0.5)+rectpuls(t-6.25,0.5)); yn=subs(y); %标记替换 figure(1),plot(t,f,“linewidth”,2); axis([0,7,-0.2,2.2]),xlabel(“Time(sec)”,“fontsize”,8),title(“输入信号”,“fontweight”,“bold”) text(3.0,1.0,“f(t)”) figure(2),plot(t,yn,“linewidth”,2); axis([0,7,-0.1,0.5]),xlabel(“Time(sec)”,“fontsize”,8),title(“输出信号”,“fontweight”,“bold”) text(3.0,0.3,“y(t)”) figure(3),plot(t,yn,“linewidth”,2); axis([8,12,-0.1,0.5]),xlabel(“Time(sec)”,“fontsize”,8),title(“输出信号稳态响应”,“fontweight”,“bold”) text(10.0,0.2,“ys(t)”) t=8:0.5:8.5; %取t=8s和t=8.5两个稳态值 ys=subs(y,t,“t”); disp(“输入为周期信号的响应的第一个周期”); y0=ilaplace(Y0); pretty(y0); %标记输出类似数值形式 disp(“输出稳态周期信号的两个值”); ys (2)波形 0 1 2 3 4 5 6 700.511.52Time(s ec)输 入 信 号f(t) 0 1 2 3 4 5 6 7-0.100.10.20.30.4Time(sec)输 出 信 号y(t)8 8.5 9 9.5 10 10.5 11 11.5 12-0.100.10.20.30.4Time(sec)输 出 信 号 稳 态 响 应ys(t) 命令窗口显示： 输入为周期信号的响应的第一个周期 heaviside(t-1/2)(exp(1/2-t)-1)-exp(-t)+ 1 输出稳态周期信号的两个值 ys = 0.1015 实验一 MATLAB 中矩阵与多项式的基本运算 实验任务 1.了解 MATLAB 命令窗口和程序文件的调用。•熟悉如下 MATLAB 的基本运算： ① 矩阵的产生、数据的输入、相关元素的显示； ② 矩阵的加法、乘法、左除、右除； ③ 特殊矩阵：单位矩阵、“ 1 ”矩阵、“ 0 ”矩阵、对角阵、随机矩阵的产生和 运算； ④ 多项式的运算：多项式求根、多项式之间的乘除。 基本命令训练 1、>> eye(2)ans = 1 0 0 1 >> eye(4)ans = 1 0 0 0 0 1 0 0 0 0 1 0 0 0 0 1 2、>> ones(2)1 1 ans =1 >> ones(4)ans = 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 >> ones(2,2)ans = 1 1 1 1 >> ones(2,3)ans = 1 1 1 1 1 1 >> ones(4,3)ans = 1 1 1 1 1 1 1 1 1 1 1 1 3、>> zeros(2) 0 0 0 0 >> zeros(4)ans = 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 >> zeros(2,2)ans = 0 0 0 0 >> zeros(2,3)ans = 0 0 0 0 0 0 >> zeros(3,2)ans = 0 0 0 0 00 4、随机阵 >> rand(2,3)ans = 0.2785 0.9575 0.1576 0.5469 0.9649 0.9706 >> rand(2,3) ans = 0.9572 0.8003 0.4218 0.4854 0.1419 0.9157 5、>> diag(5)ans = 5 >> diag(5,5)ans = 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 5 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 >> diag(2,3)ans = 0 0 0 2 0 0 0 0 0 0 0 0 0 0 0 0 6、(inv(A)为求 A 的逆矩阵)>> B=[5 3 1;2 3 8;1 1 1],inv(B)3 1 1 ans = 0.6250 0.2500-2.6250-0.7500-0.5000 4.7500 0.1250 0.2500-1.1250 >> A=[2 3;4 4],B=[5 3;3 8],inv(A),inv(B);AB,A/B,inv(A)*B,B*inv(A)A = 2 3 4 4 B = 5 3 3 8 ans =-1.0000 0.7500 1.0000-0.5000 ans =-2.7500 3.0000 3.5000-1.0000 ans = 0.2258 0.2903 0.6452 0.2581 ans =-2.7500 3.0000 3.5000 -1.0000 ans =-2.0000 2.2500 5.0000-1.7500 7、>> p =[1,-6,-72,-27], roots(p)p = 1-6-72-27 ans = 12.1229-5.7345-0.3884 >> p=[2,3,6],roots(p)p = 2 3 6 ans =-0.7500 + 1.5612i-0.7500-1.5612i 8、(A 为 n*n 的方阵)>> A=[0 1 0;-4 4 0;-2 1 2],poly(A),B=sym(A),poly(B)A = 0 1 0-4 4 0-2 1 2 ans =-6 12-8 B = [ 0, 1, 0] [-4, 4, 0] [-2, 1, 2] ans = x A 3-6*x A 2+12*x-8 9,、(conv 是多项式相乘，deconv 是多项式相除)>> u=[1 2 4 6 ],v=[5 0 0-6 7],conv(u,v)u = 1 2 4 6 v =0 0-6 7 ans =10 20 24-5-10-8 42 >> v=[1 2 4 6 ],u=[5 0 0-6 7],deconv(u,v)v = 1 2 4 6 u = 5 0 0-6 7 ans = 5-10 10、（点乘是数组的运算，没有点的乘是矩阵运算） >> a = [2 5;3 4], b =[3 1;4 7],a.*b,a*b a = 2 5 3 4 b = 3 1 4 7 ans = 6 5 12 28 ans = 26 37 25 31 >> a = [2 3];b = [4 7];a.*b = [8 21];a*b %错误 a*b“ = 29;11、（who 可以看到你用过的一些变量，来了） >> who Your variables are: A B a ans b p u >> whos Name Size Bytes whos 是把该变量及所存储的大小等信息都显示出 Class Attributes 2x2 32 double B 2x2 32 double a 1x2 16 double ans 1x2 16 double b 1x2 16 double p 1x3 24 double u 1x5 40 double v 1x4 32 double12、>> A=[2 5 3;6 5 4],disp(A),size(A),length(A)A = 2 5 3 6 5 4 2 5 3 6 5 4 ans = 2 3 ans = 3 实验二 MATLAB 绘图命令 实验任务 熟悉 MATLAB 基本绘图命令，掌握如下绘图方法：．坐标系的选择、图形的绘制；•图形注解（题目、标号、说明、分格线）的加入;3 •图形线型、符号、颜色的选取 基本命令训练 1、>>t=[0:pi/360:2*pi];x=cos(t)+ cos(t*4);y=si n(t)+ sin(t*4);xlabel(”x 轴“);ylabel(”y 轴“);plot(y,x),grid;2、>> t=0:0.1:100;x=3*t;y=4*t;z=si n(2*t);plot3(x,y,z, ”g:“) ■15 i 0 5 0 05 1 1 5 2 3、>>x = linspace(-2*pi,2*pi,40);y=si n(x);stairs(x,y)4、>> t=[0:pi/360:2*pi];x=cos(t)+ cos(t*4)+ sin(t*4);y=si n(t)+ si n(t*4);plot(y,x, ”r:“); xlabel(”x 轴“);ylabel(”y 轴“); 6、>>th=[0:pi/20:2*pi];x=exp(j*th);plot(real(x),imag(x),”r-.“);grid;text(0,0,”中心“); 5、>> th=[0:pi/1000:2*pi]”;r=cos(2*th);polar(th,r);title(“四叶草图”)270 四叶草图 7、>>x=-2:0.01:2;8、y=-2:0.01:2;9、[X,Y] = meshgrid(x,y);Z = Y.*exp(-X.A 2-Y.A 2);[C,h] = con tour(X, YZ);set(h,“ShowText”,“o n”,“TextStep”,get(h,“LevelStep”)*2)_ 1OS •I ,5 2 ■n.s o o.s 8、>>x = 0:0.2:10;y = 2*x+3;subplot(411);plot(x,y);grid;title(“y 的原函数”);subplot(412);semilogy(x,y);grid;title(“对 y 取对数”);丫的原画数 40----------1------------------------------------------1--------------------1--------------------1-----------------------------------------1--------------------1-------------------1-------- l| I p | il ■ | i| I 九 _____ 1-___ — ___ I ____ L ___ :…… ： ___ J _ _______ L ___ u i| I |l I , il _-■」 “ j I ■I __ h-_____________ I I ■ Q 」 【 I 1 F I I I I II I I IT 10 1 □ 1 2 3 4 5 6 r 6 9 10 , 对 y 取对数 对弋观对数 subplot(413);semilogx(x,y);10 10 id 10 1lZ 10 w 10 10 40 20 0 对好对数 •掌握循环、分支语句的编写，学会使用 look for、help 命令 grid;title(”对 x 取对数“);subplot(414);loglog(x,y);grid;title(”对 xy 均取对数“);9、>>x =-3:0.3:3;bar(x,exp(-x.*x),”g“)实验三 MATLAB 程序设计 实验任务 1 •熟悉 MATLAB 程序设计的方法和思路； 程序举例 1、>> f=[1,1];i=1;while f(i)+f(i+1)<1000 f(i+2)=f(i)+f(i+1);i=i+1;end f,i f = Columns 1 through 14 1 1 2 3 5 8 13 21 34 55 89 144 233 377 Columns 15 through 16 610 987 i = 15 2、>> m=3;n=4;for i=1:m for j=1:n a(i,j)=1/(i+j-1);end end format rat1/2 1/3 1/2 1/3 1/4 1/3 1/4 1/5(分数格式形式。用有理数逼近显示数据)>>m=5;n=4;for i=1:m for j=1:n a(i,j)=1/(i+j-1);end end format rat a a = 1 1/2 1/3 1/4 1/2 1/3 1/4 1/5 1/3 1/4 1/5 1/6 1/4 1/5 1/6 1/7 1/5 1/6 1/7 1/8 3、程序中没有 format rat 命令时，如果上次运行结果没有清除，输出的结果就是上次运行的 结果！但是运用 clear 命令清楚之前的运行结果之后就会正常运行。 4、>> x=input（”请输入 x 的值:’）;if x==10 y=cos(x+1)+sqrt(x*x+1);else y=x*sqrt(x+sqrt(x));1/4 1/5 1/6 end y 请输入 x 的值 :2 y = 2391/647 x=input(“请输入 x 的值:’);if x==10 y= fprintf(” 不在定义域内，请重新输入： “);return else y=1/(x-10);end y 请输入 x 的值 :2-1/8 5、>> p=[0 0 0 1 3 0 2 0 0 9];for i=1:length(p), if p(1)==0,p=p(2:length(p));end;end;p p = Columns 1 through 5 3 0 0 Columns 6 through 7 0 9 >>p=[0 0 0 1 3 0 2 0 0 9];p(p==0)=[];p p =3 6、>> e2(500) ans =1 2 3 5 8 13 21 34 55 89 144 233 377 >> lookfor ffibno e2-ffibno 计算斐波那契亚数列的函数文件 >> help e2 ffibno 计算斐波那契亚数列的函数文件 n 可取任意自然数 程序如下（用法： lookfor 关键词 在所有 M 文件中找 “关键词 ”，比如： lookfor max（即寻找关键词 “ max”） 其实就和我们平时用 CTRL+F 来查找关键词”是一样的 而 help 是显示 matlab 内置的帮助信息 用法： help 命令，比如 help inv，作用就是调用 inv 这个命令的帮助） 程序设计题 用一个 MATLAB 语言编写一个程序：输入一个自然数，判断它是否是素数，如果是，输出“ It is one prime ”，如果不是，输出“ It is n ot o ne prime.”。要求通 过调用子函数实现。 最好能具有如下功能： ①设计较好的人机对话界面，程序中 含有提示 性的输入输出语句。 ②能实现循环操作，由操作者输入相关命令来控制 是否继续进行素数的判断。 如果操作者希望停止这种判断，则可以退出程序。 ③ 如果所输入的自然数是一个合数，除了给出其不是素数的结论外，还应给出至少 一种其因数分解形式。例：输入 6，因为 6 不是素数。则程序中除了有“ It is not one prime ”的结论外，还应有：“ 6=2*3 ”的说明。 function sushu while 1 x=input（” 请输入一个自然数 “）;if x<2 disp（” 既不是质数又不是合数 “）;Else if isprime(x)==1 disp（” 这是一个素数 “）;Else disp(” 这是一个合数，可以因式分解为： “)for n=2:sqrt(x)if rem(x,n)==0 num3=x;num1=n;num2=x/n;disp([num2str(num3),”=“,num2str(num1),”x“,num2str(num2)])End End End end y=input(” 是否继续判断？继续请按 1 ,按任意键退出 ： “)； if y~=1 break end end 实验四 MATLAB 的符号计算与 SIMULINK 的使用 实验任务 1.掌握 MATLAB 符号计算的特点和常用基本命令;2.掌握 SIMULINK 的使用。 程序举例 1.求矩阵对应的行列式和特征根 >> a=sym(”[a11 a12;a21 a22]“);da=det(a)ea=eig(a)da = a11*a22-a12*a21 ea = 1/2*a11+1/2*a22+1/2*(a11 A 2-2*a11*a22+a22 A 2+4*a12*a21)A(1 /2)1/2*a11+1/2*a22-1/2*(a11 A 2-2*a11*a22+a22 A 2+4*a12*a21)A(1 /2)a=sym(”[2 3;1 5]“); da=det(a)ea=eig(a)da = 7 ea = 7/2+1/2*21A(1/2)7/2-1/2*21A(1/2)2.求方程的解(包括精确解和一定精度的解)>> r1=solve(”xA2+x-1“)rv=vpa(r1)rv4=vpa(r1,4) rv20=vpa(r1,20)r1 = 1/2*5 八(1/2)-1/2-1/2*5 八(1/2)-1/2 rv =.6894848260-1.689484826 rv4 =.6180-1.618 rv20 =.689484820-1.68948482 3 ． a=sym(”a“);b=sym(”b“);c=sym(w=10;x=5;y=-8;z=11;A=[a,b;c,d] B=[w,x;y,z] det(A)det(B)A = [ a, b] ”c“);d=sym(”d“);%定义 4 个符号变量 %定义 4 个数值变量 %建立符号矩阵 A %建立数值矩阵 B %计算符号矩阵 A 的行列式 %计算数值矩阵 B 的行列式 [ c, d] X 10-8 11 ans = a*d-b*c ans = 150 4.» syms x y;s=(-7*x A 2-8*y A 2)*(-x A 2+3*y A 2);expand(s)% 对 s 展开 collect(s,x)% 对 s 按变量 x 合并同类项（无同类项） factor(ans)% 对 ans 分解因式 ans = 7*xA4-13*xA2*yA2-24*yA4 ans = 7*xA4-13*xA2*yA2-24*yA4 ans =(8*yA2+7*xA2)*(xA2-3*yA2)5.对方程 AX=b 求解 >> A=[34,8,4;3,34,3;3,6,8];b=[4;6;2];X=linsolve(A,b)% 调用 linsolve 函数求解 Ab % 用另一种方法求解 0.0675 0.1614 diff(f)%未指定求导变量和阶数，按缺省规则处理 0.1037 ans = 0.0675 0.1614 0.1037 6 ． 对方程组求解 a11*x1+a12*x2+a13*x3=b1 a21*x1+a22*x2+a23*x3=b2 a31*x1+a32*x2+a33*x3=b3 >>syms a11 a12 a13 a21 a22 a23 a31 a32 a33 b1 b2 b3;A=[a11,a12,a13;a21,a22,a23;a31,a32,a33];b=[b1;b2;b3];XX=Ab %用左除运算求解(X=linsolve(A,b)%调用 linsolve 函数求的解)XX =(a12*a23*b3-a12*b2*a33+a13*a32*b2-a13*a22*b3+b1*a22*a33-b1*a32*a23)/(a11*a22*a33-a11*a32*a23-a12*a21*a33+a32*a21*a13-a22*a31*a13+a31*a12*a 23)-(a11*a23*b3-a11*b2*a33-a21*a13*b3-a23*a31*b1+b2*a31*a13+a21*b1*a33)/(a11*a22*a33-a11*a32*a23-a12*a21*a33+a32*a21*a13-a22*a31*a13+a31*a12* a23)(a32*a21*b1-a11*a32*b2+a11*a22*b3-a22*a31*b1-a12*a21*b3+a31*a12*b2)/(a11*a22*a33-a11*a32*a23-a12*a21*a33+a32*a21*a13-a22*a31*a13+a31*a12*a 23)7 ． syms a b t x y z;f=sqrt(1+exp(x)); %求 f 对 x 的二阶导数 %求 f 对 x 的三阶导数 %按参数方程求导公式求 y 对 x 的导数 ans = 1/2/(1+exp(x))A(1/2)*exp(x)ans =-2*si n(x)-x*cos(x)ans =-3*cos(x)+x*si n(x)ans =-b*cos(t)/a/si n(t)三、SIMULINK 的使用 f=x*cos(x);diff(f,x,2)diff(f,x,3)f1=a*cos(t);f2=b*si n(t);diff(f2)/diff(f1) diff(f)%未指定求导变量和阶数，按缺省规则处理 仿真图: 波形图: 其中: G i(R 0.14 实验五 MATLAB 在控制系统分析中的应用 实验任务 1.掌握 MATLAB 在控制系统时间响应分析中的应用;2.掌握 MATLAB 在系统根轨迹分析中的应用； 3.掌握 MATLAB 控制系统频率分析中的应用； 4.掌握 MATLAB 在控制系统稳定性分析中的应用 基本命令 1.step 2.impulse 3.in itial 4.Isim 5.rlocfi nd 6.bode 7.margin 8.nyquist 9.Nichols 10.cloop 程序举例 1.求下面系统的单位阶跃响应 1.5 G(s) t >> num=[4];den=[1 , 1 , 4];step(num , den)[y , x , t]=step(num , den);tp=spli ne(y , t , max(y))% 计算峰值时间 max(y)% 计算峰值 tp = Step Response 0.5 Time(sec)p m 1.6062 ans = 1.4441 0.18 0.16 2.求如下系统的单位阶跃响应 X i 0 1 X i 0 u X 2 6 5 X 2 1 >> a=[0,1;-6,-5];b=[0;1];c=[1,0];d=0;[y,x]=step(a,b,c,d);Plot(y) 1,0 X i X 2 3.求下面系统的单位脉冲响应: G(s)儿 >> num=[4];den=[1 , 1 ,4];impulse(nu m,de n)Time(sec)Response to Initial Conditions 4.已知二阶系统的状态方程为: e u p m A 0.25 0.2 0.15 0.1 0.05 0-0.05 L I 1- r -厂、、/ Impulse Response c=[1 , 0] ； d=[0];xO=[1 ,0] ； subplot(1 , 2,1);in itial(a , b , c ,d,x0)subplot(1 , 2,2);impulse(a , b , c , d):系统传递函数为: G(s)占 输入正弦信号时，观察输出 信号的相位差。 >> num=[1];den=[1 ,1];t=0 : 0.01 : 10;u=s in(2*t);hold on plot(t,u, ”r“)lsim(nu m,de n,u,t) Real Axis 6.有一二阶系统，求出周期为 4 秒的方波的输出响应 2s 2 5s 1 ~2 s 2s 3 >>num=[2 5 1];den=[1 2 3];t=(0:.1:10);period=4;u=(rem(t,period)>=period./2);% 看 rem 函数功能 lsim(nu m,de n,u,t);7.已知开环系统传递函数，绘制系统的根轨迹，并分析其稳定性 G(s)k(s 2)(s 2 4s 3)2 >>num=[1 2];den 仁 [1 4 3];den=conv(de n1,de n1);figure(1)-6 rlocus(num,den)[k,p]= rlocfi nd(nu m,de n)l L l ■ 匚一 L — C-------------------- --------- -r t r r -10-8-6-4-2 4 6 Root Locus 0 2 2 2 0 0 2 2--CPXA vyanma n-8 2.5 G(s)Lin ear Simulati on Results1.5 e 0.5-0.5-1-1.5-2 Time(sec)p m impulse response(k=55) 5 figure(2)k=55;num 仁 k*[1 2];den=[1 4 3];den 1=c onv(de n,den);[nu m,de n]=cloop(nu m1,de n1,-1);impulse(nu m,de n)title(”impulse resp onse(k=55)“)-1.5-1 0 200 400 600 Time(sec)800 1000 1200 5 5 0 05 0 0 o --x 10 8!------- impulse resp on se(k=56)figure(3)k=56;num 仁 k*[1 2];den=[1 4 3];den 1=c onv(de n,den);[nu m,de n]=cloop(nu m1,de n1,-1);impulse(nu m,de n)-2-4-6-8 L 0 500 1000 Time(sec)title(”impulse resp on se(k=56)“)1500 2000 2500 Select a point in the graphics win dow selected_po int =-2.5924-0.0248i 0.7133-3.4160-2.5918-0.9961 + 0.4306i-0.9961-0.4306i Bode Diagram 8.作如下系统的 bode 图 G(s)>> n=[1 , 1];d=[1 , 4,11 , 7];bode(n , d),grid on Frequency(rad/sec)s 1 s 3 4s 2 11s 7 9.系统传函如下 G(s)s 1 0.5s e(s 2)3 求有理传函的频率响应，然后在同一张图上绘出以四阶伯德近似表示的系统频 率响应 >> num=[1];de n=co nv([1 2],co nv([1 2],[1 2])); w=logspace(-1,2);t=0.5; ylabel(”gai n“);subplot(2,1,2);semilogx(w,p1,w,p2,”g--“);grid on;xlabel(”freque ncy“);ylabel(”phase“);10.已知系统模型为 求它的幅值裕度和相角裕度 >> n=[3.5];d=[1 2 3 2];[Gm,Pm,Wcg,Wcp]=margi n(n,d) G(s)3.5 s 3 2s 2 3s 2 bode plot [m1,p1]=bode(nu m,de n, 2);p1=p1-t*w”*180/pi;[n 2,d2]=pade(t,4);nu mt=c onv(n2,nu m);den t=(c onv(de n, d2));2 freque xlabel(“freque Gm = 1.1433 Pm = 7.1688 Wcg = 1.7323 Wcp = 1.6541 nyq uist(n, d1);hold on nyq uist(n, d2);nyq uist(n, d3);nyq uist(n, d4);Nyquist Diagram 80 |-------------------- --------------------[--------------------[--------- 11.二阶系统为: G(s)n s 2n S S 令 wn=1, 分别作出 E =2,1 , 0.707 , C A 5 y a 时的 nyquist 曲线。 m >> n=[1];d 仁 [1 , 4,1];d2=[1 , 2 , 1] ； d3=[1 , 1.414,1];d4=[1,1,1];Nyquist Diagram Real Axis 12.已知系统的开环传递函数为 _一-”—一 一-一_一_ ~-60 005 S S 2) S S32 S S v vr ra an nk ky ya a 14.一多环系统，其结构图如下，使用 Nyquist 频率曲线判断系统的稳定性。 16.7s(0.85s 1)(0.25s 1)(0.0625s 1)绘制系统的 Nyqusit 图, 并讨论系统的稳定性.>> G=tf(1000,co nv([1,3,2],[1,5]));nyquist(G);axis(”square“)13.分别由 w 的自动变量和人工变量作下列系统的 nyquisF 曲线： m G(s)1 s(s 1)>> n=[1];d=[1 , 1 ,0];nyquist(n ,d);% 自动变量 n=[1];d=[1 , 1 ,0];w=[0.5 : 0.1 : 3];nyq uist(n , d , w);% 人工变量-2-1-1.5-0.9-0.8-0.7-0.6-0.5 Real Axis-0.4-0.3-0.2-0.1 5 5 --0 0 S SI IX XA A u u^ ^a an nT T9 9a a 卩 G(s) 30 20 10 0 figure(2)[nu m2,de n2]=cloop(nu m,de n);impulse(nu m2,de n2); >> k1=16.70.0125;z1=[0];p1=[-1.25-4-16];[nu m1,de n1]=zp2tf(z1,p1,k1);[nu m,de n]=cloop(nu m1,de n1);[z,p,k]=tf2zp(nu m,de n);p figure(1)nyq uist(nu m,de n)-2-1 Nyquist Diagram-0.5 0 0.5 Real Axis 1 1.5 1 0.5 s a 0 n g m-0.5-1-1.5 1.5 20 Impulse Resp onse-10.5969 +36.2148i-10.5969-36.2148i-0.0562 15.已知系统为: eanML—pm -5-10-15 0 0.1 0.2 0.3 Time(sec)Nichols Chart 0.4 0.5 0.6-2 10 0 10 2 Open-Loop Phase(deg)Frequency(rad/sec) G(s)s(s 1)作该系统的 nichols 曲线。 >> n=[1];d=[1 , 1 , 0];ni chols(n , d);16.已知系统的开环传递函数为: G(s)k s(s 1)(s 2)当 k=2 时，分别作 nichols 曲线和波特图 >> num=1;den=conv(co nv([1 0],[1 1]),[0.5 1]);subplot(1,2,1);ni chols(nu m,de n);grid;B % n ichols 曲线 G subplot(1,2,2);2 g=tf(nu m,de n);bode(feedback(g,1,-1));grid;Nichols Chart koflr ea M^hnaaM msuvesaB— 50 -- 00 -270 Bode Diagram 90-- % 波特图 17.系统的开环传递函数为: 分别确定 k=2 和 k=10 时闭环系统的稳定性 >> d 仁 [1 , 3,2,0];n 仁 [2];[nc1 , dc1]=cloop(n1 , d1 ,-1);roots(dc1)d2=d1;n2=[10];[nc2 , dc2]=cloop(n2 , d2,-1);roots(dc2)ans =-2.5214-0.2393 + 0.8579i-0.2393-0.8579i ans =-3.3089 0.1545 + 1.7316i 0.1545-1.7316i 18.系统的状态方程为： X 1 4 3 0 X 1 1 X 2 1 0 0 X 2 0 u X 3 0 1 0 X 3 0 X 1 y 0 1 2 x 2 X 3 G(s)k s(s 1)(s 2) 试确定系统的稳定性。 >> a=[-4,-3,0;1,0,0;0,1,0];b=[1;0;0];c=[0,1,2];d=0;eig(a)% 求特征根 ran k(ctrb(a,b))ans = 0-1-3 ans = 3 实验六连续系统数字仿真的基本算法 实验任务 1.理解欧拉法和龙格-库塔法的基本思想； 2 •理解数值积分算法的计算精度、速度、稳定性与步长的关系;程序举例 1.取 h=0.2 ,试分别用欧拉法、RK2 法和 RK4 法求解微分方程的数值解，并 比较计算精度。 注：解析解：y(t).1 2t >> clear t(1)=0 ； y(1)=1;y_euler(1)=1;y_rk2(1)=1;y_rk4(1)=1;h=0.001;% 步长修改为 0.001 for k=1:5 y(t)y(t)y(0)2t y(t)1 t(k+1)=t(k)+h;y(k+1)=sqrt(1+2*t(k+1));end for k=1:5 y_euler(k+1)=y_euler(k)+h*(y_euler(k)-2*t(k)/y_euler(k));end for k=1:5 k1= y_rk2(k)-2*t(k)/y_rk2(k);k2=(y_rk2(k)+h*k1)-2*(t(k)+h)/(y_rk2(k)+h*k1);y_rk2(k+1)=y_rk2(k)+h*(k1+k2)/2;end for k=1:5 k1= y_rk4(k)-2*t(k)/y_rk4(k);k2=(y_rk4(k)+h*k1/2)-2*(t(k)+h/2)/(y_rk4(k)+h*k1/2);k3=(y_rk4(k)+h*k22)-2*(t(k)+h/2)/(y_rk4(k)+h*k2/2);k4=(y_rk4(k)+h*k3)-2*(t(k)+h)/(y_rk4(k)+h*k3); y_rk4(k+1)=y_rk4(k)+h*(k1+2*k2+2*k3+k4)/6;end disp(” 时间 解析解 yt=[t“, y”, y_euler“, y_rk2”, y_rk4“];disp(yt)欧拉法 RK2 法 RK4 法 ”)时间 解析解 欧拉法 RK2 法 RK4法0 1.0000 1.0000 1.0000 1.0000 0.0010 1.0010 1.0010 1.0010 1.0010 0.0020 1.0020 1.0020 1.0020 1.0020 0.0030 1.0030 1.0030 1.0030 1.0030 0.0040 1.0040 1.0040 1.0040 1.0040 0.0050 1.0050 1.0050 1.0050 1.0050 y(t)2Ry(t)y(t)0 在 0 t 10 上的数字仿真解(已知： y(0)y(0)0)，并将不同步长下的仿真结果与解析解进行精度比较。 说明： 已知该微分方程的解析解分别为： 100，yt y t 100cost(当 R 0)100e 2t c°s 仝 t 10 ^ e*s in 仝 t 2 3 2(当 R 0.5)采用 RK4 法进行计算，选择状态变量: 2.考虑如下二阶系统： x 1 y x 2 y 则有如下状态空间模型及初值条件 x 1 x 2 x 1(0)100 x 2 x 1 2Rx 2 x 2(0)0 y x 1 采用 RK4 法进行计算。 >> clear h=input(“ 请输入步长 h=”);M=round(10/h);t(1)=0;y_0(1)=100;y_05(1)=100;和 y_05 分别对应于为 R=0 和 R=0.5)x1(1)=100;x2(1)=0;y_rk4_0(1)=x1(1);y_rk4_05(1)=x1(1);% 求解析解 for k=1:M t(k+1)=t(k)+h;y_0(k+1)=100*cos(t(k+1));y_05(k+1)=100*exp(-t(k+1)/2).*cos(sqrt(3)/2*t(k+1))+100*sqrt(3)/3*exp(-t(k+1)/2).*si n(sqrt(3)/2*t(k+1));end% 输入步长 % 置总计算步数 % 置自变量初值 % 置解析解的初始值(y_0 % 置状态向量初值 % 置数值解的初值 % 利用 RK4 法求解 % R=0 for k=1:M k11=x2(k);k12=-x1(k);k21=x2(k)+h*k12/2;k22=-(x1(k)+h*k11/2);k31=x2(k)+h*k22/2;k32=-(x1(k)+h*k21/2);k41=x2(k)+h*k32;k42=-(x1(k)+h*k31);x1(k+1)=x1(k)+h*(k11+2*k21+2*k31+k41)/6;x2(k+1)=x2(k)+h*(k12+2*k22+2*k32+k42)/6;y_rk4_0(k+1)=x1(k+1);end % R=0.5 for k=1:M k11=x2(k);k12=-x1(k)-x2(k);k21=x2(k)+h*k12/2;k22=-(x1(k)+h*k11/2)-(x2(k)+h*k12/2);k31=x2(k)+h*k22/2;k32=-(x1(k)+h*k21/2)-(x2(k)+h*k22/2);k41=x2(k)+h*k32;k42=-(x1(k)+h*k31)-(x2(k)+h*k32);x1(k+1)=x1(k)+h*(k11+2*k21+2*k31+k41)/6;x2(k+1)=x2(k)+h*(k12+2*k22+2*k32+k42)/6;y_rk4_05(k+1)=x1(k+1);end % 求出误差最大值 err_0=max(abs(y_0-y_rk4_0));err_05=max(abs(y_05-y_rk4_05));% 输出结果 disp（“ 最大误差（R=0） 最大误差（R=0.5）”） err_max=[err_0,err_05];disp（err_max） 请输入步长 h=0.5 最大误差（R=0） 最大误差（R=0.5） 0.4299 0.0460 没运行一次程序，输入一个步长，记得可到相应的最大误差，将结果可列表如 下： 步长 h 0.0001 0.0005 0.001 0.005 0.01 0.05 0.1 0.5 F =0 5.4330 1.6969 1.0574 4.1107 6.6029 4.1439 6.6602 4.2988 差 最大误 X 10-10 X 10-10 X 10-10 X 10-9 X 10-8 X 10-5 X 10-4 X 10-1 13工业工程2班 李伟航 13工业工程2班 实验10 一、实验目的: 1.学习库存系统查库与订货处理的结构建模方法 2.学习用Equation模块、Equation（I）模块读写数据库的方法 3.学习用Equation模块、Equation（I）模块进行编程计算的方法 二、实验问题 1.打开上次实验你保存的文件（这个文件要保存好，下次实验还要使用），然后根据以上视频，进行操作实验。 2.简述用Equation模块计算订货量的程序逻辑。 3.简述用Equation（I）模块计算并累加总订货成本到数据库中的程序逻辑。4.Equation模块与Equation（I）模块有何不同？ 5.在本案例的假设前提下，在一笔订货的在途货物运输期间（即提前期期间），会不会再次发出订货指令？或者换句话说，每次查库时，会不会有已订但未到的货？为什么？ 三.实验过程 1.生成查库员（查库信号） 用Create模块每天生成一个库存检查员实体（实际代表一个查库信号）Create模块具体设置如下图： 2.判断是否需要订货 利用Select Item Out模块、Equation模块和Simulation Variable模块检查库存，并判断是否需要订货。若需要，就将库存检查员实体发送到Select Item Out模块的上端口输出进行后续处理；若不需要订货，就将库存检查员实体发送到下端口输出，简单地离开系统。其中，Equation模块的设置如下图。 3.无需订货的处理 由上一步Equation中设置可知，当s=1时，即无需订货的情况下，直接将库存检查员实体从Select Item Out模块下端口输出，通过Exit模块离开系统。Select Item Out模块设置如下，当s=0时从上端口输出，否则从下端输出。 4.订货处理-建立模型 当s=0时，即需要订货的情况下，库存检查员实体有Select Item Out模块上端口输出，后续订货处理模型如下图： 5.计算订货量 用一个Equation模块获取数据库中的当前库存（kc）和最大库存(ds)，计算订货量dh。Equation模块设置如下： 6.设置订货量属性 将第5步Equation模块计算得到的订货量（dh）输出给Set模块，将订货量赋值给检查实体dh属性。Set属性设置如下图。 7.用Equation（I）模块计算并累加总订货成本到数据库中 因此用Equation（I）模块获取数据库中固定费（gdf）和可变订货费（kbf），再根据订货量算出本次订货成本，把它累加到总订货费用（toc）并写入数据库中。具体Equation（I）模块设置如下图。 8.设置提前期 从下订单到所订货物入库的这段时间为订货提前期，这里用Activity模块表示，设置该提前期服从0.5-1天的均匀分布。订货提前期的设置如下。 9.到货后更新数据库中的库存数据 库存检查实体离开Activity，表示提前期结束，货物到达，到货后用Get模块获取订货量，用Write模块将订货量累加到数据库中的sc表的当前库存（kc）中，根据订货动态增加当前库存。Write模块设置如下。 10.观察运行结果 运行模型，观察数据库sc表中的库存（kc）由变为13，累计订货费用达到了11125。 四.实验问题的分析解答 1.上述过程即为实验步骤。本实验实验模型如下： 2.答：用Equation模块计算订货量的程序逻辑：是当检测到库存小于订货点xs时，用最大库存ds减去目前库存，即为订货量，即dh=ds-kc；过程如实验步骤5所示。 3.答：用Equation（I）模块计算并累加总订货成本到数据库中的程序逻辑是：因为每订货一次的成本=固定费（gdf）+可变订货费（kbf）*订货量（dh），将每次的订货成本累加到总订货成本当中去，即toc=toc+gdf+kbf*dh，即更新了总订货成本。具体步骤如步骤7所示。 4.答：Equation模块是由value库中添加的值模块，是用来读取传递系统产生的值（如产生的随机数）并进行相应的逻辑运算，而对于实体传来的属性则不可以读取；而相反的，Equation（I）模块是由Item库中添加的实物模块，用来读取和传递实体属性并进行相应的逻辑运算。两者都可以读取数据可中的数据。 5.答：一笔订货的在途货物运输期间，不会再次发出订货指令。因为在Create模块中设置了检查实体到来的间隔为1天，而订货提前期服从0.5-1的均匀分布，所以不会出现每次查库有已订但未到的货的情况。 实验11 一、实验目的 1.通过实验理解各种库存性能指标的含义 2.学习用Read模块读取数据库数据的方法 3.学习用Max&Min模块、Mean模块、Equation模块以及各种计算模块计算和采集库存系统性能指标的方法 二.实验问题 1.打开上次实验你保存的文件（这个文件要保存好，下次实验还要使用），然后根据以上视频，进行操作实验。 2.本实验中，用read模块读数据库数据时，在其option页要做何设置？ 3.用Mean Variance模块计算平均每天总成本的均值和置信区间时，模块对话框要做何设置？ 4.利用本次实验建立的模型，实验比较以下各组（s，S）下的平均每天总成本：(20, 40)(20, 60)(20, 80)(20, 100)(40, 60)(40, 80)(40, 100)(60,80)(60,100)(80,100)对每种情况重复运行20次，写出各项配置下的平均每天总成本的均值和置信区间。并写出以上最优的（使得平均每天总成本最低）的（s，S）。 5.添加必要的模块，分别计算平均每天缺货成本、平均每天储存成本、平均每天订货成本的均值和置信区间。6.用Plotter，Discrete Event模块绘制当前库存水平的波动曲线，并同时绘制一条高度为20（即订购点）的水平直线，和一条高度为0的水平直线（连接constant模块到plotter DE模块），观察并库存曲线的波动情况以及和两条直线的关系，根据你的观察，缺货情况经常发生吗？ 三.实验过程 1.在库存发生变化时读取数据库中当前库存水平 每当库存数据发生变化时，用Read模块读取库存（kc）的值，Read设置如下： 2.计算平均每天储存成本 当库存大于0时，通过max模块与0比较得到当前库存值，并通过Mean&Varience模块（运行10次）计算平均每天的库存，再用Equation模块计算平均每天储存成本。计算平均每天储存成本的模型与设置如下所示： 3.计算平均每天缺货成本 当库存小于0时，通过max模块与0比较得到当前缺货数，并通过Mean&Varience模块（运行10次）计算平均每天的缺货数，用Equation模块（缺货数取反）计算平均每天缺货成本。计算平均每天缺货成本的模型与设置如下所示： 4.计算平均每天订货成本 利用Read模块读取运行完后数据库中的总订货成本，再用Math模块除以系统运行当前时间，得到平均每天订货成本。平均每天订货成本模型和设置如下： 5.计算平均每天总成本 利用Math模块将第2、3、4步所计算出的平均每天储存成本，平均每天缺货成本和平均每天订货成本相加即得到平均每天总成本。 6.计算平均每天总成本的均值和置信区间 在Math模块后用Mean&Varience模块计算运行10次之后的平均每天总成本的均值和置信区间。 四.实验问题的分析解答 1.上述过程即为实验步骤。本实验实验模型如下： 2.答：由于本实验要求在库存发生变化时读取数据库的库存水平，因此在Read模块中的Option选项中，选择Discrete event-read data during when 以及勾中data sources change，意思即为数据发生变化时读取数据库。设置如下： 3.答：用Mean&Varience模块计算平均每天总成本的均值和置信区间时，要勾中Calculate for multiplte simulations，意思是计算出运行多次情况（这里为10）次时平均每天总成本的均值和置信区间。 4.答： 10种情况的运行结果如下，其中最优的（s，S）为（20,60）。 （20，40）（20，60）（20，80） 均值：125.8870元 均值：120.3812元 均值：121.4371元 置信区间（125.8870±2.2602）置信区间（120.3812±1.8750）置信区间（121.4371±1.0904） （20，100）（40,60）（40，80） 均值：128.1696元 均值：127.3233元 均值：127.5943元 置信区间（128.1696±1.1935）置信区间（127.3233±1.2890）置信区间（127.5943±0.8358） （40,100）（60,80）（60，100） 均值：134.0826元 均值：146.2447元 均值：147.6999元 置信区间（134.0826±1.2611）置信区间（146.2447±1.3544）置信区间（147.6999±1.1276） （80,100） 均值：167.5711元 置信区间（167.5711±1.2291） 5.答：用3个Mean&Varience模块计算平均每天存储成本、缺货成本和订货成本的均值和置信区间，Mean&Varience模块勾中Calculate for multiplte simulations。添加模块后的模型以及运行结果如下： 平均每天存储成本：平均每天缺货成本 均值约为9.5140元 均值约为16.2608元 置信区间为（9.5140±0.2234）置信区间为（16.2608±0.9608） 平均每天订货成本： 均值约为98.7903元 置信区间为（98.7903±1.3140） 6.答： Plotter Discrete Event模块设置和运行后的库存水平波动曲线如下图所示，由曲线图可以看出库存水平在0和20之间以及上下变动。由图可以看出，位于0以下的曲线分布比较多，这意味着缺货的情况发生的比较频繁；曲线在0和20之间的分布也较多，这意味着需要订货的情况也发生的比较多。 实验12 一.实验目的: 1.学习库存仿真优化方法 2.学习如何设置模块表格的内容为决策变量 二.实验问题 1.打开上次实验你保存的文件，然后根据以上视频，建立优化模型。 2.在优化模块（Optimizer）中，整数型决策变量和连续型（实数型）的输入方法有何不同？ 3.请在你的模型上，使用Extendsim软件提供的优化器来寻找（s，S）的最佳设置（总成本最小）。令s在1和99之间取值（步长为1，即为整数），S在2和100之间取值（步长为1，即为整数）。要注意s和S必须是整数而且满足s < S。a.优化器参数由先选择Quicker Defaults，写出优化结果（s，S）和平均每天总成本。 b.再将优化器参数由选择Better Defaults，写出优化结果（s，S）和平均每天总成本。 4.在上题（第2题）的基础上，通过将库存检查间隔（Evaluation Interval，目前为1天）作为变量加入到优化变量集中，来研究在每天开始工作时查看并补充（当需要时）库存是否是最佳方案，让该值在半天到5天之间连续取值，s和S的取值情况与上题相同。应用优化器求取最优设置。 a.先优化器参数选择Quicker Defaults，写出优化结果（s，S）、Evaluation Interval和平均每天总成本。 b.再将优化器参数选择Better Defaults，写出优化结果（s，S）、和平均每天总成本。 三.实验过程 1.用Data Init模块初始化最小最大库存(s，S) 2.设置Optimizer模块目标函数和决策变量 从value库中将Optimizer模块放进模型中，在查库与订货处理模型中将决策变量订货点(xs，1~99)和最大库存(ds，2~100)和输出变量平均每天总成本（toc）克隆拖放到Optimizer模块上，然后建立目标函数方程和约束条件。Optimizer模块设置如下： 3.设置Optimizer模块约束方程 在Optimizer模块中的Constraints中添加约束方程，约束条件为s和S必须是整数而且满足s < S，约束条件的代码设置如下： 4.设置Optimizer模块运行参数 在Optimizer模块下的Run Parameters中，由于模型为随机模型，所以单机Random mode下的Quicker Defaults按钮快速设置所有优化参数（速度快精度低），然后点击New Run。运行结束后选择Better Default按钮设置规模更大的优化参数（耗时长精度高），并比较两次的结果。 5.运行优化,查看结果 通过Optimizer模块下的Results可以看MinCost数值的变化，运行结束时最顶行会给出最优解。分别用Quicker Defaults参数和Better Defaults参数运行两次比较结果。结果见实验结论3。 四.实验问题的分析解答 1.实验过程如上述实验步骤所示。 2.答：整数型决策变量的输入方法为输入范围的时候不要输入小数点，如1；而实数型的输入方法为输入范围时输入带有小数点的数字，如1.0。 3.答：选择Quicker Defaults参数，（s，S）优化结果为（27,43），平均每天总成本约为122.37627元。 Quicker Defaults 选择Better Defaults参数，（s，S）优化结果为（20,58），平均每天总成本约为130.8615元。 Better Defaults 4.答：将库存检查间隔ei克隆拖入Optimizer模块上，在订货模型create模块后加队列防止其堵塞，ei的范围设置如下。分别运行Quicker Defaults参数Better Defaults进行优化。 Quicker Defaults做法 第一步：选择Quicker Defaults参数,并开始运行 第二步：分析结果。 （s，S）优化结果为（9,74），Evaluation Interval（ei）为1.0186，平均每天总成本约为123.0515元。 Quicker Defaults Better Defaults做法 第一步：选择Better Defaults参数 第二步：分析结果。（s，S）优化结果为（38,54），Evaluation Interval（ei）为1.0714，平均每天总成本约为,120.19126元。 关键词：交通仿真,交通信号控制,VC++程序设计,ADO技术 0引言 城市道路交叉口经常成为路网容量的瓶颈,人们出行的大部分交通延误是由于交叉口的通行能力不足造成的。为此,就需要建立一个交叉口信号控制仿真系统来研究改善这种现状的方法。交叉口信号控制仿真系统可以模拟交通中出现的各种问题,帮助人们了解交通系统特性,通过优化系统设计方案来改善交叉口的通行能力。本文以多个交叉口作为研究对象,通过VC技术构建一个多交叉口信号控制实时动态仿真系统。 1基于VC技术的道路交叉口仿真模型的建立 1.1 模型的整体构造 在建模过程中,把整个路网看成是由若干车道、交通标线、交叉口、信号灯、车辆等实体组成的。这些实体按性质可抽象成3种不同类型的对象:动态对象、静态对象和数据管理对象。动态对象在系统中受到其它因素的影响和制约,随时发生变化,如车辆和控制信号;静态对象在一次仿真运行开始后不再发生变化,如道路和交通规则等;数据管理对象,如路网参数的调入和仿真结果的保存等。因此,仿真模型就必须由以下几个类组成:静态模型类、动态模型类、数据管理类。其中,静态模型类包括交叉口类(CSingleIntersection)、车道类(CSingleLane)、 信号灯类(CSingleLight)等;动态模型类即车辆类(CVehicle);数据管理类则采用ADO类型库实现与数据库的交互。ADO(Active Data Object)是一种基于OLE DB的访问接口,它继承了OLE DB的优点,对OLE DB接口作了封装,定义ADO对象可以使程序开发得到简化。 1.2 仿真模型的建立 1.2.1 交叉口模型 交叉口模型由交叉口、车道、交通标线、信号灯等静态对象组成,其中交叉口类型、车道数、路段的连接等参数均由数据库中的信息决定,并通过一个路口参数设置的对话框来进行数据调入,此模型为静态模型,从仿真初始化以后就不再变化。 1.2.2 车辆产生模型 车辆的到达在某种程度上具有随机性,因此,车辆的产生应该用概率分布来描述。由于本模型是在不超车和车流量较低的情况下进行的仿真,故选用移位负指数分布来描述发车规律,其车头时距h大于每个计数间隔持续时间t概率的基本计算公式为: P(h>t)=exp[-λ(t-τ)] (t≥τ) 。 (1) 式中:λ——车流的平均到达率,辆/s; τ ——车辆之间保持的最小车头时距,s。 1.2.3 车辆跟驰模型 跟驰模型是交通系统仿真中最重要的动态模型,跟驰模型的实质是一种“刺激-反应”形式,它是根据驾驶员接受某种刺激后作出的反应来分析和研究车辆列队行驶时发生的跟驰现象,其表现形式为:反应=灵敏度×刺激。其基本表达式为: undefined。 (2) 其中:Xn(t)为前车在t 时刻的位置;Xn+1(t)为跟驰车在t 时刻的位置;undefined为t时刻前车速度;undefined为t时刻跟驰车速度;undefined为t1时刻跟驰车速度;undefined为t1时刻跟驰车加速度;α为比例常数,一般为路段畅行速度。 2基于VC的道路交叉口仿真模型的实现 2.1 路网数据获得 在开始仿真之前,把存储在数据库中的相关路网参数一一读入程序,这就需要把数据库和仿真程序连接起来。本程序以SQL Server 2000为数据源,利用ADO类型库中Connection对象和数据库建立连接,用Recordset对象获取数据库中某表数据。例如要从数据库中获取某十字路口的车道宽度,把它赋值到参数设置对话框相对应的变量中,故可写为: 在加入以上程序段之前,要用import指令来引入ADO类型库,在stdafx.h中加入如下语句: 2.2 动态模型仿真程序表达 动态仿真作为整个仿真过程的主体,不仅可再现路网的整体交通状况,还可以记录路网中的每个参与者每一时刻的信息。因此,对动态模型进行仿真编程就显得尤为重要。 2.2.1 信号灯状态程序表达式 在仿真过程中,由于车辆始终根据信号灯的颜色状态而启停,所以在信号灯类中应定义一个结构体LIGHTSTATE,记录每个小灯的当前状态。此结构体包括成员变量{COLORREF clrRGB;UINT time;}。现以单路口二相位配时方案为例来分析东、西方向信号灯在一个仿真周期内的状态变化。首先定义一个结构体数组LIGHTSTATE m_PhaseState[3],用来存储信号灯一周期内的所有变化状态,即控制方案的设置;通过这个控制方案可以判断信号灯在一个仿真周期内每一时刻的状态,随着时间段的不同信号灯的颜色也会随之变化。以下是东、西方向信号灯在一周期内状态变化的部分代码: 2.2.2 车辆程序表达式 将路网中的所有车流分成车道上的车流和交叉口处的车流两种类型,用链表形式表示。车辆程序表达式实现过程如下:类CVehicle定义为双向链表,其中包括{CVehicle *pprev;CVehicle *pnext},这是为了灵活地遍历每条车流中的汽车。由于车流被分为车道上的和交叉口处的,故应定义两个链表{CTypedPtrArray 2.3 程序仿真流程 程序仿真流程如下: (1)从存有路网参数(包括交叉口参数、交通量数据、相位配时、信号灯参数等)的数据库中读取数据,传到程序相应的参数对话框中,实现对仿真模型的初始化,十字交叉口参数设置对话框见图1,其中编辑框、列表框的数据均从数据库中读取。 (2)在给定的交通需求情况下,按照一定的随机分布使车辆在其各自的产生点一辆一辆地进入路网,然后用时间扫描法描述每一时段内每一辆车在路网上行驶的情况,直至车辆到达其目的地后从路网上消失,具体仿真界面见图2。 (3)查看仿真结果并将其写回到数据库的交通状况表中,以备将来作各种方案的对比。 3结束语 相同条件多次运行,均可获得分布比较理想的交通流画面,其性能评价指标如交叉口通行能力、饱和度、平均速度、平均延误、停车率等都与实际情况有较好的吻合。但由于仿真模型的随机性和理想化,使得仿真与实际交通环境还有相当的差距,所以本系统还有待于进一步的完善和改进。 参考文献 [1]万忠.双周期信号配时法在单路口实时控制中的应用[J].计算机工程与应用,2006,42(26):183-185. [2]钟邦秀.面向对象微观交通仿真系统的研究与实现[J].系统仿真学报,2002,14(4):418-421. 摘要： 分析了双回路制动阀的动态工作过程，探讨了双回路制动阀上、下阀芯对其前、后桥输出口的遮盖量对动态工作特性的影响.基于AMESim液压/机械多场耦合仿真平台建立了双回路制动阀的仿真模型，研究了遮盖量变化对制动压力输出特性的影响规律以及单回路制动安全性能.搭建了全液压制动系统的实验台架，对具有不同遮盖量的制动阀样品的制动性能进行了实验对比测试.实测结果表明：遮盖量与制动空行程成正相关，与前、后桥的最大制动压力成负相关；双回路制动阀的前、后回路相互独立，当其中一条回路失效时，另一条回路仍能正常工作.实验结果与仿真结果具有良好的一致性，验证了该仿真模型的有效性. 关键词： 全液压制动系统； 双回路制动阀； 遮盖量； AMESim仿真； 制动压力输出特性 中图分类号： TH 137.52+1文献标志码： A Abstract： The dynamic movement process of dual circuit braking valve was analyzed.The influence of the overlapping lengths between the upper core and the rear output port as well as the lower valve core and the front output port on its dynamic characteristics was investigated.According to the simulation model of dual circuit braking valve built by AMESim software，the effects of overlapping length on the braking output pressure characteristics and the safety of single circuit brake valve were studied.The experimental bench with a full hydraulic braking system was established to determine the braking valves output characteristics with different overlapping lengths.The experimental results showed that there were positive correlations between overlapping lengths and the length of braking spare travelling process.And the overlapping lengths were negatively correlated with the maximum output pressure.Each circuit could work independently.Therefore，when one circuit lost efficacy，the other could work well.The experimental results were consistent with the simulation results，which verified the great validity of the simulation model and its simulation method. Keywords： full hydraulic braking system； dual circuit hydraulic braking valve； overlapping length； AMESim simulation；brake pressure output characteristics 全液压制动系统具有制动力矩大、制动灵敏、稳定性强的特点，与气液制动系统相比具有显著的性能优势，已成为工程机械制动系统的主流发展趋势.由于目前国内并未完全掌握全液压制动系统的核心技术，仍依靠国外进口的系统元件，造成系统成本昂贵、升级困难.因此，实现全液压制动系统元件的国产化非常重要.双回路制动阀作为全液压制动系统压力输出的控制元件，其性能的优劣直接影响车辆制动的安全性.目前国内学者只对国外双回路制动阀样品的动、静态特性进行了仿真与实验研究[1-4]，尚未深入开展不同关键结构参数对双回路制动阀性能影响规律的研究. 本文建立了双回路制动阀工作过程的数学模型，分析了阀芯受力平衡以及动态运动过程，并基于AMESim液压/机械多场耦合仿真平台，建立了制动阀的液压仿真模型，通过变参数仿真和实验研究了关键结构参数对制动阀性能的影响规律.研究对于双回路制动阀的国产化开发具有一定的参考价值. 1双回路制动阀工作原理 1.1工作过程分析 双回路制动阀的液压原理图如图1所示，其中：A1、A2分别为前、后桥输出口；T1、T2分别为前、后桥回油口；P1、P2分别为前、后桥进油口；F为踏板力，N.双回路制动阀内部结构图如图2所示，其中：Xu为上阀芯对P1的遮盖量；Xd为下阀芯对P2的遮盖量. 双回路制动阀具有上、下两个阀芯（近制动踏板端为上阀芯，远制动踏板端为下阀芯），其工作过程可分为三个阶段：制动空行程阶段、比例制动输出阶段、制动撤销阶段. 制动空行程阶段是指上、下阀芯在F作用下左移，制动阀的前、后桥输出口A1、A2分别与回油口T1、T2由接通到逐渐断开，并与进油口P1、P2由断开到临界接通的过程. 比例制动输出阶段是指随着阀芯位移逐渐增大，A1、A2分别开始与P1、P2接通，双回路制动阀开始输出制动压力，制动输出压力与阀芯位移成比例关系.当阀芯位移达到最大值时，制动输出压力达到最大值并保持不变. 当F撤销时，阀芯在复位弹簧的作用下右移，双回路制动阀A1、A2分别与T1、T2接通，前、后桥制动器中的液压油回流至油箱，制动过程结束. 根据阀芯在不同工作阶段的受力情况，可得到其受力平衡方程.当制动阀处于空行程阶段时，阀芯受力平衡情况及动态运动方程可表示为 综上可得，上、下阀芯分别对进油口P1、P2的遮盖量Xu、Xd决定了制动过程中空行程的长短.此外，由式（4）可知，在双回路制动阀其他结构参数不变的情况下，Xu、Xd与双回路制动阀输出压力成比例关系.由于在双回路制动阀加工中主要是通过改变Xu、Xd调节不同的双回路制动阀的输出压力，因此，Xu、Xd是双回路制动阀的关键结构参数. 2双回路制动阀性能仿真 基于AMESim液压/机械多场耦合仿真平台建立的双回路制动阀仿真模型[5-6]如图3所示.为了保证该仿真模型能正确进行运算，在模型中增加了蓄能器及前、后制动器[7].主要仿真参数如表1所示. 2.1制动压力输出特性 双回路制动阀仿真模型中阀芯位移输入信号如图4（a）所示，阀芯位移变化范围为0～9 mm.双回路制动阀前、后桥输出压力仿真结果如图4（b）所示，双回路制动阀的空行程为2 mm，满足制动输出压力响应速度的要求.此后，经过一个短暂的压力脉动阶段，前、后桥制动输出压力与阀芯位移成正相关，随着阀芯位移的增大其输出压力相应增大.由于上、下阀芯运动存在不同步性，因此在输出压力比例增大阶段，后桥输出压力略小于前桥.当阀芯位移达到最大值时，双回路制动阀前、后桥的制动输出压力也达到最大值（10.4 MPa）并保持不变. 2.2制动阀特性变参数仿真 利用AMESim仿真平台中的批处理模块对制动阀上、下阀芯的Xu、Xd进行了变参数仿真分析（阀芯位移输入信号同2.1节），以确定Xu、Xd与前、后桥制动输出压力之间的关系.Xu、Xd仿真参数如表2所示. 图5为不同Xu、Xd时的仿真结果.随着Xu不断减小，双回路制动阀前桥输出口A1的最大输出压力逐渐增大；随着Xd不断减小，双回路制动阀后桥输出口A2的最大输出压力也逐渐增大.随着Xu、Xd逐渐减小，双回路制动阀前、后桥的制动空行程变小.因此，可通过改变Xu、Xd的大小改变双回路制动阀的空行程长短和前、后桥的最大输出压力. 2.3单回路安全性仿真 将双回路制动阀仿真模型的其中一条回路断开，对另一条回路的输出压力特性进行了仿真.仿真结果如图6所示.当一条回路失效时，另一条回路的制动空行程以及前、后桥最大输出压力与未失效时基本一致，即能正常实施制动，两条回路相互独立，从而提高了制动安全性. 3实验结果分析 为了对双回路制动阀的工作特性进行分析，搭建了全液压制动系统的实验台架，如图7所示.实验台架主要元件参数如表3所示. 3.1制动压力输出特性 利用由舵机、摇臂及连杆等组成的制动踏板控制机构对踏板进行匀速控制.制动踏板控制机构如图8所示.双回路制动阀的前、后桥制动输出压力的实验结果如图9所示，制动空行程约为2 mm.由图9 4结论 （1） 通过对双回路制动阀阀芯的动态运动方程进行分析，得出双回路制动阀上、下阀芯对其前、后桥输出口的遮盖量为其关键结构参数. （2） 对双回路制动阀上、下阀芯对其前、后桥输出口的遮盖量进行了变参数仿真，发现该参数与双回路制动阀的空行程成正相关，与最大输出压力成负相关；此外，对双回路制动阀的单回路安全性进行了仿真，结果显示，双回路制动阀在其中一条回路失效的情况下，另一条回路仍能正常输出压力，表明双回路制动阀具有较高的安全性. （3） 通过搭建全液压制动系统的实验台架，对双回路制动阀的输出压力特性、单回路安全性进行了实验，并对比了不同遮盖量时的制动输出压力特性. 实验结果与仿真结果具有较好的一致性，说明建立的仿真模型能对双回路制动阀的制动特性进行较准确的模拟计算.研究为双回路制动阀的产品开发尤其是结构优化提供了一定的理论参考. 参考文献： [1]程振东，田晋跃，刘刚.工程车辆全液压制动系统性能分析及仿真研究[J].工程机械，2005，36（12）：50-53. [2]陈晋市，刘昕晖，王同建，等.全液压制动系统液压制动阀的动态特性[J].哈尔滨工业大学学报，2013，45（5）：75-79. [3]张锦.全动力液压制动系统动态响应特性研究[D].太原：太原科技大学，2008. [4]ZHAN C，LIU X H，WANG Z，et al.Analysis of the regenerative brake system parameters for concrete mixing truck basded on AMESim[J].Electronic and Mechanical Engineering and Information Technology，2011：2124-2127. [5]余佑官，龚国芳，胡国良.AMESim仿真技术及其在液压系统中的应用[J].液压气动与密封，2005（3）：28-31. [6]程振东，田晋跃，刘刚.车辆全液压制动系统执行机构建模及仿真[J].系统仿真学报，2006，18（3）：778-780. [7]王磊.皮囊式蓄能器在液压系统中的应用[J].能源研究与信息，2005，21（2）：100-105. [8]王展.全液压制动系统仿真分析与实验研究[D].长春：吉林大学，2012. 【信号与系统仿真实验报告】推荐阅读： 实验报告信号与系统06-29 信号与系统实验二05-18 信号与系统综合实验05-28 民大信号与系统实验一04-24 信号与系统第二次实验06-16 信号与系统实验总结及心得体会07-24 信号与系统08-19 信号与系统作业合集08-13 信号与系统系列课程06-02信号与系统仿真实验报告 篇2
信号与系统仿真实验报告 篇3
信号与系统仿真实验报告 篇4
中南大学系统仿真实验报告 篇5
仿真与建模实验报告 篇6
信号与系统仿真实验报告 篇7
信号与系统仿真实验报告 篇8