octave-nkf: scripts/optimization/qp.m annotate

annotate scripts/optimization/qp.m @ 19898:4197fc428c7d

maint: Update copyright notices for 2015.

author	John W. Eaton <jwe@octave.org>
date	Wed, 11 Feb 2015 14:19:08 -0500
parents	a404853d2073
children	9fc020886ae9

rev	line source
19898 4197fc428c7d maint: Update copyright notices for 2015. John W. Eaton <jwe@octave.org> parents: 17931 diff changeset	1 ## Copyright (C) 2013-2015 Julien Bect
4197fc428c7d maint: Update copyright notices for 2015. John W. Eaton <jwe@octave.org> parents: 17931 diff changeset	2 ## Copyright (C) 2000-2015 Gabriele Pannocchia.
5289 b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	3 ##
b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	4 ## This file is part of Octave.
b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	5 ##
b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	6 ## Octave is free software; you can redistribute it and/or modify it
b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	7 ## under the terms of the GNU General Public License as published by
7016 93c65f2a5668 [project @ 2007-10-12 06:40:56 by jwe] jwe parents: 6848 diff changeset	8 ## the Free Software Foundation; either version 3 of the License, or (at
93c65f2a5668 [project @ 2007-10-12 06:40:56 by jwe] jwe parents: 6848 diff changeset	9 ## your option) any later version.
5289 b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	10 ##
b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	11 ## Octave is distributed in the hope that it will be useful, but
b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	12 ## WITHOUT ANY WARRANTY; without even the implied warranty of
b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	13 ## MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	14 ## General Public License for more details.
b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	15 ##
b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	16 ## You should have received a copy of the GNU General Public License
7016 93c65f2a5668 [project @ 2007-10-12 06:40:56 by jwe] jwe parents: 6848 diff changeset	17 ## along with Octave; see the file COPYING. If not, see
93c65f2a5668 [project @ 2007-10-12 06:40:56 by jwe] jwe parents: 6848 diff changeset	18 ## <http://www.gnu.org/licenses/>.
5289 b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	19
b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	20 ## -- texinfo --
10793 be55736a0783 Grammarcheck the documentation from m-files. Rik <octave@nomad.inbox5.com> parents: 10549 diff changeset	21 ## @deftypefn {Function File} {[@var{x}, @var{obj}, @var{info}, @var{lambda}] =} qp (@var{x0}, @var{H})
10059 665ad34efeed qp.m: handle optimset options Joshua Redstone <redstone@gmail.com>, John W. Eaton <jwe@octave.org> parents: 9872 diff changeset	22 ## @deftypefnx {Function File} {[@var{x}, @var{obj}, @var{info}, @var{lambda}] =} qp (@var{x0}, @var{H}, @var{q})
665ad34efeed qp.m: handle optimset options Joshua Redstone <redstone@gmail.com>, John W. Eaton <jwe@octave.org> parents: 9872 diff changeset	23 ## @deftypefnx {Function File} {[@var{x}, @var{obj}, @var{info}, @var{lambda}] =} qp (@var{x0}, @var{H}, @var{q}, @var{A}, @var{b})
665ad34efeed qp.m: handle optimset options Joshua Redstone <redstone@gmail.com>, John W. Eaton <jwe@octave.org> parents: 9872 diff changeset	24 ## @deftypefnx {Function File} {[@var{x}, @var{obj}, @var{info}, @var{lambda}] =} qp (@var{x0}, @var{H}, @var{q}, @var{A}, @var{b}, @var{lb}, @var{ub})
665ad34efeed qp.m: handle optimset options Joshua Redstone <redstone@gmail.com>, John W. Eaton <jwe@octave.org> parents: 9872 diff changeset	25 ## @deftypefnx {Function File} {[@var{x}, @var{obj}, @var{info}, @var{lambda}] =} qp (@var{x0}, @var{H}, @var{q}, @var{A}, @var{b}, @var{lb}, @var{ub}, @var{A_lb}, @var{A_in}, @var{A_ub})
665ad34efeed qp.m: handle optimset options Joshua Redstone <redstone@gmail.com>, John W. Eaton <jwe@octave.org> parents: 9872 diff changeset	26 ## @deftypefnx {Function File} {[@var{x}, @var{obj}, @var{info}, @var{lambda}] =} qp (@dots{}, @var{options})
5289 b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	27 ## Solve the quadratic program
6741 00116015904d [project @ 2007-06-18 16:07:14 by jwe] jwe parents: 6526 diff changeset	28 ## @tex
00116015904d [project @ 2007-06-18 16:07:14 by jwe] jwe parents: 6526 diff changeset	29 ## $$
00116015904d [project @ 2007-06-18 16:07:14 by jwe] jwe parents: 6526 diff changeset	30 ## \min_x {1 \over 2} x^T H x + x^T q
00116015904d [project @ 2007-06-18 16:07:14 by jwe] jwe parents: 6526 diff changeset	31 ## $$
00116015904d [project @ 2007-06-18 16:07:14 by jwe] jwe parents: 6526 diff changeset	32 ## @end tex
00116015904d [project @ 2007-06-18 16:07:14 by jwe] jwe parents: 6526 diff changeset	33 ## @ifnottex
5289 b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	34 ##
b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	35 ## @example
9051 1bf0ce0930be Grammar check TexInfo in all .m files Rik <rdrider0-list@yahoo.com> parents: 8920 diff changeset	36 ## @group
14327 4d917a6a858b doc: Use Octave coding conventions in @example blocks of docstrings. Rik <octave@nomad.inbox5.com> parents: 14138 diff changeset	37 ## min 0.5 x'Hx + x'*q
4d917a6a858b doc: Use Octave coding conventions in @example blocks of docstrings. Rik <octave@nomad.inbox5.com> parents: 14138 diff changeset	38 ## x
9051 1bf0ce0930be Grammar check TexInfo in all .m files Rik <rdrider0-list@yahoo.com> parents: 8920 diff changeset	39 ## @end group
5289 b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	40 ## @end example
b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	41 ##
6741 00116015904d [project @ 2007-06-18 16:07:14 by jwe] jwe parents: 6526 diff changeset	42 ## @end ifnottex
00116015904d [project @ 2007-06-18 16:07:14 by jwe] jwe parents: 6526 diff changeset	43 ## subject to
5289 b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	44 ## @tex
6741 00116015904d [project @ 2007-06-18 16:07:14 by jwe] jwe parents: 6526 diff changeset	45 ## $$
00116015904d [project @ 2007-06-18 16:07:14 by jwe] jwe parents: 6526 diff changeset	46 ## Ax = b \qquad lb \leq x \leq ub \qquad A_{lb} \leq A_{in} \leq A_{ub}
00116015904d [project @ 2007-06-18 16:07:14 by jwe] jwe parents: 6526 diff changeset	47 ## $$
5289 b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	48 ## @end tex
6741 00116015904d [project @ 2007-06-18 16:07:14 by jwe] jwe parents: 6526 diff changeset	49 ## @ifnottex
5289 b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	50 ##
b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	51 ## @example
9051 1bf0ce0930be Grammar check TexInfo in all .m files Rik <rdrider0-list@yahoo.com> parents: 8920 diff changeset	52 ## @group
14327 4d917a6a858b doc: Use Octave coding conventions in @example blocks of docstrings. Rik <octave@nomad.inbox5.com> parents: 14138 diff changeset	53 ## A*x = b
4d917a6a858b doc: Use Octave coding conventions in @example blocks of docstrings. Rik <octave@nomad.inbox5.com> parents: 14138 diff changeset	54 ## lb <= x <= ub
4d917a6a858b doc: Use Octave coding conventions in @example blocks of docstrings. Rik <octave@nomad.inbox5.com> parents: 14138 diff changeset	55 ## A_lb <= A_in*x <= A_ub
9051 1bf0ce0930be Grammar check TexInfo in all .m files Rik <rdrider0-list@yahoo.com> parents: 8920 diff changeset	56 ## @end group
5289 b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	57 ## @end example
10821 693e22af08ae Grammarcheck documentation of m-files Rik <octave@nomad.inbox5.com> parents: 10793 diff changeset	58 ##
6741 00116015904d [project @ 2007-06-18 16:07:14 by jwe] jwe parents: 6526 diff changeset	59 ## @end ifnottex
5289 b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	60 ## @noindent
b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	61 ## using a null-space active-set method.
b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	62 ##
b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	63 ## Any bound (@var{A}, @var{b}, @var{lb}, @var{ub}, @var{A_lb},
b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	64 ## @var{A_ub}) may be set to the empty matrix (@code{[]}) if not
b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	65 ## present. If the initial guess is feasible the algorithm is faster.
b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	66 ##
10059 665ad34efeed qp.m: handle optimset options Joshua Redstone <redstone@gmail.com>, John W. Eaton <jwe@octave.org> parents: 9872 diff changeset	67 ## @table @var
665ad34efeed qp.m: handle optimset options Joshua Redstone <redstone@gmail.com>, John W. Eaton <jwe@octave.org> parents: 9872 diff changeset	68 ## @item options
665ad34efeed qp.m: handle optimset options Joshua Redstone <redstone@gmail.com>, John W. Eaton <jwe@octave.org> parents: 9872 diff changeset	69 ## An optional structure containing the following
665ad34efeed qp.m: handle optimset options Joshua Redstone <redstone@gmail.com>, John W. Eaton <jwe@octave.org> parents: 9872 diff changeset	70 ## parameter(s) used to define the behavior of the solver. Missing elements
665ad34efeed qp.m: handle optimset options Joshua Redstone <redstone@gmail.com>, John W. Eaton <jwe@octave.org> parents: 9872 diff changeset	71 ## in the structure take on default values, so you only need to set the
665ad34efeed qp.m: handle optimset options Joshua Redstone <redstone@gmail.com>, John W. Eaton <jwe@octave.org> parents: 9872 diff changeset	72 ## elements that you wish to change from the default.
665ad34efeed qp.m: handle optimset options Joshua Redstone <redstone@gmail.com>, John W. Eaton <jwe@octave.org> parents: 9872 diff changeset	73 ##
665ad34efeed qp.m: handle optimset options Joshua Redstone <redstone@gmail.com>, John W. Eaton <jwe@octave.org> parents: 9872 diff changeset	74 ## @table @code
665ad34efeed qp.m: handle optimset options Joshua Redstone <redstone@gmail.com>, John W. Eaton <jwe@octave.org> parents: 9872 diff changeset	75 ## @item MaxIter (default: 200)
665ad34efeed qp.m: handle optimset options Joshua Redstone <redstone@gmail.com>, John W. Eaton <jwe@octave.org> parents: 9872 diff changeset	76 ## Maximum number of iterations.
665ad34efeed qp.m: handle optimset options Joshua Redstone <redstone@gmail.com>, John W. Eaton <jwe@octave.org> parents: 9872 diff changeset	77 ## @end table
665ad34efeed qp.m: handle optimset options Joshua Redstone <redstone@gmail.com>, John W. Eaton <jwe@octave.org> parents: 9872 diff changeset	78 ## @end table
665ad34efeed qp.m: handle optimset options Joshua Redstone <redstone@gmail.com>, John W. Eaton <jwe@octave.org> parents: 9872 diff changeset	79 ##
10821 693e22af08ae Grammarcheck documentation of m-files Rik <octave@nomad.inbox5.com> parents: 10793 diff changeset	80 ## @table @var
693e22af08ae Grammarcheck documentation of m-files Rik <octave@nomad.inbox5.com> parents: 10793 diff changeset	81 ## @item info
693e22af08ae Grammarcheck documentation of m-files Rik <octave@nomad.inbox5.com> parents: 10793 diff changeset	82 ## Structure containing run-time information about the algorithm. The
693e22af08ae Grammarcheck documentation of m-files Rik <octave@nomad.inbox5.com> parents: 10793 diff changeset	83 ## following fields are defined:
693e22af08ae Grammarcheck documentation of m-files Rik <octave@nomad.inbox5.com> parents: 10793 diff changeset	84 ##
5289 b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	85 ## @table @code
b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	86 ## @item solveiter
b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	87 ## The number of iterations required to find the solution.
10821 693e22af08ae Grammarcheck documentation of m-files Rik <octave@nomad.inbox5.com> parents: 10793 diff changeset	88 ##
5289 b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	89 ## @item info
10821 693e22af08ae Grammarcheck documentation of m-files Rik <octave@nomad.inbox5.com> parents: 10793 diff changeset	90 ## An integer indicating the status of the solution.
11587 c792872f8942 all script files: untabify and strip trailing whitespace John W. Eaton <jwe@octave.org> parents: 11523 diff changeset	91 ##
5289 b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	92 ## @table @asis
b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	93 ## @item 0
b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	94 ## The problem is feasible and convex. Global solution found.
10821 693e22af08ae Grammarcheck documentation of m-files Rik <octave@nomad.inbox5.com> parents: 10793 diff changeset	95 ##
5289 b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	96 ## @item 1
b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	97 ## The problem is not convex. Local solution found.
10821 693e22af08ae Grammarcheck documentation of m-files Rik <octave@nomad.inbox5.com> parents: 10793 diff changeset	98 ##
5289 b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	99 ## @item 2
b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	100 ## The problem is not convex and unbounded.
10821 693e22af08ae Grammarcheck documentation of m-files Rik <octave@nomad.inbox5.com> parents: 10793 diff changeset	101 ##
5289 b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	102 ## @item 3
b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	103 ## Maximum number of iterations reached.
10821 693e22af08ae Grammarcheck documentation of m-files Rik <octave@nomad.inbox5.com> parents: 10793 diff changeset	104 ##
5289 b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	105 ## @item 6
b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	106 ## The problem is infeasible.
b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	107 ## @end table
b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	108 ## @end table
10821 693e22af08ae Grammarcheck documentation of m-files Rik <octave@nomad.inbox5.com> parents: 10793 diff changeset	109 ## @end table
5289 b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	110 ## @end deftypefn
b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	111
13027 b9a89ca0fb75 prevent optimization functions from setting ans in workspace at startup John W. Eaton <jwe@octave.org> parents: 11589 diff changeset	112 ## PKG_ADD: ## Discard result to avoid polluting workspace with ans at startup.
b9a89ca0fb75 prevent optimization functions from setting ans in workspace at startup John W. Eaton <jwe@octave.org> parents: 11589 diff changeset	113 ## PKG_ADD: [~] = __all_opts__ ("qp");
10060 8f51a90eb8d1 implement default opts query and register opts for qp Jaroslav Hajek <highegg@gmail.com> parents: 10059 diff changeset	114
10059 665ad34efeed qp.m: handle optimset options Joshua Redstone <redstone@gmail.com>, John W. Eaton <jwe@octave.org> parents: 9872 diff changeset	115 function [x, obj, INFO, lambda] = qp (x0, H, varargin)
665ad34efeed qp.m: handle optimset options Joshua Redstone <redstone@gmail.com>, John W. Eaton <jwe@octave.org> parents: 9872 diff changeset	116
665ad34efeed qp.m: handle optimset options Joshua Redstone <redstone@gmail.com>, John W. Eaton <jwe@octave.org> parents: 9872 diff changeset	117 nargs = nargin;
665ad34efeed qp.m: handle optimset options Joshua Redstone <redstone@gmail.com>, John W. Eaton <jwe@octave.org> parents: 9872 diff changeset	118
10060 8f51a90eb8d1 implement default opts query and register opts for qp Jaroslav Hajek <highegg@gmail.com> parents: 10059 diff changeset	119 if (nargin == 1 && ischar (x0) && strcmp (x0, 'defaults'))
8f51a90eb8d1 implement default opts query and register opts for qp Jaroslav Hajek <highegg@gmail.com> parents: 10059 diff changeset	120 x = optimset ("MaxIter", 200);
8f51a90eb8d1 implement default opts query and register opts for qp Jaroslav Hajek <highegg@gmail.com> parents: 10059 diff changeset	121 return;
8f51a90eb8d1 implement default opts query and register opts for qp Jaroslav Hajek <highegg@gmail.com> parents: 10059 diff changeset	122 endif
8f51a90eb8d1 implement default opts query and register opts for qp Jaroslav Hajek <highegg@gmail.com> parents: 10059 diff changeset	123
10059 665ad34efeed qp.m: handle optimset options Joshua Redstone <redstone@gmail.com>, John W. Eaton <jwe@octave.org> parents: 9872 diff changeset	124 if (nargs > 2 && isstruct (varargin{end}))
665ad34efeed qp.m: handle optimset options Joshua Redstone <redstone@gmail.com>, John W. Eaton <jwe@octave.org> parents: 9872 diff changeset	125 options = varargin{end};
665ad34efeed qp.m: handle optimset options Joshua Redstone <redstone@gmail.com>, John W. Eaton <jwe@octave.org> parents: 9872 diff changeset	126 nargs--;
665ad34efeed qp.m: handle optimset options Joshua Redstone <redstone@gmail.com>, John W. Eaton <jwe@octave.org> parents: 9872 diff changeset	127 else
665ad34efeed qp.m: handle optimset options Joshua Redstone <redstone@gmail.com>, John W. Eaton <jwe@octave.org> parents: 9872 diff changeset	128 options = struct ();
665ad34efeed qp.m: handle optimset options Joshua Redstone <redstone@gmail.com>, John W. Eaton <jwe@octave.org> parents: 9872 diff changeset	129 endif
665ad34efeed qp.m: handle optimset options Joshua Redstone <redstone@gmail.com>, John W. Eaton <jwe@octave.org> parents: 9872 diff changeset	130
665ad34efeed qp.m: handle optimset options Joshua Redstone <redstone@gmail.com>, John W. Eaton <jwe@octave.org> parents: 9872 diff changeset	131 if (nargs >= 3)
665ad34efeed qp.m: handle optimset options Joshua Redstone <redstone@gmail.com>, John W. Eaton <jwe@octave.org> parents: 9872 diff changeset	132 q = varargin{1};
665ad34efeed qp.m: handle optimset options Joshua Redstone <redstone@gmail.com>, John W. Eaton <jwe@octave.org> parents: 9872 diff changeset	133 else
665ad34efeed qp.m: handle optimset options Joshua Redstone <redstone@gmail.com>, John W. Eaton <jwe@octave.org> parents: 9872 diff changeset	134 q = [];
665ad34efeed qp.m: handle optimset options Joshua Redstone <redstone@gmail.com>, John W. Eaton <jwe@octave.org> parents: 9872 diff changeset	135 endif
5289 b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	136
10059 665ad34efeed qp.m: handle optimset options Joshua Redstone <redstone@gmail.com>, John W. Eaton <jwe@octave.org> parents: 9872 diff changeset	137 if (nargs >= 5)
665ad34efeed qp.m: handle optimset options Joshua Redstone <redstone@gmail.com>, John W. Eaton <jwe@octave.org> parents: 9872 diff changeset	138 A = varargin{2};
665ad34efeed qp.m: handle optimset options Joshua Redstone <redstone@gmail.com>, John W. Eaton <jwe@octave.org> parents: 9872 diff changeset	139 b = varargin{3};
665ad34efeed qp.m: handle optimset options Joshua Redstone <redstone@gmail.com>, John W. Eaton <jwe@octave.org> parents: 9872 diff changeset	140 else
665ad34efeed qp.m: handle optimset options Joshua Redstone <redstone@gmail.com>, John W. Eaton <jwe@octave.org> parents: 9872 diff changeset	141 A = [];
665ad34efeed qp.m: handle optimset options Joshua Redstone <redstone@gmail.com>, John W. Eaton <jwe@octave.org> parents: 9872 diff changeset	142 b = [];
665ad34efeed qp.m: handle optimset options Joshua Redstone <redstone@gmail.com>, John W. Eaton <jwe@octave.org> parents: 9872 diff changeset	143 endif
665ad34efeed qp.m: handle optimset options Joshua Redstone <redstone@gmail.com>, John W. Eaton <jwe@octave.org> parents: 9872 diff changeset	144
665ad34efeed qp.m: handle optimset options Joshua Redstone <redstone@gmail.com>, John W. Eaton <jwe@octave.org> parents: 9872 diff changeset	145 if (nargs >= 7)
665ad34efeed qp.m: handle optimset options Joshua Redstone <redstone@gmail.com>, John W. Eaton <jwe@octave.org> parents: 9872 diff changeset	146 lb = varargin{4};
665ad34efeed qp.m: handle optimset options Joshua Redstone <redstone@gmail.com>, John W. Eaton <jwe@octave.org> parents: 9872 diff changeset	147 ub = varargin{5};
665ad34efeed qp.m: handle optimset options Joshua Redstone <redstone@gmail.com>, John W. Eaton <jwe@octave.org> parents: 9872 diff changeset	148 else
665ad34efeed qp.m: handle optimset options Joshua Redstone <redstone@gmail.com>, John W. Eaton <jwe@octave.org> parents: 9872 diff changeset	149 lb = [];
665ad34efeed qp.m: handle optimset options Joshua Redstone <redstone@gmail.com>, John W. Eaton <jwe@octave.org> parents: 9872 diff changeset	150 ub = [];
665ad34efeed qp.m: handle optimset options Joshua Redstone <redstone@gmail.com>, John W. Eaton <jwe@octave.org> parents: 9872 diff changeset	151 endif
665ad34efeed qp.m: handle optimset options Joshua Redstone <redstone@gmail.com>, John W. Eaton <jwe@octave.org> parents: 9872 diff changeset	152
665ad34efeed qp.m: handle optimset options Joshua Redstone <redstone@gmail.com>, John W. Eaton <jwe@octave.org> parents: 9872 diff changeset	153 if (nargs == 10)
665ad34efeed qp.m: handle optimset options Joshua Redstone <redstone@gmail.com>, John W. Eaton <jwe@octave.org> parents: 9872 diff changeset	154 A_lb = varargin{6};
665ad34efeed qp.m: handle optimset options Joshua Redstone <redstone@gmail.com>, John W. Eaton <jwe@octave.org> parents: 9872 diff changeset	155 A_in = varargin{7};
665ad34efeed qp.m: handle optimset options Joshua Redstone <redstone@gmail.com>, John W. Eaton <jwe@octave.org> parents: 9872 diff changeset	156 A_ub = varargin{8};
665ad34efeed qp.m: handle optimset options Joshua Redstone <redstone@gmail.com>, John W. Eaton <jwe@octave.org> parents: 9872 diff changeset	157 else
665ad34efeed qp.m: handle optimset options Joshua Redstone <redstone@gmail.com>, John W. Eaton <jwe@octave.org> parents: 9872 diff changeset	158 A_lb = [];
665ad34efeed qp.m: handle optimset options Joshua Redstone <redstone@gmail.com>, John W. Eaton <jwe@octave.org> parents: 9872 diff changeset	159 A_in = [];
665ad34efeed qp.m: handle optimset options Joshua Redstone <redstone@gmail.com>, John W. Eaton <jwe@octave.org> parents: 9872 diff changeset	160 A_ub = [];
665ad34efeed qp.m: handle optimset options Joshua Redstone <redstone@gmail.com>, John W. Eaton <jwe@octave.org> parents: 9872 diff changeset	161 endif
665ad34efeed qp.m: handle optimset options Joshua Redstone <redstone@gmail.com>, John W. Eaton <jwe@octave.org> parents: 9872 diff changeset	162
665ad34efeed qp.m: handle optimset options Joshua Redstone <redstone@gmail.com>, John W. Eaton <jwe@octave.org> parents: 9872 diff changeset	163 if (nargs == 2 \|\| nargs == 3 \|\| nargs == 5 \|\| nargs == 7 \|\| nargs == 10)
665ad34efeed qp.m: handle optimset options Joshua Redstone <redstone@gmail.com>, John W. Eaton <jwe@octave.org> parents: 9872 diff changeset	164
10064 17ce2a700a97 qp.m: Add missing semicolon. Ben Abbott <bpabbott@mac.com> parents: 10060 diff changeset	165 maxit = optimget (options, "MaxIter", 200);
5289 b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	166
b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	167 ## Checking the quadratic penalty
9872 72d6e0de76c7 fix qp, condest and krylov Jaroslav Hajek <highegg@gmail.com> parents: 9211 diff changeset	168 if (! issquare (H))
5289 b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	169 error ("qp: quadratic penalty matrix not square");
9872 72d6e0de76c7 fix qp, condest and krylov Jaroslav Hajek <highegg@gmail.com> parents: 9211 diff changeset	170 elseif (! ishermitian (H))
72d6e0de76c7 fix qp, condest and krylov Jaroslav Hajek <highegg@gmail.com> parents: 9211 diff changeset	171 ## warning ("qp: quadratic penalty matrix not hermitian");
5289 b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	172 H = (H + H')/2;
b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	173 endif
9872 72d6e0de76c7 fix qp, condest and krylov Jaroslav Hajek <highegg@gmail.com> parents: 9211 diff changeset	174 n = rows (H);
5289 b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	175
b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	176 ## Checking the initial guess (if empty it is resized to the
b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	177 ## right dimension and filled with 0)
b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	178 if (isempty (x0))
b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	179 x0 = zeros (n, 1);
10059 665ad34efeed qp.m: handle optimset options Joshua Redstone <redstone@gmail.com>, John W. Eaton <jwe@octave.org> parents: 9872 diff changeset	180 elseif (numel (x0) != n)
5289 b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	181 error ("qp: the initial guess has incorrect length");
b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	182 endif
b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	183
b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	184 ## Linear penalty.
10059 665ad34efeed qp.m: handle optimset options Joshua Redstone <redstone@gmail.com>, John W. Eaton <jwe@octave.org> parents: 9872 diff changeset	185 if (isempty (q))
665ad34efeed qp.m: handle optimset options Joshua Redstone <redstone@gmail.com>, John W. Eaton <jwe@octave.org> parents: 9872 diff changeset	186 q = zeros (n, 1);
665ad34efeed qp.m: handle optimset options Joshua Redstone <redstone@gmail.com>, John W. Eaton <jwe@octave.org> parents: 9872 diff changeset	187 elseif (numel (q) != n)
5289 b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	188 error ("qp: the linear term has incorrect length");
b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	189 endif
b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	190
b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	191 ## Equality constraint matrices
10059 665ad34efeed qp.m: handle optimset options Joshua Redstone <redstone@gmail.com>, John W. Eaton <jwe@octave.org> parents: 9872 diff changeset	192 if (isempty (A) \|\| isempty (b))
665ad34efeed qp.m: handle optimset options Joshua Redstone <redstone@gmail.com>, John W. Eaton <jwe@octave.org> parents: 9872 diff changeset	193 A = zeros (0, n);
665ad34efeed qp.m: handle optimset options Joshua Redstone <redstone@gmail.com>, John W. Eaton <jwe@octave.org> parents: 9872 diff changeset	194 b = zeros (0, 1);
5289 b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	195 n_eq = 0;
b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	196 else
b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	197 [n_eq, n1] = size (A);
b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	198 if (n1 != n)
10549 95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	199 error ("qp: equality constraint matrix has incorrect column dimension");
5289 b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	200 endif
10059 665ad34efeed qp.m: handle optimset options Joshua Redstone <redstone@gmail.com>, John W. Eaton <jwe@octave.org> parents: 9872 diff changeset	201 if (numel (b) != n_eq)
10549 95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	202 error ("qp: equality constraint matrix and vector have inconsistent dimension");
5289 b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	203 endif
b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	204 endif
b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	205
b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	206 ## Bound constraints
b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	207 Ain = zeros (0, n);
b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	208 bin = zeros (0, 1);
b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	209 n_in = 0;
10059 665ad34efeed qp.m: handle optimset options Joshua Redstone <redstone@gmail.com>, John W. Eaton <jwe@octave.org> parents: 9872 diff changeset	210 if (nargs > 5)
5289 b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	211 if (! isempty (lb))
10549 95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	212 if (numel (lb) != n)
95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	213 error ("qp: lower bound has incorrect length");
95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	214 elseif (isempty (ub))
95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	215 Ain = [Ain; eye(n)];
95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	216 bin = [bin; lb];
95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	217 endif
5289 b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	218 endif
b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	219
b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	220 if (! isempty (ub))
10549 95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	221 if (numel (ub) != n)
95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	222 error ("qp: upper bound has incorrect length");
95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	223 elseif (isempty (lb))
95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	224 Ain = [Ain; -eye(n)];
95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	225 bin = [bin; -ub];
95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	226 endif
5289 b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	227 endif
11587 c792872f8942 all script files: untabify and strip trailing whitespace John W. Eaton <jwe@octave.org> parents: 11523 diff changeset	228
8280 5ee11a81688e qp.m: convert b <= x <= b and b <= Ax <= b to equality constraints Gabriele Pannocchia <g.pannocchia@ing.unipi.it>* parents: 7795 diff changeset	229 if (! isempty (lb) && ! isempty (ub))
10549 95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	230 rtol = sqrt (eps);
95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	231 for i = 1:n
14868 5d3a684236b0 maint: Use Octave coding conventions for cuddling parentheses in scripts directory Rik <octave@nomad.inbox5.com> parents: 14327 diff changeset	232 if (abs (lb (i) - ub(i)) < rtol*(1 + max (abs (lb(i) + ub(i)))))
8280 5ee11a81688e qp.m: convert b <= x <= b and b <= Ax <= b to equality constraints Gabriele Pannocchia <g.pannocchia@ing.unipi.it>* parents: 7795 diff changeset	233 ## These are actually an equality constraint
14868 5d3a684236b0 maint: Use Octave coding conventions for cuddling parentheses in scripts directory Rik <octave@nomad.inbox5.com> parents: 14327 diff changeset	234 tmprow = zeros (1,n);
10549 95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	235 tmprow(i) = 1;
8280 5ee11a81688e qp.m: convert b <= x <= b and b <= Ax <= b to equality constraints Gabriele Pannocchia <g.pannocchia@ing.unipi.it>* parents: 7795 diff changeset	236 A = [A;tmprow];
5ee11a81688e qp.m: convert b <= x <= b and b <= Ax <= b to equality constraints Gabriele Pannocchia <g.pannocchia@ing.unipi.it>* parents: 7795 diff changeset	237 b = [b; 0.5*(lb(i) + ub(i))];
10549 95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	238 n_eq = n_eq + 1;
95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	239 else
14868 5d3a684236b0 maint: Use Octave coding conventions for cuddling parentheses in scripts directory Rik <octave@nomad.inbox5.com> parents: 14327 diff changeset	240 tmprow = zeros (1,n);
10549 95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	241 tmprow(i) = 1;
95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	242 Ain = [Ain; tmprow; -tmprow];
95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	243 bin = [bin; lb(i); -ub(i)];
95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	244 n_in = n_in + 2;
95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	245 endif
95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	246 endfor
8280 5ee11a81688e qp.m: convert b <= x <= b and b <= Ax <= b to equality constraints Gabriele Pannocchia <g.pannocchia@ing.unipi.it>* parents: 7795 diff changeset	247 endif
5289 b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	248 endif
b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	249
b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	250 ## Inequality constraints
10059 665ad34efeed qp.m: handle optimset options Joshua Redstone <redstone@gmail.com>, John W. Eaton <jwe@octave.org> parents: 9872 diff changeset	251 if (nargs > 7)
5289 b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	252 [dimA_in, n1] = size (A_in);
b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	253 if (n1 != n)
10549 95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	254 error ("qp: inequality constraint matrix has incorrect column dimension");
5289 b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	255 else
10549 95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	256 if (! isempty (A_lb))
95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	257 if (numel (A_lb) != dimA_in)
95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	258 error ("qp: inequality constraint matrix and lower bound vector inconsistent");
95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	259 elseif (isempty (A_ub))
95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	260 Ain = [Ain; A_in];
95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	261 bin = [bin; A_lb];
95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	262 endif
95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	263 endif
95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	264 if (! isempty (A_ub))
95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	265 if (numel (A_ub) != dimA_in)
95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	266 error ("qp: inequality constraint matrix and upper bound vector inconsistent");
95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	267 elseif (isempty (A_lb))
95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	268 Ain = [Ain; -A_in];
95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	269 bin = [bin; -A_ub];
95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	270 endif
95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	271 endif
11587 c792872f8942 all script files: untabify and strip trailing whitespace John W. Eaton <jwe@octave.org> parents: 11523 diff changeset	272
10549 95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	273 if (! isempty (A_lb) && ! isempty (A_ub))
95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	274 rtol = sqrt (eps);
95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	275 for i = 1:dimA_in
95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	276 if (abs (A_lb(i) - A_ub(i)) < rtol*(1 + max (abs (A_lb(i) + A_ub(i)))))
8280 5ee11a81688e qp.m: convert b <= x <= b and b <= Ax <= b to equality constraints Gabriele Pannocchia <g.pannocchia@ing.unipi.it>* parents: 7795 diff changeset	277 ## These are actually an equality constraint
10549 95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	278 tmprow = A_in(i,:);
8280 5ee11a81688e qp.m: convert b <= x <= b and b <= Ax <= b to equality constraints Gabriele Pannocchia <g.pannocchia@ing.unipi.it>* parents: 7795 diff changeset	279 A = [A;tmprow];
5ee11a81688e qp.m: convert b <= x <= b and b <= Ax <= b to equality constraints Gabriele Pannocchia <g.pannocchia@ing.unipi.it>* parents: 7795 diff changeset	280 b = [b; 0.5*(A_lb(i) + A_ub(i))];
10549 95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	281 n_eq = n_eq + 1;
95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	282 else
95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	283 tmprow = A_in(i,:);
95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	284 Ain = [Ain; tmprow; -tmprow];
95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	285 bin = [bin; A_lb(i); -A_ub(i)];
95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	286 n_in = n_in + 2;
95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	287 endif
95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	288 endfor
95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	289 endif
5289 b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	290 endif
b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	291 endif
b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	292
b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	293 ## Now we should have the following QP:
b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	294 ##
b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	295 ## min_x 0.5x'Hx + x'q
b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	296 ## s.t. A*x = b
b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	297 ## Ain*x >= bin
b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	298
b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	299 ## Discard inequality constraints that have -Inf bounds since those
b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	300 ## will never be active.
b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	301 idx = isinf (bin) & bin < 0;
6523 589c8dbba916 [project @ 2007-04-12 19:06:44 by jwe] jwe parents: 6417 diff changeset	302
6526 2a715c6409a5 [project @ 2007-04-13 20:31:50 by jwe] jwe parents: 6523 diff changeset	303 bin(idx) = [];
2a715c6409a5 [project @ 2007-04-13 20:31:50 by jwe] jwe parents: 6523 diff changeset	304 Ain(idx,:) = [];
5289 b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	305
10059 665ad34efeed qp.m: handle optimset options Joshua Redstone <redstone@gmail.com>, John W. Eaton <jwe@octave.org> parents: 9872 diff changeset	306 n_in = numel (bin);
5310 2fbcdc356fc7 [project @ 2005-04-27 15:38:53 by jwe] jwe parents: 5307 diff changeset	307
5289 b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	308 ## Check if the initial guess is feasible.
7795 df9519e9990c Handle single precision eps values David Bateman <dbateman@free.fr> parents: 7017 diff changeset	309 if (isa (x0, "single") \|\| isa (H, "single") \|\| isa (q, "single") \|\| isa (A, "single")
10549 95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	310 \|\| isa (b, "single"))
7795 df9519e9990c Handle single precision eps values David Bateman <dbateman@free.fr> parents: 7017 diff changeset	311 rtol = sqrt (eps ("single"));
df9519e9990c Handle single precision eps values David Bateman <dbateman@free.fr> parents: 7017 diff changeset	312 else
df9519e9990c Handle single precision eps values David Bateman <dbateman@free.fr> parents: 7017 diff changeset	313 rtol = sqrt (eps);
df9519e9990c Handle single precision eps values David Bateman <dbateman@free.fr> parents: 7017 diff changeset	314 endif
5289 b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	315
8280 5ee11a81688e qp.m: convert b <= x <= b and b <= Ax <= b to equality constraints Gabriele Pannocchia <g.pannocchia@ing.unipi.it>* parents: 7795 diff changeset	316 eq_infeasible = (n_eq > 0 && norm (Ax0-b) > rtol(1+abs (b)));
5ee11a81688e qp.m: convert b <= x <= b and b <= Ax <= b to equality constraints Gabriele Pannocchia <g.pannocchia@ing.unipi.it>* parents: 7795 diff changeset	317 in_infeasible = (n_in > 0 && any (Ainx0-bin < -rtol(1+abs (bin))));
5289 b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	318
b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	319 info = 0;
b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	320 if (eq_infeasible \|\| in_infeasible)
b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	321 ## The initial guess is not feasible.
b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	322 ## First define xbar that is feasible with respect to the equality
b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	323 ## constraints.
b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	324 if (eq_infeasible)
10549 95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	325 if (rank (A) < n_eq)
11589 b0084095098e missing semicolons in script files John W. Eaton <jwe@octave.org> parents: 11587 diff changeset	326 error ("qp: equality constraint matrix must be full row rank");
10549 95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	327 endif
95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	328 xbar = pinv (A) * b;
5289 b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	329 else
10549 95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	330 xbar = x0;
5289 b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	331 endif
b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	332
b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	333 ## Check if xbar is feasible with respect to the inequality
b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	334 ## constraints also.
b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	335 if (n_in > 0)
10549 95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	336 res = Ain * xbar - bin;
95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	337 if (any (res < -rtol * (1 + abs (bin))))
95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	338 ## xbar is not feasible with respect to the inequality
95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	339 ## constraints. Compute a step in the null space of the
95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	340 ## equality constraints, by solving a QP. If the slack is
95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	341 ## small, we have a feasible initial guess. Otherwise, the
95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	342 ## problem is infeasible.
95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	343 if (n_eq > 0)
95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	344 Z = null (A);
95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	345 if (isempty (Z))
95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	346 ## The problem is infeasible because A is square and full
95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	347 ## rank, but xbar is not feasible.
95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	348 info = 6;
95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	349 endif
95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	350 endif
5289 b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	351
10549 95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	352 if (info != 6)
5289 b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	353 ## Solve an LP with additional slack variables to find
10549 95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	354 ## a feasible starting point.
95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	355 gamma = eye (n_in);
95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	356 if (n_eq > 0)
95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	357 Atmp = [Ain*Z, gamma];
95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	358 btmp = -res;
95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	359 else
95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	360 Atmp = [Ain, gamma];
95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	361 btmp = bin;
95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	362 endif
95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	363 ctmp = [zeros(n-n_eq, 1); ones(n_in, 1)];
95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	364 lb = [-Inf(n-n_eq,1); zeros(n_in,1)];
95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	365 ub = [];
95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	366 ctype = repmat ("L", n_in, 1);
95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	367 [P, dummy, status] = glpk (ctmp, Atmp, btmp, lb, ub, ctype);
17897 185038fe7a16 qp.m: Fix test on GLPK's return status (bug #40536) Julien Bect <julien.bect@supelec.fr> parents: 17744 diff changeset	368 if ((status == 0)
10549 95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	369 && all (abs (P(n-n_eq+1:end)) < rtol * (1 + norm (btmp))))
95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	370 ## We found a feasible starting point
95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	371 if (n_eq > 0)
95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	372 x0 = xbar + Z*P(1:n-n_eq);
95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	373 else
95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	374 x0 = P(1:n);
5289 b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	375 endif
10549 95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	376 else
95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	377 ## The problem is infeasible
95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	378 info = 6;
95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	379 endif
95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	380 endif
95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	381 else
95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	382 ## xbar is feasible. We use it a starting point.
95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	383 x0 = xbar;
95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	384 endif
5289 b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	385 else
10549 95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	386 ## xbar is feasible. We use it a starting point.
95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10540 diff changeset	387 x0 = xbar;
5289 b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	388 endif
b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	389 endif
b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	390
b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	391 if (info == 0)
b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	392 ## The initial (or computed) guess is feasible.
b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	393 ## We call the solver.
b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	394 [x, lambda, info, iter] = __qp__ (x0, H, q, A, b, Ain, bin, maxit);
b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	395 else
b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	396 iter = 0;
b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	397 x = x0;
b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	398 lambda = [];
b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	399 endif
b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	400 obj = 0.5 * x' * H * x + q' * x;
b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	401 INFO.solveiter = iter;
b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	402 INFO.info = info;
b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	403
b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	404 else
6046 34f96dd5441b [project @ 2006-10-10 16:10:25 by jwe] jwe parents: 5310 diff changeset	405 print_usage ();
5289 b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	406 endif
b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	407
b98bf1d70a0a [project @ 2005-04-21 14:43:14 by jwe] jwe parents: diff changeset	408 endfunction
17338 1c89599167a6 maint: End m-files with 1 blank line. Rik <rik@octave.org> parents: 14868 diff changeset	409
17897 185038fe7a16 qp.m: Fix test on GLPK's return status (bug #40536) Julien Bect <julien.bect@supelec.fr> parents: 17744 diff changeset	410
17931 a404853d2073 qp.m: Condition test on HAVE_GLPK Mike Miller <mtmiller@ieee.org> parents: 17897 diff changeset	411 %% Test infeasible initial guess (bug #40536)
a404853d2073 qp.m: Condition test on HAVE_GLPK Mike Miller <mtmiller@ieee.org> parents: 17897 diff changeset	412 %!testif HAVE_GLPK
17897 185038fe7a16 qp.m: Fix test on GLPK's return status (bug #40536) Julien Bect <julien.bect@supelec.fr> parents: 17744 diff changeset	413 %!
185038fe7a16 qp.m: Fix test on GLPK's return status (bug #40536) Julien Bect <julien.bect@supelec.fr> parents: 17744 diff changeset	414 %! H = 1; q = 0; # objective: x -> 0.5 x^2
185038fe7a16 qp.m: Fix test on GLPK's return status (bug #40536) Julien Bect <julien.bect@supelec.fr> parents: 17744 diff changeset	415 %! A = 1; lb = 1; ub = +inf; # constraint: x >= 1
185038fe7a16 qp.m: Fix test on GLPK's return status (bug #40536) Julien Bect <julien.bect@supelec.fr> parents: 17744 diff changeset	416 %! x0 = 0; # infeasible initial guess
185038fe7a16 qp.m: Fix test on GLPK's return status (bug #40536) Julien Bect <julien.bect@supelec.fr> parents: 17744 diff changeset	417 %!
185038fe7a16 qp.m: Fix test on GLPK's return status (bug #40536) Julien Bect <julien.bect@supelec.fr> parents: 17744 diff changeset	418 %! [x, obj_qp, INFO, lambda] = qp (x0, H, q, [], [], [], [], lb, A, ub);
185038fe7a16 qp.m: Fix test on GLPK's return status (bug #40536) Julien Bect <julien.bect@supelec.fr> parents: 17744 diff changeset	419 %!
185038fe7a16 qp.m: Fix test on GLPK's return status (bug #40536) Julien Bect <julien.bect@supelec.fr> parents: 17744 diff changeset	420 %! assert (isstruct (INFO) && isfield (INFO, "info") && (INFO.info == 0));
185038fe7a16 qp.m: Fix test on GLPK's return status (bug #40536) Julien Bect <julien.bect@supelec.fr> parents: 17744 diff changeset	421 %! assert ([x obj_qp], [1.0 0.5], eps);
185038fe7a16 qp.m: Fix test on GLPK's return status (bug #40536) Julien Bect <julien.bect@supelec.fr> parents: 17744 diff changeset	422

Mercurial > hg > octave-nkf

annotate scripts/optimization/qp.m @ 19898:4197fc428c7d