rmathematical-optimizationlevenberg-marquardt

Optimization in R: Levenberg-Marquardt using nls.lm in minpack.lm: resetting `maxiter' to 1024


I am trying to learn how to work with nls.lm in the R library minpack.lm by using the Rosenbrock function to see if the algorithm converges to the global minimum at f(x,y) = (1,1). I do so both with and without the analytic Jacobian. In both instances, I get a warning telling me that the algorithm has decided to revert the maximum number of iterations specified in the call to nls.lm to 1024:

Warning messages:
1: In nls.lm(par = initpar, fn = objective_rosenbrock, jac = gradient_rosenbrock,  :
  resetting `maxiter' to 1024!
2: In nls.lm(par = initpar, fn = objective_rosenbrock, jac = gradient_rosenbrock,  :
  lmder: info = -1. Number of iterations has reached `maxiter' == 1024.

The algorithm never quite reaches (1,1) as a result given my initial guess of (-1.2, 1.0). I found the source code for the library on GitHub and the following lines of code are pertinent here:

https://github.com/cran/minpack.lm/blob/master/src/nls_lm.c

OS->maxiter = INTEGER_VALUE(getListElement(control, "maxiter"));
if(OS->maxiter > 1024) {
  OS->maxiter = 1024;
  warning("resetting `maxiter' to 1024!");
}

Is there any logic to why the maximum number of iterations is capped to 1024? Something with bits and 2^10? I would like to use the library for a different application, but this cap on iterations might prevent that. Any insight would be appreciated.


Solution

  • Git blame says that this code limiting the max iterations was introduced in version 1.1-0, in 2008. The NEWS file for the package only goes back as far as version 1.1-6. I can't find the code in any public repo other than the one you point to (which is only a CRAN mirror; it doesn't contain any comments/commit messages/etc. from developers that might give us clues.)

    Other than contacting the maintainer I think it's going to be hard to figure out what the rationale is for this limit.

    I do have some guesses though.

    The only places that maxiter is actually used in the code are here and here - in R code, not Fortran or C code, so it seems extremely unlikely that we are dealing with something like a 10-bit unsigned integer type (which seems an unlikely choice in any case). I think the limitation is there because we also have a buffer defined for holding trace information here:

      double rsstrace[1024];
    

    which, as you can see, is hard-coded to a length of 1024. Presumably bad things would happen if we tried to stuff 1025 iterations'-worth of tracing information into this array ...

    My suggestions:

    $ find . -type f -exec grep -Hn 1024 {} \;
    ./src/nls_lm.c:141:    if(OS->maxiter > 1024) {
    ./src/nls_lm.c:142:      OS->maxiter = 1024;
    ./src/nls_lm.c:143:      warning("resetting `maxiter' to 1024!");
    ./src/minpack_lm.h:20:  double rsstrace[1024];