Multiple initial values in optimize()

[Christian Labat]

I am attempting to minimize an objective function with multiple local minima, and would like to ask Mata to try multiple initial values so that I can identify which local minimum is the global minimum. The code below shows that obj() has local minimums around x = -6 and x = 4, and which minimum is realized depends on the initial value. How might I do this for say 30 random values in a given range?

void obj(real scalar todo, real scalar x, fval, grad, hess)
    {
        fval = x * sin(x+1)
    }
    
    // minimize
    S  = optimize_init()
    optimize_init_evaluator(S, &obj())
    optimize_init_which(S, "min")
    optimize_init_params(S, 4)
    xstar1 = optimize(S)
    xstar1
    fval1 = optimize_result_value(S)
    fval1
    
    S  = optimize_init()
    optimize_init_evaluator(S, &obj())
    optimize_init_which(S, "min")
    optimize_init_params(S, -6)
    xstar2 = optimize(S)
    xstar2
    fval2 = optimize_result_value(S)
    fval2
end

This is the output.

. mata: mata clear

. mata
------------------------------------------------- mata (type end to exit) ----------------------------------
: void obj(real scalar todo, real scalar x, fval, grad, hess)
>         {
>                 fval = x * sin(x+1)
>         }
note: argument todo unused
note: argument grad unused
note: argument hess unused

:         
:         // minimize
:         S  = optimize_init()

:         optimize_init_evaluator(S, &obj())

:         optimize_init_which(S, "min")

:         optimize_init_params(S, 4)

:         xstar1 = optimize(S)
Iteration 0:   f(p) = -3.8356971  
Iteration 1:   f(p) = -3.8392227  
Iteration 2:   f(p) = -3.8392229  

:         xstar1
  3.9597575

:         fval1 = optimize_result_value(S)

:         fval1
  -3.839222897

:         
:         S  = optimize_init()

:         optimize_init_evaluator(S, &obj())

:         optimize_init_which(S, "min")

:         optimize_init_params(S, -6)

:         xstar2 = optimize(S)
Iteration 0:   f(p) = -5.7535456  
Iteration 1:   f(p) = -5.7975924  
Iteration 2:   f(p) = -5.7976005  
Iteration 3:   f(p) = -5.7976005  

:         xstar2
  -5.880822144

:         fval2 = optimize_result_value(S)

:         fval2
  -5.797600501

: end