# 213 - A15 - Gram-Schmidt

## 836 days ago by Professor213

# Do this first cell in worksheet in Octave mode # If copying to the sagecell, then \n needs to be changed to \\n
A = floor(rand(4,4)*10)-5; det(A) # A =[0 2 2 -1; -2 -1 -5 1; 2 -1 2 4; 2 3 -5 1]; # Comment this out if you want to get another random matrix. Just check to see that det(a) is nonzero v1 = A(:,1); v2 = A(:,2); v3 = A(:,3); v4 = A(:,4); disp("Here are the original vectors") V = [v1 v2 v3 v4]; printf ("% 4.0f % 4.0f % 4.0f % 4.0f \n", V'); disp("Now apply Gram-Schmidt") u1 = v1; u2 = v2 - (v2'*u1)/(u1'*u1)*u1; u3 = v3 - (v3'*u1)/(u1'*u1)*u1 - (v3'*u2)/(u2'*u2)*u2; u4 = v4 - (v4'*u1)/(u1'*u1)*u1 - (v4'*u2)/(u2'*u2)*u2 - (v4'*u3)/(u3'*u3)*u3; Q = [u1 u2 u3 u4]; printf ("% 5.3f % 5.3f % 5.3f % 5.3f \n", Q'); # if desired, make the vectors all unit vectors too for orthonormal disp(" and after making vectors all unit vectors...") Q = [u1/norm(u1) u2/norm(u2) u3/norm(u3) u4/norm(u4)]; printf ("% 5.3f % 5.3f % 5.3f % 5.3f \n", Q'); disp(' ') disp(' Check for orthogonality') printf ("u1 with u2: % 9.8f \n",u1'*u2) printf ("u1 with u3: % 9.8f \n",u1'*u3) printf ("u1 with u4: % 9.8f \n",u1'*u4) printf ("u2 with u3: % 9.8f \n",u2'*u3) printf ("u2 with u4: % 9.8f \n",u2'*u4) printf ("u3 with u4: % 9.8f \n",u3'*u4)
 ans = 81 Here are the original vectors 4 1 -2 1 4 1 -5 -1 -4 3 0 2 -3 -1 -4 -3 Now apply Gram-Schmidt 4.000 1.070 -0.584 0.386 4.000 1.070 -3.584 -0.269 -4.000 2.930 -0.321 0.009 -3.000 -1.053 -5.130 0.144 and after making vectors all unit vectors... 0.530 0.309 -0.093 0.784 0.530 0.309 -0.569 -0.547 -0.530 0.846 -0.051 0.018 -0.397 -0.304 -0.815 0.292 Check for orthogonality u1 with u2: 0.00000000 u1 with u3: -0.00000000 u1 with u4: 0.00000000 u2 with u3: -0.00000000 u2 with u4: 0.00000000 u3 with u4: 0.00000000  ans = 81 Here are the original vectors 4 1 -2 1 4 1 -5 -1 -4 3 0 2 -3 -1 -4 -3 Now apply Gram-Schmidt 4.000 1.070 -0.584 0.386 4.000 1.070 -3.584 -0.269 -4.000 2.930 -0.321 0.009 -3.000 -1.053 -5.130 0.144 and after making vectors all unit vectors... 0.530 0.309 -0.093 0.784 0.530 0.309 -0.569 -0.547 -0.530 0.846 -0.051 0.018 -0.397 -0.304 -0.815 0.292 Check for orthogonality u1 with u2: 0.00000000 u1 with u3: -0.00000000 u1 with u4: 0.00000000 u2 with u3: -0.00000000 u2 with u4: 0.00000000 u3 with u4: 0.00000000 
