File: Nonlinear-Least_002dSquares-Iteration.html

package info (click to toggle)
gsl-ref-html 2.3-1
  • links: PTS
  • area: non-free
  • in suites: bullseye, buster, sid
  • size: 6,876 kB
  • ctags: 4,574
  • sloc: makefile: 35
file content (184 lines) | stat: -rw-r--r-- 10,675 bytes parent folder | download
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
<!DOCTYPE html PUBLIC "-//W3C//DTD HTML 4.01 Transitional//EN" "http://www.w3.org/TR/html4/loose.dtd">
<html>
<!-- Copyright (C) 1996, 1997, 1998, 1999, 2000, 2001, 2002, 2003, 2004, 2005, 2006, 2007, 2008, 2009, 2010, 2011, 2012, 2013, 2014, 2015, 2016 The GSL Team.

Permission is granted to copy, distribute and/or modify this document
under the terms of the GNU Free Documentation License, Version 1.3 or
any later version published by the Free Software Foundation; with the
Invariant Sections being "GNU General Public License" and "Free Software
Needs Free Documentation", the Front-Cover text being "A GNU Manual",
and with the Back-Cover Text being (a) (see below). A copy of the
license is included in the section entitled "GNU Free Documentation
License".

(a) The Back-Cover Text is: "You have the freedom to copy and modify this
GNU Manual." -->
<!-- Created by GNU Texinfo 5.1, http://www.gnu.org/software/texinfo/ -->
<head>
<title>GNU Scientific Library &ndash; Reference Manual: Nonlinear Least-Squares Iteration</title>

<meta name="description" content="GNU Scientific Library &ndash; Reference Manual: Nonlinear Least-Squares Iteration">
<meta name="keywords" content="GNU Scientific Library &ndash; Reference Manual: Nonlinear Least-Squares Iteration">
<meta name="resource-type" content="document">
<meta name="distribution" content="global">
<meta name="Generator" content="makeinfo">
<meta http-equiv="Content-Type" content="text/html; charset=utf-8">
<link href="index.html#Top" rel="start" title="Top">
<link href="Function-Index.html#Function-Index" rel="index" title="Function Index">
<link href="Nonlinear-Least_002dSquares-Fitting.html#Nonlinear-Least_002dSquares-Fitting" rel="up" title="Nonlinear Least-Squares Fitting">
<link href="Nonlinear-Least_002dSquares-Testing-for-Convergence.html#Nonlinear-Least_002dSquares-Testing-for-Convergence" rel="next" title="Nonlinear Least-Squares Testing for Convergence">
<link href="Nonlinear-Least_002dSquares-Function-Definition.html#Nonlinear-Least_002dSquares-Function-Definition" rel="previous" title="Nonlinear Least-Squares Function Definition">
<style type="text/css">
<!--
a.summary-letter {text-decoration: none}
blockquote.smallquotation {font-size: smaller}
div.display {margin-left: 3.2em}
div.example {margin-left: 3.2em}
div.indentedblock {margin-left: 3.2em}
div.lisp {margin-left: 3.2em}
div.smalldisplay {margin-left: 3.2em}
div.smallexample {margin-left: 3.2em}
div.smallindentedblock {margin-left: 3.2em; font-size: smaller}
div.smalllisp {margin-left: 3.2em}
kbd {font-style:oblique}
pre.display {font-family: inherit}
pre.format {font-family: inherit}
pre.menu-comment {font-family: serif}
pre.menu-preformatted {font-family: serif}
pre.smalldisplay {font-family: inherit; font-size: smaller}
pre.smallexample {font-size: smaller}
pre.smallformat {font-family: inherit; font-size: smaller}
pre.smalllisp {font-size: smaller}
span.nocodebreak {white-space:nowrap}
span.nolinebreak {white-space:nowrap}
span.roman {font-family:serif; font-weight:normal}
span.sansserif {font-family:sans-serif; font-weight:normal}
ul.no-bullet {list-style: none}
-->
</style>


</head>

