Error in the function jprod!

[amontoison]

using VecchiaMLE
using SparseArrays
n = 20
k = 10
Number_of_Samples = 100
params = [5.0, 0.2, 2.25, 0.25]
MatCov = VecchiaMLE.generate_MatCov(n, params)
samples = VecchiaMLE.generate_Samples(MatCov, n, Number_of_Samples)
iVecchiaMLE = VecchiaMLE.VecchiaMLEInput(n, k, samples, Number_of_Samples, 5, 1)
pres_chol = Matrix{eltype(iVecchiaMLE.samples)}(undef, iVecchiaMLE.n^2, iVecchiaMLE.n^2)
fill!(pres_chol, zero(eltype(iVecchiaMLE.samples)))
diagnostics = VecchiaMLE.Diagnostics(0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0, 0)
xyGrid = generate_xyGrid(iVecchiaMLE.n)
model = VecchiaMLE.VecchiaModelCPU(iVecchiaMLE.samples, iVecchiaMLE.k, xyGrid)

using NLPModelsTest
print_nlp_allocations(model)

ERROR: MethodError: no method matching (::Colon)(::Int64, ::Vector{Int64})
The function `Colon()` exists, but no method is defined for this combination of argument types.

Closest candidates are:
  (::Colon)(::T, ::Any, ::T) where T<:Real
   @ Base range.jl:50
  (::Colon)(::A, ::Any, ::C) where {A<:Real, C<:Real}
   @ Base range.jl:10
  (::Colon)(::T, ::Any, ::T) where T
   @ Base range.jl:49
  ...

Stacktrace:
 [1] jprod!(nlp::VecchiaMLE.VecchiaModel{Float64, Vector{…}}, x::Vector{Float64}, v::Vector{Float64}, Jv::Vector{Float64})
   @ VecchiaMLE ~/Bureau/git/VecchiaMLE.jl/src/models/VecchiaMLE_NLPModel.jl:242
 [2] test_allocs_nlpmodels(nlp::VecchiaMLE.VecchiaModel{Float64, Vector{Float64}}; linear_api::Bool, exclude::Vector{Any})
   @ NLPModelsTest ~/.julia/packages/NLPModelsTest/fKFni/src/allocs_model.jl:95
 [3] test_allocs_nlpmodels
   @ ~/.julia/packages/NLPModelsTest/fKFni/src/allocs_model.jl:11 [inlined]
 [4] #print_nlp_allocations#209
   @ ~/.julia/packages/NLPModelsTest/fKFni/src/allocs_model.jl:356 [inlined]
 [5] print_nlp_allocations(io::Base.TTY, nlp::VecchiaMLE.VecchiaModel{Float64, Vector{Float64}})
   @ NLPModelsTest ~/.julia/packages/NLPModelsTest/fKFni/src/allocs_model.jl:349
 [6] print_nlp_allocations(nlp::VecchiaMLE.VecchiaModel{Float64, Vector{Float64}}; kwargs::@Kwargs{})
   @ NLPModelsTest ~/.julia/packages/NLPModelsTest/fKFni/src/allocs_model.jl:346
 [7] print_nlp_allocations(nlp::VecchiaMLE.VecchiaModel{Float64, Vector{Float64}})
   @ NLPModelsTest ~/.julia/packages/NLPModelsTest/fKFni/src/allocs_model.jl:345
 [8] top-level scope
   @ REPL[29]:1
Some type information was truncated. Use `show(err)` to see complete types.

The function print_nlp_allocations is useful to detect some memory leaks

[amontoison]

We should use the following function is the tests:
https://jso.dev/NLPModelsTest.jl/dev/reference/#NLPModelsTest.test_allocs_nlpmodels-Tuple{NLPModels.AbstractNLPModel}

[amontoison]

julia> print_nlp_allocations(model)

