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)