08-Oct-2025 20:02:38 svd_test MATLAB/Octave version 6.4.0 svd() computes the singular value decomposition. svd_demo(): Test the SVD for a 3 by 3 matrix A real MxN matrix A can be factored as: A = U * S * V' where U = MxM orthogonal, S = MxN zero except for diagonal, V = NxN orthogonal. The diagonal of S contains only nonnegative numbers and these are arranged in descending order. Matrix row order M = 5 Matrix column order N = 3 Choose a "random" matrix A, with integral values. The matrix A: 3 -4 5 -8 1 0 -11 -10 1 -16 -4 -17 13 3 -15 Orthogonal factor U: -0.127077 0.192620 -0.522197 0.791238 -0.219092 0.250691 0.112051 0.434275 0.493843 0.701515 0.452779 0.293150 -0.657383 -0.315991 0.420774 0.811162 -0.415108 0.055397 0.164438 -0.373622 -0.240858 -0.831920 -0.321704 0.056318 0.378459 Diagonal factor S: Diagonal Matrix 26.8018 0 0 0 22.6148 0 0 0 8.4400 0 0 0 0 0 0 Orthogonal factor V: -0.875951 -0.341212 -0.341006 -0.288638 -0.195680 0.937229 -0.386522 0.919393 0.072919 svd_product_test(): Test the value of |A-U*S*V'| The product U * X * V': 3.0000e+00 -4.0000e+00 5.0000e+00 -8.0000e+00 1.0000e+00 -7.3353e-16 -1.1000e+01 -1.0000e+01 1.0000e+00 -1.6000e+01 -4.0000e+00 -1.7000e+01 1.3000e+01 3.0000e+00 -1.5000e+01 Frobenius Norm of A, A_NORM = 36.069378 ABSOLUTE ERROR for A = U*S*V': Frobenius norm of difference A-U*S*V' = 0.000000 RELATIVE ERROR for A = U*S*V': Ratio of DIF_NORM / A_NORM = 0.000000 rank_one_test(): Compare A to the sum of R rank one matrices. R Absolute Relative Error Error 0 36.069378 1.000000 1 24.138390 0.669221 2 8.439984 0.233993 3 0.000000 0.000000 rank_one_print_test(): Print the sums of R rank one matrices. Rank R = 0: 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 Rank R = 1: 2.9834 0.9831 1.3165 -5.8855 -1.9393 -2.5970 -10.6299 -3.5027 -4.6906 -19.0437 -6.2752 -8.4032 5.6547 1.8633 2.4952 Rank R = 2: 1.4971 0.1307 5.3214 -6.7501 -2.4352 -0.2673 -12.8920 -4.8000 1.4046 -15.8406 -4.4382 -17.0341 12.0741 5.5447 -14.8020 Rank R = 3: 3.0000e+00 -4.0000e+00 5.0000e+00 -8.0000e+00 1.0000e+00 -7.3353e-16 -1.1000e+01 -1.0000e+01 1.0000e+00 -1.6000e+01 -4.0000e+00 -1.7000e+01 1.3000e+01 3.0000e+00 -1.5000e+01 Original matrix A: 3 -4 5 -8 1 0 -11 -10 1 -16 -4 -17 13 3 -15 Pseudoinverse of A: 2.2346e-02 -2.7430e-02 7.3397e-03 -2.2486e-02 3.3422e-02 -5.8286e-02 4.4555e-02 -8.0413e-02 1.0078e-03 -2.5932e-02 5.1519e-03 4.6921e-03 -2.9148e-04 -2.8096e-02 -3.3127e-02 pseudo_product_test() The following relations MUST hold: A * A+ * A = A A+ * A * A+ = A+ ( A * A+ ) is MxM symmetric; ( A+ * A ) is NxN symmetric Here are the Frobenius norms of the errors in these relationships: A * A+ * A = A 0.000000 A+ * A * A+ = A+ 0.000000 ( A * A+ ) is MxM symmetric; 0.000000 ( A+ * A ) is NxN symmetric; 0.000000 In some cases, the matrix A * A+ may be interesting (if M <= N, then it MIGHT look like the identity.) A * A+: 0.325941 -0.237051 0.342212 -0.211967 0.038357 -0.237051 0.263996 -0.139130 0.180895 -0.293307 0.342212 -0.139130 0.723099 0.209171 -0.141450 -0.211967 0.180895 0.209171 0.833367 0.132140 0.038357 -0.293307 -0.141450 0.132140 0.853597 In some cases, the matrix A+ * A may be interesting (if N <= M, then it MIGHT look like the identity.) A+ * A: 1.0000e+00 -1.4572e-16 -1.0408e-16 -5.0654e-16 1.0000e+00 1.8735e-16 -1.7347e-16 7.6328e-17 1.0000e+00 pseudo_linear_solve_test(): Given: b = A * x1 so that b is in the range of A, solve A * x = b using the pseudoinverse: x2 = A+ * b. Norm of x1 = 20.904545 Norm of x2 = 20.904545 Norm of residual = 0.000000 For N < M, most systems A*x=b will not be exactly and uniquely solvable, except in the least squares sense. Here is an example: Given b is NOT in the range of A, solve A * x = b using the pseudoinverse: x2 = A+ * b. Norm of x2 = 20.904545 Norm of residual = 0.000000 Given: b = A' * x1 so that b is in the range of A', solve A' * x = b using the pseudoinverse: x2 = A+' * b. Norm of x1 = 329.241553 Norm of x2 = 268.792693 Norm of residual = 0.000000 svd_demo(): Test the SVD for a 5 by 5 matrix A real MxN matrix A can be factored as: A = U * S * V' where U = MxM orthogonal, S = MxN zero except for diagonal, V = NxN orthogonal. The diagonal of S contains only nonnegative numbers and these are arranged in descending order. Matrix row order M = 5 Matrix column order N = 5 Choose a "random" matrix A, with integral values. The matrix A: 7 2 13 -2 2 -9 -2 2 -5 7 7 3 -2 9 -7 4 -4 0 1 -1 8 -24 21 3 -7 Orthogonal factor U: -2.5336e-01 -4.4674e-02 9.6022e-01 -3.0595e-02 -1.0415e-01 5.4258e-02 -6.6594e-01 -6.6833e-02 4.4431e-01 -5.9305e-01 -2.8271e-02 7.3068e-01 -1.9417e-03 5.3722e-01 -4.2036e-01 -1.1674e-01 1.2906e-01 -1.2038e-01 -7.0872e-01 -6.7301e-01 -9.5835e-01 -6.3168e-02 -2.4291e-01 1.0373e-01 8.8339e-02 Diagonal factor S: Diagonal Matrix 35.0420 0 0 0 0 0 18.6677 0 0 0 0 0 12.7733 0 0 0 0 0 3.9132 0 0 0 0 0 1.7612 Orthogonal factor V: -0.302308 0.578878 0.382413 -0.627988 -0.181387 0.649716 0.237541 0.654471 0.257404 0.163874 -0.663602 -0.251799 0.567743 0.407518 0.088468 -0.085920 0.532185 -0.192032 0.581897 -0.577855 0.196798 -0.511713 0.257333 -0.186277 -0.773628 svd_product_test(): Test the value of |A-U*S*V'| The product U * X * V': 7.0000e+00 2.0000e+00 1.3000e+01 -2.0000e+00 2.0000e+00 -9.0000e+00 -2.0000e+00 2.0000e+00 -5.0000e+00 7.0000e+00 7.0000e+00 3.0000e+00 -2.0000e+00 9.0000e+00 -7.0000e+00 4.0000e+00 -4.0000e+00 2.5619e-15 1.0000e+00 -1.0000e+00 8.0000e+00 -2.4000e+01 2.1000e+01 3.0000e+00 -7.0000e+00 Frobenius Norm of A, A_NORM = 41.928511 ABSOLUTE ERROR for A = U*S*V': Frobenius norm of difference A-U*S*V' = 0.000000 RELATIVE ERROR for A = U*S*V': Ratio of DIF_NORM / A_NORM = 0.000000 rank_one_test(): Compare A to the sum of R rank one matrices. R Absolute Relative Error Error 0 41.928511 1.000000 1 23.022920 0.549099 2 13.474857 0.321377 3 4.291281 0.102348 4 1.761157 0.042004 5 0.000000 0.000000 rank_one_print_test(): Print the sums of R rank one matrices. Rank R = 0: 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 Rank R = 1: 2.6839e+00 -5.7683e+00 5.8915e+00 7.6281e-01 -1.7472e+00 -5.7478e-01 1.2353e+00 -1.2617e+00 -1.6336e-01 3.7418e-01 2.9949e-01 -6.4365e-01 6.5741e-01 8.5118e-02 -1.9496e-01 1.2367e+00 -2.6579e+00 2.7147e+00 3.5148e-01 -8.0506e-01 1.0152e+01 -2.1819e+01 2.2285e+01 2.8854e+00 -6.6090e+00 Rank R = 2: 2.2012 -5.9664 6.1015 0.3190 -1.3204 -7.7711 -1.7177 1.8685 -6.7792 6.7355 8.1954 2.5964 -2.7771 7.3441 -7.1747 2.6313 -2.0856 2.1080 1.6336 -2.0379 9.4697 -22.0993 22.5824 2.2579 -6.0056 Rank R = 3: 6.8915 2.0609 13.0650 -2.0363 1.8358 -8.0976 -2.2764 1.3839 -6.6153 6.5159 8.1859 2.5802 -2.7912 7.3489 -7.1811 2.0434 -3.0919 1.2351 1.9289 -2.4336 8.2831 -24.1300 20.8208 2.8537 -6.8040 Rank R = 4: 6.9667 2.0301 13.0162 -2.1060 1.8581 -9.1895 -1.8288 2.0924 -5.6035 6.1920 6.8657 3.1213 -1.9345 8.5722 -7.5727 3.7850 -3.8058 0.1049 0.3151 -1.9170 8.0282 -24.0255 20.9862 3.0899 -6.8796 Rank R = 5: 7.0000e+00 2.0000e+00 1.3000e+01 -2.0000e+00 2.0000e+00 -9.0000e+00 -2.0000e+00 2.0000e+00 -5.0000e+00 7.0000e+00 7.0000e+00 3.0000e+00 -2.0000e+00 9.0000e+00 -7.0000e+00 4.0000e+00 -4.0000e+00 2.5619e-15 1.0000e+00 -1.0000e+00 8.0000e+00 -2.4000e+01 2.1000e+01 3.0000e+00 -7.0000e+00 Original matrix A: 7 2 13 -2 2 -9 -2 2 -5 7 7 3 -2 9 -7 4 -4 0 1 -1 8 -24 21 3 -7 Pseudoinverse of A: 4.5185e-02 -3.3341e-02 -2.0074e-02 1.8445e-01 -2.6708e-02 3.2230e-02 -3.6849e-02 4.8970e-03 -1.1593e-01 -1.5976e-02 3.9662e-02 2.1464e-02 2.5423e-02 -1.1249e-01 2.3443e-02 1.4535e-02 2.4254e-01 2.3874e-01 1.2121e-01 -9.3599e-03 6.6353e-02 2.5657e-01 1.3885e-01 3.2275e-01 -5.2287e-02 pseudo_product_test() The following relations MUST hold: A * A+ * A = A A+ * A * A+ = A+ ( A * A+ ) is MxM symmetric; ( A+ * A ) is NxN symmetric Here are the Frobenius norms of the errors in these relationships: A * A+ * A = A 0.000000 A+ * A * A+ = A+ 0.000000 ( A * A+ ) is MxM symmetric; 0.000000 ( A+ * A ) is NxN symmetric; 0.000000 In some cases, the matrix A * A+ may be interesting (if M <= N, then it MIGHT look like the identity.) A * A+: 1.0000 0.0000 0.0000 -0.0000 0.0000 -0.0000 1.0000 0.0000 0.0000 -0.0000 0.0000 -0.0000 1.0000 -0.0000 0.0000 0.0000 0 -0.0000 1.0000 0.0000 0.0000 -0.0000 0.0000 -0.0000 1.0000 In some cases, the matrix A+ * A may be interesting (if N <= M, then it MIGHT look like the identity.) A+ * A: 1.0000e+00 -1.1102e-16 4.4409e-16 -4.1633e-16 2.7756e-16 8.3267e-17 1.0000e+00 -4.3368e-16 -1.2143e-16 1.9082e-16 4.7184e-16 -6.9389e-16 1.0000e+00 3.0878e-16 -2.8103e-16 -3.6082e-16 -2.7756e-17 1.2837e-16 1.0000e+00 6.2797e-16 -1.0547e-15 2.7756e-16 -2.9837e-16 -2.5674e-16 1.0000e+00 pseudo_linear_solve_test(): Given: b = A * x1 so that b is in the range of A, solve A * x = b using the pseudoinverse: x2 = A+ * b. Norm of x1 = 20.928450 Norm of x2 = 20.928450 Norm of residual = 0.000000 Given: b = A' * x1 so that b is in the range of A', solve A' * x = b using the pseudoinverse: x2 = A+' * b. Norm of x1 = 137.477271 Norm of x2 = 137.477271 Norm of residual = 0.000000 svd_demo(): Test the SVD for a 8 by 8 matrix A real MxN matrix A can be factored as: A = U * S * V' where U = MxM orthogonal, S = MxN zero except for diagonal, V = NxN orthogonal. The diagonal of S contains only nonnegative numbers and these are arranged in descending order. Matrix row order M = 5 Matrix column order N = 8 Choose a "random" matrix A, with integral values. The matrix A: 7 -3 -9 11 4 5 -2 -3 11 8 35 7 6 -2 8 -1 9 -14 -2 13 13 14 22 -6 16 -11 -9 -2 -3 -14 13 0 4 -2 7 -6 3 15 8 -28 Orthogonal factor U: -0.096585 -0.283689 -0.011593 -0.344169 0.889722 -0.523542 0.737849 -0.392175 -0.021593 0.164967 -0.634823 -0.464011 -0.165857 -0.447636 -0.392183 -0.089040 -0.396274 -0.615615 0.665718 0.113478 -0.552853 -0.052539 0.663004 0.487364 0.120397 Diagonal factor S: Diagonal Matrix 46.129 0 0 0 0 0 0 0 0 39.114 0 0 0 0 0 0 0 0 30.577 0 0 0 0 0 0 0 0 22.419 0 0 0 0 0 0 0 0 11.861 0 0 0 Orthogonal factor V: Columns 1 through 7: -0.342182 -0.117507 -0.427961 0.264303 0.574181 0.482596 -0.125674 0.153354 0.452888 0.152571 -0.052230 0.223599 0.374733 0.724901 -0.417388 0.831030 -0.101663 0.029310 -0.137242 -0.154284 -0.255088 -0.205613 -0.073631 -0.254301 -0.624990 0.412617 -0.506034 0.252230 -0.285542 -0.043683 -0.023538 -0.350614 -0.044594 0.271405 -0.250739 -0.333188 -0.118387 0.554935 -0.443998 -0.097348 0.385620 -0.065820 -0.510344 -0.238023 -0.309452 0.143663 -0.560612 0.016369 0.501512 0.435781 0.111684 -0.560624 -0.441861 -0.324778 0.352378 -0.106390 Column 8: 0.194216 -0.171658 0.127562 0.084123 -0.809045 0.455382 0.020678 0.217050 svd_product_test(): Test the value of |A-U*S*V'| The product U * X * V': Columns 1 through 6: 7.0000e+00 -3.0000e+00 -9.0000e+00 1.1000e+01 4.0000e+00 5.0000e+00 1.1000e+01 8.0000e+00 3.5000e+01 7.0000e+00 6.0000e+00 -2.0000e+00 9.0000e+00 -1.4000e+01 -2.0000e+00 1.3000e+01 1.3000e+01 1.4000e+01 1.6000e+01 -1.1000e+01 -9.0000e+00 -2.0000e+00 -3.0000e+00 -1.4000e+01 4.0000e+00 -2.0000e+00 7.0000e+00 -6.0000e+00 3.0000e+00 1.5000e+01 Columns 7 and 8: -2.0000e+00 -3.0000e+00 8.0000e+00 -1.0000e+00 2.2000e+01 -6.0000e+00 1.3000e+01 1.1056e-14 8.0000e+00 -2.8000e+01 Frobenius Norm of A, A_NORM = 72.360210 ABSOLUTE ERROR for A = U*S*V': Frobenius norm of difference A-U*S*V' = 0.000000 RELATIVE ERROR for A = U*S*V': Ratio of DIF_NORM / A_NORM = 0.000000 rank_one_test(): Compare A to the sum of R rank one matrices. R Absolute Relative Error Error 0 72.360210 1.000000 1 55.750515 0.770458 2 39.727093 0.549018 3 25.363550 0.350518 4 11.860869 0.163914 5 0.000000 0.000000 rank_one_print_test(): Print the sums of R rank one matrices. Rank R = 0: 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 Rank R = 1: 1.5245 -0.6832 1.8596 0.9161 1.2722 1.4845 2.2738 -1.9416 8.2639 -3.7036 10.0801 4.9657 6.8960 8.0466 12.3250 -10.5243 10.0204 -4.4908 12.2227 6.0211 8.3617 9.7570 14.9448 -12.7613 1.4054 -0.6299 1.7143 0.8445 1.1728 1.3685 2.0961 -1.7899 8.7265 -3.9109 10.6445 5.2437 7.2820 8.4971 13.0151 -11.1135 Rank R = 2: 2.8284 -5.7085 -7.3616 1.7331 1.7569 2.7981 4.9149 -3.1808 4.8726 9.3668 34.0636 2.8407 5.6353 4.6300 5.4557 -7.3011 12.1530 -12.7103 -2.8598 7.3575 9.1545 11.9056 19.2647 -14.7883 3.2268 -7.6495 -11.1664 1.9858 1.8499 3.2035 5.7854 -3.5210 8.9680 -4.8416 8.9367 5.3950 7.3718 8.7404 13.5042 -11.3430 Rank R = 3: 2.9801 -5.7626 -7.3255 1.8232 1.7652 2.6014 5.0246 -2.9821 10.0045 7.5372 35.2827 5.8901 5.9175 -2.0245 9.1665 -0.5784 14.3234 -13.4841 -2.3443 8.6471 9.2739 9.0913 20.8340 -11.9451 11.2825 -10.5214 -9.2527 6.7726 2.2929 -7.2423 11.6104 7.0319 0.2922 -1.7486 6.8757 0.2397 6.8946 19.9903 7.2308 -22.7083 Rank R = 4: 0.9407 -5.3596 -7.5517 6.6457 4.4706 6.0273 3.9161 0.4273 9.8765 7.5625 35.2685 6.1926 6.0873 -1.8095 9.0969 -0.3645 11.6709 -12.9599 -2.6384 14.9193 12.7926 13.5472 19.3922 -7.5107 15.2272 -11.3010 -8.8153 -2.5554 -2.9400 -13.8690 13.7546 0.4371 3.1801 -2.3193 7.1960 -6.5892 3.0637 15.1390 8.8006 -27.5362 Rank R = 5: Columns 1 through 6: 7.0000e+00 -3.0000e+00 -9.0000e+00 1.1000e+01 4.0000e+00 5.0000e+00 1.1000e+01 8.0000e+00 3.5000e+01 7.0000e+00 6.0000e+00 -2.0000e+00 9.0000e+00 -1.4000e+01 -2.0000e+00 1.3000e+01 1.3000e+01 1.4000e+01 1.6000e+01 -1.1000e+01 -9.0000e+00 -2.0000e+00 -3.0000e+00 -1.4000e+01 4.0000e+00 -2.0000e+00 7.0000e+00 -6.0000e+00 3.0000e+00 1.5000e+01 Columns 7 and 8: -2.0000e+00 -3.0000e+00 8.0000e+00 -1.0000e+00 2.2000e+01 -6.0000e+00 1.3000e+01 1.1056e-14 8.0000e+00 -2.8000e+01 Original matrix A: 7 -3 -9 11 4 5 -2 -3 11 8 35 7 6 -2 8 -1 9 -14 -2 13 13 14 22 -6 16 -11 -9 -2 -3 -14 13 0 4 -2 7 -6 3 15 8 -28 Pseudoinverse of A: 4.0745e-02 1.4887e-02 -1.5838e-02 2.3809e-02 6.5532e-03 1.3911e-02 8.0062e-03 -1.4661e-02 -7.3678e-03 1.9963e-03 -1.5860e-02 1.9781e-02 3.8961e-04 -6.0097e-03 9.2575e-04 4.1607e-02 1.0547e-02 3.9181e-03 -8.3478e-03 -1.2349e-02 2.9609e-03 2.4361e-03 1.3051e-02 -9.3701e-03 -5.1040e-03 8.5949e-04 -6.4957e-03 1.5064e-02 -2.3446e-02 5.5451e-03 -4.1346e-02 -2.6646e-03 2.7194e-02 8.5292e-03 -2.8414e-03 -1.9089e-02 2.5987e-04 1.5280e-02 -6.9132e-03 -3.0431e-02 pseudo_product_test() The following relations MUST hold: A * A+ * A = A A+ * A * A+ = A+ ( A * A+ ) is MxM symmetric; ( A+ * A ) is NxN symmetric Here are the Frobenius norms of the errors in these relationships: A * A+ * A = A 0.000000 A+ * A * A+ = A+ 0.000000 ( A * A+ ) is MxM symmetric; 0.000000 ( A+ * A ) is NxN symmetric; 0.000000 In some cases, the matrix A * A+ may be interesting (if M <= N, then it MIGHT look like the identity.) A * A+: 1.0000e+00 -1.0734e-16 7.6328e-17 -4.9440e-16 -1.2490e-16 9.0206e-17 1.0000e+00 -2.0817e-16 -4.6838e-17 -1.8041e-16 7.4940e-16 2.6455e-17 1.0000e+00 -3.4694e-18 2.9837e-16 -3.8164e-16 6.8522e-17 -5.3083e-16 1.0000e+00 -1.7000e-16 -1.3878e-16 -6.0715e-18 2.7756e-16 -7.0083e-16 1.0000e+00 In some cases, the matrix A+ * A may be interesting (if N <= M, then it MIGHT look like the identity.) A+ * A: Columns 1 through 6: 7.1359e-01 -5.6405e-02 1.7624e-02 2.5957e-01 -5.3602e-03 -2.8281e-01 -5.6405e-02 3.0463e-01 2.6463e-01 2.1226e-02 -5.8822e-02 -1.8622e-02 1.7624e-02 2.6463e-01 8.9485e-01 -2.4463e-02 8.1117e-02 -1.5385e-02 2.5957e-01 2.1226e-02 -2.4463e-02 6.7323e-01 2.6864e-01 1.7343e-01 -5.3602e-03 -5.8822e-02 8.1117e-02 2.6864e-01 2.0892e-01 2.4726e-01 -2.8281e-01 -1.8622e-02 -1.5385e-02 1.7343e-01 2.4726e-01 6.3959e-01 5.1111e-02 -3.6613e-01 1.2782e-01 -1.1995e-01 1.3804e-01 1.7281e-02 -2.2558e-01 -1.7667e-02 -4.5988e-04 1.8689e-01 5.3290e-02 -2.4173e-01 Columns 7 and 8: 5.1111e-02 -2.2558e-01 -3.6613e-01 -1.7667e-02 1.2782e-01 -4.5988e-04 -1.1995e-01 1.8689e-01 1.3804e-01 5.3290e-02 1.7281e-02 -2.4173e-01 7.4779e-01 4.3100e-02 4.3100e-02 8.1740e-01 pseudo_linear_solve_test(): Given: b = A * x1 so that b is in the range of A, solve A * x = b using the pseudoinverse: x2 = A+ * b. Norm of x1 = 23.727621 Norm of x2 = 20.542065 Norm of residual = 0.000000 Given: b = A' * x1 so that b is in the range of A', solve A' * x = b using the pseudoinverse: x2 = A+' * b. Norm of x1 = 117.898261 Norm of x2 = 117.898261 Norm of residual = 0.000000 For M < N, most systems A'*x=b will not be exactly and uniquely solvable, except in the least squares sense. Here is an example: Given b is NOT in the range of A', solve A' * x1 = b using the pseudoinverse: x2 = A+' * b. Norm of x2 = 0.055899 Norm of residual = 0.899198 svd_test(): Normal end of execution. 08-Oct-2025 20:02:38