Skip to content

Commit

Permalink
Merge pull request #16 from exanauts/fp/type_inference
Browse files Browse the repository at this point in the history 
Fix type inference
  • Loading branch information
@youngdae
youngdae authored Feb 24, 2021
2 parents c264f4b + cc4d452 commit 557408e
Show file tree
Hide file tree
Showing 5 changed files with 79 additions and 56 deletions.
118 changes: 70 additions & 48 deletions src/driver.jl
View file
Original file line number Diff line number Diff line change
@@ -1,10 +1,10 @@
mutable struct ExaTronProblem{VI, VD}
mutable struct ExaTronProblem{VI, VD, TM <: AbstractTronMatrix}
n::Cint # number of variables
nnz::Integer # number of Hessian entries
nnz_a::Integer # number of Hessian entries in the strict lower
A::AbstractTronMatrix
B::AbstractTronMatrix
precond::AbstractPreconditioner
nnz::Int # number of Hessian entries
nnz_a::Int # number of Hessian entries in the strict lower
A::TM
B::TM
L::TM
indfree::VI # a working array of dimension n
iwa::VI # a working array of dimension 3*n
g::VD # gradient
Expand All @@ -30,30 +30,34 @@ mutable struct ExaTronProblem{VI, VD}
f::Cdouble
gnorm_two::Float64
gnorm_inf::Float64
nfev::Integer
ngev::Integer
nhev::Integer
minor_iter::Integer
nfev::Int
ngev::Int
nhev::Int
minor_iter::Int
status::Symbol

ExaTronProblem{VI, VD}() where {VI, VD} = new{VI, VD}()
ExaTronProblem{VI, VD, TM}() where {VI, VD, TM} = new{VI, VD, TM}()
end


function set_default_options!(prob::ExaTronProblem)
prob.options = Dict{String,Any}()
set_default_options!(prob.options)
end

prob.options["max_feval"] = 500
prob.options["max_minor"] = 200
prob.options["p"] = 5
prob.options["verbose"] = 0
prob.options["tol"] = 1e-6
prob.options["fatol"] = 0
prob.options["frtol"] = 1e-12
prob.options["fmin"] = -1e32
prob.options["cgtol"] = 0.1
prob.options["tron_code"] = :Julia
prob.options["matrix_type"] = :Sparse
function set_default_options!(options::Dict{String, Any})
options["max_feval"] = 500
options["max_minor"] = 200
options["p"] = 5
options["verbose"] = 0
options["tol"] = 1e-6
options["fatol"] = 0.0
options["frtol"] = 1e-12
options["fmin"] = -1e32
options["cgtol"] = 0.1
options["tron_code"] = :Julia
options["matrix_type"] = :Sparse
options["cg_precond"] = :ICFS

