Set min/max constraint for each technology per country

[Heronimonimo]

What would be the advisable way to set a minimum/maximum constraint for technologies per country?

I was thinking about using a .csv with the same structure as the costs.csv. This could than be loaded at the add_electricity phase. Would this be the right way?

[anyexingluzk]

As far as I understand, the add_electricity rule produces the non-clustered network. I am not sure if you can set the p_min/max on the country level. I am also trying to implement the same idea as you do, and my initial plan is loading it in prepare_network phase. Let's see if it works :)

[Heronimonimo]

True, it's probably something that has to be put in as an extra optimization constraint at that state like the line cost limit.

Let me know if you have something for me to test with or have a look at.

[Heronimonimo]

@anyexingluzk: any news on this?

[anyexingluzk]

@Heronimonimo I have implemented the minimum capacity for a certain technology per country by adding p_nom_min = value, and p_nom_extendable = True in the prepare_network.py script.

[Heronimonimo]

@anyexingluzk Nice! Can you make a pull request of your changes?

[Heronimonimo]

@anyexingluzk Would you mind sharing the code you used to implement this?

[anyexingluzk]

@Heronimonimo The way I have implemented is quite straightforward as below:

def replace_pmin_by_pnom(n):
    index = (n.generators.p_nom != 0) & (n.generators.p_nom_extendable == True) 
    n.generators.loc[index, 'p_nom_min'] = n.generators.loc[index, 'p_nom']
    return n

Basically, assign the value of p_nom to p_nom_min for extendable generators. Regarding the value of p_nom, it depends on your model, where I use existing generators in my case. It is also possible to load a csv file which includes all the limits.

[Heronimonimo]

If I understand your code correctly it sets the nominal power of a generator as the minimum installed capacity. Does this also work for wind and solar?

[anyexingluzk]

It works if the p_nom_extendable is True.

[Heronimonimo]

Ah. And how would I be able to define the p_nom value for wind/solar per country?

(Sorry for al these questions, I'm quite new to PyPSA and Python)

[anyexingluzk]

No worries. I have been using PyPSA for two years, though I feel like I am still a rookie :D It is a bit complicated to add the existing renewable generators to the model, and I have attempted to make a pull request but not succeeded yet.

In a nutshell, I used the existing renewable power stations dataset from OPSD, which consists of 6 European countries. Then I create a csv file similar to powerplants.csv with the same columns except the bus, which is assigned to each generator using the build_powerplants.py script afterwards. In the end, exsiting wind and solar are attached in the add_electricity.py script similar to add conventional generatos.

[Heronimonimo]

Ah, that seems like a nice approach.

I think I'm aiming at a slightly different goal. What I would like to test is how far the national goals for 2030 and 2040 of different EU countries are from the 'optimal' EU wide coordination installed capacities. Therefore I would like to set the minimum/maximum capacity per renewable technology on a country level. The exact location of the RES installations in each country are somewhat less relevant.

P.S. Are you attending the Aarhus OpenMod Workshop?

[anyexingluzk]

If you aim to set the limits on a country level, it would be easier to implement it by reading a CSV file in the prepare_network stage after clustering into one-node-per-country.

I am one of the organisers as I work in the Aarhus Sustainable Energy System group :)

[Heronimonimo]

That's a possible approach. But that would not work when running the model with more than one cluster per country. Maybe we can set a 'country' constraint at the custer_network stage?

Maybe we can set up a session to work together on implementing new ideas and fixing outstanding bugs in PyPSA-eur during the Aarhus workshop?

[coroa]

indeed, the proposed approach will not allow you to set per country constraints for networks with more than one cluster.

PyPSA has no knowledge about constraints aggregated over more than one bus, with the notable exception of per-carrier global constraints.

Thus, you will have to modify the pyomo model after it is generated by PyPSA in solve_network.py:

The main task is to

add a function similar to add_opts_constraints from the run_lopf function (just as with add_opts_constraints). Perhaps call it add_country_carrier_generation_constraints.

This function will have to add a constraint for each (country, carrier). if we say, you already parsed your config into a pandas series mapping from country, carrier to minimal capacity, like in pseudo from

agg_p_nom_min =
(DE, "onwind") | 100 000
(FR, "onwind") | 50 000

then your function will look something like

def add_agg_p_nom_constraint(net):

    # read in agg_p_nom_min somehow from csv or config file (exercise
    # left to the reader)

    gen_in_country_b = net.generators.bus.map(net.buses.country)

    def agg_p_nom_min_rule(model, country_carrier):
        country, carrier = country_carrier
        return sum(model.generator_p_nom[gen]
            for gen in net.generators.index[(gen_in_country_b == country) & (net.generators.carrier == carrier)])
        <= agg_p_nom_min.at[country, carrier]

    n.model.agg_p_nom_min = po.Constraint(n.model, list(agg_p_nom_min.index), rule=agg_p_nom_min_rule)