# Notice, this also give us the A = QR factorization format short R = Q'*A Q*R - A
 R = 7.54983 -0.13245 -2.11925 0.13245 0.00000 3.46157 -0.94775 2.60505 0.00000 0.00000 6.29369 2.82023 0.00000 0.00000 -0.00000 0.49246 ans = 1.7764e-15 1.1102e-15 -8.8818e-16 6.6613e-16 -1.3323e-15 -3.3307e-16 8.8818e-16 -2.2204e-16 0.0000e+00 4.4409e-16 3.9311e-17 4.4409e-16 8.8818e-16 4.4409e-16 4.4409e-16 8.8818e-16 R = 7.54983 -0.13245 -2.11925 0.13245 0.00000 3.46157 -0.94775 2.60505 0.00000 0.00000 6.29369 2.82023 0.00000 0.00000 -0.00000 0.49246 ans = 1.7764e-15 1.1102e-15 -8.8818e-16 6.6613e-16 -1.3323e-15 -3.3307e-16 8.8818e-16 -2.2204e-16 0.0000e+00 4.4409e-16 3.9311e-17 4.4409e-16 8.8818e-16 4.4409e-16 4.4409e-16 8.8818e-16
# Let's do it in 3D in Sage mode # a = random_matrix(QQ,3,3)] @interact def _(show_orthogonals = checkbox(default=false)): a = matrix(QQ,[[1,-1,0],[-1,1/2,1],[0,-2,0]]) v1 =vector(QQ,[1,-1,1]) v2 =vector(QQ,[-1,1/2,-2]) v3 =vector(QQ,[0,1,0]) G = plot(v1,color='red',opacity=0.2,width=20)+plot(v2,color='blue',opacity=0.2)+plot(v3,color='orange',opacity=0.2) if show_orthogonals: u1 = v1 u2 = v2 - v2.inner_product(u1)/u1.inner_product(u1)*u1 u3 = v3 - v3.inner_product(u1)/u1.inner_product(u1)*u1 - v3.inner_product(u2)/u2.inner_product(u2)*u2 print(u1) print(u2) print(u3) G += plot(u1,color='red')+plot(u2,color='blue')+plot(u3,color='orange') G += plot(u1/norm(u1),color='black')+plot(u2/norm(u2),color='black')+plot(u3/norm(u3),color='black') G.show(spin=true)

 show_orthogonals

## Click to the left again to hide and once more to show the dynamic interactive window

u1.inner_product(u2)
 0 0
# Let's do it in 3D in Sage mode v1 =vector(QQbar,[1+I,-1,0]) v2 =vector(QQbar,[-1,1/2,-2-I]) v3 =vector(QQbar,[I,1,0]) u1 = v1 u2 = v2 - v2.inner_product(u1)/u1.inner_product(u1)*u1 u3 = v3 - v3.inner_product(u1)/u1.inner_product(u1)*u1 - v3.inner_product(u2)/u2.inner_product(u2)*u2 show(u1) show(u2) show(u3)
u1.inner_product(u2)
 0 0
