メインコンテンツへスキップ

Optimizing integral of exponential function

進行中

コメント

5件のコメント

  • Jaromił Najman
    • Gurobi Staff Gurobi Staff

    Hi Elina,

    If you can compute a closed form for \(\int_{0}^{t_i}r_i \cdot e^{-\gamma\cdot x} dx\) for every \(i\), then you should be able to use Gurobi for your model.

    If your \(r_i\) is a constant which you can compute for every \(i\), then it should hold that

    \[\begin{align*}
    \int_{0}^{t_i}r_i \cdot e^{-\gamma\cdot x} dx = \frac{r_i-r_i \cdot e^{-\gamma t_i}}{r_i}
    \end{align*}\]
    Of course, the form of the above strongly depends on whether any of \(\gamma, \mu, \Theta\) are just constant values or optimization variables. Nevertheless, if you are able to compute a closed form, you can use the addGenConstrExp method to model the exponential function.

    If you cannot compute a closed form of the integral for every \(i\), then I think that it is best to try a different approach or a different solver which can handle integration. You should be able to find some in the field of "dynamic optimization".

    Best regards, 
    Jaromił

    0
  • Elina Joksch
    • Gurobi-versary
    • Conversationalist
    • Curious

    Hello Jaromil,

    thank you very much for your quick reply!

    I tried to model this problem, but the optimal result is currently 0, which cannot be optimal. What mistake am I making? γ,μ are constants, but Θ and t I want to optimize. 

     

    from gurobipy import Model, GRB, quicksum

    c = 1
    k = 2
    a = 0.005
    b = 1/600
    T = 8*60
    e = 5
    num_work = k+1
    num_bre = k
    m = Model("Optimization")
    m.params.NonConvex = 2

    t = m.addVars(num_work, lb=0, ub=T, name="work")
    u = m.addVars(num_bre, lb=0, ub=T, name="break")
    r = m.addVars(num_work, lb=0, ub=c, name="work rate")
    l = m.addVars(num_work, name="exponent")
    y = m.addVars(num_work, name="e-Funktion")
    w = m.addVars(num_bre, name="exponent2")
    p = m.addVars(num_bre, name="e-Funktion2")
    q = m.addVars(num_work, name="Objective Function")
    z = m.addVars(num_bre, name="work rate2")

    #constraints for fatigue e-function
    m.addConstr((-b)*t[0] == l[0])
    m.addGenConstrExp(l[0],y[0])
    m.addConstr((-b)*t[1] == l[1])
    m.addGenConstrExp(l[1],y[1])
    m.addConstr((-b)*t[2] == l[2])
    m.addGenConstrExp(l[2],y[2])

    #constraints for recovery e-function
    m.addConstr((-a)*u[0] == w[0])
    m.addGenConstrExp(w[0],p[0])
    m.addConstr((-a)*u[1] == w[1])
    m.addGenConstrExp(w[1],p[1])

    # r constraint
    m.addConstr(1 == r[0])
    m.addConstr(r[0]*y[0] == z[0])
    m.addConstr((c+(z[0]-c) * p[0]) == r[1])
    m.addConstr(r[1]*y[1] == z[1])
    m.addConstr(c+(z[1]-c) * p[1] == r[2])

    m.addConstr((quicksum(t[i] for i in range(num_work)) + quicksum(u[i] for i in range(num_bre))) <= T)

    #objective function
    m.addConstr((quicksum(((r[i] - r[i]*y[i])/b) for i in range(num_work)) == quicksum(q[i] for i in range(num_work))))

    m.setObjective(quicksum(q[i] for i in range(num_work)), GRB.MAXIMIZE)

    m.optimize()
    0
  • Jaromił Najman
    • Gurobi Staff Gurobi Staff

    Your exponents can attain negative values. However, Gurobi default value for variable lower bounds is \(0\).

    Allowing for negative lower bounds for your exponents seems to do the trick.

    l = m.addVars(num_work, lb=-GRB.INFINITY, name="exponent")
    w = m.addVars(num_bre, lb=-GRB.INFINITY, name="exponent2")

    Best regards, 
    Jaromił

    0
  • Elina Joksch
    • Gurobi-versary
    • Conversationalist
    • Curious

    Thank you very much! This helped a lot!

    0
  • Elina Joksch
    • Gurobi-versary
    • Conversationalist
    • Curious

    I am once again facing a problem..I am now trying to reformulate the parts formulated with numbers into more general expressions, but it is not working. Is there another way to approach this with exponential functions in Gurobi?

    m.addConstr((-b)*t[0] == l[0])
    m.addGenConstrExp(l[0],y[0])
    m.addConstr((-b)*t[1] == l[1])
    m.addGenConstrExp(l[1],y[1])
    m.addConstr((-b)*t[2] == l[2])
    m.addGenConstrExp(l[2],y[2])

    into

    m.addConstr(((-b)*t[i] for i in range(num_work)) == (l[i] for i in range(num_work)))
    m.addGenConstrExp((l[i] for i in range(num_work)),(y[i] for i in range(num_work)))
    0

サインインしてコメントを残してください。