Slow array access during labeling algorithm

Specific Domains Optimization (Mathematical)
, ,
1

Hi all,

I switched from Python to Julia and are currently implementing a column generation algorithm for a vehicle routing problem. To generate a column (route) with negative reduced cost, I use a labeling algorithm.
The simplified version has the following structure:

struct LabelInstance 
    ...
end 

struct Label
    cost             :: Float64 
    setReps          :: SBitSet{1,UInt128}
    node             :: Int64
    time             :: Float64 
    load             :: Int64
end

function label_dominance(label_1::Label, label_2::Label, inst::LabelInstance)
    # check whether label_1 dominates label_2
    ...
end 

function check_dom_and_add_label!(all_labels::Matrix{Vector{Label}}, new_label::Label, inst::LabelInstance)
    # check whether new_label is dominated by existing labels 
    node= new_label.node
    
    lb= inst.lbTW[node]
    time_label= label.time
    for t=lb:time_label
        labels= all_labels[t,node]
        if isempty(labels)
            continue 
        end 
        for i in eachindex(labels)
            exis_label= labels[i] 
            if label_dominance(exis_label, new_label, inst)
                return false     
            end 
        end
    end 

    # check whether new_label dominates existing labels 
    ub= inst.ubTW[node]
    for t=time_label:ub
        labels= all_labels[t,node]
        if !(isempty(labels))
            deleteat!(labels, (i for i in eachindex(labels) if label_dominance(new_label, labels[i], inst))) 
        end
    end

    # add new_label
    push!(all_labels[time_label, node], new_label)                   
    
    return true
end

function extend_labels_node!(node::Int64, time::Int64, all_labels::Matrix{Vector{Label}}, inst::LabelInstance)
    for label in all_labels[time, node] 
        for node_ext=1:inst.numNodes 
            label_ext, feas_ext= extend_label(label,node_ext)     # extend label 
            if feas_ext                                           # label_ext feasible 
                label_add= check_dom_and_add_label!(all_labels, new_label, inst)
            end 
        end 
    end
    
    return nothing
end 

function solve_pricing_problem(oInst::LabelInstance)
    # initialize matrix with labels
    all_labels::Matrix{Vector{Label}}(undef, inst.numNodes, inst.ubTW[end]+1) 
    for i=1:inst.numNodes
        for t=inst.lbTW[i]:inst.ubTW[i]
            all_labels[t,i]= Vector{Label}() # stores labels for node i and time t
        end 
    end 
    
    # initialize start label
    start_node= 1             
    start_label= Label(0.0,SBitSet{1,UInt128}(start_node), start_node, inst.lbTW[start_node], inst.load[start_node])
    push!(all_labels[start_label.time,start_node], start_label)

    # extend labels
    for t=inst.lbTW[start_node]:inst.ubTW[end]
        for node=1:inst.numNodes 
            if t >= inst.lbTW[node] && t <= oInst.ubTW[node]
                extend_labels_node!(node, time, all_labels, inst) # extend labels in all_labels[t,node]
            end
        end
    end  

    return all_labels
end

Labels representing partial routes are stored in all_labels= Matrix{Vector{Label}} for each node and time point. The function check_dom_and_add_label!(all_labels::Matrix{Vector{Label}}, new_label::Label, inst::LabelInstance) checks if a new label is dominated by existing labels, removes dominated labels, and adds the new label to all_labels. This function is called many times during the algorithm, and each time the function loops trough all labels in a column of all_labels.

After profiling the code, the profiler indicated that almost all time is spend within check_dom_and_add_label! . However, about 80% of that time is spend on get_index (e.g., 60% of the time is spend on exis_label= labels[i] ). This seems counterintuitive for me, as I would expect more time to be spend inside the function label_dominance.

Can anyone explain why accessing the elements of an array is so time consuming? Could this potentially indicate other problems in my code?

Thanks!

2

Are you talking about column generation in optimization?

3

Yes indeed. I use the labeling algorithm to find a column with negative reduced cost.

4

I think you might get more and better answers if you provide a complete, runnable MWE as described in item 4 of this post.

1 Like
5

Also tagging with the specific domain tag for optimization Optimization (Mathematical) - Julia Programming Language would help route it to more people who might be able to help.

1 Like
6

Hi @Walter, nothing stands out as being the cause of the problem, and I wouldn’t expect the labels[i] to be problematic.

How did you profile the code?

Can you share a reproducible example?

I notice your check_dom_and_add_label! function has:

    time_label= label.time

but label doesn’t appear to be a local variable in the function.

