How to use current_constraint_limits and define c_rating_a? Testing and Documentation enhancement?

[friederikemeier]

Hi,

I am relatively new to PowerModels and am a former Matpower/Pandapower/Pypower user. I am struggling with the definitions of current_contraint_limits. I found the test for this in test/opf.jl in Lines 47ff and extended it a little bit to understand, what's happening. In the test case "case_5clm.m", the c_rating_a for e.g. branch 7 is limited at 240 (and then 2.4 in the PowerModels branch Dict). What unit is this? Is this in MVA at 1 p.u. as in Matpower? This would equal a current of 0.602452 kA.

Below you can see, how I calculated the resulting current on this branch. It results in 0.7206 kA. It is dar too high...

What I also found out, is that the apparent power at the to bus of this branch equals 2.64 and before I got the warning

this code only supports positive rate_a values, changing the value on branch 7 to 264.0

... But why do I get this 264.0 here in the warning? Obviously, the optimzation considers only this constraint for this branch and not the more restricitve current constraint... The warning comes from the call of calc_thermal_limits!. I get that it was constrained to the max. current at the max. voltage. But because the voltage in the results is quite low, the current is getting too high at this power value. Also, the resulting current nearly stays the same, regardless of the call of calc_thermal_limits!.

data = PowerModels.parse_file("case5_clm.m")
# in MPC: limit for branch 1: 240
for n in range(1,stop=7)
        n = string(n)
        data["branch"][n]["c_rating_a"] *= 1 
end

ipopt_solver = JuMP.optimizer_with_attributes(Ipopt.Optimizer,
                                        "print_level" => 0,
                                        "max_cpu_time" => 1e3,
                                        "tol" => 1e-8)

calc_thermal_limits!(data)
result = run_ac_opf(data, ipopt_solver)
println(result["termination_status"] == LOCALLY_SOLVED)
println(isapprox(result["objective"], 16513.6; atol = 1e0))

for n in range(1,stop=7)

        n = string(n)
        println("initial c_rating_a: " , data["branch"][n]["c_rating_a"])
        i_max = data["branch"][n]["c_rating_a"]/230/sqrt(3)*100 # 4.0
        println("according current in kA: ", i_max)
        sf = sqrt(result["solution"]["branch"][n]["qf"]^2+
                        result["solution"]["branch"][n]["pf"]^2)
        st = sqrt(result["solution"]["branch"][n]["qt"]^2+
                        result["solution"]["branch"][n]["pt"]^2)

        # voltage at from bus and to bus
        fb=data["branch"][n]["f_bus"]
        tb=data["branch"][n]["t_bus"]
        vm_f = result["solution"]["bus"][string(fb)]["vm"]
        vm_t = result["solution"]["bus"][string(tb)]["vm"]

        i_f = sf / vm_f / sqrt(3) / 230 * 100 
        i_t = st / vm_t / sqrt(3) / 230 * 100 
        println("resulting current:  ", i_t)

        f_loading = i_f/i_max
        t_loading = i_t/i_max
        println("resulting loading:  ", t_loading)
        println("reulting apparent power", st)
end

How should the c_rating_a be defined? What kind of test could I add for this? How can I get the Optimization to respect the current constraint?

I appreciate any kind of help on this and am ready to contribute if needed :)

Kind regards,

Friederike

[ccoffrin]

Hi @friederikemeier! Thanks for giving PowerModels a try.

Following the Matpower and PSSE formats, MVA flow constraints are the standard in PowerModels. There is some support for current based flow constraints, but these are a little experimental and not well documented; as you have noticed.

In the case of case5_clm.m the current values should be provided in the units of MA, as these will be rescaled by the MVA base to be compatible with the p.u. scale.

The helper functions calc_thermal_limits! and calc_current_limits! can be used to make a course conversion between the power and current ratings in a model. Because you don't know a precise voltage a priori, there is no exact one-to-one conversion between these models; the translation instead looks at the voltage bounds and makes a conservative approximation that does not remove any valid operating point within these voltage bounds.

To solve a problem with current flow limits, rather than the default MVA limits, you need to use a different model definition. For example see, _run_opf_cl, which us used for testing this feature.

[friederikemeier]

Thank you very much for the hint! I just changed everything to the call of '_run_opf_cl' and it works smoothly :) In pandapower, the line loading is always defined with respect to the current... So I would like to make this available for the pandapower-PowerModels interface. My next step is to integrate this feature into the pandapower-PowerModels interface and then I can further test it.

By the way, I completely forgot to tell you how great this package is! I am so used to working with the Pypower OPF or "homemade python optimizers" e.g. using cvxopt / linsolve. Getting OPF results in under a second still seems weird :)

I will keep you posted when I update the pandapower interface!

[ccoffrin]

That's great! I am glad to hear that PowerModels is working well for you. If the _run_opf_cl suits your needs I will consider making it a standard OPF problem specification for use in pandapower.

Let me know one your integration is working and we can close this issue.

[friederikemeier]

I was able to reproduce the test with case5_clm.m from pandapower today completely. What I stumbled upon was, that the baseMVA is not flexible with _run_opf_cl? For the other pandapower-PowerModels conversions, it was possible to give different baseMVAs. But maybe it is better to fix this to 100MVA in the pandapower-PowerModels interface? I think it shouldn't matter for our users what baseMVA is used internally and I feel like you mostly use 100 MVA for your testing? The updated interface will follow next week or so, after some more testing ;)

