x_old = zeros(_ndofs)
for i in 1:10
func(x) = f(x, x_old)
model = Model(Ipopt.Optimizer)
@variable(model, x[1:_ndofs])
@objective(model, Min, func(x) )
for i in dofs_y_top
@constraint(model, x[i] == - i/10)
end
for i in dofs_x_top
@constraint(model, x[i] == i/10)
end
optimize!(model)
solution = [value(x[i]) for i = 1:_ndofs]
x_old = solution
end

I wonder if there is a better way to implement above problem? Especially, can I use previous solution for the next optimization step ? Do I have to generate the model every time or is there a way to say to the solver just to increase some constraints value? Is it actually a good way of imposing values to some unknowns?

It depends a little on what f is. Can you share a reproducible example?

Is this the full problem you are trying to solve? Have you benchmarked this to see if problem creation is a bottleneck?

Note also that your code uses i for two different loops. Is the model building what you expect?

You could try something like this:

model = Model(Ipopt.Optimizer)
@variable(model, x[1:_ndofs])
for i in dofs_y_top
fix(x[i], -i / 10)
end
for i in dofs_x_top
fix(x[i], i / 10)
end
x_old = zeros(_ndofs)
for _ in 1:10
@objective(model, Min, f(x, x_old))
optimize!(model)
x_old = value.(x)
end

x_old = zeros(_ndofs)
for k in 1:10
func(x) = f(x, x_old)
model = Model(Ipopt.Optimizer)
@variable(model, x[1:_ndofs])
@objective(model, Min, func(x) )
for i in dofs_y_top
@constraint(model, x[i] == - k/10)
end
for i in dofs_x_top
@constraint(model, x[i] == k/10)
end
optimize!(model)
solution = [value(x[i]) for i = 1:_ndofs]
x_old = solution
end

Generally I want to use JuMP (with Ipopt) with finite element method (FEM) and my function f is just energy of a body so it is the sum over elements and integration points and I am doing this as follows:

function f(x, x_old)
acc = zero(QuadExpr)
for n in 1:n_elements
energy = element_energy(x, x_old)
JuMP.add_to_expression!(acc, energy)
end
return acc
end

so, having e.g. 1 mln elements it is painful but probably can be parallelized - I didn’t think about it yet. My constraints denote imposed deformation (i.e. displacements on the boundary of the body) but what is important for me, the energy is history dependent (and non-convex) so I can not compute whole deformation in one shot so I wanted to use loop for it and incrementally increase displacement what I do by applying constraints. Minimization after each loop should start from the previous solution (is it possible?), and the general question was: is it possible to generate the model once and then freely impose increased constraints after each loop ?

See Nonlinear Modeling · JuMP (We had some scattered docs for this, but I recently made this much more prominent)

You probably want something like this:

model = Model(Ipopt.Optimizer)
@variable(model, x[1:_ndofs])
@variable(model, x_old[1:_ndofs] in Parameter(0))
@objective(model, Min, f(x, x_old))
solution = zeros(_ndofs)
for k in 1:10
set_start_value.(x, solution)
set_parameter_value.(x_old, solution)
for i in dofs_y_top
fix(x[i], -k / 10)
end
for i in dofs_x_top
fix(x[i], k / 10)
end
optimize!(model)
solution .= value.(x)
end