第四章+MATLAB的数值计算功能

第四章MATLAB的数值计算功能Chapter4:NumericalcomputationofMATLAB数值计算是MATLAB最基本、最重要的功能，是MATLAB最具代表性的特点。MATLAB在数值计算过程中以数组和矩阵为基础。数组是MATLAB运算中的重要数据组织形式。前面章节对数组、矩阵的特征及其创建与基本运算规则等相关知识已作了较详尽的介绍，本章重点介绍常用的数值计算方法。一、多项式(Polynomial)`多项式在众多学科的计算中具有重要的作用，许多方程和定理都是多项式的形式。MATLAB提供了标准多项式运算的函数，如多项式的求根、求值和微分，还提供了一些用于更高级运算的函数，如曲线拟合和多项式展开等。1．多项式的表达与创建(ExpressionandCreatingofpolynomial)(1)多项式的表达（expressionofpolynomial）_Matlab用行矢量表达多项式系数(Coefficient)，各元素按变量的降幂顺序排列，如多项式为：P(x)=a0xn+a1xn-1+a2xn-2…an-1x+an则其系数矢量(Vectorofcoefficient)为：P=[a0a1…an-1an]如将根矢量(Vectorofroot)表示为：ar=[ar1ar2…arn]则根矢量与系数矢量之间关系为：（x-ar1）(x-ar2)…(x-arn)=a0xn+a1xn-1+a2xn-2…an-1x+an(2)多项式的创建（polynomialcreating）a)系数矢量的直接输入法利用poly2sym函数直接输入多项式的系数矢量，就可方便的建立符号形式的多项式。例1：创建多项式x3-4x2+3x+2poly2sym([1-432])ans=x^3-4*x^2+3*x+2POLYConvertrootstopolynomial.POLY(A),whenAisanNbyNmatrix,isarowvectorwithN+1elementswhicharethecoefficientsofthecharacteristicpolynomial,DET(lambda*EYE(SIZE(A))-A).POLY(V),whenVisavector,isavectorwhoseelementsarethecoefficientsofthepolynomialwhoserootsaretheelementsofV.Forvectors,ROOTSandPOLYareinversefunctionsofeachother,uptoordering,scaling,androundofferror.b)由根矢量创建多项式已知根矢量ar,通过调用函数p=poly(ar)产生多项式的系数矢量,再利用poly2sym函数就可方便的建立符号形式的多项式。例2：由根矢量创建多项式。将多项式(x-6)(x-3)(x-8)表示为系数形式的多项式。a=[638]%根矢量pa=poly(a)%求系数矢量ppa=poly2sym(pa)%以符号形式表示原多项式ezplot(ppa,[-50,50])pa=1-1790-144ppa=x^3-17*x^2+90*x-144说明：（1）根矢量元素为n，则多项式系数矢量元素为n+1；（2）函数poly2sym(pa)把多项式系数矢量表达成符号形式的多项式，缺省情况下自变量符号为x，可以指定自变量。（3）使用简单绘图函数ezplot可以直接绘制符号形式多项式的曲线。例3：由给定复数根矢量求多项式系数矢量。r=[-0.5-0.3+0.4i-0.3-0.4i];p=poly(r)pr=real(p)ppr=poly2sym(pr)p=1.00001.10000.55000.1250pr=1.00001.10000.55000.1250ppr=x^3+11/10*x^2+11/20*x+1/8说明：含复数根的根矢量所创建的多项式要注意：（1）要形成实系数多项式，根矢量中的复数根必须共轭成对；（2）含复数根的根矢量所创建的多项式系数矢量中，可能带有很小的虚部，此时可采用取实部的命令(real)把虚部滤掉。如果需要进行系数表示形式的多项式的求根运算，有两种方法可以实现，一是直接调用求根函数roots，poly和roots互为逆函数。另一种是先把多项式转化为伴随矩阵，然后再求其特征值，该特征值即是多项式的根。c)特征多项式输入法用poly函数还可实现由矩阵的特征多项式系数创建多项式。条件：特征多项式系数矢量的第一个元素必须为1。例2：求三阶方阵A的特征多项式系数，并转换为多项式形式。a=[638;756;135]Pa=poly(a)%求矩阵的特征多项式系数矢量Ppa=poly2sym(pa)Pa=1.0000-16.000038.0000-83.0000Ppa=x^3-17*x^2+90*x-144注：n阶方阵的特征多项式系数矢量一定是n+1阶的。例4：将多项式的系数表示形式转换为根表现形式。求x3-6x2-72x-27的根a=[1-6-72-27]r=roots(a)r=12.1229-5.7345-0.3884MATLAB约定，多项式系数矢量用行矢量表示，根矢量用列矢量表示。2.多项式的乘除运算(Multiplicationanddivisionofpolynomial)向量的卷积与解卷积对应着多项式的乘除法，多项式乘法（卷积）用函数conv(a,b)实现，除法（解卷积）用函数deconv(a,b)实现。长度为m的向量a和长度为n的向量b的卷积定义为：C(k)=C向量的长度为：m+n-1解卷积是卷积的逆运算，向量a对向量c进行解卷积将得到商向量q和余量r，并且满足：例1：a(s)=s2+2s+3,b(s)=4s2+5s+6,计算a(s)与b(s)的乘积。a=[123];b=[456];%建立系数矢量c=conv(a,b)cs=poly2sym(c,’s’)%建立指定变量为s的符号形式多项式c=413282718cs=4*s^4+13*s^3+28*s^2+27*s+18例2：展开(s2+2s+2)(s+4)(s+1)（多个多项式相乘）c=conv([1,2,2],conv([1,4],[1,1]))cs=poly2sym(c,’s’)%（指定变量为s）c=1716188cs=s^4+7*s^3+16*s^2+18*s+8例3：求多项式s^4+7*s^3+16*s^2+18*s+8分别被(s+4),(s+3)除后的结果。c=[1716188];[q1,r1]=deconv(c,[1,4])%q—商矢量，r—余数矢量[q2,r2]=deconv(c,[1,3])cc=conv(q2,[1,3])%对除(s+3)结果检验test=((c-r2)==cc)q1=1342r1=00000q2=1446r2=0000-10cc=17161818test=111113.其他常用的多项式运算命令(Othercomputationcommandofpolynomial)pa=polyval(p,s)按数组运算规则计算给定s时多项式p的值。pm=polyvalm(p,s)按矩阵运算规则计算给定s时多项式p的值。[r,p,k]=residue(b,a)部分分式展开，b,a分别是分子（numerator）分母（denominator）多项式系数矢量，r,p,k分别是留数、极点和直项矢量p=polyfit(x,y,n)用n阶多项式拟合x，y矢量给定的数据。polyder(p)多项式微分。注：对于多项式b(s)与不重根的n阶多项式a(s)之比，其部分分式展开为：)()()(2211skpsrLpsrpsrsasbnn式中：p1,p2,…,pn称为极点(poles)，r1,r2,…,rn称为留数(residues)，k(s)称为直项(directterms)，假如a(s)含有m重根pj，则相应部分应写成：mjmjjjjjpsrLpsrpsr)()(121RESIDUEPartial-fractionexpansion(residues).[R,P,K]=RESIDUE(B,A)findstheresidues,polesanddirecttermofapartialfractionexpansionoftheratiooftwopolynomialsB(s)/A(s).Iftherearenomultipleroots,B(s)R(1)R(2)R(n)----=--------+--------+...+--------+K(s)A(s)s-P(1)s-P(2)s-P(n)VectorsBandAspecifythecoefficientsofthenumeratoranddenominatorpolynomialsindescendingpowersofs.TheresiduesarereturnedinthecolumnvectorR,thepolelocationsincolumnvectorP,andthedirecttermsinrowvectorK.Thenumberofpolesisn=length(A)-1=length(R)=length(P).Thedirecttermcoefficientvectorisemptyiflength(B)length(A),otherwiselength(K)=length(B)-length(A)+1.IfP(j)=...=P(j+m-1)isapoleofmultplicitym,thentheexpansionincludestermsoftheformR(j)R(j+1)R(j+m-1)--------+------------+...+------------s-P(j)(s-P(j))^2(s-P(j))^m[B,A]=RESIDUE(R,P,K),with3inputargumentsand2outputarguments,convertsthepartialfractionexpansionbacktothepolynomialswithcoefficientsinBandA.例3：对(3x4+2x3+5x2+4x+6)/(x5+3x4+4x3+2x2+7x+2)做部分分式展开a=[134272];b=[32546];[r,s,k]=residue(b,a)r=1.1274+1.1513i1.1274-1.1513i-0.0232-0.0722i-0.0232+0.0722i0.7916s=-1.7680+1.2673i-1.7680-1.2673i0.4176+1.1130i0.4176-1.1130i-0.2991k=[](分母阶数高于分子阶数时，k将是空矩阵，表示无此项)例5：对一组实验数据进行多项式最小二乘拟合(leastsquarefit)x=[12345];%实验数据y=[5.543.1128290.7498.4];p=polyfit(x,y,3)%做三阶多项式拟合x2=1:.1:5;y2=polyval(p,x2);%根据给定值计算多项式结果plot(x,y,’o’,x2,y2)一．线性代数(LinearAlgebra)解线性方程(Linearequation)就是找出是否存在一个唯一的矩阵x，使得a,b满足关系：ax=b或xa=bMALAB中x=a\b是方程ax=b的解，x=b/a是方程式xa=b的解。通常线性方程多写成ax=b，“\”较多用，两者的关系为：(b/a)’=(a’\b’)系数矩阵a可能是m行n列的，有三种情况：*方阵系统:(Squarematrix)m=n可求出精确解（a必须是非奇异(nonsingular)，即满秩(fullrank)。*超定系统：(Overdetermindsystem)mn可求出最小二乘解*欠定系统：(Underdetermindsystem)mn可尝试找出含有