Skip to content

Commit

Permalink
Simplify API of SparseCallback and DenseCallback (#285)
Browse files Browse the repository at this point in the history 
* make explicit the options passed to `get_index_constraints`
* make explicit the options passed to `create_callback` and `initialize`
* add docstring for the structures in `src/nlpmodels.jl`
  • Loading branch information
@frapac
frapac authored Jan 8, 2024
1 parent 3685435 commit 14d16c2
Show file tree
Hide file tree
Showing 4 changed files with 208 additions and 134 deletions.
51 changes: 28 additions & 23 deletions src/IPM/IPM.jl
View file
Original file line number Diff line number Diff line change
Expand Up @@ -17,7 +17,7 @@ mutable struct MadNLPSolver{
IC <: AbstractInertiaCorrector,
KKTVec <: AbstractKKTVector{T, VT}
} <: AbstractMadNLPSolver{T}

nlp::Model
cb::CB
kkt::KKTSystem
Expand Down Expand Up @@ -103,29 +103,34 @@ mutable struct MadNLPSolver{
end

function MadNLPSolver(nlp::AbstractNLPModel{T,VT}; kwargs...) where {T, VT}

opt, opt_linear_solver, logger = load_options(nlp; kwargs...)
@assert is_supported(opt.linear_solver, T)

cnt = MadNLPCounters(start_time=time())
cb = create_callback(opt.callback, nlp, opt)

cb = create_callback(
opt.callback,
nlp;
fixed_variable_treatment=opt.fixed_variable_treatment,
equality_treatment=opt.equality_treatment,
)

# generic options
opt.disable_garbage_collector &&
(GC.enable(false); @warn(logger,"Julia garbage collector is temporarily disabled"))
set_blas_num_threads(opt.blas_num_threads; permanent=true)
@trace(logger,"Initializing variables.")

ind_cons = get_index_constraints(
get_lvar(nlp), get_uvar(nlp),
get_lcon(nlp), get_ucon(nlp),
opt.fixed_variable_treatment,
opt.equality_treatment
get_lcon(nlp), get_ucon(nlp);
fixed_variable_treatment=opt.fixed_variable_treatment,
equality_treatment=opt.equality_treatment
)

ind_lb = ind_cons.ind_lb
ind_ub = ind_cons.ind_ub

ns = length(ind_cons.ind_ineq)
nx = get_nvar(nlp)
n = nx+ns
Expand All @@ -148,20 +153,20 @@ function MadNLPSolver(nlp::AbstractNLPModel{T,VT}; kwargs...) where {T, VT}

x = PrimalVector(VT, nx, ns, ind_lb, ind_ub)
xl = PrimalVector(VT, nx, ns, ind_lb, ind_ub)
xu = PrimalVector(VT, nx, ns, ind_lb, ind_ub)
xu = PrimalVector(VT, nx, ns, ind_lb, ind_ub)
zl = PrimalVector(VT, nx, ns, ind_lb, ind_ub)
zu = PrimalVector(VT, nx, ns, ind_lb, ind_ub)
f = PrimalVector(VT, nx, ns, ind_lb, ind_ub)
x_trial = PrimalVector(VT, nx, ns, ind_lb, ind_ub)

d = UnreducedKKTVector(VT, n, m, nlb, nub, ind_lb, ind_ub)
p = UnreducedKKTVector(VT, n, m, nlb, nub, ind_lb, ind_ub)
_w1 = UnreducedKKTVector(VT, n, m, nlb, nub, ind_lb, ind_ub)
_w2 = UnreducedKKTVector(VT, n, m, nlb, nub, ind_lb, ind_ub)
_w3 = UnreducedKKTVector(VT, n, m, nlb, nub, ind_lb, ind_ub)
_w4 = UnreducedKKTVector(VT, n, m, nlb, nub, ind_lb, ind_ub)

jacl = VT(undef,n)
jacl = VT(undef,n)
c_trial = VT(undef, m)
y = VT(undef, m)
c = VT(undef, m)
Expand Down Expand Up @@ -190,28 +195,28 @@ function MadNLPSolver(nlp::AbstractNLPModel{T,VT}; kwargs...) where {T, VT}
VT,
n, m, nlb, nub, ind_lb, ind_ub
)

cnt.init_time = time() - cnt.start_time

return MadNLPSolver(
nlp, cb, kkt,
opt, cnt, logger,
opt, cnt, logger,
n, m, nlb, nub,
x, y, zl, zu, xl, xu,
zero(T), f, c,
jacl,
d, p,
_w1, _w2, _w3, _w4,
x_trial, c_trial, zero(T), c_slk, rhs,
ind_cons.ind_ineq, ind_cons.ind_fixed, ind_cons.ind_llb, ind_cons.ind_uub,
x_lr, x_ur, xl_r, xu_r, zl_r, zu_r, dx_lr, dx_ur, x_trial_lr, x_trial_ur,
iterator,
zero(T), f, c,
jacl,
d, p,
_w1, _w2, _w3, _w4,
x_trial, c_trial, zero(T), c_slk, rhs,
ind_cons.ind_ineq, ind_cons.ind_fixed, ind_cons.ind_llb, ind_cons.ind_uub,
x_lr, x_ur, xl_r, xu_r, zl_r, zu_r, dx_lr, dx_ur, x_trial_lr, x_trial_ur,
iterator,
zero(T), zero(T), zero(T), zero(T), zero(T), zero(T), zero(T), zero(T), zero(T),
" ",
zero(T), zero(T), zero(T),
Tuple{T, T}[],
inertia_corrector, nothing,
INITIAL, Dict(),
INITIAL, Dict(),
)

end
Expand Down
68 changes: 35 additions & 33 deletions src/IPM/solver.jl
View file
Original file line number Diff line number Diff line change
Expand Up @@ -15,8 +15,8 @@ function initialize!(solver::AbstractMadNLPSolver{T}) where T

nlp = solver.nlp
opt = solver.opt
# Initializing variables

# Initializing variables
@trace(solver.logger,"Initializing variables.")
initialize!(
solver.cb,
Expand All @@ -25,13 +25,15 @@ function initialize!(solver::AbstractMadNLPSolver{T}) where T
solver.xu,
solver.y,
solver.rhs,
solver.ind_ineq,
opt
solver.ind_ineq;
tol=opt.tol,
bound_push=opt.bound_push,
bound_fac=opt.bound_fac,
)
fill!(solver.jacl, zero(T))
fill!(solver.zl_r, one(T))
fill!(solver.zu_r, one(T))

# Initializing scaling factors
set_scaling!(
solver.cb,
Expand All @@ -50,7 +52,7 @@ function initialize!(solver::AbstractMadNLPSolver{T}) where T
# Initializing jacobian and gradient
eval_jac_wrapper!(solver, solver.kkt, solver.x)
eval_grad_f_wrapper!(solver, solver.f,solver.x)


@trace(solver.logger,"Initializing constraint duals.")
if !solver.opt.dual_initialized
Expand All @@ -65,7 +67,7 @@ function initialize!(solver::AbstractMadNLPSolver{T}) where T
copyto!(solver.y, dual(solver.d))
end
end

# Initializing
solver.obj_val = eval_f_wrapper(solver, solver.x)
eval_cons_wrapper!(solver, solver.c, solver.x)
Expand Down Expand Up @@ -207,7 +209,7 @@ function regular!(solver::AbstractMadNLPSolver{T}) where T
)
solver.inf_compl = get_inf_compl(solver.x_lr,solver.xl_r,solver.zl_r,solver.xu_r,solver.x_ur,solver.zu_r,zero(T),sc)
inf_compl_mu = get_inf_compl(solver.x_lr,solver.xl_r,solver.zl_r,solver.xu_r,solver.x_ur,solver.zu_r,solver.mu,sc)

print_iter(solver)

# evaluate termination criteria
Expand Down Expand Up @@ -244,7 +246,7 @@ function regular!(solver::AbstractMadNLPSolver{T}) where T
dual_inf_perturbation!(primal(solver.p),solver.ind_llb,solver.ind_uub,solver.mu,solver.opt.kappa_d)

inertia_correction!(solver.inertia_corrector, solver) || return ROBUST

# filter start
@trace(solver.logger,"Backtracking line search initiated.")
theta = get_theta(solver.c)
Expand Down Expand Up @@ -278,7 +280,7 @@ function regular!(solver::AbstractMadNLPSolver{T}) where T
unsuccessful_iterate = false

while true

copyto!(full(solver.x_trial), full(solver.x))
axpy!(solver.alpha, primal(solver.d), primal(solver.x_trial))
solver.obj_val_trial = eval_f_wrapper(solver, solver.x_trial)
Expand All @@ -294,7 +296,7 @@ function regular!(solver::AbstractMadNLPSolver{T}) where T
solver.filter,theta,theta_trial,varphi,varphi_trial,switching_condition,armijo_condition,
solver.theta_min,solver.opt.obj_max_inc,solver.opt.gamma_theta,solver.opt.gamma_phi,
has_constraints(solver))

