This file is indexed.

/usr/share/doc/libdogleg-doc/libdogleg.html is in libdogleg-doc 0.08-3.

This file is owned by root:root, with mode 0o644.

The actual contents of the file can be viewed below.

  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
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
<?xml version="1.0" ?>
<!DOCTYPE html PUBLIC "-//W3C//DTD XHTML 1.0 Strict//EN" "http://www.w3.org/TR/xhtml1/DTD/xhtml1-strict.dtd">
<html xmlns="http://www.w3.org/1999/xhtml">
<head>
<title>libdogleg</title>
<meta http-equiv="content-type" content="text/html; charset=utf-8" />
<link rev="made" href="mailto:root@localhost" />
</head>

<body style="background-color: white">



<ul id="index">
  <li><a href="#NAME">NAME</a></li>
  <li><a href="#DESCRIPTION">DESCRIPTION</a></li>
  <li><a href="#FUNCTIONS-AND-TYPES">FUNCTIONS AND TYPES</a>
    <ul>
      <li><a href="#Main-API">Main API</a>
        <ul>
          <li><a href="#dogleg_optimize">dogleg_optimize</a></li>
          <li><a href="#dogleg_freeContext">dogleg_freeContext</a></li>
          <li><a href="#dogleg_testGradient">dogleg_testGradient</a></li>
          <li><a href="#dogleg_callback_t">dogleg_callback_t</a></li>
          <li><a href="#dogleg_solverContext_t">dogleg_solverContext_t</a></li>
          <li><a href="#dogleg_operatingPoint_t">dogleg_operatingPoint_t</a></li>
        </ul>
      </li>
      <li><a href="#Parameters">Parameters</a>
        <ul>
          <li><a href="#dogleg_setMaxIterations">dogleg_setMaxIterations</a></li>
          <li><a href="#dogleg_setDebug">dogleg_setDebug</a></li>
          <li><a href="#dogleg_setInitialTrustregion">dogleg_setInitialTrustregion</a></li>
          <li><a href="#dogleg_setThresholds">dogleg_setThresholds</a></li>
          <li><a href="#dogleg_setTrustregionUpdateParameters">dogleg_setTrustregionUpdateParameters</a></li>
        </ul>
      </li>
    </ul>
  </li>
  <li><a href="#BUGS">BUGS</a></li>
  <li><a href="#AUTHOR">AUTHOR</a></li>
  <li><a href="#LICENSE-AND-COPYRIGHT">LICENSE AND COPYRIGHT</a></li>
</ul>

<h1 id="NAME">NAME</h1>

<p>libdogleg - A general purpose sparse optimizer to solve data fitting problems, such as sparse bundle adjustment.</p>

<h1 id="DESCRIPTION">DESCRIPTION</h1>

<p>This is a library for solving large-scale nonlinear optimization problems. By employing sparse linear algebra, it is taylored for problems that have weak coupling between the optimization variables. For appropriately sparse problems this results in <i>massive</i> performance gains.</p>

<p>The main task of this library is to find the vector <b>p</b> that minimizes</p>

<p>norm2( <b>x</b> )</p>

<p>where <b>x</b> = <i>f</i>(<b>p</b>) is a vector that has higher dimensionality than <b>p</b>. The user passes in a callback function (of type <code>dogleg_callback_t</code>) that takes in the vector <b>p</b> and returns the vector <b>x</b> and a matrix of derivatives <b>J</b> = d<b>f</b>/d<b>p</b>. <b>J</b> is a matrix with a row for each element of <i>f</i> and a column for each element of <b>p</b>. <b>J</b> is a sparse matrix, which results in substantial increases in computational efficiency if most entries of <b>J</b> are 0. <b>J</b> is stored row-first in the callback routine. libdogleg uses a column-first data representation so it references the transpose of <b>J</b> (called <b>Jt</b>). <b>J</b> stored row-first is identical to <b>Jt</b> stored column-first; this is purely a naming choice.</p>

<p>This library implements Powell&#39;s dog-leg algorithm to solve the problem. Like the more-widely-known Levenberg-Marquardt algorithm, Powell&#39;s dog-leg algorithm solves a nonlinear optimization problem by interpolating between a Gauss-Newton step and a gradient descent step. Improvements over LM are</p>

<ul>

<li><p>a more natural representation of the linearity of the operating point (trust region size vs a vague lambda term).</p>

</li>
<li><p>significant efficiency gains, since a matrix inversion isn&#39;t needed to retry a rejected step</p>

</li>
</ul>

<p>The algorithm is described in many places, originally in</p>

<p>M. Powell. A Hybrid Method for Nonlinear Equations. In P. Rabinowitz, editor, Numerical Methods for Nonlinear Algebraic Equations, pages 87-144. Gordon and Breach Science, London, 1970.</p>

