On January 7 Rodger Rosenbaum initiated the thread "Ill conditioning and calculator precision" to carry forward discussion of an earlier thread initiated by Valentin Albillo. Near the end of the thread Rodger wrote "But, back to the original purpose I had when I started this thread. I noticed that a slightly different perturbation of the column matrix has even greater effect. Change the third element of the column matrix in Valentin's original problem from 45 to 45.1--now the exact solution is:
x1 = -2128901625This is a tremendous change in the solution for an apparently insignificant change in one element."
x2 = 268387
x3 = -1898169427
x4 = -1909014362
x5 = -1317226
x6 = 1908985003
x7 = 3994038147
For reference Valentin's original problem was
1.3 x1 + 7.2 x2 + 5.7 x3 + 9.4 x4 + 9.0 x5 + 9.2 x6 + 3.5 x7 = 45.3where the system has the exact solution:
4.0 x1 + 9.3 x2 + 9.0 x3 + 9.9 x4 + 0.1 x5 + 9.5 x6 + 6.6 x7 = 48.4
4.8 x1 + 9.1 x2 + 7.1 x3 + 4.8 x4 + 9.3 x5 + 3.2 x6 + 6.7 x7 = 45.0
0.7 x1 + 9.3 x2 + 2.9 x3 + 0.2 x4 + 2.4 x5 + 2.4 x6 + 0.7 x7 = 18.6
4.1 x1 + 8.4 x2 + 4.4 x3 + 4.0 x4 + 8.2 x5 + 2.7 x6 + 4.9 x7 = 36.7
0.3 x1 + 7.2 x2 + 0.6 x3 + 3.3 x4 + 9.7 x5 + 3.4 x6 + 0.4 x7 = 24.9
4.3 x1 + 8.2 x2 + 6.6 x3 + 4.3 x4 + 8.3 x5 + 2.9 x6 + 6.1 x7 = 40.7
x1 = x2 = x3 = x4 = x5 = x6 = x7 = 1.0
In the earlier thread I provided solutions to the original problem on several different machines and to a modied problem which changed the B1 value from 45.3 to 45.4 . Solutions to the modified problem with B3 changed from 45.0 to 45.1 follow.
With my HP-28S:
B/A INV A * Bwhere the signs of the elements of the solution are correct but the magnitudes are off by an order of magnitude yielding relative errors of 14.2 . The calculations of A*X do indicate that there are difficulties with the solution.X A*X X A*X
-30272808391.1 45.04 -30272808407.9 -87.494
3816428.18496 48.302 3816428.21599 -181.961
-26991815255.9 46.374 -26991815236.6 -162.843
-27146029351.6 18.506 -27146029342.2 -7.895
-18730860.0936 36.536 -18780860.0977 -144.096
27145611882.5 24.694 27145611823 -18.702
56794898378.3 40.67 56794898341.2 -146.409
My Model 100 and TI-85 yield the following solutions
Model 100 TI-85where the signs and magnitudes are correct. At least the first three digits of each element of the Model 100 answer are correct yielding relative errors of 3.57e-04. The relative errors for the TI-85 are 9.15e-03.X A*X X A*X
-2128141186.5423 45.3 -2148392852.67 45.3014
268291.13337654 48.399 270844.2167 49.3991
-1897491405.5188 45.1 -1915548178.57 45.1
-1908332466.734 18.5999 -1926492404.72 18.5999
-1316755.4898248 36.698 -1329285.91435 36.6983
1908303118.2216 24.9 1926462776.91 24.90025
3992611486.3345 40.699 4030605692.38 40.7002
The MATRIX routine of the HP-41 Math Pac and ML-02 routine of the TI-59 Master Library yield the following results
HP-41 TI-59where the magnitudes and even most of the signs of the result from the HP-41 are wrong. The magnitudes and signs of the result from the TI-59 are correct and of the correct order of magnitude yielding relative errors of 2.4e-2. Curiously, with the HP-41 the product of Matrix A and the solution vector yield the exact values of the elements of vector B. This suggests that comparing the product A*X with B to measure the acceptability of the solution may not be a good idea if the user is interested in a real world solution as opposed to a solution which is of value only in the calculator's world. The determinants for Valentin's Matrix A wereX A*X X A*X
24666733.49 45.3 -2077198547.178 45.32
3108.726 48.4 261868.9227801 48.44
21993331.58 45.1 -1852069973.364 45.11
22118987.39 18.59 -1862651525.312 18.6
15263.2 36.7 -1285235.478686 36.7
-22118645.18 24.84 1862622879.371 24.901
-46277322.87 40.7 3897037862.598 40.71
Exact 1.e-07Rodger also wrote "Try something I described in an earlier post; add .1 to the {3,7} element of the A matrix, and subtract .1 from the {7,7} element. This tiny change reduces the condition number of the A matrix from about 3E12 to about 400, and using this slightly perturbed A matrix with the perturbed column matrix mentioned just above gives a HP48G solution:
HP-41 -86.31e-7
HP-28S 0.0703238892937e-07
TI-59 1.024890772918e-07
TI-85 0.990927530944e-07
x1 = 1.20999247798This isn't as close to [1 1 1 1 1 1 1] as we would like, but given what the solution was with just the third element in the column matrix perturbed and using the original A matrix, this inprovement is an astounding result considering that only two tiny perturbations were made to the A matrix.
x2 = 1.01325966115
x3 = 1.26616078717
x4 = 1.29157926088
x5 = .997263593179
x6 = .690061410742
x7 = .499897531935
My purpose here is to show that with proper techniques we can find reasonable solutions to very ill-conditioned systems. Systems as badly conditioned as Valentin's example are unlikely to occur in practice; they must be artificially constructed. But, if we can handle them, then we can certainly handle the ones we encounter in real life problems.
Palmer, it would be good if you tried solving the system on your various calculators, consisting of the column matrix perturbed in the third element, and the {3,7} and {7,7} elements of the A matrix perturbed as I described. I think you will find that you don't get wildly varying results with these changes. For example, on my HP48S, which uses the same matrix routines as the HP28, I get about 9 or 10 digit agreement with the solution posted above, whether I use the divide key or invert the A matrix and multiply times the column matrix."
My HP-28S yielded
B/A INV A *BMy HP-28S solutions are consistent with Rodger's results from his HP-48S. Note that when using the B/A method the product of matrix A and the solution vector yields exactly the B vector.X A*X X A*X
1.20999247808 45.3 1.20999247808 45.2999999988
1.01325466114 48.4 1.01325966123 48.3998999983
1.26616078724 45.1 1.26616078741 45.0999999989
1.29157926094 18.6 1.29157926101 18.6000000005
0.99726359318 36.7 0.997263593181 36.6999999992
0.690061410687 24.9 0.690061410466 24.9000000001
0.499897531784 40.7 0.499897531481 40.6999999991
The simultaneous equation solution with my TI-85 is
1.2099924779883where the product of the A matrix and the solution matrix yields exactly the B matrix.
1.0132596611477
1.2661607871723
1.2915792608817
0.99726359317885
0.69006141074113
0.49989753192964
The solution from my Model 100 is
X B - A*Xwhere in the second column I listed B - A*X rather than A*X to avoid typing in all those zeroes. I wasn't pleased with the TI-85 and Model 100 results. We were hoping to get a solution closer to all ones using the modification to Valentin's problem, but the TI-85 and Model 100 solutions were further away, not closer than with Valentin's original problem. It occurred to me that perhaps the exact solution to the perturbed problem wasn't all ones. I went to my HP-49G emulator and solved the problem in the exact mode with the following results1.2099924779846 -1e-12
1.0132596611482 2e-12
1.2661607871704 1e-12
1.291579260881 1e-12
0.99726359317867 0
0.69006141074154 0
0.49989753193417 1e-12
x1 = 9667493448/7989713675 = 1.2099924779845which confirmed my suspicions.
x2 = 8095654571/7989713675 = 1.013259661148
x3 = 10116362156/7989713675 = 1.2661607871699
x4 = 10319348483/7989713675 = 1.2915792608801
x5 = 7967850568/7989713675 = 0.99726359317876
x6 = 1102678618/1597942735 = 0.6900614107426
x7 = 3994038147/7989713675 = 0.49989753193502
The results obtained with the MATRIX routine in the HP-41 Math Pac module and the ML-02 routine in the TI-59 Master Library module:
HP-41 TI-59where the differences from the exact answer are only seen in the lower order digits of each machine. The determinants for the modified matrix A are1.209992482 1.209992478108
1.013259663 1.013259661129
1.266160800 1.266160787202
1.291579278 1.29157926087
0.997263593 0.9972635931871
0.690061393 0.6900614107704
0.499897511 0.4998975318243
Exact 798.9713675Some playing with the numbers on my HP-28S:
HP-41 798.9713752
HP-28S 798.971367542
TI-59 798.9713674893
TI-85 798.9713675
If I multiply Valentin's original A matrix times the X vector from from the third modification with A(3,7), A(7,7) and B3 changed I get
45.3where the RMS relative error with respect to the B vector of the third modification is 2.672e-02. Furthermore, If I multiply the original A matrix by a vector obtained by rounding the X vector to the nearest 0.25, i.e., a vector made up of the elements 1.25, 1, 1.25, 1.25, 1, 0.75, 0.5 I get
48.4
45.0500102468
18.6
36.7
24.9
40.7499897532
45.35which has n RMS relative error with respect to the B vector of the third modification of 2.988e-02
48.45
45.025
18.6
36.7
24.9
40.725
As a novice to the idea of perturbing problems to relieve ill-conditioning I am not sure what all of this means.
Finally, Rodger suggested that someone should translate a Fortran program for singular value decomposition (SVD) for use with BASIC. I took a quick look. Someone else will have to do it. I haven't done any Fortran programming since about 1970.