diff --git a/.JuliaFormatter.toml b/.JuliaFormatter.toml index a700a07..f196f10 100644 --- a/.JuliaFormatter.toml +++ b/.JuliaFormatter.toml @@ -5,4 +5,5 @@ always_for_in = true always_use_return = true margin = 80 remove_extra_newlines = true +separate_kwargs_with_semicolon = true short_to_long_function_def = true diff --git a/.github/workflows/ci.yml b/.github/workflows/ci.yml index f9c8abb..8b93ca3 100644 --- a/.github/workflows/ci.yml +++ b/.github/workflows/ci.yml @@ -49,6 +49,8 @@ jobs: with: depwarn: error - uses: julia-actions/julia-processcoverage@v1 - - uses: codecov/codecov-action@v3 + - uses: codecov/codecov-action@v5 with: - file: lcov.info + files: lcov.info + token: ${{ secrets.CODECOV_TOKEN }} + diff --git a/.github/workflows/codecov.yml b/.github/workflows/codecov.yml deleted file mode 100644 index 3217852..0000000 --- a/.github/workflows/codecov.yml +++ /dev/null @@ -1,2 +0,0 @@ -codecov: - token: 3c305759-ecc3-419d-aee0-20ab9a53d1e9 \ No newline at end of file diff --git a/.github/workflows/documentation.yml b/.github/workflows/documentation.yml index 76aff3f..31c3146 100644 --- a/.github/workflows/documentation.yml +++ b/.github/workflows/documentation.yml @@ -9,7 +9,7 @@ jobs: build: runs-on: ubuntu-latest steps: - - uses: actions/checkout@v2 + - uses: actions/checkout@v4 - uses: julia-actions/setup-julia@latest with: version: '1' diff --git a/.github/workflows/format_check.yml b/.github/workflows/format_check.yml index 88045c8..4f91f39 100644 --- a/.github/workflows/format_check.yml +++ b/.github/workflows/format_check.yml @@ -13,13 +13,13 @@ jobs: - uses: julia-actions/setup-julia@latest with: version: '1' - - uses: actions/checkout@v1 + - uses: actions/checkout@v4 - name: Format check shell: julia --color=yes {0} run: | using Pkg # If you update the version, also update the style guide docs. - Pkg.add(PackageSpec(name="JuliaFormatter", version="1.0.13")) + Pkg.add(PackageSpec(name="JuliaFormatter", version="2")) using JuliaFormatter format("."; verbose = true) out = String(read(Cmd(`git diff`))) diff --git a/JuliaFormatter.toml b/JuliaFormatter.toml deleted file mode 100644 index 74a0a71..0000000 --- a/JuliaFormatter.toml +++ /dev/null @@ -1,8 +0,0 @@ -# Configuration file for JuliaFormatter.jl -# For more information, see: https://domluna.github.io/JuliaFormatter.jl/stable/config/ - -always_for_in = true -always_use_return = true -margin = 80 -remove_extra_newlines = true -short_to_long_function_def = true \ No newline at end of file diff --git a/benchmark/benchmarks.jl b/benchmark/benchmarks.jl index 1263277..afccf0e 100644 --- a/benchmark/benchmarks.jl +++ b/benchmark/benchmarks.jl @@ -296,7 +296,7 @@ function model_sparse(F, C, P, T, D, U, V, W) # Variable creation @variable( m, - x[factory = F, customer = C, product = P, period = T], + x[factory=F, customer=C, product=P, period=T], container = IndexedVarArray ) @@ -338,10 +338,10 @@ function model_sparse_aa(F, C, P, T, D, U, V, W) @variable( m, x[ - factory = F, - customer = C, - product = P, - period = T; + factory=F, + customer=C, + product=P, + period=T; W[factory, product] == 1 && (factory, product, period) in keys(D), ] >= 0, ) @@ -455,7 +455,7 @@ REPS = 5 # ╔═╡ 04570ea7-885c-4d0e-be88-eb2a5f77da90 begin - res = DataFrame(Method = Symbol[], NC = Int[], Time = Float64[]) + res = DataFrame(; Method = Symbol[], NC = Int[], Time = Float64[]) @progress for nc in 5:10:100 for method in [ model_standard, @@ -476,7 +476,7 @@ end # ╔═╡ cc097148-23b1-4584-a150-c7f22376b65c begin - sparsity = DataFrame(Method = Symbol[], DP = Float64[], Time = Float64[]) + sparsity = DataFrame(; Method = Symbol[], DP = Float64[], Time = Float64[]) @progress for dp in 0.05:0.05:1.0 for method in [ model_standard, @@ -511,10 +511,10 @@ end # ╔═╡ b0aa0499-e920-4014-b0b2-ce8ea3da7c95 function plot(df, x = :NC, y = :Time) - CairoMakie.activate!(type = "svg") + CairoMakie.activate!(; type = "svg") return draw( data(df) * - mapping(x, y => "Time (s)", color = :Method, marker = :Method) * + mapping(x, y => "Time (s)"; color = :Method, marker = :Method) * (visual(Lines) + visual(Scatter)), ) end diff --git a/benchmark/transport.jl b/benchmark/transport.jl index 4f19faf..47b326f 100644 --- a/benchmark/transport.jl +++ b/benchmark/transport.jl @@ -360,7 +360,7 @@ function create_vars_indexedtable(m, pp) C, P, T, - V, + V; names = [:factory, :customer, :product, :period, :var], pkey = [:factory, :customer, :product, :period], ) @@ -371,7 +371,7 @@ function create_constraints_indexedtable(m, pp) flow = m[:flow] # Production capacity - pc_table = groupby(collect, flow, (:factory, :product), select = :var) + pc_table = groupby(collect, flow, (:factory, :product); select = :var) for r in rows(pc_table) if (r.factory, r.product) in keys(pp.prodcap) @constraint( @@ -384,7 +384,7 @@ function create_constraints_indexedtable(m, pp) # Customer demand cpp_table = - groupby(collect, flow, (:customer, :product, :period), select = :var) + groupby(collect, flow, (:customer, :product, :period); select = :var) for r in rows(cpp_table) if (r.customer, r.product, r.period) in keys(pp.demand) @constraint( @@ -397,7 +397,7 @@ function create_constraints_indexedtable(m, pp) # Transport capacity fc_table = - groupby(collect, flow, (:factory, :customer, :period), select = :var) + groupby(collect, flow, (:factory, :customer, :period); select = :var) for r in rows(fc_table) if (r.factory, r.customer) in keys(pp.flowcap) @constraint( diff --git a/docs/make.jl b/docs/make.jl index 8314890..ec7879a 100644 --- a/docs/make.jl +++ b/docs/make.jl @@ -9,7 +9,7 @@ pages = [ "API reference" => "reference/api.md", ] -Documenter.makedocs( +Documenter.makedocs(; sitename = "SparseVariables", format = Documenter.HTML(; prettyurls = get(ENV, "CI", "false") == "true", diff --git a/docs/notebook_juliacon2022.jl b/docs/notebook_juliacon2022.jl index 8ea55b6..3cf6ac0 100644 --- a/docs/notebook_juliacon2022.jl +++ b/docs/notebook_juliacon2022.jl @@ -408,7 +408,7 @@ REPS = 5 # ╔═╡ 04570ea7-885c-4d0e-be88-eb2a5f77da90 begin - res = DataFrame(Method = Symbol[], NC = Int[], Time = Float64[]) + res = DataFrame(; Method = Symbol[], NC = Int[], Time = Float64[]) @progress for nc in 5:10:100 for method in [ model_standard, @@ -428,7 +428,7 @@ end # ╔═╡ cc097148-23b1-4584-a150-c7f22376b65c begin - sparsity = DataFrame(Method = Symbol[], DP = Float64[], Time = Float64[]) + sparsity = DataFrame(; Method = Symbol[], DP = Float64[], Time = Float64[]) @progress for dp in 0.05:0.05:1.0 for method in [ model_standard, @@ -462,10 +462,10 @@ end # ╔═╡ b0aa0499-e920-4014-b0b2-ce8ea3da7c95 function plot(df, x = :NC, y = :Time) - CairoMakie.activate!(type = "svg") + CairoMakie.activate!(; type = "svg") return draw( data(df) * - mapping(x, y => "Time (s)", color = :Method, marker = :Method) * + mapping(x, y => "Time (s)"; color = :Method, marker = :Method) * (visual(Lines) + visual(Scatter)), ) end diff --git a/src/SparseVariables.jl b/src/SparseVariables.jl index 96b2e63..5414888 100644 --- a/src/SparseVariables.jl +++ b/src/SparseVariables.jl @@ -6,15 +6,20 @@ using LinearAlgebra using PrecompileTools include("sparsearray.jl") +include("slice.jl") +include("broadcast.jl") include("dictionaries.jl") include("indexedarray.jl") include("tables.jl") export SparseArray export IndexedVarArray +export SparseArraySlice +export slice export insertvar! export unsafe_insertvar! export SafeInsert, UnsafeInsert +export set_cache_cutoff! @setup_workload begin # Putting some things in `setup` can reduce the size of the @@ -29,11 +34,7 @@ export SafeInsert, UnsafeInsert # all calls in this block will be precompiled, regardless of whether # they belong to your package or not (on Julia 1.8 and higher) - @variable( - m, - x[r = rs, i = is, st = sts, sy = sys]; - container = IndexedVarArray - ) + @variable(m, x[r=rs, i=is, st=sts, sy=sys]; container = IndexedVarArray) for r in rs, i in is, st in sts, sy in sys insertvar!(x, r, i, st, sy) unsafe_insertvar!(x, r, i, st, sy) @@ -41,7 +42,7 @@ export SafeInsert, UnsafeInsert x[:, 1, :, :] x[10, :, :, :] x[1, :, :, :a] - @variable(m, y[i = rs, j = rs, k = rs]; container = IndexedVarArray) + @variable(m, y[i=rs, j=rs, k=rs]; container = IndexedVarArray) for i in rs, j in rs, k in rs insertvar!(y, i, j, k) end diff --git a/src/broadcast.jl b/src/broadcast.jl new file mode 100644 index 0000000..d4e4e12 --- /dev/null +++ b/src/broadcast.jl @@ -0,0 +1,129 @@ +# ------------------------------------------------------------------------------ +# Broadcasting over AbstractSparseArray +# Follows the pattern of JuMP.Containers.SparseAxisArray. +# The result of any broadcast is always a plain SparseArray. +# ------------------------------------------------------------------------------ + +""" + SparseBroadcastStyle{N,K} <: Broadcast.BroadcastStyle + +Broadcasting style for all `AbstractSparseArray` subtypes. `K` is the key tuple type. +All broadcast results are materialised as `SparseArray`. +""" +struct SparseBroadcastStyle{K} <: Broadcast.BroadcastStyle end + +function Base.BroadcastStyle(::Type{SA}) where {SA<:AbstractSparseArray} + return SparseBroadcastStyle{_keytype(SA)}() +end + +# Disallow mixing with other array types. +function Base.BroadcastStyle(::SparseBroadcastStyle, ::Base.BroadcastStyle) + return throw( + ArgumentError( + "Cannot broadcast a SparseArray with incompatible key types", + ), + ) +end + +# Scalar (0-d) broadcasting is allowed. +function Base.BroadcastStyle( + style::SparseBroadcastStyle, + ::Base.Broadcast.DefaultArrayStyle{0}, +) + return style +end + +# Fix ambiguity with Unknown. +function Base.BroadcastStyle(::SparseBroadcastStyle, ::Base.Broadcast.Unknown) + return throw( + ArgumentError( + "Cannot broadcast a SparseArray with an unknown broadcast style", + ), + ) +end + +# Bypass the default instantiate which calls axes(). +function Base.Broadcast.instantiate( + bc::Base.Broadcast.Broadcasted{<:SparseBroadcastStyle}, +) + return bc +end + +# Internal helpers +_sparse_getindex(x::AbstractSparseArray, key) = x[key] +_sparse_getindex(x::Any, ::Any) = x +_sparse_getindex(x::Ref, ::Any) = x[] + +function _sparse_getindex( + bc::Base.Broadcast.Broadcasted{<:SparseBroadcastStyle}, + key, +) + return bc.f(_sparse_get_args(bc.args, key)...) +end + +function _sparse_get_args(args::Tuple, key) + return ( + _sparse_getindex(first(args), key), + _sparse_get_args(Base.tail(args), key)..., + ) +end +_sparse_get_args(::Tuple{}, ::Any) = () + +function _sparse_check_same_keys(ref_keys, x::AbstractSparseArray, args...) + if length(ref_keys) != length(x) || any(k -> !haskey(x, k), ref_keys) + throw( + ArgumentError( + "Cannot broadcast SparseArrays with different indices", + ), + ) + end + return _sparse_check_same_keys(ref_keys, args...) +end + +function _sparse_check_same_keys(ref_keys, ::Any, args...) + return _sparse_check_same_keys(ref_keys, args...) +end +_sparse_check_same_keys(::Any) = nothing + +function _sparse_indices( + bc::Base.Broadcast.Broadcasted{<:SparseBroadcastStyle}, + rest..., +) + return _sparse_indices(bc.args..., rest...) +end + +function _sparse_indices(x::AbstractSparseArray, rest...) + ks = collect(keys(x)) + _sparse_check_same_keys(ks, rest...) + return ks +end + +_sparse_indices(::Any, rest...) = _sparse_indices(rest...) + +# Materialise + +function Base.copy( + bc::Base.Broadcast.Broadcasted{SparseBroadcastStyle{K}}, +) where {K} + indices = _sparse_indices(bc) + T = Base.Broadcast.combine_eltypes(bc.f, bc.args) + isempty(indices) && return SparseArray(Dictionary{K,T}()) + vals = [_sparse_getindex(bc, k) for k in indices] + return SparseArray(Dictionary(indices, vals)) +end + +function Base.Broadcast.broadcast_preserving_zero_d( + f, + A::AbstractSparseArray, + As..., +) + return broadcast(f, A, As...) +end +function Base.Broadcast.broadcast_preserving_zero_d( + f, + x, + A::AbstractSparseArray, + As..., +) + return broadcast(f, x, A, As...) +end diff --git a/src/dictionaries.jl b/src/dictionaries.jl index d065710..1556dad 100644 --- a/src/dictionaries.jl +++ b/src/dictionaries.jl @@ -1,7 +1,7 @@ """ make_filter_fun(c, pos) -Return function to use for filtering depending on the type and value of `c` +Return function to use for filtering depending on the type and value of `c` to apply at position `pos` """ make_filter_fun(c, pos) = x -> x[pos] == c @@ -12,7 +12,8 @@ make_filter_fun(c) = x -> x == c make_filter_fun(c::Base.Fix2) = x -> c(x) make_filter_fun(c::Function) = x -> c(x) make_filter_fun(c::Colon) = x -> true -make_filter_fun(c::UnitRange) = x -> (x ≥ c.start && x ≤ c.stop) +make_filter_fun(c::AbstractRange, pos) = x -> x[pos] in c +make_filter_fun(c::AbstractRange) = x -> x in c """ recursive_filter(fs, data) @@ -38,7 +39,7 @@ function indices_fun(some_tuple) end """ - _select_rowwise(a, pattern) + _select_rowwise(a, pattern) Filter iterable data a by tuple `pattern` by row (slow) """ @@ -47,7 +48,7 @@ function _select_rowwise(a, pattern) end """ - _select_colwise(a, pattern) + _select_colwise(a, pattern) Filter iterable data a by tuple `pattern` by column (recursively) """ @@ -59,7 +60,7 @@ end _select_gen(a, pattern) Filter iterable data `a` by tuple `pattern` by row, using generated function for speed. -See more straight-forward implementations `_select_rowwise` and `_select_colwise` for reference. +See more straight-forward implementations `_select_rowwise` and `_select_colwise` for reference. """ function _select_gen(a, pattern) return filter(x -> _select_generated(pattern, x), a) @@ -68,7 +69,7 @@ end """ _select_gen_perm(a, pattern, perm) Filter iterable `data` byt tuple `pattern` by row using generated function that permutes the sequence of -evaluation by the permutation tuple `perm` for improved control as this can give performance advantages, +evaluation by the permutation tuple `perm` for improved control as this can give performance advantages, depending on the uniqueness of the search pattern and cost of function evaluation. ## Example @@ -226,33 +227,23 @@ Works on types because it is used in generated function """ isfixed(t) = true isfixed(::Type{T} where {T<:Function}) = false -isfixed(::Type{T} where {T<:UnitRange}) = false -iscolon(t) = false -iscolon(::Type{T} where {T<:Colon}) = true +isfixed(::Type{T} where {T<:AbstractRange}) = false @generated function _getindex( sa::AbstractSparseArray{T,N}, tpl::Tuple, ) where {T,N} lookup = true - slice = true for t in fieldtypes(tpl) if !isfixed(t) lookup = false - if !iscolon(t) - slice = false - end end end if lookup return :(get(_data(sa), tpl, zero(T))) - elseif !slice - return :(retval = select(_data(sa), tpl); - length(retval) > 0 ? retval : zero(T)) - else # Return selection or zero if empty to avoid reduction of empty iterate - return :(retval = _select_var(sa, tpl); - length(retval) > 0 ? retval : zero(T)) + else + return :(_make_slice(sa, tpl)) end end diff --git a/src/indexedarray.jl b/src/indexedarray.jl index 9f52440..edc708d 100644 --- a/src/indexedarray.jl +++ b/src/indexedarray.jl @@ -14,6 +14,7 @@ struct SafeInsert end struct UnsafeInsert end _data(sa::IndexedVarArray) = sa.data +_keytype(::Type{<:IndexedVarArray{V,N,T}}) where {V,N,T} = T already_defined(var, index) = haskey(_data(var), index) @@ -59,43 +60,24 @@ function insertvar!( ::UnsafeInsert, index..., ) where {V,N,T} + clear_cache!(var) return var[index] = var.f(index...) end """ unsafe_insertvar!(var::indexedVarArray{V,N,T}, index...) -Insert a new variable with the given index withouth checking if the index is valid or +Insert a new variable with the given index without checking if the index is valid or already assigned. """ function unsafe_insertvar!(var::IndexedVarArray{V,N,T}, index...) where {V,N,T} return insertvar!(var, UnsafeInsert(), index...) end -joinex(ex1, ex2) = :($ex1..., $ex2...) -@generated function _active(idx::I, pat::P) where {I,P} - ids = fieldtypes(I) - ps = fieldtypes(P) - exs = [] - for i in 1:length(ids) - if ps[i] != Colon - if i > 2 - push!(exs, :(a1 = idx[$i],)) - else - push!(exs, :(idx[$i],)) - end - end - end - for i in 1:length(exs)-1 - exs[i+1] = joinex(exs[i], exs[i+1]) - end - return :(tuple($(exs[end])...)) -end - function build_cache!(cache, pat, sa::IndexedVarArray{V,N,T}) where {V,N,T} if isempty(cache) for v in keys(sa) - vred = _active(v, pat) + vred = _project_fixed(v, typeof(pat)) nv = get!(cache, vred, T[]) push!(nv, v) end @@ -103,13 +85,36 @@ function build_cache!(cache, pat, sa::IndexedVarArray{V,N,T}) where {V,N,T} return cache end -function _select_cached(sa::IndexedVarArray{V,N,T}, pat) where {V,N,T} - # TODO: Benchmark to find good cutoff-value for caching - # TODO: Return same type for type stability - length(_data(sa)) < 100 && return _select_gen(keys(_data(sa)), pat) +# Minimum number of entries before the index cache is used; below this a +# linear scan is assumed cheaper. Tune with set_cache_cutoff! or calibrate with +# benchmark/cutoff_benchmark.jl. +_CACHE_CUTOFF::Int = 100 + +""" + set_cache_cutoff!(n::Int) + +Set the minimum number of entries in an `IndexedVarArray` at which +selection switches from a linear scan to the pre-built +index cache. Smaller values favour caching; larger values favour the linear +scan for small arrays. Default: `100`. +""" +set_cache_cutoff!(n::Int) = (global _CACHE_CUTOFF = n; nothing) + +@generated function _is_cacheable_pattern(::Type{P}) where {P<:Tuple} + return :($(all(t == Colon || isfixed(t) for t in fieldtypes(P)))) +end + +function _select_cached( + sa::IndexedVarArray{V,N,T}, + pat, +)::Vector{T} where {V,N,T} + length(_data(sa)) < _CACHE_CUTOFF && + return collect(T, _select_gen(keys(_data(sa)), pat)) + _is_cacheable_pattern(typeof(pat)) || + return collect(T, _select_gen(keys(_data(sa)), pat)) cache = _getcache(sa, pat)::Dictionary{_decode_nonslices(sa, pat),Vector{T}} build_cache!(cache, pat, sa) - vals = _dropslices_gen(pat) + vals = _project_fixed(pat, typeof(pat)) return get!(cache, vals, T[]) end @@ -119,28 +124,6 @@ bin2int(v) = bin2int(v, Dim{length(v)}()) w = reverse([2^(i - 1) for i in 1:N]) return :(dot($w, v)) end - -function _dropslices(t::P) where {P} - return Tuple(ti for ti in t if ti != Colon()) -end - -@generated function _dropslices_gen(pat::P) where {P} - ps = fieldtypes(P) - exs = [] - for i in 1:length(ps) - if ps[i] != Colon - if i > 2 # Workaround for slurping of iterables (like strings) when passing to joinex - push!(exs, :(a2 = pat[$i],)) - else - push!(exs, :(pat[$i],)) - end - end - end - for i in 1:length(exs)-1 - exs[i+1] = joinex(exs[i], exs[i+1]) - end - return exs[end] -end """ _get_cache_index(::P) @@ -210,3 +193,21 @@ end function Base.lastindex(sa::IndexedVarArray, d) return last(sort(sa.index_names[d])) end + +# Override _view_matching_keys for IndexedVarArray parent: use index cache. +function _view_matching_keys( + v::SparseArraySlice{P,V,NF,MT}, +) where {P<:IndexedVarArray,V,NF,MT} + return _select_cached(v.parent, v.mask) +end + +# JuMP-efficient sum: build AffExpr directly via add_to_expression! for the +# standard VariableRef type. Custom AbstractVariableRef subtypes fall back to +# the generic slice sum implementation. +function Base.sum(v::SparseArraySlice{<:IndexedVarArray,VariableRef}) + result = zero(AffExpr) + for k in _view_matching_keys(v) + JuMP.add_to_expression!(result, v.parent[k]) + end + return result +end diff --git a/src/slice.jl b/src/slice.jl new file mode 100644 index 0000000..9492ef0 --- /dev/null +++ b/src/slice.jl @@ -0,0 +1,260 @@ +# Generated functions in this file cannot rely on helpers defined later in the +# include order, so selector classification stays local. +_is_exact_selector_type(::Type) = true +_is_exact_selector_type(::Type{<:Function}) = false +_is_exact_selector_type(::Type{<:AbstractRange}) = false +_is_exact_selector_type(::Type{<:Colon}) = false + +# Project tpl to the exact-match positions of MT. +@generated function _project_fixed(tpl::T, ::Type{MT}) where {T,MT} + inds = [ + i for i in 1:fieldcount(T) if _is_exact_selector_type(fieldtypes(MT)[i]) + ] + return :($(Expr(:tuple, [:(tpl[$i]) for i in inds]...))) +end + +# Project tpl to the non-exact (free) positions of MT. +@generated function _project_free(tpl::T, ::Type{MT}) where {T,MT} + inds = [ + i for + i in 1:fieldcount(T) if !_is_exact_selector_type(fieldtypes(MT)[i]) + ] + FT = Tuple{[fieldtypes(T)[i] for i in inds]...} + return :($(Expr(:tuple, [:(tpl[$i]) for i in inds]...))::$FT) +end + +# Reconstruct a full key T from the exact values in mask and the free key FT. +@generated function _reconstruct_key( + mask::MT, + free_key::FT, + ::Type{T}, +) where {MT,FT,T} + parts = Vector{Expr}(undef, fieldcount(T)) + fi = 1 + for i in 1:fieldcount(T) + if !_is_exact_selector_type(fieldtypes(MT)[i]) + parts[i] = :(free_key[$fi]) + fi += 1 + else + parts[i] = :(mask[$i]) + end + end + return :($(Expr(:tuple, parts...))::$T) +end + +""" + SparseArraySlice{P,V,NF,MT} + +A lazy, mask-filtered view of any `AbstractSparseArray`. `P` is the concrete +parent type, `NF` is the number of non-exact dimensions, and `MT` is the +mask type. Implements +`AbstractSparseArray{V,NF}`. + +Create via `slice(sa, mask...)`. +""" +struct SparseArraySlice{P<:AbstractSparseArray,V,NF,MT<:Tuple} <: + AbstractSparseArray{V,NF} + parent::P + mask::MT + _cache::Ref{Any} + function SparseArraySlice{P,V,NF,MT}(parent::P, mask::MT) where {P,V,NF,MT} + return new{P,V,NF,MT}(parent, mask, Ref{Any}(nothing)) + end +end + +function _keytype(::Type{<:SparseArraySlice{P,V,NF,MT}}) where {P,V,NF,MT} + return _free_keytype(MT, _keytype(P)) +end + +@generated function _free_keytype(::Type{MT}, ::Type{T}) where {MT,T} + inds = [ + i for + i in 1:fieldcount(T) if !_is_exact_selector_type(fieldtypes(MT)[i]) + ] + FT = Tuple{[fieldtypes(T)[i] for i in inds]...} + return :($FT) +end + +function _matches_free_key( + v::SparseArraySlice{P,V,NF,MT}, + free_key, +) where {P,V,NF,MT} + return _select_generated(_project_free(v.mask, MT), free_key) +end + +""" + slice(sa::AbstractSparseArray, mask...) + +Return a lazy `SparseArraySlice` over entries of `sa` matching `mask`. Use `:` +for wildcard dimensions, exact values for fixed dimensions, and predicates or +ranges for filtered dimensions. The result is an `AbstractSparseArray{V,NF}` +where `NF` is the number of non-exact dimensions. + +!!! note + + Slices materialise their matching keys lazily on first access and cache + them. Mutating through the slice invalidates the cache, but edits made + directly to the parent array afterwards are not tracked — re-create the + slice if the parent changes. + +# Example +```julia +v = slice(sa, :, "foo", :) # NF=2, two free dimensions +sum(v) +for (k, val) in pairs(v); ...; end +``` +""" +function slice(sa::AbstractSparseArray, mask...) + return _make_slice(sa, tuple(mask...)) +end + +@generated function _make_slice( + sa::P, + mask::MT, +) where {P<:AbstractSparseArray,MT<:Tuple} + N = ndims(P) + V = eltype(P) + fieldcount(MT) != N && return :(throw(BoundsError(sa, mask))) + NF = count(i -> !_is_exact_selector_type(fieldtypes(MT)[i]), 1:N) + return :(SparseArraySlice{$P,$V,$NF,$MT}(sa, mask)) +end + +# Default: linear scan with memoization. Subtypes may override for cached lookup. +function _view_matching_keys(v::SparseArraySlice{P,V,NF,MT}) where {P,V,NF,MT} + cached = v._cache[] + cached !== nothing && return cached::Vector{_keytype(P)} + T = _keytype(P) + keys_vec = collect(T, _select_gen(keys(_data(v.parent)), v.mask)) + v._cache[] = keys_vec + return keys_vec +end + +# Iterator traits +Base.IteratorSize(::Type{<:SparseArraySlice}) = Base.HasLength() +Base.IteratorEltype(::Type{<:SparseArraySlice}) = Base.HasEltype() +Base.eltype(::Type{<:SparseArraySlice{P,V}}) where {P,V} = V + +# Iteration: values only (AbstractArray semantics) +function Base.iterate(v::SparseArraySlice) + matching = _view_matching_keys(v) + isempty(matching) && return nothing + return (v.parent[matching[1]], (matching, 2)) +end + +function Base.iterate(v::SparseArraySlice, state::Tuple{Vector,Int}) + matching, pos = state + pos > length(matching) && return nothing + return (v.parent[matching[pos]], (matching, pos + 1)) +end + +# getindex by FT tuple: v[(f, c)] +function Base.getindex( + v::SparseArraySlice{P,V,NF,MT}, + free_key::Tuple, +) where {P,V,NF,MT} + length(free_key) == NF || throw(BoundsError(v, free_key)) + T = _keytype(P) + _matches_free_key(v, free_key) || return zero(V) + return v.parent[_reconstruct_key(v.mask, free_key, T)] +end + +# Splatted: v[f, c] or v[f] (NF==1) +function Base.getindex(v::SparseArraySlice{P,V,NF,MT}, idx...) where {P,V,NF,MT} + length(idx) == NF || throw(BoundsError(v, idx)) + return v[idx] +end + +# Forward mutation to parent array +function Base.setindex!( + v::SparseArraySlice{P,V,NF,MT}, + val, + free_key::Tuple, +) where {P,V,NF,MT} + length(free_key) == NF || throw(BoundsError(v, free_key)) + T = _keytype(P) + v.parent[_reconstruct_key(v.mask, free_key, T)] = val + v._cache[] = nothing # invalidate: a new matching key may have been added + return val +end + +# Splatted version +function Base.setindex!( + v::SparseArraySlice{P,V,NF,MT}, + val, + idx..., +) where {P,V,NF,MT} + length(idx) == NF || throw(BoundsError(v, idx)) + return setindex!(v, val, idx) +end + +function Base.size(::SparseArraySlice) + return error( + "`Base.size` is not implemented for `SparseArraySlice`. " * + "Use `length` for the number of entries.", + ) +end + +function Base.haskey( + v::SparseArraySlice{P,V,NF,MT}, + free_key::Tuple, +) where {P,V,NF,MT} + T = _keytype(P) + _matches_free_key(v, free_key) || return false + return haskey(_data(v.parent), _reconstruct_key(v.mask, free_key, T)) +end + +Base.length(v::SparseArraySlice) = length(_view_matching_keys(v)) + +function Base.keys(v::SparseArraySlice{P,V,NF,MT}) where {P,V,NF,MT} + return [_project_free(k, MT) for k in _view_matching_keys(v)] +end + +Base.values(v::SparseArraySlice) = [v.parent[k] for k in _view_matching_keys(v)] + +Base.eachindex(v::SparseArraySlice) = keys(v) + +function Base.pairs(v::SparseArraySlice{P,V,NF,MT}) where {P,V,NF,MT} + return [_project_free(k, MT) => v.parent[k] for k in _view_matching_keys(v)] +end + +function Base.firstindex(v::SparseArraySlice{P,V,NF,MT}, d) where {P,V,NF,MT} + return minimum(_project_free(k, MT)[d] for k in _view_matching_keys(v)) +end +function Base.lastindex(v::SparseArraySlice{P,V,NF,MT}, d) where {P,V,NF,MT} + return maximum(_project_free(k, MT)[d] for k in _view_matching_keys(v)) +end + +function Base.sum(v::SparseArraySlice{P,V}) where {P,V} + ks = _view_matching_keys(v) + isempty(ks) && return zero(V) + return sum(v.parent[k] for k in ks) +end + +function Base.summary(io::IO, v::SparseArraySlice) + num_entries = length(v) + return print( + io, + "SparseArraySlice with ", + num_entries, + isone(num_entries) ? " entry" : " entries", + " matching ", + v.mask, + ) +end + +function Base.show(io::IO, ::MIME"text/plain", v::SparseArraySlice) + summary(io, v) + if !iszero(length(v)) + println(io, ":") + entries = pairs(v) + if length(entries) > 20 + show(io, first(entries, 10)) + print(io, "\n ⋮\n") + show(io, last(entries, 10)) + else + show(io, entries) + end + end +end + +Base.show(io::IO, v::SparseArraySlice) = summary(io, v) diff --git a/src/sparsearray.jl b/src/sparsearray.jl index 5f6d9a7..b9c4d0f 100644 --- a/src/sparsearray.jl +++ b/src/sparsearray.jl @@ -115,3 +115,24 @@ function SparseArray{T,N,K}() where {T,N,K<:NTuple{N,Any}} end _data(sa::SparseArray) = sa.data + +# ------------------------------------------------------------------------------ +# _keytype interface +# Each concrete AbstractSparseArray subtype must implement _keytype(::Type{<:SA}). +# ------------------------------------------------------------------------------ + +""" + _keytype(sa) / _keytype(::Type{<:AbstractSparseArray}) + +Return the key tuple type used by the sparse array. Required by `slice` and +broadcasting. Implement `_keytype(::Type{MySA})` for every concrete subtype. +""" +function _keytype(::Type{<:AbstractSparseArray}) + return error( + "_keytype not implemented for this AbstractSparseArray subtype", + ) +end +_keytype(sa::AbstractSparseArray) = _keytype(typeof(sa)) +_keytype(::Type{<:SparseArray{T,N,K}}) where {T,N,K} = K + +Base.haskey(sa::AbstractSparseArray, k) = haskey(_data(sa), k) diff --git a/test/Project.toml b/test/Project.toml index 7f33635..7988d16 100644 --- a/test/Project.toml +++ b/test/Project.toml @@ -2,4 +2,5 @@ Dictionaries = "85a47980-9c8c-11e8-2b9f-f7ca1fa99fb4" HiGHS = "87dc4568-4c63-4d18-b0c0-bb2238e4078b" JuMP = "4076af6c-e467-56ae-b986-b466b2749572" +SparseVariables = "2749762c-80ed-4b14-8f33-f0736679b02b" Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40" diff --git a/test/runtests.jl b/test/runtests.jl index 742a538..2873510 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -49,7 +49,7 @@ end @variable( m, - car_vars[maker = cars, year = years, color = colors, kms = kms]; + car_vars[maker=cars, year=years, color=colors, kms=kms]; container = IndexedVarArray ) @test typeof(car_vars) == @@ -77,6 +77,10 @@ end @test car_cost["bmw", 2001] == 200 @test car_cost["bmw", 2003] == 0 + @test car_cost[endswith("s"), <(2000)] isa SparseArraySlice + @test length(car_cost[endswith("s"), <(2000)]) == 1 + @test car_cost[endswith("s"), <(2000)]["lotus", 1957] == 500 + @test car_cost[endswith("s"), <(2000)]["bmw", 2001] == 0 @test length(car_cost) == 5 @test car_cost["lotus", 1957] == 500 @@ -114,14 +118,14 @@ end (; car_cost) = testdata1() m = Model() - @variable(m, y[c = cars, i = years]; container = IndexedVarArray) + @variable(m, y[c=cars, i=years]; container = IndexedVarArray) for (c, i) in collect(keys(car_cost)) insertvar!(y, c, i) end @test typeof(y) == IndexedVarArray{VariableRef,2,Tuple{String,Int}} - @variable(m, w[c = cars, i = years], Bin; container = IndexedVarArray) + @variable(m, w[c=cars, i=years], Bin; container = IndexedVarArray) for (c, i) in collect(keys(car_cost)) insertvar!(w, c, i) end @@ -164,7 +168,7 @@ end m = Model() (; cars, year, car_cost) = testdata1(false) - @variable(m, z[cars = cars, year = year]; container = IndexedVarArray) + @variable(m, z[cars=cars, year=year]; container = IndexedVarArray) for (cr, yr) in keys(car_cost) insertvar!(z, cr, yr) @@ -181,6 +185,10 @@ end # Slicing and lookup @test length(z["bmw", :]) == 2 @test length(z[:, 2001]) == 2 + @test z[endswith("w"), isodd] isa SparseArraySlice + @test length(z[endswith("w"), isodd]) == 1 + @test haskey(z[endswith("w"), isodd], ("bmw", 2001)) + @test !haskey(z[endswith("w"), isodd], ("bmw", 2002)) @test typeof(z["bmw", 2001]) == VariableRef @test z["bmw", 20] == 0 @@ -190,7 +198,7 @@ end @test length(z) == 5 # Alternative constructor - @variable(m, z2[cars = cars, year = year], container = IndexedVarArray) + @variable(m, z2[cars=cars, year=year], container = IndexedVarArray) for k in keys(car_cost) insertvar!(z2, k...) end @@ -201,7 +209,7 @@ end @variable( m, - z3[cars = cars, year = years, color = colors, km = kms]; + z3[cars=cars, year=years, color=colors, km=kms]; container = IndexedVarArray ) for k in indices @@ -231,7 +239,7 @@ end (; cars, year, car_cost) = testdata1(false) m = Model() - @variable(m, y[car = cars, year = year] >= 0; container = IndexedVarArray) + @variable(m, y[car=cars, year=year] >= 0; container = IndexedVarArray) for c in cars insertvar!(y, c, 2002) end @@ -271,7 +279,7 @@ end # Test JuMP Extension m = Model() - @variable(m, x[i = 1:3, j = 100:102] >= 0, container = IndexedVarArray) + @variable(m, x[i=1:3, j=100:102] >= 0, container = IndexedVarArray) @test length(x) == 0 insertvar!(x, 1, 100) @test length(x) == 1 @@ -316,7 +324,7 @@ end m = Model() @variable( m, - x[i = 1:3, j = 100:102] >= 0, + x[i=1:3, j=100:102] >= 0, Mocking(), container = IndexedVarArray ) @@ -329,3 +337,278 @@ end @test sum(x) == sum(x[:, :]) @test typeof(sum(x)) <: GenericAffExpr{Float64,MockVariableRef} end + +const _test_sa = testdata_sa() + +@testset "SparseArraySlice cache invalidation" begin + sa = SparseArray( + Dict( + ("ford", 2000) => 100, + ("ford", 2001) => 150, + ("bmw", 2001) => 200, + ), + ) + + v = slice(sa, "ford", :) + @test length(v) == 2 + + # Add a matching entry through the slice + v[(2002,)] = 175 + @test length(v) == 3 # cache invalidated + recomputed + @test v[(2002,)] == 175 + @test Set(keys(v)) == Set([(2000,), (2001,), (2002,)]) + @test sort(values(v)) == [100, 150, 175] + @test sum(v) == 100 + 150 + 175 + + # Add another matching entry using splatted setindex! + v[2003] = 50 + @test length(v) == 4 + @test v[2003] == 50 + @test sa["ford", 2003] == 50 # parent updated + @test length(sa) == 5 + + # Update an existing free key: count unchanged, value refreshed + v[(2000,)] = 999 + @test length(v) == 4 + @test v[(2000,)] == 999 + + # Mutating through the slice must not affect a non-matching dim + @test v[(2001,)] == 150 + @test sa["bmw", 2001] == 200 +end + +@testset "SparseArraySlice on SparseArray" begin + sa = _test_sa + + # _keytype + @test SV._keytype(sa) == Tuple{String,Int} + + # Types + v = slice(sa, "ford", :) + @test v isa SparseArraySlice + @test v isa SV.AbstractSparseArray + @test SV._keytype(v) == Tuple{Int} # projected key type + @test eltype(v) == Int + + # length + @test length(slice(sa, "ford", :)) == 2 + @test length(slice(sa, :, 2001)) == 2 + @test length(slice(sa, :, :)) == 5 + @test length(slice(sa, "xxx", :)) == 0 + + # keys / values / eachindex / pairs + ks = sort(keys(v)) + @test ks == [(2000,), (2001,)] + @test eachindex(v) == keys(v) + @test sort(values(v)) == [100, 150] + ps = Dict(pairs(v)) + @test ps[(2000,)] == 100 + @test ps[(2001,)] == 150 + + # getindex + @test v[(2000,)] == 100 # FT-tuple + @test v[NTuple{1,Any}((2001,))] == 150 # NTuple{NF,Any} + @test v[2000] == 100 # splatted (NF==1) + @test v[2001] == 150 + + v2 = slice(sa, :, :) # NF==2 + @test v2["ford", 2000] == 100 + @test v2["bmw", 2002] == 300 + + # haskey + @test haskey(v, (2000,)) + @test !haskey(v, (1999,)) + + # sum + @test sum(v) == 250 + @test sum(slice(sa, :, 2001)) == 350 + @test sum(slice(sa, "xxx", :)) == 0 + + # firstindex / lastindex (d = free-dimension index, 1:NF) + @test SV.firstindex(v, 1) == 2000 + @test SV.lastindex(v, 1) == 2001 + + # iteration (values only) + @test sum(val for val in v) == 250 + @test Base.IteratorSize(typeof(v)) == Base.HasLength() + @test Base.IteratorEltype(typeof(v)) == Base.HasEltype() + + # setindex + v2["bmw", 2002] = 200 + @test v2["bmw", 2002] == 200 + @test sa["bmw", 2002] == 200 + v2[("bmw", 2002)] = 300 + @test v2[("bmw", 2002)] == 300 + @test sa[("bmw", 2002)] == 300 + + # show / summary + @test occursin("SparseArraySlice", sprint(summary, v)) + @test occursin("matching (\"ford\", Colon())", sprint(summary, v)) + @test occursin("(2000,) => 100", sprint(show, MIME("text/plain"), v)) + @test occursin("(2001,) => 150", sprint(show, MIME("text/plain"), v)) + + # wrong mask length + @test_throws BoundsError slice(sa, "ford", :, :) + + # empty slice + ve = slice(sa, "xxx", :) + @test length(ve) == 0 + @test isempty(keys(ve)) + @test isempty(values(ve)) + @test sum(ve) == 0 +end + +@testset "SparseArraySlice on IndexedVarArray" begin + (; cars, year, car_cost) = testdata1(false) + m = Model() + @variable(m, x[c=cars, y=year]; container = IndexedVarArray) + for k in keys(car_cost) + insertvar!(x, k...) + end + + v = slice(x, "ford", :) + + # type and _keytype + @test v isa SparseArraySlice + @test SV._keytype(x) == Tuple{String,Int} + @test SV._keytype(v) == Tuple{Int} + + # length / keys + @test length(v) == 2 + @test sort(keys(v)) == [(2000,), (2001,)] + + # JuMP sum returns AffExpr + @test sum(v) isa AffExpr + @test length(sum(v).terms) == 2 + + @test sum(slice(x, :, 2001)) isa AffExpr + @test length(sum(slice(x, :, 2001)).terms) == 2 + + # empty JuMP slice sum returns zero(AffExpr) + @test sum(slice(x, "xxx", :)) == zero(AffExpr) +end + +@testset "Broadcasting SparseArray" begin + sa = _test_sa + + # scalar broadcast + r = sa .* 2 + @test r isa SparseArray + @test length(r) == 5 + @test r["ford", 2000] == 200 + @test r["lotus", 1957] == 1000 + + r2 = 2 .* sa + @test r2["bmw", 2001] == 400 + + r3 = sa .+ 10 + @test r3["ford", 2001] == 160 + + # element-wise binary + r4 = sa .+ sa + @test r4["ford", 2000] == 200 + @test r4["bmw", 2002] == 600 + + # function broadcast + r5 = sqrt.(sa .* 1.0) + @test r5 isa SparseArray + @test r5["ford", 2000] ≈ sqrt(100.0) + + # result type + @test Base.BroadcastStyle(typeof(sa)) isa SV.SparseBroadcastStyle + + # key mismatch error + sa2 = SparseArray(Dict(("a", 1) => 1)) + @test_throws ArgumentError sa .+ sa2 + + # empty array broadcast + empty_sa = SparseArray(Dictionary{Tuple{String,Int},Int}()) + r_empty = empty_sa .* 2 + @test r_empty isa SparseArray + @test length(r_empty) == 0 +end + +@testset "Broadcasting SparseArray with JuMP scalar" begin + sa = _test_sa # SparseArray{Int} keyed by (String, Int) + m = Model() + @variable(m, t) + + # sa .* scalar-variable + ra = sa .* t + @test ra isa SparseArray + @test eltype(ra) <: JuMP.AbstractJuMPScalar + @test JuMP.isequal_canonical(ra["ford", 2000], 100 * t) + @test length(ra) == length(sa) + + # scalar-variable on the left + ra2 = t .* sa + @test JuMP.isequal_canonical(ra2["bmw", 2001], 200 * t) + + # sa .+ scalar-expression + rb = sa .+ (2t + 1) + @test rb isa SparseArray + @test JuMP.isequal_canonical(rb["ford", 2000], 100 + 2t + 1) +end + +@testset "Broadcasting SparseArraySlice" begin + sa = _test_sa + + v = slice(sa, "ford", :) + r = v .* 2 + @test r isa SparseArray + @test length(r) == 2 + @test r[(2000,)] == 200 + @test r[(2001,)] == 300 + + # slice .+ slice (same keys) + r2 = v .+ v + @test r2[(2000,)] == 200 + + # NF=2 slice broadcast + v2 = slice(sa, :, :) + r3 = v2 .* 3 + @test r3["ford", 2000] == 300 + @test r3["lotus", 1957] == 1500 + + # key mismatch between two slices + vbmw = slice(sa, "bmw", :) + @test_throws ArgumentError v .+ vbmw +end + +@testset "Broadcasting IndexedVarArray" begin + (; cars, year, car_cost) = testdata1(false) + m = Model() + @variable(m, x[c=cars, y=year] >= 0; container = IndexedVarArray) + for k in keys(car_cost) + insertvar!(x, k...) + end + @objective(m, Min, sum(x[c, y] for (c, y) in keys(car_cost))) + @constraint(m, sum(x[:, :]) == 1) + set_optimizer(m, HiGHS.Optimizer) + set_optimizer_attribute(m, MOI.Silent(), true) + optimize!(m) + + # value.(iva) → SparseArray + vals = value.(x) + @test vals isa SparseArray + @test length(vals) == length(x) + @test isapprox(sum(values(vals)), 1.0; atol = 1e-6) + + # value.(slice) → SparseArray with projected keys + vslice = value.(slice(x, "ford", :)) + @test vslice isa SparseArray + @test length(vslice) == 2 + @test eltype(vslice) == Float64 + + # iva .+ iva → SparseArray{AffExpr} + aff = x .+ x + @test aff isa SparseArray + @test length(aff) == length(x) + @test first(values(aff)) isa AffExpr + + # key mismatch error + m2 = Model() + @variable(m2, y[c=["lotus"], yr=[1957]]; container = IndexedVarArray) + insertvar!(y, "lotus", 1957) + @test_throws ArgumentError x .+ y +end diff --git a/test/testdata.jl b/test/testdata.jl index 5da7620..bb59522 100644 --- a/test/testdata.jl +++ b/test/testdata.jl @@ -26,6 +26,18 @@ function testdata(N = 998) ) end +function testdata_sa() + return SparseVariables.SparseArray( + Dict( + ("ford", 2000) => 100, + ("ford", 2001) => 150, + ("bmw", 2001) => 200, + ("bmw", 2002) => 300, + ("lotus", 1957) => 500, + ), + ) +end + function testdata1(addlotus = true) cars = ["ford", "bmw", "opel"] year = [2000, 2001, 2002, 2003] diff --git a/tutorial/sparse_tutorial.jl b/tutorial/sparse_tutorial.jl index bc2c05e..d20c8c3 100644 --- a/tutorial/sparse_tutorial.jl +++ b/tutorial/sparse_tutorial.jl @@ -6,7 +6,7 @@ using InteractiveUtils # This Pluto notebook uses @bind for interactivity. When running this notebook outside of Pluto, the following 'mock version' of @bind gives bound variables a default value (instead of an error). macro bind(def, element) - quote + return quote local iv = try Base.loaded_modules[Base.PkgId( Base.UUID("6e696c72-6542-2067-7265-42206c756150"), @@ -186,8 +186,8 @@ begin @constraint( m, sum( - x[f, c, p, t] for (f, c, p, t) in - filter(i -> i[1] == f̄ && i[4] == t̄, indices) + x[f, c, p, t] for + (f, c, p, t) in filter(i -> i[1] == f̄ && i[4] == t̄, indices) ) ≤ U[f̄, t̄] ) end @@ -244,8 +244,8 @@ begin @constraint( m, sum( - x[(f, c, p, t)] for (f, c, p, t) in - filter(i -> i[1] == f̄ && i[4] == t̄, indices) + x[(f, c, p, t)] for + (f, c, p, t) in filter(i -> i[1] == f̄ && i[4] == t̄, indices) ) ≤ U[f̄, t̄] ) end @@ -377,7 +377,7 @@ md" # ╔═╡ c49a3599-65fd-442b-b5c9-625b87e05efa begin - res = DataFrame(Method = Symbol[], NC = Int[], Time = Float64[]) + res = DataFrame(; Method = Symbol[], NC = Int[], Time = Float64[]) @progress for nc in 5:5:50 for method in [ model_standard, @@ -398,10 +398,10 @@ res # ╔═╡ 5ca68304-0c21-4344-a92b-0594c04674a4 function plot(df, x = :NC, y = :Time) - CairoMakie.activate!(type = "svg") + CairoMakie.activate!(; type = "svg") return draw( data(df) * - mapping(x, y, color = :Method, marker = :Method) * + mapping(x, y; color = :Method, marker = :Method) * (visual(Lines) + visual(Scatter)), ) end @@ -416,7 +416,7 @@ md" # ╔═╡ 1d68c8c0-1dbc-4d8b-97ae-3db1d2b06a4f begin - sparsity = DataFrame(Method = Symbol[], DP = Float64[], Time = Float64[]) + sparsity = DataFrame(; Method = Symbol[], DP = Float64[], Time = Float64[]) @progress for dp in 0.05:0.05:1.0 for method in [ model_standard, @@ -442,8 +442,12 @@ plot(sparsity, :DP, :Time) # ╔═╡ deefc30f-4846-49db-8ce8-69b4aec14924 begin - large = - DataFrame(Method = Symbol[], nc = Int[], vars = Int[], Time = Float64[]) + large = DataFrame(; + Method = Symbol[], + nc = Int[], + vars = Int[], + Time = Float64[], + ) @progress for nc in 500:500:5000 for method in [model_incremental, model_sparse] GC.gc()