use CG from iterative solvers

JuliaImageRecon · Aug 15, 2019 · 1d81120 · 1d81120
1 parent 1bddd87
commit 1d81120
Show file tree

Hide file tree

Showing 3 changed files with 7 additions and 21 deletions.
diff --git a/src/ADMM.jl b/src/ADMM.jl
@@ -85,7 +85,7 @@ function admm(A, b::Vector, reg::Regularization
               , iterationsInner::Int64=10
               , absTol::Float64=1.e-8
               , relTol::Float64=1.e-6
-              , tolInner::Float64=1.e-3
+              , tolInner::Float64=1.e-5
               , adaptRho::Bool=false
               , solverInfo = nothing
               , kargs...)
@@ -118,12 +118,7 @@ function admm(A, b::Vector, reg::Regularization
     # 1. solve arg min_x 1/2|| Ax-b ||² + ρ/2 ||x+u-z||²
     # <=> (A'A+ρ)*x = A'b+ρ(z-u)
     xᵒˡᵈ[:] = x[:]
-    # cg!(x,A'*A+ρ*opEye(length(x)),β+ρ*(z-u),Pl=precon,maxiter=iterationsInner,tol=tolInner)
-    if isa(precon, Identity)
-      x[:] = cg(A'*A+ρ*opEye(length(x)),x,β+ρ*(z-u),iterations=iterationsInner,relTol=tolInner,solverInfo=solverInfo,storeIterations=true)
-    else
-      x[:] = cg(A'*A+ρ*opEye(length(x)),x,β+ρ*(z-u),precon,iterations=iterationsInner,relTol=tolInner,solverInfo=solverInfo,storeIterations=true)
-    end
+    cg!(x,A'*A+ρ*opEye(length(x)),β+ρ*(z-u),Pl=precon,maxiter=iterationsInner,tol=tolInner)
 
     # 2. update z using the proximal map of 1/ρ*g(x)
     zᵒˡᵈ[:] = z
@@ -171,7 +166,7 @@ function fadmm(A, b::Vector, reg::Regularization
               , iterationsInner::Int64=10
               , absTol::Float64=1.e-8
               , relTol::Float64=1.e-6
-              , tolInner::Float64=1.e-3
+              , tolInner::Float64=1.e-5
               , adaptRho::Bool=false
               , solverInfo = nothing
               , kargs...)
@@ -210,12 +205,7 @@ function fadmm(A, b::Vector, reg::Regularization
     # 1. solve arg min_x 1/2|| Ax-b ||² + ρ/2 ||x+û-ẑ||²
     # <=> (A'A+ρ)*x = A'b+ρ(z-u)
     xᵒˡᵈ[:] = x[:]
-    # cg!(x,A'*A+ρ*opEye(length(x)),β+ρ*(ẑ-û),Pl=precon,maxiter=iterationsInner,tol=tolInner)
-    if isa(precon, Identity)
-      x[:] = cg(A'*A+ρ*opEye(length(x)),x,β+ρ*(z-u),iterations=iterationsInner,relTol=tolInner,solverInfo=solverInfo)
-    else
-      x[:] = cg(A'*A+ρ*opEye(length(x)),x,β+ρ*(z-u),precon,iterations=iterationsInner,relTol=tolInner,solverInfo=solverInfo)
-    end
+    cg!(x,A'*A+ρ*opEye(length(x)),β+ρ*(ẑ-û),Pl=precon,maxiter=iterationsInner,tol=tolInner)
 
     # 2. update z using the proximal map of 1/ρ*g(x)
     zᵒˡᵈ[:] = z

diff --git a/src/RegularizedLeastSquares.jl b/src/RegularizedLeastSquares.jl
@@ -36,7 +36,7 @@ include("Utils.jl")
 include("Kaczmarz.jl")
 include("DAX.jl")
 include("CGNR.jl")
-include("CG.jl")
+# include("CG.jl")
 include("Direct.jl")
 include("FusedLasso.jl")
 include("FISTA.jl")

diff --git a/src/SplitBregman.jl b/src/SplitBregman.jl
@@ -84,7 +84,7 @@ function splitBregman(A, y::Vector, reg::Vector{Regularization}
               , iterationsCG::Int64=10
               , absTol::Float64=1.e-8
               , relTol::Float64=1.e-6
-              , tolInner::Float64=1.e-3
+              , tolInner::Float64=1.e-6
               , solverInfo = nothing
               , kargs...)
 
@@ -117,11 +117,7 @@ function splitBregman(A, y::Vector, reg::Vector{Regularization}
     for k = 1:iterationsInner
       # update u analytically
       β[:] .= β_yj .+ λ*(w-bw) .+ ρ*(v-bv)
-      if isa(precon, Identity)
-        u[:] = cg(op,u,β,iterations=iterationsCG,relTol=tolInner,solverInfo=solverInfo,storeIterations=true)
-      else
-        u[:] = cg(op,u,β,precon,iterations=iterationsCG,relTol=tolInner,solverInfo=solverInfo,storeIterations=true)
-      end
+      cg!(u,op,β,Pl=precon,maxiter=iterationsCG,tol=tolInner)
 
       #  proximal map for L1 regularization
       wᵒˡᵈ[:] .= w