I’m new to JuMP and I’d like to minimize a function that has an integral in it. Here’s an example f(y) = \int_0^y y dy - y. Of course, this is just f(y)=(1/2)y^2-y which is minimized at y=1, but this is just a minimal example as a step towards something more complicated.
I don’t know how to pass this objective f into the JuMP stuff. I believe I need to register f as a “user-defined nonlinear function”, but I don’t know how to do that. According to the docs, “in addition to this list of functions, it is possible to register custom user-defined nonlinear functions.” But I don’t see how to do that and it doesn’t appear to say anything else.
Here’s the code I tried so far. The error message was: Unexpected object #89 (of type var"#89#90" in nonlinear expression.
Thanks!
using QuadGK
f(y) = quadgk(x -> x, 0, y, rtol=1e-8)[1] - y
using JuMP
function example()
model = Model(Ipopt.Optimizer)
set_silent(model)
@NLobjective(model, Min, f)
optimize!(model)
return
end
example()
# -> Unexpected object #89 (of type var"#89#90" in nonlinear expression.
You need to define a variable @variable(model, y), you need to call the function f(y), and you’ll need to register the user-defined function register(model, :f, 1, f; autodiff = true).
I haven’t run the code so there might be a typo etc, but this should point you in the right direction:
function example()
model = Model(Ipopt.Optimizer)
register(model, :f, 1, f; autodiff = true)
set_silent(model)
@variable(model, y)
@NLobjective(model, Min, f(y))
optimize!(model)
return
end
Thanks. From the docs here, I adjusted your code to replace the y’s in your second line with foo and f. (using f instead of foo works too). This worked when I use the analytically defined function such as f(y) = 0.5*y^2 - y and gives a solution of y=1. But when I used the equivalent function (up to numerical approx) f(y) = quadgk(x -> x, 0, y, rtol=1e-8)[1] - y I get the following error saying something “repeats 79984 times”:
Here’s my code. Any thoughts on why it’s not working when I use the numerical integral version of the objective?
f(y) = quadgk(x -> x, 0, y, rtol=1e-8)[1] - y
g(y) = 0.5*y^2 - y # g = f up to integration approx error
y = 0:0.1:10
function example()
model = Model(Ipopt.Optimizer)
@variable(model, y)
register(model, :foo, 1, f; autodiff = true)
set_silent(model)
@NLobjective(model, Min, foo(y))
optimize!(model)
@show value(y)
return
end
example()