python线代可视化

_光风霁月_

已于 2024-09-19 18:36:41 修改

阅读量754

点赞数 5

分类专栏：线性代数文章标签： python 开发语言线性代数

于 2024-09-19 18:35:13 首次发布

本文链接：https://blog.youkuaiyun.com/qq_45651049/article/details/142366753

版权

线性代数专栏收录该内容

1 篇文章

订阅专栏

前排提醒：本博客仅用于存档，供自己和平台上有需要的人学习使用。不可贩卖等作为他用。
如有问题请自行解决，不负责答疑。如有勘误，请联系本文作者。

Assignment 01

在这里插入图片描述

import numpy as np
import matplotlib.pyplot as plt

u = np.array([4,7])
v = np.array([8,4])

w = u + v
z = v - u

plt.quiver(0,0,*u,angles='xy',scale_units='xy',scale=1,color='r',label='u')
plt.quiver(0,0,*v,angles='xy',scale_units='xy',scale=1,color='b',label='v')
plt.quiver(u[0],u[1],*v,angles='xy',scale_units='xy',scale=1,color='b')
plt.quiver(v[0],v[1],*u,angles='xy',scale_units='xy',scale=1,color='r')
plt.quiver(0,0,*w,angles='xy',scale_units='xy',scale=1,color='g',label='u+v')
plt.text(u[0],u[1],'u',color='black')
plt.text(v[0],v[1],'v',color='black')
plt.text(w[0],w[1],'u+v',color='black')
plt.xlim(0,15)
plt.ylim(0,15)
plt.xlabel('x-axis')
plt.ylabel('y-axis')
plt.title('Vector Addition')
plt.legend()
plt.grid()
plt.show()

print("u + v = ",w)

Assignment 02

在这里插入图片描述

import numpy as np
import matplotlib.pyplot as plt

A = np.array([[1,2,3],
              [4,5,6],
              [7,8,9]])
B = np.array([[9,8,7],
              [6,5,4],
              [3,2,1]])
C = np.dot(A,B)

fig,axs = plt.subplots(1,3,figsize=(12,4))

im1 = axs[0].imshow(A,cmap='OrRd',interpolation='nearest',aspect='auto')
axs[0].set_title('Matrix A')
axs[0].axis('off')
cbar1 = plt.colorbar(im1,ax=axs[0],pad=0.02)

im2 = axs[1].imshow(B,cmap='summer',interpolation='nearest',aspect='auto')
axs[1].set_title('Matrix B')
axs[1].axis('off')
cbar2 = plt.colorbar(im2,ax=axs[1],pad=0.02)

im3 = axs[2].imshow(C,cmap='PuBu',interpolation='nearest',aspect='auto')
axs[2].set_title('Matrix C(A × B)')
axs[2].axis('off')
cbar2 = plt.colorbar(im2,ax=axs[2],pad=0.02)

for i in range(3):
    for j in range(3):
        axs[0].text(j,i,str(A[i,j]),ha='center',va='center',color='black')
        axs[1].text(j, i, str(B[i, j]), ha='center', va='center', color='black')
        axs[2].text(j, i, str(C[i, j]), ha='center', va='center', color='black')

plt.tight_layout()
plt.show()

Assignment 03

在这里插入图片描述

import numpy as np

a = np.array([[2.0,-1.0,0.0],[-1.0,3.0,-1.0],[0.0,-1.0,2.0]])
b = np.array([1.0,8.0,-5.0])
x0 = np.zeros(len(b))
n = 1000

for i in range(n):
    x = np.zeros_like(x0)
    for j in range(len(b)):
        s = np.dot(a[j,:],x0) - a[j,j] * x0[j]
        x[j] = (b[j] - s) / a[j,j]
    x0 = x

print(x.astype((int)))

Assignment 04

在这里插入图片描述

import numpy as np
import matplotlib.pyplot as plt

