Why the SPAI0-BiCGSTAB outperforms AMG-BiCGSTAB for large scale problem?

[QuanchengP]

Dear developer, the AMGCL is a wonderful choice as the solver of elastic frictional contact problem (the SPAI0-BiCGSTAB in AMGCL outperforms the default BiCGSTAB in Eigen). As the left-handed matrix is square but not symmetric, the BiCGSTAB but not CG is chosen as solver (GMRES is compared and not chosen). It is unexpected that the SPAI0-BiCGSTAB outperforms AMG-BiCGSTAB for large scale problem.

Here, SPAI0-BiCGSTAB means "the solver is BiCGSTAB, the preconditioner is SPAI0"; the AMG-BiCGSTAB means "the solver is BiCGSTAB, the preconditioner is AMG". For AMG-BiCGSTAB, different relaxation schemes are already compared and the more efficient chebyshev relaxation is chosen.

The tested matrix, RHS and code are at: https://github.com/QuanchengP/DEHW/releases/download/v1.0.0/demo.zip (100Mb). The tested results are as follows:

(1) SPAI0-BiCGSTAB

Solver
======
Type: BiCGStab
Unknowns: 135368
Memory footprint: 7.23 M

Preconditioner
==============
Relaxation as preconditioner
Unknowns: 135368
Nonzeros: 11818846
Memory: 182.41 M

Iters: 202
Error: 6.37313e-15

[Frictional contact: 2.407 s] (100.00%)
[ setup: 0.031 s] ( 1.28%)
[ solve: 2.376 s] ( 98.71%)

(2) AMG-BiCGSTAB

Solver
======
Type: BiCGStab
Unknowns: 135368
Memory footprint: 7.23 M

Preconditioner
==============
Number of levels: 2
Operator complexity: 1.02
Grid complexity: 1.02
Memory footprint: 222.20 M

level unknowns nonzeros memory
---------------------------------------------
0 135368 11818846 209.54 M (98.21%)
1 2631 215441 12.65 M ( 1.79%)

Iters: 56
Error: 7.05585e-15

[Frictional contact: 9.376 s] (100.00%)
[ setup: 0.678 s] ( 7.23%)
[ solve: 8.698 s] ( 92.77%)

For smaller problems with unknows less than 10 thousand, AMG-BiCGSTAB slightly outperforms SPAI0-BiCGSTAB. For larger problems with unknowns more than one million, the advantages of SPAI0-BiCGSTAB becomes much more obvious.

[ddemidov]

Hello, for elasticity problems it is often beneficial to use near-nullspace vectors, see https://amgcl.readthedocs.io/en/latest/tutorial/Nullspace.html

[QuanchengP]

Thank you for the reply. We belive the better implementation of AMG-BiCGSTAB will outperform SPAI0-BiCGSTAB. However, at the present, there are two obstacles for us to adopt AMG-BiCGSTAB. (1) the system matrix in realistic problem may not have block structure (3 \times 3 block in three dimension): the Dirichlet boundary may not always simutaneously constraint three displacement components of one node, the contact problem may only constraint normal direction displacement and omit other two tangential displacement components. (2) lose of block structure makes that the current calculation of nullspace vectors either difficult or uneconomic.

In this case, it would be better for us to learn more about the AMG preconditioners that is specially designed for elastic problem (hierachical block structure, strong connection definition), also better to learn more about the implementation detail of AMGCL. The tested results does not come soon and thanks greatly again.

[ddemidov]

the system matrix in realistic problem may not have block structure

if the non-regular part of the matrix is small, you could use something like schur_pressure_correction preconditioner to split it into regular and non-regular parts, and then solve the regular one using block solver

[QuanchengP]

Take the multiplier of Dirichlet boundary and contact constraint as pressure, so the composite preconditioners can be used. This way may cause serious increase of condition number ($\lambda _{max} / \lambda _{min}$) for system matrix in elastic problem. This way is also mentioned in some classical textbooks as block diagonal preconditioner/smoother. These different ways are beyond our judgement at now.

[QuanchengP]

For the above-mentioned test problem, the GMG-BiCGSTAB only needs 8 iterations to converge (less than the 56 iterations of AMG-BiCGSTAB, the detailed test results are not given because our implementation is bad and the time cost is still larger than SPAI0-BiCGSTAB in AMGCL). Here, the GMG means geometric multigrid preconditioner (special attention is given to the hierarchical nesting of fem space with Dirichlet boundary and contact constraint).