octave-nkf: scripts/optimization/fsolve.m annotate

annotate scripts/optimization/fsolve.m @ 10687:a8ce6bdecce5

Improve documentation strings.

author	Rik <octave@nomad.inbox5.com>
date	Tue, 08 Jun 2010 20:22:38 -0700
parents	95c3e38098bf
children	fbd7843974fa

rev	line source
8590 c136d313206a improvements to fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	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 4d008d9f0ccf fsolve.m: style fixes John W. Eaton <jwe@octave.org> parents: 8395 diff changeset	21 ## -- texinfo --
9143 74d5c1a4ca96 Eliminate 'unbalanced parentheses in @def...' error during texi2pdf. Rik <rdrider0-list@yahoo.com> parents: 9075 diff changeset	22 ## @deftypefn {Function File} {} fsolve (@var{fcn}, @var{x0}, @var{options})
9724 f22bbc5d56e9 Fix various incorrect usages of TeXinfo deffn and deftypefn macros Rik <rdrider0-list@yahoo.com> parents: 9628 diff changeset	23 ## @deftypefnx {Function File} {[@var{x}, @var{fvec}, @var{info}, @var{output}, @var{fjac}] =} fsolve (@var{fcn}, @dots{})
8513 352d3245d4c1 fsolve.m: doc fix John W. Eaton <jwe@octave.org> parents: 8507 diff changeset	24 ## Solve a system of nonlinear equations defined by the function @var{fcn}.
8466 4d008d9f0ccf fsolve.m: style fixes John W. Eaton <jwe@octave.org> parents: 8395 diff changeset	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 4d008d9f0ccf fsolve.m: style fixes John W. Eaton <jwe@octave.org> parents: 8395 diff changeset	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 4d008d9f0ccf fsolve.m: style fixes John W. Eaton <jwe@octave.org> parents: 8395 diff changeset	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 c136d313206a improvements to fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	35 ## @code{"TolFun"}, @code{"MaxIter"}, @code{"MaxFunEvals"},
10201 5c66978f3fdf support TypicalX and AutoScaling in fsolve/fminunc, don't autoscale by default Jaroslav Hajek <highegg@gmail.com> parents: 10050 diff changeset	36 ## @code{"Jacobian"}, @code{"Updating"}, @code{"ComplexEqn"}
5c66978f3fdf support TypicalX and AutoScaling in fsolve/fminunc, don't autoscale by default Jaroslav Hajek <highegg@gmail.com> parents: 10050 diff changeset	37 ## @code{"TypicalX"}, @code{"AutoScaling"} and @code{"FinDiffType"}.
8466 4d008d9f0ccf fsolve.m: style fixes John W. Eaton <jwe@octave.org> parents: 8395 diff changeset	38 ##
8513 352d3245d4c1 fsolve.m: doc fix John W. Eaton <jwe@octave.org> parents: 8507 diff changeset	39 ## If @code{"Jacobian"} is @code{"on"}, it specifies that @var{fcn},
352d3245d4c1 fsolve.m: doc fix John W. Eaton <jwe@octave.org> parents: 8507 diff changeset	40 ## called with 2 output arguments, also returns the Jacobian matrix
352d3245d4c1 fsolve.m: doc fix John W. Eaton <jwe@octave.org> parents: 8507 diff changeset	41 ## of right-hand sides at the requested point. @code{"TolX"} specifies
352d3245d4c1 fsolve.m: doc fix John W. Eaton <jwe@octave.org> parents: 8507 diff changeset	42 ## 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	43 ## @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	44 ## for both @code{"TolX"} and @code{"TolFun"}.
10201 5c66978f3fdf support TypicalX and AutoScaling in fsolve/fminunc, don't autoscale by default Jaroslav Hajek <highegg@gmail.com> parents: 10050 diff changeset	45 ##
5c66978f3fdf support TypicalX and AutoScaling in fsolve/fminunc, don't autoscale by default Jaroslav Hajek <highegg@gmail.com> parents: 10050 diff changeset	46 ## If @code{"AutoScaling"} is on, the variables will be automatically scaled
5c66978f3fdf support TypicalX and AutoScaling in fsolve/fminunc, don't autoscale by default Jaroslav Hajek <highegg@gmail.com> parents: 10050 diff changeset	47 ## according to the column norms of the (estimated) Jacobian. As a result,
5c66978f3fdf support TypicalX and AutoScaling in fsolve/fminunc, don't autoscale by default Jaroslav Hajek <highegg@gmail.com> parents: 10050 diff changeset	48 ## TolF becomes scaling-independent. By default, this option is off, because
5c66978f3fdf support TypicalX and AutoScaling in fsolve/fminunc, don't autoscale by default Jaroslav Hajek <highegg@gmail.com> parents: 10050 diff changeset	49 ## it may sometimes deliver unexpected (though mathematically correct) results.
5c66978f3fdf support TypicalX and AutoScaling in fsolve/fminunc, don't autoscale by default Jaroslav Hajek <highegg@gmail.com> parents: 10050 diff changeset	50 ##
8596 8833c0b18eb2 enable default settings queries in optim funcs Jaroslav Hajek <highegg@gmail.com> parents: 8592 diff changeset	51 ## If @code{"Updating"} is "on", the function will attempt to use Broyden
8590 c136d313206a improvements to fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	52 ## updates to update the Jacobian, in order to reduce the amount of jacobian
c136d313206a improvements to fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	53 ## calculations.
c136d313206a improvements to fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	54 ## If your user function always calculates the Jacobian (regardless of number
c136d313206a improvements to fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	55 ## of output arguments), this option provides no advantage and should be set to
c136d313206a improvements to fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	56 ## false.
c136d313206a improvements to fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	57 ##
8693 e5ffb52c9c61 improve fsolve and add ComplexEqn option Jaroslav Hajek <highegg@gmail.com> parents: 8661 diff changeset	58 ## @code{"ComplexEqn"} is @code{"on"}, @code{fsolve} will attempt to solve
9065 8207b833557f Cleanup documentation for arith.texi, linalg.texi, nonlin.texi Rik <rdrider0-list@yahoo.com> parents: 9051 diff changeset	59 ## complex equations in complex variables, assuming that the equations possess a
9051 1bf0ce0930be Grammar check TexInfo in all .m files Rik <rdrider0-list@yahoo.com> parents: 8997 diff changeset	60 ## 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	61 ## 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	62 ## system.
e5ffb52c9c61 improve fsolve and add ComplexEqn option Jaroslav Hajek <highegg@gmail.com> parents: 8661 diff changeset	63 ##
8466 4d008d9f0ccf fsolve.m: style fixes John W. Eaton <jwe@octave.org> parents: 8395 diff changeset	64 ## For description of the other options, see @code{optimset}.
4d008d9f0ccf fsolve.m: style fixes John W. Eaton <jwe@octave.org> parents: 8395 diff changeset	65 ##
4d008d9f0ccf fsolve.m: style fixes John W. Eaton <jwe@octave.org> parents: 8395 diff changeset	66 ## On return, @var{fval} contains the value of the function @var{fcn}
4d008d9f0ccf fsolve.m: style fixes John W. Eaton <jwe@octave.org> parents: 8395 diff changeset	67 ## evaluated at @var{x}, and @var{info} may be one of the following values:
4d008d9f0ccf fsolve.m: style fixes John W. Eaton <jwe@octave.org> parents: 8395 diff changeset	68 ##
4d008d9f0ccf fsolve.m: style fixes John W. Eaton <jwe@octave.org> parents: 8395 diff changeset	69 ## @table @asis
4d008d9f0ccf fsolve.m: style fixes John W. Eaton <jwe@octave.org> parents: 8395 diff changeset	70 ## @item 1
9051 1bf0ce0930be Grammar check TexInfo in all .m files Rik <rdrider0-list@yahoo.com> parents: 8997 diff changeset	71 ## Converged to a solution point. Relative residual error is less than specified
8466 4d008d9f0ccf fsolve.m: style fixes John W. Eaton <jwe@octave.org> parents: 8395 diff changeset	72 ## by TolFun.
4d008d9f0ccf fsolve.m: style fixes John W. Eaton <jwe@octave.org> parents: 8395 diff changeset	73 ## @item 2
4d008d9f0ccf fsolve.m: style fixes John W. Eaton <jwe@octave.org> parents: 8395 diff changeset	74 ## Last relative step size was less that TolX.
4d008d9f0ccf fsolve.m: style fixes John W. Eaton <jwe@octave.org> parents: 8395 diff changeset	75 ## @item 3
4d008d9f0ccf fsolve.m: style fixes John W. Eaton <jwe@octave.org> parents: 8395 diff changeset	76 ## Last relative decrease in residual was less than TolF.
4d008d9f0ccf fsolve.m: style fixes John W. Eaton <jwe@octave.org> parents: 8395 diff changeset	77 ## @item 0
4d008d9f0ccf fsolve.m: style fixes John W. Eaton <jwe@octave.org> parents: 8395 diff changeset	78 ## Iteration limit exceeded.
4d008d9f0ccf fsolve.m: style fixes John W. Eaton <jwe@octave.org> parents: 8395 diff changeset	79 ## @item -3
4d008d9f0ccf fsolve.m: style fixes John W. Eaton <jwe@octave.org> parents: 8395 diff changeset	80 ## The trust region radius became excessively small.
4d008d9f0ccf fsolve.m: style fixes John W. Eaton <jwe@octave.org> parents: 8395 diff changeset	81 ## @end table
8986 22c8272af34b improvements to fsolve & family Jaroslav Hajek <highegg@gmail.com> parents: 8857 diff changeset	82 ##
8466 4d008d9f0ccf fsolve.m: style fixes John W. Eaton <jwe@octave.org> parents: 8395 diff changeset	83 ## Note: If you only have a single nonlinear equation of one variable, using
4d008d9f0ccf fsolve.m: style fixes John W. Eaton <jwe@octave.org> parents: 8395 diff changeset	84 ## @code{fzero} is usually a much better idea.
8513 352d3245d4c1 fsolve.m: doc fix John W. Eaton <jwe@octave.org> parents: 8507 diff changeset	85 ## @seealso{fzero, optimset}
8986 22c8272af34b improvements to fsolve & family Jaroslav Hajek <highegg@gmail.com> parents: 8857 diff changeset	86 ##
22c8272af34b improvements to fsolve & family Jaroslav Hajek <highegg@gmail.com> parents: 8857 diff changeset	87 ## Note about user-supplied jacobians:
22c8272af34b improvements to fsolve & family Jaroslav Hajek <highegg@gmail.com> parents: 8857 diff changeset	88 ## 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	89 ## 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	90 ## 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	91 ## 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	92 ## 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	93 ## 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	94 ## 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	95 ## 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	96 ## 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	97 ## 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	98 ## example how this can be achieved follows:
8988 315828058e0d minor doc fixes in fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8986 diff changeset	99 ##
8986 22c8272af34b improvements to fsolve & family Jaroslav Hajek <highegg@gmail.com> parents: 8857 diff changeset	100 ## @example
9153 5247e89688e1 Eliminate most overfull errors when running texi2pdf for generating pdf documentation Rik <rdrider0-list@yahoo.com> parents: 9143 diff changeset	101 ## function [fvec, fjac] = user_func (x, optimvalues, state)
8986 22c8272af34b improvements to fsolve & family Jaroslav Hajek <highegg@gmail.com> parents: 8857 diff changeset	102 ## persistent sav = [], sav0 = [];
22c8272af34b improvements to fsolve & family Jaroslav Hajek <highegg@gmail.com> parents: 8857 diff changeset	103 ## if (nargin == 1)
22c8272af34b improvements to fsolve & family Jaroslav Hajek <highegg@gmail.com> parents: 8857 diff changeset	104 ## ## evaluation call
22c8272af34b improvements to fsolve & family Jaroslav Hajek <highegg@gmail.com> parents: 8857 diff changeset	105 ## if (nargout == 1)
22c8272af34b improvements to fsolve & family Jaroslav Hajek <highegg@gmail.com> parents: 8857 diff changeset	106 ## sav0.x = x; # mark saved vector
22c8272af34b improvements to fsolve & family Jaroslav Hajek <highegg@gmail.com> parents: 8857 diff changeset	107 ## ## calculate fvec, save results to sav0.
22c8272af34b improvements to fsolve & family Jaroslav Hajek <highegg@gmail.com> parents: 8857 diff changeset	108 ## elseif (nargout == 2)
22c8272af34b improvements to fsolve & family Jaroslav Hajek <highegg@gmail.com> parents: 8857 diff changeset	109 ## ## calculate fjac using sav.
22c8272af34b improvements to fsolve & family Jaroslav Hajek <highegg@gmail.com> parents: 8857 diff changeset	110 ## endif
22c8272af34b improvements to fsolve & family Jaroslav Hajek <highegg@gmail.com> parents: 8857 diff changeset	111 ## else
22c8272af34b improvements to fsolve & family Jaroslav Hajek <highegg@gmail.com> parents: 8857 diff changeset	112 ## ## outputfcn call.
22c8272af34b improvements to fsolve & family Jaroslav Hajek <highegg@gmail.com> parents: 8857 diff changeset	113 ## if (all (x == sav0.x))
22c8272af34b improvements to fsolve & family Jaroslav Hajek <highegg@gmail.com> parents: 8857 diff changeset	114 ## sav = sav0;
22c8272af34b improvements to fsolve & family Jaroslav Hajek <highegg@gmail.com> parents: 8857 diff changeset	115 ## endif
9209 923c7cb7f13f Simplify TeXinfo files by eliminating redundant @iftex followed by @tex construction. Rik <rdrider0-list@yahoo.com> parents: 9207 diff changeset	116 ## ## maybe output iteration status, etc.
8986 22c8272af34b improvements to fsolve & family Jaroslav Hajek <highegg@gmail.com> parents: 8857 diff changeset	117 ## endif
22c8272af34b improvements to fsolve & family Jaroslav Hajek <highegg@gmail.com> parents: 8857 diff changeset	118 ## endfunction
22c8272af34b improvements to fsolve & family Jaroslav Hajek <highegg@gmail.com> parents: 8857 diff changeset	119 ##
9051 1bf0ce0930be Grammar check TexInfo in all .m files Rik <rdrider0-list@yahoo.com> parents: 8997 diff changeset	120 ## ## @dots{}.
8986 22c8272af34b improvements to fsolve & family Jaroslav Hajek <highegg@gmail.com> parents: 8857 diff changeset	121 ##
9153 5247e89688e1 Eliminate most overfull errors when running texi2pdf for generating pdf documentation Rik <rdrider0-list@yahoo.com> parents: 9143 diff changeset	122 ## fsolve (@@user_func, x0, optimset ("OutputFcn", @@user_func, @dots{}))
8986 22c8272af34b improvements to fsolve & family Jaroslav Hajek <highegg@gmail.com> parents: 8857 diff changeset	123 ## @end example
22c8272af34b improvements to fsolve & family Jaroslav Hajek <highegg@gmail.com> parents: 8857 diff changeset	124 ###
8466 4d008d9f0ccf fsolve.m: style fixes John W. Eaton <jwe@octave.org> parents: 8395 diff changeset	125 ## @end deftypefn
8306 43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	126
8648 ff61b53eb294 optimization: use PKG_ADD: comments instead of PKG_ADD file John W. Eaton <jwe@octave.org> parents: 8647 diff changeset	127 ## 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	128
8596 8833c0b18eb2 enable default settings queries in optim funcs Jaroslav Hajek <highegg@gmail.com> parents: 8592 diff changeset	129 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	130
8596 8833c0b18eb2 enable default settings queries in optim funcs Jaroslav Hajek <highegg@gmail.com> parents: 8592 diff changeset	131 ## Get default options if requested.
8833c0b18eb2 enable default settings queries in optim funcs Jaroslav Hajek <highegg@gmail.com> parents: 8592 diff changeset	132 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	133 x = optimset ("MaxIter", 400, "MaxFunEvals", Inf, \
10201 5c66978f3fdf support TypicalX and AutoScaling in fsolve/fminunc, don't autoscale by default Jaroslav Hajek <highegg@gmail.com> parents: 10050 diff changeset	134 "Jacobian", "off", "TolX", 1e-7, "TolFun", 1e-7,
8693 e5ffb52c9c61 improve fsolve and add ComplexEqn option Jaroslav Hajek <highegg@gmail.com> parents: 8661 diff changeset	135 "OutputFcn", [], "Updating", "on", "FunValCheck", "off",
10201 5c66978f3fdf support TypicalX and AutoScaling in fsolve/fminunc, don't autoscale by default Jaroslav Hajek <highegg@gmail.com> parents: 10050 diff changeset	136 "ComplexEqn", "off", "FinDiffType", "central",
5c66978f3fdf support TypicalX and AutoScaling in fsolve/fminunc, don't autoscale by default Jaroslav Hajek <highegg@gmail.com> parents: 10050 diff changeset	137 "TypicalX", [], "AutoScaling", "off");
8596 8833c0b18eb2 enable default settings queries in optim funcs Jaroslav Hajek <highegg@gmail.com> parents: 8592 diff changeset	138 return;
8833c0b18eb2 enable default settings queries in optim funcs Jaroslav Hajek <highegg@gmail.com> parents: 8592 diff changeset	139 endif
8833c0b18eb2 enable default settings queries in optim funcs Jaroslav Hajek <highegg@gmail.com> parents: 8592 diff changeset	140
8833c0b18eb2 enable default settings queries in optim funcs Jaroslav Hajek <highegg@gmail.com> parents: 8592 diff changeset	141 if (nargin < 2 \|\| nargin > 3 \|\| ! ismatrix (x0))
8833c0b18eb2 enable default settings queries in optim funcs Jaroslav Hajek <highegg@gmail.com> parents: 8592 diff changeset	142 print_usage ();
8833c0b18eb2 enable default settings queries in optim funcs Jaroslav Hajek <highegg@gmail.com> parents: 8592 diff changeset	143 endif
8833c0b18eb2 enable default settings queries in optim funcs Jaroslav Hajek <highegg@gmail.com> parents: 8592 diff changeset	144
8833c0b18eb2 enable default settings queries in optim funcs Jaroslav Hajek <highegg@gmail.com> parents: 8592 diff changeset	145 if (ischar (fcn))
9464 e598248a060d safer str2func use in optim functions Jaroslav Hajek <highegg@gmail.com> parents: 9212 diff changeset	146 fcn = str2func (fcn, "global");
10050 dc88a0b6472c support old style jacobian for fsolve Jaroslav Hajek <highegg@gmail.com> parents: 9899 diff changeset	147 elseif (iscell (fcn))
dc88a0b6472c support old style jacobian for fsolve Jaroslav Hajek <highegg@gmail.com> parents: 9899 diff changeset	148 fcn = @(x) make_fcn_jac (x, fcn{1}, fcn{2});
8306 43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	149 endif
43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	150
43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	151 xsiz = size (x0);
43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	152 n = numel (x0);
43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	153
8467 77b8d4aa2743 fsolve.m, fzero.m: style fixes; use strcmpi to compare options John W. Eaton <jwe@octave.org> parents: 8466 diff changeset	154 has_jac = strcmpi (optimget (options, "Jacobian", "off"), "on");
9212 6feb27c38da1 support central differences in fminunc and fsolve Jaroslav Hajek <highegg@gmail.com> parents: 9209 diff changeset	155 cdif = strcmpi (optimget (options, "FinDiffType", "central"), "central");
8596 8833c0b18eb2 enable default settings queries in optim funcs Jaroslav Hajek <highegg@gmail.com> parents: 8592 diff changeset	156 maxiter = optimget (options, "MaxIter", 400);
8306 43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	157 maxfev = optimget (options, "MaxFunEvals", Inf);
43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	158 outfcn = optimget (options, "OutputFcn");
8596 8833c0b18eb2 enable default settings queries in optim funcs Jaroslav Hajek <highegg@gmail.com> parents: 8592 diff changeset	159 updating = strcmpi (optimget (options, "Updating", "on"), "on");
8693 e5ffb52c9c61 improve fsolve and add ComplexEqn option Jaroslav Hajek <highegg@gmail.com> parents: 8661 diff changeset	160 complexeqn = strcmpi (optimget (options, "ComplexEqn", "off"), "on");
8590 c136d313206a improvements to fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	161
10201 5c66978f3fdf support TypicalX and AutoScaling in fsolve/fminunc, don't autoscale by default Jaroslav Hajek <highegg@gmail.com> parents: 10050 diff changeset	162 ## Get scaling matrix using the TypicalX option. If set to "auto", the
5c66978f3fdf support TypicalX and AutoScaling in fsolve/fminunc, don't autoscale by default Jaroslav Hajek <highegg@gmail.com> parents: 10050 diff changeset	163 ## scaling matrix is estimated using the jacobian.
5c66978f3fdf support TypicalX and AutoScaling in fsolve/fminunc, don't autoscale by default Jaroslav Hajek <highegg@gmail.com> parents: 10050 diff changeset	164 typicalx = optimget (options, "TypicalX");
5c66978f3fdf support TypicalX and AutoScaling in fsolve/fminunc, don't autoscale by default Jaroslav Hajek <highegg@gmail.com> parents: 10050 diff changeset	165 if (isempty (typicalx))
5c66978f3fdf support TypicalX and AutoScaling in fsolve/fminunc, don't autoscale by default Jaroslav Hajek <highegg@gmail.com> parents: 10050 diff changeset	166 typicalx = ones (n, 1);
5c66978f3fdf support TypicalX and AutoScaling in fsolve/fminunc, don't autoscale by default Jaroslav Hajek <highegg@gmail.com> parents: 10050 diff changeset	167 endif
5c66978f3fdf support TypicalX and AutoScaling in fsolve/fminunc, don't autoscale by default Jaroslav Hajek <highegg@gmail.com> parents: 10050 diff changeset	168 autoscale = strcmpi (optimget (options, "AutoScaling", "off"), "on");
5c66978f3fdf support TypicalX and AutoScaling in fsolve/fminunc, don't autoscale by default Jaroslav Hajek <highegg@gmail.com> parents: 10050 diff changeset	169 if (! autoscale)
5c66978f3fdf support TypicalX and AutoScaling in fsolve/fminunc, don't autoscale by default Jaroslav Hajek <highegg@gmail.com> parents: 10050 diff changeset	170 dg = 1 ./ typicalx;
5c66978f3fdf support TypicalX and AutoScaling in fsolve/fminunc, don't autoscale by default Jaroslav Hajek <highegg@gmail.com> parents: 10050 diff changeset	171 endif
5c66978f3fdf support TypicalX and AutoScaling in fsolve/fminunc, don't autoscale by default Jaroslav Hajek <highegg@gmail.com> parents: 10050 diff changeset	172
8467 77b8d4aa2743 fsolve.m, fzero.m: style fixes; use strcmpi to compare options John W. Eaton <jwe@octave.org> parents: 8466 diff changeset	173 funvalchk = strcmpi (optimget (options, "FunValCheck", "off"), "on");
8306 43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	174
43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	175 if (funvalchk)
8604 43f831758d42 fix typo in fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8596 diff changeset	176 ## Replace fcn with a guarded version.
8693 e5ffb52c9c61 improve fsolve and add ComplexEqn option Jaroslav Hajek <highegg@gmail.com> parents: 8661 diff changeset	177 fcn = @(x) guarded_eval (fcn, x, complexeqn);
8306 43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	178 endif
43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	179
8466 4d008d9f0ccf fsolve.m: style fixes John W. Eaton <jwe@octave.org> parents: 8395 diff changeset	180 ## These defaults are rather stringent. I think that normally, user
4d008d9f0ccf fsolve.m: style fixes John W. Eaton <jwe@octave.org> parents: 8395 diff changeset	181 ## prefers accuracy to performance.
8306 43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	182
43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	183 macheps = eps (class (x0));
43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	184
10201 5c66978f3fdf support TypicalX and AutoScaling in fsolve/fminunc, don't autoscale by default Jaroslav Hajek <highegg@gmail.com> parents: 10050 diff changeset	185 tolx = optimget (options, "TolX", 1e-7);
5c66978f3fdf support TypicalX and AutoScaling in fsolve/fminunc, don't autoscale by default Jaroslav Hajek <highegg@gmail.com> parents: 10050 diff changeset	186 tolf = optimget (options, "TolFun", 1e-7);
8306 43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	187
9623 bc0739d02724 update initial TR step for fsolve and fminunc Jaroslav Hajek <highegg@gmail.com> parents: 9464 diff changeset	188 factor = 1;
8306 43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	189
8467 77b8d4aa2743 fsolve.m, fzero.m: style fixes; use strcmpi to compare options John W. Eaton <jwe@octave.org> parents: 8466 diff changeset	190 niter = 1;
8986 22c8272af34b improvements to fsolve & family Jaroslav Hajek <highegg@gmail.com> parents: 8857 diff changeset	191 nfev = 1;
8306 43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	192
43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	193 x = x0(:);
43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	194 info = 0;
43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	195
8986 22c8272af34b improvements to fsolve & family Jaroslav Hajek <highegg@gmail.com> parents: 8857 diff changeset	196 ## Initial evaluation.
8990 349d75161672 more cosmetic adjustments to fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8988 diff changeset	197 ## 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	198 fvec = fcn (reshape (x, xsiz));
22c8272af34b improvements to fsolve & family Jaroslav Hajek <highegg@gmail.com> parents: 8857 diff changeset	199 fsiz = size (fvec);
22c8272af34b improvements to fsolve & family Jaroslav Hajek <highegg@gmail.com> parents: 8857 diff changeset	200 fvec = fvec(:);
22c8272af34b improvements to fsolve & family Jaroslav Hajek <highegg@gmail.com> parents: 8857 diff changeset	201 fn = norm (fvec);
22c8272af34b improvements to fsolve & family Jaroslav Hajek <highegg@gmail.com> parents: 8857 diff changeset	202 m = length (fvec);
22c8272af34b improvements to fsolve & family Jaroslav Hajek <highegg@gmail.com> parents: 8857 diff changeset	203 n = length (x);
22c8272af34b improvements to fsolve & family Jaroslav Hajek <highegg@gmail.com> parents: 8857 diff changeset	204
22c8272af34b improvements to fsolve & family Jaroslav Hajek <highegg@gmail.com> parents: 8857 diff changeset	205 if (! isempty (outfcn))
22c8272af34b improvements to fsolve & family Jaroslav Hajek <highegg@gmail.com> parents: 8857 diff changeset	206 optimvalues.iter = niter;
22c8272af34b improvements to fsolve & family Jaroslav Hajek <highegg@gmail.com> parents: 8857 diff changeset	207 optimvalues.funccount = nfev;
22c8272af34b improvements to fsolve & family Jaroslav Hajek <highegg@gmail.com> parents: 8857 diff changeset	208 optimvalues.fval = fn;
22c8272af34b improvements to fsolve & family Jaroslav Hajek <highegg@gmail.com> parents: 8857 diff changeset	209 optimvalues.searchdirection = zeros (n, 1);
22c8272af34b improvements to fsolve & family Jaroslav Hajek <highegg@gmail.com> parents: 8857 diff changeset	210 state = 'init';
22c8272af34b improvements to fsolve & family Jaroslav Hajek <highegg@gmail.com> parents: 8857 diff changeset	211 stop = outfcn (x, optimvalues, state);
22c8272af34b improvements to fsolve & family Jaroslav Hajek <highegg@gmail.com> parents: 8857 diff changeset	212 if (stop)
22c8272af34b improvements to fsolve & family Jaroslav Hajek <highegg@gmail.com> parents: 8857 diff changeset	213 info = -1;
22c8272af34b improvements to fsolve & family Jaroslav Hajek <highegg@gmail.com> parents: 8857 diff changeset	214 break;
22c8272af34b improvements to fsolve & family Jaroslav Hajek <highegg@gmail.com> parents: 8857 diff changeset	215 endif
22c8272af34b improvements to fsolve & family Jaroslav Hajek <highegg@gmail.com> parents: 8857 diff changeset	216 endif
22c8272af34b improvements to fsolve & family Jaroslav Hajek <highegg@gmail.com> parents: 8857 diff changeset	217
22c8272af34b improvements to fsolve & family Jaroslav Hajek <highegg@gmail.com> parents: 8857 diff changeset	218 nsuciter = 0;
22c8272af34b improvements to fsolve & family Jaroslav Hajek <highegg@gmail.com> parents: 8857 diff changeset	219
8507 cadc73247d65 style fixes John W. Eaton <jwe@octave.org> parents: 8468 diff changeset	220 ## Outer loop.
8306 43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	221 while (niter < maxiter && nfev < maxfev && ! info)
43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	222
8507 cadc73247d65 style fixes John W. Eaton <jwe@octave.org> parents: 8468 diff changeset	223 ## Calculate function value and Jacobian (possibly via FD).
8306 43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	224 if (has_jac)
43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	225 [fvec, fjac] = fcn (reshape (x, xsiz));
8592 dacfd030633a handle sparse jacobians in fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8590 diff changeset	226 ## If the jacobian is sparse, disable Broyden updating.
dacfd030633a handle sparse jacobians in fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8590 diff changeset	227 if (issparse (fjac))
dacfd030633a handle sparse jacobians in fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8590 diff changeset	228 updating = false;
dacfd030633a handle sparse jacobians in fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8590 diff changeset	229 endif
8986 22c8272af34b improvements to fsolve & family Jaroslav Hajek <highegg@gmail.com> parents: 8857 diff changeset	230 fvec = fvec(:);
8306 43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	231 nfev ++;
43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	232 else
10201 5c66978f3fdf support TypicalX and AutoScaling in fsolve/fminunc, don't autoscale by default Jaroslav Hajek <highegg@gmail.com> parents: 10050 diff changeset	233 fjac = __fdjac__ (fcn, reshape (x, xsiz), fvec, typicalx, cdif);
9212 6feb27c38da1 support central differences in fminunc and fsolve Jaroslav Hajek <highegg@gmail.com> parents: 9209 diff changeset	234 nfev += (1 + cdif) * length (x);
8306 43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	235 endif
8590 c136d313206a improvements to fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	236
8763 5ce12bca4c51 update comments in fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8693 diff changeset	237 ## For square and overdetermined systems, we update a QR
8590 c136d313206a improvements to fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	238 ## factorization of the jacobian to avoid solving a full system in each
c136d313206a improvements to fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	239 ## step. In this case, we pass a triangular matrix to __dogleg__.
c136d313206a improvements to fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	240 useqr = updating && m >= n && n > 10;
c136d313206a improvements to fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	241
c136d313206a improvements to fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	242 if (useqr)
8616 3d75d717cbe0 do not pivot by default in fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8615 diff changeset	243 ## 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	244 ## 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	245 ## 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	246 ## 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	247 ## when qr updating was used.
3d75d717cbe0 do not pivot by default in fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8615 diff changeset	248 [q, r] = qr (fjac, 0);
8383 a762d9daa700 allow overdetermined systems in fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8306 diff changeset	249 endif
8306 43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	250
10201 5c66978f3fdf support TypicalX and AutoScaling in fsolve/fminunc, don't autoscale by default Jaroslav Hajek <highegg@gmail.com> parents: 10050 diff changeset	251 if (autoscale)
5c66978f3fdf support TypicalX and AutoScaling in fsolve/fminunc, don't autoscale by default Jaroslav Hajek <highegg@gmail.com> parents: 10050 diff changeset	252 ## Get column norms, use them as scaling factors.
5c66978f3fdf support TypicalX and AutoScaling in fsolve/fminunc, don't autoscale by default Jaroslav Hajek <highegg@gmail.com> parents: 10050 diff changeset	253 jcn = norm (fjac, 'columns').';
5c66978f3fdf support TypicalX and AutoScaling in fsolve/fminunc, don't autoscale by default Jaroslav Hajek <highegg@gmail.com> parents: 10050 diff changeset	254 if (niter == 1)
5c66978f3fdf support TypicalX and AutoScaling in fsolve/fminunc, don't autoscale by default Jaroslav Hajek <highegg@gmail.com> parents: 10050 diff changeset	255 dg = jcn;
5c66978f3fdf support TypicalX and AutoScaling in fsolve/fminunc, don't autoscale by default Jaroslav Hajek <highegg@gmail.com> parents: 10050 diff changeset	256 dg(dg == 0) = 1;
5c66978f3fdf support TypicalX and AutoScaling in fsolve/fminunc, don't autoscale by default Jaroslav Hajek <highegg@gmail.com> parents: 10050 diff changeset	257 else
5c66978f3fdf support TypicalX and AutoScaling in fsolve/fminunc, don't autoscale by default Jaroslav Hajek <highegg@gmail.com> parents: 10050 diff changeset	258 ## Rescale adaptively.
5c66978f3fdf support TypicalX and AutoScaling in fsolve/fminunc, don't autoscale by default Jaroslav Hajek <highegg@gmail.com> parents: 10050 diff changeset	259 ## FIXME: the original minpack used the following rescaling strategy:
5c66978f3fdf support TypicalX and AutoScaling in fsolve/fminunc, don't autoscale by default Jaroslav Hajek <highegg@gmail.com> parents: 10050 diff changeset	260 ## dg = max (dg, jcn);
5c66978f3fdf support TypicalX and AutoScaling in fsolve/fminunc, don't autoscale by default Jaroslav Hajek <highegg@gmail.com> parents: 10050 diff changeset	261 ## but it seems not good if we start with a bad guess yielding jacobian
5c66978f3fdf support TypicalX and AutoScaling in fsolve/fminunc, don't autoscale by default Jaroslav Hajek <highegg@gmail.com> parents: 10050 diff changeset	262 ## columns with large norms that later decrease, because the corresponding
5c66978f3fdf support TypicalX and AutoScaling in fsolve/fminunc, don't autoscale by default Jaroslav Hajek <highegg@gmail.com> parents: 10050 diff changeset	263 ## variable will still be overscaled. So instead, we only give the old
5c66978f3fdf support TypicalX and AutoScaling in fsolve/fminunc, don't autoscale by default Jaroslav Hajek <highegg@gmail.com> parents: 10050 diff changeset	264 ## scaling a small momentum, but do not honor it.
5c66978f3fdf support TypicalX and AutoScaling in fsolve/fminunc, don't autoscale by default Jaroslav Hajek <highegg@gmail.com> parents: 10050 diff changeset	265
5c66978f3fdf support TypicalX and AutoScaling in fsolve/fminunc, don't autoscale by default Jaroslav Hajek <highegg@gmail.com> parents: 10050 diff changeset	266 dg = max (0.1*dg, jcn);
5c66978f3fdf support TypicalX and AutoScaling in fsolve/fminunc, don't autoscale by default Jaroslav Hajek <highegg@gmail.com> parents: 10050 diff changeset	267 endif
5c66978f3fdf support TypicalX and AutoScaling in fsolve/fminunc, don't autoscale by default Jaroslav Hajek <highegg@gmail.com> parents: 10050 diff changeset	268 endif
5c66978f3fdf support TypicalX and AutoScaling in fsolve/fminunc, don't autoscale by default Jaroslav Hajek <highegg@gmail.com> parents: 10050 diff changeset	269
8306 43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	270 if (niter == 1)
43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	271 xn = norm (dg .* x);
8590 c136d313206a improvements to fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	272 ## FIXME: something better?
9628 73e6ad869f08 further correct initial TR step strategy Jaroslav Hajek <highegg@gmail.com> parents: 9623 diff changeset	273 delta = factor * max (xn, 1);
8306 43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	274 endif
43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	275
8990 349d75161672 more cosmetic adjustments to fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8988 diff changeset	276 ## 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	277 ## 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	278 ## 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	279 ## of course, how to define the above terms :)
8306 43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	280
8986 22c8272af34b improvements to fsolve & family Jaroslav Hajek <highegg@gmail.com> parents: 8857 diff changeset	281 lastratio = 0;
8306 43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	282 nfail = 0;
43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	283 nsuc = 0;
9207 25f50d2d76b3 improve TR updating strategy for fminunc and fsolve Jaroslav Hajek <highegg@gmail.com> parents: 9153 diff changeset	284 decfac = 0.5;
25f50d2d76b3 improve TR updating strategy for fminunc and fsolve Jaroslav Hajek <highegg@gmail.com> parents: 9153 diff changeset	285
8507 cadc73247d65 style fixes John W. Eaton <jwe@octave.org> parents: 8468 diff changeset	286 ## Inner loop.
8306 43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	287 while (niter <= maxiter && nfev < maxfev && ! info)
43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	288
8590 c136d313206a improvements to fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	289 ## Get trust-region model (dogleg) minimizer.
c136d313206a improvements to fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	290 if (useqr)
c136d313206a improvements to fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	291 qtf = q'*fvec;
c136d313206a improvements to fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	292 s = - __dogleg__ (r, qtf, dg, delta);
c136d313206a improvements to fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	293 w = qtf + r * s;
c136d313206a improvements to fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	294 else
c136d313206a improvements to fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	295 s = - __dogleg__ (fjac, fvec, dg, delta);
c136d313206a improvements to fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	296 w = fvec + fjac * s;
8395 88fd356b0d95 optionally allow pivoting in fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8383 diff changeset	297 endif
8590 c136d313206a improvements to fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	298
8395 88fd356b0d95 optionally allow pivoting in fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8383 diff changeset	299 sn = norm (dg .* s);
8306 43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	300 if (niter == 1)
8395 88fd356b0d95 optionally allow pivoting in fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8383 diff changeset	301 delta = min (delta, sn);
8306 43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	302 endif
43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	303
8395 88fd356b0d95 optionally allow pivoting in fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8383 diff changeset	304 fvec1 = fcn (reshape (x + s, xsiz)) (:);
8306 43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	305 fn1 = norm (fvec1);
8590 c136d313206a improvements to fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	306 nfev ++;
8306 43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	307
43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	308 if (fn1 < fn)
8507 cadc73247d65 style fixes John W. Eaton <jwe@octave.org> parents: 8468 diff changeset	309 ## Scaled actual reduction.
8306 43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	310 actred = 1 - (fn1/fn)^2;
43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	311 else
43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	312 actred = -1;
43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	313 endif
43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	314
8507 cadc73247d65 style fixes John W. Eaton <jwe@octave.org> parents: 8468 diff changeset	315 ## Scaled predicted reduction, and ratio.
8306 43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	316 t = norm (w);
43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	317 if (t < fn)
43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	318 prered = 1 - (t/fn)^2;
43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	319 ratio = actred / prered;
43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	320 else
43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	321 prered = 0;
43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	322 ratio = 0;
43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	323 endif
43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	324
8507 cadc73247d65 style fixes John W. Eaton <jwe@octave.org> parents: 8468 diff changeset	325 ## Update delta.
9207 25f50d2d76b3 improve TR updating strategy for fminunc and fsolve Jaroslav Hajek <highegg@gmail.com> parents: 9153 diff changeset	326 if (ratio < min(max(0.1, 0.8*lastratio), 0.9))
8306 43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	327 nsuc = 0;
43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	328 nfail ++;
9207 25f50d2d76b3 improve TR updating strategy for fminunc and fsolve Jaroslav Hajek <highegg@gmail.com> parents: 9153 diff changeset	329 delta *= decfac;
25f50d2d76b3 improve TR updating strategy for fminunc and fsolve Jaroslav Hajek <highegg@gmail.com> parents: 9153 diff changeset	330 decfac ^= 1.4142;
8590 c136d313206a improvements to fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	331 if (delta <= 1e1machepsxn)
8507 cadc73247d65 style fixes John W. Eaton <jwe@octave.org> parents: 8468 diff changeset	332 ## Trust region became uselessly small.
8306 43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	333 info = -3;
43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	334 break;
43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	335 endif
43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	336 else
8986 22c8272af34b improvements to fsolve & family Jaroslav Hajek <highegg@gmail.com> parents: 8857 diff changeset	337 lastratio = ratio;
9207 25f50d2d76b3 improve TR updating strategy for fminunc and fsolve Jaroslav Hajek <highegg@gmail.com> parents: 9153 diff changeset	338 decfac = 0.5;
8306 43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	339 nfail = 0;
43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	340 nsuc ++;
43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	341 if (abs (1-ratio) <= 0.1)
9207 25f50d2d76b3 improve TR updating strategy for fminunc and fsolve Jaroslav Hajek <highegg@gmail.com> parents: 9153 diff changeset	342 delta = 1.4142*sn;
8306 43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	343 elseif (ratio >= 0.5 \|\| nsuc > 1)
9207 25f50d2d76b3 improve TR updating strategy for fminunc and fsolve Jaroslav Hajek <highegg@gmail.com> parents: 9153 diff changeset	344 delta = max (delta, 1.4142*sn);
8306 43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	345 endif
43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	346 endif
43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	347
43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	348 if (ratio >= 1e-4)
8507 cadc73247d65 style fixes John W. Eaton <jwe@octave.org> parents: 8468 diff changeset	349 ## Successful iteration.
8395 88fd356b0d95 optionally allow pivoting in fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8383 diff changeset	350 x += s;
8306 43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	351 xn = norm (dg .* x);
43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	352 fvec = fvec1;
43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	353 fn = fn1;
8986 22c8272af34b improvements to fsolve & family Jaroslav Hajek <highegg@gmail.com> parents: 8857 diff changeset	354 nsuciter ++;
8306 43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	355 endif
43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	356
8857 20589a8f1a33 fix fsolve bug Jaroslav Hajek <highegg@gmail.com> parents: 8763 diff changeset	357 niter ++;
20589a8f1a33 fix fsolve bug Jaroslav Hajek <highegg@gmail.com> parents: 8763 diff changeset	358
8590 c136d313206a improvements to fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	359 ## FIXME: should outputfcn be only called after a successful iteration?
8615 7baacb6a8a65 fix typo in fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8613 diff changeset	360 if (! isempty (outfcn))
8590 c136d313206a improvements to fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	361 optimvalues.iter = niter;
c136d313206a improvements to fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	362 optimvalues.funccount = nfev;
c136d313206a improvements to fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	363 optimvalues.fval = fn;
c136d313206a improvements to fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	364 optimvalues.searchdirection = s;
c136d313206a improvements to fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	365 state = 'iter';
c136d313206a improvements to fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	366 stop = outfcn (x, optimvalues, state);
c136d313206a improvements to fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	367 if (stop)
c136d313206a improvements to fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	368 info = -1;
c136d313206a improvements to fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	369 break;
c136d313206a improvements to fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	370 endif
c136d313206a improvements to fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	371 endif
c136d313206a improvements to fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	372
8466 4d008d9f0ccf fsolve.m: style fixes John W. Eaton <jwe@octave.org> parents: 8395 diff changeset	373 ## Tests for termination conditions. A mysterious place, anything
4d008d9f0ccf fsolve.m: style fixes John W. Eaton <jwe@octave.org> parents: 8395 diff changeset	374 ## can happen if you change something here...
8306 43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	375
8466 4d008d9f0ccf fsolve.m: style fixes John W. Eaton <jwe@octave.org> parents: 8395 diff changeset	376 ## The rule of thumb (which I'm not sure M*b is quite following)
4d008d9f0ccf fsolve.m: style fixes John W. Eaton <jwe@octave.org> parents: 8395 diff changeset	377 ## is that for a tolerance that depends on scaling, only 0 makes
4d008d9f0ccf fsolve.m: style fixes John W. Eaton <jwe@octave.org> parents: 8395 diff changeset	378 ## sense as a default value. But 0 usually means uselessly long
4d008d9f0ccf fsolve.m: style fixes John W. Eaton <jwe@octave.org> parents: 8395 diff changeset	379 ## iterations, so we need scaling-independent tolerances wherever
4d008d9f0ccf fsolve.m: style fixes John W. Eaton <jwe@octave.org> parents: 8395 diff changeset	380 ## possible.
8306 43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	381
8466 4d008d9f0ccf fsolve.m: style fixes John W. Eaton <jwe@octave.org> parents: 8395 diff changeset	382 ## FIXME -- why tolfnxn? If abs (e) ~ abs(x) * eps is a vector
4d008d9f0ccf fsolve.m: style fixes John W. Eaton <jwe@octave.org> parents: 8395 diff changeset	383 ## of perturbations of x, then norm (fjace) <= epsn*xn, i.e. by
4d008d9f0ccf fsolve.m: style fixes John W. Eaton <jwe@octave.org> parents: 8395 diff changeset	384 ## 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	385 if (fn <= tolfnxn)
43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	386 info = 1;
10549 95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10201 diff changeset	387 ## The following tests done only after successful step.
9075 7b35d3b21c1c fix test in fsolve Jaroslav Hajek <highegg@gmail.com> parents: 9065 diff changeset	388 elseif (ratio >= 1e-4)
8466 4d008d9f0ccf fsolve.m: style fixes John W. Eaton <jwe@octave.org> parents: 8395 diff changeset	389 ## This one is classic. Note that we use scaled variables again,
10549 95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10201 diff changeset	390 ## but compare to scaled step, so nothing bad.
8395 88fd356b0d95 optionally allow pivoting in fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8383 diff changeset	391 if (sn <= tolx*xn)
8306 43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	392 info = 2;
8466 4d008d9f0ccf fsolve.m: style fixes John W. Eaton <jwe@octave.org> parents: 8395 diff changeset	393 ## Again a classic one. It seems weird to use the same tolf
10549 95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10201 diff changeset	394 ## for two different tests, but that's what M*b manual appears
95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10201 diff changeset	395 ## to say.
8306 43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	396 elseif (actred < tolf)
8383 a762d9daa700 allow overdetermined systems in fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8306 diff changeset	397 info = 3;
8306 43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	398 endif
43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	399 endif
43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	400
8507 cadc73247d65 style fixes John W. Eaton <jwe@octave.org> parents: 8468 diff changeset	401 ## Criterion for recalculating jacobian.
8986 22c8272af34b improvements to fsolve & family Jaroslav Hajek <highegg@gmail.com> parents: 8857 diff changeset	402 if (! updating \|\| nfail == 2 \|\| nsuciter < 2)
8306 43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	403 break;
43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	404 endif
43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	405
8507 cadc73247d65 style fixes John W. Eaton <jwe@octave.org> parents: 8468 diff changeset	406 ## Compute the scaled Broyden update.
8590 c136d313206a improvements to fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	407 if (useqr)
c136d313206a improvements to fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	408 u = (fvec1 - q*w) / sn;
c136d313206a improvements to fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	409 v = dg .* ((dg .* s) / sn);
c136d313206a improvements to fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	410
c136d313206a improvements to fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	411 ## Update the QR factorization.
c136d313206a improvements to fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	412 [q, r] = qrupdate (q, r, u, v);
c136d313206a improvements to fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	413 else
c136d313206a improvements to fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	414 u = (fvec1 - w);
c136d313206a improvements to fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	415 v = dg .* ((dg .* s) / sn);
c136d313206a improvements to fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	416
c136d313206a improvements to fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	417 ## update the jacobian
c136d313206a improvements to fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	418 fjac += u * v';
8395 88fd356b0d95 optionally allow pivoting in fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8383 diff changeset	419 endif
8306 43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	420 endwhile
43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	421 endwhile
43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	422
8507 cadc73247d65 style fixes John W. Eaton <jwe@octave.org> parents: 8468 diff changeset	423 ## Restore original shapes.
8306 43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	424 x = reshape (x, xsiz);
43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	425 fvec = reshape (fvec, fsiz);
43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	426
43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	427 output.iterations = niter;
8986 22c8272af34b improvements to fsolve & family Jaroslav Hajek <highegg@gmail.com> parents: 8857 diff changeset	428 output.successful = nsuciter;
8590 c136d313206a improvements to fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	429 output.funcCount = nfev;
8306 43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	430
43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	431 endfunction
43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	432
8466 4d008d9f0ccf fsolve.m: style fixes John W. Eaton <jwe@octave.org> parents: 8395 diff changeset	433 ## An assistant function that evaluates a function handle and checks for
4d008d9f0ccf fsolve.m: style fixes John W. Eaton <jwe@octave.org> parents: 8395 diff changeset	434 ## bad results.
8693 e5ffb52c9c61 improve fsolve and add ComplexEqn option Jaroslav Hajek <highegg@gmail.com> parents: 8661 diff changeset	435 function [fx, jx] = guarded_eval (fun, x, complexeqn)
e5ffb52c9c61 improve fsolve and add ComplexEqn option Jaroslav Hajek <highegg@gmail.com> parents: 8661 diff changeset	436 if (nargout > 1)
e5ffb52c9c61 improve fsolve and add ComplexEqn option Jaroslav Hajek <highegg@gmail.com> parents: 8661 diff changeset	437 [fx, jx] = fun (x);
e5ffb52c9c61 improve fsolve and add ComplexEqn option Jaroslav Hajek <highegg@gmail.com> parents: 8661 diff changeset	438 else
e5ffb52c9c61 improve fsolve and add ComplexEqn option Jaroslav Hajek <highegg@gmail.com> parents: 8661 diff changeset	439 fx = fun (x);
8997 187a9d9c2f04 simplifications to fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8991 diff changeset	440 jx = [];
8693 e5ffb52c9c61 improve fsolve and add ComplexEqn option Jaroslav Hajek <highegg@gmail.com> parents: 8661 diff changeset	441 endif
e5ffb52c9c61 improve fsolve and add ComplexEqn option Jaroslav Hajek <highegg@gmail.com> parents: 8661 diff changeset	442
8997 187a9d9c2f04 simplifications to fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8991 diff changeset	443 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	444 error ("fsolve:notreal", "fsolve: non-real value encountered");
8997 187a9d9c2f04 simplifications to fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8991 diff changeset	445 elseif (complexeqn && ! (isnumeric (fx) && isnumeric(jx)))
8693 e5ffb52c9c61 improve fsolve and add ComplexEqn option Jaroslav Hajek <highegg@gmail.com> parents: 8661 diff changeset	446 error ("fsolve:notnum", "fsolve: non-numeric value encountered");
e5ffb52c9c61 improve fsolve and add ComplexEqn option Jaroslav Hajek <highegg@gmail.com> parents: 8661 diff changeset	447 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	448 error ("fsolve:isnan", "fsolve: NaN value encountered");
8306 43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	449 endif
43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	450 endfunction
43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	451
10050 dc88a0b6472c support old style jacobian for fsolve Jaroslav Hajek <highegg@gmail.com> parents: 9899 diff changeset	452 function [fx, jx] = make_fcn_jac (x, fcn, fjac)
dc88a0b6472c support old style jacobian for fsolve Jaroslav Hajek <highegg@gmail.com> parents: 9899 diff changeset	453 fx = fcn (x);
dc88a0b6472c support old style jacobian for fsolve Jaroslav Hajek <highegg@gmail.com> parents: 9899 diff changeset	454 if (nargout == 2)
dc88a0b6472c support old style jacobian for fsolve Jaroslav Hajek <highegg@gmail.com> parents: 9899 diff changeset	455 jx = fjac (x);
dc88a0b6472c support old style jacobian for fsolve Jaroslav Hajek <highegg@gmail.com> parents: 9899 diff changeset	456 endif
dc88a0b6472c support old style jacobian for fsolve Jaroslav Hajek <highegg@gmail.com> parents: 9899 diff changeset	457 endfunction
dc88a0b6472c support old style jacobian for fsolve Jaroslav Hajek <highegg@gmail.com> parents: 9899 diff changeset	458
8306 43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	459 %!function retval = f (p)
43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	460 %! x = p(1);
43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	461 %! y = p(2);
43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	462 %! z = p(3);
43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	463 %! retval = zeros (3, 1);
43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	464 %! retval(1) = sin(x) + y**2 + log(z) - 7;
43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	465 %! retval(2) = 3x + 2y -z*3 + 1;
43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	466 %! retval(3) = x + y + z - 5;
43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	467 %!test
43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	468 %! x_opt = [ 0.599054;
43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	469 %! 2.395931;
43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	470 %! 2.005014 ];
43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	471 %! tol = 1.0e-5;
43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	472 %! [x, fval, info] = fsolve (@f, [ 0.5; 2.0; 2.5 ]);
43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	473 %! assert (info > 0);
8383 a762d9daa700 allow overdetermined systems in fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8306 diff changeset	474 %! assert (norm (x - x_opt, Inf) < tol);
8306 43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	475 %! assert (norm (fval) < tol);
43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	476
43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	477 %!function retval = f (p)
43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	478 %! x = p(1);
43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	479 %! y = p(2);
43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	480 %! z = p(3);
43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	481 %! w = p(4);
43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	482 %! retval = zeros (4, 1);
43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	483 %! retval(1) = 3x + 4y + exp (z + w) - 1.007;
43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	484 %! retval(2) = 6x - 4y + exp (3*z + w) - 11;
43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	485 %! retval(3) = x^4 - 4y^2 + 6z - 8*w - 20;
43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	486 %! retval(4) = x^2 + 2*y^3 + z - w - 4;
43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	487 %!test
43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	488 %! x_opt = [ -0.767297326653401, 0.590671081117440, 1.47190018629642, -1.52719341133957 ];
43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	489 %! tol = 1.0e-5;
43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	490 %! [x, fval, info] = fsolve (@f, [-1, 1, 2, -1]);
43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	491 %! assert (info > 0);
8383 a762d9daa700 allow overdetermined systems in fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8306 diff changeset	492 %! assert (norm (x - x_opt, Inf) < tol);
8306 43795cf108d0 initial implementation of fsolve Jaroslav Hajek <highegg@gmail.com> parents: diff changeset	493 %! assert (norm (fval) < tol);
8383 a762d9daa700 allow overdetermined systems in fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8306 diff changeset	494
a762d9daa700 allow overdetermined systems in fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8306 diff changeset	495 %!function retval = f (p)
a762d9daa700 allow overdetermined systems in fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8306 diff changeset	496 %! x = p(1);
a762d9daa700 allow overdetermined systems in fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8306 diff changeset	497 %! y = p(2);
a762d9daa700 allow overdetermined systems in fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8306 diff changeset	498 %! z = p(3);
a762d9daa700 allow overdetermined systems in fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8306 diff changeset	499 %! retval = zeros (3, 1);
a762d9daa700 allow overdetermined systems in fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8306 diff changeset	500 %! retval(1) = sin(x) + y**2 + log(z) - 7;
a762d9daa700 allow overdetermined systems in fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8306 diff changeset	501 %! retval(2) = 3x + 2y -z*3 + 1;
a762d9daa700 allow overdetermined systems in fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8306 diff changeset	502 %! retval(3) = x + y + z - 5;
a762d9daa700 allow overdetermined systems in fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8306 diff changeset	503 %! retval(4) = xx + y - zlog(z) - 1.36;
a762d9daa700 allow overdetermined systems in fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8306 diff changeset	504 %!test
a762d9daa700 allow overdetermined systems in fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8306 diff changeset	505 %! x_opt = [ 0.599054;
a762d9daa700 allow overdetermined systems in fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8306 diff changeset	506 %! 2.395931;
a762d9daa700 allow overdetermined systems in fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8306 diff changeset	507 %! 2.005014 ];
a762d9daa700 allow overdetermined systems in fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8306 diff changeset	508 %! tol = 1.0e-5;
8590 c136d313206a improvements to fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	509 %! [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	510 %! assert (info > 0);
a762d9daa700 allow overdetermined systems in fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8306 diff changeset	511 %! assert (norm (x - x_opt, Inf) < tol);
a762d9daa700 allow overdetermined systems in fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8306 diff changeset	512 %! assert (norm (fval) < tol);
a762d9daa700 allow overdetermined systems in fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8306 diff changeset	513
8395 88fd356b0d95 optionally allow pivoting in fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8383 diff changeset	514 %!function retval = f (p)
88fd356b0d95 optionally allow pivoting in fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8383 diff changeset	515 %! x = p(1);
88fd356b0d95 optionally allow pivoting in fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8383 diff changeset	516 %! y = p(2);
88fd356b0d95 optionally allow pivoting in fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8383 diff changeset	517 %! z = p(3);
88fd356b0d95 optionally allow pivoting in fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8383 diff changeset	518 %! retval = zeros (3, 1);
88fd356b0d95 optionally allow pivoting in fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8383 diff changeset	519 %! retval(1) = sin(x) + y**2 + log(z) - 7;
88fd356b0d95 optionally allow pivoting in fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8383 diff changeset	520 %! retval(2) = 3x + 2y -z*3 + 1;
88fd356b0d95 optionally allow pivoting in fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8383 diff changeset	521 %! retval(3) = x + y + z - 5;
88fd356b0d95 optionally allow pivoting in fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8383 diff changeset	522 %!test
88fd356b0d95 optionally allow pivoting in fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8383 diff changeset	523 %! x_opt = [ 0.599054;
88fd356b0d95 optionally allow pivoting in fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8383 diff changeset	524 %! 2.395931;
88fd356b0d95 optionally allow pivoting in fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8383 diff changeset	525 %! 2.005014 ];
88fd356b0d95 optionally allow pivoting in fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8383 diff changeset	526 %! tol = 1.0e-5;
8590 c136d313206a improvements to fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	527 %! opt = optimset ("Updating", "qrp");
8395 88fd356b0d95 optionally allow pivoting in fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8383 diff changeset	528 %! [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	529 %! assert (info > 0);
88fd356b0d95 optionally allow pivoting in fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8383 diff changeset	530 %! assert (norm (x - x_opt, Inf) < tol);
88fd356b0d95 optionally allow pivoting in fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8383 diff changeset	531 %! assert (norm (fval) < tol);
8590 c136d313206a improvements to fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	532
c136d313206a improvements to fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	533 %!test
c136d313206a improvements to fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	534 %! b0 = 3;
c136d313206a improvements to fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	535 %! a0 = 0.2;
c136d313206a improvements to fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	536 %! x = 0:.5:5;
c136d313206a improvements to fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	537 %! noise = 1e-5 * sin (100*x);
c136d313206a improvements to fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	538 %! y = exp (-a0*x) + b0 + noise;
c136d313206a improvements to fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	539 %! c_opt = [a0, b0];
c136d313206a improvements to fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	540 %! tol = 1e-5;
c136d313206a improvements to fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	541 %!
c136d313206a improvements to fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	542 %! [c, fval, info, output] = fsolve (@(c) (exp(-c(1)*x) + c(2) - y), [0, 0]);
c136d313206a improvements to fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	543 %! assert (info > 0);
c136d313206a improvements to fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	544 %! assert (norm (c - c_opt, Inf) < tol);
c136d313206a improvements to fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	545 %! assert (norm (fval) < norm (noise));
c136d313206a improvements to fsolve Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	546
8693 e5ffb52c9c61 improve fsolve and add ComplexEqn option Jaroslav Hajek <highegg@gmail.com> parents: 8661 diff changeset	547
e5ffb52c9c61 improve fsolve and add ComplexEqn option Jaroslav Hajek <highegg@gmail.com> parents: 8661 diff changeset	548 %!function y = cfun (x)
e5ffb52c9c61 improve fsolve and add ComplexEqn option Jaroslav Hajek <highegg@gmail.com> parents: 8661 diff changeset	549 %! 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	550 %! 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	551 %! 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	552
e5ffb52c9c61 improve fsolve and add ComplexEqn option Jaroslav Hajek <highegg@gmail.com> parents: 8661 diff changeset	553 %!test
e5ffb52c9c61 improve fsolve and add ComplexEqn option Jaroslav Hajek <highegg@gmail.com> parents: 8661 diff changeset	554 %! x_opt = [-1+i, 1-i, 2+i];
e5ffb52c9c61 improve fsolve and add ComplexEqn option Jaroslav Hajek <highegg@gmail.com> parents: 8661 diff changeset	555 %! x = [i, 1, 1+i];
e5ffb52c9c61 improve fsolve and add ComplexEqn option Jaroslav Hajek <highegg@gmail.com> parents: 8661 diff changeset	556 %!
e5ffb52c9c61 improve fsolve and add ComplexEqn option Jaroslav Hajek <highegg@gmail.com> parents: 8661 diff changeset	557 %! [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	558 %! tol = 1e-5;
e5ffb52c9c61 improve fsolve and add ComplexEqn option Jaroslav Hajek <highegg@gmail.com> parents: 8661 diff changeset	559 %! assert (norm (f) < tol);
e5ffb52c9c61 improve fsolve and add ComplexEqn option Jaroslav Hajek <highegg@gmail.com> parents: 8661 diff changeset	560 %! assert (norm (x - x_opt, Inf) < tol);
e5ffb52c9c61 improve fsolve and add ComplexEqn option Jaroslav Hajek <highegg@gmail.com> parents: 8661 diff changeset	561
9899 9f25290a35e8 more private function and subfunction changes John W. Eaton <jwe@octave.org> parents: 9724 diff changeset	562 ## Solve the double dogleg trust-region least-squares problem:
9f25290a35e8 more private function and subfunction changes John W. Eaton <jwe@octave.org> parents: 9724 diff changeset	563 ## Minimize norm(rx-b) subject to the constraint norm(d.x) <= delta,
9f25290a35e8 more private function and subfunction changes John W. Eaton <jwe@octave.org> parents: 9724 diff changeset	564 ## x being a convex combination of the gauss-newton and scaled gradient.
9f25290a35e8 more private function and subfunction changes John W. Eaton <jwe@octave.org> parents: 9724 diff changeset	565
9f25290a35e8 more private function and subfunction changes John W. Eaton <jwe@octave.org> parents: 9724 diff changeset	566 ## TODO: error checks
9f25290a35e8 more private function and subfunction changes John W. Eaton <jwe@octave.org> parents: 9724 diff changeset	567 ## TODO: handle singularity, or leave it up to mldivide?
9f25290a35e8 more private function and subfunction changes John W. Eaton <jwe@octave.org> parents: 9724 diff changeset	568
9f25290a35e8 more private function and subfunction changes John W. Eaton <jwe@octave.org> parents: 9724 diff changeset	569 function x = __dogleg__ (r, b, d, delta)
9f25290a35e8 more private function and subfunction changes John W. Eaton <jwe@octave.org> parents: 9724 diff changeset	570 ## Get Gauss-Newton direction.
9f25290a35e8 more private function and subfunction changes John W. Eaton <jwe@octave.org> parents: 9724 diff changeset	571 x = r \ b;
9f25290a35e8 more private function and subfunction changes John W. Eaton <jwe@octave.org> parents: 9724 diff changeset	572 xn = norm (d .* x);
9f25290a35e8 more private function and subfunction changes John W. Eaton <jwe@octave.org> parents: 9724 diff changeset	573 if (xn > delta)
9f25290a35e8 more private function and subfunction changes John W. Eaton <jwe@octave.org> parents: 9724 diff changeset	574 ## GN is too big, get scaled gradient.
9f25290a35e8 more private function and subfunction changes John W. Eaton <jwe@octave.org> parents: 9724 diff changeset	575 s = (r' * b) ./ d;
9f25290a35e8 more private function and subfunction changes John W. Eaton <jwe@octave.org> parents: 9724 diff changeset	576 sn = norm (s);
9f25290a35e8 more private function and subfunction changes John W. Eaton <jwe@octave.org> parents: 9724 diff changeset	577 if (sn > 0)
9f25290a35e8 more private function and subfunction changes John W. Eaton <jwe@octave.org> parents: 9724 diff changeset	578 ## Normalize and rescale.
9f25290a35e8 more private function and subfunction changes John W. Eaton <jwe@octave.org> parents: 9724 diff changeset	579 s = (s / sn) ./ d;
9f25290a35e8 more private function and subfunction changes John W. Eaton <jwe@octave.org> parents: 9724 diff changeset	580 ## Get the line minimizer in s direction.
9f25290a35e8 more private function and subfunction changes John W. Eaton <jwe@octave.org> parents: 9724 diff changeset	581 tn = norm (r*s);
9f25290a35e8 more private function and subfunction changes John W. Eaton <jwe@octave.org> parents: 9724 diff changeset	582 snm = (sn / tn) / tn;
9f25290a35e8 more private function and subfunction changes John W. Eaton <jwe@octave.org> parents: 9724 diff changeset	583 if (snm < delta)
10549 95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 10201 diff changeset	584 ## Get the dogleg path minimizer.
9899 9f25290a35e8 more private function and subfunction changes John W. Eaton <jwe@octave.org> parents: 9724 diff changeset	585 bn = norm (b);
9f25290a35e8 more private function and subfunction changes John W. Eaton <jwe@octave.org> parents: 9724 diff changeset	586 dxn = delta/xn; snmd = snm/delta;
9f25290a35e8 more private function and subfunction changes John W. Eaton <jwe@octave.org> parents: 9724 diff changeset	587 t = (bn/sn) * (bn/xn) * snmd;
9f25290a35e8 more private function and subfunction changes John W. Eaton <jwe@octave.org> parents: 9724 diff changeset	588 t -= dxn * snmd^2 - sqrt ((t-dxn)^2 + (1-dxn^2)*(1-snmd^2));
9f25290a35e8 more private function and subfunction changes John W. Eaton <jwe@octave.org> parents: 9724 diff changeset	589 alpha = dxn*(1-snmd^2) / t;
9f25290a35e8 more private function and subfunction changes John W. Eaton <jwe@octave.org> parents: 9724 diff changeset	590 else
9f25290a35e8 more private function and subfunction changes John W. Eaton <jwe@octave.org> parents: 9724 diff changeset	591 alpha = 0;
9f25290a35e8 more private function and subfunction changes John W. Eaton <jwe@octave.org> parents: 9724 diff changeset	592 endif
9f25290a35e8 more private function and subfunction changes John W. Eaton <jwe@octave.org> parents: 9724 diff changeset	593 else
9f25290a35e8 more private function and subfunction changes John W. Eaton <jwe@octave.org> parents: 9724 diff changeset	594 alpha = delta / xn;
9f25290a35e8 more private function and subfunction changes John W. Eaton <jwe@octave.org> parents: 9724 diff changeset	595 snm = 0;
9f25290a35e8 more private function and subfunction changes John W. Eaton <jwe@octave.org> parents: 9724 diff changeset	596 endif
9f25290a35e8 more private function and subfunction changes John W. Eaton <jwe@octave.org> parents: 9724 diff changeset	597 ## Form the appropriate convex combination.
9f25290a35e8 more private function and subfunction changes John W. Eaton <jwe@octave.org> parents: 9724 diff changeset	598 x = alpha * x + ((1-alpha) * min (snm, delta)) * s;
9f25290a35e8 more private function and subfunction changes John W. Eaton <jwe@octave.org> parents: 9724 diff changeset	599 endif
9f25290a35e8 more private function and subfunction changes John W. Eaton <jwe@octave.org> parents: 9724 diff changeset	600 endfunction
9f25290a35e8 more private function and subfunction changes John W. Eaton <jwe@octave.org> parents: 9724 diff changeset	601

Mercurial > hg > octave-nkf

annotate scripts/optimization/fsolve.m @ 10687:a8ce6bdecce5