"""Owned reduced SVD of a real two-column full-rank matrix; stdlib only."""
import math,json
def svd_two_columns(A):
    if len(A)<2 or any(len(r)!=2 for r in A):raise ValueError('only m-by-2, m>=2')
    if not all(math.isfinite(x) for r in A for x in r):raise ValueError('finite only')
    a=sum(r[0]**2 for r in A);b=sum(r[0]*r[1] for r in A);d=sum(r[1]**2 for r in A)
    radius=math.hypot((a-d)/2,b);eig=[(a+d)/2+radius,(a+d)/2-radius]
    if eig[1]<=1e-10 or radius<=1e-10:raise ValueError('rank-deficient or repeated values outside educational subset')
    theta=.5*math.atan2(2*b,a-d);c,s=math.cos(theta),math.sin(theta)
    V=[[c,-s],[s,c]];S=[math.sqrt(x) for x in eig]
    U=[[sum(row[k]*V[k][j] for k in range(2))/S[j] for j in range(2)] for row in A]
    return U,S,V
def reconstruct(U,S,V,keep):
    return [[sum(U[i][k]*S[k]*V[j][k] for k in range(keep)) for j in range(2)] for i in range(len(U))]
def gram(M):return [[sum(r[i]*r[j] for r in M) for j in range(2)] for i in range(2)]
def demo():
    A=[[1.92,.44],[1.56,1.92],[0.,0.]]
    U,S,V=svd_two_columns(A);full=reconstruct(U,S,V,2);one=reconstruct(U,S,V,1)
    error=lambda B:sum((A[i][j]-B[i][j])**2 for i in range(3) for j in range(2))
    for G in [gram(U),gram(V)]:
        assert max(abs(G[i][j]-(i==j)) for i in range(2) for j in range(2))<1e-12
    assert abs(S[0]-3)<1e-12 and abs(S[1]-1)<1e-12 and error(full)<1e-24
    assert abs(error(one)-1)<1e-12
    flippedU=[[-r[0],r[1]] for r in U];flippedV=[[-r[0],r[1]] for r in V]
    assert max(abs(reconstruct(flippedU,S,flippedV,2)[i][j]-full[i][j]) for i in range(3) for j in range(2))<1e-12
    rejected=[]
    for name,B in [('rank_one',[[1.,2.],[2.,4.]]),('equal_values',[[1.,0.],[0.,1.]])]:
        try:svd_two_columns(B)
        except ValueError:rejected.append(name)
        else:raise AssertionError(name)
    return {'input_uncentered_A':A,'reduced_U_3_by_2':U,'singular_values':S,'V_columns_2_by_2':V,'U_transpose_U':gram(U),'V_transpose_V':gram(V),'full_reconstruction':full,'rank_one_reconstruction':one,'full_squared_Frobenius_error':error(full),'rank_one_squared_Frobenius_error':error(one),'rank_one_Frobenius_error':math.sqrt(error(one)),'both_sides_first_column_sign_flip_preserves_reconstruction':True,'rejected_by_own_limited_solver_not_by_SVD_definition':rejected,'NumPy_or_training_executed':False}
if __name__=='__main__':print(json.dumps(demo(),indent=2))
