AlgebraicMultigrid.jl
AlgebraicMultigrid.jl copied to clipboard
Extend Ruge-Stuben solver to complex type
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.
Could you please add a test for this? Happy to merge after that.
@ranjanan I suggest we merge this since it only generalizes the code and put a new release out.