特征多项式法
- 格式:doc
- 大小:449.50 KB
- 文档页数:11
《MATLAB 程序设计实践》课程考核
1、
特征多项式法
算法说明:
在MA TLAB 中编程实现的特征多项式法的函数为:chapoly 。
功能:通过求矩阵特征多项式的根来求其特征值。
调用格式:l=chapoly (A ) 其中,A 为已知矩阵;
L 为求得的矩阵特征值。
Matlab 程序代码
function l=Chapoly(A)
%特征多项式法求矩阵特征值 %已知矩阵:A
%求得的矩阵特征值:l syms t; N=size(A); n=N(1,1);
y=det(A-t*eye (n,n)); l=solve(y);
l=vpa(l,5); %结果取五位精度
例题:采用特征多项式法,求矩阵A 的特征值
A=⎥⎥⎦
⎤⎢⎢⎢⎣⎡--43210658511 解:在MA TLAB 窗口输入: >> A=[11 5 8;-5 6 10;-2 3 4] >> l=Chapoly(A)
输出计算结果为:
流程图为:
2、幂法
算法说明:
在MA TLAB中编程实现的幂法的函数为:pmethod。
功能:幂法求矩阵的主特征值及主特征向量。
调用格式:[l,v,s]=pmethod(A,x0,eps)
其中,A为已知矩阵;
x0为迭代初始向量;
eps为迭代的精度;
l为求得的矩阵主特征值;
v为求得的矩阵主特征向量;
s为迭代步数。
Matlab程序代码
function[l,v,s]=pmethod(A,x0,eps)
%幂法求矩阵的主特征值及主特征向量
%已知矩阵:A
%迭代初始向量:x0
%迭代的精度:eps
%求得的矩阵特征值:l
%求得的矩阵主特征向量:v
%迭代步数:s
If nargin==2
eps=1.0e-6;
end
v=x0; %v为主特征向量
M=5000; %迭代步数限制
m=0;
l=0;
for(k=1:M)
y=A*v;
m=max(y); %m为按模最大的分量
v=y/m;
if(abs(m-1)<eps)
l=m; %到所需精度,退出,l 为主特征值 s=k; %s 为迭代步数 return; else
if(k==m)
disp(‘迭代步数太多,收敛速度太慢!’); l=m; s=m; else
l=m;
end end end
例题:幂法求取特征值实例:采用幂法求矩阵A 的主特征值和主特征向量。
A=⎥⎥⎥⎦
⎤⎢⎢⎢⎣⎡-131716251 解:在MA TLAB 窗口输入: >>A=[1 5 2;6 -1 7;1 3 1]; >>x0=[1 1 1]’;
>>[l,v,s]=pmethod(A,x0)
输出计算结果为:
流程图为:
第2题
1.因为体积一栏的数据是以粒度为变量的概率密度函数,所以先画出以粒度为自变量,体
积为因变量的条形图,从图可以看出,它们之间可能服从指数分布形式,因此给出指数
分布的函数形式F(X)=1-exp(-x/a),因为是指数形式,所以用非线性拟合得到a的值为
8.1451,并画出拟合图形,跟用原始数据画出的图形类似,并求出了置信度为95%的图
形区间。
M文件如下:x=[0 0.1 0.2 0.5 1.0 1.5 2.0 2.5 3.0 4.0 5.0 6.0 8.0 10.0 15.0];
y=[0
0.0114517
0.0249242
0.0710677
0.1316940
0.1754800 0.1923210
0.2297070
0.2317280
0.3583700
0.3853150
0.3950830
0.6645340
0.8342880
1.0000000
];
y=y';
bar(x,y); %大致判断服从什么分布
beta0=0.1;
[beta r j]=nlinfit(x,y,@funex,beta0);
Beta
ci=nlparci(beta,r,j)
nlintool(x,y,@funex,beta0,0.05,'粒度x,体积y ') 其中funex为构建的可能的函数形式
运行结果如下:beta =
8.1451
流程图如下:
0.4
对于不同的h值我们用switch开关语句来运行。
M文件入下:x=[0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8];
y=[0 2.1220 3.0244 3.2568 3.1399 2.8579 2.5140 2.1639 1.8358];
h=input('请输入h的值,h为0.1或0.2或0.4');
switch h
case 0.4
s=h/2*(y(1)+2*y(5)+y(9));
case 0.2
s=h/2*(y(1)+2*(y(3)+y(5)+y(7))+y(9));
case 0.1
s=h/2*(y(1)+2*(y(2)+y(3)+y(4)+y(5)+y(6)+y(7)+y(8))+y(9));
otherwise
disp('输入的h有误 ');
end
s
流程图如下:
a.h=0.1的时候s=1.9997
b.h=0.2的时候s= 1.9192
c.h=0.4的时候s= 1.6231。