I am trying to find the polar decomposition of a matrix.
I tried studying it from a video lecture : https://www.youtube.com/watch?v=VjBEFOwHJoo
I tried implementing the way they have done, trusting their calculations of Eigen values and Eigen vectors as shown below and it gives correct value of U
and P
matrix.
A = np.array([[4,-14,1],[22,23,18],[15,10,20]])
A_T_multiply_A = np.array(A.T @ A)
print("A transpose into A")
print(A_T_multiply_A)
x = np.array([[-1,1,1],[0,-2,1],[1,1,1]]) # Eigen vectors
y = np.array([[25,0,0],[0,225,0],[0,0,2025]]) # Eigen values
z = np.linalg.inv(np.array([[-1,1,1],[0,-2,1],[1,1,1]])) # Inverse Eigen vectors
A_transpose_A = x @ y @ z
print("A transpose A by multiplying Eigen vectors, Eigen values and Inv of Eigen vectors")
print(A_transpose_A)
P = x @ np.sqrt(y) @ z
print("Matrix P")
print(P)
U = A @ np.linalg.inv(P)
print("Matrix U")
print(U)
But when I try to implement the Eigen vectors and Eigen values using the API from numpy, the thing do not match:
import numpy as np
# Define the matrix
A = np.array([[4,-14,1],
[22,23,18],
[15,10,20]])
A_T_multiply_A = np.array(A.T @ A)
# Compute eigenvalues and eigenvectors
eigenvalues, eigenvectors = np.linalg.eig(A_T_multiply_A)
# Print the results
print("Eigenvalues:")
print(eigenvalues)
print("\nEigenvectors:")
print(eigenvectors)
Not sure why I am getting this difference? It would be great if someone can explain me this or provide me some good link.
The values of Eigen values and Eigen vectors I got using numpy :
Eigenvalues:
[2025. 25. 225.]
Eigenvectors:
[[-5.77350269e-01 -7.07106781e-01 4.08248290e-01]
[-5.77350269e-01 -3.91901695e-16 -8.16496581e-01]
[-5.77350269e-01 7.07106781e-01 4.08248290e-01]]
Also, providing the python code for calculating polar decomposition using a numpy API.
import numpy as np
from scipy.linalg import polar
A = np.array([[4,-14,1],[22,23,18],[15,10,20]])
U, P = polar(A)
print("Matrix U=")
print(U)
print("Matrix P=")
print(P)
The result is:
Matrix U=
[[ 6.00000000e-01 -8.00000000e-01 2.40023768e-16]
[ 8.00000000e-01 6.00000000e-01 3.21346710e-16]
[ 2.32191141e-16 6.61562711e-17 1.00000000e+00]]
Matrix P=
[[20. 10. 15.]
[10. 25. 10.]
[15. 10. 20.]]
Here is a solution using svd
:
import numpy as np
vecU, vals, vecV = np.linalg.svd(A)
P = vecV.T @ np.diag(vals) @ vecV
U = vecU @ vecV
P
array([[20., 10., 15.],
[10., 25., 10.],
[15., 10., 20.]])
U
array([[ 6.00000000e-01, -8.00000000e-01, 5.73090676e-16],
[ 8.00000000e-01, 6.00000000e-01, 3.76857861e-16],
[ 2.32191141e-16, -5.92653245e-17, 1.00000000e+00]])
Compare the results to what you have and you will notice they are the same
Solution using Eigen decomposition:
Evals, Evec = np.linalg.eig(A.T @ A)
P = Evec @ np.diag(np.sqrt(Evals)) @ Evec.T
U = A @ Evec @ np.diag(1/np.sqrt(Evals)) @ Evec.T
P
array([[20., 10., 15.],
[10., 25., 10.],
[15., 10., 20.]])
U
array([[ 6.00000000e-01, -8.00000000e-01, 3.33066907e-16],
[ 8.00000000e-01, 6.00000000e-01, -2.63677968e-16],
[-9.71445147e-16, -1.66533454e-16, 1.00000000e+00]])