=Paper=
{{Paper
|id=Vol-2837/paper3
|storemode=property
|title=Algebraic multigrid method with skew-Hermitian smoothers
|pdfUrl=https://ceur-ws.org/Vol-2837/paper3.pdf
|volume=Vol-2837
|authors=Tatiana S. Martynova,Galina V. Muratova,Zeng-Qi Wang
}}
==Algebraic multigrid method with skew-Hermitian smoothers==
https://ceur-ws.org/Vol-2837/paper3.pdf
Algebraic multigrid method with skew-Hermitian
smoothers
Tatiana S. Martynovaa , Galina V. Muratovaa and Zeng-Qi Wangb
a
Vorovich Institute of Mathematics, Mechanics and Computer Science, Southern Federal University, 200/1 Stachki Ave., Bld. 2,
Rostov-on-Don, 344090, Russian Federation
b
School of Mathematical Sciences, Shanghai Jiao Tong University, Shanghai 200240, P.R. China
Abstract
An algebraic multigrid method with a special type of smoothers is used for solving large sparse systems of linear
algebraic equations with a strongly non-Hermitian matrix. Namely, skew-Hermitian triangular and product-
triangular iterative methods have been used as the smoothers. When implementing the algebraic multigrid
method, the Parallel Modified Independent Set (PMIS) algorithm for coarsening grid has been used. The two-
dimensional unsteady Navier-Stokes problem for viscous incompressible fluid, written in primitive variables
"velocity-pressure", is considered as a model one. Numerical experiments have shown high efficiency of the
algebraic multigrid method with the skew-Hermitian smoothers for this problem.
Keywords
algebraic multigrid method, skew-Hermitian smoothers, unsteady Navier-Stokes problem
1. Introduction
Multigrid methods (MGMs) are called scalable or optimal because they can solve a linear system
with π unknowns with only π(π ) work. Because this work can be effectively distributed across a
parallel machine, MGMs can solve ever larger problems on proportionally larger parallel computers
in essentially constant time, making them an ideal solver for large-scale scientific simulation [1].
Nowadays MGMs are the most efficient tools for solving a large class of problems including bound-
ary value problems for partial differential equations (PDEs) or systems of PDEs. The MGM efficiency
depends on the adjustment of its components to the problem to be solved. Important parts of the MGM
are a smoothing procedure and the coarse-grid correction. We suggest to use special iterative meth-
ods as smoothers for the MGM and non-standard course-grid correction to a good approximation of
the smooth error components. A successful use of MGMs for solving hard practical problems depends
also on their robustness against bad or small parameters (e.g. geometrical parameters, singularities,
anisotropies etc.) arising in many applications.
It is known that some iteration methods have a smoothing effect on the error of approximation
[2, 3]. This property is fundamental for the multigrid idea and is connected with fast damping high-
frequency Fourier components of an initial error in decomposition on the basis from eigenvectors.
The smoothing method is the central component of the multigrid algorithm; it is the most dependent
part of the MGM on the type of the problem being solved. The role of smoothing methods is that
Far Eastern Workshop on Computational Technologies and Intelligent Systems, March 2β3, 2021, Khabarovsk, Russia
" martynova@sfedu.ru (T.S. Martynova); muratova@sfedu.ru (G.V. Muratova); wangzengqi@sjtu.edu.cn (Z. Wang)
~ https://sfedu.ru/en/person/martynova (T.S. Martynova); https://sfedu.ru/en/person/muratova (G.V. Muratova);
http://math.sjtu.edu.cn/faculty/zqwang1/zqwang.html (Z. Wang)
0000-0002-7073-1162 (T.S. Martynova); 0000-0002-9816-6759 (G.V. Muratova); 0000-0001-6087-4023 (Z. Wang)
Β© 2021 Copyright for this paper by its authors.
Use permitted under Creative Commons License Attribution 4.0 International (CC BY 4.0).
CEUR
Workshop
Proceedings
http://ceur-ws.org
ISSN 1613-0073 CEUR Workshop Proceedings (CEUR-WS.org)
they should not so much reduce the total error as smooth it (namely, suppress the high-frequency
harmonics of the error) so that the error can be well approximated on a coarse grid [2, 3].
There exist two approaches in the MGM: geometric multigrid and algebraic multigrid methods
(AMG) [1, 3]. AMG is one of the most important issues from a practical point of view. In contrast
to the geometrical multigrid method, its algebraic version recovers algebraically the ingredients for a
complete multigrid algorithm from the fine grid data only. AMG determines coarse grids, inter-grid
transfer operators, and coarse-grid equations based solely on the matrix entries.
There are two coarsening approaches in the AMG: RS (Ruge-Stuben) and PMIS (parallel changes
independent set) algorithms. The RS algorithm [4] is a traditional coarsening approach, it is based on
two heuristic criteria that achieve optimal convergence and minimal computational cost. The PMIS
is based on the same principles as the RS algorithm except that a heuristic criterion is not strictly
observed [5]. An advantage of the PMIS is the possibility of natural parallelization, which makes it
applicable for 3-D problems, as well as for problems with a huge number of unknowns.
We are interested in numerical solution of the SLAEs that arise from finite difference approximation
of the 2-D incompressible unsteady Navier-Stokes equations governing the flow of viscous Newtonian
fluids. We use two iterative methods based on the Hermitian and skew-Hermitian splitting of the ma-
trix π΄ of the SLAE as smoothers for the MGM. These are the skew-Hermitian triangular splitting (STS)
[6, 7] and the product-type skew-Hermitian triangular splitting (PSTS) [8] iteration methods. Based
on the smoothing and approximation properties convergence of the MGM with the PSTS smoothers
has been proved. Numerical experiments were carried out using the AMG (PMIS algorithm) with
STS-, PSTS- and Gauss-Seidel based smoothers.
2. Smoothers based on the STS and the PSTS iteration methods
Thus, we consider iterative solution of the large sparse SLAE
π΄π’ = π, π’, π β βπ , (1)
where π΄ β βπΓπ is a non-Hermitian and positive definite matrix.
Naturally, the matrix π΄ can be split into its Hermitian and skew-Hermitian parts as
π΄ = π΄0 + π΄1 , (2)
where
1 1
π΄0 = (π΄ + π΄β ), π΄1 = (π΄ β π΄β ), (3)
2 2
and π΄β denotes the conjugate transpose of the matrix π΄. Positive definiteness of the matrix π΄ means
that for all π₯ β βπ ⧡ {0}, π₯ β π΄0 π₯ > 0. Here π₯ β denotes the conjugate transpose of the complex vector π₯.
When π΄ is positive definite then its Hermitian part π΄0 is Hermitian positive definite, and the diagonal
of its skew-Hermitian part π΄1 is purely imaginary.
Let in some matrix norm ||| β
|||, |||π΄0 ||| << |||π΄1 |||, then the matrix π΄ is called strongly non-Hermitian
one. This situation occurs in many real applications, such as the discretization of the Navier-Stokes
equations. Notice, that ππππ(π΄1 ) = 0 when the coefficient matrix π΄ of the SLAE (1) is real.
Then, we can split the skew-Hermitian part π΄1 of the matrix π΄ β βπΓπ into
π΄1 = πΎπΏ + πΎπ , (4)
where πΎπΏ and πΎπ are the strictly lower and the strictly upper triangular parts of π΄1 , respectively.
Obviously, that πΎπΏ = βπΎπβ .
Based on the splitting (2)-(4) in [6, 7] classes of the STS iteration methods for solving SLAE (1) have
been proposed. The triangular operator of the STS uses only skew-Hermitian part of the coefficient
matrix π΄. These methods have been further developed in [9]. Moreover, based on the methodology
given in [9], the PSTS iteration method has been investigated in [8].
A particular problem when using the MGM is the choice of smoothers. There are a number of
iteration methods that can be used as smoothers, but not all of them are effective for solving strongly
non-Hermitian SLAEs. The behavior of the STS and the PSTS iteration methods is similar to the
behavior of the Gauss-Seidel method, which quickly damps the high-frequency harmonics of the
error, slowing down in the future. We give the formulas of these iteration methods.
The STS iteration method [6, 7]: Given an initial guess π’ (0) , and two positive parameters π and
π . For π = 0, 1, 2, ... until {π’ (π) } convergence, compute
π’ (π+1) = πΊ(π, π )π’ (π) + π π΅(π)β1 π,
where
πΊ(π, π ) = π΅(π)β1 (π΅(π) β π π΄),
π and π are two acceleration parameters, and π΅(π) is defined by
π΅(π) = π΅π + π((1 + π)πΎπΏ + (1 β π)πΎπ ), π = Β±1.
with π΅π β βπΓπ a prescribed Hermitian matrix.
The PSTS iteration method [8]: Given an initial guess π’ (0) , and two positive parameters π and
π . For π = 0, 1, 2, ... until {π’ (π) } convergence, compute
π’ (π+1) = πΊ(π, π )π’ (π) + π π΅(π)β1 π, (5)
where
πΊ(π, π ) = π΅(π)β1 (π΅(π) β π π΄),
π and π are two acceleration parameters, and π΅(π) is defined by
π Μ β1 π
π΅(π) = (π΅π + πΎπΏ )π΅π (π΅π + πΎΜπ ), (6)
2 2
where πΎΜπΏ = πΎπΏ + πΉ0 , πΎΜπ = πΎπ β πΉ0 , πΉ0 β βπΓπ is a Hermitian matrix, π΅π β βπΓπ is a prescribed Hermitian
positive definite matrix. Obviously that πΎΜπΏ = βπΎΜπβ , π΄1 = (πΎπΏ + πΉ0 ) + (πΎπ β πΉ0 ) = πΎΜπΏ + πΎΜπ .
It was shown in [10] that the skew-Hermitian iteration methods have smoothing effect on the error
of approximation. Therefore, these methods can be used as smoothers in the MGMs. Let us give a
proof of the MGM convergence with the PSTS smoothers.
3. MGM convergence with the PSTS smoothers
Let π΅0 (π) and π΅1 (π) be, respectively, the Hermitian and the skew-Hermitian parts of the matrix π΅(π)
defined by (6). That is to say,
π΅(π) = π΅0 (π) + π΅1 (π), (7)
where {
π΅0 (π) = 21 (π΅(π) + π΅(π)π ) = π΅π + ( π2 )2 πΎΜπΏ π΅πβ1 πΎΜπ ,
π΅1 (π) = 12 (π΅(π) β π΅(π)π ) = π2 (πΎΜπΏ + πΎΜπ ) = π2 π΄1 .
Suppose that
0 < πΌβ πΌ β€ π΄0 β€ π½β πΌ , 0 < πΌπ πΌ β€ π΅π β€ π½π πΌ , (8)
πΌπ πΌ β€ πΎΜπΏ π΅πβ1 πΎΜπ β€ π½π πΌ , πΌπ πΌ β€ π΅0 (π) β€ π½π πΌ , (9)
where πΌ is an identity matrix. Since π΄0 and π΅π are Hermitian positive definite matrices, π΅0 (π) and
πΎΜπΏ π΅πβ1 πΎΜπ are Hermitian matrices, and π΄1 is a skew-Hermitian matrix, the bounds πΌπ and π½π , π =
β, π, π, π , can be easily expressed with respect to the smallest and the largest eigenvalues or singular
values of the corresponding matrix, respectively [8].
We need to require positive definiteness of the matrix π΅(π) from (7). Since the βmatrix πΎπΏ π΅π πΎπ is
Μ β1 Μ
negative definite, then π΅0 (π) > 0, if the parameter π β (0, ππππ₯ ), where ππππ₯ = 2 (β πΌπΌππ ).
Notice, that πΌπ and π½π satisfy the following inequalities:
π 2 π 2
πΌπ β₯ πΌπ + ( ) πΌπ , π½π β€ π½π + ( ) π½π . (10)
2 2
Theorem 1. [8] Let the matrices π΄ and π΅(π) be positive definite. Then the PSTS iteration method is
convergent, i.e., the spectral radius π(πΊ(π, π )) of its iteration matrix πΊ(π, π ) is less than 1, provided that
the acceleration parameters π and π satisfy
0 < π < π, 0 < π < ππππ₯ ,
and
2 π½β
0<π < , Ξ= . (11)
Ξ πΌπ + ( π2 )2 πΌπ
Let us study the convergence of the MGMs with the PSTS smoothers.
We use some denotations and theoretical results from [11, 12].
Let π»1 β π»2 β ... be a family of nested finite-dimensional linear spaces. The dimension of π»π is ππ
and the inner product is denoted by (.,.)π with β₯ . β₯π corresponding norm in π»π , π = 1, 2...
We are interested in solution of the problem (1) in π»π
π΄π π’ = π.
Let this problem have a unique solution for some π β π»π and
Μπ + π΄Μπ ,
π΄π = π΄
where π΄ Μπ is a symmetric positive definite operator in π»π . Let ππ is the other symmetric positive
definite operator in π»π and the condition for the spectral radius of the operator πΊπ = ππ Μπ is
β1 π΄
satisfied π (πΊπ ) = 1.
Using the operators π΄Μπ and ππ we define the energy norm
Μπ π’, π’)1/2 ,
βπ’βπ΄Μπ = (π΄ βπ’ β π»π ,
and the inner product
(π’, π£)ππ = (ππ π’, π£), βπ’, π£ β π»π .
Then we may define the following norms on π»π :
β₯ π’ β₯π ,π =β₯ π’ β₯πΊππ = (πΊππ π’, π’)1/2 , π = 0, 1, ...
ππ
β₯ π’ β₯1,π = (π΄Μπ π’, π’)1/2 =β₯ π’ β₯ Μ ,
π΄π
β₯ π’ β₯0,π = (ππ π’, π’)1/2 = (π’, π’)1/2
ππ .
Define the subspace πΉπ β π»π , πΉπ = {π’ β π»π βΆ (π΄π π’, π£) = 0, βπ£ β π»πβ1 }
Now we make three basic assumptions [11, 12]:
Assumption 1. There exists π£ β π»πβ1 , πΎ , 0 < πΎ β€ 1 and πΏ < β such that for β π’ β π»π
β₯ π’ β π£ β₯21,π β€ πΏ πΎ β₯ π’ β₯21+πΎ ,π
Assumption 2. There exists ππ , ππ β 0 (π β β) such that
β£ (π΄Μπ π’, π£)π β£< ππ β₯ π’ β₯1,π β₯ π£ β₯1,π
for β π’ β πΉπ , βπ£ β π»π , π is a positive integer.
Assumption 3. There exists ππ , ππ β 0 (π β β) such that
β£ (π΄Μπ π’, π£ ) β£β€ ππ β₯ π’ β₯1,π β₯ π£ β₯0,π
π
for βπ’, π£ β π»π , π is a positive integer.
It is assumed that ππ , ππ are sufficiently small.
In [12] it was shown that assumptions 1-3 were satisfied for sufficiently wide class of boundary
problems in two dimension-limited domains with different boundary conditions. Assumption 1 means
a generalization of the usual approximation property known for the symmetric case. Assumptions 2
and 3 restrict the non-symmetric part π΄Μπ of the operator π΄π .
Let us consider two-grid algorithm. We denote the exact solution of the problem (1) by π’ β . Let
π’ 0 is an initial guess, π’ 1 is a problem solution after smoothing and π’ 2 is a problem solution after
MGM-iteration. We need to estimate the energy-norm of the contraction number defined by
β₯ π’ 2 β π’ β β₯1
π = sup , π’0 β π’β.
β₯ π’ 0 β π’ β β₯1
Denote the error by π π = π’ π β π’ β , π = 0, 1, 2.
Theorem 2. [12] Let the three basic assumptions for the two-grid method be satisfied. Furthermore, let
the following smoothing proposition also be satisfied: there exists Ξ (0 < Ξ < β) and π > 0 such that
β₯ π 1 β₯21 +π β₯ π 1 β₯22 β€ (1 + πΞ) β₯ π 0 β₯21 (12)
with π=ππ from Assumption 3. Then π β€ πΜ for the two-grid contraction number where
β§
βͺ 1/2 β«
βͺ
βͺ
βͺ π 2 +π 2 π 2 +2πππ π βͺ
βͺ
1βπ (1 + πΞ) βΆ
πΜ β‘ πΜ (π) β‘ sup β¨ ( 1+πΏ β1 π ( 1+π )2/πΎ π 2/πΎ ) β¬
βͺ
βͺ βͺ
βͺ π 2 + π 2 β 2ππ π β€ 1, π , π β₯ 0 βͺ βͺ
β© β
β₯ π 1 + π’ β β₯1 β₯ π’ β β₯1
π = , π= ,
β₯ π 1 β₯1 β₯ π 1 β₯1
where constants π and Ξ depend on properties of smoothing system and constants πΎ , πΏ, π, π taken from
assumptions 1, 2, 3 respectively, while π is calculation accuracy.
In [12] the proof of the same theorem for the multigrid method is given. For further considerations
it is necessary to have one more theorem from [11].
Theorem 3. [11] Let the problem (1) be solved by iterative method (5) - (6) written in the form of
π’ (π+1) = π’ (π) β π΅Μβ1 (π΄π’ (π) β π ) , Μ = π΄0 ,
π΄
and let the three basic assumptions be satisfied. If there exists the constant π > 0 so that inequality is
satisfied
π΅Μ + π΅Μβ β π΄ Μ β π΅Μ)β π β1 (π΄
Μ β₯ π (π΄ Μ β π΅Μ) , (13)
then the smoothing assumption (12) is satisfied with
Ξ = (1 + π ) π½ (ππ½ + 2) , (14)
where
Μ (π΅Μβ1 )β π΄
π½ = (π (π΅Μβ1 π΄ Μ))1/2 .
Now we can prove the convergence of the method.
Theorem 4. For the method (5) - (6) there exists the constant π > 0 so that inequality (13) is satisfied.
Proof. In case of using method (5)-(6) π΅Μ = 1/π π΅, π = π΄0 . Consider inequality (13) and transform
left and right parts of it. Using (7) we obtain
π΅Μ + π΅Μβ β π΄
Μ = 2/π π΅0 β π΄0 , (15)
Μ β π΅Μ)β π΄
(π΄ Μβ1 (π΄ Μ β π΅Μ) = (π΄0 β 1 π΅β ) π΄β1 1
0 (π΄0 β π π΅ ) =
π (16)
= π΄0 β π1 π΅β β π1 π΅ + π12 π΅β π΄β1 2 1 β β1
0 π΅ = π΄0 β π π΅0 + π 2 π΅ π΄0 π΅
Denote
2
π= π΅0 β π΄0 , (17)
π
then from (8) - (11) itβs easy to obtain that
π½β 2
π 2 < ,
πΌπ + ( 2 ) πΌπ π
thus, 2/π π΅0 > π΄0 for acceleration parameters π, π satisfying the conditions of the Theorem 1, and
π = π β > 0.
From (15) - (17) we obtain that for inequality (13):
1 β β1
π β₯ π (βπ + π΅ π΄0 π΅ ) .
π2
Multiplying the left and the right parts of this inequality on π β1/2 we obtain
1 β1/2 β β1 β1/2
πΌ β₯ π (βπΌ + π π΅ π΄0 π΅π ) . (18)
π2
Denote πΏ = π1 π΄β1/2
0 π΅π
β1/2 . Then (18) is transformed to
πΌ β₯ π (πΏβ πΏ β πΌ ) .
If we take
1
π= , (19)
β₯ πΏβ πΏ β πΌ β₯
then inequality (13) is satisfied. The Theorem 4 is proved.
Thus, if the conditions of the Theorem 4 are satisfied then two-grid method converges with the
PSTS-smoothers (6), π is calculated by (19) and Ξ is calculated by (14) correspondingly. The results
obtained for the two-grid method can be extended to the multigrid method.
4. Numerical experiments
Let us consider the 2-D unsteady Navier-Stokes equations of viscous incompressible fluids. The prim-
itive variables mathematical formulation of this problem is: find V and π (up to constant) so that
πV
β πΞV + (V β
β)V + βπ = f ππ Ξ© Γ (0, π ],
ππ‘
πππ£V = 0 ππ Ξ© Γ [0, π ],
V=g ππ πΞ© Γ [0, π ],
V(x, 0) = V0 (x) ππ Ξ©
where π is the kinematic viscosity coefficient; Ξ is the Laplacian, β is the gradient, πππ£ is the diver-
gence; V = (π€(π₯, π¦, π‘), π£(π₯, π¦, π‘)) is the velocity vector, π is the pressure, g is the velocity prescribed on
the boundary πΞ© and V0 is the divergence-free initial velocity field. We assume that the fluid motion
occurs in the time interval [0, π ].
One of the main difficulties in using the AMG is the constructing smoothers that are robust over
a wide range of the viscosity coefficient, especially for small values of π. Classical MGMs were de-
veloped for elliptic PDEs. When applied to nonelliptic and singular perturbation problems such as
high-Reynolds flows, performance seems to deteriorate significantly.
As a model problem, we consider the standard lid-driven cavity problem [13] in the square domain
Ξ© = (0, 1) Γ (0, 1). Let us write this problem in scalar dimensionless form
ππ€ ππ€ ππ€ ππ 1 π2π€ π2π€
+π€ +π£ + β + = 0, (20)
ππ‘ ππ₯ ππ¦ ππ₯ π
π ( ππ₯ 2 ππ¦ 2 )
ππ£ ππ£ ππ£ ππ 1 π2π£ π2π£
+π€ +π£ + β + = 0, (21)
ππ‘ ππ₯ ππ¦ ππ¦ π
π ( ππ₯ 2 ππ¦ 2 )
ππ€ ππ£
+ = 0, (22)
ππ₯ ππ¦
π€ = 0, π£ = 0, (π₯ = 0, π₯ = 1, π¦ = 0),
π€ = 1, π£ = 0, (π¦ = 1), (23)
π€ (π₯, π¦, 0) = 0, π£ (π₯, π¦, 0) = 0.
where π
π is the Reynolds number (π
π = π πΏ/π, π is a characteristic velocity of the flow, πΏ is a char-
acteristic length scale). Homogeneous Dirichlet boundary conditions are prescribed for all velocity
components with the exception of a positive unit horizontal velocity along the top face. The initial
condition assumes the fluid to be at rest. We consider this problem as a test for study the effect of the
suggested smoothers for the MGM.
The fully implicit scheme has been used for the time discretization. The use of implicit schemes re-
moves restrictions on the time integration step, which is selected based on the required computational
accuracy. We fix the time step πΏπ‘ and introduce a discrete time grid π‘ π = ππΏπ‘, π β₯ 0:
1 (π+1) 1
(V β V(π) ) + βπ (π+1) = β(V(π+1) β
β)V(π+1) + ΞV(π+1) , (24)
πΏπ‘ π
π
πππ£V(π+1) = 0, (25)
The boundary conditions are taken at π‘ = π‘ π+1 :
π€ (π+1) = 0, π£ (π+1) = 0, (π₯ = 0, π₯ = 1, π¦ = 0),
π€ (π+1) = 1, π£ (π+1) = 0, (π¦ = 1).
The uniform grid Ξ© is introduced in the domain Ξ© with steps β1 and β2 ; β1 = 1/π1 , β2 = 1/π2 , where
π1 , π2 are the number of cells in each direction. The grid cells are positioned such that the cell faces
coincide with the boundary πΞ© of Ξ©. The spatial discretization of the Navier-Stokes equations is
performed on MAC (Marker-and-Cell) [14] (staggered) grids when pressure and velocities in two-
dimensional problems are determined on three grids shifted relative to each other. So, the pressure
π is located in the center of each cell, the π₯-component velocity π€ is on the middle points of vertical
faces, the π¦-component velocity π£ is on the middle points of horizontal faces. The continuity equation
(25) is discretized at cell centers using central difference schemes.
Thus we introduce the grid sets and the corresponding spaces:
π· 1 = {π₯ππ = ((π + 1/2)β1 , πβ2 ) βΆ π = 0, ..., π1 β 1, π = 0, ..., π2 },
π· 2 = {π₯ππ = (πβ1 , (π + 1/2)β2 ) βΆ π = 0, ..., π1 , π = 0, ..., π2 β 1},
π·3 = {π₯ππ = (πβ1 , πβ2 ) βΆ π = 1, ..., π1 β 1, π = 1, ..., π2 β 1}.
Let πβ = π1,β Γ π2,β be the linear space of vector functions defined on π· 1 Γ π· 2 and vanishing at the
corresponding grid boundaries, and πβ is the space of functions defined on π·3 and orthogonal to unity.
Thus,
π1,β = {π€ππ = π€(π₯ππ ) βΆ π₯ππ β π· 1 , π€0,π = π€π1 β1,π = π€π,0 = π€π,π2 = 0},
π2,β = {π£ππ = π£(π₯ππ ) βΆ π₯ππ β π· 2 , π£0,π = π£π1 ,π = π£π,0 = π£π,π2 β1 = 0},
πβ = {πππ = π(π₯ππ ) βΆ π₯ππ β π·3 , β β1 β2 πππ = 0}.
ππ
Variables are denoted by a single set of indices, despite the fact that different variables are calculated
at different grid nodes. As a result, the indices π, π refer to a set of three mismatched points.
The equations contain the discrete nonlinear terms. For treating this non-linearity Newton lin-
earization around the old time level is used. If we want to linearize a nonlinear term π€ (π+1) ππ₯(π+1) ,
then
π€ (π+1) ππ₯(π+1) = π€ (π) ππ₯(π+1) + π€ (π+1) ππ₯(π) β π€ (π) ππ₯(π) + π(πΏπ‘ 2 ). (26)
Table 1
AMG+STS Iterations and CPU with Different π.
Grid π = 10β1 π = 10β2 π = 10β3 π = 10β4 π = 10β5
60 Γ 60 36 (38.18) 34 (36.42) 31 (32.20) 28 (30.51) 20 (22.11)
120 Γ 120 47 (62.32) 40 (55.54) 37 (53.60) 33 (50.20) 22 (29.31)
180 Γ 180 49 (87.61) 44 (84.52) 38 (79.62) 34 (76.85) 24 (48.41)
260 Γ 260 54 (122.71) 47 (129.33) 43 (132.42) 35 (110.26) 27 (83.33)
520 Γ 520 59 (251.51) 56 (248.58) 54 (235.20) 44 (215.64) 36 (151.19)
Table 2
AMG+PSTS Iterations and CPU with Different π.
Grid π = 10β1 π = 10β2 π = 10β3 π = 10β4 π = 10β5
60 Γ 60 42 (46.82) 38 (45.70) 37 (44.24) 32 (42.36) 22 (26.82)
120 Γ 120 49 (64.38) 42 (57.22) 39 (60.52) 35 (54.81) 24 (36.55)
180 Γ 180 51 (97.51) 46 (94.32) 40 (87.26) 36 (82.64) 26 (54.81)
260 Γ 260 56 (157.84) 49 (150.57) 45 (142.81) 37 (124.83) 29 (91.29)
520 Γ 520 62 (287.86) 59 (283.35) 57 (280.55) 49 (232.36) 40 (187.65)
The expression in the right-hand side of (26) is linear in the variables at the new time level and
possesses a discretization error π(πΏπ‘ 2 ).
After discretization and linearization of the problem (20)-(23) we need to solve large sparse strongly
nonsymmetric SLAEs at each time step. The choice of an appropriate iterative method and its imple-
mentation largely determine the overall efficiency of the computational algorithm. In practice, clas-
sical iterative methods are used, such as BiCG, GMRES, GMRES(m) [15] and others. We use the AMG
with the special STS- and PSTS- smoothers. Both of the smoothers are compared with the Gauss-
Seidel (GS) one which is the classic smoothing method [3]. The number of the pre-smoothing and
post-smoothing steps [2] equals 15 for all methods. Calculations are carried out using the V-cycle
[2, 3] as the least demanding from the computational point of view. For the PSTS method, π΅π = πΌ and
πΉ0 was chosen such that the matrix (πΉ0 + πΎπΏ ) was unitary [8]. In the numerical experiments, the time
step is equal to 0.01 and the Reynolds number is π
π = π β1 .
In Tables 1, 2, 3 we compare the AMG-iterations with the STS-, PSTS- and GS-based smoothers on
the different grids for small viscosity values. In our implementations, all iterations are started from
the zero vector, and terminated either when
βπ (π) β2
β€ 10β6 ,
βπ (0) β2
where π (π) = π β π΄π’ (π) is the residual vector of the SLAE (1) at the current iterate π’ (π) , and π (0) is the
initial residual. Our comparisons are done for the number of iteration steps (denoted by βITβ) and the
elapsed CPU time (in seconds and denoted by βCPUβ). The abbreviation βn.c.β in the tables means βno
convergenceβ. The experiments are run in MATLAB (version R2018b) with a machine precision 10β16 .
From Tables 1, 2, 3 an increase in the number of iterations follows with an increase in the mesh size
for all tested methods. This is due to some features of the algebraic approach in the multigrid method
(namely, the PMIS algorithm in the AMG). The traditional scalable approach (RS algorithm) works
well for many 2-D problems. In this case, the solver can be obtained with the number of iterations
that does not depend on the size of the problem. But when using AMG interpolation in combination
Table 3
AMG+GS Iterations and CPU with Different π.
Grid π = 10β1 π = 10β2 π = 10β3 π = 10β4 π = 10β5
60 Γ 60 39 (39.57) 46 (42.51) 54 (114.85) n.c. n.c.
120 Γ 120 57 (83.82) 64 (122.61) 83 (160.50) n.c. n.c.
180 Γ 180 59 (162.36) 82 (185.38) 85 (192.20) n.c. n.c.
260 Γ 260 n.c. n.c. n.c. n.c. n.c.
520 Γ 520 n.c. n.c. n.c. n.c. n.c.
with the PMIS, the AMG convergence deteriorates depending on the problem size [5]. However, the
RS algorithm has obvious drawbacks, such as poor efficiency for 3-D problems [5].
Moreover, from Tables 1, 2, 3 it follows that the AMG method with the STS- and the PSTS-
smoothers has fast convergence speed for all tested values of the viscosity coefficient (π = 10β1 Γ· 10β5 )
on all used grids, while the AMG with the GaussβSeidel smoother does not converge for π = 10β4 , 10β5
on all considered grids, and does not converge on the grids 260 Γ 260 and 520 Γ 520 nodes for all val-
ues of the viscosity coefficient. For all tests, the AMG+STS and AMG+PSTS methods outperform the
AMG+GS method with respect to both number of iteration steps and CPU time.
The AMG+STS method slightly ahead of the AMG+PSTS method, apparently due to the fact that
the first method belongs to the class of triangular, and the second one relates to the class of two-step
iterative methods. In the numerical tests the iterative parameters for all methods was chosen exper-
imentally. The decrease in the number of the AMG+STS and AMG+PSTS iterations with a decrease
in the viscosity coefficient values is explained that these methods converge the faster, the more pro-
nounced the nonsymmetry of the SLAE is. The number of iterations and CPU time of the AMG+GS
method increases with decreasing π.
5. Conclusions
Numerical experiments carried out for the 2-D model incompressible unsteady Navier-Stokes problem
have shown efficiency of the STS- and the PSTS- smoothers for the AMG and their advantage over
the classical Gauss-Seidel smoother. The use of these smoothers for the AMG allows solving this
problem for small values of the viscosity coefficient. The lack of convergence of the AMG+STS and
the AMG+PSTS methods for the smallest values of π is not observed, which could be expected for the
unsteady problem. The obtained SLAE is simpler than for the stationary problem, since the spectrum
of this system lies in the right half-plane. This is due to a shift in the spectrum when the time derivative
is approximated. It is supposed to solve the corresponding 3-D problem in the future using these
methods.
Acknowledgments
This research was funded by the Government of the Russian Federation (Project No.075-15-2019-1928),
Russia.
References
[1] R. D. Falgout, An introduction to algebraic multigrid, Computing in Science & Engineering 8
(2006) 24β33. doi:10.1109/MCSE.2006.105.
[2] W. Hackbusch, Multi-grid methods and applications, Springer-Verlag, Berlin, 1985.
[3] U. Trottenberg, C. Oosterlee, A. Schuller, Multigrid, 1st. ed., Academic Press, 2001.
[4] K. Stuben, Algebraic multigrid (amg): experiences and comparisons, Applied Mathematics and
Computation 13 (1983) 419β451. doi:10.1016/0096-3003(83)90023-1.
[5] H. de Sterck, U. M. Yang, J. J. Heys, Reducing complexity in parallel algebraic multigrid
preconditioners, SIAM Journal on Matrix Analysis and Applications 27 (2006) 1019β1039.
doi:10.1137/040615729.
[6] L. A. Krukier, Convergence acceleration of triangular iterative methods based on the skew-
symmetric part of the matrix, Applied Numerical Mathematics 30 (1999) 281β290.
[7] L. A. Krukier, L. G. Chikina, T. V. Belokon, Triangular skew-symmetric iterative solvers for
strongly nonsymmetric positive real linear system of equations, Applied Numerical Mathematics
41 (2002) 89β105.
[8] L. A. Krukier, T. S. Martynova, Z.-Z. Bai, Product-type skew-hermitian triangular splitting iter-
ation methods for strongly non-hermitian positive definite linear systems, Journal of Computa-
tional and Applied Mathematics 232 (2009) 3β16. doi:10.1016/j.cam.2008.10.033.
[9] L. Wang, Z.-Z. Bai, Skew-hermitian triangular splitting iteration methods for non-hermitian
positive definite linear systems of strong skew-hermitian parts, BIT Numerical Mathematics 44
(2004) 363β386. doi:10.1023/B:BITN.0000039428.54019.15.
[10] G. V. Muratova, L. A. Krukier, E. M. Andreeva, Fourier analysis of multigrid method with trian-
gular skew-symmetric smoothers, Communication on Applied Mathematics and Computation
27 (2013) 355β362. doi:10.3969/j.issn.1006-6330.2013.03.007.
[11] Z.-H. Cao, Convergence of multigrid methods for nonsymmetric, indefinite problems, Applied
Mathematics and Computation 28 (1988) 269β288. doi:10.1016/0096-3003(88)90076-8.
[12] J. Mandel, Multigrid convergence for nonsymmetric, indefinite variational problems and one
smoothing step, Applied Mathematics and Computation 19 (1986) 201β216. doi:10.1016/
0096-3003(86)90104-9.
[13] J.-L. Guermond, C. Migeon, G. Pineau, L. Quartapelle, Start-up flows in a three-dimensional
rectangular driven cavity of aspect ratio 1:1:2 at re=1000, Journal of Fluid Mechanics 450 (2002)
169β199. doi:10.1017/S0022112001006383.
[14] S. McKee, M. F. Tome, V. G. Ferreira, J. A. Cuminato, A. Castelo, F. Sousa, N. Mangiavacchi, The
mac method, Computers & Fluids 37 (2008) 907β930. doi:10.1016/j.compfluid.2007.10.
006.
[15] Y. Saad, Iterative methods for sparse linear systems, 2nd. ed., SIAM, 2003.