Mercurial > forge
annotate main/optim/inst/nonlin_residmin.m @ 9930:d30cfca46e8a octave-forge
optim: upgrade license to GPLv3+ and mention on DESCRIPTION the other package licenses
author | carandraug |
---|---|
date | Fri, 30 Mar 2012 15:14:48 +0000 |
parents | f80846665d3b |
children | fba8cdd5f9ad |
rev | line source |
---|---|
7951
b0fc599dcdb1
Update copyright notices of functions for 'subpackage' nonlin_residmin.
i7tiol
parents:
7871
diff
changeset
|
1 ## Copyright (C) 2010, 2011 Olaf Till <olaf.till@uni-jena.de> |
7684
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
2 ## |
9930
d30cfca46e8a
optim: upgrade license to GPLv3+ and mention on DESCRIPTION the other package licenses
carandraug
parents:
8720
diff
changeset
|
3 ## This program is free software; you can redistribute it and/or modify it under |
d30cfca46e8a
optim: upgrade license to GPLv3+ and mention on DESCRIPTION the other package licenses
carandraug
parents:
8720
diff
changeset
|
4 ## the terms of the GNU General Public License as published by the Free Software |
d30cfca46e8a
optim: upgrade license to GPLv3+ and mention on DESCRIPTION the other package licenses
carandraug
parents:
8720
diff
changeset
|
5 ## Foundation; either version 3 of the License, or (at your option) any later |
d30cfca46e8a
optim: upgrade license to GPLv3+ and mention on DESCRIPTION the other package licenses
carandraug
parents:
8720
diff
changeset
|
6 ## version. |
7684
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
7 ## |
9930
d30cfca46e8a
optim: upgrade license to GPLv3+ and mention on DESCRIPTION the other package licenses
carandraug
parents:
8720
diff
changeset
|
8 ## This program is distributed in the hope that it will be useful, but WITHOUT |
d30cfca46e8a
optim: upgrade license to GPLv3+ and mention on DESCRIPTION the other package licenses
carandraug
parents:
8720
diff
changeset
|
9 ## ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or |
d30cfca46e8a
optim: upgrade license to GPLv3+ and mention on DESCRIPTION the other package licenses
carandraug
parents:
8720
diff
changeset
|
10 ## FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License for more |
d30cfca46e8a
optim: upgrade license to GPLv3+ and mention on DESCRIPTION the other package licenses
carandraug
parents:
8720
diff
changeset
|
11 ## details. |
7684
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
12 ## |
9930
d30cfca46e8a
optim: upgrade license to GPLv3+ and mention on DESCRIPTION the other package licenses
carandraug
parents:
8720
diff
changeset
|
13 ## You should have received a copy of the GNU General Public License along with |
d30cfca46e8a
optim: upgrade license to GPLv3+ and mention on DESCRIPTION the other package licenses
carandraug
parents:
8720
diff
changeset
|
14 ## this program; if not, see <http://www.gnu.org/licenses/>. |
7684
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
15 |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
16 ## -*- texinfo -*- |
7688 | 17 ## @deftypefn {Function File} {[@var{p}, @var{resid}, @var{cvg}, @var{outp}] =} nonlin_residmin (@var{f}, @var{pin}) |
18 ## @deftypefnx {Function File} {[@var{p}, @var{resid}, @var{cvg}, @var{outp}] =} nonlin_residmin (@var{f}, @var{pin}, @var{settings}) | |
9930
d30cfca46e8a
optim: upgrade license to GPLv3+ and mention on DESCRIPTION the other package licenses
carandraug
parents:
8720
diff
changeset
|
19 ## Frontend for nonlinear minimization of residuals returned by a model |
d30cfca46e8a
optim: upgrade license to GPLv3+ and mention on DESCRIPTION the other package licenses
carandraug
parents:
8720
diff
changeset
|
20 ## function. |
7684
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
21 ## |
9930
d30cfca46e8a
optim: upgrade license to GPLv3+ and mention on DESCRIPTION the other package licenses
carandraug
parents:
8720
diff
changeset
|
22 ## The functions supplied by the user have a minimal |
7684
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
23 ## interface; any additionally needed constants (e.g. observed values) |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
24 ## can be supplied by wrapping the user functions into anonymous |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
25 ## functions. |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
26 ## |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
27 ## The following description applies to usage with vector-based |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
28 ## parameter handling. Differences in usage for structure-based |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
29 ## parameter handling will be explained in a separate section below. |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
30 ## |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
31 ## @var{f}: function returning the array of residuals. It gets a column |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
32 ## vector of real parameters as argument. In gradient determination, |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
33 ## this function may be called with an informational second argument, |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
34 ## whose content depends on the function for gradient determination. |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
35 ## |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
36 ## @var{pin}: real column vector of initial parameters. |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
37 ## |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
38 ## @var{settings}: structure whose fields stand for optional settings |
7871
d2c97f797b1d
Base on optimset mechanism if Octave version >= 3.3.55
i7tiol
parents:
7842
diff
changeset
|
39 ## referred to below. The fields can be set by @code{optimset()} with |
d2c97f797b1d
Base on optimset mechanism if Octave version >= 3.3.55
i7tiol
parents:
7842
diff
changeset
|
40 ## Octave versions 3.3.55 or greater; with older Octave versions, the |
d2c97f797b1d
Base on optimset mechanism if Octave version >= 3.3.55
i7tiol
parents:
7842
diff
changeset
|
41 ## fields must be set directly as structure-fields in the correct case. |
7684
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
42 ## |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
43 ## The returned values are the column vector of final parameters |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
44 ## @var{p}, the final array of residuals @var{resid}, an integer |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
45 ## @var{cvg} indicating if and how optimization succeeded or failed, and |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
46 ## a structure @var{outp} with additional information, curently with |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
47 ## only one field: @var{niter}, the number of iterations. @var{cvg} is |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
48 ## greater than zero for success and less than or equal to zero for |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
49 ## failure; its possible values depend on the used backend and currently |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
50 ## can be @code{0} (maximum number of iterations exceeded), @code{2} |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
51 ## (parameter change less than specified precision in two consecutive |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
52 ## iterations), or @code{3} (improvement in objective function --- e.g. |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
53 ## sum of squares --- less than specified). |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
54 ## |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
55 ## @var{settings}: |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
56 ## |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
57 ## @code{Algorithm}: String specifying the backend. Default: |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
58 ## @code{"lm_svd_feasible"}. The latter is currently the only backend |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
59 ## distributed with this package. It is described in a separate section |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
60 ## below. |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
61 ## |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
62 ## @code{dfdp}: Function computing the jacobian of the residuals with |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
63 ## respect to the parameters, assuming residuals are reshaped to a |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
64 ## vector. Default: finite differences. Will be called with the column |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
65 ## vector of parameters and an informational structure as arguments. The |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
66 ## structure has the fields @code{f}: value of residuals for current |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
67 ## parameters, reshaped to a column vector, @code{fixed}: logical vector |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
68 ## indicating which parameters are not optimized, so these partial |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
69 ## derivatives need not be computed and can be set to zero, |
7842
4f1f5a78b8ed
Allow named parameters to be arrays. Split bounds setting to lower and upper bounds.
i7tiol
parents:
7788
diff
changeset
|
70 ## @code{diffp}, @code{diff_onesided}, @code{lbound}, @code{ubound}: |
4f1f5a78b8ed
Allow named parameters to be arrays. Split bounds setting to lower and upper bounds.
i7tiol
parents:
7788
diff
changeset
|
71 ## identical to the user settings of this name, @code{plabels}: |
4f1f5a78b8ed
Allow named parameters to be arrays. Split bounds setting to lower and upper bounds.
i7tiol
parents:
7788
diff
changeset
|
72 ## 1-dimensional cell-array of column-cell-arrays, each column with |
4f1f5a78b8ed
Allow named parameters to be arrays. Split bounds setting to lower and upper bounds.
i7tiol
parents:
7788
diff
changeset
|
73 ## labels for all parameters, the first column contains the numerical |
4f1f5a78b8ed
Allow named parameters to be arrays. Split bounds setting to lower and upper bounds.
i7tiol
parents:
7788
diff
changeset
|
74 ## indices of the parameters. The default jacobian function will call |
4f1f5a78b8ed
Allow named parameters to be arrays. Split bounds setting to lower and upper bounds.
i7tiol
parents:
7788
diff
changeset
|
75 ## the model function with the second argument set with fields @code{f}: |
4f1f5a78b8ed
Allow named parameters to be arrays. Split bounds setting to lower and upper bounds.
i7tiol
parents:
7788
diff
changeset
|
76 ## as the @code{f} passed to the jacobian function, @code{plabels}: |
4f1f5a78b8ed
Allow named parameters to be arrays. Split bounds setting to lower and upper bounds.
i7tiol
parents:
7788
diff
changeset
|
77 ## cell-array of 1x1 cell-arrays with the entries of the |
4f1f5a78b8ed
Allow named parameters to be arrays. Split bounds setting to lower and upper bounds.
i7tiol
parents:
7788
diff
changeset
|
78 ## column-cell-arrays of @code{plabels} as passed to the jacobian |
4f1f5a78b8ed
Allow named parameters to be arrays. Split bounds setting to lower and upper bounds.
i7tiol
parents:
7788
diff
changeset
|
79 ## function corresponding to current parameter, @code{side}: @code{0} |
4f1f5a78b8ed
Allow named parameters to be arrays. Split bounds setting to lower and upper bounds.
i7tiol
parents:
7788
diff
changeset
|
80 ## for one-sided interval, @code{1} or @code{2}, respectively, for the |
4f1f5a78b8ed
Allow named parameters to be arrays. Split bounds setting to lower and upper bounds.
i7tiol
parents:
7788
diff
changeset
|
81 ## sides of a two-sided interval, and @code{parallel}: logical scalar |
4f1f5a78b8ed
Allow named parameters to be arrays. Split bounds setting to lower and upper bounds.
i7tiol
parents:
7788
diff
changeset
|
82 ## indicating parallel computation of partial derivatives. |
7684
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
83 ## |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
84 ## @code{diffp}: column vector of fractional intervals (doubled for |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
85 ## central intervals) supposed to be used by jacobian functions |
7688 | 86 ## performing finite differencing. Default: @code{.001 * ones (size (parameters))}. The default jacobian function will use these as |
7684
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
87 ## absolute intervals for parameters with value zero. |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
88 ## |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
89 ## @code{diff_onesided}: logical column vector indicating that one-sided |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
90 ## intervals should be used by jacobian functions performing finite |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
91 ## differencing. Default: @code{false (size (parameters))}. |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
92 ## |
8720
f80846665d3b
Use Fotios' complex step derivative approximation. And a few small fixes.
i7tiol
parents:
8088
diff
changeset
|
93 ## @code{complex_step_derivative}, @code{complex_step_derivative_inequc}, |
f80846665d3b
Use Fotios' complex step derivative approximation. And a few small fixes.
i7tiol
parents:
8088
diff
changeset
|
94 ## @code{complex_step_derivative_equc}: logical scalars, default: false. |
f80846665d3b
Use Fotios' complex step derivative approximation. And a few small fixes.
i7tiol
parents:
8088
diff
changeset
|
95 ## Estimate Jacobian of model function, general inequality constraints, |
f80846665d3b
Use Fotios' complex step derivative approximation. And a few small fixes.
i7tiol
parents:
8088
diff
changeset
|
96 ## and general equality constraints, respectively, with complex step |
f80846665d3b
Use Fotios' complex step derivative approximation. And a few small fixes.
i7tiol
parents:
8088
diff
changeset
|
97 ## derivative approximation. Use only if you know that your model |
f80846665d3b
Use Fotios' complex step derivative approximation. And a few small fixes.
i7tiol
parents:
8088
diff
changeset
|
98 ## function, function of general inequality constraints, or function of |
f80846665d3b
Use Fotios' complex step derivative approximation. And a few small fixes.
i7tiol
parents:
8088
diff
changeset
|
99 ## general equality constraints, respectively, is suitable for this. No |
f80846665d3b
Use Fotios' complex step derivative approximation. And a few small fixes.
i7tiol
parents:
8088
diff
changeset
|
100 ## user function for the respective Jacobian must be specified. |
f80846665d3b
Use Fotios' complex step derivative approximation. And a few small fixes.
i7tiol
parents:
8088
diff
changeset
|
101 ## |
f80846665d3b
Use Fotios' complex step derivative approximation. And a few small fixes.
i7tiol
parents:
8088
diff
changeset
|
102 ## @code{cstep}: scalar step size for complex step derivative |
f80846665d3b
Use Fotios' complex step derivative approximation. And a few small fixes.
i7tiol
parents:
8088
diff
changeset
|
103 ## approximation. Default: 1e-20. |
f80846665d3b
Use Fotios' complex step derivative approximation. And a few small fixes.
i7tiol
parents:
8088
diff
changeset
|
104 ## |
7684
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
105 ## @code{fixed}: logical column vector indicating which parameters |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
106 ## should not be optimized, but kept to their inital value. Fixing is |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
107 ## done independently of the backend, but the backend may choose to fix |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
108 ## additional parameters under certain conditions. |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
109 ## |
7842
4f1f5a78b8ed
Allow named parameters to be arrays. Split bounds setting to lower and upper bounds.
i7tiol
parents:
7788
diff
changeset
|
110 ## @code{lbound}, @code{ubound}: column vectors of lower and upper |
4f1f5a78b8ed
Allow named parameters to be arrays. Split bounds setting to lower and upper bounds.
i7tiol
parents:
7788
diff
changeset
|
111 ## bounds for parameters. Default: @code{-Inf} and @code{+Inf}, |
7684
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
112 ## respectively. The bounds are non-strict, i.e. parameters are allowed |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
113 ## to be exactly equal to a bound. The default jacobian function will |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
114 ## respect bounds (but no further inequality constraints) in finite |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
115 ## differencing. |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
116 ## |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
117 ## @code{inequc}: Further inequality constraints. Cell-array containing |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
118 ## up to four entries, two entries for linear inequality constraints |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
119 ## and/or one or two entries for general inequality constraints. Either |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
120 ## linear or general constraints may be the first entries, but the two |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
121 ## entries for linear constraints must be adjacent and, if two entries |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
122 ## are given for general constraints, they also must be adjacent. The |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
123 ## two entries for linear constraints are a matrix (say @code{m}) and a |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
124 ## vector (say @code{v}), specifying linear inequality constraints of |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
125 ## the form @code{m.' * parameters + v >= 0}. The first entry for |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
126 ## general constraints must be a differentiable vector valued function |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
127 ## (say @code{h}), specifying general inequality constraints of the form |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
128 ## @code{h (p[, idx]) >= 0}; @code{p} is the column vector of optimized |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
129 ## paraters and the optional argument @code{idx} is a logical index. |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
130 ## @code{h} has to return the values of all constraints if @code{idx} is |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
131 ## not given, and has to return only the indexed constraints if |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
132 ## @code{idx} is given (so computation of the other constraints can be |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
133 ## spared). In gradient determination, this function may be called with |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
134 ## an informational third argument, whose content depends on the |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
135 ## function for gradient determination. If a second entry for general |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
136 ## inequality constraints is given, it must be a function computing the |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
137 ## jacobian of the constraints with respect to the parameters. For this |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
138 ## function, the description of @code{dfdp} above applies, except that |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
139 ## it is called with 3 arguments since it has an additional argument |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
140 ## @code{idx} --- a logical index --- at second position, indicating |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
141 ## which rows of the jacobian must be returned, and except that the |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
142 ## default function calls @code{h} with 3 arguments, since the argument |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
143 ## @code{idx} is also supplied. Note that specifying linear constraints |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
144 ## as general constraints will generally waste performance, even if |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
145 ## further, non-linear, general constraints are also specified. |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
146 ## |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
147 ## @code{equc}: Equality constraints. Specified the same way as |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
148 ## inequality constraints (see @code{inequc}). |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
149 ## |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
150 ## @code{cpiv}: Function for complementary pivoting, usable in |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
151 ## algorithms for constraints. Default: @ cpiv_bard. Only the default |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
152 ## function is supplied with the package. |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
153 ## |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
154 ## @code{weights}: Array of weights for the residuals. Dimensions must |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
155 ## match. |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
156 ## |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
157 ## @code{TolFun}: Minimum fractional improvement in objective function |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
158 ## (e.g. sum of squares) in an iteration (abortion criterium). Default: |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
159 ## .0001. |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
160 ## |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
161 ## @code{MaxIter}: Maximum number of iterations (abortion criterium). |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
162 ## Default: backend-specific. |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
163 ## |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
164 ## @code{fract_prec}: Column Vector, minimum fractional change of |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
165 ## parameters in an iteration (abortion criterium if violated in two |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
166 ## consecutive iterations). Default: backend-specific. |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
167 ## |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
168 ## @code{max_fract_change}: Column Vector, enforced maximum fractional |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
169 ## change in parameters in an iteration. Default: backend-specific. |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
170 ## |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
171 ## @code{Display}: String indicating the degree of verbosity. Default: |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
172 ## @code{"off"}. Possible values are currently @code{"off"} (no |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
173 ## messages) and @code{"iter"} (some messages after each iteration). |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
174 ## Support of this setting and its exact interpretation are |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
175 ## backend-specific. |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
176 ## |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
177 ## @code{plot_cmd}: Function enabling backend to plot results or |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
178 ## intermediate results. Will be called with current computed |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
179 ## residualse. Default: plot nothing. |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
180 ## |
7782
ef1dc1539dcf
Slight changes to equality constraints. Options for testing.
i7tiol
parents:
7734
diff
changeset
|
181 ## @code{debug}: Logical scalar, default: @code{false}. Will be passed |
ef1dc1539dcf
Slight changes to equality constraints. Options for testing.
i7tiol
parents:
7734
diff
changeset
|
182 ## to the backend, which might print debugging information if true. |
ef1dc1539dcf
Slight changes to equality constraints. Options for testing.
i7tiol
parents:
7734
diff
changeset
|
183 ## |
7684
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
184 ## Structure-based parameter handling |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
185 ## |
7842
4f1f5a78b8ed
Allow named parameters to be arrays. Split bounds setting to lower and upper bounds.
i7tiol
parents:
7788
diff
changeset
|
186 ## The setting @code{param_order} is a cell-array with names of the |
4f1f5a78b8ed
Allow named parameters to be arrays. Split bounds setting to lower and upper bounds.
i7tiol
parents:
7788
diff
changeset
|
187 ## optimized parameters. If not given, and initial parameters are a |
8088
da5ad60e2221
Replace function calls: cell2cell -> num2cell, partarray -> mat2cell.
i7tiol
parents:
7951
diff
changeset
|
188 ## structure, all parameters in the structure are optimized. If initial |
da5ad60e2221
Replace function calls: cell2cell -> num2cell, partarray -> mat2cell.
i7tiol
parents:
7951
diff
changeset
|
189 ## parameters are a structure, it is an error if @code{param_order} is |
da5ad60e2221
Replace function calls: cell2cell -> num2cell, partarray -> mat2cell.
i7tiol
parents:
7951
diff
changeset
|
190 ## not given and there are any non-structure-based configuration items |
da5ad60e2221
Replace function calls: cell2cell -> num2cell, partarray -> mat2cell.
i7tiol
parents:
7951
diff
changeset
|
191 ## or functions. |
7684
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
192 ## |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
193 ## The initial parameters @var{pin} can be given as a structure |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
194 ## containing at least all fields named in @code{param_order}. In this |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
195 ## case the returned parameters @var{p} will also be a structure. |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
196 ## |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
197 ## Each user-supplied function can be called with the argument |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
198 ## containing the current parameters being a structure instead of a |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
199 ## column vector. For this, a corresponding setting must be set to |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
200 ## @code{true}: @code{f_pstruct} (model function), @code{dfdp_pstruct} |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
201 ## (jacobian of model function), @code{f_inequc_pstruct} (general |
7688 | 202 ## inequality constraints), @code{df_inequc_pstruct} (jacobian of |
203 ## general inequality constraints), @code{f_equc_pstruct} (general | |
204 ## equality constraints), and @code{df_equc_pstruct} (jacobian of | |
205 ## general equality constraints). If a jacobian-function is configured | |
206 ## in such a way, it must return the columns of the jacobian as fields | |
207 ## of a structure under the respective parameter names. | |
7684
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
208 ## |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
209 ## Similarly, for specifying linear constraints, instead of the matrix |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
210 ## (called @code{m} above), a structure containing the rows of the |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
211 ## matrix in fields under the respective parameter names can be given. |
7722 | 212 ## In this case, rows containing only zeros need not be given. |
7684
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
213 ## |
7842
4f1f5a78b8ed
Allow named parameters to be arrays. Split bounds setting to lower and upper bounds.
i7tiol
parents:
7788
diff
changeset
|
214 ## The vector-based settings @code{lbound}, @code{ubound}, |
4f1f5a78b8ed
Allow named parameters to be arrays. Split bounds setting to lower and upper bounds.
i7tiol
parents:
7788
diff
changeset
|
215 ## @code{fixed}, @code{diffp}, @code{diff_onesided}, @code{fract_prec}, |
4f1f5a78b8ed
Allow named parameters to be arrays. Split bounds setting to lower and upper bounds.
i7tiol
parents:
7788
diff
changeset
|
216 ## and @code{max_fract_change} can be replaced by the setting |
7788
ad5f34e55a1d
param_config need only contain parameters with non-defaults.
i7tiol
parents:
7782
diff
changeset
|
217 ## @code{param_config}. It is a structure that can contain fields named |
ad5f34e55a1d
param_config need only contain parameters with non-defaults.
i7tiol
parents:
7782
diff
changeset
|
218 ## in @code{param_order}. For each such field, there may be subfields |
7842
4f1f5a78b8ed
Allow named parameters to be arrays. Split bounds setting to lower and upper bounds.
i7tiol
parents:
7788
diff
changeset
|
219 ## with the same names as the above vector-based settings, but |
4f1f5a78b8ed
Allow named parameters to be arrays. Split bounds setting to lower and upper bounds.
i7tiol
parents:
7788
diff
changeset
|
220 ## containing a scalar value for the respective parameter. If |
4f1f5a78b8ed
Allow named parameters to be arrays. Split bounds setting to lower and upper bounds.
i7tiol
parents:
7788
diff
changeset
|
221 ## @code{param_config} is specified, none of the above |
4f1f5a78b8ed
Allow named parameters to be arrays. Split bounds setting to lower and upper bounds.
i7tiol
parents:
7788
diff
changeset
|
222 ## vector/matrix-based settings may be used. |
4f1f5a78b8ed
Allow named parameters to be arrays. Split bounds setting to lower and upper bounds.
i7tiol
parents:
7788
diff
changeset
|
223 ## |
4f1f5a78b8ed
Allow named parameters to be arrays. Split bounds setting to lower and upper bounds.
i7tiol
parents:
7788
diff
changeset
|
224 ## Additionally, named parameters are allowed to be non-scalar real |
4f1f5a78b8ed
Allow named parameters to be arrays. Split bounds setting to lower and upper bounds.
i7tiol
parents:
7788
diff
changeset
|
225 ## arrays. In this case, their dimensions are given by the setting |
4f1f5a78b8ed
Allow named parameters to be arrays. Split bounds setting to lower and upper bounds.
i7tiol
parents:
7788
diff
changeset
|
226 ## @code{param_dims}, a cell-array of dimension vectors, each containing |
4f1f5a78b8ed
Allow named parameters to be arrays. Split bounds setting to lower and upper bounds.
i7tiol
parents:
7788
diff
changeset
|
227 ## at least two dimensions; if not given, dimensions are taken from the |
4f1f5a78b8ed
Allow named parameters to be arrays. Split bounds setting to lower and upper bounds.
i7tiol
parents:
7788
diff
changeset
|
228 ## initial parameters, if these are given in a structure. Any |
4f1f5a78b8ed
Allow named parameters to be arrays. Split bounds setting to lower and upper bounds.
i7tiol
parents:
7788
diff
changeset
|
229 ## vector-based settings or not structure-based linear constraints then |
4f1f5a78b8ed
Allow named parameters to be arrays. Split bounds setting to lower and upper bounds.
i7tiol
parents:
7788
diff
changeset
|
230 ## must correspond to an order of parameters with all parameters |
4f1f5a78b8ed
Allow named parameters to be arrays. Split bounds setting to lower and upper bounds.
i7tiol
parents:
7788
diff
changeset
|
231 ## reshaped to vectors and concatenated in the user-given order of |
4f1f5a78b8ed
Allow named parameters to be arrays. Split bounds setting to lower and upper bounds.
i7tiol
parents:
7788
diff
changeset
|
232 ## parameter names. Structure-based settings or structure-based initial |
4f1f5a78b8ed
Allow named parameters to be arrays. Split bounds setting to lower and upper bounds.
i7tiol
parents:
7788
diff
changeset
|
233 ## parameters must contain arrays with dimensions reshapable to those of |
4f1f5a78b8ed
Allow named parameters to be arrays. Split bounds setting to lower and upper bounds.
i7tiol
parents:
7788
diff
changeset
|
234 ## the respective parameters. |
7684
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
235 ## |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
236 ## Description of backends (currently only one) |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
237 ## |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
238 ## "lm_svd_feasible" |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
239 ## |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
240 ## A Levenberg/Marquardt algorithm using singular value decomposition |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
241 ## and featuring constraints which must be met by the initial parameters |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
242 ## and are attempted to be kept met throughout the optimization. |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
243 ## |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
244 ## Parameters with identical lower and upper bounds will be fixed. |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
245 ## |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
246 ## Returned value @var{cvg} will be @code{0}, @code{2}, or @code{3}. |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
247 ## |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
248 ## Backend-specific defaults are: @code{MaxIter}: 20, @code{fract_prec}: |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
249 ## @code{zeros (size (parameters))}, @code{max_fract_change}: @code{Inf} |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
250 ## for all parameters. |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
251 ## |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
252 ## Interpretation of @code{Display}: if set to @code{"iter"}, currently |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
253 ## @code{plot_cmd} is evaluated for each iteration, and some further |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
254 ## diagnostics may be printed. |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
255 ## |
7782
ef1dc1539dcf
Slight changes to equality constraints. Options for testing.
i7tiol
parents:
7734
diff
changeset
|
256 ## Specific option: @code{lm_svd_feasible_alt_s}: if falling back to |
ef1dc1539dcf
Slight changes to equality constraints. Options for testing.
i7tiol
parents:
7734
diff
changeset
|
257 ## nearly gradient descent, do it more like original Levenberg/Marquardt |
ef1dc1539dcf
Slight changes to equality constraints. Options for testing.
i7tiol
parents:
7734
diff
changeset
|
258 ## method, with descent in each gradient component; for testing only. |
ef1dc1539dcf
Slight changes to equality constraints. Options for testing.
i7tiol
parents:
7734
diff
changeset
|
259 ## |
7684
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
260 ## @seealso {nonlin_curvefit} |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
261 ## @end deftypefn |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
262 |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
263 function [p, resid, cvg, outp] = nonlin_residmin (varargin) |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
264 |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
265 if (nargin == 1) |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
266 p = __nonlin_residmin__ (varargin{1}); |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
267 return; |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
268 endif |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
269 |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
270 if (nargin < 2 || nargin > 3) |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
271 print_usage (); |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
272 endif |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
273 |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
274 if (nargin == 2) |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
275 varargin{3} = struct (); |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
276 endif |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
277 |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
278 varargin{4} = struct (); |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
279 |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
280 [p, resid, cvg, outp] = __nonlin_residmin__ (varargin{:}); |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
281 |
d4f8cb9c7a93
Added nonlin_residmin and nonlin_curvefit with a common backend with leasqr.
i7tiol
parents:
diff
changeset
|
282 endfunction |