Best, Jonas

P.S.: I'd attend and support an openmod PyPSA-Eur-hack-a-ton, but I will not lead it!

[coroa]

(Wasn't there markdown support in the issue tracker before, confused)

[coroa]

If you get a minimal version of this running and come up with a good csv format or yaml config file encoding, feel free to open a pull request for adding it as a feature.

[Heronimonimo]

indeed, the proposed approach will not allow you to set per country constraints for networks with more than one cluster. PyPSA has no knowledge about constraints aggregated over more than one bus, with the notable exception of per-carrier global constraints. Thus, you will have to modify the pyomo model after it is generated by PyPSA in solve_network.py: The main task is to add a function similar to add_opts_constraints from the run_lopf function (just as with add_opts_constraints). Perhaps call it add_country_carrier_generation_constraints. This function will have to add a constraint for each (country, carrier). if we say, you already parsed your config into a pandas series mapping from country, carrier to minimal capacity, like in pseudo from agg_p_nom_min = (DE, "onwind") | 100 000 (FR, "onwind") | 50 000 then your function will look something like python def add_agg_p_nom_constraint(net): # read in agg_p_nom_min somehow from csv or config file (exercise # left to the reader) gen_in_country_b = net.generators.bus.map(net.buses.country) def agg_p_nom_min_rule(model, country_carrier): country, carrier = country_carrier return sum(model.generator_p_nom[gen] for gen in net.generators.index[(gen_in_country_b == country) & (net.generators.carrier == carrier)]) <= agg_p_nom_min.at[country, carrier] n.model.agg_p_nom_min = po.Constraint(n.model, list(agg_p_nom_min.index), rule=agg_p_nom_min_rule) Best, Jonas P.S.: I'd attend and support an openmod PyPSA-Eur-hack-a-ton, but I will not lead it!

Unfortunately I missed the Aarhus workshop due to a funeral.

Based on the code above I made this, but I can't get it to work:

def add_country_carrier_generation_constraints(n, opts=None):
    agg_p_nom_min = pd.read_csv("data/agg_p_nom_min.csv", index_col=list(range(2))).sort_index()

    gen_in_country_b = n.generators.bus.map(n.buses.country)

    def agg_p_nom_min_rule(model, country, carrier):
        country, carrier = country_carrier
        return sum(model.generator_p_nom[gen]
            for gen in n.generators.index[(gen_in_country_b == country) & (n.generators.carrier == carrier)]) >= agg_p_nom_min.at[country, carrier]
    n.model.agg_p_nom_min = pypsa.opt.Constraint(n.model, list(agg_p_nom_min.index), rule=agg_p_nom_min_rule)

Tried a lot of different things, but just can't get it to work. But I think it's close...

[coroa]

Hi,

i think the main problems are figuring out how Pyomo calls the rule function:

1. Comparing with the simple example from the Pyomo docs, shows that the Constraint call does not take any reference to the model.

2. I THINK pyomo will pass in the list(agg_p_nom_min.index) as tuples

so the following MAY work:

def add_country_carrier_generation_constraints(n, opts=None):
    agg_p_nom_min = pd.read_csv("data/agg_p_nom_min.csv", index_col=list(range(2))).sort_index()

    gen_country = n.generators.bus.map(n.buses.country)

    def agg_p_nom_min_rule(model, (country, carrier)):
        return sum(model.generator_p_nom[gen]
                   for gen in n.generators.index[(gen_country == country)
        & (n.generators.carrier == carrier)]) >= agg_p_nom_min.at[country, carrier])
    n.model.agg_p_nom_min = pypsa.opt.Constraint(list(agg_p_nom_min.index), rule=agg_p_nom_min_rule)

But, I would have to try myself. That means, I would -- in the jupyter notebook or ipython -- run the following:

import pypsa
n = pypsa.Network("networks/elec_s_45.nc")
model = pypsa.opf.network_lopf_build_model(n, n.snapshots[:2], formulation="kirchhoff")

# model is now a 2-hour pyomo model to which you can try to add your constraint:

agg_p_nom_min = pd.read_csv("data/agg_p_nom_min.csv", index_col=list(range(2))).sort_index()

gen_country = n.generators.bus.map(n.buses.country)

