Add more checks

docs/src/background/algebraic_modeling_languages.md

@@ -138,6 +138,9 @@ julia> function algebraic_knapsack(c, w, b)
            @objective(model, Max, sum(c[i] * x[i] for i = 1:n))
            @constraint(model, sum(w[i] * x[i] for i = 1:n) <= b)
            optimize!(model)
+           if termination_status(model) != OPTIMAL
+               error("Not solved correctly")
+           end
            return value.(x)
        end
 algebraic_knapsack (generic function with 1 method)
@@ -179,6 +182,9 @@ julia> function nonalgebraic_knapsack(c, w, b)
            con = build_constraint(error, lhs, MOI.LessThan(b))
            add_constraint(model, con)
            optimize!(model)
+           if termination_status(model) != OPTIMAL
+               error("Not solved correctly")
+           end
            return value.(x)
        end
 nonalgebraic_knapsack (generic function with 1 method)
@@ -219,6 +225,9 @@ julia> function mathoptinterface_knapsack(optimizer, c, w, b)
                MOI.LessThan(b),
            )
            MOI.optimize!(model)
+           if MOI.get(model, MOI.TerminationStatus()) != MOI.OPTIMAL
+               error("Not solved correctly")
+           end
            return MOI.get.(model, MOI.VariablePrimal(), x)
        end
 mathoptinterface_knapsack (generic function with 1 method)
@@ -257,6 +266,9 @@ julia> function highs_knapsack(c, w, b)
                w,
            )
            Highs_run(model)
+           if Highs_getModelStatus(model) != kHighsModelStatusOptimal
+               error("Not solved correctly")
+           end
            x = fill(NaN, 2)
            Highs_getSolution(model, x, C_NULL, C_NULL, C_NULL)
            Highs_destroy(model)

docs/src/tutorials/applications/optimal_power_flow.jl

@@ -420,8 +420,8 @@ optimize!(model)
 
 #-
 
-sdp_relaxation_lower_bound = round(objective_value(model); digits = 2)
 Test.@test has_optimal_solution(model; allow_almost = true)
