Skip to content

Commit

Permalink
Merge pull request #4 from rafaelmartinelli/expknap
Browse files Browse the repository at this point in the history
Pisinger's Expanding Core Algorithm
  • Loading branch information
rafaelmartinelli authored May 7, 2021
2 parents 63a1467 + de9a898 commit 1469c10
Show file tree
Hide file tree
Showing 9 changed files with 466 additions and 88 deletions.
58 changes: 44 additions & 14 deletions README.md
Original file line number Diff line number Diff line change
Expand Up @@ -3,29 +3,58 @@
<!-- [![Stable](https://img.shields.io/badge/docs-stable-blue.svg)](https://rafaelmartinelli.github.io/KnapsackLib.jl/stable)
[![Dev](https://img.shields.io/badge/docs-dev-blue.svg)](https://rafaelmartinelli.github.io/KnapsackLib.jl/dev) -->
[![Build Status](https://github.com/rafaelmartinelli/KnapsackLib.jl/workflows/CI/badge.svg)](https://github.com/rafaelmartinelli/KnapsackLib.jl/actions)
[![Coverage](https://codecov.io/gh/rafaelmartinelli/KnapsackLib.jl/branch/master/graph/badge.svg)](https://codecov.io/gh/rafaelmartinelli/KnapsackLib.jl)
[![Coverage](https://codecov.io/gh/rafaelmartinelli/KnapsackLib.jl/branch/main/graph/badge.svg)](https://codecov.io/gh/rafaelmartinelli/KnapsackLib.jl)
[![Project Status: Active – The project has reached a stable, usable state and is being actively developed.](https://www.repostatus.org/badges/latest/active.svg)](https://www.repostatus.org/#active)

This package solves Knapsack Problems (KP), by both Binary Programming and Dynamic Programming. First, it defines the `KnapsackData` type:
This package solves Knapsack Problems (KP) using different algorithms. First, it defines the `KnapItem` and `KnapData` types:
```julia
struct KnapsackData
capacity ::Int64 # Knapsack capacity
weights ::Vector{Int64} # Items' weights
profits ::Vector{Float64} # Items' profits
struct KnapItem
weight::Int64 # Item weight
profit::Int64 # Item profit
end

struct KnapData
capacity::Int64 # Knapsack capacity
items::Vector{KnapItem} # Items
end
```

Then, there are two available functions, which take a `KnapsackData` and return an `Array` with the selected items:
Then, there are three available functions, which take a `KnapData` and return the optimal value and an `Array` with the selected items:
```julia
function solveKnapsackDinProg(data::KnapsackData)
function solveKnapsackModel(data::KnapsackData, optimizer)
# Solves KP using a naïve dynamic programming
function solveKnapNaive(data::KnapsackData)
# Solves KP using a binary programming model
function solveKnapModel(data::KnapsackData, optimizer)
# Solves KP using Pisinger's expanding core algorithm
function solveKnapExpCore(data::KnapsackData)
```
Function `solveKnapsackModel` uses [JuMP](https://jump.dev/), and the user must pass the optimizer.
Function `solveKnapModel` uses [JuMP](https://jump.dev/), and the user must pass the optimizer.

For example, given a `KnapsackData` instance `data`:
For example, given a `KnapData` instance `data`:
```julia
selected = solveKnapsackDinProg(data)
selected = solveKnapsackModel(data, GLPK.Optimizer)
optimal, selected = solveKnapNaive(data)
optimal, selected = solveKnapModel(data, GLPK.Optimizer)
optimal, selected = solveKnapExpCore(data)
```

Benchmark results (time in seconds) for different maximum values for weights and profits, number of items and algorithms. Average times for 10 runs and using `@timed` (Model with GLPK).

```
------------------------------------------------------------------------------
MaxV\Items 10 50 100 500 1000 Algorithm
------------------------------------------------------------------------------
0.0000066 0.0000842 0.0003162 0.0115152 0.0585165 Naïve
10 0.0005138 0.0006685 0.0011759 0.0063126 0.0173382 Model
0.0000141 0.0000306 0.0000515 0.0002397 0.0003492 Exp Core
------------------------------------------------------------------------------
0.0000136 0.0003581 0.0014877 0.0632500 0.3000128 Naïve
50 0.0004836 0.0007185 0.0018582 0.0129172 0.0106744 Model
0.0000135 0.0000391 0.0000665 0.0002251 0.0003272 Exp Core
------------------------------------------------------------------------------
0.0000259 0.0004414 0.0023123 0.1315063 0.7479475 Naïve
100 0.0005052 0.0010954 0.0015549 0.0113286 0.0343576 Model
0.0000180 0.0000486 0.0000714 0.0002987 0.0004923 Exp Core
------------------------------------------------------------------------------
```

To install:
Expand All @@ -34,4 +63,5 @@ To install:
```

Related links:
- ...
- [David Pisinger's optimization codes](http://hjemmesider.diku.dk/~pisinger/codes.html)
- [Bin Packing and Cutting Stock Lib](https://github.com/rafaelmartinelli/BPPLib.jl)
18 changes: 18 additions & 0 deletions src/Data.jl
Original file line number Diff line number Diff line change
@@ -0,0 +1,18 @@
struct KnapItem
weight::Int64
profit::Int64
end

struct KnapData
capacity::Int64
items::Vector{KnapItem}

KnapData(capacity::Int64, items::Vector{KnapItem}) = new(capacity, items)
KnapData(capacity::Int64, weights::Vector{Int64}, profits::Vector{Int64}) = new(capacity, [ KnapItem(weights[i], profits[i]) for i in 1:length(weights) ])
end

function Base.show(io::IO, data::KnapData)
print(io, "Knapsack Data")
print(io, " ($(length(data.items)) items,")
print(io, " capacity = $(data.capacity))")
end
263 changes: 263 additions & 0 deletions src/ExpandingCore.jl
Original file line number Diff line number Diff line change
@@ -0,0 +1,263 @@
mutable struct Item
index::Int64
weight::Int64
profit::Int64
solution::Bool
end

struct Interval
first::Int64
last::Int64
weight::Int64
end

mutable struct GlobalData
items::Vector{Item}
capacity::Int64

br::Int64
wsb::Int64
psb::Int64
z::Int64
fsort::Int64
lsort::Int64

stack1::Vector{Interval}
stack2::Vector{Interval}

estack::Vector{Int64}

GlobalData() = new([], 0, 0, 0, 0, 0, 0, 0, [], [], [])
end

function solveKnapExpCore(data::KnapData)
n = length(data.items)
if n == 0 return 0, Int64[] end

gd = GlobalData()
gd.items = [ Item(i, data.items[i].weight, data.items[i].profit, false) for i in 1:n ]
gd.capacity = data.capacity

interval = partsort(gd, 1, n, 0)
gd.fsort = interval.first
gd.lsort = interval.last

gd.z = heuristic(gd, 1, n)
heur = gd.z + gd.psb

elebranch(gd, 0, gd.wsb - gd.capacity, gd.br - 1, gd.br)

return gd.z + gd.psb, definesolution(gd)
end

function det(a::Matrix{Int64})
return a[1, 1] * a[2, 2] - a[1, 2] * a[2, 1]
end

function partsort(gd::GlobalData, f::Int64, l::Int64, ws::Int64)
d = l - f + 1
if d > 1
m = f + d ÷ 2

if det([ gd.items[f].profit gd.items[f].weight; gd.items[m].profit gd.items[m].weight ]) < 0
gd.items[f], gd.items[m] = gd.items[m], gd.items[f]
end

if d > 2
if det([ gd.items[m].profit gd.items[m].weight; gd.items[l].profit gd.items[l].weight ]) < 0
gd.items[m], gd.items[l] = gd.items[l], gd.items[m]

if det([ gd.items[f].profit gd.items[f].weight; gd.items[m].profit gd.items[m].weight ]) < 0
gd.items[f], gd.items[m] = gd.items[m], gd.items[f]
end
end
end
end

if d <= 3
return Interval(f, l, ws)
else
mp = gd.items[m].profit
mw = gd.items[m].weight
i = f
j = l
wi = ws
while true
while true
wi += gd.items[i].weight
i += 1
if det([ gd.items[i].profit gd.items[i].weight; mp mw ]) <= 0 break end
end
while true
j -= 1
if det([ gd.items[j].profit gd.items[j].weight; mp mw ]) >= 0 break end
end
if i > j break end
gd.items[i], gd.items[j] = gd.items[j], gd.items[i]
end

if wi > gd.capacity
push!(gd.stack2, Interval(i, l, wi))
return partsort(gd, f, i - 1, ws)
else
push!(gd.stack1, Interval(f, i - 1, ws))
return partsort(gd, i, l, wi)
end
end
end

function heuristic(gd::GlobalData, f::Int64, l::Int64)
ps = 0
ws = gd.capacity

i = f
while gd.items[i].weight <= ws
ws -= gd.items[i].weight
ps += gd.items[i].profit
gd.items[i].solution = true
i += 1
end

gd.br = i
gd.wsb = gd.capacity - ws
gd.psb = ps

dz = (gd.capacity - gd.wsb) * gd.items[gd.br].profit ÷ gd.items[gd.br].weight

empty!(gd.estack)
gd.z = 0
if gd.z == dz return gd.z end

r = gd.capacity - gd.wsb
for i in gd.br:l
if (gd.items[i].weight <= r) && (gd.items[i].profit > gd.z)
empty!(gd.estack)
push!(gd.estack, gd.items[i].index)
gd.z = gd.items[i].profit
if gd.z == dz return gd.z end
end
end

r = gd.wsb + gd.items[gd.br].weight - gd.capacity
pb = gd.items[gd.br].profit
for i in gd.br - 1:-1:f
if (gd.items[i].weight >= r) && (pb - gd.items[i].profit > gd.z)
empty!(gd.estack)
push!(gd.estack, gd.items[gd.br].index)
push!(gd.estack, gd.items[i].index)
gd.z = pb - gd.items[i].profit
if gd.z == dz return gd.z end
end
end

return gd.z
end

function reduce(gd::GlobalData, first::Int64, last::Int64)
pb = gd.items[gd.br].profit
wb = gd.items[gd.br].weight

q = det([ gd.z + 1 gd.capacity - gd.wsb; pb wb ])
i = first
j = last

if i <= gd.br
k = gd.fsort - 1
while i <= j
if det([ -gd.items[j].profit -gd.items[j].weight; pb wb ]) < q
gd.items[i], gd.items[j] = gd.items[j], gd.items[i]
i += 1
else
gd.items[j], gd.items[k] = gd.items[k], gd.items[j]
j -= 1
k -= 1
end
end
if k == gd.fsort - 1
gd.items[first], gd.items[k] = gd.items[k], gd.items[first]
k -= 1
end
last = gd.fsort - 1
first = k + 1
else
k = gd.lsort + 1
while i <= j
if det([ gd.items[i].profit gd.items[i].weight; pb wb ]) < q
gd.items[i], gd.items[j] = gd.items[j], gd.items[i]
j -= 1
else
gd.items[i], gd.items[k] = gd.items[k], gd.items[i]
i += 1
k += 1
end
end
if k == gd.lsort + 1
gd.items[last], gd.items[k] = gd.items[k], gd.items[last]
k += 1
end
first = gd.lsort + 1
last = k - 1
end
return first, last
end

function sorti(gd::GlobalData, stack::Vector{Interval})
if isempty(stack) return false end

interval = pop!(stack)
first, last = reduce(gd, interval.first, interval.last)
interval = partsort(gd, first, last, interval.weight)
if interval.first < gd.fsort gd.fsort = interval.first end
if interval.last > gd.lsort gd.lsort = interval.last end

return true
end

function elebranch(gd::GlobalData, ps::Int64, ws::Int64, s::Int64, t::Int64)
improved = false

if ws <= 0
if ps > gd.z
improved = true
gd.z = ps
empty!(gd.estack)
end

while true
if t > gd.lsort && !sorti(gd, gd.stack2) break end
if det([ ps - (gd.z + 1) ws; gd.items[t].profit gd.items[t].weight ]) < 0 break end
if elebranch(gd, ps + gd.items[t].profit, ws + gd.items[t].weight, s, t + 1)
improved = true
push!(gd.estack, gd.items[t].index)
end
t += 1
end
else
while true
if s < gd.fsort && !sorti(gd, gd.stack1) break end
if det([ ps - (gd.z + 1) ws; gd.items[s].profit gd.items[s].weight ]) < 0 break end
if elebranch(gd, ps - gd.items[s].profit, ws - gd.items[s].weight, s - 1, t)
improved = true
push!(gd.estack, gd.items[s].index)
end
s -= 1
end
end

return improved
end

function definesolution(gd::GlobalData)
sort!(gd.items, by = x -> x.index)
for index in gd.estack
gd.items[index].solution = !gd.items[index].solution
end

solution = Int64[]
for item in gd.items
if item.solution
push!(solution, item.index)
end
end
return solution
end
11 changes: 0 additions & 11 deletions src/KnapsackData.jl

This file was deleted.

Loading

0 comments on commit 1469c10

Please sign in to comment.