return
end
Expand Down Expand Up @@ -112,14 +116,24 @@ function createProblem(n::Integer, x_l::AbstractVector{Float64}, x_u::AbstractVe
@assert n == length(x_l) == length(x_u)
@assert typeof(x_l) == typeof(x_u)

# Load options first, as we need to know :matrix_type for type inferrence
options_dict = Dict{String, Any}()
set_default_options!(options_dict)
for (name, value) in options
options_dict[string(name)] = value
end

VI = Vector{Cint}
VD = typeof(x_l)

tron = ExaTronProblem{VI, VD}()
set_default_options!(tron)
for (name, value) in options
setOption(tron, string(name), value)
if options_dict["matrix_type"] == :Sparse
TM = TronSparseMatrixCSC{VI, VD}
else
AT = isa(x_l, Array) ? Array{Float64, 2} : CuArray{Float64, 2}
TM = TronDenseMatrix{AT}
end

tron = ExaTronProblem{VI, VD, TM}()
tron.options = options_dict
instantiate_memory!(tron, n, Int64(nele_hess))
copyto!(tron.x_l, 1, x_l, 1, n)
copyto!(tron.x_u, 1, x_u, 1, n)
Expand All @@ -135,19 +149,17 @@ function createProblem(n::Integer, x_l::AbstractVector{Float64}, x_u::AbstractVe
if tron.options["matrix_type"] == :Sparse
tron.A = TronSparseMatrixCSC(tron.rows, tron.cols, tron.values, n)
tron.B = TronSparseMatrixCSC{VI, VD}(n, Int64(nele_hess))
L = TronSparseMatrixCSC{VI, VD}(n, Int64(nele_hess + n*p))
tron.L = TronSparseMatrixCSC{VI, VD}(n, Int64(nele_hess + n*p))
tron.nnz_a = tron.A.nnz
tron.precond = IncompleteCholesky(L, p, 5)
else
tron.A = TronDenseMatrix(tron.rows, tron.cols, tron.values, n)
if isa(x_l, Array)
tron.B = TronDenseMatrix{Array}(n)
L = TronDenseMatrix{Array}(n)
tron.B = TronDenseMatrix{Array{Float64, 2}}(n)
tron.L = TronDenseMatrix{Array{Float64, 2}}(n)
else
tron.B = TronDenseMatrix{CuArray}(n)
L = TronDenseMatrix{CuArray}(n)
tron.B = TronDenseMatrix{CuArray{Float64, 2}}(n)
tron.L = TronDenseMatrix{CuArray{Float64, 2}}(n)
end
tron.precond = IncompleteCholesky(L, p, 5)
tron.nnz_a = n*n
end
tron.status = :NotSolved
Expand All @@ -166,22 +178,32 @@ function solveProblem(tron::ExaTronProblem)

isave = VI(undef, 3)
dsave = tron_zeros(VD, 3)
fatol = tron.options["fatol"]
frtol = tron.options["frtol"]
fmin = tron.options["fmin"]
cgtol = tron.options["cgtol"]
gtol = tron.options["tol"]
max_feval = tron.options["max_feval"]
tcode = tron.options["tron_code"]
max_minor = tron.options["max_minor"]
fatol = tron.options["fatol"]::Float64
frtol = tron.options["frtol"]::Float64
fmin = tron.options["fmin"]::Float64
cgtol = tron.options["cgtol"]::Float64
gtol = tron.options["tol"]::Float64
max_feval = tron.options["max_feval"]::Int
tcode = tron.options["tron_code"]::Symbol
max_minor = tron.options["max_minor"]::Int

# Build preconditioner
if tron.options["cg_precond"] == :ICFS
p = tron.options["p"]::Int
precond = IncompleteCholesky(tron.L, p, 5)
elseif tron.options["cg_precond"] == :Eye
precond = EyePreconditioner()
else
error("Wrong preconditioner. Currently only `:ICFS` and `:Eye` are supported.")
end


# Sanity check
if (tcode == :Fortran) && !isa(tron.precond, IncompleteCholesky)
if (tcode == :Fortran) && !isa(precond, IncompleteCholesky)
error("Legacy Tron supports only IncompleteCholesky preconditioner." *
"Please change preconditioner or set `tron_code` option to `:Julia`")
end


# Project x into its bound. Tron requires x to be feasible.
tron.x .= max.(tron.x_l, min.(tron.x_u, tron.x))

Expand All @@ -196,7 +218,7 @@ function solveProblem(tron::ExaTronProblem)
# Evaluate the function.

if Char(task[1]) == 'F' || unsafe_string(pointer(task), 5) == "START"
tron.f = tron.eval_f(tron.x)
tron.f = tron.eval_f(tron.x)::Cdouble
tron.nfev += 1
if tron.nfev >= max_feval
tron.status = :Maximum_Iterations_Exceeded
Expand Down Expand Up @@ -246,16 +268,16 @@ function solveProblem(tron::ExaTronProblem)
tron.A.colptr, tron.A.rowval, frtol, fatol,
fmin, cgtol, itermax, delta,
task, tron.B.tril_vals, tron.B.diag_vals, tron.B.colptr,
tron.B.rowval, tron.precond.L.tril_vals, tron.precond.L.diag_vals, tron.precond.L.colptr,
tron.precond.L.rowval, tron.xc, tron.s, tron.indfree,
tron.B.rowval, tron.L.tril_vals, tron.L.diag_vals, tron.L.colptr,
tron.L.rowval, tron.xc, tron.s, tron.indfree,
isave, dsave, tron.wa, tron.iwa)
tron.delta = delta[]
else
tron.delta = ExaTron.dtron(tron.n, tron.x, tron.x_l, tron.x_u,
tron.f, tron.g, tron.A,
frtol, fatol,
fmin, cgtol, itermax, tron.delta,
task, tron.B, tron.precond,
task, tron.B, precond,
tron.xc, tron.s, tron.indfree,
isave, dsave, tron.wa, tron.iwa)
end
Expand Down
7 changes: 4 additions & 3 deletions src/dssyax.jl
View file
Original file line number Diff line number Diff line change
Expand Up @@ -20,7 +20,8 @@ end

function TronDenseMatrix{MD}(n::Int) where {MD}
@assert n >= 1
return TronDenseMatrix{MD}(n, n, tron_zeros(MD{Float64}, (n, n)))
vals = tron_zeros(MD, (n, n))
return TronDenseMatrix{MD}(n, n, vals)
end

function TronDenseMatrix(A::TronDenseMatrix)
Expand All @@ -32,9 +33,9 @@ function TronDenseMatrix(I::VI, J::VI, V::VD, n) where {VI, VD}
@assert length(I) == length(J) == length(V)

if isa(V, Array)
A = TronDenseMatrix{Array}(n, n, tron_zeros(Array{eltype(V)}, (n, n)))
A = TronDenseMatrix{Array{Float64, 2}}(n, n, tron_zeros(Array{eltype(V)}, (n, n)))
else
A = TronDenseMatrix{CuArray}(n, n, tron_zeros(CuArray{eltype(V)}, (n, n)))
A = TronDenseMatrix{CuArray{Float64, 2}}(n, n, tron_zeros(CuArray{eltype(V)}, (n, n)))
end
for i=1:length(I)
@assert 1 <= I[i] <= n && 1 <= J[i] <= n && I[i] >= J[i]
Expand Down
4 changes: 2 additions & 2 deletions src/preconditioners.jl
View file
Original file line number Diff line number Diff line change
Expand Up @@ -49,8 +49,8 @@ The incomplete factorization is computed when calling the function
`update!`, and the factorization stored inside the object `IncompleteCholesky`.
"""
struct IncompleteCholesky <: AbstractPreconditioner
L::AbstractTronMatrix
struct IncompleteCholesky{TM <: AbstractTronMatrix} <: AbstractPreconditioner
L::TM
memory::Int
nv::Int
end
Expand Down
4 changes: 2 additions & 2 deletions test/densetest.jl
View file
Original file line number Diff line number Diff line change
Expand Up @@ -5,7 +5,7 @@ import ExaTron: ICFS
Random.seed!(0)
for k=1:10
n = 10
A = ExaTron.TronDenseMatrix{Array}(n)
A = ExaTron.TronDenseMatrix{Array{Float64, 2}}(n)
B = zeros(n,n)
for j=1:n, i=j:n
v = rand(1)[1]
Expand Down Expand Up @@ -36,7 +36,7 @@ import ExaTron: ICFS
for j=1:nfree
iwa[indfree[j]] = j
end
C = ExaTron.TronDenseMatrix{Array}(n)
C = ExaTron.TronDenseMatrix{Array{Float64, 2}}(n)
ExaTron.reorder!(C, A, indfree, nfree, iwa)
@test C.vals[1:nfree,1:nfree] == A.vals[indfree,indfree]

Expand Down
2 changes: 1 addition & 1 deletion test/qptest.jl
View file
Original file line number Diff line number Diff line change
Expand Up @@ -66,7 +66,7 @@ end
prob = build_problem(; n=n)
@testset "Tron: Julia" begin
prob.x .= 0.5 .* (prob.x_l .+ prob.x_u)
prob.precond = ExaTron.EyePreconditioner()
prob.options["cg_precond"] = :Eye
ExaTron.solveProblem(prob)
@test prob.f ≈ obj♯ atol=1e-8
end
Expand Down

0 comments on commit 557408e

Please sign in to comment.