+sdp_relaxation_lower_bound = round(objective_value(model); digits = 2)
 Test.@test isapprox(sdp_relaxation_lower_bound, 2753.04; rtol = 1e-3)     #src
 println(
     "Objective value (W & V relax. lower bound): $sdp_relaxation_lower_bound",

docs/src/tutorials/conic/arbitrary_precision.jl

@@ -76,6 +76,7 @@ print(model)
 # Let's solve and inspect the solution:
 
 optimize!(model)
+@assert has_optimal_solution(model; dual = true)
 solution_summary(model)
 
 # The value of each decision variable is a `BigFloat`:
@@ -100,6 +101,7 @@ value.(x) .- [3 // 7, 3 // 14]
 set_attribute(model, "tol_gap_abs", 1e-32)
 set_attribute(model, "tol_gap_rel", 1e-32)
 optimize!(model)
+@assert has_optimal_solution(model)
 value.(x) .- [3 // 7, 3 // 14]
 
 # ## Rational arithmetic
@@ -141,6 +143,7 @@ print(model)
 # Let's solve and inspect the solution:
 
 optimize!(model)
+@assert has_optimal_solution(model)
 solution_summary(model)
 
 # The optimal values are given in exact rational arithmetic:

docs/src/tutorials/conic/dualization.jl

@@ -130,6 +130,7 @@ print(model_dual)
 
 set_optimizer(model_primal, SCS.Optimizer)
 optimize!(model_primal)
+@assert has_optimal_solution(model_primal; dual = true)
 
 # (There are five rows in the constraint matrix because SCS expects problems in
 # geometric conic form, and so JuMP has reformulated the `X, PSD` variable
@@ -153,6 +154,7 @@ objective_value(model_primal)
 
 set_optimizer(model_dual, SCS.Optimizer)
 optimize!(model_dual)
+@assert has_optimal_solution(model_dual; dual = true)
 
 # and the solution we obtain is:
 
@@ -182,6 +184,7 @@ objective_value(model_dual)
 
 set_optimizer(model_primal, Dualization.dual_optimizer(SCS.Optimizer))
 optimize!(model_primal)
+@assert has_optimal_solution(model_primal; dual = true)
 
 # The performance is the same as if we solved `model_dual`, and the correct
 # solution is returned to `X`:
@@ -197,6 +200,7 @@ dual.(primal_c)
 
 set_optimizer(model_dual, Dualization.dual_optimizer(SCS.Optimizer))
 optimize!(model_dual)
+@assert has_optimal_solution(model_dual; dual = true)
 
 #-
 

docs/src/tutorials/conic/experiment_design.jl

@@ -141,6 +141,7 @@ for i in 1:q
 end
 @objective(aOpt, Min, sum(u))
 optimize!(aOpt)
+@assert has_optimal_solution(aOpt)
 objective_value(aOpt)
 
 #-
@@ -182,6 +183,7 @@ set_silent(eOpt)
 @constraint(eOpt, sum(np) <= n)
 @objective(eOpt, Max, t)
 optimize!(eOpt)
+@assert has_optimal_solution(eOpt)
 objective_value(eOpt)
 #-
 value.(np)
@@ -212,6 +214,7 @@ set_silent(dOpt)
 E = V * LinearAlgebra.diagm(0 => np ./ n) * V'
 @constraint(dOpt, [t; 1; triangle_vec(E)] in MOI.LogDetConeTriangle(q))
 optimize!(dOpt)
+@assert has_optimal_solution(dOpt)
 objective_value(dOpt)
 #-
 value.(np)

docs/src/tutorials/conic/logistic_regression.jl

@@ -195,11 +195,12 @@ X, y = generate_dataset(n, p; shift = 10.0);
 model = build_logit_model(X, y, λ)
 set_optimizer(model, SCS.Optimizer)
 set_silent(model)
-JuMP.optimize!(model)
+optimize!(model)
+@assert has_optimal_solution(model)
 
 #-
 
-θ♯ = JuMP.value.(model[:θ])
+θ♯ = value.(model[:θ])
 
 # It appears that the speed of convergence is not that impacted by the correlation
 # of the dataset, nor by the penalty $\lambda$.
@@ -237,11 +238,12 @@ count_nonzero(v::Vector; tol = 1e-6) = sum(abs.(v) .>= tol)
 sparse_model = build_sparse_logit_model(X, y, λ)
 set_optimizer(sparse_model, SCS.Optimizer)
 set_silent(sparse_model)
-JuMP.optimize!(sparse_model)
+optimize!(sparse_model)
+@assert has_optimal_solution(sparse_model)
 
 #-
 
-θ♯ = JuMP.value.(sparse_model[:θ])
+θ♯ = value.(sparse_model[:θ])
 println(
     "Number of non-zero components: ",
     count_nonzero(θ♯),

docs/src/tutorials/conic/quantum_discrimination.jl

@@ -98,6 +98,7 @@ E = [@variable(model, [1:d, 1:d] in HermitianPSDCone()) for i in 1:N]
 # Now we optimize:
 
 optimize!(model)
+@assert has_optimal_solution(model)
 solution_summary(model)
 
 # The probability of guessing correctly is:
@@ -140,6 +141,7 @@ push!(E, E_N)
 # Then we can check that we get the same solution:
 
 optimize!(model)
+@assert has_optimal_solution(model)
 solution_summary(model)
 
 #-

docs/src/tutorials/conic/simple_examples.jl

@@ -64,6 +64,7 @@ function solve_max_cut_sdp(weights)
     @objective(model, Max, 0.25 * LinearAlgebra.dot(L, X))
     @constraint(model, LinearAlgebra.diag(X) .== 1)
     optimize!(model)
+    @assert has_optimal_solution(model)
     V = svd_cholesky(value(X))
     Random.seed!(N)
     r = rand(N)
@@ -133,6 +134,7 @@ function example_k_means_clustering()
     @constraint(model, [i = 1:m], sum(Z[i, :]) .== 1)
     @constraint(model, LinearAlgebra.tr(Z) == num_clusters)
     optimize!(model)
+    @assert has_optimal_solution(model)
     Z_val = value.(Z)
     current_cluster, visited = 0, Set{Int}()
     solution = [1, 1, 2, 1, 2, 2]  #src
@@ -185,10 +187,12 @@ function example_correlation_problem()
     @constraint(model, 0.4 <= ρ["B", "C"] <= 0.5)
     @objective(model, Max, ρ["A", "C"])
     optimize!(model)
+    @assert has_optimal_solution(model)
     println("An upper bound for ρ_AC is $(value(ρ["A", "C"]))")
     Test.@test value(ρ["A", "C"]) ≈ 0.87195 atol = 1e-4  #src
     @objective(model, Min, ρ["A", "C"])
     optimize!(model)
+    @assert has_optimal_solution(model)
     println("A lower bound for ρ_AC is $(value(ρ["A", "C"]))")
     Test.@test value(ρ["A", "C"]) ≈ -0.978 atol = 1e-3  #src
     return
@@ -380,6 +384,7 @@ function example_robust_uncertainty_sets()
     @constraint(model, [((1-ɛ)/ɛ) (u - μ)'; (u-μ) Σ] >= 0, PSDCone())
     @objective(model, Max, c' * u)
     optimize!(model)
+    @assert has_optimal_solution(model)
     exact =
         μhat' * c +
         Γ1(𝛿 / 2, N) * LinearAlgebra.norm(c) +

docs/src/tutorials/conic/start_values.jl

@@ -67,6 +67,7 @@ model = Model(SCS.Optimizer)
 @constraint(model, sum(x) <= 1)
 @objective(model, Max, sum(i * x[i] for i in 1:3))
 optimize!(model)
+@assert has_optimal_solution(model)
 
 # By looking at the log, we can see that SCS took 75 iterations to find the optimal
 # solution. Now we set the optimal solution as our starting point:

docs/src/tutorials/conic/tips_and_tricks.jl

@@ -98,6 +98,7 @@ set_silent(model)
 @constraint(model, [t; x] in SecondOrderCone())
 @objective(model, Min, t)
 optimize!(model)
+@assert has_optimal_solution(model)
 value(t), value.(x)
 
 # ## Rotated Second-Order Cone
@@ -120,6 +121,7 @@ set_silent(model)
 @constraint(model, [t; 0.5; residuals] in RotatedSecondOrderCone())
 @objective(model, Min, t)
 optimize!(model)
+@assert has_optimal_solution(model)
 value(θ), value(t)
 
 # ## Exponential Cone
@@ -142,6 +144,7 @@ set_silent(model)
 @objective(model, Min, z)
 @constraint(model, [x, 1, z] in MOI.ExponentialCone())
 optimize!(model)
+@assert has_optimal_solution(model)
 value(z), exp(1.5)
 
 # ### Logarithm
@@ -155,6 +158,7 @@ set_silent(model)
 @objective(model, Max, x)
 @constraint(model, [x, 1, z] in MOI.ExponentialCone())
 optimize!(model)
+@assert has_optimal_solution(model)
 value(x), log(1.5)
 
 # ### Log-sum-exp
@@ -216,6 +220,7 @@ set_silent(model)
 @constraint(model, A * x .<= b)
 @constraint(model, [i = 1:n], [t[i], x[i], 1] in MOI.ExponentialCone())
 optimize!(model)
+@assert has_optimal_solution(model)
 objective_value(model)
 
 # The [`MOI.ExponentialCone`](@ref) has a dual, the [`MOI.DualExponentialCone`](@ref),
@@ -234,6 +239,7 @@ set_silent(model)
 @constraint(model, A * x .<= b)
 @constraint(model, [t; ones(n); x] in MOI.RelativeEntropyCone(2n + 1))
 optimize!(model)
+@assert has_optimal_solution(model)
 objective_value(model)
 
 # ## PowerCone
@@ -255,6 +261,7 @@ set_silent(model)
 @constraint(model, [t, 1, x] in MOI.PowerCone(1 / 3))
 @objective(model, Min, t)
 optimize!(model)
+@assert has_optimal_solution(model)
 value(t), value(x)
 
 # The [`MOI.PowerCone`](@ref) has a dual, the [`MOI.DualPowerCone`](@ref),
@@ -277,6 +284,7 @@ function p_norm(x::Vector, p)
     @constraint(model, sum(r) == t)
     @objective(model, Min, t)
     optimize!(model)
+    @assert has_optimal_solution(model)
     return value(t)
 end
 
@@ -320,6 +328,7 @@ set_silent(model)
 @objective(model, Min, t)
 @constraint(model, t .* I - A in PSDCone())
 optimize!(model)
+@assert has_optimal_solution(model)
 objective_value(model)
 
 # ## GeometricMeanCone
@@ -353,6 +362,7 @@ set_silent(model)
 @constraint(model, [t; vec(X)] in MOI.RootDetConeSquare(2))
 @constraint(model, X .== [2 1; 1 3])
 optimize!(model)
+@assert has_optimal_solution(model)
 value(t), sqrt(LinearAlgebra.det(value.(X)))
 
 # If `X` is symmetric, then you can use [`MOI.RootDetConeTriangle`](@ref)
@@ -390,6 +400,7 @@ set_silent(model)
 @constraint(model, X .== [2 1; 1 3])
 @constraint(model, u == 0.5)
 optimize!(model)
+@assert has_optimal_solution(model)
 value(t), 0.5 * log(LinearAlgebra.det(value.(X) ./ 0.5))
 
 # If `X` is symmetric, then you can use [`MOI.LogDetConeTriangle`](@ref)
@@ -410,6 +421,7 @@ set_silent(model)
 @constraint(model, X .== [2 1; 1 3])
 @constraint(model, u == 0.5)
 optimize!(model)
+@assert has_optimal_solution(model)
 value(t), 0.5 * log(LinearAlgebra.det(value.(X) ./ 0.5))
 
 # ## Other Cones and Functions