Tutorials
Tutorial 1: A basic JuDGE model
Problem description
For our tutorial, we will consider the following optimization problem: Our goal is minimize the cost of a stochastic sequence of knapsack problems. We represent the stochastic process with a discrete scenario tree. At each node in the scenario tree, we solve a knapsack problem. However at any point in the tree, we have the ability to expand the capacity of our knapsack at a certain cost. Once the capacity of our bag has been expanded, we are able use the extra volume for future knapsack problems from this node of the tree forward.
In this optimization problem, we are trading off against the cost of expanding our knapsack, versus the ability to fit more into our knapsack. Deciding when to perform the knapsack expansion is the difficult part of this optimization problem.
Solving our problem using JuDGE
Let us first load the packages that we need to create and solve some simple JuDGE
models.
using JuDGE, JuMP, GLPK
The lifecycle of a JuDGEModel
is the following:
- The definition of a
Tree
; - defining the subproblems of the
JuDGEModel
; - building the
JuDGEModel
; - solving the
JuDGEModel
.
The user's job is to complete Steps 1 and 2, while JuDGE will automatically perform Steps 3 and 4.
A Tree
can be built in many different ways. A Tree
simply consists of the root node of the tree, and a list of all the nodes in the tree. This is defined as a nested set of subtrees, with the final nodes being Leaf
nodes. Each subtree simply defines its parent and children, and there are functions that facilitate the probability of arriving at the node, and the data that corresponds to the node, can be referenced through dictionaries.
For now, we will build a tree of depth 2, where each node has 2 children with uniform probabilities using narytree
:
mytree = narytree(2,2)
Subtree rooted at node 1 containing 7 nodes
mytree
is a tree which contains 7 nodes, with depth 2, and degree 2. (A depth of 0, gives only a single leaf node.) We can visualise the tree using
JuDGE.print_tree(mytree)
--1
--11
--111
--112
--12
--121
--122
Problem data
Let us now associate our tree with the problem data.
For our instance of the problem, we will use the following data: for each node, we will have a knapsack problem with five items to choose from, each with different rewards, and different volumes. The structure of this data is arbitrary; JuDGE just needs to be able to access the relevant data, based on the node being processed (dictionaries or functions are recommended).
invest_cost = Dict( zip( collect(mytree,order=:breadth), [15, 8, 8, 4, 4, 4, 4]) )
item_volume = Dict( zip( collect(mytree,order=:breadth), [ [4, 3, 3, 1, 2],
[5, 3, 4, 2, 1],
[5, 4, 2, 7, 2],
[5, 4, 1, 8, 2],
[3, 1, 5, 6, 3],
[2, 5, 8, 4, 6],
[7, 5, 4, 2, 3] ]) )
item_reward = Dict( zip( collect(mytree,order=:breadth), [ [32, 9, 9, 4, 8],
[30, 12, 40, 10, 9],
[20, 28, 12, 42, 12],
[40, 28, 9, 24, 10],
[15, 7, 20, 48, 12],
[10, 30, 54, 32, 30],
[32, 25, 24, 14, 24] ]) )
Dict{AbstractTree,Array{Int64,1}} with 7 entries:
Subtree rooted at node 11 containing 3 nodes => [30, 12, 40, 10, 9]
Leaf node 121 => [10, 30, 54, 32, 30]
Leaf node 111 => [40, 28, 9, 24, 10]
Leaf node 122 => [32, 25, 24, 14, 24]
Subtree rooted at node 12 containing 3 nodes => [20, 28, 12, 42, 12]
Subtree rooted at node 1 containing 7 nodes => [32, 9, 9, 4, 8]
Leaf node 112 => [15, 7, 20, 48, 12]
We can print the tree again and display a parameter for each node:
JuDGE.print_tree(mytree,item_reward)
--1 ([32, 9, 9, 4, 8])
--11 ([30, 12, 40, 10, 9])
--111 ([40, 28, 9, 24, 10])
--112 ([15, 7, 20, 48, 12])
--12 ([20, 28, 12, 42, 12])
--121 ([10, 30, 54, 32, 30])
--122 ([32, 25, 24, 14, 24])
We also define some other parameters that apply to all the nodes.
num_items = 5
num_invest = 6
initial_volume = 6
invest_volume = [2,2,2,3,3,3]
6-element Array{Int64,1}:
2
2
2
3
3
3
Subproblems
We need to now define the subproblems. These are JuMP models with some JuDGE- specific features. For our knapsack problem:
JuDGE_SP_Solver = optimizer_with_attributes(GLPK.Optimizer, "msg_lev" => 0, "mip_gap" => 0.0)
function sub_problems(node)
sp = Model(JuDGE_SP_Solver)
@expansion(sp, invest[1:num_invest], Bin)
@capitalcosts(sp, sum(invest[i]*invest_volume[i] for i=1:num_invest)*invest_cost[node])
@variable(sp, y[1:num_items], Bin)
@constraint(sp, BagExtension, sum(y[i]*item_volume[node][i] for i in 1:num_items) <=
initial_volume + sum(invest_volume[i] * invest[i] for i in 1:num_invest))
@objective(sp, Min, sum(-item_reward[node][i] * y[i] for i in 1:num_items))
return sp
end
sub_problems (generic function with 1 method)
The two elements of this that make it a JuDGE subproblem are:
@expansion(model, ...)
This defines the expansion variables, and supports standard JuMP vectorized variable declaration. In this case these are binary.
@capitalcosts
This declares an expression for the costs of investment; this must be linear (an AffExpr
).
The overall optimization problem at each node problem is a classical knapsack problem. We have specified the initial volume of the knapsack is initial_volume
, and the each expansion increases the volume of the knapsack by the corresponding value in invest_volume
.
Solving the JuDGE Model
Now with our tree built and the problem data referenced, we can initialize the JuDGEModel
based on our tree, subproblems, and solver.
JuDGE_MP_Solver = optimizer_with_attributes((method=GLPK.INTERIOR) -> GLPK.Optimizer(),
"msg_lev" => 0, "mip_gap" => 0.0)
judy = JuDGEModel(mytree, ConditionallyUniformProbabilities, sub_problems, JuDGE_MP_Solver)
JuDGE Model with:
Tree: Subtree rooted at node 1 containing 7 nodes
Expansion variables: invest
ConditionallyUniformProbabilities
simply applies a uniform conditional probability distribution for child nodes. Either a function or a dictionary, which maps nodes to absolute probabilities, can be used here.
At this point, we have constructed a valid JuDGEModel
.
There are a number of optional stopping criteria that can be set: abstol
, reltol
, rlx_abstol
, rlx_reltol
, time_limit
, max_iter
.
These are grouped within a Termination
struct
that can be defined with as many termination conditions as required (if any condition is met, the solve is stopped). There are defaults for all stopping criteria, so it is not necessary to provide the Termination
object to the JuDGE.solve
function.
We can now solve our model by making a call to JuDGE.solve
:
JuDGE.solve(judy,termination=Termination(rlx_abstol=10^-7),verbose=1)
[ Info: Solving JuDGE model for tree: Subtree rooted at node 1 containing 7 nodes
Relaxed ObjVal | Upper Bound Lower Bound | Absolute Diff Relative Diff | Fractional | Time Iter
2.250000e+02 | 2.250000e+02 -Inf | Inf NaN | 0 | 11.829 1
-1.550000e+02 | -1.550000e+02 -2.410000e+02 | 8.600000e+01 3.568465e-01 | 0 | 12.208 2
-1.625000e+02 | -1.625000e+02 -2.410000e+02 | 7.850000e+01 3.257261e-01 | 0 | 12.211 3
-1.630000e+02 | -1.630000e+02 -2.410000e+02 | 7.800000e+01 3.236515e-01 | 0 | 12.212 4
-1.635000e+02 | -1.630000e+02 -1.965000e+02 | 3.300000e+01 1.679389e-01 | 7 | 12.214 5
-1.640000e+02 | -1.640000e+02 -1.812500e+02 | 1.725000e+01 9.517241e-02 | 0 | 12.217 6
-1.640000e+02 | -1.640000e+02 -1.762500e+02 | 1.225000e+01 6.950355e-02 | 0 | 12.219 7
-1.640000e+02 | -1.640000e+02 -1.680000e+02 | 4.000000e+00 2.380952e-02 | 0 | 12.221 8
-1.640000e+02 | -1.640000e+02 -1.650000e+02 | 1.000000e+00 6.060606e-03 | 0 | 12.223 9
-1.640000e+02 | -1.640000e+02 -1.650000e+02 | 1.000000e+00 6.060606e-03 | 0 | 12.224 10
-1.640000e+02 | -1.640000e+02 -1.640000e+02 | 0.000000e+00 0.000000e+00 | 0 | 12.226 11
Convergence criteria met.
Currently, we recommend using JuDGE with Gurobi as the subproblem and master problem solvers. Any solvers can be specified, but the master problem must return duals, and an interior point method is recommended for reliable convergence. The subproblems can be solved with any method. (If you do not solve the subproblems to a zero bound-gap, the upper and lower bounds will not fully converge; it is therefore also necessary to set appropriate convergence tolerances in the master problem.)
We can view the optimal solution to our problem by calling
println("Objective: "*string(JuDGE.get_objval(judy))
JuDGE.print_expansions(judy)
Finally, if we want to recover the optimal solutions for the nodes, we must fix the investments and resolve each subproblem, after which we can write the solution to a CSV file.
println("Re-solved Objective: " * string(resolve_subproblems(judy)))
Re-solved Objective: -164.0
JuDGE.write_solution_to_file(judy,joinpath(@__DIR__,"knapsack_solution.csv"))
Tutorial 2: Formatting output
Since the expansion variables are all binary, the print_expansions
function doesn't directly convey how much capacity is built. In order to output more information about the capacity it's possible to write a custom function that is supplied as an optional argument format
. First, let us simply provide the capacity associated with each decision.
function format_output(s::Symbol,value)
if s==:invest
output=Dict{Int,Float64}()
for i in 1:num_invest
output[i]=invest_volume[i]*value[i]
end
return output
end
return nothing
end
JuDGE.print_expansions(judy, format=format_output)
JuDGE Expansions
Node 11: "invest[4]" 3.0
Node 11: "invest[5]" 3.0
Node 112: "invest[6]" 3.0
Node 12: "invest[4]" 3.0
Node 12: "invest[5]" 3.0
Node 12: "invest[6]" 3.0
Node 121: "invest[2]" 2.0
Node 121: "invest[3]" 2.0
Node 122: "invest[2]" 2.0
Node 122: "invest[3]" 2.0
Node 122: "invest[1]" 2.0
We can also aggregate the capacities of all the expansions.
function format_output(s::Symbol,value)
if s==:invest
return sum(invest_volume[i]*value[i] for i in 1:num_invest)
end
return nothing
end
JuDGE.print_expansions(judy, format=format_output)
JuDGE Expansions
Node 11: "invest" 6.0
Node 112: "invest" 3.0
Node 12: "invest" 9.0
Node 121: "invest" 4.0
Node 122: "invest" 6.0
Tutorial 3: Ongoing costs
Depending on the capacity planning application that JuDGE is being applied to there may be ongoing upkeep / maintenance costs for the expansions. This is modelled within JuDGE by using the @ongoingcosts
macro to specify the cost of the expansions being available at each node in the scenario tree. (That is, the corresponding expansion variable has been set to 1 in the master problem.) If there are costs that can be avoided, by not making use of capacity that is granted, then this should be modelled within the @objective
.
function sub_problems_ongoing(node)
sp = Model(JuDGE_SP_Solver)
@expansion(sp, invest[1:num_invest], Bin)
@capitalcosts(sp, sum(invest[i]*invest_volume[i] for i=1:num_invest)*invest_cost[node])
@ongoingcosts(sp, sum(invest[i]*invest_volume[i] for i=1:num_invest)*2)
@variable(sp, y[1:num_items], Bin)
@constraint(sp, BagExtension, sum(y[i]*item_volume[node][i] for i in 1:num_items) <=
initial_volume + sum(invest_volume[i] * invest[i] for i in 1:num_invest))
@objective(sp, Min, sum(-item_reward[node][i] * y[i] for i in 1:num_items))
return sp
end
judy = JuDGEModel(mytree, ConditionallyUniformProbabilities, sub_problems_ongoing,
JuDGE_MP_Solver)
JuDGE.solve(judy,verbose=1)
println("Objective: "*string(judy.bounds.UB))
JuDGE.print_expansions(judy, format=format_output)
[ Info: Establishing JuDGE model for tree: Subtree rooted at node 1 containing 7 nodes
[ Info: Checking sub-problem format
[ Info: Sub-problem format is valid
[ Info: Building master problem
[ Info: Master problem built
[ Info: Solving JuDGE model for tree: Subtree rooted at node 1 containing 7 nodes
Relaxed ObjVal | Upper Bound Lower Bound | Absolute Diff Relative Diff | Fractional | Time Iter
3.150000e+02 | 3.150000e+02 -Inf | Inf NaN | 0 | 0.003 1
-1.017500e+02 | -1.017500e+02 -2.397500e+02 | 1.380000e+02 5.755996e-01 | 0 | 0.005 2
-1.327500e+02 | -1.327500e+02 -2.397500e+02 | 1.070000e+02 4.462982e-01 | 0 | 0.007 3
-1.327500e+02 | -1.327500e+02 -2.117500e+02 | 7.900000e+01 3.730815e-01 | 0 | 0.009 4
-1.327500e+02 | -1.327500e+02 -1.935000e+02 | 6.075000e+01 3.139535e-01 | 0 | 0.011 5
-1.327500e+02 | -1.327500e+02 -1.833750e+02 | 5.062500e+01 2.760736e-01 | 0 | 0.013 6
-1.342500e+02 | -1.327500e+02 -1.762500e+02 | 4.200000e+01 2.382979e-01 | 5 | 0.016 7
-1.342500e+02 | -1.327500e+02 -1.517500e+02 | 1.750000e+01 1.153213e-01 | 5 | 0.018 8
-1.357500e+02 | -1.357500e+02 -1.455000e+02 | 9.750000e+00 6.701031e-02 | 0 | 0.019 9
-1.357500e+02 | -1.357500e+02 -1.427500e+02 | 7.000000e+00 4.903678e-02 | 0 | 0.021 10
-1.357500e+02 | -1.357500e+02 -1.392500e+02 | 3.500000e+00 2.513465e-02 | 0 | 0.023 11
-1.357500e+02 | -1.357500e+02 -1.373750e+02 | 1.625000e+00 1.182894e-02 | 0 | 0.025 12
-1.357500e+02 | -1.357500e+02 -1.357500e+02 | 0.000000e+00 0.000000e+00 | 0 | 0.027 13
Convergence criteria met.
Objective: -135.75
JuDGE Expansions
Node 11: "invest" 4.0
Node 121: "invest" 6.0
Node 122: "invest" 3.0
Tutorial 4: Deterministic equivalent
JuDGE can use the tree, and subproblems to automatically construct the deterministic equivalent of the stochastic capacity expansion problem. This is created by defining a DetEqModel
with the same arguments as a JuDGEModel
, as follows:
JuDGE_DE_Solver = optimizer_with_attributes(GLPK.Optimizer, "msg_lev" => 2, "mip_gap" => 0)
deteq = DetEqModel(mytree, ConditionallyUniformProbabilities, sub_problems, JuDGE_DE_Solver)
JuDGE.solve(deteq)
Establishing deterministic equivalent model for tree: Subtree rooted at node 1 containing 7 nodes
Checking sub-problem format...Passed
Building deterministic equivalent problem...Complete
Solving deterministic equivalent formulation...* 0: obj = 0.000000000e+00 inf = 0.000e+00 (4)
* 109: obj = -1.645000000e+02 inf = 2.220e-16 (0)
+ 109: mip = not found yet >= -inf (1; 0)
+ 144: >>>>> -1.620000000e+02 >= -1.640000000e+02 1.2% (17; 0)
+ 148: >>>>> -1.640000000e+02 >= -1.640000000e+02 0.0% (11; 11)
+ 148: mip = -1.640000000e+02 >= tree is empty 0.0% (0; 33)
Solved.
The solution can be printed to the REPL or a CSV in the same way as a JUDGEModel
.
println("Deterministic Equivalent Objective: " * string(objective_value(deteq.problem)))
JuDGE.print_expansions(deteq, format=format_output)
Deterministic Equivalent Objective: -164.0
Deterministic Equivalent Expansions
Node 11: "invest" 6.0
Node 112: "invest" 3.0
Node 12: "invest" 9.0
Node 121: "invest" 4.0
Node 122: "invest" 6.0
JuDGE.write_solution_to_file(judy,joinpath(@__DIR__,"knapsack_solution.csv"))
Tutorial 5: Lag and duration
Depending on the particular expansion problem being modelled, there may be some delay (lag) between when the expansion decision is made, and when the capacity becomes available. This can be modelled in JuDGE be specifying a lag when defining the expansion variables in the subproblems, the @expansion
macro allows additional named arguments lag
and duration
. For a lag of 1 you would define the expansion variables as follows:
@expansion(sp, invest[1:num_invest], Bin, lag=1)
The duration of an expansion is set to ∞ by default. However, if an expansion is temporary or otherwise has a limited lifespan, we can set the duration
argument when defining the expansion variable. For example if the lag is 0 and the duration is 2, we can define it as follows:
@expansion(sp, invest[1:num_invest], Bin, duration=2)
We will now redefine our subproblems and re-solve our model. (Note that the invest_cost
has been divided by 2.)
function sub_problems_lag(node)
sp = Model(JuDGE_SP_Solver)
@expansion(sp, invest[1:num_invest], Bin, lag=1)
@capitalcosts(sp, sum(invest[i]*invest_volume[i] for i=1:num_invest)*invest_cost[node]/2)
@variable(sp, y[1:num_items], Bin)
@constraint(sp, BagExtension, sum(y[i]*item_volume[node][i] for i in 1:num_items) <=
initial_volume + sum(invest_volume[i] * invest[i] for i in 1:num_invest))
@objective(sp, Min, sum(-item_reward[node][i] * y[i] for i in 1:num_items))
return sp
end
judy = JuDGEModel(mytree, ConditionallyUniformProbabilities, sub_problems_lag, JuDGE_MP_Solver)
JuDGE.solve(judy,verbose=1)
println("Objective: "*string(judy.bounds.UB))
println("Lower Bound: "*string(judy.bounds.LB))
JuDGE.print_expansions(judy, format=format_output)
[ Info: Establishing JuDGE model for tree: Subtree rooted at node 1 containing 7 nodes
[ Info: Checking sub-problem format
[ Info: Sub-problem format is valid
[ Info: Building master problem
[ Info: Master problem built
[ Info: Solving JuDGE model for tree: Subtree rooted at node 1 containing 7 nodes
Relaxed ObjVal | Upper Bound Lower Bound | Absolute Diff Relative Diff | Fractional | Time Iter
Inf | Inf -Inf | Inf NaN | 0 | 0.002 1
0.000000e+00 | 0.000000e+00 -Inf | Inf NaN | 0 | 0.004 2
-8.175000e+01 | -8.175000e+01 -3.687500e+02 | 2.870000e+02 7.783051e-01 | 0 | 0.006 3
-1.072500e+02 | -1.072500e+02 -2.947500e+02 | 1.875000e+02 6.361323e-01 | 0 | 0.008 4
-1.222500e+02 | -1.222500e+02 -2.315000e+02 | 1.092500e+02 4.719222e-01 | 0 | 0.010 5
-1.222500e+02 | -1.222500e+02 -2.257500e+02 | 1.035000e+02 4.584718e-01 | 0 | 0.012 6
-1.224167e+02 | -1.222500e+02 -2.248750e+02 | 1.024583e+02 4.556235e-01 | 3 | 0.014 7
-1.235833e+02 | -1.222500e+02 -2.093333e+02 | 8.575000e+01 4.096338e-01 | 3 | 0.017 8
-1.635833e+02 | -1.222500e+02 -1.782778e+02 | 1.469444e+01 8.242443e-02 | 3 | 0.019 9
-1.644167e+02 | -1.222500e+02 -1.782778e+02 | 1.386111e+01 7.775008e-02 | 9 | 0.022 10
-1.655833e+02 | -1.222500e+02 -1.782778e+02 | 1.269444e+01 7.120598e-02 | 10 | 0.024 11
-1.668333e+02 | -1.222500e+02 -1.754167e+02 | 8.583333e+00 4.893112e-02 | 10 | 0.027 12
-1.668333e+02 | -1.222500e+02 -1.749444e+02 | 8.111111e+00 4.636393e-02 | 10 | 0.029 13
-1.668333e+02 | -1.222500e+02 -1.749444e+02 | 8.111111e+00 4.636393e-02 | 10 | 0.032 14
-1.668333e+02 | -1.222500e+02 -1.708333e+02 | 4.000000e+00 2.341463e-02 | 10 | 0.034 15
-1.668333e+02 | -1.222500e+02 -1.688333e+02 | 2.000000e+00 1.184600e-02 | 10 | 0.037 16
-1.668333e+02 | -1.222500e+02 -1.668333e+02 | 0.000000e+00 0.000000e+00 | 10 | 0.039 17
-1.668333e+02 | -1.640000e+02 -1.668333e+02 | 2.833333e+00 1.698302e-02 | 0 | 0.070 18*
Convergence criteria met.
Objective: -164.0
Lower Bound: -166.83333333333337
JuDGE Expansions
Node 1: "invest" 9.0
Node 11: "invest" 3.0
Node 12: "invest" 6.0
We see that although this solution has ostensibly converged, the best integer is greater than the lower bound. The * in the final line of the output means that the integer solution has been found using a MIP solve for the generated columns. In order to find a provably optimal solution we must use branch-and-price.
Tutorial 6: Branch and Price
JuDGE implements a branch-and-price algorithm for problems which are not naturally integer. It can be run with default settings as follows:
judy = JuDGEModel(mytree, ConditionallyUniformProbabilities, sub_problems_lag, JuDGE_MP_Solver)
best = JuDGE.branch_and_price(judy,search=:lowestLB,branch_method=JuDGE.variable_branch,verbose=1)
JuDGE.print_expansions(best, format=format_output)
[ Info: Establishing JuDGE model for tree: Subtree rooted at node 1 containing 7 nodes
[ Info: Checking sub-problem format
[ Info: Sub-problem format is valid
[ Info: Building master problem
[ Info: Master problem built
Model 1 of 1. UB: Inf, LB:-Inf, Time: 0.0s
[ Info: Solving JuDGE model for tree: Subtree rooted at node 1 containing 7 nodes
Relaxed ObjVal | Upper Bound Lower Bound | Absolute Diff Relative Diff | Fractional | Time Iter
Inf | Inf -Inf | Inf NaN | 0 | 0.003 1
0.000000e+00 | 0.000000e+00 -Inf | Inf NaN | 0 | 0.004 2
-8.175000e+01 | -8.175000e+01 -3.687500e+02 | 2.870000e+02 7.783051e-01 | 0 | 0.006 3
-1.072500e+02 | -1.072500e+02 -2.947500e+02 | 1.875000e+02 6.361323e-01 | 0 | 0.008 4
-1.222500e+02 | -1.222500e+02 -2.315000e+02 | 1.092500e+02 4.719222e-01 | 0 | 0.010 5
-1.222500e+02 | -1.222500e+02 -2.257500e+02 | 1.035000e+02 4.584718e-01 | 0 | 0.013 6
-1.224167e+02 | -1.222500e+02 -2.248750e+02 | 1.024583e+02 4.556235e-01 | 3 | 0.015 7
-1.235833e+02 | -1.222500e+02 -2.093333e+02 | 8.575000e+01 4.096338e-01 | 3 | 0.017 8
-1.635833e+02 | -1.222500e+02 -1.782778e+02 | 1.469444e+01 8.242443e-02 | 3 | 0.034 9
-1.644167e+02 | -1.222500e+02 -1.782778e+02 | 1.386111e+01 7.775008e-02 | 9 | 0.036 10
-1.655833e+02 | -1.222500e+02 -1.782778e+02 | 1.269444e+01 7.120598e-02 | 10 | 0.039 11
-1.668333e+02 | -1.222500e+02 -1.754167e+02 | 8.583333e+00 4.893112e-02 | 10 | 0.042 12
-1.668333e+02 | -1.222500e+02 -1.749444e+02 | 8.111111e+00 4.636393e-02 | 10 | 0.044 13
-1.668333e+02 | -1.222500e+02 -1.749444e+02 | 8.111111e+00 4.636393e-02 | 10 | 0.047 14
-1.668333e+02 | -1.222500e+02 -1.708333e+02 | 4.000000e+00 2.341463e-02 | 10 | 0.049 15
-1.668333e+02 | -1.222500e+02 -1.688333e+02 | 2.000000e+00 1.184600e-02 | 10 | 0.052 16
-1.668333e+02 | -1.222500e+02 -1.668333e+02 | 0.000000e+00 0.000000e+00 | 10 | 0.054 17
-1.668333e+02 | -1.640000e+02 -1.668333e+02 | 2.833333e+00 1.698302e-02 | 0 | 0.056 18*
Convergence criteria met.
Attempting to branch.
Adding 2 new nodes to B&P tree.
Model 2 of 3. UB: -164.0, LB:-166.83333333333337, Time: 0.32s
[ Info: Solving JuDGE model for tree: Subtree rooted at node 1 containing 7 nodes
Relaxed ObjVal | Upper Bound Lower Bound | Absolute Diff Relative Diff | Fractional | Time Iter
-1.650000e+02 | Inf -1.668333e+02 | 1.833333e+00 1.098901e-02 | 5 | 0.578 1
-1.650500e+02 | Inf -1.668333e+02 | 1.783333e+00 1.068931e-02 | 8 | 0.581 2
-1.655000e+02 | Inf -1.668000e+02 | 1.300000e+00 7.793765e-03 | 6 | 0.584 3
-1.655000e+02 | Inf -1.657500e+02 | 2.500000e-01 1.508296e-03 | 6 | 0.588 4
-1.655000e+02 | Inf -1.655000e+02 | 0.000000e+00 0.000000e+00 | 6 | 0.591 5
-1.655000e+02 | -1.640000e+02 -1.655000e+02 | 1.500000e+00 9.063444e-03 | 0 | 0.592 6*
Convergence criteria met.
Attempting to branch.
Adding 2 new nodes to B&P tree.
Model 3 of 5. UB: -164.0, LB:-166.83333333333337, Time: 0.913s
[ Info: Solving JuDGE model for tree: Subtree rooted at node 1 containing 7 nodes
Relaxed ObjVal | Upper Bound Lower Bound | Absolute Diff Relative Diff | Fractional | Time Iter
-1.649583e+02 | Inf -1.668333e+02 | 1.875000e+00 1.123876e-02 | 12 | 0.225 1
-1.649583e+02 | Inf -1.668333e+02 | 1.875000e+00 1.123876e-02 | 12 | 0.227 2
-1.655000e+02 | Inf -1.668333e+02 | 1.333333e+00 7.992008e-03 | 10 | 0.230 3
-1.655000e+02 | Inf -1.668333e+02 | 1.333333e+00 7.992008e-03 | 10 | 0.233 4
-1.655000e+02 | Inf -1.655000e+02 | 0.000000e+00 0.000000e+00 | 10 | 0.235 5
-1.655000e+02 | -1.635000e+02 -1.655000e+02 | 2.000000e+00 1.208459e-02 | 0 | 0.237 6*
Convergence criteria met.
Attempting to branch.
Adding 2 new nodes to B&P tree.
Model 4 of 7. UB: -164.0, LB:-165.5, Time: 1.151s
[ Info: Solving JuDGE model for tree: Subtree rooted at node 1 containing 7 nodes
Relaxed ObjVal | Upper Bound Lower Bound | Absolute Diff Relative Diff | Fractional | Time Iter
-1.640000e+02 | -1.640000e+02 -1.655000e+02 | 1.500000e+00 9.063444e-03 | 0 | 0.130 1
-1.640000e+02 | -1.640000e+02 -1.655000e+02 | 1.500000e+00 9.063444e-03 | 0 | 0.133 2
-1.640000e+02 | -1.640000e+02 -1.655000e+02 | 1.500000e+00 9.063444e-03 | 0 | 0.135 3
-1.640000e+02 | -1.640000e+02 -1.647500e+02 | 7.500000e-01 4.552352e-03 | 0 | 0.138 4
-1.640000e+02 | -1.640000e+02 -1.640000e+02 | 0.000000e+00 0.000000e+00 | 0 | 0.141 5
Convergence criteria met.
Attempting to branch.
Model 5 of 7. UB: -164.0, LB:-165.5, Time: 1.292s
[ Info: Solving JuDGE model for tree: Subtree rooted at node 1 containing 7 nodes
Relaxed ObjVal | Upper Bound Lower Bound | Absolute Diff Relative Diff | Fractional | Time Iter
-1.652500e+02 | Inf -1.655000e+02 | 2.500000e-01 1.510574e-03 | 4 | 0.004 1
-1.652500e+02 | Inf -1.652500e+02 | 0.000000e+00 0.000000e+00 | 4 | 0.008 2
-1.652500e+02 | -1.640000e+02 -1.652500e+02 | 1.250000e+00 7.564297e-03 | 0 | 0.008 3*
Convergence criteria met.
Attempting to branch.
Adding 2 new nodes to B&P tree.
Model 6 of 9. UB: -164.0, LB:-165.5, Time: 1.302s
[ Info: Solving JuDGE model for tree: Subtree rooted at node 1 containing 7 nodes
Relaxed ObjVal | Upper Bound Lower Bound | Absolute Diff Relative Diff | Fractional | Time Iter
-1.618750e+02 | Inf -1.655000e+02 | 3.625000e+00 2.190332e-02 | 7 | 0.004 1
-1.621250e+02 | Inf -1.628750e+02 | 7.500000e-01 4.604758e-03 | 2 | 0.008 2
Dominated by incumbent.
Model 7 of 9. UB: -164.0, LB:-165.5, Time: 1.31s
[ Info: Solving JuDGE model for tree: Subtree rooted at node 1 containing 7 nodes
Relaxed ObjVal | Upper Bound Lower Bound | Absolute Diff Relative Diff | Fractional | Time Iter
-1.648333e+02 | Inf -1.655000e+02 | 6.666667e-01 4.028197e-03 | 3 | 0.004 1
-1.648333e+02 | Inf -1.655000e+02 | 6.666667e-01 4.028197e-03 | 3 | 0.007 2
-1.648333e+02 | Inf -1.653333e+02 | 5.000000e-01 3.024194e-03 | 3 | 0.009 3
-1.648333e+02 | Inf -1.650833e+02 | 2.500000e-01 1.514387e-03 | 3 | 0.012 4
-1.648333e+02 | Inf -1.650000e+02 | 1.666667e-01 1.010101e-03 | 3 | 0.015 5
-1.648333e+02 | Inf -1.648333e+02 | 0.000000e+00 0.000000e+00 | 3 | 0.017 6
-1.648333e+02 | -1.635000e+02 -1.648333e+02 | 1.333333e+00 8.088979e-03 | 0 | 0.019 7*
Convergence criteria met.
Attempting to branch.
Adding 2 new nodes to B&P tree.
Model 8 of 11. UB: -164.0, LB:-165.25, Time: 1.33s
[ Info: Solving JuDGE model for tree: Subtree rooted at node 1 containing 7 nodes
Relaxed ObjVal | Upper Bound Lower Bound | Absolute Diff Relative Diff | Fractional | Time Iter
-1.640000e+02 | -1.640000e+02 -1.652500e+02 | 1.250000e+00 7.564297e-03 | 0 | 0.004 1
-1.640000e+02 | -1.640000e+02 -1.652500e+02 | 1.250000e+00 7.564297e-03 | 0 | 0.020 2
-1.640000e+02 | -1.640000e+02 -1.645000e+02 | 5.000000e-01 3.039514e-03 | 0 | 0.023 3
-1.640000e+02 | -1.640000e+02 -1.645000e+02 | 5.000000e-01 3.039514e-03 | 0 | 0.026 4
-1.640000e+02 | -1.640000e+02 -1.640000e+02 | 0.000000e+00 0.000000e+00 | 0 | 0.029 5
Convergence criteria met.
Attempting to branch.
Model 9 of 11. UB: -164.0, LB:-165.25, Time: 1.359s
[ Info: Solving JuDGE model for tree: Subtree rooted at node 1 containing 7 nodes
Relaxed ObjVal | Upper Bound Lower Bound | Absolute Diff Relative Diff | Fractional | Time Iter
-1.640000e+02 | -1.640000e+02 -1.652500e+02 | 1.250000e+00 7.564297e-03 | 0 | 0.004 1
-1.640000e+02 | -1.640000e+02 -1.652500e+02 | 1.250000e+00 7.564297e-03 | 0 | 0.006 2
-1.640000e+02 | -1.640000e+02 -1.645000e+02 | 5.000000e-01 3.039514e-03 | 0 | 0.009 3
-1.640000e+02 | -1.640000e+02 -1.640000e+02 | 0.000000e+00 0.000000e+00 | 0 | 0.011 4
Convergence criteria met.
Attempting to branch.
Model 10 of 11. UB: -164.0, LB:-164.83333333333331, Time: 1.371s
[ Info: Solving JuDGE model for tree: Subtree rooted at node 1 containing 7 nodes
Relaxed ObjVal | Upper Bound Lower Bound | Absolute Diff Relative Diff | Fractional | Time Iter
-1.645000e+02 | Inf -1.648333e+02 | 3.333333e-01 2.022245e-03 | 2 | 0.005 1
-1.645000e+02 | Inf -1.645000e+02 | 0.000000e+00 0.000000e+00 | 2 | 0.007 2
-1.645000e+02 | -1.635000e+02 -1.645000e+02 | 1.000000e+00 6.079027e-03 | 0 | 0.010 3*
Convergence criteria met.
Attempting to branch.
Adding 2 new nodes to B&P tree.
Model 11 of 13. UB: -164.0, LB:-164.83333333333331, Time: 1.382s
[ Info: Solving JuDGE model for tree: Subtree rooted at node 1 containing 7 nodes
Relaxed ObjVal | Upper Bound Lower Bound | Absolute Diff Relative Diff | Fractional | Time Iter
-1.640000e+02 | -1.640000e+02 -1.648333e+02 | 8.333333e-01 5.055612e-03 | 0 | 0.004 1
-1.640000e+02 | -1.640000e+02 -1.648333e+02 | 8.333333e-01 5.055612e-03 | 0 | 0.006 2
-1.640000e+02 | -1.640000e+02 -1.642500e+02 | 2.500000e-01 1.522070e-03 | 0 | 0.009 3
-1.640000e+02 | -1.640000e+02 -1.640000e+02 | 0.000000e+00 0.000000e+00 | 0 | 0.011 4
Convergence criteria met.
Attempting to branch.
Model 12 of 13. UB: -164.0, LB:-164.5, Time: 1.394s
[ Info: Solving JuDGE model for tree: Subtree rooted at node 1 containing 7 nodes
Relaxed ObjVal | Upper Bound Lower Bound | Absolute Diff Relative Diff | Fractional | Time Iter
-1.635000e+02 | -1.635000e+02 -1.645000e+02 | 1.000000e+00 6.079027e-03 | 0 | 0.005 1
-1.640000e+02 | -1.640000e+02 -1.640000e+02 | 0.000000e+00 0.000000e+00 | 0 | 0.008 2
Convergence criteria met.
Attempting to branch.
Model 13 of 13. UB: -164.0, LB:-164.5, Time: 1.402s
[ Info: Solving JuDGE model for tree: Subtree rooted at node 1 containing 7 nodes
Relaxed ObjVal | Upper Bound Lower Bound | Absolute Diff Relative Diff | Fractional | Time Iter
-1.640000e+02 | -1.640000e+02 -1.645000e+02 | 5.000000e-01 3.039514e-03 | 0 | 0.005 1
-1.640000e+02 | -1.640000e+02 -1.645000e+02 | 5.000000e-01 3.039514e-03 | 0 | 0.008 2
-1.640000e+02 | -1.640000e+02 -1.645000e+02 | 5.000000e-01 3.039514e-03 | 0 | 0.010 3
-1.640000e+02 | -1.640000e+02 -1.644167e+02 | 4.166667e-01 2.534212e-03 | 0 | 0.013 4
-1.640000e+02 | -1.640000e+02 -1.644167e+02 | 4.166667e-01 2.534212e-03 | 0 | 0.016 5
-1.640000e+02 | -1.640000e+02 -1.644167e+02 | 4.166667e-01 2.534212e-03 | 0 | 0.019 6
-1.640000e+02 | -1.640000e+02 -1.644167e+02 | 4.166667e-01 2.534212e-03 | 0 | 0.021 7
-1.640000e+02 | -1.640000e+02 -1.640000e+02 | 0.000000e+00 0.000000e+00 | 0 | 0.024 8
Convergence criteria met.
Objective value of best integer-feasible solution: -164.0
Objective value of lower bound: -164.0
Solve time: 1.433s
JuDGE Expansions
Node 1: "invest" 9.0
Node 11: "invest" 3.0
Node 12: "invest" 6.0
We now see that we have found a better solution, and proved it is optimal.
There are several options for the search: :lowestLB
always chooses to branch on the node with the lowest lower bound; :depth_first_dive
performs a depth-first search of the branch and bound tree, once it find a node with an integer relaxation it searches adjacent nodes within the tree; :depth_first_resurface
performs a depth-first search of the branch and bound tree, but once it finds a node with an integer relaxation it returns to the root node and explores the other branch; :breadth_first
performs a breadth-first search of the tree.
There is a default branching method: JuDGE.variable_branch
, but it is also possible to write custom branching methods; see the API for more details.
Tutorial 7: Risk aversion
JuDGE implements risk aversion using the risk measure CVaR over the accumulated profits up to each of leaf nodes in the scenario tree. The objective function minimized is a convex combination of expectation and CVaR, with the parameter λ=1 meaning at all the weight is placed on CVaR. In our implementation CVaR represents the expected cost of the 100α% worst scenarios. In order to implement CVaR, we supply the optional argument risk=Risk(λ,α)
when we construct the JuDGEModel
.
judy = JuDGEModel(mytree, ConditionallyUniformProbabilities, sub_problems, JuDGE_MP_Solver,
risk=Risk(0.5,0.1))
best = JuDGE.branch_and_price(judy,verbose=0)
println("Objective: "*string(best.bounds.UB))
println("Lower Bound: "*string(best.bounds.LB))
println("Expected Costs: "*string(JuDGE.get_objval(best,risk=JuDGE.RiskNeutral()))
JuDGE.print_expansions(best, format=format_output)
[ Info: Establishing JuDGE model for tree: Subtree rooted at node 1 containing 7 nodes
[ Info: Checking sub-problem format
[ Info: Sub-problem format is valid
[ Info: Building master problem
[ Info: Master problem built
Model 1 of 1. UB: Inf, LB:-Inf, Time: 0.0s
Objective value of best integer-feasible solution: -160.5
Objective value of lower bound: -160.5
Solve time: 0.502s
Objective: -160.5
Lower Bound: -160.5
Tutorial 8: Shutdown variables
JuDGE supports shutdown decisions. These are variables that remove capacity when they are activated. The @capitalcosts
of these decisions may be negative, reflecting some salvage value; moreover, the @ongoingcosts
may also be negative, reflecting avoided maintenance costs. Given that this is a shutdown variable, it is important to remember that the capacity being removed should be accounted for elsewhere within the subproblem.
function sub_problems_shutdown(node)
model = Model(JuDGE_SP_Solver)
@shutdown(model, bag, Bin)
@capitalcosts(model, -bag*invest_cost[node])
@variable(model, y[1:num_items], Bin)
@constraint(model, BagExtension, sum(y[i]*item_volume[node][i] for i in 1:num_items) <=
7 - bag)
@objective(model, Min, sum(-item_reward[node][i] * y[i] for i in 1:num_items))
return model
end
judy = JuDGEModel(mytree, ConditionallyUniformProbabilities, sub_problems_shutdown,
JuDGE_MP_Solver)
JuDGE.solve(judy,verbose=1)
println("Objective: "*string(objective_value(judy.master_problem)))
JuDGE.print_expansions(judy)
[ Info: Establishing JuDGE model for tree: Subtree rooted at node 1 containing 7 nodes
[ Info: Checking sub-problem format
[ Info: Sub-problem format is valid
[ Info: Building master problem
[ Info: Master problem built
[ Info: Solving JuDGE model for tree: Subtree rooted at node 1 containing 7 nodes
Relaxed ObjVal | Upper Bound Lower Bound | Absolute Diff Relative Diff | Fractional | Time Iter
-1.500000e+01 | -1.500000e+01 -Inf | Inf NaN | 0 | 0.138 1
-1.427500e+02 | -1.427500e+02 -1.597500e+02 | 1.700000e+01 1.064163e-01 | 0 | 0.187 2
-1.450000e+02 | -1.450000e+02 -1.597500e+02 | 1.475000e+01 9.233177e-02 | 0 | 0.188 3
-1.450000e+02 | -1.450000e+02 -1.597500e+02 | 1.475000e+01 9.233177e-02 | 0 | 0.190 4
-1.450000e+02 | -1.450000e+02 -1.450000e+02 | 0.000000e+00 0.000000e+00 | 0 | 0.191 5
Convergence criteria met.
Objective: -145.0
JuDGE Expansions
Node 111: "bag" 1.0
Node 121: "bag" 1.0
Tutorial 9: Side-constraints
JuDGE supports side-constraints being added to the master problem. These can be constraints across expansion variables at a single node, or can be constraints on variables corresponding to different nodes. The JuDGE.history
function can be useful if applying logical constraints on expansion variables. In order to apply a budget constraint at each node we can define the following function:
function budget(model,tree)
for node in collect(tree)
@constraint(model,sum(invest_cost[node]*invest_volume[i]*invest[node][i]
for i in 1:num_invest)<=40)
end
end
We now can define a JuDGEModel
with these side-constraints, and solve it using branch and price.
judy = JuDGEModel(mytree, ConditionallyUniformProbabilities, sub_problems, JuDGE_MP_Solver,
sideconstraints=budget)
judy = JuDGE.branch_and_price(judy, search=:lowestLB, branch_method=JuDGE.constraint_branch)
JuDGE.print_expansions(judy, format=format_output)
Tutorial 10: State variables
JuDGE has experimental support for state variables. These are variables that are part of the master problem, and are either increased or decreased by actions within the subproblems. See the example inventory.jl for details of how these variables can be implemented.
Tutorial 11: Set partitioning / packing models
JuDGE has experimental support for set partitioning / packing models. The documentation for this hasn't been written yet, but examples of a vehicle routing model (vrp.jl) and a cutting stock model (cutting_stock.jl) are provided in the examples directory.