vectors_set_1 = np.array([[2,4,5],[-4,-8,-6]])
vectors_set_2 = np.array([[3,5,7],[5,8,4]])

dot_product_set_1 = np.dot(vectors_set_1[0],vectors_set_1[1])
dot_product_set_2 = np.dot(vectors_set_2[0],vectors_set_2[1])

if (dot_product_set_1 == 0):
    isDependent_1 = 1
else:
    isDependent_1 = 0

if(dot_product_set_2 == 0):
    isDependent_2 = 1
else:
    isDependent_2 = 0

fig = plt.figure()
ax = fig.add_subplot(111,projection = "3d")

origin = np.zeros(3)

for vector in vectors_set_1:
    ax.quiver(*origin,*vector,color="blue",label="set1")
for vector in vectors_set_2:
    ax.quiver(*origin,*vector,color="red",label="set2")

ax.set_xlim([-10,10])
ax.set_ylim([-10,10])
ax.set_zlim([-10,10])

if (isDependent_1):
    print("Vectors set1 is dependent.")
else:
    print("Vectors set1 is independent.")

if (isDependent_2):
    print("Vectors set2 is dependent.")
else:
    print("Vectors set2 is independent.")

plt.legend()
plt.show()

Assignment 05

在这里插入图片描述
在文件目录下新建一个txt文件：matrix_xishu.txt，内容如下：

import numpy as np
import matplotlib.pyplot as plt
import networkx as nx

def SpMV(rptA, cidA, x, y):
    y[:] = np.zeros_like(y)
    for i in range(len(rptA) - 1):
        y[i] = np.sum(x[cidA[rptA[i]:rptA[i+1]]])
    y[y != 0] = 1

def BFS(rpt, cid, start_node):
    x = np.zeros(m)
    y = np.zeros(m)
    x[start_node] = 1

    node_colors[start_node] = 'lightgreen'
    edge_colors = ['black'] * len(G.edges)

    plt.clf()
    nx.draw(G, pos, with_labels=True, edge_color=edge_colors, node_color=node_colors, node_size=500)
    plt.pause(1)

    loop_num = 0
    while True:
        print('level', loop_num, ':')
        loop_num += 1
        SpMV(rpt, cid, x, y)

        xid = np.where(x == 1)[0]
        yid = np.where(y == 1)[0]
        print('vector x:', x, '->input:', xid)
        print('vector y:', y, '->output:', yid)

        if np.sum(y) == 0:
            edge_colors = ['black'] * len(G.edges)
            plt.clf()
            nx.draw(G, pos, with_labels=True, edge_color=edge_colors, node_color=node_colors, node_size=500)
            plt.pause(1)
            break
        else:
            for i in xid:
                temp = 0
                for j, edge in enumerate(G.edges):
                    if edge[0] == i:
                        edge_colors[temp] = 'red'
                    temp += 1
            for i in yid:
                node_colors[i] = 'lightgreen'
            plt.clf()
            nx.draw(G, pos, with_labels=True, edge_color=edge_colors, node_color=node_colors, node_size=500)
            plt.pause(1)

            x = y
            y = np.zeros(m)

def readGraphfile():
    with open(filename, 'r') as f:
        m, n, num = map(int, f.readline().split())
        adj_matrix = np.zeros((m, n))
        for line in f.readlines():
            src, dst = map(int, line.split())
            adj_matrix[dst][src] = 1
    return m, n, num, adj_matrix

