數(shù)值線性代數(shù)第二版徐樹(shù)方高立張平文上機(jī)習(xí)題第三章實(shí)驗(yàn)報(bào)告

1、QR 分解求解線性方程組和線性最小二乘問(wèn)題的y=at +bt+c-11x0ayHvx1v規(guī)格化為 vx/|x|Am2-0.75 0.8125a1x11Ix |x|20時(shí)，令 vv/v1，相應(yīng)的，代替n,，使得在殘向量的-0.5 0.75a x2,a2,2wwTe1,-(x21變?yōu)閙 n轉(zhuǎn)化為上三角矩陣2 范數(shù)下最小的意義下擬合表0 12,a11分別表示稅、浴室數(shù)目、占地面積、車(chē)庫(kù)數(shù)目、房屋數(shù)目、居室數(shù)目、房IQR 分解求解線性方程組和線性最小二乘問(wèn)題的y=at +bt+c-11x0ayHvx1v規(guī)格化為 vx/|x|Am2-0.75 0.8125a1x11Ix |x|20時(shí)，令 vv/v1，相

2、應(yīng)的，代替n,，使得在殘向量的-0.5 0.75a x2,a2,2wwTe1,-(x21變?yōu)閙 n轉(zhuǎn)化為上三角矩陣2 范數(shù)下最小的意義下擬合表0 12,a11分別表示稅、浴室數(shù)目、占地面積、車(chē)庫(kù)數(shù)目、房屋數(shù)目、居室數(shù)目、房I2/(2x12v2 /(HnHn3.20.251.3125a11x11vvTvTv)x2n)|x|2110.51.75，則計(jì)算一個(gè) Householder變換 H 等價(jià)于計(jì)算相應(yīng)的。；vTv)H1A，則有0.752.3125、v。第三章上機(jī)習(xí)題用你所熟悉的的計(jì)算機(jī)語(yǔ)言編制利用通用子程序，并用你編制的子程序完成下面的計(jì)算任務(wù)：（1）求解第一章上機(jī)習(xí)題中的三個(gè)線性方程組， 并將

3、所得的計(jì)算結(jié)果與前面的結(jié)果相比較，說(shuō)明各方法的優(yōu)劣；（2）求一個(gè)二次多項(xiàng)式中的數(shù)據(jù)；表 3.2 tiyi（3）在房產(chǎn)估價(jià)的線性模型y中，齡、建筑類型、 戶型及壁爐數(shù)目， 代表房屋價(jià)格。 現(xiàn)根據(jù)表 3.3 和表 3.4 給出的 28 組數(shù)據(jù)，求出模型中參數(shù)的最小二乘結(jié)果。（表 3.3 和表 3.4 見(jiàn)課本 P99-100）解 分析：（1）計(jì)算一個(gè) Householder變換 H：由于其中在實(shí)際計(jì)算中，為避免出現(xiàn)兩個(gè)相近的數(shù)出現(xiàn)的情形，當(dāng)為便于儲(chǔ)存，將為防止溢出現(xiàn)象，用（2）QR 分解：利用 Householder變換逐步將- 1 - R0jvj 1QHAxAc1RxR 的結(jié)果最后一行為零，故使用

4、前代法時(shí)不計(jì)最后一x84。v,belta=house(x)Q,R=QRfenjie(A)- ，其中1:n jj0Hjb或最小二乘問(wèn)題的步驟為QT1 b，其中 Qc1H1H2，若，j1Hn，Rm xjQ(:,1:n)(1:n,:)。，依次計(jì)算，將，( j :m, j)vjQ對(duì)應(yīng)的(2:mH1H2R0jvj 1QHAxAc1RxR 的結(jié)果最后一行為零，故使用前代法時(shí)不計(jì)最后一x84。v,belta=house(x)Q,R=QRfenjie(A)- ，其中1:n jj0Hjb或最小二乘問(wèn)題的步驟為QT1 b，其中 Qc1H1H2，若，j1Hn，Rm xjQ(:,1:n)(1:n,:)。，依次計(jì)算，將

5、，( j :m, j)vjQ對(duì)應(yīng)的(2:mH1H2(Hj 1)H nj儲(chǔ)存到n（)( jm Q(1:時(shí)m k 1)(m, j)，H H2mk 1)j1儲(chǔ)存到Hn-1）d( j)，迭代結(jié)束后再次計(jì)A Q Q在實(shí)際計(jì)算中，從即對(duì)應(yīng)的I算 ，有0（3）求解線性方程組i 計(jì)算 的 QR分解；ii 計(jì)算iii 利用回代法求解上三角方程組（4）對(duì)第一章第一個(gè)線性方程組，由于行，而用運(yùn)行結(jié)果計(jì)算運(yùn)算 matlab 程序?yàn)? 計(jì)算 Householder 變換function v,belta=house(x) n=length(x); x=x/norm(x,inf); sigma=x(2:n)*x(2:n);

