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