annotate scripts/optimization/fminunc.m @ 27918:b442ec6dda5c

use centralized file for copyright info for individual contributors * COPYRIGHT.md: New file. * In most other files, use "Copyright (C) YYYY-YYYY The Octave Project Developers" instead of tracking individual names in separate source files. The motivation is to reduce the effort required to update the notices each year. Until now, the Octave source files contained copyright notices that list individual contributors. I adopted these file-scope copyright notices because that is what everyone was doing 30 years ago in the days before distributed version control systems. But now, with many contributors and modern version control systems, having these file-scope copyright notices causes trouble when we update copyright years or refactor code. Over time, the file-scope copyright notices may become outdated as new contributions are made or code is moved from one file to another. Sometimes people contribute significant patches but do not add a line claiming copyright. Other times, people add a copyright notice for their contribution but then a later refactoring moves part or all of their contribution to another file and the notice is not moved with the code. As a practical matter, moving such notices is difficult -- determining what parts are due to a particular contributor requires a time-consuming search through the project history. Even managing the yearly update of copyright years is problematic. We have some contributors who are no longer living. Should we update the copyright dates for their contributions when we release new versions? Probably not, but we do still want to claim copyright for the project as a whole. To minimize the difficulty of maintaining the copyright notices, I would like to change Octave's sources to use what is described here: https://softwarefreedom.org/resources/2012/ManagingCopyrightInformation.html in the section "Maintaining centralized copyright notices": The centralized notice approach consolidates all copyright notices in a single location, usually a top-level file. This file should contain all of the copyright notices provided project contributors, unless the contribution was clearly insignificant. It may also credit -- without a copyright notice -- anyone who helped with the project but did not contribute code or other copyrighted material. This approach captures less information about contributions within individual files, recognizing that the DVCS is better equipped to record those details. As we mentioned before, it does have one disadvantage as compared to the file-scope approach: if a single file is separated from the distribution, the recipient won't see the contributors' copyright notices. But this can be easily remedied by including a single copyright notice in each file's header, pointing to the top-level file: Copyright YYYY-YYYY The Octave Project Developers See the COPYRIGHT file at the top-level directory of this distribution or at https://octave.org/COPYRIGHT.html. followed by the usual GPL copyright statement. For more background, see the discussion here: https://lists.gnu.org/archive/html/octave-maintainers/2020-01/msg00009.html Most files in the following directories have been skipped intentinally in this changeset: doc libgui/qterminal liboctave/external m4
author John W. Eaton <jwe@octave.org>
date Mon, 06 Jan 2020 15:38:17 -0500
parents b0626d075ad7
children 1891570abac8
Ignore whitespace changes - Everywhere: Within whitespace: At end of lines:
rev   line source
27918
b442ec6dda5c use centralized file for copyright info for individual contributors
John W. Eaton <jwe@octave.org>
parents: 27586
diff changeset
1 ## Copyright (C) 2008-2019 The Octave Project Developers
b442ec6dda5c use centralized file for copyright info for individual contributors
John W. Eaton <jwe@octave.org>
parents: 27586
diff changeset
2 ##
b442ec6dda5c use centralized file for copyright info for individual contributors
John W. Eaton <jwe@octave.org>
parents: 27586
diff changeset
3 ## See the file COPYRIGHT.md in the top-level directory of this distribution
b442ec6dda5c use centralized file for copyright info for individual contributors
John W. Eaton <jwe@octave.org>
parents: 27586
diff changeset
4 ## or <https://octave.org/COPYRIGHT.html/>.
b442ec6dda5c use centralized file for copyright info for individual contributors
John W. Eaton <jwe@octave.org>
parents: 27586
diff changeset
5 ##
9084
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
6 ##
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
7 ## This file is part of Octave.
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
8 ##
24534
194eb4bd202b maint: Update punctuation for GPL v3 license text.
Rik <rik@octave.org>
parents: 23220
diff changeset
9 ## 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
10 ## 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
11 ## 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
12 ## (at your option) any later version.
9084
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
13 ##
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
14 ## 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
15 ## WITHOUT ANY WARRANTY; without even the implied warranty of
22755
3a2b891d0b33 maint: Standardize Copyright formatting.
Rik <rik@octave.org>
parents: 22323
diff changeset
16 ## MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
3a2b891d0b33 maint: Standardize Copyright formatting.
Rik <rik@octave.org>
parents: 22323
diff changeset
17 ## GNU General Public License for more details.
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 ## 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
20 ## 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
21 ## <https://www.gnu.org/licenses/>.
9084
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
22 ##
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
23 ## Author: Jaroslav Hajek <highegg@gmail.com>
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
24
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
25 ## -*- texinfo -*-
20852
516bb87ea72e 2015 Code Sprint: remove class of function from docstring for all m-files.
Rik <rik@octave.org>
parents: 20735
diff changeset
26 ## @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
27 ## @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
28 ## @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
29 ## 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
30 ## @var{fcn}.
14895
e0525ecf156e Add new function fminsearch.m
Andy Adler <andy@analyti.ca>
parents: 14868
diff changeset
31 ##
18609
923614060f1d fminunc.m: Improve documentation.
Rik <rik@octave.org>
parents: 18608
diff changeset
32 ## @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
33 ## @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
34 ##
27390
baeed03c3766 doc: Use common verbiage to describe input FUN in optimization functions.
Rik <rik@octave.org>
parents: 27069
diff changeset
35 ## @var{fun} is a function handle, inline function, or string containing the
baeed03c3766 doc: Use common verbiage to describe input FUN in optimization functions.
Rik <rik@octave.org>
parents: 27069
diff changeset
36 ## name of the function to evaluate. @var{fcn} should accept a vector (array)
baeed03c3766 doc: Use common verbiage to describe input FUN in optimization functions.
Rik <rik@octave.org>
parents: 27069
diff changeset
37 ## defining the unknown variables, and return the objective function value,
baeed03c3766 doc: Use common verbiage to describe input FUN in optimization functions.
Rik <rik@octave.org>
parents: 27069
diff changeset
38 ## optionally with gradient.
baeed03c3766 doc: Use common verbiage to describe input FUN in optimization functions.
Rik <rik@octave.org>
parents: 27069
diff changeset
39 ##
20165
f1d0f506ee78 doc: Update more docstrings to have one sentence summary as first line.
Rik <rik@octave.org>
parents: 19833
diff changeset
40 ## @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
41 ## 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
42 ##
26137
1ae11ca7dceb fminunc.m: Change algorithm defaults to match Matlab.
Rik <rik@octave.org>
parents: 26118
diff changeset
43 ## @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
44 ## 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
45 ## @qcode{"AutoScaling"}, @qcode{"FinDiffType"}, @qcode{"FunValCheck"},
1ae11ca7dceb fminunc.m: Change algorithm defaults to match Matlab.
Rik <rik@octave.org>
parents: 26118
diff changeset
46 ## @qcode{"GradObj"}, @qcode{"MaxFunEvals"}, @qcode{"MaxIter"},
1ae11ca7dceb fminunc.m: Change algorithm defaults to match Matlab.
Rik <rik@octave.org>
parents: 26118
diff changeset
47 ## @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
48 ##
1ae11ca7dceb fminunc.m: Change algorithm defaults to match Matlab.
Rik <rik@octave.org>
parents: 26118
diff changeset
49 ## 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
50 ## 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
51 ## 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
52 ## 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
53 ## unexpected (though mathematically correct) results.
9084
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
54 ##
27585
757a7119e319 doc: Create en-dashes and em-dashes correctly in documentation.
Rik <rik@octave.org>
parents: 26376
diff changeset
55 ## If @qcode{"GradObj"} is @qcode{"on"}, it specifies that @var{fcn}---when
26137
1ae11ca7dceb fminunc.m: Change algorithm defaults to match Matlab.
Rik <rik@octave.org>
parents: 26118
diff changeset
56 ## 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
57 ## partial first derivatives at the requested point.
1ae11ca7dceb fminunc.m: Change algorithm defaults to match Matlab.
Rik <rik@octave.org>
parents: 26118
diff changeset
58 ##
1ae11ca7dceb fminunc.m: Change algorithm defaults to match Matlab.
Rik <rik@octave.org>
parents: 26118
diff changeset
59 ## @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
60 ## 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
61 ## @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
62 ## The value must be a positive integer.
1ae11ca7dceb fminunc.m: Change algorithm defaults to match Matlab.
Rik <rik@octave.org>
parents: 26118
diff changeset
63 ##
1ae11ca7dceb fminunc.m: Change algorithm defaults to match Matlab.
Rik <rik@octave.org>
parents: 26118
diff changeset
64 ## @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
65 ## 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
66 ## The value must be a positive integer.
1ae11ca7dceb fminunc.m: Change algorithm defaults to match Matlab.
Rik <rik@octave.org>
parents: 26118
diff changeset
67 ##
1ae11ca7dceb fminunc.m: Change algorithm defaults to match Matlab.
Rik <rik@octave.org>
parents: 26118
diff changeset
68 ## @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
69 ## @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
70 ## 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
71 ##
18609
923614060f1d fminunc.m: Improve documentation.
Rik <rik@octave.org>
parents: 18608
diff changeset
72 ## 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
73 ##
18609
923614060f1d fminunc.m: Improve documentation.
Rik <rik@octave.org>
parents: 18608
diff changeset
74 ## 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
75 ## 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
76 ##
f1d0f506ee78 doc: Update more docstrings to have one sentence summary as first line.
Rik <rik@octave.org>
parents: 19833
diff changeset
77 ## @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
78 ##
9084
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
79 ## @table @asis
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
80 ## @item 1
10793
be55736a0783 Grammarcheck the documentation from m-files.
Rik <octave@nomad.inbox5.com>
parents: 10687
diff changeset
81 ## Converged to a solution point. Relative gradient error is less than
18609
923614060f1d fminunc.m: Improve documentation.
Rik <rik@octave.org>
parents: 18608
diff changeset
82 ## specified by @code{TolFun}.
10821
693e22af08ae Grammarcheck documentation of m-files
Rik <octave@nomad.inbox5.com>
parents: 10793
diff changeset
83 ##
9084
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
84 ## @item 2
18609
923614060f1d fminunc.m: Improve documentation.
Rik <rik@octave.org>
parents: 18608
diff changeset
85 ## 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
86 ##
9084
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
87 ## @item 3
18609
923614060f1d fminunc.m: Improve documentation.
Rik <rik@octave.org>
parents: 18608
diff changeset
88 ## 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
89 ##
9084
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
90 ## @item 0
18816
e275d15c27b5 doc: Periodic spellcheck of documentation.
Rik <rik@octave.org>
parents: 18609
diff changeset
91 ## Iteration limit exceeded---either maximum number of algorithm iterations
18609
923614060f1d fminunc.m: Improve documentation.
Rik <rik@octave.org>
parents: 18608
diff changeset
92 ## @code{MaxIter} or maximum number of function evaluations @code{MaxFunEvals}.
923614060f1d fminunc.m: Improve documentation.
Rik <rik@octave.org>
parents: 18608
diff changeset
93 ##
923614060f1d fminunc.m: Improve documentation.
Rik <rik@octave.org>
parents: 18608
diff changeset
94 ## @item -1
18816
e275d15c27b5 doc: Periodic spellcheck of documentation.
Rik <rik@octave.org>
parents: 18609
diff changeset
95 ## Algorithm terminated by @code{OutputFcn}.
10821
693e22af08ae Grammarcheck documentation of m-files
Rik <octave@nomad.inbox5.com>
parents: 10793
diff changeset
96 ##
9084
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
97 ## @item -3
11587
c792872f8942 all script files: untabify and strip trailing whitespace
John W. Eaton <jwe@octave.org>
parents: 11523
diff changeset
98 ## The trust region radius became excessively small.
9084
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
99 ## @end table
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
100 ##
21546
f7f97d7e9294 doc: Wrap m-file docstrings to 79 characters + newline (80 total).
Rik <rik@octave.org>
parents: 21178
diff changeset
101 ## 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
102 ## 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
103 ## 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
104 ## @var{x}.
9627
5bcfa0b346e8 fix extra outputs in fminunc
Jaroslav Hajek <highegg@gmail.com>
parents: 9623
diff changeset
105 ##
20695
6faaab833605 fminunc.m: Clean up function to meet Octave coding standards.
Rik <rik@octave.org>
parents: 20694
diff changeset
106 ## 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
107 ## 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
108 ##
22092
d3cef63f79ac * fminunc.m: Fix incorrect reference to fminunc as fminsearch (bug #48479)
Mike Miller <mtmiller@octave.org>
parents: 21758
diff changeset
109 ## 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
110 ## 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
111 ## 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
112 ## @code{fminsearch}.
14895
e0525ecf156e Add new function fminsearch.m
Andy Adler <andy@analyti.ca>
parents: 14868
diff changeset
113 ## @seealso{fminbnd, fminsearch, optimset}
9084
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
114 ## @end deftypefn
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
115
13027
b9a89ca0fb75 prevent optimization functions from setting ans in workspace at startup
John W. Eaton <jwe@octave.org>
parents: 12578
diff changeset
116 ## 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
117 ## PKG_ADD: [~] = __all_opts__ ("fminunc");
9084
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
118
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
119 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
120
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
121 ## 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
122 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
123 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
124 "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
125 "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
126 "TolFun", 1e-6, "TolX", 1e-6, "TypicalX", []);
9084
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
127 return;
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
128 endif
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
129
21171
2935d56203a4 Fix regressions caused by ismatrix definition change (partial fix bug #47036).
Rik <rik@octave.org>
parents: 20165
diff changeset
130 if (nargin < 2 || nargin > 3 || ! isnumeric (x0))
9084
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
131 print_usage ();
11587
c792872f8942 all script files: untabify and strip trailing whitespace
John W. Eaton <jwe@octave.org>
parents: 11523
diff changeset
132 endif
9084
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 if (ischar (fcn))
26118
7502fce4cd3a str2func: eliminate optional second "global" argument
John W. Eaton <jwe@octave.org>
parents: 25054
diff changeset
135 fcn = str2func (fcn);
9084
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
136 endif
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
137
20695
6faaab833605 fminunc.m: Clean up function to meet Octave coding standards.
Rik <rik@octave.org>
parents: 20694
diff changeset
138 xsz = size (x0);
9084
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
139 n = numel (x0);
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
140
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
141 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
142 cdif = strcmpi (optimget (options, "FinDiffType", "forward"), "central");
9084
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
143 maxiter = optimget (options, "MaxIter", 400);
26137
1ae11ca7dceb fminunc.m: Change algorithm defaults to match Matlab.
Rik <rik@octave.org>
parents: 26118
diff changeset
144 maxfev = optimget (options, "MaxFunEvals", 100*n);
9084
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
145 outfcn = optimget (options, "OutputFcn");
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
146
20695
6faaab833605 fminunc.m: Clean up function to meet Octave coding standards.
Rik <rik@octave.org>
parents: 20694
diff changeset
147 ## 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
148 ## 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
149 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
150 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
151 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
152 endif
5c66978f3fdf support TypicalX and AutoScaling in fsolve/fminunc, don't autoscale by default
Jaroslav Hajek <highegg@gmail.com>
parents: 9899
diff changeset
153 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
154 if (! autoscale)
5c66978f3fdf support TypicalX and AutoScaling in fsolve/fminunc, don't autoscale by default
Jaroslav Hajek <highegg@gmail.com>
parents: 9899
diff changeset
155 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
156 endif
5c66978f3fdf support TypicalX and AutoScaling in fsolve/fminunc, don't autoscale by default
Jaroslav Hajek <highegg@gmail.com>
parents: 9899
diff changeset
157
9084
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
158 funvalchk = strcmpi (optimget (options, "FunValCheck", "off"), "on");
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
159
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
160 if (funvalchk)
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
161 ## Replace fcn with a guarded version.
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
162 fcn = @(x) guarded_eval (fcn, x);
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
163 endif
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
164
20695
6faaab833605 fminunc.m: Clean up function to meet Octave coding standards.
Rik <rik@octave.org>
parents: 20694
diff changeset
165 ## 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
166 ## prefers accuracy to performance.
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
167
10201
5c66978f3fdf support TypicalX and AutoScaling in fsolve/fminunc, don't autoscale by default
Jaroslav Hajek <highegg@gmail.com>
parents: 9899
diff changeset
168 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
169 tolf = optimget (options, "TolFun", 1e-7);
9084
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
170
9623
bc0739d02724 update initial TR step for fsolve and fminunc
Jaroslav Hajek <highegg@gmail.com>
parents: 9464
diff changeset
171 factor = 0.1;
9084
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
172 ## FIXME: TypicalX corresponds to user scaling (???)
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
173 autodg = true;
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
174
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
175 niter = 1;
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
176 nfev = 0;
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 x = x0(:);
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
179 info = 0;
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
180
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
181 ## Initial evaluation.
20695
6faaab833605 fminunc.m: Clean up function to meet Octave coding standards.
Rik <rik@octave.org>
parents: 20694
diff changeset
182 fval = fcn (reshape (x, xsz));
9084
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
183 n = length (x);
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
184
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
185 if (! isempty (outfcn))
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
186 optimvalues.iter = niter;
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
187 optimvalues.funccount = nfev;
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
188 optimvalues.fval = fval;
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
189 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
190 state = "init";
9084
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
191 stop = outfcn (x, optimvalues, state);
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
192 if (stop)
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
193 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
194 return;
9084
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
195 endif
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
196 endif
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
197
21099
52af4092f863 For optimization scripts, correctly choose tolerance (eps) based on class of fun and X0.
Rik <rik@octave.org>
parents: 20852
diff changeset
198 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
199 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
200 else
52af4092f863 For optimization scripts, correctly choose tolerance (eps) based on class of fun and X0.
Rik <rik@octave.org>
parents: 20852
diff changeset
201 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
202 endif
52af4092f863 For optimization scripts, correctly choose tolerance (eps) based on class of fun and X0.
Rik <rik@octave.org>
parents: 20852
diff changeset
203
9084
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
204 nsuciter = 0;
9199
399884c9d4a1 import the step adaptation strategy from fsolve to fminunc
Jaroslav Hajek <highegg@gmail.com>
parents: 9084
diff changeset
205 lastratio = 0;
9084
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 grad = [];
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
208
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
209 ## Outer loop.
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
210 while (niter < maxiter && nfev < maxfev && ! info)
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
211
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
212 grad0 = grad;
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
213
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
214 ## Calculate function value and gradient (possibly via FD).
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
215 if (has_grad)
20695
6faaab833605 fminunc.m: Clean up function to meet Octave coding standards.
Rik <rik@octave.org>
parents: 20694
diff changeset
216 [fval, grad] = fcn (reshape (x, xsz));
9084
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
217 grad = grad(:);
20695
6faaab833605 fminunc.m: Clean up function to meet Octave coding standards.
Rik <rik@octave.org>
parents: 20694
diff changeset
218 nfev += 1;
9084
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
219 else
20695
6faaab833605 fminunc.m: Clean up function to meet Octave coding standards.
Rik <rik@octave.org>
parents: 20694
diff changeset
220 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
221 nfev += (1 + cdif) * length (x);
9084
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
222 endif
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
223
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
224 if (niter == 1)
11587
c792872f8942 all script files: untabify and strip trailing whitespace
John W. Eaton <jwe@octave.org>
parents: 11523
diff changeset
225 ## Initialize by identity matrix.
9084
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
226 hesr = eye (n);
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
227 else
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
228 ## Use the damped BFGS formula.
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
229 y = grad - grad0;
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
230 sBs = sumsq (w);
20695
6faaab833605 fminunc.m: Clean up function to meet Octave coding standards.
Rik <rik@octave.org>
parents: 20694
diff changeset
231 Bs = hesr' * w;
6faaab833605 fminunc.m: Clean up function to meet Octave coding standards.
Rik <rik@octave.org>
parents: 20694
diff changeset
232 sy = y' * s;
9633
ecc2c556f844 simplify damped BFGS formula
Jaroslav Hajek <highegg@gmail.com>
parents: 9631
diff changeset
233 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
234 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
235 hesr = cholupdate (hesr, r / sqrt (s' * r), "+");
9084
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
236 [hesr, info] = cholupdate (hesr, Bs / sqrt (sBs), "-");
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
237 if (info)
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
238 hesr = eye (n);
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
239 endif
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
240 endif
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
241
10201
5c66978f3fdf support TypicalX and AutoScaling in fsolve/fminunc, don't autoscale by default
Jaroslav Hajek <highegg@gmail.com>
parents: 9899
diff changeset
242 if (autoscale)
20695
6faaab833605 fminunc.m: Clean up function to meet Octave coding standards.
Rik <rik@octave.org>
parents: 20694
diff changeset
243 ## Second derivatives approximate the Hessian.
6faaab833605 fminunc.m: Clean up function to meet Octave coding standards.
Rik <rik@octave.org>
parents: 20694
diff changeset
244 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
245 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
246 dg = d2f;
5c66978f3fdf support TypicalX and AutoScaling in fsolve/fminunc, don't autoscale by default
Jaroslav Hajek <highegg@gmail.com>
parents: 9899
diff changeset
247 else
5c66978f3fdf support TypicalX and AutoScaling in fsolve/fminunc, don't autoscale by default
Jaroslav Hajek <highegg@gmail.com>
parents: 9899
diff changeset
248 ## 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
249 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
250 endif
5c66978f3fdf support TypicalX and AutoScaling in fsolve/fminunc, don't autoscale by default
Jaroslav Hajek <highegg@gmail.com>
parents: 9899
diff changeset
251 endif
5c66978f3fdf support TypicalX and AutoScaling in fsolve/fminunc, don't autoscale by default
Jaroslav Hajek <highegg@gmail.com>
parents: 9899
diff changeset
252
9084
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
253 if (niter == 1)
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
254 xn = norm (dg .* x);
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
255 ## FIXME: something better?
9628
73e6ad869f08 further correct initial TR step strategy
Jaroslav Hajek <highegg@gmail.com>
parents: 9627
diff changeset
256 delta = factor * max (xn, 1);
9084
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
257 endif
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
258
18609
923614060f1d fminunc.m: Improve documentation.
Rik <rik@octave.org>
parents: 18608
diff changeset
259 ## 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
260 ## 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
261 ## 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
262 if (norm (grad) <= tolf*n*xn)
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
263 info = 1;
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
264 break;
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
265 endif
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
266
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
267 suc = false;
9207
25f50d2d76b3 improve TR updating strategy for fminunc and fsolve
Jaroslav Hajek <highegg@gmail.com>
parents: 9199
diff changeset
268 decfac = 0.5;
25f50d2d76b3 improve TR updating strategy for fminunc and fsolve
Jaroslav Hajek <highegg@gmail.com>
parents: 9199
diff changeset
269
9084
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
270 ## Inner loop.
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
271 while (! suc && niter <= maxiter && nfev < maxfev && ! info)
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
272
9631
00958d0c4e3c split __dogleg__ > __doglegm__
Jaroslav Hajek <highegg@gmail.com>
parents: 9628
diff changeset
273 s = - __doglegm__ (hesr, grad, dg, delta);
9084
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
274
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
275 sn = norm (dg .* s);
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
276 if (niter == 1)
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
277 delta = min (delta, sn);
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
278 endif
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
279
20695
6faaab833605 fminunc.m: Clean up function to meet Octave coding standards.
Rik <rik@octave.org>
parents: 20694
diff changeset
280 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
281 nfev += 1;
9084
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
282
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
283 if (fval1 < fval)
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
284 ## Scaled actual reduction.
19833
9fc020886ae9 maint: Clean up m-files to follow Octave coding conventions.
Rik <rik@octave.org>
parents: 19697
diff changeset
285 actred = (fval - fval1) / (abs (fval1) + abs (fval));
9084
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
286 else
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
287 actred = -1;
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
288 endif
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
289
20695
6faaab833605 fminunc.m: Clean up function to meet Octave coding standards.
Rik <rik@octave.org>
parents: 20694
diff changeset
290 w = hesr * s;
9084
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
291 ## Scaled predicted reduction, and ratio.
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
292 t = 1/2 * sumsq (w) + grad'*s;
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
293 if (t < 0)
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
294 prered = -t/(abs (fval) + abs (fval + t));
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
295 ratio = actred / prered;
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
296 else
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
297 prered = 0;
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
298 ratio = 0;
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
299 endif
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
300
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
301 ## Update delta.
14868
5d3a684236b0 maint: Use Octave coding conventions for cuddling parentheses in scripts directory
Rik <octave@nomad.inbox5.com>
parents: 14552
diff changeset
302 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
303 delta *= decfac;
25f50d2d76b3 improve TR updating strategy for fminunc and fsolve
Jaroslav Hajek <highegg@gmail.com>
parents: 9199
diff changeset
304 decfac ^= 1.4142;
20695
6faaab833605 fminunc.m: Clean up function to meet Octave coding standards.
Rik <rik@octave.org>
parents: 20694
diff changeset
305 if (delta <= 10*macheps*xn)
9084
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
306 ## Trust region became uselessly small.
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
307 info = -3;
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
308 break;
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
309 endif
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
310 else
9199
399884c9d4a1 import the step adaptation strategy from fsolve to fminunc
Jaroslav Hajek <highegg@gmail.com>
parents: 9084
diff changeset
311 lastratio = ratio;
9207
25f50d2d76b3 improve TR updating strategy for fminunc and fsolve
Jaroslav Hajek <highegg@gmail.com>
parents: 9199
diff changeset
312 decfac = 0.5;
9084
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
313 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
314 delta = 1.4142*sn;
9084
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
315 elseif (ratio >= 0.5)
9207
25f50d2d76b3 improve TR updating strategy for fminunc and fsolve
Jaroslav Hajek <highegg@gmail.com>
parents: 9199
diff changeset
316 delta = max (delta, 1.4142*sn);
9084
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
317 endif
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
318 endif
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
319
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
320 if (ratio >= 1e-4)
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
321 ## Successful iteration.
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
322 x += s;
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
323 xn = norm (dg .* x);
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
324 fval = fval1;
20735
418ae0cb752f Replace ++,-- with in-place operators for performance.
Rik <rik@octave.org>
parents: 20715
diff changeset
325 nsuciter += 1;
9084
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
326 suc = true;
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
327 endif
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
328
20695
6faaab833605 fminunc.m: Clean up function to meet Octave coding standards.
Rik <rik@octave.org>
parents: 20694
diff changeset
329 niter += 1;
9084
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
330
18609
923614060f1d fminunc.m: Improve documentation.
Rik <rik@octave.org>
parents: 18608
diff changeset
331 ## 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
332 if (! isempty (outfcn))
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
333 optimvalues.iter = niter;
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
334 optimvalues.funccount = nfev;
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
335 optimvalues.fval = fval;
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
336 optimvalues.searchdirection = s;
20695
6faaab833605 fminunc.m: Clean up function to meet Octave coding standards.
Rik <rik@octave.org>
parents: 20694
diff changeset
337 state = "iter";
9084
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
338 stop = outfcn (x, optimvalues, state);
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
339 if (stop)
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
340 info = -1;
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
341 break;
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
342 endif
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
343 endif
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
344
20695
6faaab833605 fminunc.m: Clean up function to meet Octave coding standards.
Rik <rik@octave.org>
parents: 20694
diff changeset
345 ## 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
346 ## 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
347
9084
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
348 ## 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
349 ## 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
350 ## 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
351 ## iterations, so we need scaling-independent tolerances wherever
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
352 ## possible.
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
353
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
354 ## The following tests done only after successful step.
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
355 if (ratio >= 1e-4)
20695
6faaab833605 fminunc.m: Clean up function to meet Octave coding standards.
Rik <rik@octave.org>
parents: 20694
diff changeset
356 ## 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
357 ## but compare to scaled step, so nothing bad.
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
358 if (sn <= tolx*xn)
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
359 info = 2;
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
360 ## Again a classic one.
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
361 elseif (actred < tolf)
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
362 info = 3;
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
363 endif
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
364 endif
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
365
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
366 endwhile
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
367 endwhile
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
368
20695
6faaab833605 fminunc.m: Clean up function to meet Octave coding standards.
Rik <rik@octave.org>
parents: 20694
diff changeset
369 ## 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
370 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
371 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
372 if (has_grad)
20695
6faaab833605 fminunc.m: Clean up function to meet Octave coding standards.
Rik <rik@octave.org>
parents: 20694
diff changeset
373 [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
374 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
375 else
20695
6faaab833605 fminunc.m: Clean up function to meet Octave coding standards.
Rik <rik@octave.org>
parents: 20694
diff changeset
376 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
377 endif
9084
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
378
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
379 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
380 ## 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
381 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
382 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
383 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
384 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
385 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
386 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
387 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
388 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
389 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
390 ## 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
391 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
392 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
393
9084
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
394 ## Restore original shapes.
20695
6faaab833605 fminunc.m: Clean up function to meet Octave coding standards.
Rik <rik@octave.org>
parents: 20694
diff changeset
395 x = reshape (x, xsz);
9084
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
396
20695
6faaab833605 fminunc.m: Clean up function to meet Octave coding standards.
Rik <rik@octave.org>
parents: 20694
diff changeset
397 if (nargout > 3)
6faaab833605 fminunc.m: Clean up function to meet Octave coding standards.
Rik <rik@octave.org>
parents: 20694
diff changeset
398 output.iterations = niter;
6faaab833605 fminunc.m: Clean up function to meet Octave coding standards.
Rik <rik@octave.org>
parents: 20694
diff changeset
399 output.successful = nsuciter;
6faaab833605 fminunc.m: Clean up function to meet Octave coding standards.
Rik <rik@octave.org>
parents: 20694
diff changeset
400 output.funcCount = nfev;
6faaab833605 fminunc.m: Clean up function to meet Octave coding standards.
Rik <rik@octave.org>
parents: 20694
diff changeset
401 endif
9084
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
402
9627
5bcfa0b346e8 fix extra outputs in fminunc
Jaroslav Hajek <highegg@gmail.com>
parents: 9623
diff changeset
403 if (nargout > 5)
5bcfa0b346e8 fix extra outputs in fminunc
Jaroslav Hajek <highegg@gmail.com>
parents: 9623
diff changeset
404 hess = hesr'*hesr;
5bcfa0b346e8 fix extra outputs in fminunc
Jaroslav Hajek <highegg@gmail.com>
parents: 9623
diff changeset
405 endif
5bcfa0b346e8 fix extra outputs in fminunc
Jaroslav Hajek <highegg@gmail.com>
parents: 9623
diff changeset
406
9084
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
407 endfunction
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
408
18609
923614060f1d fminunc.m: Improve documentation.
Rik <rik@octave.org>
parents: 18608
diff changeset
409 ## 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
410 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
411
9084
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
412 if (nargout > 1)
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
413 [fx, gx] = fun (x);
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
414 else
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
415 fx = fun (x);
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
416 gx = [];
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
417 endif
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
418
14383
07c55bceca23 Fix guarded_eval() subfunction in fminunc (bug #35534).
Olaf Till <olaf.till@uni-jena.de>
parents: 14138
diff changeset
419 if (! (isreal (fx) && isreal (gx)))
27069
0a62d9a6aa2d Place Octave's warning and error IDs in to the "Octave" namespace (bug #56213).
Rik <rik@octave.org>
parents: 26376
diff changeset
420 error ("Octave:fminunc:notreal", "fminunc: non-real value encountered");
9084
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
421 elseif (any (isnan (fx(:))))
27069
0a62d9a6aa2d Place Octave's warning and error IDs in to the "Octave" namespace (bug #56213).
Rik <rik@octave.org>
parents: 26376
diff changeset
422 error ("Octave: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
423 elseif (any (isinf (fx(:))))
27069
0a62d9a6aa2d Place Octave's warning and error IDs in to the "Octave" namespace (bug #56213).
Rik <rik@octave.org>
parents: 26376
diff changeset
424 error ("Octave:fminunc:isinf", "fminunc: Inf value encountered");
9084
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
425 endif
21758
ffad2baa90f7 maint: Use newlines to make code more readable.
Rik <rik@octave.org>
parents: 21751
diff changeset
426
9084
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
427 endfunction
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
428
14386
59aab666f2bf Extend "FunValCheck" option to optimization routines to detect Inf values.
Rik <octave@nomad.inbox5.com>
parents: 14383
diff changeset
429
20695
6faaab833605 fminunc.m: Clean up function to meet Octave coding standards.
Rik <rik@octave.org>
parents: 20694
diff changeset
430 %!function f = __rosenb__ (x)
9084
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
431 %! n = length (x);
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
432 %! 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
433 %!endfunction
19833
9fc020886ae9 maint: Clean up m-files to follow Octave coding conventions.
Rik <rik@octave.org>
parents: 19697
diff changeset
434 %!
9084
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
435 %!test
20695
6faaab833605 fminunc.m: Clean up function to meet Octave coding standards.
Rik <rik@octave.org>
parents: 20694
diff changeset
436 %! [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
437 %! tol = 2e-5;
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
438 %! assert (info > 0);
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
439 %! assert (x, ones (1, 2), tol);
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
440 %! assert (fval, 0, tol);
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
441 %!test
20695
6faaab833605 fminunc.m: Clean up function to meet Octave coding standards.
Rik <rik@octave.org>
parents: 20694
diff changeset
442 %! [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
443 %! tol = 2e-5;
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
444 %! assert (info > 0);
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
445 %! assert (x, ones (1, 4), tol);
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
446 %! assert (fval, 0, tol);
19833
9fc020886ae9 maint: Clean up m-files to follow Octave coding conventions.
Rik <rik@octave.org>
parents: 19697
diff changeset
447
9fc020886ae9 maint: Clean up m-files to follow Octave coding conventions.
Rik <rik@octave.org>
parents: 19697
diff changeset
448 ## Test FunValCheck works correctly
26137
1ae11ca7dceb fminunc.m: Change algorithm defaults to match Matlab.
Rik <rik@octave.org>
parents: 26118
diff changeset
449 %!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
450 %!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
451 %!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
452 %!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
453
9084
b7210faa3ed0 implement fminunc using factored trust-region BFGS
Jaroslav Hajek <highegg@gmail.com>
parents:
diff changeset
454
9899
9f25290a35e8 more private function and subfunction changes
John W. Eaton <jwe@octave.org>
parents: 9758
diff changeset
455 ## 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
456 ## 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
457 ## 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
458
21578
683a1beee538 maint: Use "FIXME:" for all code blocks needing further attention.
Rik <rik@octave.org>
parents: 21546
diff changeset
459 ## FIXME: error checks
683a1beee538 maint: Use "FIXME:" for all code blocks needing further attention.
Rik <rik@octave.org>
parents: 21546
diff changeset
460 ## 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
461
9f25290a35e8 more private function and subfunction changes
John W. Eaton <jwe@octave.org>
parents: 9758
diff changeset
462 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
463
9899
9f25290a35e8 more private function and subfunction changes
John W. Eaton <jwe@octave.org>
parents: 9758
diff changeset
464 ## Get Gauss-Newton direction.
9f25290a35e8 more private function and subfunction changes
John W. Eaton <jwe@octave.org>
parents: 9758
diff changeset
465 b = r' \ g;
9f25290a35e8 more private function and subfunction changes
John W. Eaton <jwe@octave.org>
parents: 9758
diff changeset
466 x = r \ b;
9f25290a35e8 more private function and subfunction changes
John W. Eaton <jwe@octave.org>
parents: 9758
diff changeset
467 xn = norm (d .* x);
9f25290a35e8 more private function and subfunction changes
John W. Eaton <jwe@octave.org>
parents: 9758
diff changeset
468 if (xn > delta)
9f25290a35e8 more private function and subfunction changes
John W. Eaton <jwe@octave.org>
parents: 9758
diff changeset
469 ## GN is too big, get scaled gradient.
9f25290a35e8 more private function and subfunction changes
John W. Eaton <jwe@octave.org>
parents: 9758
diff changeset
470 s = g ./ d;
9f25290a35e8 more private function and subfunction changes
John W. Eaton <jwe@octave.org>
parents: 9758
diff changeset
471 sn = norm (s);
9f25290a35e8 more private function and subfunction changes
John W. Eaton <jwe@octave.org>
parents: 9758
diff changeset
472 if (sn > 0)
9f25290a35e8 more private function and subfunction changes
John W. Eaton <jwe@octave.org>
parents: 9758
diff changeset
473 ## Normalize and rescale.
9f25290a35e8 more private function and subfunction changes
John W. Eaton <jwe@octave.org>
parents: 9758
diff changeset
474 s = (s / sn) ./ d;
9f25290a35e8 more private function and subfunction changes
John W. Eaton <jwe@octave.org>
parents: 9758
diff changeset
475 ## Get the line minimizer in s direction.
9f25290a35e8 more private function and subfunction changes
John W. Eaton <jwe@octave.org>
parents: 9758
diff changeset
476 tn = norm (r*s);
9f25290a35e8 more private function and subfunction changes
John W. Eaton <jwe@octave.org>
parents: 9758
diff changeset
477 snm = (sn / tn) / tn;
9f25290a35e8 more private function and subfunction changes
John W. Eaton <jwe@octave.org>
parents: 9758
diff changeset
478 if (snm < delta)
10549
95c3e38098bf Untabify .m scripts
Rik <code@nomad.inbox5.com>
parents: 10201
diff changeset
479 ## Get the dogleg path minimizer.
9899
9f25290a35e8 more private function and subfunction changes
John W. Eaton <jwe@octave.org>
parents: 9758
diff changeset
480 bn = norm (b);
9f25290a35e8 more private function and subfunction changes
John W. Eaton <jwe@octave.org>
parents: 9758
diff changeset
481 dxn = delta/xn; snmd = snm/delta;
9f25290a35e8 more private function and subfunction changes
John W. Eaton <jwe@octave.org>
parents: 9758
diff changeset
482 t = (bn/sn) * (bn/xn) * snmd;
9f25290a35e8 more private function and subfunction changes
John W. Eaton <jwe@octave.org>
parents: 9758
diff changeset
483 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
484 alpha = dxn*(1-snmd^2) / t;
9f25290a35e8 more private function and subfunction changes
John W. Eaton <jwe@octave.org>
parents: 9758
diff changeset
485 else
9f25290a35e8 more private function and subfunction changes
John W. Eaton <jwe@octave.org>
parents: 9758
diff changeset
486 alpha = 0;
9f25290a35e8 more private function and subfunction changes
John W. Eaton <jwe@octave.org>
parents: 9758
diff changeset
487 endif
9f25290a35e8 more private function and subfunction changes
John W. Eaton <jwe@octave.org>
parents: 9758
diff changeset
488 else
9f25290a35e8 more private function and subfunction changes
John W. Eaton <jwe@octave.org>
parents: 9758
diff changeset
489 alpha = delta / xn;
9f25290a35e8 more private function and subfunction changes
John W. Eaton <jwe@octave.org>
parents: 9758
diff changeset
490 snm = 0;
9f25290a35e8 more private function and subfunction changes
John W. Eaton <jwe@octave.org>
parents: 9758
diff changeset
491 endif
9f25290a35e8 more private function and subfunction changes
John W. Eaton <jwe@octave.org>
parents: 9758
diff changeset
492 ## Form the appropriate convex combination.
9f25290a35e8 more private function and subfunction changes
John W. Eaton <jwe@octave.org>
parents: 9758
diff changeset
493 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
494 endif
21758
ffad2baa90f7 maint: Use newlines to make code more readable.
Rik <rik@octave.org>
parents: 21751
diff changeset
495
9899
9f25290a35e8 more private function and subfunction changes
John W. Eaton <jwe@octave.org>
parents: 9758
diff changeset
496 endfunction