# Here is another octave-only version of this in one cell A = [2 -3 5;1 4 0;-1 3 2;1 2 3] # Let's make the columns of A orthogonal to each other A1 = A(:,1); A2 = A(:,2); A3 = A(:,3); disp("First we check for existing orthogonality") A1'*A2 A1'*A3 A2'*A3 disp("Pick one to start") v1 = A1 disp("Then, remove the part of the second column that is parallel to v1") v2 = A2 - (A2'*v1)/(v1'*v1)*v1 disp("Now, remove the parts of the third column that are parallel to v1 and v2") v3 = A3 - (A3'*v1)/(v1'*v1)*v1 - (A3'*v2)/(v2'*v2)*v2 disp("And check for orthogonality of all three vectors") v1'*v2 v1'*v3 v2'*v3 disp("Finally, let's make all the vectors 'unit vectors'") v1 = v1/norm(v1) v2 = v2/norm(v2) v3 = v3/norm(v3) disp("And create a matrix of these three normalized and orthogonalized vectors that came from A") Q = [v1 v2 v3] disp("Notice that after this process, the 'A = QR' factorization can be obtained by Q'*A = R") format short R = Q'*A disp("And we can check our factorization") A - Q*R disp("Let's solve a system using Least Squares") b = [4;-2;-3;1] mat = [A b] rref(mat) ans =rref([A'*A A'*b]) x = ans(:,4) r = b- A*x norm(r)
 WARNING: Output truncated! full_output.txt A = 2 -3 5 1 4 0 -1 3 2 1 2 3 First we check for existing orthogonality ans = -3 ans = 11 ans = -3 Pick one to start v1 = 2 1 -1 1 Then, remove the part of the second column that is parallel to v1 v2 = -2.1429 4.4286 2.5714 2.4286 Now, remove the parts of the third column that are parallel to v1 and v2 v3 = 1.9572 -1.7782 3.4514 1.3152 And check for orthogonality of all three vectors ans = 0 ans = 6.6613e-16 ans = 4.4409e-16 Finally, let's make all the vectors 'unit vectors' v1 = 0.75593 0.37796 -0.37796 0.37796 v2 = -0.35365 0.73088 0.42438 0.40081 v3 = 0.43086 -0.39146 0.75979 ... -0.37796 0.42438 0.75979 0.37796 0.40081 0.28953 Notice that after this process, the 'A = QR' factorization can be obtained by Q'*A = R R = 2.64575 -1.13389 4.15761 0.00000 6.05923 0.28292 0.00000 -0.00000 4.54249 And we can check our factorization ans = 0.0000e+00 0.0000e+00 0.0000e+00 0.0000e+00 -1.7764e-15 4.4409e-16 -1.1102e-16 -4.4409e-16 -8.8818e-16 0.0000e+00 -4.4409e-16 0.0000e+00 Let's solve a system using Least Squares b = 4 -2 -3 1 mat = 2 -3 5 4 1 4 0 -2 -1 3 2 -3 1 2 3 1 ans = 1 0 0 0 0 1 0 0 0 0 1 0 0 0 0 1 ans = 1.00000 0.00000 0.00000 0.98246 0.00000 1.00000 0.00000 -0.62399 0.00000 0.00000 1.00000 0.11371 x = 0.98246 -0.62399 0.11371 r = -0.40543 -0.48652 -0.37300 0.92438 ans = 1.1810 WARNING: Output truncated! full_output.txt A = 2 -3 5 1 4 0 -1 3 2 1 2 3 First we check for existing orthogonality ans = -3 ans = 11 ans = -3 Pick one to start v1 = 2 1 -1 1 Then, remove the part of the second column that is parallel to v1 v2 = -2.1429 4.4286 2.5714 2.4286 Now, remove the parts of the third column that are parallel to v1 and v2 v3 = 1.9572 -1.7782 3.4514 1.3152 And check for orthogonality of all three vectors ans = 0 ans = 6.6613e-16 ans = 4.4409e-16 Finally, let's make all the vectors 'unit vectors' v1 = 0.75593 0.37796 -0.37796 0.37796 v2 = -0.35365 0.73088 0.42438 0.40081 v3 = 0.43086 -0.39146 0.75979 ... -0.37796 0.42438 0.75979 0.37796 0.40081 0.28953 Notice that after this process, the 'A = QR' factorization can be obtained by Q'*A = R R = 2.64575 -1.13389 4.15761 0.00000 6.05923 0.28292 0.00000 -0.00000 4.54249 And we can check our factorization ans = 0.0000e+00 0.0000e+00 0.0000e+00 0.0000e+00 -1.7764e-15 4.4409e-16 -1.1102e-16 -4.4409e-16 -8.8818e-16 0.0000e+00 -4.4409e-16 0.0000e+00 Let's solve a system using Least Squares b = 4 -2 -3 1 mat = 2 -3 5 4 1 4 0 -2 -1 3 2 -3 1 2 3 1 ans = 1 0 0 0 0 1 0 0 0 0 1 0 0 0 0 1 ans = 1.00000 0.00000 0.00000 0.98246 0.00000 1.00000 0.00000 -0.62399 0.00000 0.00000 1.00000 0.11371 x = 0.98246 -0.62399 0.11371 r = -0.40543 -0.48652 -0.37300 0.92438 ans = 1.1810