Skip to main content

Problem in optimization of MIP model

Answered
atefeh Behzadi
Gurobi-versary
Conversationalist
First Question

Hello,

I have an optimization model in Gurobi optimization in python which has two binary variables and one non convex constraint, I have a problem with finding the initial point. I think The BestBd of Objective Bounds is low and model can not find the first initial point after times. the direction_AB and direction_BA are binary variables and the rest of variable are continuous.

This is important constraints and objective function in this code:

I have an optimization model in Gurobi optimization in Python, which includes two binary variables and one non-convex constraint. I have difficulty in finding the initial point. I believe the BestBd of Objective Bounds is low, causing the model to struggle in identifying the initial point after multiple attempts. The variables 'direction_AB' and 'direction_BA' are binary, while the remaining variables are continuous.

The following are crucial constraints and the objective function in this code:

model.addConstrs(
    (direction_AB[edge] + direction_BA[edge] <= 1) for edge in df_edges.index
)

for node in df_nodes.index:
    if (
        df_nodes.at[node, "node_type"] != "resource"
        and df_nodes.at[node, "node_type"] != "building"
    ):
        in_edges = df_edges[df_edges["B"] == node]
        out_edges = df_edges[df_edges["A"] == node]

        for edge_idx in in_edges.index:
            model.addConstr(
                power_inflow[edge_idx]
                == (power_edge[edge_idx]) * direction_AB[edge_idx]
                - (power_edge[edge_idx]) * (direction_BA[edge_idx])
            )
            # model.addConstr(power_inflow[edge_idx] == power_flow[edge_idx] * Bp[edge_idx])
        for edge_idx in out_edges.index:
            model.addConstr(
                power_outflow[edge_idx]
                == (power_edge[edge_idx]) * direction_AB[edge_idx]
                - (power_edge[edge_idx]) * (direction_BA[edge_idx])
            )
            # model.addConstr(power_outflow[edge_idx] == power_flow[edge_idx] * Bp[edge_idx])
        # Add energy balance constraint for the node
        # model.addConstr(power_inflow == power_outflow, name= 'Power balance')
        model.addConstr(
            quicksum(power_inflow[edge_idx] for edge_idx in in_edges.index)
            == quicksum(power_outflow[edge_idx] for edge_idx in out_edges.index)
        )

# For all pipes
for edge in df_edges.index:
    model.addConstr(
        mass_edge[edge] == (power_edge[edge]) / (cp * dT * (1 / (3.6e6)))
    )  ### kg/hr

fix = (4 * 1000 * 1000) / (3600 * 3.14 * ro * V_max)

for edge in df_edges.index:
    model.addConstr(D_edge[edge] * D_edge[edge] == fix * mass_edge[edge])

My objective function is such as follow:

model.addConstr(
    obj
    == quicksum(
        (5.36 * D_edge[edge] + 139.15)
        * df_edges.at[edge, "length"]
        * (direction_AB[edge] + direction_BA[edge])
        for edge in df_edges.index
    )
)
model.setObjective(obj, GRB.MINIMIZE)

 

Please provide guidance on optimizing my code. Is it possible for you to modify the constraints and objective function for greater efficiency, ensuring that Gurobi can solve it effectively?

0

Comments

9 comments

Date Votes
  • Riley Clement
    Gurobi Staff Gurobi Staff

    Hi Atefeh,

    I suggest replacing this:

    # For all pipes
for edge in df_edges.index:
    model.addConstr(
        mass_edge[edge] == (power_edge[edge]) / (cp * dT * (1 / (3.6e6)))
    )  ### kg/hr

fix = (4 * 1000 * 1000) / (3600 * 3.14 * ro * V_max)

for edge in df_edges.index:
    model.addConstr(D_edge[edge] * D_edge[edge] == fix * mass_edge[edge])

    with this:

    # For all pipes
for edge in df_edges.index:
    model.addConstr(
        cp * dT * mass_edge[edge] == 3.6e6 power_edge[edge]
    )  ### kg/hr

