计算方法大作业(西电机电院)

计算方法大作业第一部分、计算方法学习回顾1.主要知识点：误差分析，插值法与拟合，数值积分，数值微分，线性方程组的直接解法和迭代解法，非线性方程求根，矩阵特征值问题计算、常微分方程初值问题数值解法2.解决的问题：（1）最小二乘拟合方法的实际应用曲线拟合就是拟合测量数据曲线。所选择的曲线有时通过数据点，但在其他点上，曲线接近它们而不必通过它们13,41~在大多数情况下，选择曲线使得数据点的平方误差和最小。这种选择就是最小二乘曲线拟合。下面介绍一下最小二乘法拟合的基本原理。设已知个数据点)(i=0，1，„，一1)，求(m一1)次最小二乘拟合多项式：其中设拟合多项式为各正交多项式：的线性组合：则继续往向下推导得：继续推导最后可得一般形式的m一1次多项式：即为最小二乘拟合多项式其拟合精度由下式来评定：应用实例：某建筑物176d水平位移测量数据如下表所示，在程序编制过程中，为了防止运算溢出，𝑥𝑖用来代替，其中𝑥̅=(∑)𝑥𝑖𝑛−1𝑖=0/𝑛。此时，拟合多项式的形式为：此时，拟合多项式的形式为：运用最小二乘多项式拟合时，拟合多项式的次数越高，其拟合精度未必越高。以拟合最高次数l9次为例，拟合系数如表2，拟合的精度评定见表3。根据水平位移的观测数据，实现了累计观测时间与水平位移的曲线拟合，在有限的测量数据条件下，表述了时间与该建筑物水平位移之间的函数关系。曲线拟合的最小二乘法在解决这类问题的数据处理和误差分析中应用非常广泛，提高了数据处理的效率和精确度，最d-乘曲线拟合实现方法简明、适用，可应用于类似的测量数据处理和实验研究。举例：在万能拉拨机中有一个园柱形凸轮，其底园半径R=300mm，凸轮的上端面不在同一平面上，而要根据动杆位移变化的需要进行设计制造。按设计要求，将底园周18等分，旋转一周。第个分点对应柱高，数据见下表。为了数控加工，需要计算出园周上任一点的柱高。凸轮高度的数据（单位：mm）分点i0和1812345柱高502.8525.0514.3451.0326.5188.6分点i67891011柱高92.259.662.2102.7147.1191.6分点i121314151617柱高236.0280.5324.9369.4413.8458.3我们将园周展开，借助MATLAB软件画出对应的柱高曲线散点图(左下图)。clear;close;x=linspace(0,2*pi*300,19);y=[502.8,525.0,514.3,451.0,326.5,188.6,92.2,59.6,62.2,102.7,147.1,191.6,236.0,280.5,324.9,369.4,413.8,458.3,502.8];plot(x,y,’o’);axis([0,2000,0,550]);可见，可以用三次多项式插值，下面给出借助MATLAB软件画出的柱高插值曲线图(右上图)。xi=0:2*pi*300;yi=interp1(x,y,xi,’cubic’);plot(xi,yi);（2）数据插值方法及应用多项式插值设有m次多项式mmmmaxaxaxaxP1110)(通过所有1n个点),(,),,(),,(1100nnyxyxyx，那么就有niyaxaxaxaimimmimi,,1,0,1110可以证明当nm且nxxx10时，这样的多项式存在且唯一。若要求得到函数表达式，可直接解上面方程组。若只要求得函数在插值点处数值，可用下列Lagrange插值公式)()(,00nijjjijniinxxxxyxP多项式插值光滑但不具有收敛性，一般不宜采用高次多项式（如7m）插值。第二部分、水槽的液位控制1问题的提出#2泵#1泵#1水槽#2水槽出口1出口2探针#1探针#2挡板耦合水槽如图3.7所示，对耦合水槽建模如下[i]：11111312dHAQHHHdt22222312dHAQHHHdt(3-39)其中，21236.52AAcm为1号和2号水槽的截面积，1H，2H分别是两个水槽的液面高度，1Q,2Q分别是两个水泵向水槽中注水的体积流速（3cms），1,2,3分别是与1H,2H,12HH对应的比例常数，通过升高或者降低挡板可以控制两个水槽之间的流速，关闭2号水泵，通过1号水泵向1号水槽中注水以控制2号水槽液面高度2用MATLAB语言编程的方法实现水槽液位模糊控制系统的设计与仿真（1）建立模糊推理系统(FIS)要想设计整个模糊控制系统,首先需要建立模糊推理系统,这也是建模的主要工作。采用MATLAB语言编程的方法建立模糊推理系统的过程如下所述。1)建立新的FIS结构设所建的FIS结构为FC.fis,则命令行为:FC=newfis(’FC’);这时便建立了一个新的名叫FC的mamdani型系统。此时它的输入、输出及规则都还没有定义。2)添加输入输出变量由于本例的模糊推理系统有两个输入变量时间t和注入1号水的速度v,一个输出变量为2号水槽液面高度h2,它们所对应的模糊语言变量T的论域为[2，200],V的论域为[0，75],H2的论域为[0，15]。则用addvar函数添加这三个变量的命令行分别为如下:FC=addvar(FC,’input’,’T’,[2，200]);FC=addvar(FC,’input’,’V’,[0，75]);FC=addvar(FC,’output’,’H’,[0，15]);3)添加隶属度函数已经加入模糊语言变量之后,下一步就是给每个语言变量添加隶属度函数。给变量T添加隶属度函数的语句如下:FC=addmf(FC,’input’,1,’NB’,’trapmf’,[-6,-6,-5,-3])FC=addmf(FC,’input’,1,’NS’,’trapmf’,[-5,-3,-2,0]);FC=addmf(FC,’input’,1,’ZR’,’trimf’,[-2,0,2]);FC=addmf(FC,’input’,1,’PS’,’trapmf’,[0,2,3,5]);FC=addmf(FC,’input’,1,’PB’,’trapmf’,[3,5,6,6]);变量V和H2隶属度函数添加的方法与T类似。4)设计并添加模糊推理规则表用编号的格式给本系统建立模糊列表。根据编制模糊列表的要求,给变量T、V及U的五个语言值NB、NS、ZR、PS、PB依次编号为1、2、3、4、5,则根据前面定义好的规则表rulelist=1111122…33…44…555551234512…12…12…123455544354…44…43…33221然后可用addrule函数给模糊推理系统添加规则列表,添加规则列表的命令行如下:FC=addrule(FC,rulelist);到此便基本上完成了模糊推理系统的建立。对于建立好的FIS结构FC.fis,需要将它保存起来,以便在模糊控制系统中直接调用它。可以通过writefis(FC)指令将它保存在MATLAB选取的当前路径下。（2）对整个模糊系统进行设计和仿真结合被控对象和已经建立的模糊推理系统FC.fis,对整个模糊控制系统进行建模,并设置相应的仿真参数:给定输入R=1.0;采样周期T=1s,仿真时间长度t为200秒;一型模糊控制系统选取参数如下：0.07,0.75,80erKKH二型模糊控制系统选取参数如下：0.07,0.75,80erKKH，121.2,0.12采用一定的算法（一型模糊控制系统,KM迭代降型的区间二型模糊PI控制系统，Wu-Mendel近似法降型的区间二型模糊PI控制系统），计算出相应的响应输出。最后可得出模糊控制系统的阶跃响应曲线如下图所示。其中,横坐标代表时间t,纵坐标代表输出响应y。3分析三个控制系统对应的性能指标如表3.1中所示，使用ITAE（时间乘绝对误差的积分）作为衡量控制系统的性能标准。从图3.8及表3.1的性能指标中，可以看出，使用KM迭代降型法的IT2模糊PI控制系统性能最优，然后为使用Wu-Mendel近似降型法的IT2模糊PI控制系统，一型模糊PI控制系统的性能最差。与一型模糊PI控制系统相比，使用KM迭代降型的IT2模糊PI控制系统在保持上升时间（从稳态值10%上升至90%所需时间）不增加的情况下，减小了超调量05010015020005101520时间/s液面高度目标液面T1PIFS实际液面IT2PIFS1实际液面IT2PIFS2实际液面0501001502001414.51515.516时间/s液面高度目标液面T1PIFS实际液面IT2PIFS1实际液面IT2PIFS2实际液面和调节时间，而使用Wu-Mendel近似降型的IT2模糊PI控制系统在上升时间不明显增加时，超调量和调节时间只是略微下降，如表3.3所示。不同控制系统的性能对比系统类型T1模糊系统IT2模糊系统KMWu-Mendel上升时间(秒)303031超调量4.6289%2.4333%4.2686%调节时间(秒)=2%937381性能指标(ITAE)6747.65927.96384.4注：0()TITAEtetdt，T为仿真截止时间程序：MATLAB程序：m1=min([P(:,2),P(:,3)]);m2=max([P(:,3),P(:,2)]);L=length(x);u=[];l=[];fori=1:Lifx(i)=m1&x(i)=(m1+m2)/2mu1=1;mu2=gaussmf(x(i),[sigma,m2]);elseifx(i)(m1+m2)/2&x(i)=m2mu1=1;mu2=gaussmf(x(i),[sigma,m1]);elsemu1=gaussmf(x(i),[sigma,m1]);mu2=gaussmf(x(i),[sigma,m2]);endu=[u,max([mu1,mu2])];l=[l,min([mu1,mu2])];endfunction[CA,CA1,CAr]=centroid(A)left=min(A(1:8));right=max(A(1:8));dx=(right-left)/500;Forj=1:501X=left+(j-1)*dx;Ifx=A(1)|x=A(4);upperA(j)=0;ElseifxA(1)&xA(2);upperA(j)=(x-A(1))/(A(2)-A(1));Elseifx=A(2)&x=A(3);upperA(j)=1;ElseupperA(j)=(A(4)-x)/(A(4)-A(3));EndIfx=A(5)|x=A(8);lowerA(j)=0;ElseifxA(6)lowerA(j)=A(9)*(x-A(5))/(A(6)-A(5));ElseifxA(7)lowerA(j)=A(9);elselowerA(j)=A(9)*(A(8)-x)/(A(8)-A(7));end仿真结果：心得体会：本次水槽得液位控制仿真有利用到数学建模、MATLAB仿真，由于个人知识有限，可能理解偏差有失误，不过还是收货颇多。