第六讲：Matlab多项式和代数方程求解器
合集下载
- 1、下载文档前请自行甄别文档内容的完整性，平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益，请联系客服反馈,我们会尽快为您处理(人工客服工作时间：9:00-18:30)。
x= -0.1819 -1.6630 2.2172 -0.4467
线性方程组的解析解
linsolve(A,b)：解线性方程组 Ax b
例：解方程组
x
2yz xz3
2
x 3y 8
>> A=[1 2 –1; 1 0 1; 1 3 0]; >> b=[2;3;8]; >> x=linsolve(A,b)
Matlab 没有提供专门进行多项式加减运算的函数，事实 上，多项式的加减就是其所对应的系数向量的加减运算
对于次数相同的多项式，可以直接对其系数向量进行 加减运算；
如果两个多项式次数不同，则应该把低次多项式中系 数不足的高次项用 0 补足，然后进行加减运算。
例： p1 2x3 x2 3
p2 2x 1 p1 p2 2x3 x2 2x 4
Matlab 多项式运算 与代数方程求解器
Matlab 多项式运算
Matlab 中多项式的表示方法
在 Matlab 中，n 次多项式是用一个长度为 n+1的向量来 表示，缺少的幂次项系数为 0。
p( x) an x n an1x n1 a1x a0 在 Matlab中表示为向量：[an ,an1, ,a1,a0 ]
1.恰定方程组的解
方程ax=b(a为非奇异) x=a-1 b 矩阵逆
两种解: x=inv(a)b — 采用求逆运算解方程 x=a\b — 采用左除运算解方程
例: x1+2x2=8 2x1+3x2=13
方程ax=b
1 2
2 3
x1
x
2
=
8 13
a x=b
a=[1 2;2 3];b=[8;13];
1.系数向量的直接输入法
例： 2x3 x2 3
[2, 1, 0, 3]
注：系数中的零不能省！ 按降幂顺序
多项式的符号形式：poly2sym， p1=poly2str(p,‘x’) >> poly2sym([2,-1,0,3])
例：输入多项式： x3 5 x2 6x 33 解：p=[1 -5 6 -33]; Poly2sym(p)
[2, 1, 0, 3] [ 0, 0[, 22,,11]] [2, 1, 2, 4]
多项式四则运算
多项式乘法运算： k = conv(p,q)
例：计算多项式 2x3 x2 3 和 2x 1 的乘积 >> p=[2,-1,0,3]; >> q=[2,1]; >> k=conv(p,q);
多项式除法运算： [k,r] = deconv(p,q) 其中 k 返回的是多项式 p 除以 q 的商，r 是余式。 [k,r]=deconv(p,q) <==> p=conv(q,k)+r
计算两多项式的乘除法
p=[2 -5 6 -1 9]; poly2sym(p) d=[3 -90 -18]; poly2sym(d) pd=conv(p,d) poly2sym(pd) p1=deconv(pd,d)
多项式的求导
polyder
k=polyder(p) : 多项式 p 的导数； k=polyder(p,q): p*q 的导数； [k,d]=polyder(p,q): p/q 的导数，k 是分子，d 是分母
x=inv(a)*b
x=a\b
x=
x=
2.00
2.00
3.00
3.00
2.超定方程组的解
方程 ax=b ,m<n时此时不存在唯一解。 方程解 (a ' a)x=a ' b
x=(a' a)-1 a ' b —— 求逆法 x=a\b —— matlab用最小二乘法找一
个准确地基本解。
例: x1+2x2=1 2x1+3x2=2 3x1+4x2=3
是基于伪逆pinv求得的。
x1+2x2+3x3=1 2x1+3x2+4x3=2
x1
1 2
2 3
3 4
x2
x3
= 12
a=[1 2 3;2 3 4];b=[1;2];
a x =b
x=a\b
x=pinv(a)b
x=
x=
1.00
0.83
0
0.33
0
-0.17
内部包含许多自适应算法：对超定方程用最小二乘 法，对欠定方程它将给出范数最小的一个解。
b是列向量！
非线性方程的根
非线性方程的解法：
1.二分法
function y=erfen(fun,a,b,esp) if nargin<4 esp=1e-4;end if feval(fun,a)*feval(fun,b)<0
n=1; c=(a+b)/2; while c>esp
if feval(fun,a)*feval(fun,c)<0 b=c;c=(a+b)/2;
多项式的值
计算多项式在给定点的值,有两种形式，对应两种算法。
代数多项式求值
y = polyval(p,x): 计算多项式 p 在 x 点的值
注：若 x 是向量或矩阵，则采用数组运算 (点运算)！
例：已知 p(x) 2x3 x2 3，分别取 x=2 和一个 22 矩阵， 求 p(x) 在 x 处的值
elseif feval(fun,c)*feval(fun,b)<0 a=c;c=(a+b)/2;
else y=c;esp=10000; end n=n+1; end y=c; elseif feval(fun,a)==0 y=a; elseif feval(fun,b)==0 y=b; else disp('These may not be a root on the intercal'); end n
求解下列方程组 0.4096x1+0.1234x2+0.3678x3+0.2943x4=0.4043 0.2246x1+0.3872x2+0.4015x3+0.1129x4=0.1550 0.3645x1+0.1920x2+0.3781x3+0.0643x4=0.4240 0.1784x1+0.4002x2+0.2786x3+0.3927x4=-0.2557
例：已知 p(x) 2x3 x2 3，则
>> p=[2,-1,0,3]; >> x=[-1, 2;-2,1];polyval(p,x) >> polyvalm(p,x)
多项式的零点（根）
x=roots(p)：若 p 是 n 次多项式，则输出是 p=0 的 n 个根组成的 n 维向量。
例：已知 p(x) 2x3 x2 3，求 p(x) 的零点。
不涉及符号计算！
x = roots(p)
polyfit(X,Y,n)
线性方程组求解
对于方程ax=b，a 为an×m矩阵，有三种情况：
当n=m时，此方程成为“恰定”方程,寻求精确 解 当n>m时，此方程成为“超定”方程，寻求最小 二乘解 当n<m时，此方程成为“欠定”方程，寻求基本 解
matlab定义的除运算可以很方便地解上 述三种方程
2.由矩阵求其特征多项式，由函数poly实现。
a=[1 2 3;2 3 4; 3 4 5]; p1=poly(a) poly2sym(p1)
由特征多项式生成的多项式的首项系数一定为1。 N阶矩阵一般生成N次多项式。
3.由给定的根求其对应的多项式，也由poly函 数实现。
root=[-5 -3+4i -3-4i] p=poly(root)
1 2
1
2 3
3 4
x1 x2
=2
3
a x=b
a=[1 2;2 3;3 4];b=[1;2;3];
解1 x=a\b 解2 x=inv(a'a) a' b
x=
x=
1.0000
0000
-0.0000
0.0000
3.欠定方程组的解
当方程数少于未知量个数时,即不定 情况,有无穷多个解存在。 matlab可求出两个解： 用除法求的解x是具有最多零元素的解 是具有最小长度或范数的解，这个解
用二分法求方程 x2 x 1 0 的正根，要求误差小于0.05。
先编制函数文件fc.m： function y=fc(x) y=x^2-x-1;
在命令窗口输入： erfen('fc',0,10,0.05)
n= 56
ans = 1.6180
Matlab 非线性方程的数值求解
fzero(f,x0)：求方程 f=0 在 x0 附近的根。
因式分解和展开
MATLAB提供了对符号表达式、符号矩阵等 进行因式分解、同类项合并和展开等操作。
因式分解的主要函数有：factor、 collet、expand
多项式的简化
简化的主要函数有： pretty(s)%将代数式s转化为手写形式。 simplify%单一化简命令 simple
多项式运算小结
若要生成实系数多项式，则根中的复数必定对 应共轭。
4.多项式的符号表示 若矩阵A是符号表示其调用格式为：
Poly(A)%返回A的用符号’x’或’t’表示的特征多项式
Poly(A,v) )%返回A的用符号’v’表示的特征多项式
Poly2sym(C)%将向量C表示的多项式转化为用符号表示。 其变量用符号’x’表示。
>> p=[2,-1,0,3]; >> x=2; y=polyval(p,x) >> x=[-1, 2;-2,1]; y=polyval(p,x)
多项式的值
矩阵多项式求值
Y=polyvalm(p,X)
与代数多项式求值不同之处： 采用的是普通矩阵运算； X 必须是方阵
polyvalm(p,A) = 2*A*A*A - A*A + 3*eye(size(A)) polyval(P,A) = 2*A.*A.*A - A.*A + 3*ones(size(A))
例：已知 p(x) 2x3 x2 3，q(x) 2x 1 ， 求 p',( p q)',( p / q)'
>> k1=polyder([2,-1,0,3]); >> k2=polyder([2,-1,0,3],[2,1]); >> [k2,d]=polyder([2,-1,0,3],[2,1]);
>> p=[2,-1,0,3]; >> x=roots(p)
多项式的拟合
拟合函数：polyfit(X,Y,n) polyfit(X,Y,n) %X,Y为拟合的数据，n为拟合多项式的阶 数。 [p,s]=polyfit(X,Y,n) % p为拟合多项式系数向量，s为拟 合多项式系数向量的结构信息。
例：用5阶多项式对[0，/2]上的正弦函数值进行 最小二乘法拟合
poly2sym(p)
k = conv(p,q) [k,r] = deconv(p,q)
k = polyder(p) k = polyder(p,q) [k,d] = polyder(p,q)
y = polyval(p,x) Y = polyvalm(p,X)
多项式运算中， 使用的是多项式
系数向量，
4
原曲线
3.5
拟合曲线
3
2.5
2
1.5
1
0.5
0
-0.5
-1
0
1
2
3
4
5
6
7
多项式的部分分式展开
调用格式： [r,p,k]=residue(b,a) %多项式b(s)/a(s) 的部分分式展开，r是留数，p是极点,k是直 向量 [b,a]=residue(r,p,k)
a=[2 45 82];b=[3 5 8 3]; [r,p,k]=residue(b,a) [bb,aa]=residue(r,p,k)
Poly2sym(C,’v’) )%将向量C表示的多项式转化为用符号 表示。其变量用符号’v’表示。
syms a b c d v A=[a b;c d];
C=[1 0 -2 -5];
poly(A) %poly(A,v)
%poly2sym(C) poly2sym(C，’t’)
多项式四则运算
多项式加减运算
x=0:pi/20:pi/2;y=sin(x); a=polyfit(x,y,5); x1=0:pi/30:2*pi;y1=sin(x1); y2=a(1)*x1.^5+a(2)*x1.^4+a(3)*x1.^3+a(4)*x1.^2 +a(5)*x1+a(6); plot(x1,y1,'b-',x1,y2,'r*') legend(‘原曲线’,‘拟合曲线') axis([0, 7,-1.2,4])
在MATLAB命令窗口输入
a=[0.4096 0.1234 0.3678 0.2943 0.2246 0.3872 0.4015 0.1129 0.3645 0.1920 0.3781 0.0643 0.1784 0.4002 0.2786 0.3927]; b=[0.4043 0.1550 0.4240 -0.2557]; x=a\b
线性方程组的解析解
linsolve(A,b)：解线性方程组 Ax b
例：解方程组
x
2yz xz3
2
x 3y 8
>> A=[1 2 –1; 1 0 1; 1 3 0]; >> b=[2;3;8]; >> x=linsolve(A,b)
Matlab 没有提供专门进行多项式加减运算的函数，事实 上，多项式的加减就是其所对应的系数向量的加减运算
对于次数相同的多项式，可以直接对其系数向量进行 加减运算；
如果两个多项式次数不同，则应该把低次多项式中系 数不足的高次项用 0 补足，然后进行加减运算。
例： p1 2x3 x2 3
p2 2x 1 p1 p2 2x3 x2 2x 4
Matlab 多项式运算 与代数方程求解器
Matlab 多项式运算
Matlab 中多项式的表示方法
在 Matlab 中，n 次多项式是用一个长度为 n+1的向量来 表示，缺少的幂次项系数为 0。
p( x) an x n an1x n1 a1x a0 在 Matlab中表示为向量：[an ,an1, ,a1,a0 ]
1.恰定方程组的解
方程ax=b(a为非奇异) x=a-1 b 矩阵逆
两种解: x=inv(a)b — 采用求逆运算解方程 x=a\b — 采用左除运算解方程
例: x1+2x2=8 2x1+3x2=13
方程ax=b
1 2
2 3
x1
x
2
=
8 13
a x=b
a=[1 2;2 3];b=[8;13];
1.系数向量的直接输入法
例： 2x3 x2 3
[2, 1, 0, 3]
注：系数中的零不能省！ 按降幂顺序
多项式的符号形式：poly2sym， p1=poly2str(p,‘x’) >> poly2sym([2,-1,0,3])
例：输入多项式： x3 5 x2 6x 33 解：p=[1 -5 6 -33]; Poly2sym(p)
[2, 1, 0, 3] [ 0, 0[, 22,,11]] [2, 1, 2, 4]
多项式四则运算
多项式乘法运算： k = conv(p,q)
例：计算多项式 2x3 x2 3 和 2x 1 的乘积 >> p=[2,-1,0,3]; >> q=[2,1]; >> k=conv(p,q);
多项式除法运算： [k,r] = deconv(p,q) 其中 k 返回的是多项式 p 除以 q 的商，r 是余式。 [k,r]=deconv(p,q) <==> p=conv(q,k)+r
计算两多项式的乘除法
p=[2 -5 6 -1 9]; poly2sym(p) d=[3 -90 -18]; poly2sym(d) pd=conv(p,d) poly2sym(pd) p1=deconv(pd,d)
多项式的求导
polyder
k=polyder(p) : 多项式 p 的导数； k=polyder(p,q): p*q 的导数； [k,d]=polyder(p,q): p/q 的导数，k 是分子，d 是分母
x=inv(a)*b
x=a\b
x=
x=
2.00
2.00
3.00
3.00
2.超定方程组的解
方程 ax=b ,m<n时此时不存在唯一解。 方程解 (a ' a)x=a ' b
x=(a' a)-1 a ' b —— 求逆法 x=a\b —— matlab用最小二乘法找一
个准确地基本解。
例: x1+2x2=1 2x1+3x2=2 3x1+4x2=3
是基于伪逆pinv求得的。
x1+2x2+3x3=1 2x1+3x2+4x3=2
x1
1 2
2 3
3 4
x2
x3
= 12
a=[1 2 3;2 3 4];b=[1;2];
a x =b
x=a\b
x=pinv(a)b
x=
x=
1.00
0.83
0
0.33
0
-0.17
内部包含许多自适应算法：对超定方程用最小二乘 法，对欠定方程它将给出范数最小的一个解。
b是列向量！
非线性方程的根
非线性方程的解法：
1.二分法
function y=erfen(fun,a,b,esp) if nargin<4 esp=1e-4;end if feval(fun,a)*feval(fun,b)<0
n=1; c=(a+b)/2; while c>esp
if feval(fun,a)*feval(fun,c)<0 b=c;c=(a+b)/2;
多项式的值
计算多项式在给定点的值,有两种形式，对应两种算法。
代数多项式求值
y = polyval(p,x): 计算多项式 p 在 x 点的值
注：若 x 是向量或矩阵，则采用数组运算 (点运算)！
例：已知 p(x) 2x3 x2 3，分别取 x=2 和一个 22 矩阵， 求 p(x) 在 x 处的值
elseif feval(fun,c)*feval(fun,b)<0 a=c;c=(a+b)/2;
else y=c;esp=10000; end n=n+1; end y=c; elseif feval(fun,a)==0 y=a; elseif feval(fun,b)==0 y=b; else disp('These may not be a root on the intercal'); end n
求解下列方程组 0.4096x1+0.1234x2+0.3678x3+0.2943x4=0.4043 0.2246x1+0.3872x2+0.4015x3+0.1129x4=0.1550 0.3645x1+0.1920x2+0.3781x3+0.0643x4=0.4240 0.1784x1+0.4002x2+0.2786x3+0.3927x4=-0.2557
例：已知 p(x) 2x3 x2 3，则
>> p=[2,-1,0,3]; >> x=[-1, 2;-2,1];polyval(p,x) >> polyvalm(p,x)
多项式的零点（根）
x=roots(p)：若 p 是 n 次多项式，则输出是 p=0 的 n 个根组成的 n 维向量。
例：已知 p(x) 2x3 x2 3，求 p(x) 的零点。
不涉及符号计算！
x = roots(p)
polyfit(X,Y,n)
线性方程组求解
对于方程ax=b，a 为an×m矩阵，有三种情况：
当n=m时，此方程成为“恰定”方程,寻求精确 解 当n>m时，此方程成为“超定”方程，寻求最小 二乘解 当n<m时，此方程成为“欠定”方程，寻求基本 解
matlab定义的除运算可以很方便地解上 述三种方程
2.由矩阵求其特征多项式，由函数poly实现。
a=[1 2 3;2 3 4; 3 4 5]; p1=poly(a) poly2sym(p1)
由特征多项式生成的多项式的首项系数一定为1。 N阶矩阵一般生成N次多项式。
3.由给定的根求其对应的多项式，也由poly函 数实现。
root=[-5 -3+4i -3-4i] p=poly(root)
1 2
1
2 3
3 4
x1 x2
=2
3
a x=b
a=[1 2;2 3;3 4];b=[1;2;3];
解1 x=a\b 解2 x=inv(a'a) a' b
x=
x=
1.0000
0000
-0.0000
0.0000
3.欠定方程组的解
当方程数少于未知量个数时,即不定 情况,有无穷多个解存在。 matlab可求出两个解： 用除法求的解x是具有最多零元素的解 是具有最小长度或范数的解，这个解
用二分法求方程 x2 x 1 0 的正根，要求误差小于0.05。
先编制函数文件fc.m： function y=fc(x) y=x^2-x-1;
在命令窗口输入： erfen('fc',0,10,0.05)
n= 56
ans = 1.6180
Matlab 非线性方程的数值求解
fzero(f,x0)：求方程 f=0 在 x0 附近的根。
因式分解和展开
MATLAB提供了对符号表达式、符号矩阵等 进行因式分解、同类项合并和展开等操作。
因式分解的主要函数有：factor、 collet、expand
多项式的简化
简化的主要函数有： pretty(s)%将代数式s转化为手写形式。 simplify%单一化简命令 simple
多项式运算小结
若要生成实系数多项式，则根中的复数必定对 应共轭。
4.多项式的符号表示 若矩阵A是符号表示其调用格式为：
Poly(A)%返回A的用符号’x’或’t’表示的特征多项式
Poly(A,v) )%返回A的用符号’v’表示的特征多项式
Poly2sym(C)%将向量C表示的多项式转化为用符号表示。 其变量用符号’x’表示。
>> p=[2,-1,0,3]; >> x=2; y=polyval(p,x) >> x=[-1, 2;-2,1]; y=polyval(p,x)
多项式的值
矩阵多项式求值
Y=polyvalm(p,X)
与代数多项式求值不同之处： 采用的是普通矩阵运算； X 必须是方阵
polyvalm(p,A) = 2*A*A*A - A*A + 3*eye(size(A)) polyval(P,A) = 2*A.*A.*A - A.*A + 3*ones(size(A))
例：已知 p(x) 2x3 x2 3，q(x) 2x 1 ， 求 p',( p q)',( p / q)'
>> k1=polyder([2,-1,0,3]); >> k2=polyder([2,-1,0,3],[2,1]); >> [k2,d]=polyder([2,-1,0,3],[2,1]);
>> p=[2,-1,0,3]; >> x=roots(p)
多项式的拟合
拟合函数：polyfit(X,Y,n) polyfit(X,Y,n) %X,Y为拟合的数据，n为拟合多项式的阶 数。 [p,s]=polyfit(X,Y,n) % p为拟合多项式系数向量，s为拟 合多项式系数向量的结构信息。
例：用5阶多项式对[0，/2]上的正弦函数值进行 最小二乘法拟合
poly2sym(p)
k = conv(p,q) [k,r] = deconv(p,q)
k = polyder(p) k = polyder(p,q) [k,d] = polyder(p,q)
y = polyval(p,x) Y = polyvalm(p,X)
多项式运算中， 使用的是多项式
系数向量，
4
原曲线
3.5
拟合曲线
3
2.5
2
1.5
1
0.5
0
-0.5
-1
0
1
2
3
4
5
6
7
多项式的部分分式展开
调用格式： [r,p,k]=residue(b,a) %多项式b(s)/a(s) 的部分分式展开，r是留数，p是极点,k是直 向量 [b,a]=residue(r,p,k)
a=[2 45 82];b=[3 5 8 3]; [r,p,k]=residue(b,a) [bb,aa]=residue(r,p,k)
Poly2sym(C,’v’) )%将向量C表示的多项式转化为用符号 表示。其变量用符号’v’表示。
syms a b c d v A=[a b;c d];
C=[1 0 -2 -5];
poly(A) %poly(A,v)
%poly2sym(C) poly2sym(C，’t’)
多项式四则运算
多项式加减运算
x=0:pi/20:pi/2;y=sin(x); a=polyfit(x,y,5); x1=0:pi/30:2*pi;y1=sin(x1); y2=a(1)*x1.^5+a(2)*x1.^4+a(3)*x1.^3+a(4)*x1.^2 +a(5)*x1+a(6); plot(x1,y1,'b-',x1,y2,'r*') legend(‘原曲线’,‘拟合曲线') axis([0, 7,-1.2,4])
在MATLAB命令窗口输入
a=[0.4096 0.1234 0.3678 0.2943 0.2246 0.3872 0.4015 0.1129 0.3645 0.1920 0.3781 0.0643 0.1784 0.4002 0.2786 0.3927]; b=[0.4043 0.1550 0.4240 -0.2557]; x=a\b