QuadraticModels.jl icon indicating copy to clipboard operation
QuadraticModels.jl copied to clipboard

Track allocs

Open tmigot opened this issue 3 years ago • 9 comments

There seems to be two functions of the API allocating

  • [ ] obj
  • [x] hess_coord

I used the following script to track allocations in QuadraticModels

using Pkg
Pkg.activate(".")

# stdlib
using LinearAlgebra, Printf, SparseArrays, Test

# our packages
using ADNLPModels,
  LinearOperators,
  NLPModels,
  NLPModelsModifiers,
  NLPModelsTest, # main version of NLPModelsTest
  QPSReader,
  QuadraticModels,
  SparseMatricesCOO

function only_nonzeros(table)
  for k in keys(table)
    if table[k] == 0
      pop!(table, k)
    end
  end
  return table
end

# Definition of quadratic problems
qp_problems_Matrix = ["bndqp", "eqconqp"]
qp_problems_COO = ["uncqp", "ineqconqp"]
for qp in [qp_problems_Matrix; qp_problems_COO]
  include(joinpath("problems", "$qp.jl"))
end

for problem in qp_problems_Matrix
  @info "Checking allocs of dense problem $(problem)_QPSData"
  nlp_qps = eval(Symbol(problem * "_QPSData"))()
  print_nlp_allocations(nlp_qps, only_nonzeros(test_allocs_nlpmodels(nlp_qps)))
end

for problem in qp_problems_Matrix
  @info "Checking allocs of dense problem $(problem)_QP_dense"
  nlp_qm_dense = eval(Symbol(problem * "_QP_dense"))()
  print_nlp_allocations(nlp_qm_dense, only_nonzeros(test_allocs_nlpmodels(nlp_qm_dense)))
end

for problem in qp_problems_Matrix
  @info "Checking allocs of dense problem $(problem)_QP_sparse"
  nlp_qm_sparse = eval(Symbol(problem * "_QP_sparse"))()
  print_nlp_allocations(nlp_qm_sparse, only_nonzeros(test_allocs_nlpmodels(nlp_qm_sparse)))
end

for problem in qp_problems_Matrix
  @info "Checking allocs of dense problem $(problem)_QP_symmetric"
  nlp_qm_symmetric = eval(Symbol(problem * "_QP_symmetric"))()
  print_nlp_allocations(nlp_qm_symmetric, only_nonzeros(test_allocs_nlpmodels(nlp_qm_symmetric)))
end

for problem in qp_problems_COO
  @info "Checking allocs of COO problem $(problem)_QPSData"
  nlp_qps = eval(Symbol(problem * "_QPSData"))()
  print_nlp_allocations(nlp_qps, only_nonzeros(test_allocs_nlpmodels(nlp_qps)))
end

for problem in qp_problems_COO
  @info "Checking allocs of COO problem $(problem)_QP"
  nlp_qm_dense = eval(Symbol(problem * "_QP"))()
  print_nlp_allocations(nlp_qm_dense, only_nonzeros(test_allocs_nlpmodels(nlp_qm_dense)))
end

for problem in NLPModelsTest.nlp_problems
  @info "Testing allocs of quadratic approximation of problem $problem"
  nlp = eval(Symbol(problem))()
  x = nlp.meta.x0
  nlp_qm = QuadraticModel(nlp, x)
  print_nlp_allocations(nlp_qm, only_nonzeros(test_allocs_nlpmodels(nlp_qm)))
end

and the results

