Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Jump to bottom

very slow transfer from JuMP to SCS? #266

Closed
kalmarek opened this issue Apr 6, 2023 · 11 comments
Closed

very slow transfer from JuMP to SCS? #266

kalmarek opened this issue Apr 6, 2023 · 11 comments

Comments

@kalmarek
Copy link
Collaborator

kalmarek commented Apr 6, 2023

I'm actually not sure if this is the correct place to report this, but here it is:

My problem has a single psd constraint of size 7_381 with 22_043_332 linear constraints (veeeery sparse ones). After preprocessing (i.e. symmetry reduction:) I get 10 psd constraints of sizes

[334, 332, 182, 167, 165, 153, 459, 447, 449, 439]

and 459_849 (somewhat denser, but still on average 1e-4 to 1e-5 density) constraints.

The latter problem I can create easily in JuMP (it takes about ~75 minutes). However calling optimize! with SCS.Optimizer attached results in

  • long waiting times (>36h),
  • very slowly growing memory consumption (from ~60GB to 128GB, until julia process is killed)
  • uninterruptable julia process (sending SIGINT has no effect so I was not able for now to inspect what is really happening).

However I believe all of this happens before the problem is handled to libscs (e.g. I get no printout of SCS header).

Are these reasonable sizes for JuMP/MOI/bridges (? I'm shooting in the dark here) to behave well? Where should I start debugging?

@odow, @blegat
@kalmarek
Copy link
Collaborator Author

kalmarek commented Apr 6, 2023

Upon killing the process I got this:

QDLDL_factor at /workspace/srcdir/scs/linsys/external/qdldl/qdldl.c:200                                                                                                                      
ldl_factor at /workspace/srcdir/scs/linsys/cpu/direct/private.c:66 [inlined]                                                                                                                 
scs_init_lin_sys_work at /workspace/srcdir/scs/linsys/cpu/direct/private.c:238                                                                                                               
init_work at /workspace/srcdir/scs/src/scs.c:869 [inlined]                                                                                                                                   
scs_init at /workspace/srcdir/scs/src/scs.c:1198                                                                                                                                             
scs_init at /users/kalmar/.julia/packages/SCS/h63sM/src/linear_solvers/direct.jl:25 [inlined]                                                                                                
_unsafe_scs_solve at /users/kalmar/.julia/packages/SCS/h63sM/src/c_wrapper.jl:390                                                                                                            
#scs_solve#13 at /users/kalmar/.julia/packages/SCS/h63sM/src/c_wrapper.jl:349                                                                                                                
scs_solve##kw at /users/kalmar/.julia/packages/SCS/h63sM/src/c_wrapper.jl:278                                                                                                                
unknown function (ip: 0x7f3544ae8be6)                                                                                                                                                        
_jl_invoke at /cache/build/default-amdci4-2/julialang/julia-release-1-dot-8/src/gf.c:2377 [inlined]                                                                                          
ijl_apply_generic at /cache/build/default-amdci4-2/julialang/julia-release-1-dot-8/src/gf.c:2559                                                                                             
optimize! at /users/kalmar/.julia/packages/SCS/h63sM/src/MOI_wrapper/MOI_wrapper.jl:348       
unknown function (ip: 0x7f3544ae41db)                                                                                                                                                        
_jl_invoke at /cache/build/default-amdci4-2/julialang/julia-release-1-dot-8/src/gf.c:2377 [inlined]                                                                                          
ijl_apply_generic at /cache/build/default-amdci4-2/julialang/julia-release-1-dot-8/src/gf.c:2559                                                                                             
optimize! at /users/kalmar/.julia/packages/SCS/h63sM/src/MOI_wrapper/MOI_wrapper.jl:422       
optimize! at /users/kalmar/.julia/packages/MathOptInterface/wx5Ea/src/Utilities/cachingoptimizer.jl:316

so the qdldl factorisation seems to struggle with problems of such sizes:

ok, I think this is more of a scs issue, and no printout wass related to buffered io; I managed to capture the following header when directly logging to a file:

------------------------------------------------------------------
               SCS v3.2.1 - Splitting Conic Solver
        (c) Brendan O'Donoghue, Stanford University, 2012
------------------------------------------------------------------
problem:  variables n: 570684, constraints m: 1030531
cones:    z: primal zero / dual free vars: 459848
          s: psd vars: 570683, ssize: 10
settings: eps_abs: 1.0e-10, eps_rel: 1.0e-10, eps_infeas: 1.0e-07
          alpha: 1.95, scale: 1.00e-01, adaptive_scale: 1
          max_iters: 50000, normalize: 1, rho_x: 1.00e-06
          acceleration_lookback: 50, acceleration_interval: 10
lin-sys:  sparse-direct-amd-qdldl
          nnz(A): 23505312, nnz(P): 0

@bodono
Copy link
Contributor

bodono commented Apr 6, 2023

Try the MKL linear solver backend, it tends to be faster for fractorizing very large matrices.

@kalmarek
Copy link
Collaborator Author

kalmarek commented Apr 6, 2023

@bodono that's what I did!
MKL took just 6min for the initial factorisation :wow: using just a fraction of memory (<9GB). Overall MKL solver is a great success, thanks for wrapping it!

@bodono
Copy link
Contributor

bodono commented Apr 7, 2023

Great, glad to hear it worked well! Was SCS able to solve the problem in the end?

@kalmarek
Copy link
Collaborator Author

kalmarek commented Apr 7, 2023

well, it will take scs some time to finish :)

 iter | pri res | dua res |   gap   |   obj   |  scale  | time (s) 
