Simple electricity market examples

This Jupyter notebook is meant for teaching purposes. To use it, you need to install a Python environment with Jupyter notebooks, and the Python for Power System Analysis (PyPSA) library. See

https://pypsa.org/doc/installation.html

for tips on installation.

It gradually builds up more and more complicated energy-only electricity markets in PyPSA, starting from a single bidding zone, going up to multiple bidding zones connected with transmission (NTCs) along with variable renewables and storage.

Available as a Jupyter notebook at http://www.pypsa.org/examples/simple-electricity-market-examples.ipynb.

Preliminaries

Here libraries are imported and data is defined.

In [1]:
import pypsa, numpy as np
In [2]:
#marginal costs in EUR/MWh
marginal_costs = {"Wind" : 0,
                  "Hydro" : 0,
                  "Coal" : 30,
                  "Gas" : 60,
                  "Oil" : 80}

#power plant capacities (nominal powers in MW) in each country (not necessarily realistic)
power_plant_p_nom = {"South Africa" : {"Coal" : 35000,
                                       "Wind" : 3000,
                                       "Gas" : 8000,
                                       "Oil" : 2000
                                      },
                     "Mozambique" : {"Hydro" : 1200,
                                    },
                     "Swaziland" : {"Hydro" : 600,
                                    },
                    }

#transmission capacities in MW (not necessarily realistic)
transmission = {"South Africa" : {"Mozambique" : 500,
                                  "Swaziland" : 250},
                "Mozambique" : {"Swaziland" : 100}}

#country electrical loads in MW (not necessarily realistic)
loads = {"South Africa" : 42000,
         "Mozambique" : 650,
         "Swaziland" : 250}

Single bidding zone with fixed load, one period

In this example we consider a single market bidding zone, South Africa.

The inelastic load has essentially infinite marginal utility (or higher than the marginal cost of any generator).

In [3]:
country = "South Africa"

network = pypsa.Network()

network.add("Bus",country)

for tech in power_plant_p_nom[country]:
    network.add("Generator",
                "{} {}".format(country,tech),
                bus=country,
                p_nom=power_plant_p_nom[country][tech],
                marginal_cost=marginal_costs[tech])


network.add("Load",
            "{} load".format(country),
            bus=country,
            p_set=loads[country])
In [4]:
#Run optimisation to determine market dispatch
network.lopf()
INFO:pypsa.pf:Slack bus for sub-network 0 is South Africa
INFO:pypsa.opf:Performed preliminary steps
INFO:pypsa.opf:Building pyomo model using `angles` formulation
INFO:pypsa.opf:Solving model using glpk
INFO:pypsa.opf:Optimization successful
Out[4]:
('ok', 'optimal')
In [5]:
#print the load active power (P) consumption
network.loads_t.p
Out[5]:
South Africa load
now 42000.0
In [6]:
#print the generator active power (P) dispatch
network.generators_t.p
Out[6]:
South Africa Oil South Africa Wind South Africa Gas South Africa Coal
now 0.0 3000.0 4000.0 35000.0
In [7]:
#print the clearing price (corresponding to gas)
network.buses_t.marginal_price
Out[7]:
South Africa
now 60.0

Two bidding zones connected by transmission, one period

In this example we have bidirectional transmission capacity between two bidding zones. The power transfer is treated as controllable (like an A/NTC (Available/Net Transfer Capacity) or HVDC line). Note that in the physical grid, power flows passively according to the network impedances.

In [8]:
network = pypsa.Network()

countries = ["Mozambique", "South Africa"]

for country in countries:
    
    network.add("Bus",country)

    for tech in power_plant_p_nom[country]:
        network.add("Generator",
                    "{} {}".format(country,tech),
                    bus=country,
                    p_nom=power_plant_p_nom[country][tech],
                    marginal_cost=marginal_costs[tech])


    network.add("Load",
                "{} load".format(country),
                bus=country,
                p_set=loads[country])
    
    #add transmission as controllable Link
    if country not in transmission:
        continue
    
    for other_country in countries:
        if other_country not in transmission[country]:
            continue
        
        #NB: Link is by default unidirectional, so have to set p_min_pu = -1
        #to allow bidirectional (i.e. also negative) flow
        network.add("Link",
                   "{} - {} link".format(country, other_country),
                   bus0=country,
                   bus1=other_country,
                   p_nom=transmission[country][other_country],
                   p_min_pu=-1)
