1 /*
2 ////////////////////////////////////////////////////////////////////////////////////
3 //
4 //  Prototypes and definitions for the Levenberg - Marquardt minimization algorithm
5 //  Copyright (C) 2004  Manolis Lourakis (lourakis at ics forth gr)
6 //  Institute of Computer Science, Foundation for Research & Technology - Hellas
7 //  Heraklion, Crete, Greece.
8 //
9 //  This program is free software; you can redistribute it and/or modify
10 //  it under the terms of the GNU General Public License as published by
11 //  the Free Software Foundation; either version 2 of the License, or
12 //  (at your option) any later version.
13 //
14 //  This program is distributed in the hope that it will be useful,
15 //  but WITHOUT ANY WARRANTY; without even the implied warranty of
16 //  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
17 //  GNU General Public License for more details.
18 //
19 ////////////////////////////////////////////////////////////////////////////////////
20 */
21 
22 #ifndef _LEVMAR_H_
23 #define _LEVMAR_H_
24 
25 /************************************* Start of configuration options *************************************/
26 /* Note that when compiling with CMake, this configuration section is automatically generated
27  * based on the user's input, see levmar.h.in
28  */
29 
30 /* specifies whether to use LAPACK or not. Using LAPACK is strongly recommended */
31 #define HAVE_LAPACK
32 
33 /* specifies whether the PLASMA parallel library for multicore CPUs is available */
34 /* #undef HAVE_PLASMA */
35 
36 /* to avoid the overhead of repeated mallocs(), routines in Axb.c can be instructed to
37  * retain working memory between calls. Such a choice, however, renders these routines
38  * non-reentrant and is not safe in a shared memory multiprocessing environment.
39  * Bellow, an attempt is made to issue a warning if this option is turned on and OpenMP
40  * is being used (note that this will work only if omp.h is included before levmar.h)
41  */
42 #define LINSOLVERS_RETAIN_MEMORY
43 #if (defined(_OPENMP))
44 # ifdef LINSOLVERS_RETAIN_MEMORY
45 #  ifdef _MSC_VER
46 #  pragma message("LINSOLVERS_RETAIN_MEMORY is not safe in a multithreaded environment and should be turned off!")
47 #  else
48 #  warning LINSOLVERS_RETAIN_MEMORY is not safe in a multithreaded environment and should be turned off!
49 #  endif /* _MSC_VER */
50 # endif /* LINSOLVERS_RETAIN_MEMORY */
51 #endif /* _OPENMP */
52 
53 /* specifies whether double precision routines will be compiled or not */
54 #define LM_DBL_PREC
55 /* specifies whether single precision routines will be compiled or not */
56 #define LM_SNGL_PREC
57 
58 /****************** End of configuration options, no changes necessary beyond this point ******************/
59 
60 
61 #ifdef __cplusplus
62 extern "C" {
63 #endif
64 
65 /* work arrays size for ?levmar_der and ?levmar_dif functions.
66  * should be multiplied by sizeof(double) or sizeof(float) to be converted to bytes
67  */
68 #define LM_DER_WORKSZ(npar, nmeas) (2*(nmeas) + 4*(npar) + (nmeas)*(npar) + (npar)*(npar))
69 #define LM_DIF_WORKSZ(npar, nmeas) (4*(nmeas) + 4*(npar) + (nmeas)*(npar) + (npar)*(npar))
70 
71 /* work arrays size for ?levmar_bc_der and ?levmar_bc_dif functions.
72  * should be multiplied by sizeof(double) or sizeof(float) to be converted to bytes
73  */
74 #define LM_BC_DER_WORKSZ(npar, nmeas) (2*(nmeas) + 4*(npar) + (nmeas)*(npar) + (npar)*(npar))
75 #define LM_BC_DIF_WORKSZ(npar, nmeas) LM_BC_DER_WORKSZ((npar), (nmeas)) /* LEVMAR_BC_DIF currently implemented using LEVMAR_BC_DER()! */
76 
77 /* work arrays size for ?levmar_lec_der and ?levmar_lec_dif functions.
78  * should be multiplied by sizeof(double) or sizeof(float) to be converted to bytes
79  */
80 #define LM_LEC_DER_WORKSZ(npar, nmeas, nconstr) LM_DER_WORKSZ((npar)-(nconstr), (nmeas))
81 #define LM_LEC_DIF_WORKSZ(npar, nmeas, nconstr) LM_DIF_WORKSZ((npar)-(nconstr), (nmeas))
82 
83 /* work arrays size for ?levmar_blec_der and ?levmar_blec_dif functions.
84  * should be multiplied by sizeof(double) or sizeof(float) to be converted to bytes
85  */
86 #define LM_BLEC_DER_WORKSZ(npar, nmeas, nconstr) LM_LEC_DER_WORKSZ((npar), (nmeas)+(npar), (nconstr))
87 #define LM_BLEC_DIF_WORKSZ(npar, nmeas, nconstr) LM_LEC_DIF_WORKSZ((npar), (nmeas)+(npar), (nconstr))
88 
89 /* work arrays size for ?levmar_bleic_der and ?levmar_bleic_dif functions.
90  * should be multiplied by sizeof(double) or sizeof(float) to be converted to bytes
91  */
92 #define LM_BLEIC_DER_WORKSZ(npar, nmeas, nconstr1, nconstr2) LM_BLEC_DER_WORKSZ((npar)+(nconstr2), (nmeas)+(nconstr2), (nconstr1)+(nconstr2))
93 #define LM_BLEIC_DIF_WORKSZ(npar, nmeas, nconstr1, nconstr2) LM_BLEC_DIF_WORKSZ((npar)+(nconstr2), (nmeas)+(nconstr2), (nconstr1)+(nconstr2))
94 
95 #define LM_OPTS_SZ    	 5 /* max(4, 5) */
96 #define LM_INFO_SZ    	 10
97 #define LM_ERROR         -1
98 #define LM_INIT_MU    	 1E-03
99 #define LM_STOP_THRESH	 1E-17
100 #define LM_DIFF_DELTA    1E-06
101 #define LM_VERSION       "2.6 (November 2011)"
102 
103 #ifdef LM_DBL_PREC
104 /* double precision LM, with & without Jacobian */
105 /* unconstrained minimization */
106 extern int dlevmar_der(
107       void (*func)(double *p, double *hx, int m, int n, void *adata),
108       void (*jacf)(double *p, double *j, int m, int n, void *adata),
109       double *p, double *x, int m, int n, int itmax, double *opts,
110       double *info, double *work, double *covar, void *adata);
111 
112 extern int dlevmar_dif(
113       void (*func)(double *p, double *hx, int m, int n, void *adata),
114       double *p, double *x, int m, int n, int itmax, double *opts,
115       double *info, double *work, double *covar, void *adata);
116 
117 /* box-constrained minimization */
118 extern int dlevmar_bc_der(
119        void (*func)(double *p, double *hx, int m, int n, void *adata),
120        void (*jacf)(double *p, double *j, int m, int n, void *adata),
121        double *p, double *x, int m, int n, double *lb, double *ub, double *dscl,
122        int itmax, double *opts, double *info, double *work, double *covar, void *adata);
123 
124 extern int dlevmar_bc_dif(
125        void (*func)(double *p, double *hx, int m, int n, void *adata),
126        double *p, double *x, int m, int n, double *lb, double *ub, double *dscl,
127        int itmax, double *opts, double *info, double *work, double *covar, void *adata);
128 
129 #ifdef HAVE_LAPACK
130 /* linear equation constrained minimization */
131 extern int dlevmar_lec_der(
132       void (*func)(double *p, double *hx, int m, int n, void *adata),
133       void (*jacf)(double *p, double *j, int m, int n, void *adata),
134       double *p, double *x, int m, int n, double *A, double *b, int k,
135       int itmax, double *opts, double *info, double *work, double *covar, void *adata);
136 
137 extern int dlevmar_lec_dif(
138       void (*func)(double *p, double *hx, int m, int n, void *adata),
139       double *p, double *x, int m, int n, double *A, double *b, int k,
140       int itmax, double *opts, double *info, double *work, double *covar, void *adata);
141 
142 /* box & linear equation constrained minimization */
143 extern int dlevmar_blec_der(
144       void (*func)(double *p, double *hx, int m, int n, void *adata),
145       void (*jacf)(double *p, double *j, int m, int n, void *adata),
146       double *p, double *x, int m, int n, double *lb, double *ub, double *A, double *b, int k, double *wghts,
147       int itmax, double *opts, double *info, double *work, double *covar, void *adata);
148 
149 extern int dlevmar_blec_dif(
150       void (*func)(double *p, double *hx, int m, int n, void *adata),
151       double *p, double *x, int m, int n, double *lb, double *ub, double *A, double *b, int k, double *wghts,
152       int itmax, double *opts, double *info, double *work, double *covar, void *adata);
153 
154 /* box, linear equations & inequalities constrained minimization */
155 extern int dlevmar_bleic_der(
156       void (*func)(double *p, double *hx, int m, int n, void *adata),
157       void (*jacf)(double *p, double *j, int m, int n, void *adata),
158       double *p, double *x, int m, int n, double *lb, double *ub,
159       double *A, double *b, int k1, double *C, double *d, int k2,
160       int itmax, double *opts, double *info, double *work, double *covar, void *adata);
161 
162 extern int dlevmar_bleic_dif(
163       void (*func)(double *p, double *hx, int m, int n, void *adata),
164       double *p, double *x, int m, int n, double *lb, double *ub,
165       double *A, double *b, int k1, double *C, double *d, int k2,
166       int itmax, double *opts, double *info, double *work, double *covar, void *adata);
167 
168 /* box & linear inequality constraints */
169 extern int dlevmar_blic_der(
170       void (*func)(double *p, double *hx, int m, int n, void *adata),
171       void (*jacf)(double *p, double *j, int m, int n, void *adata),
172       double *p, double *x, int m, int n, double *lb, double *ub, double *C, double *d, int k2,
173       int itmax, double opts[4], double info[LM_INFO_SZ], double *work, double *covar, void *adata);
174 
175 extern int dlevmar_blic_dif(
176       void (*func)(double *p, double *hx, int m, int n, void *adata),
177       double *p, double *x, int m, int n, double *lb, double *ub, double *C, double *d, int k2,
178       int itmax, double opts[5], double info[LM_INFO_SZ], double *work, double *covar, void *adata);
179 
180 /* linear equation & inequality constraints */
181 extern int dlevmar_leic_der(
182       void (*func)(double *p, double *hx, int m, int n, void *adata),
183       void (*jacf)(double *p, double *j, int m, int n, void *adata),
184       double *p, double *x, int m, int n, double *A, double *b, int k1, double *C, double *d, int k2,
185       int itmax, double opts[4], double info[LM_INFO_SZ], double *work, double *covar, void *adata);
186 
187 extern int dlevmar_leic_dif(
188       void (*func)(double *p, double *hx, int m, int n, void *adata),
189       double *p, double *x, int m, int n, double *A, double *b, int k1, double *C, double *d, int k2,
190       int itmax, double opts[5], double info[LM_INFO_SZ], double *work, double *covar, void *adata);
191 
192 /* linear inequality constraints */
193 extern int dlevmar_lic_der(
194       void (*func)(double *p, double *hx, int m, int n, void *adata),
195       void (*jacf)(double *p, double *j, int m, int n, void *adata),
196       double *p, double *x, int m, int n, double *C, double *d, int k2,
197       int itmax, double opts[4], double info[LM_INFO_SZ], double *work, double *covar, void *adata);
198 
199 extern int dlevmar_lic_dif(
200       void (*func)(double *p, double *hx, int m, int n, void *adata),
201       double *p, double *x, int m, int n, double *C, double *d, int k2,
202       int itmax, double opts[5], double info[LM_INFO_SZ], double *work, double *covar, void *adata);
203 #endif /* HAVE_LAPACK */
204 
205 #endif /* LM_DBL_PREC */
206 
207 
208 #ifdef LM_SNGL_PREC
209 /* single precision LM, with & without Jacobian */
210 /* unconstrained minimization */
211 extern int slevmar_der(
212       void (*func)(float *p, float *hx, int m, int n, void *adata),
213       void (*jacf)(float *p, float *j, int m, int n, void *adata),
214       float *p, float *x, int m, int n, int itmax, float *opts,
215       float *info, float *work, float *covar, void *adata);
216 
217 extern int slevmar_dif(
218       void (*func)(float *p, float *hx, int m, int n, void *adata),
219       float *p, float *x, int m, int n, int itmax, float *opts,
220       float *info, float *work, float *covar, void *adata);
221 
222 /* box-constrained minimization */
223 extern int slevmar_bc_der(
224        void (*func)(float *p, float *hx, int m, int n, void *adata),
225        void (*jacf)(float *p, float *j, int m, int n, void *adata),
226        float *p, float *x, int m, int n, float *lb, float *ub, float *dscl,
227        int itmax, float *opts, float *info, float *work, float *covar, void *adata);
228 
229 extern int slevmar_bc_dif(
230        void (*func)(float *p, float *hx, int m, int n, void *adata),
231        float *p, float *x, int m, int n, float *lb, float *ub, float *dscl,
232        int itmax, float *opts, float *info, float *work, float *covar, void *adata);
233 
234 #ifdef HAVE_LAPACK
235 /* linear equation constrained minimization */
236 extern int slevmar_lec_der(
237       void (*func)(float *p, float *hx, int m, int n, void *adata),
238       void (*jacf)(float *p, float *j, int m, int n, void *adata),
239       float *p, float *x, int m, int n, float *A, float *b, int k,
240       int itmax, float *opts, float *info, float *work, float *covar, void *adata);
241 
242 extern int slevmar_lec_dif(
243       void (*func)(float *p, float *hx, int m, int n, void *adata),
244       float *p, float *x, int m, int n, float *A, float *b, int k,
245       int itmax, float *opts, float *info, float *work, float *covar, void *adata);
246 
247 /* box & linear equation constrained minimization */
248 extern int slevmar_blec_der(
249       void (*func)(float *p, float *hx, int m, int n, void *adata),
250       void (*jacf)(float *p, float *j, int m, int n, void *adata),
251       float *p, float *x, int m, int n, float *lb, float *ub, float *A, float *b, int k, float *wghts,
252       int itmax, float *opts, float *info, float *work, float *covar, void *adata);
253 
254 extern int slevmar_blec_dif(
255       void (*func)(float *p, float *hx, int m, int n, void *adata),
256       float *p, float *x, int m, int n, float *lb, float *ub, float *A, float *b, int k, float *wghts,
257       int itmax, float *opts, float *info, float *work, float *covar, void *adata);
258 
259 /* box, linear equations & inequalities constrained minimization */
260 extern int slevmar_bleic_der(
261       void (*func)(float *p, float *hx, int m, int n, void *adata),
262       void (*jacf)(float *p, float *j, int m, int n, void *adata),
263       float *p, float *x, int m, int n, float *lb, float *ub,
264       float *A, float *b, int k1, float *C, float *d, int k2,
265       int itmax, float *opts, float *info, float *work, float *covar, void *adata);
266 
267 extern int slevmar_bleic_dif(
268       void (*func)(float *p, float *hx, int m, int n, void *adata),
269       float *p, float *x, int m, int n, float *lb, float *ub,
270       float *A, float *b, int k1, float *C, float *d, int k2,
271       int itmax, float *opts, float *info, float *work, float *covar, void *adata);
272 
273 /* box & linear inequality constraints */
274 extern int slevmar_blic_der(
275       void (*func)(float *p, float *hx, int m, int n, void *adata),
276       void (*jacf)(float *p, float *j, int m, int n, void *adata),
277       float *p, float *x, int m, int n, float *lb, float *ub, float *C, float *d, int k2,
278       int itmax, float opts[4], float info[LM_INFO_SZ], float *work, float *covar, void *adata);
279 
280 extern int slevmar_blic_dif(
281       void (*func)(float *p, float *hx, int m, int n, void *adata),
282       float *p, float *x, int m, int n, float *lb, float *ub, float *C, float *d, int k2,
283       int itmax, float opts[5], float info[LM_INFO_SZ], float *work, float *covar, void *adata);
284 
285 /* linear equality & inequality constraints */
286 extern int slevmar_leic_der(
287       void (*func)(float *p, float *hx, int m, int n, void *adata),
288       void (*jacf)(float *p, float *j, int m, int n, void *adata),
289       float *p, float *x, int m, int n, float *A, float *b, int k1, float *C, float *d, int k2,
290       int itmax, float opts[4], float info[LM_INFO_SZ], float *work, float *covar, void *adata);
291 
292 extern int slevmar_leic_dif(
293       void (*func)(float *p, float *hx, int m, int n, void *adata),
294       float *p, float *x, int m, int n, float *A, float *b, int k1, float *C, float *d, int k2,
295       int itmax, float opts[5], float info[LM_INFO_SZ], float *work, float *covar, void *adata);
296 
297 /* linear inequality constraints */
298 extern int slevmar_lic_der(
299       void (*func)(float *p, float *hx, int m, int n, void *adata),
300       void (*jacf)(float *p, float *j, int m, int n, void *adata),
301       float *p, float *x, int m, int n, float *C, float *d, int k2,
302       int itmax, float opts[4], float info[LM_INFO_SZ], float *work, float *covar, void *adata);
303 
304 extern int slevmar_lic_dif(
305       void (*func)(float *p, float *hx, int m, int n, void *adata),
306       float *p, float *x, int m, int n, float *C, float *d, int k2,
307       int itmax, float opts[5], float info[LM_INFO_SZ], float *work, float *covar, void *adata);
308 #endif /* HAVE_LAPACK */
309 
310 #endif /* LM_SNGL_PREC */
311 
312 /* linear system solvers */
313 #ifdef HAVE_LAPACK
314 
315 #ifdef LM_DBL_PREC
316 extern int dAx_eq_b_QR(double *A, double *B, double *x, int m);
317 extern int dAx_eq_b_QRLS(double *A, double *B, double *x, int m, int n);
318 extern int dAx_eq_b_Chol(double *A, double *B, double *x, int m);
319 extern int dAx_eq_b_LU(double *A, double *B, double *x, int m);
320 extern int dAx_eq_b_SVD(double *A, double *B, double *x, int m);
321 extern int dAx_eq_b_BK(double *A, double *B, double *x, int m);
322 #endif /* LM_DBL_PREC */
323 
324 #ifdef LM_SNGL_PREC
325 extern int sAx_eq_b_QR(float *A, float *B, float *x, int m);
326 extern int sAx_eq_b_QRLS(float *A, float *B, float *x, int m, int n);
327 extern int sAx_eq_b_Chol(float *A, float *B, float *x, int m);
328 extern int sAx_eq_b_LU(float *A, float *B, float *x, int m);
329 extern int sAx_eq_b_SVD(float *A, float *B, float *x, int m);
330 extern int sAx_eq_b_BK(float *A, float *B, float *x, int m);
331 #endif /* LM_SNGL_PREC */
332 
333 #else /* no LAPACK */
334 
335 #ifdef LM_DBL_PREC
336 extern int dAx_eq_b_LU_noLapack(double *A, double *B, double *x, int n);
337 #endif /* LM_DBL_PREC */
338 
339 #ifdef LM_SNGL_PREC
340 extern int sAx_eq_b_LU_noLapack(float *A, float *B, float *x, int n);
341 #endif /* LM_SNGL_PREC */
342 
343 #endif /* HAVE_LAPACK */
344 
345 #ifdef HAVE_PLASMA
346 #ifdef LM_DBL_PREC
347 extern int dAx_eq_b_PLASMA_Chol(double *A, double *B, double *x, int m);
348 #endif
349 #ifdef LM_SNGL_PREC
350 extern int sAx_eq_b_PLASMA_Chol(float *A, float *B, float *x, int m);
351 #endif
352 extern void levmar_PLASMA_setnbcores(int cores);
353 #endif /* HAVE_PLASMA */
354 
355 /* Jacobian verification, double & single precision */
356 #ifdef LM_DBL_PREC
357 extern void dlevmar_chkjac(
358     void (*func)(double *p, double *hx, int m, int n, void *adata),
359     void (*jacf)(double *p, double *j, int m, int n, void *adata),
360     double *p, int m, int n, void *adata, double *err);
361 #endif /* LM_DBL_PREC */
362 
363 #ifdef LM_SNGL_PREC
364 extern void slevmar_chkjac(
365     void (*func)(float *p, float *hx, int m, int n, void *adata),
366     void (*jacf)(float *p, float *j, int m, int n, void *adata),
367     float *p, int m, int n, void *adata, float *err);
368 #endif /* LM_SNGL_PREC */
369 
370 /* miscellaneous: standard deviation, coefficient of determination (R2),
371  *                Pearson's correlation coefficient for best-fit parameters
372  */
373 #ifdef LM_DBL_PREC
374 extern double dlevmar_stddev( double *covar, int m, int i);
375 extern double dlevmar_corcoef(double *covar, int m, int i, int j);
376 extern double dlevmar_R2(void (*func)(double *p, double *hx, int m, int n, void *adata), double *p, double *x, int m, int n, void *adata);
377 
378 #endif /* LM_DBL_PREC */
379 
380 #ifdef LM_SNGL_PREC
381 extern float slevmar_stddev( float *covar, int m, int i);
382 extern float slevmar_corcoef(float *covar, int m, int i, int j);
383 extern float slevmar_R2(void (*func)(float *p, float *hx, int m, int n, void *adata), float *p, float *x, int m, int n, void *adata);
384 
385 extern void slevmar_locscale(
386         void (*func)(float *p, float *hx, int m, int n, void *adata),
387         float *p, float *x, int m, int n, void *adata,
388         int howto, float locscl[2], float **residptr);
389 
390 extern int slevmar_outlid(float *r, int n, float thresh, float ls[2], char *outlmap);
391 
392 #endif /* LM_SNGL_PREC */
393 
394 #ifdef __cplusplus
395 }
396 #endif
397 
398 #endif /* _LEVMAR_H_ */
399