Include existing MOMIBB implementation

.gitignore

@@ -0,0 +1,6 @@
+## GitIgnore for HyPaD
+### MATLAB temporary files
+*.aux
+
+### Gurobi problem file 
+*.lp

README.md

@@ -1 +1,20 @@
 # MOMIBB
+This repository contains a MATLAB implementation of the Multiobjective Mixed-Integer Branch-and-Bound Algorithm (MOMIBB).
+It allows to solve both convex and nonconvex multiobjective mixed-integer optimization problems.
+
+## Getting Started
+You can start using MOMIBB by downloading or cloning this repository using the green button near the top of the [GitHub page](https://github.com/LeoWarnow/MOMIBB).
+Then, just open the [UserFile.m](https://github.com/LeoWarnow/MOMIBB/blob/main/UserFile.m).
+It serves as an interface to the MOMIBB algorithm and provides all instructions that you need to get started.
+
+## References
+As this implementation is based on the scientific work by [Gabriele Eichfelder](https://www.tu-ilmenau.de/mmor/team/gabriele-eichfelder/), [Oliver Stein](https://kop.ior.kit.edu/stein.php), and [Leo Warnow](https://www.tu-ilmenau.de/mmor/team/leo-warnow/), please cite the corresponding preprint when using this code:
+````
+@Misc{ESW2022,
+  author       = {Gabriele Eichfelder and Oliver Stein and Leo Warnow},
+  howpublished = {\url{https://optimization-online.org/?p=18696}},
+  title        = {A deterministic solver for multiobjective mixed-integer convex and nonconvex optimization},
+  year         = {2022},
+  journal      = {Optimization Online},
+}
+````

UserFile.m

@@ -0,0 +1,35 @@
+%%% UserFile
+% This is the main interface file for MOMIBB. Please make sure that you 
+% have installed and started Gurobi and INTLAB before calling AdEnA.
+% Call initSession.m prior to first use
+
+%% Some clean-up first
+clear;
+close all;
+clc;
+
+%% Please enter your parameters below
+% Your Problem
+problem = 'T4';
+param = [2;2];
+
+% Provide initial bounds yourself or leave empty to auto-compute (INTLAB is
+% required for auto-compute)
+L = [];
+U = [];
+
+% Set quality epsilon and offset
+EPSILON = 0.1;
+OFFSET = EPSILON*1e-3;
+
+% Specify cutting strategy [0 == none, 1 == continuous, 2 == continuous+mixed integer]
+CUT_MODE = 1;
+
+% Specify handling of quadratic instances [1 == normal, 2 == directly using Gurobi]
+SOL_MODE = 1;
+
+% Should the result be plotted (m = 2 and m = 3 only) [1 == yes, 0 == no]
+plot_result = 1;
+
+%% Call solver
+[L,U,N,box_struct,it,exitflag,time] = callSolver(problem,param,L,U,CUT_MODE,SOL_MODE,EPSILON,OFFSET,plot_result);

callSolver.m

@@ -0,0 +1,96 @@
+function [L,U,N,box_struct,it,exitflag,time] = callSolver(problem_name,problem_param,L,U,CUT_MODE,SOL_MODE,EPSILON,OFFSET,plot_result)
+
+%% Load problem
+problem = str2func(problem_name);
+problem_quadr = [];
+
+%% Handling of incorrect or missing input
+[~,~,p,~,f,~,~,~,~,~,~,~,lb,ub,~,~,is_quadratic] = problem(problem_param);
+
+% Applying standard tolerances
+if isempty(EPSILON)
+    EPSILON = 0.1;
+end
+if isempty(OFFSET)
+    OFFSET = EPSILON*1e-3;
+end
+
+% Checking MOMIBB parameters
+if isempty(SOL_MODE)
+    SOL_MODE = 1;
+    disp('Selected SOL_MODE = 1');
+elseif SOL_MODE == 2
+	if ~is_quadratic
+        SOL_MODE = 1;
+        warning('No quadratic formulation available. Changing to SOL_MODE = 1.');
+    else
+        problem_quadr = str2func(problem_name+"_quadr");
+        if ~isempty(CUT_MODE) && (CUT_MODE > 0)
+            CUT_MODE = 0;
+            warning('Cuts are not supported for SOL_MODE = 2. Changing to CUT_MODE = 0.');
+        end
+	end
+end
+if isempty(CUT_MODE)
+    CUT_MODE = 1;
+    disp('Selected cutting strategy using continuous cuts only (CUT_MODE == 1)');
+elseif CUT_MODE == 2
+    if ~is_quadratic
+        CUT_MODE = 1;
+        warning('No quadratic formulation available. Changing to CUT_MODE = 1.');
+    else
+        problem_quadr = str2func(problem_name+"_quadr");
+    end
+end
+
+% Initialization of L,U
+if isempty(L) || isempty(U)
+    startbox = infsup(lb, ub);
+    B = intval;
+    for i=1:p
+        B(i) = (1:p==i)*f(startbox);
+    end
+    if isempty(L)
+        L = B.inf';
+    end
+    if isempty(U)
+        U = B.sup';
+    end
+end
+L = L-OFFSET;
+U = U+OFFSET;
+
+%% Call MOMIBB
+if SOL_MODE == 1
+    tic;
+    [L,U,N,box_struct,it,exitflag] = MOMIBB(problem,problem_quadr,problem_param,CUT_MODE,L,U,EPSILON,OFFSET);
+    time = toc;
+elseif SOL_MODE == 2
+	tic;
+    [L,U,N,box_struct,it,exitflag] = MOMIBB_direct(problem,problem_quadr,problem_param,L,U,EPSILON,OFFSET);
+    time = toc; 
+end
+
+%% Plot if wanted
+if plot_result > 0
+    plotBoxes(L,U,p);
+    if p < 3
+        figure;
+        hold on;
+        plot(L(1,:),L(2,:),'LineStyle','none','Marker','.','Color','blue');
+        plot(U(1,:),U(2,:),'LineStyle','none','Marker','.','Color',[1, 0.4745, 0]);
+        grid on;
+        xlabel('f_1');
+        ylabel('f_2');
+    elseif p < 4
+        figure;
+        hold on;
+        plot3(L(1,:),L(2,:),L(3,:),'LineStyle','none','Marker','.','Color','blue');
+        plot3(U(1,:),U(2,:),U(3,:),'LineStyle','none','Marker','.','Color',[1, 0.4745, 0]);
+        grid on;
+        xlabel('f_1');
+        ylabel('f_2');
+        zlabel('f_3');
+    end
+end
+end

initSession.m

@@ -0,0 +1,18 @@
+%initSession initializes a session by loading all needed solvers and
+% paths
+
+%% Add paths for test problems and solvers
+addpath(genpath('problems'));
+addpath(genpath('solver'));
+
+%% Start subsolvers
+% Gurobi
+cd 'your_gurobi_path\matlab'
+gurobi_setup
+
+% Intlab
+cd 'your_intlab_path'
+startintlab
+
+%% Switch back to main path
+cd 'path_to_this_file'

problems/P2.m

@@ -0,0 +1,36 @@
+function [n,m,p,q,f,g,Df,Dg,Aineq,bineq,Aeq,beq,lb,ub,x0,is_convex,is_quadratic] = P2(params)
+%P2 A simple non convex example
+
+% Dimension of decision and criterion space
+n = params(1); % Continuous variables
+m = params(2); % Integer variables
+p = 2; % Dimension criterion space
+q = 2; % Number of constraints
+assert(mod(n,2)==0,'Number of continuous variables has to be even.')
+assert(mod(m,2)==0,'Number of integer variables has to be even.')
+
+% Problem type
+is_convex = false;
+is_quadratic = true;
+
+% Objective function
+f = @(x) [sum(x(1:n/2))+sum(x(n+1:n+m/2));sum(x(n/2+1:n))+sum(x(n+m/2+1:n+m))];
+Df = @(x) [ones(1,n/2),zeros(1,n/2),ones(1,m/2),zeros(1,m/2);zeros(1,n/2),ones(1,n/2),zeros(1,m/2),ones(1,m/2)];
+
+% Linear constraints (Aineq*x <= bineq, Aeq*x = beq)
+Aineq = [];
+bineq = [];
+Aeq = [];
+beq = [];
+
+% Lower and upper bounds (lb <= x <= ub)
+lb = [zeros(n,1);-3.*ones(m,1)];
+ub = 3.*ones(n+m,1);
+
+% Start point x0
+x0 = ceil((lb+ub)/2);
+
+% Non-linear constraints (g(x) <= 0)
+g = @(x) [-sum(x(1:n).^2)+1;sum(x(n+1:n+m).^2)-9];
+Dg = @(x) [-2.*x(1:n)',zeros(1,m);zeros(1,n),2.*x(n+1:n+m)'];
+end

problems/P2_quadr.m

@@ -0,0 +1,25 @@
+function [Qfun,qfun,cfun,Qcon,qcon,ccon] = P2_quadr(params)
+
+% Dimension of decision and criterion space
+n = params(1); % Continuous variables
+m = params(2); % Integer variables
+r = n+m;
+assert(mod(n,2)==0,'Number of continuous variables has to be even.')
+assert(mod(n,2)==0,'Number of continuous variables has to be even.')
+
+% Objective function
+Qfun{1} = zeros(r);
+qfun{1} = [ones(1,n/2),zeros(1,n/2),ones(1,m/2),zeros(1,m/2)]';
+cfun{1} = 0;
+Qfun{2} = zeros(r);
+qfun{2} = [zeros(1,n/2),ones(1,n/2),zeros(1,m/2),ones(1,m/2)]';
+cfun{2} = 0;
+
+% Constraints
+Qcon{1} = diag([-ones(1,n),zeros(1,m)]);
+qcon{1} = zeros(r,1);
+ccon{1} = -1;
+Qcon{2} = diag([zeros(1,n),ones(1,m)]);
+qcon{2} = zeros(r,1);
+ccon{2} = 9;
+end

problems/P3.m

@@ -0,0 +1,35 @@
+function [n,m,p,q,f,g,Df,Dg,Aineq,bineq,Aeq,beq,lb,ub,x0,is_convex,is_quadratic] = P3(params)
+%P3 A simple non convex example
+
+% Dimension of decision and criterion space
+n = params(1); % Continuous variables
+m = 1; % Integer variables
+p = 2; % Dimension criterion space
+q = 1; % Number of constraints
+assert(mod(n,2)==0,'Number of continuous variables has to be even.')
+
+% Problem type
+is_convex = false;
+is_quadratic = false;
+
+% Objective function
+f = @(x) [sum(x(1:n/2))+x(n+m);sum(x(n/2+1:n))-exp(x(n+m))];
+Df = @(x) [ones(1,n/2),zeros(1,n/2),1;zeros(1,n/2),ones(1,n/2),-exp(x(n+m))];
+
+% Linear constraints (Aineq*x <= bineq, Aeq*x = beq)
+Aineq = [];
+bineq = [];
+Aeq = [];
+beq = [];
+
+% Lower and upper bounds (lb <= x <= ub)
+lb = [zeros(n,1);-4.*ones(m,1)];
+ub = 1.*ones(n+m,1);
+
+% Start point x0
+x0 = ceil((lb+ub)/2);
+
+% Non-linear constraints (g(x) <= 0)
+g = @(x) [-sum(x(1:n).^2)+1];
+Dg = @(x) [-2.*x(1:n)',0];
+end

problems/P8.m

@@ -0,0 +1,34 @@
+function [n,m,p,q,f,g,Df,Dg,Aineq,bineq,Aeq,beq,lb,ub,x0,is_convex,is_quadratic] = P8(~)
+%P8 A non-convex non-quadratic tri-objective test instance
+
+% Dimension of decision and criterion space
+n = 3; % Continuous variables
+m = 1; % Integer variables
+p = 3; % Dimension criterion space
+q = 3; % Number of constraints
+
+% Problem type
+is_convex = false;
+is_quadratic = false;
+
+% Objective function
+f = @(x) [x(1)+x(4);x(2)-x(4);x(3)-exp(x(4))-3];
+Df = @(x) [];
+
+% Linear constraints (Aineq*x <= bineq, Aeq*x = beq)
+Aineq = [];
+bineq = [];
+Aeq = [];
+beq = [];
+
+% Lower and upper bounds (lb <= x <= ub)
+lb = -2.*ones(n+m,1);
+ub = 2.*ones(n+m,1);
+
+% Start point x0
+x0 = ceil((lb+ub)/2);
+
+% Non-linear constraints (g(x) <= 0)
+g = @(x) [x(1)^2+x(2)^2-1;x(1)*x(2)*(1-x(3))-1;exp(x(3))-1];
+Dg = @(x) [];
+end

problems/T4.m

@@ -0,0 +1,41 @@
+function [n,m,p,q,f,g,Df,Dg,Aineq,bineq,Aeq,beq,lb,ub,x0,is_convex,is_quadratic] = T4(params)
+%T4 A scalable test instance
+%   This example was taken from:
+%   Marianna De Santis, Gabriele Eichfelder, Julia Niebling, Stefan
+%   Rocktäschel. Solving Multiobjective Mixed Integer Convex Optimization
+%   Problems, SIAM Journal on Optimization (2020)
+
+% Dimension of decision and criterion space
+n = params(1); % Continuous variables
+m = params(2); % Integer variables
+p = 2; % Dimension criterion space
+q = 1; % Number of constraints
+assert(mod(n,2)==0,'Number of continuous variables has to be even.')
+
+% Problem type
+is_convex = true;
+is_quadratic = true;
+
+% Objective function
+f = @(x) [sum(x(1:(n/2)))+sum(x(n+1:n+m));sum(x((n/2+1):n))-sum(x(n+1:n+m))];
+Df = @(x) [[ones(1,n/2),zeros(1,n/2),ones(1,m)];[zeros(1,n/2),ones(1,n/2),-ones(1,m)]];
+
+% Linear constraints (Aineq*x <= bineq, Aeq*x = beq)
+Aineq = [];
+bineq = [];
+Aeq = [];
+beq = [];
+
+% Lower and upper bounds (lb <= x <= ub)
+lb = -2.*ones(n+m,1);
+ub = 2.*ones(n+m,1);
+% z = [-2*m-n/2;-2*m-n/2];
+% Z = [2*m+n/2;2*m+n/2];
+
+% Start point x0
+x0 = ceil((lb+ub)/2);
+
+% Non-linear constraints (g(x) <= 0)
+g = @(x) [sum(x(1:n).^2)-1];
+Dg = @(x) [2.*x(1:n)',zeros(1,m)];
+end

problems/T4_quadr.m

@@ -0,0 +1,22 @@
+function [Qfun,qfun,cfun,Qcon,qcon,ccon] = T4_quadr(params)
+
+% Dimension of decision and criterion space
+n = params(1); % Continuous variables
+m = params(2); % Integer variables
+r = n+m;
+assert(mod(n,2)==0,'Number of continuous variables has to be even.')
+
+% Objective function
+Qfun{1} = zeros(r);
+qfun{1} = [ones(n/2,1); zeros(n/2,1); ones(m,1)];
+cfun{1} = 0;
+
+Qfun{2} = zeros(r);
+qfun{2} = [zeros(n/2,1); ones(n/2,1); -ones(m,1)];
+cfun{2} = 0;
+
+% Constraint
+Qcon{1} = diag([ones(1,n), zeros(1,m)]);
+qcon{1} = zeros(r,1);
+ccon{1} = 1;
+end