diff --git a/README.md b/README.md
index 842383736e0cb1dabd8c9d46ee578b7c282c200b..aeba261e0a96db8103a227e872fcadb9f6c7ea8b 100644
--- a/README.md
+++ b/README.md
@@ -4,19 +4,15 @@ This repository contains numerical procedures linked to the paper
**> ADD CITATION**
-If you encounter a bug or something unexpected, please let me know by raising an issue on the project page or by contacting me by email at `theo.guyard@insa-rennes.fr`.
+If you encounter a bug or something unexpected, please let me know by raising an issue on the project page or by contacting me by [mail](mailto:theo.guyard@insa-rennes.fr).
## Requirements
This repository works with Julia v1.5+. Please refer to the [Julia download page](https://julialang.org/downloads/) for install instructions.
-Dependencies :
-* [JuMP.jl](https://github.com/jump-dev/JuMP.jl)
-* [CPLEX.jl](https://github.com/jump-dev/CPLEX.jl)
-* [Distributions.jl](https://github.com/JuliaStats/Distributions.jl)
-* [StatsBase](https://github.com/JuliaStats/StatsBase.jl)
+**Dependencies :** [JuMP.jl](https://github.com/jump-dev/JuMP.jl), [CPLEX.jl](https://github.com/jump-dev/CPLEX.jl), [Distributions.jl](https://github.com/JuliaStats/Distributions.jl), [StatsBase.jl](https://github.com/JuliaStats/StatsBase.jl)
-Our code works with the CPLEX solver through [CPLEX.jl](https://github.com/jump-dev/CPLEX.jl). You cannot use CPLEX without having purchased and installed a copy of CPLEX Optimization Studio from [IBM](https://www.ibm.com). However, CPLEX is available for free to [academics and students](https://community.ibm.com/community/user/datascience/blogs/xavier-nodet1/2020/07/09/cplex-free-for-students). We recommend to follow carefully the install instructions described at [CPLEX.jl](https://github.com/jump-dev/CPLEX.jl).
+Our code works with the [CPLEX](https://www.ibm.com/fr-fr/analytics/cplex-optimizer) solver through [CPLEX.jl](https://github.com/jump-dev/CPLEX.jl). You cannot use CPLEX without having purchased and installed a copy of CPLEX Optimization Studio from [IBM](https://www.ibm.com). However, CPLEX is available for free to [academics and students](https://community.ibm.com/community/user/datascience/blogs/xavier-nodet1/2020/07/09/cplex-free-for-students). We recommend you to follow carefully the install instructions described at [CPLEX.jl](https://github.com/jump-dev/CPLEX.jl) to install the package.
## Install instructions
@@ -25,7 +21,7 @@ Our code works with the CPLEX solver through [CPLEX.jl](https://github.com/jump-
git clone https://gitlab.insa-rennes.fr/Theo.Guyard/bnb-screening.git
```
-2. Enter in the project folder
+2. Enter in the project root folder
```bash
cd bnb-screening
```
@@ -39,27 +35,32 @@ julia -v
4. Install dependencies using [Julia's Pkg](https://docs.julialang.org/en/v1/stdlib/Pkg/) in REPL mode
```julia
(@v1.5) pkg> activate .
-(@v1.5) pkg> add JuMP
-(@v1.5) pkg> add CPLEX
-(@v1.5) pkg> add Distributions
-(@v1.5) pkg> add StatsBase
+(bnb-screening) pkg> instantiate
+```
+It is possible that you also have to build [CPLEX.jl](https://github.com/jump-dev/CPLEX.jl) manually with the command
+```julia
+(bnb-screening) pkg> build CPLEX
```
## ICASSP experiments
-Experiments presented in the paper corresponding to this repository are located in the folder `bnb-screening/exp/ICASSP/`. The available experiment are
-* `performances.jl` corresponding to table 1
+Experiments presented in the paper submitted to [ICASSP 2022](https://2022.ieeeicassp.org) are located in the folder `bnb-screening/exp/ICASSP/`. The available experiments are
+* `performances.jl` corresponding to results of table 1
-Run experiments directly from the project root folder :
+Run experiments directly from the project root folder as follows.
```bash
julia --project=. -t 1 exp/ICASSP/performances.jl
```
-The `--project=.` argument installs package dependencies if needed. The flag `-t 1` restricts computations to only 1 CPU core in order to avoid bias due to parallelization capabilities. Logs are saved in the same folder as `performances.jl`. You can edit the experimental setup directly in the script, which is well commented.
+The `--project=.` flag installs and build dependencies if needed. The `-t 1` flag restricts computations to only 1 CPU core in order to avoid bias due to parallelization capabilities. Logs of the experiment are saved in the same folder as `performances.jl`. You can edit the experimental setup directly in the script, which is well commented.
## Running demo examples
-The `example` folder contains portions of code explaining how to use our package.
+The `bnb-screening/examples/` folder contains portions of code explaining how to use our package. You can run these scripts as follows.
+```bash
+julia --project=. examples/solve_bnb.jl
+julia --project=. examples/solve_mip.jl
+```
## Licence
@@ -69,4 +70,4 @@ This software is distributed under the [MIT Licence](https://mit-license.org).
If you use this package for your own work, please consider citing it as :
-**ADD CITATION**
+**> ADD CITATION**
diff --git a/examples/solve_bnb.jl b/examples/solve_bnb.jl
index d08ade7c26090f19a3cbbfbd3309b2ca4a9bff3d..4f04cb56d527206221e54d5f0325aa800242c076 100644
--- a/examples/solve_bnb.jl
+++ b/examples/solve_bnb.jl
@@ -28,12 +28,5 @@ bnbparams = BnbParams(
showevery=1, # If verbosity=true, number of nodes between two displays
)
-# Solve the problem using the BnB algorithm
-result_bnb = solve_bnb(A,y,λ,M,bnbparams=BnbParams(verbosity=false))
-
-# Display BnB results and statistics
-println("BnB")
-println(" Status : $(result_bnb.termination_status)")
-println(" Objective : $(result_bnb.objective_value)")
-println(" Nodes : $(result_bnb.node_count)")
-println(" Timer : $(result_bnb.solve_time)")
\ No newline at end of file
+# Solve the problem using the BnB algorithm with appropriate parameters
+result_bnb = solve_bnb(A,y,λ,M,bnbparams=bnbparams)
diff --git a/examples/solve_mip.jl b/examples/solve_mip.jl
index 4630f50e81c704c6472a5b73f0d907964fad585d..5af1f4e23974c747cf1404902365a6e6a912a493 100644
--- a/examples/solve_mip.jl
+++ b/examples/solve_mip.jl
@@ -16,11 +16,11 @@ mipparams = MipParams(
silent = false, # Toogle solver verbosity
)
-# Solve the problem using a MIP solver
-result_mip = solve_mip(A,y,λ,M)
+# Solve the problem using a MIP solver with appropriate parameters
+result_mip = solve_mip(A,y,λ,M,mipparams=mipparams)
# Display solver results and statistics
-println("MIP")
+println("\nMIP")
println(" Status : $(result_mip.termination_status)")
println(" Objective : $(result_mip.objective_value)")
println(" Nodes : $(result_mip.node_count)")
diff --git a/exp/ICASSP/performances.jl b/exp/ICASSP/performances.jl
index a668323127ebfc72914f268073a2a38a95fd2071..771abcbccd1dedb95562b4b9f7998c74cbae3c95 100644
--- a/exp/ICASSP/performances.jl
+++ b/exp/ICASSP/performances.jl
@@ -8,17 +8,22 @@ using StatsBase
# This is the only part of code you may have to edit
setup = "toeplitz" # Setup to use : "gaussian" or "toeplitz"
-k = 7 # Sparsity level
+k = 5 # Sparsity level
m, n = 500, 300 # Dictionnary dimension
snr = 10. * log10(10.) # SNR of the observation noise
maxtime = 1000 # Maximum solution time (in second) allowed
maxiter = 1000 # Maximum Active-Set iterations allowed per relaxation resolution
tol = 1.e-8 # Integer tolerance (ie, |x|<tol <=> x=0)
-repeats = 100 # Number of repeats in the experiment
+repeats = 10 # Number of repeats in the experiment
################################## Experiment ##################################
+function print_and_log(s::String, f::IOStream)
+ print(s)
+ write(f,s)
+end
+
function precompilation(
setup::String,
k::Int,
@@ -29,7 +34,7 @@ function precompilation(
maxiter::Int,
tol::Float64,
)
- println("Precomiling methods ...")
+ print_and_log("Precompiling methods ...\n",logfile)
try
if setup == "gaussian"
A, y, λ, M = data_gaussian(k,m,n,snr)
@@ -42,10 +47,12 @@ function precompilation(
result_bnb = solve_bnb(A,y,λ,M,bnbparams=BnbParams(maxtime=maxtime,maxiter=maxiter,tol=tol,verbosity=false))
result_bnbs = solve_bnb(A,y,λ,M,bnbparams=BnbParams(maxtime=maxtime,maxiter=maxiter,tol=tol,verbosity=false,screening=true))
catch e
- e == InterruptException() && rethrow(e)
- println("Precompilation failed with the following error")
- throw(e)
+ print_and_log("Precompilation failed with the following error\n",logfile)
+ print_and_log("$e\n",logfile)
+ print_and_log("Stopping experiment\n",logfile)
+ rethrow(e)
end
+ print_and_log("\n",logfile)
return nothing
end
@@ -57,29 +64,18 @@ logfile_name = string(
)
logfile = open(logfile_name, "w")
-println("--------------------------------")
-println("ICASSP experiment - Performances\n")
-println("Configuration")
-println(" setup : $(setup)")
-println(" k : $(k)")
-println(" m, n : $(m), $(n)")
-println(" snr : $(snr)")
-println(" maxtime : $(maxtime)")
-println(" maxiter : $(maxiter)")
-println(" tol : $(tol)")
-println(" repeats : $(repeats)")
-println()
-
-write(logfile,"ICASSP experiment - Performances\n\n")
-write(logfile,"Configuration\n")
-write(logfile," setup : $(setup)\n")
-write(logfile," k : $(k)\n")
-write(logfile," m, n : $(m), $(n)\n")
-write(logfile," snr : $(snr)\n")
-write(logfile," maxtime : $(maxtime)\n")
-write(logfile," maxiter : $(maxiter)\n")
-write(logfile," tol : $(tol)\n")
-write(logfile," repeats : $(repeats)\n")
+print_and_log("--------------------------------\n",logfile)
+print_and_log("ICASSP experiment - Performances\n\n",logfile)
+print_and_log("Configuration\n",logfile)
+print_and_log(" setup : $(setup)\n",logfile)
+print_and_log(" k : $(k)\n",logfile)
+print_and_log(" m, n : $(m), $(n)\n",logfile)
+print_and_log(" snr : $(snr)\n",logfile)
+print_and_log(" maxtime : $(maxtime)\n",logfile)
+print_and_log(" maxiter : $(maxiter)\n",logfile)
+print_and_log(" tol : $(tol)\n",logfile)
+print_and_log(" repeats : $(repeats)\n",logfile)
+print_and_log("\n",logfile)
fail_mip = Vector{Union{Bool,Missing}}(missing,repeats)
fail_bnb = Vector{Union{Bool,Missing}}(missing,repeats)
@@ -94,7 +90,7 @@ time_bnbs = Vector{Union{Float64,Missing}}(missing,repeats)
precompilation(setup,k,m,n,snr,maxtime,maxiter,tol)
for i in 1:repeats
- println("Repeat $i/$repeats")
+ print_and_log("Repeat $i/$repeats\n",logfile)
try
# Generate data
if setup == "gaussian"
@@ -130,42 +126,48 @@ for i in 1:repeats
else
fail_bnbs[i] = true
end
+
+ # Log results
+ print_and_log(" Fails\n",logfile)
+ print_and_log(" MIP : $(fail_mip[i])\n",logfile)
+ print_and_log(" BnB : $(fail_bnb[i])\n",logfile)
+ print_and_log(" SBnB : $(fail_bnbs[i])\n",logfile)
+ print_and_log(" Node count\n",logfile)
+ print_and_log(" MIP : $(node_mip[i])\n",logfile)
+ print_and_log(" BnB : $(node_bnb[i])\n",logfile)
+ print_and_log(" SBnB : $(node_bnbs[i])\n",logfile)
+ print_and_log(" Solve time\n",logfile)
+ print_and_log(" MIP : $(round(time_mip[i],digits=3)) seconds\n",logfile)
+ print_and_log(" BnB : $(round(time_bnb[i],digits=3)) seconds\n",logfile)
+ print_and_log(" SBnB : $(round(time_bnbs[i],digits=3)) seconds\n",logfile)
+
catch e
- e == InterruptException() && rethrow(e)
- println("warning : Repeat $i failed with the following error")
- println(e)
+ if e == InterruptException()
+ close(logfile)
+ rethrow(e)
+ end
+ print_and_log("warning : Repeat $i failed with the following error\n",logfile)
+ print_and_log("$e\n",logfile)
+ print_and_log("Skipping this repeat\n",logfile)
end
end
################################### Results ####################################
-println()
-println("Fails")
-println(" MIP : $(100 * mean(skipmissing(fail_mip)))%")
-println(" BnB : $(100 * mean(skipmissing(fail_bnb)))%")
-println(" SBnB : $(100 * mean(skipmissing(fail_bnbs)))%")
-println("Node count")
-println(" MIP : $(mean(skipmissing(node_mip)))")
-println(" BnB : $(mean(skipmissing(node_bnb)))")
-println(" SBnB : $(mean(skipmissing(node_bnbs)))")
-println("Solve time")
-println(" MIP : $(round(mean(skipmissing(time_mip)),digits=3)) seconds")
-println(" BnB : $(round(mean(skipmissing(time_bnb)),digits=3)) seconds")
-println(" SBnB : $(round(mean(skipmissing(time_bnbs)),digits=3)) seconds")
-println("\n--------------------------------")
-
-write(logfile,"\n")
-write(logfile,"Fails\n")
-write(logfile," MIP : $(100 * mean(skipmissing(fail_mip)))%\n")
-write(logfile," BnB : $(100 * mean(skipmissing(fail_bnb)))%\n")
-write(logfile," SBnB : $(100 * mean(skipmissing(fail_bnbs)))%\n")
-write(logfile,"Node count\n")
-write(logfile," MIP : $(mean(skipmissing(node_mip)))\n")
-write(logfile," BnB : $(mean(skipmissing(node_bnb)))\n")
-write(logfile," SBnB : $(mean(skipmissing(node_bnbs)))\n")
-write(logfile,"Solve time\n")
-write(logfile," MIP : $(round(mean(skipmissing(time_mip)),digits=3)) seconds\n")
-write(logfile," BnB : $(round(mean(skipmissing(time_bnb)),digits=3)) seconds\n")
-write(logfile," SBnB : $(round(mean(skipmissing(time_bnbs)),digits=3)) seconds\n")
-close(logfile)
\ No newline at end of file
+print_and_log("\n",logfile)
+print_and_log("Results\n",logfile)
+print_and_log(" Fails (mean over $repeats repeats)\n",logfile)
+print_and_log(" MIP : $(100 * mean(skipmissing(fail_mip)))%\n",logfile)
+print_and_log(" BnB : $(100 * mean(skipmissing(fail_bnb)))%\n",logfile)
+print_and_log(" SBnB : $(100 * mean(skipmissing(fail_bnbs)))%\n",logfile)
+print_and_log(" Node count (mean over $repeats repeats)\n",logfile)
+print_and_log(" MIP : $(mean(skipmissing(node_mip)))\n",logfile)
+print_and_log(" BnB : $(mean(skipmissing(node_bnb)))\n",logfile)
+print_and_log(" SBnB : $(mean(skipmissing(node_bnbs)))\n",logfile)
+print_and_log(" Solve time (mean over $repeats repeats)\n",logfile)
+print_and_log(" MIP : $(round(mean(skipmissing(time_mip)),digits=3)) seconds\n",logfile)
+print_and_log(" BnB : $(round(mean(skipmissing(time_bnb)),digits=3)) seconds\n",logfile)
+print_and_log(" SBnB : $(round(mean(skipmissing(time_bnbs)),digits=3)) seconds\n",logfile)
+print_and_log("--------------------------------",logfile)
+close(logfile)
diff --git a/src/data.jl b/src/data.jl
index d1d1f55a3957e4a0cdfcab86d61134aecdee8505..f3d730788600ebec214d5f62c60d8ad3880225fe 100644
--- a/src/data.jl
+++ b/src/data.jl
@@ -9,6 +9,8 @@ function validate_data(
M > 0 || throw(ArgumentError("M must be positive"))
end
+# ============================================================================ #
+
"""
data_gaussian(
k::Int=5,