Mercurial > octave-nkf
changeset 5237:652e8aa49fa7
[project @ 2005-03-23 21:28:45 by jwe]
| author | jwe |
|---|---|
| date | Wed, 23 Mar 2005 21:28:46 +0000 |
| parents | 6879f10db3a4 |
| children | d432b7809fe5 |
| files | scripts/ChangeLog scripts/configure.in scripts/optimization/glpk.m scripts/optimization/glpkmex.m scripts/optimization/glpkparams.m scripts/optimization/glpktest1 scripts/optimization/glpktest2 src/ChangeLog src/DLD-FUNCTIONS/__glpk__.cc src/toplev.cc |
| diffstat | 10 files changed, 276 insertions(+), 142 deletions(-) [+] |
line wrap: on
line diff
--- a/scripts/ChangeLog Wed Mar 23 17:11:08 2005 +0000 +++ b/scripts/ChangeLog Wed Mar 23 21:28:46 2005 +0000 @@ -1,5 +1,14 @@ +2005-03-23 John W. Eaton <jwe@octave.org> + + * optimization/glpk.m: Simplify interface. By default, solve + standard LP min C'*x s.t. A*x = b, x >= 0. + * optimization/glpkmex.m: New file. + 2005-03-22 John W. Eaton <jwe@octave.org> + * configure.in (AC_CONFIG_FILES): Add optimization/Makefile to the + list. + * optimization/glpk.m: Adapt to Octave coding style. No need for varargout or varargin. Print usage message if nargin > 11.
--- a/scripts/configure.in Wed Mar 23 17:11:08 2005 +0000 +++ b/scripts/configure.in Wed Mar 23 21:28:46 2005 +0000 @@ -32,7 +32,8 @@ control/obsolete/Makefile control/system/Makefile \ control/util/Makefile deprecated/Makefile elfun/Makefile \ finance/Makefile general/Makefile image/Makefile io/Makefile \ - linear-algebra/Makefile miscellaneous/Makefile plot/Makefile \ + linear-algebra/Makefile miscellaneous/Makefile \ + optimization/Makefile plot/Makefile \ polynomial/Makefile quaternion/Makefile set/Makefile \ signal/Makefile sparse/Makefile specfun/Makefile \ special-matrix/Makefile startup/Makefile statistics/Makefile \
--- a/scripts/optimization/glpk.m Wed Mar 23 17:11:08 2005 +0000 +++ b/scripts/optimization/glpk.m Wed Mar 23 21:28:46 2005 +0000 @@ -19,19 +19,28 @@ ## GLPK - An Octave Interface for the GNU GLPK library ## -## This routine calls the glpk library to solve an LP/MIP problem. A typical -## LP problem has following structure: +## This routine calls the glpk library to solve an LP/MIP problem. By +## default, glpk solves the following standard LP: +## +## min C'*x +## +## subject to +## +## A*x = b +## x >= 0 ## -## [min|max] C'x -## s.t. -## Ax ["="|"<="|">="] b -## {x <= UB} -## {x >= LB} +## but may also solve problems of the form +## +## [min|max] C'x +## +## subject to +## +## Ax ["="|"<="|">="] b +## {x <= UB} +## {x >= LB} ## ## The calling syntax is: -## [XOPT,FOPT,STATUS,EXTRA]=glpk(SENSE,C,... -## A,B,CTYPE,LB,UB,... -## VARTYPE,PARAM,LPSOLVER,SAVE) +## [XOPT,FOPT,STATUS,EXTRA]=glpk(C,A,B,LB,UB,CTYPE,VARTYPE,SENSE,PARAM) ## ## For a quick reference to the syntax just type glpk at command prompt. ## @@ -40,11 +49,6 @@ ## ## --- INPUTS --- ## -## SENSE: indicates whether the problem is a minimization -## or maximization problem. -## SENSE = 1 minimize -## SENSE = -1 maximize. -## ## C: A column array containing the objective function ## coefficients. ## @@ -53,6 +57,12 @@ ## B: A column array containing the right-hand side value for ## each constraint in the constraint matrix. ## +## LB: An array of at least length numcols containing the lower +## bound on each of the variables. +## +## UB: An array of at least length numcols containing the upper +## bound on each of the variables. +## ## CTYPE: A column array containing the sense of each constraint ## in the constraint matrix. ## CTYPE(i) = 'F' Free (unbounded) variable @@ -62,30 +72,20 @@ ## CTYPE(i) = 'D' Double-bounded variable ## (This is case sensitive). ## -## LB: An array of at least length numcols containing the lower -## bound on each of the variables. -## -## UB: An array of at least length numcols containing the upper -## bound on each of the variables. -## ## VARTYPE: A column array containing the types of the variables. ## VARTYPE(i) = 'C' continuous variable ## VARTYPE(i) = 'I' Integer variable ## (This is case sensitive). ## +## SENSE: indicates whether the problem is a minimization +## or maximization problem. +## SENSE = 1 minimize +## SENSE = -1 maximize. +## ## PARAM: A structure containing some parameters used to define ## the behavior of solver. For more details type ## HELP GLPKPARAMS. ## -## LPSOLVER: Selects which solver using to solve LP problems. -## LPSOLVER=1 Revised Simplex Method -## LPSOLVER=2 Interior Point Method -## If the problem is a MIP problem this flag will be ignored. -## -## SAVE: Saves a copy of the problem if SAVE<>0. -## The file name can not be specified and defaults to "outpb.lp". -## The output file is CPLEX LP format. -## ## --- OUTPUTS --- ## ## XOPT: The optimizer. @@ -145,11 +145,11 @@ ## Author: Nicolo' Giorgetti <giorgetti@dii.unisi.it> ## Adapted-by: jwe -function [xopt, fmin, status, extra] = glpk (sense, c, a, b, ctype, lb, ub, vartype, param, lpsolver, savepb) +function [xopt, fmin, status, extra] = glpk (c, a, b, lb, ub, ctype, vartype, sense, param) ## If there is no input output the version and syntax - if (nargin < 4 || nargin > 11) - usage ("[xopt, fopt, status, extra] = glpk (sense, c, a, b, ctype, lb, ub, vartype, param, lpsolver, savepb"); + if (nargin < 3 || nargin > 9) + usage ("[xopt, fopt, status, extra] = glpk (c, a, b, lb, ub, ctype, vartype, sense, param)"); return; endif @@ -161,7 +161,7 @@ ## Force column vector. c = c(:); - ## 3) Matrix constraint + ## 2) Matrix constraint if (isempty (a)) error ("A cannot be an empty matrix"); @@ -173,7 +173,7 @@ return; endif - ## 4) RHS + ## 3) RHS if (isempty (b)) error ("B cannot be an empty vector"); @@ -184,39 +184,22 @@ return; endif - ## 5) Sense of each constraint + ## 4) Vector with the lower bound of each variable - if (nargin > 4) - if (isempty (ctype)) - ctype = repmat ("U", nc, 1); - elseif (! ischar (ctype) || all (size (ctype) > 1) || length (ctype) != nc) - error ("CTYPE must be a char valued %d by 1 column vector", nc); - return; - elseif (! all (ctype == "F" | ctype == "U" | ctype == "S" - | ctype == "L" | ctype == "D")) - error ("CTYPE must contain only F, U, S, L, or D"); - return; - endif - else - ctype = repmat ("U", nc, 1); - end - - ## 6) Vector with the lower bound of each variable - - if (nargin > 5) + if (nargin > 3) if (isempty (lb)) - lb = repmat (-Inf, nx, 1); + lb = zeros (0, nx, 1); elseif (! isreal (lb) || all (size (lb) > 1) || length (lb) != nx) error ("LB must be a real valued %d by 1 column vector", nx); return; endif else - lb = repmat (-Inf, nx, 1); + lb = zeros (nx, 1); end - ## 7) Vector with the upper bound of each variable + ## 5) Vector with the upper bound of each variable - if (nargin > 6) + if (nargin > 4) if (isempty (ub)) ub = repmat (Inf, nx, 1); elseif (! isreal (ub) || all (size (ub) > 1) || length (ub) != nx) @@ -227,9 +210,26 @@ ub = repmat (Inf, nx, 1); end - ## 8) Vector with the type of variables + ## 6) Sense of each constraint - if (nargin > 7) + if (nargin > 5) + if (isempty (ctype)) + ctype = repmat ("S", nc, 1); + elseif (! ischar (ctype) || all (size (ctype) > 1) || length (ctype) != nc) + error ("CTYPE must be a char valued %d by 1 column vector", nc); + return; + elseif (! all (ctype == "F" | ctype == "U" | ctype == "S" + | ctype == "L" | ctype == "D")) + error ("CTYPE must contain only F, U, S, L, or D"); + return; + endif + else + ctype = repmat ("S", nc, 1); + end + + ## 7) Vector with the type of variables + + if (nargin > 6) if isempty (vartype) vartype = repmat ("C", nx, 1); elseif (! ischar (vartype) || all (size (vartype) > 1) @@ -245,7 +245,7 @@ vartype = repmat ("C", nx, 1); endif - ## 9) Parameters vector + ## 8) Parameters vector if (nargin > 8) if (! isstruct (param)) @@ -256,29 +256,7 @@ param = struct (); endif - ## 10) Select solver method: simplex or interior point - - if (nargin > 9) - if (! (isreal (lpsolver) && isscalar (lpsolver))) - error ("LPSOLVER must be a real scalar value"); - return; - endif - else - lpsolver = 1; - endif - - ## 11) Save the problem - - if (nargin > 10) - if (! (isreal (savepb) && isscalar (lpsolver))) - error ("LPSOLVER must be a real scalar value"); - return; - endif - else - savepb = 0; - end - [xopt, fmin, status, extra] = ... - __glpk__ (sense, c, a, b, ctype, lb, ub, vartype, param, lpsolver, savepb); + __glpk__ (c, a, b, lb, ub, ctype, vartype, sense, param); endfunction
--- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/scripts/optimization/glpkmex.m Wed Mar 23 21:28:46 2005 +0000 @@ -0,0 +1,108 @@ +## Copyright (C) 2005 Nicolo' Giorgetti +## +## This file is part of Octave. +## +## Octave is free software; you can redistribute it and/or modify it +## under the terms of the GNU General Public License as published by +## the Free Software Foundation; either version 2, or (at your option) +## any later version. +## +## Octave is distributed in the hope that it will be useful, but +## WITHOUT ANY WARRANTY; without even the implied warranty of +## MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU +## General Public License for more details. +## +## You should have received a copy of the GNU General Public License +## along with Octave; see the file COPYING. If not, write to the Free +## Software Foundation, 59 Temple Place - Suite 330, Boston, MA +## 02111-1307, USA. + +## GLPKMEX - An Octave Interface for the GNU GLPK library +## +## This function is provided for compatibility with the old Matlab +## interface to GLPK. + +function [xopt, fopt, status, extra] = glpkmex (varargin) + + ## If there is no input output the version and syntax + if (nargin < 4 || nargin > 11) + usage ("[xopt, fopt, status, extra] = glpkmex (sense, c, a, b, ctype, lb, ub, vartype, param, lpsolver, savepb"); + return; + endif + + ## reorder args: + ## + ## glpkmex glpk + ## + ## 1 sense c + ## 2 c a + ## 3 a b + ## 4 b lb + ## 5 ctype ub + ## 6 lb ctype + ## 7 ub vartype + ## 8 vartype sense + ## 9 param param + ## 10 lpsolver + ## 11 savepb + + sense = varargin{1}; + c = varargin{2}; + a = varargin{3}; + b = varargin{4}; + + nx = length (c); + + if (nargin > 4) + ctype = varargin{5}; + else + ctype = repmat ("U", nx, 1); + endif + + if (nargin > 5) + lb = varargin{6}; + else + lb = repmat (-Inf, nx, 1); + endif + + if (nargin > 6) + ub = varargin{7}; + else + ub = repmat (Inf, nx, 1); + endif + + if (nargin > 7) + vartype = varargin{8}; + else + vartype = repmat ("C", nx, 1); + endif + + if (nargin > 8) + param = varargin{9}; + else + param = struct (); + endif + + if (nargin > 9 && ! isfield (param, "lpsolver")) + param.lpsolver = varargin{10}; + endif + + if (nargin > 10 && ! isfield (param, "save")) + param.lpsolver = varargin{11}; + endif + + if (nargout == 0) + glpk (c, a, b, lb, ub, ctype, vartype, sense, param); + elseif (nargout == 1) + xopt = glpk (c, a, b, lb, ub, ctype, vartype, sense, param); + elseif (nargout == 2) + [xopt, fopt] = glpk (c, a, b, lb, ub, ctype, vartype, sense, param); + elseif (nargout == 3) + [xopt, fopt, status] = ... + glpk (c, a, b, lb, ub, ctype, vartype, sense, param); + else + [xopt, fopt, status, extra] = ... + glpk (c, a, b, lb, ub, ctype, vartype, sense, param); + endif + +endfunction
--- a/scripts/optimization/glpkparams.m Wed Mar 23 17:11:08 2005 +0000 +++ b/scripts/optimization/glpkparams.m Wed Mar 23 21:28:46 2005 +0000 @@ -70,6 +70,16 @@ % problem using the built-in LP presolver. Otherwise the LP % presolver is not used. % +% lpsolver type: int, default: 1 +% Select which solver to use: +% 1 revised simplex method +% 2 interior point method +% If the problem is a MIP problem this flag will be ignored. +% +% save type: int, default 0 +% If this parameter is nonzero, save a copy of the problem +% problem in CPLEX LP format to the file "outpb.lp". There +% is currently no way to change the name of the output file. % % ----------------------- % 2 Real parameters
--- a/scripts/optimization/glpktest1 Wed Mar 23 17:11:08 2005 +0000 +++ b/scripts/optimization/glpktest1 Wed Mar 23 21:28:46 2005 +0000 @@ -1,17 +1,20 @@ -clear; - -disp('1st problem'); -s=-1; -c=[10,6,4]'; -a=[1,1,1;... - 10,4,5;... - 2,2,6]; -b=[100,600,300]'; -ctype=['U','U','U']'; -lb=[0,0,0]'; -ub=[]'; -vartype=['C','C','C']'; -param.msglev=3; -lpsolver=1; -save_pb=1; -[xmin,fmin,status,extra]=glpk(s,c,a,b,ctype,lb,ub,vartype,param,lpsolver,save_pb) +clear; + +disp('1st problem'); +s=-1; +c=[10,6,4]'; +a=[1,1,1;... + 10,4,5;... + 2,2,6]; +b=[100,600,300]'; +ctype=['U','U','U']'; +lb=[0,0,0]'; +ub=[]'; +vartype=['C','C','C']'; +param.msglev=3; +param.lpsolver=1; +param.save=1; +[xmin,fmin,status,extra]=glpk(c,a,b,lb,ub,ctype,vartype,s,param) +lpsolver = param.lpsolver; +save_pb = param.save; +[xmin,fmin,status,extra]=glpkmex(s,c,a,b,ctype,lb,ub,vartype,param,lpsolver,save_pb)
--- a/scripts/optimization/glpktest2 Wed Mar 23 17:11:08 2005 +0000 +++ b/scripts/optimization/glpktest2 Wed Mar 23 21:28:46 2005 +0000 @@ -9,7 +9,8 @@ lb=[0;0]; ub=[]; vartype=['I';'I']; param.msglev=1; -[xmin,fmin,status,extra]=glpk(s,c,a,b,ctype,lb,ub,vartype,param) +[xmin,fmin,status,extra]=glpk(c,a,b,lb,ub,ctype,vartype,s,param) +[xmin,fmin,status,extra]=glpkmex(s,c,a,b,ctype,lb,ub,vartype,param) pause; disp('3rd problem'); @@ -22,4 +23,5 @@ ctype=['S','S','S']'; lb=[0,0,0,0,0]'; ub=[]; vartype=['C','C','C','C','C']'; -[xmin,fmin,status,extra]=glpk(s,c,a,b,ctype,lb,ub,vartype) +[xmin,fmin,status,extra]=glpk(c,a,b,lb,ub,ctype,vartype,s) +[xmin,fmin,status,extra]=glpkmex(s,c,a,b,ctype,lb,ub,vartype)
--- a/src/ChangeLog Wed Mar 23 17:11:08 2005 +0000 +++ b/src/ChangeLog Wed Mar 23 21:28:46 2005 +0000 @@ -1,3 +1,12 @@ +2005-03-23 John W. Eaton <jwe@octave.org> + + * DLD-FUNCTIONS/__glpk__.cc (F__glpk__): Accept lpsolver and + save_pb in param arg instead of as separate args. + Arg list now matches new interface for glpk.m. + + * toplev.cc (do_octave_atexit): Call reset_error_handler before + each call to feval. + 2005-03-22 John W. Eaton <jwe@octave.org> * Makefile.in: Add special rule for __glpk__.oct.
--- a/src/DLD-FUNCTIONS/__glpk__.cc Wed Mar 23 17:11:08 2005 +0000 +++ b/src/DLD-FUNCTIONS/__glpk__.cc Wed Mar 23 21:28:46 2005 +0000 @@ -378,26 +378,18 @@ return retval; } - //-- 1st Input. Sense of optimization. - volatile int sense; - double SENSE = args(0).scalar_value (); - if (SENSE >= 0) - sense = 1; - else - sense = -1; + //-- 1nd Input. A column array containing the objective function + //-- coefficients. + int mrowsc = args(0).rows(); - //-- 2nd Input. A column array containing the objective function - //-- coefficients. - int mrowsc = args(1).rows(); - - Matrix C (args(1).matrix_value ()); + Matrix C (args(0).matrix_value ()); double *c = C.fortran_vec (); - //-- 3rd Input. A matrix containing the constraints coefficients. + //-- 2nd Input. A matrix containing the constraints coefficients. // If matrix A is NOT a sparse matrix // if(!mxIsSparse(A_IN)){ - int mrowsA = args(2).rows(); - Matrix A (args(2).matrix_value ()); // get the matrix + int mrowsA = args(1).rows(); + Matrix A (args(1).matrix_value ()); // get the matrix Array<int> rn (mrowsA*mrowsc+1); Array<int> cn (mrowsA*mrowsc+1); ColumnVector a (mrowsA*mrowsc+1, 0.0); @@ -449,19 +441,14 @@ // } // } - //-- 4th Input. A column array containing the right-hand side value + //-- 3rd Input. A column array containing the right-hand side value // for each constraint in the constraint matrix. - Matrix B (args(3).matrix_value ()); + Matrix B (args(2).matrix_value ()); double *b = B.fortran_vec (); - //-- 5th Input. A column array containing the sense of each constraint - //-- in the constraint matrix. - charMatrix CTYPE (args(4).char_matrix_value ()); - char *ctype = CTYPE.fortran_vec (); - - //-- 6th Input. An array of length mrowsc containing the lower + //-- 4th Input. An array of length mrowsc containing the lower //-- bound on each of the variables. - Matrix LB (args(5).matrix_value ()); + Matrix LB (args(3).matrix_value ()); double *lb = LB.fortran_vec (); //-- LB argument, default: Free @@ -477,9 +464,9 @@ freeLB(i) = 0; } - //-- 7th Input. An array of at least length numcols containing the upper + //-- 5th Input. An array of at least length numcols containing the upper //-- bound on each of the variables. - Matrix UB (args(6).matrix_value ()); + Matrix UB (args(4).matrix_value ()); double *ub = UB.fortran_vec (); @@ -495,8 +482,13 @@ freeUB(i) = 0; } - //-- 8th Input. A column array containing the types of the variables. - charMatrix VTYPE (args(7).char_matrix_value ()); + //-- 6th Input. A column array containing the sense of each constraint + //-- in the constraint matrix. + charMatrix CTYPE (args(5).char_matrix_value ()); + char *ctype = CTYPE.fortran_vec (); + + //-- 7th Input. A column array containing the types of the variables. + charMatrix VTYPE (args(6).char_matrix_value ()); Array<int> vartype (mrowsc); volatile int isMIP = 0; @@ -511,6 +503,14 @@ vartype(i) = LPX_CV; } + //-- 8th Input. Sense of optimization. + volatile int sense; + double SENSE = args(7).scalar_value (); + if (SENSE >= 0) + sense = 1; + else + sense = -1; + //-- 9th Input. A structure containing the control parameters. Octave_map PARAM = args(8).map_value (); @@ -645,6 +645,28 @@ lpxIntParam[16] = static_cast<int> (numtmp); } + //-- LPsolver option + volatile int lpsolver = 1; + if (PARAM.contains ("lpsolver")) + { + octave_value tmp = PARAM.contents (PARAM.seek ("lpsolver"))(0); + double numtmp = tmp.scalar_value (); + if (numtmp != 1 && numtmp != 2) + { + OCTOUT << "'lpsolver' parameter must be only:\n\t1 - simplex method,\n\t2 - interior point method\n"; + return retval; + } + lpsolver = static_cast<int> (numtmp); + } + + //-- Save option + volatile int save_pb = 0; + if (PARAM.contains ("save")) + { + octave_value tmp = PARAM.contents (PARAM.seek ("save"))(0); + save_pb = (tmp.scalar_value () != 0); + } + //-- ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ //-- Real parameters //-- ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ @@ -716,16 +738,6 @@ lpxRealParam[9] = tmp.scalar_value (); } - //-- 10th Input. If the problem is a LP problem you may select which solver - //-- use: RSM (Revised Simplex Method) or IPM (Interior Point Method). - //-- If the problem is a MIP problem this field will be ignored. - octave_value tmp = args(9).scalar_value (); - int lpsolver = static_cast<int> (tmp.scalar_value ()); - - //-- 11th Input. Saves a copy of the problem if SAVE<>0. - tmp = args(10).scalar_value(); - int save_pb = (tmp.scalar_value() != 0); - //-- Assign pointers to the output parameters ColumnVector xmin (mrowsc); ColumnVector fmin (1);