if __name__ == '__main__':
    filename = input("Please input the graph file name:")

    print('Read file and create the adjacency matrix...')
    m, n, num, adj_matrix = readGraphfile()
    print('The size (n) of the adjacency matrix:', n)
    print('The number of nonzeros (nnz) in the adjacency matrix:', num)
    print('')

    csrRowPtr = np.zeros(m + 1, dtype=int)
    csrColIdx = np.zeros(num, dtype=int)
    cooRowIdx = np.zeros(num, dtype=int)
    cooColIdx = np.zeros(num, dtype=int)
    nnz_ptr = 0

    for i in range(m):
        csrRowPtr[i + 1] = csrRowPtr[i] + np.sum(adj_matrix[i] == 1)
        indices = np.where(adj_matrix[i] == 1)[0]
        csrColIdx[nnz_ptr:nnz_ptr+len(indices)] = indices
        cooRowIdx[nnz_ptr:nnz_ptr+len(indices)] = i
        cooColIdx[nnz_ptr:nnz_ptr+len(indices)] = indices
        nnz_ptr += len(indices)
    print('===================================')
    print('Coordinate (COO) format:')
    print('RowIdx:', cooRowIdx)
    print('ColIdx:', cooColIdx)
    print('The total length of COO format data: len(RowIdx) + len(ColIdx) -> nnz + nnz ->', len(cooRowIdx) + len(cooColIdx))
    print('')

    print('===================================')
    print('Compressed Sparse Row (CSR) format:')
    print('RowPtr:', csrRowPtr)
    print('ColIdx:', csrColIdx)
    print('The total length of CSR format data: len(RowPtr) + len(ColIdx) -> (n + 1) + nnz ->', len(csrRowPtr) + len(csrColIdx))

    G = nx.DiGraph()
    G.add_nodes_from(range(m))
    edges = [(csrColIdx[j], i) for i in range(m) for j in range(csrRowPtr[i], csrRowPtr[i + 1])]
    G.add_edges_from(edges)
    pos = nx.shell_layout(G)

    node_colors = ['grey'] * len(G.nodes)
    edge_colors = ['black'] * len(G.edges)

    plt.figure()
    nx.draw(G, pos, with_labels=True, edge_color=edge_colors, node_color=node_colors, node_size=500)

    print('===================================')
    print('Breadth First Search (BFS) Process:')
    start_node = 0
    BFS(csrRowPtr, csrColIdx, start_node)

    plt.show()

Assignment 06

在这里插入图片描述

import numpy as np
import matplotlib.pyplot as plt

u1 = np.array([2,5,-1])
u2 = np.array([-2,1,1])
y = np.array([1,2,3])

projection_u1 = ((np.dot(y,u1)/np.dot(u1,u1))*u1)
projection_u2 = ((np.dot(y,u2)/np.dot(u2,u2))*u2)
projection = projection_u1 + projection_u2

fig = plt.figure()
ax = fig.add_subplot(111,projection='3d')

ax.quiver(0,0,0,u1[0],u1[1],u1[2],color='r',label='u1')
ax.quiver(0,0,0,u2[0],u2[1],u2[2],color='b',label='u2')
ax.quiver(0,0,0,y[0],y[1],y[2],color='g',label='y')
ax.quiver(0,0,0,projection[0],projection[1],projection[2],color='b',label='projection')

for vec in [projection_u1,projection_u2]:
    ax.plot([y[0],vec[0]],[y[1],vec[1]],[y[2],vec[2]],linestyle='dashed',color='gray')

ax.set_xlim([-3,3])
ax.set_ylim([-3,3])
ax.set_zlim([-3,3])

ax.set_xlabel('X')
ax.set_ylabel('Y')
ax.set_zlabel('Z')
ax.legend()

plt.show()

Assignment 07

在这里插入图片描述

