原创 MATLAB多项式的表达方式及其操作

2009-11-16 13:05 9466 6 4 分类: 软件与OS
3.8 多项式的表达方式及其操作


3.8.1.2 多项式行向量的创建方法

【 * 例 3.8.1 .2-1 】求 3 阶方阵 A 的特征多项式。
A=[11 12 13;14 15 16;17 18 19];
PA=poly(A) %A 的特征多项式
PPA=poly2str(PA,'s') % 以较习惯的方式显示多项式
PA =
1.0000 -45.0000 -18.0000 -0.0000
PPA =
s^3 - 45 s^2 - 18 s - 2.8387e-015


【 * 例 3.8.1 .2-2 】由给定根向量求多项式系数向量。
R=[-0.5,-0.3+0.4*i,-0.3-0.4*i]; % 根向量
P=poly(R) %R 的特征多项式
PR=real(P) % 求 PR 的实部
PPR=poly2str(PR,'x')
P =
1.0000 1.1000 0.5500 0.1250
PR =
1.0000 1.1000 0.5500 0.1250
PPR =
x^3 + 1.1 x^2 + 0.55 x + 0.125


3.8.2 多项式运算函数

【 * 例 3.8.2 -1 】求 点击查看完整图片的“商”及“余”多项式。
p1=conv([1,0,2],conv([1,4],[1,1])); % 计算分子多项式
p2=[1 0 1 1]; % 注意缺项补零
[q,r]=deconv(p1,p2);
cq=' 商多项式为 '; cr=' 余多项式为 ';
disp([cq,poly2str(q,'s')]),disp([cr,poly2str(r,'s')])
商多项式为 s + 5
余多项式为 5 s^2 + 4 s + 3


【 * 例 3.8.2 -2 】两种多项式求值指令的差别。
S=pascal(4) % 生成一个 4 阶方阵
P=poly(S);PP=poly2str(P,'s')
PA=polyval(P,S) % 独立变量取数组 S 元素时的多项式值
PM=polyvalm(P,S) % 独立变量取矩阵 S 时的多项式值
S =
1 1 1 1
1 2 3 4
1 3 6 10
1 4 10 20
PP =
s^4 - 29 s^3 + 72 s^2 - 29 s + 1
PA =
1.0e+004 *
0.0016 0.0016 0.0016 0.0016
0.0016 0.0015 -0.0140 -0.0563
0.0016 -0.0140 -0.2549 -1.2089
0.0016 -0.0563 -1.2089 -4.3779
PM =
1.0e-011 *
-0.0077 0.0053 -0.0096 0.0430
-0.0068 0.0481 -0.0110 0.1222
0.0075 0.1400 -0.0095 0.2608
0.0430 0.2920 -0.0007 0.4737


【 * 例 3.8.2 -3 】部分分式展开。
a=[1,3,4,2,7,2]; % 分母多项式系数向量
b=[3,2,5,4,6]; % 分子多项式系数向量
[r,s,k]=residue(b,a)
r =
1.1274 + 1.1513i
1.1274 - 1.1513i
-0.0232 - 0.0722i
-0.0232 + 0.0722i
0.7916
s =
-1.7680 + 1.2673i
-1.7680 - 1.2673i
0.4176 + 1.1130i
0.4176 - 1.1130i
-0.2991
k =
[]



10.1 根


找出多项式的根,即多项式为零的值,可能是许多学科共同的问题,。MATLAB求解这个问题,并提供其它的多项式操作工具。MATLAB里,多项式由一个行向量表示,它的系数是按降序排列。例如,输入多项式x412x30x225x116


? p=[1-12025116]


p =


1-12025116


注意,必须包括具有零系数的项。除非特别地辨认,MATLAB无法知道哪一项为零。给出这种形式,用函数roots找出一个多项式的根。


? r="roots"(p)


r =


11.7473


2.7028


-1.2251 + 1.4672i


-1.2251 - 1.4672i


因为在MATLAB中,无论是一个多项式,还是它的根,都是向量,MATLAB按惯例规定,多项式是行向量,根是列向量。给出一


? pp="poly"(r)


pp =


1.0e+002 *


Columns 1 through 4


0.0100-0.12000.00000.2500


Column 5


1.1600 + 0.0000i



? pp="real"(pp) %throw away spurious imaginary part


pp =


1.0000-12.00000.000025.0000116.0000


因为MATLAB无隙地处理复数,当用根重组多项式时,如果一些根有虚部,由于截断误差,则poly的结果有一些小的虚部,这是很普通的。消除虚假的虚部,如上所示,只要使用函数real抽取实部。


10.2 乘法


函数conv支持多项式乘法(执行两个数组的卷积)。考虑两个多项式a(x)=x32x23x4b(x)= x34x29x16的乘积:


? a=[1234] ;b=[14916];


