Merge pull request #12 from ranjanan/moreopt

Fix preconditioning

Author: ranjanan

src/multilevel.jl

@@ -58,17 +58,18 @@ abstract type Cycle end
 struct V <: Cycle
 end
 
-function solve{T}(ml::MultiLevel, b::Vector{T}; maxiter = 100,
-                                                cycle = V(),
-                                                tol = 1e-5)
+function solve{T}(ml::MultiLevel, b::Vector{T},
+                                    maxiter = 100,
+                                    cycle = V(),
+                                    tol = 1e-5)
     x = zeros(T, size(b))
     residuals = Vector{T}()
     A = ml.levels[1].A
     normb = norm(b)
     if normb != 0
         tol *= normb
     end
-    push!(residuals, norm(b - A*x))
+    push!(residuals, T(norm(b - A*x)))
 
     lvl = 1
     while length(residuals) <= maxiter && residuals[end] > tol
@@ -77,7 +78,7 @@ function solve{T}(ml::MultiLevel, b::Vector{T}; maxiter = 100,
         else
             x = __solve(cycle, ml, x, b, lvl)
         end
-        push!(residuals, norm(b - A * x))
+        push!(residuals, T(norm(b - A * x)))
     end
     x
 end
@@ -103,4 +104,5 @@ function __solve{T}(v::V, ml, x::Vector{T}, b::Vector{T}, lvl)
     x
 end
 
-coarse_solver(::Pinv, A, b) = pinv(full(A)) * b
+coarse_solver{Tv,Ti}(::Pinv, A::SparseMatrixCSC{Tv,Ti}, b::Vector{Tv}) =
+                                                        pinv(full(A)) * b

src/preconditioner.jl

@@ -7,7 +7,7 @@ end
 aspreconditioner(ml::MultiLevel) = Preconditioner(ml)
 
 \(p::Preconditioner, b) = p * b
-*(p::Preconditioner, b) = solve(p.ml, b; cycle = V(), maxiter = 1, tol = 1e-12)
+*(p::Preconditioner, b) = solve(p.ml, b, 1, V(), 1e-12)
 
 A_ldiv_B!(x, p::Preconditioner, b) = copy!(x, p \ b)
 A_mul_B!(b, p::Preconditioner, x) = A_mul_B!(b, p.ml.levels[1].A, x)

src/splitting.jl

@@ -4,18 +4,6 @@ const U_NODE = 2
 
 struct RS
 end
-#=function split_nodes(::RS, S)
-	n = size(S, 1)
-	for i = 1:n
-		for j in nzrange(S, i)
-			row = S.rowval[j]
-			if row == i
-				S.nzval[j] = 0
-			end
-		end
-	end
-	RS_CF_splitting(sparse(i,j,v,n,n), sparse(j,i,v,n,n))
-end=#
 
 function remove_diag!(a)
 	n = size(a, 1)

test/runtests.jl

@@ -2,6 +2,7 @@ using AMG
 using Base.Test
 using JLD
 using IterativeSolvers
+import AMG: V, coarse_solver, Pinv, Classical
 
 @testset "AMG Tests" begin
 
@@ -35,7 +36,7 @@ S = sprand(10,10,0.1); S = S + S'
 @test split_nodes(RS(), S) ==  [0, 1, 1, 0, 0, 0, 0, 0, 1, 1]
 
 a = load("thing.jld")["G"]
-S, T = AMG.strength_of_connection(AMG.Classical(0.25), a)
+S, T = AMG.strength_of_connection(Classical(0.25), a)
 @test split_nodes(RS(), S) == [0, 0, 1, 0, 1, 1, 0, 1, 0, 0, 1, 0, 1, 1, 0, 0, 1, 0,
 0, 0, 1, 0, 1, 0, 1, 0, 0, 1, 0, 0, 1, 0, 1, 0, 1, 0, 1, 0, 0, 0, 0, 1, 1, 0,
 1, 0]
@@ -62,9 +63,9 @@ ml = ruge_stuben(A)
 end
 
 @testset "Coarse Solver" begin
-A = poisson(10)
+A = float.(poisson(10))
 b = A * ones(10)
-@test sum(abs2, AMG.coarse_solver(AMG.Pinv(), A, b) - ones(10)) < 1e-6
+@test sum(abs2, coarse_solver(Pinv(), A, b) - ones(10)) < 1e-6
 end
 
 @testset "Multilevel" begin
@@ -120,7 +121,7 @@ p = aspreconditioner(ml)
 b = zeros(n)
 b[1] = 1
 b[2] = -1
-x = solve(p.ml, A * ones(n); maxiter = 1, tol = 1e-12)
+x = solve(p.ml, A * ones(n), 1, V(), 1e-12)
 diff = x - [  1.88664780e-16,   2.34982727e-16,   2.33917697e-16,
          8.77869044e-17,   7.16783490e-17,   1.43415460e-16,
          3.69199021e-17,   9.70950385e-17,   4.77034895e-17,
@@ -138,7 +139,7 @@ diff = x - [  1.88664780e-16,   2.34982727e-16,   2.33917697e-16,
         -6.76965535e-16,  -7.00643227e-16,  -6.23581397e-16,
         -7.03016682e-16]
 @test sum(abs2, diff) < 1e-8
-x = solve(p.ml, b; maxiter = 1, tol = 1e-12)
+x = solve(p.ml, b, 1, V(), 1e-12)
 diff = x - [ 0.76347046, -0.5498286 , -0.2705487 , -0.15047352, -0.10248021,
         0.60292674, -0.11497073, -0.08460548, -0.06931461,  0.38230708,
        -0.055664  , -0.04854558, -0.04577031,  0.09964325,  0.01825624,