-1--1-第三章上机习题用你所熟悉的的计算机语言编制利用QR分解求解线性方程组和线性最小二乘问题的通用子程序,并用你编制的子程序完成下面的计算任务:(1)求解第一章上机习题中的三个线性方程组,并将所得的计算结果与前面的结果相比较,说明各方法的优劣;(2)求一个二次多项式+bt+cy=at2,使得在残向量的2范数下最小的意义下拟合表3.2中的数据;表3.2ti-1-0.75-0.500.250.50.75yi10.81250.7511.31251.752.3125(3)在房产估价的线性模型111122110xaxaxaxy中,1121,,,aaa分别表示税、浴室数目、占地面积、车库数目、房屋数目、居室数目、房龄、建筑类型、户型及壁炉数目,y代表房屋价格。现根据表3.3和表3.4给出的28组数据,求出模型中参数的最小二乘结果。(表3.3和表3.4见课本P99-100)解分析:(1)计算一个Householder变换H:由于TTvvIwwIH2,则计算一个Householder变换H等价于计算相应的v、。其中)/(2,||||12vvexxvT。在实际计算中,为避免出现两个相近的数出现的情形,当01x时,令212221||||)(-xxxxvn;为便于储存,将v规格化为1/vvv,相应的,变为)/(221vvvT为防止溢出现象,用||||/xx代替(2)QR分解:利用Householder变换逐步将nmAnm,转化为上三角矩阵AHHHnn11,则有-2-0RQA,其中nHHHQ21,:),:1(nR。在实际计算中,从nj:1,若mj,依次计算)),:((jmjAx对应的)1()1()~(kmkmjH即对应的jv,j,将)1:2(jmvj储存到),:1(jmjA,j储存到)(jd,迭代结束后再次计算Q,有~001jjjHIH,nHHHQ21(mn时1-21nHHHQ)(3)求解线性方程组bAx或最小二乘问题的步骤为i计算A的QR分解;ii计算bQcT11,其中):1(:,1nQQiii利用回代法求解上三角方程组1cRx(4)对第一章第一个线性方程组,由于R的结果最后一行为零,故使用前代法时不计最后一行,而用运行结果计算84x。运算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);ifsigma==0belta=0;elsealpha=sqrt(x(1)^2+sigma);ifx(1)=0v(1)=x(1)-alpha;elsev(1)=-sigma/(x(1)+alpha);endbelta=2*v(1)^2/(sigma+v(1)^2);v=v/v(1,1);endend-3-2计算A的QR分解[Q,R]=QRfenjie(A)function[Q,R]=QRfenjie(A)[m,n]=size(A);Q=eye(m);forj=1:nifjm[v,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,:));forj=1:nifjmH=eye(m);temp=[1;A(j+1:m,j)];H(j:m,j:m)=H(j:m,j:m)-d(j)*temp*temp';Q=Q*H;endendend3解下三角形方程组的前代法x=qiandaifa(L,b)functionx=qiandaifa(L,b)n=length(b);forj=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*ones(82,1);14];n=length(A);-4-%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(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]=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: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');%改进的平方根法[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.75-0.500.250.50.75];y=[10.81250.7511.31251.752.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;A=xlsread('E:\temporary\专业课\数值代数\cha3_3_4.xls','A2:L29');y=xlsread('E:\temporary\专业课\数值代数\cha3_3_4.xls','M2:M29');[Q,R]=QRfenjie(A);-6-Q1=Q(:,1:12);c=Q1'*y;x=huidaifa(R,c);x=x'计算结果为(1)第一章上机习题中的三个线性方程组结果对比图依次为05010000.20.40.60.811.21.41.61.8QR分解050100-6-4-20246x108Gauss05010011111111PGauss050100-20246810QR分解050100-20246810Gauss050100-20246810PGauss050100-20246810平方根法050100-20246810改进的平方根法-7-02040-2000-1500-1000-50005001000150020002500QR分解02040-200-150-100-50050100150200Gauss02040-300-200-1000100200300400PGauss02040-5-4-3-2-1012345x107平方根法02040-80-60-40-200204060改进的平方根法以第二个线性方程组为例,比较各方法的运行速度。依次为QR分解,不选主元的Gauss消去法,列主元Gauss消去法,平方根法,改进的平方根法。Elapsedtimeis0.034588seconds.Elapsedtimeis0.006237seconds.Elapsedtimeis0.009689seconds.Elapsedtimeis0.030862seconds.Elapsedtimeis0.007622seconds.(2)二次多项式的系数为x=1.00001.00001.0000(3)房产估价的线性模型的系数为x=Columns1through62.07750.71899.68020.153513.67961.9868Columns7through12-0.9582-0.4840-0.07361.01871.44352.9028结果分析对第一章上机习题中的第二个线性方程组利用五种求解方法求解所需时间可知,不选主元的Gauss消去法,列主元Gauss消去法,改进的平方根法较快,所需时间大致在一个数量级,QR分解,平方根法,所需时间较慢,所需时间在一个数量级上。