Mercurial > octave
annotate scripts/optimization/fminunc.m @ 26376:00f796120a6d stable
maint: Update copyright dates in all source files.
author | John W. Eaton <jwe@octave.org> |
---|---|
date | Wed, 02 Jan 2019 16:32:43 -0500 |
parents | 804e18e3e320 |
children | 0a62d9a6aa2d 757a7119e319 |
rev | line source |
---|---|
26376
00f796120a6d
maint: Update copyright dates in all source files.
John W. Eaton <jwe@octave.org>
parents:
26138
diff
changeset
|
1 ## Copyright (C) 2008-2019 VZLU Prague, a.s. |
9084
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
2 ## |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
3 ## This file is part of Octave. |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
4 ## |
24534
194eb4bd202b
maint: Update punctuation for GPL v3 license text.
Rik <rik@octave.org>
parents:
23220
diff
changeset
|
5 ## Octave is free software: you can redistribute it and/or modify it |
9084
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
6 ## under the terms of the GNU General Public License as published by |
24534
194eb4bd202b
maint: Update punctuation for GPL v3 license text.
Rik <rik@octave.org>
parents:
23220
diff
changeset
|
7 ## the Free Software Foundation, either version 3 of the License, or |
22755
3a2b891d0b33
maint: Standardize Copyright formatting.
Rik <rik@octave.org>
parents:
22323
diff
changeset
|
8 ## (at your option) any later version. |
9084
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
9 ## |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
10 ## Octave is distributed in the hope that it will be useful, but |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
11 ## WITHOUT ANY WARRANTY; without even the implied warranty of |
22755
3a2b891d0b33
maint: Standardize Copyright formatting.
Rik <rik@octave.org>
parents:
22323
diff
changeset
|
12 ## MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the |
3a2b891d0b33
maint: Standardize Copyright formatting.
Rik <rik@octave.org>
parents:
22323
diff
changeset
|
13 ## GNU General Public License for more details. |
9084
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
14 ## |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
15 ## You should have received a copy of the GNU General Public License |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
16 ## along with Octave; see the file COPYING. If not, see |
24534
194eb4bd202b
maint: Update punctuation for GPL v3 license text.
Rik <rik@octave.org>
parents:
23220
diff
changeset
|
17 ## <https://www.gnu.org/licenses/>. |
9084
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
18 ## |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
19 ## Author: Jaroslav Hajek <highegg@gmail.com> |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
20 |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
21 ## -*- texinfo -*- |
20852
516bb87ea72e
2015 Code Sprint: remove class of function from docstring for all m-files.
Rik <rik@octave.org>
parents:
20735
diff
changeset
|
22 ## @deftypefn {} {} fminunc (@var{fcn}, @var{x0}) |
516bb87ea72e
2015 Code Sprint: remove class of function from docstring for all m-files.
Rik <rik@octave.org>
parents:
20735
diff
changeset
|
23 ## @deftypefnx {} {} fminunc (@var{fcn}, @var{x0}, @var{options}) |
516bb87ea72e
2015 Code Sprint: remove class of function from docstring for all m-files.
Rik <rik@octave.org>
parents:
20735
diff
changeset
|
24 ## @deftypefnx {} {[@var{x}, @var{fval}, @var{info}, @var{output}, @var{grad}, @var{hess}] =} fminunc (@var{fcn}, @dots{}) |
10793
be55736a0783
Grammarcheck the documentation from m-files.
Rik <octave@nomad.inbox5.com>
parents:
10687
diff
changeset
|
25 ## Solve an unconstrained optimization problem defined by the function |
be55736a0783
Grammarcheck the documentation from m-files.
Rik <octave@nomad.inbox5.com>
parents:
10687
diff
changeset
|
26 ## @var{fcn}. |
14895
e0525ecf156e
Add new function fminsearch.m
Andy Adler <andy@analyti.ca>
parents:
14868
diff
changeset
|
27 ## |
20165
f1d0f506ee78
doc: Update more docstrings to have one sentence summary as first line.
Rik <rik@octave.org>
parents:
19833
diff
changeset
|
28 ## @var{fcn} should accept a vector (array) defining the unknown variables, and |
f1d0f506ee78
doc: Update more docstrings to have one sentence summary as first line.
Rik <rik@octave.org>
parents:
19833
diff
changeset
|
29 ## return the objective function value, optionally with gradient. |
18609 | 30 ## @code{fminunc} attempts to determine a vector @var{x} such that |
20165
f1d0f506ee78
doc: Update more docstrings to have one sentence summary as first line.
Rik <rik@octave.org>
parents:
19833
diff
changeset
|
31 ## @code{@var{fcn} (@var{x})} is a local minimum. |
f1d0f506ee78
doc: Update more docstrings to have one sentence summary as first line.
Rik <rik@octave.org>
parents:
19833
diff
changeset
|
32 ## |
f1d0f506ee78
doc: Update more docstrings to have one sentence summary as first line.
Rik <rik@octave.org>
parents:
19833
diff
changeset
|
33 ## @var{x0} determines a starting guess. The shape of @var{x0} is preserved in |
f1d0f506ee78
doc: Update more docstrings to have one sentence summary as first line.
Rik <rik@octave.org>
parents:
19833
diff
changeset
|
34 ## all calls to @var{fcn}, but otherwise is treated as a column vector. |
f1d0f506ee78
doc: Update more docstrings to have one sentence summary as first line.
Rik <rik@octave.org>
parents:
19833
diff
changeset
|
35 ## |
26137
1ae11ca7dceb
fminunc.m: Change algorithm defaults to match Matlab.
Rik <rik@octave.org>
parents:
26118
diff
changeset
|
36 ## @var{options} is a structure specifying additional parameters which |
1ae11ca7dceb
fminunc.m: Change algorithm defaults to match Matlab.
Rik <rik@octave.org>
parents:
26118
diff
changeset
|
37 ## control the algorithm. Currently, @code{fminunc} recognizes these options: |
1ae11ca7dceb
fminunc.m: Change algorithm defaults to match Matlab.
Rik <rik@octave.org>
parents:
26118
diff
changeset
|
38 ## @qcode{"AutoScaling"}, @qcode{"FinDiffType"}, @qcode{"FunValCheck"}, |
1ae11ca7dceb
fminunc.m: Change algorithm defaults to match Matlab.
Rik <rik@octave.org>
parents:
26118
diff
changeset
|
39 ## @qcode{"GradObj"}, @qcode{"MaxFunEvals"}, @qcode{"MaxIter"}, |
1ae11ca7dceb
fminunc.m: Change algorithm defaults to match Matlab.
Rik <rik@octave.org>
parents:
26118
diff
changeset
|
40 ## @qcode{"OutputFcn"}, @qcode{"TolFun"}, @qcode{"TolX"}, @qcode{"TypicalX"}. |
1ae11ca7dceb
fminunc.m: Change algorithm defaults to match Matlab.
Rik <rik@octave.org>
parents:
26118
diff
changeset
|
41 ## |
1ae11ca7dceb
fminunc.m: Change algorithm defaults to match Matlab.
Rik <rik@octave.org>
parents:
26118
diff
changeset
|
42 ## If @qcode{"AutoScaling"} is @qcode{"on"}, the variables will be |
1ae11ca7dceb
fminunc.m: Change algorithm defaults to match Matlab.
Rik <rik@octave.org>
parents:
26118
diff
changeset
|
43 ## automatically scaled according to the column norms of the (estimated) |
1ae11ca7dceb
fminunc.m: Change algorithm defaults to match Matlab.
Rik <rik@octave.org>
parents:
26118
diff
changeset
|
44 ## Jacobian. As a result, @qcode{"TolFun"} becomes scaling-independent. By |
1ae11ca7dceb
fminunc.m: Change algorithm defaults to match Matlab.
Rik <rik@octave.org>
parents:
26118
diff
changeset
|
45 ## default, this option is @qcode{"off"} because it may sometimes deliver |
1ae11ca7dceb
fminunc.m: Change algorithm defaults to match Matlab.
Rik <rik@octave.org>
parents:
26118
diff
changeset
|
46 ## unexpected (though mathematically correct) results. |
9084
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
47 ## |
26137
1ae11ca7dceb
fminunc.m: Change algorithm defaults to match Matlab.
Rik <rik@octave.org>
parents:
26118
diff
changeset
|
48 ## If @qcode{"GradObj"} is @qcode{"on"}, it specifies that @var{fcn}--when |
1ae11ca7dceb
fminunc.m: Change algorithm defaults to match Matlab.
Rik <rik@octave.org>
parents:
26118
diff
changeset
|
49 ## called with two output arguments---also returns the Jacobian matrix of |
1ae11ca7dceb
fminunc.m: Change algorithm defaults to match Matlab.
Rik <rik@octave.org>
parents:
26118
diff
changeset
|
50 ## partial first derivatives at the requested point. |
1ae11ca7dceb
fminunc.m: Change algorithm defaults to match Matlab.
Rik <rik@octave.org>
parents:
26118
diff
changeset
|
51 ## |
1ae11ca7dceb
fminunc.m: Change algorithm defaults to match Matlab.
Rik <rik@octave.org>
parents:
26118
diff
changeset
|
52 ## @qcode{"MaxFunEvals"} proscribes the maximum number of function evaluations |
1ae11ca7dceb
fminunc.m: Change algorithm defaults to match Matlab.
Rik <rik@octave.org>
parents:
26118
diff
changeset
|
53 ## before optimization is halted. The default value is |
1ae11ca7dceb
fminunc.m: Change algorithm defaults to match Matlab.
Rik <rik@octave.org>
parents:
26118
diff
changeset
|
54 ## @code{100 * number_of_variables}, i.e., @code{100 * length (@var{x0})}. |
1ae11ca7dceb
fminunc.m: Change algorithm defaults to match Matlab.
Rik <rik@octave.org>
parents:
26118
diff
changeset
|
55 ## The value must be a positive integer. |
1ae11ca7dceb
fminunc.m: Change algorithm defaults to match Matlab.
Rik <rik@octave.org>
parents:
26118
diff
changeset
|
56 ## |
1ae11ca7dceb
fminunc.m: Change algorithm defaults to match Matlab.
Rik <rik@octave.org>
parents:
26118
diff
changeset
|
57 ## @qcode{"MaxIter"} proscribes the maximum number of algorithm iterations |
1ae11ca7dceb
fminunc.m: Change algorithm defaults to match Matlab.
Rik <rik@octave.org>
parents:
26118
diff
changeset
|
58 ## before optimization is halted. The default value is 400. |
1ae11ca7dceb
fminunc.m: Change algorithm defaults to match Matlab.
Rik <rik@octave.org>
parents:
26118
diff
changeset
|
59 ## The value must be a positive integer. |
1ae11ca7dceb
fminunc.m: Change algorithm defaults to match Matlab.
Rik <rik@octave.org>
parents:
26118
diff
changeset
|
60 ## |
1ae11ca7dceb
fminunc.m: Change algorithm defaults to match Matlab.
Rik <rik@octave.org>
parents:
26118
diff
changeset
|
61 ## @qcode{"TolX"} specifies the termination tolerance for the unknown variables |
1ae11ca7dceb
fminunc.m: Change algorithm defaults to match Matlab.
Rik <rik@octave.org>
parents:
26118
diff
changeset
|
62 ## @var{x}, while @qcode{"TolFun"} is a tolerance for the objective function |
1ae11ca7dceb
fminunc.m: Change algorithm defaults to match Matlab.
Rik <rik@octave.org>
parents:
26118
diff
changeset
|
63 ## value @var{fval}. The default is @code{1e-6} for both options. |
11587
c792872f8942
all script files: untabify and strip trailing whitespace
John W. Eaton <jwe@octave.org>
parents:
11523
diff
changeset
|
64 ## |
18609 | 65 ## For a description of the other options, see @code{optimset}. |
9084
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
66 ## |
18609 | 67 ## On return, @var{x} is the location of the minimum and @var{fval} contains |
20165
f1d0f506ee78
doc: Update more docstrings to have one sentence summary as first line.
Rik <rik@octave.org>
parents:
19833
diff
changeset
|
68 ## the value of the objective function at @var{x}. |
f1d0f506ee78
doc: Update more docstrings to have one sentence summary as first line.
Rik <rik@octave.org>
parents:
19833
diff
changeset
|
69 ## |
f1d0f506ee78
doc: Update more docstrings to have one sentence summary as first line.
Rik <rik@octave.org>
parents:
19833
diff
changeset
|
70 ## @var{info} may be one of the following values: |
11587
c792872f8942
all script files: untabify and strip trailing whitespace
John W. Eaton <jwe@octave.org>
parents:
11523
diff
changeset
|
71 ## |
9084
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
72 ## @table @asis |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
73 ## @item 1 |
10793
be55736a0783
Grammarcheck the documentation from m-files.
Rik <octave@nomad.inbox5.com>
parents:
10687
diff
changeset
|
74 ## Converged to a solution point. Relative gradient error is less than |
18609 | 75 ## specified by @code{TolFun}. |
10821
693e22af08ae
Grammarcheck documentation of m-files
Rik <octave@nomad.inbox5.com>
parents:
10793
diff
changeset
|
76 ## |
9084
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
77 ## @item 2 |
18609 | 78 ## Last relative step size was less than @code{TolX}. |
10821
693e22af08ae
Grammarcheck documentation of m-files
Rik <octave@nomad.inbox5.com>
parents:
10793
diff
changeset
|
79 ## |
9084
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
80 ## @item 3 |
18609 | 81 ## Last relative change in function value was less than @code{TolFun}. |
10821
693e22af08ae
Grammarcheck documentation of m-files
Rik <octave@nomad.inbox5.com>
parents:
10793
diff
changeset
|
82 ## |
9084
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
83 ## @item 0 |
18816
e275d15c27b5
doc: Periodic spellcheck of documentation.
Rik <rik@octave.org>
parents:
18609
diff
changeset
|
84 ## Iteration limit exceeded---either maximum number of algorithm iterations |
18609 | 85 ## @code{MaxIter} or maximum number of function evaluations @code{MaxFunEvals}. |
86 ## | |
87 ## @item -1 | |
18816
e275d15c27b5
doc: Periodic spellcheck of documentation.
Rik <rik@octave.org>
parents:
18609
diff
changeset
|
88 ## Algorithm terminated by @code{OutputFcn}. |
10821
693e22af08ae
Grammarcheck documentation of m-files
Rik <octave@nomad.inbox5.com>
parents:
10793
diff
changeset
|
89 ## |
9084
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
90 ## @item -3 |
11587
c792872f8942
all script files: untabify and strip trailing whitespace
John W. Eaton <jwe@octave.org>
parents:
11523
diff
changeset
|
91 ## The trust region radius became excessively small. |
9084
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
92 ## @end table |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
93 ## |
21546
f7f97d7e9294
doc: Wrap m-file docstrings to 79 characters + newline (80 total).
Rik <rik@octave.org>
parents:
21178
diff
changeset
|
94 ## Optionally, @code{fminunc} can return a structure with convergence |
f7f97d7e9294
doc: Wrap m-file docstrings to 79 characters + newline (80 total).
Rik <rik@octave.org>
parents:
21178
diff
changeset
|
95 ## statistics (@var{output}), the output gradient (@var{grad}) at the |
f7f97d7e9294
doc: Wrap m-file docstrings to 79 characters + newline (80 total).
Rik <rik@octave.org>
parents:
21178
diff
changeset
|
96 ## solution @var{x}, and approximate Hessian (@var{hess}) at the solution |
f7f97d7e9294
doc: Wrap m-file docstrings to 79 characters + newline (80 total).
Rik <rik@octave.org>
parents:
21178
diff
changeset
|
97 ## @var{x}. |
9627
5bcfa0b346e8
fix extra outputs in fminunc
Jaroslav Hajek <highegg@gmail.com>
parents:
9623
diff
changeset
|
98 ## |
20695
6faaab833605
fminunc.m: Clean up function to meet Octave coding standards.
Rik <rik@octave.org>
parents:
20694
diff
changeset
|
99 ## Application Notes: If the objective function is a single nonlinear equation |
6faaab833605
fminunc.m: Clean up function to meet Octave coding standards.
Rik <rik@octave.org>
parents:
20694
diff
changeset
|
100 ## of one variable then using @code{fminbnd} is usually a better choice. |
20165
f1d0f506ee78
doc: Update more docstrings to have one sentence summary as first line.
Rik <rik@octave.org>
parents:
19833
diff
changeset
|
101 ## |
22092
d3cef63f79ac
* fminunc.m: Fix incorrect reference to fminunc as fminsearch (bug #48479)
Mike Miller <mtmiller@octave.org>
parents:
21758
diff
changeset
|
102 ## The algorithm used by @code{fminunc} is a gradient search which depends |
20165
f1d0f506ee78
doc: Update more docstrings to have one sentence summary as first line.
Rik <rik@octave.org>
parents:
19833
diff
changeset
|
103 ## on the objective function being differentiable. If the function has |
f1d0f506ee78
doc: Update more docstrings to have one sentence summary as first line.
Rik <rik@octave.org>
parents:
19833
diff
changeset
|
104 ## discontinuities it may be better to use a derivative-free algorithm such as |
f1d0f506ee78
doc: Update more docstrings to have one sentence summary as first line.
Rik <rik@octave.org>
parents:
19833
diff
changeset
|
105 ## @code{fminsearch}. |
14895
e0525ecf156e
Add new function fminsearch.m
Andy Adler <andy@analyti.ca>
parents:
14868
diff
changeset
|
106 ## @seealso{fminbnd, fminsearch, optimset} |
9084
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
107 ## @end deftypefn |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
108 |
13027
b9a89ca0fb75
prevent optimization functions from setting ans in workspace at startup
John W. Eaton <jwe@octave.org>
parents:
12578
diff
changeset
|
109 ## 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:
12578
diff
changeset
|
110 ## PKG_ADD: [~] = __all_opts__ ("fminunc"); |
9084
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
111 |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
112 function [x, fval, info, output, grad, hess] = fminunc (fcn, x0, options = struct ()) |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
113 |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
114 ## Get default options if requested. |
20695
6faaab833605
fminunc.m: Clean up function to meet Octave coding standards.
Rik <rik@octave.org>
parents:
20694
diff
changeset
|
115 if (nargin == 1 && strcmp (fcn, "defaults")) |
26138
804e18e3e320
Reenable query of optimization options (bugs #54952 and #55089).
Kai T. Ohlhus <k.ohlhus@gmail.com>
parents:
26137
diff
changeset
|
116 x = struct ("AutoScaling", "off", "FunValCheck", "off", |
804e18e3e320
Reenable query of optimization options (bugs #54952 and #55089).
Kai T. Ohlhus <k.ohlhus@gmail.com>
parents:
26137
diff
changeset
|
117 "FinDiffType", "forward", "GradObj", "off", |
804e18e3e320
Reenable query of optimization options (bugs #54952 and #55089).
Kai T. Ohlhus <k.ohlhus@gmail.com>
parents:
26137
diff
changeset
|
118 "MaxFunEvals", [], "MaxIter", 400, "OutputFcn", [], |
804e18e3e320
Reenable query of optimization options (bugs #54952 and #55089).
Kai T. Ohlhus <k.ohlhus@gmail.com>
parents:
26137
diff
changeset
|
119 "TolFun", 1e-6, "TolX", 1e-6, "TypicalX", []); |
9084
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
120 return; |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
121 endif |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
122 |
21171
2935d56203a4
Fix regressions caused by ismatrix definition change (partial fix bug #47036).
Rik <rik@octave.org>
parents:
20165
diff
changeset
|
123 if (nargin < 2 || nargin > 3 || ! isnumeric (x0)) |
9084
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
124 print_usage (); |
11587
c792872f8942
all script files: untabify and strip trailing whitespace
John W. Eaton <jwe@octave.org>
parents:
11523
diff
changeset
|
125 endif |
9084
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
126 |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
127 if (ischar (fcn)) |
26118
7502fce4cd3a
str2func: eliminate optional second "global" argument
John W. Eaton <jwe@octave.org>
parents:
25054
diff
changeset
|
128 fcn = str2func (fcn); |
9084
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
129 endif |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
130 |
20695
6faaab833605
fminunc.m: Clean up function to meet Octave coding standards.
Rik <rik@octave.org>
parents:
20694
diff
changeset
|
131 xsz = size (x0); |
9084
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
132 n = numel (x0); |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
133 |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
134 has_grad = strcmpi (optimget (options, "GradObj", "off"), "on"); |
26137
1ae11ca7dceb
fminunc.m: Change algorithm defaults to match Matlab.
Rik <rik@octave.org>
parents:
26118
diff
changeset
|
135 cdif = strcmpi (optimget (options, "FinDiffType", "forward"), "central"); |
9084
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
136 maxiter = optimget (options, "MaxIter", 400); |
26137
1ae11ca7dceb
fminunc.m: Change algorithm defaults to match Matlab.
Rik <rik@octave.org>
parents:
26118
diff
changeset
|
137 maxfev = optimget (options, "MaxFunEvals", 100*n); |
9084
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
138 outfcn = optimget (options, "OutputFcn"); |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
139 |
20695
6faaab833605
fminunc.m: Clean up function to meet Octave coding standards.
Rik <rik@octave.org>
parents:
20694
diff
changeset
|
140 ## Get scaling matrix using the TypicalX option. If set to "auto", the |
6faaab833605
fminunc.m: Clean up function to meet Octave coding standards.
Rik <rik@octave.org>
parents:
20694
diff
changeset
|
141 ## scaling matrix is estimated using the Jacobian. |
10201
5c66978f3fdf
support TypicalX and AutoScaling in fsolve/fminunc, don't autoscale by default
Jaroslav Hajek <highegg@gmail.com>
parents:
9899
diff
changeset
|
142 typicalx = optimget (options, "TypicalX"); |
5c66978f3fdf
support TypicalX and AutoScaling in fsolve/fminunc, don't autoscale by default
Jaroslav Hajek <highegg@gmail.com>
parents:
9899
diff
changeset
|
143 if (isempty (typicalx)) |
5c66978f3fdf
support TypicalX and AutoScaling in fsolve/fminunc, don't autoscale by default
Jaroslav Hajek <highegg@gmail.com>
parents:
9899
diff
changeset
|
144 typicalx = ones (n, 1); |
5c66978f3fdf
support TypicalX and AutoScaling in fsolve/fminunc, don't autoscale by default
Jaroslav Hajek <highegg@gmail.com>
parents:
9899
diff
changeset
|
145 endif |
5c66978f3fdf
support TypicalX and AutoScaling in fsolve/fminunc, don't autoscale by default
Jaroslav Hajek <highegg@gmail.com>
parents:
9899
diff
changeset
|
146 autoscale = strcmpi (optimget (options, "AutoScaling", "off"), "on"); |
5c66978f3fdf
support TypicalX and AutoScaling in fsolve/fminunc, don't autoscale by default
Jaroslav Hajek <highegg@gmail.com>
parents:
9899
diff
changeset
|
147 if (! autoscale) |
5c66978f3fdf
support TypicalX and AutoScaling in fsolve/fminunc, don't autoscale by default
Jaroslav Hajek <highegg@gmail.com>
parents:
9899
diff
changeset
|
148 dg = 1 ./ typicalx; |
5c66978f3fdf
support TypicalX and AutoScaling in fsolve/fminunc, don't autoscale by default
Jaroslav Hajek <highegg@gmail.com>
parents:
9899
diff
changeset
|
149 endif |
5c66978f3fdf
support TypicalX and AutoScaling in fsolve/fminunc, don't autoscale by default
Jaroslav Hajek <highegg@gmail.com>
parents:
9899
diff
changeset
|
150 |
9084
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
151 funvalchk = strcmpi (optimget (options, "FunValCheck", "off"), "on"); |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
152 |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
153 if (funvalchk) |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
154 ## Replace fcn with a guarded version. |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
155 fcn = @(x) guarded_eval (fcn, x); |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
156 endif |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
157 |
20695
6faaab833605
fminunc.m: Clean up function to meet Octave coding standards.
Rik <rik@octave.org>
parents:
20694
diff
changeset
|
158 ## These defaults are rather stringent. I think that normally, user |
9084
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
159 ## prefers accuracy to performance. |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
160 |
10201
5c66978f3fdf
support TypicalX and AutoScaling in fsolve/fminunc, don't autoscale by default
Jaroslav Hajek <highegg@gmail.com>
parents:
9899
diff
changeset
|
161 tolx = optimget (options, "TolX", 1e-7); |
5c66978f3fdf
support TypicalX and AutoScaling in fsolve/fminunc, don't autoscale by default
Jaroslav Hajek <highegg@gmail.com>
parents:
9899
diff
changeset
|
162 tolf = optimget (options, "TolFun", 1e-7); |
9084
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
163 |
9623
bc0739d02724
update initial TR step for fsolve and fminunc
Jaroslav Hajek <highegg@gmail.com>
parents:
9464
diff
changeset
|
164 factor = 0.1; |
9084
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
165 ## FIXME: TypicalX corresponds to user scaling (???) |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
166 autodg = true; |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
167 |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
168 niter = 1; |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
169 nfev = 0; |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
170 |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
171 x = x0(:); |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
172 info = 0; |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
173 |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
174 ## Initial evaluation. |
20695
6faaab833605
fminunc.m: Clean up function to meet Octave coding standards.
Rik <rik@octave.org>
parents:
20694
diff
changeset
|
175 fval = fcn (reshape (x, xsz)); |
9084
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
176 n = length (x); |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
177 |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
178 if (! isempty (outfcn)) |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
179 optimvalues.iter = niter; |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
180 optimvalues.funccount = nfev; |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
181 optimvalues.fval = fval; |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
182 optimvalues.searchdirection = zeros (n, 1); |
20695
6faaab833605
fminunc.m: Clean up function to meet Octave coding standards.
Rik <rik@octave.org>
parents:
20694
diff
changeset
|
183 state = "init"; |
9084
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
184 stop = outfcn (x, optimvalues, state); |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
185 if (stop) |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
186 info = -1; |
22785
9c6661004167
error if break statement is in script file separate from loop (bug #39168)
John W. Eaton <jwe@octave.org>
parents:
22755
diff
changeset
|
187 return; |
9084
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
188 endif |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
189 endif |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
190 |
21099
52af4092f863
For optimization scripts, correctly choose tolerance (eps) based on class of fun and X0.
Rik <rik@octave.org>
parents:
20852
diff
changeset
|
191 if (isa (x0, "single") || isa (fval, "single")) |
52af4092f863
For optimization scripts, correctly choose tolerance (eps) based on class of fun and X0.
Rik <rik@octave.org>
parents:
20852
diff
changeset
|
192 macheps = eps ("single"); |
52af4092f863
For optimization scripts, correctly choose tolerance (eps) based on class of fun and X0.
Rik <rik@octave.org>
parents:
20852
diff
changeset
|
193 else |
52af4092f863
For optimization scripts, correctly choose tolerance (eps) based on class of fun and X0.
Rik <rik@octave.org>
parents:
20852
diff
changeset
|
194 macheps = eps ("double"); |
52af4092f863
For optimization scripts, correctly choose tolerance (eps) based on class of fun and X0.
Rik <rik@octave.org>
parents:
20852
diff
changeset
|
195 endif |
52af4092f863
For optimization scripts, correctly choose tolerance (eps) based on class of fun and X0.
Rik <rik@octave.org>
parents:
20852
diff
changeset
|
196 |
9084
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
197 nsuciter = 0; |
9199
399884c9d4a1
import the step adaptation strategy from fsolve to fminunc
Jaroslav Hajek <highegg@gmail.com>
parents:
9084
diff
changeset
|
198 lastratio = 0; |
9084
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
199 |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
200 grad = []; |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
201 |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
202 ## Outer loop. |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
203 while (niter < maxiter && nfev < maxfev && ! info) |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
204 |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
205 grad0 = grad; |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
206 |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
207 ## Calculate function value and gradient (possibly via FD). |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
208 if (has_grad) |
20695
6faaab833605
fminunc.m: Clean up function to meet Octave coding standards.
Rik <rik@octave.org>
parents:
20694
diff
changeset
|
209 [fval, grad] = fcn (reshape (x, xsz)); |
9084
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
210 grad = grad(:); |
20695
6faaab833605
fminunc.m: Clean up function to meet Octave coding standards.
Rik <rik@octave.org>
parents:
20694
diff
changeset
|
211 nfev += 1; |
9084
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
212 else |
20695
6faaab833605
fminunc.m: Clean up function to meet Octave coding standards.
Rik <rik@octave.org>
parents:
20694
diff
changeset
|
213 grad = __fdjac__ (fcn, reshape (x, xsz), fval, typicalx, cdif)(:); |
9212
6feb27c38da1
support central differences in fminunc and fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
9207
diff
changeset
|
214 nfev += (1 + cdif) * length (x); |
9084
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
215 endif |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
216 |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
217 if (niter == 1) |
11587
c792872f8942
all script files: untabify and strip trailing whitespace
John W. Eaton <jwe@octave.org>
parents:
11523
diff
changeset
|
218 ## Initialize by identity matrix. |
9084
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
219 hesr = eye (n); |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
220 else |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
221 ## Use the damped BFGS formula. |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
222 y = grad - grad0; |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
223 sBs = sumsq (w); |
20695
6faaab833605
fminunc.m: Clean up function to meet Octave coding standards.
Rik <rik@octave.org>
parents:
20694
diff
changeset
|
224 Bs = hesr' * w; |
6faaab833605
fminunc.m: Clean up function to meet Octave coding standards.
Rik <rik@octave.org>
parents:
20694
diff
changeset
|
225 sy = y' * s; |
9633
ecc2c556f844
simplify damped BFGS formula
Jaroslav Hajek <highegg@gmail.com>
parents:
9631
diff
changeset
|
226 theta = 0.8 / max (1 - sy / sBs, 0.8); |
9084
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
227 r = theta * y + (1-theta) * Bs; |
20695
6faaab833605
fminunc.m: Clean up function to meet Octave coding standards.
Rik <rik@octave.org>
parents:
20694
diff
changeset
|
228 hesr = cholupdate (hesr, r / sqrt (s' * r), "+"); |
9084
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
229 [hesr, info] = cholupdate (hesr, Bs / sqrt (sBs), "-"); |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
230 if (info) |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
231 hesr = eye (n); |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
232 endif |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
233 endif |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
234 |
10201
5c66978f3fdf
support TypicalX and AutoScaling in fsolve/fminunc, don't autoscale by default
Jaroslav Hajek <highegg@gmail.com>
parents:
9899
diff
changeset
|
235 if (autoscale) |
20695
6faaab833605
fminunc.m: Clean up function to meet Octave coding standards.
Rik <rik@octave.org>
parents:
20694
diff
changeset
|
236 ## Second derivatives approximate the Hessian. |
6faaab833605
fminunc.m: Clean up function to meet Octave coding standards.
Rik <rik@octave.org>
parents:
20694
diff
changeset
|
237 d2f = norm (hesr, "columns").'; |
10201
5c66978f3fdf
support TypicalX and AutoScaling in fsolve/fminunc, don't autoscale by default
Jaroslav Hajek <highegg@gmail.com>
parents:
9899
diff
changeset
|
238 if (niter == 1) |
5c66978f3fdf
support TypicalX and AutoScaling in fsolve/fminunc, don't autoscale by default
Jaroslav Hajek <highegg@gmail.com>
parents:
9899
diff
changeset
|
239 dg = d2f; |
5c66978f3fdf
support TypicalX and AutoScaling in fsolve/fminunc, don't autoscale by default
Jaroslav Hajek <highegg@gmail.com>
parents:
9899
diff
changeset
|
240 else |
5c66978f3fdf
support TypicalX and AutoScaling in fsolve/fminunc, don't autoscale by default
Jaroslav Hajek <highegg@gmail.com>
parents:
9899
diff
changeset
|
241 ## FIXME: maybe fixed lower and upper bounds? |
5c66978f3fdf
support TypicalX and AutoScaling in fsolve/fminunc, don't autoscale by default
Jaroslav Hajek <highegg@gmail.com>
parents:
9899
diff
changeset
|
242 dg = max (0.1*dg, d2f); |
5c66978f3fdf
support TypicalX and AutoScaling in fsolve/fminunc, don't autoscale by default
Jaroslav Hajek <highegg@gmail.com>
parents:
9899
diff
changeset
|
243 endif |
5c66978f3fdf
support TypicalX and AutoScaling in fsolve/fminunc, don't autoscale by default
Jaroslav Hajek <highegg@gmail.com>
parents:
9899
diff
changeset
|
244 endif |
5c66978f3fdf
support TypicalX and AutoScaling in fsolve/fminunc, don't autoscale by default
Jaroslav Hajek <highegg@gmail.com>
parents:
9899
diff
changeset
|
245 |
9084
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
246 if (niter == 1) |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
247 xn = norm (dg .* x); |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
248 ## FIXME: something better? |
9628
73e6ad869f08
further correct initial TR step strategy
Jaroslav Hajek <highegg@gmail.com>
parents:
9627
diff
changeset
|
249 delta = factor * max (xn, 1); |
9084
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
250 endif |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
251 |
18609 | 252 ## FIXME: why tolf*n*xn? If abs (e) ~ abs(x) * eps is a vector |
21751
b571fc85953f
maint: Use two spaces after period to indicate sentence break.
Rik <rik@octave.org>
parents:
21578
diff
changeset
|
253 ## of perturbations of x, then norm (hesr*e) <= eps*xn, i.e., by |
9084
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
254 ## tolf ~ eps we demand as much accuracy as we can expect. |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
255 if (norm (grad) <= tolf*n*xn) |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
256 info = 1; |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
257 break; |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
258 endif |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
259 |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
260 suc = false; |
9207
25f50d2d76b3
improve TR updating strategy for fminunc and fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
9199
diff
changeset
|
261 decfac = 0.5; |
25f50d2d76b3
improve TR updating strategy for fminunc and fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
9199
diff
changeset
|
262 |
9084
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
263 ## Inner loop. |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
264 while (! suc && niter <= maxiter && nfev < maxfev && ! info) |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
265 |
9631
00958d0c4e3c
split __dogleg__ > __doglegm__
Jaroslav Hajek <highegg@gmail.com>
parents:
9628
diff
changeset
|
266 s = - __doglegm__ (hesr, grad, dg, delta); |
9084
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
267 |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
268 sn = norm (dg .* s); |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
269 if (niter == 1) |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
270 delta = min (delta, sn); |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
271 endif |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
272 |
20695
6faaab833605
fminunc.m: Clean up function to meet Octave coding standards.
Rik <rik@octave.org>
parents:
20694
diff
changeset
|
273 fval1 = fcn (reshape (x + s, xsz))(:); |
6faaab833605
fminunc.m: Clean up function to meet Octave coding standards.
Rik <rik@octave.org>
parents:
20694
diff
changeset
|
274 nfev += 1; |
9084
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
275 |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
276 if (fval1 < fval) |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
277 ## Scaled actual reduction. |
19833
9fc020886ae9
maint: Clean up m-files to follow Octave coding conventions.
Rik <rik@octave.org>
parents:
19697
diff
changeset
|
278 actred = (fval - fval1) / (abs (fval1) + abs (fval)); |
9084
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
279 else |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
280 actred = -1; |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
281 endif |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
282 |
20695
6faaab833605
fminunc.m: Clean up function to meet Octave coding standards.
Rik <rik@octave.org>
parents:
20694
diff
changeset
|
283 w = hesr * s; |
9084
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
284 ## Scaled predicted reduction, and ratio. |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
285 t = 1/2 * sumsq (w) + grad'*s; |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
286 if (t < 0) |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
287 prered = -t/(abs (fval) + abs (fval + t)); |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
288 ratio = actred / prered; |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
289 else |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
290 prered = 0; |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
291 ratio = 0; |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
292 endif |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
293 |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
294 ## Update delta. |
14868
5d3a684236b0
maint: Use Octave coding conventions for cuddling parentheses in scripts directory
Rik <octave@nomad.inbox5.com>
parents:
14552
diff
changeset
|
295 if (ratio < min (max (0.1, 0.8*lastratio), 0.9)) |
9207
25f50d2d76b3
improve TR updating strategy for fminunc and fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
9199
diff
changeset
|
296 delta *= decfac; |
25f50d2d76b3
improve TR updating strategy for fminunc and fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
9199
diff
changeset
|
297 decfac ^= 1.4142; |
20695
6faaab833605
fminunc.m: Clean up function to meet Octave coding standards.
Rik <rik@octave.org>
parents:
20694
diff
changeset
|
298 if (delta <= 10*macheps*xn) |
9084
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
299 ## Trust region became uselessly small. |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
300 info = -3; |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
301 break; |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
302 endif |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
303 else |
9199
399884c9d4a1
import the step adaptation strategy from fsolve to fminunc
Jaroslav Hajek <highegg@gmail.com>
parents:
9084
diff
changeset
|
304 lastratio = ratio; |
9207
25f50d2d76b3
improve TR updating strategy for fminunc and fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
9199
diff
changeset
|
305 decfac = 0.5; |
9084
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
306 if (abs (1-ratio) <= 0.1) |
9207
25f50d2d76b3
improve TR updating strategy for fminunc and fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
9199
diff
changeset
|
307 delta = 1.4142*sn; |
9084
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
308 elseif (ratio >= 0.5) |
9207
25f50d2d76b3
improve TR updating strategy for fminunc and fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
9199
diff
changeset
|
309 delta = max (delta, 1.4142*sn); |
9084
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
310 endif |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
311 endif |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
312 |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
313 if (ratio >= 1e-4) |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
314 ## Successful iteration. |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
315 x += s; |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
316 xn = norm (dg .* x); |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
317 fval = fval1; |
20735
418ae0cb752f
Replace ++,-- with in-place operators for performance.
Rik <rik@octave.org>
parents:
20715
diff
changeset
|
318 nsuciter += 1; |
9084
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
319 suc = true; |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
320 endif |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
321 |
20695
6faaab833605
fminunc.m: Clean up function to meet Octave coding standards.
Rik <rik@octave.org>
parents:
20694
diff
changeset
|
322 niter += 1; |
9084
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
323 |
18609 | 324 ## FIXME: should outputfcn be called only after a successful iteration? |
9084
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
325 if (! isempty (outfcn)) |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
326 optimvalues.iter = niter; |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
327 optimvalues.funccount = nfev; |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
328 optimvalues.fval = fval; |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
329 optimvalues.searchdirection = s; |
20695
6faaab833605
fminunc.m: Clean up function to meet Octave coding standards.
Rik <rik@octave.org>
parents:
20694
diff
changeset
|
330 state = "iter"; |
9084
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
331 stop = outfcn (x, optimvalues, state); |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
332 if (stop) |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
333 info = -1; |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
334 break; |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
335 endif |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
336 endif |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
337 |
20695
6faaab833605
fminunc.m: Clean up function to meet Octave coding standards.
Rik <rik@octave.org>
parents:
20694
diff
changeset
|
338 ## Tests for termination conditions. A mysterious place, anything |
9084
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
339 ## can happen if you change something here... |
11587
c792872f8942
all script files: untabify and strip trailing whitespace
John W. Eaton <jwe@octave.org>
parents:
11523
diff
changeset
|
340 |
9084
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
341 ## The rule of thumb (which I'm not sure M*b is quite following) |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
342 ## is that for a tolerance that depends on scaling, only 0 makes |
20695
6faaab833605
fminunc.m: Clean up function to meet Octave coding standards.
Rik <rik@octave.org>
parents:
20694
diff
changeset
|
343 ## sense as a default value. But 0 usually means uselessly long |
9084
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
344 ## iterations, so we need scaling-independent tolerances wherever |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
345 ## possible. |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
346 |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
347 ## The following tests done only after successful step. |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
348 if (ratio >= 1e-4) |
20695
6faaab833605
fminunc.m: Clean up function to meet Octave coding standards.
Rik <rik@octave.org>
parents:
20694
diff
changeset
|
349 ## This one is classic. Note that we use scaled variables again, |
9084
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
350 ## but compare to scaled step, so nothing bad. |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
351 if (sn <= tolx*xn) |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
352 info = 2; |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
353 ## Again a classic one. |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
354 elseif (actred < tolf) |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
355 info = 3; |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
356 endif |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
357 endif |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
358 |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
359 endwhile |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
360 endwhile |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
361 |
20695
6faaab833605
fminunc.m: Clean up function to meet Octave coding standards.
Rik <rik@octave.org>
parents:
20694
diff
changeset
|
362 ## When info != 1, recalculate the gradient and Hessian using the final x. |
20694
a7dbc4fc3732
fminunc.m: Correctly calculate gradient and hessian even when fminunc fails (bug #40101)
José Vallet <jose.vallet@aalto.fi>
parents:
20165
diff
changeset
|
363 if (nargout > 4 && (info == -1 || info == 2 || info == 3)) |
20715
5b7643257978
Remove trailing whitespace at end of lines.
Rik <rik@octave.org>
parents:
20695
diff
changeset
|
364 grad0 = grad; |
20694
a7dbc4fc3732
fminunc.m: Correctly calculate gradient and hessian even when fminunc fails (bug #40101)
José Vallet <jose.vallet@aalto.fi>
parents:
20165
diff
changeset
|
365 if (has_grad) |
20695
6faaab833605
fminunc.m: Clean up function to meet Octave coding standards.
Rik <rik@octave.org>
parents:
20694
diff
changeset
|
366 [fval, grad] = fcn (reshape (x, xsz)); |
20694
a7dbc4fc3732
fminunc.m: Correctly calculate gradient and hessian even when fminunc fails (bug #40101)
José Vallet <jose.vallet@aalto.fi>
parents:
20165
diff
changeset
|
367 grad = grad(:); |
a7dbc4fc3732
fminunc.m: Correctly calculate gradient and hessian even when fminunc fails (bug #40101)
José Vallet <jose.vallet@aalto.fi>
parents:
20165
diff
changeset
|
368 else |
20695
6faaab833605
fminunc.m: Clean up function to meet Octave coding standards.
Rik <rik@octave.org>
parents:
20694
diff
changeset
|
369 grad = __fdjac__ (fcn, reshape (x, xsz), fval, typicalx, cdif)(:); |
20694
a7dbc4fc3732
fminunc.m: Correctly calculate gradient and hessian even when fminunc fails (bug #40101)
José Vallet <jose.vallet@aalto.fi>
parents:
20165
diff
changeset
|
370 endif |
9084
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
371 |
20694
a7dbc4fc3732
fminunc.m: Correctly calculate gradient and hessian even when fminunc fails (bug #40101)
José Vallet <jose.vallet@aalto.fi>
parents:
20165
diff
changeset
|
372 if (nargout > 5) |
a7dbc4fc3732
fminunc.m: Correctly calculate gradient and hessian even when fminunc fails (bug #40101)
José Vallet <jose.vallet@aalto.fi>
parents:
20165
diff
changeset
|
373 ## Use the damped BFGS formula. |
a7dbc4fc3732
fminunc.m: Correctly calculate gradient and hessian even when fminunc fails (bug #40101)
José Vallet <jose.vallet@aalto.fi>
parents:
20165
diff
changeset
|
374 y = grad - grad0; |
a7dbc4fc3732
fminunc.m: Correctly calculate gradient and hessian even when fminunc fails (bug #40101)
José Vallet <jose.vallet@aalto.fi>
parents:
20165
diff
changeset
|
375 sBs = sumsq (w); |
a7dbc4fc3732
fminunc.m: Correctly calculate gradient and hessian even when fminunc fails (bug #40101)
José Vallet <jose.vallet@aalto.fi>
parents:
20165
diff
changeset
|
376 Bs = hesr' * w; |
a7dbc4fc3732
fminunc.m: Correctly calculate gradient and hessian even when fminunc fails (bug #40101)
José Vallet <jose.vallet@aalto.fi>
parents:
20165
diff
changeset
|
377 sy = y' * s; |
a7dbc4fc3732
fminunc.m: Correctly calculate gradient and hessian even when fminunc fails (bug #40101)
José Vallet <jose.vallet@aalto.fi>
parents:
20165
diff
changeset
|
378 theta = 0.8 / max (1 - sy / sBs, 0.8); |
a7dbc4fc3732
fminunc.m: Correctly calculate gradient and hessian even when fminunc fails (bug #40101)
José Vallet <jose.vallet@aalto.fi>
parents:
20165
diff
changeset
|
379 r = theta * y + (1-theta) * Bs; |
a7dbc4fc3732
fminunc.m: Correctly calculate gradient and hessian even when fminunc fails (bug #40101)
José Vallet <jose.vallet@aalto.fi>
parents:
20165
diff
changeset
|
380 hesr = cholupdate (hesr, r / sqrt (s' * r), "+"); |
a7dbc4fc3732
fminunc.m: Correctly calculate gradient and hessian even when fminunc fails (bug #40101)
José Vallet <jose.vallet@aalto.fi>
parents:
20165
diff
changeset
|
381 hesr = cholupdate (hesr, Bs / sqrt (sBs), "-"); |
a7dbc4fc3732
fminunc.m: Correctly calculate gradient and hessian even when fminunc fails (bug #40101)
José Vallet <jose.vallet@aalto.fi>
parents:
20165
diff
changeset
|
382 endif |
a7dbc4fc3732
fminunc.m: Correctly calculate gradient and hessian even when fminunc fails (bug #40101)
José Vallet <jose.vallet@aalto.fi>
parents:
20165
diff
changeset
|
383 ## Return the gradient in the same shape as x |
20695
6faaab833605
fminunc.m: Clean up function to meet Octave coding standards.
Rik <rik@octave.org>
parents:
20694
diff
changeset
|
384 grad = reshape (grad, xsz); |
20694
a7dbc4fc3732
fminunc.m: Correctly calculate gradient and hessian even when fminunc fails (bug #40101)
José Vallet <jose.vallet@aalto.fi>
parents:
20165
diff
changeset
|
385 endif |
a7dbc4fc3732
fminunc.m: Correctly calculate gradient and hessian even when fminunc fails (bug #40101)
José Vallet <jose.vallet@aalto.fi>
parents:
20165
diff
changeset
|
386 |
9084
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
387 ## Restore original shapes. |
20695
6faaab833605
fminunc.m: Clean up function to meet Octave coding standards.
Rik <rik@octave.org>
parents:
20694
diff
changeset
|
388 x = reshape (x, xsz); |
9084
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
389 |
20695
6faaab833605
fminunc.m: Clean up function to meet Octave coding standards.
Rik <rik@octave.org>
parents:
20694
diff
changeset
|
390 if (nargout > 3) |
6faaab833605
fminunc.m: Clean up function to meet Octave coding standards.
Rik <rik@octave.org>
parents:
20694
diff
changeset
|
391 output.iterations = niter; |
6faaab833605
fminunc.m: Clean up function to meet Octave coding standards.
Rik <rik@octave.org>
parents:
20694
diff
changeset
|
392 output.successful = nsuciter; |
6faaab833605
fminunc.m: Clean up function to meet Octave coding standards.
Rik <rik@octave.org>
parents:
20694
diff
changeset
|
393 output.funcCount = nfev; |
6faaab833605
fminunc.m: Clean up function to meet Octave coding standards.
Rik <rik@octave.org>
parents:
20694
diff
changeset
|
394 endif |
9084
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
395 |
9627
5bcfa0b346e8
fix extra outputs in fminunc
Jaroslav Hajek <highegg@gmail.com>
parents:
9623
diff
changeset
|
396 if (nargout > 5) |
5bcfa0b346e8
fix extra outputs in fminunc
Jaroslav Hajek <highegg@gmail.com>
parents:
9623
diff
changeset
|
397 hess = hesr'*hesr; |
5bcfa0b346e8
fix extra outputs in fminunc
Jaroslav Hajek <highegg@gmail.com>
parents:
9623
diff
changeset
|
398 endif |
5bcfa0b346e8
fix extra outputs in fminunc
Jaroslav Hajek <highegg@gmail.com>
parents:
9623
diff
changeset
|
399 |
9084
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
400 endfunction |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
401 |
18609 | 402 ## A helper function that evaluates a function and checks for bad results. |
9084
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
403 function [fx, gx] = guarded_eval (fun, x) |
21758
ffad2baa90f7
maint: Use newlines to make code more readable.
Rik <rik@octave.org>
parents:
21751
diff
changeset
|
404 |
9084
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
405 if (nargout > 1) |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
406 [fx, gx] = fun (x); |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
407 else |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
408 fx = fun (x); |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
409 gx = []; |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
410 endif |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
411 |
14383
07c55bceca23
Fix guarded_eval() subfunction in fminunc (bug #35534).
Olaf Till <olaf.till@uni-jena.de>
parents:
14138
diff
changeset
|
412 if (! (isreal (fx) && isreal (gx))) |
11587
c792872f8942
all script files: untabify and strip trailing whitespace
John W. Eaton <jwe@octave.org>
parents:
11523
diff
changeset
|
413 error ("fminunc:notreal", "fminunc: non-real value encountered"); |
9084
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
414 elseif (any (isnan (fx(:)))) |
11587
c792872f8942
all script files: untabify and strip trailing whitespace
John W. Eaton <jwe@octave.org>
parents:
11523
diff
changeset
|
415 error ("fminunc:isnan", "fminunc: NaN value encountered"); |
14386
59aab666f2bf
Extend "FunValCheck" option to optimization routines to detect Inf values.
Rik <octave@nomad.inbox5.com>
parents:
14383
diff
changeset
|
416 elseif (any (isinf (fx(:)))) |
59aab666f2bf
Extend "FunValCheck" option to optimization routines to detect Inf values.
Rik <octave@nomad.inbox5.com>
parents:
14383
diff
changeset
|
417 error ("fminunc:isinf", "fminunc: Inf value encountered"); |
9084
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
418 endif |
21758
ffad2baa90f7
maint: Use newlines to make code more readable.
Rik <rik@octave.org>
parents:
21751
diff
changeset
|
419 |
9084
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
420 endfunction |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
421 |
14386
59aab666f2bf
Extend "FunValCheck" option to optimization routines to detect Inf values.
Rik <octave@nomad.inbox5.com>
parents:
14383
diff
changeset
|
422 |
20695
6faaab833605
fminunc.m: Clean up function to meet Octave coding standards.
Rik <rik@octave.org>
parents:
20694
diff
changeset
|
423 %!function f = __rosenb__ (x) |
9084
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
424 %! n = length (x); |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
425 %! f = sumsq (1 - x(1:n-1)) + 100 * sumsq (x(2:n) - x(1:n-1).^2); |
13305
63463570d9fe
Add %!endfunction block keyword to test.m
Rik <octave@nomad.inbox5.com>
parents:
13027
diff
changeset
|
426 %!endfunction |
19833
9fc020886ae9
maint: Clean up m-files to follow Octave coding conventions.
Rik <rik@octave.org>
parents:
19697
diff
changeset
|
427 %! |
9084
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
428 %!test |
20695
6faaab833605
fminunc.m: Clean up function to meet Octave coding standards.
Rik <rik@octave.org>
parents:
20694
diff
changeset
|
429 %! [x, fval, info, out] = fminunc (@__rosenb__, [5, -5]); |
9084
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
430 %! tol = 2e-5; |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
431 %! assert (info > 0); |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
432 %! assert (x, ones (1, 2), tol); |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
433 %! assert (fval, 0, tol); |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
434 %!test |
20695
6faaab833605
fminunc.m: Clean up function to meet Octave coding standards.
Rik <rik@octave.org>
parents:
20694
diff
changeset
|
435 %! [x, fval, info, out] = fminunc (@__rosenb__, zeros (1, 4)); |
9084
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
436 %! tol = 2e-5; |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
437 %! assert (info > 0); |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
438 %! assert (x, ones (1, 4), tol); |
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
439 %! assert (fval, 0, tol); |
19833
9fc020886ae9
maint: Clean up m-files to follow Octave coding conventions.
Rik <rik@octave.org>
parents:
19697
diff
changeset
|
440 |
9fc020886ae9
maint: Clean up m-files to follow Octave coding conventions.
Rik <rik@octave.org>
parents:
19697
diff
changeset
|
441 ## Test FunValCheck works correctly |
26137
1ae11ca7dceb
fminunc.m: Change algorithm defaults to match Matlab.
Rik <rik@octave.org>
parents:
26118
diff
changeset
|
442 %!assert (fminunc (@(x) x^2, 1, optimset ("FunValCheck", "on")), 0, 1e-6) |
14386
59aab666f2bf
Extend "FunValCheck" option to optimization routines to detect Inf values.
Rik <octave@nomad.inbox5.com>
parents:
14383
diff
changeset
|
443 %!error <non-real value> fminunc (@(x) x + i, 1, optimset ("FunValCheck", "on")) |
59aab666f2bf
Extend "FunValCheck" option to optimization routines to detect Inf values.
Rik <octave@nomad.inbox5.com>
parents:
14383
diff
changeset
|
444 %!error <NaN value> fminunc (@(x) x + NaN, 1, optimset ("FunValCheck", "on")) |
59aab666f2bf
Extend "FunValCheck" option to optimization routines to detect Inf values.
Rik <octave@nomad.inbox5.com>
parents:
14383
diff
changeset
|
445 %!error <Inf value> fminunc (@(x) x + Inf, 1, optimset ("FunValCheck", "on")) |
59aab666f2bf
Extend "FunValCheck" option to optimization routines to detect Inf values.
Rik <octave@nomad.inbox5.com>
parents:
14383
diff
changeset
|
446 |
9084
b7210faa3ed0
implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
447 |
9899
9f25290a35e8
more private function and subfunction changes
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
448 ## Solve the double dogleg trust-region minimization problem: |
9f25290a35e8
more private function and subfunction changes
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
449 ## Minimize 1/2*norm(r*x)^2 subject to the constraint norm(d.*x) <= delta, |
9f25290a35e8
more private function and subfunction changes
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
450 ## x being a convex combination of the gauss-newton and scaled gradient. |
9f25290a35e8
more private function and subfunction changes
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
451 |
21578
683a1beee538
maint: Use "FIXME:" for all code blocks needing further attention.
Rik <rik@octave.org>
parents:
21546
diff
changeset
|
452 ## FIXME: error checks |
683a1beee538
maint: Use "FIXME:" for all code blocks needing further attention.
Rik <rik@octave.org>
parents:
21546
diff
changeset
|
453 ## FIXME: handle singularity, or leave it up to mldivide? |
9899
9f25290a35e8
more private function and subfunction changes
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
454 |
9f25290a35e8
more private function and subfunction changes
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
455 function x = __doglegm__ (r, g, d, delta) |
21758
ffad2baa90f7
maint: Use newlines to make code more readable.
Rik <rik@octave.org>
parents:
21751
diff
changeset
|
456 |
9899
9f25290a35e8
more private function and subfunction changes
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
457 ## Get Gauss-Newton direction. |
9f25290a35e8
more private function and subfunction changes
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
458 b = r' \ g; |
9f25290a35e8
more private function and subfunction changes
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
459 x = r \ b; |
9f25290a35e8
more private function and subfunction changes
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
460 xn = norm (d .* x); |
9f25290a35e8
more private function and subfunction changes
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
461 if (xn > delta) |
9f25290a35e8
more private function and subfunction changes
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
462 ## GN is too big, get scaled gradient. |
9f25290a35e8
more private function and subfunction changes
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
463 s = g ./ d; |
9f25290a35e8
more private function and subfunction changes
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
464 sn = norm (s); |
9f25290a35e8
more private function and subfunction changes
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
465 if (sn > 0) |
9f25290a35e8
more private function and subfunction changes
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
466 ## Normalize and rescale. |
9f25290a35e8
more private function and subfunction changes
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
467 s = (s / sn) ./ d; |
9f25290a35e8
more private function and subfunction changes
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
468 ## Get the line minimizer in s direction. |
9f25290a35e8
more private function and subfunction changes
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
469 tn = norm (r*s); |
9f25290a35e8
more private function and subfunction changes
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
470 snm = (sn / tn) / tn; |
9f25290a35e8
more private function and subfunction changes
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
471 if (snm < delta) |
10549 | 472 ## Get the dogleg path minimizer. |
9899
9f25290a35e8
more private function and subfunction changes
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
473 bn = norm (b); |
9f25290a35e8
more private function and subfunction changes
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
474 dxn = delta/xn; snmd = snm/delta; |
9f25290a35e8
more private function and subfunction changes
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
475 t = (bn/sn) * (bn/xn) * snmd; |
9f25290a35e8
more private function and subfunction changes
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
476 t -= dxn * snmd^2 - sqrt ((t-dxn)^2 + (1-dxn^2)*(1-snmd^2)); |
9f25290a35e8
more private function and subfunction changes
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
477 alpha = dxn*(1-snmd^2) / t; |
9f25290a35e8
more private function and subfunction changes
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
478 else |
9f25290a35e8
more private function and subfunction changes
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
479 alpha = 0; |
9f25290a35e8
more private function and subfunction changes
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
480 endif |
9f25290a35e8
more private function and subfunction changes
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
481 else |
9f25290a35e8
more private function and subfunction changes
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
482 alpha = delta / xn; |
9f25290a35e8
more private function and subfunction changes
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
483 snm = 0; |
9f25290a35e8
more private function and subfunction changes
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
484 endif |
9f25290a35e8
more private function and subfunction changes
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
485 ## Form the appropriate convex combination. |
9f25290a35e8
more private function and subfunction changes
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
486 x = alpha * x + ((1-alpha) * min (snm, delta)) * s; |
9f25290a35e8
more private function and subfunction changes
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
487 endif |
21758
ffad2baa90f7
maint: Use newlines to make code more readable.
Rik <rik@octave.org>
parents:
21751
diff
changeset
|
488 |
9899
9f25290a35e8
more private function and subfunction changes
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
489 endfunction |