grimme-lab / nlopt-f

Fortran bindings for the NLopt library
Apache License 2.0
28 stars 7 forks source link

Alternative to associating objective functions #13

Open awvwgk opened 2 years ago

awvwgk commented 2 years ago

Currently we implement objective functions by

  1. creating callback wrappers with bind(c) and pass those as procedure pointers to NLopt
  2. safe the actual user-defined objective function as data pointer (void*)
  3. resolve the user-defined objective function from our callback wrapper
  4. invoke the actual objective function

The drawback of this approach is that we only have a pointer to the objective function and the function object must at least exist until we call the optimize method. Therefore, using the following (intuitive) approach will fail:

call opt%set_min_objective(nlopt_func(myfunc))

call opt%add_inequality_constraint(nlopt_func(myconstraint, d1), 1.0e-8_wp)
call opt%add_inequality_constraint(nlopt_func(myconstraint, d2), 1.0e-8_wp)

call opt%set_xtol_rel(xtol)

x = [1.234_wp, 5.678_wp]
call opt%optimize(x, minf, stat)

Since the function objects contain only a procedure pointer and a class polymorphic pointer, which don't own any data, we could in principle copy those objects into our nlopt_opt instance on the Fortran side and ensure they have the same lifetime as the instance itself.

This would be done in one of the hooks to register an objective function like

https://github.com/grimme-lab/nlopt-f/blob/696b0e225fe393f1b77006912f09882c25599c48/src/nlopt_wrap.f90#L281-L290

However, any mechanism to own the function object in the nlopt_opt instance must ensure that the address of the function object is not changed by any operation on the wrapper, e.g. adding another objective function and resizing an array of function objects.

ivan-pi commented 2 years ago

I've been playing around with the idea of placing the function and constraints inside of a nonlinear optimization problem derived type. The interface becomes closer to the Julia NLopt.jl wrapper. Here's a sample:

program nlopt_alternative_demo

   use example_funcs
   use nlopt_alternative
   implicit none

   type(nlopt_prob(n = 2)) :: prob

   type :: constraint_data
      real(wp) :: a, b
   end type

   type(nlopt_func) :: cnstr(2)
   type(constraint_data), target :: d1, d2

   call create(prob, "LD_MMA")  ! initializes with default settings

   prob%lower_bounds(2) = 0.0_wp
   prob%xtol_rel = 1.e-4_wp

   prob%min_objective = nlopt_func(myfunc)

   d1 = constraint_data(2.0_wp, 0.0_wp)
   d2 = constraint_data(-1.0_wp, 1.0_wp)

   cnstr(1) = nlopt_func(myconstraint, d1)
   cnstr(2) = nlopt_func(myconstraint, d2)

   call add_inequality_constraints(opt, cnstr)

   x = [1.234_wp, 5.678_wp]

   call optimize(prob, x)

   print *, "Found minimum ", prob%optf, "at ", x
   print *, "Number of function evaluations ", prob%numevals

end program

A few other choices I've been playing with:

The actual settings get passed through the C-API directly within the call to optimize.

awvwgk commented 2 years ago

One problem is the assignment(=) which will require to use nlopt_copy to duplicate the handle and all the associated data, like the procedure/data pointer. The solution to actually defer all this to the optimize call is useful, however in case of an error, at which point do you become aware about passing an incorrect value or missing an upper bound?

I like the idea of having a nlopt_problem object/module, this could be built on top of the existing wrapper, especially if the creation of the handle and the setup is completely deferred to the optimize call. I think we can support something along this lines in nlopt-f.