diff --git a/DantzigWolfeDecomposition.jl b/DantzigWolfeDecomposition.jl new file mode 100644 index 0000000..b0844f8 --- /dev/null +++ b/DantzigWolfeDecomposition.jl @@ -0,0 +1,61 @@ +# Dantzig-Wolfe Reformulation and Column Generation +# Edward J. Xu +# 2019.5.27 +######################################################################################################################## +module DantzigWolfeDecomposition +using JuMP +using Gurobi +using CPLEX +using GLPKMathProgInterface +using LinearAlgebra +using MathProgBase +using Random +using SparseArrays +include("FuncSub.jl") +include("FuncMas.jl") +include("FuncStab.jl") +include("FuncOptim.jl") +export doDWDecomp, Sense, leq, geq, eq +######################################################################################################################## + + +@enum Sense begin + leq = 1 + geq = 2 + eq = 3 +end + + +function doDWDecomp(mat_a, vec_b, vec_c, vecSenseAll, indexMas, blocks, indexSub, numXPerSub, + dualPen, dualPenMult, dualPenThreshold, epsilon, whePrint::Bool, whiSolver) + if whiSolver == 1 + gurobi_env = Gurobi.Env() + elseif whiSolver == 2 + gurobi_env = [] + else + gurobi_env = [] + end + println("#### 1/4, Set vecModelSub #####################################################") ######################## + vecModelSub = setModelSub(mat_a, vec_b, vec_c, vecSenseAll, indexMas, blocks, indexSub, gurobi_env, whiSolver) + numSub = length(vecModelSub) + numQ = size(vecModelSub[1].mat_e)[1] # [num of constraints in matA_0] + vecStrNameVar = getStrNameVar(numSub, numXPerSub) # define variable names + println("#### 2/4, Set modelMas ########################################################") ######################## + vecSenseP = deepcopy(vecSenseAll[collect(indexMas)]) + vecP = deepcopy(vec_b[collect(indexMas)]) + (modMas, vecConsRef, vecConsConvex, vecLambda, vecMuMinus, vecMuPlus, vecMuMinusConv, vecMuPlusConv) = + setModelMas(numQ, vecP, numSub, vecSenseP, dualPen, gurobi_env, whiSolver) + ## 5, Optimization + println("#### 3/4, Begin Optim #########################################################") ######################## + (vecLambdaResult, extremePointForSub, extremePoints) = doOptim( + vecModelSub, modMas, vecConsRef, vecConsConvex, vecLambda, # Para for ColGen + dualPen, dualPenMult, dualPenThreshold, vecMuMinus, vecMuPlus, vecMuMinusConv, vecMuPlusConv, # Para for Stable + epsilon, whePrint + ) + println("#### 4/4, Print Result ########################################################") ######################## + printResult(vecStrNameVar, indexSub, vecLambda, vecLambdaResult, extremePointForSub, extremePoints) + println("################################################################################") +end + + +end diff --git a/FuncMas.jl b/FuncMas.jl index 253a6ce..079e1d5 100644 --- a/FuncMas.jl +++ b/FuncMas.jl @@ -5,11 +5,17 @@ ######################################################################################################################## -function setModelMas(numQ, vecP, numSub, vecSenseP, dualPen) +function setModelMas(numQ, vecP, numSub, vecSenseP, dualPen, gurobi_env, whiSolver) vecDualGuessPi = 100 * ones(Float64, numQ) vecDualGuessKappa = 0 * ones(Float64, numSub) ## X1: initial matrix of extreme points - modMas = Model(solver = GurobiSolver(OutputFlag = 0, gurobi_env)) + if whiSolver == 1 + modMas = Model(solver = GurobiSolver(OutputFlag = 0, gurobi_env)) + elseif whiSolver == 2 + modMas = Model(solver = CplexSolver(CPX_PARAM_SCRIND=0)) + else + modMas = Model(solver = GLPKSolverLP()) + end K = 1 # In this case we do not use a starting set of extreme points. # we just use a dummy starting column diff --git a/FuncDW.jl b/FuncOptim.jl similarity index 71% rename from FuncDW.jl rename to FuncOptim.jl index f9e7477..429b7d1 100644 --- a/FuncDW.jl +++ b/FuncOptim.jl @@ -111,28 +111,3 @@ function printResult(vecStrNameVar, indexSub, vecLambda, vecLambdaResult, extrem end end end - - -function doDWDecomp(mat_a, vec_b, vec_c, vecSenseAll, indexMas, blocks, indexSub, numXPerSub, - dualPen, dualPenMult, dualPenThreshold, epsilon, whePrint::Bool) - println("#### 2/5, Set vecModelSub #####################################################") ######################## - vecModelSub = setModelSub(mat_a, vec_b, vec_c, vecSenseAll, indexMas, blocks, indexSub) - numSub = length(vecModelSub) - numQ = size(vecModelSub[1].mat_e)[1] # [num of constraints in matA_0] - vecStrNameVar = getStrNameVar(numSub, numXPerSub) # define variable names - println("#### 3/5, Set modelMas ########################################################") ######################## - vecSenseP = deepcopy(vecSenseAll[collect(indexMas)]) - vecP = deepcopy(vec_b[collect(indexMas)]) - (modMas, vecConsRef, vecConsConvex, vecLambda, vecMuMinus, vecMuPlus, vecMuMinusConv, vecMuPlusConv) = - setModelMas(numQ, vecP, numSub, vecSenseP, dualPen) - ## 5, Optimization - println("#### 4/5, Begin Optim #########################################################") ######################## - (vecLambdaResult, extremePointForSub, extremePoints) = doOptim( - vecModelSub, modMas, vecConsRef, vecConsConvex, vecLambda, # Para for ColGen - dualPen, dualPenMult, dualPenThreshold, vecMuMinus, vecMuPlus, vecMuMinusConv, vecMuPlusConv, # Para for Stable - epsilon, whePrint - ) - println("#### 5/5, Print Result ########################################################") ######################## - printResult(vecStrNameVar, indexSub, vecLambda, vecLambdaResult, extremePointForSub, extremePoints) - println("################################################################################") -end diff --git a/FuncSub.jl b/FuncSub.jl index bd65d9d..940ad88 100644 --- a/FuncSub.jl +++ b/FuncSub.jl @@ -34,26 +34,33 @@ function setVariableSub(modSub::JuMP.Model, mat_d, vec_q, vecSense) end -function setBranch(vecModelSub, i, j) - vec_new = zeros(1, 15) - vec_new[j] = -1 - vecModelSub[i].mat_d = vcat(vecModelSub[i].mat_d, vec_new) - vecModelSub[i].vec_q = vcat(vecModelSub[i].vec_q, [-1]) - println(vecModelSub[i].mat_d) - println(vecModelSub[i].vec_q) - return vecModelSub -end +# function setBranch(vecModelSub, i, j) +# vec_new = zeros(1, 15) +# vec_new[j] = -1 +# vecModelSub[i].mat_d = vcat(vecModelSub[i].mat_d, vec_new) +# vecModelSub[i].vec_q = vcat(vecModelSub[i].vec_q, [-1]) +# println(vecModelSub[i].mat_d) +# println(vecModelSub[i].vec_q) +# return vecModelSub +# end -function setModelSub(mat_a, vec_b, vec_c, vecSense, indexMas, blocks, indexSub) +function setModelSub(mat_a, vec_b, vec_c, vecSense, indexMas, blocks, indexSub, gurobi_env, whiSolver) numSub = length(blocks) vecModelSub = Vector{ModelSub}(undef, numSub) + if whiSolver == 1 + modWhiSolver = Model(solver = GurobiSolver(OutputFlag = 0, gurobi_env)) + elseif whiSolver == 2 + modWhiSolver = Model(solver = CplexSolver(CPX_PARAM_SCRIND=0)) + else + modWhiSolver = Model(solver = GLPKSolverLP()) + end for k = 1:numSub vecModelSub[k] = ModelSub( - Model(solver = GurobiSolver(OutputFlag = 0, gurobi_env)), # mod ! - deepcopy(mat_a[collect(indexMas), indexSub[k]]), # mat_e <- A0 + modWhiSolver, # mod + deepcopy(mat_a[collect(indexMas), indexSub[k]]), # mat_e deepcopy(vec_c[indexSub[k]]), # vec_l - deepcopy(mat_a[blocks[k], indexSub[k]]), # mat_d <- ASub + deepcopy(mat_a[blocks[k], indexSub[k]]), # mat_d deepcopy(vec_b[blocks[k]]), # vec_q 0, # vec_x deepcopy(vecSense[blocks[k]]) # vecSense diff --git a/README.md b/README.md index f14f5c8..6d1fd0d 100644 --- a/README.md +++ b/README.md @@ -19,36 +19,62 @@ More detailed explanation can be found in [edxu96/DantzigWolfeDecomposition/wiki ## 2, How to use -The code is in four files: +The code is in five files, which should be stored in your working directory: ``` -FuncDW.jl +DantzigWolfeDecomposition.jl + FuncDW.jl FuncSub.jl FuncMas.jl FuncStab.jl ``` -To use the function `doDWDecomp`: +Load other necessary packages: ```Julia +using JuMP +using CPLEX +using Gurobi +using GLPKMathProgInterface +using LinearAlgebra +using MathProgBase +using SparseArrays +``` + +There is only `doDWDecomp` in the module: + +```Julia +using DantzigWolfeDecomposition ## Set parameter and start DW-Decomposition dualPen = 10 dualPenMult = 0.1 dualPenThreshold = 0.01 - 1e-5 -epsilon = 0.00001 -whePrint = false +epsilon = 0.00001 # [maximum difference between two bounds] +whePrint = false # [whether to print detailed info] +whiSolver = 2 # 1: Gurobi, 2: CPLEX, 3: GLPK doDWDecomp( mat_a, vec_b, vec_c, # Data in LP Problem vecSenseAll, indexMas, blocks, indexSub, numXPerSub, # Data for DW-Decomp dualPen, dualPenMult, dualPenThreshold, # Para for Stable - epsilon, whePrint # Control Para + epsilon, whePrint, whiSolver # Control Para ) ``` -In `main.jl`, you can see the example to use `doDWDecomp` in solving Generalized Assignment Problem. +Attention: The support for `GLPK` is not very well. The `Gurobi` solver does better than `CPLEX`. If you don't have `Gurobi` or `CPLEX`, you may need to comment out the two lines in the `DantzigWolfeDecomposition.jl` file. + +```Julia +using CPLEX +using Gurobi +``` + +In `Test/test.jl`, you can see the example to use `doDWDecomp` to solve Generalized Assignment Problem. ## 3, Contributions Edward J. Xu is maintaining the project. It's originally inspired by Professor Stefan Røpke, DTU Management. + +1. Tebboth, J. R., 2001. A computational study of Dantzig-Wolfe decomposition. University of Buckingham. + +[1]: http://eaton.math.rpi.edu/CourseMaterials/PreviousSemesters/PreviousSemesters/Spring08/JM6640/tebboth.pdf diff --git a/Data_GAP_1.jl b/Test/Data_GAP_1.jl similarity index 100% rename from Data_GAP_1.jl rename to Test/Data_GAP_1.jl diff --git a/Data_GAP_2.jl b/Test/Data_GAP_2.jl similarity index 100% rename from Data_GAP_2.jl rename to Test/Data_GAP_2.jl diff --git a/FuncData.jl b/Test/FuncData.jl similarity index 100% rename from FuncData.jl rename to Test/FuncData.jl diff --git a/main.jl b/Test/test.jl similarity index 60% rename from main.jl rename to Test/test.jl index a0f26f9..dcf22b1 100644 --- a/main.jl +++ b/Test/test.jl @@ -2,47 +2,44 @@ # Edward J. Xu # 2019.5.26 ######################################################################################################################## -# include("$(homedir())/Documents/Github/DantzigWolfeDecomposition/main.jl") +# include("$(homedir())/Documents/Github/DantzigWolfeDecomposition/Test/test.jl") push!(LOAD_PATH, "$(homedir())/Documents/Github/DantzigWolfeDecomposition") cd("$(homedir())/Documents/Github/DantzigWolfeDecomposition") using JuMP -# using CPLEX +using CPLEX +using GLPKMathProgInterface using Gurobi using LinearAlgebra using MathProgBase -using Random using SparseArrays -include("FuncSub.jl") -include("FuncMas.jl") -include("FuncDW.jl") +######################################################################################################################## include("FuncData.jl") -include("FuncProcess.jl") -include("FuncStab.jl") +using DantzigWolfeDecomposition ######################################################################################################################## function main(strFileName::String, series::Int64) - timeStart = time() - global gurobi_env = Gurobi.Env() - println("#### 1/5, Prepare Data ########################################################") ######################## + println("#### 1/2, Prepare Data ########################################################") ######################## (mat_a, vec_b, vec_c, vecSenseAll, indexMas, blocks, indexSub, numXPerSub) = setGeneralAssignProbDate( strFileName, series) + println("#### 2/2, Begin Optim #########################################################") ######################## + timeStart = time() ## Set parameter and start DW-Decomposition dualPen = 10 dualPenMult = 0.1 dualPenThreshold = 0.01 - 1e-5 - epsilon = 0.00001 - whePrint = false + epsilon = 0.00001 # [maximum difference between two bounds] + whePrint = false # [whether to print detailed info] + whiSolver = 3 # 1: Gurobi, 2: CPLEX, 3: GLPK doDWDecomp( mat_a, vec_b, vec_c, # Data in LP Problem vecSenseAll, indexMas, blocks, indexSub, numXPerSub, # Data for DW-Decomp dualPen, dualPenMult, dualPenThreshold, # Para for Stable - epsilon, whePrint # Control Para + epsilon, whePrint, whiSolver # Control Para ) - println("Elapsed time is $(time() - timeStart) seconds.\n", - "################################################################################") + println("Elapsed time is $(time() - timeStart) seconds.") + println("#### End #######################################################################") ######################## end -@enum Sense leq = 1 geq = 2 eq = 3 main("Data/gapa.txt", 1)