fix_num = 4 * 1000  # numerator
    fix_den = 3.6 * 3.14 * ro * V_max  # denominator

    for edge in df_edges.index:
    model.addConstr(fix_den * D_edge[edge] * D_edge[edge] == fix_num * mass_edge[edge])

    The modifications are done to avoid creating coefficients in your model which cannot be accurately represented by floating point arithmetic, which leads to numerical issues (which can contribute to difficulty finding solutions). We see this frequently when users convert from seconds to hours and divide by 3600.  Please our Guidelines for Numerical Issues in our reference manual, or our YouTube video Mastering Numerical Challenges for more information.

    If you see any warnings in your log relating to the NumericFocus parameter then I would follow the suggestions provided.

    You may also find experimenting with the following strategies useful in finding an initial solution:

    • Using our No Relaxation heuristic by setting NoRelHeurTime (or NoRelHeurWork) to a positive value
    • Setting MIPFocus=1.
    • Setting your objective to 0, so the problem becomes a feasibility problem, set SolutionLimit=1, then once a solution has been found, replace the objective function with the original objective function, call resetParams(), then optimize() to resume the optimization.

    - Riley

     

     

    0
  • atefeh Behzadi
    Gurobi-versary
    Conversationalist
    First Question

    Thank you very much for your assistance; I have implemented all of your suggestions. However, I have a challenge in finding the Incumbent of the objective function. I believe this issue may be attributed to a non-convex constraint. When I removed this particular constraint, the program successfully identified the Incumbent. Here is an image illustrating my result

    Do you have any suggestion for finding the Incumbent in my code? thanks again.

    0
  • Riley Clement
    Gurobi Staff Gurobi Staff

    What does your non-convex constraint look like?

    0
  • atefeh Behzadi
    Gurobi-versary
    Conversationalist
    First Question 
    fix_num = 4 * 1000  # numerator
    fix_den = 3.6 * 3.14 * ro * V_max  # denominator

    for edge in df_edges.index:
         model.addConstr(fix_den * D_edge_2[edge] * D_edge_2[edge] == fix_num * mass_edge[edge])
    0
  • Riley Clement
    Gurobi Staff Gurobi Staff

    Can you try fixing D_edge_2[edge] to a reasonable value, through lower and upper bounds, and solving the resulting convex problem to obtain a feasible solution?

    0
  • atefeh Behzadi
    Gurobi-versary
    Conversationalist
    First Question

    Yes, I have tried it, and I have set bounds for it. Also, These are the parameters that I have set: 

    model.Params.ScaleFlag = 2
    model.Params.MIPGap = 0.03
    # Set NonConvex parameter to 2
    model.setParam('NonConvex', 2)
    model.setParam("Threads", 4)
    model.setParam(GRB.Param.Method, 1)
    model.setParam(GRB.Param.Heuristics, 0.05)  # Determines the amount of time spent in MIP heuristics (%)
    #model.setParam(GRB.Param.LogToConsole, 1) 
    model.setParam(GRB.Param.MIPFocus, 1)
    model.setParam(GRB.Param.IntFeasTol, 1e-06)
    model.Params.OptimalityTol = 1e-06
    model.Params.IntFeasTol = 1e-06
    0
  • Matthias Miltenberger
    Gurobi Staff Gurobi Staff

    Please make sure to use the latest version of Gurobi (11.0). This might already improve the performance compared to previous versions.

    Independent of that, it can quickly happen that nonconvex models get very hard to solve, and you might want to think about a reformulation to make this constraint easier to handle.

    0
  • atefeh Behzadi
    Gurobi-versary
    Conversationalist
    First Question

    Thanks to your recommendation, I upgraded my Gurobi from version 10 to 11.0.0. However, I have issues with the updated version, and it is not functioning as expected, even for problems that worked seamlessly in the previous version. Specifically, there seems to be a problem with convex problems.

    Could you kindly provide guidance on potential issues with this new version of Gurobi? I am unable to execute the model that ran successfully in Gurobi 10. Your assistance would be greatly appreciated. Thank you.

    0
  • Riley Clement
    Gurobi Staff Gurobi Staff

    Hi Atefah,

    I am creating a ticket through our Help Centre to carry this discussion further. You will receive an email to the effect shortly.

    - Riley

    0

Please sign in to leave a comment.

Didn't find what you were looking for?

New post