------------------------------------------------------------------
[...]
63000| 8.41e-05  3.84e-06  1.44e-04 -1.76e+00  2.21e-01  8.92e+04

(to certify the solution I need eps=1e-9 or so)

@kalmarek
Copy link
Collaborator Author

so the convergence of scs is pretty slow, but when I had added (an artificial) bound on the objective <=1.0 (to get a feasible solution) it solved the problem immediately!

------------------------------------------------------------------
               SCS v3.2.1 - Splitting Conic Solver
        (c) Brendan O'Donoghue, Stanford University, 2012
------------------------------------------------------------------
problem:  variables n: 570684, constraints m: 1030532
cones:    z: primal zero / dual free vars: 459848
          l: linear vars: 1
          s: psd vars: 570683, ssize: 10
settings: eps_abs: 1.0e-09, eps_rel: 1.0e-09, eps_infeas: 1.0e-07
          alpha: 1.95, scale: 1.00e-01, adaptive_scale: 1
          max_iters: 100000, normalize: 1, rho_x: 1.00e-06
          acceleration_lookback: 50, acceleration_interval: 10
lin-sys:  sparse-direct-mkl-pardiso
          nnz(A): 23505313, nnz(P): 0
------------------------------------------------------------------
 iter | pri res | dua res |   gap   |   obj   |  scale  | time (s)
------------------------------------------------------------------
     0| 4.00e+01  1.00e+00  2.03e+01 -1.19e+01  1.00e-01  4.78e+02 
   250| 4.53e-03  1.26e-05  4.59e-03 -1.00e+00  2.13e+00  1.75e+03 
   500| 2.57e-07  3.12e-07  2.17e-07 -1.00e+00  2.13e+00  2.13e+03 
   600| 3.08e-09  8.33e-10  6.20e-10 -1.00e+00  6.64e-01  2.73e+03 
------------------------------------------------------------------
status:  solved
timings: total: 2.73e+03s = setup: 4.75e+02s + solve: 2.25e+03s
         lin-sys: 8.48e+02s, cones: 6.93e+01s, accel: 1.79e+00s
------------------------------------------------------------------
objective = -1.000000
------------------------------------------------------------------

@odow
Copy link
Member

