split every algorithm in *Init() + while(running) { *OneStep() }

This commit is contained in:
Thomas Capricelli 2009-08-25 23:43:33 +02:00
parent bbd44ef0ad
commit eac9293449
3 changed files with 1149 additions and 984 deletions

View File

@ -37,6 +37,17 @@ public:
Matrix< Scalar, Dynamic, 1 > &x, Matrix< Scalar, Dynamic, 1 > &x,
const Scalar tol = ei_sqrt(epsilon<Scalar>()) const Scalar tol = ei_sqrt(epsilon<Scalar>())
); );
Status solveInit(
Matrix< Scalar, Dynamic, 1 > &x,
const Parameters &parameters,
const int mode=1
);
Status solveOneStep(
Matrix< Scalar, Dynamic, 1 > &x,
const Parameters &parameters,
const int mode=1
);
Status solve( Status solve(
Matrix< Scalar, Dynamic, 1 > &x, Matrix< Scalar, Dynamic, 1 > &x,
const Parameters &parameters, const Parameters &parameters,
@ -47,6 +58,17 @@ public:
Matrix< Scalar, Dynamic, 1 > &x, Matrix< Scalar, Dynamic, 1 > &x,
const Scalar tol = ei_sqrt(epsilon<Scalar>()) const Scalar tol = ei_sqrt(epsilon<Scalar>())
); );
Status solveNumericalDiffInit(
Matrix< Scalar, Dynamic, 1 > &x,
const Parameters &parameters,
const int mode=1
);
Status solveNumericalDiffOneStep(
Matrix< Scalar, Dynamic, 1 > &x,
const Parameters &parameters,
const int mode=1
);
Status solveNumericalDiff( Status solveNumericalDiff(
Matrix< Scalar, Dynamic, 1 > &x, Matrix< Scalar, Dynamic, 1 > &x,
const Parameters &parameters, const Parameters &parameters,
@ -107,11 +129,9 @@ HybridNonLinearSolver<FunctorType,Scalar>::solve(
); );
} }
template<typename FunctorType, typename Scalar> template<typename FunctorType, typename Scalar>
typename HybridNonLinearSolver<FunctorType,Scalar>::Status typename HybridNonLinearSolver<FunctorType,Scalar>::Status
HybridNonLinearSolver<FunctorType,Scalar>::solve( HybridNonLinearSolver<FunctorType,Scalar>::solveInit(
Matrix< Scalar, Dynamic, 1 > &x, Matrix< Scalar, Dynamic, 1 > &x,
const Parameters &parameters, const Parameters &parameters,
const int mode const int mode
@ -124,7 +144,6 @@ HybridNonLinearSolver<FunctorType,Scalar>::solve(
qtf.resize(n); qtf.resize(n);
R.resize( (n*(n+1))/2); R.resize( (n*(n+1))/2);
fjac.resize(n, n); fjac.resize(n, n);
fvec.resize(n);
if (mode != 2) if (mode != 2)
diag.resize(n); diag.resize(n);
assert( (mode!=2 || diag.size()==n) || "When using mode==2, the caller must provide a valid 'diag'"); assert( (mode!=2 || diag.size()==n) || "When using mode==2, the caller must provide a valid 'diag'");
@ -158,9 +177,17 @@ HybridNonLinearSolver<FunctorType,Scalar>::solve(
nslow1 = 0; nslow1 = 0;
nslow2 = 0; nslow2 = 0;
/* beginning of the outer loop. */ return Running;
}
while (true) { template<typename FunctorType, typename Scalar>
typename HybridNonLinearSolver<FunctorType,Scalar>::Status
HybridNonLinearSolver<FunctorType,Scalar>::solveOneStep(
Matrix< Scalar, Dynamic, 1 > &x,
const Parameters &parameters,
const int mode
)
{
int i, j, l, iwa[1]; int i, j, l, iwa[1];
jeval = true; jeval = true;
@ -370,8 +397,23 @@ HybridNonLinearSolver<FunctorType,Scalar>::solve(
jeval = false; jeval = false;
} }
/* end of the outer loop. */ /* end of the outer loop. */
return Running;
} }
assert(false); // should never be reached
template<typename FunctorType, typename Scalar>
typename HybridNonLinearSolver<FunctorType,Scalar>::Status
HybridNonLinearSolver<FunctorType,Scalar>::solve(
Matrix< Scalar, Dynamic, 1 > &x,
const Parameters &parameters,
const int mode
)
{
Status status = solveInit(x, parameters, mode);
while (status==Running)
status = solveOneStep(x, parameters, mode);
return status;
} }
@ -403,10 +445,9 @@ HybridNonLinearSolver<FunctorType,Scalar>::solveNumericalDiff(
); );
} }
template<typename FunctorType, typename Scalar> template<typename FunctorType, typename Scalar>
typename HybridNonLinearSolver<FunctorType,Scalar>::Status typename HybridNonLinearSolver<FunctorType,Scalar>::Status
HybridNonLinearSolver<FunctorType,Scalar>::solveNumericalDiff( HybridNonLinearSolver<FunctorType,Scalar>::solveNumericalDiffInit(
Matrix< Scalar, Dynamic, 1 > &x, Matrix< Scalar, Dynamic, 1 > &x,
const Parameters &parameters, const Parameters &parameters,
const int mode const int mode
@ -428,6 +469,7 @@ HybridNonLinearSolver<FunctorType,Scalar>::solveNumericalDiff(
diag.resize(n); diag.resize(n);
assert( (mode!=2 || diag.size()==n) || "When using mode==2, the caller must provide a valid 'diag'"); assert( (mode!=2 || diag.size()==n) || "When using mode==2, the caller must provide a valid 'diag'");
/* Function Body */ /* Function Body */
nfev = 0; nfev = 0;
@ -457,11 +499,23 @@ HybridNonLinearSolver<FunctorType,Scalar>::solveNumericalDiff(
nslow1 = 0; nslow1 = 0;
nslow2 = 0; nslow2 = 0;
/* beginning of the outer loop. */ return Running;
}
while (true) { template<typename FunctorType, typename Scalar>
typename HybridNonLinearSolver<FunctorType,Scalar>::Status
HybridNonLinearSolver<FunctorType,Scalar>::solveNumericalDiffOneStep(
Matrix< Scalar, Dynamic, 1 > &x,
const Parameters &parameters,
const int mode
)
{
int i, j, l, iwa[1]; int i, j, l, iwa[1];
jeval = true; jeval = true;
int nsub = parameters.nb_of_subdiagonals;
int nsup = parameters.nb_of_superdiagonals;
if (nsub<0) nsub= n-1;
if (nsup<0) nsup= n-1;
/* calculate the jacobian matrix. */ /* calculate the jacobian matrix. */
@ -670,7 +724,21 @@ HybridNonLinearSolver<FunctorType,Scalar>::solveNumericalDiff(
jeval = false; jeval = false;
} }
/* end of the outer loop. */ /* end of the outer loop. */
}
assert(false); // should never be reached return Running;
}
template<typename FunctorType, typename Scalar>
typename HybridNonLinearSolver<FunctorType,Scalar>::Status
HybridNonLinearSolver<FunctorType,Scalar>::solveNumericalDiff(
Matrix< Scalar, Dynamic, 1 > &x,
const Parameters &parameters,
const int mode
)
{
Status status = solveNumericalDiffInit(x, parameters, mode);
while (status==Running)
status = solveNumericalDiffOneStep(x, parameters, mode);
return status;
} }

View File

@ -46,6 +46,16 @@ public:
const Parameters &parameters, const Parameters &parameters,
const int mode=1 const int mode=1
); );
Status minimizeInit(
Matrix< Scalar, Dynamic, 1 > &x,
const Parameters &parameters,
const int mode=1
);
Status minimizeOneStep(
Matrix< Scalar, Dynamic, 1 > &x,
const Parameters &parameters,
const int mode=1
);
Status minimizeNumericalDiff( Status minimizeNumericalDiff(
Matrix< Scalar, Dynamic, 1 > &x, Matrix< Scalar, Dynamic, 1 > &x,
@ -57,6 +67,16 @@ public:
const Parameters &parameters, const Parameters &parameters,
const int mode=1 const int mode=1
); );
Status minimizeNumericalDiffInit(
Matrix< Scalar, Dynamic, 1 > &x,
const Parameters &parameters,
const int mode=1
);
Status minimizeNumericalDiffOneStep(
Matrix< Scalar, Dynamic, 1 > &x,
const Parameters &parameters,
const int mode=1
);
Status minimizeOptimumStorage( Status minimizeOptimumStorage(
Matrix< Scalar, Dynamic, 1 > &x, Matrix< Scalar, Dynamic, 1 > &x,
@ -68,6 +88,16 @@ public:
const Parameters &parameters, const Parameters &parameters,
const int mode=1 const int mode=1
); );
Status minimizeOptimumStorageInit(
Matrix< Scalar, Dynamic, 1 > &x,
const Parameters &parameters,
const int mode=1
);
Status minimizeOptimumStorageOneStep(
Matrix< Scalar, Dynamic, 1 > &x,
const Parameters &parameters,
const int mode=1
);
Matrix< Scalar, Dynamic, 1 > fvec; Matrix< Scalar, Dynamic, 1 > fvec;
Matrix< Scalar, Dynamic, Dynamic > fjac; Matrix< Scalar, Dynamic, Dynamic > fjac;
@ -126,6 +156,20 @@ LevenbergMarquardt<FunctorType,Scalar>::minimize(
const Parameters &parameters, const Parameters &parameters,
const int mode const int mode
) )
{
Status status = minimizeInit(x, parameters, mode);
while (status==Running)
status = minimizeOneStep(x, parameters, mode);
return status;
}
template<typename FunctorType, typename Scalar>
typename LevenbergMarquardt<FunctorType,Scalar>::Status
LevenbergMarquardt<FunctorType,Scalar>::minimizeInit(
Matrix< Scalar, Dynamic, 1 > &x,
const Parameters &parameters,
const int mode
)
{ {
n = x.size(); n = x.size();
m = functor.nbOfFunctions(); m = functor.nbOfFunctions();
@ -167,9 +211,17 @@ LevenbergMarquardt<FunctorType,Scalar>::minimize(
par = 0.; par = 0.;
iter = 1; iter = 1;
/* beginning of the outer loop. */ return Running;
}
while (true) { template<typename FunctorType, typename Scalar>
typename LevenbergMarquardt<FunctorType,Scalar>::Status
LevenbergMarquardt<FunctorType,Scalar>::minimizeOneStep(
Matrix< Scalar, Dynamic, 1 > &x,
const Parameters &parameters,
const int mode
)
{
int i, j, l; int i, j, l;
/* calculate the jacobian matrix. */ /* calculate the jacobian matrix. */
@ -352,9 +404,9 @@ LevenbergMarquardt<FunctorType,Scalar>::minimize(
/* end of the inner loop. repeat if iteration unsuccessful. */ /* end of the inner loop. repeat if iteration unsuccessful. */
} while (ratio < Scalar(1e-4)); } while (ratio < Scalar(1e-4));
/* end of the outer loop. */ /* end of the outer loop. */
}
}
return Running;
}
template<typename FunctorType, typename Scalar> template<typename FunctorType, typename Scalar>
typename LevenbergMarquardt<FunctorType,Scalar>::Status typename LevenbergMarquardt<FunctorType,Scalar>::Status
@ -385,7 +437,7 @@ LevenbergMarquardt<FunctorType,Scalar>::minimizeNumericalDiff(
template<typename FunctorType, typename Scalar> template<typename FunctorType, typename Scalar>
typename LevenbergMarquardt<FunctorType,Scalar>::Status typename LevenbergMarquardt<FunctorType,Scalar>::Status
LevenbergMarquardt<FunctorType,Scalar>::minimizeNumericalDiff( LevenbergMarquardt<FunctorType,Scalar>::minimizeNumericalDiffInit(
Matrix< Scalar, Dynamic, 1 > &x, Matrix< Scalar, Dynamic, 1 > &x,
const Parameters &parameters, const Parameters &parameters,
const int mode const int mode
@ -429,9 +481,17 @@ LevenbergMarquardt<FunctorType,Scalar>::minimizeNumericalDiff(
par = 0.; par = 0.;
iter = 1; iter = 1;
/* beginning of the outer loop. */ return Running;
}
while (true) { template<typename FunctorType, typename Scalar>
typename LevenbergMarquardt<FunctorType,Scalar>::Status
LevenbergMarquardt<FunctorType,Scalar>::minimizeNumericalDiffOneStep(
Matrix< Scalar, Dynamic, 1 > &x,
const Parameters &parameters,
const int mode
)
{
int i, j, l; int i, j, l;
/* calculate the jacobian matrix. */ /* calculate the jacobian matrix. */
@ -615,8 +675,22 @@ LevenbergMarquardt<FunctorType,Scalar>::minimizeNumericalDiff(
/* end of the inner loop. repeat if iteration unsuccessful. */ /* end of the inner loop. repeat if iteration unsuccessful. */
} while (ratio < Scalar(1e-4)); } while (ratio < Scalar(1e-4));
/* end of the outer loop. */ /* end of the outer loop. */
return Running;
} }
assert(false); // should never be reached
template<typename FunctorType, typename Scalar>
typename LevenbergMarquardt<FunctorType,Scalar>::Status
LevenbergMarquardt<FunctorType,Scalar>::minimizeNumericalDiff(
Matrix< Scalar, Dynamic, 1 > &x,
const Parameters &parameters,
const int mode
)
{
Status status = minimizeNumericalDiffInit(x, parameters, mode);
while (status==Running)
status = minimizeNumericalDiffOneStep(x, parameters, mode);
return status;
} }
@ -651,7 +725,7 @@ LevenbergMarquardt<FunctorType,Scalar>::minimizeOptimumStorage(
template<typename FunctorType, typename Scalar> template<typename FunctorType, typename Scalar>
typename LevenbergMarquardt<FunctorType,Scalar>::Status typename LevenbergMarquardt<FunctorType,Scalar>::Status
LevenbergMarquardt<FunctorType,Scalar>::minimizeOptimumStorage( LevenbergMarquardt<FunctorType,Scalar>::minimizeOptimumStorageInit(
Matrix< Scalar, Dynamic, 1 > &x, Matrix< Scalar, Dynamic, 1 > &x,
const Parameters &parameters, const Parameters &parameters,
const int mode const int mode
@ -670,9 +744,6 @@ LevenbergMarquardt<FunctorType,Scalar>::minimizeOptimumStorage(
assert( (mode!=2 || diag.size()==n) || "When using mode==2, the caller must provide a valid 'diag'"); assert( (mode!=2 || diag.size()==n) || "When using mode==2, the caller must provide a valid 'diag'");
qtf.resize(n); qtf.resize(n);
/* Local variables */
bool sing;
/* Function Body */ /* Function Body */
nfev = 0; nfev = 0;
njev = 0; njev = 0;
@ -700,10 +771,20 @@ LevenbergMarquardt<FunctorType,Scalar>::minimizeOptimumStorage(
par = 0.; par = 0.;
iter = 1; iter = 1;
/* beginning of the outer loop. */ return Running;
}
while (true) {
template<typename FunctorType, typename Scalar>
typename LevenbergMarquardt<FunctorType,Scalar>::Status
LevenbergMarquardt<FunctorType,Scalar>::minimizeOptimumStorageOneStep(
Matrix< Scalar, Dynamic, 1 > &x,
const Parameters &parameters,
const int mode
)
{
int i, j, l; int i, j, l;
bool sing;
/* compute the qr factorization of the jacobian matrix */ /* compute the qr factorization of the jacobian matrix */
/* calculated one row at a time, while simultaneously */ /* calculated one row at a time, while simultaneously */
@ -902,8 +983,23 @@ LevenbergMarquardt<FunctorType,Scalar>::minimizeOptimumStorage(
/* end of the inner loop. repeat if iteration unsuccessful. */ /* end of the inner loop. repeat if iteration unsuccessful. */
} while (ratio < Scalar(1e-4)); } while (ratio < Scalar(1e-4));
/* end of the outer loop. */ /* end of the outer loop. */
}
assert(false); // should never be reached return Running;
}
template<typename FunctorType, typename Scalar>
typename LevenbergMarquardt<FunctorType,Scalar>::Status
LevenbergMarquardt<FunctorType,Scalar>::minimizeOptimumStorage(
Matrix< Scalar, Dynamic, 1 > &x,
const Parameters &parameters,
const int mode
)
{
Status status = minimizeOptimumStorageInit(x, parameters, mode);
while (status==Running)
status = minimizeOptimumStorageOneStep(x, parameters, mode);
return status;
} }

View File

@ -1836,6 +1836,7 @@ void test_NonLinear()
printf("x[1] : %.32g\n", x[1]); printf("x[1] : %.32g\n", x[1]);
printf("x[2] : %.32g\n", x[2]); printf("x[2] : %.32g\n", x[2]);
printf("x[3] : %.32g\n", x[3]); printf("x[3] : %.32g\n", x[3]);
printf("fvec.squaredNorm() : %.32g\n", fvec.squaredNorm()); printf("fvec.blueNorm() : %.32g\n", solver.fvec.blueNorm());
printf("fvec.blueNorm() : %.32g\n", lm.fvec.blueNorm());
*/ */