数值线性代数课程设计

《数值线性代数》课程设计
题目
学生
指导教师
徐州师范大学
课程设计任务书
院班学生
课程设计课题：
矩阵的QR 分解
课程设计任务要求（包括课题来源、类型、目的和意义、基本要求、参考资料等）：
✧ 来源与意义:
本课题来源于教材第三章第三节——正交化算法，目的是通过QR 分解求解方程Ax=b ，以求解最小二乘问题。

✧ 基本要求: 1. 要求自编程序；
2. 掌握编程思想，学会一门编程语言；
3. 报告要有较强的理论分析；有较强说服力的数据表或图像；
4. 对结果进行分析；
5. 给出相应结论。

✧ 参考资料:
《数值线性代数》：徐树芳 高立 张平文
指导教师签字：
一、 问题提出：
在探究最小二乘问题时，我们得出结论，求解最小二乘问题等价于求
2min ||||T T Q Ax Q b 的问题，其中Q 为正交矩阵，因此，我们希望通过选取适当的正交矩阵Q ，使原问题转化为较易求解的最小二乘问题。

这里，我们将讨论对任意A 的QR 分解的方法。

二、理论基础
1、QR 分解定理：设()m n A R m n ⨯∈≥，则A 有QR 分解：
0R A Q ⎡⎤=⎢⎥⎣⎦
，
其中m m Q R ⨯∈是正交矩阵，n n R R ⨯∈是具有非负对角元的上三角阵；而且当m n =且A 非奇异时，上述分解还是唯一的。

3）说明：在这个定理的基础上，我们就可以放心的研究任意矩阵的QR 分解，因为它总是存在的。

2、Householder 变换：
1) 定义：设n R ω∈满足21ω=，定义n n H R ⨯∈为
2T H I ωω=-，
则称H 为Householder 变换。

2）作用：H x ⨯可以将12
n x x x x ⎡⎤
⎢⎥
⎢⎥=⎢⎥
⎢⎥
⎣⎦ 化成(1)100x x ⎡⎤⎢⎥⎢⎥=⎢⎥⎢⎥⎢⎥⎣⎦
的形式。

3）说明：利用Householder 变换，我们可以将A 的对角线下元素逐列的消为0。

3、Givens 变换
1）定义：形如G=11i c s
j s c i j
⎛⎫
⎪
⎪ ⎪ ⎪ ⎪
⎪- ⎪
⎪ ⎪⎝⎭
的矩阵称为Givens 矩阵。

2）作用：G x ⨯可以将12
n x x x x ⎡⎤
⎢⎥
⎢⎥=⎢⎥
⎢⎥
⎣⎦
化成0
0i x j ⎡
⎤⎢⎥⎢⎥⎢⎥=⎢⎥
⎢⎥⎢⎥
⎢⎥
⎢⎥⎣⎦
的形式。

3）说明：利用Givens 变换，我们可以将A 的对角线下元素逐个消为0。

三、实验内容
模型一：利用householder变换对矩阵A进行QR分解。