odow commented Apr 11, 2023

it solved the problem immediately

time (s)
2.73e+03

I guess "immediately" is relative 😆

@bodono
Copy link
Contributor

bodono commented Apr 12, 2023

so the convergence of scs is pretty slow, but when I had added (an artificial) bound on the objective <=1.0 (to get a feasible solution) it solved the problem immediately!

------------------------------------------------------------------
               SCS v3.2.1 - Splitting Conic Solver
        (c) Brendan O'Donoghue, Stanford University, 2012
------------------------------------------------------------------
problem:  variables n: 570684, constraints m: 1030532
cones:    z: primal zero / dual free vars: 459848
          l: linear vars: 1
          s: psd vars: 570683, ssize: 10
settings: eps_abs: 1.0e-09, eps_rel: 1.0e-09, eps_infeas: 1.0e-07
          alpha: 1.95, scale: 1.00e-01, adaptive_scale: 1
          max_iters: 100000, normalize: 1, rho_x: 1.00e-06
          acceleration_lookback: 50, acceleration_interval: 10
lin-sys:  sparse-direct-mkl-pardiso
          nnz(A): 23505313, nnz(P): 0
------------------------------------------------------------------
 iter | pri res | dua res |   gap   |   obj   |  scale  | time (s)
------------------------------------------------------------------
     0| 4.00e+01  1.00e+00  2.03e+01 -1.19e+01  1.00e-01  4.78e+02 
   250| 4.53e-03  1.26e-05  4.59e-03 -1.00e+00  2.13e+00  1.75e+03 
   500| 2.57e-07  3.12e-07  2.17e-07 -1.00e+00  2.13e+00  2.13e+03 
   600| 3.08e-09  8.33e-10  6.20e-10 -1.00e+00  6.64e-01  2.73e+03 
------------------------------------------------------------------
status:  solved
timings: total: 2.73e+03s = setup: 4.75e+02s + solve: 2.25e+03s
         lin-sys: 8.48e+02s, cones: 6.93e+01s, accel: 1.79e+00s
------------------------------------------------------------------
objective = -1.000000
------------------------------------------------------------------

Woah! So adding a redundant constraint actually improved performance? That's wild. This could open up a whole new line of research.

@kalmarek
Copy link
Collaborator Author

I think the constraint is not really redundant:

previously I had

max t
s.t. [...]

now instead of asking for max t I asked for a single feasible solution with t=1.0:

max t
s.t. t ≤ 1.0
     [...] # same constraints as above

Is this what what you meant by "redundant"?

I've already noted much earlier (2017? :) that approaching this way the maximum (from below) is sometimes much faster than running scs directly with t "unconstrained"...

@bodono
Copy link
Contributor

bodono commented Apr 13, 2023

Oh I see, I thought the constraint was redundant because the objective is -1 and the constraint was t <= 1, but I see you're doing a maximization problem.

Does approaching from below involve solving a series of these problems and using warm-starting etc? Is it faster than just solving the problem directly?

@kalmarek
Copy link
Collaborator Author

Yes, I think JuMP reformulates the problem and flips some signs so that maximization becomes minimization etc. In particular SCS reports negative objective when JuMP.objective_value(model)` is positive. Sorry for the confusion. The problem SCS sees looks is

min t
s.t t ≥ -1.0
    [ ... ] # # same constraints as above

Yes, I know that my problem will have a (trivial) solution when t=0; I'm interested in solutions with t << 0 (and sufficiently large, so that I can certify the solution).

So sometimes I start with t ≥ -0.1 and then get warmstart, decrement the lower bound in by 0.05 and solve again. when SCS can't get the optimal solution quickly this method works reasonably well.

@kalmarek kalmarek mentioned this issue Jun 28, 2023
Open
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Assignees
No one assigned
Labels
None yet
Milestone
No milestone
Development

No branches or pull requests
3 participants
@bodono @kalmarek @odow