import matplotlib.pyplot as plt
import numpy as np
def draw_basis_vector(u1,u2,u3,X,Y,Z):
    fig=plt.figure(figsize=(12,12))
    ax=fig.add_subplot(projection='3d')
    ax.plot_wireframe(X,Y,Z,linewidth=1.5,alpha=0.3)

    vectors = [
        {'data':u1,'label':r'$\mathbf{u}_1$','color':'red'},
        {'data':u2,'label':r'$\mathbf{u}_2$','color':'blue'},
        {'data':u3,'label':r'$\mathbf{u}_3$','color':'green'}
    ]

    for vector in vectors:
        ax.quiver(0,0,0,vector['data'][0],vector['data'][1],vector['data'][2],length=1,normalize=False,
                  color=vector['color'],alpha=0.6,arrow_length_ratio=0.08,pivot='tail',linewidth=3)
        ax.text(vector['data'][0],vector['data'][1],vector['data'][2],vector['label'],size=15)

    ax.set_xlabel('x-axis')
    ax.set_ylabel('y-axis')
    ax.set_zlabel('z-axis')

    ax.set_xlim([-10,10])
    ax.set_ylim([-10,10])
    ax.set_zlim([-10,10])
    plt.show()

def draw_orthogonal_basis_vector_v1_v2(u1,u2,u3,X,Y,Z):
    v1=u1
    v2=u2-(u2@v1)/(v1@v1)
    alpha=(u2@v1)/(v1@v1)
    projW_u2=[alpha*v1[0],alpha*v1[1],alpha*v1[2]]
    fig=plt.figure(figsize=(12,12))
    ax=fig.add_subplot(projection='3d')
    ax.plot_wireframe(X,Y,Z,linewidth=1.5,alpha=0.3)

    vectors = [
        {'data':u1,'label':r'$\mathbf{u}_1=\mathbf{v}_1$','color':'red'},
        {'data':u2,'label':r'$\mathbf{u}_2$','color':'blue'},
        {'data':u3,'label':r'%\mathbf{u}_3$','color':'green'},
        {'data':v2,'label':r'$\mathbf{v}_2$','color':'purple'},
        {'data':projW_u2,'label':r'$\mathbf{\hat{u}}_2$','color':'blue'}
    ]
    for vector in vectors:
        ax.quiver(0,0,0,vector['data'][0],vector['data'][1],vector['data'][2],length=1,normalize=False,
                  color=vector['color'],alpha=0.6,arrow_length_ratio=0.08,pivot='tail',linewidth=3)
        ax.text(vector['data'][0],vector['data'][1],vector['data'][2],vector['label'],size=15)
    ax.plot([projW_u2[0],u2[0]],[projW_u2[1],u2[1]],[projW_u2[2],u2[2]],c='b',lw=3.5,alpha=0.5,ls='--')
    ax.plot([v2[0],u2[0]],[v2[1],u2[1]],[v2[2],u2[2]],c='b',lw=3.5,alpha=0.5,ls='--')
    ax.set_xlabel('x-axis')
    ax.set_ylabel('y-axis')
    ax.set_zlabel('z-axis')

    ax.set_xlim([-10, 10])
    ax.set_ylim([-10, 10])
    ax.set_zlim([-10, 10])
    plt.show()
    return v1,v2,projW_u2
def draw_orthogonal_basis_vector_v1_v3(u1,u2,u3,v1,v2,X,Y,Z,projW_u2):
    projW_u3=(u3@v1)/(v1@v1)*v1+(u3@v2)/(v2@v2)*v2
    v3=u3-projW_u3
    fig=plt.figure(figsize=(12,12))
    ax=fig.add_subplot(projection='3d')
    ax.plot_wireframe(X,Y,Z,linewidth=1.5,alpha=0.3)
    vectors = [
        {'data':u1,'label':r'$\mathbf{u}_1=\mathbf{v}_1$','color':'red'},
        {'data':u2,'label':r'$\mathbf{u}_2$','color':'red'},
        {'data':u3,'label':r'$\mathbf{u}_3$','color':'red'},
        {'data':v2,'label':r'$$\mathf{v}_2$','color':'purple'},
        {'data':projW_u3,'label':r'$\hat{\mathbf{u}}_3$','color':'black'},
        {'data':projW_u2,'label':r'$\mathbf{\hat{u}}_2$','color':'blue'},
        {'data':v3,'label':r'$\mathbf{v}_3$','color':'purple'}
    ]
    for vector in vectors:
        ax.quiver(0,0,0,vector['data'][0],vector['data'][1],vector['data'][2],length=1,normalize=False,
                  color=vector['color'],alpha=0.6,arrow_length_ratio=0.08,pivot='tail',linewidth=3)
        ax.text(vector['data'][0],vector['data'][1],vector['data'][2],vector['label'],size=15)

    ax.plot([projW_u2[0],u2[0]],[projW_u2[1],u2[1]],[projW_u2[2],u2[2]],c='b',lw=3.5,alpha=0.5,ls='--')
    ax.plot([v2[0],u2[0]],[v2[1],u2[1]],[v2[2],u2[2]],c='b',lw=3.5,alpha=0.5,ls='--')
    ax.plot([projW_u3[0],u3[0]],[projW_u3[1],u3[1]],[projW_u3[2],u3[2]],c='b',lw=3.5,alpha=0.5,ls='--')
    ax.set_xlabel('x-axis')
    ax.set_ylabel('y-axis')
    ax.set_zlabel('z-axis')

    ax.set_xlim([-10, 10])
    ax.set_ylim([-10, 10])
    ax.set_zlim([-10, 10])
    plt.show()

