"Как сделать суммирование к определенному периоду у Юлии?" - PullRequest
Новая
       15

"Как сделать суммирование к определенному периоду у Юлии?"

0 голосов
/ 28 марта 2019

Я формулирую математическое программирование, и теперь у меня есть проблема со знаком суммирования в ограничении. Я хочу подвести итоги за определенный период в горизонте планирования, я пробовал некоторые из них, но Джулия сообщает об ошибке.

математическая формула выглядит примерно так: constraint1

и это: ограничение2

вот код, который я пробовал: 

Horizon = 12
Section = 5
TBetPM = [4 6 9 8 5]
LPM = [1 4 5 4 4]
MaxPM = [9 8 7 10 6]
PrevPM = [3 3 2 5 2]

tam=zeros(Float64,1,5)
for i=1:Section
tam[i] = TBetPM[i]-LPM[i]
end

tar = zeros(Float64,1,5)
for i=1:Section
tar[i] = Pi[i]*(MaxPM[i]-PrevPM[i])-LPM[i]
end

@constraint(mod, [i=1:Section],
sum(m[i,t] for t=1:Horizon if t<=tam[i]) >= 1
)
@constraint(mod, [i=1:Section],
sum(r[i,t] for t=1:Horizon if t<=tar[i]) >= 1
)

Я тоже пробовал это, но оба не работают 

@constraint(mod, [i=1:Section],
sum(m[i,t] for t=1:tam[i]) >= 1
)
@constraint(mod, [i=1:Section],
sum(r[i,t] for t=1:tar[i]) >= 1
)

Заранее спасибо за все ответы:)

Ответы [ 2 ]

0 голосов
/ 28 марта 2019

Спасибо за ответ, Дж. Хамфоусон.

Вот полный код, который я пробовал: 

using JuMP, CPLEX, Gurobi, GLPKMathProgInterface
sex = Model(solver = GurobiSolver())
Horizon = 12
Section = 5
TBetPM = [4 6 9 8 5]
TlastPM = [0 0 0 0 0]
MaxPM = [9 8 7 10 6]
PrevPM = [3 3 2 5 2]
taf1=zeros(Float64,1,5)
for i=1:Section
taf1[i] = TBetPM[i]-TlastPM[i]
end
tafr1 = zeros(Float64,1,5)
for i=1:Section
tafr1[i] = TBetPM[i]*(MaxPM[i]-PrevPM[i])-TlastPM[i]
end
tafr = zeros(Float64,1,5)
for i=1:Section
tafr[i] = TBetPM[i]*MaxPM[i]
end
mdur = [9 6 8 10 3]
rdur = [18 16 23 16 12]
Maxdur = 24
mcost = [2 6 5.5 4 4]
rcost = [6 15 20 18 25]
scost = [0.6 1.17 0.81 0.66 1.4]
pcost = 2
ccost = 0.001
Ncus = 100
@variable(sex, m[1:Section,1:Horizon]>=0, Bin) # 1 if section s maintain in week h
@variable(sex, r[1:Section,1:Horizon]>=0, Bin) # 1 if section s renew in week h
@variable(sex, p[1:Horizon]>=0,Bin) #1 if possession in week h

@objective(sex, Min, sum(pcost*p[t] for t=1:Horizon)+
sum(mcost[i]*m[i,t] for t=1:Horizon, i=1:Section)+
sum(rcost[i]*r[i,t] for t=1:Horizon, i=1:Section)
)

#select first maintenance time (correct)
@constraint(sex, [i=1:Section],
sum(m[i,t] for t=1:Horizon if t==taf1[i]) >= 1
)
#select next maintenance time (correct)
@constraint(sex, [i=1:Section, k=1:(Horizon-TBetPM[i])],
sum(m[i,t] for t=1+k:TBetPM[i]+k) >= 1
)
#select first renewal time
@constraint(sex, [i=1:Section],
sum(r[i,t] for t in tafr1[i]) >= 1
)
#select next renewal time
@constraint(sex, [i=1:Section, k=1:(Horizon-tafr[i])],
sum(r[i,t] for t=1+k:tafr[i]+k) >= 1
)
# if there is maintenance, there is possession
@constraint(sex, [i=1:Section, h=1:Horizon],
m[i,h] <= p[h]
)
# if there is renewal, there is possession
@constraint(sex, [i=1:Section, h=1:Horizon],
r[i,h] <= p[h]
)


solve(sex)
#print(sex)
mVal = getvalue(m)
for i=1:Section, t=1:Horizon
    if (mVal[i,t] ==1)
        println("section[$i] repair in period [$t]", )
    end
end
rVal = getvalue(r)
for i=1:Section, t=1:Horizon
    if (rVal[i,t] ==1)
        println("section[$i] renew in period [$t]", )
    end
end
0 голосов
/ 28 марта 2019

Вы использовали модель mod в макросе @constraint. Это не было определено. Вы должны создать его, если используете JuMP. 

using JuMP
mod = Model()

Вы предварительно инициализировали tam и tar как два размерных массива размером 1x5 (матрицы). Я думаю, что вам нужен одномерный массив, поскольку вы обращаетесь к ним как к векторам: 

tam = zeros(Float64,5)
tar = zeros(Float64,5)

Вы не определили переменные m и r в вашей модели: 

 @variable(mod, m[i=1:Section, t=1:Horizon])
 @variable(mod, r[i=1:Section, t=1:Horizon])

Наконец, вы можете захотеть, чтобы JuMP решил вашу модель, это можно сделать с помощью: 

using GLPK # feel free here to use the solver you prefer 
optimize!(mod, with_optimizer(GLPK.Optimizer))

И печать решений: 

 if termination_status(mod) == MOI.OPTIMAL
     optimal_solution = value.(m), value.(r)
     optimal_objective = objective_value(mod)

     @show optimal_solution
     @show optimal_objective
 else
     error("The model was not solved correctly.")
 end

Весь рабочий код (Julia v1.1.0, JuMP v0.19.0, GPLK v0.9.1): 

Horizon = 12
Section = 5
TBetPM = [4 6 9 8 5]
LPM = [1 4 5 4 4]
MaxPM = [9 8 7 10 6]
PrevPM = [3 3 2 5 2]


using JuMP
using GLPK # feel free here to use the solver you prefer

mod = Model()

tam = zeros(Float64,5)
for i in 1:Section
    tam[i] = TBetPM[i]-LPM[i]
end

tar = zeros(Float64,5)
for i=1:Section
    tar[i] = pi*(MaxPM[i] - PrevPM[i]) - LPM[i]
end


@variable(mod, m[i=1:Section, t=1:Horizon])
@variable(mod, r[i=1:Section, t=1:Horizon])

@constraint(mod, [i=1:Section],
    sum(m[i,t] for t in 1:Horizon if t <= tam[i]) >= 1)
@constraint(mod, [i=1:Section],
    sum(r[i,t] for t in 1:Horizon if t <= tar[i]) >= 1)


# Solve model and printing solution
optimize!(mod, with_optimizer(GLPK.Optimizer))
if termination_status(mod) == MOI.OPTIMAL
    optimal_solution = value.(m), value.(r)
    optimal_objective = objective_value(mod)

    @show optimal_solution
    @show optimal_objective
else
    error("The model was not solved correctly.")
end

Я не знаю, происходит ли это из-за копирования / вставки или чего-то еще, но вы должны сделать отступ в своем коде, даже если его компиляция не обязательна: p

...