<p>Various enhancements to Powell&#39;s original method are described in the literature; at this time only the original algorithm is implemented here.</p>

<p>The sparse matrix algebra is handled by the CHOLMOD library, written by Tim Davis. Parts of CHOLMOD are licensed under the GPL and parts under the LGPL. Only the LGPL pieces are used here, allowing libdogleg to be licensed under the LGPL as well. Due to this I lose some convenience (all simple sparse matrix arithmetic in CHOLMOD is GPL-ed) and some performance (the fancier computational methods, such as supernodal analysis are GPL-ed). For my current applications the performance losses are minor.</p>

<h1 id="FUNCTIONS-AND-TYPES">FUNCTIONS AND TYPES</h1>

<h2 id="Main-API">Main API</h2>

<h3 id="dogleg_optimize">dogleg_optimize</h3>

<p>This is the main call to the library. Its declared as</p>

<pre><code> double dogleg_optimize(double* p, unsigned int Nstate,
                        unsigned int Nmeas, unsigned int NJnnz,
                        dogleg_callback_t* f, void* cookie,
                        dogleg_solverContext_t** returnContext);</code></pre>

<ul>

<li><p><b>p</b> is the initial estimate of the state vector (and holds the final result)</p>

</li>
<li><p><code>Nstate</code> specifies the number of optimization variables (elements of <b>p</b>)</p>

</li>
<li><p><code>Nmeas</code> specifies the number of measurements (elements of <b>x</b>). <code>Nmeas &gt;= Nstate</code> is a requirement</p>

</li>
<li><p><code>NJnnz</code> specifies the number of non-zero elements of the jacobian matrix d<b>f</b>/d<b>p</b>. In a dense matrix <code>Jnnz = Nstate*Nmeas</code>. We are dealing with sparse jacobians, so <code>NJnnz</code> should be far less. If not, libdogleg is not an appropriate routine to solve this problem.</p>

</li>
<li><p><code>f</code> specifies the callback function that the optimization routine calls to sample the problem being solved</p>

</li>
<li><p><code>cookie</code> is an arbitrary data pointer passed untouched to <code>f</code></p>

</li>
<li><p>If not <code>NULL</code>, <code>returnContext</code> can be used to retrieve the full context structure from the solver. This can be useful since this structure contains the latest operating point values. It also has an active <code>cholmod_common</code> structure, which can be reused if more CHOLMOD routines need to be called externally. <i>If this data is requested, the user is required to free it with <code>dogleg_freeContext</code> when done</i>.</p>

</li>
</ul>

<p><code>dogleg_optimize</code> returns norm2( <b>x</b> ) at the minimum, or a negative value if an error occurred.</p>

<h3 id="dogleg_freeContext">dogleg_freeContext</h3>

<p>Used to deallocate memory used for an optimization cycle. Defined as:</p>

<pre><code> void dogleg_freeContext(dogleg_solverContext_t** ctx);</code></pre>

<p>If a pointer to a context is not requested (by passing <code>returnContext = NULL</code> to <code>dogleg_optimize</code>), libdogleg calls this routine automatically. If the user <i>did</i> retrieve this pointer, though, it must be freed with <code>dogleg_freeContext</code> when the user is finished.</p>

<h3 id="dogleg_testGradient">dogleg_testGradient</h3>

<p>libdogleg requires the user to compute the jacobian matrix <b>J</b>. This is a performance optimization, since <b>J</b> could be computed by differences of <b>x</b>. This optimization is often worth the extra effort, but it creates a possibility that <b>J</b> will have a mistake and <b>J</b> = d<b>f</b>/d<b>p</b> would not be true. To find these types of issues, the user can call</p>

<pre><code> void dogleg_testGradient(unsigned int var, const double* p0,
                          unsigned int Nstate, unsigned int Nmeas, unsigned int NJnnz,
                          dogleg_callback_t* f, void* cookie);</code></pre>

<p>This function computes the jacobian with center differences and compares the result with the jacobian computed by the callback function. It is recommended to do this for every variable while developing the program that uses libdogleg.</p>

<ul>

<li><p><code>var</code> is the index of the variable being tested</p>

</li>
<li><p><code>p0</code> is the state vector <b>p</b> where we&#39;re evaluating the jacobian</p>

</li>
<li><p><code>Nstate</code>, <code>Nmeas</code>, <code>NJnnz</code> are the number of state variables, measurements and non-zero jacobian elements, as before</p>

</li>
<li><p><code>f</code> is the callback function, as before</p>

</li>
<li><p><code>cookie</code> is the user data, as before</p>

</li>
</ul>