if solver.ftype in ["f","h"]
@trace(solver.logger,"Step accepted with type $(solver.ftype)")
break
Expand All @@ -308,7 +310,7 @@ function regular!(solver::AbstractMadNLPSolver{T}) where T
end
end

unsuccessful_iterate = true
unsuccessful_iterate = true
solver.alpha /= 2
solver.cnt.l += 1
if solver.alpha < alpha_min
Expand All @@ -333,7 +335,7 @@ function regular!(solver::AbstractMadNLPSolver{T}) where T
empty!(solver.filter)
push!(solver.filter,(solver.theta_max,-Inf))
solver.cnt.k+=1

return REGULAR
end
end
Expand Down Expand Up @@ -378,7 +380,7 @@ function regular!(solver::AbstractMadNLPSolver{T}) where T
primal(solver.x),
solver.mu,solver.opt.kappa_sigma,
)

eval_grad_f_wrapper!(solver, solver.f,solver.x)

if !switching_condition || !armijo_condition
Expand Down Expand Up @@ -462,7 +464,7 @@ function restore!(solver::AbstractMadNLPSolver{T}) where T
end

adjust_boundary!(solver.x_lr,solver.xl_r,solver.x_ur,solver.xu_r,solver.mu)

F = F_trial

theta = get_theta(solver.c)
Expand Down Expand Up @@ -561,20 +563,20 @@ function robust!(solver::MadNLPSolver{T}) where T
eval_lag_hess_wrapper!(solver, solver.kkt, solver.x, solver.y; is_resto=true)
end
set_aug_RR!(solver.kkt, solver, RR)