u1 = np.array([5, 6, 4], dtype=np.float64)
u2 = np.array([2, 4, 8], dtype=np.float64)
u3 = np.array([6, -6, 4], dtype=np.float64)
s = np.linspace(-1,1,10)
t = np.linspace(-1,1,10)
S,T = np.meshgrid(s,t)
X = u1[0]*S+u2[0]*T
Y = u1[1]*S+u2[1]*T
Z = u1[2]*S+u2[2]*T
draw_basis_vector(u1,u2,u3,X,Y,Z)
v1,v2,projW_u2 = draw_orthogonal_basis_vector_v1_v2(u1,u2,u3,X,Y,Z)
draw_orthogonal_basis_vector_v1_v3(u1,u2,u3,v1,v2,X,Y,Z,projW_u2)

Assignment 08

在这里插入图片描述

import numpy as np

def compute_least_squares(A,b):
    Q,R = np.linalg.qr(A)
    x = np.linalg.solve(R,np.dot(Q.T,b))

    print("Least-squares solutions:",x)
A = np.array([[1,-1],
              [1,4],
              [1,-1],
              [1,4]])
b = np.array([-1,6,5,7])

compute_least_squares(A,b)

Assignment 09

在这里插入图片描述

import matplotlib.pyplot as plt
import numpy as np
from mpl_toolkits.mplot3d.art3d import Poly3DCollection
vertices=[]
def plot_arrow(ax,start,scale,color):
    end=[start[0]+scale[0],start[1]+scale[1],start[2]+scale[2]]

    ax.quiver(start[0],start[1],start[2],scale[0],scale[1],scale[2],color=color,arrow_length_ratio=0.08,linewidths=3)
    ax.text(end[0],end[1],end[2],f"({end[0]},{end[1]},{end[2]})")

    if end not in vertices:
        vertices.append(end)

def plot_3D_det(matrix):
    fig=plt.figure(figsize=(8,8))
    ax=fig.add_subplot(111,projection='3d')

    det=np.linalg.det(matrix)

    ax.text(0,0,0,f"({0,0},{0,0},{0,0})")
    vertices.append([0,0,0])
    color=['green','blue','red']

    for i in range(3):
        start=[[0,0,0],matrix[(i+1)%3],matrix[(i+2)%3],[matrix[(i+1)%3][j]+matrix[(i+2)%3][j] for j in range(3)]]
        for j in range(4):
            plot_arrow(ax,start[j],matrix[i],color[i])

    faces=[[vertices[1],vertices[3],vertices[4],vertices[2]],
           [vertices[5],vertices[6],vertices[4],vertices[2]],
           [vertices[0],vertices[7],vertices[6],vertices[5]],
           [vertices[0],vertices[1],vertices[3],vertices[7]],
           [vertices[0],vertices[7],vertices[3],vertices[1]],
           [vertices[7],vertices[6],vertices[4],vertices[3]]]
    for face in faces:
        ax.add_collection3d(Poly3DCollection([face],facecolors='gray',linewidths=1,edgecolors='k',alpha=0.1))

    ax.text(matrix[0][0]/2,matrix[1][1]/2,matrix[2][2]/2,f'determinant={det:.2f}',size=12)

    ax.set_xlabel('X')
    ax.set_ylabel('Y')
    ax.set_zlabel('Z')
    ax.grid()
    plt.show()