In [9]:
network.lopf()
INFO:pypsa.pf:Slack bus for sub-network 0 is Mozambique
INFO:pypsa.pf:Slack bus for sub-network 1 is South Africa
INFO:pypsa.opf:Performed preliminary steps
INFO:pypsa.opf:Building pyomo model using `angles` formulation
INFO:pypsa.opf:Solving model using glpk
INFO:pypsa.opf:Optimization successful
Out[9]:
('ok', 'optimal')
In [10]:
network.loads_t.p
Out[10]:
Mozambique load South Africa load
now 650.0 42000.0
In [11]:
network.generators_t.p
Out[11]:
Mozambique Hydro South Africa Oil South Africa Wind South Africa Gas South Africa Coal
now 1150.0 0.0 3000.0 3500.0 35000.0
In [12]:
network.links_t.p0
Out[12]:
South Africa - Mozambique link
now -500.0
In [13]:
#print the clearing price (corresponding to water in Mozambique and gas in SA)
network.buses_t.marginal_price
Out[13]:
Mozambique South Africa
now 0.0 60.0
In [15]:
#link shadow prices
network.links_t.mu_lower
Out[15]:
South Africa - Mozambique link
now 60.0

Three bidding zones connected by transmission, one period

In this example we have bidirectional transmission capacity between three bidding zones. The power transfer is treated as controllable (like an A/NTC (Available/Net Transfer Capacity) or HVDC line). Note that in the physical grid, power flows passively according to the network impedances.

In [16]:
network = pypsa.Network()

countries = ["Swaziland", "Mozambique", "South Africa"]

for country in countries:
    
    network.add("Bus",country)

    for tech in power_plant_p_nom[country]:
        network.add("Generator",
                    "{} {}".format(country,tech),
                    bus=country,
                    p_nom=power_plant_p_nom[country][tech],
                    marginal_cost=marginal_costs[tech])


    network.add("Load",
                "{} load".format(country),
                bus=country,
                p_set=loads[country])
    
    #add transmission as controllable Link
    if country not in transmission:
        continue
    
    for other_country in countries:
        if other_country not in transmission[country]:
            continue
        
        #NB: Link is by default unidirectional, so have to set p_min_pu = -1
        #to allow bidirectional (i.e. also negative) flow
        network.add("Link",
                   "{} - {} link".format(country, other_country),
                   bus0=country,
                   bus1=other_country,
                   p_nom=transmission[country][other_country],
                   p_min_pu=-1)
In [17]:
network.lopf()
INFO:pypsa.pf:Slack bus for sub-network 0 is Swaziland
INFO:pypsa.pf:Slack bus for sub-network 1 is Mozambique
INFO:pypsa.pf:Slack bus for sub-network 2 is South Africa
INFO:pypsa.opf:Performed preliminary steps
INFO:pypsa.opf:Building pyomo model using `angles` formulation
INFO:pypsa.opf:Solving model using glpk
INFO:pypsa.opf:Optimization successful
Out[17]:
('ok', 'optimal')
In [18]:
network.loads_t.p
Out[18]:
Swaziland load Mozambique load South Africa load
now 250.0 650.0 42000.0
In [19]:
network.generators_t.p
Out[19]:
Swaziland Hydro Mozambique Hydro South Africa Oil South Africa Wind South Africa Gas South Africa Coal
now 600.0 1050.0 0.0 3000.0 3250.0 35000.0
In [20]:
network.links_t.p0
Out[20]:
Mozambique - Swaziland link South Africa - Swaziland link South Africa - Mozambique link
now -100.0 -250.0 -500.0
In [21]:
#print the clearing price (corresponding to hydro in S and M, and gas in SA)
network.buses_t.marginal_price
Out[21]:
Swaziland Mozambique South Africa
now 0.0 0.0 60.0
In [22]:
#link shadow prices
network.links_t.mu_lower
Out[22]:
Mozambique - Swaziland link South Africa - Swaziland link South Africa - Mozambique link
now 0.0 60.0 60.0

Single bidding zone with price-sensitive industrial load, one period

In this example we consider a single market bidding zone, South Africa.

Now there is a large industrial load with a marginal utility which is low enough to interact with the generation marginal cost.

In [17]:
country = "South Africa"

network = pypsa.Network()

network.add("Bus",country)

for tech in power_plant_p_nom[country]:
    network.add("Generator",
                "{} {}".format(country,tech),
                bus=country,
                p_nom=power_plant_p_nom[country][tech],
                marginal_cost=marginal_costs[tech])

#standard high marginal utility consumers
network.add("Load",
            "{} load".format(country),
            bus=country,
            p_set=loads[country])

#add an industrial load as a dummy negative-dispatch generator with marginal utility of 70 EUR/MWh for 8000 MW
network.add("Generator",
            "{} industrial load".format(country),
            bus=country,
            p_max_pu=0,
            p_min_pu=-1,
            p_nom=8000,
            marginal_cost=70)
In [18]:
network.lopf()
INFO:pypsa.pf:Slack bus for sub-network 0 is South Africa
INFO:pypsa.opf:Performed preliminary steps
INFO:pypsa.opf:Building pyomo model using `angles` formulation
INFO:pypsa.opf:Solving model using glpk
INFO:pypsa.opf:Optimization successful
Out[18]:
('ok', 'optimal')
In [19]:
network.loads_t.p
Out[19]:
South Africa load
now 42000.0
In [20]:
#NB only half of industrial load is served, because this maxes out 
#Gas. Oil is too expensive with a marginal cost of 80 EUR/MWh
network.generators_t.p
Out[20]:
South Africa Wind South Africa Oil South Africa Gas South Africa Coal South Africa industrial load
now 3000.0 0.0 8000.0 35000.0 -4000.0
In [22]:
network.buses_t.marginal_price
Out[22]:
Swaziland Mozambique South Africa
now 0.0 0.0 60.0

