Summary 2071 lab 6 M. M. Sussman $Id: summary.txt,v 1.12 2017/03/18 20:30:09 mike Exp $ {@(#) Thu Jan 21 16:50:15 2016 } EXERCISE 1 1. exer1.m: Ainv=???, xSolution=???, name,date 2. using "matlab -singleCompThread" to get 1 thread on my laptop Time to compute inverse matrices n time ratio 161 0.0037 2.1920 321 0.0080 6.6466 641 0.0535 6.7805 1281 0.3627 7.4577 2561 2.7046 7.4786 5121 20.2263 7.6735 10241 155.2077 3. yes, consistent with O(n^3) EXERCISE 2 1. Time to compute solutions n time ratio 161 0.0021 1.8561 321 0.0038 6.1150 641 0.0235 6.0935 1281 0.1432 6.9328 2561 0.9930 7.0661 5121 7.0167 7.5437 10241 52.9325 2. yes 3. Comparison of times n (time for inverse)/(time for solution) 161 1.7755 =1/0.5632 321 2.1332 =1/0.4688 641 2.2576 =1/0.4430 1281 2.5140 =1/0.3978 2561 2.7141 =1/0.3684 5121 2.8708 =1/0.3483 10241 2.9334 =1/0.3409 4. yes, roughly 3 times faster EXERCISE 3 1. norm((A\b)-(AA\b))/norm(A\b)=roundoff: yes, same 2. Time to compute solutions Size time ratio 10240 0.0011 102400 0.0027 1024000 0.0307 10240000 0.3017 3. p=1 4. sparse time 10240=.0011, full time= 52 EXERCISE 4 1. dif2 matrix size error determinant cond. no. 10 1.6522e-15 11.0000 48.3742 40 1.2154e-14 41.0000 680.6171 160 1.3467e-13 161.0000 1.0505e+04 2. Hilbert matrix size error determinant cond. no. 10 9.4097e-05 2.1644e-53 1.6025e+13 15 1.6972 -1.8701e-120 4.4333e+17 20 21.4973 1.9979e-194 2.0383e+18 3. Frank matrix size error determinant cond. no. 10 3.3283e-11 1 2.8543e+07 15 1.3641e-04 1 1.3710e+13 20 0.7172 18.1766 2.4122e+18 4. Pascal matrix size error determinant cond. no. 10 1.2946e-09 1 4.1552e+09 15 0.0034 1.0011 2.8389e+15 20 550.3993 -27.5964 1.3747e+20 EXERCISE 5 1. L = 1.0000 0 0 0 0 0.5000 1.0000 0 0 0 0.3333 0 1.0000 0 0 0.2500 0 0 1.0000 0 0.2000 0 0 0 1.0000 U = 1.0000 0.5000 0.3333 0.2500 0.2000 0 0.0833 0.0833 0.0750 0.0667 0 0.0833 0.0889 0.0833 0.0762 0 0.0750 0.0833 0.0804 0.0750 0 0.0667 0.0762 0.0750 0.0711 2. norm(L*U-A,'fro')=0 3. gauss_lu.m, wrap a loop on Jcol around given code 4. L = 1.0000 0 0 0 0 0.5000 1.0000 0 0 0 0.3333 1.0000 1.0000 0 0 0.2500 0.9000 1.5000 1.0000 0 0.2000 0.8000 1.7143 2.0000 1.0000 U = 1.0000 0.5000 0.3333 0.2500 0.2000 0 0.0833 0.0833 0.0750 0.0667 0 -0.0000 0.0056 0.0083 0.0095 0 0 0 0.0004 0.0007 0 0 0 0 0.00002 5. U(5,5)=2.2676e-05 6. visually, L is lower-triangular, U is upper-triangular, 7. norm(L*U-A,'fro')=roundoff 8. norm(LR*UR-R,'fro')=roundoff norm(tril(UR,-1),'fro')=roundoff norm(triu(LR,1),'fro')=roundoff also: norm(LR-tril(LR))=roundoff EXERCISE 6. 1. Fails 2. Irow=4, Jcol=3: 0 on diag=>0/0 on diag, n/0 off 3. det(A1)=-8 /= 0; cond(A1)=15.4, neither singular nor ill-conditioned. EXERCISE 7. 1. A2 = 16 2 3 13 9 7 6 12 5 11 10 8 4 14 15 1 2. A3 = 16 3 2 13 5 10 11 8 9 6 7 12 4 15 14 1 3. P3 interchange first, fourth, 0 0 0 1 0 1 0 0 0 0 1 0 1 0 0 0 4. P=P1*P2 = 1 0 0 0 0 0 1 0 0 0 0 1 0 1 0 0 5. A4=P*A = 16 2 3 13 9 7 6 12 4 14 15 1 5 11 10 8 EXERCISE 8 1,2,3,4 gauss_plu.m: P returned, pivot code, debug code 5. A=pascal(5) A = 1 1 1 1 1 1 2 3 4 5 1 3 6 10 15 1 4 10 20 35 1 5 15 35 70 P= 1 0 0 0 0 0 0 0 0 1 0 0 1 0 0 0 0 0 1 0 0 1 0 0 0 L= 1. 0. 0. 0. 0. 1. 1. 0. 0. 0. 1. 0.5 1. 0. 0. 1. 0.75 0.75 1. 0. 1. 0.25 0.75 -1. 1. U = 1. 1. 1. 1. 1. 0. 4. 14. 34. 69. 0. 0. -2. -8. -20.5 0. 0. 0. -0.5 -2.37500 0. 0. 0. 0. -0.25 6. A1 from before P1 = 1 0 0 0 0 0 0 0 1 0 0 0 0 0 1 0 1 0 0 0 0 0 1 0 0 L1 = 1 0 0 0 0 -0.5000 1 0 0 0 0 -0.6667 1 0 0 -0.5000 -0.3333 -1 1 0 0 0 0 -1 1 U1 = -2 1 0 0 0 0 -1.5000 1 0 0 0 0 -1.3333 1 0 0 0 0 -1 0 0 0 0 0 -2 7. remove debug code. EXERCISE 9. 1. P*L*U= [ 2 6 12 1 3 8 4 4 8] 2. Step 0 b = [ 28; 18; 16] 1 z = [ 16; 28; 18] 2 y = [ 16; 20; 4] 3 x = [ -1; 1; 2] EXERCISE 10 1. p_solve.m, z=P'*b 2. check: z = [ 16; 28; 18] 3. l_solve.m: rhs = rhs - L(Irow,Jcol)*y(Jcol); y(Irow) = rhs/L(Irow,Irow); or just y(Irow)=rhs; 4. check: y = [ 16; 20; 4] 5. u_solve.m: rhs = rhs - U(Irow,Jcol)*x(Jcol); x(Irow) = rhs/U(Irow,Irow); 6. x = [ -1; 1; 2] 7. x=plu_solve(P,L,U,b)=[-1;1;2] relErr = norm(P*L*U*x-b)/norm(b) = 0 9. Random 100x100 matrix. Rel Difference=roundoff 7.0961e-14 EXERCISE 11 1. check "debugging check" is gone 2. bels.m, Factoring each time yields about sec. (104.1853 on oz) 3. bels1.m, check soln AND TIME against bels 1 t.s. It is an error to leave the factorization inside the loop! 4. Factoring once yields about ____ sec. (3.8601 on oz) 5. Two solutions are identical EXERCISE 12 (Extra) 1. P = 1 0 0 0 0 0 interchanges rows 4 and 6 0 1 0 0 0 0 0 0 1 0 0 0 0 0 0 0 0 1 0 0 0 0 1 0 0 0 0 1 0 0 L = 1. 0. 0. 0. 0. 0. 0.16667 1. 0. 0. 0. 0. 0.16667 0.14286 1. 0. 0. 0. 0. 0.17143 0.15000 1. 0. 0. 0.16667 -0.02857 0.15000 0.15789 1. 0. 0.16667 0.14286 -0.05000 0.15789 0.13636 1. U = 6. 1. 1. 0. 1. 1. 0. 5.83333 0.83333 1. -0.16667 0.83333 0. 0. 5.71429 0.85714 0.85714 -0.28571 0. 0. -0. 5.7 0.9 0.9 0. 0. 0. 0. 5.55789 0.75789 0. 0. 0. 0. 0. 5.45455 2. P = 1 0 0 0 0 0 0 1 0 0 0 0 0 0 1 0 0 0 0 0 0 1 0 0 0 0 0 0 1 0 0 0 0 0 0 1 L = 1. 0. 0. 0. 0. 0. 0.16667 1. 0. 0. 0. 0. 0.16667 0.14286 1. 0. 0. 0. 0.16667 0.14286 -0.05000 1. 0. 0. 0.16667 -0.02857 0.15000 0. 1. 0. 0. 0.17143 0.15000 0. 0. 1. U = 6. 1. 1. 0. 1. 1. 0. 5.83333 0.83333 1. -0.16667 0.83333 0. 0. 5.71429 0.85714 0.85714 -0.28571 0. 0. 0. 0.9 0.9 5.7 0. 0. 0. 0.9 5.7 0.9 0. 0. -0. 5.7 0.9 0.9 3. Pi = 1 0 0 0 0 0 0 1 0 0 0 0 0 0 1 0 0 0 0 0 0 0 0 1 0 0 0 0 1 0 0 0 0 1 0 0 Lambda = 0. 0. 0. 0. 0. 0. 0.16667 0. 0. 0. 0. 0. 0.16667 0.14286 0. 0. 0. 0. 0.16667 0.14286 -0.05000 0. 0. 0. 0.16667 -0.02857 0.15000 0. 0. 0. 0. 0.17143 0.15000 0. 0. 0. norm(Lambda-Lambda*Pi') Lambda*Pi'=Lambda because on step k, Lambda has non-zero values only in columns 1-(k-1) and Pi' interchanges column k with a later column. 4. new_plu.m + check tests: test suite is graded.