I use the following code for a manual single value decomposition using numpy. Depending on the array I choose it sometimes works out well and I can verify the svd and sometimes it does not work out straight aways and requires a sign flip.
Code
import numpy as np
array = np.array([[-3,3,6],
[3,8,7]]) # here you can choose any matrix values. Some matrices work, some don't.
# left singular vectors
AAT = np.matmul(array, array.T)
eAAT_values, eAAT_vectors = np.linalg.eig(AAT)
idx = eAAT_values.argsort()[::-1] #sorting: largest singular values first
eAAT_values = eAAT_values[idx]
eAAT_vectors = eAAT_vectors[:,idx]
SL = eAAT_vectors
# Some arrays require an additional line:
# SL[:,0] = -SL[:,0]
# right singular vectors
ATA = np.matmul(array.T,array)
eATA_values, eATA_vectors = np.linalg.eig(ATA)
idx = eATA_values.argsort()[::-1] #sorting: largest singular values first
eATA_values = eATA_values[idx]
eATA_vectors = eATA_vectors[:,idx]
SR = eATA_vectors.T
# singular values
S0 = np.zeros(np.shape(array))
np.fill_diagonal(S0, np.sqrt(eATA_values), wrap=True)
# verifying. Expected proof == array
Proof = np.matmul(SL,np.matmul(S0,SR)) # works out with some arrays with others it does not
# Python linalg.svd works consistently
U, S, Vt = np.linalg.svd(array)
Sm = np.zeros(np.shape(array))
np.fill_diagonal(Sm, S, wrap=True)
proof2 = np.matmul(U, np.matmul(Sm,Vt)) # always works out
The array in the code above works well. Contrary examples where the svd does not work out without an additional sign flip are presented here.
array_not_working1 = np.array([[4,5,9],
[3,2,6]])
array_not_working2 = np.array([[3,-1,4],
[1,5,9]])
How can I stabilize the outcome or determine when a signflip of the first eigenvector from eAAT_vectors is necessary?
Another working array without an additional sign flip:
array = np.array([[7,-12,2], [-4,3,1]]).
This looks the same as I implemented an SVD, but sometimes it is not possible to restore the original matrix
You may look at my answer and the comments below.