# Help on model a problem/function to use with Optimization.jl
**URL:** https://discourse.julialang.org/t/help-on-model-a-problem-function-to-use-with-optimization-jl/95942
**Category:** Optimization (Mathematical)
**Tags:** question
**Created:** [March 11, 2023, 10:53pm UTC](https://discourse.julialang.org/t/help-on-model-a-problem-function-to-use-with-optimization-jl/95942 "2023-03-11T22:53:38Z")
**Posts on this page:** 18
**Page:** 1
### Author:  [@rodolforbcoutinho](https://discourse.julialang.org/u/rodolforbcoutinho)
#### Post date: [March 11, 2023, 10:53pm UTC](https://discourse.julialang.org/t/help-on-model-a-problem-function-to-use-with-optimization-jl/95942/1 "2023-03-11T22:53:38Z")
Hello to all.
I’m trying to use Optimization.jl to solve a problem, but a I dont know how to define the problem in a way the package needs it. In my problem I have an Array of a Struct with several parameters and for each element of the Array I need the proper index of a field that is a StepRange type.
I’ll write here a smaller version of my struct and equation. Any help or hint is really appreciated, because I dont know even how to start it. Thanks in advance!!!
```julia
#Struct with problem values - needed an Array of it
struct L
a::Float64
b::Float64
s::StepRange{Time,Minute}
t1::Time
t2::Time
t3::Time
L() = new()
function L(a,b,t1,t2,t3)
s = t1:t2:t3
new(a,b,s,t1,t2,t3)
end
end
#GOAL is find all values of x that minimize this function
function Eval_EQ_SUM(H::Array{L,1},x::Array{Int64,1})
if length(H) ≠ length(x)
@error "H and x should have the same length"
return nothing
end
S = Vector{Float64}(undef,length(H))
for i in eachindex(H)
S[i] = H[i].a*EQ(H[i].s[x[i]]), H[i].b)
end
return sum(S)
end
#auxiliar function for better reading
function EQ(t::time, b::Float64)
if t < Time(2)
F = 2*b
elseif t > Time(4) && t < Time(12)
F = sqrt(b)
else
F = b^3
end
return F
end
```
---
### Author:  [@ChrisRackauckas](https://discourse.julialang.org/u/ChrisRackauckas)
#### Post date: [March 12, 2023, 11:32am UTC](https://discourse.julialang.org/t/help-on-model-a-problem-function-to-use-with-optimization-jl/95942/2 "2023-03-12T11:32:53Z")
It looks like this is an integer programming problem?
---
### Author:  [@rodolforbcoutinho](https://discourse.julialang.org/u/rodolforbcoutinho)
#### Post date: [March 12, 2023, 12:18pm UTC](https://discourse.julialang.org/t/help-on-model-a-problem-function-to-use-with-optimization-jl/95942/3 "2023-03-12T12:18:16Z")
Yes, values of x[i] are index to a range H[i].s and its values should be restricted to bounds of the vector H[i].s.
I’m not an expert in optimization at all but i think using integer variables should be easy to solve the problem. Past colleagues used Genetic Algorithm and PSO to solve MIN the equation, but the programing interface was Matlab/Python. Thats why Im looking into Optimization.jl
---
### Author:  [@gdalle](https://discourse.julialang.org/u/gdalle)
#### Post date: [March 12, 2023, 12:22pm UTC](https://discourse.julialang.org/t/help-on-model-a-problem-function-to-use-with-optimization-jl/95942/4 "2023-03-12T12:22:43Z")
To me it looks more like a constraint programming problem, cause the syntax `s[x[i]]` is not natively supported by integer programming solvers
---
### Author:  [@rodolforbcoutinho](https://discourse.julialang.org/u/rodolforbcoutinho)
#### Post date: [March 12, 2023, 12:33pm UTC](https://discourse.julialang.org/t/help-on-model-a-problem-function-to-use-with-optimization-jl/95942/5 "2023-03-12T12:33:54Z")
Well, as I said, I’m not an expert in optimization.
I’m just thinking, ‘x’ is integer so it should be an integer programming problem.
If it could help, I can replace `H[i].s[x[i]]` for a variable `t[i], typeof(t) = Vector{Time}`
My original problem is a schedule type. I should position the task within an window of time (t1:t3) in a way it minimize a not continuous function based on the schedule point
---
### Author:  [@ChrisRackauckas](https://discourse.julialang.org/u/ChrisRackauckas)
#### Post date: [March 12, 2023, 5:06pm UTC](https://discourse.julialang.org/t/help-on-model-a-problem-function-to-use-with-optimization-jl/95942/6 "2023-03-12T17:06:31Z")
This is probably more in the domain of JuMP
---
### Author:  [@gdalle](https://discourse.julialang.org/u/gdalle)
#### Post date: [March 12, 2023, 5:23pm UTC](https://discourse.julialang.org/t/help-on-model-a-problem-function-to-use-with-optimization-jl/95942/7 "2023-03-12T17:23:24Z")
@rodolforbcoutinho could you point me to the mathematical formulation of your optimization problem?
---
### Author:  [@rodolforbcoutinho](https://discourse.julialang.org/u/rodolforbcoutinho)
#### Post date: [March 12, 2023, 6:33pm UTC](https://discourse.julialang.org/t/help-on-model-a-problem-function-to-use-with-optimization-jl/95942/8 "2023-03-12T18:33:47Z")
Here are some link for the papers I’m trying to reproduce the base model.
> **[Proposed Architecture for Energy Efficiency and Comfort Optimization in Smart...](https://link.springer.com/article/10.1007/s40313-018-0410-y)**
>
> In this paper a smart home controller proposal is formalized as a multi-objective integer linear programming problem that minimizes energy consumption and maximizes comfort. A comfort objective function is tested for several tariff scenarios...
> **[An Improved Optimization Function to Integrate the User’s Comfort Perception...](https://www.mdpi.com/1424-8220/23/6/3021)**
>
> Scheduling residential loads for financial savings and user comfort may be performed by smart home controllers (SHCs). For this purpose, the electricity utility’s tariff variation costs, the lowest tariff cost schedules, the user’s preferences, and...
The image above also have the main equation I’m trying to work with. As you could see its a bit more complicated system than the one I wrote in the initial post. But I’m sure if I could model the function **`Eval_EQ_SUM`** , I could model the complete system.
I found a way to eliminate the inner Sum, but to do that I should call a function with some if/else logic. Dont see how could i put it in a single statemente like the one in examples (2x + 3y)
I’ll do some tests now with JuMP.jl. Its a way more clear in its tutorial “_who goes where_”.
In/with JuMP.jl, can I do something like this ???
```julia
function f(a)
if a > 10
return 2*a
elseif a < 4
return 3*a
else
return a^2
end
end
@objective(model, Min, f(x) + 20y)
```
---
### Author:  [@odow](https://discourse.julialang.org/u/odow)
#### Post date: [March 12, 2023, 6:47pm UTC](https://discourse.julialang.org/t/help-on-model-a-problem-function-to-use-with-optimization-jl/95942/9 "2023-03-12T18:47:31Z")
> In/with JuMP.jl, can I do something like this ???
It depends on the rest of the problem structure. But that objective function probably doesn’t make sense, because it is discontinuous.
If you mean something like your `EQ` function, is `t` a decision variable? If it is not, then yes.
I’d start by providing a reproducible example of what you’re trying to do. Use the `f(a)` syntax for now, even if JuMP complains and gives an error. Once things are started, people can probably point you in the right direction.
---
### Author:  [@rodolforbcoutinho](https://discourse.julialang.org/u/rodolforbcoutinho)
#### Post date: [March 12, 2023, 7:09pm UTC](https://discourse.julialang.org/t/help-on-model-a-problem-function-to-use-with-optimization-jl/95942/10 "2023-03-12T19:09:28Z")
The function I must handle its **discontinuous**. As I said before its something realy like f(a) in the last post. It allow me to cut of about 80% memory allocation and its about 10x faster than the original formulation. So I’ll keep trying on.
Also **`vector t`** is my decision variable.
But I’ll try to do like you said, and bring on the real problem.
---
### Author:  [@gdalle](https://discourse.julialang.org/u/gdalle)
#### Post date: [March 13, 2023, 9:07am UTC](https://discourse.julialang.org/t/help-on-model-a-problem-function-to-use-with-optimization-jl/95942/11 "2023-03-13T09:07:50Z")
I think you might want to try actually modeling the variable `u` in your JuMP problem, and using an integer programming solver like HiGHS
---
### Author:  [@rodolforbcoutinho](https://discourse.julialang.org/u/rodolforbcoutinho)
#### Post date: [March 26, 2023, 9:21pm UTC](https://discourse.julialang.org/t/help-on-model-a-problem-function-to-use-with-optimization-jl/95942/12 "2023-03-26T21:21:21Z")
Hello again!
Finally I have some real working code to share!
**But First what didn’t worked at all**
```julia
function JuMP_Cost_alfa(L::BasicLoad)
model = Model(HiGHS.Optimizer)
u = length(L.s)
@variable(model, 1 < x < u, Int)
@objective(model, Min, Cost(L,x))
optimize!(model)
return value.(x)
end
```
**Now the code that works**
```julia
function JuMP_Cost(L::BasicLoad)
model = Model(HiGHS.Optimizer)
u = length(L.s)
@variable(model, x[1:u], Bin)
@constraint(model,sum(x)==1)
@objective(model, Min, sum(Cost(L,i)*x[i] for i in 1:u) )
optimize!(model)
A = value.(x)
ix = findfirst(x->x==1.0,A)
return ix
end
```
**Some Complementar code**
```julia
#this is discontinuous function with lots of logic constraints but is working fine
#just to point the use of index i inside L.s[i]
Cost(L,i) = Eval_BasicLoadCost(L, L.s[i])
```
```julia
#this is the struct i'm working on
struct BasicLoad
# (_release time_)
r::Time
#(_expected ON time_)
e::Time
# (_dead line_)
d::Time
#(_start time_)
s::StepRange{Time, Minute}
# discrete time step
Δt::Time
#Load duration time
C::Time
#discrete time duration
CΔ::Int64
#Power over time model
# OBS: typeof(P) == Array{Number,1} OR Function t -> f(t)
# OBS: length(P{Array}) = CΔ
P
#Load Average Power
P_av::Number
#Load Peak Power
P_pk::Number
BasicLoad() = new()
function BasicLoad(r,e,d,Δt,C,P; allow_rand_peak = true)
CΔ = div(C.instant,Δt.instant)
s = r:Minute(Δt):(Time(d-C))
if typeof(P) <: Function
V = P.(1:CΔ)
elseif typeof(P) <: Array
if length(P) == CΔ
V = P
else
@error "Unable to create new Load. P::Array{Number,1} must be compatible with Δt size"
return nothing
end
else
@error "Unable to create new Load. P must be a Function t -> f(t) or Array{Number,1} type"
return nothing
end
P_av = mean(V)
P_pk = maximum(V)
if (P_pk == P_av)&&(allow_rand_peak)
P_pk = P_av*rand(1.1:0.1:1.5)
end
new(r,e,d,s,Δt,C,CΔ,P,P_av,P_pk)
end
end
```
---
### Author:  [@odow](https://discourse.julialang.org/u/odow)
#### Post date: [March 26, 2023, 9:25pm UTC](https://discourse.julialang.org/t/help-on-model-a-problem-function-to-use-with-optimization-jl/95942/13 "2023-03-26T21:25:29Z")
> [@rodolforbcoutinho](#):
>
> ```julia
> function JuMP_Cost(L::BasicLoad)
>
> model = Model(HiGHS.Optimizer)
> u = length(L.s)
> @variable(model, x[1:u], Bin)
> @constraint(model,sum(x)==1)
> @objective(model, Min, sum(Cost(L,i)*x[i] for i in 1:u) )
> optimize!(model)
>
> A = value.(x)
> ix = findfirst(x->x==1.0,A)
> return ix
> end
>
> ```
Yes, this approach looks good. You can’t use user-defined functions like `Cost(L, x)` with HiGHS. Everything must be explicitly written as a MILP.
You should be careful with:
```julia
A = value.(x)
ix = findfirst(x->x==1.0,A)
```
HiGHS uses numerical tolerances, so you’re not guaranteed to return a value `==1.0`.
Do instead:
```julia
A = value.(x)
ix = findfirst(x -> x > 0.5, A)
# or
ix = findfirst(xi -> round(Int, value(xi)) == 1, x)
```
---
### Author:  [@rodolforbcoutinho](https://discourse.julialang.org/u/rodolforbcoutinho)
#### Post date: [March 26, 2023, 9:27pm UTC](https://discourse.julialang.org/t/help-on-model-a-problem-function-to-use-with-optimization-jl/95942/14 "2023-03-26T21:27:15Z")
Now I’m trying to solve the same problem, but for an Array of the struct BasicLoad
Here the last code, it has a bug i couldn’t find.
```julia
function JuMP_VectorCost(H::Vector{BasicLoad})
model = Model(HiGHS.Optimizer)
n = length(H)
u = Vector{Int64}(undef,n)
for i in 1:n
u[i] = length(H[i].s)
end
@variable(model, x[i = n, j = u[i]], Bin)
@constraint(model, [i = n, j = u[i]], sum(x[i,:])==1)
@objective(model, Min, sum(Cost(H[i],j)*x[i,j] for i in 1:n, j in 1:u[i]) )
optimize!(model)
return value.(x)
end
```
_This is the error message_
```julia
KeyError: key (1, 1) not found
1. **getindex** @*dict.jl:498* [inlined]
2. **getindex** (::JuMP.Containers.SparseAxisArray{JuMP.VariableRef, 2, Tuple{Int64, Int64}}, ::Int64, ::Int64)@*SparseAxisArray.jl:128*
3. **macro expansion** @*[Other: 322](http://localhost:1234/edit?id=7f4cfdf0-cbfd-11ed-3860-9df06b1a7405#)* [inlined]
4. **macro expansion** @*[Other: 1297](http://localhost:1234/edit?id=7f4cfdf0-cbfd-11ed-3860-9df06b1a7405#)* [inlined]
5. **JuMP_VectorCost** (::Vector{Main.var"RAW_IDEAS.jl".BasicLoad})@*[Other: 13](http://localhost:1234/edit?id=7f4cfdf0-cbfd-11ed-3860-9df06b1a7405#)*
6. **top-level scope** @*[Local: 1](http://localhost:1234/edit?id=7f4cfdf0-cbfd-11ed-3860-9df06b1a7405#)* [inlined]
```
Any help on the `JuMP_VectorCost function` would be nice! Thanks in advance!!!
---
### Author:  [@rodolforbcoutinho](https://discourse.julialang.org/u/rodolforbcoutinho)
#### Post date: [March 26, 2023, 9:32pm UTC](https://discourse.julialang.org/t/help-on-model-a-problem-function-to-use-with-optimization-jl/95942/15 "2023-03-26T21:32:43Z")
> [@odow](#):
>
> ```julia
> A = value.(x)
> ix = findfirst(x -> x > 0.5, A)
> # or
> ix = findfirst(xi -> round(Int, value(xi)) == 1, x)
>
> ```
Yeah thats a better idea indeed! The other code there was a “late night” solution kind.
---
### Author:  [@odow](https://discourse.julialang.org/u/odow)
#### Post date: [March 26, 2023, 9:55pm UTC](https://discourse.julialang.org/t/help-on-model-a-problem-function-to-use-with-optimization-jl/95942/16 "2023-03-26T21:55:14Z")
> `@variable(model, x[i = n, j = u[i]], Bin)`
You probably need `@variable(model, x[i = 1:n, j = 1:u[i]], Bin)`
---
### Author:  [@rodolforbcoutinho](https://discourse.julialang.org/u/rodolforbcoutinho)
#### Post date: [March 27, 2023, 12:56am UTC](https://discourse.julialang.org/t/help-on-model-a-problem-function-to-use-with-optimization-jl/95942/17 "2023-03-27T00:56:03Z")
It works nice now!!! I could even do it multiobjective with JuMP and HiGHS.
In full problem I should minimize COST and maximize COMFORT
Bellow the final solution for the record
```julia
function JuMP_VectorMulti(H::Vector{Raw.BasicLoad})
model = Model()
set_optimizer(model, () -> MOA.Optimizer(HiGHS.Optimizer))
n = length(H)
u = Vector{Int64}(undef,n)
for i in 1:n
u[i] = length(H[i].s)
end
@variable(model, x[i = 1:n, j = 1:u[i]], Bin)
@constraint(model, [i = 1:n, j = 1:u[i]], sum(x[i,:])==1)
#MIN
@expression(model, Cost_expr, sum(Cost(H[i],j)*x[i,j] for i in 1:n, j in 1:u[i]))
#MAX
@expression(model, Comf_expr, sum(Comf(H[i],j)*x[i,j] for i in 1:n, j in 1:u[i]))
@objective(model, Max, [-Cost_expr, Comf_expr])
optimize!(model)
solution_summary(model)
return value.(x)
end
```
Should I change the title or add some other mark? maybe this solution could be usefull to others.
---
### Author:  [@odow](https://discourse.julialang.org/u/odow)
#### Post date: [March 27, 2023, 2:25am UTC](https://discourse.julialang.org/t/help-on-model-a-problem-function-to-use-with-optimization-jl/95942/18 "2023-03-27T02:25:55Z")
> [@rodolforbcoutinho](#):
>
> `@constraint(model, [i = 1:n, j = 1:u[i]], sum(x[i,:])==1)`
I think this should be `@constraint(model, [i = 1:n], sum(x[i, :]) == 1)`.
You could also try different algorithms in MOA that return a set of solutions, not just a single point. See [Multi-objective knapsack · JuMP](https://jump.dev/JuMP.jl/stable/tutorials/linear/multi_objective_knapsack/).