JuliaLinearAlgebra / AlgebraicMultigrid.jl

Algebraic Multigrid in Julia
Other
116 stars 23 forks source link

Extend Ruge-Stuben solver to complex type #96

Open learning-chip opened 2 years ago

learning-chip commented 2 years ago

This PR shows minimum changes to allow ruge_stuben to take complex-valued matrix.

It seems to converge well for simple problems:

using LinearAlgebra
using AlgebraicMultigrid
A = poisson((32, 32, 32)) + 1.0im * I
ml = ruge_stuben(A)

n = size(A)[1]
b = ones(n) * (1.0 + 1.0im)
x = AlgebraicMultigrid._solve(ml, b; verbose=true, maxiter=10)
@assert norm(b - A*x) < 1e-6
Norm of residual at iteration      1 is 2.5600e+02
Norm of residual at iteration      2 is 4.5756e+00
Norm of residual at iteration      3 is 5.9946e-02
Norm of residual at iteration      4 is 9.5474e-04
Norm of residual at iteration      5 is 1.7687e-05

From literature on complex-valued AMG (Algebraic Multigrid Solvers for Complex-Valued Matrices and A multigrid-based shifted Laplacian preconditioner for a fourth-order Helmholtz discretization; see also pyamg/pyamg#327), the AMG interpolation formula carries straightly from real to complex cases. So this PR just tweaks max(), > operators to support complex numbers.

ranjanan commented 1 year ago

Could you please add a test for this? Happy to merge after that.

ViralBShah commented 6 months ago

@ranjanan I suggest we merge this since it only generalizes the code and put a new release out.