Skip to content

Add option for Bk in sparse format #313

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

Draft
MaxenceGollier wants to merge 3 commits into
base: master
Choose a base branch
from
Draft

Add option for Bk in sparse format #313

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
3 commits
Select commit Hold shift + click to select a range
33a583d
add sparse hessian
MaxenceGollier Feb 13, 2026
af57c36
add sparse matrix option
MaxenceGollier Mar 11, 2026
53aa80f
fix ci
MaxenceGollier Mar 11, 2026
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
6 changes: 5 additions & 1 deletion src/R2N.jl
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -35,6 +35,7 @@ function R2NSolver(
reg_nlp::AbstractRegularizedNLPModel{T, V};
subsolver = R2Solver,
m_monotone::Int = 1,
sparse = false
) where {T, V}
x0 = reg_nlp.model.meta.x0
l_bound = reg_nlp.model.meta.lvar
Expand Down Expand Up @@ -69,7 +70,8 @@ function R2NSolver(
has_bnds ? shifted(reg_nlp.h, xk, l_bound_m_x, u_bound_m_x, reg_nlp.selected) :
shifted(reg_nlp.h, xk)

Bk = hess_op(reg_nlp, xk)
sparse = isa(reg_nlp.model, QuasiNewtonModel) ? false : sparse
Bk = sparse ? hess(reg_nlp, xk) : hess_op(reg_nlp, xk)
Copy link
Member

@dpo dpo Mar 11, 2026

Choose a reason for hiding this comment

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

Why would we want to evaluate the Hessian, even if it is available and sparse?

Copy link
Collaborator Author

@MaxenceGollier MaxenceGollier Mar 11, 2026

Choose a reason for hiding this comment

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

For my penalty method, for example :)

Copy link
Member

@dpo dpo Mar 11, 2026

Choose a reason for hiding this comment

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

Ok, but generally speaking, we don't want to evaluate the Hessian explicitly. Yours is a special case. It feels like we need a level of abstraction here.

Copy link
Collaborator Author

@MaxenceGollier MaxenceGollier Mar 11, 2026

Choose a reason for hiding this comment

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

I could dispatch hess_op to hess for my penalized subproblem type, though i don't see how i could justify this

sub_nlp = QuadraticModel(∇fk, Bk, σ = T(1), x0 = x0)
Comment on lines +73 to 75
Copy link

Copilot AI Mar 11, 2026

Choose a reason for hiding this comment

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

When sparse=true, Bk = hess(reg_nlp, xk) is evaluated while xk has only been allocated (via similar(x0)) and hasn’t been initialized yet. This will compute the initial Hessian at undefined values, and subsequent computations (e.g. initial operator norm / subproblem solve) may use incorrect data.

To fix the sparse-matrix path, initialize xk to a valid point before calling hess, or compute/refresh the Hessian at the start of solve! after xk is set from the x keyword, and before any use of solver.subpb.model.data.H.

Copilot uses AI. Check for mistakes.
Comment on lines +73 to 75
Copy link

Copilot AI Mar 11, 2026

Choose a reason for hiding this comment

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

This introduces a new sparse-Hessian path (sparse ? hess(...) : hess_op(...)) and an in-place Hessian refresh using hess_coord!, but the existing R2N tests don’t appear to exercise sparse=true.

Given that the repo already has automated tests covering R2N behavior, please add a test that runs R2N(...; sparse=true) (with a model that yields a sparse Hessian) to ensure the option is forwarded and the Hessian update path doesn’t error (e.g. on .nzval) and produces consistent results.

Copilot uses AI. Check for mistakes.
subpb = RegularizedNLPModel(sub_nlp, ψ)
substats = RegularizedExecutionStats(subpb)
Expand Down Expand Up @@ -464,6 +466,8 @@ function SolverCore.solve!(
qn_update_y!(nlp, solver, stats)
push!(nlp, s, solver.y)
qn_copy!(nlp, solver, stats)
elseif isa(solver.subpb.model.data.H, AbstractMatrix)
hess_coord!(reg_nlp, xk, solver.subpb.model.data.H.nzval)
Comment on lines +469 to +470
Copy link

Copilot AI Mar 11, 2026

Choose a reason for hiding this comment

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

isa(solver.subpb.model.data.H, AbstractMatrix) is not sufficient for safely using .nzval: dense matrices and many AbstractMatrix wrappers do not have an nzval field. If hess(reg_nlp, xk) returns a dense matrix, this branch will throw.

If you only intend to support in-place updates for sparse CSC Hessians, narrow the type check to the supported sparse type(s) and use a safe accessor for the stored values (or handle dense matrices separately by updating the full matrix).

Copilot uses AI. Check for mistakes.
end

if opnorm_maxiter ≤ 0
Expand Down
6 changes: 5 additions & 1 deletion src/TR_alg.jl
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -36,6 +36,7 @@ function TRSolver(
χ::X = NormLinf(one(T)),
subsolver = R2Solver,
m_monotone::Int = 1,
sparse::Bool = false,
) where {T, V, X}
x0 = reg_nlp.model.meta.x0
l_bound = reg_nlp.model.meta.lvar
Expand Down Expand Up @@ -69,7 +70,8 @@ function TRSolver(
shifted(reg_nlp.h, xk, l_bound_m_x, u_bound_m_x, reg_nlp.selected) :
shifted(reg_nlp.h, xk, T(1), χ)

Bk = hess_op(reg_nlp, xk)
sparse = isa(reg_nlp.model, QuasiNewtonModel) ? false : sparse
Bk = sparse ? hess(reg_nlp, xk) : hess_op(reg_nlp, xk)
sub_nlp = QuadraticModel(∇fk, Bk, x0 = x0)
subpb = RegularizedNLPModel(sub_nlp, ψ)
substats = RegularizedExecutionStats(subpb)
Expand Down Expand Up @@ -468,6 +470,8 @@ function SolverCore.solve!(
if quasiNewtTest
@. ∇fk⁻ = ∇fk - ∇fk⁻
push!(nlp, s, ∇fk⁻) # update QN operator
elseif isa(solver.subpb.model.data.H, AbstractMatrix)
hess_coord!(reg_nlp, xk, solver.subpb.model.data.H.nzval)
Comment on lines +473 to +474
Copy link

Copilot AI Mar 11, 2026

Choose a reason for hiding this comment

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

The guard isa(solver.subpb.model.data.H, AbstractMatrix) is too broad for accessing .nzval: dense matrices (and many AbstractMatrix wrappers like Symmetric) do not have an nzval field, so this branch will throw if hess(reg_nlp, xk) returns a dense matrix.

If the intent is to update a sparse CSC Hessian in-place, restrict the check to the specific sparse type(s) you support (e.g. SparseMatrixCSC / AbstractSparseMatrix) and use a safe way to access values (e.g. nonzeros(H)), and handle the dense-matrix case by updating via a dense Hessian fill/copy instead of .nzval.

Copilot uses AI. Check for mistakes.
end

if opnorm_maxiter ≤ 0
Expand Down
Loading