<p>This function returns nothing, but prints out the test results.</p>

<h3 id="dogleg_callback_t">dogleg_callback_t</h3>

<p>The main user callback that specifies the optimization problem has type</p>

<pre><code> typedef void (dogleg_callback_t)(const double*   p,
                                  double*         x,
                                  cholmod_sparse* Jt,
                                  void*           cookie);</code></pre>

<ul>

<li><p><b>p</b> is the current state vector</p>

</li>
<li><p><b>x</b> is the resulting <i>f</i>(<b>p</b>)</p>

</li>
<li><p><b>Jt</b> is the transpose of d<b>f</b>/d<b>p</b> at <b>p</b>. As mentioned previously, <b>Jt</b> is stored column-first by CHOLMOD, which can be interpreted as storing <b>J</b> row-first by the user-defined callback routine</p>

</li>
<li><p>The <code>cookie</code> is the user-defined arbitrary data passed into <code>dogleg_optimize</code>.</p>

</li>
</ul>

<h3 id="dogleg_solverContext_t">dogleg_solverContext_t</h3>

<p>This is the solver context that can be retrieved through the <code>returnContext</code> parameter of the <code>dogleg_optimize</code> call. This structure contains <i>all</i> the internal state used by the solver. If requested, the user is responsible for calling <code>dogleg_freeContext</code> when done. This structure is defined as:</p>

<pre><code> typedef struct
 {
   cholmod_common  common;
 
   dogleg_callback_t* f;
   void*              cookie;
 
   // between steps, beforeStep contains the operating point of the last step.
   // afterStep is ONLY used while making the step. Externally, use beforeStep
   // unless you really know what you&#39;re doing
   dogleg_operatingPoint_t* beforeStep;
   dogleg_operatingPoint_t* afterStep;

   // The result of the last JtJ factorization performed. Note that JtJ is not
   // necessarily factorized at every step, so this is NOT guaranteed to contain
   // the factorization of the most recent JtJ
   cholmod_factor*          factorization;
 
   // Have I ever seen a singular JtJ? If so, I add a small constant to the
   // diagonal from that point on. This is a simple and fast way to deal with
   // singularities. This is suboptimal but works for me for now.
   int               wasPositiveSemidefinite;
 } dogleg_solverContext_t;</code></pre>

<p>Some of the members are copies of the data passed into <code>dogleg_optimize</code>; some others are internal state. Of potential interest are</p>

<ul>

<li><p><code>common</code> is a cholmod_common structure used by all CHOLMOD calls. This can be used for any extra CHOLMOD work the user may want to do</p>

</li>
<li><p><code>beforeStep</code> contains the operating point of the optimum solution. The user can analyze this data without the need to re-call the callback routine.</p>

</li>
</ul>

<h3 id="dogleg_operatingPoint_t">dogleg_operatingPoint_t</h3>

<p>An operating point of the solver. This is a part of <code>dogleg_solverContext_t</code>. Various variables describing the operating point such as <b>p</b>, <b>J</b>, <b>x</b>, <b>norm2(x)</b> and <b>Jt x</b> are available. All of the just-mentioned variables are computed at every step and are thus always valid.</p>

<pre><code> // an operating point of the solver
 typedef struct
 {
   double*         p;
   double*         x;
   double          norm2_x;
   cholmod_sparse* Jt;
   double*         Jt_x;
 
   // the cached update vectors. It&#39;s useful to cache these so that when a step is rejected, we can
   // reuse these when we retry
   double*        updateCauchy;
   cholmod_dense* updateGN_cholmoddense;
   double         updateCauchy_lensq, updateGN_lensq; // update vector lengths
 
   // whether the current update vectors are correct or not
   int updateCauchy_valid, updateGN_valid;
 
   int didStepToEdgeOfTrustRegion;
 } dogleg_operatingPoint_t;</code></pre>

<h2 id="Parameters">Parameters</h2>

<p>It is not required to call any of these, but it&#39;s highly recommended to set the initial trust-region size and the termination thresholds to match the problem being solved. Furthermore, it&#39;s highly recommended for the problem being solved to be scaled so that every state variable affects the objective norm2( <b>x</b> ) equally. This makes this method&#39;s concept of &quot;trust region&quot; much more well-defined and makes the termination criteria work correctly.</p>

<h3 id="dogleg_setMaxIterations">dogleg_setMaxIterations</h3>

<p>To set the maximum number of solver iterations, call</p>

<pre><code> void dogleg_setMaxIterations(int n);</code></pre>

<h3 id="dogleg_setDebug">dogleg_setDebug</h3>

<p>To turn on debug output, call</p>

<pre><code> void dogleg_setDebug(int debug);</code></pre>