Have you verified that your code is working as expected (excepting the performance issue) in a new Julia instance? I think you might be unintentionally re-using global variables in places.

1 Like
7

Hi @odow,

Thank you for your response. I now also recognise the mistake in the example (I subtracted it from my original code). However, my code does not use global variables. Below is an reproducible example:

using StaticBitSets
using Distances 

struct LabelInstance 
    numNodes         :: Int64               # number of nodes
    dem              :: Vector{Int64}       # demand 
    dur              :: Vector{Int64}       # duration 
    lbTW             :: Vector{Int64}       # lower bound of time window
    ubTW             :: Vector{Int64}       # upper bound of time window 
    travelTimes      :: Matrix{Int64}       # travel times 
    cost             :: Matrix{Float64}     # cost matrix       
    cap              :: Int64               # capacity of vehicle
end 

struct Label
    cost             :: Float64             # reduced cost of the label
    setReps          :: SBitSet{1,UInt128}  # set representation partial route 
    node             :: Int64               # node of the label
    time             :: Int64               # time of the label
    load             :: Int64               # load of the label
    complete         :: Bool                # complete label
    timePrevLab      :: Int64               # time of the previous label (used for backtracking)
    nodePrevLab      :: Int64               # node of the previous label (used for backtracking)
    idxPrevLabel     :: Int64               # index of the previous label (used for backtracking)
    # other attributes of the label, currently not used 
    ##=
    a1               :: Int64 
    a2               :: Int64
    a3               :: Int64
    a4               :: Int64 
    a5               :: Int64
    a6               :: Bool 
    a7               :: Bool  
    ##=#
end

function label_dominance(label_1::Label, label_2::Label, inst::LabelInstance; eps::Float64=0.000001)
    # check whether label_1 dominates label_2
    if label_1.time > label_2.time && label_2.complete == false # for current structure current condution always true
        return false 
    elseif !(label_1.setReps ⊆ label_2.setReps) 
        return false
    elseif label_1.cost > label_2.cost + eps
        return false  
    end 
    # some other checks (using inst and some attributes from a1-a7)
    # ...
    return true 
end 

@inline function single_node_bitset_1xUInt128(n::Integer)
    @assert 1 ≤ n ≤ 128
    SBitSet{1,UInt128}((UInt128(1) << (n - 1),))
end

function dummy_label()
    #Label(0.0, SBitSet{1,UInt128}(1), 0, 0, 0, false, 0,0,0)
    Label(0.0, SBitSet{1,UInt128}(1), 0, 0, 0, false, 0,0,0,0,0,0,0,0,false, false)
end 
function extend_label(label::Label, node_ext::Integer, inst::LabelInstance, idx_label::Int64)
    # extend label to node_ext 
    arr_time= label.time + inst.dur[label.node] + inst.travelTimes[label.node, node_ext] 
    start_time= max(arr_time, inst.lbTW[node_ext])                                       
    if start_time > inst.ubTW[node_ext] # no arrival possible within time-window
        return false, dummy_label()
    end
    
    load= label.load + inst.dem[node_ext]
    if load > inst.cap                    # load exceeds vehicle capacity
        return false, dummy_label() 
    end

    if node_ext in label.setReps    # node_ext already in patial route 
        return false, dummy_label()
    end

    set_repr= label.setReps ∪ single_node_bitset_1xUInt128(node_ext) 
    cost= label.cost + inst.cost[label.node, node_ext]  
    complete = node_ext == inst.numNodes ? true : false 

    #label_ext= Label(cost, set_repr, node_ext, start_time, load, complete, label.time, label.node, idx_label)
    label_ext= Label(cost, set_repr, node_ext, start_time, load, complete, label.time, label.node, idx_label, 0,0,0,0,0,false,false)
    
    return true, label_ext 
end 
            
function check_dom_and_add_label!(all_labels::Matrix{Vector{Label}}, new_label::Label, inst::LabelInstance)
    # check whether new_label is dominated by existing labels 
    node= new_label.node
    
    lb= inst.lbTW[node]
    time_label= new_label.time
    for t=lb:time_label
        labels= all_labels[t+1,node]
        if isempty(labels)
            continue 
        end 
        for i in eachindex(labels)
            exis_label= labels[i] 
            if label_dominance(exis_label, new_label, inst)
                return false     
            end 
        end
    end 

    # check whether new_label dominates existing labels 
    ub= inst.ubTW[node]
    for t=time_label:ub
        labels= all_labels[t+1,node]
        if !(isempty(labels))
            deleteat!(labels, (i for i in eachindex(labels) if label_dominance(new_label, labels[i], inst))) 
        end
    end

    # add new_label
    push!(all_labels[time_label+1, node], new_label)                   
    
    return true
