数值分析报告Word下载.docx
《数值分析报告Word下载.docx》由会员分享，可在线阅读，更多相关《数值分析报告Word下载.docx（10页珍藏版）》请在冰豆网上搜索。
A=[12123
34521
54264
22590];
v=[1111]'
;
u（:
1）=v（:
1）;
fori=1:
100
v（:
i+1）=A*u（:
i）;
ifabs（max（v（:
i+1））-max（v（:
i）））<
10^-4
break
end
u（:
i+1）=v（:
i+1）/max（v（:
i+1））;
disp（u（:
i））;
disp（max（v（:
i）））;
k=u（:
i）'
*A*u（:
i）/（u（:
*u（:
disp（k）;
[x,c]=eig（A）;
disp（c）;
disp（x）;
disp（i）;
3.结果比较和结论
初始矩阵
A=[12123
34521
54264
22590];
初始向量
v=[1111]'
幂法求得特征向量
12.5147
14.9140
25.6935
39.6330
归一化后特征向量
0.3158
0.3763
0.6483
1.0000
列向量最大值近似主特征值
39.6330
Rayleigh商求出主特征值
用eig（）函数算出的特征值和特征向量
39.6331000
0-18.2401+7.4985i00
00-18.2401-7.4985i0
0008.8471
-0.2450-0.0471+0.0965i-0.0471-0.0965i-0.0574
-0.29190.1147-0.0939i0.1147+0.0939i-0.1747
-0.50290.2862-0.1187i0.2862+0.1187i0.1535
-0.7758-0.9330-0.93300.9709
达到精度要求所需次数：
17
结论：
可以看出初始列向量经过多次迭代后，用幂法求出的特征值和用eig（）函数求出的A的特征值，满足计算精度在，并且特征向量也具有数乘关系。
另外发现取特征向量v的列最大值近似的特征值和用Rayleigh商求出的特征值相同。
部分二
（QR分解与特征值）
1.
幂法求主特征值思路
2.源程序
Household函数代码
function[QR]=hshd（a）
n=length（a）;
h=eye（size（a））;
n-1
if（a（i,i）>
0）
sgn=1;
else
sgn=-1;
x=sgn*（sum（a（i:
n,i）.^2））^0.5;
u=a（i:
n,i）+x*eye（n-i+1,1）;
ro=norm（u,2）^2/2;
H=eye（n-i+1）-u*u'
/ro;
H1=eye（size（a））;
H1（i:
n,i:
n）=H;
a=H1*a;
h=H1*h;
Q=inv（h）;
R=a;
return
QR分解代码
function[]=qrfenjie（a）
50
[Q,R]=hshd（a）;
a=R*Q;
disp（a）
3.结果与结论
假设a矩阵为如下3*3矩阵
a=[221
014
302]
计算出a的特征值为5，+2.2361i，-2.2361i三个
用QR分解逆序相乘法后得到矩阵
5.0000-1.06901.2344
-0.00000.35712.7630
-0.0000-1.8558-0.3571
三对角线上要么是实特征值，要么是两两共轭的虚特征值。
对于上述矩阵，有实特征5和2*2矩阵构成的虚特征值，求得虚特征值正好为2.2361i。
如果最后求得的矩阵为上三角矩阵，那么特征值就是对角线上元素，而且特征值均为实数。
用eig（）函数求得矩阵特征值如下
5.0000
0.0000+2.2361i
0.0000-2.2361i
结论
可以看到两种算法在10^-4次精度上特征值是相等的，这也验证了QR算法的正确性。
部分三
（矩阵特征值与飞机横侧向稳定性分析）
根据飞机平飞时横向动力学方程
上述方程可以看做状态方程X’=AX+B的形式，A为状态空间，如果能求得A的特征值，根据飞行力学的知识可以得到飞机横向飞行的模态。
所以求出A的特征值是解决这一问题的基本步骤。
在前两个部分已经提到，幂法和QR分解法都能求出A的某一特征值及全部特征值，对于本问题，不能采用幂法因为1.根据经验，A矩阵有复数特征值，并且特征值可能为负。
2.要求出所有飞行模态必须求出所有特征值。
基于上述两点，决定采用QR分解来计算特征值。
对于某型号运输机得知如下A矩阵
-0.11500-1.00000.04150
-7.9283-1.08300.531100
1.6327-0.0108-0.161800
01.0000000
001.000000
经过QR分解50步迭代之后得到
-0.00641.0071-0.29290.3413-0.1398
-1.7685-0.4813-5.96105.18330.4769
0.1014-0.0485-0.86420.5776-1.2995
-0.00000.00000.0000-0.0079-0.1203
0000.0000-0.0000
对于3阶顺序主子式
-0.00641.0071-0.2929
-1.7685-0.4813-5.9610
0.1014-0.0485-0.8642
继续做QR分解得到
-0.43171.39354.7202
-1.37350.29973.5888
-0.00120.0014-1.2199
可以得到A矩阵的5个特征值近似为
-1.2199
-0.0079
-0.0660+1.3343i
-0.0660-1.3343i
事实上，用eig（）函数求出的特征值为
0
-1.2233
-0.0079
-0.0643+1.3359i
-0.0643-1.3359i
可见QR算法具有一定的精度
那么根据飞行力学知识可以知道特征值分别对应的飞行模态
0航向中立
-1.2233滚转模态
-0.0079螺旋模态
-0.0643+1.3359i荷兰滚
-0.0643-1.3359i荷兰滚
得到飞行模态之后就可以分析它的稳定性，阻尼，超调等等更为细致的问题。
为分析影响飞机横侧向稳定因素提供了基础。
小结
矩阵特征值这一章节主要讲的就是如何避开解高阶的特征方程，寻求另一种可机器运算的算法来求解特征值。
通过分析矩阵特性最终找到QR分解->
逆序相乘->
迭代的算法。
这种算法运算方法单一，重复度高，计算量大非常适合计算机运算。
避免了求解高次方程的困难。
数值分析主要的思路就是用数值解代替精确解，从而寻找在一定精度内可以接受的数值算法。
这是一种重要的简化问题的思路。
在工程上，往往不要求精确解，所以不需要那些复杂但是精确的算法，更多时候可能也根本没有精确算法。
当计算机能力越来越强，数值解的应用更为广泛。
数值分析这门课提到最多的思路就是如何设计算法寻找一种近似。
它看似与数学的思想背道而驰，却恰恰很好地体现了哲学的精神。
世界上没有绝对的相同，也就是说求出精确解实际上却是不可能达到的。
精确解和实际环境脱节，应用价值不大，而数值分析的思想不仅可以简化问题而且也能保证精度要求。
升级会员 