Improving the accuracy of solving ill-conditioned SLAEs using the Gaussian method

Most problems in computational mathematics ultimately come down to solving systems of linear equations. At the moment, there are a huge number of algorithms for solving such systems. They are divided into two large groups: iterative and direct. Direct methods allow you to get the exact values ​​of the unknowns if the calculations are carried out accurately. Next, we will consider the Gauss method. This method can be used to solve systems of algebraic equations with so-called general matrices. At the same time, it allows one to determine whether the system is consistent and, if it is consistent, is the solution unique? However, for all methods there is a problem associated with the difficulties of ill-conditioned systems. This is the most common way to solve SLAE, which is based on the idea of ​​successive elimination of unknowns (this method will be described in more detail below).

Algorithm

The algorithm of the method uses sequential elimination of unknowns. With the help of elementary transformations, the system of equations, in the case of the existence of a unique solution, is reduced to an equivalent system of a triangular form, from which all unknown unknowns are found sequentially, starting from the last (by number).

The algorithm for solving SLAE by the Gaussian method consists of two stages

At the first stage, which is usually called a direct move, by elementary transformations over the equations of the system, it is brought to a stepped or upper triangular form. The forward move algorithm with a partial choice of the leading element in the column assumes that at the beginning, among the elements of the first column of the matrix, a nonzero one is chosen, it is moved to the uppermost position by permuting the rows, and the first row obtained after the permutation is subtracted from the remaining rows, multiplying it by a value equal to the ratio of the first element of each of these rows to the first element of the first row, thus nullifying the column containing the first component of the vector of unknowns below it. After these transformations have been completed, the first row and the first column are mentally crossed out, and similar actions are continued until a system with an upper triangular matrix is ​​obtained.

At the second stage, which is called the reverse move, the components of the vector of unknowns are sequentially determined, starting with the last one. The number of operations required to implement the forward stroke in the Gaussian method can be calculated using the formula:

Q_1=2/3 n^3+n^2/2-n/6

In the process of implementing the reverse move, the algorithm performs Q2 operations:

The paper investigates the influence of the conditionality of the SLAE on the error of the resulting solution, including for systems with a Hilbert matrix.

Software implementation

The program that implements the Gaussian algorithm in Python is shown below.

For the convenience of entering the system matrix and unknown members, the split() method is used

This method allows you to enter the coefficients of the equation and the free term separated by a space. Entering the next equation requires pressing the ENTER key.
The full text of the program with input and output information is as follows:

Program testing

To test the program’s performance, a number of tests were carried out to solve various SLAEs.

First test

First test

Second test

Second test

Third test

Third test

Tests were also carried out with ill-conditioned systems. For them, we took the Hilbert matrix with dimensions n×n as a matrix of coefficients. The Hilbert matrix is ​​a square matrix A with entries

A_ij=1/(i+j-1);i,j=1,2,3,…,n

SLAE with this matrix is ​​a good example of an ill-conditioned system, since the condition number for order 5 y is equal to μ(A)= 4.8* 105
The condition number of a matrix shows how close the matrix is ​​to a matrix of incomplete rank (for square matrices, to degeneracy). Those. if A is almost degenerate, then we can expect that small changes in A And b will cause significant changes in x. The condition number can be found using the following formulas:

cond(A)≡‖A‖*‖A^{-1}‖

Where A-1 is the matrix inverse for a nonsingular matrix A, ‖A‖,‖A-1 are matrix norms A And A-1 dimensions n×n respectively. Each of the norms can be calculated by the formula:

If cond(A)≥103then the matrix can be considered ill-conditioned.

Test 4. SLAE with Hilbert matrix H. Let the vector b = (b1,b2,b3… ,bn)T is performed under the assumption of the already known vector x =(1, 2, 3, …, n)T according to the formula:

b_i=∑_{k=1}^nh_{ik} x_k

As a result, the system takes the form: Hx=b

5th order matrix

5th order matrix

10 order matrix

10 order matrix

20 order matrix

20 order matrix

As can be seen from the results, the error in the calculations is already noticeable at n=5. The critical difference between the approximate and exact solutions can be observed already at n = 10. This is due to the approximate representation of numbers in floating point form in a computer under conditions of ill-conditioned Hilbert matrix. In order to exclude perturbations of input and intermediate data caused by rounding errors, we use the Fractions library, which allows us to represent data, including calculation results, in the form of numerical algebraic fractions

5th Order Matrix and Fractions

5th Order Matrix and Fractions

10th Order Matrix and Fractions

10th Order Matrix and Fractions

20th Order Matrix and Fractions

20th Order Matrix and Fractions

As you can see, the solutions obtained as a result of calculations using the Fraction library coincide with their “exact” representation as a vector x = (1, 2, 3, …, n)T.

All numbers until the final result is obtained are represented as a pair of numbers (numerator and denominator). Naturally, the time required to calculate the solution increases significantly, since the program spends a large amount of time finding a common denominator when adding two fractions. Below is the dependence of the system solution time on the order of the SLAE. The blue color shows the time of work using the Fractions method, and the red color shows the time using the standard float data type in which all calculations take place. On the vertical axis – time in seconds, on the horizontal – the order of the Hilbert matrix.

Solution time from order

Solution time from order

Summing up

When solving SLAE, it is important to take into account the conditionality of the system matrix. In the case of a large value of the condition number, one should take into account the possibility of a significant and sometimes catastrophic effect of rounding errors and input data assignment errors on the result (solution). In this case, to eliminate the influence of rounding errors, you can use the representation of numbers in the form of numerical algebraic fractions. In this case, the time to solve the problem can increase for systems of order n, n>100 many dozens of times.

Similar Posts

Leave a Reply

Your email address will not be published. Required fields are marked *