matrix=np.array([[2,0,0],[0,3,0],[0,3,4]])
plot_3D_det(matrix)

Assignment 10

在这里插入图片描述

import random
import math
import numpy as np
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D

def matrix_vector_product(matrix,vector):
    return np.dot(matrix,vector)
def normalize(vector):
    norm=np.linalg.norm(vector)
    if norm == 0:
        return vector
    return vector / norm
def dot_product(vector1,vector2):
    return np.dot(vector1,vector2)
def print_vector(vector):
    n = len(vector)
    for i in range(n):
        if i != n-1:
            print("{:.4f},".format(vector[i]),end=" ")
        else:
            print("{:.4f}".format(vector[i]),end="")
    print("]")

def power_method(matrix,max_iterations=2000,tolerance=1e-5):
    initial_vector=[1.0,1.0,1.0]
    right_vector=[0.41,0.41,0.82]

    fig=plt.figure()
    ax=fig.add_subplot(111,projection='3d')

    for iteration in range(max_iterations):
        eigenvector = matrix_vector_product(matrix,initial_vector)
        eigenvector=normalize(eigenvector)
        eigenvalue=dot_product(eigenvector,matrix_vector_product(matrix,eigenvector))
        diff_squared_sum=np.sum((eigenvector-initial_vector) ** 2)
        if math.sqrt(diff_squared_sum)<tolerance:
            break
        print("Iter",iteration,end=" ")
        print(",eigenvector=[",end="")
        print_vector(eigenvector)

        ax.cla()

        ax.quiver(0,0,0,eigenvector[0],eigenvector[1],eigenvector[2],color='blue',label='eigen vector',
                  alpha=0.8,
                  lw=1.5,arrow_length_ratio=0.2)
        ax.quiver(0,0,0,right_vector[0],right_vector[1],right_vector[2],color='red',label='right eigen vector',
                  alpha=0.8,lw=1.5,arrow_length_ratio=0.2)

        ax.text(eigenvector[0],eigenvector[1],eigenvector[2],
                f'({eigenvector[0]:.2f},{eigenvector[1]:.2f},{eigenvector[2]:.2f})')
        ax.text(right_vector[0],right_vector[1],right_vector[2],
                f'({right_vector[0]:.2f},{right_vector[1]:.2f},{right_vector[2]:.2f})')

        ax.set_xlim([-1,1])
        ax.set_ylim([-1,1])
        ax.set_zlim([-1,1])

        ax.set_xlabel('X')
        ax.set_ylabel('Y')
        ax.set_zlabel('Z')

        ax.legend()

        ax.set_title(f'Iteration{iteration+1}')

        ax.grid(True)

        plt.pause(2)

        initial_vector = eigenvector

    plt.close()
    return eigenvalue,eigenvector
matrix=[[1,-3,3],
        [3,-5,3],
        [6,-6,4]]
max_eigenvalue,eigenvector=power_method(matrix)
print("The Maximum Eigenvalue = ","{:.4f}".format(max_eigenvalue))
print("The Corresponding Eigenvector = [",end="")
print_vector(eigenvector)

Assignment 11

在这里插入图片描述

from math import sqrt
import matplotlib.pyplot as plt
import seaborn as sns

def dot_product(vector1,vector2):
    result=0
    for i in range(len(vector1)):
        result += vector1[i]*vector2[i]
    return result
