In [1]:
# Code attribution: Yiyin Shen, Tyler Caraza-Harter
# Imports
import scipy
import numpy
In [2]:
# Create a Hilbert matrix
h = scipy.linalg.hilbert(5)
h
Out[2]:
array([[1. , 0.5 , 0.33333333, 0.25 , 0.2 ],
[0.5 , 0.33333333, 0.25 , 0.2 , 0.16666667],
[0.33333333, 0.25 , 0.2 , 0.16666667, 0.14285714],
[0.25 , 0.2 , 0.16666667, 0.14285714, 0.125 ],
[0.2 , 0.16666667, 0.14285714, 0.125 , 0.11111111]])
In [3]:
# Solve h x = h @ [1, 1, 1, 1, 1] for x -> x should be [1, 1, 1, 1, 1]
x = scipy.linalg.solve(h, h @ ([1] * 5))
x
Out[3]:
array([1., 1., 1., 1., 1.])
In [4]:
# Creat a larger Hilbert matrix
h = scipy.linalg.hilbert(50)
h
Out[4]:
array([[1. , 0.5 , 0.33333333, ..., 0.02083333, 0.02040816,
0.02 ],
[0.5 , 0.33333333, 0.25 , ..., 0.02040816, 0.02 ,
0.01960784],
[0.33333333, 0.25 , 0.2 , ..., 0.02 , 0.01960784,
0.01923077],
...,
[0.02083333, 0.02040816, 0.02 , ..., 0.01052632, 0.01041667,
0.01030928],
[0.02040816, 0.02 , 0.01960784, ..., 0.01041667, 0.01030928,
0.01020408],
[0.02 , 0.01960784, 0.01923077, ..., 0.01030928, 0.01020408,
0.01010101]])
In [5]:
# Solve h x = h @ [1, 1, 1, 1, 1] for x -> x should be [1, 1, 1, 1, 1]
x = scipy.linalg.solve(h, h @ ([1] * 50))
x
C:\Users\young\AppData\Local\Temp\ipykernel_28360\2369993951.py:2: LinAlgWarning: Ill-conditioned matrix (rcond=3.66446e-20): result may not be accurate. x = scipy.linalg.solve(h, h @ ([1] * 50))
Out[5]:
array([ 1.00000007, 0.99999158, 1.00016222, 1.0011918 ,
0.94023228, 1.67868573, -2.87605487, 13.57667066,
-21.99608675, 21.1753999 , -1.86476998, 1.57826089,
-14.76788412, 7.73021038, 14.18600245, -1.0017191 ,
-10.69517354, 8.14104972, 5.30813439, 2.11098027,
-43.2867539 , 24.86775979, 52.62730875, -12.4174761 ,
-39.58786349, -10.76829926, 49.22531073, -16.91015691,
-20.70245705, 27.2848732 , -14.92041968, 28.31102569,
-35.54797011, 25.68298935, -13.49669524, 17.23147797,
-0.80052785, -6.66775122, -4.00359584, 20.38809736,
-23.42301574, 0.71219121, 26.89047453, -45.89097618,
52.06348848, -1.64428045, -32.81983181, 32.06026043,
-20.13780275, 8.45533223])
In [6]:
# Compute the coondition number
numpy.linalg.cond(h)
Out[6]:
1.0993264246156683e+19
In [7]:
# Creat a even larger Hilbert matrix
h = scipy.linalg.hilbert(500)
h
Out[7]:
array([[1. , 0.5 , 0.33333333, ..., 0.00200803, 0.00200401,
0.002 ],
[0.5 , 0.33333333, 0.25 , ..., 0.00200401, 0.002 ,
0.00199601],
[0.33333333, 0.25 , 0.2 , ..., 0.002 , 0.00199601,
0.00199203],
...,
[0.00200803, 0.00200401, 0.002 , ..., 0.00100503, 0.00100402,
0.00100301],
[0.00200401, 0.002 , 0.00199601, ..., 0.00100402, 0.00100301,
0.001002 ],
[0.002 , 0.00199601, 0.00199203, ..., 0.00100301, 0.001002 ,
0.001001 ]])
In [8]:
# Solve h x = h @ [1, 1, 1, 1, 1] for x -> x should be [1, 1, 1, 1, 1]
x = scipy.linalg.solve(h, h @ ([1] * 500))
[min(x), max(x)]
C:\Users\young\AppData\Local\Temp\ipykernel_28360\4270393687.py:2: LinAlgWarning: Ill-conditioned matrix (rcond=7.68731e-22): result may not be accurate. x = scipy.linalg.solve(h, h @ ([1] * 500))
Out[8]:
[-346.8311692965331, 325.42502787608424]
In [9]:
# Compute the coondition number
numpy.linalg.cond(h)
Out[9]:
1.8115962657688702e+20