19 #ifndef EIGEN_LEVENBERGMARQUARDT_H 20 #define EIGEN_LEVENBERGMARQUARDT_H 24 namespace LevenbergMarquardtSpace {
28 ImproperInputParameters = 0,
29 RelativeReductionTooSmall = 1,
30 RelativeErrorTooSmall = 2,
31 RelativeErrorAndReductionTooSmall = 3,
33 TooManyFunctionEvaluation = 5,
41 template <
typename _Scalar,
int NX=Dynamic,
int NY=Dynamic>
44 typedef _Scalar Scalar;
46 InputsAtCompileTime = NX,
47 ValuesAtCompileTime = NY
49 typedef Matrix<Scalar,InputsAtCompileTime,1> InputType;
50 typedef Matrix<Scalar,ValuesAtCompileTime,1> ValueType;
51 typedef Matrix<Scalar,ValuesAtCompileTime,InputsAtCompileTime> JacobianType;
52 typedef ColPivHouseholderQR<JacobianType> QRSolver;
53 const int m_inputs, m_values;
55 DenseFunctor() : m_inputs(InputsAtCompileTime), m_values(ValuesAtCompileTime) {}
56 DenseFunctor(
int inputs,
int values) : m_inputs(inputs), m_values(values) {}
58 int inputs()
const {
return m_inputs; }
59 int values()
const {
return m_values; }
68 template <
typename _Scalar,
typename _Index>
71 typedef _Scalar Scalar;
73 typedef Matrix<Scalar,Dynamic,1> InputType;
74 typedef Matrix<Scalar,Dynamic,1> ValueType;
75 typedef SparseMatrix<Scalar, ColMajor, Index> JacobianType;
76 typedef SparseQR<JacobianType, COLAMDOrdering<int> > QRSolver;
78 InputsAtCompileTime = Dynamic,
79 ValuesAtCompileTime = Dynamic
82 SparseFunctor(
int inputs,
int values) : m_inputs(inputs), m_values(values) {}
84 int inputs()
const {
return m_inputs; }
85 int values()
const {
return m_values; }
87 const int m_inputs, m_values;
96 template <
typename QRSolver,
typename VectorType>
97 void lmpar2(
const QRSolver &qr,
const VectorType &diag,
const VectorType &qtb,
98 typename VectorType::Scalar m_delta,
typename VectorType::Scalar &par,
109 template<
typename _FunctorType>
113 typedef _FunctorType FunctorType;
114 typedef typename FunctorType::QRSolver QRSolver;
115 typedef typename FunctorType::JacobianType JacobianType;
116 typedef typename JacobianType::Scalar Scalar;
117 typedef typename JacobianType::RealScalar RealScalar;
118 typedef typename QRSolver::StorageIndex PermIndex;
119 typedef Matrix<Scalar,Dynamic,1> FVectorType;
120 typedef PermutationMatrix<Dynamic,Dynamic> PermutationType;
123 : m_functor(functor),m_nfev(0),m_njev(0),m_fnorm(0.0),m_gnorm(0),
124 m_isInitialized(
false),m_info(InvalidInput)
127 m_useExternalScaling=
false;
130 LevenbergMarquardtSpace::Status minimize(FVectorType &x);
131 LevenbergMarquardtSpace::Status minimizeInit(FVectorType &x);
132 LevenbergMarquardtSpace::Status minimizeOneStep(FVectorType &x);
133 LevenbergMarquardtSpace::Status lmder1(
135 const Scalar tol = std::sqrt(NumTraits<Scalar>::epsilon())
137 static LevenbergMarquardtSpace::Status lmdif1(
138 FunctorType &functor,
141 const Scalar tol = std::sqrt(NumTraits<Scalar>::epsilon())
151 m_ftol = sqrt(NumTraits<RealScalar>::epsilon());
152 m_xtol = sqrt(NumTraits<RealScalar>::epsilon());
158 void setXtol(RealScalar xtol) { m_xtol = xtol; }
161 void setFtol(RealScalar ftol) { m_ftol = ftol; }
164 void setGtol(RealScalar gtol) { m_gtol = gtol; }
167 void setFactor(RealScalar factor) { m_factor = factor; }
179 RealScalar
xtol()
const {
return m_xtol; }
182 RealScalar
ftol()
const {
return m_ftol; }
185 RealScalar
gtol()
const {
return m_gtol; }
188 RealScalar
factor()
const {
return m_factor; }
191 RealScalar
epsilon()
const {
return m_epsfcn; }
197 FVectorType&
diag() {
return m_diag; }
203 Index
nfev() {
return m_nfev; }
206 Index
njev() {
return m_njev; }
209 RealScalar
fnorm() {
return m_fnorm; }
212 RealScalar
gnorm() {
return m_gnorm; }
219 FVectorType&
fvec() {
return m_fvec; }
250 JacobianType m_rfactor;
251 FunctorType &m_functor;
252 FVectorType m_fvec, m_qtf, m_diag;
267 bool m_useExternalScaling;
268 PermutationType m_permutation;
269 FVectorType m_wa1, m_wa2, m_wa3, m_wa4;
271 bool m_isInitialized;
272 ComputationInfo m_info;
275 template<
typename FunctorType>
276 LevenbergMarquardtSpace::Status
279 LevenbergMarquardtSpace::Status status = minimizeInit(x);
280 if (status==LevenbergMarquardtSpace::ImproperInputParameters) {
281 m_isInitialized =
true;
286 status = minimizeOneStep(x);
287 }
while (status==LevenbergMarquardtSpace::Running);
288 m_isInitialized =
true;
292 template<
typename FunctorType>
293 LevenbergMarquardtSpace::Status
297 m = m_functor.values();
299 m_wa1.resize(n); m_wa2.resize(n); m_wa3.resize(n);
305 if (!m_useExternalScaling)
307 eigen_assert( (!m_useExternalScaling || m_diag.size()==n) &&
"When m_useExternalScaling is set, the caller must provide a valid 'm_diag'");
315 if (n <= 0 || m < n || m_ftol < 0. || m_xtol < 0. || m_gtol < 0. || m_maxfev <= 0 || m_factor <= 0.){
316 m_info = InvalidInput;
317 return LevenbergMarquardtSpace::ImproperInputParameters;
320 if (m_useExternalScaling)
321 for (Index j = 0; j < n; ++j)
324 m_info = InvalidInput;
325 return LevenbergMarquardtSpace::ImproperInputParameters;
331 if ( m_functor(x, m_fvec) < 0)
332 return LevenbergMarquardtSpace::UserAsked;
333 m_fnorm = m_fvec.stableNorm();
339 return LevenbergMarquardtSpace::NotStarted;
342 template<
typename FunctorType>
343 LevenbergMarquardtSpace::Status
350 m = m_functor.values();
353 if (n <= 0 || m < n || tol < 0.)
354 return LevenbergMarquardtSpace::ImproperInputParameters;
359 m_maxfev = 100*(n+1);
365 template<
typename FunctorType>
366 LevenbergMarquardtSpace::Status
368 FunctorType &functor,
375 Index m = functor.values();
378 if (n <= 0 || m < n || tol < 0.)
379 return LevenbergMarquardtSpace::ImproperInputParameters;
388 LevenbergMarquardtSpace::Status info = LevenbergMarquardtSpace::Status(lm.minimize(x));
396 #endif // EIGEN_LEVENBERGMARQUARDT_H RealScalar fnorm()
Definition: LevenbergMarquardt/LevenbergMarquardt.h:209
Performs non linear optimization over a non-linear function, using a variant of the Levenberg Marquar...
Definition: LevenbergMarquardt/LevenbergMarquardt.h:110
RealScalar ftol() const
Definition: LevenbergMarquardt/LevenbergMarquardt.h:182
JacobianType & jacobian()
Definition: LevenbergMarquardt/LevenbergMarquardt.h:223
Namespace containing all symbols from the Eigen library.
Definition: AdolcForward:45
FVectorType & fvec()
Definition: LevenbergMarquardt/LevenbergMarquardt.h:219
RealScalar gtol() const
Definition: LevenbergMarquardt/LevenbergMarquardt.h:185
PermutationType permutation()
Definition: LevenbergMarquardt/LevenbergMarquardt.h:232
void setFtol(RealScalar ftol)
Definition: LevenbergMarquardt/LevenbergMarquardt.h:161
Index maxfev() const
Definition: LevenbergMarquardt/LevenbergMarquardt.h:194
ComputationInfo info() const
Reports whether the minimization was successful.
Definition: LevenbergMarquardt/LevenbergMarquardt.h:243
void setEpsilon(RealScalar epsfcn)
Definition: LevenbergMarquardt/LevenbergMarquardt.h:170
RealScalar xtol() const
Definition: LevenbergMarquardt/LevenbergMarquardt.h:179
void setMaxfev(Index maxfev)
Definition: LevenbergMarquardt/LevenbergMarquardt.h:173
Index nfev()
Definition: LevenbergMarquardt/LevenbergMarquardt.h:203
Index njev()
Definition: LevenbergMarquardt/LevenbergMarquardt.h:206
RealScalar lm_param(void)
Definition: LevenbergMarquardt/LevenbergMarquardt.h:215
void resetParameters()
Definition: LevenbergMarquardt/LevenbergMarquardt.h:145
void setFactor(RealScalar factor)
Definition: LevenbergMarquardt/LevenbergMarquardt.h:167
RealScalar factor() const
Definition: LevenbergMarquardt/LevenbergMarquardt.h:188
FVectorType & diag()
Definition: LevenbergMarquardt/LevenbergMarquardt.h:197
void setXtol(RealScalar xtol)
Definition: LevenbergMarquardt/LevenbergMarquardt.h:158
JacobianType & matrixR()
Definition: LevenbergMarquardt/LevenbergMarquardt.h:228
Index iterations()
Definition: LevenbergMarquardt/LevenbergMarquardt.h:200
RealScalar epsilon() const
Definition: LevenbergMarquardt/LevenbergMarquardt.h:191
Definition: NumericalDiff.h:36
void setGtol(RealScalar gtol)
Definition: LevenbergMarquardt/LevenbergMarquardt.h:164
void setExternalScaling(bool value)
Definition: LevenbergMarquardt/LevenbergMarquardt.h:176
RealScalar gnorm()
Definition: LevenbergMarquardt/LevenbergMarquardt.h:212