《矩阵论上机作业》由会员分享，可在线阅读，更多相关《矩阵论上机作业（13页珍藏版）》请在金锄头文库上搜索。
1、 矩阵论上机作业班级：070921学号：07092002姓名：黄涛完成时间： 时间：2011 年 10 月 21 号一、 doolittle 方法定义： 设A Cn n . 如果A 可分解成A = LR,其中L 是对角元素为1 的下三角矩阵(称为单位下三角矩阵) , R 是上三角矩阵, 则称之为A 的Doolit tle分解5 .根据定义获得的计算公式r1 j = a1 j ( j = 1 ,2 , , n) ,li1 = ai1 / r11 ( i = 2 ,3 , , n) ,rkj = akj - k- 1t = 1lkt3 rtj ( j = k , k + 1 , , n; k =2
2、，3 , , n) ,lik = 1/ rkk3 aik - k- 1t = 1lit3 rtk ( i = k + 1 , ,n; k = 2 ,3 , , n) .由定义知,Doolit tle 分解是在n 行n 列的矩阵上进行的,由公式知道当i , j , k 较大时, lik , rkj 的计算量是比较大的,当矩阵的行数与列数不相同时是无法进行Doolittle 分解的, 为克服这些缺陷, 从分块的角度出发研究Doolit tle 分解. 使Doolit tle 分解能够适合可分任何矩阵。 程序代码#include #include #define N 20 int main() in
3、t n;coutn;static double aNN;int i,j;for(i=0;iaij; /input Acout“A=“endl;for(i=0;in;i+)for(j=0;jn;j+)printf(“%-8.2f“,aij);coutendl;/output Astatic double lNN,uNN;int k; /step numberfor(i=0;in;i+) lii=1;for(k=0;kn;k+) /用 k 步实现，第 k 步计算：U 的第 k 行，L 的第 k 列for(j=k;jn;j+) /U 的第 k 行ukj=akj;for(i=0;i=k-1;i+)uk
4、j-=lki*uij;for(i=k+1;in;i+) /L 的第 k 列lik=aik;for(j=0;j=k-1;j+)lik-=lij*ujk;lik/=ukk;cout“L=“endl; for(i=0;in;i+)for(j=0;jn;j+)printf(“%-8.2f“,lij);coutendl; /output Lcout“U=“endl; for(i=0;in;i+)for(j=0;jn;j+)printf(“%-8.2f“,uij);coutendl; /output Ureturn 0;案例1、 crout分解%本函数将一个满秩方阵按Crout方式分解function L
5、,U=Crout(A)b=size(A);%b(1)行%b(2)列n=b(1);%这里只处理n*n的非奇异矩阵%错误检查if b(1)=b(2)%非方阵错误error(MATLAB:Crout:Input Matrix should be a Square matrix. See Crout.);endif n=rank(A)%非满秩矩阵错误error(MATLAB:Crout:Input Matrix should be FULL RANK. See Crout.);end %初始化L=zeros(n,n);U=zeros(n,n); %U中的主对角线元素均为1for i=1:nU(i,i)
6、=1;end%开始计算for k=1:nfor i=k:n %L中计算的方式是行为外循环，列为内循环temp_sum=0; %临时和for t=1:k-1temp_sum=temp_sum+L(i,t)*U(t,k);endL(i,k)=A(i,k)-temp_sum; %计算L的第k列元素endfor j=k+1:n %U中计算的方式是列为外循环，行为内循环temp_sum=0; %临时和for t=1:k-1temp_sum=temp_sum+L(k,t)*U(t,j); endU(k,j)=(A(k,j)-temp_sum)/L(k,k);%计算U的第k行元素endendend %end
7、 of function三、cholesky 分解A=input(A=) m,n=size(A); if m=n %判断输入的矩阵是不是方阵disp(输入的矩阵不是方阵，请重新输入);return; endd=eig(A); %根据方阵的特征值判定是不是正定矩阵for i=1:nif d(i)=0disp(输入的矩阵不是正定矩阵，请重新输入);return;elsebreak;end end disp(输入的矩阵可以进行 Cholesky 分解); %如果是正定矩阵，可以进行下面的分解操作 S(1,1)=sqrt(A(1,1); %确定第 1 列元素for i=2:nS(i,1)=A(i,1)
8、/S(1,1); end sum1=0; for j=2:n-1 %确定第 j 列元素for k=1:j-1sum1=sum1+S(j,k)2;endS(j,j)=sqrt(A(j,j)-sum1);sum1=0;for i=j+1:nfor k=1:j-1sum1=sum1+S(i,k)*S(j,k);endS(i,j)=(A(i,j)-sum1)/S(j,j);endsum1=0; end for k=1:n-1 %确定第 n 列元素sum1=sum1+S(n,k)2; endS(n,n)=sqrt(A(n,n)-sum1); S P=S*S %验证分解是否正确 disp(请验证 Chol
9、esky 分解正确（是/否）？);end任取 A=，则其仿真结果为：A =1 6 9 6 5 33 5 6 4 2 43 7 9 0 4 35 4 6 5 4 93 8 5 7 6 89 3 6 4 9 6输入的矩阵可以进行 Cholesky 分解ans = Columns 1 through 4 1.0000 3.0000 3.0000 5.0000 0 0 - 2.0000i 0 - 1.0000i 0 -10.0000i 0 0 1.0000 1.0000 0 0 0 8.8882 0 0 0 0 0 0 0 0 Column 5 through 6 3.0000 9.0000 0-12
10、.5000 i 0 -28.5000i3.5000 6.0000 12.7697 38.7593 0-4.6795i 0 -37.8279i 0 25.0981 P =1.0e+003 *0.0010 0.0030 0.0030 0.0050 0.0030 0.00900.0030 0.0130 0.0110 0.0350 0.0340 0.08400.0030 0.0110 0.0110 0.0260 0.0250 0.06150.0050 0.0350 0.0260 0.2050 0.2570 0.68050.0030 0.0340 0.0250 0.2570 0.3626 1.07690
11、.0090 0.0840 0.0615 0.6805 1.0769 4.4924请验证 Cholesky 分解正确（是/否）？4、吉文斯方法的 QR 分解function Q,R=qrgv(A)% 基于 Givens 变换，将方阵 A 分解为 A=QR，其中 Q 为正交矩阵，R 为上三角阵% 参数说明% A：需要进行 QR 分解的方阵% Q：分解得到的正交矩阵% R：分解得到的上三角阵% Q,R=qr(A) % 调用 MATLAB 自带的 QR 分解函数进行验证% q,r=qrgv(A) % 调用本函数进行 QR 分解% q*r-A % 验证 A=QR% q*q % 验证 q 的正交性% no
12、rm(q) % 验证 q 的标准化，即二范数等于 1% 线性代数基础知识% 1.B=P*A*inv(P)，称 A 与 B 相似，相似矩阵具有相同的特征值% 2.Q*Q=I，称 Q 为正交矩阵，正交矩阵的乘积仍为正交矩阵% by dynamic of Matlab 技术论坛% see also % contact me % 2010-01-17 22:51:18%n=size(A,1); R=A; Q=eye(n); for i=1:n-1for j=2:n-i+1x=R(i:n,i);rt=givens(x,1,j);r=blkdiag(eye(i-1),rt);Q=Q*r;R=r*R;end
13、 end function R,y=givens(x,i,j) % 求解标准正交的 Given 变换矩阵 R，使用 Rx=y，其中 y(j)=0，y(i)=sqrt(x(i)2+x(j)2) % % 参数说明 % x：需要进行 Givens 变换的列向量 % i：变为 sqrt(x(i)2+x(j)2)的元素下标 % j：变为0的元素的下标 % R：Givens 变换矩阵 % y：Givens 变换结果 % % 实例说明 % x=1 3 5 9 6; % 将3等效到9上 % R,y=givens(x,4,2) % 注意3的下标为2，9的下标为4 % R*x-y % 验证 Rx=y % R*R % 验证正交性 % norm(R) % 验证标准性，就是范数为1 % % 关于 Givens 变换说明 % 1.Givens 矩阵是标准正交矩阵，也叫平面旋转矩阵，它是通过坐标旋转的原理将元素 j 的数值等效到元素 i 上 % 2.Givens 变换每次只能将一个元素变为0，而 Householder 变
