Skip to content

cortner/PyAMG.jl

Folders and files

NameName
Last commit message
Last commit date

Latest commit

 

History

68 Commits
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 

Repository files navigation

PyAMG.jl

Build Status

Convenience wrapper module for the PyAMG library. Requires a Python installation with scipy and pyamg. If an Anconda distribution is used (including the Conda.jl package manager) then pyamg should be automatically installed on first use. Otherwise, follow the straightforward instructions.

Basic Usage

In all examples it is assumed that A is a sparse matrix and b a vector and amg is an AMGSolver instance constructed from A. The classical example would be the Dirichlet problem on a square,

using SparseArrays, LinearAlgebra
N = 100
L1 = spdiagm( -1 =>  -ones(N-1),
               0 => 2*ones(N),
               1 =>  -ones(N-1) ) * N^2
Id = sparse(1.0*I, N, N)
A = kron(Id, L1) + kron(L1, Id)
b = ones(size(A,1))

Blackbox solver

using PyAMG
x = PyAMG.solve(A, b)

Multiple solves

To initialise, call

using PyAMG
amg = RugeStubenSolver(A)

Then the system Ax=b can be solved using

x = amg \ b
x = solve(amg, b; tol=1e-6)

or, one can specify a different 'outer solver'

x = solve(amg, b; tol=1e-6, accel="cg")

see ?solve for more options. In particular, note the that default keyword arguments can be set via set_kwargs! or in the construction of the AMG solver, which will then be used by both \ and solve. E.g.,

amg = RugeStubenSolver(A, tol=1e-6, accel="cg")
x = amg \ b

As Preconditioner

After initialising, we can construct a preconditioner via

M = aspreconditioner(amg)

The line M \ b then performes a single MG cycle. This is e.g. compatible with the IterativeSolvers package:

using PyAMG, IterativeSolvers
TOL = 1e-3
M = aspreconditioner(RugeStubenSolver(A))
IterativeSolvers.cg(A, b, M; tol=TOL)

Solver history

To extract the solver history as a vector of residuals, use

amg = RugeStubenSolver(A)
r = Float64[]
x = PyAMG.solve(amg, b, residuals=r)
@show r

Since version 3.2.1.dev0+2227b77 the residuals can also be returned for the blackbox solver variant.

(NOTE: although pyamg needs residuals to be a list, PyAMG.jl will detect if residuals is a numpy vector and replace it with a list, then convert back to a types Julia vector.)

List of Types and Methods

(this section is incomplete)

Types

  • AMGSolver{T} : encapsulates the pyamg solver PyObject
  • RugeStubenSolver : typealias for AMGSolver{RugeStuben}
  • SmoothedAggregationSolver : typealias for AMGSolver{SmoothedAggregation}
  • AMGPreconditioner : encapsulates the output of aspreconditioner to use PyAMG as a preconditioner for iterative linear algebra.