def outer_product(vector1,vector2):
    result=[[0]*len(vector2) for _ in range(len(vector1))]
    for i in range(len(vector1)):
        for j in range(len(vector2)):
            result[i][j]=vector1[i]*vector2[j]
    return result
def matrix_vector_multiply(matrix,vector):
    rows,cols=len(matrix),len(matrix[0])
    result=[0]*rows
    for i in range(rows):
        for j in range(cols):
            result[i]+=matrix[i][j]*vector[j]
    return result
def matrix_matrix_multiply(matrix1,matrix2):
    rows1=len(matrix1)
    cols1=len(matrix1[0])
    cols2=len(matrix2[0])
    result=[[0]*cols2 for _ in range(rows1)]
    for i in range(rows1):
        for j in range(cols2):
            for k in range(cols1):
                result[i][j] += matrix1[i][k]*matrix2[k][j]
    return result
def matrix_transpose(matrix):
    rows,cols=len(matrix),len(matrix[0])
    transposed=[[0]*rows for _ in range(cols)]
    for i in range(rows):
        for j in range(cols):
            transposed[j][i]=matrix[i][j]
    return transposed
def norm(vector):
    result=0
    for x in vector:
        result +=x ** 2
    return sqrt(result)
def normalize(vector):
    result_norm=norm(vector)
    return [element / result_norm for element in vector]
def print_vector(vector):
    n = len(vector)
    for i in range(n):
        if i != n-1:
            print("{:.4f},".format(vector[i]),end=" ")
        else:
            print("{:.4f}".format(vector[i]),end="")
    print("]")

def power_method(matrix,max_iterations=2000,tolerance=1e-8):
    initial_vector=[1,1,0]
    diff=[0,0,0]
    eigenvector=initial_vector
    eigenvalue=0
    final_iteration=0
    for iteration in range(max_iterations):
        eigenvector=matrix_vector_multiply(matrix,initial_vector)
        eigenvector=normalize(eigenvector)
        eigenvalue=dot_product(eigenvector,matrix_vector_multiply(matrix,eigenvector))

        for i in range(len(eigenvector)):
            diff[i]=eigenvector[i]-initial_vector[i]
        if norm(diff) < tolerance:
            final_iteration = iteration
            break
        initial_vector = eigenvector
    print("Iter",final_iteration,end=" ")
    print(",eigenvalue = ","{:.4f},".format(eigenvalue),end=" ")
    print(",eigenvector = [",end="")
    print_vector(eigenvector)
    return eigenvector,eigenvalue
def svd(A):
    n,m = len(A),len(A[0])
    k=min(n,m)
    singularValues = [0] * k
    us = [[0]*k for _ in range(n)]
    vs = [[0]*m for _ in range(k)]
    matrixForsvd = [[0] * m for _ in range(n)]
    for i in range(n):
        for j in range(m):
            matrixForsvd[i][j] = A[i][j]
    for i in range(k):
        v,eigenvalue=power_method(
            matrix_matrix_multiply(matrix_transpose(matrixForsvd),matrixForsvd))
        sigma=sqrt(eigenvalue)
        u = [element / sigma for element in matrix_vector_multiply(A,v)]
        singularValues[i] = sigma
        for row in range(len(us)):
            us[row][i] = u[row]
        for col in range(len(vs[0])):
            vs[i][col] = v[col]
        outer_u_v = outer_product(u,v)
        for p in range(n):
            for q in range(m):
                matrixForsvd[p][q] -= sigma * outer_u_v[p][q]
    return us, singularValues, vs
