Fortran 77 implementation of NEWUOA/BOBYQA is buggy and challenging to maintain. Switch to the PRIMA implementation?

[zaikunzhang]

Dear maintainers of Elmer FEM,

This is Dr. Zaikun Zhang from the Hong Kong Polytechnic University. Together with Professor N.I.M. Gould, I am responsible for maintaining the derivative-free optimization solvers of the late Professor M.J.D. Powell. I have developed PRIMA, which is a package for solving general nonlinear optimization problems without using derivatives. It provides the reference implementation for Powell's derivative-free optimization methods, i.e., COBYLA, UOBYQA, NEWUOA, BOBYQA, and LINCOA.

Thank you for making NEWUOA and BOBYQA available in Elmer FEM. I note that the current version is based on the original Fortran 77 implementation (refactored to Fortran 90 without essential changes). However, The Fortran 77 implementation is buggy and not maintained anymore.

Although the Fortran 77 code is truly a masterpiece, it contains many bugs, most of which are due to the language itself. For example, see Section 4.4 of our recent paper and the GitHub issues / requests listed below.

• The Fortran 77 solvers may get stuck in infinite loops.

  • Gradients example doesn't make progress on the GPU runners

  • optimize: COBYLA hangs / infinite loop #8998

  • BUG: Scipy.optimize / COBYLA hangs on some CPUs #15527

  • COBYLA freezes (though maxeval and maxtime are given) #370

  • COBYLA hangs #118

  • NEWUOA_BOUND stuck in infinite loop inside MMA #117

  • Cobyla freezes in 0T1.16rc1 #1651

  • Optimization freezes #25

  • BOBYQA gets stuck in infinite loop. #7

  • Algorithm turns into infinite loop and never finishes #3

• The Fortran 77 solvers may crash with segmentation faults due to uninitialized variables that are used as indices.

  • Fix all uninitialized variable warnings #134

  • BOBYQA uninitialised variables in rare cases #133

  • Use of uninitialized variable in BOBYQA altmov #36

• Fortran 77 COBYLA may not return the best point that is evaluated; sometimes, the returned point can have a large constraint violation even though the starting point is feasible.

  • nlopt COBYLA optimizer gives unexpected output #182

  • Last Result Returned Not Optimized Result #110

  • COBYLA returns last evaluated function which might not be minimum #57

  • Successful termination when constraints violated #1

To avoid the problems originating from the Fortran 77 code, I suggest you use the PRIMA implementation of Powell's solvers. PRIMA provides the modern implementation of the solvers in Fortran 2008. It fixes bugs in the original Fortran 77 code. In addition, it introduces improvements that boost the performance in terms of the number of function evaluations, which is the standard measure of computational costs in derivative-free optimization.

See the GitHub repo of PRIMA for more information. I will be glad to provide assistance if help is needed.

Thanks.

Best regards, Zaikun Zhang, Ph.D. and Assistant Professor Department of Applied Mathematics The Hong Kong Polytechnic University

[raback]

Thank you for making us aware of these recent developments! It indeed seems that it would make sense to replace the code with the improved one. Currently this is not maybe the highest priority since there are no issues with this currently.

Many people do not know that these routines exist in Elmer. They are intended for optimization problems (minimization of some cost function) within a single Elmer run. This may be optimal in that the simulation does not need to be restarted all the time. Even more useful is the fact that the optimization may be used inside the simulation in ways that go beyond standard optimization.

NEWUOA/BOBYQA were selected at the time partly because they were provided as F90 routines making the interplay with Elmer code straight-forward. The only test cases in public repo are: OptimizeBobyqaFourHeaters OptimizeNewuoaFourHeaters.