diff options
| -rw-r--r-- | macros/Checkdims.bin | bin | 0 -> 8188 bytes | |||
| -rw-r--r-- | macros/Checklhs.bin | bin | 0 -> 10312 bytes | |||
| -rw-r--r-- | macros/Checkrhs.bin | bin | 0 -> 12852 bytes | |||
| -rw-r--r-- | macros/Checktype.bin | bin | 0 -> 9156 bytes | |||
| -rw-r--r-- | macros/Checkvector.bin | bin | 0 -> 9696 bytes | |||
| -rw-r--r-- | macros/cbcintlinprog.bin | bin | 0 -> 40320 bytes | |||
| -rw-r--r-- | macros/cbcintlinprog.sci | 234 | ||||
| -rw-r--r-- | macros/cbcmatrixintlinprog.bin | bin | 0 -> 38928 bytes | |||
| -rw-r--r-- | macros/cbcmatrixintlinprog.sci | 303 | ||||
| -rw-r--r-- | macros/cbcmpsintlinprog.bin | bin | 0 -> 7528 bytes | |||
| -rw-r--r-- | macros/cbcmpsintlinprog.sci | 86 | ||||
| -rw-r--r-- | macros/ecos.bin | bin | 0 -> 67656 bytes | |||
| -rw-r--r-- | macros/ecos.sci | 453 | ||||
| -rw-r--r-- | macros/fgoalattain.bin | bin | 0 -> 76660 bytes | |||
| -rw-r--r-- | macros/fminbnd.bin | bin | 0 -> 64156 bytes | |||
| -rw-r--r-- | macros/fmincon.bin | bin | 0 -> 156136 bytes | |||
| -rw-r--r-- | macros/fminimax.bin | bin | 0 -> 86532 bytes | |||
| -rw-r--r-- | macros/fminunc.bin | bin | 0 -> 69644 bytes | |||
| -rw-r--r-- | macros/intfminbnd.bin | bin | 0 -> 52768 bytes | |||
| -rw-r--r-- | macros/intfminbnd.sci | 360 | ||||
| -rw-r--r-- | macros/intfmincon.bin | bin | 0 -> 84944 bytes | |||
| -rw-r--r-- | macros/intfmincon.sci | 589 | ||||
| -rw-r--r-- | macros/intfminimax.bin | bin | 0 -> 73232 bytes | |||
| -rw-r--r-- | macros/intfminimax.sci | 456 | ||||
| -rw-r--r-- | macros/intfminunc.bin | bin | 0 -> 50884 bytes | |||
| -rw-r--r-- | macros/intfminunc.sci | 350 | ||||
| -rw-r--r-- | macros/intqpipopt.bin | bin | 0 -> 63384 bytes | |||
| -rw-r--r-- | macros/intqpipopt.sci | 461 | ||||
| -rw-r--r-- | macros/lib | bin | 0 -> 1120 bytes | |||
| -rw-r--r-- | macros/linprog.bin | bin | 0 -> 30012 bytes | |||
| -rw-r--r-- | macros/lsqlin.bin | bin | 0 -> 64528 bytes | |||
| -rw-r--r-- | macros/lsqnonlin.bin | bin | 0 -> 48680 bytes | |||
| -rw-r--r-- | macros/lsqnonneg.bin | bin | 0 -> 33808 bytes | |||
| -rw-r--r-- | macros/matrix_linprog.bin | bin | 0 -> 30584 bytes | |||
| -rw-r--r-- | macros/mps_linprog.bin | bin | 0 -> 10064 bytes | |||
| -rw-r--r-- | macros/names | 31 | ||||
| -rw-r--r-- | macros/qpipopt.bin | bin | 0 -> 63432 bytes | |||
| -rw-r--r-- | macros/qpipoptmat.bin | bin | 0 -> 65684 bytes | |||
| -rw-r--r-- | macros/setOptions.bin | bin | 0 -> 3040 bytes | |||
| -rw-r--r-- | macros/symphony.bin | bin | 0 -> 62048 bytes | |||
| -rw-r--r-- | macros/symphony_call.bin | bin | 0 -> 4580 bytes | |||
| -rw-r--r-- | macros/symphonymat.bin | bin | 0 -> 65440 bytes |
42 files changed, 3323 insertions, 0 deletions
diff --git a/macros/Checkdims.bin b/macros/Checkdims.bin Binary files differnew file mode 100644 index 0000000..40e385a --- /dev/null +++ b/macros/Checkdims.bin diff --git a/macros/Checklhs.bin b/macros/Checklhs.bin Binary files differnew file mode 100644 index 0000000..7156107 --- /dev/null +++ b/macros/Checklhs.bin diff --git a/macros/Checkrhs.bin b/macros/Checkrhs.bin Binary files differnew file mode 100644 index 0000000..2c45876 --- /dev/null +++ b/macros/Checkrhs.bin diff --git a/macros/Checktype.bin b/macros/Checktype.bin Binary files differnew file mode 100644 index 0000000..bbe1585 --- /dev/null +++ b/macros/Checktype.bin diff --git a/macros/Checkvector.bin b/macros/Checkvector.bin Binary files differnew file mode 100644 index 0000000..dfe03ab --- /dev/null +++ b/macros/Checkvector.bin diff --git a/macros/cbcintlinprog.bin b/macros/cbcintlinprog.bin Binary files differnew file mode 100644 index 0000000..fd6f801 --- /dev/null +++ b/macros/cbcintlinprog.bin diff --git a/macros/cbcintlinprog.sci b/macros/cbcintlinprog.sci new file mode 100644 index 0000000..00e129f --- /dev/null +++ b/macros/cbcintlinprog.sci @@ -0,0 +1,234 @@ +// Copyright (C) 2016 - IIT Bombay - FOSSEE +// +// This file must be used under the terms of the CeCILL. +// This source file is licensed as described in the file COPYING, which +// you should have received as part of this distribution. The terms +// are also available at +// http://www.cecill.info/licences/Licence_CeCILL_V2-en.txt +// Author: Pranav Deshpande and Akshay Miterani +// Organization: FOSSEE, IIT Bombay +// Email: toolbox@scilab.in + +function [xopt,fopt,exitflag,output] = cbcintlinprog (varargin) + // Solves a mixed integer linear programming constrained optimization problem in intlinprog format. + // + // Calling Sequence + // xopt = cbcintlinprog(c,intcon,A,b) + // xopt = cbcintlinprog(c,intcon,A,b,Aeq,beq) + // xopt = cbcintlinprog(c,intcon,A,b,Aeq,beq,lb,ub) + // xopt = cbcintlinprog(c,intcon,A,b,Aeq,beq,lb,ub,options) + // xopt = cbcintlinprog('path_to_mps_file') + // xopt = cbcintlinprog('path_to_mps_file',options) + // [xopt,fopt,status,output] = cbcintlinprog( ... ) + // + // Parameters + // c : a vector of double, contains coefficients of the variables in the objective + // intcon : Vector of integer constraints, specified as a vector of positive integers. The values in intcon indicate the // components of the decision variable x that are integer-valued. intcon has values from 1 through number of variable. + // A : a matrix of double, represents the linear coefficients in the inequality constraints A⋅x ≤ b. + // b : a vector of double, represents the linear coefficients in the inequality constraints A⋅x ≤ b. + // Aeq : a matrix of double, represents the linear coefficients in the equality constraints Aeq⋅x = beq. + // beq : a vector of double, represents the linear coefficients in the equality constraints Aeq⋅x = beq. + // lb : Lower bounds, specified as a vector or array of double. lb represents the lower bounds elementwise in lb ≤ x ≤ ub. + // ub : Upper bounds, specified as a vector or array of double. ub represents the upper bounds elementwise in lb ≤ x ≤ ub. + // options : a list containing the parameters to be set. + // xopt : a vector of double, the computed solution of the optimization problem. + // fopt : a double, the value of the function at x. + // status : status flag returned from symphony. See below for details. + // output : The output data structure contains detailed information about the optimization process. See below for details. + // + // Description + // Search the minimum or maximum of a constrained mixed integer linear programming optimization problem specified by : + // + // <latex> + // \begin{eqnarray} + // &\mbox{min}_{x} + // & C^T⋅x \\ + // & \text{subject to} & A⋅x \leq b \\ + // & & Aeq⋅x = beq \\ + // & & lb \leq x \leq ub \\ + // & & x_i \in \!\, \mathbb{Z}, i \in \!\, intcon\\ + // \end{eqnarray} + // </latex> + // + // Examples + // // Objective function + // // Reference: Westerberg, Carl-Henrik, Bengt Bjorklund, and Eskil Hultman. "An application of mixed integer programming in a Swedish steel mill." Interfaces 7, no. 2 (1977): 39-43. + // c = [350*5,330*3,310*4,280*6,500,450,400,100]'; + // // Lower Bound of variable + // lb = repmat(0,1,8); + // // Upper Bound of variables + // ub = [repmat(1,1,4) repmat(%inf,1,4)]; + // // Constraint Matrix + // Aeq = [5,3,4,6,1,1,1,1; + // 5*0.05,3*0.04,4*0.05,6*0.03,0.08,0.07,0.06,0.03; + // 5*0.03,3*0.03,4*0.04,6*0.04,0.06,0.07,0.08,0.09;] + // beq = [ 25, 1.25, 1.25] + // intcon = [1 2 3 4]; + // // Calling Symphony + // [x,f,status,output] = cbcintlinprog(c,intcon,[],[],Aeq,beq,lb,ub) + // // Press ENTER to continue + // + // Examples + // // An advanced case where we set some options in symphony + // // This problem is taken from + // // P.C.Chu and J.E.Beasley + // // "A genetic algorithm for the multidimensional knapsack problem", + // // Journal of Heuristics, vol. 4, 1998, pp63-86. + // // The problem to be solved is: + // // Max sum{j=1,...,n} p(j)x(j) + // // st sum{j=1,...,n} r(i,j)x(j) <= b(i) i=1,...,m + // // x(j)=0 or 1 + // // The function to be maximize i.e. P(j) + // c = -1*[ 504 803 667 1103 834 585 811 856 690 832 846 813 868 793 .. + // 825 1002 860 615 540 797 616 660 707 866 647 746 1006 608 .. + // 877 900 573 788 484 853 942 630 591 630 640 1169 932 1034 .. + // 957 798 669 625 467 1051 552 717 654 388 559 555 1104 783 .. + // 959 668 507 855 986 831 821 825 868 852 832 828 799 686 .. + // 510 671 575 740 510 675 996 636 826 1022 1140 654 909 799 .. + // 1162 653 814 625 599 476 767 954 906 904 649 873 565 853 1008 632]'; + // // Constraint Matrix + // A = [ //Constraint 1 + // 42 41 523 215 819 551 69 193 582 375 367 478 162 898 .. + // 550 553 298 577 493 183 260 224 852 394 958 282 402 604 .. + // 164 308 218 61 273 772 191 117 276 877 415 873 902 465 .. + // 320 870 244 781 86 622 665 155 680 101 665 227 597 354 .. + // 597 79 162 998 849 136 112 751 735 884 71 449 266 420 .. + // 797 945 746 46 44 545 882 72 383 714 987 183 731 301 .. + // 718 91 109 567 708 507 983 808 766 615 554 282 995 946 651 298; + // //Constraint 2 + // 509 883 229 569 706 639 114 727 491 481 681 948 687 941 .. + // 350 253 573 40 124 384 660 951 739 329 146 593 658 816 .. + // 638 717 779 289 430 851 937 289 159 260 930 248 656 833 .. + // 892 60 278 741 297 967 86 249 354 614 836 290 893 857 .. + // 158 869 206 504 799 758 431 580 780 788 583 641 32 653 .. + // 252 709 129 368 440 314 287 854 460 594 512 239 719 751 .. + // 708 670 269 832 137 356 960 651 398 893 407 477 552 805 881 850; + // //Constraint 3 + // 806 361 199 781 596 669 957 358 259 888 319 751 275 177 .. + // 883 749 229 265 282 694 819 77 190 551 140 442 867 283 .. + // 137 359 445 58 440 192 485 744 844 969 50 833 57 877 .. + // 482 732 968 113 486 710 439 747 174 260 877 474 841 422 .. + // 280 684 330 910 791 322 404 403 519 148 948 414 894 147 .. + // 73 297 97 651 380 67 582 973 143 732 624 518 847 113 .. + // 382 97 905 398 859 4 142 110 11 213 398 173 106 331 254 447 ; + // //Constraint 4 + // 404 197 817 1000 44 307 39 659 46 334 448 599 931 776 .. + // 263 980 807 378 278 841 700 210 542 636 388 129 203 110 .. + // 817 502 657 804 662 989 585 645 113 436 610 948 919 115 .. + // 967 13 445 449 740 592 327 167 368 335 179 909 825 614 .. + // 987 350 179 415 821 525 774 283 427 275 659 392 73 896 .. + // 68 982 697 421 246 672 649 731 191 514 983 886 95 846 .. + // 689 206 417 14 735 267 822 977 302 687 118 990 323 993 525 322; + // //Constraint 5 + // 475 36 287 577 45 700 803 654 196 844 657 387 518 143 .. + // 515 335 942 701 332 803 265 922 908 139 995 845 487 100 .. + // 447 653 649 738 424 475 425 926 795 47 136 801 904 740 .. + // 768 460 76 660 500 915 897 25 716 557 72 696 653 933 .. + // 420 582 810 861 758 647 237 631 271 91 75 756 409 440 .. + // 483 336 765 637 981 980 202 35 594 689 602 76 767 693 .. + // 893 160 785 311 417 748 375 362 617 553 474 915 457 261 350 635 ; + // ]; + // nbVar = size(c,1); + // b=[11927 13727 11551 13056 13460 ]; + // // Lower Bound of variables + // lb = repmat(0,1,nbVar); + // // Upper Bound of variables + // ub = repmat(1,1,nbVar); + // // Lower Bound of constrains + // intcon = []; + // for i = 1:nbVar + // intcon = [intcon i]; + // end + // options = list('MaxTime', 25); + // // The expected solution : + // // Output variables + // xopt = [0 1 1 0 0 1 0 1 0 1 0 0 0 0 0 0 0 1 0 0 0 0 1 0 1 1 0 1 1 0 1 .. + // 0 0 0 0 0 0 0 0 0 0 0 1 0 0 0 0 0 1 0 0 0 0 0 0 1 0 0 0 0 1 1 .. + // 0 0 1 0 0 1 0 1 0 0 1 0 0 1 0 1 0 0 0 0 0 1 1 0 0 0 0 0 1 1 0 0 1 0 0 1 0]; + // // Optimal value + // fopt = [ 24381 ] + // // Calling cbc + // [x,f,status,output] = cbcintlinprog(c,intcon,A,b,[],[],lb,ub,options); + // Authors + // Akshay Miterani and Pranav Deshpande + + if(type(varargin(1))==1) then + + //To check the number of input and output argument + [lhs , rhs] = argn(); + + //To check the number of argument given by user + if ( rhs < 4 | rhs == 5 | rhs == 7 | rhs > 9 ) then + errmsg = msprintf(gettext("%s: Unexpected number of input arguments : %d provided while should be in the set [4 6 8 9]"), "cbcintlinprog", rhs); + error(errmsg); + end + + c = []; + intcon = []; + A = []; + b = []; + Aeq = []; + beq = []; + lb = []; + ub = []; + options = list(); + + c = varargin(1) + intcon = varargin(2) + A = varargin(3) + b = varargin(4) + + if(size(c,2) == 0) then + errmsg = msprintf(gettext("%s: Cannot determine the number of variables because input objective coefficients is empty"),"cbcintlinprog"); + error(errmsg); + end + + if (size(c,2)~=1) then + errmsg = msprintf(gettext("%s: Objective Coefficients should be a column matrix"), "cbcintlinprog"); + error(errmsg); + end + + nbVar = size(c,1); + + if ( rhs<5 ) then + Aeq = [] + beq = [] + else + Aeq = varargin(5); + beq = varargin(6); + end + + if ( rhs<7 ) then + lb = repmat(-%inf,1,nbVar); + ub = repmat(%inf,1,nbVar); + else + lb = varargin(7); + ub = varargin(8); + end + + if (rhs<9|size(varargin(9))==0) then + options = list(); + else + options = varargin(9); + end + [xopt,fopt,exitflag,output]=cbcmatrixintlinprog(c,intcon,A,b,Aeq,beq,lb,ub,options); + elseif(type(varargin(1))==10) then + + [lhs , rhs] = argn(); + + //To check the number of argument given by user + if ( rhs < 1 | rhs > 2) then + errmsg = msprintf(gettext("%s: Unexpected number of input arguments : %d provided while should be in the set of [1 2]"),"cbcintlinprog",rhs); + error(errmsg) + end + mpsFile = varargin(1); + if ( rhs<2 | size(varargin(2)) ==0 ) then + param = list(); + else + param =varargin(2); + end + + [xopt,fopt,exitflag,output]=cbcmpsintlinprog(mpsFile,param); + end + +endfunction diff --git a/macros/cbcmatrixintlinprog.bin b/macros/cbcmatrixintlinprog.bin Binary files differnew file mode 100644 index 0000000..189ce16 --- /dev/null +++ b/macros/cbcmatrixintlinprog.bin diff --git a/macros/cbcmatrixintlinprog.sci b/macros/cbcmatrixintlinprog.sci new file mode 100644 index 0000000..b37f360 --- /dev/null +++ b/macros/cbcmatrixintlinprog.sci @@ -0,0 +1,303 @@ +// Copyright (C) 2016 - IIT Bombay - FOSSEE +// +// This file must be used under the terms of the CeCILL. +// This source file is licensed as described in the file COPYING, which +// you should have received as part of this distribution. The terms +// are also available at +// http://www.cecill.info/licences/Licence_CeCILL_V2-en.txt +// Author: Pranav Deshpande and Akshay Miterani +// Organization: FOSSEE, IIT Bombay +// Email: toolbox@scilab.in + +function [xopt,fopt,status,output] = cbcmatrixintlinprog (varargin) + // Sci file wrapper for the matrix_cbcintlinprog.cpp file + + // To check the number of input and output argument + [lhs , rhs] = argn(); + + // To check the number of argument given by user + if ( rhs < 4 | rhs == 5 | rhs == 7 | rhs > 9 ) then + errmsg = msprintf(gettext("%s: Unexpected number of input arguments : %d provided while should be in the set [4 6 8 9]"), "cbcintlinprog", rhs); + error(errmsg); + end + + c = []; + intcon = []; + A = []; + b = []; + Aeq = []; + beq = []; + lb = []; + ub = []; + options = list(); + + c = varargin(1) + intcon = varargin(2) + A = varargin(3) + b = varargin(4) + + if(size(c,2) == 0) then + errmsg = msprintf(gettext("%s: Cannot determine the number of variables because input objective coefficients is empty"),"cbcintlinprog"); + error(errmsg); + end + + if (size(c,2)~=1) then + errmsg = msprintf(gettext("%s: Objective Coefficients should be a column matrix"), "cbcintlinprog"); + error(errmsg); + end + + nbVar = size(c,1); + + if ( rhs<5 ) then + Aeq = [] + beq = [] + else + Aeq = varargin(5); + beq = varargin(6); + end + + if ( rhs<7 ) then + lb = repmat(-%inf,1,nbVar); + ub = repmat(%inf,1,nbVar); + else + lb = varargin(7); + ub = varargin(8); + end + + if (rhs<9|size(varargin(9))==0) then + options = list(); + else + options = varargin(9); + end + + //Check type of variables + Checktype("cbcintlinprog", c, "c", 1, "constant") + Checktype("cbcintlinprog", intcon, "intcon", 2, "constant") + Checktype("cbcintlinprog", A, "A", 3, "constant") + Checktype("cbcintlinprog", b, "b", 4, "constant") + Checktype("cbcintlinprog", Aeq, "Aeq", 5, "constant") + Checktype("cbcintlinprog", beq, "beq", 6, "constant") + Checktype("cbcintlinprog", lb, "lb", 7, "constant") + Checktype("cbcintlinprog", ub, "ub", 8, "constant") + + // Check if the user gives empty matrix + if (size(lb,2)==0) then + lb = repmat(-%inf,nbVar,1); + end + + if (size(intcon,2)==0) then + intcon = []; + end + + if (size(ub,2)==0) then + ub = repmat(%inf,nbVar,1); + end + + // Calculating the size of equality and inequality constraints + nbConInEq = size(A,1); + nbConEq = size(Aeq,1); + + // Check if the user gives row vector + // and Changing it to a column matrix + + if (size(lb,2)== [nbVar]) then + lb = lb'; + end + + if (size(ub,2)== [nbVar]) then + ub = ub'; + end + + if (size(b,2)== [nbConInEq]) then + b = b'; + end + + if (size(beq,2)== [nbConEq]) then + beq = beq'; + end + + for i=1:size(intcon,2) + if(intcon(i)>nbVar) then + errmsg = msprintf(gettext("%s: The values inside intcon should be less than the number of variables"), "cbcintlinprog"); + error(errmsg); + end + + if (intcon(i)<0) then + errmsg = msprintf(gettext("%s: The values inside intcon should be greater than 0 "), "cbcintlinprog"); + error(errmsg); + end + + if(modulo(intcon(i),1)) then + errmsg = msprintf(gettext("%s: The values inside intcon should be an integer "), "cbcintlinprog"); + error(errmsg); + end + end + + //Check the size of inequality constraint which should equal to the number of inequality constraints + if ( size(A,2) ~= nbVar & size(A,2) ~= 0) then + errmsg = msprintf(gettext("%s: The size of inequality constraint is not equal to the number of variables"), "cbcintlinprog"); + error(errmsg); + end + + + //Check the size of lower bound of inequality constraint which should equal to the number of constraints + if ( size(b,1) ~= size(A,1)) then + errmsg = msprintf(gettext("%s: The Lower Bound of inequality constraint is not equal to the number of constraint"), "cbcintlinprog"); + error(errmsg); + end + + //Check the size of equality constraint which should equal to the number of inequality constraints + if ( size(Aeq,2) ~= nbVar & size(Aeq,2) ~= 0) then + errmsg = msprintf(gettext("%s: The size of equality constraint is not equal to the number of variables"), "cbcintlinprog"); + error(errmsg); + end + + //Check the size of upper bound of equality constraint which should equal to the number of constraints + if ( size(beq,1) ~= size(Aeq,1)) then + errmsg = msprintf(gettext("%s: The equality constraint upper bound is not equal to the number of equality constraint"), "cbcintlinprog"); + error(errmsg); + end + + //Check the size of Lower Bound which should equal to the number of variables + if ( size(lb,1) ~= nbVar) then + errmsg = msprintf(gettext("%s: The Lower Bound is not equal to the number of variables"), "cbcintlinprog"); + error(errmsg); + end + + //Check the size of Upper Bound which should equal to the number of variables + if ( size(ub,1) ~= nbVar) then + errmsg = msprintf(gettext("%s: The Upper Bound is not equal to the number of variables"), "cbcintlinprog"); + error(errmsg); + end + + if (type(options) ~= 15) then + errmsg = msprintf(gettext("%s: Options should be a list "), "cbcintlinprog"); + error(errmsg); + end + + + if (modulo(size(options),2)) then + errmsg = msprintf(gettext("%s: Size of parameters should be even"), "cbcintlinprog"); + error(errmsg); + end + + //Check if the user gives a matrix instead of a vector + + if (((size(intcon,1)~=1)& (size(intcon,2)~=1))&(size(intcon,2)~=0)) then + errmsg = msprintf(gettext("%s: intcon should be a vector"), "cbcintlinprog"); + error(errmsg); + end + + if (size(lb,1)~=1)& (size(lb,2)~=1) then + errmsg = msprintf(gettext("%s: Lower Bound should be a vector"), "cbcintlinprog"); + error(errmsg); + end + + if (size(ub,1)~=1)& (size(ub,2)~=1) then + errmsg = msprintf(gettext("%s: Upper Bound should be a vector"), "cbcintlinprog"); + error(errmsg); + end + + if (nbConInEq) then + if ((size(b,1)~=1)& (size(b,2)~=1)) then + errmsg = msprintf(gettext("%s: Constraint Lower Bound should be a vector"), "cbcintlinprog"); + error(errmsg); + end + end + + if (nbConEq) then + if (size(beq,1)~=1)& (size(beq,2)~=1) then + errmsg = msprintf(gettext("%s: Constraint Upper Bound should be a vector"), "cbcintlinprog"); + error(errmsg); + end + end + + + //Changing the inputs in symphony's format + conMatrix = [A;Aeq] + nbCon = size(conMatrix,1); + conLB = [repmat(-%inf,size(A,1),1);beq]; + conUB = [b;beq] ; + + isInt = repmat(%f,1,nbVar); + // Changing intcon into column vector + intcon = intcon(:); + for i=1:size(intcon,1) + isInt(intcon(i)) = %t + end + + objSense = 1.0; + + //Changing into row vector + lb = lb'; + ub = ub'; + c = c'; + + //Pusing options as required to a double array + optval = []; + if length(options) == 0 then + optval = [0 0 0 0]; + else + optval = [0 0 0 0]; + for i=1:2:length(options) + select options(i) + case 'IntegerTolerance' then + optval(1) = options(i+1); + case 'MaxNodes' then + optval(2) = options(i+1); + case 'MaxTime' then + optval(3) = options(i+1); + case 'AllowableGap' then + optval(4) = options(i+1); + else + error(999, 'Unknown string argument passed.'); + end + end + end + + [xopt,fopt,status,nodes,nfpoints,L,U,niter] = sci_matrix_intlinprog(nbVar,nbCon,c,intcon,conMatrix,conLB,conUB,lb,ub,objSense,optval); + + //Debugging Prints + //disp(c);disp(intcon);disp(conMatrix);disp(conLB);disp(conUB);disp(lb);disp(ub);disp(Aeq);disp(beq);disp(xopt); + //disp(L);disp(U); + //disp(options); + + output = struct("relativegap" , [],.. + "absolutegap" , [],.. + "numnodes" , [],.. + "numfeaspoints" , [],.. + "numiterations" , [],.. + "constrviolation" , [],.. + "message" , ''); + + output.numnodes=[nodes]; + output.numfeaspoints=[nfpoints]; + output.numiterations=[niter]; + output.relativegap=(U-L)/(abs(U)+1); + output.absolutegap=(U-L); + output.constrviolation = max([0;norm(Aeq*xopt-beq, 'inf');(lb'-xopt);(xopt-ub');(A*xopt-b)]); + + select status + + case 0 then + output.message="Optimal Solution" + case 1 then + output.message="Primal Infeasible" + case 2 then + output.message="Solution Limit is reached" + case 3 then + output.message="Node Limit is reached" + case 4 then + output.message="Numerical Difficulties" + case 5 then + output.message="Time Limit Reached" + case 6 then + output.message="Continuous Solution Unbounded" + case 7 then + output.message="Dual Infeasible" + else + output.message="Invalid status returned. Notify the Toolbox authors" + break; + end + +endfunction diff --git a/macros/cbcmpsintlinprog.bin b/macros/cbcmpsintlinprog.bin Binary files differnew file mode 100644 index 0000000..ee042a0 --- /dev/null +++ b/macros/cbcmpsintlinprog.bin diff --git a/macros/cbcmpsintlinprog.sci b/macros/cbcmpsintlinprog.sci new file mode 100644 index 0000000..1d46b5a --- /dev/null +++ b/macros/cbcmpsintlinprog.sci @@ -0,0 +1,86 @@ +// Copyright (C) 2016 - IIT Bombay - FOSSEE +// +// This file must be used under the terms of the CeCILL. +// This source file is licensed as described in the file COPYING, which +// you should have received as part of this distribution. The terms +// are also available at +// http://www.cecill.info/licences/Licence_CeCILL_V2-en.txt +// Author: Pranav Deshpande and Akshay Miterani +// Organization: FOSSEE, IIT Bombay +// Email: toolbox@scilab.in + +function [xopt,fopt,status,output] = cbcmpsintlinprog(varargin) + // Sci File Wrapper for the mps_cbcintlinprog.cpp file + + // Number of input and output arguments + [nOutput, nInput] = argn(); + + // To check the number of arguments given by the user + if (nInput<1 | nInput>2) then + error(999, 'Check the number of input arguments!'); + end + + mpsFile = varargin(1); + optval = [0,0,0,0]; + if(nInput==2) then + options=varargin(2); + if length(options) == 0 then + optval = [0 0 0 0]; + else + optval = [0 0 0 0]; + for i=1:2:length(options) + select options(i) + case 'IntegerTolerance' then + optval(1) = options(i+1); + case 'MaxNodes' then + optval(2) = options(i+1); + case 'MaxTime' then + optval(3) = options(i+1); + case 'AllowableGap' then + optval(4) = options(i+1); + else + error(999, 'Unknown string argument passed.'); + end + end + end + end + + [xopt,fopt,status,nodes,nfpoints,L,U,niter] = sci_mps_intlinprog(mpsFile, optval) + + output = struct("relativegap" , [],.. + "absolutegap" , [],.. + "numnodes" , [],.. + "numfeaspoints" , [],.. + "numiterations" , [],.. + "message" , ''); + + output.numnodes = [nodes]; + output.numfeaspoints = [nfpoints]; + output.numiterations = [niter]; + output.relativegap = (U-L)/(abs(U)+1); + output.absolutegap = (U-L); + + select status + + case 0 then + output.message="Optimal Solution" + case 1 then + output.message="Primal Infeasible" + case 2 then + output.message="Solution Limit is reached" + case 3 then + output.message="Node Limit is reached" + case 4 then + output.message="Numerical Difficulties" + case 5 then + output.message="Time Limit Reached" + case 6 then + output.message="Continuous Solution Unbounded" + case 7 then + output.message="Dual Infeasible" + else + output.message="Invalid status returned. Notify the Toolbox authors" + break; + end + +endfunction diff --git a/macros/ecos.bin b/macros/ecos.bin Binary files differnew file mode 100644 index 0000000..85c0c2f --- /dev/null +++ b/macros/ecos.bin diff --git a/macros/ecos.sci b/macros/ecos.sci new file mode 100644 index 0000000..225a986 --- /dev/null +++ b/macros/ecos.sci @@ -0,0 +1,453 @@ +function [x,y,s,z,info,status] = ecos(varargin) + // Solves conic optimization problems with integer and boolean constraints. + // + // Calling Sequence + // x = ecos(c,G,h,dims) + // x = ecos(c,G,h,dims,A,b) + // x = ecos(c,G,h,dims,A,b,options) + // [x,y,s,z,info,status] = ecos( ... ) + // + // Parameters + // c : a vector of double, contains coefficients of the variables in the objective. + // G : a matrix of double, represents the linear coefficients in the inequality constraints G⋅x ≤ h with respect to cone K. + // h : a vector of double, represents the linear coefficients in the inequality constraints G⋅x ≤ h with respect to cone K. + // dims : a list containing the number of n-dimensional positive orthant(l),the second order cones(q) and exponential cones(e). + // A : a matrix of double, represents the linear coefficients in the equality constraints A⋅x = b. + // b : a vector of double, represents the linear coefficients in the equality constraints A⋅x = b. + // options : a list containing the parameters to be set. + // x : a vector of double, the primal solution variable of the optimization problem. + // y : a vector of double, the dual variables fo equality constraints of the optimization problem. + // s : a vector of double, the slack for inequality constraints of the optimization problem. + // z : a vector of double, the dula varibale for inequality constraints of the optimization problem. + // info : a struct defining different parameters of the ecos solver. See below for details. + // status : The status returns the solver status after the optimization process. See below for details. + // + // Description + // Search the minimum of a conic constrained mixed integer programming optimization problem specified by : + // + // <latex> + // \begin{eqnarray} + // &\mbox{min}_{x} + // & C^T⋅x \\ + // & \text{subject to}& A⋅x = b \\ + // & & G⋅x \preceq_K h \\ + // & & x_i \in \{0,1\}, i \in \!\, bool_vars\ + // & & x_j \in \!\, \mathbb{Z}, i \in \!\, int_vars\\ + // \end{eqnarray} + // </latex> + // + // Examples + // // Objective function + // c = [-750 -1000]; + // // Constraint inequality matrix LHS + // G = [ + // 1 1; + // 1 2; + // 4 3; + // ]; + // // Constraint inequality matrix RHS + // h = [10 15 25]'; + // // Constraint equality matrix LHS + // A = [ + // 0.5 1 + // ]; + // // Constraint equality matrix RHS + // b=[7.5]; + // // Dimension of positive orthant + // l = [3]; + // q = []; + // e = []; + // dims=list("l",l,"q",q,"e",e) + // //Calling ecos + // [x,y,s,z,info,status] =ecos(c,G,h,dims,A,b); + + // Examples + // + // // Objective function + // c = [0 0 0 0 1]; + // //Constraint inequality matrix LHS + // G = [ + // 0.4167578 0.0562668 0. 0. 0. + // 2.1361961 -1.6402708 0. 0. 0. + // 1.7934356 0.8417474 0. 0. 0. + // 0. 0. 0.4167578 0.0562668 0. + // 0. 0. 2.1361961 -1.6402708 0. + // 0. 0. 1.7934356 0.8417474 0. + // 0. 0. 0. 0. -1. + // -1. 0. 0. 0. 0. + // 0. -1. 0. 0. 0. + // 0. 0. -1. 0. 0. + // 0. 0. 0. -1. 0. + // ]; + // //Constraint inequality matrix RHS + // h = [0 0 0 0 0 0 0 0 0 0 0]'; + // // Dimension of positive orthant + // l = [6]; + // q = [5]; + // e = [0] + // dims=list("l",l,"q",q,"e",e) + // [x,y,s,z,info,status] =ecos(c,G,h,dims); + // Author + // Georgey John + + function [A1,b1,s0]=linconcheck(A,b,inputs_name) + //Function to check the linear inputs A,b and Aeq,beq + + if(size(b,2)>1) then + errmsg = msprintf(gettext("%s: Expected Column vector of size (Number of constraints) for %s"), "ecos",inputs_name(2)); + error(errmsg); + end + + s0=size(A); + //To check for correct size of A + if(s0(2)==0) then + if(size(b,2)~=0) then + errmsg = msprintf(gettext("%s: As Linear Inequality Constraint coefficient Matrix %s is empty, %s should also be empty"), "ecos",inputs_name(1),inputs_name(2)); + error(errmsg); + end + else + if((size(b,1)~=1) & (size(b,2)~=1)) then + errmsg = msprintf(gettext("%s: Expected Non empty Row/Column Vector for %s for your Inputs "), "ecos",inputs_name(2)); + error(errmsg); + elseif(size(b,1)~=s0(1) & size(b,2)==1) then + errmsg = msprintf(gettext("%s: Expected Column Vector (number of linear inequality constraints X 1) for %s"), "ecos",inputs_name(2)); + error(errmsg); + // elseif(size(b,1)==1 & size(b,2)~=s0(1)) then + // errmsg = msprintf(gettext("%s: Expected Row Vector (1 X number of linear inequality constraints) for %s"), "ecos",inputs_name(2)); + // error(errmsg); + end + end + b1=b; + + //To check for corrcet size of A + if(size(A,1)~=size(b,1) & size(A,2)~=0) then + errmsg = msprintf(gettext("%s: Expected Matrix of size (No of linear inequality constraints X No of Variables) or an Empty Matrix for Linear Inequality Constraint coefficient Matrix %s"), "ecos",inputs_name(1)); + error(errmsg); + end + A1=A; + + endfunction + + [lhs ,rhs]=argn() + + if ( rhs<4 | rhs>9 ) then + errmsg = msprintf(gettext("%s: Unexpected number of input arguments : %d provided while it should be 4,6,7,8 or 9"), "ecos", rhs); + error(errmsg); + end + + if (rhs==5) then + errmsg = msprintf(gettext("%s: Unexpected number of input arguments : %d provided while it should be 4,6,8,10 or 11"), "ecos", rhs); + error(errmsg); + end + + c = varargin(1); + G = varargin(2); + h = varargin(3); + dims = varargin(4); + + A=[];Apr=[];Ajc=[];Air=[]; + b=[]; + param=list(); + m=0,n=0,p=0; + + if (rhs>4 & rhs<7) then + A = varargin(5); + b = varargin(6); + else + A=[]; + b=[]; + end + + if (rhs ==7) then + param = varargin(7); + end + + Checktype("ecos", c, "c", 1, "constant"); + + // Error check for objective matrix(c) + if (size(c,1)~=1) then + errmsg = msprintf(gettext("%s: c should be a row vector"), "ecos"); + error(errmsg); + elseif (size(c,1)==0) then + errmsg = msprintf(gettext("%s: c cannot be empty"), "ecos"); + error(errmsg); + end + + s=size(c,2); + + + Checktype("ecos", G, "G", 2, "constant"); + Checktype("ecos", h, "h", 3, "constant"); + + [m,n]=size(G); + + // Error check and converting inequaity matrix(G) to CCS format required by ecos + if (n~=size(c,2)) then + errmsg = msprintf(gettext("%s: Number of columns of G and c do not match"), "ecos"); + error(errmsg); + end + + [G,h,s1]=linconcheck(G,h,["G","h"]); + [Gjc,Gir,Gpr] = sp2adj(sparse(G)); + + // Extracting values and Error checks for dims + Checktype("ecos", dims, "dims", 4, "list"); + + if (isempty(dims) | modulo(size(dims),2)) then + errmsg = msprintf(gettext("%s: dims cannot be empty and should be list of even size"), "ecos"); + error(errmsg); + end + + for i = 1:(size(dims))/2 + select convstr(dims(2*i-1),'l') + case "l" then + Checktype("ecos", dims(2*i), "l", 4, "constant"); + if (isempty(dims(2*i))) then + l=0; + elseif (dims(2*i)<0 | modulo(dims(2*i),1)) then + errmsg = msprintf(gettext("%s: l in dims should be a positive integer"), "ecos"); + error(errmsg); + else + l=dims(2*i); + // if (prod(size(l))>1)) then + // errmsg = msprintf(gettext("%s: l in dims should be a positive integer"), "ecos"); + // error(errmsg); + // end + end + case "q" then + Checktype("ecos", dims(2*i), "q", 4, "constant"); + if (isempty(dims(2*i))) then + q=[]; + elseif (dims(2*i)<0 | modulo(dims(2*i),1)) then + errmsg = msprintf(gettext("%s: q in dims should be a positive vector"), "ecos"); + error(errmsg); + else + q=dims(2*i); + end + case "e" then + Checktype("ecos", dims(2*i), "e", 4, "constant"); + if (isempty(dims(2*i))) then + e=0; + elseif (dims(2*i)<0 | modulo(dims(2*i),1)) then + errmsg = msprintf(gettext("%s: e in dims should be a positive integer"), "ecos"); + error(errmsg); + else + e=dims(2*i); + // if (prod(size(e))>1)) then + // errmsg = msprintf(gettext("%s: e in dims should be a positive integer"), "ecos"); + // error(errmsg); + // end + end + else + errmsg = msprintf(gettext("%s: Unrecognized parameter name %s."), "ecos", dims(2*i-1)); + error(errmsg); + end + end + + // Error check and converting equaity matrix(A) to CCS format required by ecos + if ((size(A,1)*size(b,1)==0) & (size(A,1)+size(b,1)~=0)) then + errmsg = msprintf(gettext("%s: One of %s and %s is an empty matrix"), "ecos","A","b"); + error(errmsg); + end + + if (size(A,1)~=0) then + Checktype("ecos", A, "A", 5, "constant"); + Checktype("ecos", b, "b", 6, "constant"); + + [A,b,s2]=linconcheck(A,b,["A","b"]); + [Ajc,Air,Apr] = sp2adj(sparse(A)); + [m2,p]= size(A); + + if (n~=p) then + errmsg = msprintf(gettext("%s: Number of columns of G and A do not match"), "ecos"); + error(errmsg); + end + + if (p~=size(c,2)) then + errmsg = msprintf(gettext("%s: Number of columns of A and c do not match"), "ecos"); + error(errmsg); + end + end + + if (size(G,1)==0 & (size(h,1)==0)) then + if (size(A,1)==0 & (size(b,1)==0)) then + errmsg = msprintf(gettext("%s: At most one of the pair (G, h) or (A, b) is allowed to be absent"), "ecos"); + error(errmsg); + end + end + + // Extracting values and Error checks for options(param) + Checktype("ecos", param, "param", 1, "list"); + if (modulo(size(param),2)) then + errmsg = msprintf(gettext("%s: Size of Options (list) should be even"), "ecos"); + error(errmsg); + end + + option = list("maxiter", [100], "feastol", [1e-8],"reltol",[1e-8],"abstol",[1e-8],"feastol_inacc",[1e-4],"abstol_inacc",[5e-5],"reltol_inacc",[5e-5],"verbose",[0],"mi_max_iters",[1000],"mi_int_tol",[1e-4],"mi_abs_eps",[1e-6],"mi_rel_eps",[1e-6]); + + for i = 1:(size(param))/2 + select convstr(param(2*i-1),'l') + case "maxit" then + if (type(option(2*i))~=1 | modulo(option(2*i),1)) then + errmsg = msprintf(gettext("%s: Value for Maximum Iteration should be a Constant integer"), "ecos"); + error(errmsg); + else + option(2) = param(2*i); + end + + case "feastol" then + if (type(option(2*i))~=1) then + errmsg = msprintf(gettext("%s: Value for tolerance should be a Constant"), "ecos"); + error(errmsg); + else + option(4) = param(2*i); + end + + case "reltol" then + if (type(option(2*i))~=1) then + errmsg = msprintf(gettext("%s: Value for relative tolerance should be a Constant"), "ecos"); + error(errmsg); + else + option(6) = param(2*i); + end + + case "abstol" then + if (type(option(2*i))~=1) then + errmsg = msprintf(gettext("%s: Value for absolute tolerance should be a Constant"), "ecos"); + error(errmsg); + else + option(8) = param(2*i); + end + + case "feastol_inacc" then + if (type(option(2*i))~=1) then + errmsg = msprintf(gettext("%s: Value for tolerance with reduced precision should be a Constant"), "ecos"); + error(errmsg); + else + option(10) = param(2*i); + end + + case "reltol_inacc" then + if (type(option(2*i))~=1) then + errmsg = msprintf(gettext("%s: Value for relative tolerance with reduced precision should be a Constant"), "ecos"); + error(errmsg); + else + option(12) = param(2*i); + end + + case "abstol_inacc" then + if (type(option(2*i))~=1) then + errmsg = msprintf(gettext("%s: Value for absolute tolerance with reduced precision should be a Constant"), "ecos"); + error(errmsg); + else + option(14) = param(2*i); + end + + case "verbose" then + if (type(option(2*i))~=1 | modulo(option(2*i),1)) then + errmsg = msprintf(gettext("%s: Value for verbose level should be a Constant integer"), "ecos"); + error(errmsg); + else + option(16) = param(2*i); + end + + case "mi_max_iters" then + if (type(option(2*i))~=1 | modulo(option(2*i),1)) then + errmsg = msprintf(gettext("%s: Value for maximum branch and bound iterations should be a Constant integer"), "ecos"); + error(errmsg); + else + option(18) = param(2*i); + end + + case "mi_int_tol" then + if (type(option(2*i))~=1) then + errmsg = msprintf(gettext("%s: Value for integer tolerance should be a Constant"), "ecos"); + error(errmsg); + else + option(20) = param(2*i); + end + + case "mi_abs_eps" then + if (type(option(2*i))~=1) then + errmsg = msprintf(gettext("%s: Value for absolute tolerance between upper and lower bounds should be a Constant"), "ecos"); + error(errmsg); + else + option(22) = param(2*i); + end + + case "mi_rel_eps" then + if (type(option(2*i))~=1) then + errmsg = msprintf(gettext("%s: Value for relative tolerance between upper and lower bounds should be a Constant"), "ecos"); + error(errmsg); + else + option(24) = param(2*i); + end + + else + errmsg = msprintf(gettext("%s: Unrecognized parameter name %s."), "ecos", param(2*i-1)); + error(errmsg); + end + end + + // number of second order cones is size of q + ncones = size(q,2); + + //Converting form 1-based(Scilab) indexing to 0-based indexing(C) + Gjc=Gjc-1; + Gir=Gir-1; + if (p ~=0) + Ajc=Ajc-1; + Air=Air-1; + end + + // Calling ecos + [x,y,info1,s,z]=solveecos(c,Gpr,int32(Gjc),int32(Gir),h,Apr,int32(Ajc),int32(Air),b,l,int32(q),e,option,m,n,p,ncones) + + // Assigning output parameters + info=struct(); + info.setup_time=info1(1); + info.solve_time = info1(2); + info.primal_objective_cost =info1(3); + info.dual_objective_cost =info1(4); + info.Normalized_primal_residual = info1(5); + info.Normalized_dual_residual = info1(6); + info.primal_infeasibile = info1(7); + info.dual_infeasibile = info1(8); + info.primal_infeasibility_measure = info1(9); + info.dual_infeasibility_measure = info1(10); + info.Complementarity_gap = info1(11); + info.Normalized_complementarity_gap = info1(12); + info.Iterations = info1(13); + info.exitflag = info1(14); + + // Status + select info.exitflag + case 0 then + status='Optimal solution found'; + case 1 then + status='Certificate of primal infeasibility found'; + case 2 then + status='Certificate of dual infeasibility found'; + case 10 then + status='Optimal solution found subject to reduced tolerances'; + case 11 then + status='Certificate of primal infeasibility found subject to reduced tolerances'; + case 12 then + status='Certificate of dual infeasibility found subject to reduced tolerances'; + case -1 then + status='Maximum number of iterations reached'; + case -2 then + status='Numerical problems (unreliable search direction)'; + case -3 then + status='Numerical problems (slacks or multipliers outside cone)'; + case -4 then + status='Interrupted by signal or CTRL-C'; + case -7 then + status='Unknown problem in solver'; + case -8 then + status='ecos setup error'; + else + status='Unknown problem in toolbox,contact toolbox authors'; + end + disp(status); +endfunction
\ No newline at end of file diff --git a/macros/fgoalattain.bin b/macros/fgoalattain.bin Binary files differnew file mode 100644 index 0000000..faa0821 --- /dev/null +++ b/macros/fgoalattain.bin diff --git a/macros/fminbnd.bin b/macros/fminbnd.bin Binary files differnew file mode 100644 index 0000000..24cd387 --- /dev/null +++ b/macros/fminbnd.bin diff --git a/macros/fmincon.bin b/macros/fmincon.bin Binary files differnew file mode 100644 index 0000000..77a2852 --- /dev/null +++ b/macros/fmincon.bin diff --git a/macros/fminimax.bin b/macros/fminimax.bin Binary files differnew file mode 100644 index 0000000..c023720 --- /dev/null +++ b/macros/fminimax.bin diff --git a/macros/fminunc.bin b/macros/fminunc.bin Binary files differnew file mode 100644 index 0000000..b5dee4a --- /dev/null +++ b/macros/fminunc.bin diff --git a/macros/intfminbnd.bin b/macros/intfminbnd.bin Binary files differnew file mode 100644 index 0000000..dc0f380 --- /dev/null +++ b/macros/intfminbnd.bin diff --git a/macros/intfminbnd.sci b/macros/intfminbnd.sci new file mode 100644 index 0000000..2304bbf --- /dev/null +++ b/macros/intfminbnd.sci @@ -0,0 +1,360 @@ +// Copyright (C) 2015 - IIT Bombay - FOSSEE +// +// This file must be used under the terms of the CeCILL. +// This source file is licensed as described in the file COPYING, which +// you should have received as part of this distribution. The terms +// are also available at +// http://www.cecill.info/licences/Licence_CeCILL_V2-en.txt +// Author: Harpreet Singh, Pranav Deshpande and Akshay Miterani +// Organization: FOSSEE, IIT Bombay +// Email: toolbox@scilab.in + +function [xopt,fopt,exitflag,gradient,hessian] = intfminbnd (varargin) + // Solves a multi-variable mixed integer non linear programming optimization problem on a bounded interval + // + // Calling Sequence + // xopt = intfminbnd(f,intcon,x1,x2) + // xopt = intfminbnd(f,intcon,x1,x2,options) + // [xopt,fopt] = intfminbnd(.....) + // [xopt,fopt,exitflag]= intfminbnd(.....) + // [xopt,fopt,exitflag,output]=intfminbnd(.....) + // [xopt,fopt,exitflag,gradient,hessian]=intfminbnd(.....) + // + // Parameters + // f : a function, representing the objective function of the problem + // x1 : a vector, containing the lower bound of the variables. + // x2 : a vector, containing the upper bound of the variables. + // intcon : a vector of integers, represents which variables are constrained to be integers + // options : a list, containing the option for user to specify. See below for details. + // xopt : a vector of doubles, containing the the computed solution of the optimization problem. + // fopt : a scalar of double, containing the the function value at x. + // exitflag : a scalar of integer, containing the flag which denotes the reason for termination of algorithm. See below for details. + // gradient : a vector of doubles, containing the Objective's gradient of the solution. + // hessian : a matrix of doubles, containing the Objective's hessian of the solution. + // + // Description + // Search the minimum of a multi-variable mixed integer non linear programming optimization on bounded interval specified by : + // Find the minimum of f(x) such that + // + // <latex> + // \begin{eqnarray} + // &\mbox{min}_{x} + // & f(x)\\ + // & \text{subject to} & x1 \ < x \ < x2 \\ + // & x_i \in \!\, \mathbb{Z}, i \in \!\, I + // \end{eqnarray} + // </latex> + // + // The routine calls Bonmin for solving the Bounded Optimization problem, Bonmin is a library written in C++. + // + // The options allows the user to set various parameters of the Optimization problem. + // It should be defined as type "list" and contains the following fields. + // <itemizedlist> + // <listitem>Syntax : options= list("IntegerTolerance", [---], "MaxNodes",[---], "MaxIter", [---], "AllowableGap",[---] "CpuTime", [---],"gradobj", "off", "hessian", "off" );</listitem> + // <listitem>IntegerTolerance : a Scalar, a number with that value of an integer is considered integer..</listitem> + // <listitem>MaxNodes : a Scalar, containing the Maximum Number of Nodes that the solver should search.</listitem> + // <listitem>CpuTime : a Scalar, containing the Maximum amount of CPU Time that the solver should take.</listitem> + // <listitem>AllowableGap : a Scalar, to stop the tree search when the gap between the objective value of the best known solution is reached.</listitem> + // <listitem>MaxIter : a Scalar, containing the Maximum Number of Iteration that the solver should take.</listitem> + // <listitem>gradobj : a string, to turn on or off the user supplied objective gradient.</listitem> + // <listitem>hessian : a Scalar, to turn on or off the user supplied objective hessian.</listitem> + // <listitem>Default Values : options = list('integertolerance',1d-06,'maxnodes',2147483647,'cputime',1d10,'allowablegap',0,'maxiter',2147483647,'gradobj',"off",'hessian',"off")</listitem> + // </itemizedlist> + // + // The exitflag allows to know the status of the optimization which is given back by Ipopt. + // <itemizedlist> + // <listitem>exitflag=0 : Optimal Solution Found </listitem> + // <listitem>exitflag=1 : Maximum Number of Iterations Exceeded. Output may not be optimal.</listitem> + // <listitem>exitflag=2 : Maximum CPU Time exceeded. Output may not be optimal.</listitem> + // <listitem>exitflag=3 : Stop at Tiny Step.</listitem> + // <listitem>exitflag=4 : Solved To Acceptable Level.</listitem> + // <listitem>exitflag=5 : Converged to a point of local infeasibility.</listitem> + // </itemizedlist> + // + // For more details on exitflag see the Bonmin documentation, go to http://www.coin-or.org/Bonmin + // + // Examples + // //Find x in R^6 such that it minimizes: + // //f(x)= sin(x1) + sin(x2) + sin(x3) + sin(x4) + sin(x5) + sin(x6) + // //-2 <= x1,x2,x3,x4,x5,x6 <= 2 + // //Objective function to be minimised + // function y=f(x) + // y=0 + // for i =1:6 + // y=y+sin(x(i)); + // end + // endfunction + // //Variable bounds + // x1 = [-2, -2, -2, -2, -2, -2]; + // x2 = [2, 2, 2, 2, 2, 2]; + // intcon = [2 3 4] + // //Options + // options=list("MaxIter",[1500],"CpuTime", [100]) + // [x,fval] =intfminbnd(f ,intcon, x1, x2, options) + // // Press ENTER to continue + // + // Examples + // //Find x in R such that it minimizes: + // //f(x)= 1/x^2 + // //0 <= x <= 1000 + // //Objective function to be minimised + // function y=f(x) + // y=1/x^2; + // endfunction + // //Variable bounds + // x1 = [0]; + // x2 = [1000]; + // intcon = [1]; + // [x,fval,exitflag,output,lambda] =intfminbnd(f,intcon , x1, x2) + // // Press ENTER to continue + // + // Examples + // //The below problem is an unbounded problem: + // //Find x in R^2 such that it minimizes: + // //f(x)= -[(x1-1)^2 + (x2-1)^2] + // //-inf <= x1,x2 <= inf + // //Objective function to be minimised + // function y=f(x) + // y=-((x(1)-1)^2+(x(2)-1)^2); + // endfunction + // //Variable bounds + // x1 = [-%inf , -%inf]; + // x2 = [ %inf , %inf]; + // //Options + // options=list("MaxIter",[1500],"CpuTime", [100]) + // [x,fval,exitflag,output,lambda] =intfminbnd(f,intcon, x1, x2, options) + // Authors + // Harpreet Singh + + //To check the number of input and output arguments + [lhs , rhs] = argn(); + + //To check the number of arguments given by the user + if ( rhs<4 | rhs>5 ) then + errmsg = msprintf(gettext("%s: Unexpected number of input arguments : %d provided while should be int [4 5] "), "intfminbnd", rhs); + error(errmsg); + end + + //Storing the Input Parameters + fun = varargin(1); + intcon = varargin(2); + x1 = varargin(3); + x2 = varargin(4); + nbvar = size(x1,"*"); + + param = list(); + //To check whether options has been entered by user + if ( rhs>=5 ) then + param =varargin(5); + end + + //To check whether the Input arguments + Checktype("intfminbnd", fun, "fun", 1, "function"); + Checktype("intfminbnd", intcon, "intcon", 2, "constant"); + Checktype("intfminbnd", x1, "x1", 3, "constant"); + Checktype("intfminbnd", x2, "x2", 4, "constant"); + Checktype("intfminbnd", param, "options", 5, "list"); + + + if(nbvar==0) then + errmsg = msprintf(gettext("%s: x1 cannot be an empty"), "intfminbnd"); + error(errmsg); + end + + ///////////////// To check vectors ///////////////// + Checkvector("intfminbnd", x1, "x1", 3, nbvar) + x1 = x1(:); + Checkvector("intfminbnd", x2, "x2", 4, nbvar) + x2 = x2(:); + Checkvector("intfminbnd", intcon, "intcon", 2, size(intcon,"*")) + intcon = intcon(:); + + if(~isequal(size(x1),size(x2))) then + errmsg = msprintf(gettext("%s: x1 and x2 should be of same size"), "intfminbnd"); + error(errmsg); + end + + for i=1:size(intcon,1) + if(intcon(i)>nbvar) then + errmsg = msprintf(gettext("%s: The values inside intcon should be less than the number of variables"), "intfminbnd"); + error(errmsg); + end + + if (intcon(i)<0) then + errmsg = msprintf(gettext("%s: The values inside intcon should be greater than 0 "), "intfminbnd"); + error(errmsg); + end + + if(modulo(intcon(i),1)) then + errmsg = msprintf(gettext("%s: The values inside intcon should be an integer "), "intfminbnd"); + error(errmsg); + end + end + +options = list('integertolerance',1d-06,'maxnodes',2147483647,'cputime',1d10,'allowablegap',0,'maxiter',2147483647,'gradobj',"off",'hessian',"off") + + //Pushing param into default value + + for i = 1:(size(param))/2 + select convstr(param(2*i-1),'l') + case 'integertolerance' then + Checktype("intfminbnd_options", param(2*i), param(2*i-1), 2*i, "constant"); + options(2) = param(2*i); + case 'maxnodes' then + Checktype("intfminbnd_options", param(2*i), param(2*i-1), 2*i, "constant"); + options(4) = options(2*i); + case 'cputime' then + Checktype("intfminbnd_options", param(2*i), param(2*i-1), 2*i, "constant"); + options(6) = options(2*i); + case 'allowablegap' then + Checktype("intfminbnd_options", param(2*i), param(2*i-1), 2*i, "constant"); + options(8) = options(2*i); + case 'maxiter' then + Checktype("intfminbnd_options", param(2*i), param(2*i-1), 2*i, "constant"); + options(10) = options(2*i); + case 'gradobj' then + Checktype("intfminbnd_options", param(2*i), param(2*i-1), 2*i, "string"); + if(convstr(options(2*i),'l') == "on") then + options(12) = "on" + elseif(convstr(options(2*i),'l') == "off") then + options(12) = "off" + else + error(999, 'Unknown string passed in gradobj.'); + end + case 'hessian' then + Checktype("intfminbnd_options", param(2*i), param(2*i-1), 2*i, "string"); + if(convstr(options(2*i),'l') == "on") then + options(14) = "on"; + elseif(convstr(options(2*i),'l') == "off") then + options(14) = "off"; + else + error(999, 'Unknown string passed in hessian.'); + end + else + error(999, 'Unknown string argument passed.'); + end + end + + ///////////////// Functions Check ///////////////// + + //To check the match between f (1st Parameter) and x1 (2nd Parameter) + if(execstr('init=fun(x1)','errcatch')==21) then + errmsg = msprintf(gettext("%s: Objective function and x1 did not match"), "intfminbnd"); + error(errmsg); + end + + if(options(12) == "on") then + if(execstr('[grad_y,grad_dy]=fun(x1)','errcatch')==59) then + errmsg = msprintf(gettext("%s: Gradient of objective function is not provided"), "intfminbnd"); + error(errmsg); + end + Checkvector("intfminbnd_options", grad_dy, "dy", 12, nbvar); + end + + if(options(14) == "on") then + if(execstr('[hessian_y,hessian_dy,hessian]=fun(x1)','errcatch')==59) then + errmsg = msprintf(gettext("%s: Gradient of objective function is not provided"), "intfminbnd"); + error(errmsg); + end + + if ( ~isequal(size(hessian) == [nbvar nbvar]) ) then + errmsg = msprintf(gettext("%s: Size of hessian should be nbvar X nbvar"), "intfminbnd"); + error(errmsg); + end + end + + //Converting the User defined Objective function into Required form (Error Detectable) + function [y,check] = _f(x) + try + y=fun(x) + [y,check] = checkIsreal(y) + catch + y=0; + check=1; + end + endfunction + + //Defining a function to calculate Hessian if the respective user entry is OFF + function [hessy,check]=_gradhess(x) + if (options(14) == "on") then + try + [obj,dy,hessy] = fun(x) + [hessy,check] = checkIsreal(hessy) + catch + hessy = 0; + check=1; + end + else + try + [dy,hessy]=numderivative(fun,x) + [hessy,check] = checkIsreal(hessy) + catch + hessy=0; + check=1; + end + end + endfunction + + //Defining an inbuilt Objective gradient function + function [dy,check] = _gradf(x) + if (options(12) =="on") then + try + [y,dy]=fun(x) + [dy,check] = checkIsreal(dy) + catch + dy = 0; + check=1; + end + else + try + dy=numderivative(fun,x) + [dy,check] = checkIsreal(dy) + catch + dy=0; + check=1; + end + end + endfunction + + intconsize = size(intcon,"*"); + + [xopt,fopt,exitflag] = inter_fminbnd(_f,_gradf,_gradhess,x1,x2,intcon,options,nbvar); + + //In the cases of the problem not being solved, return NULL to the output matrices + if( exitflag~=0 & exitflag~=3 ) then + gradient = []; + hessian = []; + else + [ gradient, hessian] = numderivative(_f, xopt, [], [], "blockmat"); + end + + //To print output message + select exitflag + + case 0 then + printf("\nOptimal Solution Found.\n"); + case 1 then + printf("\nInFeasible Solution.\n"); + case 2 then + printf("\nnObjective Function is Continuous Unbounded.\n"); + case 3 then + printf("\Limit Exceeded.\n"); + case 4 then + printf("\nUser Interrupt.\n"); + case 5 then + printf("\nMINLP Error.\n"); + else + printf("\nInvalid status returned. Notify the Toolbox authors\n"); + break; + end +endfunction + +function [y, check] = checkIsreal(x) + if ((~isreal(x))) then + y = 0 + check=1; + else + y = x; + check=0; + end +endfunction diff --git a/macros/intfmincon.bin b/macros/intfmincon.bin Binary files differnew file mode 100644 index 0000000..b925f44 --- /dev/null +++ b/macros/intfmincon.bin diff --git a/macros/intfmincon.sci b/macros/intfmincon.sci new file mode 100644 index 0000000..cbdd2ef --- /dev/null +++ b/macros/intfmincon.sci @@ -0,0 +1,589 @@ +// Copyright (C) 2015 - IIT Bombay - FOSSEE +// +// This file must be used under the terms of the CeCILL. +// This source file is licensed as described in the file COPYING, which +// you should have received as part of this distribution. The terms +// are also available at +// http://www.cecill.info/licences/Licence_CeCILL_V2-en.txt +// Author: Harpreet Singh, Pranav Deshpande and Akshay Miterani +// Organization: FOSSEE, IIT Bombay +// Email: toolbox@scilab.in + +function [xopt,fopt,exitflag,gradient,hessian] = intfmincon (varargin) + // Solves a constrainted multi-variable mixed integer non linear programming problem + // + // Calling Sequence + // xopt = intfmincon(f,x0,intcon,A,b) + // xopt = intfmincon(f,x0,intcon,A,b,Aeq,beq) + // xopt = intfmincon(f,x0,intcon,A,b,Aeq,beq,lb,ub) + // xopt = intfmincon(f,x0,intcon,A,b,Aeq,beq,lb,ub,nlc) + // xopt = intfmincon(f,x0,intcon,A,b,Aeq,beq,lb,ub,nlc,options) + // [xopt,fopt] = intfmincon(.....) + // [xopt,fopt,exitflag]= intfmincon(.....) + // [xopt,fopt,exitflag,gradient]=intfmincon(.....) + // [xopt,fopt,exitflag,gradient,hessian]=intfmincon(.....) + // + // Parameters + // f : a function, representing the objective function of the problem + // x0 : a vector of doubles, containing the starting values of variables. + // intcon : a vector of integers, represents which variables are constrained to be integers + // A : a matrix of double, represents the linear coefficients in the inequality constraints A⋅x ≤ b. + // b : a vector of double, represents the linear coefficients in the inequality constraints A⋅x ≤ b. + // Aeq : a matrix of double, represents the linear coefficients in the equality constraints Aeq⋅x = beq. + // beq : a vector of double, represents the linear coefficients in the equality constraints Aeq⋅x = beq. + // lb : Lower bounds, specified as a vector or array of double. lb represents the lower bounds elementwise in lb ≤ x ≤ ub. + // ub : Upper bounds, specified as a vector or array of double. ub represents the upper bounds elementwise in lb ≤ x ≤ ub. + // nlc : a function, representing the Non-linear Constraints functions(both Equality and Inequality) of the problem. It is declared in such a way that non-linear inequality constraints are defined first as a single row vector (c), followed by non-linear equality constraints as another single row vector (ceq). Refer Example for definition of Constraint function. + // options : a list, containing the option for user to specify. See below for details. + // xopt : a vector of doubles, containing the the computed solution of the optimization problem. + // fopt : a scalar of double, containing the the function value at x. + // exitflag : a scalar of integer, containing the flag which denotes the reason for termination of algorithm. See below for details. + // gradient : a vector of doubles, containing the Objective's gradient of the solution. + // hessian : a matrix of doubles, containing the Objective's hessian of the solution. + // + // Description + // Search the minimum of a mixed integer constrained optimization problem specified by : + // Find the minimum of f(x) such that + // + // <latex> + // \begin{eqnarray} + // &\mbox{min}_{x} + // & f(x) \\ + // & \text{subject to} & A*x \leq b \\ + // & & Aeq*x \ = beq\\ + // & & c(x) \leq 0\\ + // & & ceq(x) \ = 0\\ + // & & lb \leq x \leq ub \\ + // & & x_i \in \!\, \mathbb{Z}, i \in \!\, I + // \end{eqnarray} + // </latex> + // + // The routine calls Bonmin for solving the Bounded Optimization problem, Bonmin is a library written in C++. + // + // The options allows the user to set various parameters of the Optimization problem. + // It should be defined as type "list" and contains the following fields. + // <itemizedlist> + // <listitem>Syntax : options= list("IntegerTolerance", [---], "MaxNodes",[---], "MaxIter", [---], "AllowableGap",[---] "CpuTime", [---],"gradobj", "off", "hessian", "off" );</listitem> + // <listitem>IntegerTolerance : a Scalar, a number with that value of an integer is considered integer..</listitem> + // <listitem>MaxNodes : a Scalar, containing the Maximum Number of Nodes that the solver should search.</listitem> + // <listitem>CpuTime : a Scalar, containing the Maximum amount of CPU Time that the solver should take.</listitem> + // <listitem>AllowableGap : a Scalar, to stop the tree search when the gap between the objective value of the best known solution is reached.</listitem> + // <listitem>MaxIter : a Scalar, containing the Maximum Number of Iteration that the solver should take.</listitem> + // <listitem>gradobj : a string, to turn on or off the user supplied objective gradient.</listitem> + // <listitem>hessian : a Scalar, to turn on or off the user supplied objective hessian.</listitem> + // <listitem>Default Values : options = list('integertolerance',1d-06,'maxnodes',2147483647,'cputime',1d10,'allowablegap',0,'maxiter',2147483647,'gradobj',"off",'hessian',"off")</listitem> + // </itemizedlist> + // + // The exitflag allows to know the status of the optimization which is given back by Ipopt. + // <itemizedlist> + // <listitem>exitflag=0 : Optimal Solution Found </listitem> + // <listitem>exitflag=1 : InFeasible Solution.</listitem> + // <listitem>exitflag=2 : Objective Function is Continuous Unbounded.</listitem> + // <listitem>exitflag=3 : Limit Exceeded.</listitem> + // <listitem>exitflag=4 : User Interrupt.</listitem> + // <listitem>exitflag=5 : MINLP Error.</listitem> + // </itemizedlist> + // + // For more details on exitflag see the Bonmin documentation, go to http://www.coin-or.org/Bonmin + // + // Examples + // //Find x in R^2 such that it minimizes: + // //f(x)= -x1 -x2/3 + // //x0=[0,0] + // //constraint-1 (c1): x1 + x2 <= 2 + // //constraint-2 (c2): x1 + x2/4 <= 1 + // //constraint-3 (c3): x1 - x2 <= 2 + // //constraint-4 (c4): -x1/4 - x2 <= 1 + // //constraint-5 (c5): -x1 - x2 <= -1 + // //constraint-6 (c6): -x1 + x2 <= 2 + // //constraint-7 (c7): x1 + x2 = 2 + // //Objective function to be minimised + // function [y,dy]=f(x) + // y=-x(1)-x(2)/3; + // dy= [-1,-1/3]; + // endfunction + // //Starting point, linear constraints and variable bounds + // x0=[0 , 0]; + // intcon = [1] + // A=[1,1 ; 1,1/4 ; 1,-1 ; -1/4,-1 ; -1,-1 ; -1,1]; + // b=[2;1;2;1;-1;2]; + // Aeq=[1,1]; + // beq=[2]; + // lb=[]; + // ub=[]; + // nlc=[]; + // //Options + // options=list("GradObj", "on"); + // //Calling Ipopt + // [x,fval,exitflag,grad,hessian] =intfmincon(f, x0,intcon,A,b,Aeq,beq,lb,ub,nlc,options) + // // Press ENTER to continue + // + // Examples + // //Find x in R^3 such that it minimizes: + // //f(x)= x1*x2 + x2*x3 + // //x0=[0.1 , 0.1 , 0.1] + // //constraint-1 (c1): x1^2 - x2^2 + x3^2 <= 2 + // //constraint-2 (c2): x1^2 + x2^2 + x3^2 <= 10 + // //Objective function to be minimised + // function [y,dy]=f(x) + // y=x(1)*x(2)+x(2)*x(3); + // dy= [x(2),x(1)+x(3),x(2)]; + // endfunction + // //Starting point, linear constraints and variable bounds + // x0=[0.1 , 0.1 , 0.1]; + // intcon = [2] + // A=[]; + // b=[]; + // Aeq=[]; + // beq=[]; + // lb=[]; + // ub=[]; + // //Nonlinear constraints + // function [c,ceq,cg,cgeq]=nlc(x) + // c = [x(1)^2 - x(2)^2 + x(3)^2 - 2 , x(1)^2 + x(2)^2 + x(3)^2 - 10]; + // ceq = []; + // cg=[2*x(1) , -2*x(2) , 2*x(3) ; 2*x(1) , 2*x(2) , 2*x(3)]; + // cgeq=[]; + // endfunction + // //Options + // options=list("MaxIter", [1500], "CpuTime", [500], "GradObj", "on","GradCon", "on"); + // //Calling Ipopt + // [x,fval,exitflag,output] =intfmincon(f, x0,intcon,A,b,Aeq,beq,lb,ub,nlc,options) + // // Press ENTER to continue + // + // Examples + // //The below problem is an unbounded problem: + // //Find x in R^3 such that it minimizes: + // //f(x)= -(x1^2 + x2^2 + x3^2) + // //x0=[0.1 , 0.1 , 0.1] + // // x1 <= 0 + // // x2 <= 0 + // // x3 <= 0 + // //Objective function to be minimised + // function y=f(x) + // y=-(x(1)^2+x(2)^2+x(3)^2); + // endfunction + // //Starting point, linear constraints and variable bounds + // x0=[0.1 , 0.1 , 0.1]; + // intcon = [3] + // A=[]; + // b=[]; + // Aeq=[]; + // beq=[]; + // lb=[]; + // ub=[0,0,0]; + // //Options + // options=list("MaxIter", [1500], "CpuTime", [500]); + // //Calling Ipopt + // [x,fval,exitflag,grad,hessian] =intfmincon(f, x0,intcon,A,b,Aeq,beq,lb,ub,[],options) + // // Press ENTER to continue + // + // Examples + // //The below problem is an infeasible problem: + // //Find x in R^3 such that in minimizes: + // //f(x)=x1*x2 + x2*x3 + // //x0=[1,1,1] + // //constraint-1 (c1): x1^2 <= 1 + // //constraint-2 (c2): x1^2 + x2^2 <= 1 + // //constraint-3 (c3): x3^2 <= 1 + // //constraint-4 (c4): x1^3 = 0.5 + // //constraint-5 (c5): x2^2 + x3^2 = 0.75 + // // 0 <= x1 <=0.6 + // // 0.2 <= x2 <= inf + // // -inf <= x3 <= 1 + // //Objective function to be minimised + // function [y,dy]=f(x) + // y=x(1)*x(2)+x(2)*x(3); + // dy= [x(2),x(1)+x(3),x(2)]; + // endfunction + // //Starting point, linear constraints and variable bounds + // x0=[1,1,1]; + // intcon = [2] + // A=[]; + // b=[]; + // Aeq=[]; + // beq=[]; + // lb=[0 0.2,-%inf]; + // ub=[0.6 %inf,1]; + // //Nonlinear constraints + // function [c,ceq,cg,cgeq]=nlc(x) + // c=[x(1)^2-1,x(1)^2+x(2)^2-1,x(3)^2-1]; + // ceq=[x(1)^3-0.5,x(2)^2+x(3)^2-0.75]; + // cg = [2*x(1),0,0;2*x(1),2*x(2),0;0,0,2*x(3)]; + // cgeq = [3*x(1)^2,0,0;0,2*x(2),2*x(3)]; + // endfunction + // //Options + // options=list("MaxIter", [1500], "CpuTime", [500], "GradObj", "on","GradCon", "on"); + // //Calling Ipopt + // [x,fval,exitflag,grad,hessian] =intfmincon(f, x0,intcon,A,b,Aeq,beq,lb,ub,nlc,options) + // // Press ENTER to continue + // Authors + // Harpreet Singh + + //To check the number of input and output arguments + [lhs , rhs] = argn(); + + //To check the number of arguments given by the user + if ( rhs<4 | rhs>11 ) then + errmsg = msprintf(gettext("%s: Unexpected number of input arguments : %d provided while should be int [4 5] "), "intfmincon", rhs); + error(errmsg); + end + + //Storing the Input Parameters + fun = varargin(1); + x0 = varargin(2); + intcon = varargin(3); + A = varargin(4); + b = varargin(5); + Aeq = []; + beq = []; + lb = []; + ub = []; + nlc = []; + + if (rhs>5) then + Aeq = varargin(6); + beq = varargin(7); + end + + if (rhs>7) then + lb = varargin(8); + ub = varargin(9); + end + + if (rhs>9) then + nlc = varargin(10); + end + + param = list(); + //To check whether options has been entered by user + if ( rhs> 10) then + param =varargin(11); + end + + //To check whether the Input arguments + Checktype("intfmincon", fun, "fun", 1, "function"); + Checktype("intfmincon", x0, "x0", 2, "constant"); + Checktype("intfmincon", intcon, "intcon", 3, "constant"); + Checktype("intfmincon", A, "A", 4, "constant"); + Checktype("intfmincon", b, "b", 5, "constant"); + Checktype("intfmincon", Aeq, "Aeq", 6, "constant"); + Checktype("intfmincon", beq, "beq", 7, "constant"); + Checktype("intfmincon", lb, "lb", 8, "constant"); + Checktype("intfmincon", ub, "ub", 9, "constant"); + Checktype("intfmincon", nlc, "nlc", 10, ["constant","function"]); + Checktype("intfmincon", param, "options", 11, "list"); + + + nbVar = size(x0,"*"); + if(nbVar==0) then + errmsg = msprintf(gettext("%s: x0 cannot be an empty"), "intfmincon"); + error(errmsg); + end + + if(size(lb,"*")==0) then + lb = repmat(-%inf,nbVar,1); + end + + if(size(ub,"*")==0) then + ub = repmat(%inf,nbVar,1); + end + + //////////////// To Check linear constraints ///////// + + //To check for correct size of A(3rd paramter) + if(size(A,2)~=nbVar & size(A,2)~=0) then + errmsg = msprintf(gettext("%s: Expected Matrix of size (No of linear inequality constraints X No of Variables) or an Empty Matrix for Linear Inequality Constraint coefficient Matrix A"), "intfmincon"); + error(errmsg); + end + nbConInEq=size(A,"r"); + + //To check for the correct size of Aeq (5th paramter) + if(size(Aeq,2)~=nbVar & size(Aeq,2)~=0) then + errmsg = msprintf(gettext("%s: Expected Matrix of size (No of linear equality constraints X No of Variables) or an Empty Matrix for Linear Equality Constraint coefficient Matrix Aeq"), "intfmincon"); + error(errmsg); + end + nbConEq=size(Aeq,"r"); + + ///////////////// To check vectors ///////////////// + + Checkvector("intfmincon", x0, "x0", 2, nbVar); + x0 = x0(:); + if(size(intcon,"*")) then + Checkvector("intfmincon", intcon, "intcon", 3, size(intcon,"*")) + intcon = intcon(:); + end + if(nbConInEq) then + Checkvector("intfmincon", b, "b", 5, nbConInEq); + b = b(:); + end + if(nbConEq) then + Checkvector("intfmincon", beq, "beq", 7, nbConEq); + beq = beq(:); + end + Checkvector("intfmincon", lb, "lb", 8, nbVar); + lb = lb(:); + + Checkvector("intfmincon", ub, "ub", 9, nbVar); + ub = ub(:); + + /////////////// To check integer ////////////////////// + for i=1:size(intcon,1) + if(intcon(i)>nbVar) then + errmsg = msprintf(gettext("%s: The values inside intcon should be less than the number of variables"), "intfmincon"); + error(errmsg); + end + + if (intcon(i)<0) then + errmsg = msprintf(gettext("%s: The values inside intcon should be greater than 0 "), "intfmincon"); + error(errmsg); + end + + if(modulo(intcon(i),1)) then + errmsg = msprintf(gettext("%s: The values inside intcon should be an integer "), "intfmincon"); + error(errmsg); + end + end + +options = list('integertolerance',1d-06,'maxnodes',2147483647,'cputime',1d10,'allowablegap',0,'maxiter',2147483647,'gradobj',"off",'hessian',"off",'gradcon',"off") + + //Pushing param into default value + + for i = 1:(size(param))/2 + select convstr(param(2*i-1),'l') + case 'integertolerance' then + Checktype("intfmincon_options", param(2*i), param(2*i-1), 2*i, "constant"); + options(2) = param(2*i); + case 'maxnodes' then + Checktype("intfmincon_options", param(2*i), param(2*i-1), 2*i, "constant"); + options(4) = param(2*i); + case 'cputime' then + Checktype("intfmincon_options", param(2*i), param(2*i-1), 2*i, "constant"); + options(6) = param(2*i); + case 'allowablegap' then + Checktype("intfmincon_options", param(2*i), param(2*i-1), 2*i, "constant"); + options(8) = param(2*i); + case 'maxiter' then + Checktype("intfmincon_options", param(2*i), param(2*i-1), 2*i, "constant"); + options(10) = param(2*i); + case 'gradobj' then + Checktype("intfmincon_options", param(2*i), param(2*i-1), 2*i, "string"); + if(convstr(param(2*i),'l') == "on") then + options(12) = "on" + elseif(convstr(param(2*i),'l') == "off") then + options(12) = "off" + else + error(999, 'Unknown string passed in gradobj.'); + end + case 'hessian' then + Checktype("intfmincon_options", param(2*i), param(2*i-1), 2*i, "function"); + options(14) = param(2*i); + case 'gradcon' then + Checktype("intfmincon_options", param(2*i), param(2*i-1), 2*i, "string"); + if(convstr(param(2*i),'l') == "on") then + options(16) = "on" + elseif(convstr(param(2*i),'l') == "off") then + options(16) = "off" + else + error(999, 'Unknown string passed in gradcon.'); + end + else + error(999, 'Unknown string argument passed.'); + end + end + + ///////////////// Functions Check ///////////////// + + //To check the match between f (1st Parameter) and x0 (2nd Parameter) + if(execstr('init=fun(x0)','errcatch')==21) then + errmsg = msprintf(gettext("%s: Objective function and x0 did not match"), "intfmincon"); + error(errmsg); + end + + if(options(12) == "on") then + if(execstr('[grad_y,grad_dy]=fun(x0)','errcatch')==59) then + errmsg = msprintf(gettext("%s: Gradient of objective function is not provided"), "intfmincon"); + error(errmsg); + end + if(grad_dy<>[]) then + Checkvector("intfmincon_options", grad_dy, "dy", 12, nbVar); + end + end + + if(options(14) == "on") then + if(execstr('[hessian_y,hessian_dy,hessian]=fun(x0)','errcatch')==59) then + errmsg = msprintf(gettext("%s: Gradient of objective function is not provided"), "intfmincon"); + error(errmsg); + end + if ( ~isequal(size(hessian) == [nbVar nbVar]) ) then + errmsg = msprintf(gettext("%s: Size of hessian should be nbVar X nbVar"), "intfmincon"); + error(errmsg); + end + end + + numNlic = 0; + numNlec = 0; + numNlc = 0; + + if (type(nlc) == 13 | type(nlc) == 11) then + [sample_c,sample_ceq] = nlc(x0); + if(execstr('[sample_c,sample_ceq] = nlc(x0)','errcatch')==21) then + errmsg = msprintf(gettext("%s: Non-Linear Constraint function and x0 did not match"), "intfmincon"); + error(errmsg); + end + numNlic = size(sample_c,"*"); + numNlec = size(sample_ceq,"*"); + numNlc = numNlic + numNlec; + end + + /////////////// Creating conLb and conUb //////////////////////// + + conLb = [repmat(-%inf,numNlic,1);repmat(0,numNlec,1);repmat(-%inf,nbConInEq,1);beq;] + conUb = [repmat(0,numNlic,1);repmat(0,numNlec,1);b;beq;] + + //Converting the User defined Objective function into Required form (Error Detectable) + function [y,check] = _f(x) + try + y=fun(x) + [y,check] = checkIsreal(y) + catch + y=0; + check=1; + end + endfunction + + //Defining an inbuilt Objective gradient function + function [dy,check] = _gradf(x) + if (options(12) =="on") then + try + [y,dy]=fun(x) + [dy,check] = checkIsreal(dy) + catch + dy = 0; + check=1; + end + else + try + dy=numderivative(fun,x) + [dy,check] = checkIsreal(dy) + catch + dy=0; + check=1; + end + end + endfunction + + function [y,check] = _addnlc(x) + x= x(:) + c = [] + ceq = [] + try + if((type(nlc) == 13 | type(nlc) == 11) & numNlc~=0) then + [c,ceq]=nlc(x) + end + ylin = [A*x;Aeq*x]; + y = [c(:);ceq(:);ylin(:);]; + [y,check] = checkIsreal(y) + catch + y=0; + check=1; + end + endfunction + + //Defining an inbuilt jacobian of constraints function + function [dy,check] = _gradnlc(x) + if (options(16) =="on") then + try + [y1,y2,dy1,dy2]=nlc(x) + //Adding derivative of Linear Constraint + dylin = [A;Aeq] + dy = [dy1;dy2;dylin]; + [dy,check] = checkIsreal(dy) + catch + dy = 0; + check=1; + end + else + try + dy=numderivative(_addnlc,x) + [dy,check] = checkIsreal(dy) + catch + dy=0; + check=1; + end + end + endfunction + + //Defining a function to calculate Hessian if the respective user entry is OFF + function [hessy,check]=_gradhess(x,obj_factor,lambda) + x=x(:); + if (type(options(14)) == "function") then + try + [obj,dy,hessy] = fun(x,obj_factor,lambda) + [hessy,check] = checkIsreal(hessy) + catch + hessy = 0; + check=1; + end + else + try + [dy,hessfy]=numderivative(_f,x) + hessfy = matrix(hessfy,nbVar,nbVar) + if((type(nlc) == 13 | type(nlc) == 11) & numNlc~=0) then + [dy,hessny]=numderivative(nlc,x) + end + hessianc = [] + for i = 1:numNlc + hessianc = hessianc + lambda(i)*matrix(hessny(i,:),nbVar,nbVar) + end + hessy = obj_factor*hessfy + hessianc; + [hessy,check] = checkIsreal(hessy) + catch + hessy=0; + check=1; + end + end + endfunction + + intconsize = size(intcon,"*") + + [xopt,fopt,exitflag] = inter_fmincon(_f,_gradf,_addnlc,_gradnlc,_gradhess,x0,lb,ub,conLb,conUb,intcon,options,nbConInEq+nbConEq); + + //In the cases of the problem not being solved, return NULL to the output matrices + if( exitflag~=0 & exitflag~=3 ) then + gradient = []; + hessian = []; + else + [ gradient, hessian] = numderivative(_f, xopt) + end + + //To print output message + select exitflag + + case 0 then + printf("\nOptimal Solution Found.\n"); + case 1 then + printf("\nInFeasible Solution.\n"); + case 2 then + printf("\nObjective Function is Continuous Unbounded.\n"); + case 3 then + printf("\Limit Exceeded.\n"); + case 4 then + printf("\nUser Interrupt.\n"); + case 5 then + printf("\nMINLP Error.\n"); + else + printf("\nInvalid status returned. Notify the Toolbox authors\n"); + break; + end +endfunction + +function [y, check] = checkIsreal(x) + if ((~isreal(x))) then + y = 0 + check=1; + else + y = x; + check=0; + end +endfunction diff --git a/macros/intfminimax.bin b/macros/intfminimax.bin Binary files differnew file mode 100644 index 0000000..8d7af05 --- /dev/null +++ b/macros/intfminimax.bin diff --git a/macros/intfminimax.sci b/macros/intfminimax.sci new file mode 100644 index 0000000..bf3ec97 --- /dev/null +++ b/macros/intfminimax.sci @@ -0,0 +1,456 @@ +// Copyright (C) 2015 - IIT Bombay - FOSSEE +// +// This file must be used under the terms of the CeCILL. +// This source file is licensed as described in the file COPYING, which +// you should have received as part of this distribution. The terms +// are also available at +// http://www.cecill.info/licences/Licence_CeCILL_V2-en.txt +// Authors: Animesh Baranawal +// Organization: FOSSEE, IIT Bombay +// Email: toolbox@scilab.in + +function [x,fval,maxfval,exitflag] = intfminimax(varargin) + // Solves minimax constraint problem + // + // Calling Sequence + // xopt = intfminimax(fun,x0,intcon) + // xopt = intfminimax(fun,x0,intcon,A,b) + // xopt = intfminimax(fun,x0,intcon,A,b,Aeq,beq) + // xopt = intfminimax(fun,x0,intcon,A,b,Aeq,beq,lb,ub) + // xopt = intfminimax(fun,x0,intcon,A,b,Aeq,beq,lb,ub,nonlinfun) + // xopt = intfminimax(fun,x0,intcon,A,b,Aeq,beq,lb,ub,nonlinfun,options) + // [xopt, fval] = intfminimax(.....) + // [xopt, fval, maxfval]= intfminimax(.....) + // [xopt, fval, maxfval, exitflag]= intfminimax(.....) + // + // Parameters + // fun: The function to be minimized. fun is a function that accepts a vector x and returns a vector F, the objective functions evaluated at x. + // x0 : a vector of double, contains initial guess of variables. + // A : a matrix of double, represents the linear coefficients in the inequality constraints A⋅x ≤ b. + // intcon : a vector of integers, represents which variables are constrained to be integers + // b : a vector of double, represents the linear coefficients in the inequality constraints A⋅x ≤ b. + // Aeq : a matrix of double, represents the linear coefficients in the equality constraints Aeq⋅x = beq. + // beq : a vector of double, represents the linear coefficients in the equality constraints Aeq⋅x = beq. + // lb : a vector of double, contains lower bounds of the variables. + // ub : a vector of double, contains upper bounds of the variables. + // nonlinfun: function that computes the nonlinear inequality constraints c⋅x ≤ 0 and nonlinear equality constraints c⋅x = 0. + // xopt : a vector of double, the computed solution of the optimization problem. + // fopt : a double, the value of the function at x. + // maxfval: a 1x1 matrix of doubles, the maximum value in vector fval + // exitflag : The exit status. See below for details. + // output : The structure consist of statistics about the optimization. See below for details. + // lambda : The structure consist of the Lagrange multipliers at the solution of problem. See below for details. + // + // Description + // intfminimax minimizes the worst-case (largest) value of a set of multivariable functions, starting at an initial estimate. This is generally referred to as the minimax problem. + // + // <latex> + // \min_{x} \max_{i} F_{i}(x)\: \textrm{such that} \:\begin{cases} + // & c(x) \leq 0 \\ + // & ceq(x) = 0 \\ + // & A.x \leq b \\ + // & Aeq.x = beq \\ + // & lb \leq x \leq ub + // & x_i \in \!\, \mathbb{Z}, i \in \!\, I + // \end{cases} + // </latex> + // + // Currently, intfminimax calls intfmincon which uses the bonmin algorithm. + // + // max-min problems can also be solved with intfminimax, using the identity + // + // <latex> + // \max_{x} \min_{i} F_{i}(x) = -\min_{x} \max_{i} \left( -F_{i}(x) \right) + // </latex> + // + // The options allows the user to set various parameters of the Optimization problem. + // It should be defined as type "list" and contains the following fields. + // <itemizedlist> + // <listitem>Syntax : options= list("IntegerTolerance", [---], "MaxNodes",[---], "MaxIter", [---], "AllowableGap",[---] "CpuTime", [---],"gradobj", "off", "hessian", "off" );</listitem> + // <listitem>IntegerTolerance : a Scalar, a number with that value of an integer is considered integer..</listitem> + // <listitem>MaxNodes : a Scalar, containing the Maximum Number of Nodes that the solver should search.</listitem> + // <listitem>CpuTime : a Scalar, containing the Maximum amount of CPU Time that the solver should take.</listitem> + // <listitem>AllowableGap : a Scalar, to stop the tree search when the gap between the objective value of the best known solution is reached.</listitem> + // <listitem>MaxIter : a Scalar, containing the Maximum Number of Iteration that the solver should take.</listitem> + // <listitem>gradobj : a string, to turn on or off the user supplied objective gradient.</listitem> + // <listitem>hessian : a Scalar, to turn on or off the user supplied objective hessian.</listitem> + // <listitem>Default Values : options = list('integertolerance',1d-06,'maxnodes',2147483647,'cputime',1d10,'allowablegap',0,'maxiter',2147483647,'gradobj',"off",'hessian',"off")</listitem> + // </itemizedlist> + // The objective function must have header : + // <programlisting> + // F = fun(x) + // </programlisting> + // where x is a n x 1 matrix of doubles and F is a m x 1 matrix of doubles where m is the total number of objective functions inside F. + // On input, the variable x contains the current point and, on output, the variable F must contain the objective function values. + // + // By default, the gradient options for intfminimax are turned off and and intfmincon does the gradient opproximation of objective function. In case the GradObj option is off and GradConstr option is on, intfminimax approximates Objective function gradient using numderivative toolbox. + // + // If we can provide exact gradients, we should do so since it improves the convergence speed of the optimization algorithm. + // + // + // The exitflag allows to know the status of the optimization which is given back by Ipopt. + // <itemizedlist> + // <listitem>exitflag=0 : Optimal Solution Found </listitem> + // <listitem>exitflag=1 : InFeasible Solution.</listitem> + // <listitem>exitflag=2 : Objective Function is Continuous Unbounded.</listitem> + // <listitem>exitflag=3 : Limit Exceeded.</listitem> + // <listitem>exitflag=4 : User Interrupt.</listitem> + // <listitem>exitflag=5 : MINLP Error.</listitem> + // </itemizedlist> + // + // For more details on exitflag see the ipopt documentation, go to http://www.coin-or.org/bonmin/ + // + // Examples + // // A basic case : + // // we provide only the objective function and the nonlinear constraint + // // function + // function f = myfun(x) + // f(1)= 2*x(1)^2 + x(2)^2 - 48*x(1) - 40*x(2) + 304; //Objectives + // f(2)= -x(1)^2 - 3*x(2)^2; + // f(3)= x(1) + 3*x(2) -18; + // f(4)= -x(1) - x(2); + // f(5)= x(1) + x(2) - 8; + // endfunction + // // The initial guess + // x0 = [0.1,0.1]; + // // The expected solution : only 4 digits are guaranteed + // xopt = [4 4] + // fopt = [0 -64 -2 -8 0] + // intcon = [1] + // maxfopt = 0 + // // Run fminimax + // [x,fval,maxfval,exitflag] = intfminimax(myfun, x0,intcon) + // // Press ENTER to continue + // + // Examples + // // A case where we provide the gradient of the objective + // // functions and the Jacobian matrix of the constraints. + // // The objective function and its gradient + // function [f,G] = myfun(x) + // f(1)= 2*x(1)^2 + x(2)^2 - 48*x(1) - 40*x(2) + 304; + // f(2)= -x(1)^2 - 3*x(2)^2; + // f(3)= x(1) + 3*x(2) -18; + // f(4)= -x(1) - x(2); + // f(5)= x(1) + x(2) - 8; + // G = [ 4*x(1) - 48, -2*x(1), 1, -1, 1; + // 2*x(2) - 40, -6*x(2), 3, -1, 1; ]' + // endfunction + // // The nonlinear constraints + // function [c,ceq,DC,DCeq] = confun(x) + // // Inequality constraints + // c = [1.5 + x(1)*x(2) - x(1) - x(2), -x(1)*x(2) - 10] + // // No nonlinear equality constraints + // ceq=[] + // DC= [x(2)-1, -x(2); + // x(1)-1, -x(1)]' + // DCeq = []' + // endfunction + // // Test with both gradient of objective and gradient of constraints + // minimaxOptions = list("GradObj","on","GradCon","on"); + // // The initial guess + // x0 = [0,10]; + // intcon = [2] + // // Run intfminimax + // [x,fval,maxfval,exitflag] = intfminimax(myfun,x0,intcon,[],[],[],[],[],[], confun, minimaxOptions) + // Authors + // Harpreet Singh + + // Check number of input and output arguments + [minmaxLhs,minmaxRhs] = argn() + Checkrhs("fminimax", minmaxRhs, [2 3 5 7 9 10 11]) + Checklhs("fminimax", minmaxLhs, 1:7) + + // Proper initialisation of objective function + minmaxObjfun = varargin(1) + Checktype("fminimax", minmaxObjfun, "minmaxObjfun", 1, "function") + + // Proper initialisation of starting point + minmaxStartpoint = varargin(2) + Checktype("fminimax", minmaxStartpoint, "minmaxStartpoint", 2, "constant") + + minmaxNumvar = size(minmaxStartpoint,"*") + Checkvector("fminimax", minmaxStartpoint, "minmaxStartpoint", 2, minmaxNumvar) + minmaxStartpoint = minmaxStartpoint(:) + + if(minmaxRhs < 3) then // if A and b are not provided, declare as empty + intcon = 0; + else + intcon = varargin(3); + end + + // Proper initialisation of A and b + if(minmaxRhs < 4) then // if A and b are not provided, declare as empty + minmaxA = [] + minmaxB = [] + else + minmaxA = varargin(4) + minmaxB = varargin(5) + end + + Checktype("fminimax", minmaxA, "A", 4, "constant") + Checktype("fminimax", minmaxB, "b", 5, "constant") + + // Check if A and b of proper dimensions + if(minmaxA <> [] & minmaxB == []) then + errmsg = msprintf(gettext("%s: Incompatible input arguments #%d and #%d: matrix A is empty, but the column vector b is not empty"), "fminimax", 4, 5) + error(errmsg) + end + + if(minmaxA == [] & minmaxB <> []) then + errmsg = msprintf(gettext("%s: Incompatible input arguments #%d and #%d: matrix A is not empty, but the column vector b is empty"), "fminimax", 4, 5) + error(errmsg) + end + + minmaxNumrowA = size(minmaxA,"r") + if(minmaxA <> []) then + Checkdims("fminimax", minmaxA, "A", 4, [minmaxNumrowA minmaxNumvar]) + Checkvector("fminimax", minmaxB, "b", 5, minmaxNumrowA) + minmaxB = minmaxB(:) + end + + // Proper initialisation of Aeq and beq + if(minmaxRhs < 6) then // if Aeq and beq are not provided, declare as empty + minmaxAeq = [] + minmaxBeq = [] + else + minmaxAeq = varargin(6) + minmaxBeq = varargin(7) + end + + Checktype("fminimax", minmaxAeq, "Aeq", 6, "constant") + Checktype("fminimax", minmaxBeq, "beq", 7, "constant") + + // Check if Aeq and beq of proper dimensions + if(minmaxAeq <> [] & minmaxBeq == []) then + errmsg = msprintf(gettext("%s: Incompatible input arguments #%d and #%d: matrix Aeq is empty, but the column vector beq is not empty"), "fminimax", 6, 7) + error(errmsg) + end + + if(minmaxAeq == [] & minmaxBeq <> []) then + errmsg = msprintf(gettext("%s: Incompatible input arguments #%d and #%d: matrix Aeq is not empty, but the column vector beq is empty"), "fminimax", 6, 7) + error(errmsg) + end + + minmaxNumrowAeq = size(minmaxAeq,"r") + if(minmaxAeq <> []) then + Checkdims("fminimax", minmaxAeq, "Aeq", 6, [minmaxNumrowAeq minmaxNumvar]) + Checkvector("fminimax", minmaxBeq, "beq", 7, minmaxNumrowAeq) + minmaxBeq = minmaxBeq(:) + end + + // Proper initialisation of minmaxLb and minmaxUb + if(minmaxRhs < 6) then // if minmaxLb and minmaxUb are not provided, declare as empty + minmaxLb = [] + minmaxUb = [] + else + minmaxLb = varargin(6) + minmaxUb = varargin(7) + end + + Checktype("fminimax", minmaxLb, "lb", 6, "constant") + Checktype("fminimax", minmaxUb, "ub", 7, "constant") + + // Check dimensions of minmaxLb and minmaxUb + if(minmaxLb <> []) then + Checkvector("fminimax", minmaxLb, "lb", 8, minmaxNumvar) + minmaxLb = minmaxLb(:) + end + + if(minmaxUb <> []) then + Checkvector("fminimax", minmaxUb, "ub", 9, minmaxNumvar) + minmaxUb = minmaxUb(:) + end + + // Proper Initialisation of minmaxNonlinfun + if(minmaxRhs < 10) then // if minmaxNonlinfun is not provided, declare as empty + minmaxNonlinfun = [] + else + minmaxNonlinfun = varargin(10) + end + if(minmaxNonlinfun<>[]) then + Checktype("fminimax", minmaxNonlinfun, "nonlinfun", 10, "function") + end + + //To check, Whether minimaxOptions is been entered by user + if ( minmaxRhs<11 ) then + minmaxUserOptions = list(); + else + minmaxUserOptions = varargin(11); //Storing the 3rd Input minmaxUserOptionseter in intermediate list named 'minmaxUserOptions' + end + + //If minimaxOptions is entered then checking its type for 'list' + if (type(minmaxUserOptions) ~= 15) then + errmsg = msprintf(gettext("%s: minimaxOptions (10th parameter) should be a list"), "fminimax"); + error(errmsg); + end + + //If minimaxOptions is entered then checking whether even number of entires are entered + if (modulo(size(minmaxUserOptions),2)) then + errmsg = msprintf(gettext("%s: Size of minimaxOptions (list) should be even"), "fminimax"); + error(errmsg); + end + + /////////////// To check integer ////////////////////// + for i=1:size(intcon,1) + if(intcon(i)>minmaxNumvar) then + errmsg = msprintf(gettext("%s: The values inside intcon should be less than the number of variables"), "intfminimax"); + error(errmsg); + end + + if (intcon(i)<0) then + errmsg = msprintf(gettext("%s: The values inside intcon should be greater than 0 "), "intfminimax"); + error(errmsg); + end + + if(modulo(intcon(i),1)) then + errmsg = msprintf(gettext("%s: The values inside intcon should be an integer "), "intfminimax"); + error(errmsg); + end + end + + //If minimaxOptions is entered then checking its type for 'list' + if (type(minmaxUserOptions) ~= 15) then + errmsg = msprintf(gettext("%s: minimaxOptions (10th parameter) should be a list"), "intfminimax"); + error(errmsg); + end + + //If minimaxOptions is entered then checking whether even number of entires are entered + if (modulo(size(minmaxUserOptions),2)) then + errmsg = msprintf(gettext("%s: Size of minimaxOptions (list) should be even"), "intfminimax"); + error(errmsg); + end + +minmaxoptions = list('integertolerance',1d-06,'maxnodes',2147483647,'cputime',1d10,'allowablegap',0,'maxiter',2147483647,'gradobj',"off",'gradcon',"off") + + //Pushing minmaxUserOptions into default value + + for i = 1:(size(minmaxUserOptions))/2 + select convstr(minmaxUserOptions(2*i-1),'l') + case 'integertolerance' then + Checktype("intfminimax_options", minmaxUserOptions(2*i), minmaxUserOptions(2*i-1), 2*i, "constant"); + minmaxoptions(2) = minmaxUserOptions(2*i); + case 'maxnodes' then + Checktype("intfminimax_options", minmaxUserOptions(2*i), minmaxUserOptions(2*i-1), 2*i, "constant"); + minmaxoptions(4) = minmaxUserOptions(2*i); + case 'cputime' then + Checktype("intfminimax_options", minmaxUserOptions(2*i), minmaxUserOptions(2*i-1), 2*i, "constant"); + minmaxoptions(6) = minmaxUserOptions(2*i); + case 'allowablegap' then + Checktype("intfminimax_options", minmaxUserOptions(2*i), minmaxUserOptions(2*i-1), 2*i, "constant"); + minmaxoptions(8) = minmaxUserOptions(2*i); + case 'maxiter' then + Checktype("intfminimax_options", minmaxUserOptions(2*i), minmaxUserOptions(2*i-1), 2*i, "constant"); + minmaxoptions(10) = minmaxUserOptions(2*i); + case 'gradobj' then + Checktype("intfminimax_options", minmaxUserOptions(2*i), minmaxUserOptions(2*i-1), 2*i, "string"); + if(convstr(minmaxUserOptions(2*i),'l') == "on") then + minmaxoptions(12) = "on" + elseif(convstr(minmaxUserOptions(2*i),'l') == "off") then + minmaxoptions(12) = "off" + else + error(999, 'Unknown string passed in gradobj.'); + end + case 'gradcon' then + Checktype("intfminimax_options", minmaxUserOptions(2*i), minmaxUserOptions(2*i-1), 2*i, "string"); + if(convstr(minmaxUserOptions(2*i),'l') == "on") then + minmaxoptions(14) = "on" + elseif(convstr(minmaxUserOptions(2*i),'l') == "off") then + minmaxoptions(14) = "off" + else + error(999, 'Unknown string passed in gradcon.'); + end + else + error(999, 'Unknown string argument passed.'); + end + end + + // Reformulating the problem fminimax to fmincon + minmaxObjfunval = minmaxObjfun(minmaxStartpoint) + minmaxStartpoint(minmaxNumvar+1) = max(minmaxObjfunval) + + if(minmaxA <> []) then + minmaxA = [minmaxA, zeros(minmaxNumrowA,1)] + end + if(minmaxAeq <> []) then + minmaxAeq = [minmaxAeq, zeros(minmaxNumrowAeq,1)] + end + if(minmaxLb <> []) then + minmaxLb(minmaxNumvar+1) = -%inf + end + if(minmaxUb <> []) then + minmaxUb(minmaxNumvar+1) = +%inf + end + + // function handle defining the additional inequalities + function temp = minmaxAddIneq(z) + temp = minmaxObjfun(z) - z(minmaxNumvar+1) + endfunction + + // function handle defining minmaxNonlinfun derivative using numderivative + function [dc,dceq] = minmaxNonlinDer(z) + // function handle extracting c and ceq components from minmaxNonlinfun + function foo = minmaxC(z) + [foo,tmp1] = minmaxNonlinfun(z) + foo = foo' + endfunction + + function foo = minmaxCEQ(z) + [tmp1,foo] = minmaxNonlinfun(z) + foo = foo' + endfunction + + dc = numderivative(minmaxC,z) + dceq = numderivative(minmaxCEQ,z) + endfunction + + // function handle defining new objective function + function newfunc = newObjfun(z) + newfunc = z(minmaxNumvar+1) + endfunction + + // function handle defining new minmaxNonlinfun function + function [nc,nceq,dnc,dnceq] = newNonlinfun(z) + dnc = []; + dnceq = []; + nc = []; + nceq= []; + if (minmaxNonlinfun<>[]) then + [nc,nceq] = minmaxNonlinfun(z) + end + // add inequalities of the form Fi(x) - y <= 0 + tmp = [minmaxObjfun(z) - z(minmaxNumvar+1)]' + nc = [nc, tmp] + if(options(14) =="on") then + [temp1,temp2,dnc, dnceq] = minmaxNonlinfun(z) + dnc = [dnc, zeros(size(dnc,'r'),1)] + dnceq = [dnceq, zeros(size(dnceq,'r'),1)] + else + // else use numderivative method to calculate gradient of constraints + if (minmaxNonlinfun<>[]) then + [dnc, dnceq] = minmaxNonlinDer(z) + end + end + + if(options(12) =="on") then + [temp,derObjfun] = minmaxObjfun(z); + mderObjfun = [derObjfun, -1*ones(size(derObjfun,'r'),1)]; + dnc = [dnc; mderObjfun]; + else + // else use numderivative to calculate gradient of set of obj functions + derObjfun = numderivative(minmaxAddIneq,z) + dnc = [dnc; derObjfun] + end + endfunction + + if( minmaxoptions(12)=="on"| minmaxoptions(12)=="on") then + options(14)="on"; + end + + minmaxoptions(12)="off"; + [x,fval,exitflag,gradient,hessian] = ... + intfmincon(newObjfun,minmaxStartpoint,intcon,minmaxA,minmaxB,minmaxAeq,minmaxBeq,minmaxLb,minmaxUb,newNonlinfun,minmaxoptions) + + x = x(1:minmaxNumvar) + fval = minmaxObjfun(x) + maxfval = max(fval) +endfunction diff --git a/macros/intfminunc.bin b/macros/intfminunc.bin Binary files differnew file mode 100644 index 0000000..3802be4 --- /dev/null +++ b/macros/intfminunc.bin diff --git a/macros/intfminunc.sci b/macros/intfminunc.sci new file mode 100644 index 0000000..2fe73fb --- /dev/null +++ b/macros/intfminunc.sci @@ -0,0 +1,350 @@ +// Copyright (C) 2015 - IIT Bombay - FOSSEE +// +// This file must be used under the terms of the CeCILL. +// This source file is licensed as described in the file COPYING, which +// you should have received as part of this distribution. The terms +// are also available at +// http://www.cecill.info/licences/Licence_CeCILL_V2-en.txt +// Author: Harpreet Singh, Pranav Deshpande and Akshay Miterani +// Organization: FOSSEE, IIT Bombay +// Email: toolbox@scilab.in + +function [xopt,fopt,exitflag,gradient,hessian] = intfminunc (varargin) + // Solves an unconstrainted multi-variable mixed integer non linear programming optimization problem + // + // Calling Sequence + // xopt = intfminunc(f,x0) + // xopt = intfminunc(f,x0,intcon) + // xopt = intfminunc(f,x0,intcon,options) + // [xopt,fopt] = intfminunc(.....) + // [xopt,fopt,exitflag]= intfminunc(.....) + // [xopt,fopt,exitflag,gradient,hessian]= intfminunc(.....) + // + // Parameters + // f : a function, representing the objective function of the problem + // x0 : a vector of doubles, containing the starting of variables. + // intcon : a vector of integers, represents which variables are constrained to be integers + // options: a list, containing the option for user to specify. See below for details. + // xopt : a vector of doubles, the computed solution of the optimization problem. + // fopt : a scalar of double, the function value at x. + // exitflag : a scalar of integer, containing the flag which denotes the reason for termination of algorithm. See below for details. + // gradient : a vector of doubles, containing the Objective's gradient of the solution. + // hessian : a matrix of doubles, containing the Objective's hessian of the solution. + // + // Description + // Search the minimum of a multi-variable mixed integer non linear programming unconstrained optimization problem specified by : + // Find the minimum of f(x) such that + // + // <latex> + // \begin{eqnarray} + // &\mbox{min}_{x} + // & f(x) + // & x_i \in \!\, \mathbb{Z}, i \in \!\, I + // \end{eqnarray} + // </latex> + // + // The routine calls Bonmin for solving the Un-constrained Optimization problem, Bonmin is a library written in C++. + // + // The options allows the user to set various parameters of the Optimization problem. + // It should be defined as type "list" and contains the following fields. + // <itemizedlist> + // <listitem>Syntax : options= list("IntegerTolerance", [---], "MaxNodes", [---], "CpuTime", [---], "AllowableGap", [---], "MaxIter", [---]);</listitem> + // <listitem>IntegerTolerance : a Scalar, containing the Integer tolerance value that the solver should take.</listitem> + // <listitem>MaxNodes : a Scalar, containing the maximum nodes that the solver should make.</listitem> + // <listitem>MaxIter : a Scalar, containing the Maximum Number of Iteration that the solver should take.</listitem> + // <listitem>AllowableGap : a Scalar, containing the allowable gap value that the solver should take.</listitem> + // <listitem>CpuTime : a Scalar, containing the Maximum amount of CPU Time that the solver should take.</listitem> + // <listitem>gradobj : a string, to turn on or off the user supplied objective gradient.</listitem> + // <listitem>hessian : a Scalar, to turn on or off the user supplied objective hessian.</listitem> + // <listitem>Default Values : options = list('integertolerance',1d-06,'maxnodes',2147483647,'cputime',1d10,'allowablegap',0,'maxiter',2147483647,'gradobj',"off",'hessian',"off")</listitem> + // </itemizedlist> + // + // The exitflag allows to know the status of the optimization which is given back by Bonmin. + // <itemizedlist> + // <listitem>exitflag=0 : Optimal Solution Found. </listitem> + // <listitem>exitflag=1 : InFeasible Solution.</listitem> + // <listitem>exitflag=2 : Output is Continuous Unbounded.</listitem> + // <listitem>exitflag=3 : Limit Exceeded.</listitem> + // <listitem>exitflag=4 : User Interrupt.</listitem> + // <listitem>exitflag=5 : MINLP Error.</listitem> + // </itemizedlist> + // + // For more details on exitflag see the Bonmin page, go to http://www.coin-or.org/Bonmin + // + // Examples + // //Find x in R^2 such that it minimizes the Rosenbrock function + // //f = 100*(x2 - x1^2)^2 + (1-x1)^2 + // //Objective function to be minimised + // function y= f(x) + // y= 100*(x(2) - x(1)^2)^2 + (1-x(1))^2; + // endfunction + // //Starting point + // x0=[-1,2]; + // intcon = [2] + // //Options + // options=list("MaxIter", [1500], "CpuTime", [500]); + // //Calling + // [xopt,fopt,exitflag,gradient,hessian]=intfminunc(f,x0,intcon,options) + // // Press ENTER to continue + // + // Examples + // //Find x in R^2 such that the below function is minimum + // //f = x1^2 + x2^2 + // //Objective function to be minimised + // function y= f(x) + // y= x(1)^2 + x(2)^2; + // endfunction + // //Starting point + // x0=[2,1]; + // intcon = [1]; + // [xopt,fopt]=intfminunc(f,x0,intcon) + // // Press ENTER to continue + // + // Examples + // //The below problem is an unbounded problem: + // //Find x in R^2 such that the below function is minimum + // //f = - x1^2 - x2^2 + // //Objective function to be minimised + // function [y,g,h] = f(x) + // y = -x(1)^2 - x(2)^2; + // g = [-2*x(1),-2*x(2)]; + // h = [-2,0;0,-2]; + // endfunction + // //Starting point + // x0=[2,1]; + // intcon = [1] + // options = list("gradobj","ON","hessian","on"); + // [xopt,fopt,exitflag,gradient,hessian]=intfminunc(f,x0,intcon,options) + + //To check the number of input and output arguments + [lhs , rhs] = argn(); + + //To check the number of arguments given by the user + if ( rhs<2 | rhs>4 ) then + errmsg = msprintf(gettext("%s: Unexpected number of input arguments : %d provided while should be int [2 3 4] "), "intfminunc", rhs); + error(errmsg); + end + + //Storing the 1st and 2nd Input Parameters + fun = varargin(1); + x0 = varargin(2); + //To add intcon + intcon=[]; + if ( rhs >=3 ) then + intcon = varargin(3); + end + + param = list(); + //To check whether options has been entered by user + if ( rhs>=4 ) then + param =varargin(4); + end + + nbvar = size(x0,"*"); + + ///////////////// To check whether the Input arguments ///////////////// + Checktype("intfminunc", fun, "fun", 1, "function"); + Checktype("intfminunc", x0, "x0", 2, "constant"); + Checktype("intfminunc", intcon, "intcon", 3, "constant"); + Checktype("intfminunc", param, "options", 4, "list"); + + ///////////////// To check x0 ///////////////// + Checkvector("intfminunc", x0, "x0", 2, nbvar) + x0 = x0(:); + Checkvector("intfminbnd", intcon, "intcon", 3, size(intcon,"*")) + intcon = intcon(:); + + //Error Checks for intcon + for i=1:size(intcon,1) + if(intcon(i)>nbvar) then + errmsg = msprintf(gettext("%s: The values inside intcon should be less than the number of variables"), "intfminunc"); + error(errmsg); + end + + if (intcon(i)<0) then + errmsg = msprintf(gettext("%s: The values inside intcon should be greater than 0 "), "intfminunc"); + error(errmsg); + end + + if(modulo(intcon(i),1)) then + errmsg = msprintf(gettext("%s: The values inside intcon should be an integer "), "intfminunc"); + error(errmsg); + end + end + + //If options has been entered, then check whether an even number of entires has been entered + if (modulo(size(param),2)) then + errmsg = msprintf(gettext("%s: Size of parameters should be even"), "intfminunc"); + error(errmsg); + end + + intconSize = length(intcon); + + options = list('integertolerance',1d-06,'maxnodes',2147483647,'cputime',1d10,'allowablegap',0,'maxiter',2147483647,'gradobj',"off",'hessian', "off") + + //Pushing param into default value + + for i = 1:(size(param))/2 + select convstr(param(2*i-1),'l') + case 'integertolerance' then + Checktype("intfminbnd_options", param(2*i), param(2*i-1), 2*i, "constant"); + options(2) = param(2*i); + case 'maxnodes' then + Checktype("intfminbnd_options", param(2*i), param(2*i-1), 2*i, "constant"); + options(4) = options(2*i); + case 'cputime' then + Checktype("intfminbnd_options", param(2*i), param(2*i-1), 2*i, "constant"); + options(6) = options(2*i); + case 'allowablegap' then + Checktype("intfminbnd_options", param(2*i), param(2*i-1), 2*i, "constant"); + options(8) = options(2*i); + case 'maxiter' then + Checktype("intfminbnd_options", param(2*i), param(2*i-1), 2*i, "constant"); + options(10) = options(2*i); + case 'gradobj' then + Checktype("intfminbnd_options", param(2*i), param(2*i-1), 2*i, "string"); + if(convstr(param(2*i),'l') == "on") then + options(12) = "on" + elseif(convstr(param(2*i),'l') == "off") then + options(12) = "off" + else + error(999, 'Unknown string passed in gradobj.'); + end + case 'hessian' then + Checktype("intfminbnd_options", param(2*i), param(2*i-1), 2*i, "string"); + if(convstr(param(2*i),'l') == "on") then + options(14) = "on"; + elseif(convstr(param(2*i),'l') == "off") then + options(14) = "off"; + else + error(999, 'Unknown string passed in hessian.'); + end + else + error(999, 'Unknown string argument passed.'); + end + end + + ///////////////// Functions Check ///////////////// + + //To check the match between fun (1st Parameter) and x0 (2nd Parameter) + if(execstr('init=fun(x0)','errcatch')==21) then + errmsg = msprintf(gettext("%s: Objective function and x0 did not match"), "intfminunc"); + error(errmsg); + end + + if(options(12) == "on") then + + if(execstr('[grad_y,grad_dy]=fun(x0)','errcatch')==59) then + errmsg = msprintf(gettext("%s: Gradient of objective function is not provided"), "intfminunc"); + error(errmsg); + end + + Checkvector("intfminunc_options", grad_dy, "dy", 12, nbvar); + end + + if(options(14) == "on") then + + if(execstr('[hessian_y,hessian_dy,hessian]=fun(x0)','errcatch')==59) then + errmsg = msprintf(gettext("%s: Gradient of objective function is not provided"), "intfminunc"); + error(errmsg); + end + + if ( ~isequal(size(hessian),[nbvar nbvar]) ) then + errmsg = msprintf(gettext("%s: Size of hessian should be nbvar X nbvar"), "intfminunc"); + error(errmsg); + end + end + + //Converting the User defined Objective function into Required form (Error Detectable) + function [y,check] = _f(x) + try + y=fun(x) + [y,check] = checkIsreal(y) + catch + y=0; + check=1; + end + endfunction + + //Defining an inbuilt Objective gradient function + function [dy,check] = _gradf(x) + if (options(12) =="on") then + try + [y,dy]=fun(x); + [dy,check] = checkIsreal(dy); + catch + dy = 0; + check=1; + end + else + try + dy=numderivative(fun,x); + [dy,check] = checkIsreal(dy); + catch + dy=0; + check=1; + end + end + endfunction + + //Defining a function to calculate Hessian if the respective user entry is OFF + function [hessy,check]=_gradhess(x) + if (options(14) == "on") then + try + [obj,dy,hessy] = fun(x) + [hessy,check] = checkIsreal(hessy) + catch + hessy = 0; + check=1; + end + else + try + [dy,hessy]=numderivative(fun,x) + [hessy,check] = checkIsreal(hessy) + catch + hessy=0; + check=1; + end + end + endfunction + + //Calling the bonmin function for solving the above problem + [xopt,fopt,exitflag] = inter_fminunc(_f,_gradf,_gradhess,x0,intcon,options,intconSize,nbvar); + + //In the cases of the problem not being solved, return NULL to the output matrices + if( exitflag~=0 & exitflag~=3 ) then + gradient = []; + hessian = []; + else + [ gradient, hessian] = numderivative(_f, xopt, [], [], "blockmat"); + end + + //To print output message + select exitflag + case 0 then + printf("\nOptimal Solution Found.\n"); + case 1 then + printf("\nInFeasible Solution.\n"); + case 2 then + printf("\nObjective Function is Continuous Unbounded.\n"); + case 3 then + printf("\Limit Exceeded.\n"); + case 4 then + printf("\nUser Interrupt.\n"); + case 5 then + printf("\nMINLP Error.\n"); + else + printf("\nInvalid status returned. Notify the Toolbox authors\n"); + break; + end +endfunction + +function [y, check] = checkIsreal(x) + if ((~isreal(x))) then + y = 0 + check=1; + else + y = x; + check=0; + end +endfunction diff --git a/macros/intqpipopt.bin b/macros/intqpipopt.bin Binary files differnew file mode 100644 index 0000000..c0ef46f --- /dev/null +++ b/macros/intqpipopt.bin diff --git a/macros/intqpipopt.sci b/macros/intqpipopt.sci new file mode 100644 index 0000000..f45d9cb --- /dev/null +++ b/macros/intqpipopt.sci @@ -0,0 +1,461 @@ +// Copyright (C) 2016 - IIT Bombay - FOSSEE +// +// This file must be used under the terms of the CeCILL. +// This source file is licensed as described in the file COPYING, which +// you should have received as part of this distribution. The terms +// are also available at +// http://www.cecill.info/licences/Licence_CeCILL_V2-en.txt +// Author: Pranav Deshpande and Akshay Miterani +// Organization: FOSSEE, IIT Bombay +// Email: toolbox@scilab.in + +function varargout = intqpipopt (varargin) + // Solves a linear quadratic problem. + // + // Calling Sequence + // xopt = intqpipopt(H,f) + // xopt = intqpipopt(H,f,intcon) + // xopt = intqpipopt(H,f,intcon,A,b) + // xopt = intqpipopt(H,f,intcon,A,b,Aeq,beq) + // xopt = intqpipopt(H,f,intcon,A,b,Aeq,beq,lb,ub) + // xopt = intqpipopt(H,f,intcon,A,b,Aeq,beq,lb,ub,x0) + // xopt = intqpipopt(H,f,intcon,A,b,Aeq,beq,lb,ub,x0,options) + // xopt = intqpipopt(H,f,intcon,A,b,Aeq,beq,lb,ub,x0,options,"file_path") + // [xopt,fopt,exitflag,output] = intqpipopt( ... ) + // + // Parameters + // H : a symmetric matrix of double, represents coefficients of quadratic in the quadratic problem. + // f : a vector of double, represents coefficients of linear in the quadratic problem + // intcon : a vector of integers, represents which variables are constrained to be integers + // A : a matrix of double, represents the linear coefficients in the inequality constraints A⋅x ≤ b. + // b : a vector of double, represents the linear coefficients in the inequality constraints A⋅x ≤ b. + // Aeq : a matrix of double, represents the linear coefficients in the equality constraints Aeq⋅x = beq. + // beq : a vector of double, represents the linear coefficients in the equality constraints Aeq⋅x = beq. + // lb : a vector of double, contains lower bounds of the variables. + // ub : a vector of double, contains upper bounds of the variables. + // x0 : a vector of double, contains initial guess of variables. + // options : a list containing the parameters to be set. + // file_path : path to bonmin opt file if used. + // xopt : a vector of double, the computed solution of the optimization problem. + // fopt : a double, the value of the function at x. + // exitflag : The exit status. See below for details. + // output : The structure consist of statistics about the optimization. See below for details. + // + // Description + // Search the minimum of a constrained linear quadratic optimization problem specified by : + // + // <latex> + // \begin{eqnarray} + // &\mbox{min}_{x} + // & 1/2⋅x^T⋅H⋅x + f^T⋅x \\ + // & \text{subject to} & A⋅x \leq b \\ + // & & Aeq⋅x = beq \\ + // & & lb \leq x \leq ub \\ + // & & x_i \in \!\, \mathbb{Z}, i \in \!\, intcon\\ + // \end{eqnarray} + // </latex> + // + // The routine calls Bonmin for solving the quadratic problem, Bonmin is a library written in C++. + // + // The exitflag allows to know the status of the optimization which is given back by Bonmin. + // <itemizedlist> + // <listitem>exitflag=0 : Optimal Solution Found. </listitem> + // <listitem>exitflag=1 : InFeasible Solution.</listitem> + // <listitem>exitflag=2 : Output is Continuous Unbounded.</listitem> + // <listitem>exitflag=3 : Limit Exceeded.</listitem> + // <listitem>exitflag=4 : User Interrupt.</listitem> + // <listitem>exitflag=5 : MINLP Error.</listitem> + // </itemizedlist> + // + // For more details on exitflag see the Bonmin page, go to http://www.coin-or.org/Bonmin + // + // The output data structure contains detailed informations about the optimization process. + // It has type "struct" and contains the following fields. + // <itemizedlist> + // <listitem>output.constrviolation: The max-norm of the constraint violation.</listitem> + // </itemizedlist> + // + // Examples + // H = [1 -1; -1 2]; + // f = [-2; -6]; + // A = [1 1; -1 2; 2 1]; + // b = [2; 2; 3]; + // lb=[0,0]; + // ub=[%inf, %inf]; + // intcon = [1 2]; + //[xopt,fopt,status,output]=intqpipopt(H,f,intcon,A,b,[],[],lb,ub) + // + // Examples + // //Find x in R^6 such that: + // Aeq= [1,-1,1,0,3,1; + // -1,0,-3,-4,5,6; + // 2,5,3,0,1,0]; + // beq=[1; 2; 3]; + // A= [0,1,0,1,2,-1; + // -1,0,2,1,1,0]; + // b = [-1; 2.5]; + // lb=[-1000; -10000; 0; -1000; -1000; -1000]; + // ub=[10000; 100; 1.5; 100; 100; 1000]; + // x0 = repmat(0,6,1); + // param = list("MaxIter", 300, "CpuTime", 100); + // //and minimize 0.5*x'*H*x + f'*x with + // f=[1; 2; 3; 4; 5; 6]; H=eye(6,6); + // intcon = [2 4]; + // [xopt,fopt,exitflag,output]=intqpipopt(H,f,intcon,A,b,Aeq,beq,lb,ub,x0,param) + // Authors + // Akshay Miterani and Pranav Deshpande + + //To check the number of input and output argument + [lhs , rhs] = argn(); + + //To check the number of argument given by user + if ( rhs < 2 | rhs == 4 | rhs == 6 | rhs == 8 | rhs > 12 ) then + errmsg = msprintf(gettext("%s: Unexpected number of input arguments : %d provided while should be in the set of [2 3 5 7 9 10 11 12]"), "intqpipopt", rhs); + error(errmsg) + end + + //To check the number of output arguments + if lhs > 4 then + errmsg = msprintf(gettext("%s: Unexpected number of output arguments: %d provided while should be in the set of [1 2 3 4]"), "intqpipopt", lhs); + end + + H = []; + f = []; + intcon = []; + A = []; + b = []; + Aeq = []; + beq = []; + lb = []; + ub = []; + options = list(); + bonmin_options_file = ''; + + H = varargin(1); + f = varargin(2); + nbVar = size(H,1); + + if(nbVar == 0) then + errmsg = msprintf(gettext("%s: Cannot determine the number of variables because input objective coefficients is empty"), "intqpipopt"); + error(errmsg); + end + + if ( rhs>=3 ) then + intcon=varargin(3); + end + + if ( rhs<=3 ) then + A = [] + b = [] + else + A = varargin(4); + b = varargin(5); + end + + if ( rhs<6 ) then + Aeq = [] + beq = [] + else + Aeq = varargin(6); + beq = varargin(7); + end + + if ( rhs<8 ) then + lb = repmat(-%inf,nbVar,1); + ub = repmat(%inf,nbVar,1); + else + lb = varargin(8); + ub = varargin(9); + end + + + if ( rhs<10 | size(varargin(10)) ==0 ) then + x0 = repmat(0,nbVar,1) + else + x0 = varargin(10); + end + + if (rhs<11|size(varargin(11))==0) then + options = list(); + else + options = varargin(11); + end + + if(rhs==12) then + bonmin_options_file=varargin(12); + end + + if (size(lb,2)==0) then + lb = repmat(-%inf,nbVar,1); + end + + if (size(ub,2)==0) then + ub = repmat(%inf,nbVar,1); + end + + if (size(f,2)==0) then + f = repmat(0,nbVar,1); + end + + //Check type of variables + Checktype("intqpipopt", H, "H", 1, "constant") + Checktype("intqpipopt", f, "f", 2, "constant") + Checktype("intqpipopt", intcon, "intcon", 3, "constant") + Checktype("intqpipopt", A, "A", 4, "constant") + Checktype("intqpipopt", b, "b", 5, "constant") + Checktype("intqpipopt", Aeq, "Aeq", 6, "constant") + Checktype("intqpipopt", beq, "beq", 7, "constant") + Checktype("intqpipopt", lb, "lb", 8, "constant") + Checktype("intqpipopt", ub, "ub", 9, "constant") + Checktype("intqpipopt", x0, "x0", 10, "constant") + Checktype("intqpipopt", options, "options", 11, "list") + Checktype("intqpipopt", bonmin_options_file, "bonmin_options_file", 12, "string") + + nbConInEq = size(A,1); + nbConEq = size(Aeq,1); + + // Check if the user gives row vector + // and Changing it to a column matrix + + if (size(f,2)== [nbVar]) then + f=f'; + end + + if (size(lb,2)== [nbVar]) then + lb = lb'; + end + + if (size(ub,2)== [nbVar]) then + ub = ub'; + end + + if (size(b,2)==nbConInEq) then + b = b'; + end + + if (size(beq,2)== nbConEq) then + beq = beq'; + end + + if (size(x0,2)== [nbVar]) then + x0=x0'; + end + + //Checking the H matrix which needs to be a symmetric matrix + if ( ~isequal(H,H')) then + errmsg = msprintf(gettext("%s: H is not a symmetric matrix"), "intqpipopt"); + error(errmsg); + end + + //Check the size of f which should equal to the number of variable + if ( size(f,1) ~= [nbVar]) then + errmsg = msprintf(gettext("%s: The number of rows and columns in H must be equal the number of elements of f"), "intqpipopt"); + error(errmsg); + end + + //Error Checks for intcon + + for i=1:size(intcon,2) + if(intcon(i)>nbVar) then + errmsg = msprintf(gettext("%s: The values inside intcon should be less than the number of variables"), "intqpipopt"); + error(errmsg); + end + + if (intcon(i)<0) then + errmsg = msprintf(gettext("%s: The values inside intcon should be greater than 0 "), "intqpipopt"); + error(errmsg); + end + + if(modulo(intcon(i),1)) then + errmsg = msprintf(gettext("%s: The values inside intcon should be an integer "), "intqpipopt"); + error(errmsg); + end + end + + //Check the size of inequality constraint which should be equal to the number of variables + if ( size(A,2) ~= nbVar & size(A,2) ~= 0) then + errmsg = msprintf(gettext("%s: The number of columns in A must be the same as the number of elements of f"), "intqpipopt"); + error(errmsg); + end + + //Check the size of equality constraint which should be equal to the number of variables + if ( size(Aeq,2) ~= nbVar & size(Aeq,2) ~= 0 ) then + errmsg = msprintf(gettext("%s: The number of columns in Aeq must be the same as the number of elements of f"), "intqpipopt"); + error(errmsg); + end + + //Check the size of Lower Bound which should be equal to the number of variables + if ( size(lb,1) ~= nbVar) then + errmsg = msprintf(gettext("%s: The Lower Bound is not equal to the number of variables"), "intqpipopt"); + error(errmsg); + end + + //Check the size of Upper Bound which should equal to the number of variables + if ( size(ub,1) ~= nbVar) then + errmsg = msprintf(gettext("%s: The Upper Bound is not equal to the number of variables"), "intqpipopt"); + error(errmsg); + end + + //Check the size of constraints of Lower Bound which should equal to the number of constraints + if ( size(b,1) ~= nbConInEq & size(b,1) ~= 0) then + errmsg = msprintf(gettext("%s: The number of rows in A must be the same as the number of elements of b"), "intqpipopt"); + error(errmsg); + end + + //Check the size of constraints of Upper Bound which should equal to the number of constraints + if ( size(beq,1) ~= nbConEq & size(beq,1) ~= 0) then + errmsg = msprintf(gettext("%s: The number of rows in Aeq must be the same as the number of elements of beq"), "intqpipopt"); + error(errmsg); + end + + //Check the size of initial of variables which should equal to the number of variables + if ( size(x0,1) ~= nbVar) then + warnmsg = msprintf(gettext("%s: Ignoring initial guess of variables as it is not equal to the number of variables"), "intqpipopt"); + warning(warnmsg); + x0 = repmat(0,nbVar,1); + end + + //Check if the user gives a matrix instead of a vector + + if ((size(f,1)~=1)& (size(f,2)~=1)) then + errmsg = msprintf(gettext("%s: f should be a vector"), "intqpipopt"); + error(errmsg); + end + + if (size(lb,1)~=1)& (size(ub,2)~=1) then + errmsg = msprintf(gettext("%s: Lower Bound should be a vector"), "intqpipopt"); + error(errmsg); + end + + if (size(ub,1)~=1)& (size(ub,2)~=1) then + errmsg = msprintf(gettext("%s: Upper Bound should be a vector"), "intqpipopt"); + error(errmsg); + end + + if (nbConInEq) then + if ((size(b,1)~=1)& (size(b,2)~=1)) then + errmsg = msprintf(gettext("%s: Constraint Lower Bound should be a vector"), "intqpipopt"); + error(errmsg); + end + end + + if (nbConEq) then + if (size(beq,1)~=1)& (size(beq,2)~=1) then + errmsg = msprintf(gettext("%s: Constraint should be a vector"), "intqpipopt"); + error(errmsg); + end + end + + for i = 1:nbConInEq + if (b(i) == -%inf) then + errmsg = msprintf(gettext("%s: Value of b can not be negative infinity"), "intqpipopt"); + error(errmsg); + end + end + + for i = 1:nbConEq + if (beq(i) == -%inf) then + errmsg = msprintf(gettext("%s: Value of beq can not be negative infinity"), "intqpipopt"); + error(errmsg); + end + end + + for i = 1:nbVar + if(lb(i)>ub(i)) then + errmsg = msprintf(gettext("%s: Problem has inconsistent variable bounds"), "intqpipopt"); + error(errmsg); + end + end + + // Checking if the specified options file exists or not + if (~isfile(bonmin_options_file) & ~bonmin_options_file=='') then + error(999, 'The specified options file does not exist!'); + end + + // Validating Options + if (type(options) ~= 15) then + errmsg = msprintf(gettext("%s: Options should be a list "), "cbcintlinprog"); + error(errmsg); + end + + + if (modulo(size(options),2)) then + errmsg = msprintf(gettext("%s: Size of parameters should be even"), "cbcintlinprog"); + error(errmsg); + end + + //Pusing options as required to a double array + optval = []; + ifval =[]; + if length(options) == 0 then + optval = [0 0 0 0 0]; + ifval = [0 0 0 0 0]; + else + optval = [0 0 0 0 0]; + ifval = [0 0 0 0 0]; + for i=1:2:length(options) + select convstr(options(i),'l') + case 'integertolerance' then + ifval(1) = 1; + optval(1) = options(i+1); + case 'maxnodes' then + ifval(2) = 1; + optval(2) = options(i+1); + case 'cputime' then + ifval(3) = 1; + optval(3) = options(i+1); + case 'allowablegap' then + ifval(4) = 1; + optval(4) = options(i+1); + case 'maxiter' then + ifval(5) = 1; + optval(5) = options(i+1); + else + error(999, 'Unknown string argument passed.'); + end + end + end + + //Converting it into bonmin format + f = f'; + lb = lb'; + ub = ub'; + x0 = x0'; + conMatrix = [Aeq;A]; + nbCon = size(conMatrix,1); + conLB = [beq; repmat(-%inf,nbConInEq,1)]'; + conUB = [beq;b]'; + intcon = intcon' + intconSize = length(intcon); + xopt=[] + fopt=[] + status=[] + [xopt,fopt,status] = sci_intqpipopt(nbVar,nbCon,intconSize,H,f,intcon,conMatrix,conLB,conUB,lb,ub,x0,optval,ifval,bonmin_options_file); + xopt = xopt'; + + output.ConstrViolation = max([0;norm(Aeq*xopt-beq, 'inf');(lb'-xopt);(xopt-ub');(A*xopt-b)]); + + varargout(1) = xopt + varargout(2) = fopt + varargout(3) = status + varargout(4) = output + + select status + case 0 then + printf("\nOptimal Solution Found.\n"); + case 1 then + printf("\nInFeasible Solution.\n"); + case 2 then + printf("\nOutput is Continuous Unbounded.s\n"); + case 3 then + printf("\nLimit Exceeded.\n"); + case 4 then + printf("\nUser Interrupt.\n"); + case 5 then + printf("\nMINLP Error.\n"); + else + printf("\nInvalid status returned. Notify the Toolbox authors\n"); + break; + end + +endfunction diff --git a/macros/lib b/macros/lib Binary files differnew file mode 100644 index 0000000..37d0bcd --- /dev/null +++ b/macros/lib diff --git a/macros/linprog.bin b/macros/linprog.bin Binary files differnew file mode 100644 index 0000000..1d3a5aa --- /dev/null +++ b/macros/linprog.bin diff --git a/macros/lsqlin.bin b/macros/lsqlin.bin Binary files differnew file mode 100644 index 0000000..7baf7be --- /dev/null +++ b/macros/lsqlin.bin diff --git a/macros/lsqnonlin.bin b/macros/lsqnonlin.bin Binary files differnew file mode 100644 index 0000000..c90bc0d --- /dev/null +++ b/macros/lsqnonlin.bin diff --git a/macros/lsqnonneg.bin b/macros/lsqnonneg.bin Binary files differnew file mode 100644 index 0000000..32620d0 --- /dev/null +++ b/macros/lsqnonneg.bin diff --git a/macros/matrix_linprog.bin b/macros/matrix_linprog.bin Binary files differnew file mode 100644 index 0000000..a51bc14 --- /dev/null +++ b/macros/matrix_linprog.bin diff --git a/macros/mps_linprog.bin b/macros/mps_linprog.bin Binary files differnew file mode 100644 index 0000000..c0d2e3e --- /dev/null +++ b/macros/mps_linprog.bin diff --git a/macros/names b/macros/names new file mode 100644 index 0000000..17cb430 --- /dev/null +++ b/macros/names @@ -0,0 +1,31 @@ +Checkdims +Checklhs +Checkrhs +Checktype +Checkvector +cbcintlinprog +cbcmatrixintlinprog +cbcmpsintlinprog +ecos +fgoalattain +fminbnd +fmincon +fminimax +fminunc +intfminbnd +intfmincon +intfminimax +intfminunc +intqpipopt +linprog +lsqlin +lsqnonlin +lsqnonneg +matrix_linprog +mps_linprog +qpipopt +qpipoptmat +setOptions +symphony +symphony_call +symphonymat diff --git a/macros/qpipopt.bin b/macros/qpipopt.bin Binary files differnew file mode 100644 index 0000000..e2ba3de --- /dev/null +++ b/macros/qpipopt.bin diff --git a/macros/qpipoptmat.bin b/macros/qpipoptmat.bin Binary files differnew file mode 100644 index 0000000..b9c741f --- /dev/null +++ b/macros/qpipoptmat.bin diff --git a/macros/setOptions.bin b/macros/setOptions.bin Binary files differnew file mode 100644 index 0000000..8d23e73 --- /dev/null +++ b/macros/setOptions.bin diff --git a/macros/symphony.bin b/macros/symphony.bin Binary files differnew file mode 100644 index 0000000..3dacf06 --- /dev/null +++ b/macros/symphony.bin diff --git a/macros/symphony_call.bin b/macros/symphony_call.bin Binary files differnew file mode 100644 index 0000000..a883f05 --- /dev/null +++ b/macros/symphony_call.bin diff --git a/macros/symphonymat.bin b/macros/symphonymat.bin Binary files differnew file mode 100644 index 0000000..9d7271a --- /dev/null +++ b/macros/symphonymat.bin |