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:

  1. The definition of a Tree;
  2. defining the subproblems of the JuDGEModel;
  3. building the JuDGEModel;
  4. 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.