1.用MATLAB建立m文件house.m,程序如下：function [v,b]=house(x)
n=length(x);
v=zeros(n,1);
eta=norm(x,inf);
x=x/eta;
v(1)=1;
v(2:n)=x(2:n);
sig=x(2:n)'*x(2:n);
if (sig==0)
b==0;
else
a=(x(1)*x(1)+sig)^(1/2);
if (x(1)<=0)
v(1)=x(1)-a;
else
v(1)=-sig/(x(1)+a);
end
b=2*v(1)*v(1)/(sig+v(1)*v(1));
v=v/v(1);
end
2.用MATLAB建立m文件：qr1.m
A=input('A=');
[m,n]=size(A);
for j=1:n
if j<=m
[v,b]=house(A(j:m,j))
A(j:m,j:n)=(eye(m-j+1)-b*v*v')*A(j:m,j:n); d(j)=b;
A(j+1:m,j)=v(2:m-j+1);
end
end
R=A(1:n,1:n);
R=triu(R);
Q=eye(m);
for j=1:n
A(j,j)=1;
end
for j=1:n
H=eye(m);
H(j:m,j:m)=eye(m-j+1)-d(j)*A(j:m,j)*A(j:m,j)';
Q=Q*H;
end
Q
R
Q*Q' %这里是为了检验Q是否为正交矩阵
3.结果：利用Householder变换，我们可以将A的对角线下元素逐列的消为0，从而得到R，每一步中的Householder变换矩阵H的乘积就是我们所求的Q。

（程序运行的结果数据见附录一）。

模型二：利用givens变换对矩阵A进行QR分解。

1.用MATLAB建立m文件givens.m,程序如下：
function[c,s]=givens (a,b)
if (b==0)
c=1;
s=0;
else if(abs(b)>abs(a))
r=a/b;
s=1/sqrt(1+r*r);
c=s*r;
else
r=b/a;
c=1/sqrt(1+r*r);
s=c*r;
end
end
2. 用MATLAB建立m文件：qr2.m
A=input('A=');
[m,n]=size(A);
Q=eye(m);
for j=1:n
for i=j+1:m
[c,s]=givens(A(j,j),A(i,j));
G=eye(m);
G(j,j)=c;
G(i,i)=c;
G(j,i)=s;
G(i,j)=-s;
A=G*A;
Q=Q*G;
end
end
R=A
QR=A
3.结果：利用Givens变换，我们可以将A的对角线下元素逐个消为0，从而得到R，每一步中的Givens变换矩阵G的乘积就是我们所求的Q。

（程序运行的结果数据见附录二）。

四、结果分析
1.两种方法都可以求得A的QR分解，但从程序的运行代价上来看，Givens 变换需要的存储更高，速度更慢。

2.两种方法得到的Q和R可能不同，但是R一定是上三角矩阵，Q一定是正交矩阵。

3.我们也可以试着探寻基于其他方法的QR分解。

五、结论
Householder变换和Givens变换都能够较好的实现QR分解，在这个基础上，我们可以更好更方便的解决最小二乘问题。

但是它们所消耗的内存空间都比较大，我们也可以试着探寻更好的方法，以节省内存。

六、参考文献
《数值线性代数》：徐树芳高立张平文
附录一：
qr1.m的运行结果：
任意键入一矩阵A=[1,5,6;1,3,5;5,8,9;8,5,6]，按enter键，则有结果v =
1.0000
-0.1171
-0.5855
-0.9368
b =
0.8952
v =
1.0000
-1.5136
-0.2854
b =
0.5931
v =
1.0000
-5.6333
b =
0.0611
Q =
0.1048 0.6551 -0.0553 -0.7462 0.1048 0.3302 0.9075 0.2374 0.5241 0.5141 -0.3925 0.5540 0.8386 -0.4445 0.1388 -0.2827
R =
9.5394 9.2249 10.9022
0 6.1564 7.5415
0 0 1.5061
ans =
1.0000 -0.0000 0.0000 -0.0000 -0.0000 1.0000 -0.0000 0.0000 0.0000 -0.0000 1.0000 -0.0000 -0.0000 0.0000 -0.0000 1.0000
附录二：
qr1.m的运行结果：
c =
0.7071
s =
0.7071
c =
0.2722 s =
0.9623
c =
0.5447 s =
0.8386
c =
0.3974 s =
0.9177
c =
0.5781 s =
0.8160
c =
0.8976
s =
-0.4408
R =
9.5394 9.2249 10.9022
0 -6.1564 -7.5415
0 0 -1.5061
0 0 -0.0000
Q =
0.1048 -0.3302 0.7427 -0.5730
-0.1048 0.6551 0.6243 0.4125
-0.5241 0.5141 -0.1984 -0.6493
-0.8386 -0.4445 0.1388 0.2827
ans =
1.0000 0.0000 0.0000 -0.0000
0.0000 1.0000 0 0.0000
0.0000 0 1.0000 0.0000
-0.0000 0.0000 0.0000 1.0000
- 10 -。