Hilbert Matrices and Their Inverses
This example shows how to compute the inverse of a Hilbert matrix using Symbolic Math Toolbox™.
Definition : A Hilbert matrix is a square matrix with entries being the unit fraction. . For example, the 3x3
Hilbert matrix is
Symbolic computations give accurate results for these ill-conditioned matrices, while purely numerical methods fail.
Create a 20-by-20 numeric Hilbert matrix.
H = hilb(20);
Find the condition number of this matrix. Hilbert matrices are ill-conditioned, meaning that they have large condition numbers indicating that such matrices are nearly singular. Note that computing condition numbers is also prone to numeric errors.
cond(H)
ans = 5.1944e+19
Therefore, inverting Hilbert matrices is numerically unstable. When you compute a matrix inverse, H*inv(H)
must return an identity matrix or a matrix close to the identity matrix within some error margin.
First, compute the inverse of H
by using the inv
function. A warning is thrown due to the numerical instability.
H*inv(H)
Warning: Matrix is close to singular or badly scaled. Results may be inaccurate. RCOND = 1.309732e-19.
ans = 20×20
1.0000 -0.0000 -0.0004 -0.0133 -0.0078 0.1494 0.4224 -0.1053 -1.5628 1.5000 0.4620 0.9110 2.8195 -0.3235 -0.6250 -5.6448 -0.9671 -4.5672 3.7431 -0.3125
-0.0558 1.0000 0.0000 -0.0050 0.0078 0.1464 0.2495 0.3592 -1.6037 -0.1875 -0.0668 0.9002 2.9993 -1.2306 -1.6875 -6.1419 -1.9599 -2.3883 0.7089 0.3750
-0.2370 0.5524 1.0000 -0.0034 0.0195 0.2554 0.0272 0.5033 -1.1266 0.1875 -0.3678 0.9255 0.9066 0.7549 -0.8125 -6.0203 -1.3259 0.5149 -0.5031 0.1875
-0.3796 1.1272 -0.2467 0.9963 0.0273 -0.0371 0.4644 -0.5493 -0.1693 -0.9375 -0.2952 0.1800 2.9372 -1.2208 -2.0000 -1.5512 -1.3187 0.1597 0.2638 0
-0.4789 1.5782 -0.4757 -0.0590 0.9922 0.1828 -0.0083 0.3378 -0.1781 -0.5000 -0.1041 0.6158 2.4213 -1.8971 -1.2500 -1.5460 -0.2774 -2.0273 1.0890 -0.1875
-0.5455 1.8979 -0.5939 -0.2425 0.0664 1.0698 0.2080 0.4684 -0.9442 -0.9375 -0.5483 0.7483 1.0083 0.8748 1.0625 -4.4184 -0.7901 1.6298 -0.7021 0.1250
-0.5889 2.1108 -0.5976 -0.5488 0.2051 0.1245 1.0389 -0.3502 -0.7717 -0.3750 -0.2460 0.5758 1.9188 -0.3853 -0.3125 -3.8710 -1.0324 -0.5902 0.1486 0.0625
-0.6160 2.2432 -0.5093 -0.9648 0.4648 -0.1813 0.7114 -0.1442 0.2887 -1.0000 -0.2212 -1.5915 1.7356 -1.7679 -1.4375 -1.7909 -4.9411 5.2257 -2.2542 0.6250
-0.6317 2.3172 -0.3569 -1.4416 0.6758 0.0942 -0.2973 0.4413 -0.9635 1.0000 -0.9617 1.9564 -3.2818 2.5536 1.9062 -7.7002 2.1550 -1.9163 -0.7808 0.2500
-0.6393 2.3492 -0.1620 -1.9724 1.0527 -0.1114 0.2407 0.0095 -0.2193 1.5000 -0.7571 1.0671 -0.7490 -0.6288 0.1250 -4.0656 0.5029 -0.3419 0.5529 0.0938
⋮
Now, use the MATLAB® invhilb
function that offers better accuracy for Hilbert matrices. This function finds exact inverses of Hilbert matrices up to 15-by-15. For a 20-by-20 Hilbert matrix, invhilb
finds the approximation of the matrix inverse.
H*invhilb(20)
ans = 20×20
1010 ×
0.0000 -0.0000 0.0000 -0.0000 -0.0000 -0.0004 0.0013 -0.0037 0.0047 -0.1057 0.4019 0.2620 -0.4443 -7.5099 1.5435 3.2884 -1.1618 0.2301 0.0437 0.0029
-0.0000 0.0000 -0.0000 0.0000 -0.0000 -0.0001 -0.0009 0.0172 -0.0628 0.2441 -0.8975 2.7711 -2.8924 -1.4888 -4.3218 6.4897 -3.0581 0.7925 -0.0037 0.0029
0.0000 0.0000 0.0000 -0.0000 -0.0000 -0.0001 0.0009 0.0042 -0.0303 -0.1460 -0.1028 0.2862 -0.9267 -4.0597 -4.8050 4.7727 -1.3255 0.4667 -0.0211 -0.0033
0.0000 0.0000 -0.0000 -0.0000 -0.0000 -0.0001 0.0004 0.0002 -0.0056 -0.1653 -0.1136 0.3616 -1.6920 -3.7214 0.9328 4.4169 -1.1602 0.5541 -0.0024 -0.0050
-0.0000 0.0000 -0.0000 0.0000 -0.0000 -0.0000 -0.0006 0.0068 -0.0394 0.1946 -0.7818 1.9314 -2.1273 -1.0745 1.2885 4.6349 -1.8923 0.3793 -0.0435 0.0045
0.0000 -0.0000 0.0000 -0.0000 0.0000 -0.0002 0.0014 -0.0028 0.0304 0.0562 -0.0724 -0.2073 -0.3769 -5.0729 4.9325 2.3488 -0.5046 0.2240 0.0675 -0.0138
-0.0000 0.0000 -0.0000 0.0000 -0.0000 0.0000 -0.0013 0.0101 -0.0520 0.0721 -0.7715 2.9297 -3.3374 1.4084 0.5570 7.1053 -2.7578 0.8815 -0.1648 -0.0038
0.0000 -0.0000 0.0000 -0.0000 0.0000 -0.0002 0.0018 -0.0040 0.0231 -0.2278 0.2000 0.2766 -0.8262 -4.7229 3.8319 2.1857 -0.4749 -0.1747 0.1257 -0.0141
-0.0000 0.0000 -0.0000 0.0000 -0.0000 0.0000 0.0004 0.0144 -0.0513 0.1850 -0.7063 1.9776 -2.0602 -1.0763 -0.0201 3.8539 -2.4580 0.5675 -0.0189 0.0036
0.0000 0.0000 0.0000 -0.0000 -0.0000 -0.0002 0.0008 -0.0019 0.0027 -0.0679 -0.0497 0.8865 -0.6720 -1.2729 0.8925 2.3625 -1.2616 0.2443 0.0184 -0.0091
⋮
To avoid round-off errors, use exact symbolic computations. For this, create the symbolic Hilbert matrix.
Hsym = sym(H)
Hsym =
Get the value of the condition number. It has been derived by symbolic methods and is free of numerical errors.
vpa(cond(Hsym))
ans =
Although its condition number is large, you can compute the exact inverse of the matrix.
Hsym*inv(Hsym)
ans =