Single bidding zone with fixed load, several periods

In this example we consider a single market bidding zone, South Africa.

We consider multiple time periods (labelled [0,1,2,3]) to represent variable wind generation.

In [23]:
country = "South Africa"

network = pypsa.Network()

#snapshots labelled by [0,1,2,3]
network.set_snapshots(range(4))

network.add("Bus",country)

#p_max_pu is variable for wind
for tech in power_plant_p_nom[country]:
    network.add("Generator",
                "{} {}".format(country,tech),
                bus=country,
                p_nom=power_plant_p_nom[country][tech],
                marginal_cost=marginal_costs[tech],
                p_max_pu=([0.3,0.6,0.4,0.5] if tech == "Wind" else 1),
                )

#load which varies over the snapshots
network.add("Load",
            "{} load".format(country),
            bus=country,
            p_set=loads[country] + np.array([0,1000,3000,4000]))
In [24]:
#specify that we consider all snapshots
network.lopf(network.snapshots)
INFO:pypsa.pf:Slack bus for sub-network 0 is South Africa
INFO:pypsa.opf:Performed preliminary steps
INFO:pypsa.opf:Building pyomo model using `angles` formulation
INFO:pypsa.opf:Solving model using glpk
INFO:pypsa.opf:Optimization successful
Out[24]:
('ok', 'optimal')
In [25]:
network.loads_t.p
Out[25]:
South Africa load
0 42000.0
1 43000.0
2 45000.0
3 46000.0
In [26]:
network.generators_t.p
Out[26]:
South Africa Oil South Africa Wind South Africa Gas South Africa Coal
0 0.0 900.0 6100.0 35000.0
1 0.0 1800.0 6200.0 35000.0
2 800.0 1200.0 8000.0 35000.0
3 1500.0 1500.0 8000.0 35000.0
In [27]:
network.buses_t.marginal_price
Out[27]:
South Africa
0 60.0
1 60.0
2 80.0
3 80.0

Single bidding zone with fixed load and storage, several periods

In this example we consider a single market bidding zone, South Africa.

We consider multiple time periods (labelled [0,1,2,3]) to represent variable wind generation. Storage is allowed to do price arbitrage to reduce oil consumption.

In [28]:
country = "South Africa"

network = pypsa.Network()

#snapshots labelled by [0,1,2,3]
network.set_snapshots(range(4))

network.add("Bus",country)

#p_max_pu is variable for wind
for tech in power_plant_p_nom[country]:
    network.add("Generator",
                "{} {}".format(country,tech),
                bus=country,
                p_nom=power_plant_p_nom[country][tech],
                marginal_cost=marginal_costs[tech],
                p_max_pu=([0.3,0.6,0.4,0.5] if tech == "Wind" else 1),
                )

#load which varies over the snapshots
network.add("Load",
            "{} load".format(country),
            bus=country,
            p_set=loads[country] + np.array([0,1000,3000,4000]))

#storage unit to do price arbitrage
network.add("StorageUnit",
            "{} pumped hydro".format(country),
            bus=country,
            p_nom=1000,
            max_hours=6, #energy storage in terms of hours at full power
           )
In [29]:
network.lopf(network.snapshots)
INFO:pypsa.pf:Slack bus for sub-network 0 is South Africa
INFO:pypsa.opf:Performed preliminary steps
INFO:pypsa.opf:Building pyomo model using `angles` formulation
INFO:pypsa.opf:Solving model using glpk
INFO:pypsa.opf:Optimization successful
Out[29]:
('ok', 'optimal')
In [30]:
network.loads_t.p
Out[30]:
South Africa load
0 42000.0
1 43000.0
2 45000.0
3 46000.0
In [31]:
network.generators_t.p
Out[31]:
South Africa Oil South Africa Wind South Africa Gas South Africa Coal
0 0.0 900.0 6900.0 35000.0
1 0.0 1800.0 7200.0 35000.0
2 0.0 1200.0 8000.0 35000.0
3 500.0 1500.0 8000.0 35000.0
In [32]:
network.storage_units_t.p
Out[32]:
South Africa pumped hydro
0 -800.0
1 -1000.0
2 800.0
3 1000.0
In [33]:
network.storage_units_t.state_of_charge
Out[33]:
South Africa pumped hydro
0 800.0
1 1800.0
2 1000.0
3 0.0
In [34]:
network.buses_t.marginal_price
Out[34]:
South Africa
0 60.0
1 60.0
2 60.0
3 80.0