Note
You can download this example as a Jupyter notebook or start it in interactive mode.
Screening curve analysis
Compute the long-term equilibrium power plant investment for a given load duration curve (1000-1000z for z \(\in\) [0,1]) and a given set of generator investment options.
[1]:
import pypsa
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
%matplotlib inline
Generator marginal (m) and capital (c) costs in EUR/MWh - numbers chosen for simple answer.
[2]:
generators = {"coal" : {"m" : 2, "c" : 15},
"gas" : {"m" : 12, "c": 10},
"load-shedding" : {"m" : 1012, "c" : 0}}
The screening curve intersections are at 0.01 and 0.5.
[3]:
x = np.linspace(0,1,101)
df = pd.DataFrame({key : pd.Series(item["c"] + x*item["m"],x) for key,item in generators.items()})
df.plot(ylim=[0,50],title="Screening Curve", figsize = (9,5))
plt.tight_layout()
[4]:
n = pypsa.Network()
num_snapshots = 1001
n.snapshots = np.linspace(0,1,num_snapshots)
n.snapshot_weightings = n.snapshot_weightings/num_snapshots
n.add("Bus",name="bus")
n.add("Load",name="load",bus="bus",
p_set=1000-1000*n.snapshots.values)
for gen in generators:
n.add("Generator",name=gen,bus="bus",
p_nom_extendable=True,
marginal_cost=float(generators[gen]["m"]),
capital_cost=float(generators[gen]["c"]))
[5]:
n.loads_t.p_set.plot.area(title="Load Duration Curve", figsize = (9,5), ylabel="MW")
plt.tight_layout()
[6]:
n.lopf(solver_name="cbc")
n.objective
INFO:pypsa.opf:Performed preliminary steps
INFO:pypsa.opf:Building pyomo model using `kirchhoff` formulation
INFO:pypsa.opf:Solving model using cbc
WARNING:pyomo.solvers:Could not locate the 'cbc' executable, which is required for solver cbc
---------------------------------------------------------------------------
ApplicationError Traceback (most recent call last)
/tmp/ipykernel_391/4187261325.py in <module>
----> 1 n.lopf(solver_name="cbc")
2 n.objective
~/checkouts/readthedocs.org/user_builds/pypsa-docs-staging/envs/latest/lib/python3.7/site-packages/pypsa/components.py in lopf(self, snapshots, pyomo, solver_name, solver_options, solver_logfile, formulation, keep_files, extra_functionality, multi_investment_periods, **kwargs)
646
647 if pyomo:
--> 648 return network_lopf(self, **args)
649 else:
650 return network_lopf_lowmem(self, **args)
~/checkouts/readthedocs.org/user_builds/pypsa-docs-staging/envs/latest/lib/python3.7/site-packages/pypsa/opf.py in network_lopf(network, snapshots, solver_name, solver_io, skip_pre, extra_functionality, multi_investment_periods, solver_logfile, solver_options, keep_files, formulation, ptdf_tolerance, free_memory, extra_postprocessing)
1663 solver_logfile=solver_logfile, solver_options=solver_options,
1664 keep_files=keep_files, free_memory=free_memory,
-> 1665 extra_postprocessing=extra_postprocessing)
~/checkouts/readthedocs.org/user_builds/pypsa-docs-staging/envs/latest/lib/python3.7/site-packages/pypsa/opf.py in network_lopf_solve(network, snapshots, formulation, solver_options, solver_logfile, keep_files, free_memory, extra_postprocessing)
1563 network.results = network.opt.solve(*args, suffixes=["dual"], keepfiles=keep_files, logfile=solver_logfile, options=solver_options)
1564 else:
-> 1565 network.results = network.opt.solve(*args, suffixes=["dual"], keepfiles=keep_files, logfile=solver_logfile, options=solver_options)
1566
1567 if logger.isEnabledFor(logging.INFO):
~/checkouts/readthedocs.org/user_builds/pypsa-docs-staging/envs/latest/lib/python3.7/site-packages/pyomo/opt/base/solvers.py in solve(self, *args, **kwds)
510 """ Solve the problem """
511
--> 512 self.available(exception_flag=True)
513 #
514 # If the inputs are models, then validate that they have been
~/checkouts/readthedocs.org/user_builds/pypsa-docs-staging/envs/latest/lib/python3.7/site-packages/pyomo/opt/solver/shellcmd.py in available(self, exception_flag)
123 if exception_flag:
124 msg = "No executable found for solver '%s'"
--> 125 raise ApplicationError(msg % self.name)
126 return False
127 return True
ApplicationError: No executable found for solver 'cbc'
The capacity is set by total electricity required.
NB: No load shedding since all prices are below 10 000.
[7]:
n.generators.p_nom_opt.round(2)
[7]:
coal 0.0
gas 0.0
load-shedding 0.0
Name: p_nom_opt, dtype: float64
[8]:
n.buses_t.marginal_price.plot(title="Price Duration Curve", figsize = (9,4))
plt.tight_layout()
---------------------------------------------------------------------------
TypeError Traceback (most recent call last)
/tmp/ipykernel_391/2851249190.py in <module>
----> 1 n.buses_t.marginal_price.plot(title="Price Duration Curve", figsize = (9,4))
2 plt.tight_layout()
~/checkouts/readthedocs.org/user_builds/pypsa-docs-staging/envs/latest/lib/python3.7/site-packages/pandas/plotting/_core.py in __call__(self, *args, **kwargs)
970 data.columns = label_name
971
--> 972 return plot_backend.plot(data, kind=kind, **kwargs)
973
974 __call__.__doc__ = __doc__
~/checkouts/readthedocs.org/user_builds/pypsa-docs-staging/envs/latest/lib/python3.7/site-packages/pandas/plotting/_matplotlib/__init__.py in plot(data, kind, **kwargs)
69 kwargs["ax"] = getattr(ax, "left_ax", ax)
70 plot_obj = PLOT_CLASSES[kind](data, **kwargs)
---> 71 plot_obj.generate()
72 plot_obj.draw()
73 return plot_obj.result
~/checkouts/readthedocs.org/user_builds/pypsa-docs-staging/envs/latest/lib/python3.7/site-packages/pandas/plotting/_matplotlib/core.py in generate(self)
284 def generate(self):
285 self._args_adjust()
--> 286 self._compute_plot_data()
287 self._setup_subplots()
288 self._make_plot()
~/checkouts/readthedocs.org/user_builds/pypsa-docs-staging/envs/latest/lib/python3.7/site-packages/pandas/plotting/_matplotlib/core.py in _compute_plot_data(self)
451 # no non-numeric frames or series allowed
452 if is_empty:
--> 453 raise TypeError("no numeric data to plot")
454
455 self.data = numeric_data.apply(self._convert_to_ndarray)
TypeError: no numeric data to plot
The prices correspond either to VOLL (1012) for first 0.01 or the marginal costs (12 for 0.49 and 2 for 0.5)
Except for (infinitesimally small) points at the screening curve intersections, which correspond to changing the load duration near the intersection, so that capacity changes. This explains 7 = (12+10 - 15) (replacing coal with gas) and 22 = (12+10) (replacing load-shedding with gas).
Note: What remains unclear is what is causing :nbsphinx-math:`l `= 0… it should be 2.
[9]:
n.buses_t.marginal_price.round(2).sum(axis=1).value_counts()
[9]:
0.0 1001
dtype: int64
[10]:
n.generators_t.p.plot(ylim=[0,600],title="Generation Dispatch", figsize = (9,5))
plt.tight_layout()
---------------------------------------------------------------------------
TypeError Traceback (most recent call last)
/tmp/ipykernel_391/2429873976.py in <module>
----> 1 n.generators_t.p.plot(ylim=[0,600],title="Generation Dispatch", figsize = (9,5))
2 plt.tight_layout()
~/checkouts/readthedocs.org/user_builds/pypsa-docs-staging/envs/latest/lib/python3.7/site-packages/pandas/plotting/_core.py in __call__(self, *args, **kwargs)
970 data.columns = label_name
971
--> 972 return plot_backend.plot(data, kind=kind, **kwargs)
973
974 __call__.__doc__ = __doc__
~/checkouts/readthedocs.org/user_builds/pypsa-docs-staging/envs/latest/lib/python3.7/site-packages/pandas/plotting/_matplotlib/__init__.py in plot(data, kind, **kwargs)
69 kwargs["ax"] = getattr(ax, "left_ax", ax)
70 plot_obj = PLOT_CLASSES[kind](data, **kwargs)
---> 71 plot_obj.generate()
72 plot_obj.draw()
73 return plot_obj.result
~/checkouts/readthedocs.org/user_builds/pypsa-docs-staging/envs/latest/lib/python3.7/site-packages/pandas/plotting/_matplotlib/core.py in generate(self)
284 def generate(self):
285 self._args_adjust()
--> 286 self._compute_plot_data()
287 self._setup_subplots()
288 self._make_plot()
~/checkouts/readthedocs.org/user_builds/pypsa-docs-staging/envs/latest/lib/python3.7/site-packages/pandas/plotting/_matplotlib/core.py in _compute_plot_data(self)
451 # no non-numeric frames or series allowed
452 if is_empty:
--> 453 raise TypeError("no numeric data to plot")
454
455 self.data = numeric_data.apply(self._convert_to_ndarray)
TypeError: no numeric data to plot
Demonstrate zero-profit condition.
The total cost is given by
[11]:
(n.generators.p_nom_opt*n.generators.capital_cost
+ n.generators_t.p.multiply(n.snapshot_weightings.generators,axis=0).sum()*n.generators.marginal_cost)
[11]:
coal NaN
gas NaN
load-shedding NaN
dtype: float64
The total revenue by
[12]:
(n.generators_t.p.multiply(n.snapshot_weightings.generators, axis=0)
.multiply(n.buses_t.marginal_price["bus"], axis=0).sum(0))
---------------------------------------------------------------------------
KeyError Traceback (most recent call last)
~/checkouts/readthedocs.org/user_builds/pypsa-docs-staging/envs/latest/lib/python3.7/site-packages/pandas/core/indexes/base.py in get_loc(self, key, method, tolerance)
3360 try:
-> 3361 return self._engine.get_loc(casted_key)
3362 except KeyError as err:
~/checkouts/readthedocs.org/user_builds/pypsa-docs-staging/envs/latest/lib/python3.7/site-packages/pandas/_libs/index.pyx in pandas._libs.index.IndexEngine.get_loc()
~/checkouts/readthedocs.org/user_builds/pypsa-docs-staging/envs/latest/lib/python3.7/site-packages/pandas/_libs/index.pyx in pandas._libs.index.IndexEngine.get_loc()
pandas/_libs/hashtable_class_helper.pxi in pandas._libs.hashtable.PyObjectHashTable.get_item()
pandas/_libs/hashtable_class_helper.pxi in pandas._libs.hashtable.PyObjectHashTable.get_item()
KeyError: 'bus'
The above exception was the direct cause of the following exception:
KeyError Traceback (most recent call last)
/tmp/ipykernel_391/3942873135.py in <module>
1 (n.generators_t.p.multiply(n.snapshot_weightings.generators, axis=0)
----> 2 .multiply(n.buses_t.marginal_price["bus"], axis=0).sum(0))
~/checkouts/readthedocs.org/user_builds/pypsa-docs-staging/envs/latest/lib/python3.7/site-packages/pandas/core/frame.py in __getitem__(self, key)
3456 if self.columns.nlevels > 1:
3457 return self._getitem_multilevel(key)
-> 3458 indexer = self.columns.get_loc(key)
3459 if is_integer(indexer):
3460 indexer = [indexer]
~/checkouts/readthedocs.org/user_builds/pypsa-docs-staging/envs/latest/lib/python3.7/site-packages/pandas/core/indexes/base.py in get_loc(self, key, method, tolerance)
3361 return self._engine.get_loc(casted_key)
3362 except KeyError as err:
-> 3363 raise KeyError(key) from err
3364
3365 if is_scalar(key) and isna(key) and not self.hasnans:
KeyError: 'bus'
Now, take the capacities from the above long-term equilibrium, then disallow expansion.
Show that the resulting market prices are identical.
This holds in this example, but does NOT necessarily hold and breaks down in some circumstances (for example, when there is a lot of storage and inter-temporal shifting).
[13]:
n.generators.p_nom_extendable = False
n.generators.p_nom = n.generators.p_nom_opt
[14]:
n.lopf();
INFO:pypsa.opf:Performed preliminary steps
INFO:pypsa.opf:Building pyomo model using `kirchhoff` formulation
INFO:pypsa.opf:Solving model using glpk
WARNING:pyomo.solvers:Could not locate the 'glpsol' executable, which is required for solver 'glpk'
---------------------------------------------------------------------------
ApplicationError Traceback (most recent call last)
/tmp/ipykernel_391/1812968345.py in <module>
----> 1 n.lopf();
~/checkouts/readthedocs.org/user_builds/pypsa-docs-staging/envs/latest/lib/python3.7/site-packages/pypsa/components.py in lopf(self, snapshots, pyomo, solver_name, solver_options, solver_logfile, formulation, keep_files, extra_functionality, multi_investment_periods, **kwargs)
646
647 if pyomo:
--> 648 return network_lopf(self, **args)
649 else:
650 return network_lopf_lowmem(self, **args)
~/checkouts/readthedocs.org/user_builds/pypsa-docs-staging/envs/latest/lib/python3.7/site-packages/pypsa/opf.py in network_lopf(network, snapshots, solver_name, solver_io, skip_pre, extra_functionality, multi_investment_periods, solver_logfile, solver_options, keep_files, formulation, ptdf_tolerance, free_memory, extra_postprocessing)
1663 solver_logfile=solver_logfile, solver_options=solver_options,
1664 keep_files=keep_files, free_memory=free_memory,
-> 1665 extra_postprocessing=extra_postprocessing)
~/checkouts/readthedocs.org/user_builds/pypsa-docs-staging/envs/latest/lib/python3.7/site-packages/pypsa/opf.py in network_lopf_solve(network, snapshots, formulation, solver_options, solver_logfile, keep_files, free_memory, extra_postprocessing)
1563 network.results = network.opt.solve(*args, suffixes=["dual"], keepfiles=keep_files, logfile=solver_logfile, options=solver_options)
1564 else:
-> 1565 network.results = network.opt.solve(*args, suffixes=["dual"], keepfiles=keep_files, logfile=solver_logfile, options=solver_options)
1566
1567 if logger.isEnabledFor(logging.INFO):
~/checkouts/readthedocs.org/user_builds/pypsa-docs-staging/envs/latest/lib/python3.7/site-packages/pyomo/opt/base/solvers.py in solve(self, *args, **kwds)
510 """ Solve the problem """
511
--> 512 self.available(exception_flag=True)
513 #
514 # If the inputs are models, then validate that they have been
~/checkouts/readthedocs.org/user_builds/pypsa-docs-staging/envs/latest/lib/python3.7/site-packages/pyomo/opt/solver/shellcmd.py in available(self, exception_flag)
123 if exception_flag:
124 msg = "No executable found for solver '%s'"
--> 125 raise ApplicationError(msg % self.name)
126 return False
127 return True
ApplicationError: No executable found for solver 'glpk'
[15]:
n.buses_t.marginal_price.plot(title="Price Duration Curve", figsize = (9,5))
plt.tight_layout()
---------------------------------------------------------------------------
TypeError Traceback (most recent call last)
/tmp/ipykernel_391/3992284789.py in <module>
----> 1 n.buses_t.marginal_price.plot(title="Price Duration Curve", figsize = (9,5))
2 plt.tight_layout()
~/checkouts/readthedocs.org/user_builds/pypsa-docs-staging/envs/latest/lib/python3.7/site-packages/pandas/plotting/_core.py in __call__(self, *args, **kwargs)
970 data.columns = label_name
971
--> 972 return plot_backend.plot(data, kind=kind, **kwargs)
973
974 __call__.__doc__ = __doc__
~/checkouts/readthedocs.org/user_builds/pypsa-docs-staging/envs/latest/lib/python3.7/site-packages/pandas/plotting/_matplotlib/__init__.py in plot(data, kind, **kwargs)
69 kwargs["ax"] = getattr(ax, "left_ax", ax)
70 plot_obj = PLOT_CLASSES[kind](data, **kwargs)
---> 71 plot_obj.generate()
72 plot_obj.draw()
73 return plot_obj.result
~/checkouts/readthedocs.org/user_builds/pypsa-docs-staging/envs/latest/lib/python3.7/site-packages/pandas/plotting/_matplotlib/core.py in generate(self)
284 def generate(self):
285 self._args_adjust()
--> 286 self._compute_plot_data()
287 self._setup_subplots()
288 self._make_plot()
~/checkouts/readthedocs.org/user_builds/pypsa-docs-staging/envs/latest/lib/python3.7/site-packages/pandas/plotting/_matplotlib/core.py in _compute_plot_data(self)
451 # no non-numeric frames or series allowed
452 if is_empty:
--> 453 raise TypeError("no numeric data to plot")
454
455 self.data = numeric_data.apply(self._convert_to_ndarray)
TypeError: no numeric data to plot
[16]:
n.buses_t.marginal_price.sum(axis=1).value_counts()
[16]:
0.0 1001
dtype: int64
Demonstrate zero-profit condition. Differences are due to singular times, see above, not a problem
Total costs
[17]:
(n.generators.p_nom * n.generators.capital_cost
+ n.generators_t.p.multiply(n.snapshot_weightings.generators, axis=0).sum() * n.generators.marginal_cost)
[17]:
coal NaN
gas NaN
load-shedding NaN
dtype: float64
Total revenue
[18]:
(n.generators_t.p.multiply(n.snapshot_weightings.generators, axis=0)
.multiply(n.buses_t.marginal_price["bus"],axis=0).sum())
---------------------------------------------------------------------------
KeyError Traceback (most recent call last)
~/checkouts/readthedocs.org/user_builds/pypsa-docs-staging/envs/latest/lib/python3.7/site-packages/pandas/core/indexes/base.py in get_loc(self, key, method, tolerance)
3360 try:
-> 3361 return self._engine.get_loc(casted_key)
3362 except KeyError as err:
~/checkouts/readthedocs.org/user_builds/pypsa-docs-staging/envs/latest/lib/python3.7/site-packages/pandas/_libs/index.pyx in pandas._libs.index.IndexEngine.get_loc()
~/checkouts/readthedocs.org/user_builds/pypsa-docs-staging/envs/latest/lib/python3.7/site-packages/pandas/_libs/index.pyx in pandas._libs.index.IndexEngine.get_loc()
pandas/_libs/hashtable_class_helper.pxi in pandas._libs.hashtable.PyObjectHashTable.get_item()
pandas/_libs/hashtable_class_helper.pxi in pandas._libs.hashtable.PyObjectHashTable.get_item()
KeyError: 'bus'
The above exception was the direct cause of the following exception:
KeyError Traceback (most recent call last)
/tmp/ipykernel_391/4015668699.py in <module>
1 (n.generators_t.p.multiply(n.snapshot_weightings.generators, axis=0)
----> 2 .multiply(n.buses_t.marginal_price["bus"],axis=0).sum())
~/checkouts/readthedocs.org/user_builds/pypsa-docs-staging/envs/latest/lib/python3.7/site-packages/pandas/core/frame.py in __getitem__(self, key)
3456 if self.columns.nlevels > 1:
3457 return self._getitem_multilevel(key)
-> 3458 indexer = self.columns.get_loc(key)
3459 if is_integer(indexer):
3460 indexer = [indexer]
~/checkouts/readthedocs.org/user_builds/pypsa-docs-staging/envs/latest/lib/python3.7/site-packages/pandas/core/indexes/base.py in get_loc(self, key, method, tolerance)
3361 return self._engine.get_loc(casted_key)
3362 except KeyError as err:
-> 3363 raise KeyError(key) from err
3364
3365 if is_scalar(key) and isna(key) and not self.hasnans:
KeyError: 'bus'