I am trying to solve a simple function with two variables (a profit function). I tried two ways: using Optim and hence function, gradient and hessian, as well as NLsolve, and hence only gradient and hessian. However, both either give me “DomainError: Exponentiation yielding a complex result requires a complex argument.” but mostly “MethodError: no method matching”. My first idea was that there is a type conversion within the iteration that let’s the solvers crash. Do you have any ideas???

Could you post a MWE with your error? The above gives me μ not defined, so you might have defined a global somewhere earlier in your code?

Your errors seem to indicate two problems:

(1) DomainError: this likely comes from exponentiating negative numbers with exponents between -1 and 1 (e.g. in x[2]^β2); you could get around this by constraining your search to positive x's (if you’re looking at a profit function I’d assume these are factor inputs and therefore positive in any case?)

(2) MethodError: this is harder to diagnose without seeing all your actual function calls, but given that your only function seems to be f! this likely comes from writing the function definition in a way that is incompatible with the way the function is called in the solver. Optim.optimize for example would expect a function which takes one input array (in your case x), and returns the function value to minimise; depending on the method of optimize you’re calling you might also need to include the gradient (as a separate function!). Your function seems to return the [2,2] element of the Jacobian matrix, which is probably not what you’re actually trying to minimise?

first of all, one mistake I figured out is that the gradient and the hessian are not defined for zero values since I divide by the inputs and have negative exponents in gradient and hessian. This indeed would be solved by restricting the input values >0.
Let me give a different way I coded it using optim (so that I premultiplied function, gradient and hessian with -1 since I want to maximize)

Thanks, that makes more sense. What you probably want to do then is use box constrained optimization to set lower bounds, described in the Optim.jl manual here. The function call would then be:

The problem with this is that for some reasons when I tried it just now, Fminbox(inner_optimizer)) throws an error, which seems unrelated to your code (the error also shows up when simply running the example from the Optim.jl docs - I’ll file an issue)

One other option that I’ve used in my own thesis is the NLopt package which wraps the widely used NLopt library. For this you would do something like:

(Note that I haven’t thought at all about your problem so the choice of bounds and/or optimisation algorithm here is likely suboptimal but it should hopefully get you going!)

Hope that helps!

P.S. send my regards to Sevi - he was one of the examiners in my viva (small world)!