Dual infeasible problem when using Inf as upper limit

[danigiro]

R version: 4.3.1 (2023-06-16)
OSQP version: 0.6.3.2

I noticed that osqp has problems reaching convergence when the upper limit is considered to be Inf.
I am not sure what the problem is due to, and it does not always occur. It occurs mostly when the variable x in my problem contains very high values.
I leave an example to reproduce it.

library("osqp")
x <- c(24347.5410, 7743.3586, 4957.3505, 5588.0308, 1445.2660, 3002.4793, 643.5833, 240.3561)
A <- matrix(c(1,-1,-1,-1,-1,-1,-1,-1,
              0, 1, 0, 0, 0, 0, 0, 0), 2, byrow = TRUE)
# constraints:  Ax  =>  x[1] = x[2] + x[3] + x[4] + x[5] + x[6] + x[7] + x[8]
#                       x[2] >= 0
W <- diag(c(2280147, 428471, 184644, 472915, 45151, 124546, 25157, 30286))
P <- solve(W)
q <- (-1) * t(P) %*% as.vector(x)
rec1 <- solve_osqp(P, q, A, c(0, rep(0, NROW(A)-1)), c(0, rep(Inf, NROW(A)-1)))
#> -----------------------------------------------------------------
#>   OSQP v0.6.3  -  Operator Splitting QP Solver
#> (c) Bartolomeo Stellato,  Goran Banjac
#> University of Oxford  -  Stanford University 2021
#> -----------------------------------------------------------------
#>   problem:  variables n = 8, constraints m = 2
#> nnz(P) + nnz(A) = 17
#> settings: linear system solver = qdldl,
#> eps_abs = 1.0e-03, eps_rel = 1.0e-03,
#> eps_prim_inf = 1.0e-04, eps_dual_inf = 1.0e-04,
#> rho = 1.00e-01 (adaptive),
#> sigma = 1.00e-06, alpha = 1.60, max_iter = 4000
#> check_termination: on (interval 25),
#> scaling: on, scaled_termination: off
#> warm start: on, polish: off, time_limit: off
#> 
#> iter   objective    pri res    dua res    rho        time
#> 1  -1.9594e+02   1.16e-03   2.17e-02   1.00e-01   6.53e-05s
#> 25  -1.0000e+30   9.16e-07   2.01e-02   1.00e-01   1.17e-04s
#> 
#> status:               dual infeasible
#> number of iterations: 25
#> run time:             1.61e-04s
#> optimal rho estimate: 7.87e-06

rec2 <- solve_osqp(P, q, A, c(0, rep(0, NROW(A)-1)), c(0, rep(9e25, NROW(A)-1)))
#> -----------------------------------------------------------------
#>   OSQP v0.6.3  -  Operator Splitting QP Solver
#> (c) Bartolomeo Stellato,  Goran Banjac
#> University of Oxford  -  Stanford University 2021
#> -----------------------------------------------------------------
#>   problem:  variables n = 8, constraints m = 2
#> nnz(P) + nnz(A) = 17
#> settings: linear system solver = qdldl,
#> eps_abs = 1.0e-03, eps_rel = 1.0e-03,
#> eps_prim_inf = 1.0e-04, eps_dual_inf = 1.0e-04,
#> rho = 1.00e-01 (adaptive),
#> sigma = 1.00e-06, alpha = 1.60, max_iter = 4000
#> check_termination: on (interval 25),
#> scaling: on, scaled_termination: off
#> warm start: on, polish: off, time_limit: off
#> 
#> iter     objective    pri res    dua res         rho          time
#> 1      -1.9594e+02   1.16e-03   2.17e-02     1.00e-01    6.67e-05s
#> 50     -3.6796e+02   7.30e-04   1.70e-09     7.87e-06    1.25e-04s
#> 
#> status:               solved
#> number of iterations: 50
#> optimal objective:    -367.9584
#> run time:             1.58e-04s
#> optimal rho estimate: 1.05e-05

[ben-schwen]

I think this is not a R problem but an OSQP problem, since all big values are replaced by OSQP_INFTY which is 1e30.