Search code examples
pythonoptimizationconstraintspyomo

Pyomo Constraint iteration for Optimal Power Flow


I am solving a DC Optimal Power Flow Problem and I am trying to brainstorm the most efficient way to iterate over a constraint in Pyomo.

The following is the data structure, where i and k are the nodes connected through a branch, and X is the reactance, a property of the branch.

Sample Branch Data

The constraint I am having trouble with is the following:

constraint

Where the symbol "delta" and "p" is a variable in the constraint, each node has a single value of delta and p. What this constraint basically does, is it insures that all the powers flowing into the node i, from all the connected nodes k, are equal to the existing power value in the same node.

Here is an example for i=1 and i=2, iterations of the constraint.

Sample constraint

So I am trying to find the most efficient way to state this constraint in pyomo, where instead of having several constraint iterations written like this:

def P1_rule(modelo):
    return modelo.p[0]-L[0]== ((modelo.d[0]-modelo.d[1])/0.1)+((modelo.d[0]-modelo.d[2])/0.1)
model.P1 = Constraint(rule=P1_rule)

def P2_rule(modelo):
    return modelo.p[1]-L[1]==((modelo.d[1]-modelo.d[0])/0.1)+((modelo.d[1]-modelo.d[2])/0.1)
model.P2 = Constraint(rule=P2_rule)

def P3_rule(modelo):
    return modelo.p[2]-L[2] ==((modelo.d[2]-modelo.d[0])/0.1)+((modelo.d[2]-modelo.d[1])/0.1)
model.P3 = Constraint(rule=P3_rule)

I want a single line like this, so that it easily generalizes over a huge network:

def P3_rule(modelo):
    return modelo.p[i] ==((modelo.d[i]-modelo.d[k])/X[k])
model.P3 = Constraint(rule=P3_rule)

I have came up with a way that includes restructuring the data and creating new arrays of indices etc... I would like to see if i can apply the constraint using the data keeping the same structure and more directly.


Solution

  • Ok so I figured out how to do it. The best way, which I didn't know was possible, is through making an if statement within the summation, so basically make a full conditional iteration inside of the summation. In the code bellow, G is the list of nodes, while "From" and "To" are the branch number columns in the FullBranch data table.

    def Try_rule(mod,g):
            return mod.p[g] - L[g] == sum((mod.d[i-1]-mod.d[k-1])/FullBranch.loc[x,"X"] for x,(i,k) in enumerate(zip(FullBranch["From"], FullBranch["To"])) if i == g+1)
    model.Try = Constraint(G,rule=Try_rule)