Hello everyone,

I have an optimization (minimization) problem* in Julia to solve and I’m a looking for a package with which, most importantly, I can do it fast.

Specifically:

- the problem is in two variables (but it might reach even about 10 in the future, as I will develop it),
- the problem for these variables is on an a priori defined domain,
- the problem has no constraints,
- the objective function is (most likely) bounded,
- the objective function is nonlinear,
- the objective function has a subdomain (not precisely defined a priori) for which it is (only a posteriori) not defined. I was thinking of evading this problem in the minimization by setting on this subdomain an objective function to be very high compared to the values on the rest of the domain, e.g. =1000. For this reason, the function would not only be not differentiable, but also not continuous. If that helps, a priori I only know that in a rectangle domain of the two variables, the objective function is not defined on a bounded subdomain; if the rectangle was roughly split into two right-triangle-like subdomains, on the one of them the objective function is not defined; for certain in the one corner of the rectangle it is not, but the exact place of the hypotenuse (more generally: the curve rather than a straight line) is not known a priori. Also, the maximum of the objective function is most likely very close, just next to the line (curve) that is splitting the domain into these two subdomains (this “hypotenuse”), i.e. something similar to the tangent function close to \pi/2, but bounded. I’m however not sure if such a structure of the domain would also be the case if I extend the number of variables from two to, say, 5 or 10. Please, forgive my poor mathemetical description, I would do better if I knew how.

Now, I realize there are many packages and solution methods within them to solve the problem, and ideally I should experiment with all of them, but because just one evaluation of the objective function for a candidate pair of two variables takes about a minute, it would take ages to conduct such a meaningful experiment, in a situation in which I’m pretty time constrained. Moreover, I am not sure if the experiments performed with two variables, would be suggestive for the case when I would have to use more than two of them, say 5 or 10. For this reason, I’m looking for recommendations of a package as well as the method within it to be able to solve such a problem as fast as possible. I was thinking of something like, first, evaluating the objective function on (either a regular or a random) grid on the rectangle, and second, using Nelder-Mead with the initial guesses for the couple of such best gridpoints on the rectangle.

*Actually, the problem is the simulated method of moments estimation of a structural economic model. First, I try to estimate the two parameters, but later I will have to extend the set of estimated parameters.

Many thanks in advance!

Edit 1:

Also,

- My objective is a function of some statistics of a simulated economy, for which I don’t have any closed-form solutions.
- I have a 6-physical-core laptop, so I can take advantage of parallelism.

Edit 2:

My question is very general and I don’t think it hinges on a particular code, but below is the general structure of the objective function that I would minimize in some range on θ₁ and θ₂.

```
function ObjectiveFunction(θ₁,θ₂)
...
return (mₛ₁-mₕ₁)^2 + (mₛ₂-mₕ₂)^2 + (mₛ₃-mₕ₃)^2 # the number of terms here is at least two, also, they may come with some weights, ω
end
```

where the symbols come from the simulations that I perform for some parameters \Theta=(\theta_1,\theta_2), with the following structure:

- solve for the agents’ state variables’ policy functions,
- simulate the stochastic shocks for N agents over T periods,
- using 1. & 2., simulate the state variables for the N agents over T periods,
- compute some moments (statistics) based on these simulated state variables, then I minimize the quadratic form being the weighted combination of the differences of values of these simulated (s) moments with some empirical targets from the historical data (h) that I have, i.e. the objective function is of the form f(\Theta)=(m_{s,1}(\Theta)-m_{h,1})^2+(m_{s,2}(\Theta)-m_{h,2})^2.

When I say that I have no closed-form solution, I mean that I don’t have analytical/algebraic expression for m_{s,1}(\Theta) and m_{s,2}(\Theta), which I just take from the simulated values.