基于MATLAB的数值分析编程上机作业1

基于MATLAB的数值分析编程上机作业（1）

二、用Householder变换法求矩阵A的正交分解A＝QR。 1、基本原理： 任一实列满秩的m×n矩阵A，可以分解成两个矩阵的乘积，即A＝QR，其中Q是具有法正交列向量的m×n矩阵，R是非奇异的n阶上三角阵。 算法：

1、 输入n阶矩阵A

2、 对A?k:n,k?，求Househoulder初等反射阵的?k,u(k)。k?1,2,?,n?1 3、 计算上三角阵R，仍然存储在A

HkA(k)?A(k?1)

j?k,k?1,?,n

1?n?k?(k)?tj???uiaij??k?i?k? (k?1)(k)aij?aij?tjui?k?

4、计算正交阵Q

Q(1)?IQHk?Qi?1,2,?,n

(k)(k?1)

ti?

1?k?qj?kn(k)iju?jk?

j?k,k?1,?,n(k?1)(k)qij?qij?tiu?jk?5、按要求输出 结束

2、程序代码： function [Q,A]=qrhh(A) %QRHH 用Householder变换法对n阶矩阵A作正交分解A=QR； %程序功能：函数qrhh用Householder变换法对矩阵A作正交分解A=QR； %输入：n阶矩阵A； %输出：[Q,A]。Q是具有法正交列向量的n阶矩阵， % A（即R）是非奇异的n阶上三角阵，仍用输入的矩阵A存储。 %引用函数： % holder2；示例 [p,u]=holder2(x); %使用举例： % [Q,R]=qrhh(A) %变量说明： % A － 输入的n阶矩阵，同时用于存储上三角阵R； % n － 矩（方）阵A的阶数；

4

基于MATLAB的数值分析编程上机作业（1）

% Q － Q是具有法正交列向量的n阶矩阵； % p,u － 向量A（k:n,k），对应初等反射阵的ρ,u % k，jj，ii － 循环变量； % t1 － 计算上三角阵R的系数tj； % t2 － 计算正交矩阵Q的系数ti； [n,n]=size(A); %求矩（方）阵A的阶数； Q=eye(n); %构造正交矩阵Q(1)=I； for k=1:n-1 [p,u]=holder2(A(k:n,k)); %向量A（k:n,k），对应初等反射阵的ρ,u for jj=k:n %计算上三角阵R（仍存贮于A） t1=dot(u,A(k:n,jj))/p; %利用向量内积求和 A(k:n,jj)=A(k:n,jj)-t1*u; end for ii=1:n %计算正交矩阵Q t2=dot(u,Q(ii,k:n))/p; %利用向量内积求和 Q(ii,k:n)=Q(ii,k:n)-t2*u'; end end

3、计算实例： >> A=[2 1 0; 1 3 1; 0 1 4] A = 2 1 0 1 3 1 0 1 4 >> [q,r]=qrhh(A) q = -0.8944 0.4082 0.1826 -0.4472 -0.8165 -0.3651 0 -0.4082 0.9129 r = -2.2361 -2.2361 -0.4472 0 -2.4495 -2.4495 0 -0.0000 3.2863 >> q*r ans = 2.0000 1.0000 -0.0000 1.0000 3.0000 1.0000 0 1.0000 4.0000

5