1)
a)
A = matrix([[3, 7], [1, 3]])
print "A = ",
A
A = \newcommand{\Bold}[1]{\mathbf{#1}}\left(\begin{array}{rr}
3 & 7 \\
1 & 3
\end{array}\right)
A = \newcommand{\Bold}[1]{\mathbf{#1}}\left(\begin{array}{rr}
3 & 7 \\
1 & 3
\end{array}\right)
|
x1 = matrix([[0.3536],[0.9354]])
y1 = A*x1
x2 = matrix([[0.9354], [0.3536]])
y2 = A*x2
x3 = matrix([[-0.3536], [0.9354]])
y3 = A*x3
x4 = matrix([[-0.9354], [0.3536]])
y4 = A*x4
myPlot = list_plot([(0, 0), (x1[0,0], x1[1,0])], plotjoined=true, legend_label ='i) x', color='blue', marker='o', markersize=4, linestyle='-.') + list_plot([(0, 0), (y1[0,0], y1[1,0])], plotjoined=true, legend_label ='i) y', color='blue', marker='o', markersize=4) + list_plot([(0, 0), (x2[0,0], x2[1,0])], plotjoined=true, legend_label ='ii) x', color='purple', marker='o', markersize=4, linestyle='-.') + list_plot([(0, 0), (y2[0,0], y2[1,0])], plotjoined=true, legend_label ='ii) y', color='purple', marker='o', markersize=4) + list_plot([(0, 0), (x3[0,0], x3[1,0])], plotjoined=true, legend_label ='iii) x', color='orange', marker='o', markersize=4, linestyle='-.') + list_plot([(0, 0), (y3[0,0], y3[1,0])], plotjoined=true, legend_label ='iii) y', color='orange', marker='o', markersize=4) + list_plot([(0, 0), (x4[0,0], x4[1,0])], plotjoined=true, legend_label ='iv) x', color='green', marker='o', markersize=4, linestyle='-.') + list_plot([(0, 0), (y4[0,0], y4[1,0])], plotjoined=true, legend_label ='iv) y', color='green', marker='o', markersize=4)
myPlot.set_legend_options(handletextpad=1.0, labelspacing=0.3, handlelength=2, loc=(0.14, 0.40))
myPlot
|
|
The vector x in ii) and iv) are eigenvectors because the x and y are along the same line.
b)
lam1 = ( ((A*x1).transpose()*x1) * ((x1.transpose()*x1).inverse()) )[0][0]
print "i)\n\nAx =\n", A*x1, "\n\n(lamda)x =\n", lam1*x1
lam2 = ( ((A*x2).transpose()*x2) * ((x2.transpose()*x2).inverse()) )[0][0]
print "\nii)\n\nAx =\n", A*x2, "\n\n(lamda)x =\n", lam2*x2
lam3 = ( ((A*x3).transpose()*x3) * ((x3.transpose()*x3).inverse()) )[0][0]
print "\niii)\n\nAx =\n", A*x3, "\n\n(lamda)x =\n", lam3*x3
lam4 = ( ((A*x4).transpose()*x4) * ((x4.transpose()*x4).inverse()) )[0][0]
print "\niv)\n\nAx =\n", A*x4, "\n\n(lamda)x =\n", lam4*x4
i) Ax = [7.60860000000000] [3.15980000000000] (lamda)x = [1.99644092014957] [5.28130892734136] ii) Ax = [5.28140000000000] [1.99620000000000] (lamda)x = [5.28130892734136] [1.99644092014957] iii) Ax = [5.48700000000000] [2.45260000000000] (lamda)x = [-0.125159079850431] [ 0.331091072658635] iv) Ax = [-0.331000000000000] [ 0.125400000000000] (lamda)x = [-0.331091072658635] [ 0.125159079850431] i) Ax = [7.60860000000000] [3.15980000000000] (lamda)x = [1.99644092014957] [5.28130892734136] ii) Ax = [5.28140000000000] [1.99620000000000] (lamda)x = [5.28130892734136] [1.99644092014957] iii) Ax = [5.48700000000000] [2.45260000000000] (lamda)x = [-0.125159079850431] [ 0.331091072658635] iv) Ax = [-0.331000000000000] [ 0.125400000000000] (lamda)x = [-0.331091072658635] [ 0.125159079850431] |
c)
x1 = matrix([[0.8660],[0.5000]])
y1 = A*x1
x2 = matrix([[0.5000], [0.8660]])
y2 = A*x2
x3 = matrix([[-0.5000], [0.8660]])
y3 = A*x3
x4 = matrix([[-0.8660], [0.5000]])
y4 = A*x4
myPlot = list_plot([(0, 0), (x1[0,0], x1[1,0])], plotjoined=true, legend_label ='i) x', color='blue', marker='o', markersize=4, linestyle='-.') + list_plot([(0, 0), (y1[0,0], y1[1,0])], plotjoined=true, legend_label ='i) y', color='blue', marker='o', markersize=4) + list_plot([(0, 0), (x2[0,0], x2[1,0])], plotjoined=true, legend_label ='ii) x', color='purple', marker='o', markersize=4, linestyle='-.') + list_plot([(0, 0), (y2[0,0], y2[1,0])], plotjoined=true, legend_label ='ii) y', color='purple', marker='o', markersize=4) + list_plot([(0, 0), (x3[0,0], x3[1,0])], plotjoined=true, legend_label ='iii) x', color='orange', marker='o', markersize=4, linestyle='-.') + list_plot([(0, 0), (y3[0,0], y3[1,0])], plotjoined=true, legend_label ='iii) y', color='orange', marker='o', markersize=4) + list_plot([(0, 0), (x4[0,0], x4[1,0])], plotjoined=true, legend_label ='iv) x', color='green', marker='o', markersize=4, linestyle='-.') + list_plot([(0, 0), (y4[0,0], y4[1,0])], plotjoined=true, legend_label ='iv) y', color='green', marker='o', markersize=4)
myPlot.set_legend_options(handletextpad=1.0, labelspacing=0.3, handlelength=2, loc=(0.14, 0.40))
myPlot
|
|
c)
None of the vectors x is an eigenvector.
lam1 = ( ((A*x1).transpose()*x1) * ((x1.transpose()*x1).inverse()) )[0][0]
print "i)\n\nAx =\n", A*x1, "\n\n(lamda)x =\n", lam1*x1
lam2 = ( ((A*x2).transpose()*x2) * ((x2.transpose()*x2).inverse()) )[0][0]
print "\nii)\n\nAx =\n", A*x2, "\n\n(lamda)x =\n", lam2*x2
lam3 = ( ((A*x3).transpose()*x3) * ((x3.transpose()*x3).inverse()) )[0][0]
print "\niii)\n\nAx =\n", A*x3, "\n\n(lamda)x =\n", lam3*x3
lam4 = ( ((A*x4).transpose()*x4) * ((x4.transpose()*x4).inverse()) )[0][0]
print "\niv)\n\nAx =\n", A*x4, "\n\n(lamda)x =\n", lam4*x4
i) Ax = [6.09800000000000] [2.36600000000000] (lamda)x = [5.59795599806391] [3.23207621135330] ii) Ax = [7.56200000000000] [3.09800000000000] (lamda)x = [3.23207621135330] [5.59795599806391] iii) Ax = [4.56200000000000] [2.09800000000000] (lamda)x = [ 0.232076211353300] [-0.401955998063915] iv) Ax = [0.902000000000000] [0.634000000000000] (lamda)x = [ 0.401955998063915] [-0.232076211353300] i) Ax = [6.09800000000000] [2.36600000000000] (lamda)x = [5.59795599806391] [3.23207621135330] ii) Ax = [7.56200000000000] [3.09800000000000] (lamda)x = [3.23207621135330] [5.59795599806391] iii) Ax = [4.56200000000000] [2.09800000000000] (lamda)x = [ 0.232076211353300] [-0.401955998063915] iv) Ax = [0.902000000000000] [0.634000000000000] (lamda)x = [ 0.401955998063915] [-0.232076211353300] |
2)
a)
# Generate random 6x6 matrices
A = random_matrix(ZZ, 6, x=0, y=10)
# Divide matrix A into 4 3x3 matrix blocks using matrix surgery
A1 = A.submatrix(0, 0, 3, 3)
A2 = A.submatrix(0, 3, 3, 3)
A3 = A.submatrix(3, 0, 3, 3)
A4 = A.submatrix(3, 3, 3, 3)
print "det(A) = ", A.det()
print "det(A1)det(A4) - det(A3)det(A2) = ", A1.det()*A4.det() - A3.det()*A2.det()
det(A) = 111153 det(A1)det(A4) - det(A3)det(A2) = 12408 det(A) = 111153 det(A1)det(A4) - det(A3)det(A2) = 12408 |
It appears that the formula det(A) = det(A1)det(A4) - det(A3)det(A2) is not valid in general.
b)
A3 = matrix(3, 3, 0)
blockTop = (A1).augment(A2)
blockBottom = (A3).augment(A4)
A = (blockTop).stack(blockBottom)
print "det(A) = ", A.det()
print "det(A1)det(A4) = ", A1.det()*A4.det()
det(A) = 2688 det(A1)det(A4) = 2688 det(A) = 2688 det(A1)det(A4) = 2688 |
c)
# Generate random 12x12 matrices
A = random_matrix(ZZ, 12, x=0, y=10)
# Divide matrix A into 4 3x3 matrix blocks using matrix surgery
A11 = A.submatrix(0, 0, 4, 4)
A12 = A.submatrix(0, 4, 4, 4)
A13 = A.submatrix(0, 8, 4, 4)
A21 = matrix(4, 4, 0)
A22 = A.submatrix(4, 4, 4, 4)
A23 = A.submatrix(4, 8, 4, 4)
A31 = matrix(4, 4, 0)
A32 = matrix(4, 4, 0)
A33 = A.submatrix(8, 8, 4, 4)
blockTop = A11.augment(A12).augment(A13)
blockMiddle = A21.augment(A22).augment(A23)
blockBottom = A31.augment(A32).augment(A33)
A = blockTop.stack(blockMiddle).stack(blockBottom)
print "det(A) = ", A.det()
print "det(A11)det(A22)det(33) = ", A11.det()*A22.det()*A33.det()
det(A) = -247665600 det(A11)det(A22)det(33) = -247665600 det(A) = -247665600 det(A11)det(A22)det(33) = -247665600 |
3)
a)
A = matrix(RR, [[3, -1,-1], [-12, 0, 5], [4, -2, -1]])
x = matrix(RR, [[1], [1], [1]])
|
|
print "x0 =\n", x, "\n"
for i in range(10):
x = A*x
print "x\r", i+1, " =\n", x, "\n"
x0 = [1.00000000000000] [1.00000000000000] [1.00000000000000] x1 = [ 1.00000000000000] [-7.00000000000000] [ 1.00000000000000] x2 = [ 9.00000000000000] [-7.00000000000000] [ 17.0000000000000] x3 = [ 17.0000000000000] [-23.0000000000000] [ 33.0000000000000] x4 = [ 41.0000000000000] [-39.0000000000000] [ 81.0000000000000] x5 = [ 81.0000000000000] [-87.0000000000000] [ 161.000000000000] x6 = [ 169.000000000000] [-167.000000000000] [ 337.000000000000] x7 = [ 337.000000000000] [-343.000000000000] [ 673.000000000000] x8 = [ 681.000000000000] [-679.000000000000] [ 1361.00000000000] x9 = [ 1361.00000000000] [-1367.00000000000] [ 2721.00000000000] x10 = [ 2729.00000000000] [-2727.00000000000] [ 5457.00000000000] x0 = [1.00000000000000] [1.00000000000000] [1.00000000000000] x1 = [ 1.00000000000000] [-7.00000000000000] [ 1.00000000000000] x2 = [ 9.00000000000000] [-7.00000000000000] [ 17.0000000000000] x3 = [ 17.0000000000000] [-23.0000000000000] [ 33.0000000000000] x4 = [ 41.0000000000000] [-39.0000000000000] [ 81.0000000000000] x5 = [ 81.0000000000000] [-87.0000000000000] [ 161.000000000000] x6 = [ 169.000000000000] [-167.000000000000] [ 337.000000000000] x7 = [ 337.000000000000] [-343.000000000000] [ 673.000000000000] x8 = [ 681.000000000000] [-679.000000000000] [ 1361.00000000000] x9 = [ 1361.00000000000] [-1367.00000000000] [ 2721.00000000000] x10 = [ 2729.00000000000] [-2727.00000000000] [ 5457.00000000000] |
g = x[0,0]
x[2,0] = x[2,0] / g
h = x[2,0] / g
x[1,0] = x[1,0] / g
h = x[2,0] / g
x[0,0] = x[0,0] / g
print "x10 =",
x
x10 = \newcommand{\Bold}[1]{\mathbf{#1}}\left(\begin{array}{r}
1.00000000000000 \\
-0.999267130817149 \\
1.99963356540857
\end{array}\right)
x10 = \newcommand{\Bold}[1]{\mathbf{#1}}\left(\begin{array}{r}
1.00000000000000 \\
-0.999267130817149 \\
1.99963356540857
\end{array}\right)
|
x = A*x
print "x11 = \n", x
g = x[0,0]
x[2,0] = x[2,0] / g
h = x[2,0] / g
x[1,0] = x[1,0] / g
h = x[2,0] / g
x[0,0] = x[0,0] / g
print "=\n", x, "\n"
x = A*x
print "x12 = \n", x
g = x[0,0]
x[2,0] = x[2,0] / g
h = x[2,0] / g
x[1,0] = x[1,0] / g
h = x[2,0] / g
x[0,0] = x[0,0] / g
print "=\n", x, "\n"
x = A*x
print "x13 = \n", x
g = x[0,0]
x[2,0] = x[2,0] / g
h = x[2,0] / g
x[1,0] = x[1,0] / g
h = x[2,0] / g
x[0,0] = x[0,0] / g
print "=\n", x
x11 = [ 1.99963356540857] [-2.00183217295713] [ 3.99890069622572] = [ 1.00000000000000] [-1.00109950522265] [ 1.99981674912956] x12 = [ 2.00128275609309] [-2.00091625435221] [ 4.00238226131574] = [ 1.00000000000000] [-0.999816866587309] [ 1.99990843329365] x13 = [ 1.99990843329365] [-2.00045783353173] [ 3.99972529988096] = [ 1.00000000000000] [-1.00027471269630] [ 1.99995421455062] x11 = [ 1.99963356540857] [-2.00183217295713] [ 3.99890069622572] = [ 1.00000000000000] [-1.00109950522265] [ 1.99981674912956] x12 = [ 2.00128275609309] [-2.00091625435221] [ 4.00238226131574] = [ 1.00000000000000] [-0.999816866587309] [ 1.99990843329365] x13 = [ 1.99990843329365] [-2.00045783353173] [ 3.99972529988096] = [ 1.00000000000000] [-1.00027471269630] [ 1.99995421455062] |
b)
Since xk+1 ≈ λ1xk and xk+1 ≈ 2xk, then λ1 ≈ 2. Thus the dominant eigenvalue is ≈ 2.
c)
x = matrix([[1], [2], [3]])
print "x0 =\n", x, "\n"
for i in range(10):
x = A*x
x = x / x.norm()
print "x\r", i+1, "=\n", x, "\n"
x0 = [1] [2] [3] x1 = [-0.426401432711221] [ 0.639602149066831] [-0.639602149066831] x2 = [-0.388922234131299] [ 0.583383351196948] [-0.713024095907381] x3 = [-0.412294457624059] [ 0.438062861225562] [-0.798820511646614] x4 = [-0.405971957998941] [ 0.441793013116494] [-0.800003564292032] x5 = [-0.409183387378042] [ 0.414866489980515] [-0.812683672153612] x6 = [-0.407767145362268] [ 0.416145922321766] [-0.812741365071370] x7 = [-0.408477576428033] [ 0.409857568240290] [-0.815575161043810] x8 = [-0.408132827478516] [ 0.410194104384974] [-0.815578562654880] x9 = [-0.408305336227165] [ 0.408647874260912] [-0.816268134420583] x10 = [-0.408219716780686] [ 0.408732985494158] [-0.816268343990216] x11 = [-0.408262534746855] [ 0.408348016852120] [-0.816439587388445] x12 = [-0.408241165154219] [ 0.408369354760192] [-0.816439600439781] x13 = [-0.408251850463421] [ 0.408273211485317] [-0.816482339904947] x0 = [1] [2] [3] x1 = [-0.426401432711221] [ 0.639602149066831] [-0.639602149066831] x2 = [-0.388922234131299] [ 0.583383351196948] [-0.713024095907381] x3 = [-0.412294457624059] [ 0.438062861225562] [-0.798820511646614] x4 = [-0.405971957998941] [ 0.441793013116494] [-0.800003564292032] x5 = [-0.409183387378042] [ 0.414866489980515] [-0.812683672153612] x6 = [-0.407767145362268] [ 0.416145922321766] [-0.812741365071370] x7 = [-0.408477576428033] [ 0.409857568240290] [-0.815575161043810] x8 = [-0.408132827478516] [ 0.410194104384974] [-0.815578562654880] x9 = [-0.408305336227165] [ 0.408647874260912] [-0.816268134420583] x10 = [-0.408219716780686] [ 0.408732985494158] [-0.816268343990216] x11 = [-0.408262534746855] [ 0.408348016852120] [-0.816439587388445] x12 = [-0.408241165154219] [ 0.408369354760192] [-0.816439600439781] x13 = [-0.408251850463421] [ 0.408273211485317] [-0.816482339904947] |
A = matrix([[2.8, -1.6], [-1.6, 5.2]])
x = matrix([[1/sqrt(2)], [1/sqrt(2)]])
p = list_plot([(0, 0), (x[0,0], x[1,0])], plotjoined=true, marker='o', markersize=4)
for i in range(10):
x = A*x
x = x / x.norm()
p = p + list_plot([(0, 0), (x[0,0], x[1,0])], plotjoined=true, marker='o', markersize=4)
p
|
|
for i in range(5):
Ax = A*x
Beta = x.transpose()*(Ax)*(x.transpose()*x).inverse()
print "Beta\r", i, "=", Beta
x = Ax
Beta0 = [5.99999999999983] Beta1 = [5.99999999999998] Beta2 = [6.00000000000000] Beta3 = [6.00000000000000] Beta4 = [6.00000000000000] Beta0 = [5.99999999999983] Beta1 = [5.99999999999998] Beta2 = [6.00000000000000] Beta3 = [6.00000000000000] Beta4 = [6.00000000000000] |