<p>with a non-zero value for <code>debug</code>. By default, debug output is disabled.</p>

<h3 id="dogleg_setInitialTrustregion">dogleg_setInitialTrustregion</h3>

<p>The optimization method keeps track of a trust region size. Here, the trust region is a ball in R^Nstate. When the method takes a step <b>p</b> -&gt; <b>p + delta_p</b>, it makes sure that</p>

<p><span style="white-space: nowrap;">sqrt( norm2( <b>delta_p</b> ) ) &lt; trust region size</span>.</p>

<p>The initial value of the trust region size can be set with</p>

<pre><code> void dogleg_setInitialTrustregion(double t);</code></pre>

<p>The dogleg algorithm is efficient when recomputing a rejected step for a smaller trust region, so set the initial trust region size to a value larger to a reasonable estimate; the method will quickly shrink the trust region to the correct size.</p>

<h3 id="dogleg_setThresholds">dogleg_setThresholds</h3>

<p>The routine exits when the maximum number of iterations is exceeded, or a termination threshold is hit, whichever happens first. The termination thresholds are all designed to trigger when very slow progress is being made. If all went well, this slow progress is due to us finding the optimum. There are 3 termination thresholds:</p>

<ul>

<li><p>The function being minimized is E = norm2( <b>x</b> ) where <b>x</b> = <i>f</i>(<b>p</b>).</p>

<p>dE/d<b>p</b> = 2*<b>Jt</b>*<b>x</b> where <b>Jt</b> is transpose(d<b>x</b>/d<b>p</b>).</p>

<pre><code> if( for every i  fabs(Jt_x[i]) &lt; JT_X_THRESHOLD )
 { we are done; }</code></pre>

</li>
<li><p>The method takes discrete steps: <b>p</b> -&gt; <b>p + delta_p</b></p>

<pre><code> if( for every i  fabs(delta_p[i]) &lt; UPDATE_THRESHOLD)
 { we are done; }</code></pre>

</li>
<li><p>The method dynamically controls the trust region.</p>

<pre><code> if(trustregion &lt; TRUSTREGION_THRESHOLD)
 { we are done; }</code></pre>

</li>
</ul>

<p>To set these threholds, call</p>

<pre><code> void dogleg_setThresholds(double Jt_x, double update, double trustregion);</code></pre>

<p>To leave a particular threshold alone, specify a negative value.</p>

<h3 id="dogleg_setTrustregionUpdateParameters">dogleg_setTrustregionUpdateParameters</h3>

<p>This function sets the parameters that control when and how the trust region is updated. The default values should work well in most cases, and shouldn&#39;t need to be tweaked.</p>

<p>Declaration looks like</p>

<pre><code> void dogleg_setTrustregionUpdateParameters(double downFactor, double downThreshold,
                                            double upFactor,   double upThreshold);</code></pre>

<p>To see what the parameters do, look at <code>evaluateStep_adjustTrustRegion</code> in the source. Again, these should just work as is.</p>

<h1 id="BUGS">BUGS</h1>

<p>The current implementation doesn&#39;t handle a singular <b>JtJ</b> gracefully (<b>JtJ</b> = <b>Jt</b> * <b>J</b>). Analytically, <b>JtJ</b> is at worst positive semi-definite (has 0 eigenvalues). If a singular <b>JtJ</b> is ever encountered, from that point on, <b>JtJ</b> + lambda*<b>I</b> is inverted instead for some positive constant lambda. This makes the matrix strictly positive definite, but is sloppy. At least I should vary lambda. In my current applications, a singular <b>JtJ</b> only occurs if at a particular operating point the vector <b>x</b> has no dependence at all on some elements of <b>p</b>. In the general case other causes could exist, though.</p>

<p>There&#39;s an inefficiency in that the callback always returns <b>x</b> and <b>J</b>. When I evaluate and reject a step, I do not end up using <b>J</b> at all. Dependng on the callback function, it may be better to ask for <b>x</b> and then, if the step is accepted, to ask for <b>J</b>.</p>

<h1 id="AUTHOR">AUTHOR</h1>

<p>Dima Kogan, <code>&lt;dima@secretsauce.net&gt;</code></p>

<h1 id="LICENSE-AND-COPYRIGHT">LICENSE AND COPYRIGHT</h1>

<p>Copyright 2011 Oblong Industries</p>

<p>This program is free software: you can redistribute it and/or modify it under the terms of the GNU Lesser General Public License as published by the Free Software Foundation, either version 3 of the License, or (at your option) any later version.</p>

<p>This program is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License for more details.</p>

<p>The full text of the license is available at <a href="http://www.gnu.org/licenses">http://www.gnu.org/licenses</a></p>


</body>

</html>