# without inertia correction,
@trace(solver.logger,"Solving restoration phase primal-dual system.")
set_aug_rhs_RR!(solver, solver.kkt, RR, solver.opt.rho)

inertia_correction!(solver.inertia_corrector, solver) || return RESTORATION_FAILED


finish_aug_solve_RR!(
RR.dpp,RR.dnn,RR.dzp,RR.dzn,solver.y,dual(solver.d),
RR.pp,RR.nn,RR.zp,RR.zn,RR.mu_R,solver.opt.rho
)


theta_R = get_theta_R(solver.c,RR.pp,RR.nn)
varphi_R = get_varphi_R(RR.obj_val_R,solver.x_lr,solver.xl_r,solver.xu_r,solver.x_ur,RR.pp,RR.nn,RR.mu_R)
varphi_d_R = get_varphi_d_R(
Expand Down Expand Up @@ -623,7 +625,7 @@ function robust!(solver::MadNLPSolver{T}) where T
varphi_R_trial = get_varphi_R(
RR.obj_val_R_trial,solver.x_trial_lr,solver.xl_r,solver.xu_r,solver.x_trial_ur,RR.pp_trial,RR.nn_trial,RR.mu_R)

armijo_condition = is_armijo(varphi_R_trial,varphi_R,solver.opt.eta_phi,solver.alpha,varphi_d_R)
armijo_condition = is_armijo(varphi_R_trial,varphi_R,solver.opt.eta_phi,solver.alpha,varphi_d_R)

small_search_norm && break
solver.ftype = get_ftype(
Expand All @@ -643,7 +645,7 @@ function robust!(solver::MadNLPSolver{T}) where T
# (experimental) while giving up directly
# we give MadNLP.jl second chance to explore
# some possibility at the current iterate

fill!(solver.y, zero(T))
fill!(solver.zl_r, one(T))
fill!(solver.zu_r, one(T))
Expand Down Expand Up @@ -722,7 +724,7 @@ function robust!(solver::MadNLPSolver{T}) where T
else
copyto!(solver.y, dual(solver.d))
end

solver.cnt.k+=1
solver.cnt.t+=1

Expand Down Expand Up @@ -806,7 +808,7 @@ function inertia_correction!(
inertia_corrector::InertiaBased,
solver::MadNLPSolver{T}
) where {T}

n_trial = 0
solver.del_w = del_w_prev = zero(T)

Expand All @@ -815,14 +817,14 @@ function inertia_correction!(
factorize_wrapper!(solver)

num_pos,num_zero,num_neg = inertia(solver.kkt.linear_solver)


solve_status = !is_inertia_correct(solver.kkt, num_pos, num_zero, num_neg) ?
false : solve_refine_wrapper!(
solver.d, solver, solver.p, solver._w4,
)


while !solve_status
@debug(solver.logger,"Primal-dual perturbed.")

Expand All @@ -837,7 +839,7 @@ function inertia_correction!(
return false
end
end
solver.del_c = num_neg == 0 ? zero(T) : solver.opt.jacobian_regularization_value * solver.mu^(solver.opt.jacobian_regularization_exponent)
solver.del_c = num_neg == 0 ? zero(T) : solver.opt.jacobian_regularization_value * solver.mu^(solver.opt.jacobian_regularization_exponent)
regularize_diagonal!(solver.kkt, solver.del_w - del_w_prev, solver.del_c)
del_w_prev = solver.del_w

Expand All @@ -850,15 +852,15 @@ function inertia_correction!(
)
n_trial += 1
end

solver.del_w != 0 && (solver.del_w_last = solver.del_w)
return true
end

function inertia_correction!(
inertia_corrector::InertiaFree,
solver::MadNLPSolver{T}
) where T
) where T

n_trial = 0
solver.del_w = del_w_prev = zero(T)
Expand Down Expand Up @@ -922,7 +924,7 @@ function inertia_correction!(
inertia_corrector::InertiaIgnore,
solver::MadNLPSolver{T}
) where T

n_trial = 0
solver.del_w = del_w_prev = zero(T)

Expand All @@ -946,7 +948,7 @@ function inertia_correction!(
return false
end
end
solver.del_c = solver.opt.jacobian_regularization_value * solver.mu^(solver.opt.jacobian_regularization_exponent)
solver.del_c = solver.opt.jacobian_regularization_value * solver.mu^(solver.opt.jacobian_regularization_exponent)
regularize_diagonal!(solver.kkt, solver.del_w - del_w_prev, solver.del_c)
del_w_prev = solver.del_w

Expand Down
Loading

0 comments on commit 14d16c2

Please sign in to comment.