end

function extend_labels_node!(node::Int64, time::Int64, all_labels::Matrix{Vector{Label}}, inst::LabelInstance, neg_cost_labels::Vector{Label}; eps::Float64=0.000001)
    for (idx,label) in enumerate(all_labels[time+1, node]) 
        for node_ext=2:inst.numNodes 
            feas_label, new_label= extend_label(label, node_ext, inst, idx)  # extend label 
            if feas_label == true                                            # label_ext feasible 
                label_add= check_dom_and_add_label!(all_labels, new_label, inst)
                (label_add == true && node_ext == inst.numNodes && new_label.cost < -eps) ? push!(neg_cost_labels, new_label) : nothing
            end 
        end 
    end
    
    return nothing
end 

function solve_pricing_problem(inst::LabelInstance, time_limit::Float64=60.0)
    # initialize matrix with labels
    all_labels= Matrix{Vector{Label}}(undef,inst.ubTW[end]+1, inst.numNodes) 
    for i=1:inst.numNodes
        for t=inst.lbTW[i]:inst.ubTW[i]
            all_labels[t+1,i]= Vector{Label}() # stores labels for node i and time t
        end 
    end 
    neg_cost_labels= Vector{Label}() # stores labels with negative reduced cost
    
    # initialize start label
    start_node= 1             
    #start_label= Label(0.0,SBitSet{1,UInt128}(start_node), start_node, inst.lbTW[start_node], inst.dem[start_node], false, 0,0,0)
    start_label= Label(0.0,SBitSet{1,UInt128}(start_node), start_node, inst.lbTW[start_node], inst.dem[start_node], false, 0,0,0,0,0,0,0,0,false,false)
    push!(all_labels[start_label.time+1,start_node], start_label)

    # extend labels
    time_start_alg= time()
    for t=inst.lbTW[start_node]:inst.ubTW[end]
        if time()- time_start_alg < time_limit && length(neg_cost_labels) < 100
            for node=1:inst.numNodes 
                if t >= inst.lbTW[node] && t <= inst.ubTW[node] 
                    extend_labels_node!(node, t, all_labels, inst, neg_cost_labels) # extend labels in all_labels[t,node]
                end
            end
        end 
    end  

    return all_labels, neg_cost_labels
end 

