Simple and practical! The linear iteration method created by three Germans surpassed an era

December 26th this year marks the 200th anniversary of the invention of the first linear iteration method by German mathematician Gauss. In the article "How to Iterate to Solve Linear Equations?" published not long ago, we started from the basic concepts and clarified the preparation knowledge and specific process of iterative methods for solving linear equations. This time we will explore the basic ideas and convergence of the two most famous and simplest first-order constant iteration methods in history: the Jacobi iteration method and the Gauss-Seidel iteration method.

Written by Ding Jiu (Professor of Mathematics at the University of Southern Mississippi)

In my previous article "How to Iterate to Solve Linear Equations?", in order to link it with the contraction mapping theorem applicable to general complete distance spaces, I only gave a simple convergence proof for iterative solutions to linear equations.

All roads lead to Rome

Fortunately, we have other norms available. As a direct generalization of the usual notion of line segment length, the n-dimensional Euclidean space

Eigenvalues ​​and spectral radius

However, there is a necessary and sufficient condition for convergence waiting for us to explore. This condition no longer uses the norm of the matrix, but the spectrum of the matrix. The spectrum gives an intrinsic property of the matrix that does not depend on the norm of the matrix. The spectrum of a square matrix is ​​a finite set of complex numbers consisting of all its eigenvalues. Given a square matrix M of order n, we call the complex number λ an eigenvalue of M if the complex matrix M - λI is not invertible, that is, there exists a complex non-zero column vector x such that Mx = λx. This x is called an eigenvector of M corresponding to the eigenvalue λ. We have mentioned that a matrix is ​​singular if and only if its determinant is equal to zero, so λ is an eigenvalue of M if and only if det(M - λI) = 0, where the symbol det represents the determinant. If λ on the left side of this equation is regarded as a variable, according to the definition of the determinant, the expansion result of det(M - λI) is an n-order polynomial in λ, so an n-order square matrix M has at most n different eigenvalues. We call the maximum absolute value of all eigenvalues ​​of M the spectral radius of M, denoted by ρ(M).

The main purpose of studying linear iteration is to numerically solve the linear equation system Ax = b, where the coefficient matrix A is non-singular.

Conversely, if ρ(M) < 1, then there exists a sufficiently small positive number ε such that ρ(M) + ε < 1. This seemingly simple fact is a fundamental property of real numbers. On the other hand, for every eigenvalue λ of M and any norm ||•|| on Rn, let x be an eigenvector corresponding to λ. From the inequalities ||M|| ||x|| ≥ ||Mx|| = ||λx|| = |λ| ||x|| and x ≠ 0, we can see that both sides of the inequality just obtained can be divided by the positive number ||x||, and we get the inequality ||M|| ≥ |λ|, which shows that: the spectral radius ρ(M) of the matrix M is a lower bound of the matrix norm induced by all vector norms in M ​​||M||. In fact, this lower bound is the "lowest bound" of the set of real numbers corresponding to all norm values ​​||M|| of the fixed matrix M, that is, the "greatest lower bound" among all lower bounds of this set. The concept of "upper and lower bounds" is the most basic and important concept in calculus. Once you really master it, the limit theory and the concepts that follow, such as continuity, derivatives, integrals, series, etc., can even be learned easily.

Jacobi Iteration Method and Gauss-Seidel Iteration Method

Now we can focus on the Jacobi and Gauss-Seidel methods for solving the linear system Ax = b, where A is an n-order reversible matrix. Both of these classic German iterative methods are based on the special selection of N and P in the matrix decomposition A = N - P, and their iterative formats can be written as

Let's compare the significant differences between the Jacobi method and the Gauss-Seidel method in the process of calculating iterative point components. From the iterative formula of n components in the Jacobi iteration method, it can be seen that after all components of the k-1th iteration vector are calculated, the components of the k-th iteration vector are calculated one by one; in other words, all components obtained in the k-1th iteration are used to calculate all components of the k-th iteration. For the Gauss-Seidel iteration method, whenever the i-th component of the k-th iteration vector needs to be calculated, not only the components from the i+1th to the nth components of the k-1th iteration vector that have been calculated are needed to help, but also the 1st to i-1th components of the k-th iteration vector that have been obtained in the current iteration are needed to help; in other words, the Gauss-Seidel method is "more eager for success" than the Jacobi method, and orders the iteration point component that has just been captured in the current iteration step to replace the same subscript component that was "captured" in the previous iteration step, and "go into battle". This shows that Gauss and Seidl were better at "using troops" than Jacobi and were unwilling to waste any morale.

