yakubin’s notes

Floating-point Gaussian elimination errors for Hilbert matrices

I’ve written a program which calculates the solution to a system of n affine equations of n variables using arbitrary-precision fractions, as well as IEEE754 64-bit binary floating-point numbers with the Gaussian elimination algorithm with partial and complete pivoting, and calculates the numerical error of the latter 2 methods based on the exact solution. The input matrices I’ve used are known as Hilbert matrices, i.e. H[i][j] = 1/(i + j - 1). The results show that without a thorough analysis of the numerical conditioning of a problem, even using a numerically-stable method (like Gaussian elimination with complete pivoting), floating-point numbers lead to catastrophic errors.

Results: hilbert.csv. R program used to generate the plots: hilbert-error-plots.R.

Gaussian elimination relative error for Hilbert matrices Gaussian elimination absolute error for Hilbert matrices