30 #define omp_get_max_threads() 1 56 static void lm_initialize_control(lm_control_type * control);
59 static double lm_enorm(
int,
double *);
62 static void lm_minimize(
int m_dat,
int n_par,
double *par,
63 void (*evaluate) (
double *par,
int m_dat,
double *fvec,
64 void *data,
int *info),
65 void (*printout) (
const char *func,
int n_par,
66 double *par,
int m_dat,
67 double *fvec,
void *data,
int iflag,
69 void *data, lm_control_type * control);
76 static void lm_lmdif(
int m,
int n,
double *x,
double *fvec,
double ftol,
77 double xtol,
double gtol,
int maxfev,
double epsfcn,
78 double *diag,
int mode,
double factor,
int *info,
79 int *nfev,
double *fjac,
int *ipvt,
double *qtf,
80 double *wa1,
double *wa2,
double *wa3,
double *wa4,
81 void (*evaluate) (
double *par,
int m_dat,
double *fvec,
82 void *data,
int *info),
83 void (*printout) (
const char *func,
int n_par,
84 double *par,
int m_dat,
85 double *fvec,
void *data,
int iflag,
99 #include "muse_optimize.h" 100 #include "muse_utils.h" 110 cpl_error_code retval;
111 } muse_cpl_optimize_eval_struct;
120 lm_evaluate_cpl (
double *par,
int m_dat,
double *fvec,
void *data,
int *info) {
121 muse_cpl_optimize_eval_struct *s = data;
123 cpl_array *cpl_par = cpl_array_wrap_double(par, s->npar);
124 cpl_array *cpl_fvec = cpl_array_wrap_double(fvec, m_dat);
126 s->retval = (*(s->func))(s->data, cpl_par, cpl_fvec);
128 cpl_array_unwrap(cpl_par);
129 cpl_array_unwrap(cpl_fvec);
131 if (s->retval != CPL_ERROR_NONE) {
136 static const cpl_error_code error_code_tbl[] = {
137 CPL_ERROR_ILLEGAL_INPUT,
141 CPL_ERROR_SINGULAR_MATRIX,
146 CPL_ERROR_ILLEGAL_OUTPUT,
150 static const char *lm_infmsg[] = {
151 "fatal coding error (improper input parameters)",
152 "success (the relative error in the sum of squares is at most tol)",
153 "success (the relative error between x and the solution is at most tol)",
154 "success (both errors are at most tol)",
155 "trapped by degeneracy (fvec is orthogonal to the columns of the jacobian)" 156 "timeout (number of calls to fcn has reached maxcall*(n+1))",
157 "failure (ftol<tol: cannot reduce sum of squares any further)",
158 "failure (xtol<tol: cannot improve approximate solution any further)",
159 "failure (gtol<tol: cannot improve approximate solution any further)",
160 "exception (not enough memory)",
161 "exception (break requested within function evaluation)" 178 int m_dat,
double *fvec,
179 void *data,
int iflag,
int iter,
int nfev) {
181 cpl_msg_debug(func,
"trying step in gradient direction");
182 }
else if (iflag == 1) {
183 cpl_msg_debug(func,
"determining gradient (iteration %d)", iter);
184 }
else if (iflag == 0) {
185 cpl_msg_debug(func,
"starting minimization");
186 }
else if (iflag == -1) {
187 cpl_msg_debug(func,
"terminated after %d evaluations", nfev);
190 char *parstring = cpl_calloc(n_par * 15 + 30,
sizeof(
char));
191 snprintf(parstring, 5,
"par:");
193 for (i = 0; i < n_par; ++i) {
194 snprintf(parstring + strlen(parstring), 15,
" %7.3g", par[i]);
196 snprintf(parstring + strlen(parstring), 25,
" => norm: %7g",
197 lm_enorm(m_dat, fvec));
198 cpl_msg_debug(func,
"%s", parstring);
203 muse_cpl_optimize_eval_struct *mydata
204 = (muse_cpl_optimize_eval_struct *) data;
207 cpl_msg_debug(func,
" fitting data as follows:", mydata);
208 for (i = 0; i < m_dat; ++i) {
209 t = (mydata->tvec)[i];
210 y = (mydata->yvec)[i];
211 f = mydata->f(t, par);
212 cpl_msg_debug(func,
" t[%2d]=%12g y=%12g fit=%12g residue=%12g",
217 UNUSED_ARGUMENT(data);
255 cpl_ensure_code(aPar != NULL, CPL_ERROR_NULL_INPUT);
256 int npars = cpl_array_get_size(aPar);
257 cpl_ensure_code(npars > 0, CPL_ERROR_ILLEGAL_INPUT);
258 cpl_ensure_code(aSize > 0, CPL_ERROR_ILLEGAL_INPUT);
259 cpl_ensure_code(aFunction != NULL, CPL_ERROR_NULL_INPUT);
261 muse_cpl_optimize_eval_struct s = {
268 void (*debug_func) (
const char *, int,
double *, int,
double *,
void *,
269 int, int, int) = NULL;
270 lm_control_type control;
271 lm_initialize_control(&control);
274 control.maxcall = aCtrl->
maxcall;
276 if (aCtrl->
ftol > 0) {
277 control.ftol = aCtrl->
ftol;
279 if (aCtrl->
xtol > 0) {
280 control.xtol = aCtrl->
xtol;
282 if (aCtrl->
gtol > 0) {
283 control.gtol = aCtrl->
gtol;
285 if (aCtrl->
debug == CPL_TRUE) {
290 double timeinit = cpl_test_get_walltime();
291 double cpuinit = cpl_test_get_cputime();
293 cpl_array_get_size(aPar),
294 cpl_array_get_data_double(aPar),
299 double timefini = cpl_test_get_walltime();
300 double cpufini = cpl_test_get_cputime();
302 cpl_msg_debug(__func__,
"Minimizing finished after %i steps: %s",
303 control.nfev / ((
int)cpl_array_get_size(aPar)+1),
304 lm_infmsg[control.info]);
306 cpl_msg_debug(__func__,
"processing time %.3fs (%.3fs CPU, %d CPUs)",
307 timefini - timeinit, cpufini - cpuinit, omp_get_max_threads());
310 cpl_error_code res = error_code_tbl[control.info];
312 cpl_ensure_code(res == CPL_ERROR_NONE, res);
347 #define LM_MACHEP DBL_EPSILON 348 #define LM_DWARF DBL_MIN 350 #define LM_SQRT_DWARF sqrt(DBL_MIN) 351 #define LM_SQRT_GIANT sqrt(DBL_MAX) 353 #define LM_SQRT_DWARF 1.e-150 354 #define LM_SQRT_GIANT 1.e150 356 #define LM_USERTOL 30*LM_MACHEP 373 static void lm_initialize_control( lm_control_type * control ) {
374 control->maxcall = 100;
375 control->epsilon = LM_USERTOL;
376 control->stepbound = 100.;
377 control->ftol = LM_USERTOL;
378 control->xtol = LM_USERTOL;
379 control->gtol = LM_USERTOL;
382 static void lm_minimize(
int m_dat,
int n_par,
double *par,
383 void (*evaluate) (
double *par,
int m_dat,
double *fvec,
384 void *data,
int *info),
385 void (*printout) (
const char *func,
386 int n_par,
double *par,
int m_dat,
387 double *fvec,
void *data,
int iflag,
389 void *data, lm_control_type * control ) {
393 double *fvec, *diag, *fjac, *qtf, *wa1, *wa2, *wa3, *wa4;
399 if ( (fvec = (
double *) cpl_malloc(m *
sizeof(
double))) == NULL ||
400 (diag = (
double *) cpl_malloc(n *
sizeof(
double))) == NULL ||
401 (qtf = (
double *) cpl_malloc(n *
sizeof(
double))) == NULL ||
402 (fjac = (
double *) cpl_malloc(n*m*
sizeof(
double))) == NULL ||
403 (wa1 = (
double *) cpl_malloc(n *
sizeof(
double))) == NULL ||
404 (wa2 = (
double *) cpl_malloc(n *
sizeof(
double))) == NULL ||
405 (wa3 = (
double *) cpl_malloc(n *
sizeof(
double))) == NULL ||
406 (wa4 = (
double *) cpl_malloc(m *
sizeof(
double))) == NULL ||
407 (ipvt = (
int *) cpl_malloc(n *
sizeof(
int) )) == NULL ) {
418 lm_lmdif( m, n, par, fvec, control->ftol, control->xtol, control->gtol,
419 control->maxcall * (n + 1), control->epsilon, diag, 1,
420 control->stepbound, &(control->info),
421 &(control->nfev), fjac, ipvt, qtf, wa1, wa2, wa3, wa4,
422 evaluate, printout, data );
425 (*printout) (__func__, n, par, m, fvec, data, -1, 0, control->nfev);
426 control->fnorm = lm_enorm(m, fvec);
427 if ( control->info < 0 )
447 static const char *lm_infmsg[] = {
448 "fatal coding error (improper input parameters)",
449 "success (the relative error in the sum of squares is at most tol)",
450 "success (the relative error between x and the solution is at most tol)",
451 "success (both errors are at most tol)",
452 "trapped by degeneracy (fvec is orthogonal to the columns of the jacobian)" 453 "timeout (number of calls to fcn has reached maxcall*(n+1))",
454 "failure (ftol<tol: cannot reduce sum of squares any further)",
455 "failure (xtol<tol: cannot improve approximate solution any further)",
456 "failure (gtol<tol: cannot improve approximate solution any further)",
457 "exception (not enough memory)",
458 "exception (break requested within function evaluation)" 461 static const char *lm_shortmsg[] = {
484 static void lm_qrfac(
int m,
int n,
double *a,
int pivot,
int *ipvt,
485 double *rdiag,
double *acnorm,
double *wa);
486 static void lm_qrsolv(
int n,
double *r,
int ldr,
int *ipvt,
double *diag,
487 double *qtb,
double *x,
double *sdiag,
double *wa);
488 static void lm_lmpar(
int n,
double *r,
int ldr,
int *ipvt,
double *diag,
489 double *qtb,
double delta,
double *par,
double *x,
490 double *sdiag,
double *wa1,
double *wa2);
492 #define MIN(a,b) (((a)<=(b)) ? (a) : (b)) 493 #define MAX(a,b) (((a)>=(b)) ? (a) : (b)) 494 #define SQR(x) (x)*(x) 673 static void lm_lmdif(
int m,
int n,
double *x,
double *fvec,
double ftol,
674 double xtol,
double gtol,
int maxfev,
double epsfcn,
675 double *diag,
int mode,
double factor,
int *info,
676 int *nfev,
double *fjac,
int *ipvt,
double *qtf,
677 double *wa1,
double *wa2,
double *wa3,
double *wa4,
678 void (*evaluate) (
double *par,
int m_dat,
double *fvec,
679 void *data,
int *info),
680 void (*printout) (
const char *func,
int n_par,
double *par,
int m_dat,
681 double *fvec,
void *data,
int iflag,
686 double actred, delta, dirder, eps, fnorm, fnorm1, gnorm, par, pnorm,
687 prered, ratio, sum, temp, temp1, temp2, temp3, xnorm;
688 static double p1 = 0.1;
689 static double p0001 = 1.0e-4;
696 temp = MAX(epsfcn, LM_MACHEP);
701 if ((n <= 0) || (m < n) || (ftol < 0.)
702 || (xtol < 0.) || (gtol < 0.) || (maxfev <= 0) || (factor <= 0.)) {
707 for (j = 0; j < n; j++) {
708 if (diag[j] <= 0.0) {
721 (*evaluate) (x, m, fvec, data, info); ++(*nfev);
723 (*printout) (__func__, n, x, m, fvec, data, 0, 0, *nfev);
726 fnorm = lm_enorm(m, fvec);
732 cpl_msg_debug(__func__,
"iter=%d nfev=%d fnorm=%g\n",
739 #pragma omp parallel for default(none) \ 740 shared(n, info, m, x, eps, evaluate, data, printout, nfev, iter, \ 743 for (j = 0; j < n; j++) {
746 double *wa5 = (
double *) cpl_malloc(m *
sizeof(
double));
747 double *xd = (
double *) cpl_malloc(n *
sizeof(
double));
748 memcpy(xd, x, n *
sizeof(
double));
749 double step = eps * fabs(xd[j]);
754 (*evaluate) (xd, m, wa5, data, &infod); ++(*nfev);
756 (*printout) (__func__, n, xd, m, wa5, data, 1, iter, *nfev);
764 double diff = xd[j] - x[j];
765 double *fjacm = fjac + j * m;
766 for (i = 0; i < m; i++) {
767 fjacm[i] = (wa5[i] - fvec[i]) / diff;
779 for (i = 0; i < m; i++) {
780 for (j = 0; j < n; j++)
781 printf(
"%.5e ", fjac[j * m + i]);
789 lm_qrfac(m, n, fjac, 1, ipvt, wa1, wa2, wa3);
794 for (j = 0; j < n; j++) {
801 for (j = 0; j < n; j++)
802 wa3[j] = diag[j] * x[j];
803 xnorm = lm_enorm(n, wa3);
805 delta = factor * xnorm;
813 for (i = 0; i < m; i++)
816 for (j = 0; j < n; j++) {
817 temp3 = fjac[j * m + j];
820 for (i = j; i < m; i++)
821 sum += fjac[j * m + i] * wa4[i];
823 for (i = j; i < m; i++)
824 wa4[i] += fjac[j * m + i] * temp;
826 fjac[j * m + j] = wa1[j];
834 for (j = 0; j < n; j++) {
835 if (wa2[ipvt[j]] == 0)
839 for (i = 0; i <= j; i++)
840 sum += fjac[j * m + i] * qtf[i] / fnorm;
841 gnorm = MAX(gnorm, fabs(sum / wa2[ipvt[j]]));
853 for (j = 0; j < n; j++)
854 diag[j] = MAX(diag[j], wa2[j]);
860 printf(
"lmdif/ inner loop iter=%d nfev=%d\n", iter, *nfev);
865 lm_lmpar(n, fjac, m, ipvt, diag, qtf, delta, &par,
870 for (j = 0; j < n; j++) {
872 wa2[j] = x[j] + wa1[j];
873 wa3[j] = diag[j] * wa1[j];
875 pnorm = lm_enorm(n, wa3);
880 delta = MIN(delta, pnorm);
885 (*evaluate) (wa2, m, wa4, data, info); ++(*nfev);
887 (*printout) (__func__, n, x, m, wa4, data, 2, iter, *nfev);
891 fnorm1 = lm_enorm(m, wa4);
893 printf(
"lmdif/ pnorm %.10e fnorm1 %.10e fnorm %.10e" 894 " delta=%.10e par=%.10e\n",
895 pnorm, fnorm1, fnorm, delta, par);
900 if (p1 * fnorm1 < fnorm)
901 actred = 1 - SQR(fnorm1 / fnorm);
908 for (j = 0; j < n; j++) {
910 for (i = 0; i <= j; i++)
911 wa3[i] += fjac[j * m + i] * wa1[ipvt[j]];
913 temp1 = lm_enorm(n, wa3) / fnorm;
914 temp2 = sqrt(par) * pnorm / fnorm;
915 prered = SQR(temp1) + 2 * SQR(temp2);
916 dirder = -(SQR(temp1) + SQR(temp2));
920 ratio = prered != 0 ? actred / prered : 0;
922 printf(
"lmdif/ actred=%.10e prered=%.10e ratio=%.10e" 923 " sq(1)=%.10e sq(2)=%.10e dd=%.10e\n",
924 actred, prered, prered != 0 ? ratio : 0.,
925 SQR(temp1), SQR(temp2), dirder);
934 temp = 0.5 * dirder / (dirder + 0.55 * actred);
935 if (p1 * fnorm1 >= fnorm || temp < p1)
937 delta = temp * MIN(delta, pnorm / p1);
939 }
else if (par == 0. || ratio >= 0.75) {
946 if (ratio >= p0001) {
948 for (j = 0; j < n; j++) {
950 wa2[j] = diag[j] * x[j];
952 for (i = 0; i < m; i++)
954 xnorm = lm_enorm(n, wa2);
960 printf(
"ATTN: iteration considered unsuccessful\n");
967 if (fabs(actred) <= ftol && prered <= ftol && 0.5 * ratio <= 1)
969 if (delta <= xtol * xnorm)
978 if (fabs(actred) <= LM_MACHEP &&
979 prered <= LM_MACHEP && 0.5 * ratio <= 1)
981 if (delta <= LM_MACHEP * xnorm)
983 if (gnorm <= LM_MACHEP)
990 }
while (ratio < p0001);
1068 static void lm_lmpar(
int n,
double *r,
int ldr,
int *ipvt,
double *diag,
1069 double *qtb,
double delta,
double *par,
double *x,
1070 double *sdiag,
double *wa1,
double *wa2)
1072 int i, iter, j, nsing;
1073 double dxnorm, fp, fp_old, gnorm, parc, parl, paru;
1075 static double p1 = 0.1;
1085 for (j = 0; j < n; j++) {
1087 if (r[j * ldr + j] == 0 && nsing == n)
1093 printf(
"nsing %d ", nsing);
1095 for (j = nsing - 1; j >= 0; j--) {
1096 wa1[j] = wa1[j] / r[j + ldr * j];
1098 for (i = 0; i < j; i++)
1099 wa1[i] -= r[j * ldr + i] * temp;
1102 for (j = 0; j < n; j++)
1103 x[ipvt[j]] = wa1[j];
1109 for (j = 0; j < n; j++)
1110 wa2[j] = diag[j] * x[j];
1111 dxnorm = lm_enorm(n, wa2);
1112 fp = dxnorm - delta;
1113 if (fp <= p1 * delta) {
1115 printf(
"lmpar/ terminate (fp<p1*delta)\n");
1127 for (j = 0; j < n; j++)
1128 wa1[j] = diag[ipvt[j]] * wa2[ipvt[j]] / dxnorm;
1130 for (j = 0; j < n; j++) {
1132 for (i = 0; i < j; i++)
1133 sum += r[j * ldr + i] * wa1[i];
1134 wa1[j] = (wa1[j] - sum) / r[j + ldr * j];
1136 temp = lm_enorm(n, wa1);
1137 parl = fp / delta / temp / temp;
1142 for (j = 0; j < n; j++) {
1144 for (i = 0; i <= j; i++)
1145 sum += r[j * ldr + i] * qtb[i];
1146 wa1[j] = sum / diag[ipvt[j]];
1148 gnorm = lm_enorm(n, wa1);
1149 paru = gnorm / delta;
1151 paru = LM_DWARF / MIN(delta, p1);
1156 *par = MAX(*par, parl);
1157 *par = MIN(*par, paru);
1159 *par = gnorm / dxnorm;
1161 printf(
"lmpar/ parl %.4e par %.4e paru %.4e\n", parl, *par, paru);
1171 *par = MAX(LM_DWARF, 0.001 * paru);
1173 for (j = 0; j < n; j++)
1174 wa1[j] = temp * diag[j];
1175 lm_qrsolv(n, r, ldr, ipvt, wa1, qtb, x, sdiag, wa2);
1176 for (j = 0; j < n; j++)
1177 wa2[j] = diag[j] * x[j];
1178 dxnorm = lm_enorm(n, wa2);
1180 fp = dxnorm - delta;
1186 if (fabs(fp) <= p1 * delta
1187 || (parl == 0. && fp <= fp_old && fp_old < 0.)
1193 for (j = 0; j < n; j++)
1194 wa1[j] = diag[ipvt[j]] * wa2[ipvt[j]] / dxnorm;
1196 for (j = 0; j < n; j++) {
1197 wa1[j] = wa1[j] / sdiag[j];
1198 for (i = j + 1; i < n; i++)
1199 wa1[i] -= r[j * ldr + i] * wa1[j];
1201 temp = lm_enorm(n, wa1);
1202 parc = fp / delta / temp / temp;
1207 parl = MAX(parl, *par);
1209 paru = MIN(paru, *par);
1214 *par = MAX(parl, *par + parc);
1270 static void lm_qrfac(
int m,
int n,
double *a,
int pivot,
int *ipvt,
1271 double *rdiag,
double *acnorm,
double *wa)
1275 static double p05 = 0.05;
1279 for (j = 0; j < n; j++) {
1280 acnorm[j] = lm_enorm(m, &a[j * m]);
1281 rdiag[j] = acnorm[j];
1293 for (j = 0; j < minmn; j++) {
1294 double *am = a + j *m;
1300 for (k = j + 1; k < n; k++)
1301 if (rdiag[k] > rdiag[kmax])
1304 double *akmaxm = a + kmax * m;
1305 for (i = 0; i < m; i++) {
1306 double temp = am[i];
1310 rdiag[kmax] = rdiag[j];
1313 ipvt[j] = ipvt[kmax];
1320 ajnorm = lm_enorm(m - j, &am[j]);
1328 for (i = j; i < m; i++)
1335 for (k = j + 1; k < n; k++) {
1337 double *akm = a + k * m;
1339 for (i = j; i < m; i++)
1340 sum += am[i] * akm[i];
1342 double temp = sum / am[j];
1344 for (i = j; i < m; i++)
1345 akm[i] -= temp * am[i];
1347 if (pivot && rdiag[k] != 0.) {
1348 temp = akm[j] / rdiag[k];
1349 temp = MAX(0., 1 - temp * temp);
1350 rdiag[k] *= sqrt(temp);
1351 temp = rdiag[k] / wa[k];
1352 if (p05 * SQR(temp) <= LM_MACHEP) {
1353 rdiag[k] = lm_enorm(m - j - 1, &akm[j + 1]);
1425 lm_qrsolv(
int n,
double *r,
int ldr,
int *ipvt,
double *diag,
1426 double *qtb,
double *x,
double *sdiag,
double *wa)
1428 int i, kk, j, k, nsing;
1429 double qtbpj, sum, temp;
1430 double _sin, _cos, _tan, _cot;
1435 for (j = 0; j < n; j++) {
1436 for (i = j; i < n; i++)
1437 r[j * ldr + i] = r[i * ldr + j];
1438 x[j] = r[j * ldr + j];
1447 for (j = 0; j < n; j++) {
1452 if (diag[ipvt[j]] == 0.)
1454 for (k = j; k < n; k++)
1456 sdiag[j] = diag[ipvt[j]];
1463 for (k = j; k < n; k++) {
1471 if (fabs(r[kk]) < fabs(sdiag[k])) {
1472 _cot = r[kk] / sdiag[k];
1473 _sin = 1 / sqrt(1 + SQR(_cot));
1476 _tan = sdiag[k] / r[kk];
1477 _cos = 1 / sqrt(1 + SQR(_tan));
1484 r[kk] = _cos * r[kk] + _sin * sdiag[k];
1485 temp = _cos * wa[k] + _sin * qtbpj;
1486 qtbpj = -_sin * wa[k] + _cos * qtbpj;
1491 for (i = k + 1; i < n; i++) {
1492 temp = _cos * r[k * ldr + i] + _sin * sdiag[i];
1493 sdiag[i] = -_sin * r[k * ldr + i] + _cos * sdiag[i];
1494 r[k * ldr + i] = temp;
1502 sdiag[j] = r[j * ldr + j];
1503 r[j * ldr + j] = x[j];
1510 for (j = 0; j < n; j++) {
1511 if (sdiag[j] == 0. && nsing == n)
1517 for (j = nsing - 1; j >= 0; j--) {
1519 for (i = j + 1; i < nsing; i++)
1520 sum += r[j * ldr + i] * wa[i];
1521 wa[j] = (wa[j] - sum) / sdiag[j];
1526 for (j = 0; j < n; j++)
1548 lm_enorm(
int n,
double *x) {
1550 double agiant, s1, s2, s3, xabs, x1max, x3max, temp;
1557 agiant = LM_SQRT_GIANT / ((double) n);
1560 for (i = 0; i < n; i++) {
1562 if (xabs > LM_SQRT_DWARF && xabs < agiant) {
1568 if (xabs > LM_SQRT_DWARF) {
1571 temp = x1max / xabs;
1572 s1 = 1 + s1 * SQR(temp);
1575 temp = xabs / x1max;
1582 temp = x3max / xabs;
1583 s3 = 1 + s3 * SQR(temp);
1587 temp = xabs / x3max;
1596 return x1max * sqrt(s1 + (s2 / x1max) / x1max);
1599 return sqrt(s2 * (1 + (x3max / s2) * (x3max * s3)));
1601 return sqrt(x3max * ((s2 / x3max) + (x3max * s3)));
1604 return x3max * sqrt(s3);
cpl_error_code muse_cpl_optimize_lvmq(void *aData, cpl_array *aPar, int aSize, muse_cpl_evaluate_func *aFunction, muse_cpl_optimize_control_t *aCtrl)
Minimize a function with the Levenberg-Marquardt algorithm.
int debug
Flag to switch on debugging messages. Default value: CPL_FALSE.
int maxcall
Maximum number of iterations. Default value (when set to -1): 100.
double xtol
Relative error between last two approximations. Default value (when set to -1): 30 * DBL_EPSILON...
double ftol
Relative error desired in the sum of squares. Default value (when set to -1): 30 * DBL_EPSILON...
double gtol
Orthogonality desired between fvec and its derivs. Default value (when set to -1): 30 * DBL_EPSILON...
Optimization control parameters.
cpl_error_code( muse_cpl_evaluate_func)(void *aData, cpl_array *aPar, cpl_array *aRetval)
User provided function to evaluate in muse_cpl_optimize_lvmq().
static void muse_cpl_optimize_print(const char *func, int n_par, double *par, int m_dat, double *fvec, void *data, int iflag, int iter, int nfev)
Default printout method adopted to CPL.