<body lang="en" bgcolor="#FFFFFF" text="#000000" link="#0000FF" vlink="#800080" alink="#FF0000">
<a name="Nonlinear-Least_002dSquares-Iteration"></a>
<div class="header">
<p>
Next: <a href="Nonlinear-Least_002dSquares-Testing-for-Convergence.html#Nonlinear-Least_002dSquares-Testing-for-Convergence" accesskey="n" rel="next">Nonlinear Least-Squares Testing for Convergence</a>, Previous: <a href="Nonlinear-Least_002dSquares-Function-Definition.html#Nonlinear-Least_002dSquares-Function-Definition" accesskey="p" rel="previous">Nonlinear Least-Squares Function Definition</a>, Up: <a href="Nonlinear-Least_002dSquares-Fitting.html#Nonlinear-Least_002dSquares-Fitting" accesskey="u" rel="up">Nonlinear Least-Squares Fitting</a> &nbsp; [<a href="Function-Index.html#Function-Index" title="Index" rel="index">Index</a>]</p>
</div>
<hr>
<a name="Iteration-4"></a>
<h3 class="section">39.7 Iteration</h3>

<p>The following functions drive the iteration of each algorithm.  Each
function performs one iteration of the trust region method and updates
the state of the solver.
</p>
<dl>
<dt><a name="index-gsl_005fmultifit_005fnlinear_005fiterate"></a>Function: <em>int</em> <strong>gsl_multifit_nlinear_iterate</strong> <em>(gsl_multifit_nlinear_workspace * <var>w</var>)</em></dt>
<dt><a name="index-gsl_005fmultilarge_005fnlinear_005fiterate"></a>Function: <em>int</em> <strong>gsl_multilarge_nlinear_iterate</strong> <em>(gsl_multilarge_nlinear_workspace * <var>w</var>)</em></dt>
<dd><p>These functions perform a single iteration of the solver <var>w</var>.  If
the iteration encounters an unexpected problem then an error code will
be returned.  The solver workspace maintains a current estimate of the
best-fit parameters at all times.
</p></dd></dl>

<p>The solver workspace <var>w</var> contains the following entries, which can
be used to track the progress of the solution:
</p>
<dl compact="compact">
<dt><code>gsl_vector * x</code></dt>
<dd><p>The current position, length <em>p</em>.
</p>
</dd>
<dt><code>gsl_vector * f</code></dt>
<dd><p>The function residual vector at the current position <em>f(x)</em>, length
<em>n</em>.
</p>
</dd>
<dt><code>gsl_matrix * J</code></dt>
<dd><p>The Jacobian matrix at the current position <em>J(x)</em>, size
<em>n</em>-by-<em>p</em> (only for <code>gsl_multifit_nlinear</code> interface).
</p>
</dd>
<dt><code>gsl_vector * dx</code></dt>
<dd><p>The difference between the current position and the previous position,
i.e. the last step <em>\delta</em>, taken as a vector, length <em>p</em>.
</p>
</dd>
</dl>

<p>These quantities can be accessed with the following functions,
</p>
<dl>
<dt><a name="index-gsl_005fmultifit_005fnlinear_005fposition"></a>Function: <em>gsl_vector *</em> <strong>gsl_multifit_nlinear_position</strong> <em>(const gsl_multifit_nlinear_workspace * <var>w</var>)</em></dt>
<dt><a name="index-gsl_005fmultilarge_005fnlinear_005fposition"></a>Function: <em>gsl_vector *</em> <strong>gsl_multilarge_nlinear_position</strong> <em>(const gsl_multilarge_nlinear_workspace * <var>w</var>)</em></dt>
<dd><p>These functions return the current position <em>x</em> (i.e. best-fit
parameters) of the solver <var>w</var>.
</p></dd></dl>

<dl>
<dt><a name="index-gsl_005fmultifit_005fnlinear_005fresidual"></a>Function: <em>gsl_vector *</em> <strong>gsl_multifit_nlinear_residual</strong> <em>(const gsl_multifit_nlinear_workspace * <var>w</var>)</em></dt>
<dt><a name="index-gsl_005fmultilarge_005fnlinear_005fresidual"></a>Function: <em>gsl_vector *</em> <strong>gsl_multilarge_nlinear_residual</strong> <em>(const gsl_multilarge_nlinear_workspace * <var>w</var>)</em></dt>
<dd><p>These functions return the current residual vector <em>f(x)</em> of the
solver <var>w</var>.  For weighted systems, the residual vector includes the
weighting factor <em>\sqrt{W}</em>.
</p></dd></dl>

