Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Fix type inference #16

Merged
merged 2 commits into from
Feb 24, 2021
Merged

Fix type inference #16

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
2 commits
Select commit Hold shift + click to select a range
429bb52
fix type inference for AbstractTronMatrix
frapac Feb 15, 2021
cc4d452
fix type inference with IncompletePreconditioner
frapac Feb 15, 2021
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
118 changes: 70 additions & 48 deletions src/driver.jl
View file
Open in desktop
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
Open in desktop
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)))
Copy link
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

We need to specify explicitly the signature of the Array here, if we want to compile on Array{Float64, 2} instead of the more generic type Array

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
Open in desktop
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
Open in desktop
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
Open in desktop
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