Dict{Symbol, Float64} with 8 entries:
  :obj             => 70544.0
  :cons!           => 192.0
  :hess_coord!     => 3424.0
  :grad!           => 57616.0
  :hess_lag_coord! => 160.0
  :jac_coord!      => 16.0
  :hess_structure! => 1.60547e6
  :jac_structure!  => 19520.0

@CalebDerrickson 💀

[amontoison]

The issue is related to type inference.

julia> @code_warntype NLPModels.grad!(model, x, gx)
MethodInstance for grad!(::VecchiaMLE.VecchiaModel{Float64, Vector{Float64}}, ::Vector{Float64}, ::Vector{Float64})
  from grad!(nlp::VecchiaMLE.VecchiaModel, x::AbstractVector, gx::AbstractVector) @ VecchiaMLE ~/Bureau/git/VecchiaMLE.jl/src/models/VecchiaMLE_NLPModel.jl:126
Arguments
  #self#::Core.Const(NLPModels.grad!)
  nlp::VecchiaMLE.VecchiaModel{Float64, Vector{Float64}}
  x::Vector{Float64}
  gx::Vector{Float64}
Locals
  @_5::Union{Nothing, Tuple{Int64, Int64}}
  gx_z::SubArray{Float64, 1, Vector{Float64}, Tuple{UnitRange{Int64}}, true}
  pos::Int64
  j::Int64
  gxj::SubArray{Float64, 1, Vector{Float64}, Tuple{UnitRange{Int64}}, true}
  xj::SubArray{Float64, 1, Vector{Float64}, Tuple{UnitRange{Int64}}, true}
  Bj::Any  # <-------------------------------- NOT GOOD
Body::Vector{Float64}

[amontoison]

Ok, I firxed many of the type instabilities:

julia> print_nlp_allocations(model)

  Problem name: Vecchia_manual
             obj: ⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅ 0.0   
           cons!: ⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅ 0.0   
     hess_coord!: ⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅ 0.0   
           grad!: ⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅ 0.0   
 hess_lag_coord!: ⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅ 0.0   
      jac_coord!: ⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅ 0.0   
 hess_structure!: ████████████████████ 1.605472e6
  jac_structure!: █⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅ 19424.0

Dict{Symbol, Float64} with 8 entries:
  :obj             => 0.0
  :cons!           => 0.0
  :hess_coord!     => 0.0
  :grad!           => 0.0
  :hess_lag_coord! => 0.0
  :jac_coord!      => 0.0
  :hess_structure! => 1.60547e6
  :jac_structure!  => 19424.0

[CalebDerrickson]

Bj::Any # <-------------------------------- NOT GOOD

Are you saying that the type Any is not good, or that Bj shouldn't be in the grad function?

Also, I fixed the error in jprod. See #16.

[amontoison]

Type Any is not good, the compiler was not able to determine the type of Bj and needed to allocate a lot of memory to handle that

[amontoison]

julia> print_nlp_allocations(model)
  Problem name: Vecchia_manual
             obj: ⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅ 0.0   
           cons!: ⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅ 0.0   
     hess_coord!: ⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅ 0.0   
           grad!: ⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅ 0.0   
 hess_lag_coord!: ⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅ 0.0   
      jac_coord!: ⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅ 0.0   
 hess_structure!: ⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅ 0.0   
  jac_structure!: ████████████████████ 19424.0

Dict{Symbol, Float64} with 8 entries:
  :obj             => 0.0
  :cons!           => 0.0
  :hess_coord!     => 0.0
  :grad!           => 0.0
  :hess_lag_coord! => 0.0
  :jac_coord!      => 0.0
  :hess_structure! => 0.0
  :jac_structure!  => 19424.0

[amontoison]

Good job Caleb with #17 👍

[CalebDerrickson]

@amontoison I just ran @code_warntype NLPModels.grad! and saw that Bj was Matrix{Float64}, and the issue with jprod was resolved. Should this issue be closed?

[amontoison]

Can you add a test with NLPModelsTest.test_allocs_nlpmodels before such we keep the functions allocation-free with the future modifications?