<dl>
<dt><a name="index-gsl_005fmultifit_005fnlinear_005fjac"></a>Function: <em>gsl_matrix *</em> <strong>gsl_multifit_nlinear_jac</strong> <em>(const gsl_multifit_nlinear_workspace * <var>w</var>)</em></dt>
<dd><p>This function returns a pointer to the <em>n</em>-by-<em>p</em> Jacobian matrix for the
current iteration of the solver <var>w</var>. This function is available only for the
<code>gsl_multifit_nlinear</code> interface.
</p></dd></dl>

<dl>
<dt><a name="index-gsl_005fmultifit_005fnlinear_005fniter"></a>Function: <em>size_t</em> <strong>gsl_multifit_nlinear_niter</strong> <em>(const gsl_multifit_nlinear_workspace * <var>w</var>)</em></dt>
<dt><a name="index-gsl_005fmultilarge_005fnlinear_005fniter"></a>Function: <em>size_t</em> <strong>gsl_multilarge_nlinear_niter</strong> <em>(const gsl_multilarge_nlinear_workspace * <var>w</var>)</em></dt>
<dd><p>These functions return the number of iterations performed so far.
The iteration counter is updated on each call to the
<code>_iterate</code> functions above, and reset to 0 in the
<code>_init</code> functions.
</p></dd></dl>

<dl>
<dt><a name="index-gsl_005fmultifit_005fnlinear_005frcond"></a>Function: <em>int</em> <strong>gsl_multifit_nlinear_rcond</strong> <em>(double * <var>rcond</var>, const gsl_multifit_nlinear_workspace * <var>w</var>)</em></dt>
<dt><a name="index-gsl_005fmultilarge_005fnlinear_005frcond"></a>Function: <em>int</em> <strong>gsl_multilarge_nlinear_rcond</strong> <em>(double * <var>rcond</var>, const gsl_multilarge_nlinear_workspace * <var>w</var>)</em></dt>
<dd><p>This function estimates the reciprocal condition number
of the Jacobian matrix at the current position <em>x</em> and
stores it in <var>rcond</var>. The computed value is only an estimate
to give the user a guideline as to the conditioning of their particular
problem. Its calculation is based on which factorization
method is used (Cholesky, QR, or SVD). 
</p>
<ul>
<li> For the Cholesky solver, the matrix <em>J^T J</em> is factored at each
iteration. Therefore this function will estimate the 1-norm condition number
<em>rcond^2 = 1/(||J^T J||_1 \cdot ||(J^T J)^{-1}||_1)</em>

</li><li> For the QR solver, <em>J</em> is factored as <em>J = Q R</em> at each
iteration. For simplicity, this function calculates the 1-norm conditioning of
only the <em>R</em> factor, <em>rcond = 1 / (||R||_1 \cdot ||R^{-1}||_1)</em>.
This can be computed efficiently since <em>R</em> is upper triangular.

</li><li> For the SVD solver, in order to efficiently solve the trust region
subproblem, the matrix which is factored is <em>J D^{-1}</em>, instead of
<em>J</em> itself. The resulting singular values are used to provide
the 2-norm reciprocal condition number, as <em>rcond = \sigma_{min} / \sigma_{max}</em>.
Note that when using Mor&eacute; scaling, <em>D \ne I</em> and the resulting
<var>rcond</var> estimate may be significantly different from the true
<var>rcond</var> of <em>J</em> itself.

</li></ul>
</dd></dl>

<hr>
<div class="header">
<p>
Next: <a href="Nonlinear-Least_002dSquares-Testing-for-Convergence.html#Nonlinear-Least_002dSquares-Testing-for-Convergence" accesskey="n" rel="next">Nonlinear Least-Squares Testing for Convergence</a>, Previous: <a href="Nonlinear-Least_002dSquares-Function-Definition.html#Nonlinear-Least_002dSquares-Function-Definition" accesskey="p" rel="previous">Nonlinear Least-Squares Function Definition</a>, Up: <a href="Nonlinear-Least_002dSquares-Fitting.html#Nonlinear-Least_002dSquares-Fitting" accesskey="u" rel="up">Nonlinear Least-Squares Fitting</a> &nbsp; [<a href="Function-Index.html#Function-Index" title="Index" rel="index">Index</a>]</p>
</div>



</body>
</html>