function label_instance()
    dem= [0, 10, 7, 13, 19, 26, 3, 5, 9, 16, 16, 12, 19, 23, 20, 8, 19, 2, 12, 17, 9, 11, 18, 29, 3, 6, 17, 16, 16, 9, 21, 27, 23, 11, 14, 8, 5, 8, 16, 31, 9, 5, 5, 7, 18, 16, 1, 27, 36, 30, 13, 0]
    dur= [0, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 0]
    lb= [0, 707, 143, 527, 678, 34, 415, 331, 404, 400, 577, 206, 228, 690, 33, 175, 272, 733, 377, 269, 581, 214, 409, 206, 704, 817, 588, 104, 114, 190, 259, 152, 660, 45, 529, 686, 42, 606, 302, 34, 360, 396, 26, 620, 233, 30, 515, 85, 773, 501, 547, 0]
    ub= [1000, 848, 282, 584, 801, 209, 514, 410, 481, 497, 632, 325, 345, 827, 243, 300, 373, 870, 434, 378, 666, 331, 494, 325, 847, 956, 667, 255, 255, 313, 354, 275, 777, 200, 614, 813, 208, 693, 405, 224, 437, 511, 172, 705, 340, 189, 628, 250, 906, 540, 642, 1000]
    xy_coord= [35 35; 41 49; 35 17; 55 45; 55 20; 15 30; 25 30; 20 50; 10 43; 55 60; 30 60; 20 65; 50 35; 30 25; 15 10; 30 5; 10 20; 5 30; 20 40; 15 60; 45 65; 45 20; 45 10; 55 5; 65 35; 65 20; 45 30; 35 40; 41 37; 64 42; 40 60; 31 52; 35 69; 53 52; 65 55; 63 65; 2 60; 20 20; 5 5; 60 12; 40 25; 42 7; 24 12; 23 3; 11 14; 6 38; 2 48; 8 56; 13 52; 6 68; 47 47; 35 35]
    travel_time=ceil.(Int64,pairwise(Euclidean(), xy_coord', dims=2))
    ass= [-56.1, 88.1, 13.45, 102.1, 106.1, 8.0, 80.1, 29.0, 110.1, 122.1, 108.1, 12.45, 22.0, 80.1, 39.0, 16.0, 48.0, 118.1, 10.55, 26.0, 120.1, 17.55, 110.1, 28.45, 116.1, 124.1, 80.1, 10.55, 7.0, 27.55, 29.0, 17.0, 124.1, 34.0, 130.1, 140.1, 28.55, 100.1, 22.55, 42.55, 22.0, 45.0, 12.55, 126.1, 10.0, 14.0, 128.1, 34.0, 112.1, 144.1, 90.1, 0.0]
    cost= travel_time .- ass 
    cap_vehicle= 500 
    num_nodes= length(dur)

    inst= LabelInstance(num_nodes, dem, dur, lb, ub, travel_time, cost, cap_vehicle)

    return inst 
end 

function main(time_limit::Float64)
    inst= label_instance()
    s= time()
    all_labels, neg_cost_labels= solve_pricing_problem(inst, time_limit)
    println(s-time())

    return nothing 
end 

main(5.0)

using Profile
using ProfileView
Profile.clear()
Profile.init()
VSCodeServer.@profview main(50.0)

The profiler gives the following output:

large_label
large_label967×667 28.8 KB
.
Here, the #label_dominance#11 corresponds to the fraction getindex that is indicated for my original code. Of the 25% time spend in deleteat!, 19% is spend in getindex (specifically line eval(:(getindex(A::Array, i1::Int) = arrayref($(Expr(:boundscheck)), A, i1))))

When I reduce the label by removing attributes a1,…,a7 (this changes nothing to the algorithm), surprisingly, the run time decreases dramatically from 40 to 20 seconds and the profiler looks like this:

image
image979×649 29.7 KB

8

Here, the #label_dominance#11 corresponds to the fraction getindex that is indicated for my original code

I don’t follow this part.

I haven’t used the profiler with VSCode. I prefer using PProf.jl. If you do:

using PProf
PProf.@pprof main(5.0)

then open the webpage, you can click around (view>graph, then click the check_dom_and_add_label! node, then view>source) until you find:

image
image2344×1134 227 KB

That seems to be what I expect. It’s probably something like the exis_label= labels[i] is inlined into the next line, and the Julia profiler doesn’t pull back the correct line information or something.

1 Like
9

Thank you very much @odow. This definitely helps.

This profiler seems to work much better, and now most of the time is attributed to the function label_dominance where I expect it. With the sentence “#label_dominance#11 corresponds to the fraction getindex that is indicated for my original code”, I referred to the fact that initially the profiler in VS code indicated a similar fraction of time spend on the getindex function.

When I run the code block

using StaticBitSets
using Distances 

struct LabelInstance 
    numNodes         :: Int64               # number of nodes
    dem              :: Vector{Int64}       # demand 
    dur              :: Vector{Int64}       # duration 
    lbTW             :: Vector{Int64}       # lower bound of time window
    ubTW             :: Vector{Int64}       # upper bound of time window 
    travelTimes      :: Matrix{Int64}       # travel times 
    cost             :: Matrix{Float64}     # cost matrix       
    cap              :: Int64               # capacity of vehicle
end 

struct Label
    cost             :: Float64             # reduced cost of the label
    setReps          :: SBitSet{1,UInt128}  # set representation partial route 
    node             :: Int64               # node of the label
    time             :: Int64               # time of the label
    load             :: Int64               # load of the label
    complete         :: Bool                # complete label
    timePrevLab      :: Int64               # time of the previous label (used for backtracking)
    nodePrevLab      :: Int64               # node of the previous label (used for backtracking)
    idxPrevLabel     :: Int64               # index of the previous label (used for backtracking)
    # other attributes of the label, currently not used 
    ##=
    a1               :: Int64 
    a2               :: Int64
    a3               :: Int64
    a4               :: Int64 
    a5               :: Int64
    a6               :: Bool 
    a7               :: Bool  
    ##=#
end

function label_dominance(label_1::Label, label_2::Label, inst::LabelInstance; eps::Float64=0.000001)
    # check whether label_1 dominates label_2
    if label_1.time > label_2.time && label_2.complete == false # for current structure current condution always true
        return false 
    elseif !(label_1.setReps ⊆ label_2.setReps) 
        return false
    elseif label_1.cost > label_2.cost + eps
        return false  
    end 
    # some other checks (using inst and some attributes from a1-a7)
    # ...
    return true 
end 

@inline function single_node_bitset_1xUInt128(n::Integer)
    @assert 1 ≤ n ≤ 128
    SBitSet{1,UInt128}((UInt128(1) << (n - 1),))
end

function dummy_label()
    #Label(0.0, SBitSet{1,UInt128}(1), 0, 0, 0, false, 0,0,0)
    Label(0.0, SBitSet{1,UInt128}(1), 0, 0, 0, false, 0,0,0,0,0,0,0,0,false, false)
end 
function extend_label(label::Label, node_ext::Integer, inst::LabelInstance, idx_label::Int64)
    # extend label to node_ext 
    arr_time= label.time + inst.dur[label.node] + inst.travelTimes[label.node, node_ext] 
    start_time= max(arr_time, inst.lbTW[node_ext])                                       
    if start_time > inst.ubTW[node_ext] # no arrival possible within time-window
        return false, dummy_label()
    end
    
    load= label.load + inst.dem[node_ext]
    if load > inst.cap                    # load exceeds vehicle capacity
        return false, dummy_label() 
    end

    if node_ext in label.setReps    # node_ext already in patial route 
        return false, dummy_label()
    end

    set_repr= label.setReps ∪ single_node_bitset_1xUInt128(node_ext) 
    cost= label.cost + inst.cost[label.node, node_ext]  
    complete = node_ext == inst.numNodes ? true : false 

    #label_ext= Label(cost, set_repr, node_ext, start_time, load, complete, label.time, label.node, idx_label)
    label_ext= Label(cost, set_repr, node_ext, start_time, load, complete, label.time, label.node, idx_label, 0,0,0,0,0,false,false)
    
    return true, label_ext 
end 
            
function check_dom_and_add_label!(all_labels::Matrix{Vector{Label}}, new_label::Label, inst::LabelInstance)
    # check whether new_label is dominated by existing labels 
    node= new_label.node
    
    lb= inst.lbTW[node]
    time_label= new_label.time
    for t=lb:time_label
        labels= all_labels[t+1,node]
        if isempty(labels)
            continue 
        end 
        for i in eachindex(labels)
            exis_label= labels[i] 
            if label_dominance(exis_label, new_label, inst)
                return false     
            end 
        end
    end 

    # check whether new_label dominates existing labels 
    ub= inst.ubTW[node]
    for t=time_label:ub
        labels= all_labels[t+1,node]
        if !(isempty(labels))
            deleteat!(labels, (i for i in eachindex(labels) if label_dominance(new_label, labels[i], inst))) 
        end
    end

    # add new_label
    push!(all_labels[time_label+1, node], new_label)                   
    
    return true
end

function extend_labels_node!(node::Int64, time::Int64, all_labels::Matrix{Vector{Label}}, inst::LabelInstance, neg_cost_labels::Vector{Label}; eps::Float64=0.000001)
    for (idx,label) in enumerate(all_labels[time+1, node]) 
        for node_ext=2:inst.numNodes 
            feas_label, new_label= extend_label(label, node_ext, inst, idx)  # extend label 
            if feas_label == true                                            # label_ext feasible 
                label_add= check_dom_and_add_label!(all_labels, new_label, inst)
                (label_add == true && node_ext == inst.numNodes && new_label.cost < -eps) ? push!(neg_cost_labels, new_label) : nothing
            end 
        end 
    end
    
    return nothing
end 

function solve_pricing_problem(inst::LabelInstance, time_limit::Float64=60.0)
    # initialize matrix with labels
    all_labels= Matrix{Vector{Label}}(undef,inst.ubTW[end]+1, inst.numNodes) 
    for i=1:inst.numNodes
        for t=inst.lbTW[i]:inst.ubTW[i]
            all_labels[t+1,i]= Vector{Label}() # stores labels for node i and time t
        end 
    end 
    neg_cost_labels= Vector{Label}() # stores labels with negative reduced cost
    
    # initialize start label
    start_node= 1             
    #start_label= Label(0.0,SBitSet{1,UInt128}(start_node), start_node, inst.lbTW[start_node], inst.dem[start_node], false, 0,0,0)
    start_label= Label(0.0,SBitSet{1,UInt128}(start_node), start_node, inst.lbTW[start_node], inst.dem[start_node], false, 0,0,0,0,0,0,0,0,false,false)
    push!(all_labels[start_label.time+1,start_node], start_label)

    # extend labels
    time_start_alg= time()
    for t=inst.lbTW[start_node]:inst.ubTW[end]
        if time()- time_start_alg < time_limit && length(neg_cost_labels) < 100
            for node=1:inst.numNodes 
                if t >= inst.lbTW[node] && t <= inst.ubTW[node] 
                    extend_labels_node!(node, t, all_labels, inst, neg_cost_labels) # extend labels in all_labels[t,node]
                end
            end
        end 
    end  

    return all_labels, neg_cost_labels
end 

function label_instance()
    dem= [0, 10, 7, 13, 19, 26, 3, 5, 9, 16, 16, 12, 19, 23, 20, 8, 19, 2, 12, 17, 9, 11, 18, 29, 3, 6, 17, 16, 16, 9, 21, 27, 23, 11, 14, 8, 5, 8, 16, 31, 9, 5, 5, 7, 18, 16, 1, 27, 36, 30, 13, 0]
    dur= [0, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 0]
    lb= [0, 707, 143, 527, 678, 34, 415, 331, 404, 400, 577, 206, 228, 690, 33, 175, 272, 733, 377, 269, 581, 214, 409, 206, 704, 817, 588, 104, 114, 190, 259, 152, 660, 45, 529, 686, 42, 606, 302, 34, 360, 396, 26, 620, 233, 30, 515, 85, 773, 501, 547, 0]
    ub= [1000, 848, 282, 584, 801, 209, 514, 410, 481, 497, 632, 325, 345, 827, 243, 300, 373, 870, 434, 378, 666, 331, 494, 325, 847, 956, 667, 255, 255, 313, 354, 275, 777, 200, 614, 813, 208, 693, 405, 224, 437, 511, 172, 705, 340, 189, 628, 250, 906, 540, 642, 1000]
    xy_coord= [35 35; 41 49; 35 17; 55 45; 55 20; 15 30; 25 30; 20 50; 10 43; 55 60; 30 60; 20 65; 50 35; 30 25; 15 10; 30 5; 10 20; 5 30; 20 40; 15 60; 45 65; 45 20; 45 10; 55 5; 65 35; 65 20; 45 30; 35 40; 41 37; 64 42; 40 60; 31 52; 35 69; 53 52; 65 55; 63 65; 2 60; 20 20; 5 5; 60 12; 40 25; 42 7; 24 12; 23 3; 11 14; 6 38; 2 48; 8 56; 13 52; 6 68; 47 47; 35 35]
    travel_time=ceil.(Int64,pairwise(Euclidean(), xy_coord', dims=2))
    ass= [-56.1, 88.1, 13.45, 102.1, 106.1, 8.0, 80.1, 29.0, 110.1, 122.1, 108.1, 12.45, 22.0, 80.1, 39.0, 16.0, 48.0, 118.1, 10.55, 26.0, 120.1, 17.55, 110.1, 28.45, 116.1, 124.1, 80.1, 10.55, 7.0, 27.55, 29.0, 17.0, 124.1, 34.0, 130.1, 140.1, 28.55, 100.1, 22.55, 42.55, 22.0, 45.0, 12.55, 126.1, 10.0, 14.0, 128.1, 34.0, 112.1, 144.1, 90.1, 0.0]
    cost= travel_time .- ass 
    cap_vehicle= 500 
    num_nodes= length(dur)

    inst= LabelInstance(num_nodes, dem, dur, lb, ub, travel_time, cost, cap_vehicle)

    return inst 
end 

function main(time_limit::Float64)
    inst= label_instance()
    s= time()
    all_labels, neg_cost_labels= solve_pricing_problem(inst, time_limit)
    println(s-time())

    return nothing 
end 

main(5.0)

using Profile
using PProf
PProf.@pprof main(40.0)

it takes about 30 seconds, and still 20% of the time is spend on getindex

big_label_dom
big_label_dom906×889 48.7 KB
.

When I reduce the number of elements in the label (this only changes the amount of memory that is used to store each label, but nothing about the actual computations) and run the code


using StaticBitSets
using Distances 

struct LabelInstance 
    numNodes         :: Int64               # number of nodes
    dem              :: Vector{Int64}       # demand 
    dur              :: Vector{Int64}       # duration 
    lbTW             :: Vector{Int64}       # lower bound of time window
    ubTW             :: Vector{Int64}       # upper bound of time window 
    travelTimes      :: Matrix{Int64}       # travel times 
    cost             :: Matrix{Float64}     # cost matrix       
    cap              :: Int64               # capacity of vehicle
end 

struct Label
    cost             :: Float64             # reduced cost of the label
    setReps          :: SBitSet{1,UInt128}  # set representation partial route 
    node             :: Int64               # node of the label
    time             :: Int64               # time of the label
    load             :: Int64               # load of the label
    complete         :: Bool                # complete label
    timePrevLab      :: Int64               # time of the previous label (used for backtracking)
    nodePrevLab      :: Int64               # node of the previous label (used for backtracking)
    idxPrevLabel     :: Int64               # index of the previous label (used for backtracking)
    # other attributes of the label, currently not used 
    #=
    a1               :: Int64 
    a2               :: Int64
    a3               :: Int64
    a4               :: Int64 
    a5               :: Int64
    a6               :: Bool 
    a7               :: Bool  
    =#
end

function label_dominance(label_1::Label, label_2::Label, inst::LabelInstance; eps::Float64=0.000001)
    # check whether label_1 dominates label_2
    if label_1.time > label_2.time && label_2.complete == false # for current structure current condution always true
        return false 
    elseif !(label_1.setReps ⊆ label_2.setReps) 
        return false
    elseif label_1.cost > label_2.cost + eps
        return false  
    end 
    # some other checks (using inst and some attributes from a1-a7)
    # ...
    return true 
end 

@inline function single_node_bitset_1xUInt128(n::Integer)
    @assert 1 ≤ n ≤ 128
    SBitSet{1,UInt128}((UInt128(1) << (n - 1),))
end

function dummy_label()
    Label(0.0, SBitSet{1,UInt128}(1), 0, 0, 0, false, 0,0,0)
    #Label(0.0, SBitSet{1,UInt128}(1), 0, 0, 0, false, 0,0,0,0,0,0,0,0,false, false)
end 
function extend_label(label::Label, node_ext::Integer, inst::LabelInstance, idx_label::Int64)
    # extend label to node_ext 
    arr_time= label.time + inst.dur[label.node] + inst.travelTimes[label.node, node_ext] 
    start_time= max(arr_time, inst.lbTW[node_ext])                                       
    if start_time > inst.ubTW[node_ext] # no arrival possible within time-window
        return false, dummy_label()
    end
    
    load= label.load + inst.dem[node_ext]
    if load > inst.cap                    # load exceeds vehicle capacity
        return false, dummy_label() 
    end

    if node_ext in label.setReps    # node_ext already in patial route 
        return false, dummy_label()
    end

    set_repr= label.setReps ∪ single_node_bitset_1xUInt128(node_ext) 
    cost= label.cost + inst.cost[label.node, node_ext]  
    complete = node_ext == inst.numNodes ? true : false 

    label_ext= Label(cost, set_repr, node_ext, start_time, load, complete, label.time, label.node, idx_label)
    #label_ext= Label(cost, set_repr, node_ext, start_time, load, complete, label.time, label.node, idx_label, 0,0,0,0,0,false,false)
    
    return true, label_ext 
end 
            
function check_dom_and_add_label!(all_labels::Matrix{Vector{Label}}, new_label::Label, inst::LabelInstance)
    # check whether new_label is dominated by existing labels 
    node= new_label.node
    
    lb= inst.lbTW[node]
    time_label= new_label.time
    for t=lb:time_label
        labels= all_labels[t+1,node]
        if isempty(labels)
            continue 
        end 
        for i in eachindex(labels)
            exis_label= labels[i] 
            if label_dominance(exis_label, new_label, inst)
                return false     
            end 
        end
    end 

    # check whether new_label dominates existing labels 
    ub= inst.ubTW[node]
    for t=time_label:ub
        labels= all_labels[t+1,node]
        if !(isempty(labels))
            deleteat!(labels, (i for i in eachindex(labels) if label_dominance(new_label, labels[i], inst))) 
        end
    end

    # add new_label
    push!(all_labels[time_label+1, node], new_label)                   
    
    return true
end

function extend_labels_node!(node::Int64, time::Int64, all_labels::Matrix{Vector{Label}}, inst::LabelInstance, neg_cost_labels::Vector{Label}; eps::Float64=0.000001)
    for (idx,label) in enumerate(all_labels[time+1, node]) 
        for node_ext=2:inst.numNodes 
            feas_label, new_label= extend_label(label, node_ext, inst, idx)  # extend label 
            if feas_label == true                                            # label_ext feasible 
                label_add= check_dom_and_add_label!(all_labels, new_label, inst)
                (label_add == true && node_ext == inst.numNodes && new_label.cost < -eps) ? push!(neg_cost_labels, new_label) : nothing
            end 
        end 
    end
    
    return nothing
end 

function solve_pricing_problem(inst::LabelInstance, time_limit::Float64=60.0)
    # initialize matrix with labels
    all_labels= Matrix{Vector{Label}}(undef,inst.ubTW[end]+1, inst.numNodes) 
    for i=1:inst.numNodes
        for t=inst.lbTW[i]:inst.ubTW[i]
            all_labels[t+1,i]= Vector{Label}() # stores labels for node i and time t
        end 
    end 
    neg_cost_labels= Vector{Label}() # stores labels with negative reduced cost
    
    # initialize start label
    start_node= 1             
    start_label= Label(0.0,SBitSet{1,UInt128}(start_node), start_node, inst.lbTW[start_node], inst.dem[start_node], false, 0,0,0)
    #start_label= Label(0.0,SBitSet{1,UInt128}(start_node), start_node, inst.lbTW[start_node], inst.dem[start_node], false, 0,0,0,0,0,0,0,0,false,false)
    push!(all_labels[start_label.time+1,start_node], start_label)

    # extend labels
    time_start_alg= time()
    for t=inst.lbTW[start_node]:inst.ubTW[end]
        if time()- time_start_alg < time_limit && length(neg_cost_labels) < 100
            for node=1:inst.numNodes 
                if t >= inst.lbTW[node] && t <= inst.ubTW[node] 
                    extend_labels_node!(node, t, all_labels, inst, neg_cost_labels) # extend labels in all_labels[t,node]
                end
            end
        end 
    end  

    return all_labels, neg_cost_labels
end 

function label_instance()
    dem= [0, 10, 7, 13, 19, 26, 3, 5, 9, 16, 16, 12, 19, 23, 20, 8, 19, 2, 12, 17, 9, 11, 18, 29, 3, 6, 17, 16, 16, 9, 21, 27, 23, 11, 14, 8, 5, 8, 16, 31, 9, 5, 5, 7, 18, 16, 1, 27, 36, 30, 13, 0]
    dur= [0, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 0]
    lb= [0, 707, 143, 527, 678, 34, 415, 331, 404, 400, 577, 206, 228, 690, 33, 175, 272, 733, 377, 269, 581, 214, 409, 206, 704, 817, 588, 104, 114, 190, 259, 152, 660, 45, 529, 686, 42, 606, 302, 34, 360, 396, 26, 620, 233, 30, 515, 85, 773, 501, 547, 0]
    ub= [1000, 848, 282, 584, 801, 209, 514, 410, 481, 497, 632, 325, 345, 827, 243, 300, 373, 870, 434, 378, 666, 331, 494, 325, 847, 956, 667, 255, 255, 313, 354, 275, 777, 200, 614, 813, 208, 693, 405, 224, 437, 511, 172, 705, 340, 189, 628, 250, 906, 540, 642, 1000]
    xy_coord= [35 35; 41 49; 35 17; 55 45; 55 20; 15 30; 25 30; 20 50; 10 43; 55 60; 30 60; 20 65; 50 35; 30 25; 15 10; 30 5; 10 20; 5 30; 20 40; 15 60; 45 65; 45 20; 45 10; 55 5; 65 35; 65 20; 45 30; 35 40; 41 37; 64 42; 40 60; 31 52; 35 69; 53 52; 65 55; 63 65; 2 60; 20 20; 5 5; 60 12; 40 25; 42 7; 24 12; 23 3; 11 14; 6 38; 2 48; 8 56; 13 52; 6 68; 47 47; 35 35]
    travel_time=ceil.(Int64,pairwise(Euclidean(), xy_coord', dims=2))
    ass= [-56.1, 88.1, 13.45, 102.1, 106.1, 8.0, 80.1, 29.0, 110.1, 122.1, 108.1, 12.45, 22.0, 80.1, 39.0, 16.0, 48.0, 118.1, 10.55, 26.0, 120.1, 17.55, 110.1, 28.45, 116.1, 124.1, 80.1, 10.55, 7.0, 27.55, 29.0, 17.0, 124.1, 34.0, 130.1, 140.1, 28.55, 100.1, 22.55, 42.55, 22.0, 45.0, 12.55, 126.1, 10.0, 14.0, 128.1, 34.0, 112.1, 144.1, 90.1, 0.0]
    cost= travel_time .- ass 
    cap_vehicle= 500 
    num_nodes= length(dur)

    inst= LabelInstance(num_nodes, dem, dur, lb, ub, travel_time, cost, cap_vehicle)

    return inst 
end 

function main(time_limit::Float64)
    inst= label_instance()
    s= time()
    all_labels, neg_cost_labels= solve_pricing_problem(inst, time_limit)
    println(s-time())

    return nothing 
end 

main(5.0)

using Profile
using PProf

PProf.@pprof main(40.0)

it takes only about 20 seconds, and about 12% of the time is spend on getindex.

small_label_dom
small_label_dom1180×931 58.5 KB
.

How is it possible that when I reduce the amount of memory for each label, the code is 33% faster (this is unfortunately no option for me)? And why is 20% of the time spend on getindex?

10

I would bet that it is related to the memory bandwidth varying (increasing) for smaller problem sizes (cache effect).

1 Like