def agg_p_nom_min_rule(model, (country, carrier)):
    return sum(model.generator_p_nom[gen]
                for gen in n.generators.index[(gen_country == country)
    & (n.generators.carrier == carrier)]) >= agg_p_nom_min.at[country, carrier])

model.agg_p_nom_min = pypsa.opt.Constraint(list(agg_p_nom_min.index), rule=agg_p_nom_min_rule)

[Heronimonimo]

Hi,

i think the main problems are figuring out how Pyomo calls the rule function:

1. Comparing with the simple example from the Pyomo docs, shows that the Constraint call does not take any reference to the model.
2. I THINK pyomo will pass in the list(agg_p_nom_min.index) as tuples

so the following MAY work:

def add_country_carrier_generation_constraints(n, opts=None):
    agg_p_nom_min = pd.read_csv("data/agg_p_nom_min.csv", index_col=list(range(2))).sort_index()

    gen_country = n.generators.bus.map(n.buses.country)

    def agg_p_nom_min_rule(model, (country, carrier)):
        return sum(model.generator_p_nom[gen]
                   for gen in n.generators.index[(gen_country == country)
        & (n.generators.carrier == carrier)]) >= agg_p_nom_min.at[country, carrier])
    n.model.agg_p_nom_min = pypsa.opt.Constraint(list(agg_p_nom_min.index), rule=agg_p_nom_min_rule)

But, I would have to try myself. That means, I would -- in the jupyter notebook or ipython -- run the following:

import pypsa
n = pypsa.Network("networks/elec_s_45.nc")
model = pypsa.opf.network_lopf_build_model(n, n.snapshots[:2], formulation="kirchhoff")

# model is now a 2-hour pyomo model to which you can try to add your constraint:

agg_p_nom_min = pd.read_csv("data/agg_p_nom_min.csv", index_col=list(range(2))).sort_index()

gen_country = n.generators.bus.map(n.buses.country)

def agg_p_nom_min_rule(model, (country, carrier)):
    return sum(model.generator_p_nom[gen]
                for gen in n.generators.index[(gen_country == country)
    & (n.generators.carrier == carrier)]) >= agg_p_nom_min.at[country, carrier])

model.agg_p_nom_min = pypsa.opt.Constraint(list(agg_p_nom_min.index), rule=agg_p_nom_min_rule)

Using extra () brackets in a definition gives a syntax error.

Removing the brackets gives the following error:

ERROR:pyomo.core:Rule failed when generating expression for constraint agg_p_nom_min with index ('BE', 'solar'):
KeyError: 'solar'
ERROR:pyomo.core:Constructing component 'agg_p_nom_min' from data=None failed:
KeyError: 'solar'

Not experienced enough unfortunately to think of a solution right away.

[coroa]

Using extra () brackets in a definition gives a syntax error.

Ah, it had escaped my notice, that they removed tuple parameter unpacking in Python 3. It is equivalent to:

def agg_p_nom_min_rule(model, country_carrier):
    country, carrier = country_carrier
    return sum(model.generator_p_nom[gen]
                for gen in n.generators.index[(gen_country == country)
    & (n.generators.carrier == carrier)]) >= agg_p_nom_min.at[country, carrier])

[Heronimonimo]

Testing it like this still gives me an error:

import pypsa
n = pypsa.Network("networks/elec_s_12.nc")
model = pypsa.opf.network_lopf_build_model(n, n.snapshots[:2], formulation="kirchhoff")

# model is now a 2-hour pyomo model to which you can try to add your constraint:

agg_p_nom_min = pd.read_csv("data/agg_p_nom_min.csv", index_col=list(range(2))).sort_index()

gen_country = n.generators.bus.map(n.buses.country)

def agg_p_nom_min_rule(model, country_carrier):
    country, carrier = country_carrier
    return sum(model.generator_p_nom[gen]
                for gen in n.generators.index[(gen_country == country)
    & (n.generators.carrier == carrier)]) >= agg_p_nom_min.at[country, carrier]

model.agg_p_nom_min = pypsa.opt.Constraint(list(agg_p_nom_min.index), rule=agg_p_nom_min_rule)

It seems that passing the double index is not working:

ERROR: Rule failed when generating expression for constraint agg_p_nom_min
    with index ('BE', 'solar'): TypeError: agg_p_nom_min_rule() takes 2
    positional arguments but 3 were given
ERROR: Constructing component 'agg_p_nom_min' from data=None failed:
    TypeError: agg_p_nom_min_rule() takes 2 positional arguments but 3 were
    given

[coroa]

Did you know you can easily debug your code in ipython or the notebook by running %debug just after an error happened; from there it is easy to determine that tuples are unpacked by Pyomo.