[ Info: Checking allocs of dense problem bndqp_QPSData
  Problem name: Generic
                        obj: ████████████████████ 80.0
                hess_coord!: ████████████████████ 80.0

[ Info: Checking allocs of dense problem eqconqp_QPSData
  Problem name: Generic
                        obj: ████████████████████ 496.0
                hess_coord!: ████████████████████ 496.0
            hess_lag_coord!: ████████████████████ 496.0

[ Info: Checking allocs of dense problem bndqp_QP_dense
  Problem name: bndqp_QP
                        obj: ████████████████████ 80.0

[ Info: Checking allocs of dense problem eqconqp_QP_dense
  Problem name: eqconqp_QP
                        obj: ████████████████████ 496.0

[ Info: Checking allocs of dense problem bndqp_QP_sparse
  Problem name: bndqp_QP
                        obj: ████████████████████ 80.0

[ Info: Checking allocs of dense problem eqconqp_QP_sparse
  Problem name: eqconqp_QP
                        obj: ████████████████████ 496.0

[ Info: Checking allocs of dense problem bndqp_QP_symmetric
  Problem name: bndqp_QP
                        obj: ████████████████████ 80.0

[ Info: Checking allocs of dense problem eqconqp_QP_symmetric
  Problem name: eqconqp_QP
                        obj: ████████████████████ 496.0

[ Info: Checking allocs of COO problem uncqp_QPSData
  Problem name: Generic
                        obj: ████████████████████ 80.0
                hess_coord!: ████████████████████ 80.0

[ Info: Checking allocs of COO problem ineqconqp_QPSData
  Problem name: Generic
                        obj: ████████████████████ 80.0
                hess_coord!: ████████████████████ 80.0
            hess_lag_coord!: ████████████████████ 80.0

[ Info: Checking allocs of COO problem uncqp_QP
  Problem name: uncqp_QP
                        obj: ████████████████████ 80.0
                hess_coord!: ████████████████████ 80.0

[ Info: Checking allocs of COO problem ineqconqp_QP
  Problem name: ineqconqp_QP
                        obj: ████████████████████ 80.0
                hess_coord!: ████████████████████ 80.0
            hess_lag_coord!: ████████████████████ 80.0

[ Info: Testing allocs of quadratic approximation of problem BROWNDEN
  Problem name: Generic
                        obj: ██████████████⋅⋅⋅⋅⋅⋅ 96.0
                hess_coord!: ████████████████████ 144.0

[ Info: Testing allocs of quadratic approximation of problem HS5
  Problem name: Generic
                        obj: ████████████████████ 80.0
                hess_coord!: ████████████████████ 80.0

[ Info: Testing allocs of quadratic approximation of problem HS6
  Problem name: Generic
                        obj: ████████████████████ 80.0
                hess_coord!: ████████████████⋅⋅⋅⋅ 64.0
            hess_lag_coord!: ████████████████⋅⋅⋅⋅ 64.0

[ Info: Testing allocs of quadratic approximation of problem HS10
  Problem name: Generic
                        obj: ████████████████████ 80.0
                hess_coord!: ████████████████████ 80.0
            hess_lag_coord!: ████████████████████ 80.0

[ Info: Testing allocs of quadratic approximation of problem HS11
  Problem name: Generic
                        obj: ████████████████████ 80.0
                hess_coord!: ████████████████████ 80.0
            hess_lag_coord!: ████████████████████ 80.0

[ Info: Testing allocs of quadratic approximation of problem HS13
  Problem name: Generic
                        obj: ████████████████████ 80.0
                hess_coord!: ████████████████████ 80.0
            hess_lag_coord!: ████████████████████ 80.0

[ Info: Testing allocs of quadratic approximation of problem HS14
                        obj: ████████████████████ 80.0
                hess_coord!: ████████████████████ 80.0
            hess_lag_coord!: ████████████████████ 80.0

[ Info: Testing allocs of quadratic approximation of problem LINCON
  Problem name: Generic
                        obj: ████████████████████ 176.0
                hess_coord!: ████████████████████ 176.0
            hess_lag_coord!: ████████████████████ 176.0

[ Info: Testing allocs of quadratic approximation of problem LINSV
  Problem name: Generic
                        obj: ████████████████████ 80.0
                hess_coord!: ████████████████⋅⋅⋅⋅ 64.0
            hess_lag_coord!: ████████████████⋅⋅⋅⋅ 64.0

[ Info: Testing allocs of quadratic approximation of problem MGH01Feas
  Problem name: Generic
                        obj: ████████████████████ 80.0
                hess_coord!: ████████████████⋅⋅⋅⋅ 64.0
            hess_lag_coord!: ████████████████⋅⋅⋅⋅ 64.0

tmigot avatar Sep 09 '22 13:09 tmigot

After #107 only obj allocate.

tmigot avatar Sep 09 '22 14:09 tmigot

Fixing obj would require writing a function that computes x' * H * x without allocating (for all H types).

geoffroyleconte avatar Sep 09 '22 14:09 geoffroyleconte

Isn't there a 3-arguments dot?

help?> dot
search: dot @__dot__ stdout DEPOT_PATH adjoint Adjoint adjoint! fieldcount codepoint fieldoffset bound_constrained isdisjoint lowrankdowndate lowrankdowndate! denominator VoidListofStates add_to_list! PKGMODE_PROJECT modifyproperty!

  dot(x, y)
  x ⋅ y

  Compute the dot product between two vectors. For complex vectors, the first vector is conjugated.

  dot also works on arbitrary iterable objects, including arrays of any dimension, as long as dot is defined on the elements.

  dot is semantically equivalent to sum(dot(vx,vy) for (vx,vy) in zip(x, y)), with the added restriction that the arguments must have equal lengths.

  x ⋅ y (where ⋅ can be typed by tab-completing \cdot in the REPL) is a synonym for dot(x, y).

  Examples
  ≡≡≡≡≡≡≡≡≡≡

  julia> dot([1; 1], [2; 3])
  5
  
  julia> dot([im; im], [1; 1])
  0 - 2im
  
  julia> dot(1:5, 2:6)
  70
  
  julia> x = fill(2., (5,5));
  
  julia> y = fill(3., (5,5));
  
  julia> dot(x, y)
  150.0

  ──────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────  

  dot(x, A, y)

  Compute the generalized dot product dot(x, A*y) between two vectors x and y, without storing the intermediate result of A*y. As for the two-argument dot(_,_), this acts recursively. Moreover, for complex vectors, the first vector is        
  conjugated.

  │ Julia 1.4
  │
  │  Three-argument dot requires at least Julia 1.4.

  Examples
  ≡≡≡≡≡≡≡≡≡≡

  julia> dot([1; 1], [1 2; 3 4], [2; 3])
  26
  
  julia> dot(1:5, reshape(1:25, 5, 5), 2:6)
  4850
  
  julia> ⋅(1:5, reshape(1:25, 5, 5), 2:6) == dot(1:5, reshape(1:25, 5, 5), 2:6)
  true

tmigot avatar Sep 09 '22 14:09 tmigot

Oh I didn't know about this function. It would have to be added to SparseMatricesCOO. Also I would keep the existing behaviour in case H is a LinearOperator.

geoffroyleconte avatar Sep 09 '22 15:09 geoffroyleconte

Ok, so if I understand. We should implement this 3-args dot in SparseMatricesCOO and in LinearOperators?

tmigot avatar Sep 09 '22 15:09 tmigot

Yes for SparseMatricesCOO. For LinearOperators this is complicated, I would keep the existing code even though it allocates.

geoffroyleconte avatar Sep 09 '22 17:09 geoffroyleconte

How about preallocating storage in the quadratic model to store H * x? It could be reused for the gradient.

dpo avatar Sep 09 '22 23:09 dpo

Has this been resolved since ? I'd be interested in using this repo but I need there to be no allocations. I can open a PR to implement @dpo's last suggestion, I think it is the easiest.

MaxenceGollier avatar Feb 03 '25 21:02 MaxenceGollier

Feel free to go ahead @MaxenceGollier !

dpo avatar Feb 04 '25 00:02 dpo