1、- 1 - 1 -第三章上机习题用你所熟悉的的计算机语言编制利用 QR 分解求解线性方程组和线性最小二乘问题的通用子程序，并用你编制的子程序完成下面的计算任务：（1 ） 求解第一章上机习题中的三个线性方程组，并将所得的计算结果与前面的结果相比较，说明各方法的优劣；（2 ） 求一个二次多项式 ，使得在残向量的 2 范数下最小的意义下拟合表 3.2+btcy=a2中的数据；表 3.2ti -1 -0.75 -0.5 0 0.25 0.5 0.75yi 1 0.8125 0.75 1 1.3125 1.75 2.3125（3 ） 在房产估价的线性模型 1210xaxay中， 分别表示税、浴室数目、占
2、地面积、车库数目、房屋数目、居室数目、121,a房龄、建筑类型、户型及壁炉数目， 代表房屋价格。现根据表 3.3 和表 3.4 给出的 28 组数据，求出模型中参数的最小二乘结果。（表 3.3 和表 3.4 见课本 P99-100）解 分析：（1 ） 计算一个 Householder 变换 H：由于 ，则计算一个 Householder 变换 H 等价于计算相应的 。TTvIwIH2 v、其中 。)/(,|1exv在实际计算中，为避免出现两个相近的数出现的情形，当 时，令 ；01x21|)(-xvn为便于储存，将 规格化为 ，相应的， 变为v1/v)/(21vT为防止溢出现象，用 代替|/x（
3、2 ） QR 分解：利用 Householder 变换逐步将 转化为上三角矩阵 ，则有nmA, AHn1- 2 -，其中 ， 。0RQAnH21:),1(R在实际计算中，从 ，若 ，依次计算 对应的j:mj),:(jmAx )1()(kmjH即对应的 ， ，将 储存到 ， 储存到 ，迭代结束jvj)12(vj 1jjd后再次计算 ，有 ， （ 时 ）Q0jjjHI nHQ211-21nQ（3 ） 求解线性方程组 或最小二乘问题的步骤为bAxi 计算 的 QR 分解；ii 计算 ，其中cT1):1(,niii 利用回代法求解上三角方程组 cRx（4）对第一章第一个线性方程组，由于 R 的结果最后
4、一行为零，故使用前代法时不计最后一行，而用运行结果计算 。84运算 matlab 程序为 1 计算 Householder 变换 v,belta=house(x)function v,belta=house(x)n=length(x);x=x/norm(x,inf);sigma=x(2:n)*x(2:n);v=zeros(n,1);v(2:n,1)=x(2:n);if sigma=0belta=0;else alpha=sqrt(x(1)2+sigma);if x(1)=0v(1)=x(1)-alpha;elsev(1)=-sigma/(x(1)+alpha);endbelta=2*v(1)2
5、/(sigma+v(1)2);v=v/v(1,1);endend- 3 -2 计算 的 QR 分解 Q,R=QRfenjie(A)Afunction Q,R=QRfenjie(A)m,n=size(A);Q=eye(m);for j=1:nif jmv,belta=house(A(j:m,j);H=eye(m-j+1)-belta*v*v;A(j:m,j:n)=H*A(j:m,j:n);d(j)=belta;A(j+1:m,j)=v(2:m-j+1);endendR=triu(A(1:n,:);for j=1:nif jmH=eye(m);temp=1;A(j+1:m,j);H(j:m,j:m
6、)=H(j:m,j:m)-d(j)*temp*temp;Q=Q*H;endendend3 解下三角形方程组的前代法 x=qiandaifa(L,b)function x=qiandaifa(L,b)n=length(b);for j=1:n-1b(j)=b(j)/L(j,j);b(j+1:n)=b(j+1:n)-b(j)*L(j+1:n,j);endb(n)=b(n)/L(n,n);x=b; end4 求解第一章上机习题中的三个线性方程组 ex3_1clear;clc;%第一题A=6*eye(84)+diag(8*ones(1,83),-1)+diag(ones(1,83),1);b=7;15
7、*ones(82,1);14;- 4 -n=length(A);%QR 分解Q,R=QRfenjie(A);c=Q*b;x1=huidaifa(R(1:n-1,1:n-1),c(1:n-1);x1(n)=c(n)-R(n,1:n-1)*x1;%不选主元 Gauss 消去法L,U=GaussLA(A);x1_1=Gauss(A,b,L,U);%列主元 Gauss 消去法L,U,P=GaussCol(A);x1_2=Gauss(A,b,L,U,P);%解的比较figure(1);subplot(1,3,1);plot(1:n,x1);title(QR 分解);subplot(1,3,2);plot
8、(1:84,x1_1);title(Gauss);subplot(1,3,3);plot(1:84,x1_2);title(PGauss);%第二题第一问A=10*eye(100)+diag(ones(1,99),-1)+diag(ones(1,99),1);b=round(100*rand(100,1);n=length(A);%QR 分解tic;Q,R=QRfenjie(A);c=Q*b;x2=huidaifa(R,c);toc;%不选主元 Gauss 消去法tic;L,U=GaussLA(A);x2_1=Gauss(A,b,L,U);toc;%列主元 Gauss 消去法tic;L,U,P
9、=GaussCol(A);x2_2=Gauss(A,b,L,U,P);toc;%平方根法tic;L=Cholesky(A);x2_3=Gauss(A,b,L,L);toc;%改进的平方根法tic;L,D=LDLt(A);x2_4=Gauss(A,b,L,D*L);toc;%解的比较figure(2);subplot(1,5,1);plot(1:n,x2);title(QR 分解);subplot(1,5,2);plot(1:n,x2_1);title(Gauss);subplot(1,5,3);plot(1:n,x2_2);title(PGauss);subplot(1,5,4);plot(1
10、:n,x2_3);title(平方根法);subplot(1,5,5);plot(1:n,x2_4);title(改进的平方根法);- 5 -%第二题第二问A=hilb(40);b=sum(A);b=b;n=length(A);Q,R=QRfenjie(A);c=Q*b;x3=huidaifa(R,c);%不选主元 Gauss 消去法L,U=GaussLA(A);x3_1=Gauss(A,b,L,U);%列主元 Gauss 消去法L,U,P=GaussCol(A);x3_2=Gauss(A,b,L,U,P);%平方根法L=Cholesky(A);x3_3=Gauss(A,b,L,L);%改进的
11、平方根法L,D=LDLt(A);x3_4=Gauss(A,b,L,D*L);%解的比较figure(3);subplot(1,5,1);plot(1:n,x3);title(QR 分解);subplot(1,5,2);plot(1:n,x3_1);title(Gauss);subplot(1,5,3);plot(1:n,x3_2);title(PGauss);subplot(1,5,4);plot(1:n,x3_3);title(平方根法);subplot(1,5,5);plot(1:n,x3_4);title(改进的平方根法);5 求解二次多项式 ex3_2clear;clc;t=-1 -0
12、.75 -0.5 0 0.25 0.5 0.75;y=1 0.8125 0.75 1 1.3125 1.75 2.3125;A=ones(7,3);A(:,1)=t.2;A(:,2)=t;Q,R=QRfenjie(A);Q1=Q(:,1:3);c=Q1*y;x=huidaifa(R,c)6 求解房产估价的线性模型 ex3_3clear;clc;- 6 -A=xlsread(E:temporary专业课数值代数cha3_3_4.xls,A2:L29);y=xlsread(E:temporary专业课数值代数cha3_3_4.xls,M2:M29);Q,R=QRfenjie(A);Q1=Q(:,1
13、:12);c=Q1*y;x=huidaifa(R,c);x=x计算结果为（1）第一章上机习题中的三个线性方程组结果对比图依次为0 50 10000.20.40.60.811.21.41.61.8 QR分分0 50 100-6-4-20246x 108 Gauss0 50 10011111111 PGauss0 50 100-20246810 QR分分0 50 100-20246810 Gauss0 50 100-20246810 PGauss0 50 100-20246810分分分分0 50 100-20246810分分分分分分分- 7 -0 20 40-2000-1500-1000-5000
14、5001000150020002500QR分分0 20 40-200-150-100-50050100150200 Gauss0 20 40-300-200-1000100200300400 PGauss0 20 40-5-4-3-2-1012345x 107分分分分0 20 40-80-60-40-200204060分分分分分分分以第二个线性方程组为例，比较各方法的运行速度。依次为 QR 分解，不选主元的 Gauss 消去法，列主元 Gauss 消去法，平方根法，改进的平方根法。Elapsed time is 0.034588 seconds.Elapsed time is 0.006237
15、 seconds.Elapsed time is 0.009689 seconds.Elapsed time is 0.030862 seconds.Elapsed time is 0.007622 seconds.（2）二次多项式的系数为x =1.00001.00001.0000（3）房产估价的线性模型的系数为x =Columns 1 through 62.0775 0.7189 9.6802 0.1535 13.6796 1.9868Columns 7 through 12-0.9582 -0.4840 -0.0736 1.0187 1.4435 2.9028结果分析对第一章上机习题中的第二个线性方程组利用五种求解方法求解所需时间可知，不选主元的 Gauss 消去法，列主元 Gauss 消去法，改进的平方根法较快，所需时间大致在一个数量级，QR 分解，平方根法，所需时间较慢，所需时间在一个数量- 8 -级上。
