Solving SOCOPF with Clarabel128 mostly results in ALMOST_OPTIMAL and NEARLY_FEASIBLE_POINT

[ivanightingale]

I noticed that most if not all SOCOPF solutions solved by Clarabel under Float128 precision have the ALMOST_OPTIMAL termination status and the NEARLY_FEASIBLE_POINT primal and dual status.

However, I also noticed that these solutions in fact have very small constraint violation.

Reproducible example:

julia --project=. exp/sampler.jl config_debug.toml 1 64

The content of config_debug.toml is:

pglib_case = "14_ieee"
export_dir = "data/14_ieee_128"

floating_point_type = "float64"

[sampler]

[sampler.load]
noise_type = "ScaledUniform"
l          = 0.80               # Lower bound of base load factor
u          = 1.20               # Upper bound of base load factor
sigma      = 0.20               # Relative (multiplicative) noise level.

[OPF]
# OPF formulations to solve for each sample

[OPF.SOCOPF128]
# These settings are meant to solve `SOCWRConic` problems to high precision
type = "SOCOPF"
solver.name = "Clarabel128"
solver.attributes.max_iter = 2000
solver.attributes.max_step_fraction = 0.995
# Disabling scaling can avoid violations (especially in the dual)
# /!\ This results in worse performance /!\
solver.attributes.equilibrate_enable = true
# Tight tolerances
solver.attributes.tol_gap_abs    = 1e-12
solver.attributes.tol_gap_rel    = 1e-12
solver.attributes.tol_feas       = 1e-12
solver.attributes.tol_infeas_rel = 1e-12
solver.attributes.tol_ktratio    = 1e-10
# Reduced accuracy settings
solver.attributes.reduced_tol_gap_abs    = 1e-8
solver.attributes.reduced_tol_gap_rel    = 1e-8
solver.attributes.reduced_tol_feas       = 1e-8
solver.attributes.reduced_tol_infeas_abs = 1e-8
solver.attributes.reduced_tol_infeas_rel = 1e-8
solver.attributes.reduced_tol_ktratio    = 1e-7
# Additional linear solver-related settings, to match precision of Float128
solver.attributes.static_regularization_enable = false  # typically not needed unless tolerances go below 1e-14
solver.attributes.dynamic_regularization_enable = true  # kept for safety
solver.attributes.dynamic_regularization_eps = 1e-28
solver.attributes.dynamic_regularization_delta = 1e-14
solver.attributes.iterative_refinement_reltol = 1e-18
solver.attributes.iterative_refinement_abstol = 1e-18

Clarabel: v0.11.0

[klamike]

Do you think this is something we can fix here? Or it is a Clarabel issue?

@mtanneau since he wrote the settings

[mtanneau]

The Jabr SOC formulation is known to be numerically shaky, which is why we use Float128 (it's 10x slower but avoids much of the numerical crashes we otherwise see).

The corresponding tolerances (in the snippet shared above) were designed to achieve very high accuracy (ie very low constraint violations and complementarity gaps), hence the 1e-12 feasibility tolerances. The reduced accuracy settings correspond to standard tolerances of 1e-8 when solving a problem in Float64 arithmetic.

I have no issue with relaxing the tolerances so that the return status ends up being OPTIMAL instead of ALMOST_OPTIMAL, as long as they remain in the order of 1e-8.

[ivanightingale]

@copilot Please do the following:

1. Set up the project on your machine.
2. Create config_debug.toml with the following content:

pglib_case = "14_ieee"
export_dir = "data/14_ieee"

floating_point_type = "float64"

[sampler]

[sampler.load]
noise_type = "ScaledUniform"
l          = 0.80               # Lower bound of base load factor
u          = 1.20               # Upper bound of base load factor
sigma      = 0.20               # Relative (multiplicative) noise level.

[OPF]
# OPF formulations to solve for each sample

[OPF.SOCOPF128]
# These settings are meant to solve `SOCWRConic` problems to high precision
type = "SOCOPF"
solver.name = "Clarabel128"
# Tight tolerances
solver.attributes.tol_gap_abs    = 1e-8
solver.attributes.tol_gap_rel    = 1e-8
solver.attributes.tol_feas       = 1e-8
solver.attributes.tol_infeas_rel = 1e-8
solver.attributes.tol_ktratio    = 1e-8
# Reduced accuracy settings
solver.attributes.reduced_tol_gap_abs    = 1e-8
solver.attributes.reduced_tol_gap_rel    = 1e-8
solver.attributes.reduced_tol_feas       = 1e-8
solver.attributes.reduced_tol_infeas_abs = 1e-8
solver.attributes.reduced_tol_infeas_rel = 1e-8
solver.attributes.reduced_tol_ktratio    = 1e-8

3. Execute

julia --project=. exp/sampler.jl config_debug.toml 1 64

4. You should now have an h5 file in data/14_ieee/res_h5. Load the file with HDF5.
5. Check whether the values of the "termination_status" group under the "meta" group are mostly "OPTIMAL" instead of "ALMOST_OPTIMAL".
6. Check whether the values of the "primal_status" group under the "meta" group are mostly "FEASIBLE_POINT".
7. Check whether the values of the "dual_status" group under the "meta" group are mostly "FEASIBLE_POINT".
8. Check that the values of the "primal_objective_value" group and the "dual_objective_value" group under the "meta" group differ by a magnitude no greater than 1e-8 at all indices where the "termination_status" group has value "OPTIMAL".
9. Comment in this issue your results.