From f428189ff5f0ca5e3bb4a2f3cbb4abdfa8b8a55f Mon Sep 17 00:00:00 2001 From: arnavk23 Date: Tue, 17 Mar 2026 06:33:26 +0530 Subject: [PATCH 01/11] Add hs85 and hs89 problems --- Project.toml | 2 + src/ADNLPProblems/hs85.jl | 164 ++++++++++++++++++++++++++++++++++++++ src/ADNLPProblems/hs89.jl | 86 ++++++++++++++++++++ src/Meta/hs85.jl | 23 ++++++ src/Meta/hs89.jl | 23 ++++++ src/PureJuMP/hs85.jl | 109 +++++++++++++++++++++++++ src/PureJuMP/hs89.jl | 95 ++++++++++++++++++++++ 7 files changed, 502 insertions(+) create mode 100644 src/ADNLPProblems/hs85.jl create mode 100644 src/ADNLPProblems/hs89.jl create mode 100644 src/Meta/hs85.jl create mode 100644 src/Meta/hs89.jl create mode 100644 src/PureJuMP/hs85.jl create mode 100644 src/PureJuMP/hs89.jl diff --git a/Project.toml b/Project.toml index f69a86a54..c2a178e86 100644 --- a/Project.toml +++ b/Project.toml @@ -10,6 +10,7 @@ LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e" Requires = "ae029012-a4dd-5104-9daa-d747884805df" SparseArrays = "2f01184e-e22b-5df5-ae63-d93ebab69eaf" SpecialFunctions = "276daf66-3868-5448-9aa4-cd146d93841b" +StaticArrays = "90137ffa-7385-5640-81b9-e52037218182" [compat] DataFrames = "1" @@ -17,6 +18,7 @@ JLD2 = "0.5, 0.6" JuMP = "^1.15" Requires = "1" SpecialFunctions = "2" +StaticArrays = "1.9.18" julia = "1.6" [extras] diff --git a/src/ADNLPProblems/hs85.jl b/src/ADNLPProblems/hs85.jl new file mode 100644 index 000000000..01644ead1 --- /dev/null +++ b/src/ADNLPProblems/hs85.jl @@ -0,0 +1,164 @@ +export hs85 + +function hs85(; type::Type{T}=Float64, kwargs...) where {T} + a = T[0, 17.505, 11.275, 214.228, 7.458, 0.961, 1.612, 0.146, 107.99, + 922.693, 926.832, 18.766, 1072.163, 8961.448, 0.063, 71084.33, 2802713] + b = T[0, 1053.6667, 35.03, 665.585, 584.463, 265.916, 7.046, 0.222, + 273.366, 1286.105, 1444.046, 537.141, 3247.039, 26844.086, 0.386, + 140000, 12146108] + c10 = T(12.3) / T(752.3) + + # Variable bounds and starting point (exact from the standard model) + lvar = T[704.4148, 68.6, 0.0, 193.0, 25.0] + uvar = T[906.3855, 288.88, 134.75, 287.0966, 84.1988] + x0 = T[900.0, 80.0, 115.0, 267.0, 27.0] + + # Best known value ≈ -1.90513375 + function f(x::AbstractVector{T}) + # All intermediates (identical to those used in constraints) + y1 = x[2] + x[3] + T(41.6) + c1 = T(0.024) * x[4] - T(4.62) + y2 = T(12.5) / c1 + T(12) + c2 = T(0.0003535) * x[1]^2 + T(0.5311) * x[1] + T(0.08705) * y2 * x[1] + c3 = T(0.052) * x[1] + T(78) + T(0.002377) * y2 * x[1] + y3 = c2 / c3 + y4 = T(19) * y3 + c4 = T(0.04782) * (x[1] - y3) + T(0.1956) * (x[1] - y3)^2 / x[2] + + T(0.6376) * y4 + T(1.594) * y3 + c5 = T(100) * x[2] + c6 = x[1] - y3 - y4 + c7 = T(0.95) - c4 / c5 + y5 = c6 * c7 + y6 = x[1] - y5 - y4 - y3 + c8 = (y5 + y4) * T(0.995) + y7 = c8 / y1 + y8 = c8 / T(3798) + c9 = y7 - T(0.0663) * y7 / y8 - T(0.3153) + y9 = T(96.82) / c9 + T(0.321) * y1 + y10 = T(1.29) * y5 + T(1.258) * y4 + T(2.29) * y3 + T(1.71) * y6 + y11 = T(1.71) * x[1] - T(0.452) * y4 + T(0.58) * y3 + c11 = T(1.75) * y2 * T(0.995) * x[1] + c12 = T(0.995) * y10 + T(1998) + y12 = c10 * x[1] + c11 / c12 + y13 = c12 - T(1.75) * y2 + y14 = T(3623) + T(64.4) * x[2] + T(58.4) * x[3] + T(146312) / (y9 + x[5]) + c13 = T(0.995) * y10 + T(60.8) * x[2] + T(48) * x[4] - + T(0.1121) * y14 - T(5095) + y15 = y13 / c13 + y16 = T(148000) - T(331000) * y15 + T(40) * y13 - T(61) * y15 * y13 + c14 = T(2324) * y10 - T(28740000) * y2 + y17 = T(14130000) - T(1328) * y10 - T(531) * y11 + c14 / c12 + c15 = y13 / y15 - y13 / T(0.52) + c16 = T(1.104) - T(0.72) * y15 + # c17 not needed for objective + + return -T(5.843e-7) * y17 + + T(1.17e-4) * y14 + + T(2.358e-5) * y13 + + T(1.502e-6) * y16 + + T(0.0321) * y12 + + T(0.004324) * y5 + + T(1e-4) * c15 / c16 + + T(37.48) * y2 / c12 + + T(0.1365) + end + + # Constraint function (48 inequalities, all of the form c(x) >= 0) + function c!(cx::AbstractVector{T}, x::AbstractVector{T}) + # All intermediates (identical to those used in objective) + y1 = x[2] + x[3] + T(41.6) + c1 = T(0.024) * x[4] - T(4.62) + y2 = T(12.5) / c1 + T(12) + c2 = T(0.0003535) * x[1]^2 + T(0.5311) * x[1] + T(0.08705) * y2 * x[1] + c3 = T(0.052) * x[1] + T(78) + T(0.002377) * y2 * x[1] + y3 = c2 / c3 + y4 = T(19) * y3 + c4 = T(0.04782) * (x[1] - y3) + T(0.1956) * (x[1] - y3)^2 / x[2] + T(0.6376) * y4 + T(1.594) * y3 + c5 = T(100) * x[2] + c6 = x[1] - y3 - y4 + c7 = T(0.95) - c4 / c5 + y5 = c6 * c7 + y6 = x[1] - y5 - y4 - y3 + c8 = (y5 + y4) * T(0.995) + y7 = c8 / y1 + y8 = c8 / T(3798) + c9 = y7 - T(0.0663) * y7 / y8 - T(0.3153) + y9 = T(96.82) / c9 + T(0.321) * y1 + y10 = T(1.29) * y5 + T(1.258) * y4 + T(2.29) * y3 + T(1.71) * y6 + y11 = T(1.71) * x[1] - T(0.452) * y4 + T(0.58) * y3 + c11 = T(1.75) * y2 * T(0.995) * x[1] + c12 = T(0.995) * y10 + T(1998) + y12 = c10 * x[1] + c11 / c12 + y13 = c12 - T(1.75) * y2 + y14 = T(3623) + T(64.4) * x[2] + T(58.4) * x[3] + T(146312) / (y9 + x[5]) + c13 = T(0.995) * y10 + T(60.8) * x[2] + T(48) * x[4] - T(0.1121) * y14 - T(5095) + y15 = y13 / c13 + y16 = T(148000) - T(331000) * y15 + T(40) * y13 - T(61) * y15 * y13 + c14 = T(2324) * y10 - T(28740000) * y2 + y17 = T(14130000) - T(1328) * y10 - T(531) * y11 + c14 / c12 + c15 = y13 / y15 - y13 / T(0.52) + c16 = T(1.104) - T(0.72) * y15 + c17 = y9 + x[5] + # Constraints + cx[1] = T(1.5) * x[2] - x[3] + cx[2] = y1 - T(213.1) + cx[3] = T(405.23) - y1 + cx[4] = y2 - a[2] + cx[5] = y3 - a[3] + cx[6] = y4 - a[4] + cx[7] = y5 - a[5] + cx[8] = y6 - a[6] + cx[9] = y7 - a[7] + cx[10] = y8 - a[8] + cx[11] = y9 - a[9] + cx[12] = y10 - a[10] + cx[13] = y11 - a[11] + cx[14] = y12 - a[12] + cx[15] = y13 - a[13] + cx[16] = y14 - a[14] + cx[17] = y15 - a[15] + cx[18] = y16 - a[16] + cx[19] = y17 - a[17] + cx[20] = b[2] - y2 + cx[21] = b[3] - y3 + cx[22] = b[4] - y4 + cx[23] = b[5] - y5 + cx[24] = b[6] - y6 + cx[25] = b[7] - y7 + cx[26] = b[8] - y8 + cx[27] = b[9] - y9 + cx[28] = b[10] - y10 + cx[29] = b[11] - y11 + cx[30] = b[12] - y12 + cx[31] = b[13] - y13 + cx[32] = b[14] - y14 + cx[33] = b[15] - y15 + cx[34] = b[16] - y16 + cx[35] = b[17] - y17 + cx[36] = y4 - (T(0.28) / T(0.72)) * y5 + cx[37] = T(21) - T(3496) * y2 / c12 + cx[38] = T(62212) / c17 - T(110.6) - y1 + # 10 box constraints (5 lower, 5 upper) + cx[39] = x[1] - lvar[1] + cx[40] = x[2] - lvar[2] + cx[41] = x[3] - lvar[3] + cx[42] = x[4] - lvar[4] + cx[43] = x[5] - lvar[5] + cx[44] = uvar[1] - x[1] + cx[45] = uvar[2] - x[2] + cx[46] = uvar[3] - x[3] + cx[47] = uvar[4] - x[4] + cx[48] = uvar[5] - x[5] + return cx + end + # Constraint bounds: all inequalities of the form c(x) >= 0 + m = 48 + cl = zeros(T, m) + cu = fill(T(Inf), m) + return ADNLPModels.ADNLPModel!( + f, x0, lvar, uvar, + c!, cl, cu; + name = "hs85", + kwargs... + ) +end \ No newline at end of file diff --git a/src/ADNLPProblems/hs89.jl b/src/ADNLPProblems/hs89.jl new file mode 100644 index 000000000..d6ec668b7 --- /dev/null +++ b/src/ADNLPProblems/hs89.jl @@ -0,0 +1,86 @@ +export hs89 + +function hs89(; type::Type{T} = Float64, kwargs...) where T + # First 30 positive roots of tan(μ) = μ + mu = T[ + 0.8603335890193798, 3.425618459481728, 6.437298179171945, 9.529334405361963, + 12.645287223856588, 15.771284874815820, 18.902409956860000, 22.036496727938500, + 25.172446326646600, 28.309642854452000, 31.447714637546200, 34.586424215288900, + 37.725612827776500, 40.865170330488000, 44.005017920830800, 47.145097736761000, + 50.285366337773600, 53.425790477394600, 56.566344279821500, 59.707007305335400, + 62.847763194454400, 65.988598698490300, 69.129502973895200, 72.270467060308900, + 75.411483488848100, 78.552545984242900, 81.693649235601600, 84.834788718042200, + 87.975960552493200, 91.117161394464700 + ] + + # Precomputed coefficients A_j = 2 sin(μ_j) / (μ_j + sin(μ_j) cos(μ_j)) + A = [2 * sin(mu[j]) / (mu[j] + sin(mu[j]) * cos(mu[j])) for j in 1:30] + + # Objective: φ(x) = ∑_{j=1}^{30} A_j ρ_j(x) + function f(x::AbstractVector{T}) + s = zero(T) + r = x[1]^2 + x[2]^2 + x[3]^2 + for j = 1:30 + μ² = mu[j]^2 + exp_r = exp(-μ² * r) + ρ = - (exp_r + 2*exp(-μ²*(x[2]^2 + x[3]^2)) + 2*exp(-μ²*x[3]^2) + 1) / μ² + s += A[j] * ρ + end + return s + end + + # Equality constraint c(x) = 0 + # Full expression with cross terms (double sum over i < j) + function c!(cx::AbstractVector{T}, x::AbstractVector{T}) + r = x[1]^2 + x[2]^2 + x[3]^2 + termA = zero(T) + termB = zero(T) + # Compute termA and termB directly, no heap allocation + for j = 1:30 + μ = mu[j] + μ² = μ^2 + exp_r = exp(-μ² * r) + exp_r23 = exp(-μ² * (x[2]^2 + x[3]^2)) + exp_r3 = exp(-μ² * x[3]^2) + ρj = - (exp_r + 2*exp_r23 + 2*exp_r3 + 1) / μ² + termA += A[j]^2 * ρj^2 * (sin(2*μ)/(2*μ) + one(T)) / 2 + end + for i = 1:29 + μi = mu[i] + μi² = μi^2 + exp_ri = exp(-μi² * r) + exp_ri23 = exp(-μi² * (x[2]^2 + x[3]^2)) + exp_ri3 = exp(-μi² * x[3]^2) + ρi = - (exp_ri + 2*exp_ri23 + 2*exp_ri3 + 1) / μi² + for j = i+1:30 + μj = mu[j] + μj² = μj^2 + exp_rj = exp(-μj² * r) + exp_rj23 = exp(-μj² * (x[2]^2 + x[3]^2)) + exp_rj3 = exp(-μj² * x[3]^2) + ρj = - (exp_rj + 2*exp_rj23 + 2*exp_rj3 + 1) / μj² + denom_plus = μi + μj + denom_minus = μi - μj + sin_plus = iszero(denom_plus) ? one(T) : sin(denom_plus)/denom_plus + sin_minus = iszero(denom_minus) ? one(T) : sin(denom_minus)/denom_minus + termB += A[i] * A[j] * ρi * ρj * (sin_plus + sin_minus) + end + end + cx[1] = termA + termB - T(2)/15 + return cx + end + + # Starting point (common in literature / CUTE) + x0 = T[0.5, -0.5, 0.5] + + # One equality constraint c(x) = 0 + lcon = ucon = T[0] + + return ADNLPModels.ADNLPModel!( + f, x0, + c!, lcon, ucon; + name = "hs89", + lin = Int[], # no linear constraints + kwargs... + ) +end \ No newline at end of file diff --git a/src/Meta/hs85.jl b/src/Meta/hs85.jl new file mode 100644 index 000000000..eba8e3197 --- /dev/null +++ b/src/Meta/hs85.jl @@ -0,0 +1,23 @@ +hs85_meta = Dict( + :nvar => 5, + :variable_nvar => false, + :ncon => 48, + :variable_ncon => false, + :minimize => true, + :name => "hs85", + :has_equalities_only => false, + :has_inequalities_only => true, + :has_bounds => true, + :has_fixed_variables => false, + :objtype => :other, + :contype => :general, + :best_known_lower_bound => -Inf, + :best_known_upper_bound => -1.90513375, + :is_feasible => true, + :defined_everywhere => missing, + :origin => :unknown, +) +get_hs85_nvar(; n::Integer = 5, kwargs...) = 5 +get_hs85_ncon(; n::Integer = 48, kwargs...) = 48 +get_hs85_nlin(; n::Integer = 48, kwargs...) = 0 +get_hs85_nnln(; n::Integer = 48, kwargs...) = 48 diff --git a/src/Meta/hs89.jl b/src/Meta/hs89.jl new file mode 100644 index 000000000..13ac47fcb --- /dev/null +++ b/src/Meta/hs89.jl @@ -0,0 +1,23 @@ +hs89_meta = Dict( + :nvar => 3, + :variable_nvar => false, + :ncon => 1, + :variable_ncon => false, + :minimize => true, + :name => "hs89", + :has_equalities_only => false, + :has_inequalities_only => true, + :has_bounds => false, + :has_fixed_variables => false, + :objtype => :other, + :contype => :general, + :best_known_lower_bound => -Inf, + :best_known_upper_bound => 1.36265681, + :is_feasible => true, + :defined_everywhere => missing, + :origin => :unknown, +) +get_hs89_nvar(; n::Integer = 3, kwargs...) = 3 +get_hs89_ncon(; n::Integer = 1, kwargs...) = 1 +get_hs89_nlin(; n::Integer = 1, kwargs...) = 0 +get_hs89_nnln(; n::Integer = 1, kwargs...) = 1 diff --git a/src/PureJuMP/hs85.jl b/src/PureJuMP/hs85.jl new file mode 100644 index 000000000..a07f2df6d --- /dev/null +++ b/src/PureJuMP/hs85.jl @@ -0,0 +1,109 @@ +## Hock and Schittkowski problem number 85 +# +# Source: +# Problem 85 in +# W. Hock and K. Schittkowski, +# Test examples for nonlinear programming codes, +# Lectures Notes in Economics and Mathematical Systems 187, +# Springer Verlag, Heidelberg, 1981. +# +# classification QGR-P1-(1,...,6) +# +# Implementation: AI/JSO, 03/2026 + +export hs85 + +"HS85 model" +function hs85_jump() + m = Model(Ipopt.Optimizer) + set_attribute(m, "tol", 1e-10) + set_attribute(m, "print_level", 5) + + # Decision variables + @variable(m, 704.4148 ≤ x1 ≤ 906.3855) + @variable(m, 68.6 ≤ x2 ≤ 288.88) + @variable(m, 0.0 ≤ x3 ≤ 134.75) + @variable(m, 193.0 ≤ x4 ≤ 287.0966) + @variable(m, 25.0 ≤ x5 ≤ 84.1988) + + # Intermediates (as @expression or @variable + @NLconstraint) + @expression(m, y1, x2 + x3 + 41.6) + @expression(m, c1, 0.024 * x4 - 4.62) + @expression(m, y2, 12.5 / c1 + 12) + @expression(m, c2, 0.0003535 * x1^2 + 0.5311 * x1 + 0.08705 * y2 * x1) + @expression(m, c3, 0.052 * x1 + 78 + 0.002377 * y2 * x1) + @expression(m, y3, c2 / c3) + @expression(m, y4, 19 * y3) + @expression(m, c4, 0.04782 * (x1 - y3) + 0.1956 * (x1 - y3)^2 / x2 + + 0.6376 * y4 + 1.594 * y3) + @expression(m, c5, 100 * x2) + @expression(m, c6, x1 - y3 - y4) + @expression(m, c7, 0.95 - c4 / c5) + @expression(m, y5, c6 * c7) + @expression(m, y6, x1 - y5 - y4 - y3) + @expression(m, c8, (y5 + y4) * 0.995) + @expression(m, y7, c8 / y1) + @expression(m, y8, c8 / 3798) + @expression(m, c9, y7 - 0.0663 * y7 / y8 - 0.3153) + @expression(m, y9, 96.82 / c9 + 0.321 * y1) + @expression(m, y10, 1.29 * y5 + 1.258 * y4 + 2.29 * y3 + 1.71 * y6) + @expression(m, y11, 1.71 * x1 - 0.452 * y4 + 0.58 * y3) + @expression(m, c11, 1.75 * y2 * 0.995 * x1) + @expression(m, c12, 0.995 * y10 + 1998) + @expression(m, y12, (12.3/752.3) * x1 + c11 / c12) + @expression(m, y13, c12 - 1.75 * y2) + @expression(m, y14, 3623 + 64.4 * x2 + 58.4 * x3 + 146312 / (y9 + x5)) + @expression(m, c13, 0.995 * y10 + 60.8 * x2 + 48 * x4 - + 0.1121 * y14 - 5095) + @expression(m, y15, y13 / c13) + @expression(m, y16, 148000 - 331000 * y15 + 40 * y13 - 61 * y15 * y13) + @expression(m, c14, 2324 * y10 - 28740000 * y2) + @expression(m, y17, 14130000 - 1328 * y10 - 531 * y11 + c14 / c12) + @expression(m, c15, y13 / y15 - y13 / 0.52) + @expression(m, c16, 1.104 - 0.72 * y15) + @expression(m, c17, y9 + x5) + + # Bounds on y_i + a = [0, 17.505, 11.275, 214.228, 7.458, 0.961, 1.612, 0.146, 107.99, + 922.693, 926.832, 18.766, 1072.163, 8961.448, 0.063, 71084.33, 2802713] + + b = [0, 1053.6667, 35.03, 665.585, 584.463, 265.916, 7.046, 0.222, + 273.366, 1286.105, 1444.046, 537.141, 3247.039, 26844.086, 0.386, + 140000, 12146108] + + for i in 2:17 + yi = Symbol("y$i") + @constraint(m, getfield(Main, yi) >= a[i]) + @constraint(m, getfield(Main, yi) <= b[i]) + end + + # Other inequalities + @constraint(m, 1.5 * x2 - x3 >= 0) + @constraint(m, y1 >= 213.1) + @constraint(m, y1 <= 405.23) + @constraint(m, y4 >= (0.28 / 0.72) * y5) + @constraint(m, 3496 * y2 / c12 <= 21) + @constraint(m, 62212 / c17 - 110.6 >= y1) + + # Objective + @NLobjective(m, Min, + -5.843e-7 * y17 + + 1.17e-4 * y14 + + 2.358e-5 * y13 + + 1.502e-6 * y16 + + 0.0321 * y12 + + 0.004324 * y5 + + 1e-4 * (c15 / c16) + + 37.48 * y2 / c12 + + 0.1365 + ) + + # Good starting point helps a lot + set_start_value(x1, 900.0) + set_start_value(x2, 80.0) + set_start_value(x3, 115.0) + set_start_value(x4, 267.0) + set_start_value(x5, 27.0) + + return m +end \ No newline at end of file diff --git a/src/PureJuMP/hs89.jl b/src/PureJuMP/hs89.jl new file mode 100644 index 000000000..a36876fa4 --- /dev/null +++ b/src/PureJuMP/hs89.jl @@ -0,0 +1,95 @@ +# Hock and Schittkowski problem number 89 +# +# Source: +# Problem 89 in +# W. Hock and K. Schittkowski, +# Test examples for nonlinear programming codes, +# Lectures Notes in Economics and Mathematical Systems 187, +# Springer Verlag, Heidelberg, 1981. +# +# classification QGR-P1-(1,...,6) +# +# Implementation: AI/JSO, 03/2026 + +export hs89 + +"HS89 model" +function hs89_jump( + optimizer = Ipopt.Optimizer, + optimizer_attributes = Dict( + "tol" => 1e-10, + "print_level" => 5, + "max_iter" => 10000, + "acceptable_tol" => 1e-8 + ) +) + model = Model(optimizer) + + # Apply solver-specific options + for (k, v) in optimizer_attributes + set_optimizer_attribute(model, k, v) + end + + # Variables (no simple bounds in standard HS89, but loose ones help numerics) + @variable(model, x1) + @variable(model, x2) + @variable(model, x3) + + # First 30 positive roots of tan(μ) = μ + mu = [ + 0.8603335890193798, 3.425618459481728, 6.437298179171945, 9.529334405361963, + 12.645287223856588, 15.771284874815820, 18.902409956860000, 22.036496727938500, + 25.172446326646600, 28.309642854452000, 31.447714637546200, 34.586424215288900, + 37.725612827776500, 40.865170330488000, 44.005017920830800, 47.145097736761000, + 50.285366337773600, 53.425790477394600, 56.566344279821500, 59.707007305335400, + 62.847763194454400, 65.988598698490300, 69.129502973895200, 72.270467060308900, + 75.411483488848100, 78.552545984242900, 81.693649235601600, 84.834788718042200, + 87.975960552493200, 91.117161394464700 + ] + + # Coefficients Aⱼ + A = [2 * sin(mu[j]) / (mu[j] + sin(mu[j]) * cos(mu[j])) for j in 1:30] + + # Precompute nothing — JuMP + AD will handle it + + # Helper expressions for ρⱼ (makes model more readable) + @expression(model, ρ[j=1:30], + let μ² = mu[j]^2, + r = x1^2 + x2^2 + x3^2, + r23 = x2^2 + x3^2, + r3 = x3^2; + -( exp(-μ² * r) + 2*exp(-μ² * r23) + 2*exp(-μ² * r3) + 1 ) / μ² + end + ) + + # Objective: ∑ Aⱼ ρⱼ + @NLobjective(model, Min, + sum(A[j] * ρ[j] for j in 1:30) + ) + + # Constraint: termA + termB = 2/15 + @NLconstraint(model, eq, + sum( + A[j]^2 * ρ[j]^2 * (sin(2*mu[j])/(2*mu[j]) + 1)/2 + for j in 1:30 + ) + + + sum( + sum( + A[i] * A[j] * ρ[i] * ρ[j] * ( + sin(mu[i] + mu[j]) / (mu[i] + mu[j]) + + sin(mu[i] - mu[j]) / (mu[i] - mu[j]) + ) + for j in i+1:30 + ) + for i in 1:29 + ) + == 2/15 + ) + + # Good starting point (from CUTE / literature) + set_start_value(x1, 0.5) + set_start_value(x2, -0.5) + set_start_value(x3, 0.5) + return model +end \ No newline at end of file From 763e55bb0ebaf9eb7253307315b02bfb11a939f0 Mon Sep 17 00:00:00 2001 From: arnavk23 Date: Tue, 17 Mar 2026 06:46:38 +0530 Subject: [PATCH 02/11] doc files --- docs/RULES_hs85_hs89.md | 22 ++++++++++++++ docs/check_rules.jl | 49 ++++++++++++++++++++++++++++++ docs/extract_hs85.json | 66 +++++++++++++++++++++++++++++++++++++++++ docs/extract_hs89.json | 21 +++++++++++++ 4 files changed, 158 insertions(+) create mode 100644 docs/RULES_hs85_hs89.md create mode 100644 docs/check_rules.jl create mode 100644 docs/extract_hs85.json create mode 100644 docs/extract_hs89.json diff --git a/docs/RULES_hs85_hs89.md b/docs/RULES_hs85_hs89.md new file mode 100644 index 000000000..c5da50116 --- /dev/null +++ b/docs/RULES_hs85_hs89.md @@ -0,0 +1,22 @@ +# OptimizationProblems.jl Rules — HS85 & HS89 Example + +| Rule | HS85 | HS89 | How to Check | +|------|------|------|--------------| +| File Structure | src/ADNLPProblems/hs85.jl, src/PureJuMP/hs85.jl, src/Meta/hs85.md | src/ADNLPProblems/hs89.jl, src/PureJuMP/hs89.jl, src/Meta/hs89.md | Check for all three files per problem | +| Header | HS-style header: source, classification, implementation date | Same | Parse header comments for required fields | +| Mathematical Expressions | All intermediates, constraints, objective explicit, match paper | Same | Compare expressions to paper and extracted JSON | +| Variable Bounds | Explicit bounds and x0 in ADNLP and PureJuMP | Same | Check for bounds and x0 in both files | +| Metadata | References, classification, origin, scalable, etc. in Meta file | Same | Parse Meta file for required metadata | +| Naming | hs85.jl, hs85_jump.jl, hs85.md | hs89.jl, hs89_jump.jl, hs89.md | Check file/function names | +| Allocation | Constraint function minimizes allocations (relaxed for complex) | Same | Run allocation test or allow allocation | +| Ipopt Solve | Problem solves with Ipopt (PureJuMP) | Same | Run Ipopt solver and check for solution | +| Reviewer Markdown | Summary, PDF screenshot, extraction uncertainties, test results | Same | Generate reviewer markdown file | +| Duplication | No duplicate (by name, structure, metadata) | Same | Check for similar names/metadata | +| Traceability | Origin clear and referenced in header/Meta | Same | Check Meta/header | +| Scalability | Marked if scalable (not for hs85/hs89) | Same | Check Meta/implementation | +| Multiple Problems | One PDF = one problem (for hs85/hs89) | Same | Allow multiple if needed | +| Uncertainty | Warn if extraction unclear (not for hs85/hs89) | Same | Parse extraction JSON for uncertainties | + +--- + +This table is tailored to HS85 and HS89. Use as a checklist for similar problems and PRs. diff --git a/docs/check_rules.jl b/docs/check_rules.jl new file mode 100644 index 000000000..70716b042 --- /dev/null +++ b/docs/check_rules.jl @@ -0,0 +1,49 @@ +# Script to check OptimizationProblems.jl rules for HS85/HS89 +using JSON + +function check_problem(problem::String) + result = Dict() + # File structure + result["ADNLP"] = isfile("src/ADNLPProblems/$problem.jl") + result["PureJuMP"] = isfile("src/PureJuMP/$problem.jl") + result["Meta"] = isfile("src/Meta/$problem.md") + # Header + header_ok = false + if result["ADNLP"] + header = readlines("src/ADNLPProblems/$problem.jl")[1:10] + header_ok = any(occursin("Hock and Schittkowski", h) for h in header) + end + result["Header"] = header_ok + # Bounds/x0 + bounds_ok = false + if result["ADNLP"] + file = String(read("src/ADNLPProblems/$problem.jl")) + bounds_ok = occursin("lvar", file) && occursin("uvar", file) && occursin("x0", file) + end + result["Bounds"] = bounds_ok + # Metadata + meta_ok = false + if result["Meta"] + meta = String(read("src/Meta/$problem.md")) + meta_ok = occursin("Source", meta) && occursin("classification", meta) + end + result["Metadata"] = meta_ok + # Allocation (relaxed) + result["Allocation"] = false + # Ipopt solve (not checked here) + result["IpoptSolve"] = "manual" + # Reviewer markdown + result["ReviewerMarkdown"] = isfile("docs/review_$problem.md") + return result +end + +function main() + problems = ["hs85", "hs89"] + results = Dict() + for p in problems + results[p] = check_problem(p) + end + println(JSON.json(results)) +end + +main() diff --git a/docs/extract_hs85.json b/docs/extract_hs85.json new file mode 100644 index 000000000..e45013a2b --- /dev/null +++ b/docs/extract_hs85.json @@ -0,0 +1,66 @@ +{ + "problem": "hs85", + "source": "Hock and Schittkowski, Problem 85", + "variables": ["x1", "x2", "x3", "x4", "x5"], + "bounds": { + "lvar": [704.4148, 68.6, 0.0, 193.0, 25.0], + "uvar": [906.3855, 288.88, 134.75, 287.0966, 84.1988], + "x0": [900.0, 80.0, 115.0, 267.0, 27.0] + }, + "objective": "-5.843e-7 * y17 + 1.17e-4 * y14 + 2.358e-5 * y13 + 1.502e-6 * y16 + 0.0321 * y12 + 0.004324 * y5 + 1e-4 * c15 / c16 + 37.48 * y2 / c12 + 0.1365", + "constraints": [ + "1.5 * x2 - x3 >= 0", + "y1 - 213.1 >= 0", + "405.23 - y1 >= 0", + "y2 - a2 >= 0", + "y3 - a3 >= 0", + "y4 - a4 >= 0", + "y5 - a5 >= 0", + "y6 - a6 >= 0", + "y7 - a7 >= 0", + "y8 - a8 >= 0", + "y9 - a9 >= 0", + "y10 - a10 >= 0", + "y11 - a11 >= 0", + "y12 - a12 >= 0", + "y13 - a13 >= 0", + "y14 - a14 >= 0", + "y15 - a15 >= 0", + "y16 - a16 >= 0", + "y17 - a17 >= 0", + "b2 - y2 >= 0", + "b3 - y3 >= 0", + "b4 - y4 >= 0", + "b5 - y5 >= 0", + "b6 - y6 >= 0", + "b7 - y7 >= 0", + "b8 - y8 >= 0", + "b9 - y9 >= 0", + "b10 - y10 >= 0", + "b11 - y11 >= 0", + "b12 - y12 >= 0", + "b13 - y13 >= 0", + "b14 - y14 >= 0", + "b15 - y15 >= 0", + "b16 - y16 >= 0", + "b17 - y17 >= 0", + "y4 - (0.28 / 0.72) * y5 >= 0", + "21 - 3496 * y2 / c12 >= 0", + "62212 / c17 - 110.6 - y1 >= 0", + "x1 - lvar1 >= 0", + "x2 - lvar2 >= 0", + "x3 - lvar3 >= 0", + "x4 - lvar4 >= 0", + "x5 - lvar5 >= 0", + "uvar1 - x1 >= 0", + "uvar2 - x2 >= 0", + "uvar3 - x3 >= 0", + "uvar4 - x4 >= 0", + "uvar5 - x5 >= 0" + ], + "metadata": { + "classification": "QGR-P1-(1,...,6)", + "implementation": "AI/JSO, 03/2026" + }, + "uncertainties": [] +} diff --git a/docs/extract_hs89.json b/docs/extract_hs89.json new file mode 100644 index 000000000..172f3476f --- /dev/null +++ b/docs/extract_hs89.json @@ -0,0 +1,21 @@ +{ + "problem": "hs89", + "source": "Hock and Schittkowski, Problem 89", + "variables": ["x1", "x2", "x3"], + "bounds": { + "lvar": "none", + "uvar": "none", + "x0": "none" + }, + "objective": "sum_{j=1}^{30} A_j * rho_j(x)\nwhere rho_j(x) = - (exp(-mu_j^2 * r) + 2*exp(-mu_j^2*(x2^2 + x3^2)) + 2*exp(-mu_j^2*x3^2) + 1) / mu_j^2\nr = x1^2 + x2^2 + x3^2\nA_j = 2*sin(mu_j)/(mu_j + sin(mu_j)*cos(mu_j))\nmu_j: first 30 positive roots of tan(mu) = mu", + "constraints": [ + "c(x) = termA + termB - 2/15 = 0", + "termA = sum_{j=1}^{30} A_j^2 * rho_j(x)^2 * (sin(2*mu_j)/(2*mu_j) + 1)/2", + "termB = sum_{i=1}^{29} sum_{j=i+1}^{30} A_i * A_j * rho_i(x) * rho_j(x) * (sin(mu_i+mu_j)/(mu_i+mu_j) + sin(mu_i-mu_j)/(mu_i-mu_j))" + ], + "metadata": { + "classification": "QGR-P1-(1,...,6)", + "implementation": "AI/JSO, 03/2026" + }, + "uncertainties": [] +} From 726ca3298227acc4254022a271d2d9057b858005 Mon Sep 17 00:00:00 2001 From: arnavk23 Date: Tue, 31 Mar 2026 03:01:01 +0530 Subject: [PATCH 03/11] changes --- docs/check_rules.jl | 16 ++++++++-------- 1 file changed, 8 insertions(+), 8 deletions(-) diff --git a/docs/check_rules.jl b/docs/check_rules.jl index 70716b042..5c56927e3 100644 --- a/docs/check_rules.jl +++ b/docs/check_rules.jl @@ -3,36 +3,36 @@ using JSON function check_problem(problem::String) result = Dict() - # File structure + # Check that the required files exist for the problem: ADNLPProblems, PureJuMP, and Meta. result["ADNLP"] = isfile("src/ADNLPProblems/$problem.jl") result["PureJuMP"] = isfile("src/PureJuMP/$problem.jl") result["Meta"] = isfile("src/Meta/$problem.md") - # Header + # Check that the header of the ADNLPProblems file for this problem contains a reference to "Hock and Schittkowski" in the first 10 lines, to ensure proper attribution/documentation. header_ok = false if result["ADNLP"] header = readlines("src/ADNLPProblems/$problem.jl")[1:10] header_ok = any(occursin("Hock and Schittkowski", h) for h in header) end result["Header"] = header_ok - # Bounds/x0 - bounds_ok = false + # Check that the ADNLPProblems file for problem defines lower bounds (lvar), upper bounds (uvar), and an initial point (x0).bounds_ok = false if result["ADNLP"] file = String(read("src/ADNLPProblems/$problem.jl")) bounds_ok = occursin("lvar", file) && occursin("uvar", file) && occursin("x0", file) end result["Bounds"] = bounds_ok - # Metadata + # Check that the Meta file for this problem contains both a "Source" and a "classification" entry. meta_ok = false if result["Meta"] meta = String(read("src/Meta/$problem.md")) meta_ok = occursin("Source", meta) && occursin("classification", meta) end result["Metadata"] = meta_ok - # Allocation (relaxed) + # Allocation check is currently not implemented (set to false as a placeholder). result["Allocation"] = false - # Ipopt solve (not checked here) + # Ipopt solve check is not performed by this script. + # Marked as "manual" for reviewer attention. result["IpoptSolve"] = "manual" - # Reviewer markdown + # Check that a reviewer markdown file exists for this problem in docs/review_$problem.md. result["ReviewerMarkdown"] = isfile("docs/review_$problem.md") return result end From f10898a8c4b9a538b07b05906f7e4161f286da4f Mon Sep 17 00:00:00 2001 From: arnavk23 Date: Tue, 31 Mar 2026 03:22:48 +0530 Subject: [PATCH 04/11] run JuliaFormatter --- src/ADNLPProblems/hs85.jl | 344 +++++++++++++++++++++----------------- src/ADNLPProblems/hs89.jl | 177 +++++++++++--------- src/PureJuMP/hs85.jl | 199 +++++++++++++--------- src/PureJuMP/hs89.jl | 151 +++++++++-------- 4 files changed, 486 insertions(+), 385 deletions(-) diff --git a/src/ADNLPProblems/hs85.jl b/src/ADNLPProblems/hs85.jl index 01644ead1..dc64284a4 100644 --- a/src/ADNLPProblems/hs85.jl +++ b/src/ADNLPProblems/hs85.jl @@ -1,164 +1,192 @@ export hs85 -function hs85(; type::Type{T}=Float64, kwargs...) where {T} - a = T[0, 17.505, 11.275, 214.228, 7.458, 0.961, 1.612, 0.146, 107.99, - 922.693, 926.832, 18.766, 1072.163, 8961.448, 0.063, 71084.33, 2802713] - b = T[0, 1053.6667, 35.03, 665.585, 584.463, 265.916, 7.046, 0.222, - 273.366, 1286.105, 1444.046, 537.141, 3247.039, 26844.086, 0.386, - 140000, 12146108] - c10 = T(12.3) / T(752.3) +function hs85(; type::Type{T} = Float64, kwargs...) where {T} + a = T[ + 0, + 17.505, + 11.275, + 214.228, + 7.458, + 0.961, + 1.612, + 0.146, + 107.99, + 922.693, + 926.832, + 18.766, + 1072.163, + 8961.448, + 0.063, + 71084.33, + 2802713, + ] + b = T[ + 0, + 1053.6667, + 35.03, + 665.585, + 584.463, + 265.916, + 7.046, + 0.222, + 273.366, + 1286.105, + 1444.046, + 537.141, + 3247.039, + 26844.086, + 0.386, + 140000, + 12146108, + ] + c10 = T(12.3) / T(752.3) - # Variable bounds and starting point (exact from the standard model) - lvar = T[704.4148, 68.6, 0.0, 193.0, 25.0] - uvar = T[906.3855, 288.88, 134.75, 287.0966, 84.1988] - x0 = T[900.0, 80.0, 115.0, 267.0, 27.0] + # Variable bounds and starting point (exact from the standard model) + lvar = T[704.4148, 68.6, 0.0, 193.0, 25.0] + uvar = T[906.3855, 288.88, 134.75, 287.0966, 84.1988] + x0 = T[900.0, 80.0, 115.0, 267.0, 27.0] - # Best known value ≈ -1.90513375 - function f(x::AbstractVector{T}) - # All intermediates (identical to those used in constraints) - y1 = x[2] + x[3] + T(41.6) - c1 = T(0.024) * x[4] - T(4.62) - y2 = T(12.5) / c1 + T(12) - c2 = T(0.0003535) * x[1]^2 + T(0.5311) * x[1] + T(0.08705) * y2 * x[1] - c3 = T(0.052) * x[1] + T(78) + T(0.002377) * y2 * x[1] - y3 = c2 / c3 - y4 = T(19) * y3 - c4 = T(0.04782) * (x[1] - y3) + T(0.1956) * (x[1] - y3)^2 / x[2] + - T(0.6376) * y4 + T(1.594) * y3 - c5 = T(100) * x[2] - c6 = x[1] - y3 - y4 - c7 = T(0.95) - c4 / c5 - y5 = c6 * c7 - y6 = x[1] - y5 - y4 - y3 - c8 = (y5 + y4) * T(0.995) - y7 = c8 / y1 - y8 = c8 / T(3798) - c9 = y7 - T(0.0663) * y7 / y8 - T(0.3153) - y9 = T(96.82) / c9 + T(0.321) * y1 - y10 = T(1.29) * y5 + T(1.258) * y4 + T(2.29) * y3 + T(1.71) * y6 - y11 = T(1.71) * x[1] - T(0.452) * y4 + T(0.58) * y3 - c11 = T(1.75) * y2 * T(0.995) * x[1] - c12 = T(0.995) * y10 + T(1998) - y12 = c10 * x[1] + c11 / c12 - y13 = c12 - T(1.75) * y2 - y14 = T(3623) + T(64.4) * x[2] + T(58.4) * x[3] + T(146312) / (y9 + x[5]) - c13 = T(0.995) * y10 + T(60.8) * x[2] + T(48) * x[4] - - T(0.1121) * y14 - T(5095) - y15 = y13 / c13 - y16 = T(148000) - T(331000) * y15 + T(40) * y13 - T(61) * y15 * y13 - c14 = T(2324) * y10 - T(28740000) * y2 - y17 = T(14130000) - T(1328) * y10 - T(531) * y11 + c14 / c12 - c15 = y13 / y15 - y13 / T(0.52) - c16 = T(1.104) - T(0.72) * y15 - # c17 not needed for objective + # Best known value ≈ -1.90513375 + function f(x::AbstractVector{T}) + # All intermediates (identical to those used in constraints) + y1 = x[2] + x[3] + T(41.6) + c1 = T(0.024) * x[4] - T(4.62) + y2 = T(12.5) / c1 + T(12) + c2 = T(0.0003535) * x[1]^2 + T(0.5311) * x[1] + T(0.08705) * y2 * x[1] + c3 = T(0.052) * x[1] + T(78) + T(0.002377) * y2 * x[1] + y3 = c2 / c3 + y4 = T(19) * y3 + c4 = + T(0.04782) * (x[1] - y3) + T(0.1956) * (x[1] - y3)^2 / x[2] + T(0.6376) * y4 + T(1.594) * y3 + c5 = T(100) * x[2] + c6 = x[1] - y3 - y4 + c7 = T(0.95) - c4 / c5 + y5 = c6 * c7 + y6 = x[1] - y5 - y4 - y3 + c8 = (y5 + y4) * T(0.995) + y7 = c8 / y1 + y8 = c8 / T(3798) + c9 = y7 - T(0.0663) * y7 / y8 - T(0.3153) + y9 = T(96.82) / c9 + T(0.321) * y1 + y10 = T(1.29) * y5 + T(1.258) * y4 + T(2.29) * y3 + T(1.71) * y6 + y11 = T(1.71) * x[1] - T(0.452) * y4 + T(0.58) * y3 + c11 = T(1.75) * y2 * T(0.995) * x[1] + c12 = T(0.995) * y10 + T(1998) + y12 = c10 * x[1] + c11 / c12 + y13 = c12 - T(1.75) * y2 + y14 = T(3623) + T(64.4) * x[2] + T(58.4) * x[3] + T(146312) / (y9 + x[5]) + c13 = T(0.995) * y10 + T(60.8) * x[2] + T(48) * x[4] - T(0.1121) * y14 - T(5095) + y15 = y13 / c13 + y16 = T(148000) - T(331000) * y15 + T(40) * y13 - T(61) * y15 * y13 + c14 = T(2324) * y10 - T(28740000) * y2 + y17 = T(14130000) - T(1328) * y10 - T(531) * y11 + c14 / c12 + c15 = y13 / y15 - y13 / T(0.52) + c16 = T(1.104) - T(0.72) * y15 + # c17 not needed for objective - return -T(5.843e-7) * y17 + - T(1.17e-4) * y14 + - T(2.358e-5) * y13 + - T(1.502e-6) * y16 + - T(0.0321) * y12 + - T(0.004324) * y5 + - T(1e-4) * c15 / c16 + - T(37.48) * y2 / c12 + - T(0.1365) - end + return -T(5.843e-7) * y17 + + T(1.17e-4) * y14 + + T(2.358e-5) * y13 + + T(1.502e-6) * y16 + + T(0.0321) * y12 + + T(0.004324) * y5 + + T(1e-4) * c15 / c16 + + T(37.48) * y2 / c12 + + T(0.1365) + end - # Constraint function (48 inequalities, all of the form c(x) >= 0) - function c!(cx::AbstractVector{T}, x::AbstractVector{T}) - # All intermediates (identical to those used in objective) - y1 = x[2] + x[3] + T(41.6) - c1 = T(0.024) * x[4] - T(4.62) - y2 = T(12.5) / c1 + T(12) - c2 = T(0.0003535) * x[1]^2 + T(0.5311) * x[1] + T(0.08705) * y2 * x[1] - c3 = T(0.052) * x[1] + T(78) + T(0.002377) * y2 * x[1] - y3 = c2 / c3 - y4 = T(19) * y3 - c4 = T(0.04782) * (x[1] - y3) + T(0.1956) * (x[1] - y3)^2 / x[2] + T(0.6376) * y4 + T(1.594) * y3 - c5 = T(100) * x[2] - c6 = x[1] - y3 - y4 - c7 = T(0.95) - c4 / c5 - y5 = c6 * c7 - y6 = x[1] - y5 - y4 - y3 - c8 = (y5 + y4) * T(0.995) - y7 = c8 / y1 - y8 = c8 / T(3798) - c9 = y7 - T(0.0663) * y7 / y8 - T(0.3153) - y9 = T(96.82) / c9 + T(0.321) * y1 - y10 = T(1.29) * y5 + T(1.258) * y4 + T(2.29) * y3 + T(1.71) * y6 - y11 = T(1.71) * x[1] - T(0.452) * y4 + T(0.58) * y3 - c11 = T(1.75) * y2 * T(0.995) * x[1] - c12 = T(0.995) * y10 + T(1998) - y12 = c10 * x[1] + c11 / c12 - y13 = c12 - T(1.75) * y2 - y14 = T(3623) + T(64.4) * x[2] + T(58.4) * x[3] + T(146312) / (y9 + x[5]) - c13 = T(0.995) * y10 + T(60.8) * x[2] + T(48) * x[4] - T(0.1121) * y14 - T(5095) - y15 = y13 / c13 - y16 = T(148000) - T(331000) * y15 + T(40) * y13 - T(61) * y15 * y13 - c14 = T(2324) * y10 - T(28740000) * y2 - y17 = T(14130000) - T(1328) * y10 - T(531) * y11 + c14 / c12 - c15 = y13 / y15 - y13 / T(0.52) - c16 = T(1.104) - T(0.72) * y15 - c17 = y9 + x[5] - # Constraints - cx[1] = T(1.5) * x[2] - x[3] - cx[2] = y1 - T(213.1) - cx[3] = T(405.23) - y1 - cx[4] = y2 - a[2] - cx[5] = y3 - a[3] - cx[6] = y4 - a[4] - cx[7] = y5 - a[5] - cx[8] = y6 - a[6] - cx[9] = y7 - a[7] - cx[10] = y8 - a[8] - cx[11] = y9 - a[9] - cx[12] = y10 - a[10] - cx[13] = y11 - a[11] - cx[14] = y12 - a[12] - cx[15] = y13 - a[13] - cx[16] = y14 - a[14] - cx[17] = y15 - a[15] - cx[18] = y16 - a[16] - cx[19] = y17 - a[17] - cx[20] = b[2] - y2 - cx[21] = b[3] - y3 - cx[22] = b[4] - y4 - cx[23] = b[5] - y5 - cx[24] = b[6] - y6 - cx[25] = b[7] - y7 - cx[26] = b[8] - y8 - cx[27] = b[9] - y9 - cx[28] = b[10] - y10 - cx[29] = b[11] - y11 - cx[30] = b[12] - y12 - cx[31] = b[13] - y13 - cx[32] = b[14] - y14 - cx[33] = b[15] - y15 - cx[34] = b[16] - y16 - cx[35] = b[17] - y17 - cx[36] = y4 - (T(0.28) / T(0.72)) * y5 - cx[37] = T(21) - T(3496) * y2 / c12 - cx[38] = T(62212) / c17 - T(110.6) - y1 - # 10 box constraints (5 lower, 5 upper) - cx[39] = x[1] - lvar[1] - cx[40] = x[2] - lvar[2] - cx[41] = x[3] - lvar[3] - cx[42] = x[4] - lvar[4] - cx[43] = x[5] - lvar[5] - cx[44] = uvar[1] - x[1] - cx[45] = uvar[2] - x[2] - cx[46] = uvar[3] - x[3] - cx[47] = uvar[4] - x[4] - cx[48] = uvar[5] - x[5] - return cx - end - # Constraint bounds: all inequalities of the form c(x) >= 0 - m = 48 - cl = zeros(T, m) - cu = fill(T(Inf), m) - return ADNLPModels.ADNLPModel!( - f, x0, lvar, uvar, - c!, cl, cu; - name = "hs85", - kwargs... - ) -end \ No newline at end of file + # Constraint function (48 inequalities, all of the form c(x) >= 0) + function c!(cx::AbstractVector{T}, x::AbstractVector{T}) + # All intermediates (identical to those used in objective) + y1 = x[2] + x[3] + T(41.6) + c1 = T(0.024) * x[4] - T(4.62) + y2 = T(12.5) / c1 + T(12) + c2 = T(0.0003535) * x[1]^2 + T(0.5311) * x[1] + T(0.08705) * y2 * x[1] + c3 = T(0.052) * x[1] + T(78) + T(0.002377) * y2 * x[1] + y3 = c2 / c3 + y4 = T(19) * y3 + c4 = + T(0.04782) * (x[1] - y3) + T(0.1956) * (x[1] - y3)^2 / x[2] + T(0.6376) * y4 + T(1.594) * y3 + c5 = T(100) * x[2] + c6 = x[1] - y3 - y4 + c7 = T(0.95) - c4 / c5 + y5 = c6 * c7 + y6 = x[1] - y5 - y4 - y3 + c8 = (y5 + y4) * T(0.995) + y7 = c8 / y1 + y8 = c8 / T(3798) + c9 = y7 - T(0.0663) * y7 / y8 - T(0.3153) + y9 = T(96.82) / c9 + T(0.321) * y1 + y10 = T(1.29) * y5 + T(1.258) * y4 + T(2.29) * y3 + T(1.71) * y6 + y11 = T(1.71) * x[1] - T(0.452) * y4 + T(0.58) * y3 + c11 = T(1.75) * y2 * T(0.995) * x[1] + c12 = T(0.995) * y10 + T(1998) + y12 = c10 * x[1] + c11 / c12 + y13 = c12 - T(1.75) * y2 + y14 = T(3623) + T(64.4) * x[2] + T(58.4) * x[3] + T(146312) / (y9 + x[5]) + c13 = T(0.995) * y10 + T(60.8) * x[2] + T(48) * x[4] - T(0.1121) * y14 - T(5095) + y15 = y13 / c13 + y16 = T(148000) - T(331000) * y15 + T(40) * y13 - T(61) * y15 * y13 + c14 = T(2324) * y10 - T(28740000) * y2 + y17 = T(14130000) - T(1328) * y10 - T(531) * y11 + c14 / c12 + c15 = y13 / y15 - y13 / T(0.52) + c16 = T(1.104) - T(0.72) * y15 + c17 = y9 + x[5] + # Constraints + cx[1] = T(1.5) * x[2] - x[3] + cx[2] = y1 - T(213.1) + cx[3] = T(405.23) - y1 + cx[4] = y2 - a[2] + cx[5] = y3 - a[3] + cx[6] = y4 - a[4] + cx[7] = y5 - a[5] + cx[8] = y6 - a[6] + cx[9] = y7 - a[7] + cx[10] = y8 - a[8] + cx[11] = y9 - a[9] + cx[12] = y10 - a[10] + cx[13] = y11 - a[11] + cx[14] = y12 - a[12] + cx[15] = y13 - a[13] + cx[16] = y14 - a[14] + cx[17] = y15 - a[15] + cx[18] = y16 - a[16] + cx[19] = y17 - a[17] + cx[20] = b[2] - y2 + cx[21] = b[3] - y3 + cx[22] = b[4] - y4 + cx[23] = b[5] - y5 + cx[24] = b[6] - y6 + cx[25] = b[7] - y7 + cx[26] = b[8] - y8 + cx[27] = b[9] - y9 + cx[28] = b[10] - y10 + cx[29] = b[11] - y11 + cx[30] = b[12] - y12 + cx[31] = b[13] - y13 + cx[32] = b[14] - y14 + cx[33] = b[15] - y15 + cx[34] = b[16] - y16 + cx[35] = b[17] - y17 + cx[36] = y4 - (T(0.28) / T(0.72)) * y5 + cx[37] = T(21) - T(3496) * y2 / c12 + cx[38] = T(62212) / c17 - T(110.6) - y1 + # 10 box constraints (5 lower, 5 upper) + cx[39] = x[1] - lvar[1] + cx[40] = x[2] - lvar[2] + cx[41] = x[3] - lvar[3] + cx[42] = x[4] - lvar[4] + cx[43] = x[5] - lvar[5] + cx[44] = uvar[1] - x[1] + cx[45] = uvar[2] - x[2] + cx[46] = uvar[3] - x[3] + cx[47] = uvar[4] - x[4] + cx[48] = uvar[5] - x[5] + return cx + end + # Constraint bounds: all inequalities of the form c(x) >= 0 + m = 48 + cl = zeros(T, m) + cu = fill(T(Inf), m) + return ADNLPModels.ADNLPModel!(f, x0, lvar, uvar, c!, cl, cu; name = "hs85", kwargs...) +end diff --git a/src/ADNLPProblems/hs89.jl b/src/ADNLPProblems/hs89.jl index d6ec668b7..3c360ee0d 100644 --- a/src/ADNLPProblems/hs89.jl +++ b/src/ADNLPProblems/hs89.jl @@ -1,86 +1,111 @@ export hs89 -function hs89(; type::Type{T} = Float64, kwargs...) where T - # First 30 positive roots of tan(μ) = μ - mu = T[ - 0.8603335890193798, 3.425618459481728, 6.437298179171945, 9.529334405361963, - 12.645287223856588, 15.771284874815820, 18.902409956860000, 22.036496727938500, - 25.172446326646600, 28.309642854452000, 31.447714637546200, 34.586424215288900, - 37.725612827776500, 40.865170330488000, 44.005017920830800, 47.145097736761000, - 50.285366337773600, 53.425790477394600, 56.566344279821500, 59.707007305335400, - 62.847763194454400, 65.988598698490300, 69.129502973895200, 72.270467060308900, - 75.411483488848100, 78.552545984242900, 81.693649235601600, 84.834788718042200, - 87.975960552493200, 91.117161394464700 - ] +function hs89(; type::Type{T} = Float64, kwargs...) where {T} + # First 30 positive roots of tan(μ) = μ + mu = T[ + 0.8603335890193798, + 3.425618459481728, + 6.437298179171945, + 9.529334405361963, + 12.645287223856588, + 15.771284874815820, + 18.902409956860000, + 22.036496727938500, + 25.172446326646600, + 28.309642854452000, + 31.447714637546200, + 34.586424215288900, + 37.725612827776500, + 40.865170330488000, + 44.005017920830800, + 47.145097736761000, + 50.285366337773600, + 53.425790477394600, + 56.566344279821500, + 59.707007305335400, + 62.847763194454400, + 65.988598698490300, + 69.129502973895200, + 72.270467060308900, + 75.411483488848100, + 78.552545984242900, + 81.693649235601600, + 84.834788718042200, + 87.975960552493200, + 91.117161394464700, + ] - # Precomputed coefficients A_j = 2 sin(μ_j) / (μ_j + sin(μ_j) cos(μ_j)) - A = [2 * sin(mu[j]) / (mu[j] + sin(mu[j]) * cos(mu[j])) for j in 1:30] + # Precomputed coefficients A_j = 2 sin(μ_j) / (μ_j + sin(μ_j) cos(μ_j)) + A = [2 * sin(mu[j]) / (mu[j] + sin(mu[j]) * cos(mu[j])) for j = 1:30] - # Objective: φ(x) = ∑_{j=1}^{30} A_j ρ_j(x) - function f(x::AbstractVector{T}) - s = zero(T) - r = x[1]^2 + x[2]^2 + x[3]^2 - for j = 1:30 - μ² = mu[j]^2 - exp_r = exp(-μ² * r) - ρ = - (exp_r + 2*exp(-μ²*(x[2]^2 + x[3]^2)) + 2*exp(-μ²*x[3]^2) + 1) / μ² - s += A[j] * ρ - end - return s + # Objective: φ(x) = ∑_{j=1}^{30} A_j ρ_j(x) + function f(x::AbstractVector{T}) + s = zero(T) + r = x[1]^2 + x[2]^2 + x[3]^2 + for j = 1:30 + μ² = mu[j]^2 + exp_r = exp(-μ² * r) + ρ = - (exp_r + 2*exp(-μ²*(x[2]^2 + x[3]^2)) + 2*exp(-μ²*x[3]^2) + 1) / μ² + s += A[j] * ρ end + return s + end - # Equality constraint c(x) = 0 - # Full expression with cross terms (double sum over i < j) - function c!(cx::AbstractVector{T}, x::AbstractVector{T}) - r = x[1]^2 + x[2]^2 + x[3]^2 - termA = zero(T) - termB = zero(T) - # Compute termA and termB directly, no heap allocation - for j = 1:30 - μ = mu[j] - μ² = μ^2 - exp_r = exp(-μ² * r) - exp_r23 = exp(-μ² * (x[2]^2 + x[3]^2)) - exp_r3 = exp(-μ² * x[3]^2) - ρj = - (exp_r + 2*exp_r23 + 2*exp_r3 + 1) / μ² - termA += A[j]^2 * ρj^2 * (sin(2*μ)/(2*μ) + one(T)) / 2 - end - for i = 1:29 - μi = mu[i] - μi² = μi^2 - exp_ri = exp(-μi² * r) - exp_ri23 = exp(-μi² * (x[2]^2 + x[3]^2)) - exp_ri3 = exp(-μi² * x[3]^2) - ρi = - (exp_ri + 2*exp_ri23 + 2*exp_ri3 + 1) / μi² - for j = i+1:30 - μj = mu[j] - μj² = μj^2 - exp_rj = exp(-μj² * r) - exp_rj23 = exp(-μj² * (x[2]^2 + x[3]^2)) - exp_rj3 = exp(-μj² * x[3]^2) - ρj = - (exp_rj + 2*exp_rj23 + 2*exp_rj3 + 1) / μj² - denom_plus = μi + μj - denom_minus = μi - μj - sin_plus = iszero(denom_plus) ? one(T) : sin(denom_plus)/denom_plus - sin_minus = iszero(denom_minus) ? one(T) : sin(denom_minus)/denom_minus - termB += A[i] * A[j] * ρi * ρj * (sin_plus + sin_minus) - end - end - cx[1] = termA + termB - T(2)/15 - return cx + # Equality constraint c(x) = 0 + # Full expression with cross terms (double sum over i < j) + function c!(cx::AbstractVector{T}, x::AbstractVector{T}) + r = x[1]^2 + x[2]^2 + x[3]^2 + termA = zero(T) + termB = zero(T) + # Compute termA and termB directly, no heap allocation + for j = 1:30 + μ = mu[j] + μ² = μ^2 + exp_r = exp(-μ² * r) + exp_r23 = exp(-μ² * (x[2]^2 + x[3]^2)) + exp_r3 = exp(-μ² * x[3]^2) + ρj = - (exp_r + 2*exp_r23 + 2*exp_r3 + 1) / μ² + termA += A[j]^2 * ρj^2 * (sin(2*μ)/(2*μ) + one(T)) / 2 end + for i = 1:29 + μi = mu[i] + μi² = μi^2 + exp_ri = exp(-μi² * r) + exp_ri23 = exp(-μi² * (x[2]^2 + x[3]^2)) + exp_ri3 = exp(-μi² * x[3]^2) + ρi = - (exp_ri + 2*exp_ri23 + 2*exp_ri3 + 1) / μi² + for j = (i + 1):30 + μj = mu[j] + μj² = μj^2 + exp_rj = exp(-μj² * r) + exp_rj23 = exp(-μj² * (x[2]^2 + x[3]^2)) + exp_rj3 = exp(-μj² * x[3]^2) + ρj = - (exp_rj + 2*exp_rj23 + 2*exp_rj3 + 1) / μj² + denom_plus = μi + μj + denom_minus = μi - μj + sin_plus = iszero(denom_plus) ? one(T) : sin(denom_plus)/denom_plus + sin_minus = iszero(denom_minus) ? one(T) : sin(denom_minus)/denom_minus + termB += A[i] * A[j] * ρi * ρj * (sin_plus + sin_minus) + end + end + cx[1] = termA + termB - T(2)/15 + return cx + end - # Starting point (common in literature / CUTE) - x0 = T[0.5, -0.5, 0.5] + # Starting point (common in literature / CUTE) + x0 = T[0.5, -0.5, 0.5] - # One equality constraint c(x) = 0 - lcon = ucon = T[0] + # One equality constraint c(x) = 0 + lcon = ucon = T[0] - return ADNLPModels.ADNLPModel!( - f, x0, - c!, lcon, ucon; - name = "hs89", - lin = Int[], # no linear constraints - kwargs... - ) -end \ No newline at end of file + return ADNLPModels.ADNLPModel!( + f, + x0, + c!, + lcon, + ucon; + name = "hs89", + lin = Int[], # no linear constraints + kwargs..., + ) +end diff --git a/src/PureJuMP/hs85.jl b/src/PureJuMP/hs85.jl index a07f2df6d..dd55570d7 100644 --- a/src/PureJuMP/hs85.jl +++ b/src/PureJuMP/hs85.jl @@ -15,95 +15,128 @@ export hs85 "HS85 model" function hs85_jump() - m = Model(Ipopt.Optimizer) - set_attribute(m, "tol", 1e-10) - set_attribute(m, "print_level", 5) + m = Model(Ipopt.Optimizer) + set_attribute(m, "tol", 1e-10) + set_attribute(m, "print_level", 5) - # Decision variables - @variable(m, 704.4148 ≤ x1 ≤ 906.3855) - @variable(m, 68.6 ≤ x2 ≤ 288.88) - @variable(m, 0.0 ≤ x3 ≤ 134.75) - @variable(m, 193.0 ≤ x4 ≤ 287.0966) - @variable(m, 25.0 ≤ x5 ≤ 84.1988) + # Decision variables + @variable(m, 704.4148 ≤ x1 ≤ 906.3855) + @variable(m, 68.6 ≤ x2 ≤ 288.88) + @variable(m, 0.0 ≤ x3 ≤ 134.75) + @variable(m, 193.0 ≤ x4 ≤ 287.0966) + @variable(m, 25.0 ≤ x5 ≤ 84.1988) - # Intermediates (as @expression or @variable + @NLconstraint) - @expression(m, y1, x2 + x3 + 41.6) - @expression(m, c1, 0.024 * x4 - 4.62) - @expression(m, y2, 12.5 / c1 + 12) - @expression(m, c2, 0.0003535 * x1^2 + 0.5311 * x1 + 0.08705 * y2 * x1) - @expression(m, c3, 0.052 * x1 + 78 + 0.002377 * y2 * x1) - @expression(m, y3, c2 / c3) - @expression(m, y4, 19 * y3) - @expression(m, c4, 0.04782 * (x1 - y3) + 0.1956 * (x1 - y3)^2 / x2 + - 0.6376 * y4 + 1.594 * y3) - @expression(m, c5, 100 * x2) - @expression(m, c6, x1 - y3 - y4) - @expression(m, c7, 0.95 - c4 / c5) - @expression(m, y5, c6 * c7) - @expression(m, y6, x1 - y5 - y4 - y3) - @expression(m, c8, (y5 + y4) * 0.995) - @expression(m, y7, c8 / y1) - @expression(m, y8, c8 / 3798) - @expression(m, c9, y7 - 0.0663 * y7 / y8 - 0.3153) - @expression(m, y9, 96.82 / c9 + 0.321 * y1) - @expression(m, y10, 1.29 * y5 + 1.258 * y4 + 2.29 * y3 + 1.71 * y6) - @expression(m, y11, 1.71 * x1 - 0.452 * y4 + 0.58 * y3) - @expression(m, c11, 1.75 * y2 * 0.995 * x1) - @expression(m, c12, 0.995 * y10 + 1998) - @expression(m, y12, (12.3/752.3) * x1 + c11 / c12) - @expression(m, y13, c12 - 1.75 * y2) - @expression(m, y14, 3623 + 64.4 * x2 + 58.4 * x3 + 146312 / (y9 + x5)) - @expression(m, c13, 0.995 * y10 + 60.8 * x2 + 48 * x4 - - 0.1121 * y14 - 5095) - @expression(m, y15, y13 / c13) - @expression(m, y16, 148000 - 331000 * y15 + 40 * y13 - 61 * y15 * y13) - @expression(m, c14, 2324 * y10 - 28740000 * y2) - @expression(m, y17, 14130000 - 1328 * y10 - 531 * y11 + c14 / c12) - @expression(m, c15, y13 / y15 - y13 / 0.52) - @expression(m, c16, 1.104 - 0.72 * y15) - @expression(m, c17, y9 + x5) + # Intermediates (as @expression or @variable + @NLconstraint) + @expression(m, y1, x2 + x3 + 41.6) + @expression(m, c1, 0.024 * x4 - 4.62) + @expression(m, y2, 12.5 / c1 + 12) + @expression(m, c2, 0.0003535 * x1^2 + 0.5311 * x1 + 0.08705 * y2 * x1) + @expression(m, c3, 0.052 * x1 + 78 + 0.002377 * y2 * x1) + @expression(m, y3, c2 / c3) + @expression(m, y4, 19 * y3) + @expression(m, c4, 0.04782 * (x1 - y3) + 0.1956 * (x1 - y3)^2 / x2 + 0.6376 * y4 + 1.594 * y3) + @expression(m, c5, 100 * x2) + @expression(m, c6, x1 - y3 - y4) + @expression(m, c7, 0.95 - c4 / c5) + @expression(m, y5, c6 * c7) + @expression(m, y6, x1 - y5 - y4 - y3) + @expression(m, c8, (y5 + y4) * 0.995) + @expression(m, y7, c8 / y1) + @expression(m, y8, c8 / 3798) + @expression(m, c9, y7 - 0.0663 * y7 / y8 - 0.3153) + @expression(m, y9, 96.82 / c9 + 0.321 * y1) + @expression(m, y10, 1.29 * y5 + 1.258 * y4 + 2.29 * y3 + 1.71 * y6) + @expression(m, y11, 1.71 * x1 - 0.452 * y4 + 0.58 * y3) + @expression(m, c11, 1.75 * y2 * 0.995 * x1) + @expression(m, c12, 0.995 * y10 + 1998) + @expression(m, y12, (12.3/752.3) * x1 + c11 / c12) + @expression(m, y13, c12 - 1.75 * y2) + @expression(m, y14, 3623 + 64.4 * x2 + 58.4 * x3 + 146312 / (y9 + x5)) + @expression(m, c13, 0.995 * y10 + 60.8 * x2 + 48 * x4 - 0.1121 * y14 - 5095) + @expression(m, y15, y13 / c13) + @expression(m, y16, 148000 - 331000 * y15 + 40 * y13 - 61 * y15 * y13) + @expression(m, c14, 2324 * y10 - 28740000 * y2) + @expression(m, y17, 14130000 - 1328 * y10 - 531 * y11 + c14 / c12) + @expression(m, c15, y13 / y15 - y13 / 0.52) + @expression(m, c16, 1.104 - 0.72 * y15) + @expression(m, c17, y9 + x5) - # Bounds on y_i - a = [0, 17.505, 11.275, 214.228, 7.458, 0.961, 1.612, 0.146, 107.99, - 922.693, 926.832, 18.766, 1072.163, 8961.448, 0.063, 71084.33, 2802713] + # Bounds on y_i + a = [ + 0, + 17.505, + 11.275, + 214.228, + 7.458, + 0.961, + 1.612, + 0.146, + 107.99, + 922.693, + 926.832, + 18.766, + 1072.163, + 8961.448, + 0.063, + 71084.33, + 2802713, + ] - b = [0, 1053.6667, 35.03, 665.585, 584.463, 265.916, 7.046, 0.222, - 273.366, 1286.105, 1444.046, 537.141, 3247.039, 26844.086, 0.386, - 140000, 12146108] + b = [ + 0, + 1053.6667, + 35.03, + 665.585, + 584.463, + 265.916, + 7.046, + 0.222, + 273.366, + 1286.105, + 1444.046, + 537.141, + 3247.039, + 26844.086, + 0.386, + 140000, + 12146108, + ] - for i in 2:17 - yi = Symbol("y$i") - @constraint(m, getfield(Main, yi) >= a[i]) - @constraint(m, getfield(Main, yi) <= b[i]) - end + for i = 2:17 + yi = Symbol("y$i") + @constraint(m, getfield(Main, yi) >= a[i]) + @constraint(m, getfield(Main, yi) <= b[i]) + end - # Other inequalities - @constraint(m, 1.5 * x2 - x3 >= 0) - @constraint(m, y1 >= 213.1) - @constraint(m, y1 <= 405.23) - @constraint(m, y4 >= (0.28 / 0.72) * y5) - @constraint(m, 3496 * y2 / c12 <= 21) - @constraint(m, 62212 / c17 - 110.6 >= y1) + # Other inequalities + @constraint(m, 1.5 * x2 - x3 >= 0) + @constraint(m, y1 >= 213.1) + @constraint(m, y1 <= 405.23) + @constraint(m, y4 >= (0.28 / 0.72) * y5) + @constraint(m, 3496 * y2 / c12 <= 21) + @constraint(m, 62212 / c17 - 110.6 >= y1) - # Objective - @NLobjective(m, Min, - -5.843e-7 * y17 + - 1.17e-4 * y14 + - 2.358e-5 * y13 + - 1.502e-6 * y16 + - 0.0321 * y12 + - 0.004324 * y5 + - 1e-4 * (c15 / c16) + - 37.48 * y2 / c12 + - 0.1365 - ) + # Objective + @NLobjective( + m, + Min, + -5.843e-7 * y17 + + 1.17e-4 * y14 + + 2.358e-5 * y13 + + 1.502e-6 * y16 + + 0.0321 * y12 + + 0.004324 * y5 + + 1e-4 * (c15 / c16) + + 37.48 * y2 / c12 + + 0.1365 + ) - # Good starting point helps a lot - set_start_value(x1, 900.0) - set_start_value(x2, 80.0) - set_start_value(x3, 115.0) - set_start_value(x4, 267.0) - set_start_value(x5, 27.0) + # Good starting point helps a lot + set_start_value(x1, 900.0) + set_start_value(x2, 80.0) + set_start_value(x3, 115.0) + set_start_value(x4, 267.0) + set_start_value(x5, 27.0) - return m -end \ No newline at end of file + return m +end diff --git a/src/PureJuMP/hs89.jl b/src/PureJuMP/hs89.jl index a36876fa4..158972922 100644 --- a/src/PureJuMP/hs89.jl +++ b/src/PureJuMP/hs89.jl @@ -15,81 +15,96 @@ export hs89 "HS89 model" function hs89_jump( - optimizer = Ipopt.Optimizer, - optimizer_attributes = Dict( - "tol" => 1e-10, - "print_level" => 5, - "max_iter" => 10000, - "acceptable_tol" => 1e-8 - ) + optimizer = Ipopt.Optimizer, + optimizer_attributes = Dict( + "tol" => 1e-10, + "print_level" => 5, + "max_iter" => 10000, + "acceptable_tol" => 1e-8, + ), ) - model = Model(optimizer) + model = Model(optimizer) - # Apply solver-specific options - for (k, v) in optimizer_attributes - set_optimizer_attribute(model, k, v) - end + # Apply solver-specific options + for (k, v) in optimizer_attributes + set_optimizer_attribute(model, k, v) + end - # Variables (no simple bounds in standard HS89, but loose ones help numerics) - @variable(model, x1) - @variable(model, x2) - @variable(model, x3) + # Variables (no simple bounds in standard HS89, but loose ones help numerics) + @variable(model, x1) + @variable(model, x2) + @variable(model, x3) - # First 30 positive roots of tan(μ) = μ - mu = [ - 0.8603335890193798, 3.425618459481728, 6.437298179171945, 9.529334405361963, - 12.645287223856588, 15.771284874815820, 18.902409956860000, 22.036496727938500, - 25.172446326646600, 28.309642854452000, 31.447714637546200, 34.586424215288900, - 37.725612827776500, 40.865170330488000, 44.005017920830800, 47.145097736761000, - 50.285366337773600, 53.425790477394600, 56.566344279821500, 59.707007305335400, - 62.847763194454400, 65.988598698490300, 69.129502973895200, 72.270467060308900, - 75.411483488848100, 78.552545984242900, 81.693649235601600, 84.834788718042200, - 87.975960552493200, 91.117161394464700 - ] + # First 30 positive roots of tan(μ) = μ + mu = [ + 0.8603335890193798, + 3.425618459481728, + 6.437298179171945, + 9.529334405361963, + 12.645287223856588, + 15.771284874815820, + 18.902409956860000, + 22.036496727938500, + 25.172446326646600, + 28.309642854452000, + 31.447714637546200, + 34.586424215288900, + 37.725612827776500, + 40.865170330488000, + 44.005017920830800, + 47.145097736761000, + 50.285366337773600, + 53.425790477394600, + 56.566344279821500, + 59.707007305335400, + 62.847763194454400, + 65.988598698490300, + 69.129502973895200, + 72.270467060308900, + 75.411483488848100, + 78.552545984242900, + 81.693649235601600, + 84.834788718042200, + 87.975960552493200, + 91.117161394464700, + ] - # Coefficients Aⱼ - A = [2 * sin(mu[j]) / (mu[j] + sin(mu[j]) * cos(mu[j])) for j in 1:30] + # Coefficients Aⱼ + A = [2 * sin(mu[j]) / (mu[j] + sin(mu[j]) * cos(mu[j])) for j = 1:30] - # Precompute nothing — JuMP + AD will handle it + # Precompute nothing — JuMP + AD will handle it - # Helper expressions for ρⱼ (makes model more readable) - @expression(model, ρ[j=1:30], - let μ² = mu[j]^2, - r = x1^2 + x2^2 + x3^2, - r23 = x2^2 + x3^2, - r3 = x3^2; - -( exp(-μ² * r) + 2*exp(-μ² * r23) + 2*exp(-μ² * r3) + 1 ) / μ² - end - ) + # Helper expressions for ρⱼ (makes model more readable) + @expression( + model, + ρ[j = 1:30], + let μ² = mu[j]^2, r = x1^2 + x2^2 + x3^2, r23 = x2^2 + x3^2, r3 = x3^2; + -(exp(-μ² * r) + 2*exp(-μ² * r23) + 2*exp(-μ² * r3) + 1) / μ² + end + ) - # Objective: ∑ Aⱼ ρⱼ - @NLobjective(model, Min, - sum(A[j] * ρ[j] for j in 1:30) - ) + # Objective: ∑ Aⱼ ρⱼ + @NLobjective(model, Min, sum(A[j] * ρ[j] for j = 1:30)) - # Constraint: termA + termB = 2/15 - @NLconstraint(model, eq, - sum( - A[j]^2 * ρ[j]^2 * (sin(2*mu[j])/(2*mu[j]) + 1)/2 - for j in 1:30 - ) - + - sum( - sum( - A[i] * A[j] * ρ[i] * ρ[j] * ( - sin(mu[i] + mu[j]) / (mu[i] + mu[j]) + - sin(mu[i] - mu[j]) / (mu[i] - mu[j]) - ) - for j in i+1:30 - ) - for i in 1:29 - ) - == 2/15 - ) + # Constraint: termA + termB = 2/15 + @NLconstraint( + model, + eq, + sum(A[j]^2 * ρ[j]^2 * (sin(2*mu[j])/(2*mu[j]) + 1)/2 for j = 1:30) + sum( + sum( + A[i] * + A[j] * + ρ[i] * + ρ[j] * + (sin(mu[i] + mu[j]) / (mu[i] + mu[j]) + sin(mu[i] - mu[j]) / (mu[i] - mu[j])) for + j = (i + 1):30 + ) for i = 1:29 + ) == 2/15 + ) - # Good starting point (from CUTE / literature) - set_start_value(x1, 0.5) - set_start_value(x2, -0.5) - set_start_value(x3, 0.5) - return model -end \ No newline at end of file + # Good starting point (from CUTE / literature) + set_start_value(x1, 0.5) + set_start_value(x2, -0.5) + set_start_value(x3, 0.5) + return model +end From f30253121215a1dd043955bbcfc30604274728d8 Mon Sep 17 00:00:00 2001 From: arnavk23 Date: Mon, 6 Apr 2026 04:49:48 +0530 Subject: [PATCH 05/11] add notebook --- docs/check_rules.jl | 54 +++++++++++++++++++++++++++++----------- docs/hs85_workflow.ipynb | 0 2 files changed, 39 insertions(+), 15 deletions(-) create mode 100644 docs/hs85_workflow.ipynb diff --git a/docs/check_rules.jl b/docs/check_rules.jl index 5c56927e3..59465e04d 100644 --- a/docs/check_rules.jl +++ b/docs/check_rules.jl @@ -3,28 +3,52 @@ using JSON function check_problem(problem::String) result = Dict() + adnlp_path = "src/ADNLPProblems/$problem.jl" + jump_path = "src/PureJuMP/$problem.jl" + meta_path = "src/Meta/$problem.jl" # Check that the required files exist for the problem: ADNLPProblems, PureJuMP, and Meta. - result["ADNLP"] = isfile("src/ADNLPProblems/$problem.jl") - result["PureJuMP"] = isfile("src/PureJuMP/$problem.jl") - result["Meta"] = isfile("src/Meta/$problem.md") + result["ADNLP"] = isfile(adnlp_path) + result["PureJuMP"] = isfile(jump_path) + result["Meta"] = isfile(meta_path) # Check that the header of the ADNLPProblems file for this problem contains a reference to "Hock and Schittkowski" in the first 10 lines, to ensure proper attribution/documentation. header_ok = false - if result["ADNLP"] - header = readlines("src/ADNLPProblems/$problem.jl")[1:10] - header_ok = any(occursin("Hock and Schittkowski", h) for h in header) + if result["ADNLP"] || result["PureJuMP"] + adnlp_header_ok = false + jump_header_ok = false + if result["ADNLP"] + adnlp_lines = readlines(adnlp_path) + k = min(length(adnlp_lines), 20) + adnlp_header_ok = any(occursin("Hock and Schittkowski", adnlp_lines[i]) for i in 1:k) + end + if result["PureJuMP"] + jump_lines = readlines(jump_path) + k = min(length(jump_lines), 20) + jump_header_ok = any(occursin("Hock and Schittkowski", jump_lines[i]) for i in 1:k) + end + header_ok = adnlp_header_ok || jump_header_ok end result["Header"] = header_ok - # Check that the ADNLPProblems file for problem defines lower bounds (lvar), upper bounds (uvar), and an initial point (x0).bounds_ok = false + # Check that the ADNLPProblems file for problem defines lower bounds (lvar), upper bounds (uvar), and an initial point (x0). + bounds_ok = false if result["ADNLP"] - file = String(read("src/ADNLPProblems/$problem.jl")) - bounds_ok = occursin("lvar", file) && occursin("uvar", file) && occursin("x0", file) + adnlp_file = String(read(adnlp_path)) + has_x0 = occursin("x0", adnlp_file) + has_bounds_in_code = occursin("lvar", adnlp_file) && occursin("uvar", adnlp_file) + has_bounds_required = true + if result["Meta"] + meta_file = String(read(meta_path)) + if occursin(":has_bounds => false", meta_file) + has_bounds_required = false + end + end + bounds_ok = has_x0 && (has_bounds_in_code || !has_bounds_required) end result["Bounds"] = bounds_ok # Check that the Meta file for this problem contains both a "Source" and a "classification" entry. meta_ok = false if result["Meta"] - meta = String(read("src/Meta/$problem.md")) - meta_ok = occursin("Source", meta) && occursin("classification", meta) + meta = String(read(meta_path)) + meta_ok = occursin(":name => \"$problem\"", meta) && occursin(":best_known_upper_bound", meta) end result["Metadata"] = meta_ok # Allocation check is currently not implemented (set to false as a placeholder). @@ -32,18 +56,18 @@ function check_problem(problem::String) # Ipopt solve check is not performed by this script. # Marked as "manual" for reviewer attention. result["IpoptSolve"] = "manual" - # Check that a reviewer markdown file exists for this problem in docs/review_$problem.md. - result["ReviewerMarkdown"] = isfile("docs/review_$problem.md") + # Check whether reviewer markdown exists (optional, but useful for PR traceability). + result["ReviewerMarkdown"] = isfile("docs/review_$problem.md") || isfile("docs/review/$problem.md") return result end function main() - problems = ["hs85", "hs89"] + problems = isempty(ARGS) ? ["hs85", "hs89"] : ARGS results = Dict() for p in problems results[p] = check_problem(p) end - println(JSON.json(results)) + println(JSON.json(results, 2)) end main() diff --git a/docs/hs85_workflow.ipynb b/docs/hs85_workflow.ipynb new file mode 100644 index 000000000..e69de29bb From 355f0686e236319dd1840df65a7b2fe3e661bbe8 Mon Sep 17 00:00:00 2001 From: arnavk23 Date: Mon, 6 Apr 2026 04:50:19 +0530 Subject: [PATCH 06/11] update --- docs/hs85_workflow.ipynb | 247 +++++++++++++++++++++++++++++++++++++++ 1 file changed, 247 insertions(+) diff --git a/docs/hs85_workflow.ipynb b/docs/hs85_workflow.ipynb index e69de29bb..0f4096ffe 100644 --- a/docs/hs85_workflow.ipynb +++ b/docs/hs85_workflow.ipynb @@ -0,0 +1,247 @@ +{ + "cells": [ + { + "cell_type": "markdown", + "id": "ff939ce8", + "metadata": {}, + "source": [ + "# HS85 End-to-End Workflow (PR #412 style)\n", + "\n", + "This notebook demonstrates one complete flow for adding/verifying a problem in `OptimizationProblems.jl` using **HS85** as a focused test case.\n", + "\n", + "Flow covered:\n", + "1. Load extraction artifact (`docs/extract_hs85.json`).\n", + "2. Run repository rule checks (files, header, bounds, metadata).\n", + "3. Compare extracted bounds with implementation bounds.\n", + "4. Instantiate the ADNLP problem and inspect dimensions.\n", + "5. Optionally solve the PureJuMP model with Ipopt (if available).\n", + "\n", + "The goal is to make the process reviewable and reproducible in one place." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "dd6abb58", + "metadata": { + "vscode": { + "languageId": "julia" + } + }, + "outputs": [], + "source": [ + "using Pkg\n", + "\n", + "repo_root = if isfile(joinpath(pwd(), \"Project.toml\"))\n", + " pwd()\n", + "elseif isfile(joinpath(pwd(), \"..\", \"Project.toml\"))\n", + " normpath(joinpath(pwd(), \"..\"))\n", + "else\n", + " error(\"Could not locate Project.toml from current working directory: $(pwd())\")\n", + "end\n", + "\n", + "Pkg.activate(repo_root)\n", + "\n", + "using JSON\n", + "using OptimizationProblems\n", + "using ADNLPModels\n", + "using NLPModels\n", + "\n", + "println(\"Activated project at: \", repo_root)" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "cfaf1484", + "metadata": { + "vscode": { + "languageId": "julia" + } + }, + "outputs": [], + "source": [ + "extract_path = joinpath(repo_root, \"docs\", \"extract_hs85.json\")\n", + "extract = JSON.parsefile(extract_path)\n", + "\n", + "println(\"Loaded extraction file: \", extract_path)\n", + "println(\"Problem: \", extract[\"problem\"])\n", + "println(\"Variables: \", length(extract[\"variables\"]))\n", + "println(\"Constraints listed in JSON: \", length(extract[\"constraints\"]))\n", + "println(\"Uncertainties count: \", length(extract[\"uncertainties\"]))" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "7654a919", + "metadata": { + "vscode": { + "languageId": "julia" + } + }, + "outputs": [], + "source": [ + "function first_lines_have_pattern(path::String, pattern::String; n::Int = 20)\n", + " if !isfile(path)\n", + " return false\n", + " end\n", + " lines = readlines(path)\n", + " k = min(length(lines), n)\n", + " return any(occursin(pattern, lines[i]) for i in 1:k)\n", + "end\n", + "\n", + "function check_problem_files_and_rules(repo_root::String, problem::String)\n", + " adnlp_path = joinpath(repo_root, \"src\", \"ADNLPProblems\", \"$(problem).jl\")\n", + " jump_path = joinpath(repo_root, \"src\", \"PureJuMP\", \"$(problem).jl\")\n", + " meta_path = joinpath(repo_root, \"src\", \"Meta\", \"$(problem).jl\")\n", + "\n", + " file_ok = Dict(\n", + " \"ADNLP\" => isfile(adnlp_path),\n", + " \"PureJuMP\" => isfile(jump_path),\n", + " \"Meta\" => isfile(meta_path),\n", + " )\n", + "\n", + " header_ok = first_lines_have_pattern(adnlp_path, \"Hock and Schittkowski\")\n", + "\n", + " adnlp_text = file_ok[\"ADNLP\"] ? read(adnlp_path, String) : \"\"\n", + " bounds_ok = occursin(\"lvar\", adnlp_text) && occursin(\"uvar\", adnlp_text) && occursin(\"x0\", adnlp_text)\n", + "\n", + " meta_text = file_ok[\"Meta\"] ? read(meta_path, String) : \"\"\n", + " metadata_ok = occursin(\":name => \\\"$(problem)\\\"\", meta_text) &&\n", + " occursin(\":best_known_upper_bound\", meta_text)\n", + "\n", + " return Dict(\n", + " \"problem\" => problem,\n", + " \"files\" => file_ok,\n", + " \"header\" => header_ok,\n", + " \"bounds\" => bounds_ok,\n", + " \"metadata\" => metadata_ok,\n", + " \"all_core_checks_pass\" => all(values(file_ok)) && header_ok && bounds_ok && metadata_ok,\n", + " )\n", + "end\n", + "\n", + "checks = check_problem_files_and_rules(repo_root, \"hs85\")\n", + "println(JSON.json(checks, 2))" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "78215d1d", + "metadata": { + "vscode": { + "languageId": "julia" + } + }, + "outputs": [], + "source": [ + "function parse_numeric_vector_literal(text::String, name::String)\n", + " m = match(Regex(name * raw\"\\s*=\\s*T\\[(.*?)\\]\", \"s\"), text)\n", + " m === nothing && error(\"Could not locate vector literal for $(name)\")\n", + " body = m.captures[1]\n", + " tokens = split(replace(body, '\\n' => ' '), ',')\n", + " vals = Float64[]\n", + " for t in tokens\n", + " s = strip(t)\n", + " isempty(s) && continue\n", + " push!(vals, parse(Float64, s))\n", + " end\n", + " return vals\n", + "end\n", + "\n", + "adnlp_text = read(joinpath(repo_root, \"src\", \"ADNLPProblems\", \"hs85.jl\"), String)\n", + "code_lvar = parse_numeric_vector_literal(adnlp_text, \"lvar\")\n", + "code_uvar = parse_numeric_vector_literal(adnlp_text, \"uvar\")\n", + "code_x0 = parse_numeric_vector_literal(adnlp_text, \"x0\")\n", + "\n", + "json_lvar = Float64.(extract[\"bounds\"][\"lvar\"])\n", + "json_uvar = Float64.(extract[\"bounds\"][\"uvar\"])\n", + "json_x0 = Float64.(extract[\"bounds\"][\"x0\"])\n", + "\n", + "println(\"Bounds and x0 consistency with extraction JSON:\")\n", + "println(\"lvar match: \", code_lvar == json_lvar)\n", + "println(\"uvar match: \", code_uvar == json_uvar)\n", + "println(\"x0 match: \", code_x0 == json_x0)" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "c8224124", + "metadata": { + "vscode": { + "languageId": "julia" + } + }, + "outputs": [], + "source": [ + "try\n", + " nlp = OptimizationProblems.ADNLPProblems.hs85()\n", + " x0 = nlp.meta.x0\n", + " f0 = NLPModels.obj(nlp, x0)\n", + " c0 = NLPModels.cons(nlp, x0)\n", + "\n", + " println(\"ADNLP model summary:\")\n", + " println(\"name: \", nlp.meta.name)\n", + " println(\"nvar: \", nlp.meta.nvar)\n", + " println(\"ncon: \", nlp.meta.ncon)\n", + " println(\"objective at x0: \", f0)\n", + " println(\"min constraint value at x0 (for c(x) >= 0 style constraints): \", minimum(c0))\n", + "catch err\n", + " println(\"Skipped ADNLP instantiation step. Reason: \", err)\n", + "end" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "5863d263", + "metadata": { + "vscode": { + "languageId": "julia" + } + }, + "outputs": [], + "source": [ + "try\n", + " using JuMP, Ipopt\n", + " model = OptimizationProblems.PureJuMP.hs85_jump()\n", + " optimize!(model)\n", + " println(\"Ipopt termination status: \", termination_status(model))\n", + " println(\"Objective value: \", objective_value(model))\n", + "catch err\n", + " println(\"Skipped solve step. Reason: \", err)\n", + " println(\"This is expected if Ipopt/JuMP are not available in the current notebook environment.\")\n", + "end" + ] + }, + { + "cell_type": "markdown", + "id": "46e5809a", + "metadata": {}, + "source": [ + "## Manual vs AI Intervention Notes\n", + "\n", + "What this notebook automates well:\n", + "- Structural checks (file existence, headers, basic metadata patterns).\n", + "- Consistency checks between extraction JSON and implementation bounds/x0.\n", + "- Basic model instantiation and optional solve sanity-check.\n", + "\n", + "What still needs manual review:\n", + "- Ambiguous PDF interpretation and formula disambiguation.\n", + "- Naming choices and duplicate detection decisions.\n", + "- Semantic fidelity of objective/constraints against original publication.\n", + "\n", + "This balance matches PR #412 workflow: AI accelerates drafting/checking, while humans validate mathematical intent and project policy decisions." + ] + } + ], + "metadata": { + "language_info": { + "name": "python" + } + }, + "nbformat": 4, + "nbformat_minor": 5 +} From dbd831160165d73bd87aab6348123433cd8df543 Mon Sep 17 00:00:00 2001 From: Arnav Kapoor Date: Mon, 6 Apr 2026 05:24:53 +0530 Subject: [PATCH 07/11] Apply suggestions from code review Co-authored-by: Copilot <175728472+Copilot@users.noreply.github.com> --- docs/RULES_hs85_hs89.md | 4 ++-- docs/extract_hs89.json | 2 +- src/Meta/hs89.jl | 4 ++-- src/PureJuMP/hs85.jl | 6 ++---- src/PureJuMP/hs89.jl | 5 ++++- 5 files changed, 11 insertions(+), 10 deletions(-) diff --git a/docs/RULES_hs85_hs89.md b/docs/RULES_hs85_hs89.md index c5da50116..6b8555a9f 100644 --- a/docs/RULES_hs85_hs89.md +++ b/docs/RULES_hs85_hs89.md @@ -2,12 +2,12 @@ | Rule | HS85 | HS89 | How to Check | |------|------|------|--------------| -| File Structure | src/ADNLPProblems/hs85.jl, src/PureJuMP/hs85.jl, src/Meta/hs85.md | src/ADNLPProblems/hs89.jl, src/PureJuMP/hs89.jl, src/Meta/hs89.md | Check for all three files per problem | +| File Structure | src/ADNLPProblems/hs85.jl, src/PureJuMP/hs85.jl, src/Meta/hs85.jl | src/ADNLPProblems/hs89.jl, src/PureJuMP/hs89.jl, src/Meta/hs89.jl | Check for all three files per problem | | Header | HS-style header: source, classification, implementation date | Same | Parse header comments for required fields | | Mathematical Expressions | All intermediates, constraints, objective explicit, match paper | Same | Compare expressions to paper and extracted JSON | | Variable Bounds | Explicit bounds and x0 in ADNLP and PureJuMP | Same | Check for bounds and x0 in both files | | Metadata | References, classification, origin, scalable, etc. in Meta file | Same | Parse Meta file for required metadata | -| Naming | hs85.jl, hs85_jump.jl, hs85.md | hs89.jl, hs89_jump.jl, hs89.md | Check file/function names | +| Naming | hs85.jl in ADNLPProblems, hs85.jl in PureJuMP, hs85.jl in Meta | hs89.jl in ADNLPProblems, hs89.jl in PureJuMP, hs89.jl in Meta | Check file/function names against repository conventions | | Allocation | Constraint function minimizes allocations (relaxed for complex) | Same | Run allocation test or allow allocation | | Ipopt Solve | Problem solves with Ipopt (PureJuMP) | Same | Run Ipopt solver and check for solution | | Reviewer Markdown | Summary, PDF screenshot, extraction uncertainties, test results | Same | Generate reviewer markdown file | diff --git a/docs/extract_hs89.json b/docs/extract_hs89.json index 172f3476f..12668aa5c 100644 --- a/docs/extract_hs89.json +++ b/docs/extract_hs89.json @@ -5,7 +5,7 @@ "bounds": { "lvar": "none", "uvar": "none", - "x0": "none" + "x0": [0.0, 0.0, 0.0] }, "objective": "sum_{j=1}^{30} A_j * rho_j(x)\nwhere rho_j(x) = - (exp(-mu_j^2 * r) + 2*exp(-mu_j^2*(x2^2 + x3^2)) + 2*exp(-mu_j^2*x3^2) + 1) / mu_j^2\nr = x1^2 + x2^2 + x3^2\nA_j = 2*sin(mu_j)/(mu_j + sin(mu_j)*cos(mu_j))\nmu_j: first 30 positive roots of tan(mu) = mu", "constraints": [ diff --git a/src/Meta/hs89.jl b/src/Meta/hs89.jl index 13ac47fcb..5163b9bb9 100644 --- a/src/Meta/hs89.jl +++ b/src/Meta/hs89.jl @@ -5,8 +5,8 @@ hs89_meta = Dict( :variable_ncon => false, :minimize => true, :name => "hs89", - :has_equalities_only => false, - :has_inequalities_only => true, + :has_equalities_only => true, + :has_inequalities_only => false, :has_bounds => false, :has_fixed_variables => false, :objtype => :other, diff --git a/src/PureJuMP/hs85.jl b/src/PureJuMP/hs85.jl index dd55570d7..2b364f265 100644 --- a/src/PureJuMP/hs85.jl +++ b/src/PureJuMP/hs85.jl @@ -14,10 +14,8 @@ export hs85 "HS85 model" -function hs85_jump() - m = Model(Ipopt.Optimizer) - set_attribute(m, "tol", 1e-10) - set_attribute(m, "print_level", 5) +function hs85() + m = Model() # Decision variables @variable(m, 704.4148 ≤ x1 ≤ 906.3855) diff --git a/src/PureJuMP/hs89.jl b/src/PureJuMP/hs89.jl index 158972922..c939819d3 100644 --- a/src/PureJuMP/hs89.jl +++ b/src/PureJuMP/hs89.jl @@ -21,8 +21,11 @@ function hs89_jump( "print_level" => 5, "max_iter" => 10000, "acceptable_tol" => 1e-8, - ), + ); + n::Int = 3, + kwargs..., ) + _ = n, kwargs model = Model(optimizer) # Apply solver-specific options From a4c4bdc58525eedc428306e1eaf4c95c2929a2b1 Mon Sep 17 00:00:00 2001 From: Arnav Kapoor Date: Mon, 6 Apr 2026 05:32:14 +0530 Subject: [PATCH 08/11] Update src/PureJuMP/hs89.jl Co-authored-by: Copilot <175728472+Copilot@users.noreply.github.com> --- src/PureJuMP/hs89.jl | 19 ++++++++----------- 1 file changed, 8 insertions(+), 11 deletions(-) diff --git a/src/PureJuMP/hs89.jl b/src/PureJuMP/hs89.jl index c939819d3..749884922 100644 --- a/src/PureJuMP/hs89.jl +++ b/src/PureJuMP/hs89.jl @@ -18,19 +18,16 @@ function hs89_jump( optimizer = Ipopt.Optimizer, optimizer_attributes = Dict( "tol" => 1e-10, - "print_level" => 5, - "max_iter" => 10000, - "acceptable_tol" => 1e-8, - ); - n::Int = 3, - kwargs..., + optimizer = nothing, + optimizer_attributes = nothing, ) - _ = n, kwargs - model = Model(optimizer) + model = optimizer === nothing ? Model() : Model(optimizer) - # Apply solver-specific options - for (k, v) in optimizer_attributes - set_optimizer_attribute(model, k, v) + # Apply solver-specific options only when an optimizer and attributes are provided + if optimizer_attributes !== nothing + for (k, v) in optimizer_attributes + set_optimizer_attribute(model, k, v) + end end # Variables (no simple bounds in standard HS89, but loose ones help numerics) From d30b774f52339e881b14725f317772be59848a29 Mon Sep 17 00:00:00 2001 From: arnavk23 Date: Mon, 6 Apr 2026 05:49:23 +0530 Subject: [PATCH 09/11] adding copilot suggestion --- docs/hs85_workflow.ipynb | 47 +++++++++++++-------------------------- src/ADNLPProblems/hs85.jl | 15 ++----------- src/Meta/hs85.jl | 8 +++---- src/PureJuMP/hs85.jl | 8 +++---- src/PureJuMP/hs89.jl | 14 ++++++------ 5 files changed, 32 insertions(+), 60 deletions(-) diff --git a/docs/hs85_workflow.ipynb b/docs/hs85_workflow.ipynb index 0f4096ffe..2c4671c4f 100644 --- a/docs/hs85_workflow.ipynb +++ b/docs/hs85_workflow.ipynb @@ -23,11 +23,7 @@ "cell_type": "code", "execution_count": null, "id": "dd6abb58", - "metadata": { - "vscode": { - "languageId": "julia" - } - }, + "metadata": {}, "outputs": [], "source": [ "using Pkg\n", @@ -54,11 +50,7 @@ "cell_type": "code", "execution_count": null, "id": "cfaf1484", - "metadata": { - "vscode": { - "languageId": "julia" - } - }, + "metadata": {}, "outputs": [], "source": [ "extract_path = joinpath(repo_root, \"docs\", \"extract_hs85.json\")\n", @@ -75,11 +67,7 @@ "cell_type": "code", "execution_count": null, "id": "7654a919", - "metadata": { - "vscode": { - "languageId": "julia" - } - }, + "metadata": {}, "outputs": [], "source": [ "function first_lines_have_pattern(path::String, pattern::String; n::Int = 20)\n", @@ -129,11 +117,7 @@ "cell_type": "code", "execution_count": null, "id": "78215d1d", - "metadata": { - "vscode": { - "languageId": "julia" - } - }, + "metadata": {}, "outputs": [], "source": [ "function parse_numeric_vector_literal(text::String, name::String)\n", @@ -169,11 +153,7 @@ "cell_type": "code", "execution_count": null, "id": "c8224124", - "metadata": { - "vscode": { - "languageId": "julia" - } - }, + "metadata": {}, "outputs": [], "source": [ "try\n", @@ -197,16 +177,12 @@ "cell_type": "code", "execution_count": null, "id": "5863d263", - "metadata": { - "vscode": { - "languageId": "julia" - } - }, + "metadata": {}, "outputs": [], "source": [ "try\n", " using JuMP, Ipopt\n", - " model = OptimizationProblems.PureJuMP.hs85_jump()\n", + " model = OptimizationProblems.PureJuMP.hs85()\n", " optimize!(model)\n", " println(\"Ipopt termination status: \", termination_status(model))\n", " println(\"Objective value: \", objective_value(model))\n", @@ -238,10 +214,17 @@ } ], "metadata": { + "kernelspec": { + "display_name": "julia 1.11", + "language": "julia", + "name": "julia" + }, "language_info": { - "name": "python" + "name": "julia", + "version": "3.11.9" } }, "nbformat": 4, "nbformat_minor": 5 } + diff --git a/src/ADNLPProblems/hs85.jl b/src/ADNLPProblems/hs85.jl index dc64284a4..135620684 100644 --- a/src/ADNLPProblems/hs85.jl +++ b/src/ADNLPProblems/hs85.jl @@ -95,7 +95,7 @@ function hs85(; type::Type{T} = Float64, kwargs...) where {T} T(0.1365) end - # Constraint function (48 inequalities, all of the form c(x) >= 0) + # Constraint function (38 nonlinear inequalities, all of the form c(x) >= 0) function c!(cx::AbstractVector{T}, x::AbstractVector{T}) # All intermediates (identical to those used in objective) y1 = x[2] + x[3] + T(41.6) @@ -171,21 +171,10 @@ function hs85(; type::Type{T} = Float64, kwargs...) where {T} cx[36] = y4 - (T(0.28) / T(0.72)) * y5 cx[37] = T(21) - T(3496) * y2 / c12 cx[38] = T(62212) / c17 - T(110.6) - y1 - # 10 box constraints (5 lower, 5 upper) - cx[39] = x[1] - lvar[1] - cx[40] = x[2] - lvar[2] - cx[41] = x[3] - lvar[3] - cx[42] = x[4] - lvar[4] - cx[43] = x[5] - lvar[5] - cx[44] = uvar[1] - x[1] - cx[45] = uvar[2] - x[2] - cx[46] = uvar[3] - x[3] - cx[47] = uvar[4] - x[4] - cx[48] = uvar[5] - x[5] return cx end # Constraint bounds: all inequalities of the form c(x) >= 0 - m = 48 + m = 38 cl = zeros(T, m) cu = fill(T(Inf), m) return ADNLPModels.ADNLPModel!(f, x0, lvar, uvar, c!, cl, cu; name = "hs85", kwargs...) diff --git a/src/Meta/hs85.jl b/src/Meta/hs85.jl index eba8e3197..0cbd12a77 100644 --- a/src/Meta/hs85.jl +++ b/src/Meta/hs85.jl @@ -1,7 +1,7 @@ hs85_meta = Dict( :nvar => 5, :variable_nvar => false, - :ncon => 48, + :ncon => 38, :variable_ncon => false, :minimize => true, :name => "hs85", @@ -18,6 +18,6 @@ hs85_meta = Dict( :origin => :unknown, ) get_hs85_nvar(; n::Integer = 5, kwargs...) = 5 -get_hs85_ncon(; n::Integer = 48, kwargs...) = 48 -get_hs85_nlin(; n::Integer = 48, kwargs...) = 0 -get_hs85_nnln(; n::Integer = 48, kwargs...) = 48 +get_hs85_ncon(; n::Integer = 38, kwargs...) = 38 +get_hs85_nlin(; n::Integer = 38, kwargs...) = 0 +get_hs85_nnln(; n::Integer = 38, kwargs...) = 38 diff --git a/src/PureJuMP/hs85.jl b/src/PureJuMP/hs85.jl index 2b364f265..f9db0bfce 100644 --- a/src/PureJuMP/hs85.jl +++ b/src/PureJuMP/hs85.jl @@ -14,7 +14,7 @@ export hs85 "HS85 model" -function hs85() +function hs85(args...; kwargs...) m = Model() # Decision variables @@ -100,10 +100,10 @@ function hs85() 12146108, ] + y = Any[y1, y2, y3, y4, y5, y6, y7, y8, y9, y10, y11, y12, y13, y14, y15, y16, y17] for i = 2:17 - yi = Symbol("y$i") - @constraint(m, getfield(Main, yi) >= a[i]) - @constraint(m, getfield(Main, yi) <= b[i]) + @constraint(m, y[i] >= a[i]) + @constraint(m, y[i] <= b[i]) end # Other inequalities diff --git a/src/PureJuMP/hs89.jl b/src/PureJuMP/hs89.jl index 749884922..b03773f3b 100644 --- a/src/PureJuMP/hs89.jl +++ b/src/PureJuMP/hs89.jl @@ -14,15 +14,18 @@ export hs89 "HS89 model" -function hs89_jump( - optimizer = Ipopt.Optimizer, - optimizer_attributes = Dict( - "tol" => 1e-10, +function hs89( + args...; optimizer = nothing, optimizer_attributes = nothing, + kwargs..., ) model = optimizer === nothing ? Model() : Model(optimizer) + if optimizer !== nothing && optimizer_attributes === nothing + optimizer_attributes = Dict("tol" => 1e-10) + end + # Apply solver-specific options only when an optimizer and attributes are provided if optimizer_attributes !== nothing for (k, v) in optimizer_attributes @@ -72,9 +75,6 @@ function hs89_jump( # Coefficients Aⱼ A = [2 * sin(mu[j]) / (mu[j] + sin(mu[j]) * cos(mu[j])) for j = 1:30] - # Precompute nothing — JuMP + AD will handle it - - # Helper expressions for ρⱼ (makes model more readable) @expression( model, ρ[j = 1:30], From 3684253bf3815e1783048e5472d1593a5ac0ebb6 Mon Sep 17 00:00:00 2001 From: arnavk23 Date: Fri, 4 Sep 2026 03:12:22 +0530 Subject: [PATCH 10/11] Fix hs85/hs89 PureJuMP/ADNLP compatibility and missing meta getters --- Project.toml | 2 - docs/extract_hs89.json | 2 +- src/ADNLPProblems/hs85.jl | 110 +++++++++++++++++++++++--------------- src/ADNLPProblems/hs89.jl | 4 +- src/Meta/hs85.jl | 10 ++-- src/Meta/hs89.jl | 10 ++-- src/PureJuMP/hs85.jl | 2 +- src/PureJuMP/hs89.jl | 13 ++--- 8 files changed, 87 insertions(+), 66 deletions(-) diff --git a/Project.toml b/Project.toml index c2a178e86..f69a86a54 100644 --- a/Project.toml +++ b/Project.toml @@ -10,7 +10,6 @@ LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e" Requires = "ae029012-a4dd-5104-9daa-d747884805df" SparseArrays = "2f01184e-e22b-5df5-ae63-d93ebab69eaf" SpecialFunctions = "276daf66-3868-5448-9aa4-cd146d93841b" -StaticArrays = "90137ffa-7385-5640-81b9-e52037218182" [compat] DataFrames = "1" @@ -18,7 +17,6 @@ JLD2 = "0.5, 0.6" JuMP = "^1.15" Requires = "1" SpecialFunctions = "2" -StaticArrays = "1.9.18" julia = "1.6" [extras] diff --git a/docs/extract_hs89.json b/docs/extract_hs89.json index 12668aa5c..53b5548ca 100644 --- a/docs/extract_hs89.json +++ b/docs/extract_hs89.json @@ -5,7 +5,7 @@ "bounds": { "lvar": "none", "uvar": "none", - "x0": [0.0, 0.0, 0.0] + "x0": [0.5, -0.5, 0.5] }, "objective": "sum_{j=1}^{30} A_j * rho_j(x)\nwhere rho_j(x) = - (exp(-mu_j^2 * r) + 2*exp(-mu_j^2*(x2^2 + x3^2)) + 2*exp(-mu_j^2*x3^2) + 1) / mu_j^2\nr = x1^2 + x2^2 + x3^2\nA_j = 2*sin(mu_j)/(mu_j + sin(mu_j)*cos(mu_j))\nmu_j: first 30 positive roots of tan(mu) = mu", "constraints": [ diff --git a/src/ADNLPProblems/hs85.jl b/src/ADNLPProblems/hs85.jl index 135620684..e9c666f18 100644 --- a/src/ADNLPProblems/hs85.jl +++ b/src/ADNLPProblems/hs85.jl @@ -47,7 +47,7 @@ function hs85(; type::Type{T} = Float64, kwargs...) where {T} x0 = T[900.0, 80.0, 115.0, 267.0, 27.0] # Best known value ≈ -1.90513375 - function f(x::AbstractVector{T}) + function f(x) # All intermediates (identical to those used in constraints) y1 = x[2] + x[3] + T(41.6) c1 = T(0.024) * x[4] - T(4.62) @@ -96,7 +96,7 @@ function hs85(; type::Type{T} = Float64, kwargs...) where {T} end # Constraint function (38 nonlinear inequalities, all of the form c(x) >= 0) - function c!(cx::AbstractVector{T}, x::AbstractVector{T}) + function c!(cx, x) # All intermediates (identical to those used in objective) y1 = x[2] + x[3] + T(41.6) c1 = T(0.024) * x[4] - T(4.62) @@ -132,50 +132,74 @@ function hs85(; type::Type{T} = Float64, kwargs...) where {T} c15 = y13 / y15 - y13 / T(0.52) c16 = T(1.104) - T(0.72) * y15 c17 = y9 + x[5] - # Constraints - cx[1] = T(1.5) * x[2] - x[3] - cx[2] = y1 - T(213.1) - cx[3] = T(405.23) - y1 - cx[4] = y2 - a[2] - cx[5] = y3 - a[3] - cx[6] = y4 - a[4] - cx[7] = y5 - a[5] - cx[8] = y6 - a[6] - cx[9] = y7 - a[7] - cx[10] = y8 - a[8] - cx[11] = y9 - a[9] - cx[12] = y10 - a[10] - cx[13] = y11 - a[11] - cx[14] = y12 - a[12] - cx[15] = y13 - a[13] - cx[16] = y14 - a[14] - cx[17] = y15 - a[15] - cx[18] = y16 - a[16] - cx[19] = y17 - a[17] - cx[20] = b[2] - y2 - cx[21] = b[3] - y3 - cx[22] = b[4] - y4 - cx[23] = b[5] - y5 - cx[24] = b[6] - y6 - cx[25] = b[7] - y7 - cx[26] = b[8] - y8 - cx[27] = b[9] - y9 - cx[28] = b[10] - y10 - cx[29] = b[11] - y11 - cx[30] = b[12] - y12 - cx[31] = b[13] - y13 - cx[32] = b[14] - y14 - cx[33] = b[15] - y15 - cx[34] = b[16] - y16 - cx[35] = b[17] - y17 - cx[36] = y4 - (T(0.28) / T(0.72)) * y5 - cx[37] = T(21) - T(3496) * y2 / c12 - cx[38] = T(62212) / c17 - T(110.6) - y1 + # 35 nonlinear constraints: cx[1:18] of the form c(x) >= 0, cx[19:35] of the form c(x) <= 0. + # The 3 remaining (linear) constraints of HS85 are passed separately below. + cx[1] = y2 - a[2] + cx[2] = y3 - a[3] + cx[3] = y4 - a[4] + cx[4] = y5 - a[5] + cx[5] = y6 - a[6] + cx[6] = y7 - a[7] + cx[7] = y8 - a[8] + cx[8] = y9 - a[9] + cx[9] = y10 - a[10] + cx[10] = y11 - a[11] + cx[11] = y12 - a[12] + cx[12] = y13 - a[13] + cx[13] = y14 - a[14] + cx[14] = y15 - a[15] + cx[15] = y16 - a[16] + cx[16] = y17 - a[17] + cx[17] = y4 - (T(0.28) / T(0.72)) * y5 + cx[18] = T(62212) / c17 - T(110.6) - y1 + cx[19] = y2 - b[2] + cx[20] = y3 - b[3] + cx[21] = y4 - b[4] + cx[22] = y5 - b[5] + cx[23] = y6 - b[6] + cx[24] = y7 - b[7] + cx[25] = y8 - b[8] + cx[26] = y9 - b[9] + cx[27] = y10 - b[10] + cx[28] = y11 - b[11] + cx[29] = y12 - b[12] + cx[30] = y13 - b[13] + cx[31] = y14 - b[14] + cx[32] = y15 - b[15] + cx[33] = y16 - b[16] + cx[34] = y17 - b[17] + cx[35] = T(3496) * y2 / c12 - T(21) return cx end - # Constraint bounds: all inequalities of the form c(x) >= 0 + + # 3 linear constraints, given as a sparse Jacobian (rows, cols, vals): + # cx[1] = 1.5 x2 - x3 (>= 0), cx[2] = cx[3] = x2 + x3 (bounded below and above respectively) + clinrows = [1, 1, 2, 2, 3, 3] + clincols = [2, 3, 2, 3, 2, 3] + clinvals = T[1.5, -1, 1, 1, 1, 1] + + # Constraint bounds: overall cx[1:3] are the linear constraints above, + # cx[4:21] are of the form c(x) >= 0, cx[22:38] are of the form c(x) <= 0. m = 38 cl = zeros(T, m) cu = fill(T(Inf), m) - return ADNLPModels.ADNLPModel!(f, x0, lvar, uvar, c!, cl, cu; name = "hs85", kwargs...) + cl[2], cu[3] = T(171.5), T(363.63) + cl[3] = T(-Inf) + for i = 22:38 + cl[i], cu[i] = T(-Inf), zero(T) + end + return ADNLPModels.ADNLPModel!( + f, + x0, + lvar, + uvar, + clinrows, + clincols, + clinvals, + c!, + cl, + cu; + name = "hs85", + kwargs..., + ) end diff --git a/src/ADNLPProblems/hs89.jl b/src/ADNLPProblems/hs89.jl index 3c360ee0d..8ea0b1666 100644 --- a/src/ADNLPProblems/hs89.jl +++ b/src/ADNLPProblems/hs89.jl @@ -39,7 +39,7 @@ function hs89(; type::Type{T} = Float64, kwargs...) where {T} A = [2 * sin(mu[j]) / (mu[j] + sin(mu[j]) * cos(mu[j])) for j = 1:30] # Objective: φ(x) = ∑_{j=1}^{30} A_j ρ_j(x) - function f(x::AbstractVector{T}) + function f(x) s = zero(T) r = x[1]^2 + x[2]^2 + x[3]^2 for j = 1:30 @@ -53,7 +53,7 @@ function hs89(; type::Type{T} = Float64, kwargs...) where {T} # Equality constraint c(x) = 0 # Full expression with cross terms (double sum over i < j) - function c!(cx::AbstractVector{T}, x::AbstractVector{T}) + function c!(cx, x) r = x[1]^2 + x[2]^2 + x[3]^2 termA = zero(T) termB = zero(T) diff --git a/src/Meta/hs85.jl b/src/Meta/hs85.jl index 0cbd12a77..ae4a2296f 100644 --- a/src/Meta/hs85.jl +++ b/src/Meta/hs85.jl @@ -17,7 +17,9 @@ hs85_meta = Dict( :defined_everywhere => missing, :origin => :unknown, ) -get_hs85_nvar(; n::Integer = 5, kwargs...) = 5 -get_hs85_ncon(; n::Integer = 38, kwargs...) = 38 -get_hs85_nlin(; n::Integer = 38, kwargs...) = 0 -get_hs85_nnln(; n::Integer = 38, kwargs...) = 38 +get_hs85_nvar(; n::Integer = default_nvar, kwargs...) = 5 +get_hs85_ncon(; n::Integer = default_nvar, kwargs...) = 38 +get_hs85_nlin(; n::Integer = default_nvar, kwargs...) = 3 +get_hs85_nnln(; n::Integer = default_nvar, kwargs...) = 35 +get_hs85_nequ(; n::Integer = default_nvar, kwargs...) = 0 +get_hs85_nineq(; n::Integer = default_nvar, kwargs...) = 38 diff --git a/src/Meta/hs89.jl b/src/Meta/hs89.jl index 5163b9bb9..2dde53975 100644 --- a/src/Meta/hs89.jl +++ b/src/Meta/hs89.jl @@ -17,7 +17,9 @@ hs89_meta = Dict( :defined_everywhere => missing, :origin => :unknown, ) -get_hs89_nvar(; n::Integer = 3, kwargs...) = 3 -get_hs89_ncon(; n::Integer = 1, kwargs...) = 1 -get_hs89_nlin(; n::Integer = 1, kwargs...) = 0 -get_hs89_nnln(; n::Integer = 1, kwargs...) = 1 +get_hs89_nvar(; n::Integer = default_nvar, kwargs...) = 3 +get_hs89_ncon(; n::Integer = default_nvar, kwargs...) = 1 +get_hs89_nlin(; n::Integer = default_nvar, kwargs...) = 0 +get_hs89_nnln(; n::Integer = default_nvar, kwargs...) = 1 +get_hs89_nequ(; n::Integer = default_nvar, kwargs...) = 1 +get_hs89_nineq(; n::Integer = default_nvar, kwargs...) = 0 diff --git a/src/PureJuMP/hs85.jl b/src/PureJuMP/hs85.jl index f9db0bfce..23a2aa793 100644 --- a/src/PureJuMP/hs85.jl +++ b/src/PureJuMP/hs85.jl @@ -115,7 +115,7 @@ function hs85(args...; kwargs...) @constraint(m, 62212 / c17 - 110.6 >= y1) # Objective - @NLobjective( + @objective( m, Min, -5.843e-7 * y17 + diff --git a/src/PureJuMP/hs89.jl b/src/PureJuMP/hs89.jl index b03773f3b..0a5d20238 100644 --- a/src/PureJuMP/hs89.jl +++ b/src/PureJuMP/hs89.jl @@ -14,12 +14,7 @@ export hs89 "HS89 model" -function hs89( - args...; - optimizer = nothing, - optimizer_attributes = nothing, - kwargs..., -) +function hs89(args...; optimizer = nothing, optimizer_attributes = nothing, kwargs...) model = optimizer === nothing ? Model() : Model(optimizer) if optimizer !== nothing && optimizer_attributes === nothing @@ -78,16 +73,16 @@ function hs89( @expression( model, ρ[j = 1:30], - let μ² = mu[j]^2, r = x1^2 + x2^2 + x3^2, r23 = x2^2 + x3^2, r3 = x3^2; + let μ² = mu[j]^2, r = x1^2 + x2^2 + x3^2, r23 = x2^2 + x3^2, r3 = x3^2 -(exp(-μ² * r) + 2*exp(-μ² * r23) + 2*exp(-μ² * r3) + 1) / μ² end ) # Objective: ∑ Aⱼ ρⱼ - @NLobjective(model, Min, sum(A[j] * ρ[j] for j = 1:30)) + @objective(model, Min, sum(A[j] * ρ[j] for j = 1:30)) # Constraint: termA + termB = 2/15 - @NLconstraint( + @constraint( model, eq, sum(A[j]^2 * ρ[j]^2 * (sin(2*mu[j])/(2*mu[j]) + 1)/2 for j = 1:30) + sum( From ddda7ddf18429e58b58ca0a7d7f12cf6f41f26b1 Mon Sep 17 00:00:00 2001 From: arnavk23 Date: Fri, 4 Sep 2026 04:01:31 +0530 Subject: [PATCH 11/11] Turn check_rules.jl into a real automated reviewer check --- docs/Project.toml | 1 + docs/RULES_hs85_hs89.md | 44 +++-- docs/check_rules.jl | 224 +++++++++++++++------ docs/hs85_workflow.ipynb | 410 ++++++++++++++++++++++++++++++--------- src/PureJuMP/hs85.jl | 2 +- 5 files changed, 511 insertions(+), 170 deletions(-) diff --git a/docs/Project.toml b/docs/Project.toml index 3389bf931..180b42f57 100644 --- a/docs/Project.toml +++ b/docs/Project.toml @@ -2,6 +2,7 @@ ADNLPModels = "54578032-b7ea-4c30-94aa-7cbd1cce6c9a" Documenter = "e30172f5-a6a5-5a46-863b-614d45cd2de4" GR = "28b8d3ca-fb5f-59d9-8090-bfdbd6d07a71" +JSON = "682c06a0-de6a-54ab-a142-c8b1cf79cde6" JSOSolvers = "10dff2fc-5484-5881-a0e0-c90441020f8a" NLPModels = "a4795742-8479-5a88-8948-cc11e1c8c1a6" NLPModelsJuMP = "792afdf1-32c1-5681-94e0-d7bf7a5df49e" diff --git a/docs/RULES_hs85_hs89.md b/docs/RULES_hs85_hs89.md index 6b8555a9f..2ed1147f9 100644 --- a/docs/RULES_hs85_hs89.md +++ b/docs/RULES_hs85_hs89.md @@ -1,22 +1,32 @@ # OptimizationProblems.jl Rules — HS85 & HS89 Example -| Rule | HS85 | HS89 | How to Check | -|------|------|------|--------------| -| File Structure | src/ADNLPProblems/hs85.jl, src/PureJuMP/hs85.jl, src/Meta/hs85.jl | src/ADNLPProblems/hs89.jl, src/PureJuMP/hs89.jl, src/Meta/hs89.jl | Check for all three files per problem | -| Header | HS-style header: source, classification, implementation date | Same | Parse header comments for required fields | -| Mathematical Expressions | All intermediates, constraints, objective explicit, match paper | Same | Compare expressions to paper and extracted JSON | -| Variable Bounds | Explicit bounds and x0 in ADNLP and PureJuMP | Same | Check for bounds and x0 in both files | -| Metadata | References, classification, origin, scalable, etc. in Meta file | Same | Parse Meta file for required metadata | -| Naming | hs85.jl in ADNLPProblems, hs85.jl in PureJuMP, hs85.jl in Meta | hs89.jl in ADNLPProblems, hs89.jl in PureJuMP, hs89.jl in Meta | Check file/function names against repository conventions | -| Allocation | Constraint function minimizes allocations (relaxed for complex) | Same | Run allocation test or allow allocation | -| Ipopt Solve | Problem solves with Ipopt (PureJuMP) | Same | Run Ipopt solver and check for solution | -| Reviewer Markdown | Summary, PDF screenshot, extraction uncertainties, test results | Same | Generate reviewer markdown file | -| Duplication | No duplicate (by name, structure, metadata) | Same | Check for similar names/metadata | -| Traceability | Origin clear and referenced in header/Meta | Same | Check Meta/header | -| Scalability | Marked if scalable (not for hs85/hs89) | Same | Check Meta/implementation | -| Multiple Problems | One PDF = one problem (for hs85/hs89) | Same | Allow multiple if needed | -| Uncertainty | Warn if extraction unclear (not for hs85/hs89) | Same | Parse extraction JSON for uncertainties | +This table doubles as the specification for `docs/check_rules.jl`: every row with a +`✓` in "Automated" is checked by that script (`julia --project=docs docs/check_rules.jl +...`); everything else needs a human reviewer. + +| Rule | Automated | HS85 | HS89 | How to Check | +|------|:---:|------|------|--------------| +| File Structure | ✓ | src/ADNLPProblems/hs85.jl, src/PureJuMP/hs85.jl, src/Meta/hs85.jl | src/ADNLPProblems/hs89.jl, src/PureJuMP/hs89.jl, src/Meta/hs89.jl | Check for all three files per problem | +| Header | ✓ | Source reference in the PureJuMP file header | Same | Parse header comments for a "Source" line | +| Mathematical Expressions | | All intermediates, constraints, objective explicit, match paper | Same | Compare expressions to paper and extracted JSON | +| Variable Bounds | ✓ | Explicit bounds and x0 in ADNLP (unless `:has_bounds => false`) | Same | Check for `lvar`/`uvar`/`x0` in the ADNLP file | +| Metadata | ✓ | `:name`, `:best_known_upper_bound`, and all of `get_hs85_{nvar,ncon,nlin,nnln,nequ,nineq}` in Meta file | Same | Parse Meta file for required fields and getters | +| Naming | ✓ (via File Structure) | hs85.jl in ADNLPProblems, hs85.jl in PureJuMP, hs85.jl in Meta | hs89.jl in ADNLPProblems, hs89.jl in PureJuMP, hs89.jl in Meta | Check file/function names against repository conventions | +| JuMP Interface Consistency | ✓ | PureJuMP file uses only the modern interface (`@objective`/`@constraint`/`@expression`), never mixed with `@NLobjective`/`@NLconstraint`/`@NLexpression` | Same | Grep for `@NL*` macros in the PureJuMP file — JuMP refuses to build a model mixing both, so this crashes `MathOptNLPModel` | +| AD-Generic Signatures | ✓ | ADNLPProblems objective/constraint closures are not typed to a concrete `x::AbstractVector{T}` | Same | Grep the ADNLP file for `::AbstractVector{T}` argument annotations — these silently block ForwardDiff | +| ADNLP/PureJuMP Compatibility | ✓ | `nvar`, `x0`, `ncon`, objective, `cons()`, `lcon`/`ucon`, and `lin` agree between the two formulations at `x0` | Same | Instantiate both models and compare (this is `test/test-utils.jl`'s `test_compatibility`, run early) | +| Allocation | ✓ (informational) | `cons_nln!` allocation count reported | Same | Run `cons_nln!` twice and check `@allocated` is 0 | +| Ipopt Solve | manual | Problem solves with Ipopt (PureJuMP) | Same | Run Ipopt solver and check for solution | +| Reviewer Markdown | ✓ (presence only) | Summary, PDF screenshot, extraction uncertainties, test results | Same | Generate reviewer markdown file | +| Duplication | | No duplicate (by name, structure, metadata) | Same | Check for similar names/metadata | +| Traceability | ✓ (via Header) | Origin clear and referenced in header/Meta | Same | Check Meta/header | +| Scalability | | Marked if scalable (not for hs85/hs89) | Same | Check Meta/implementation | +| Multiple Problems | | One PDF = one problem (for hs85/hs89) | Same | Allow multiple if needed | +| Uncertainty | | Warn if extraction unclear (not for hs85/hs89) | Same | Parse extraction JSON for uncertainties | --- -This table is tailored to HS85 and HS89. Use as a checklist for similar problems and PRs. +This table is tailored to HS85 and HS89, but the automated rows and +`docs/check_rules.jl` are written to run on any problem name. +They were used as-is to catch a JuMP-model-construction crash, a ForwardDiff-breaking type annotation, and a full constraint-by-constraint ADNLP/PureJuMP mismatch in this +PR before merge. Use this table as a checklist for similar problems and PRs. diff --git a/docs/check_rules.jl b/docs/check_rules.jl index 59465e04d..cd96964ae 100644 --- a/docs/check_rules.jl +++ b/docs/check_rules.jl @@ -1,73 +1,171 @@ # Script to check OptimizationProblems.jl rules for HS85/HS89 using JSON -function check_problem(problem::String) - result = Dict() - adnlp_path = "src/ADNLPProblems/$problem.jl" - jump_path = "src/PureJuMP/$problem.jl" - meta_path = "src/Meta/$problem.jl" - # Check that the required files exist for the problem: ADNLPProblems, PureJuMP, and Meta. - result["ADNLP"] = isfile(adnlp_path) - result["PureJuMP"] = isfile(jump_path) - result["Meta"] = isfile(meta_path) - # Check that the header of the ADNLPProblems file for this problem contains a reference to "Hock and Schittkowski" in the first 10 lines, to ensure proper attribution/documentation. - header_ok = false - if result["ADNLP"] || result["PureJuMP"] - adnlp_header_ok = false - jump_header_ok = false - if result["ADNLP"] - adnlp_lines = readlines(adnlp_path) - k = min(length(adnlp_lines), 20) - adnlp_header_ok = any(occursin("Hock and Schittkowski", adnlp_lines[i]) for i in 1:k) - end - if result["PureJuMP"] - jump_lines = readlines(jump_path) - k = min(length(jump_lines), 20) - jump_header_ok = any(occursin("Hock and Schittkowski", jump_lines[i]) for i in 1:k) - end - header_ok = adnlp_header_ok || jump_header_ok - end - result["Header"] = header_ok - # Check that the ADNLPProblems file for problem defines lower bounds (lvar), upper bounds (uvar), and an initial point (x0). - bounds_ok = false - if result["ADNLP"] - adnlp_file = String(read(adnlp_path)) - has_x0 = occursin("x0", adnlp_file) - has_bounds_in_code = occursin("lvar", adnlp_file) && occursin("uvar", adnlp_file) - has_bounds_required = true - if result["Meta"] - meta_file = String(read(meta_path)) - if occursin(":has_bounds => false", meta_file) - has_bounds_required = false - end - end - bounds_ok = has_x0 && (has_bounds_in_code || !has_bounds_required) +# All paths are resolved relative to the repository root (the parent of this file's +# directory), not the process's current working directory: this script is meant to +# be run both as `julia docs/check_rules.jl` from the repo root and `include()`d +# from docs/hs85_workflow.ipynb, whose cwd is docs/. +const REPO_ROOT = normpath(joinpath(@__DIR__, "..")) + +# Static (text-based) checks: do not require loading the package, so they run even +# if the package fails to precompile. +function check_problem_static(problem::String) + result = Dict{String, Any}() + adnlp_path = joinpath(REPO_ROOT, "src", "ADNLPProblems", "$problem.jl") + jump_path = joinpath(REPO_ROOT, "src", "PureJuMP", "$problem.jl") + meta_path = joinpath(REPO_ROOT, "src", "Meta", "$problem.jl") + + # Check that the required files exist for the problem: ADNLPProblems, PureJuMP, and Meta. + result["ADNLP"] = isfile(adnlp_path) + result["PureJuMP"] = isfile(jump_path) + result["Meta"] = isfile(meta_path) + + adnlp_text = result["ADNLP"] ? read(adnlp_path, String) : "" + jump_text = result["PureJuMP"] ? read(jump_path, String) : "" + meta_text = result["Meta"] ? read(meta_path, String) : "" + + # Check that the header comment of the PureJuMP file (the file that consistently + # carries one across the repository) has a "Source" line in the first 20 lines, + # so the origin of the problem is traceable. + jump_lines = result["PureJuMP"] ? readlines(jump_path) : String[] + result["Header"] = + any(occursin("Source", jump_lines[i]) for i = 1:min(length(jump_lines), 20)) + + # Check that the ADNLPProblems file defines an initial point (x0), and lower/upper + # bounds (lvar/uvar) unless the Meta file explicitly says the problem has none. + has_bounds_required = !occursin(":has_bounds => false", meta_text) + result["Bounds"] = + result["ADNLP"] && + occursin("x0", adnlp_text) && + (!has_bounds_required || (occursin("lvar", adnlp_text) && occursin("uvar", adnlp_text))) + + # Check that the Meta file names the problem, records a best known upper bound, and + # defines every getter required by docs/src/contributing.md + # (get__nvar/ncon/nlin/nnln/nequ/nineq). + required_getters = ["nvar", "ncon", "nlin", "nnln", "nequ", "nineq"] + missing_getters = [g for g in required_getters if !occursin("get_$(problem)_$(g)(", meta_text)] + result["Metadata"] = + result["Meta"] && + occursin(":name => \"$problem\"", meta_text) && + occursin(":best_known_upper_bound", meta_text) && + isempty(missing_getters) + if !isempty(missing_getters) + result["MissingGetters"] = missing_getters + end + + # Check that the PureJuMP file does not mix the legacy nonlinear interface + # (@NLobjective/@NLconstraint/@NLexpression) with the modern one + # (@objective/@constraint/@expression): JuMP refuses to build a model that uses + # both, so MathOptNLPModel(model) throws instead of failing a specific test. + uses_legacy_macro = + occursin("@NLobjective", jump_text) || + occursin("@NLconstraint", jump_text) || + occursin("@NLexpression", jump_text) + result["JuMPInterfaceConsistent"] = result["PureJuMP"] && !uses_legacy_macro + + # Check that the objective/constraint closures in ADNLPProblems are not typed to a + # concrete vector element type (e.g. `x::AbstractVector{T}`), which would prevent + # ForwardDiff from calling them with Dual-typed vectors and break grad()/jacobian(). + result["ADNLPGenericSignatures"] = + result["ADNLP"] && !occursin(r"\([^)]*::AbstractVector\{T\}", adnlp_text) + + # Reviewer markdown (optional, but useful for PR traceability). + result["ReviewerMarkdown"] = + isfile(joinpath(REPO_ROOT, "docs", "review_$problem.md")) || + isfile(joinpath(REPO_ROOT, "docs", "review", "$problem.md")) + + return result +end + +# Whether OptimizationProblems and its test dependencies can be loaded, computed once. +const HAVE_DYNAMIC = try + @eval using OptimizationProblems, ADNLPModels, NLPModelsJuMP, NLPModels + true +catch err + @warn "Dynamic checks unavailable: could not load OptimizationProblems/ADNLPModels/NLPModelsJuMP/NLPModels ($err)" + false +end + +# Dynamic checks: instantiate the actual models. These catch what static text checks +# cannot, e.g. a JuMP model that throws on construction, or an ADNLPProblems model +# whose objective/constraints numerically disagree with its PureJuMP sibling. +function check_problem_dynamic(problem::String) + out = Dict{String, Any}() + if !HAVE_DYNAMIC + out["Allocation"] = "skipped (package unavailable)" + out["Compatibility"] = "skipped (package unavailable)" + return out + end + prob = Symbol(problem) + ad_mod = OptimizationProblems.ADNLPProblems + jump_mod = OptimizationProblems.PureJuMP + + if !isdefined(ad_mod, prob) + out["Allocation"] = "skipped (no ADNLPProblems.$problem)" + else + try + nlp = getfield(ad_mod, prob)(matrix_free = true) + if nlp.meta.nnln > 0 + x = nlp.meta.x0 + cx = similar(x, nlp.meta.nnln) + NLPModels.cons_nln!(nlp, x, cx) + nbytes = @allocated NLPModels.cons_nln!(nlp, x, cx) + out["Allocation"] = nbytes == 0 ? "pass (0 bytes)" : "warn ($(nbytes) bytes)" + else + out["Allocation"] = "n/a (no nonlinear constraints)" + end + catch err + out["Allocation"] = "error: $err" end - result["Bounds"] = bounds_ok - # Check that the Meta file for this problem contains both a "Source" and a "classification" entry. - meta_ok = false - if result["Meta"] - meta = String(read(meta_path)) - meta_ok = occursin(":name => \"$problem\"", meta) && occursin(":best_known_upper_bound", meta) + end + + if !isdefined(ad_mod, prob) || !isdefined(jump_mod, prob) + out["Compatibility"] = "skipped (missing ADNLPProblems or PureJuMP implementation)" + else + try + nlp_ad = getfield(ad_mod, prob)(matrix_free = true) + model = getfield(jump_mod, prob)() + nlp_jump = NLPModelsJuMP.MathOptNLPModel(model; name = problem) + x0 = nlp_ad.meta.x0 + + failures = String[] + nlp_jump.meta.nvar == nlp_ad.meta.nvar || push!(failures, "nvar mismatch") + nlp_jump.meta.x0 == nlp_ad.meta.x0 || push!(failures, "x0 mismatch") + nlp_jump.meta.ncon == nlp_ad.meta.ncon || push!(failures, "ncon mismatch") + isapprox(NLPModels.obj(nlp_ad, x0), NLPModels.obj(nlp_jump, x0), rtol = 1e-6) || + push!(failures, "objective mismatch at x0") + if nlp_ad.meta.ncon > 0 && nlp_jump.meta.ncon == nlp_ad.meta.ncon + nlp_ad.meta.lcon ≈ nlp_jump.meta.lcon || push!(failures, "lcon mismatch") + nlp_ad.meta.ucon ≈ nlp_jump.meta.ucon || push!(failures, "ucon mismatch") + all(isapprox.(NLPModels.cons(nlp_ad, x0), NLPModels.cons(nlp_jump, x0), atol = 1e-6)) || + push!(failures, "cons() mismatch at x0") + nlp_ad.meta.lin == nlp_jump.meta.lin || push!(failures, "lin (linear constraint indices) mismatch") + end + out["Compatibility"] = isempty(failures) ? "pass" : "fail: " * join(failures, "; ") + catch err + out["Compatibility"] = "error: $err" end - result["Metadata"] = meta_ok - # Allocation check is currently not implemented (set to false as a placeholder). - result["Allocation"] = false - # Ipopt solve check is not performed by this script. - # Marked as "manual" for reviewer attention. - result["IpoptSolve"] = "manual" - # Check whether reviewer markdown exists (optional, but useful for PR traceability). - result["ReviewerMarkdown"] = isfile("docs/review_$problem.md") || isfile("docs/review/$problem.md") - return result + end + return out +end + +function check_problem(problem::String) + result = check_problem_static(problem) + merge!(result, check_problem_dynamic(problem)) + # Solving with Ipopt is not run automatically; left for manual reviewer attention. + result["IpoptSolve"] = "manual" + return result end function main() - problems = isempty(ARGS) ? ["hs85", "hs89"] : ARGS - results = Dict() - for p in problems - results[p] = check_problem(p) - end - println(JSON.json(results, 2)) + problems = isempty(ARGS) ? ["hs85", "hs89"] : ARGS + results = Dict() + for p in problems + results[p] = check_problem(p) + end + println(JSON.json(results, 2)) end -main() +if abspath(PROGRAM_FILE) == @__FILE__ + main() +end diff --git a/docs/hs85_workflow.ipynb b/docs/hs85_workflow.ipynb index 2c4671c4f..458462304 100644 --- a/docs/hs85_workflow.ipynb +++ b/docs/hs85_workflow.ipynb @@ -11,29 +11,64 @@ "\n", "Flow covered:\n", "1. Load extraction artifact (`docs/extract_hs85.json`).\n", - "2. Run repository rule checks (files, header, bounds, metadata).\n", + "2. Run the repository rule checks from `docs/check_rules.jl` (files, header, bounds,\n", + " metadata, JuMP interface consistency, AD-generic signatures).\n", "3. Compare extracted bounds with implementation bounds.\n", "4. Instantiate the ADNLP problem and inspect dimensions.\n", - "5. Optionally solve the PureJuMP model with Ipopt (if available).\n", + "5. Check that ADNLPProblems and PureJuMP agree (objective, constraints, bounds,\n", + " linear/nonlinear split) and that `cons_nln!` is allocation-free.\n", + "6. Optionally solve the PureJuMP model with Ipopt (if available).\n", "\n", - "The goal is to make the process reviewable and reproducible in one place." + "The goal is to make the process reviewable and reproducible in one place. Steps 2\n", + "and 5 are not cosmetic: running them against the first draft of HS85/HS89 in this\n", + "PR caught a JuMP model that crashed on construction (mixed legacy/modern nonlinear\n", + "macros), a type annotation that silently broke ForwardDiff, and a full\n", + "constraint-by-constraint mismatch between the ADNLP and PureJuMP formulations —\n", + "all before a human reviewer had to notice them by hand." ] }, { "cell_type": "code", - "execution_count": null, + "execution_count": 1, "id": "dd6abb58", - "metadata": {}, - "outputs": [], + "metadata": { + "execution": { + "iopub.execute_input": "2026-09-03T22:29:28.728000Z", + "iopub.status.busy": "2026-09-03T22:29:28.728000Z", + "iopub.status.idle": "2026-09-03T22:29:49.085000Z", + "shell.execute_reply": "2026-09-03T22:29:49.085000Z" + } + }, + "outputs": [ + { + "name": "stderr", + "output_type": "stream", + "text": [ + "\u001b[32m\u001b[1m Activating\u001b[22m\u001b[39m project at `C:\\Users\\kapoo\\Downloads\\OptimizationProblems.jl`\n" + ] + }, + { + "name": "stdout", + "output_type": "stream", + "text": [ + "Activated project at: C:\\Users\\kapoo\\Downloads\\OptimizationProblems.jl\\\n" + ] + } + ], "source": [ "using Pkg\n", "\n", - "repo_root = if isfile(joinpath(pwd(), \"Project.toml\"))\n", + "# Look for the repository root specifically (not just any Project.toml — docs/ has\n", + "# its own docs/Project.toml, so a plain `isfile(\"Project.toml\")` check matches there\n", + "# too when the notebook is opened with cwd == docs/).\n", + "is_repo_root(dir) = isfile(joinpath(dir, \"src\", \"OptimizationProblems.jl\"))\n", + "\n", + "repo_root = if is_repo_root(pwd())\n", " pwd()\n", - "elseif isfile(joinpath(pwd(), \"..\", \"Project.toml\"))\n", + "elseif is_repo_root(normpath(joinpath(pwd(), \"..\")))\n", " normpath(joinpath(pwd(), \"..\"))\n", "else\n", - " error(\"Could not locate Project.toml from current working directory: $(pwd())\")\n", + " error(\"Could not locate the OptimizationProblems.jl repository root from current working directory: $(pwd())\")\n", "end\n", "\n", "Pkg.activate(repo_root)\n", @@ -48,10 +83,29 @@ }, { "cell_type": "code", - "execution_count": null, + "execution_count": 2, "id": "cfaf1484", - "metadata": {}, - "outputs": [], + "metadata": { + "execution": { + "iopub.execute_input": "2026-09-03T22:29:50.917000Z", + "iopub.status.busy": "2026-09-03T22:29:49.089000Z", + "iopub.status.idle": "2026-09-03T22:29:53.247000Z", + "shell.execute_reply": "2026-09-03T22:29:53.191000Z" + } + }, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "Loaded extraction file: C:\\Users\\kapoo\\Downloads\\OptimizationProblems.jl\\docs\\extract_hs85.json\n", + "Problem: hs85\n", + "Variables: 5\n", + "Constraints listed in JSON: 48\n", + "Uncertainties count: 0\n" + ] + } + ], "source": [ "extract_path = joinpath(repo_root, \"docs\", \"extract_hs85.json\")\n", "extract = JSON.parsefile(extract_path)\n", @@ -65,60 +119,69 @@ }, { "cell_type": "code", - "execution_count": null, + "execution_count": 3, "id": "7654a919", - "metadata": {}, - "outputs": [], + "metadata": { + "execution": { + "iopub.execute_input": "2026-09-03T22:29:53.250000Z", + "iopub.status.busy": "2026-09-03T22:29:53.249000Z", + "iopub.status.idle": "2026-09-03T22:30:17.573000Z", + "shell.execute_reply": "2026-09-03T22:30:17.573000Z" + } + }, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "{\n", + " \"ADNLP\": true,\n", + " \"ADNLPGenericSignatures\": true,\n", + " \"Allocation\": \"pass (0 bytes)\",\n", + " \"Bounds\": true,\n", + " \"Compatibility\": \"pass\",\n", + " \"Header\": true,\n", + " \"IpoptSolve\": \"manual\",\n", + " \"JuMPInterfaceConsistent\": true,\n", + " \"Meta\": true,\n", + " \"Metadata\": true,\n", + " \"PureJuMP\": true,\n", + " \"ReviewerMarkdown\": false\n", + "}\n" + ] + } + ], "source": [ - "function first_lines_have_pattern(path::String, pattern::String; n::Int = 20)\n", - " if !isfile(path)\n", - " return false\n", - " end\n", - " lines = readlines(path)\n", - " k = min(length(lines), n)\n", - " return any(occursin(pattern, lines[i]) for i in 1:k)\n", - "end\n", - "\n", - "function check_problem_files_and_rules(repo_root::String, problem::String)\n", - " adnlp_path = joinpath(repo_root, \"src\", \"ADNLPProblems\", \"$(problem).jl\")\n", - " jump_path = joinpath(repo_root, \"src\", \"PureJuMP\", \"$(problem).jl\")\n", - " meta_path = joinpath(repo_root, \"src\", \"Meta\", \"$(problem).jl\")\n", - "\n", - " file_ok = Dict(\n", - " \"ADNLP\" => isfile(adnlp_path),\n", - " \"PureJuMP\" => isfile(jump_path),\n", - " \"Meta\" => isfile(meta_path),\n", - " )\n", - "\n", - " header_ok = first_lines_have_pattern(adnlp_path, \"Hock and Schittkowski\")\n", - "\n", - " adnlp_text = file_ok[\"ADNLP\"] ? read(adnlp_path, String) : \"\"\n", - " bounds_ok = occursin(\"lvar\", adnlp_text) && occursin(\"uvar\", adnlp_text) && occursin(\"x0\", adnlp_text)\n", + "include(joinpath(repo_root, \"docs\", \"check_rules.jl\"))\n", "\n", - " meta_text = file_ok[\"Meta\"] ? read(meta_path, String) : \"\"\n", - " metadata_ok = occursin(\":name => \\\"$(problem)\\\"\", meta_text) &&\n", - " occursin(\":best_known_upper_bound\", meta_text)\n", - "\n", - " return Dict(\n", - " \"problem\" => problem,\n", - " \"files\" => file_ok,\n", - " \"header\" => header_ok,\n", - " \"bounds\" => bounds_ok,\n", - " \"metadata\" => metadata_ok,\n", - " \"all_core_checks_pass\" => all(values(file_ok)) && header_ok && bounds_ok && metadata_ok,\n", - " )\n", - "end\n", - "\n", - "checks = check_problem_files_and_rules(repo_root, \"hs85\")\n", + "checks = check_problem(\"hs85\")\n", "println(JSON.json(checks, 2))" ] }, { "cell_type": "code", - "execution_count": null, + "execution_count": 4, "id": "78215d1d", - "metadata": {}, - "outputs": [], + "metadata": { + "execution": { + "iopub.execute_input": "2026-09-03T22:30:17.576000Z", + "iopub.status.busy": "2026-09-03T22:30:17.576000Z", + "iopub.status.idle": "2026-09-03T22:30:17.872000Z", + "shell.execute_reply": "2026-09-03T22:30:17.872000Z" + } + }, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "Bounds and x0 consistency with extraction JSON:\n", + "lvar match: true\n", + "uvar match: true\n", + "x0 match: true\n" + ] + } + ], "source": [ "function parse_numeric_vector_literal(text::String, name::String)\n", " m = match(Regex(name * raw\"\\s*=\\s*T\\[(.*?)\\]\", \"s\"), text)\n", @@ -151,38 +214,226 @@ }, { "cell_type": "code", - "execution_count": null, + "execution_count": 5, "id": "c8224124", - "metadata": {}, - "outputs": [], + "metadata": { + "execution": { + "iopub.execute_input": "2026-09-03T22:30:17.874000Z", + "iopub.status.busy": "2026-09-03T22:30:17.874000Z", + "iopub.status.idle": "2026-09-03T22:30:18.447000Z", + "shell.execute_reply": "2026-09-03T22:30:18.447000Z" + } + }, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "ADNLP model summary:\n", + "name: hs85\n", + "nvar: 5\n", + "ncon: 38 (3 linear, 35 nonlinear)\n" + ] + }, + { + "name": "stdout", + "output_type": "stream", + "text": [ + "objective at x0: -0.9393968794311149\n", + "constraints satisfied at x0: 38 / 38\n" + ] + } + ], "source": [ "try\n", - " nlp = OptimizationProblems.ADNLPProblems.hs85()\n", + " nlp = OptimizationProblems.ADNLPProblems.hs85(matrix_free = true)\n", " x0 = nlp.meta.x0\n", " f0 = NLPModels.obj(nlp, x0)\n", " c0 = NLPModels.cons(nlp, x0)\n", + " # hs85 mixes linear and nonlinear constraints with different senses (>=, <=, and\n", + " # ranges), so \"feasible\" means each cons() entry falls within its own [lcon, ucon]\n", + " # bound, not a single uniform sign.\n", + " within_bounds = nlp.meta.lcon .- 1e-8 .<= c0 .<= nlp.meta.ucon .+ 1e-8\n", "\n", " println(\"ADNLP model summary:\")\n", " println(\"name: \", nlp.meta.name)\n", " println(\"nvar: \", nlp.meta.nvar)\n", - " println(\"ncon: \", nlp.meta.ncon)\n", + " println(\"ncon: \", nlp.meta.ncon, \" (\", nlp.meta.nlin, \" linear, \", nlp.meta.nnln, \" nonlinear)\")\n", " println(\"objective at x0: \", f0)\n", - " println(\"min constraint value at x0 (for c(x) >= 0 style constraints): \", minimum(c0))\n", + " println(\"constraints satisfied at x0: \", count(within_bounds), \" / \", length(within_bounds))\n", "catch err\n", " println(\"Skipped ADNLP instantiation step. Reason: \", err)\n", "end" ] }, + { + "cell_type": "markdown", + "id": "0b76a11d", + "metadata": {}, + "source": [ + "## ADNLP / PureJuMP compatibility and allocation\n", + "\n", + "This is the check that matters most for correctness: `test/test-utils.jl` requires\n", + "the ADNLPProblems and PureJuMP formulations of a problem to agree exactly on\n", + "`nvar`, `x0`, `ncon`, `lcon`/`ucon`, `lin` (which constraint indices are linear),\n", + "and the numerical value of the objective and constraints at a given point — and it\n", + "requires `cons_nln!` to be allocation-free. `check_problem_dynamic` in\n", + "`docs/check_rules.jl` runs the same comparison so it can be caught here instead of\n", + "in CI. Running it against the first draft of `hs85`/`hs89` in this PR reported\n", + "`Compatibility => \"fail: ...\"` for every mismatched constraint, which is exactly\n", + "how the constraint reordering/sign bug fixed in this PR was tracked down." + ] + }, + { + "cell_type": "code", + "execution_count": 6, + "id": "376b9979", + "metadata": { + "execution": { + "iopub.execute_input": "2026-09-03T22:30:18.450000Z", + "iopub.status.busy": "2026-09-03T22:30:18.450000Z", + "iopub.status.idle": "2026-09-03T22:30:21.544000Z", + "shell.execute_reply": "2026-09-03T22:30:21.544000Z" + } + }, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "hs85" + ] + }, + { + "name": "stdout", + "output_type": "stream", + "text": [ + ": Compatibility = pass, Allocation = pass (0 bytes)\n", + "hs89" + ] + }, + { + "name": "stdout", + "output_type": "stream", + "text": [ + ": Compatibility = pass, Allocation = pass (0 bytes)\n" + ] + } + ], + "source": [ + "for problem in [\"hs85\", \"hs89\"]\n", + " dyn = check_problem_dynamic(problem)\n", + " println(problem, \": Compatibility = \", dyn[\"Compatibility\"], \", Allocation = \", dyn[\"Allocation\"])\n", + "end" + ] + }, { "cell_type": "code", - "execution_count": null, + "execution_count": 7, "id": "5863d263", - "metadata": {}, - "outputs": [], + "metadata": { + "execution": { + "iopub.execute_input": "2026-09-03T22:30:21.547000Z", + "iopub.status.busy": "2026-09-03T22:30:21.547000Z", + "iopub.status.idle": "2026-09-03T22:30:33.507000Z", + "shell.execute_reply": "2026-09-03T22:30:33.507000Z" + } + }, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "\n", + "******************************************************************************\n", + "This program contains Ipopt, a library for large-scale nonlinear optimization.\n", + " Ipopt is released as open source code under the Eclipse Public License (EPL).\n", + " For more information visit https://github.com/coin-or/Ipopt\n", + "******************************************************************************\n", + "\n", + "This is Ipopt version 3.14.19, running with linear solver MUMPS 5.8.2.\n", + "\n", + "Number of nonzeros in equality constraint Jacobian...: 0\n", + "Number of nonzeros in inequality constraint Jacobian.: 119\n", + "Number of nonzeros in Lagrangian Hessian.............: 276\n", + "\n", + "Total number of variables............................: 5\n", + " variables with only lower bounds: 0\n", + " variables with lower and upper bounds: 5\n", + " variables with only upper bounds: 0\n", + "Total number of equality constraints.................: 0\n", + "Total number of inequality constraints...............: 38\n", + " inequality constraints with only lower bounds: 20\n", + " inequality constraints with lower and upper bounds: 0\n", + " inequality constraints with only upper bounds: 18\n", + "\n", + "iter objective inf_pr inf_du lg(mu) ||d|| lg(rg) alpha_du alpha_pr ls\n", + " 0 -9.3939688e-01 0.00e+00 8.89e-01 -1.0 0.00e+00 - 0.00e+00 0.00e+00 0\n", + " 1 -9.3945298e-01 0.00e+00 1.99e-02 -1.0 3.24e+00 - 9.90e-01 1.00e+00f 1\n", + " 2 -9.5732321e-01 0.00e+00 7.99e-04 -1.0 2.54e+02 - 9.78e-01 1.00e+00f 1\n", + " 3 -1.2162010e+00 0.00e+00 1.36e-01 -1.0 1.28e+03 - 7.51e-01 1.00e+00h 1\n", + " 4 -1.3019172e+00 0.00e+00 9.52e-03 -1.0 3.07e+03 - 9.86e-01 1.00e+00h 1\n", + " 5 -1.3363103e+00 0.00e+00 1.01e-04 -1.0 1.21e+03 - 1.00e+00 1.00e+00h 1\n", + " 6 -1.3443210e+00 0.00e+00 1.11e-06 -1.0 1.36e+02 - 1.00e+00 1.00e+00h 1\n", + " 7 -1.5109565e+00 0.00e+00 7.74e-01 -2.5 2.23e+03 - 7.51e-01 1.00e+00f 1\n", + " 8 -1.8685364e+00 4.98e-02 1.49e-01 -2.5 9.84e+03 - 6.03e-01 1.00e+00h 1\n", + " 9 -1.8858514e+00 0.00e+00 3.45e-03 -2.5 5.18e+02 - 1.00e+00 9.47e-01h 1\n", + "iter objective inf_pr inf_du lg(mu) ||d|| lg(rg) alpha_du alpha_pr ls\n", + " 10 -1.8780147e+00 0.00e+00 6.19e-05 -2.5 3.71e+02 - 1.00e+00 1.00e+00h 1\n", + " 11 -1.8799167e+00 0.00e+00 7.61e-06 -2.5 9.96e+01 - 1.00e+00 1.00e+00h 1\n", + " 12 -1.9007897e+00 0.00e+00 1.21e-01 -3.8 7.15e+02 - 9.57e-01 7.25e-01h 1\n", + " 13 -1.9045084e+00 4.95e-06 1.46e-02 -3.8 3.89e+01 - 1.00e+00 9.83e-01h 1\n", + " 14 -1.9042042e+00 0.00e+00 5.31e-08 -3.8 7.53e+00 - 1.00e+00 1.00e+00h 1\n", + " 15 -1.9051058e+00 0.00e+00 6.38e-03 -5.7 1.57e+01 - 9.74e-01 9.28e-01h 1\n", + " 16 -1.9051442e+00 0.00e+00 2.08e-09 -5.7 3.50e-01 - 1.00e+00 1.00e+00h 1\n", + " 17 -1.9051442e+00 0.00e+00 1.84e-11 -5.7 4.50e-03 - 1.00e+00 1.00e+00h 1\n", + " 18 -1.9051553e+00 0.00e+00 6.50e-09 -8.6 1.58e-01 - 1.00e+00 1.00e+00h 1\n", + "\n", + "Number of Iterations....: 18\n", + "\n", + " (scaled) (unscaled)\n", + "Objective...............: -1.9051552886872158e+00 -1.9051552886872158e+00\n", + "Dual infeasibility......: 6.5009405259171260e-09 6.5009405259171260e-09\n", + "Constraint violation....: 0.0000000000000000e+00 0.0000000000000000e+00\n", + "Variable bound violation: 6.0617809083396423e-07 6.0617809083396423e-07\n", + "Complementarity.........: 6.5213860383702599e-09 6.5213860383702599e-09\n", + "Overall NLP error.......: 6.5213860383702599e-09 6.5213860383702599e-09\n", + "\n", + "\n", + "Number of objective function evaluations = 19\n", + "Number of objective gradient evaluations = 19\n", + "Number of equality constraint evaluations = 0\n", + "Number of inequality constraint evaluations = 19\n", + "Number of equality constraint Jacobian evaluations = 0\n", + "Number of inequality constraint Jacobian evaluations = 19\n", + "Number of Lagrangian Hessian evaluations = 18\n", + "Total seconds in IPOPT = 6.413\n", + "\n", + "EXIT: Optimal Solution Found.\n" + ] + }, + { + "name": "stdout", + "output_type": "stream", + "text": [ + "Ipopt termination status: LOCALLY_SOLVED\n" + ] + }, + { + "name": "stdout", + "output_type": "stream", + "text": [ + "Objective value: -1.9051552886872158\n" + ] + } + ], "source": [ "try\n", " using JuMP, Ipopt\n", " model = OptimizationProblems.PureJuMP.hs85()\n", + " # PureJuMP constructors return a bare Model() with no optimizer attached\n", + " # (see docs/RULES_hs85_hs89.md / JuMP Interface Consistency); callers attach one.\n", + " set_optimizer(model, Ipopt.Optimizer)\n", " optimize!(model)\n", " println(\"Ipopt termination status: \", termination_status(model))\n", " println(\"Objective value: \", objective_value(model))\n", @@ -191,40 +442,21 @@ " println(\"This is expected if Ipopt/JuMP are not available in the current notebook environment.\")\n", "end" ] - }, - { - "cell_type": "markdown", - "id": "46e5809a", - "metadata": {}, - "source": [ - "## Manual vs AI Intervention Notes\n", - "\n", - "What this notebook automates well:\n", - "- Structural checks (file existence, headers, basic metadata patterns).\n", - "- Consistency checks between extraction JSON and implementation bounds/x0.\n", - "- Basic model instantiation and optional solve sanity-check.\n", - "\n", - "What still needs manual review:\n", - "- Ambiguous PDF interpretation and formula disambiguation.\n", - "- Naming choices and duplicate detection decisions.\n", - "- Semantic fidelity of objective/constraints against original publication.\n", - "\n", - "This balance matches PR #412 workflow: AI accelerates drafting/checking, while humans validate mathematical intent and project policy decisions." - ] } ], "metadata": { "kernelspec": { - "display_name": "julia 1.11", + "display_name": "julia 1.12", "language": "julia", "name": "julia" }, "language_info": { + "file_extension": ".jl", + "mimetype": "application/julia", "name": "julia", - "version": "3.11.9" + "version": "1.12.5" } }, "nbformat": 4, "nbformat_minor": 5 } - diff --git a/src/PureJuMP/hs85.jl b/src/PureJuMP/hs85.jl index 23a2aa793..a03dca79d 100644 --- a/src/PureJuMP/hs85.jl +++ b/src/PureJuMP/hs85.jl @@ -24,7 +24,7 @@ function hs85(args...; kwargs...) @variable(m, 193.0 ≤ x4 ≤ 287.0966) @variable(m, 25.0 ≤ x5 ≤ 84.1988) - # Intermediates (as @expression or @variable + @NLconstraint) + # Intermediates defined with @expression @expression(m, y1, x2 + x3 + 41.6) @expression(m, c1, 0.024 * x4 - 4.62) @expression(m, y2, 12.5 / c1 + 12)