Skip to content

Commit

Permalink
New API (#95)
Browse files Browse the repository at this point in the history
  • Loading branch information
lopezm94 committed Nov 11, 2016
1 parent e5cc6d5 commit e2a7fc4
Show file tree
Hide file tree
Showing 17 changed files with 1,406 additions and 582 deletions.
2 changes: 2 additions & 0 deletions REQUIRE
Original file line number Diff line number Diff line change
@@ -1 +1,3 @@
julia 0.5-

UnicodePlots
108 changes: 98 additions & 10 deletions src/cg.jl
Original file line number Diff line number Diff line change
@@ -1,21 +1,34 @@
export cg, cg!

cg(A, b, Pl=1; kwargs...) = cg!(zerox(A,b), A, b, Pl; kwargs...)
####################
# API method calls #
####################

function cg!(x, A, b, Pl=1; tol::Real=size(A,2)*eps(), maxiter::Int=size(A,2))
history = ConvergenceHistory()
history[:tol] = tol
reserve!(history,:resnorm, maxiter)
cg(A, b; kwargs...) = cg!(zerox(A,b), A, b; kwargs...)

function cg!(x, A, b;
tol::Real=size(A,2)*eps(), maxiter::Integer=size(A,2),
plot=false, log::Bool=false, kwargs...
)
(plot & !log) && error("Can't plot when log keyword is false")
K = KrylovSubspace(A, length(b), 1, Vector{Adivtype(A,b)}[])
init!(K, x)
cg_method!(history, x, K, b, Pl; tol=tol, maxiter=maxiter)
x, history
history = ConvergenceHistory(partial=!log)
history[:tol] = tol
reserve!(history,:resnorm, maxiter)
cg_method!(history, x, K, b; tol=tol, maxiter=maxiter, kwargs...)
plot && (shrink!(history); showplot(history))
log ? (x, history) : x
end

function cg_method!(log::ConvergenceHistory, x, K::KrylovSubspace, b, Pl=1;
tol::Real=size(K.A,2)*eps(), maxiter::Integer=size(K.A,2))
#########################
# Method Implementation #
#########################

function cg_method!(log::ConvergenceHistory, x, K, b;
Pl=1,tol::Real=size(K.A,2)*eps(),maxiter::Integer=size(K.A,2), verbose::Bool=false
)
verbose && @printf("=== cg ===\n%4s\t%7s\n","iter","resnorm")
tol = tol * norm(b)
r = b - nextvec(K)
p = z = Pl\r
Expand All @@ -30,14 +43,89 @@ function cg_method!(log::ConvergenceHistory, x, K::KrylovSubspace, b, Pl=1;
r -= α*q
resnorm = norm(r)
push!(log,:resnorm,resnorm)
verbose && @printf("%3d\t%1.2e\n",iter,resnorm)
resnorm < tol && break
z = Pl\r
oldγ = γ
γ = dot(r, z)
β = γ/oldγ
p = z + β*p
end
shrink!(log)
verbose && @printf("\n")
setconv(log, 0<=norm(r)<tol)
x
end

#################
# Documentation #
#################

let
#Initialize parameters
doc_call = """ cg(A, b)
"""
doc!_call = """ cg!(x, A, b)
"""

doc_msg = "Solve A*x=b with the conjugate gradients method."
doc!_msg = "Overwrite `x`.\n\n" * doc_msg

doc_arg = ""
doc!_arg = """* `x`: initial guess, overwrite final estimation."""

doc_version = (doc_call, doc_msg, doc_arg)
doc!_version = (doc!_call, doc!_msg, doc!_arg)

i = 0
docstring = Vector(2)

#Build docs
for (call, msg, arg) in [doc_version, doc!_version] #Start
i+=1
docstring[i] = """
$call
$msg
If `log` is set to `true` is given, method will output a tuple `x, ch`. Where
`ch` is a `ConvergenceHistory` object. Otherwise it will only return `x`.
The `plot` attribute can only be used when `log` is set version.
**Arguments**
$arg
* `A`: linear operator.
* `b`: right hand side.
*Keywords*
* `Pl = 1`: left preconditioner of the method.
* `tol::Real = size(A,2)*eps()`: stopping tolerance.
* `maxiter::Integer = size(A,2)`: maximum number of iterations.
* `verbose::Bool = false`: print method information.
* `log::Bool = false`: output an extra element of type `ConvergenceHistory`
containing extra information of the method execution.
* `plot::Bool = false`: plot data. (Only when `log` is set)
**Output**
*`log` is `false`:*
* `x`: approximated solution.
*`log` is `true`:*
* `x`: approximated solution.
* `ch`: convergence history.
*ConvergenceHistory keys*
* `:tol` => `::Real`: stopping tolerance.
* `:resnom` => `::Vector`: residual norm at each iteration.
"""
end

@doc docstring[1] -> cg
@doc docstring[2] -> cg!
end
117 changes: 103 additions & 14 deletions src/chebyshev.jl
Original file line number Diff line number Diff line change
@@ -1,30 +1,43 @@
export chebyshev, chebyshev!

chebyshev(A, b, λmin::Real, λmax::Real, Pr = 1, n = size(A,2);
tol::Real = sqrt(eps(typeof(real(b[1])))), maxiter::Int = n^3) =
chebyshev!(zerox(A, b), A, b, λmin, λmax, Pr, n; tol=tol,maxiter=maxiter)
####################
# API method calls #
####################

function chebyshev!(x, A, b, λmin::Real, λmax::Real, Pr = 1, n = size(A,2);
tol::Real = sqrt(eps(typeof(real(b[1])))), maxiter::Int = n^3
)
history = ConvergenceHistory()
history[:tol] = tol
reserve!(history,:resnorm,maxiter)
chebyshev(A, b, λmin::Real, λmax::Real; kwargs...) =
chebyshev!(zerox(A, b), A, b, λmin, λmax; kwargs...)

function chebyshev!(x, A, b, λmin::Real, λmax::Real;
n::Int=size(A,2), tol::Real = sqrt(eps(typeof(real(b[1])))),
maxiter::Int = n^3, plot::Bool=false, log::Bool=false, kwargs...
)
K = KrylovSubspace(A, n, 1, Adivtype(A, b))
init!(K, x)
chebyshev_method!(history, x, K, b, λmin, λmax, Pr; tol = tol, maxiter = maxiter)
x, history

(plot & !log) && error("Can't plot when log keyword is false")
history = ConvergenceHistory(partial=!log)
history[:tol] = tol
reserve!(history,:resnorm,maxiter)
chebyshev_method!(history, x, K, b, λmin, λmax; tol=tol, maxiter=maxiter, kwargs...)
plot && (shrink!(history); showplot(history))
log ? (x, history) : x
end

#########################
# Method Implementation #
#########################

function chebyshev_method!(
log::ConvergenceHistory, x, K::KrylovSubspace, b, λmin::Real, λmax::Real,
Pr = 1; tol::Real = sqrt(eps(typeof(real(b[1])))), maxiter::Int = K.n^3
log::ConvergenceHistory, x, K::KrylovSubspace, b, λmin::Real, λmax::Real;
Pr = 1, tol::Real = sqrt(eps(typeof(real(b[1])))), maxiter::Int = K.n^3,
verbose::Bool=false
)

verbose && @printf("=== chebyshev ===\n%4s\t%7s\n","iter","resnorm")
local α, p
K.order = 1
tol = tol*norm(b)
log.mvps=1
r = b - nextvec(K)
d::eltype(b) = (λmax + λmin)/2
c::eltype(b) = (λmax - λmin)/2
Expand All @@ -45,9 +58,85 @@ function chebyshev_method!(
#Check convergence
resnorm = norm(r)
push!(log, :resnorm, resnorm)
verbose && @printf("%3d\t%1.2e\n",iter,resnorm)
resnorm < tol && break
end
shrink!(log)
verbose && @printf("\n")
setconv(log, 0<=norm(r)<tol)
x
end

#################
# Documentation #
#################

let
#Initialize parameters
doc_call = """ chebyshev!(A, b, λmin, λmax)
"""
doc!_call = """ chebyshev!(x, A, b, λmin, λmax)
"""

doc_msg = "Solve A*x=b using the chebyshev method."
doc!_msg = "Overwrite `x`.\n\n" * doc_msg

doc_arg = ""
doc!_arg = """* `x`: initial guess, overwrite final estimation."""

doc_version = (doc_call, doc_msg, doc_arg)
doc!_version = (doc!_call, doc!_msg, doc!_arg)

i=0
docstring = Vector(2)

#Build docs
for (call, msg, arg) in [doc_version, doc!_version] #Start
i+=1
docstring[i] = """
$call
$msg
If `log` is set to `true` is given, method will output a tuple `x, ch`. Where
`ch` is a `ConvergenceHistory` object. Otherwise it will only return `x`.
The `plot` attribute can only be used when `log` is set version.
**Arguments**
$arg
* `A`: linear operator.
* `b`: right hand side.
*Keywords*
* `Pr = 1`: right preconditioner of the method.
* `tol::Real = sqrt(eps())`: stopping tolerance.
* `maxiter::Integer = size(A,2)^3`: maximum number of iterations.
* `verbose::Bool = false`: print method information.
* `log::Bool = false`: output an extra element of type `ConvergenceHistory`
containing extra information of the method execution.
* `plot::Bool = false`: plot data. (Only with `Master` version)
**Output**
*`log` is `false`:*
* `x`: approximated solution.
*`log` is `true`:*
* `x`: approximated solution.
* `ch`: convergence history.
*ConvergenceHistory keys*
* `:tol` => `::Real`: stopping tolerance.
* `:resnom` => `::Vector`: residual norm at each iteration.
"""
end

@doc docstring[1] -> chebyshev
@doc docstring[2] -> chebyshev!
end

0 comments on commit e2a7fc4

Please sign in to comment.