Convergence conditions

Naturally, neither the Jacobi iteration method nor the Gauss-Seidel iteration method can guarantee convergence without certain additional assumptions about the matrix A. So, when the prerequisite that "all diagonal elements of A cannot be zero" is met, what is the sufficient condition to ensure their convergence?

Here, the strict diagonal dominance requirement to ensure the convergence of the Jacobi iteration can be changed to irreducible weak diagonal dominance, that is, all the strict inequalities above are replaced by general inequalities "≥", and at least one inequality is strictly true, but with the additional condition that the matrix is ​​irreducible. If the rows and columns of a square matrix can be rearranged in the same way (no more rearrangement is also considered "rearrangement" because the identity matrix is ​​also a permutation matrix), the result becomes a 2×2 block upper triangular matrix, that is, the lower left corner of the matrix has at least a 1×1 zero matrix, then we say that the matrix is ​​reducible, otherwise it is irreducible.

When A is a strictly diagonally dominant matrix, the ∞-norm of the Gauss-Seidel iteration matrix is ​​no greater than the ∞-norm of the Jacobi iteration matrix. The smaller this norm is, the faster the convergence speed is usually. Therefore, in terms of convergence speed, the former method is no less than the latter method.

As the final mathematical dish of this article, we present another unique convergence proposition of the Gauss-Seidel method:

If the coefficient matrix A of the linear system Ax = b is symmetric positive definite, then the Gauss-Seidel iteration method converges.

end

The torch of linear iteration method ignited by German mathematicians two hundred years ago has been passed across the ocean a hundred years later to the American continent where science and technology are booming and computers are born. In order to optimize the iteration efficiency, starting from the classic Gauss-Seidel method, in 1950, American mathematician and computer scientist David M. Young Jr. (1923-2008, "Fanpu" will introduce him separately) and American physicist and computer scientist Stanley Phillips Frankel (1919-1978) almost simultaneously introduced a relaxation factor ω for a certain affine combination, leading to the "method of successive over-relaxation (SOR)". The former proposed the SOR method in his doctoral dissertation Iterative Methods for Solving Partial Difference Equations of Elliptic Type at Harvard University; the latter published a paper Convergence rates of iterative treatments of partial differential equations in the Journal of the American Mathematical Society (renamed Mathematics of Computation ten years later), in which he conducted a comprehensive analysis of the SOR method he designed and called the "Extrapolated Liebmann method". These two pioneering studies that proposed the SOR method are both related to partial differential equations that appear in large numbers in science and engineering. When using electronic computers to solve large sparse linear equations obtained by discretizing these continuous equations, iterative methods are the first choice. This is the era in which the SOR method came into being.

When the relaxation factor ω is 1, the SOR method is reduced to the Gauss-Seidel method. We will not elaborate on this method and end this article here, because we have already known the basic idea of ​​the convergence criterion of linear iterative methods by understanding the two classical iterative methods. In short, the inventions of three German mathematicians in the 19th century are the source of the modern wave of research on iterative solutions of large, sparse, and structured linear equations.

Written on Sunday, November 19, 2023, Summer Hills, Hattiesburg, USA

This article is supported by the Science Popularization China Starry Sky Project

Produced by: China Association for Science and Technology Department of Science Popularization

Producer: China Science and Technology Press Co., Ltd., Beijing Zhongke Xinghe Culture Media Co., Ltd.

Special Tips

1. Go to the "Featured Column" at the bottom of the menu of the "Fanpu" WeChat public account to read a series of popular science articles on different topics.

2. Fanpu provides a function to search articles by month. Follow the official account and reply with the four-digit year + month, such as "1903", to get the article index for March 2019, and so on.

Copyright statement: Personal forwarding is welcome. Any form of media or organization is not allowed to reprint or excerpt without authorization. For reprint authorization, please contact the backstage of the "Fanpu" WeChat public account.