[ccoffrin]

That sounds like a bug, if different base-mva values work for run_opf it should also be working with _run_opf_cl as well. A key point to check is the PandaPower function that converts the PandaPower data model into the PowerModel data model. It needs to be setup to convert the current limits into p.u. I presume, can you check on that?

[friederikemeier]

Hi, Here's an example of what I mean concerning the change of baseMVA. I used the file for your test opf-var.jl from line 25.

using PowerModels
using Ipopt
import JuMP

pm = PowerModels.parse_file("case5_clm.m")
ipopt_solver = JuMP.optimizer_with_attributes(Ipopt.Optimizer, "print_level" => 0)

result_100 = PowerModels._run_opf_cl(pm, ACPPowerModel, ipopt_solver)

    # print assertions from original test
println(result_100["termination_status"] == LOCALLY_SOLVED)
println(isapprox(result_100["objective"], 17015.5; atol = 1e0))

### now convert to baseMVA = 1
pm["baseMVA"] = 1

for (key, value) in pm["gen"]
   value["pmax"] *= 100
   value["qmax"] *= 100
   value["qmin"] *= 100
   value["pg"] *= 100
   value["qg"] *= 100
   value["cost"] *= 0.01
end

for (key, value) in pm["branch"]
   value["c_rating_a"] *= 100
end

for (key, value) in pm["load"]
   value["pd"] *= 100
   value["qd"] *= 100
end

result_1 = PowerModels._run_opf_cl(pm, ACPPowerModel, ipopt_solver)
println(result_1["termination_status"] == LOCALLY_SOLVED)
println(isapprox(result_1["objective"], 17015.5; atol = 1e0))

After converting it to baseMVA = 1, it doesn't solve anymore.

[ccoffrin]

Ok I understand your question. At a first glance something like this could work,

using PowerModels
using Ipopt
import JuMP

pm = PowerModels.parse_file("case5_clm.m")
ipopt_solver = JuMP.optimizer_with_attributes(Ipopt.Optimizer, "print_level" => 0)

result_100 = PowerModels._run_opf_cl(pm, ACPPowerModel, ipopt_solver)
println(result_100["termination_status"] == LOCALLY_SOLVED)
println(isapprox(result_100["objective"], 17015.5; atol = 1e0))

### now convert to baseMVA = 1
make_mixed_units!(pm)
pm["baseMVA"] = 1
make_per_unit!(pm)

result_1 = PowerModels._run_opf_cl(pm, ACPPowerModel, ipopt_solver)
println(result_1["termination_status"] == LOCALLY_SOLVED)
println(isapprox(result_1["objective"], 17015.5; atol = 1e0))

The problem is that in Matpower, PSSE, and similar transmission level data models (inlc. PowerModels), the line parameters are provided already in per-unit. So if you change the given "baseMVA" value the line parameters also need to be adjusted and this is not done by default. For example, this adjustment almost works,

using PowerModels
using Ipopt
import JuMP

pm = PowerModels.parse_file("case5_clm.m")
ipopt_solver = JuMP.optimizer_with_attributes(Ipopt.Optimizer, "print_level" => 0)

result_100 = PowerModels._run_opf_cl(pm, ACPPowerModel, ipopt_solver)
println(result_100["termination_status"] == LOCALLY_SOLVED)
println(isapprox(result_100["objective"], 17015.5; atol = 1e0))

### now convert to baseMVA = 1
make_mixed_units!(pm)
pm["baseMVA"] = 1
make_per_unit!(pm)

for (i,branch) in pm["branch"]
   branch["br_r"] = branch["br_r"]/100
   branch["br_x"] = branch["br_x"]/100
   branch["b_fr"] = branch["b_fr"]*100
   branch["b_to"] = branch["b_to"]*100
end

result_1 = PowerModels._run_opf_cl(pm, ACPPowerModel, ipopt_solver)
println(result_1["termination_status"] == LOCALLY_SOLVED)
println(isapprox(result_1["objective"], 17015.5; atol = 1e0))

A careful derivation would be required to get the line parameter conversion right. One issue I recall is that not all datasets have good base_kv values, which presents some difficulties in changing the p.u. base.

If your final goal is to integrate with pandapower my recommendation would be to have the pandapower data export function put everything into p.u. before sending it to PowerModels. In this approach pandapower can focus on getting the units correct and PowerModels can focus on the optimization of models in p.u., which is its focus.

However, if you prefer to work in PowerModels exclusively we could consider adding support for changing of the base mva when all of the required data is available.

[ccoffrin]

@friederikemeier, if there are no further questions I propose to close this issue.

[friederikemeier]

@ccoffrin Thank you for all your advice, I am currently working on improving the pandapower interface to powermodels and will let you know, once I finish it. So I am closing it here.

[friederikemeier]

Thanks for all the advice, I was able to push this as a new feature to the pandapower-powermodels interface. Also, I was benchmarking the powermodels results against matpower/pypower a lot with the current constraints. We will soon make this the standard way for calling the powermodels OPF. Why is _run_opf_cl still an "experimental" feature? Works just like it is supposed to! 👍

[ccoffrin]

I am glad to hear this and that _run_opf_cl works as expected. I suppose I consider this feature still a little experimental because I don't use it, but now that you say it behaves as expected, then maybe we consider it complete once issue #625 is closed.