Mercurial > hg > octave-lyh
annotate scripts/optimization/fsolve.m @ 9051:1bf0ce0930be
Grammar check TexInfo in all .m files
Cleanup documentation sources to follow a few consistent rules.
Spellcheck was NOT done. (but will be in another changeset)
| author | Rik <rdrider0-list@yahoo.com> |
|---|---|
| date | Fri, 27 Mar 2009 22:31:03 -0700 |
| parents | 187a9d9c2f04 |
| children | 8207b833557f |
| rev | line source |
|---|---|
| 8590 | 1 ## Copyright (C) 2008, 2009 VZLU Prague, a.s. |
|
8306
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
2 ## |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
3 ## This file is part of Octave. |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
4 ## |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
5 ## Octave is free software; you can redistribute it and/or modify it |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
6 ## under the terms of the GNU General Public License as published by |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
7 ## the Free Software Foundation; either version 3 of the License, or (at |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
8 ## your option) any later version. |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
9 ## |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
10 ## Octave is distributed in the hope that it will be useful, but |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
11 ## WITHOUT ANY WARRANTY; without even the implied warranty of |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
12 ## MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
13 ## General Public License for more details. |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
14 ## |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
15 ## You should have received a copy of the GNU General Public License |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
16 ## along with Octave; see the file COPYING. If not, see |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
17 ## <http://www.gnu.org/licenses/>. |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
18 ## |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
19 ## Author: Jaroslav Hajek <highegg@gmail.com> |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
20 |
| 8466 | 21 ## -*- texinfo -*- |
| 8513 | 22 ## @deftypefn{Function File} {} fsolve (@var{fcn}, @var{x0}, @var{options}) |
| 8466 | 23 ## @deftypefnx{Function File} {[@var{x}, @var{fvec}, @var{info}, @var{output}, @var{fjac}]} = fsolve (@var{fcn}, @dots{}) |
| 8513 | 24 ## Solve a system of nonlinear equations defined by the function @var{fcn}. |
| 8466 | 25 ## @var{fcn} should accepts a vector (array) defining the unknown variables, |
|
9051
1bf0ce0930be
Grammar check TexInfo in all .m files
Rik <rdrider0-list@yahoo.com>
parents:
8997
diff
changeset
|
26 ## and return a vector of left-hand sides of the equations. Right-hand sides |
| 8466 | 27 ## are defined to be zeros. |
|
8514
39867b4aca52
fsolve.m: additional doc fix
John W. Eaton <jwe@octave.org>
parents:
8513
diff
changeset
|
28 ## In other words, this function attempts to determine a vector @var{x} such |
|
39867b4aca52
fsolve.m: additional doc fix
John W. Eaton <jwe@octave.org>
parents:
8513
diff
changeset
|
29 ## that @code{@var{fcn} (@var{x})} gives (approximately) all zeros. |
|
9051
1bf0ce0930be
Grammar check TexInfo in all .m files
Rik <rdrider0-list@yahoo.com>
parents:
8997
diff
changeset
|
30 ## @var{x0} determines a starting guess. The shape of @var{x0} is preserved |
| 8466 | 31 ## in all calls to @var{fcn}, but otherwise it is treated as a column vector. |
|
8515
ec2715c76039
fzero.m, fsolve.m: additional doc fixes
John W. Eaton <jwe@octave.org>
parents:
8514
diff
changeset
|
32 ## @var{options} is a structure specifying additional options. |
|
ec2715c76039
fzero.m, fsolve.m: additional doc fixes
John W. Eaton <jwe@octave.org>
parents:
8514
diff
changeset
|
33 ## Currently, @code{fsolve} recognizes these options: |
|
ec2715c76039
fzero.m, fsolve.m: additional doc fixes
John W. Eaton <jwe@octave.org>
parents:
8514
diff
changeset
|
34 ## @code{"FunValCheck"}, @code{"OutputFcn"}, @code{"TolX"}, |
| 8590 | 35 ## @code{"TolFun"}, @code{"MaxIter"}, @code{"MaxFunEvals"}, |
|
8693
e5ffb52c9c61
improve fsolve and add ComplexEqn option
Jaroslav Hajek <highegg@gmail.com>
parents:
8661
diff
changeset
|
36 ## @code{"Jacobian"}, @code{"Updating"} and @code{"ComplexEqn"}. |
| 8466 | 37 ## |
| 8513 | 38 ## If @code{"Jacobian"} is @code{"on"}, it specifies that @var{fcn}, |
| 39 ## called with 2 output arguments, also returns the Jacobian matrix | |
| 40 ## of right-hand sides at the requested point. @code{"TolX"} specifies | |
| 41 ## the termination tolerance in the unknown variables, while | |
|
9051
1bf0ce0930be
Grammar check TexInfo in all .m files
Rik <rdrider0-list@yahoo.com>
parents:
8997
diff
changeset
|
42 ## @code{"TolFun"} is a tolerance for equations. Default is @code{1e-7} |
|
8596
8833c0b18eb2
enable default settings queries in optim funcs
Jaroslav Hajek <highegg@gmail.com>
parents:
8592
diff
changeset
|
43 ## for both @code{"TolX"} and @code{"TolFun"}. |
|
8833c0b18eb2
enable default settings queries in optim funcs
Jaroslav Hajek <highegg@gmail.com>
parents:
8592
diff
changeset
|
44 ## If @code{"Updating"} is "on", the function will attempt to use Broyden |
| 8590 | 45 ## updates to update the Jacobian, in order to reduce the amount of jacobian |
| 46 ## calculations. | |
| 47 ## If your user function always calculates the Jacobian (regardless of number | |
| 48 ## of output arguments), this option provides no advantage and should be set to | |
| 49 ## false. | |
| 50 ## | |
|
8693
e5ffb52c9c61
improve fsolve and add ComplexEqn option
Jaroslav Hajek <highegg@gmail.com>
parents:
8661
diff
changeset
|
51 ## @code{"ComplexEqn"} is @code{"on"}, @code{fsolve} will attempt to solve |
|
e5ffb52c9c61
improve fsolve and add ComplexEqn option
Jaroslav Hajek <highegg@gmail.com>
parents:
8661
diff
changeset
|
52 ## complex equations in complex variables, assuming that the equations posess a |
|
9051
1bf0ce0930be
Grammar check TexInfo in all .m files
Rik <rdrider0-list@yahoo.com>
parents:
8997
diff
changeset
|
53 ## complex derivative (i.e., are holomorphic). If this is not what you want, |
|
8693
e5ffb52c9c61
improve fsolve and add ComplexEqn option
Jaroslav Hajek <highegg@gmail.com>
parents:
8661
diff
changeset
|
54 ## should unpack the real and imaginary parts of the system to get a real |
|
e5ffb52c9c61
improve fsolve and add ComplexEqn option
Jaroslav Hajek <highegg@gmail.com>
parents:
8661
diff
changeset
|
55 ## system. |
|
e5ffb52c9c61
improve fsolve and add ComplexEqn option
Jaroslav Hajek <highegg@gmail.com>
parents:
8661
diff
changeset
|
56 ## |
| 8466 | 57 ## For description of the other options, see @code{optimset}. |
| 58 ## | |
| 59 ## On return, @var{fval} contains the value of the function @var{fcn} | |
| 60 ## evaluated at @var{x}, and @var{info} may be one of the following values: | |
| 61 ## | |
| 62 ## @table @asis | |
| 63 ## @item 1 | |
|
9051
1bf0ce0930be
Grammar check TexInfo in all .m files
Rik <rdrider0-list@yahoo.com>
parents:
8997
diff
changeset
|
64 ## Converged to a solution point. Relative residual error is less than specified |
| 8466 | 65 ## by TolFun. |
| 66 ## @item 2 | |
| 67 ## Last relative step size was less that TolX. | |
| 68 ## @item 3 | |
| 69 ## Last relative decrease in residual was less than TolF. | |
| 70 ## @item 0 | |
| 71 ## Iteration limit exceeded. | |
| 72 ## @item -3 | |
| 73 ## The trust region radius became excessively small. | |
| 74 ## @end table | |
|
8986
22c8272af34b
improvements to fsolve & family
Jaroslav Hajek <highegg@gmail.com>
parents:
8857
diff
changeset
|
75 ## |
| 8466 | 76 ## Note: If you only have a single nonlinear equation of one variable, using |
| 77 ## @code{fzero} is usually a much better idea. | |
| 8513 | 78 ## @seealso{fzero, optimset} |
|
8986
22c8272af34b
improvements to fsolve & family
Jaroslav Hajek <highegg@gmail.com>
parents:
8857
diff
changeset
|
79 ## |
|
22c8272af34b
improvements to fsolve & family
Jaroslav Hajek <highegg@gmail.com>
parents:
8857
diff
changeset
|
80 ## Note about user-supplied jacobians: |
|
22c8272af34b
improvements to fsolve & family
Jaroslav Hajek <highegg@gmail.com>
parents:
8857
diff
changeset
|
81 ## As an inherent property of the algorithm, jacobian is always requested for a |
|
22c8272af34b
improvements to fsolve & family
Jaroslav Hajek <highegg@gmail.com>
parents:
8857
diff
changeset
|
82 ## solution vector whose residual vector is already known, and it is the last |
|
9051
1bf0ce0930be
Grammar check TexInfo in all .m files
Rik <rdrider0-list@yahoo.com>
parents:
8997
diff
changeset
|
83 ## accepted successful step. Often this will be one of the last two calls, but |
|
1bf0ce0930be
Grammar check TexInfo in all .m files
Rik <rdrider0-list@yahoo.com>
parents:
8997
diff
changeset
|
84 ## not always. If the savings by reusing intermediate results from residual |
|
8986
22c8272af34b
improvements to fsolve & family
Jaroslav Hajek <highegg@gmail.com>
parents:
8857
diff
changeset
|
85 ## calculation in jacobian calculation are significant, the best strategy is to |
|
22c8272af34b
improvements to fsolve & family
Jaroslav Hajek <highegg@gmail.com>
parents:
8857
diff
changeset
|
86 ## employ OutputFcn: After a vector is evaluated for residuals, if OutputFcn is |
|
22c8272af34b
improvements to fsolve & family
Jaroslav Hajek <highegg@gmail.com>
parents:
8857
diff
changeset
|
87 ## called with that vector, then the intermediate results should be saved for |
|
22c8272af34b
improvements to fsolve & family
Jaroslav Hajek <highegg@gmail.com>
parents:
8857
diff
changeset
|
88 ## future jacobian evaluation, and should be kept until a jacobian evaluation |
|
22c8272af34b
improvements to fsolve & family
Jaroslav Hajek <highegg@gmail.com>
parents:
8857
diff
changeset
|
89 ## is requested or until outputfcn is called with a different vector, in which |
|
9051
1bf0ce0930be
Grammar check TexInfo in all .m files
Rik <rdrider0-list@yahoo.com>
parents:
8997
diff
changeset
|
90 ## case they should be dropped in favor of this most recent vector. A short |
|
8986
22c8272af34b
improvements to fsolve & family
Jaroslav Hajek <highegg@gmail.com>
parents:
8857
diff
changeset
|
91 ## example how this can be achieved follows: |
|
8988
315828058e0d
minor doc fixes in fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
8986
diff
changeset
|
92 ## |
|
8986
22c8272af34b
improvements to fsolve & family
Jaroslav Hajek <highegg@gmail.com>
parents:
8857
diff
changeset
|
93 ## @example |
|
9051
1bf0ce0930be
Grammar check TexInfo in all .m files
Rik <rdrider0-list@yahoo.com>
parents:
8997
diff
changeset
|
94 ## @group |
|
8986
22c8272af34b
improvements to fsolve & family
Jaroslav Hajek <highegg@gmail.com>
parents:
8857
diff
changeset
|
95 ## function [fvec, fjac] = my_optim_func (x, optimvalues, state) |
|
22c8272af34b
improvements to fsolve & family
Jaroslav Hajek <highegg@gmail.com>
parents:
8857
diff
changeset
|
96 ## persistent sav = [], sav0 = []; |
|
22c8272af34b
improvements to fsolve & family
Jaroslav Hajek <highegg@gmail.com>
parents:
8857
diff
changeset
|
97 ## if (nargin == 1) |
|
22c8272af34b
improvements to fsolve & family
Jaroslav Hajek <highegg@gmail.com>
parents:
8857
diff
changeset
|
98 ## ## evaluation call |
|
22c8272af34b
improvements to fsolve & family
Jaroslav Hajek <highegg@gmail.com>
parents:
8857
diff
changeset
|
99 ## if (nargout == 1) |
|
22c8272af34b
improvements to fsolve & family
Jaroslav Hajek <highegg@gmail.com>
parents:
8857
diff
changeset
|
100 ## sav0.x = x; # mark saved vector |
|
22c8272af34b
improvements to fsolve & family
Jaroslav Hajek <highegg@gmail.com>
parents:
8857
diff
changeset
|
101 ## ## calculate fvec, save results to sav0. |
|
22c8272af34b
improvements to fsolve & family
Jaroslav Hajek <highegg@gmail.com>
parents:
8857
diff
changeset
|
102 ## elseif (nargout == 2) |
|
22c8272af34b
improvements to fsolve & family
Jaroslav Hajek <highegg@gmail.com>
parents:
8857
diff
changeset
|
103 ## ## calculate fjac using sav. |
|
22c8272af34b
improvements to fsolve & family
Jaroslav Hajek <highegg@gmail.com>
parents:
8857
diff
changeset
|
104 ## endif |
|
22c8272af34b
improvements to fsolve & family
Jaroslav Hajek <highegg@gmail.com>
parents:
8857
diff
changeset
|
105 ## else |
|
22c8272af34b
improvements to fsolve & family
Jaroslav Hajek <highegg@gmail.com>
parents:
8857
diff
changeset
|
106 ## ## outputfcn call. |
|
22c8272af34b
improvements to fsolve & family
Jaroslav Hajek <highegg@gmail.com>
parents:
8857
diff
changeset
|
107 ## if (all (x == sav0.x)) |
|
22c8272af34b
improvements to fsolve & family
Jaroslav Hajek <highegg@gmail.com>
parents:
8857
diff
changeset
|
108 ## sav = sav0; |
|
22c8272af34b
improvements to fsolve & family
Jaroslav Hajek <highegg@gmail.com>
parents:
8857
diff
changeset
|
109 ## endif |
|
22c8272af34b
improvements to fsolve & family
Jaroslav Hajek <highegg@gmail.com>
parents:
8857
diff
changeset
|
110 ## ## maybe output iteration status etc. |
|
22c8272af34b
improvements to fsolve & family
Jaroslav Hajek <highegg@gmail.com>
parents:
8857
diff
changeset
|
111 ## endif |
|
22c8272af34b
improvements to fsolve & family
Jaroslav Hajek <highegg@gmail.com>
parents:
8857
diff
changeset
|
112 ## endfunction |
|
22c8272af34b
improvements to fsolve & family
Jaroslav Hajek <highegg@gmail.com>
parents:
8857
diff
changeset
|
113 ## |
|
9051
1bf0ce0930be
Grammar check TexInfo in all .m files
Rik <rdrider0-list@yahoo.com>
parents:
8997
diff
changeset
|
114 ## ## @dots{}. |
|
8986
22c8272af34b
improvements to fsolve & family
Jaroslav Hajek <highegg@gmail.com>
parents:
8857
diff
changeset
|
115 ## |
|
8988
315828058e0d
minor doc fixes in fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
8986
diff
changeset
|
116 ## fsolve (@@my_optim_func, x0, optimset ("OutputFcn", @@my_optim_func, @dots{})) |
|
9051
1bf0ce0930be
Grammar check TexInfo in all .m files
Rik <rdrider0-list@yahoo.com>
parents:
8997
diff
changeset
|
117 ## @end group |
|
8986
22c8272af34b
improvements to fsolve & family
Jaroslav Hajek <highegg@gmail.com>
parents:
8857
diff
changeset
|
118 ## @end example |
|
22c8272af34b
improvements to fsolve & family
Jaroslav Hajek <highegg@gmail.com>
parents:
8857
diff
changeset
|
119 ### |
| 8466 | 120 ## @end deftypefn |
|
8306
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
121 |
|
8648
ff61b53eb294
optimization: use PKG_ADD: comments instead of PKG_ADD file
John W. Eaton <jwe@octave.org>
parents:
8647
diff
changeset
|
122 ## PKG_ADD: __all_opts__ ("fsolve"); |
|
ff61b53eb294
optimization: use PKG_ADD: comments instead of PKG_ADD file
John W. Eaton <jwe@octave.org>
parents:
8647
diff
changeset
|
123 |
|
8596
8833c0b18eb2
enable default settings queries in optim funcs
Jaroslav Hajek <highegg@gmail.com>
parents:
8592
diff
changeset
|
124 function [x, fvec, info, output, fjac] = fsolve (fcn, x0, options = struct ()) |
|
8306
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
125 |
|
8596
8833c0b18eb2
enable default settings queries in optim funcs
Jaroslav Hajek <highegg@gmail.com>
parents:
8592
diff
changeset
|
126 ## Get default options if requested. |
|
8833c0b18eb2
enable default settings queries in optim funcs
Jaroslav Hajek <highegg@gmail.com>
parents:
8592
diff
changeset
|
127 if (nargin == 1 && ischar (fcn) && strcmp (fcn, 'defaults')) |
|
8833c0b18eb2
enable default settings queries in optim funcs
Jaroslav Hajek <highegg@gmail.com>
parents:
8592
diff
changeset
|
128 x = optimset ("MaxIter", 400, "MaxFunEvals", Inf, \ |
|
8647
06f5dd901f30
implement registering of optimization options
Jaroslav Hajek <highegg@gmail.com>
parents:
8616
diff
changeset
|
129 "Jacobian", "off", "TolX", 1.5e-8, "TolFun", 1.5e-8, |
|
8693
e5ffb52c9c61
improve fsolve and add ComplexEqn option
Jaroslav Hajek <highegg@gmail.com>
parents:
8661
diff
changeset
|
130 "OutputFcn", [], "Updating", "on", "FunValCheck", "off", |
|
e5ffb52c9c61
improve fsolve and add ComplexEqn option
Jaroslav Hajek <highegg@gmail.com>
parents:
8661
diff
changeset
|
131 "ComplexEqn", "off"); |
|
8596
8833c0b18eb2
enable default settings queries in optim funcs
Jaroslav Hajek <highegg@gmail.com>
parents:
8592
diff
changeset
|
132 return; |
|
8833c0b18eb2
enable default settings queries in optim funcs
Jaroslav Hajek <highegg@gmail.com>
parents:
8592
diff
changeset
|
133 endif |
|
8833c0b18eb2
enable default settings queries in optim funcs
Jaroslav Hajek <highegg@gmail.com>
parents:
8592
diff
changeset
|
134 |
|
8833c0b18eb2
enable default settings queries in optim funcs
Jaroslav Hajek <highegg@gmail.com>
parents:
8592
diff
changeset
|
135 if (nargin < 2 || nargin > 3 || ! ismatrix (x0)) |
|
8833c0b18eb2
enable default settings queries in optim funcs
Jaroslav Hajek <highegg@gmail.com>
parents:
8592
diff
changeset
|
136 print_usage (); |
|
8833c0b18eb2
enable default settings queries in optim funcs
Jaroslav Hajek <highegg@gmail.com>
parents:
8592
diff
changeset
|
137 endif |
|
8833c0b18eb2
enable default settings queries in optim funcs
Jaroslav Hajek <highegg@gmail.com>
parents:
8592
diff
changeset
|
138 |
|
8833c0b18eb2
enable default settings queries in optim funcs
Jaroslav Hajek <highegg@gmail.com>
parents:
8592
diff
changeset
|
139 if (ischar (fcn)) |
|
8833c0b18eb2
enable default settings queries in optim funcs
Jaroslav Hajek <highegg@gmail.com>
parents:
8592
diff
changeset
|
140 fcn = str2func (fcn); |
|
8306
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
141 endif |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
142 |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
143 xsiz = size (x0); |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
144 n = numel (x0); |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
145 |
|
8467
77b8d4aa2743
fsolve.m, fzero.m: style fixes; use strcmpi to compare options
John W. Eaton <jwe@octave.org>
parents:
8466
diff
changeset
|
146 has_jac = strcmpi (optimget (options, "Jacobian", "off"), "on"); |
|
8596
8833c0b18eb2
enable default settings queries in optim funcs
Jaroslav Hajek <highegg@gmail.com>
parents:
8592
diff
changeset
|
147 maxiter = optimget (options, "MaxIter", 400); |
|
8306
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
148 maxfev = optimget (options, "MaxFunEvals", Inf); |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
149 outfcn = optimget (options, "OutputFcn"); |
|
8596
8833c0b18eb2
enable default settings queries in optim funcs
Jaroslav Hajek <highegg@gmail.com>
parents:
8592
diff
changeset
|
150 updating = strcmpi (optimget (options, "Updating", "on"), "on"); |
|
8693
e5ffb52c9c61
improve fsolve and add ComplexEqn option
Jaroslav Hajek <highegg@gmail.com>
parents:
8661
diff
changeset
|
151 complexeqn = strcmpi (optimget (options, "ComplexEqn", "off"), "on"); |
| 8590 | 152 |
|
8467
77b8d4aa2743
fsolve.m, fzero.m: style fixes; use strcmpi to compare options
John W. Eaton <jwe@octave.org>
parents:
8466
diff
changeset
|
153 funvalchk = strcmpi (optimget (options, "FunValCheck", "off"), "on"); |
|
8306
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
154 |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
155 if (funvalchk) |
| 8604 | 156 ## Replace fcn with a guarded version. |
|
8693
e5ffb52c9c61
improve fsolve and add ComplexEqn option
Jaroslav Hajek <highegg@gmail.com>
parents:
8661
diff
changeset
|
157 fcn = @(x) guarded_eval (fcn, x, complexeqn); |
|
8306
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
158 endif |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
159 |
| 8466 | 160 ## These defaults are rather stringent. I think that normally, user |
| 161 ## prefers accuracy to performance. | |
|
8306
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
162 |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
163 macheps = eps (class (x0)); |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
164 |
|
8613
38482007c834
relax scaling in fsolve and use class-dependent default tolerances
Jaroslav Hajek <highegg@gmail.com>
parents:
8604
diff
changeset
|
165 tolx = optimget (options, "TolX", sqrt (macheps)); |
|
38482007c834
relax scaling in fsolve and use class-dependent default tolerances
Jaroslav Hajek <highegg@gmail.com>
parents:
8604
diff
changeset
|
166 tolf = optimget (options, "TolFun", sqrt (macheps)); |
|
8306
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
167 |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
168 factor = 100; |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
169 |
|
8467
77b8d4aa2743
fsolve.m, fzero.m: style fixes; use strcmpi to compare options
John W. Eaton <jwe@octave.org>
parents:
8466
diff
changeset
|
170 niter = 1; |
|
8986
22c8272af34b
improvements to fsolve & family
Jaroslav Hajek <highegg@gmail.com>
parents:
8857
diff
changeset
|
171 nfev = 1; |
|
8306
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
172 |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
173 x = x0(:); |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
174 info = 0; |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
175 |
|
8986
22c8272af34b
improvements to fsolve & family
Jaroslav Hajek <highegg@gmail.com>
parents:
8857
diff
changeset
|
176 ## Initial evaluation. |
|
8990
349d75161672
more cosmetic adjustments to fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
8988
diff
changeset
|
177 ## Handle arbitrary shapes of x and f and remember them. |
|
8986
22c8272af34b
improvements to fsolve & family
Jaroslav Hajek <highegg@gmail.com>
parents:
8857
diff
changeset
|
178 fvec = fcn (reshape (x, xsiz)); |
|
22c8272af34b
improvements to fsolve & family
Jaroslav Hajek <highegg@gmail.com>
parents:
8857
diff
changeset
|
179 fsiz = size (fvec); |
|
22c8272af34b
improvements to fsolve & family
Jaroslav Hajek <highegg@gmail.com>
parents:
8857
diff
changeset
|
180 fvec = fvec(:); |
|
22c8272af34b
improvements to fsolve & family
Jaroslav Hajek <highegg@gmail.com>
parents:
8857
diff
changeset
|
181 fn = norm (fvec); |
|
22c8272af34b
improvements to fsolve & family
Jaroslav Hajek <highegg@gmail.com>
parents:
8857
diff
changeset
|
182 m = length (fvec); |
|
22c8272af34b
improvements to fsolve & family
Jaroslav Hajek <highegg@gmail.com>
parents:
8857
diff
changeset
|
183 n = length (x); |
|
22c8272af34b
improvements to fsolve & family
Jaroslav Hajek <highegg@gmail.com>
parents:
8857
diff
changeset
|
184 |
|
22c8272af34b
improvements to fsolve & family
Jaroslav Hajek <highegg@gmail.com>
parents:
8857
diff
changeset
|
185 if (! isempty (outfcn)) |
|
22c8272af34b
improvements to fsolve & family
Jaroslav Hajek <highegg@gmail.com>
parents:
8857
diff
changeset
|
186 optimvalues.iter = niter; |
|
22c8272af34b
improvements to fsolve & family
Jaroslav Hajek <highegg@gmail.com>
parents:
8857
diff
changeset
|
187 optimvalues.funccount = nfev; |
|
22c8272af34b
improvements to fsolve & family
Jaroslav Hajek <highegg@gmail.com>
parents:
8857
diff
changeset
|
188 optimvalues.fval = fn; |
|
22c8272af34b
improvements to fsolve & family
Jaroslav Hajek <highegg@gmail.com>
parents:
8857
diff
changeset
|
189 optimvalues.searchdirection = zeros (n, 1); |
|
22c8272af34b
improvements to fsolve & family
Jaroslav Hajek <highegg@gmail.com>
parents:
8857
diff
changeset
|
190 state = 'init'; |
|
22c8272af34b
improvements to fsolve & family
Jaroslav Hajek <highegg@gmail.com>
parents:
8857
diff
changeset
|
191 stop = outfcn (x, optimvalues, state); |
|
22c8272af34b
improvements to fsolve & family
Jaroslav Hajek <highegg@gmail.com>
parents:
8857
diff
changeset
|
192 if (stop) |
|
22c8272af34b
improvements to fsolve & family
Jaroslav Hajek <highegg@gmail.com>
parents:
8857
diff
changeset
|
193 info = -1; |
|
22c8272af34b
improvements to fsolve & family
Jaroslav Hajek <highegg@gmail.com>
parents:
8857
diff
changeset
|
194 break; |
|
22c8272af34b
improvements to fsolve & family
Jaroslav Hajek <highegg@gmail.com>
parents:
8857
diff
changeset
|
195 endif |
|
22c8272af34b
improvements to fsolve & family
Jaroslav Hajek <highegg@gmail.com>
parents:
8857
diff
changeset
|
196 endif |
|
22c8272af34b
improvements to fsolve & family
Jaroslav Hajek <highegg@gmail.com>
parents:
8857
diff
changeset
|
197 |
|
22c8272af34b
improvements to fsolve & family
Jaroslav Hajek <highegg@gmail.com>
parents:
8857
diff
changeset
|
198 nsuciter = 0; |
|
22c8272af34b
improvements to fsolve & family
Jaroslav Hajek <highegg@gmail.com>
parents:
8857
diff
changeset
|
199 |
| 8507 | 200 ## Outer loop. |
|
8306
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
201 while (niter < maxiter && nfev < maxfev && ! info) |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
202 |
| 8507 | 203 ## Calculate function value and Jacobian (possibly via FD). |
|
8306
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
204 if (has_jac) |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
205 [fvec, fjac] = fcn (reshape (x, xsiz)); |
|
8592
dacfd030633a
handle sparse jacobians in fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
8590
diff
changeset
|
206 ## If the jacobian is sparse, disable Broyden updating. |
|
dacfd030633a
handle sparse jacobians in fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
8590
diff
changeset
|
207 if (issparse (fjac)) |
|
dacfd030633a
handle sparse jacobians in fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
8590
diff
changeset
|
208 updating = false; |
|
dacfd030633a
handle sparse jacobians in fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
8590
diff
changeset
|
209 endif |
|
8986
22c8272af34b
improvements to fsolve & family
Jaroslav Hajek <highegg@gmail.com>
parents:
8857
diff
changeset
|
210 fvec = fvec(:); |
|
8306
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
211 nfev ++; |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
212 else |
|
8986
22c8272af34b
improvements to fsolve & family
Jaroslav Hajek <highegg@gmail.com>
parents:
8857
diff
changeset
|
213 fjac = __fdjac__ (fcn, reshape (x, xsiz), fvec); |
|
22c8272af34b
improvements to fsolve & family
Jaroslav Hajek <highegg@gmail.com>
parents:
8857
diff
changeset
|
214 nfev += length (x); |
|
8306
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
215 endif |
| 8590 | 216 |
|
8763
5ce12bca4c51
update comments in fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
8693
diff
changeset
|
217 ## For square and overdetermined systems, we update a QR |
| 8590 | 218 ## factorization of the jacobian to avoid solving a full system in each |
| 219 ## step. In this case, we pass a triangular matrix to __dogleg__. | |
| 220 useqr = updating && m >= n && n > 10; | |
| 221 | |
| 222 if (useqr) | |
|
8616
3d75d717cbe0
do not pivot by default in fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
8615
diff
changeset
|
223 ## FIXME: Currently, pivoting is mostly useless because the \ operator |
|
3d75d717cbe0
do not pivot by default in fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
8615
diff
changeset
|
224 ## cannot exploit the resulting props of the triangular factor. |
|
3d75d717cbe0
do not pivot by default in fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
8615
diff
changeset
|
225 ## Unpivoted QR is significantly faster so it doesn't seem right to pivot |
|
3d75d717cbe0
do not pivot by default in fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
8615
diff
changeset
|
226 ## just to get invariance. Original MINPACK didn't pivot either, at least |
|
3d75d717cbe0
do not pivot by default in fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
8615
diff
changeset
|
227 ## when qr updating was used. |
|
3d75d717cbe0
do not pivot by default in fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
8615
diff
changeset
|
228 [q, r] = qr (fjac, 0); |
|
8383
a762d9daa700
allow overdetermined systems in fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
8306
diff
changeset
|
229 endif |
|
8306
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
230 |
| 8590 | 231 ## Get column norms, use them as scaling factors. |
|
8306
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
232 jcn = norm (fjac, 'columns').'; |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
233 if (niter == 1) |
|
8990
349d75161672
more cosmetic adjustments to fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
8988
diff
changeset
|
234 dg = jcn; |
|
349d75161672
more cosmetic adjustments to fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
8988
diff
changeset
|
235 dg(dg == 0) = 1; |
|
8306
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
236 xn = norm (dg .* x); |
| 8590 | 237 ## FIXME: something better? |
| 238 delta = max (factor * xn, 1); | |
|
8306
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
239 endif |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
240 |
|
8990
349d75161672
more cosmetic adjustments to fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
8988
diff
changeset
|
241 ## Rescale adaptively. |
|
349d75161672
more cosmetic adjustments to fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
8988
diff
changeset
|
242 ## FIXME: the original minpack used the following rescaling strategy: |
|
349d75161672
more cosmetic adjustments to fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
8988
diff
changeset
|
243 ## dg = max (dg, jcn); |
|
349d75161672
more cosmetic adjustments to fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
8988
diff
changeset
|
244 ## but it seems not good if we start with a bad guess yielding jacobian |
|
349d75161672
more cosmetic adjustments to fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
8988
diff
changeset
|
245 ## columns with large norms that later decrease, because the corresponding |
|
349d75161672
more cosmetic adjustments to fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
8988
diff
changeset
|
246 ## variable will still be overscaled. So instead, we only give the old |
|
349d75161672
more cosmetic adjustments to fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
8988
diff
changeset
|
247 ## scaling a small momentum, but do not honor it. |
|
8613
38482007c834
relax scaling in fsolve and use class-dependent default tolerances
Jaroslav Hajek <highegg@gmail.com>
parents:
8604
diff
changeset
|
248 |
|
8990
349d75161672
more cosmetic adjustments to fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
8988
diff
changeset
|
249 dg = max (0.1*dg, jcn); |
|
8613
38482007c834
relax scaling in fsolve and use class-dependent default tolerances
Jaroslav Hajek <highegg@gmail.com>
parents:
8604
diff
changeset
|
250 |
|
8990
349d75161672
more cosmetic adjustments to fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
8988
diff
changeset
|
251 ## It also seems that in the case of fast (and inhomogeneously) changing |
|
349d75161672
more cosmetic adjustments to fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
8988
diff
changeset
|
252 ## jacobian, the Broyden updates are of little use, so maybe we could |
|
349d75161672
more cosmetic adjustments to fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
8988
diff
changeset
|
253 ## skip them if a big disproportional change is expected. The question is, |
|
349d75161672
more cosmetic adjustments to fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
8988
diff
changeset
|
254 ## of course, how to define the above terms :) |
|
8306
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
255 |
|
8986
22c8272af34b
improvements to fsolve & family
Jaroslav Hajek <highegg@gmail.com>
parents:
8857
diff
changeset
|
256 lastratio = 0; |
|
8306
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
257 nfail = 0; |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
258 nsuc = 0; |
| 8507 | 259 ## Inner loop. |
|
8306
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
260 while (niter <= maxiter && nfev < maxfev && ! info) |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
261 |
| 8590 | 262 ## Get trust-region model (dogleg) minimizer. |
| 263 if (useqr) | |
| 264 qtf = q'*fvec; | |
| 265 s = - __dogleg__ (r, qtf, dg, delta); | |
| 266 w = qtf + r * s; | |
| 267 else | |
| 268 s = - __dogleg__ (fjac, fvec, dg, delta); | |
| 269 w = fvec + fjac * s; | |
|
8395
88fd356b0d95
optionally allow pivoting in fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
8383
diff
changeset
|
270 endif |
| 8590 | 271 |
|
8395
88fd356b0d95
optionally allow pivoting in fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
8383
diff
changeset
|
272 sn = norm (dg .* s); |
|
8306
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
273 if (niter == 1) |
|
8395
88fd356b0d95
optionally allow pivoting in fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
8383
diff
changeset
|
274 delta = min (delta, sn); |
|
8306
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
275 endif |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
276 |
|
8395
88fd356b0d95
optionally allow pivoting in fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
8383
diff
changeset
|
277 fvec1 = fcn (reshape (x + s, xsiz)) (:); |
|
8306
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
278 fn1 = norm (fvec1); |
| 8590 | 279 nfev ++; |
|
8306
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
280 |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
281 if (fn1 < fn) |
| 8507 | 282 ## Scaled actual reduction. |
|
8306
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
283 actred = 1 - (fn1/fn)^2; |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
284 else |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
285 actred = -1; |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
286 endif |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
287 |
| 8507 | 288 ## Scaled predicted reduction, and ratio. |
|
8306
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
289 t = norm (w); |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
290 if (t < fn) |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
291 prered = 1 - (t/fn)^2; |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
292 ratio = actred / prered; |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
293 else |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
294 prered = 0; |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
295 ratio = 0; |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
296 endif |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
297 |
| 8507 | 298 ## Update delta. |
|
8986
22c8272af34b
improvements to fsolve & family
Jaroslav Hajek <highegg@gmail.com>
parents:
8857
diff
changeset
|
299 if (ratio < min(max(0.1, lastratio), 0.9)) |
|
8306
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
300 nsuc = 0; |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
301 nfail ++; |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
302 delta *= 0.5; |
| 8590 | 303 if (delta <= 1e1*macheps*xn) |
| 8507 | 304 ## Trust region became uselessly small. |
|
8306
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
305 info = -3; |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
306 break; |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
307 endif |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
308 else |
|
8986
22c8272af34b
improvements to fsolve & family
Jaroslav Hajek <highegg@gmail.com>
parents:
8857
diff
changeset
|
309 lastratio = ratio; |
|
8306
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
310 nfail = 0; |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
311 nsuc ++; |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
312 if (abs (1-ratio) <= 0.1) |
|
8395
88fd356b0d95
optionally allow pivoting in fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
8383
diff
changeset
|
313 delta = 2*sn; |
|
8306
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
314 elseif (ratio >= 0.5 || nsuc > 1) |
|
8395
88fd356b0d95
optionally allow pivoting in fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
8383
diff
changeset
|
315 delta = max (delta, 2*sn); |
|
8306
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
316 endif |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
317 endif |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
318 |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
319 if (ratio >= 1e-4) |
| 8507 | 320 ## Successful iteration. |
|
8395
88fd356b0d95
optionally allow pivoting in fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
8383
diff
changeset
|
321 x += s; |
|
8306
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
322 xn = norm (dg .* x); |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
323 fvec = fvec1; |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
324 fn = fn1; |
|
8986
22c8272af34b
improvements to fsolve & family
Jaroslav Hajek <highegg@gmail.com>
parents:
8857
diff
changeset
|
325 nsuciter ++; |
|
8306
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
326 endif |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
327 |
| 8857 | 328 niter ++; |
| 329 | |
| 8590 | 330 ## FIXME: should outputfcn be only called after a successful iteration? |
| 8615 | 331 if (! isempty (outfcn)) |
| 8590 | 332 optimvalues.iter = niter; |
| 333 optimvalues.funccount = nfev; | |
| 334 optimvalues.fval = fn; | |
| 335 optimvalues.searchdirection = s; | |
| 336 state = 'iter'; | |
| 337 stop = outfcn (x, optimvalues, state); | |
| 338 if (stop) | |
| 339 info = -1; | |
| 340 break; | |
| 341 endif | |
| 342 endif | |
| 343 | |
| 8466 | 344 ## Tests for termination conditions. A mysterious place, anything |
| 345 ## can happen if you change something here... | |
|
8306
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
346 |
| 8466 | 347 ## The rule of thumb (which I'm not sure M*b is quite following) |
| 348 ## is that for a tolerance that depends on scaling, only 0 makes | |
| 349 ## sense as a default value. But 0 usually means uselessly long | |
| 350 ## iterations, so we need scaling-independent tolerances wherever | |
| 351 ## possible. | |
|
8306
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
352 |
| 8466 | 353 ## FIXME -- why tolf*n*xn? If abs (e) ~ abs(x) * eps is a vector |
| 354 ## of perturbations of x, then norm (fjac*e) <= eps*n*xn, i.e. by | |
| 355 ## tolf ~ eps we demand as much accuracy as we can expect. | |
|
8306
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
356 if (fn <= tolf*n*xn) |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
357 info = 1; |
| 8466 | 358 ## The following tests done only after successful step. |
|
8306
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
359 elseif (actred > 0) |
| 8466 | 360 ## This one is classic. Note that we use scaled variables again, |
| 361 ## but compare to scaled step, so nothing bad. | |
|
8395
88fd356b0d95
optionally allow pivoting in fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
8383
diff
changeset
|
362 if (sn <= tolx*xn) |
|
8306
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
363 info = 2; |
| 8466 | 364 ## Again a classic one. It seems weird to use the same tolf |
| 365 ## for two different tests, but that's what M*b manual appears | |
| 366 ## to say. | |
|
8306
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
367 elseif (actred < tolf) |
|
8383
a762d9daa700
allow overdetermined systems in fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
8306
diff
changeset
|
368 info = 3; |
|
8306
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
369 endif |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
370 endif |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
371 |
| 8507 | 372 ## Criterion for recalculating jacobian. |
|
8986
22c8272af34b
improvements to fsolve & family
Jaroslav Hajek <highegg@gmail.com>
parents:
8857
diff
changeset
|
373 if (! updating || nfail == 2 || nsuciter < 2) |
|
8306
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
374 break; |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
375 endif |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
376 |
| 8507 | 377 ## Compute the scaled Broyden update. |
| 8590 | 378 if (useqr) |
| 379 u = (fvec1 - q*w) / sn; | |
| 380 v = dg .* ((dg .* s) / sn); | |
| 381 | |
| 382 ## Update the QR factorization. | |
| 383 [q, r] = qrupdate (q, r, u, v); | |
| 384 else | |
| 385 u = (fvec1 - w); | |
| 386 v = dg .* ((dg .* s) / sn); | |
| 387 | |
| 388 ## update the jacobian | |
| 389 fjac += u * v'; | |
|
8395
88fd356b0d95
optionally allow pivoting in fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
8383
diff
changeset
|
390 endif |
|
8306
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
391 endwhile |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
392 endwhile |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
393 |
| 8507 | 394 ## Restore original shapes. |
|
8306
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
395 x = reshape (x, xsiz); |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
396 fvec = reshape (fvec, fsiz); |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
397 |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
398 output.iterations = niter; |
|
8986
22c8272af34b
improvements to fsolve & family
Jaroslav Hajek <highegg@gmail.com>
parents:
8857
diff
changeset
|
399 output.successful = nsuciter; |
| 8590 | 400 output.funcCount = nfev; |
|
8306
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
401 |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
402 endfunction |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
403 |
| 8466 | 404 ## An assistant function that evaluates a function handle and checks for |
| 405 ## bad results. | |
|
8693
e5ffb52c9c61
improve fsolve and add ComplexEqn option
Jaroslav Hajek <highegg@gmail.com>
parents:
8661
diff
changeset
|
406 function [fx, jx] = guarded_eval (fun, x, complexeqn) |
|
e5ffb52c9c61
improve fsolve and add ComplexEqn option
Jaroslav Hajek <highegg@gmail.com>
parents:
8661
diff
changeset
|
407 if (nargout > 1) |
|
e5ffb52c9c61
improve fsolve and add ComplexEqn option
Jaroslav Hajek <highegg@gmail.com>
parents:
8661
diff
changeset
|
408 [fx, jx] = fun (x); |
|
e5ffb52c9c61
improve fsolve and add ComplexEqn option
Jaroslav Hajek <highegg@gmail.com>
parents:
8661
diff
changeset
|
409 else |
|
e5ffb52c9c61
improve fsolve and add ComplexEqn option
Jaroslav Hajek <highegg@gmail.com>
parents:
8661
diff
changeset
|
410 fx = fun (x); |
|
8997
187a9d9c2f04
simplifications to fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
8991
diff
changeset
|
411 jx = []; |
|
8693
e5ffb52c9c61
improve fsolve and add ComplexEqn option
Jaroslav Hajek <highegg@gmail.com>
parents:
8661
diff
changeset
|
412 endif |
|
e5ffb52c9c61
improve fsolve and add ComplexEqn option
Jaroslav Hajek <highegg@gmail.com>
parents:
8661
diff
changeset
|
413 |
|
8997
187a9d9c2f04
simplifications to fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
8991
diff
changeset
|
414 if (! complexeqn && ! (isreal (fx) && isreal (jx))) |
|
8468
866492035ecf
fsolve.m, fzero.m: undo part of previous change
John W. Eaton <jwe@octave.org>
parents:
8467
diff
changeset
|
415 error ("fsolve:notreal", "fsolve: non-real value encountered"); |
|
8997
187a9d9c2f04
simplifications to fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
8991
diff
changeset
|
416 elseif (complexeqn && ! (isnumeric (fx) && isnumeric(jx))) |
|
8693
e5ffb52c9c61
improve fsolve and add ComplexEqn option
Jaroslav Hajek <highegg@gmail.com>
parents:
8661
diff
changeset
|
417 error ("fsolve:notnum", "fsolve: non-numeric value encountered"); |
|
e5ffb52c9c61
improve fsolve and add ComplexEqn option
Jaroslav Hajek <highegg@gmail.com>
parents:
8661
diff
changeset
|
418 elseif (any (isnan (fx(:)))) |
|
8468
866492035ecf
fsolve.m, fzero.m: undo part of previous change
John W. Eaton <jwe@octave.org>
parents:
8467
diff
changeset
|
419 error ("fsolve:isnan", "fsolve: NaN value encountered"); |
|
8306
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
420 endif |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
421 endfunction |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
422 |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
423 %!function retval = f (p) |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
424 %! x = p(1); |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
425 %! y = p(2); |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
426 %! z = p(3); |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
427 %! retval = zeros (3, 1); |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
428 %! retval(1) = sin(x) + y**2 + log(z) - 7; |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
429 %! retval(2) = 3*x + 2**y -z**3 + 1; |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
430 %! retval(3) = x + y + z - 5; |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
431 %!test |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
432 %! x_opt = [ 0.599054; |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
433 %! 2.395931; |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
434 %! 2.005014 ]; |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
435 %! tol = 1.0e-5; |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
436 %! [x, fval, info] = fsolve (@f, [ 0.5; 2.0; 2.5 ]); |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
437 %! assert (info > 0); |
|
8383
a762d9daa700
allow overdetermined systems in fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
8306
diff
changeset
|
438 %! assert (norm (x - x_opt, Inf) < tol); |
|
8306
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
439 %! assert (norm (fval) < tol); |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
440 |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
441 %!function retval = f (p) |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
442 %! x = p(1); |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
443 %! y = p(2); |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
444 %! z = p(3); |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
445 %! w = p(4); |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
446 %! retval = zeros (4, 1); |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
447 %! retval(1) = 3*x + 4*y + exp (z + w) - 1.007; |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
448 %! retval(2) = 6*x - 4*y + exp (3*z + w) - 11; |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
449 %! retval(3) = x^4 - 4*y^2 + 6*z - 8*w - 20; |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
450 %! retval(4) = x^2 + 2*y^3 + z - w - 4; |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
451 %!test |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
452 %! x_opt = [ -0.767297326653401, 0.590671081117440, 1.47190018629642, -1.52719341133957 ]; |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
453 %! tol = 1.0e-5; |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
454 %! [x, fval, info] = fsolve (@f, [-1, 1, 2, -1]); |
|
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
455 %! assert (info > 0); |
|
8383
a762d9daa700
allow overdetermined systems in fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
8306
diff
changeset
|
456 %! assert (norm (x - x_opt, Inf) < tol); |
|
8306
43795cf108d0
initial implementation of fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
diff
changeset
|
457 %! assert (norm (fval) < tol); |
|
8383
a762d9daa700
allow overdetermined systems in fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
8306
diff
changeset
|
458 |
|
a762d9daa700
allow overdetermined systems in fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
8306
diff
changeset
|
459 %!function retval = f (p) |
|
a762d9daa700
allow overdetermined systems in fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
8306
diff
changeset
|
460 %! x = p(1); |
|
a762d9daa700
allow overdetermined systems in fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
8306
diff
changeset
|
461 %! y = p(2); |
|
a762d9daa700
allow overdetermined systems in fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
8306
diff
changeset
|
462 %! z = p(3); |
|
a762d9daa700
allow overdetermined systems in fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
8306
diff
changeset
|
463 %! retval = zeros (3, 1); |
|
a762d9daa700
allow overdetermined systems in fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
8306
diff
changeset
|
464 %! retval(1) = sin(x) + y**2 + log(z) - 7; |
|
a762d9daa700
allow overdetermined systems in fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
8306
diff
changeset
|
465 %! retval(2) = 3*x + 2**y -z**3 + 1; |
|
a762d9daa700
allow overdetermined systems in fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
8306
diff
changeset
|
466 %! retval(3) = x + y + z - 5; |
|
a762d9daa700
allow overdetermined systems in fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
8306
diff
changeset
|
467 %! retval(4) = x*x + y - z*log(z) - 1.36; |
|
a762d9daa700
allow overdetermined systems in fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
8306
diff
changeset
|
468 %!test |
|
a762d9daa700
allow overdetermined systems in fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
8306
diff
changeset
|
469 %! x_opt = [ 0.599054; |
|
a762d9daa700
allow overdetermined systems in fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
8306
diff
changeset
|
470 %! 2.395931; |
|
a762d9daa700
allow overdetermined systems in fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
8306
diff
changeset
|
471 %! 2.005014 ]; |
|
a762d9daa700
allow overdetermined systems in fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
8306
diff
changeset
|
472 %! tol = 1.0e-5; |
| 8590 | 473 %! [x, fval, info] = fsolve (@f, [ 0.5; 2.0; 2.5 ]); |
|
8383
a762d9daa700
allow overdetermined systems in fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
8306
diff
changeset
|
474 %! assert (info > 0); |
|
a762d9daa700
allow overdetermined systems in fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
8306
diff
changeset
|
475 %! assert (norm (x - x_opt, Inf) < tol); |
|
a762d9daa700
allow overdetermined systems in fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
8306
diff
changeset
|
476 %! assert (norm (fval) < tol); |
|
a762d9daa700
allow overdetermined systems in fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
8306
diff
changeset
|
477 |
|
8395
88fd356b0d95
optionally allow pivoting in fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
8383
diff
changeset
|
478 %!function retval = f (p) |
|
88fd356b0d95
optionally allow pivoting in fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
8383
diff
changeset
|
479 %! x = p(1); |
|
88fd356b0d95
optionally allow pivoting in fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
8383
diff
changeset
|
480 %! y = p(2); |
|
88fd356b0d95
optionally allow pivoting in fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
8383
diff
changeset
|
481 %! z = p(3); |
|
88fd356b0d95
optionally allow pivoting in fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
8383
diff
changeset
|
482 %! retval = zeros (3, 1); |
|
88fd356b0d95
optionally allow pivoting in fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
8383
diff
changeset
|
483 %! retval(1) = sin(x) + y**2 + log(z) - 7; |
|
88fd356b0d95
optionally allow pivoting in fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
8383
diff
changeset
|
484 %! retval(2) = 3*x + 2**y -z**3 + 1; |
|
88fd356b0d95
optionally allow pivoting in fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
8383
diff
changeset
|
485 %! retval(3) = x + y + z - 5; |
|
88fd356b0d95
optionally allow pivoting in fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
8383
diff
changeset
|
486 %!test |
|
88fd356b0d95
optionally allow pivoting in fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
8383
diff
changeset
|
487 %! x_opt = [ 0.599054; |
|
88fd356b0d95
optionally allow pivoting in fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
8383
diff
changeset
|
488 %! 2.395931; |
|
88fd356b0d95
optionally allow pivoting in fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
8383
diff
changeset
|
489 %! 2.005014 ]; |
|
88fd356b0d95
optionally allow pivoting in fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
8383
diff
changeset
|
490 %! tol = 1.0e-5; |
| 8590 | 491 %! opt = optimset ("Updating", "qrp"); |
|
8395
88fd356b0d95
optionally allow pivoting in fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
8383
diff
changeset
|
492 %! [x, fval, info] = fsolve (@f, [ 0.5; 2.0; 2.5 ], opt); |
|
88fd356b0d95
optionally allow pivoting in fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
8383
diff
changeset
|
493 %! assert (info > 0); |
|
88fd356b0d95
optionally allow pivoting in fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
8383
diff
changeset
|
494 %! assert (norm (x - x_opt, Inf) < tol); |
|
88fd356b0d95
optionally allow pivoting in fsolve
Jaroslav Hajek <highegg@gmail.com>
parents:
8383
diff
changeset
|
495 %! assert (norm (fval) < tol); |
| 8590 | 496 |
| 497 %!test | |
| 498 %! b0 = 3; | |
| 499 %! a0 = 0.2; | |
| 500 %! x = 0:.5:5; | |
| 501 %! noise = 1e-5 * sin (100*x); | |
| 502 %! y = exp (-a0*x) + b0 + noise; | |
| 503 %! c_opt = [a0, b0]; | |
| 504 %! tol = 1e-5; | |
| 505 %! | |
| 506 %! [c, fval, info, output] = fsolve (@(c) (exp(-c(1)*x) + c(2) - y), [0, 0]); | |
| 507 %! assert (info > 0); | |
| 508 %! assert (norm (c - c_opt, Inf) < tol); | |
| 509 %! assert (norm (fval) < norm (noise)); | |
| 510 | |
|
8693
e5ffb52c9c61
improve fsolve and add ComplexEqn option
Jaroslav Hajek <highegg@gmail.com>
parents:
8661
diff
changeset
|
511 |
|
e5ffb52c9c61
improve fsolve and add ComplexEqn option
Jaroslav Hajek <highegg@gmail.com>
parents:
8661
diff
changeset
|
512 %!function y = cfun (x) |
|
e5ffb52c9c61
improve fsolve and add ComplexEqn option
Jaroslav Hajek <highegg@gmail.com>
parents:
8661
diff
changeset
|
513 %! y(1) = (1+i)*x(1)^2 - (1-i)*x(2) - 2; |
|
e5ffb52c9c61
improve fsolve and add ComplexEqn option
Jaroslav Hajek <highegg@gmail.com>
parents:
8661
diff
changeset
|
514 %! y(2) = sqrt (x(1)*x(2)) - (1-2i)*x(3) + (3-4i); |
|
e5ffb52c9c61
improve fsolve and add ComplexEqn option
Jaroslav Hajek <highegg@gmail.com>
parents:
8661
diff
changeset
|
515 %! y(3) = x(1) * x(2) - x(3)^2 + (3+2i); |
|
e5ffb52c9c61
improve fsolve and add ComplexEqn option
Jaroslav Hajek <highegg@gmail.com>
parents:
8661
diff
changeset
|
516 |
|
e5ffb52c9c61
improve fsolve and add ComplexEqn option
Jaroslav Hajek <highegg@gmail.com>
parents:
8661
diff
changeset
|
517 %!test |
|
e5ffb52c9c61
improve fsolve and add ComplexEqn option
Jaroslav Hajek <highegg@gmail.com>
parents:
8661
diff
changeset
|
518 %! x_opt = [-1+i, 1-i, 2+i]; |
|
e5ffb52c9c61
improve fsolve and add ComplexEqn option
Jaroslav Hajek <highegg@gmail.com>
parents:
8661
diff
changeset
|
519 %! x = [i, 1, 1+i]; |
|
e5ffb52c9c61
improve fsolve and add ComplexEqn option
Jaroslav Hajek <highegg@gmail.com>
parents:
8661
diff
changeset
|
520 %! |
|
e5ffb52c9c61
improve fsolve and add ComplexEqn option
Jaroslav Hajek <highegg@gmail.com>
parents:
8661
diff
changeset
|
521 %! [x, f, info] = fsolve (@cfun, x, optimset ("ComplexEqn", "on")); |
|
e5ffb52c9c61
improve fsolve and add ComplexEqn option
Jaroslav Hajek <highegg@gmail.com>
parents:
8661
diff
changeset
|
522 %! tol = 1e-5; |
|
e5ffb52c9c61
improve fsolve and add ComplexEqn option
Jaroslav Hajek <highegg@gmail.com>
parents:
8661
diff
changeset
|
523 %! assert (norm (f) < tol); |
|
e5ffb52c9c61
improve fsolve and add ComplexEqn option
Jaroslav Hajek <highegg@gmail.com>
parents:
8661
diff
changeset
|
524 %! assert (norm (x - x_opt, Inf) < tol); |
|
e5ffb52c9c61
improve fsolve and add ComplexEqn option
Jaroslav Hajek <highegg@gmail.com>
parents:
8661
diff
changeset
|
525 |