There was one more detail wrong. It's quite unintuitive how to use .at with MultiIndices. All in all:

agg_p_nom_min = pd.read_csv("data/agg_p_nom_min.csv", index_col=list(range(2))).sort_index()

gen_country = n.generators.bus.map(n.buses.country)

def agg_p_nom_min_rule(model, country, carrier):
    return sum(model.generator_p_nom[gen]
                for gen in n.generators.index[(gen_country == country)
    & (n.generators.carrier == carrier)]) >= agg_p_nom_min.at[(country, carrier),]

model.agg_p_nom_min = pypsa.opt.Constraint(list(agg_p_nom_min.index), rule=agg_p_nom_min_rule)

[Heronimonimo]

Did you know you can easily debug your code in ipython or the notebook by running %debug just after an error happened; from there it is easy to determine that tuples are unpacked by Pyomo.

There was one more detail wrong. It's quite unintuitive how to use .at with MultiIndices. All in all:

agg_p_nom_min = pd.read_csv("data/agg_p_nom_min.csv", index_col=list(range(2))).sort_index()

gen_country = n.generators.bus.map(n.buses.country)

def agg_p_nom_min_rule(model, country, carrier):
    return sum(model.generator_p_nom[gen]
                for gen in n.generators.index[(gen_country == country)
    & (n.generators.carrier == carrier)]) >= agg_p_nom_min.at[(country, carrier),]

model.agg_p_nom_min = pypsa.opt.Constraint(list(agg_p_nom_min.index), rule=agg_p_nom_min_rule)

Still not getting it to work. Rewrote it to:

agg_p_nom_min = pd.read_csv("data/agg_p_nom_min.csv", index_col=list(range(1))).sort_index()

gen_country = n.generators.bus.map(n.buses.country)

def agg_p_nom_min_rule(model, country_carrier):
    country, carrier = country_carrier.split("_")
    gens = n.generators.index[(gen_country == country) & (n.generators.carrier == carrier)]
    return sum(model.generator_p_nom[gen] for gen in gens) <= agg_p_nom_min.loc[country_carrier,'capacity']

for country_carrier in agg_p_nom_min.index:
    model.agg_p_nom_min = pypsa.opt.Constraint(country_carrier, rule=agg_p_nom_min_rule)

Now getting:

WARNING:pyomo.core:Element f already exists in set agg_p_nom_min_index; no action taken.
ERROR:pyomo.core:Rule failed when generating expression for constraint agg_p_nom_min with index c:
ValueError: not enough values to unpack (expected 2, got 1)
ERROR:pyomo.core:Constructing component 'agg_p_nom_min' from data=None failed:
ValueError: not enough values to unpack (expected 2, got 1)

[coroa]

The version in my comment above is correct and tested (with a synthetic data series of course).

[Heronimonimo]

Rewrote it to this and now I got it running. But unfortunately the capacities in the results don't seem to change:

def add_country_carrier_generation_constraints(n, opts=None):
    agg_p_nom_min = pd.read_csv("data/agg_p_nom_min.csv", index_col=list(range(1))).sort_index()

    gen_country = n.generators.bus.map(n.buses.country)

    for country_carrier in agg_p_nom_min.index:
        country, carrier = country_carrier.split("_")
        mincap = int(agg_p_nom_min.at[country_carrier,'capacity'])
        gens = n.generators.index[(gen_country == country) & (n.generators.carrier == carrier)]
        n.model.agg_p_nom_min = pypsa.opt.Constraint(expr=sum(n.model.generator_p_nom[gen] for gen in gens) >= mincap)

The version in my comment above is correct and tested (with a synthetic data series of course).

What is the structure of your .csv?

[coroa]

You are repeatedly overwriting the same constraint, only the last country_carrier combination will stay, the others are gone.

What goes wrong with the one, I showed above?

What is the structure of your .csv?

I didn't use one, I only used a hard-coded two country example.

In [3]: pd.Series([1,2], pd.MultiIndex.from_tuples([('DE', "onwind"), ('FR', "onwind")]))
Out[3]: 
DE  onwind    1
FR  onwind    2
dtype: int64

[Heronimonimo]

Found the problem! :)

Because the csv has colnames, the '.at' becomes: agg_p_nom_min.at[(country, carrier),'capacity']

Will now test the code and make both a min and max option.

[Heronimonimo]

Well, that was a good learning experience. Thanks for all your support @coroa !

I made a pull request #46 including all the changes. Please feel free to make further improvements and/or changes!

[coroa]

Fixed by #46