int main() { int m = 3 , n = 3; int lmn = min(m,n); array double complex a[3][3] = {7, 8, 1, 3, 6, 4, 3, 5, 7} ; /* m-by-n matrix */ array complex z[3][3] = {complex(1,2), 8, 1, 3, 6, 4, 3, 5, 7} ; /* m-by-n matrix */ array double s[lmn]; array double complex u[m][m]; array double complex vt[n][n]; int status; svd(a, s, NULL, NULL); printf("s =\n%f\n", s); svd(a, s, u, vt); printf("u =\n%f\n", u); printf("s =\n%f\n", s); printf("vt =\n%f\n", vt); printf("u*s*v^t =\n%f\n", u*diagonalmatrix(s)*transpose(vt)); status = svd(z, s, u, vt); if(status ==0) { printf("u =\n%f\n", u); printf("s =\n%f\n", s); printf("vt =\n%f\n", vt); printf("u*s*v^t =\n%f\n", u*diagonalmatrix(s)*transpose(vt)); } else printf("error: numerical error in svd()\n"); }
s = 15.071167 5.471953 0.957939 u = complex(-0.660443,0.000000) complex(0.704551,0.000000) complex(0.259659,0.000000) complex(-0.511773,0.000000) complex(-0.169318,0.000000) complex(-0.842271,0.000000) complex(-0.549458,0.000000) complex(-0.689158,0.000000) complex(0.472395,0.000000) s = 15.071167 5.471953 0.957939 vt = complex(-0.517995,0.000000) complex(0.430638,0.000000) complex(0.739075,-0.000000) complex(-0.736603,0.000000) complex(0.214678,0.000000) complex(-0.641349,0.000000) complex(-0.434853,0.000000) complex(-0.876621,0.000000) complex(0.206007,-0.000000) u*s*v^t = complex(7.000000,0.000000) complex(8.000000,0.000000) complex(1.000000,0.000000) complex(3.000000,0.000000) complex(6.000000,0.000000) complex(4.000000,0.000000) complex(3.000000,0.000000) complex(5.000000,0.000000) complex(7.000000,0.000000) u = complex(-0.504577,-0.182672) complex(0.536523,-0.575026) complex(-0.247592,0.179497) complex(-0.551746,-0.131934) complex(-0.031610,0.061460) complex(0.625115,-0.531624) complex(-0.607942,-0.143559) complex(-0.372381,0.487898) complex(-0.372831,0.306730) s = 13.678758 5.090474 0.989266 vt = complex(-0.317937,0.000000) complex(-0.358611,0.000000) complex(0.877675,-0.000000) complex(-0.759338,0.217182) complex(0.440159,0.352023) complex(-0.095224,0.222508) complex(-0.509341,0.125401) complex(-0.431508,-0.606249) complex(-0.360819,-0.202282) u*s*v^t = complex(1.000000,2.000000) complex(8.000000,0.000000) complex(1.000000,-0.000000) complex(3.000000,0.000000) complex(6.000000,0.000000) complex(4.000000,0.000000) complex(3.000000,0.000000) complex(5.000000,0.000000) complex(7.000000,0.000000)