def add_annotations(matrix):
    for i in range(len(matrix)):
        for j in range(len(matrix[0])):
            if abs(matrix[i][j]) <= 1e-8 and abs(matrix[i][j]) != 0:
                plt.text(j + 0.5,i+0.5, '~0,00',
                         ha='center',va='center',color='black',fontsize=12)
            elif abs(matrix[i][j]) == 0:
                plt.text(j+0.5,i+0.5,'Not Stored',
                         ha='center',va='center',color='black',fontsize=7)
            else:
                plt.text(j+0.5,i+0.5,f'{matrix[i][j]:.2f}',
                         ha='center',va='center',color='black',fontsize=12)
def draw(A,u,sigma,v,axs,index):
    title_A = ['$A=U*S*V$','$U\' *S\'*V\'$', '$U\'\'*S\'\'*V\'\'$',
               '$U*S*V$']
    title_U = ['$U$','$U\' \ (Compression \ rate \ 66 \\%)$', '$U\'\' \ (Compression \ rate \ 33\\%)$',
               '$U$']
    title_S = ['$S$','$S\' \ (Compression \ rate \ 66 \\%)$', '$S\'\' \ (Compression \ rate \ 33\\%)$',
               '$S$']
    title_V = ['$V$','$V\' \ (Compression \ rate \ 66 \\%)$', '$V\'\' \ (Compression \ rate \ 33\\%)$',
               '$V$']

    plt.sca(axs[0])
    ax=sns.heatmap(A,cmap='Wistia',vmax=2,vmin=1)
    ax.set_aspect(0.6)
    plt.title(title_A[index],size=10)
    add_annotations(A)
    ax.tick_params(axis='both',labelsize=8)
    cbar=ax.collections[0].colorbar
    cbar.ax.tick_params(labelsize=8)

    plt.sca(axs[1])
    ax=sns.heatmap(u,cmap='Wistia',vmax=1,vmin=-1)
    ax.set_aspect(0.6)
    plt.title(title_U[index],size=10)
    add_annotations(u)
    ax.tick_params(axis='both',labelsize=8)
    cbar=ax.collections[0].colorbar
    cbar.ax.tick_params(labelsize=8)

    plt.sca(axs[2])
    ax=sns.heatmap(sigma,cmap='Wistia',vmax=7,vmin=0)
    ax.set_aspect(0.6)
    plt.title(title_S[index],size=10)
    add_annotations(sigma)
    ax.tick_params(axis='both',labelsize=8)
    cbar=ax.collections[0].colorbar
    cbar.ax.tick_params(labelsize=8)

    plt.sca(axs[3])
    ax=sns.heatmap(v,cmap='Wistia',vmax=1,vmin=-1)
    ax.set_aspect(0.6)
    plt.title(title_V[index],size=10)
    add_annotations(v)
    ax.tick_params(axis='both',labelsize=8)
    cbar=ax.collections[0].colorbar
    cbar.ax.tick_params(labelsize=8)

A=[[2.0,2.0,2.0],
   [2.0,1.0,1.0],
   [2.0,2.0,2.0],
   [1.0,1.0,2.0],
   [2.0,2.0,2.0]]

fig,axs = plt.subplots(nrows=4,ncols=4,figsize=(12,10))
m=len(A)
n=len(A[0])

u,sigma,v = svd(A)
sigma = [[sigma[i] if j == i else 0 for j in range(n)] for i in range(n)]

draw(A,u,sigma,v,axs[0],0)
u_new = [[0] * n for _ in range(m)]
sigma_new = [[0] * n for _ in range(n)]
v_new = [[0] * n for _ in range(n)]
for k in range(0,len(sigma)):
    for i in range(len(u)):
        for j in range(k+1):
            u_new[i][j] = u[i][j]
    for i in range(k+1):
        for j in range(len(sigma[0])):
            sigma_new[i][j] = sigma[i][j]
    for i in range(k+1):
        for j in range(len(v[0])):
            v_new[i][j] = v[i][j]
    A_new = matrix_matrix_multiply(matrix_matrix_multiply(u_new,sigma_new),v_new)

    draw(A_new,u_new,sigma_new,v_new,axs[k+1],k+1)

plt.tight_layout()
plt.show()