How can I find which equations are infeasible in PyPSA lopf?

[marilena]

Hello all, 

I want to optimize a simple unit commitment problem using PyPSA. The optimization is infeasible, and I am trying to find out which equations or variable bounds are causing the infeasibility. 

To do so, I am trying to call the functions that are in the pyomo.util.infeasible module. But I am not certain, where/how/if they should be called! 

In particular, I am trying something like: 

...

import logging

logging.basicConfig(level=logging.ERROR)

...

pypsa.opf.logger.setLevel(logging.INFO)

from pyomo.util.infeasible import log_infeasible_constraints 

 ...

network.lopf(network.snapshots)

log_infeasible_constraints(network)

but this does not return any log information regarding the infeasibility. 

Next, I thought that maybe I should call the infeasibility functions within the opf.py file, maybe even within the network_lopf function, but that didn't work either.  

Final trial I made was to change keep_files to True and search for the log files, but I was unable to read the files generated in raw and lp format... 

It would be great if you could provide guidance as to where/how I should "ask" pyPSA to show me the infeasible constraints. 

Thanks for taking the time to read my question!

Marilena

[Fabian Neumann]

Hello Marilena,

I had a brief look at
https://github.com/Pyomo/pyomo/blob/master/pyomo/util/infeasible.py.

First of all, make sure you call log_infeasible_constraints on the pyomo
model, not the PyPSA network, i.e.
log_infeasible_constraints(network.model). Also, from what I read, this
function might not tell you all combinations of constraints which are
infeasible.

Some solvers have the capability to find a so-called irreducible
infeasible set/subsystem (IIS), i.e. a (minimal) set of constraints
which make the problem infeasible (feasible set = {}). Solvers, which
can provide this are e.g. Gurobi or CPLEX (there are probably others as
well).

For example code with pyomo see:
*
https://github.com/Pyomo/pyomo/blob/master/examples/pyomo/suffixes/gurobi_ampl_iis.py
* https://github.com/Pyomo/pyomo/pull/596

If you already have a clue which subset of constraints might be causing
the infeasibility, you can also check manually by relaxing the
constraints in question and rerun the problem.

I hope I could help,

Fabian

> --
> You received this message because you are subscribed to the Google
> Groups "pypsa" group.
> To unsubscribe from this group and stop receiving emails from it, send
> an email to pypsa+un...@googlegroups.com
> <mailto:pypsa+un...@googlegroups.com>.
> To post to this group, send email to py...@googlegroups.com
> <mailto:py...@googlegroups.com>.
> To view this discussion on the web, visit
> https://groups.google.com/d/msgid/pypsa/280e4ac9-08a5-426f-9c9e-d0bd55f8595a%40googlegroups.com
> <https://groups.google.com/d/msgid/pypsa/280e4ac9-08a5-426f-9c9e-d0bd55f8595a%40googlegroups.com?utm_medium=email&utm_source=footer>.
> For more options, visit https://groups.google.com/d/optout.

--
Karlsruhe Institute of Technology (KIT)
Institute for Automation and Applied Informatics (IAI)

Fabian Neumann, M.Sc.
PhD Student in Energy System Modelling

Phone: +49 721 608 25707
Mobile: +49 171 2943831
Fax: +49 721 608 22602
Website: https://www.iai.kit.edu/ESM.php
Github: fneum
LinkedIn: fabian-neumann-renewables
ORCID: 0000 0001 8551 1480

Building 445, Office 304
Campus North
Hermann-von-Helmholtz-Platz 1
76344 Eggenstein-Leopoldshafen

KIT – the Research University in the Helmholtz Association

[A.Magh]

Hey Marilena,

were you able to log the infeasible constraints ? 

I am facing a similar problem, and when i try to use the gurobi IIS() function that Fabian mentioned, it doesn't work because IIS() only works with abstract models. Pypsa lopf() generates concrete models, not abstract.

Ammar.

[marilena zambara]

Hello Amar, 

Using Gurobi did help in my case.. In order to get the proper reporting I made a specific inquiry to "Gurobi-support" and they recommended following the info in this link: https://stackoverflow.com/questions/51522718/get-gurobi-iis-using-pyomo-and-gurobipy/51994135#51994135. I basically added the highlighted lines in my code in order to "ask" gurobi to generate an ilp file, named model.ilp, which contains the infeasible constraints info. 

....

...

solver_parameters = "ResultFile=model.ilp" # write an ILP file to print the IIS
model = pypsa.opf.network_lopf_build_model(network, network.snapshots)
extra_functionality(network, network.snapshots)
opt = pypsa.opf.network_lopf_prepare_solver(network, solver_name="gurobi")
network.results=opt.solve(model, options_string=solver_parameters)
network.results.write()

Note that I am not using network.lopf() directly... This is not related to using or not the IIS. If you are using network.lopf() I am not sure where to include the "options_string=solver_parameters", - maybe someone else in the forum can help? 

After you succeed in generating the ilp file, it will look something like this: 

\ Model x112194_copy
\ LP format - for model browsing. Use MPS format to capture full model detail.
Minimize
 
Subject To
 c_u_x322884_: x108479 <= 50000
 c_e_x323627_: - x106250 - x108479 + x106993 - x107736 = -416085.9319441
Bounds
 -infinity <= x106250 <= 50000
 x108479 free
 -infinity <= x107736 <= 130.3763441
End

The above is the model.ilp for a dummy example I was trying. Understanding what each x is referring to is hard. I was able to understand what is going on because it was easy to understand what the rhs values where (the 50000, -416085.9319441 , 130.3763441), my model is relatively small... 

As it gets larger it will be harder to clear it out without some kind of equation-naming ... But it did work. 

I hope this can help you too! 

Meanwhile, finding ways to locate infeasibility issues is an ongoing task for me as I was only using Gurobi through a demo license... Now I am back to glpk! 

So, I will be in touch in case I find alternate ways!

Best regards, 

Marilena

--

You received this message because you are subscribed to the Google Groups "pypsa" group.

To unsubscribe from this group and stop receiving emails from it, send an email to pypsa+un...@googlegroups.com.
To post to this group, send email to py...@googlegroups.com.
To view this discussion on the web, visit https://groups.google.com/d/msgid/pypsa/2f32bc4a-4c05-41d3-aca8-dad6b0172a29%40googlegroups.com.

[Ramiz Qussous]

If you are using network.lopf() I am not sure where to include the "options_string=solver_

parameters", - maybe someone else in the forum can help? 

This can be done using

network.lopf(solver_name='gurobi', solver_options={'ResultFile':'model.ilp'})

To unsubscribe from this group and stop receiving emails from it, send an email to py...@googlegroups.com.

[Fajar Andy Setyawan]

Hi Marilena, 

Have you figured out how to locate the infeasibility issues using GLPK?
Currently, I run around 200 buses and got infeasibility issue.

Thanks