6、 v=zeros(n,1); v(2:n,1)=x(2:n);if sigma=0 belta=0;else alpha=sqrt(x(1)2+sigma);if x(1)=0 v(1)=x(1)-alpha;else v(1)=-sigma/(x(1)+alpha);end belta=2*v(1)2/(sigma+v(1)2); v=v/v(1,1);endend 2 計(jì)算 A的 QR 分解-完整版學(xué)習(xí)資料分享x=qiandaifa(L,b)- -WORD 格式-可編輯-專業(yè)資料- x=qiandaifa(L,b)- function Q,R=QRfenjie(A)m,n=size(A);

7、Q=eye(m);for j=1:nif jm 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,:);for j=1:nif jm H=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;endendend 3 解下三角形方程組的前代法function x=qiandaifa(L,b)n=length(b

8、);for j=1:n-1 b(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; end 4 求解第一章上機(jī)習(xí)題中的三個(gè)線性方程組 ex3_1 clear;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); %QR分解Q,R=QRfenjie(A); c=Q*b; x1=huidaifa(R(1:n-1,1:n-1),c(1:n-1); -完整版學(xué)習(xí)資料

9、分享- -WORD 格式-可編輯-專業(yè)資料- - 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(

10、 PGauss ); %第二題第一問(wèn)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=Ch

11、olesky(A); x2_3=Gauss(A,b,L,L);toc; %改進(jìn)的平方根法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);

12、plot(1:n,x2_4);title( 改進(jìn)的平方根法 ); %第二題第二問(wèn)A=hilb(40); b=sum(A); b=b; n=length(A); Q,R=QRfenjie(A); c=Q*b; -完整版學(xué)習(xí)資料分享求解房產(chǎn)估價(jià)的線性模型專業(yè)課 數(shù)值代數(shù)cha3_3_4.xls專業(yè)課 數(shù)值代數(shù)cha3_3_4.xls- ex3_3 ,A2:L29求解房產(chǎn)估價(jià)的線性模型專業(yè)課 數(shù)值代數(shù)cha3_3_4.xls專業(yè)課 數(shù)值代數(shù)cha3_3_4.xls- ex3_3 ,A2:L29 );,M2:M29 );x3=huidaifa(R,c); %不選主元 Gauss消去法L,U=Gauss

13、LA(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); %改進(jìn)的平方根法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_

14、2);title( PGauss ); subplot(1,5,4);plot(1:n,x3_3);title( 平方根法 ); subplot(1,5,5);plot(1:n,x3_4);title( 改進(jìn)的平方根法 ); 5 求解二次多項(xiàng)式 ex3_2 clear;clc;t=-1 -0.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 clear;cl

15、c;A=xlsread(y=xlsread(Q,R=QRfenjie(A);Q1=Q(:,1:12);c=Q1*y;x=huidaifa(R,c);x=x計(jì)算結(jié)果為（1）第一章上機(jī)習(xí)題中的三個(gè)線性方程組結(jié)果對(duì)比圖依次為-完整版學(xué)習(xí)資料分享x 1061211-21-650Gauss1086420-250- 811100PGauss1086420-2100Gauss0平 方根法 改進(jìn)的平方 根法1086420-20PGauss501086420-2501001000050501001000501000 x 1061211-21-650Gauss1086420-250- 811100PGauss10

16、86420-2100Gauss0平 方根法 改進(jìn)的平方 根法1086420-20PGauss501086420-250100100005050100100050100050100QR分 解1.81.641.411.2100.80.610.4-40.200QR分解1086420-20-完整版學(xué)習(xí)資料分享Gauss2001503100210000-2-100-3-1500不選主元的- PGauss4004300201-1-40-200-420平 根 法 改 進(jìn)的平方 根法5 60400-20-6040 x 10方-200702040-30002040-5020Gauss20015031002100

17、00-2-100-3-1500不選主元的- PGauss4004300201-1-40-200-420平 根 法 改 進(jìn)的平方 根法5 60400-20-6040 x 10方-200702040-30002040-502040-8002040QR分解25002000150020010005050000-50-500-100-1000-1500-2000以第二個(gè)線性方程組為例， 比較各方法的運(yùn)行速度。 依次為 QR分解，不選主元的 Gauss消去法，列主元 Gauss消去法，平方根法，改進(jìn)的平方根法。Elapsed time is 0.034588 seconds. Elapsed time is 0.006237 seconds. Elapsed time is 0.009689 seconds. Elapsed time is 0.030862 seconds. Elapsed time is 0.007622 seconds. （2）二次多項(xiàng)式的系數(shù)為x = 1.0000 1.0000 1.00