? c="conv"(a , b)


c =


162050758464


结果是c(x)=x66x520x450x375x284x64。两个以上的多项式的乘法需要重复使用conv


10.3 加法


对多项式加法,MATLAB不提供一个直接的函数。如果两个多项式向量大小相同,标准的数组加法有效。把多项式a(x)与上面给出的b(x)相加。


? d="a"+b


d =


261220


结果是d(x)= 2x36x212x20。当两个多项式阶次不同,低阶的多项式必须用首零填补,使其与高阶多项式有同样的阶次。考虑上面多项式cd相加:


? e="c"+[000d]


e =


162052819684


结果是e(x)= x66x520x452x381x296x84。要求首零而不是尾零,是因为相关的系数象x幂次一样,必须整齐。


如果需要,可用一个文件编辑器创建一个函数M文件来执行一般的多项式加法。精通MATLAB工具箱包含下列实现:


function p="mmpadd"(a,b)


%MMPADD Polynomial addition.


%MMPADD(A,B) adds the polynomial A and B


%Copyright (c) 1996 by Prentice Hall,Inc.


if nargin<2


error(' Not enough input arguments ')


end


a=a(.' ;%make sure inputs are polynomial row vectors


b=b(.' ;


na=length(a) ;%find lengths of a and b


nb=length(b) ;


p=[zeros(1,nb-na) a]+[zeros(1,na-nb) b] ;%add zeros as necessary


现在,为了阐述mmpadd的使用,再考虑前一页的例子。


? f="mmpadd"(c,d)


f =


162052819684


它与上面的e相同。当然,mmpadd也用于减法。


?g=mmpadd(c , -d)


g =


162048697244


结果是g(x)= x66x520x448x369x272x44


10.4 除法


在一些特殊情况,一个多项式需要除以另一个多项式。在MATLAB中,这由函数deconv完成。用上面的多项式bc


? [q , r]=deconv(c , b)


q =


1234


r =


0000000


这个结果是bc除,给出商多项式q和余数r,在现在情况下r是零,因为bq的乘积恰好是c


10.5 导数


由于一个多项式的导数表示简单,MATLAB为多项式求导提供了函数polyder


? g


g =


162048697244


? h="polyder"(g)


h =


6308014413872


10.6估值


根据多项式系数的行向量,可对多项式进行加,减,乘,除和求导,也应该能对它们进行估值。在MATLAB中,这由函数polyval来完成。


? x="linspace"(-1, 3) ; %choose 100 data points between -1and 3.


? p=[14-7-10] ;%uses polynomial p(x) = x34x27x10


? v="polyval"(p , x) ;


计算x值上的p(x),把结果存在v里。然后用函数plot绘出结果。


? plot(x , v),title(' x^3+4x^2-7x-10 '),xlabel(' x ')


点击查看完整图片


10.1多项式估值



10.7有理多项式


在许多应用中,例如富里哀(Fourier),拉普拉斯(Laplace)Z变换,出现有理多项式或两个多项式之比。在MATLAB中,有理多项式由它们的分子多项式和分母多项式表示。对有理多项式进行运算的两个函数是residuepolyder。函数residue执行部分分式展开。


? num="10"*[12] ; %numerator polynomial


? den="poly"([-1; -3; -4]) ;%denominator polynomial


? [res, poles, k]=residue(num, den)


res =


-6.6667


5.0000


1.6667


poles =


-4.0000


-3.0000


-1.0000


k =


[]


结果是余数、极点和部分分式展开的常数项。上面的结果说明了该问题:


点击查看完整图片


这个函数也执行逆运算。


? [n, d]=residue(res, poles, k)


n =


0.000010.000020.0000


d =


1.00008.000019.000012.0000


? roots(d)


ans =


-4.0000


-3.0000


-1.0000


在截断误差内,这与我们开始时的分子和分母多项式一致。residue也能处理重极点的情况,尽管这里没有考虑。


正如前面所述,函数polyder,对多项式求导。除此之外,如果给出两个输入,则它对有理多项式求导。


? [b , a]=polyder(num , den)


b =


-20-140-320-260


a =


116102328553456144


该结果证实:


点击查看完整图片


10.8 M文件举例


本章说明了在精通MATLAB工具箱中两个函数。这些函数说明了本章所论述的多项式概念和如何写M文件函数。关于M文件的更多信息,参阅第8章。


在讨论M文件函数的内部结构之前,我们考虑这些函数做些什么。


? n%earlier data


n =


0.000010.000020.0000


? b%earlier data


b =


-20-140-320-260


? mmpsim(n)%strip away negligible leading term


ans =


10.000020.0000


? mmp2str(b)%convert polynomial to string


ans =


-20s^3 - 140s^2 - 320s^1 - 260


? mmp2str(b , ' x ')


ans =


-20x^3 - 140x^2 - 320x^1 - 260


? mmp2str(b , [] , 1)


ans =


-20*(s^3 + 7s^2 + 16s^1 + 13)


? mmp2str(b , ' x ' , 1)


ans =


-20*(x^3 + 7x^2 + 16x^1 + 13)


这里函数mmpsim删除在多项式n中近似为零的第一个系数。函数mmp2str把数值多项式变换成等价形式的字符串表达式。该两个函数的主体是:


function y="mmpsim"(x,tol)


%MMPSIM Polynomial Simplification,Strip Leading Zero Terms.


%MMPSIM(A) Delets leading zeros and small coefficients in the


%polynomial A(s). Coefficients are considered small if their


%magnitude is less than both one and norm(A)*1000*eps.


%MMPSIM(A,TOL) uses TOL for its smallness tolerance.


%Copyright (c) 1996 by Prentice-Hall,Inc.


if nargin<2, tol="norm"(x)*1000*eps; end


x=x(.' ;%make sure input is a row


i=find(abs(x)<.99&abs(x)<tol) ;%find insignificant indices


x(i)=zeros(1, length(i)) ;%set them to zero


i=find(x~=0) ;%find significant indices


if isempty(i)


y=0 ;%extreme case: nothing left!


else


y=x(i(1) : length(x)) ;%start with first term


end%and go to end of polynomial


function s="mmp2str"(p,v,ff)


%MMP2STR Polynomial Vector to String Conversion.


%MMP2STR(P) converts the polynomial vector P into a string.


%For example: P = [234] becomes the string ' 2s^2 + 3s + 4 '


%


%MMP2STR(P,V) generates the string using the variable V


%instead of s. MMP2STR([234],' z ') becomes ' 2z^2 + 3z + 4 '


%


%MMP2STR(P,V,1) factors the polynomial into the product of a


%constant and a monic polynomial.


%MMP2STR([234],[],1) becomes ' 2(s^2 + 1.5s + 2) '


%Copyright (c) 1996 by Prentice-Hall,Inc.


if nargin<3, ff="0"; end%factored form is False


if nargin <2, v=' s ' ; end%default variable is ' s '


if isempty(v), v=' s ' ; end%default variable is ' s '


v=v(1) ;%variable must be scalar


p=mmpsim(p) ;%strip insignificant terms


n=length(p) ;


if ff%put in factored form


K=p(1) ; Ka="abs"(K) ; p="p/K";


if abs(K-1)<1e-4


pp=[]; pe=[] ;


elseif abs(K+1)<1e-4


pp=' -(' ; pe= ') ' ;


elseif abs(Ka-round(Ka))<=1e-5*Ka


pp=[sprintf(' %.0f ', K) '*( ' ] ; pe= ') ' ;


else


pp=[sprintf(' %.4g ' , K) '*(' ] ; pe= ') ' ;


end


else%not factored form


K=p(1);


pp=sprintf(' %.4g ' , K) ;


pe=[];


end


if n==1%polynomial is just a constant


s=sprintf(' %.4g ',K);


return


end


s=[pp v ' ^ ' sprintf(' %.0f ',n-1)];%begin string construction


for i="2:n-1"%catenate center terms in polynomial


if p(i)<0, pm= ' -' ;else,if p(i)<0,pm= ' ;end


if p(i)= =1,pp=[] ; else,pp=sprintf(' %.4g ', abs(p(i))) ;end


if p(i)~ =0,s=[spmppv' ^ ' sprintf(' %.0f ',n-i)] ;end


end


if p(n)~ =0,pp=sprintf(' %.4g ',abs(p(n))); else,pp=[];end


if p(n)<0 , pm= ' -' ; elseifp(n)>0 , pm= ' +' ; else,pm=[];end


s=;%add final terms


10.9 小结


下列表格概括了在本章所讨论的多项式操作特性。


10.1




项式函数


conv(a, b)


乘法


[q, r]=deconv(a, b)


除法


poly(r)


用根构造多项式


polyder(a)


对多项式或有理多项式求导


polyfit(x, y, n)


多项式数据拟合


polyval(p, x)


计算x点中多项式值


[r, p, k]=residue(a, b)


部分分式展开式


[a, b]=residue(r, p, k)


部分分式组合


roots(a)


求多项式的根


10.2



精 通 MATLAB 多 项 式 操 作


mmp2str(a)


多项式向量到字符串变换,a(s)


mmp2str(a, ' x ')


多项式向量到字符串变换,a(x)


mmp2str(a, ' x ', 1)


常数和符号多项式变换


mmpadd(a, b)


多项式加法


mmpsim(a)


多项式简化

文章评论0条评论)

登录后参与讨论
我要评论
0
6
关闭 站长推荐上一条 /5 下一条