merge both c methods lmstr/lmstr1 into one class

LevenbergMarquardtOptimumStorage with two methods.
This commit is contained in:
Thomas Capricelli 2009-08-25 14:18:38 +02:00
parent 3f1b81e129
commit 201f58e528
5 changed files with 90 additions and 57 deletions

View File

@ -48,13 +48,13 @@ namespace Eigen {
#include "src/NonLinear/dogleg.h" #include "src/NonLinear/dogleg.h"
#include "src/NonLinear/covar.h" #include "src/NonLinear/covar.h"
#include "src/NonLinear/chkder.h"
#include "src/NonLinear/lmder.h" #include "src/NonLinear/lmder.h"
#include "src/NonLinear/hybrd.h" #include "src/NonLinear/hybrd.h"
#include "src/NonLinear/lmstr.h" #include "src/NonLinear/lmstr.h"
#include "src/NonLinear/lmdif.h" #include "src/NonLinear/lmdif.h"
#include "src/NonLinear/hybrj.h" #include "src/NonLinear/hybrj.h"
#include "src/NonLinear/lmstr1.h"
#include "src/NonLinear/chkder.h"
//@} //@}

View File

@ -50,7 +50,7 @@ int LevenbergMarquardt<FunctorType,Scalar>::minimize(
/* check the input parameters for errors. */ /* check the input parameters for errors. */
if (n <= 0 || m < n || tol < 0.) { if (n <= 0 || m < n || tol < 0.) {
printf("ei_lmder1 bad args : m,n,tol,..."); printf("LevenbergMarquardt::minimize() bad args : m,n,tol,...");
return 0; return 0;
} }

View File

@ -1,7 +1,73 @@
template<typename FunctorType, typename Scalar> template<typename FunctorType, typename Scalar>
int ei_lmstr( class LevenbergMarquardtOptimumStorage
const FunctorType &Functor, {
public:
LevenbergMarquardtOptimumStorage(const FunctorType &_functor)
: functor(_functor) {}
int minimize(
Matrix< Scalar, Dynamic, 1 > &x,
Matrix< Scalar, Dynamic, 1 > &fvec,
const Scalar tol = ei_sqrt(epsilon<Scalar>())
);
int minimize(
Matrix< Scalar, Dynamic, 1 > &x,
Matrix< Scalar, Dynamic, 1 > &fvec,
int &nfev,
int &njev,
Matrix< Scalar, Dynamic, Dynamic > &fjac,
VectorXi &ipvt,
Matrix< Scalar, Dynamic, 1 > &qtf,
Matrix< Scalar, Dynamic, 1 > &diag,
int mode=1,
Scalar factor = 100.,
int maxfev = 400,
Scalar ftol = ei_sqrt(epsilon<Scalar>()),
Scalar xtol = ei_sqrt(epsilon<Scalar>()),
Scalar gtol = Scalar(0.),
int nprint=0
);
private:
const FunctorType &functor;
};
template<typename FunctorType, typename Scalar>
int LevenbergMarquardtOptimumStorage<FunctorType,Scalar>::minimize(
Matrix< Scalar, Dynamic, 1 > &x,
Matrix< Scalar, Dynamic, 1 > &fvec,
Scalar tol
)
{
const int n = x.size(), m=fvec.size();
int info, nfev=0, njev=0;
Matrix< Scalar, Dynamic, Dynamic > fjac(m, n);
Matrix< Scalar, Dynamic, 1> diag, qtf;
VectorXi ipvt(n);
/* check the input parameters for errors. */
if (n <= 0 || m < n || tol < 0.) {
printf("LevenbergMarquardtOptimumStorage::minimize() bad args : m,n,tol,...");
return 0;
}
info = minimize(
x, fvec,
nfev, njev,
fjac, ipvt, qtf, diag,
1,
100.,
(n+1)*100,
tol, tol, Scalar(0.)
);
return (info==8)?4:info;
}
template<typename FunctorType, typename Scalar>
int LevenbergMarquardtOptimumStorage<FunctorType,Scalar>::minimize(
Matrix< Scalar, Dynamic, 1 > &x, Matrix< Scalar, Dynamic, 1 > &x,
Matrix< Scalar, Dynamic, 1 > &fvec, Matrix< Scalar, Dynamic, 1 > &fvec,
int &nfev, int &nfev,
@ -10,13 +76,13 @@ int ei_lmstr(
VectorXi &ipvt, VectorXi &ipvt,
Matrix< Scalar, Dynamic, 1 > &qtf, Matrix< Scalar, Dynamic, 1 > &qtf,
Matrix< Scalar, Dynamic, 1 > &diag, Matrix< Scalar, Dynamic, 1 > &diag,
int mode=1, int mode,
Scalar factor = 100., Scalar factor,
int maxfev = 400, int maxfev,
Scalar ftol = ei_sqrt(epsilon<Scalar>()), Scalar ftol,
Scalar xtol = ei_sqrt(epsilon<Scalar>()), Scalar xtol,
Scalar gtol = Scalar(0.), Scalar gtol,
int nprint=0 int nprint
) )
{ {
const int m = fvec.size(), n = x.size(); const int m = fvec.size(), n = x.size();
@ -56,7 +122,7 @@ int ei_lmstr(
/* evaluate the function at the starting point */ /* evaluate the function at the starting point */
/* and calculate its norm. */ /* and calculate its norm. */
iflag = Functor.f(x, fvec); iflag = functor.f(x, fvec);
nfev = 1; nfev = 1;
if (iflag < 0) if (iflag < 0)
goto algo_end; goto algo_end;
@ -71,12 +137,12 @@ int ei_lmstr(
while (true) { while (true) {
/* if requested, call Functor.f to enable printing of iterates. */ /* if requested, call functor.f to enable printing of iterates. */
if (nprint > 0) { if (nprint > 0) {
iflag = 0; iflag = 0;
if ((iter - 1) % nprint == 0) if ((iter - 1) % nprint == 0)
iflag = Functor.debug(x, fvec, wa3); iflag = functor.debug(x, fvec, wa3);
if (iflag < 0) if (iflag < 0)
break; break;
} }
@ -90,7 +156,7 @@ int ei_lmstr(
fjac.fill(0.); fjac.fill(0.);
iflag = 2; iflag = 2;
for (i = 0; i < m; ++i) { for (i = 0; i < m; ++i) {
if (Functor.df(x, wa3, iflag) < 0) if (functor.df(x, wa3, iflag) < 0)
break; break;
temp = fvec[i]; temp = fvec[i];
ei_rwupdt<Scalar>(n, fjac.data(), fjac.rows(), wa3.data(), qtf.data(), &temp, wa1.data(), wa2.data()); ei_rwupdt<Scalar>(n, fjac.data(), fjac.rows(), wa3.data(), qtf.data(), &temp, wa1.data(), wa2.data());
@ -195,7 +261,7 @@ int ei_lmstr(
/* evaluate the function at x + p and calculate its norm. */ /* evaluate the function at x + p and calculate its norm. */
iflag = Functor.f(wa2, wa4); iflag = functor.f(wa2, wa4);
++nfev; ++nfev;
if (iflag < 0) if (iflag < 0)
goto algo_end; goto algo_end;
@ -292,7 +358,7 @@ algo_end:
if (iflag < 0) if (iflag < 0)
info = iflag; info = iflag;
if (nprint > 0) if (nprint > 0)
iflag = Functor.debug(x, fvec, wa3); iflag = functor.debug(x, fvec, wa3);
return info; return info;
} }

View File

@ -1,35 +0,0 @@
template<typename FunctorType, typename Scalar>
int ei_lmstr1(
const FunctorType &Functor,
Matrix< Scalar, Dynamic, 1 > &x,
Matrix< Scalar, Dynamic, 1 > &fvec,
VectorXi &ipvt,
Scalar tol = ei_sqrt(epsilon<Scalar>())
)
{
const int n = x.size(), m=fvec.size();
int info, nfev=0, njev=0;
Matrix< Scalar, Dynamic, Dynamic > fjac(m, n);
Matrix< Scalar, Dynamic, 1> diag, qtf;
/* check the input parameters for errors. */
if (n <= 0 || m < n || tol < 0.) {
printf("ei_lmstr1 bad args : m,n,tol,...");
return 0;
}
ipvt.resize(n);
info = ei_lmstr(
Functor,
x, fvec,
nfev, njev,
fjac, ipvt, qtf, diag,
1,
100.,
(n+1)*100,
tol, tol, Scalar(0.)
);
return (info==8)?4:info;
}

View File

@ -453,13 +453,14 @@ void testLmstr1()
int m=15, n=3, info; int m=15, n=3, info;
VectorXd x(n), fvec(m); VectorXd x(n), fvec(m);
VectorXi ipvt;
/* the following starting values provide a rough fit. */ /* the following starting values provide a rough fit. */
x.setConstant(n, 1.); x.setConstant(n, 1.);
// do the computation // do the computation
info = ei_lmstr1(lmstr_functor(), x, fvec, ipvt); lmstr_functor functor;
LevenbergMarquardtOptimumStorage<lmstr_functor,double> lm(functor);
info = lm.minimize(x, fvec);
// check return value // check return value
VERIFY( 1 == info); VERIFY( 1 == info);
@ -486,8 +487,9 @@ void testLmstr()
x.setConstant(n, 1.); x.setConstant(n, 1.);
// do the computation // do the computation
info = ei_lmstr(lmstr_functor(), x, fvec, nfev, njev, fjac, ipvt, qtf, diag); lmstr_functor functor;
VectorXd wa(n); LevenbergMarquardtOptimumStorage<lmstr_functor,double> lm(functor);
info = lm.minimize(x, fvec, nfev, njev, fjac, ipvt, qtf, diag);
// check return values // check return values
VERIFY( 1 == info); VERIFY( 1 == info);